hydro/nudging/module_nudging_io.F

   1 !  Program Name:
   2 !  Author(s)/Contact(s):
   3 !  Abstract:
   4 !  History Log:
   5 !
   6 !  Usage:
   7 !  Parameters: <Specify typical arguments passed>
   8 !  Input Files:
   9 !        <list file names and briefly describe the data they include>
  10 !  Output Files:
  11 !        <list file names and briefly describe the information they include>
  12 !
  13 !  Condition codes:
  14 !        <list exit condition or error codes returned >
  15 !        If appropriate, descriptive troubleshooting instructions or
  16 !        likely causes for failures could be mentioned here with the
  17 !        appropriate error code
  18 !
  19 !  User controllable options: <if applicable>
  20
  21 module module_nudging_io
  22
  23 use netcdf
  24 use config_base, only: nlst
  25 use module_RT_data,  only: rt_domain
  26 use module_hydro_stop, only: HYDRO_stop
  27
  28 implicit none
  29 #include <netcdf.inc>
  30
  31 !========================
  32 ! lastObs structure, corresponding to nudgingLastObs.YYYY-mm-dd_HH:MM:ss.nc
  33 ! How observations from the past are carried forward.
  34 ! This type is extended in module_stream_nudging
  35 type lastObsStructure
  36    character(len=15)                            :: usgsId
  37    real,              allocatable, dimension(:) :: lastObsDischarge
  38    real,              allocatable, dimension(:) :: lastObsModelDischarge
  39    character(len=19), allocatable, dimension(:) :: lastObsTime
  40    real,              allocatable, dimension(:) :: lastObsQuality
  41 end type lastObsStructure
  42
  43 type lastObsStructure_SoA
  44    character(len=15), allocatable, dimension(:)   :: usgsId
  45    real,              allocatable, dimension(:,:) :: lastObsDischarge
  46    real,              allocatable, dimension(:,:) :: lastObsModelDischarge
  47    character(len=19), allocatable, dimension(:,:) :: lastObsTime
  48    real,              allocatable, dimension(:,:) :: lastObsQuality
  49 end type lastObsStructure_SoA
  50
  51 integer, parameter :: did = 1
  52
  53 contains
  54
  55 !===================================================================================================
  56 ! Program Name:
  57 !   subroutine find_timeslice_file
  58 ! Author(s)/Contact(s):
  59 !   James L McCreight <jamesmcc><ucar><edu>
  60 ! Abstract:
  61 !   Return a single file path/names of timeslice files given a single date.
  62 ! History Log:
  63 !   6/4/15 - Created,
  64 ! Usage:
  65 ! Parameters:
  66 !   q
  67 ! Input Files:
  68 ! Output Files:
  69 ! Condition codes:
  70 ! User controllable options: namelist option for nlst_rt(did)%timeSlicePath
  71 ! Notes:
  72
  73 function find_timeslice_file(date, obsResolution)
  74 implicit none
  75 character(len=256) :: find_timeslice_file      ! Output
  76 character(len=19), intent(in) :: date          ! Input
  77 character(len=2),  intent(in) :: obsResolution ! Input
  78 !Internals
  79 character(len=256) :: tmpTimeSlice
  80 logical :: fileExists
  81
  82 #ifdef HYDRO_D
  83 print*,'Ndg: start find_timeslice_file'
  84 call flush(6)
  85 #endif
  86
  87 find_timeslice_file=''
  88 ! is there a file with this name?
  89 ! note files are not resolved below minutes.
  90 tmpTimeSlice = trim(nlst(did)%timeSlicePath) // '/' //date // "." // &
  91                     obsResolution // 'min.usgsTimeSlice.ncdf'
  92 inquire(FILE=tmpTimeSlice, EXIST=fileExists)
  93 if (fileExists) find_timeslice_file = trim(tmpTimeSlice)
  94
  95 #ifdef HYDRO_D
  96 print*,'Ndg: timeSlice file: ', tmpTimeSlice
  97 print*,'Ndg: file found: ', fileExists
  98 print*,'Ndg: finish find_timeslice_file'
  99 call flush(6)
 100 #endif
 101
 102 end function find_timeslice_file
 103
 104 !===================================================================================================
 105 ! Program Name:
 106 !   subroutine read_timeslice_file
 107 ! Author(s)/Contact(s):
 108 !   James L McCreight <jamesmcc><ucar><edu>
 109 ! Abstract:
 110 !   Return the contents of a timeslice file.
 111 ! History Log:
 112 !   6/4/15 - Created,
 113 ! Usage:
 114 ! Parameters:
 115 !
 116 ! Input Files:
 117 ! Output Files:
 118 ! Condition codes:
 119 ! User controllable options: namelist option for nlst_rt(did)%timeSlicePath
 120 ! Notes:
 121
 122 subroutine read_timeslice_file(  &
 123      timeSliceFile,              &  !! IN: char, file path/name
 124      sanityQcDischarge,          &  !! IN: perform sanity check qc?
 125      sliceTime,                  &  !! OUT: the file time.
 126      updateTime,                 &  !! OUT: the time the file was updated.
 127      sliceResolution,            &  !! OUT: temporal resolution of the file in minutes.
 128      gageId,                     &  !! OUT: integer, the USGS gage ID.
 129      gageTime,                   &  !! OUT: output converted for comparision to model time??
 130      gageQC,                     &  !! OUT: quality control flag for discharge
 131      gageDischarge,              &  !! OUT: real, the m3/s observed discharge
 132      fatalErr,                   &  !! IN:  logical, are IO errors fatal?
 133      errStatus                   &  !! OUT: count of errors encountered
 134                                  )
 135
 136 implicit none
 137 character(len=256), intent(in)                :: timeSliceFile
 138 logical,            intent(in)                :: sanityQcDischarge
 139 character(len=19),  intent(out)               :: sliceTime
 140 character(len=19),  intent(out)               :: updateTime
 141 character(len=2),   intent(out)               :: sliceResolution
 142 character(len=15),  intent(out), dimension(:) :: gageId
 143 character(len=19),  intent(out), dimension(:) :: gageTime
 144 real,               intent(out), dimension(:) :: gageQC
 145 real,               intent(out), dimension(:) :: gageDischarge
 146 logical,            intent(in)                :: fatalErr
 147 integer,            intent(out)               :: errStatus
 148
 149 integer*2, allocatable, dimension(:) :: gageQCIn
 150 integer :: iRet, ncId, fileNameLen, errStatusOut
 151 character(len=19) :: caller
 152 real, parameter :: zero = 0.0000000000000
 153 real, parameter :: oneHundred = 100.0000000000000
 154
 155 #ifdef HYDRO_D
 156 print*,'Ndg: start read_timeslice_file'
 157 write(*,'("Ndg: read_timeslice_file: ''", A, "''")') trim(timeSliceFile)
 158 call flush(6)
 159 #endif
 160
 161 caller = 'read_timeslice_file'
 162 errStatus=0
 163
 164 iRet = nf90_open(trim(timeSliceFile), nf90_nowrite, ncid)
 165 if (iRet /= nf90_NoErr) then
 166    write(*,'("Problem opening timeSliceFile: ''", A, "''")') trim(timeSliceFile)
 167    if(fatalErr) call hydro_stop('read_timeslice_file')
 168    errStatus=errStatus+1
 169 endif
 170
 171 !! global atts
 172 iRet = nf90_get_att(ncid, nf90_global, 'sliceCenterTimeUTC',         sliceTime)
 173 if (iRet /= nf90_NoErr) errStatus=errStatus+1
 174 iRet = nf90_get_att(ncid, nf90_global, 'fileUpdateTimeUTC',          updateTime)
 175 if (iRet /= nf90_NoErr) errStatus=errStatus+1
 176 iRet = nf90_get_att(ncid, nf90_global, 'sliceTimeResolutionMinutes', sliceResolution)
 177 if (iRet /= nf90_NoErr) errStatus=errStatus+1
 178
 179 ! variables
 180 call get_1d_netcdf_text(ncid, 'stationId', gageId,        caller, &
 181                         fatalErr, errStatusOut)
 182 errStatus=errStatus+errStatusOut
 183 call get_1d_netcdf_text(ncid, 'time',      gageTime,      caller, &
 184                         fatalErr, errStatusOut)
 185 errStatus=errStatus+errStatusOut
 186 call get_1d_netcdf_real(ncid, 'discharge', gageDischarge, caller, &
 187                         fatalErr, errStatusOut)
 188 errStatus=errStatus+errStatusOut
 189
 190 !! the quality short integer needs scaled
 191 allocate(gageQCIn(size(gageQC)))
 192 call get_1d_netcdf_int2(ncid, 'discharge_quality', gageQCIn,        caller, &
 193                         fatalErr, errStatusOut)
 194 errStatus=errStatus+errStatusOut
 195 gageQC=gageQCIn/oneHundred
 196 deallocate(gageQCIn)
 197
 198 !! JLM TODO fix temporary hardwire QC sanity checks
 199 if (sanityQcDischarge) then
 200    !! First, quality check on the quality flag!
 201    where(gageQC .lt. 0 .or. gageQC .gt. 1) gageQC=0
 202    !! Now flow QC.
 203    where(gageDischarge .le. 0.000) gageQC=0
 204    !! peak flow on MS river *2
 205    !! http://nwis.waterdata.usgs.gov/nwis/peak?site_no=07374000&agency_cd=USGS&format=html
 206    !! baton rouge 1945: 1,473,000cfs=41,711cms, multiply it roughly by 2
 207    where(gageDischarge .ge. 90000.0) gageQC=0
 208    if(any(gageQc .eq. 0)) then
 209       write(6,*) 'Ndg: some gageQC set to zero'
 210    endif
 211 endif
 212
 213 iRet = nf90_close(ncId)
 214 if (iRet /= nf90_NoErr) then
 215    write(*,'("Problem closing timeSliceFile: ''", A, "''")') trim(timeSliceFile)
 216    if(fatalErr) call hydro_stop('read_timeslice_file')
 217    errStatus=errStatus+1
 218 endif
 219
 220 #ifdef HYDRO_D
 221 print*,'Ndg: finish read_timeslice_file'
 222 call flush(6)!
 223 #endif
 224
 225 end subroutine read_timeslice_file
 226
 227 !===================================================================================================
 228 ! Program Name:
 229 !   subroutine read_reach_gage_collocation
 230 ! Author(s)/Contact(s):
 231 !   James L McCreight <jamesmcc><ucar><edu>
 232 ! Abstract:
 233 !   Read the gages column from the "RouteLink.nc" netcdf file specifing the channel topology.
 234 ! History Log:
 235 !   7/20/15 -Created, JLM.
 236 ! Usage:
 237 ! Parameters:
 238 ! Input Files:
 239 !   netcdf file RouteLink.nc or other name ending with .nc.
 240 ! Output Files:
 241 ! Condition codes:
 242 ! User controllable options:
 243 ! Notes: Currently no csv support, but basic code is in place though commented out.
 244
 245 subroutine read_reach_gage_collocation(gageIds)
 246   use config_base, only: nlst
 247 implicit none
 248 character(len=15), dimension(:), intent(out) :: gageIds
 249 !integer, dimension(:), intent(out) :: gageIds
 250
 251 integer               :: ncid , iRet, varId
 252 logical               :: fatal_if_error
 253 character(len=256) :: route_link_f_r, route_link_f, varName
 254 integer :: lenRouteLinkFR ! so the preceeding var chan be changed without changing code
 255 logical :: routeLinkNetcdf
 256
 257 #ifdef HYDRO_D
 258 print*,'Ndg: Starting read_reach_gage_collocation'
 259 #endif
 260
 261 varName = 'gages'
 262 !! is RouteLink file netcdf (*.nc) or csv (*.csv)
 263 route_link_f   = nlst(did)%route_link_f
 264 route_link_f_r = adjustr(route_link_f)
 265 lenRouteLinkFR = len(route_link_f_r)
 266 routeLinkNetcdf = route_link_f_r( (lenRouteLinkFR-2):lenRouteLinkFR) .eq. '.nc'
 267
 268 if(routeLinkNetcdf) then
 269
 270 !! part of this could become a get_1d_netcdf
 271    iRet = nf90_open(trim(route_link_f), nf90_NoWrite, ncId)
 272    if (iRet /= nf90_NoErr) then
 273       write(6,'("read_reach_gage_collocation: Problem opening file ''", A, "''")') trim(route_link_f)
 274       call hydro_stop("read_reach_gage_collocation")
 275    endif
 276
 277    iRet = nf90_inq_varid(ncid, varName, varid)
 278    if (iRet /= nf90_NoErr) then
 279       print*,"Ndg: read_reach_gage_collocation: variable: " // trim(varName)
 280       call hydro_stop("read_reach_gage_collocation")
 281    end if
 282
 283    iRet = nf90_get_var(ncid, varid, gageIds)
 284    if (iRet /= nf90_NoErr) then
 285       print*,"Ndg: read_reach_gage_collocation: value: " // trim(varName)
 286       call hydro_stop("read_reach_gage_collocation")
 287    end if
 288
 289    iRet = nf90_close(ncid)
 290    if (iRet /= nf90_NoErr) then
 291       print*,"Ndg: read_reach_gage_collocation: closing: " // trim(varName)
 292       call hydro_stop("read_reach_gage_collocation")
 293    end if
 294
 295 else
 296    call hydro_stop('read_reach_gage_collocation: csv not currently supported.')
 297    !! here's some code to start from if/when we do want to support csv for this.
 298    !open(unit=79,file=trim(route_link_f),form='formatted',status='old')
 299    !read(79,*)  header
 300    !print *, "header ", header, "NLINKSL = ", NLINKSL, GNLINKSL
 301    !call flush(6)
 302    !do i=1,GNLINKSL
 303    !   read (79,*) tmpLINKID(i),   tmpFROM_NODE(i), tmpTO_NODE(i), tmpCHLON(i),    &
 304    !        tmpCHLAT(i),    tmpZELEV(i),     tmpTYPEL(i),   tmpORDER(i),    &
 305    !        tmpQLINK(i,1),  tmpMUSK(i),      tmpMUSX(i),    tmpCHANLEN(i),  &
 306    !        tmpMannN(i),    tmpSo(i),        tmpChSSlp(i),  tmpBw(i),       &
 307    !        tmpChannK(i),                                                   &
 308    !        tmpHRZAREA(i),  tmpLAKEMAXH(i),  tmpWEIRC(i),   tmpWEIRL(i),    &
 309    !        tmpORIFICEC(i), tmpORIFICEA(i),  tmpORIFICEE(i)
 310    !   ! if (So(i).lt.0.005) So(i) = 0.005  !-- impose a minimum slope requireement
 311    !   if (tmpORDER(i) .gt. MAXORDER) MAXORDER = tmpORDER(i)
 312    !end do
 313    !close(79)
 314 end if ! routeLinkNetcdf
 315
 316 #ifdef HYDRO_D
 317 print*,'Ndg: Finish read_reach_gage_collocation'
 318 call flush(6)
 319 #endif
 320
 321 end subroutine read_reach_gage_collocation
 322
 323 !===================================================================================================
 324 ! Program Name:
 325 !   subroutine read_gridded_nudging_frxst_gage_csv
 326 ! Author(s)/Contact(s):
 327 !   James L McCreight <jamesmcc><ucar><edu>
 328 ! Abstract:
 329 !   For gridded channel routine, the mechanism is to use forecast points to
 330 !   identify gages. This file determines the collocation of the two. See the
 331 !   file in the nudging/ directory for more info.
 332 ! History Log:
 333 !   6/22/15 -Created, JLM.
 334 ! Usage:
 335 ! Parameters:
 336 ! Input Files: Nudging_frxst_gage.csv
 337 ! Output Files:
 338 ! Condition codes:
 339 ! User controllable options:
 340 ! Notes:
 341
 342 subroutine read_gridded_nudging_frxst_gage_csv(frxstId, gageId, nGages)
 343 implicit none
 344 integer,           dimension(:), intent(out) :: frxstId
 345 character(len=15), dimension(:), intent(out) :: gageId
 346 integer,                         intent(out) :: nGages
 347
 348 character(len=5000) :: line
 349 integer :: ii, colCount, lineCount, lastCommaPos, tmpInt
 350
 351 #ifdef HYDRO_D
 352 print*,'Ndg: start read_gridded_nudging_frxst_gage_csv'
 353 #endif
 354 frxstId=-9999
 355
 356 lineCount = 0
 357
 358 #ifndef NCEP_WCOSS
 359 open(117, file = 'Nudging_frxst_gage.csv' )
 360 #else
 361 open(26)
 362 #endif
 363 do
 364 #ifndef NCEP_WCOSS
 365    read(117, '( a )', end = 200 ) line
 366 #else
 367    read(26, '( a )', end = 200 ) line
 368 #endif
 369    if( line( 1:1 ) == '!' ) cycle
 370    colCount = 1
 371    lastCommaPos = 0
 372    do ii = 1, len( line )
 373
 374       if( line( ii:ii ) == '!' ) exit
 375
 376       if (ii .eq. 1) lineCount = lineCount + 1
 377
 378       if( line( ii:ii ) == ',' ) then
 379          if(colCount .eq. 1) then
 380             read (line((lastCommaPos+1):(ii-1)),'(I5)') tmpInt
 381             frxstId(lineCount) = tmpInt
 382          end if
 383          if(colCount .eq. 2) gageId(lineCount)  =  trim(line((lastCommaPos+1):(ii-1)))
 384          colCount = colCount + 1
 385          lastCommaPos = ii
 386       end if
 387
 388       if(colCount .eq. 3) cycle
 389
 390    end do
 391 end do
 392
 393 #ifndef NCEP_WCOSS
 394 close(117)
 395 #else
 396 close(26)
 397 #endif
 398
 399 200 continue
 400
 401 #ifdef HYDRO_D
 402 nGages = lineCount
 403 print*,"Ndg: nGages: ",  nGages
 404 print*,"Ndg: frxstId: ", frxstId(1:nGages)
 405 print*,"Ndg: gageId:",   gageId(1:nGages)
 406 print*,'Ndg: finish read_gridded_nudging_frxst_gage_csv'
 407 call flush(6)
 408 #endif
 409
 410 end subroutine read_gridded_nudging_frxst_gage_csv
 411
 412 !===================================================================================================
 413 ! Program Name:
 414 !   subroutine read_nudging_param_file
 415 ! Author(s)/Contact(s):
 416 !   James L McCreight <jamesmcc><ucar><edu>
 417 ! Abstract:
 418 !   Set up the default nudging parameters with the model initialization.
 419 ! History Log:
 420 !   6/22/15 -Created, JLM.
 421 ! Usage:
 422 ! Parameters:
 423 ! Input Files: Nudging_frxst_gage.csv
 424 ! Output Files:
 425 ! Condition codes:
 426 ! User controllable options:
 427 ! Notes:
 428
 429 subroutine read_nudging_param_file(  &
 430      paramFile,                      &  !! IN: char, file path/name
 431      gageId,                         &  !! OUT: integer, the USGS gage ID.
 432      gageR,                          &  !! OUT: output converted for comparision to model time??
 433      gageG,                          &  !! OUT: quality control flag for discharge
 434      gageTau,                        &  !! OUT: real the half-window for NN/IDW interp
 435      gageQThresh,                    &  !! OUT: real the half-window for NN/IDW interp
 436      gageExpCoeff                    &  !! OUT: real, params for exponent-based temporal persistence
 437 )
 438 implicit none
 439 character(len=256), intent(in)                    :: paramFile
 440 character(len=15),  intent(out), dimension(:)     :: gageId
 441 real,               intent(out), dimension(:)     :: gageR
 442 real,               intent(out), dimension(:)     :: gageG
 443 real,               intent(out), dimension(:)     :: gageTau
 444 real,               intent(out), dimension(:,:,:) :: gageQThresh
 445 real,               intent(out), dimension(:,:,:) :: gageExpCoeff
 446
 447 integer :: iRet, ncId, varId, dimId, fileNameLen, errStatus
 448 character(len=23) :: caller
 449
 450 #ifdef HYDRO_D
 451 print*,'Ndg: start read_nudging_param_file'
 452 write(*,'("Ndg: paramFile: ''", A, "''")') trim(paramFile)
 453 call flush(6)
 454 #endif
 455
 456 caller='read_nudging_param_file'
 457
 458 iRet = nf90_open(trim(paramFile), nf90_NoWrite, ncid)
 459 if (iRet /= nf90_NoErr) then
 460    write(6,'("read_nudging_param_file: Problem opening file ''", A, "''")') trim(paramFile)
 461    call hydro_stop("read_nudging_param_file")
 462 endif
 463
 464 call get_1d_netcdf_text(ncid, 'stationId', gageId,  caller, .TRUE., errStatus)
 465 call get_1d_netcdf_real(ncid, 'R',         gageR,   caller, .TRUE., errStatus)
 466 call get_1d_netcdf_real(ncid, 'G',         gageG,   caller, .TRUE., errStatus)
 467 call get_1d_netcdf_real(ncid, 'tau',       gageTau, caller, .TRUE., errStatus)
 468
 469 if(size(gageExpCoeff,2) .ne. 0 .and. size(gageExpCoeff,3) .ne. 0) then
 470    call get_3d_netcdf_real(ncid, 'qThresh',  gageQThresh,  caller, .true., errStatus)
 471    call get_3d_netcdf_real(ncid, 'expCoeff', gageExpCoeff, caller, .true., errStatus)
 472 end if
 473
 474 iRet = nf90_close(ncId)
 475 if (iRet /= nf90_NoErr) then
 476    write(*,'("Problem closing paramFile: ''", A, "''")') trim(paramFile)
 477    call hydro_stop('read_nudging_param_file')
 478 endif
 479
 480 #ifdef HYDRO_D
 481 print*,'Ndg: finish read_nudging_param_file'
 482 call flush(6)
 483 #endif
 484
 485 end subroutine read_nudging_param_file
 486
 487
 488 !===================================================================================================
 489 ! Subroutine Name:
 490 !   subroutine write_nwis_not_in_RLAndParams
 491 ! Author(s)/Contact(s):
 492 !   James L McCreight <jamesmcc><ucar><edu>
 493 ! Abstract:
 494 !   Write a log file of gages supplied by NWIS but which are not found in
 495 !   the intersection of our Route_Link gages and our parameter file gages.
 496 !   Perhaps they were not picked up by us as available,
 497 !   or they are simply not in NHD+.
 498 ! History Log:
 499 !   11/04/15 -Created, JLM.
 500 ! Usage:
 501 ! Parameters:
 502 ! Input Files:
 503 ! Output Files:
 504 ! Condition codes:
 505 ! User controllable options: None.
 506 ! Notes:
 507
 508 subroutine write_nwis_not_in_RLAndParams( gageId, count )
 509 implicit none
 510 character(len=15), dimension(:),  intent(in) :: gageId
 511 integer,                          intent(in) :: count
 512 integer :: cc
 513
 514 #ifndef NCEP_WCOSS
 515 open (unit=919,file='NWIS_gages_not_in_RLAndParams.txt',status='unknown',position='append')
 516 #else
 517 open (unit=27,status='unknown',position='append')
 518 #endif
 519
 520 do cc=1,count
 521
 522 #ifndef NCEP_WCOSS
 523    write(919,'(A15)') gageId(cc)
 524 #else
 525    write(27,'(A15)') gageId(cc)
 526 #endif
 527
 528 end do
 529
 530 #ifndef NCEP_WCOSS
 531 close(919)
 532 #else
 533 close(27)
 534 #endif
 535
 536 end subroutine write_nwis_not_in_RLAndParams
 537
 538
 539 !===================================================================================================
 540 ! Subroutine Name:
 541 !   subroutine write_nudging_last_obs
 542 ! Author(s)/Contact(s):
 543 !   James L McCreight <jamesmcc><ucar><edu>
 544 ! Abstract:
 545 !   Write out the last observations collected over time.
 546 ! History Log:
 547 !   02/03/16 -Created, JLM.
 548 ! Usage:
 549 ! Parameters:
 550 ! Input Files:
 551 ! Output Files:
 552 ! Condition codes:
 553 ! User controllable options: None.
 554 ! Notes:  Needs better error handling...
 555
 556 subroutine write_nudging_last_obs(lastObsStr, modelTime, g_nudge)
 557
 558 implicit none
 559
 560 !Arugments
 561 type(lastObsStructure), intent(in), dimension(:) :: lastObsStr
 562 character(len=19),      intent(in)               :: modelTime
 563 real,                   intent(in), dimension(:) :: g_nudge
 564
 565 ! Local
 566 character(len=37) :: outFileName
 567 integer :: nSpace, nTime, nFeature, tt, ss
 568 integer :: ncid, iret, varid, dimIdTimeStr, dimIdStnStr, dimIdNTime, dimIdNStn, dimIdNFeature
 569 character(len=15), allocatable, dimension(:) :: charArr15
 570 character(len=19), allocatable, dimension(:,:) :: charArr19
 571 real,      allocatable, dimension(:,:) :: tmpFloat
 572 integer*2, allocatable, dimension(:,:) :: qualityInt2
 573
 574 #ifdef HYDRO_D
 575 print*,'Ndg: start write_nudging_last_obs'
 576 flush(6)
 577 #endif
 578
 579 nSpace = size(lastObsStr)
 580 nTime  = size(lastObsStr(1)%lastObsDischarge)
 581 nFeature = size(g_nudge)
 582
 583 !           !1234567891123456789212345678931234567
 584 outFileName='nudgingLastObs.' // modelTime // '.nc'
 585
 586 ! create file
 587 iret = nf90_create(outFileName, nf90_clobber, ncid)
 588
 589 ! create dimensions
 590 ! station id has length of the number of gages
 591 ! feature id has length of the number of reaches/features.
 592 iret = nf90_def_dim(ncid, "timeStrLen",      19,             dimIdTimeStr)
 593 iret = nf90_def_dim(ncid, "timeInd",         nTime,          dimIdNTime)
 594 iret = nf90_def_dim(ncid, "stationIdStrLen", 15,             dimIdStnStr)
 595 iret = nf90_def_dim(ncid, "stationIdInd",    nf90_unlimited, dimIdNStn)
 596 iret = nf90_def_dim(ncid, "feature_id",      nFeature,       dimIdNFeature)
 597
 598 !! gageId def var
 599 iret = nf90_def_var(ncid, "stationId", nf90_char, (/ dimIdStnStr, dimIdNstn /), varid)
 600 iret = nf90_put_att(ncid, varid, 'long_name', "USGS station identifer of length 15")
 601 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 602
 603 !! time def var
 604 iret = nf90_def_var(ncid, "time", nf90_char, (/ dimIdTimeStr, dimIdNTime, dimIdNStn /), varid)
 605 iret = nf90_put_att(ncid, varid, 'units', "UTC")
 606 iret = nf90_put_att(ncid, varid, 'long_name', "YYYY-mm-dd_HH:MM:SS UTC")
 607 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 608
 609 ! discharge def var
 610 iret = nf90_def_var(ncid, "discharge", nf90_float, (/ dimIdNTime, dimIdNstn /), varid)
 611 iret = nf90_put_att(ncid, varid, 'units', "m^3/s")
 612 iret = nf90_put_att(ncid, varid, 'long_name', "Discharge.cubic_meters_per_second")
 613 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 614
 615 ! model discharge def var
 616 iret = nf90_def_var(ncid, "model_discharge", nf90_float, (/ dimIdNTime, dimIdNstn /), varid)
 617 iret = nf90_put_att(ncid, varid, 'units', "m^3/s")
 618 iret = nf90_put_att(ncid, varid, 'long_name', "modelDischarge.cubic_meters_per_second")
 619 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 620
 621 ! discharge quality def var
 622 iret = nf90_def_var(ncid, "discharge_quality", nf90_short, (/ dimIdNTime, dimIdNstn /), varid)
 623 iret = nf90_put_att(ncid, varid, 'units', "-")
 624 iret = nf90_put_att(ncid, varid, 'long_name', "Discharge quality 0 to 100 to be scaled by 100.")
 625 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 626
 627 ! nudge def var
 628 iret = nf90_def_var(ncid, "nudge", nf90_float, (/ dimIdNFeature /), varid)
 629 iret = nf90_put_att(ncid, varid, 'units', "m3 s-1")
 630 iret = nf90_put_att(ncid, varid, 'long_name', "Amount of stream flow alteration")
 631 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 632
 633 !! global attributes
 634 iret = nf90_put_att(ncid, nf90_global, "modelTimeAtOutput", modelTime)
 635
 636 !! end definition
 637 iret = nf90_enddef(ncid)
 638
 639 !! gageId write var
 640 iret = nf90_inq_varid(ncid, "stationId", varid)
 641 allocate(charArr15(nSpace))
 642 charArr15=lastObsStr(:)%usgsId
 643
 644 iret = nf90_put_var(ncid, varid, charArr15)
 645 deallocate(charArr15)
 646
 647 ! time write var
 648 iret = nf90_inq_varid(ncid, "time", varid)
 649 allocate(charArr19(nTime,nSpace))
 650 do tt=1,nTime
 651    do ss=1,nSpace
 652       charArr19(tt,ss) = lastObsStr(ss)%lastObsTime(tt)
 653    end do
 654 end do
 655 iret = nf90_put_var(ncid, varid, charArr19)
 656 deallocate(charArr19)
 657
 658 ! discharge write var
 659 iret = nf90_inq_varid(ncid, "discharge", varid)
 660 allocate(tmpFloat(nTime,nSpace))
 661 do tt=1,nTime
 662    do ss=1,nSpace
 663       tmpFloat(tt,ss) = lastObsStr(ss)%lastObsDischarge(tt)
 664    end do
 665 end do
 666 iret = nf90_put_var(ncid, varid, tmpFloat)
 667 deallocate(tmpFloat)
 668
 669 ! model discharge write var
 670 iret = nf90_inq_varid(ncid, "model_discharge", varid)
 671 allocate(tmpFloat(nTime,nSpace))
 672 do tt=1,nTime
 673    do ss=1,nSpace
 674       tmpFloat(tt,ss) = lastObsStr(ss)%lastObsModelDischarge(tt)
 675    end do
 676 end do
 677 iret = nf90_put_var(ncid, varid, tmpFloat)
 678 deallocate(tmpFloat)
 679
 680 ! discharge_quality write var
 681 iret = nf90_inq_varid(ncid, "discharge_quality", varid)
 682 allocate(qualityInt2(nTime,nSpace))
 683 do tt=1,nTime
 684    do ss=1,nSpace
 685    qualityInt2(tt,ss) = nint(100*lastObsStr(ss)%lastObsQuality(tt), kind=2)
 686 end do
 687 end do
 688 iret = nf90_put_var(ncid, varid, qualityInt2)
 689 deallocate(qualityInt2)
 690
 691 ! discharge_quality write var
 692 iret = nf90_inq_varid(ncid, "nudge", varid)
 693 iret = nf90_put_var(ncid, varid, g_nudge)
 694
 695 !! close
 696 iret= nf90_close(ncid)
 697
 698 #ifdef HYDRO_D
 699 print*,'Ndg: end write_nudging_last_obs'
 700 flush(6)
 701 #endif
 702
 703 end subroutine write_nudging_last_obs
 704
 705 subroutine write_nudging_last_obs_soa(lastObsStr, modelTime, g_nudge)
 706
 707 implicit none
 708
 709 !Arguments
 710 type(lastObsStructure_SoA), intent(in) :: lastObsStr
 711 character(len=19),      intent(in)               :: modelTime
 712 real,                   intent(in), dimension(:) :: g_nudge
 713
 714 ! Local
 715 character(len=37) :: outFileName
 716 integer :: nSpace, nTime, nFeature, tt, ss
 717 integer :: ncid, iret, varid, dimIdTimeStr, dimIdStnStr, dimIdNTime, dimIdNStn, dimIdNFeature
 718 character(len=15), allocatable, dimension(:) :: charArr15
 719 character(len=19), allocatable, dimension(:,:) :: charArr19
 720 real,      allocatable, dimension(:,:) :: tmpFloat
 721 integer*2, allocatable, dimension(:,:) :: qualityInt2
 722
 723 #ifdef HYDRO_D
 724 print*,'Ndg: start write_nudging_last_obs'
 725 flush(6)
 726 #endif
 727
 728 nSpace = size(lastObsStr%usgsId)
 729 nTime  = size(lastObsStr%lastObsDischarge,1)
 730 nFeature = size(g_nudge)
 731
 732 !           !1234567891123456789212345678931234567
 733 outFileName='nudgingLastObs.' // modelTime // '.nc'
 734
 735 ! create file
 736 iret = nf90_create(outFileName, nf90_clobber, ncid)
 737
 738 ! create dimensions
 739 ! station id has length of the number of gages
 740 ! feature id has length of the number of reaches/features.
 741 iret = nf90_def_dim(ncid, "timeStrLen",      19,             dimIdTimeStr)
 742 iret = nf90_def_dim(ncid, "timeInd",         nTime,          dimIdNTime)
 743 iret = nf90_def_dim(ncid, "stationIdStrLen", 15,             dimIdStnStr)
 744 iret = nf90_def_dim(ncid, "stationIdInd",    nf90_unlimited, dimIdNStn)
 745 iret = nf90_def_dim(ncid, "feature_id",      nFeature,       dimIdNFeature)
 746
 747 ! gageId def var
 748 iret = nf90_def_var(ncid, "stationId", nf90_char, (/ dimIdStnStr, dimIdNstn /), varid)
 749 iret = nf90_put_att(ncid, varid, 'long_name', "USGS station identifer of length 15")
 750 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 751
 752 ! time def var
 753 iret = nf90_def_var(ncid, "time", nf90_char, (/ dimIdTimeStr, dimIdNTime, dimIdNStn /), varid)
 754 iret = nf90_put_att(ncid, varid, 'units', "UTC")
 755 iret = nf90_put_att(ncid, varid, 'long_name', "YYYY-mm-dd_HH:MM:SS UTC")
 756 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 757
 758 ! discharge def var
 759 iret = nf90_def_var(ncid, "discharge", nf90_float, (/ dimIdNTime, dimIdNstn /), varid)
 760 iret = nf90_put_att(ncid, varid, 'units', "m^3/s")
 761 iret = nf90_put_att(ncid, varid, 'long_name', "Discharge.cubic_meters_per_second")
 762 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 763
 764 ! model discharge def var
 765 iret = nf90_def_var(ncid, "model_discharge", nf90_float, (/ dimIdNTime, dimIdNstn /), varid)
 766 iret = nf90_put_att(ncid, varid, 'units', "m^3/s")
 767 iret = nf90_put_att(ncid, varid, 'long_name', "modelDischarge.cubic_meters_per_second")
 768 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 769
 770 ! discharge quality def var
 771 iret = nf90_def_var(ncid, "discharge_quality", nf90_short, (/ dimIdNTime, dimIdNstn /), varid)
 772 iret = nf90_put_att(ncid, varid, 'units', "-")
 773 iret = nf90_put_att(ncid, varid, 'long_name', "Discharge quality 0 to 100 to be scaled by 100.")
 774 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 775
 776 ! nudge def var
 777 iret = nf90_def_var(ncid, "nudge", nf90_float, (/ dimIdNFeature /), varid)
 778 iret = nf90_put_att(ncid, varid, 'units', "m3 s-1")
 779 iret = nf90_put_att(ncid, varid, 'long_name', "Amount of stream flow alteration")
 780 iret = nf90_def_var_deflate(ncid, varid, 0, 1, 2)
 781
 782 ! global attributes
 783 iret = nf90_put_att(ncid, nf90_global, "modelTimeAtOutput", modelTime)
 784
 785 ! end definition
 786 iret = nf90_enddef(ncid)
 787
 788 ! gageId write var
 789 iret = nf90_inq_varid(ncid, "stationId", varid)
 790 !allocate(charArr15(nSpace))
 791 !charArr15=lastObsStr%usgsId(:)
 792
 793 iret = nf90_put_var(ncid, varid, lastObsStr%usgsId)
 794 !deallocate(charArr15)
 795
 796 ! time write var
 797 iret = nf90_inq_varid(ncid, "time", varid)
 798 !allocate(charArr19(nTime,nSpace))
 799 !do tt=1,nTime
 800 !   do ss=1,nSpace
 801 !      charArr19(tt,ss) = lastObsStr%lastObsTime(tt,ss)
 802 !   end do
 803 !end do
 804 iret = nf90_put_var(ncid, varid, lastObsStr%lastObsTime)
 805 !deallocate(charArr19)
 806
 807 ! discharge write var
 808 iret = nf90_inq_varid(ncid, "discharge", varid)
 809 !allocate(tmpFloat(nTime,nSpace))
 810 !do tt=1,nTime
 811 !   do ss=1,nSpace
 812 !      tmpFloat(tt,ss) = lastObsStr%lastObsDischarge(tt,ss)
 813 !   end do
 814 !end do
 815 iret = nf90_put_var(ncid, varid, lastObsStr%lastObsDischarge)
 816 !deallocate(tmpFloat)
 817
 818 ! model discharge write var
 819 iret = nf90_inq_varid(ncid, "model_discharge", varid)
 820 !allocate(tmpFloat(nTime,nSpace))
 821 !do tt=1,nTime
 822 !   do ss=1,nSpace
 823 !      tmpFloat(tt,ss) = lastObsStr%lastObsModelDischarge(tt,ss)
 824 !   end do
 825 !end do
 826 iret = nf90_put_var(ncid, varid, lastObsStr%lastObsModelDischarge)
 827 !deallocate(tmpFloat)
 828
 829 ! discharge_quality write var
 830 iret = nf90_inq_varid(ncid, "discharge_quality", varid)
 831 allocate(qualityInt2(nTime,nSpace))
 832 do tt=1,nTime
 833    do ss=1,nSpace
 834    qualityInt2(tt,ss) = nint(100*lastObsStr%lastObsQuality(tt,ss), kind=2)
 835 end do
 836 end do
 837 iret = nf90_put_var(ncid, varid, qualityInt2)
 838 deallocate(qualityInt2)
 839
 840 ! discharge_quality write var
 841 iret = nf90_inq_varid(ncid, "nudge", varid)
 842 iret = nf90_put_var(ncid, varid, g_nudge)
 843
 844 ! close
 845 iret= nf90_close(ncid)
 846
 847 #ifdef HYDRO_D
 848 print*,'Ndg: end write_nudging_last_obs_soa'
 849 flush(6)
 850 #endif
 851
 852 end subroutine write_nudging_last_obs_soa
 853
 854
 855 !===================================================================================================
 856 ! Program Name:
 857 !   subroutine find_nudging_last_obs_file
 858 ! Author(s)/Contact(s):
 859 !   James L McCreight <jamesmcc><ucar><edu>
 860 ! Abstract:
 861 !   Return a single file path/names of nudging last obs file given a single date.
 862 ! History Log:
 863 !   2/9/16 - Created,
 864 ! Usage:
 865 ! Parameters:
 866 !   date
 867 ! Input Files:  nudgingLastObs.YYYY-mm-dd_HH:MM:SS.nc
 868 ! Output Files:
 869 ! Condition codes:
 870 ! User controllable options:
 871 ! Notes:
 872
 873 function find_nudging_last_obs_file(fileName)
 874 implicit none
 875 character(len=256) :: find_nudging_last_obs_file  ! Output
 876 character(len=256), intent(in) :: fileName
 877 !Internals
 878 logical :: fileExists
 879 #ifdef HYDRO_D
 880 print*,'Ndg: start find_nudging_last_obs_file'
 881 flush(6)
 882 #endif
 883
 884 find_nudging_last_obs_file=''
 885 ! is there a file with this name?
 886 ! note files are not resolved below minutes.
 887 inquire(FILE=fileName, EXIST=fileExists)
 888 if (fileExists) find_nudging_last_obs_file = trim(fileName)
 889
 890 #ifdef HYDRO_D
 891 print*,'Ndg: last obs file: ', trim(fileName)
 892 print*,'Ndg: file found: ', fileExists
 893 print*,'Ndg: finish find_nudging_last_obs_file'
 894 flush(6)
 895 #endif
 896
 897 end function find_nudging_last_obs_file
 898
 899
 900 !===================================================================================================
 901 ! Subroutine Name:
 902 !   subroutine read_nudging_last_obs
 903 ! Author(s)/Contact(s):
 904 !   James L McCreight <jamesmcc><ucar><edu>
 905 ! Abstract:
 906 !   read in the last observations collected over time.
 907 ! History Log:
 908 !   02/03/16 -Created, JLM.
 909 ! Usage:
 910 ! Parameters:
 911 ! Input Files:
 912 ! Output Files:
 913 ! Condition codes:
 914 ! User controllable options: None.
 915 ! Notes: Needs better error handling...
 916
 917 subroutine read_nudging_last_obs(fileName, lastObsStr, g_nudge)
 918
 919 implicit none
 920 character(len=*),       intent(in)                  :: fileName
 921 type(lastObsStructure), intent(inout), dimension(:) :: lastObsStr
 922 real,                   intent(inout), dimension(:) :: g_nudge
 923
 924 integer :: nSpace, nTime, tt, ss
 925 integer:: ncid, iret, varid
 926 real,              allocatable, dimension(:,:) :: tmpFloat
 927 integer*2,         allocatable, dimension(:,:) :: qualityInt2
 928 character(len=19), allocatable, dimension(:,:) :: charArr19
 929 character(len=15), allocatable, dimension(:)   :: charArr15
 930
 931 print*,'read_nudging_last_obs'
 932
 933 nSpace = size(lastObsStr)
 934 nTime  = size(lastObsStr(1)%lastObsDischarge)
 935
 936 iret = nf90_open(fileName, nf90_nowrite, ncid)
 937
 938 ! gageId read var
 939 allocate(charArr15(nSpace))
 940 iret = nf90_inq_varid(ncid, "stationId", varid)
 941 iret = nf90_get_var(ncid, varid, charArr15)
 942 lastObsStr(:)%usgsId=charArr15
 943 deallocate(charArr15)
 944
 945 ! time read var
 946 iret = nf90_inq_varid(ncid, "time", varid)
 947 allocate(charArr19(nTime,nSpace))
 948 iret = nf90_get_var(ncid, varid, charArr19)
 949 do tt=1,nTime
 950    do ss=1,nSpace
 951       lastObsStr(ss)%lastObsTime(tt)=charArr19(tt,ss)
 952    end do
 953 end do
 954 deallocate(charArr19)
 955
 956 ! discharge read var
 957 iret = nf90_inq_varid(ncid, "discharge", varid)
 958 allocate(tmpFloat(nTime,nSpace))
 959 iret = nf90_get_var(ncid, varid, tmpFloat)
 960 do tt=1,nTime
 961    do ss=1,nSpace
 962       lastObsStr(ss)%lastObsDischarge(tt) = tmpFloat(tt,ss)
 963    end do
 964 end do
 965 deallocate(tmpFloat)
 966
 967 ! model discharge read var
 968 iret = nf90_inq_varid(ncid, "model_discharge", varid)
 969 allocate(tmpFloat(nTime,nSpace))
 970 iret = nf90_get_var(ncid, varid, tmpFloat)
 971 do tt=1,nTime
 972    do ss=1,nSpace
 973       lastObsStr(ss)%lastObsModelDischarge(tt) = tmpFloat(tt,ss)
 974    end do
 975 end do
 976 deallocate(tmpFloat)
 977
 978 !discharge _quality read var
 979 iret = nf90_inq_varid(ncid, "discharge_quality", varid)
 980 allocate(qualityInt2(nTime,nSpace))
 981 iret = nf90_get_var(ncid, varid, qualityInt2)
 982 do tt=1,nTime
 983    do ss=1,nSpace
 984       lastObsStr(ss)%lastObsQuality(tt) = qualityInt2(tt,ss)/real(100)
 985    end do
 986 end do
 987 deallocate(qualityInt2)
 988
 989 iret = nf90_inq_varid(ncid, "nudge", varid)
 990 if(iret .eq. nf90_noerr) then
 991    iret = nf90_get_var(ncid, varid, g_nudge)
 992 else
 993    g_nudge=0.0000000000
 994 end if
 995
 996 ! close
 997 iret= nf90_close(ncid)
 998
 999 #ifdef HYDRO_D
1000 print*,'Ndg: end read_nudging_last_obs'
1001 flush(6)
1002 #endif
1003
1004 end subroutine read_nudging_last_obs
1005
1006 subroutine read_nudging_last_obs_soa(fileName, lastObsStr, g_nudge)
1007
1008 implicit none
1009 character(len=*),       intent(in)                  :: fileName
1010 type(lastObsStructure_Soa), intent(inout) :: lastObsStr
1011 real,                   intent(inout), dimension(:) :: g_nudge
1012
1013 integer :: nSpace, nTime, tt, ss
1014 integer:: ncid, iret, varid
1015 real,              allocatable, dimension(:,:) :: tmpFloat
1016 integer*2,         allocatable, dimension(:,:) :: qualityInt2
1017 character(len=19), allocatable, dimension(:,:) :: charArr19
1018 character(len=15), allocatable, dimension(:)   :: charArr15
1019
1020 print*,'read_nudging_last_obs_soa'
1021
1022 nSpace = size(lastObsStr%usgsId)
1023 nTime  = size(lastObsStr%lastObsDischarge,1)
1024
1025 iret = nf90_open(fileName, nf90_nowrite, ncid)
1026
1027 !! gageId read var
1028 !allocate(charArr15(nSpace))
1029 iret = nf90_inq_varid(ncid, "stationId", varid)
1030 iret = nf90_get_var(ncid, varid, lastObsStr%usgsId(:))
1031 !lastObsStr%usgsId(:) = charArr15
1032 !deallocate(charArr15)
1033
1034 ! time read var
1035 iret = nf90_inq_varid(ncid, "time", varid)
1036 !allocate(charArr19(nTime,nSpace))
1037 iret = nf90_get_var(ncid, varid, lastObsStr%lastObsTime)
1038 !do tt=1,nTime
1039 !   do ss=1,nSpace
1040 !      lastObsStr%lastObsTime(tt,ss)=charArr19(tt,ss)
1041 !   end do
1042 !end do
1043 !deallocate(charArr19)
1044
1045 ! discharge read var
1046 iret = nf90_inq_varid(ncid, "discharge", varid)
1047 !allocate(tmpFloat(nTime,nSpace))
1048 iret = nf90_get_var(ncid, varid, lastObsStr%lastObsDischarge)
1049 !do tt=1,nTime
1050 !   do ss=1,nSpace
1051 !      lastObsStr%lastObsDischarge(tt,ss) = tmpFloat(tt,ss)
1052 !   end do
1053 !end do
1054 !deallocate(tmpFloat)
1055
1056 ! model discharge read var
1057 iret = nf90_inq_varid(ncid, "model_discharge", varid)
1058 !allocate(tmpFloat(nTime,nSpace))
1059 iret = nf90_get_var(ncid, varid, lastObsStr%lastObsModelDischarge)
1060 !do tt=1,nTime
1061 !   do ss=1,nSpace
1062 !      lastObsStr%lastObsModelDischarge(tt,ss) = tmpFloat(tt,ss)
1063 !   end do
1064 !end do
1065 !deallocate(tmpFloat)
1066
1067 !discharge _quality read var
1068 iret = nf90_inq_varid(ncid, "discharge_quality", varid)
1069 allocate(qualityInt2(nTime,nSpace))
1070 iret = nf90_get_var(ncid, varid, qualityInt2)
1071 do tt=1,nTime
1072    do ss=1,nSpace
1073       lastObsStr%lastObsQuality(tt,ss) = qualityInt2(tt,ss)*0.01!/real(100)
1074    end do
1075 end do
1076 deallocate(qualityInt2)
1077
1078 iret = nf90_inq_varid(ncid, "nudge", varid)
1079 if(iret .eq. nf90_noerr) then
1080    iret = nf90_get_var(ncid, varid, g_nudge)
1081 else
1082    g_nudge=0.0000000000
1083 end if
1084
1085 !! close
1086 iret= nf90_close(ncid)
1087
1088 #ifdef HYDRO_D
1089 print*,'Ndg: end read_nudging_last_obs'
1090 flush(6)
1091 #endif
1092
1093 end subroutine read_nudging_last_obs_soa
1094
1095
1096 !===================================================================================================
1097 ! Subroutine Name:
1098 !   subroutine read_network_reexpression
1099 ! Author(s)/Contact(s):
1100 !   James L McCreight <jamesmcc><ucar><edu>
1101 ! Abstract:
1102 !   Read the three netcdf files which allow the stream network to be traversed with
1103 !   indexing.
1104 ! History Log:
1105 !   7/23/15 -Created, JLM.
1106 ! Usage:
1107 ! Parameters:
1108 ! Input Files:
1109 ! Output Files:
1110 ! Condition codes:
1111 ! User controllable options: None.
1112 ! Notes:
1113
1114 subroutine read_network_reexpression( &
1115                            file,      & ! file with dims of the stream ntwk
1116                            upGo,      & ! where each ind came from, upstream
1117                            upStart,   & ! where each ind's upstream links start in upGo
1118                            upEnd,     & ! where each ind's upstream links end   in upGo
1119                            downGo,    & ! where each ind goes, downstream
1120                            downStart, & ! where each ind's downstream links start in downGo
1121                            downEnd    ) ! where each ind's downstream links end   in downGo
1122
1123 implicit none
1124
1125 character(len=*), intent(in) :: file
1126 integer*4, dimension(:), intent(out) :: upGo
1127 integer*4, dimension(:), intent(out) :: upStart
1128 integer*4, dimension(:), intent(out) :: upEnd
1129 integer*4, dimension(:), intent(out) :: downGo
1130 integer*4, dimension(:), intent(out) :: downStart
1131 integer*4, dimension(:), intent(out) :: downEnd
1132
1133 character(len=25) :: subroutineName
1134 integer :: iRet, ncid, errStatus
1135 #ifdef HYDRO_D
1136 print*,"Ndg: start read_network_reexpression"
1137 flush(6)
1138 #endif
1139
1140 subroutineName = 'read_network_reexpression'
1141
1142 iRet = nf90_open(trim(file), nf90_nowrite, ncid)
1143 if (iRet /= nf90_noerr) then
1144    write(*,'("read_network_reexpression: Problem opening: ''", A, "''")') trim(file)
1145    call hydro_stop(subroutineName // ": Problem opening file: " // file)
1146 endif
1147
1148 call get_1d_netcdf_int4(ncid, 'upGo',      upGo,      subroutineName, .TRUE., errStatus)
1149 call get_1d_netcdf_int4(ncid, 'upStart',   upStart,   subroutineName, .TRUE., errStatus)
1150 call get_1d_netcdf_int4(ncid, 'upEnd',     upEnd,     subroutineName, .TRUE., errStatus)
1151 call get_1d_netcdf_int4(ncid, 'downGo',    downGo,    subroutineName, .TRUE., errStatus)
1152 call get_1d_netcdf_int4(ncid, 'downStart', downStart, subroutineName, .TRUE., errStatus)
1153 call get_1d_netcdf_int4(ncid, 'downEnd',   downEnd,   subroutineName, .TRUE., errStatus)
1154
1155 iRet = nf90_close(ncId)
1156 if (iRet /= nf90_noerr) then
1157    write(*,'("read_network_reexpression: Problem closing: ''", A, "''")') trim(file)
1158    call hydro_stop(subroutineName // ": Problem closing file: " // file)
1159 end if
1160
1161 #ifdef HYDRO_D
1162 print*,"Ndg: finish read_network_reexpression"
1163 flush(6)
1164 #endif
1165
1166 end subroutine read_network_reexpression
1167
1168 !===================================================================================================
1169 ! Subroutine Name:
1170 !   subroutine output_chan_connectivity
1171 ! Author(s)/Contact(s):
1172 !   James L McCreight <jamesmcc><ucar><edu>
1173 ! Abstract:
1174 !   For gridded channel routing, write the channel connectivity to a netcdf file
1175 !   so channel analyses can be performed offline. Do it here so any changes to the
1176 !   topology calculation are maintained without external code.
1177 ! History Log:
1178 !   5/27/15 -Created, JLM.
1179 ! Usage:
1180 ! Parameters:
1181 ! Input Files:
1182 ! Output Files:
1183 ! Condition codes:
1184 ! User controllable options:
1185 ! Notes  :
1186 !   One day it might be worth writing code to read the files output by this code,
1187 !   depending on the time spent on the calculations when restarting a domain:
1188 !   is recalculation faster than reading in the file written by this routine?
1189
1190 subroutine output_chan_connectivity( &
1191      inCHLAT,     inCHLON,             &   !! Channel grid lat, lon.
1192      inCHANLEN,                        &   !! The distance between channel grid centers in m.
1193      inFROM_NODE, inTO_NODE,           &   !! Index of a given cell and the index which it flows to.
1194      inCHANXI,    inCHANYJ,            &   !! Index on fine/routing grid of grid cells.
1195      inTYPEL,     inLAKENODE           &   !! Lake type? and node? indications.
1196 )
1197
1198
1199 #ifdef MPP_LAND
1200 use module_mpp_land
1201 #endif
1202
1203 implicit none
1204 #include <netcdf.inc>
1205
1206 !! These are the names used in module_HYDRO_io.F: SUBROUTINE READ_CHROUTING1
1207 real,    dimension(:),  intent(in) :: inCHLAT, inCHLON, inCHANLEN
1208 integer, dimension(:),  intent(in) :: inFROM_NODE, inTO_NODE, inCHANXI, inCHANYJ, inTYPEL, inLAKENODE
1209
1210 integer            :: nStreamCells, streamCellDimID
1211 integer            :: iret, projInfo_flag
1212 integer            :: ncid, ncstatic, varid
1213
1214
1215 real               :: long_cm, lat_po, fe, fn
1216 real, dimension(2) :: sp
1217 ! JLM this could go to a namelist for flexibility
1218 character(len=256), parameter :: output_flnm = "CHANNEL_CONNECTIVITY.nc"
1219
1220 real,    allocatable, dimension(:) :: CHLAT, CHLON, CHANLEN
1221 integer, allocatable, dimension(:) :: FROM_NODE, TO_NODE, CHANXI, CHANYJ, TYPEL, LAKENODE
1222
1223 !! handle the parallelization in this routine instead of in the main code.
1224
1225 #ifdef MPP_LAND
1226
1227 if(my_id .eq. io_id) then
1228    allocate( chLat(rt_domain(did)%gnlinks),    chLon(rt_domain(did)%gnlinks)     )
1229    allocate( chanLen(rt_domain(did)%gnlinks),  from_node(rt_domain(did)%gnlinks) )
1230    allocate( to_node(rt_domain(did)%gnlinks),  chanXI(rt_domain(did)%gnlinks)    )
1231    allocate( chanYJ(rt_domain(did)%gnlinks),   typeL(rt_domain(did)%gnlinks)     )
1232    allocate( lakeNode(rt_domain(did)%gnlinks) )
1233 else
1234    allocate(chLat(1),  chLon(1),  chanLen(1), from_node(1), to_node(1) )
1235    allocate(chanXI(1), chanYJ(1), typeL(1),  lakeNode(1))
1236 end if
1237
1238 call write_chanel_real(inChLat,     rt_domain(did)%map_l2g, &
1239                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, chLat)
1240 call write_chanel_real(inChLon,     rt_domain(did)%map_l2g, &
1241                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, chLon)
1242 call write_chanel_real(inChanLen,   rt_domain(did)%map_l2g, &
1243                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, chanLen)
1244 call write_chanel_int(inFrom_node, rt_domain(did)%map_l2g, &
1245                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, from_node)
1246 call write_chanel_int(inTo_node,   rt_domain(did)%map_l2g, &
1247                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, to_node)
1248 call write_chanel_int(inChanXI,    rt_domain(did)%map_l2g, &
1249                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, chanXI)
1250 call write_chanel_int(inChanYJ,    rt_domain(did)%map_l2g, &
1251                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, chanYJ)
1252 call write_chanel_int(inTypeL,     rt_domain(did)%map_l2g, &
1253                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, typeL)
1254 call write_chanel_int(inLakeNode,  rt_domain(did)%map_l2g, &
1255                        rt_domain(did)%gnlinks, rt_domain(did)%nlinks, lakeNode)
1256
1257 #else
1258
1259 allocate( chLat(rt_domain(did)%nlinks),    chLon(rt_domain(did)%nlinks)      )
1260 allocate( chanLen(rt_domain(did)%nlinks),  from_node(rt_domain(did)%nlinks)  )
1261 allocate( to_node(rt_domain(did)%nlinks),  chanXI(rt_domain(did)%nlinks)     )
1262 allocate( chanYJ(rt_domain(did)%nlinks),   typeL(rt_domain(did)%nlinks)      )
1263 allocate( lakeNode(rt_domain(did)%nlinks) )
1264
1265 chLat     = inChLat
1266 chLon     = inChLon
1267 chanlen     = inChanlen
1268 from_node = inFrom_node
1269 to_node   = inTo_node
1270 chanXI    = inChanXI
1271 chanYJ    = inChanYJ
1272 typeL     = inTypeL
1273 lakeNode  = inLakeNode
1274
1275 #endif
1276
1277
1278 #ifdef MPP_LAND
1279 if(my_id .eq. io_id) then
1280 #endif
1281
1282    !! jlm: check that all the input variables have the same dimensions? Or will that happen
1283    !! by defult when trying to write to ncdf?
1284
1285    ! Open the  finemesh static files to obtain projection information
1286    ! jlm: this seems optional, might be nice if output is used for plotting/GIS.
1287    ! jlm: set did:=1 in nlst_rt
1288 #ifdef HYDRO_D
1289    write(*,'("geo_finegrid_flnm: ''", A, "''")') trim(nlst(did)%geo_finegrid_flnm)
1290 #endif
1291
1292    iret = nf_open(trim(nlst(1)%geo_finegrid_flnm), NF_NOWRITE, ncstatic)
1293
1294    if (iret /= 0) then
1295       write(*,'("Problem opening geo_finegrid file: ''", A, "''")') &
1296            trim(nlst(1)%geo_finegrid_flnm)
1297       write(*,*) "HIRES_OUTPUT will not be georeferenced..."
1298       projInfo_flag = 0
1299    else
1300       projInfo_flag = 1
1301    endif
1302
1303    if(projInfo_flag.eq.1) then !if/then hires_georef
1304       ! Get projection information from finegrid netcdf file
1305       iret = NF_INQ_VARID(ncstatic,'lambert_conformal_conic',varid)
1306       if(iret .eq. 0) &
1307            iret = NF_GET_ATT_REAL(ncstatic, varid, 'longitude_of_central_meridian', long_cm)
1308       iret = NF_GET_ATT_REAL(ncstatic, varid, 'latitude_of_projection_origin', lat_po)
1309       iret = NF_GET_ATT_REAL(ncstatic, varid, 'false_easting', fe)
1310       iret = NF_GET_ATT_REAL(ncstatic, varid, 'false_northing', fn)
1311       iret = NF_GET_ATT_REAL(ncstatic, varid, 'standard_parallel', sp)
1312    end if  !endif hires_georef
1313    iret = nf_close(ncstatic)
1314
1315
1316    ! Create the channel connectivity file
1317 #ifdef HYDRO_D
1318    print*,'Ndg: output_flnm = "'//trim(output_flnm)//'"'
1319    flush(6)
1320 #endif
1321
1322    iret = nf90_create(trim(output_flnm), OR(NF90_CLOBBER, NF90_NETCDF4), ncid)
1323
1324    if (iret /= 0) then
1325       print*,"Ndg: Problem nf90_create"
1326       call hydro_stop("output_channel_connectivity")
1327    endif
1328
1329    nStreamCells=size(CHLON,1)
1330 #ifdef HYDRO_D
1331    print*,'Ndg: nStreamCells:', nStreamCells
1332    flush(6)
1333 #endif
1334
1335    ! Dimension definitions
1336    iret = nf_def_dim(ncid, "nStreamCells", nStreamCells, streamCellDimID)
1337
1338    ! Variable definitions
1339    ! LATITUDE - float
1340    iret = nf_def_var(ncid, "LATITUDE", NF_FLOAT, 1, (/ streamCellDimID /), varid)
1341    iret = nf_put_att_text(ncid,varid, 'long_name',     22, 'Upstream cell latitude')
1342    iret = nf_put_att_text(ncid,varid, 'standard_name',  8, 'LATITUDE')
1343    iret = nf_put_att_text(ncid,varid, 'units',          5, 'deg North')
1344
1345    ! LONGITUDE - float
1346    iret = nf_def_var(ncid, "LONGITUDE", NF_FLOAT, 1, (/ streamCellDimID /), varid)
1347    iret = nf_put_att_text(ncid, varid, 'long_name',     23, 'Upstream cell longitude')
1348    iret = nf_put_att_text(ncid, varid, 'standard_name',  9, 'LONGITUDE')
1349    iret = nf_put_att_text(ncid, varid, 'units',          8, 'deg East')
1350
1351    ! CHANLEN - float
1352    ! JLM: should check if pour points have chanLen, should they?
1353    iret = nf_def_var(ncid, "CHANLEN", NF_FLOAT, 1, (/ streamCellDimID /), varid)
1354    iret = nf_put_att_text(ncid, varid, 'units',      1,'m')
1355    iret = nf_put_att_text(ncid, varid, 'long_name', 58, &
1356         'distance between stream cell center points with downstream')
1357    iret = nf_put_att_real(ncid, varid, 'missing_value', NF_REAL, 1, -9E15)
1358
1359
1360
1361    ! FROM_NODE - integer
1362    iret = nf_def_var(ncid, "FROM_NODE", NF_INT, 1, (/ streamCellDimID /), varid)
1363    iret = nf_put_att_text(ncid, varid, 'units',      5, 'index')
1364    iret = nf_put_att_text(ncid, varid, 'long_name', 19, 'Upstream cell index')
1365    iret = nf_put_att_int(ncid, varid, 'missing_value', NF_INT, 1, -9999)
1366
1367    ! TO_NODE - integer
1368    iret = nf_def_var(ncid, "TO_NODE", NF_INT, 1, (/ streamCellDimID /), varid)
1369    iret = nf_put_att_text(ncid, varid, 'units',      5, 'index')
1370    iret = nf_put_att_text(ncid, varid, 'long_name', 21, 'Downstream cell index')
1371    iret = nf_put_att_int(ncid, varid, 'missing_value', NF_INT, 1, -9999)
1372
1373    ! CHANXI - integer
1374    iret = nf_def_var(ncid, "CHANXI", NF_INT, 1, (/ streamCellDimID /), varid)
1375    iret = nf_put_att_text(ncid, varid, 'units',      5, 'index')
1376    iret = nf_put_att_text(ncid, varid, 'long_name', 34, 'Upstream cell x index on fine grid')
1377    iret = nf_put_att_int(ncid, varid, 'missing_value', NF_INT, 1, -9999)
1378
1379    ! CHANYJ - integer
1380    iret = nf_def_var(ncid, "CHANYJ", NF_INT, 1, (/ streamCellDimID /), varid)
1381    iret = nf_put_att_text(ncid, varid, 'units',      5, 'index')
1382    iret = nf_put_att_text(ncid, varid, 'long_name', 34, 'Upstream cell y index on fine grid')
1383    iret = nf_put_att_int(ncid, varid, 'missing_value', NF_INT, 1, -9999)
1384
1385    ! TYPEL - integer
1386    iret = nf_def_var(ncid, "TYPEL", NF_INT, 1, (/ streamCellDimID /), varid)
1387    iret = nf_put_att_text(ncid, varid, 'units',      5, 'code')
1388    iret = nf_put_att_text(ncid, varid, 'long_name', 80, &
1389         'Link Type 0 is channel 1 is pour point crit depth downstream 2 is reservoir lake')
1390
1391    ! LAKENODE - integer
1392    iret = nf_def_var(ncid, "LAKENODE", NF_INT, 1, (/ streamCellDimID /), varid)
1393    iret = nf_put_att_text(ncid, varid, 'units',      5, 'index')
1394    iret = nf_put_att_text(ncid, varid, 'long_name', 32, 'Index of lake in downstream cell')
1395
1396
1397    ! Projection information
1398    if(projInfo_flag .eq. 1) then
1399       iret = nf_def_var(ncid, "lambert_conformal_conic", NF_INT, 0, 0, varid)
1400       iret = nf_put_att_text(ncid, varid, 'grid_mapping_name', 23, 'lambert_conformal_conic')
1401       iret = nf_put_att_real(ncid, varid, 'longitude_of_central_meridian', NF_FLOAT, 1, long_cm)
1402       iret = nf_put_att_real(ncid, varid, 'latitude_of_projection_origin', NF_FLOAT, 1, lat_po)
1403       iret = nf_put_att_real(ncid, varid, 'false_easting',                 NF_FLOAT, 1, fe)
1404       iret = nf_put_att_real(ncid, varid, 'false_northing',                NF_FLOAT, 1, fn)
1405       iret = nf_put_att_real(ncid, varid, 'standard_parallel',             NF_FLOAT, 2, sp)
1406    end if
1407
1408    ! End NCDF definition section
1409    iret = nf_enddef(ncid)
1410
1411    ! Put data in to the file
1412
1413    ! Data for the dim? JLM: no, seems pointless index, if not necessary
1414    !iret = nf_inq_varid(ncid,"nStreamCells", varid)
1415    !iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), 1:nStreamCells - or however you)
1416
1417    ! Reals
1418    iret = nf_inq_varid(ncid, "LATITUDE", varid)
1419    iret = nf_put_vara_real(ncid, varid, (/1/), (/ nStreamCells /), CHLAT)
1420
1421    iret = nf_inq_varid(ncid, "LONGITUDE", varid)
1422    iret = nf_put_vara_real(ncid, varid, (/1/), (/ nStreamCells /), CHLON)
1423
1424    iret = nf_inq_varid(ncid, "CHANLEN", varid)
1425    iret = nf_put_vara_real(ncid, varid, (/1/), (/ nStreamCells /), CHANLEN)
1426
1427    ! Integers
1428    iret = nf_inq_varid(ncid, "FROM_NODE", varid)
1429    iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), FROM_NODE)
1430
1431    iret = nf_inq_varid(ncid, "TO_NODE", varid)
1432    iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), TO_NODE)
1433
1434    iret = nf_inq_varid(ncid, "CHANXI", varid)
1435    iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), CHANXI)
1436
1437    iret = nf_inq_varid(ncid, "CHANYJ", varid)
1438    iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), CHANYJ)
1439
1440    iret = nf_inq_varid(ncid, "TYPEL", varid)
1441    iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), TYPEL)
1442
1443    iret = nf_inq_varid(ncid, "LAKENODE", varid)
1444    iret = nf_put_vara_int(ncid, varid, (/1/), (/ nStreamCells /), LAKENODE)
1445
1446
1447    ! Close the file
1448    iret = nf_close(ncid)
1449
1450 #ifdef MPP_LAND
1451 endif
1452
1453 deallocate(chLat, chLon, chanLen, from_node, to_node, chanXI, chanYJ, typeL, lakeNode)
1454
1455 if(my_id .eq. io_id) then
1456 #endif
1457 #ifdef HYDRO_D
1458    write(6,*) "end of output_chan_connectivity"
1459    flush(6)
1460 #endif
1461 #ifdef MPP_LAND
1462 endif
1463 #endif
1464
1465 end subroutine output_chan_connectivity
1466
1467 !===================================================================================================
1468 ! Program Names:
1469 !   get_1d_netcdf_real, get_1d_netcdf_int4
1470 ! Author(s)/Contact(s):
1471 !   James L McCreight <jamesmcc><ucar><edu>
1472 ! Abstract:
1473 !   Read a variable of real or integer type from an open netcdf file, respectively.
1474 ! History Log:
1475 !   7/17/15 -Created, JLM.
1476 ! Usage:
1477 ! Parameters:
1478 !   See definitions.
1479 ! Input Files:
1480 !   This file is refered to by it's "ncid" obtained from nc_open
1481 !   prior to calling this routine.
1482 ! Output Files:
1483 !   None.
1484 ! Condition codes:
1485 !   hydro_stop is passed "get_1d_netcdf".
1486 ! User controllable options:
1487 ! Notes:
1488 !   Could define an interface for these.
1489
1490 subroutine get_1d_netcdf_int2(ncid, varName, var, callingRoutine, fatal_if_error, errStatus)
1491 implicit none
1492 integer*4,               intent(in)  :: ncid !! the file identifier
1493 character(len=*),        intent(in)  :: varName
1494 integer*2, dimension(:), intent(out) :: var
1495 character(len=*),        intent(in)  :: callingRoutine
1496 logical,                 intent(in)  :: fatal_if_error
1497 integer,                 intent(out) :: errStatus
1498 integer :: varid, iret
1499 errStatus=0
1500 iRet = nf90_inq_varid(ncid, varName, varid)
1501 if (iret /= nf90_NoErr) then
1502    print*, trim(callingRoutine) // ": get_1d_netcdf_int2: variable: " // trim(varName)
1503    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_int2")
1504    errStatus=errStatus+1
1505 end if
1506 iRet = nf90_get_var(ncid, varid, var)
1507 if (iret /= nf90_NoErr) then
1508    print*, trim(callingRoutine) // ": get_1d_netcdf_int2: values: " // trim(varName)
1509    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_int2")
1510    errStatus=errStatus+1
1511 end if
1512 end subroutine get_1d_netcdf_int2
1513
1514
1515 subroutine get_1d_netcdf_int4(ncid, varName, var, callingRoutine, fatal_if_error, errStatus)
1516 implicit none
1517 integer*4,             intent(in)  :: ncid !! the file identifier
1518 character(len=*),      intent(in)  :: varName
1519 integer, dimension(:), intent(out) :: var
1520 character(len=*),      intent(in)  :: callingRoutine
1521 logical,               intent(in)  :: fatal_if_error
1522 integer,               intent(out) :: errStatus
1523 integer :: varid, iret
1524 errStatus=0
1525 iRet = nf90_inq_varid(ncid, varName, varid)
1526 if (iret /= nf90_NoErr) then
1527    print*, trim(callingRoutine) // ": get_1d_netcdf_int4: variable: " // trim(varName)
1528    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_int4")
1529    errStatus=errStatus+1
1530 end if
1531 iRet = nf90_get_var(ncid, varid, var)
1532 if (iret /= nf90_NoErr) then
1533    print*, trim(callingRoutine) // ": get_1d_netcdf_int4: values: " // trim(varName)
1534    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_int4")
1535    errStatus=errStatus+1
1536 end if
1537 end subroutine get_1d_netcdf_int4
1538
1539
1540 subroutine get_1d_netcdf_real(ncid, varName, var, callingRoutine, fatal_if_error, errStatus)
1541 implicit none
1542 integer,            intent(in)  :: ncid !! the file identifier
1543 character(len=*),   intent(in)  :: varName
1544 real, dimension(:), intent(out) :: var
1545 character(len=*),   intent(in)  :: callingRoutine
1546 logical,            intent(in)  :: fatal_if_error
1547 integer,            intent(out) :: errStatus
1548 integer :: varId, iRet
1549 errStatus=0
1550 iRet = nf90_inq_varid(ncid, varName, varid)
1551 if (iret /= nf90_NoErr) then
1552    print*, trim(callingRoutine) // ": get_1d_netcdf_real: variable: " // trim(varName)
1553    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_real")
1554    errStatus=errStatus+1
1555 end if
1556 iRet = nf90_get_var(ncid, varid, var)
1557 if (iret /= nf90_NoErr) then
1558    print*, trim(callingRoutine) // ": get_1d_netcdf_real: values: " // trim(varName)
1559    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_real")
1560    errStatus=errStatus+1
1561 end if
1562 end subroutine get_1d_netcdf_real
1563
1564
1565 subroutine get_1d_netcdf_text(ncid, varName, var, callingRoutine, fatal_if_error, errStatus)
1566 implicit none
1567 character(len=*), dimension(:), intent(out) :: var
1568 integer,                        intent(in)  :: ncid !! the file identifier
1569 character(len=*),               intent(in)  :: varName
1570 character(len=*),               intent(in)  :: callingRoutine
1571 logical,                        intent(in)  :: fatal_if_error
1572 integer,                        intent(out) :: errStatus
1573 integer :: varId, iRet
1574 errStatus=0
1575 iRet = nf90_inq_varid(ncid, varName, varid)
1576 if (iret /= nf90_NoErr) then
1577    print*, trim(callingRoutine) // ": get_1d_netcdf_text: variable: " // trim(varName)
1578    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_text")
1579    errStatus=errStatus+1
1580 end if
1581 iRet = nf90_get_var(ncid, varid, var)
1582 if (iret /= nf90_NoErr) then
1583    print*, trim(callingRoutine) // ": get_1d_netcdf_text: values: " // trim(varName)
1584    if (fatal_IF_ERROR) call hydro_stop("get_1d_netcdf_text")
1585    errStatus=errStatus+1
1586 end if
1587 end subroutine get_1d_netcdf_text
1588
1589
1590 !==============================================================================
1591 ! 3D real
1592 subroutine get_3d_netcdf_real(ncid, varName, var, callingRoutine, fatal_if_error, errStatus)
1593 implicit none
1594 integer,                intent(in)  :: ncid !! the file identifier
1595 character(len=*),       intent(in)  :: varName
1596 real, dimension(:,:,:), intent(out) :: var
1597 character(len=*),       intent(in)  :: callingRoutine
1598 logical,                intent(in)  :: fatal_if_error
1599 integer,                intent(out) :: errStatus
1600 integer :: varId, iRet
1601 errStatus=0
1602 iRet = nf90_inq_varid(ncid, varName, varid)
1603 if (iret /= nf90_NoErr) then
1604    print*, trim(callingRoutine) // ": get_3d_netcdf_real: variable: " // trim(varName)
1605    if (fatal_IF_ERROR) call hydro_stop("get_3d_netcdf_real")
1606    errStatus=errStatus+1
1607 end if
1608 iRet = nf90_get_var(ncid, varid, var)
1609 if (iret /= nf90_NoErr) then
1610    print*, trim(callingRoutine) // ": get_3d_netcdf_real: values: " // trim(varName)
1611    if (fatal_IF_ERROR) call hydro_stop("get_3d_netcdf_real")
1612    errStatus=errStatus+1
1613 end if
1614 end subroutine get_3d_netcdf_real
1615
1616 !===================================================================================================
1617 ! Program Names:
1618 !   get_netcdf_dim
1619 ! Author(s)/Contact(s):
1620 !   James L McCreight <jamesmcc><ucar><edu>
1621 ! Abstract:
1622 !   Get the length of a provided dimension.
1623 ! History Log:
1624 !   7/23/15 -Created, JLM.
1625 ! Usage:
1626 ! Parameters:
1627 !   file: character, the file to query
1628 !   dimName: character, the name of the dimension
1629 !   callingRoutine: character, the name of the calling routine for error messages
1630 ! Input Files:
1631 !   Specified argument.
1632 ! Output Files:
1633 ! Condition codes:
1634 !   hydro_stop is called. .
1635 ! User controllable options:
1636 ! Notes:
1637
1638 function get_netcdf_dim(file, dimName, callingRoutine, fatalErr)
1639 implicit none
1640 integer :: get_netcdf_dim  !! return value, zero if failure
1641 character(len=*), intent(in)   :: file, dimName, callingRoutine
1642 logical, optional, intent(in) :: fatalErr
1643 logical :: fatalErr_local
1644 integer :: ncId, dimId, iRet
1645
1646 fatalErr_local = .false.
1647 if(present(fatalErr)) fatalErr_local=fatalErr
1648
1649 #ifdef HYDRO_D
1650 write(*,'("getting dimension from file: ''", A, "''")') trim(file)
1651 flush(6)
1652 #endif
1653
1654 iRet = nf90_open(trim(file), nf90_NOWRITE, ncId)
1655 if (iret /= nf90_noerr) then
1656    write(*,'("Problem opening file: ''", A, "''")') trim(file)
1657    if(fatalErr_local) call hydro_stop(trim(callingRoutine) // ': get_netcdf_dim')
1658    get_netcdf_dim = 0
1659    return
1660 endif
1661
1662 iRet = nf90_inq_dimid(ncId, trim(dimName), dimId)
1663 if (iret /= nf90_noerr) then
1664    write(*,'("Problem getting the dimension ID: ''", A, "''")') trim(file)
1665    if(fatalErr_local) call hydro_stop(trim(callingRoutine) // ': get_netcdf_dim')
1666    get_netcdf_dim = 0
1667    return
1668 endif
1669
1670 iRet = nf90_inquire_dimension(ncId, dimId, len= get_netcdf_dim)
1671 if (iret /= nf90_noerr) then
1672    write(*,'("Problem getting the dimension length: ''", A, "''")') trim(file)
1673    if(fatalErr_local) call hydro_stop(trim(callingRoutine) // ': get_netcdf_dim')
1674    get_netcdf_dim = 0
1675    return
1676 endif
1677
1678 iRet = nf90_close(ncId)
1679 if (iret /= nf90_noerr) then
1680    write(*,'("Problem closing file: ''", A, "''")') trim(file)
1681    if(fatalErr_local) call hydro_stop(trim(callingRoutine) // ': get_netcdf_dim')
1682    get_netcdf_dim = 0
1683    return
1684 endif
1685 end function get_netcdf_dim
1686
1687
1688 !===================================================================================================
1689 end module module_nudging_io