hydro/nudging/module_nudging_io.F90

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