| 1 | #include "cppdefs.h"
|
|---|
| 2 | MODULE nf_fread4d_mod
|
|---|
| 3 | !
|
|---|
| 4 | !svn $Id$
|
|---|
| 5 | !================================================== Hernan G. Arango ===
|
|---|
| 6 | ! Copyright (c) 2002-2018 The ROMS/TOMS Group !
|
|---|
| 7 | ! Licensed under a MIT/X style license !
|
|---|
| 8 | ! See License_ROMS.txt !
|
|---|
| 9 | !=======================================================================
|
|---|
| 10 | ! !
|
|---|
| 11 | ! This function reads in a generic floating point 4D array from an !
|
|---|
| 12 | ! input NetCDF file. !
|
|---|
| 13 | ! !
|
|---|
| 14 | ! On Input: !
|
|---|
| 15 | ! !
|
|---|
| 16 | ! ng Nested grid number (integer) !
|
|---|
| 17 | ! model Calling model identifier (integer) !
|
|---|
| 18 | ! ncname NetCDF file name (string) !
|
|---|
| 19 | ! ncid NetCDF file ID (integer) !
|
|---|
| 20 | ! ncvname NetCDF variable name (string) !
|
|---|
| 21 | ! ncvarid NetCDF variable ID (integer) !
|
|---|
| 22 | ! tindex NetCDF time record index to read (integer) !
|
|---|
| 23 | ! gtype C-grid type (integer) !
|
|---|
| 24 | ! Vsize Variable dimensions in NetCDF file (integer 1D array) !
|
|---|
| 25 | ! LBi I-dimension Lower bound (integer) !
|
|---|
| 26 | ! UBi I-dimension Upper bound (integer) !
|
|---|
| 27 | ! LBj J-dimension Lower bound (integer) !
|
|---|
| 28 | ! UBj J-dimension Upper bound (integer) !
|
|---|
| 29 | ! LBk K-dimension Lower bound (integer) !
|
|---|
| 30 | ! UBk K-dimension Upper bound (integer) !
|
|---|
| 31 | ! LBt Time-dimension Lower bound (integer) !
|
|---|
| 32 | ! UBt Time-dimension Upper bound (integer) !
|
|---|
| 33 | ! Ascl Factor to scale field after reading (real). !
|
|---|
| 34 | ! Amask Land/Sea mask, if any (real 4D array) !
|
|---|
| 35 | ! !
|
|---|
| 36 | ! On Output: !
|
|---|
| 37 | ! !
|
|---|
| 38 | ! Amin Field minimum value (real) !
|
|---|
| 39 | ! Amax Field maximum value (real) !
|
|---|
| 40 | ! A Field to read in (real 4D array) !
|
|---|
| 41 | ! nf_fread4d Error flag (integer) !
|
|---|
| 42 | ! !
|
|---|
| 43 | !=======================================================================
|
|---|
| 44 | !
|
|---|
| 45 | implicit none
|
|---|
| 46 |
|
|---|
| 47 | CONTAINS
|
|---|
| 48 |
|
|---|
| 49 | #if defined PARALLEL_IO && defined DISTRIBUTE
|
|---|
| 50 | !
|
|---|
| 51 | !***********************************************************************
|
|---|
| 52 | FUNCTION nf_fread4d (ng, model, ncname, ncid, &
|
|---|
| 53 | & ncvname, ncvarid, &
|
|---|
| 54 | & tindex, gtype, Vsize, &
|
|---|
| 55 | & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, &
|
|---|
| 56 | & Ascl, Amin, Amax, &
|
|---|
| 57 | # ifdef MASKING
|
|---|
| 58 | & Amask, &
|
|---|
| 59 | # endif
|
|---|
| 60 | & A)
|
|---|
| 61 | !***********************************************************************
|
|---|
| 62 | !
|
|---|
| 63 | USE mod_param
|
|---|
| 64 | USE mod_parallel
|
|---|
| 65 | USE mod_iounits
|
|---|
| 66 | USE mod_ncparam
|
|---|
| 67 | USE mod_netcdf
|
|---|
| 68 | USE mod_scalars
|
|---|
| 69 | !
|
|---|
| 70 | USE distribute_mod, ONLY : mp_bcasti, mp_reduce
|
|---|
| 71 | # if defined MASKING && defined READ_WATER
|
|---|
| 72 | USE distribute_mod, ONLY : mp_collect
|
|---|
| 73 | # endif
|
|---|
| 74 | USE strings_mod, ONLY : FoundError
|
|---|
| 75 | !
|
|---|
| 76 | implicit none
|
|---|
| 77 | !
|
|---|
| 78 | ! Imported variable declarations.
|
|---|
| 79 | !
|
|---|
| 80 | integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
|
|---|
| 81 | integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt
|
|---|
| 82 | integer, intent(in) :: Vsize(4)
|
|---|
| 83 |
|
|---|
| 84 | real(r8), intent(in) :: Ascl
|
|---|
| 85 | real(r8), intent(out) :: Amin
|
|---|
| 86 | real(r8), intent(out) :: Amax
|
|---|
| 87 |
|
|---|
| 88 | character (len=*), intent(in) :: ncname
|
|---|
| 89 | character (len=*), intent(in) :: ncvname
|
|---|
| 90 |
|
|---|
| 91 | # ifdef ASSUMED_SHAPE
|
|---|
| 92 | # ifdef MASKING
|
|---|
| 93 | real(r8), intent(in) :: Amask(LBi:,LBj:)
|
|---|
| 94 | # endif
|
|---|
| 95 | real(r8), intent(out) :: A(LBi:,LBj:,LBk:,LBt:)
|
|---|
| 96 | # else
|
|---|
| 97 | # ifdef MASKING
|
|---|
| 98 | real(r8), intent(in) :: Amask(LBi:UBi,LBj:UBj)
|
|---|
| 99 | # endif
|
|---|
| 100 | real(r8), intent(out) :: A(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt)
|
|---|
| 101 | # endif
|
|---|
| 102 | !
|
|---|
| 103 | ! Local variable declarations.
|
|---|
| 104 | !
|
|---|
| 105 | logical, dimension(3) :: foundit
|
|---|
| 106 |
|
|---|
| 107 | integer :: i, ic, ij, j, jc, k, kc, l, lc, np, Npts
|
|---|
| 108 | integer :: Imin, Imax, Isize, Jmin, Jmax, Jsize, IJsize
|
|---|
| 109 | integer :: Istr, Iend
|
|---|
| 110 | integer :: Ioff, Joff, Koff, Loff
|
|---|
| 111 | integer :: Ilen, Jlen, Klen, Llen, IJlen
|
|---|
| 112 | integer :: Cgrid, MyType, ghost, status, wtype
|
|---|
| 113 |
|
|---|
| 114 | integer, dimension(5) :: start, total
|
|---|
| 115 |
|
|---|
| 116 | integer :: nf_fread4d
|
|---|
| 117 |
|
|---|
| 118 | real(r8) :: Afactor, Aoffset, Aspval
|
|---|
| 119 |
|
|---|
| 120 | real(r8), parameter :: IniVal= 0.0_r8
|
|---|
| 121 |
|
|---|
| 122 | real(r8), dimension(2) :: buffer
|
|---|
| 123 | real(r8), dimension(3) :: AttValue
|
|---|
| 124 |
|
|---|
| 125 | # if defined MASKING && defined READ_WATER
|
|---|
| 126 | real(r8), allocatable :: A2d(:)
|
|---|
| 127 | # endif
|
|---|
| 128 | real(r8), allocatable :: wrk(:)
|
|---|
| 129 |
|
|---|
| 130 | character (len= 3), dimension(2) :: op_handle
|
|---|
| 131 | character (len=12), dimension(3) :: AttName
|
|---|
| 132 | !
|
|---|
| 133 | !-----------------------------------------------------------------------
|
|---|
| 134 | ! Set starting and ending indices to process.
|
|---|
| 135 | !-----------------------------------------------------------------------
|
|---|
| 136 | !
|
|---|
| 137 | ! Set first and last grid point according to staggered C-grid
|
|---|
| 138 | ! classification. Set the offsets for variables with starting
|
|---|
| 139 | ! zero-index. Recall the NetCDF does not support a zero-index.
|
|---|
| 140 | !
|
|---|
| 141 | ! Notice that (Imin,Jmin) and (Imax,Jmax) are the corner of the
|
|---|
| 142 | ! computational tile. If ghost=0, ghost points are not processed.
|
|---|
| 143 | ! They will be processed elsewhere by the appropriate call to any
|
|---|
| 144 | ! of the routines in "mp_exchange.F". If ghost=1, the ghost points
|
|---|
| 145 | ! are read.
|
|---|
| 146 | !
|
|---|
| 147 | # ifdef NO_READ_GHOST
|
|---|
| 148 | ghost=0 ! non-overlapping, no ghost points
|
|---|
| 149 | # else
|
|---|
| 150 | IF (model.eq.iADM) THEN
|
|---|
| 151 | ghost=0 ! non-overlapping, no ghost points
|
|---|
| 152 | ELSE
|
|---|
| 153 | ghost=1 ! overlapping, read ghost points
|
|---|
| 154 | END IF
|
|---|
| 155 | # endif
|
|---|
| 156 |
|
|---|
| 157 | MyType=gtype
|
|---|
| 158 |
|
|---|
| 159 | SELECT CASE (ABS(MyType))
|
|---|
| 160 | CASE (p2dvar, p3dvar)
|
|---|
| 161 | Cgrid=1
|
|---|
| 162 | Isize=IOBOUNDS(ng)%xi_psi
|
|---|
| 163 | Jsize=IOBOUNDS(ng)%eta_psi
|
|---|
| 164 | CASE (r2dvar, r3dvar, w3dvar)
|
|---|
| 165 | Cgrid=2
|
|---|
| 166 | Isize=IOBOUNDS(ng)%xi_rho
|
|---|
| 167 | Jsize=IOBOUNDS(ng)%eta_rho
|
|---|
| 168 | CASE (u2dvar, u3dvar)
|
|---|
| 169 | Cgrid=3
|
|---|
| 170 | Isize=IOBOUNDS(ng)%xi_u
|
|---|
| 171 | Jsize=IOBOUNDS(ng)%eta_u
|
|---|
| 172 | CASE (v2dvar, v3dvar)
|
|---|
| 173 | Cgrid=4
|
|---|
| 174 | Isize=IOBOUNDS(ng)%xi_v
|
|---|
| 175 | Jsize=IOBOUNDS(ng)%eta_v
|
|---|
| 176 | CASE DEFAULT
|
|---|
| 177 | Cgrid=2
|
|---|
| 178 | Isize=IOBOUNDS(ng)%xi_rho
|
|---|
| 179 | Jsize=IOBOUNDS(ng)%eta_rho
|
|---|
| 180 | END SELECT
|
|---|
| 181 |
|
|---|
| 182 | Imin=BOUNDS(ng)%Imin(Cgrid,ghost,MyRank)
|
|---|
| 183 | Imax=BOUNDS(ng)%Imax(Cgrid,ghost,MyRank)
|
|---|
| 184 | Jmin=BOUNDS(ng)%Jmin(Cgrid,ghost,MyRank)
|
|---|
| 185 | Jmax=BOUNDS(ng)%Jmax(Cgrid,ghost,MyRank)
|
|---|
| 186 |
|
|---|
| 187 | Ilen=Imax-Imin+1
|
|---|
| 188 | Jlen=Jmax-Jmin+1
|
|---|
| 189 | Klen=UBk-LBk+1
|
|---|
| 190 | Llen=UBt-LBt+1
|
|---|
| 191 | !
|
|---|
| 192 | ! Check if the following attributes: "scale_factor", "add_offset", and
|
|---|
| 193 | ! "_FillValue" are present in the input NetCDF variable:
|
|---|
| 194 | !
|
|---|
| 195 | ! If the "scale_value" attribute is present, the data is multiplied by
|
|---|
| 196 | ! this factor after reading.
|
|---|
| 197 | ! If the "add_offset" attribute is present, this value is added to the
|
|---|
| 198 | ! data after reading.
|
|---|
| 199 | ! If both "scale_factor" and "add_offset" attributes are present, the
|
|---|
| 200 | ! data are first scaled before the offset is added.
|
|---|
| 201 | ! If the "_FillValue" attribute is present, the data having this value
|
|---|
| 202 | ! is treated as missing and it is replaced with zero. This feature it
|
|---|
| 203 | ! is usually related with the land/sea masking.
|
|---|
| 204 | !
|
|---|
| 205 | AttName(1)='scale_factor'
|
|---|
| 206 | AttName(2)='add_offset '
|
|---|
| 207 | AttName(3)='_FillValue '
|
|---|
| 208 |
|
|---|
| 209 | CALL netcdf_get_fatt (ng, model, ncname, ncvarid, AttName, &
|
|---|
| 210 | & AttValue, foundit, &
|
|---|
| 211 | & ncid = ncid)
|
|---|
| 212 | IF (FoundError(exit_flag, NoError, __LINE__, &
|
|---|
| 213 | & __FILE__)) THEN
|
|---|
| 214 | nf_fread4d=ioerror
|
|---|
| 215 | RETURN
|
|---|
| 216 | END IF
|
|---|
| 217 |
|
|---|
| 218 | IF (.not.foundit(1)) THEN
|
|---|
| 219 | Afactor=1.0_r8
|
|---|
| 220 | ELSE
|
|---|
| 221 | Afactor=AttValue(1)
|
|---|
| 222 | END IF
|
|---|
| 223 |
|
|---|
| 224 | IF (.not.foundit(2)) THEN
|
|---|
| 225 | Aoffset=0.0_r8
|
|---|
| 226 | ELSE
|
|---|
| 227 | Aoffset=AttValue(2)
|
|---|
| 228 | END IF
|
|---|
| 229 |
|
|---|
| 230 | IF (.not.foundit(3)) THEN
|
|---|
| 231 | Aspval=spval_check
|
|---|
| 232 | ELSE
|
|---|
| 233 | Aspval=AttValue(3)
|
|---|
| 234 | END IF
|
|---|
| 235 | !
|
|---|
| 236 | !-----------------------------------------------------------------------
|
|---|
| 237 | ! Parallel I/O: Read in tile data from requested field and scale it.
|
|---|
| 238 | ! Processing both water and land points.
|
|---|
| 239 | !-----------------------------------------------------------------------
|
|---|
| 240 | !
|
|---|
| 241 | IF (gtype.gt.0) THEN
|
|---|
| 242 | !
|
|---|
| 243 | ! Set offsets due the NetCDF dimensions. Recall that some output
|
|---|
| 244 | ! variables not always start at one.
|
|---|
| 245 | !
|
|---|
| 246 | SELECT CASE (ABS(MyType))
|
|---|
| 247 | CASE (p2dvar, p3dvar)
|
|---|
| 248 | Ioff=0
|
|---|
| 249 | Joff=0
|
|---|
| 250 | CASE (r2dvar, r3dvar, w3dvar)
|
|---|
| 251 | Ioff=1
|
|---|
| 252 | Joff=1
|
|---|
| 253 | CASE (u2dvar, u3dvar)
|
|---|
| 254 | Ioff=0
|
|---|
| 255 | Joff=1
|
|---|
| 256 | CASE (v2dvar, v3dvar)
|
|---|
| 257 | Ioff=1
|
|---|
| 258 | Joff=0
|
|---|
| 259 | CASE DEFAULT
|
|---|
| 260 | Ioff=1
|
|---|
| 261 | Joff=1
|
|---|
| 262 | END SELECT
|
|---|
| 263 |
|
|---|
| 264 | IF (LBk.eq.0) THEN
|
|---|
| 265 | Koff=1
|
|---|
| 266 | ELSE
|
|---|
| 267 | Koff=0
|
|---|
| 268 | END IF
|
|---|
| 269 |
|
|---|
| 270 | IF (LBt.eq.0) THEN
|
|---|
| 271 | Loff=1
|
|---|
| 272 | ELSE
|
|---|
| 273 | Loff=0
|
|---|
| 274 | END IF
|
|---|
| 275 |
|
|---|
| 276 | Npts=Ilen*Jlen*Klen*Llen
|
|---|
| 277 | !
|
|---|
| 278 | ! Allocate scratch work array.
|
|---|
| 279 | !
|
|---|
| 280 | IF (.not.allocated(wrk)) THEN
|
|---|
| 281 | allocate ( wrk(Npts) )
|
|---|
| 282 | wrk=0.0_r8
|
|---|
| 283 | END IF
|
|---|
| 284 | !
|
|---|
| 285 | ! Read in data: all parallel nodes read their own tile data.
|
|---|
| 286 | !
|
|---|
| 287 | start(1)=Imin+Ioff
|
|---|
| 288 | total(1)=Ilen
|
|---|
| 289 | start(2)=Jmin+Joff
|
|---|
| 290 | total(2)=Jlen
|
|---|
| 291 | start(3)=LBk+Koff
|
|---|
| 292 | total(3)=Klen
|
|---|
| 293 | start(4)=LBt+Loff
|
|---|
| 294 | total(4)=Llen
|
|---|
| 295 | start(5)=tindex
|
|---|
| 296 | total(5)=1
|
|---|
| 297 |
|
|---|
| 298 | status=nf90_get_var(ncid, ncvarid, wrk, start, total)
|
|---|
| 299 | nf_fread4d=status
|
|---|
| 300 | !
|
|---|
| 301 | ! Scale read data and process fill values, if any. Compute minimum
|
|---|
| 302 | ! and maximum values.
|
|---|
| 303 | !
|
|---|
| 304 | IF (status.eq.nf90_noerr) THEN
|
|---|
| 305 | Amin=spval
|
|---|
| 306 | Amax=-spval
|
|---|
| 307 | DO i=1,Npts
|
|---|
| 308 | IF (ABS(wrk(i)).ge.ABS(Aspval)) THEN
|
|---|
| 309 | wrk(i)=0.0_r8 ! masked with _FillValue
|
|---|
| 310 | ELSE
|
|---|
| 311 | wrk(i)=Ascl*(Afactor*wrk(i)+Aoffset)
|
|---|
| 312 | Amin=MIN(Amin,wrk(i))
|
|---|
| 313 | Amax=MAX(Amax,wrk(i))
|
|---|
| 314 | END IF
|
|---|
| 315 | END DO
|
|---|
| 316 | !
|
|---|
| 317 | ! Set minimum and maximum values: global reduction.
|
|---|
| 318 | !
|
|---|
| 319 | buffer(1)=Amin
|
|---|
| 320 | op_handle(1)='MIN'
|
|---|
| 321 | buffer(2)=Amax
|
|---|
| 322 | op_handle(2)='MAX'
|
|---|
| 323 | CALL mp_reduce (ng, model, 2, buffer, op_handle)
|
|---|
| 324 | Amin=buffer(1)
|
|---|
| 325 | Amax=buffer(2)
|
|---|
| 326 | !
|
|---|
| 327 | ! Unpack read data.
|
|---|
| 328 | !
|
|---|
| 329 | ic=0
|
|---|
| 330 | DO l=LBt,UBt
|
|---|
| 331 | DO k=LBk,UBk
|
|---|
| 332 | DO j=Jmin,Jmax
|
|---|
| 333 | DO i=Imin,Imax
|
|---|
| 334 | ic=ic+1
|
|---|
| 335 | A(i,j,k,l)=wrk(ic)
|
|---|
| 336 | END DO
|
|---|
| 337 | END DO
|
|---|
| 338 | END DO
|
|---|
| 339 | END DO
|
|---|
| 340 | ELSE
|
|---|
| 341 | exit_flag=2
|
|---|
| 342 | ioerror=status
|
|---|
| 343 | END IF
|
|---|
| 344 | END IF
|
|---|
| 345 |
|
|---|
| 346 | # if defined MASKING && defined READ_WATER
|
|---|
| 347 | !
|
|---|
| 348 | !-----------------------------------------------------------------------
|
|---|
| 349 | ! Parallel I/O: Read in tile data from requested field and scale it.
|
|---|
| 350 | ! Processing water points only.
|
|---|
| 351 | !-----------------------------------------------------------------------
|
|---|
| 352 | !
|
|---|
| 353 | IF (gtype.lt.0) THEN
|
|---|
| 354 | !
|
|---|
| 355 | ! Set number of points to process, grid type switch, and offsets due
|
|---|
| 356 | ! array packing into 1D array in column-major order.
|
|---|
| 357 | !
|
|---|
| 358 | SELECT CASE (ABS(MyType))
|
|---|
| 359 | CASE (p3dvar)
|
|---|
| 360 | IJlen=IOBOUNDS(ng)%xy_psi
|
|---|
| 361 | wtype=p2dvar
|
|---|
| 362 | Ioff=0
|
|---|
| 363 | Joff=1
|
|---|
| 364 | CASE (r3dvar, w3dvar)
|
|---|
| 365 | IJlen=IOBOUNDS(ng)%xy_rho
|
|---|
| 366 | wtype=r2dvar
|
|---|
| 367 | Ioff=1
|
|---|
| 368 | Joff=0
|
|---|
| 369 | CASE (u3dvar)
|
|---|
| 370 | IJlen=IOBOUNDS(ng)%xy_u
|
|---|
| 371 | wtype=u2dvar
|
|---|
| 372 | Ioff=0
|
|---|
| 373 | Joff=0
|
|---|
| 374 | CASE (v3dvar)
|
|---|
| 375 | IJlen=IOBOUNDS(ng)%xy_v
|
|---|
| 376 | wtype=v2dvar
|
|---|
| 377 | Ioff=1
|
|---|
| 378 | Joff=1
|
|---|
| 379 | CASE DEFAULT
|
|---|
| 380 | IJlen=IOBOUNDS(ng)%xy_rho
|
|---|
| 381 | wtype=r2dvar
|
|---|
| 382 | Ioff=1
|
|---|
| 383 | Joff=0
|
|---|
| 384 | END SELECT
|
|---|
| 385 |
|
|---|
| 386 | IF (LBk.eq.0) THEN
|
|---|
| 387 | Koff=0
|
|---|
| 388 | ELSE
|
|---|
| 389 | Koff=1
|
|---|
| 390 | END IF
|
|---|
| 391 |
|
|---|
| 392 | IF (LBt.eq.0) THEN
|
|---|
| 393 | Loff=1
|
|---|
| 394 | ELSE
|
|---|
| 395 | Loff=0
|
|---|
| 396 | END IF
|
|---|
| 397 |
|
|---|
| 398 | Npts=IJlen*Klen*Llen
|
|---|
| 399 | IJsize=Isize*Jsize
|
|---|
| 400 | !
|
|---|
| 401 | ! Allocate scratch work arrays.
|
|---|
| 402 | !
|
|---|
| 403 | IF (.not.allocated(A2d)) THEN
|
|---|
| 404 | allocate ( A2d(IJsize) )
|
|---|
| 405 | END IF
|
|---|
| 406 | IF (.not.allocated(wrk)) THEN
|
|---|
| 407 | allocate ( wrk(Npts) )
|
|---|
| 408 | wrk=IniVal
|
|---|
| 409 | END IF
|
|---|
| 410 | !
|
|---|
| 411 | ! Read in data: all parallel nodes read a segment of the 1D data.
|
|---|
| 412 | ! Recall that water points are pack in the NetCDF file in a single
|
|---|
| 413 | ! dimension.
|
|---|
| 414 | !
|
|---|
| 415 | CALL tile_bounds_1d (ng, MyRank, Npts, Istr, Iend)
|
|---|
| 416 |
|
|---|
| 417 | start(1)=Istr
|
|---|
| 418 | total(1)=Iend-Istr+1
|
|---|
| 419 | start(2)=1
|
|---|
| 420 | total(2)=tindex
|
|---|
| 421 |
|
|---|
| 422 | status=nf90_get_var(ncid, ncvarid, wrk(Istr:), start, total)
|
|---|
| 423 | nf_fread4d=status
|
|---|
| 424 | !
|
|---|
| 425 | ! Global reduction of work array. We need this because the packing
|
|---|
| 426 | ! of the water point only affects the model tile partition.
|
|---|
| 427 | !
|
|---|
| 428 | IF (status.eq.nf90_noerr) THEN
|
|---|
| 429 | CALL mp_collect (ng, model, Npts, IniVal, wrk)
|
|---|
| 430 | !
|
|---|
| 431 | ! Scale read data and process fill values, if any. Compute minimum
|
|---|
| 432 | ! and maximum values.
|
|---|
| 433 | !
|
|---|
| 434 | Amin=spval
|
|---|
| 435 | Amax=-spval
|
|---|
| 436 | DO i=1,Npts
|
|---|
| 437 | IF (ABS(wrk(i)).ge.ABS(Aspval)) THEN
|
|---|
| 438 | wrk(i)=0.0_r8 ! set _FillValue to zero
|
|---|
| 439 | ELSE
|
|---|
| 440 | wrk(i)=Ascl*(Afactor*wrk(i)+Aoffset)
|
|---|
| 441 | Amin=MIN(Amin,wrk(i))
|
|---|
| 442 | Amax=MAX(Amax,wrk(i))
|
|---|
| 443 | END IF
|
|---|
| 444 | END DO
|
|---|
| 445 | !
|
|---|
| 446 | ! Unpack read data. This is tricky in parallel I/O. The cheapeast
|
|---|
| 447 | ! thing to do is reconstruct a packed 2D global array and then select
|
|---|
| 448 | ! the appropriate values for the tile.
|
|---|
| 449 | !
|
|---|
| 450 | DO l=LBt,UBt
|
|---|
| 451 | lc=(l-Loff)*IJlen*Klen
|
|---|
| 452 | DO k=LBk,UBk
|
|---|
| 453 | kc=(k-Koff)*IJlen+lc
|
|---|
| 454 | A2d=IniVal
|
|---|
| 455 | DO np=1,IJlen
|
|---|
| 456 | ij=SCALARS(ng)%IJwater(np,wtype)
|
|---|
| 457 | A2d(ij)=wrk(np+kc)
|
|---|
| 458 | END DO
|
|---|
| 459 | DO j=Jmin,Jmax
|
|---|
| 460 | jc=(j-Joff)*Isize
|
|---|
| 461 | DO i=Imin,Imax
|
|---|
| 462 | ij=i+Ioff+jc
|
|---|
| 463 | A(i,j,k,l)=A2d(ij)
|
|---|
| 464 | END DO
|
|---|
| 465 | END DO
|
|---|
| 466 | END DO
|
|---|
| 467 | END DO
|
|---|
| 468 | ELSE
|
|---|
| 469 | exit_flag=2
|
|---|
| 470 | ioerror=status
|
|---|
| 471 | END IF
|
|---|
| 472 | END IF
|
|---|
| 473 | # endif
|
|---|
| 474 | !
|
|---|
| 475 | !-----------------------------------------------------------------------
|
|---|
| 476 | ! Deallocate scratch work vector.
|
|---|
| 477 | !-----------------------------------------------------------------------
|
|---|
| 478 | !
|
|---|
| 479 | # if defined MASKING && defined READ_WATER
|
|---|
| 480 | IF (allocated(A2d)) THEN
|
|---|
| 481 | deallocate (A2d)
|
|---|
| 482 | END IF
|
|---|
| 483 | # endif
|
|---|
| 484 |
|
|---|
| 485 | IF (allocated(wrk)) THEN
|
|---|
| 486 | deallocate (wrk)
|
|---|
| 487 | END IF
|
|---|
| 488 |
|
|---|
| 489 | RETURN
|
|---|
| 490 | END FUNCTION nf_fread4d
|
|---|
| 491 |
|
|---|
| 492 | #else
|
|---|
| 493 |
|
|---|
| 494 | !
|
|---|
| 495 | !***********************************************************************
|
|---|
| 496 | FUNCTION nf_fread4d (ng, model, ncname, ncid, &
|
|---|
| 497 | & ncvname, ncvarid, &
|
|---|
| 498 | & tindex, gtype, Vsize, &
|
|---|
| 499 | & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, &
|
|---|
| 500 | & Ascl, Amin, Amax, &
|
|---|
| 501 | # ifdef MASKING
|
|---|
| 502 | & Amask, &
|
|---|
| 503 | # endif
|
|---|
| 504 | & A)
|
|---|
| 505 | !***********************************************************************
|
|---|
| 506 | !
|
|---|
| 507 | USE mod_param
|
|---|
| 508 | USE mod_parallel
|
|---|
| 509 | USE mod_iounits
|
|---|
| 510 | USE mod_ncparam
|
|---|
| 511 | USE mod_netcdf
|
|---|
| 512 | USE mod_scalars
|
|---|
| 513 | !
|
|---|
| 514 | # ifdef DISTRIBUTE
|
|---|
| 515 | USE distribute_mod, ONLY : mp_bcasti
|
|---|
| 516 | # ifdef INLINE_2DIO
|
|---|
| 517 | USE distribute_mod, ONLY : mp_scatter2d
|
|---|
| 518 | # else
|
|---|
| 519 | USE distribute_mod, ONLY : mp_scatter3d
|
|---|
| 520 | # endif
|
|---|
| 521 |
|
|---|
| 522 | # endif
|
|---|
| 523 | USE strings_mod, ONLY : FoundError
|
|---|
| 524 | !
|
|---|
| 525 | implicit none
|
|---|
| 526 | !
|
|---|
| 527 | ! Imported variable declarations.
|
|---|
| 528 | !
|
|---|
| 529 | integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
|
|---|
| 530 | integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt
|
|---|
| 531 | integer, intent(in) :: Vsize(4)
|
|---|
| 532 |
|
|---|
| 533 | real(r8), intent(in) :: Ascl
|
|---|
| 534 | real(r8), intent(out) :: Amin
|
|---|
| 535 | real(r8), intent(out) :: Amax
|
|---|
| 536 |
|
|---|
| 537 | character (len=*), intent(in) :: ncname
|
|---|
| 538 | character (len=*), intent(in) :: ncvname
|
|---|
| 539 |
|
|---|
| 540 | # ifdef ASSUMED_SHAPE
|
|---|
| 541 | # ifdef MASKING
|
|---|
| 542 | real(r8), intent(in) :: Amask(LBi:,LBj:)
|
|---|
| 543 | # endif
|
|---|
| 544 | real(r8), intent(out) :: A(LBi:,LBj:,LBk:,LBt:)
|
|---|
| 545 | # else
|
|---|
| 546 | # ifdef MASKING
|
|---|
| 547 | real(r8), intent(in) :: Amask(LBi:UBi,LBj:UBj)
|
|---|
| 548 | # endif
|
|---|
| 549 | real(r8), intent(out) :: A(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt)
|
|---|
| 550 | # endif
|
|---|
| 551 | !
|
|---|
| 552 | ! Local variable declarations.
|
|---|
| 553 | !
|
|---|
| 554 | logical, dimension(3) :: foundit
|
|---|
| 555 |
|
|---|
| 556 | integer :: i, j, k, ic, fourth, Npts, NWpts, status, wtype
|
|---|
| 557 | integer :: Imin, Imax, Jmin, Jmax, Loff, Koff
|
|---|
| 558 | integer :: Ilen, Jlen, Klen, IJlen, MyType
|
|---|
| 559 | # ifdef DISTRIBUTE
|
|---|
| 560 | integer :: Nghost
|
|---|
| 561 | # endif
|
|---|
| 562 | integer, dimension(5) :: start, total
|
|---|
| 563 |
|
|---|
| 564 | integer :: nf_fread4d
|
|---|
| 565 |
|
|---|
| 566 | real(r8) :: Afactor, Aoffset, Aspval
|
|---|
| 567 |
|
|---|
| 568 | real(r8), dimension(3) :: AttValue
|
|---|
| 569 |
|
|---|
| 570 | # if defined INLINE_2DIO && defined DISTRIBUTE
|
|---|
| 571 | real(r8), dimension(2+(Lm(ng)+2)*(Mm(ng)+2)) :: wrk
|
|---|
| 572 | # else
|
|---|
| 573 | real(r8), dimension(2+(Lm(ng)+2)*(Mm(ng)+2)*(UBk-LBk+1)) :: wrk
|
|---|
| 574 | # endif
|
|---|
| 575 |
|
|---|
| 576 | character (len=12), dimension(3) :: AttName
|
|---|
| 577 | !
|
|---|
| 578 | !-----------------------------------------------------------------------
|
|---|
| 579 | ! Set starting and ending indices to process.
|
|---|
| 580 | !-----------------------------------------------------------------------
|
|---|
| 581 | !
|
|---|
| 582 | ! Set first and last grid point according to staggered C-grid
|
|---|
| 583 | ! classification. Set loops offsets.
|
|---|
| 584 | !
|
|---|
| 585 | MyType=gtype
|
|---|
| 586 |
|
|---|
| 587 | SELECT CASE (ABS(MyType))
|
|---|
| 588 | CASE (p2dvar, p3dvar)
|
|---|
| 589 | Imin=IOBOUNDS(ng)%ILB_psi
|
|---|
| 590 | Imax=IOBOUNDS(ng)%IUB_psi
|
|---|
| 591 | Jmin=IOBOUNDS(ng)%JLB_psi
|
|---|
| 592 | Jmax=IOBOUNDS(ng)%JUB_psi
|
|---|
| 593 | CASE (r2dvar, r3dvar, w3dvar)
|
|---|
| 594 | Imin=IOBOUNDS(ng)%ILB_rho
|
|---|
| 595 | Imax=IOBOUNDS(ng)%IUB_rho
|
|---|
| 596 | Jmin=IOBOUNDS(ng)%JLB_rho
|
|---|
| 597 | Jmax=IOBOUNDS(ng)%JUB_rho
|
|---|
| 598 | CASE (u2dvar, u3dvar)
|
|---|
| 599 | Imin=IOBOUNDS(ng)%ILB_u
|
|---|
| 600 | Imax=IOBOUNDS(ng)%IUB_u
|
|---|
| 601 | Jmin=IOBOUNDS(ng)%JLB_u
|
|---|
| 602 | Jmax=IOBOUNDS(ng)%JUB_u
|
|---|
| 603 | CASE (v2dvar, v3dvar)
|
|---|
| 604 | Imin=IOBOUNDS(ng)%ILB_v
|
|---|
| 605 | Imax=IOBOUNDS(ng)%IUB_v
|
|---|
| 606 | Jmin=IOBOUNDS(ng)%JLB_v
|
|---|
| 607 | Jmax=IOBOUNDS(ng)%JUB_v
|
|---|
| 608 | CASE DEFAULT
|
|---|
| 609 | Imin=IOBOUNDS(ng)%ILB_rho
|
|---|
| 610 | Imax=IOBOUNDS(ng)%IUB_rho
|
|---|
| 611 | Jmin=IOBOUNDS(ng)%JLB_rho
|
|---|
| 612 | Jmax=IOBOUNDS(ng)%JUB_rho
|
|---|
| 613 | END SELECT
|
|---|
| 614 |
|
|---|
| 615 | Ilen=Imax-Imin+1
|
|---|
| 616 | Jlen=Jmax-Jmin+1
|
|---|
| 617 | Klen=UBk-LBk+1
|
|---|
| 618 | IJlen=Ilen*Jlen
|
|---|
| 619 |
|
|---|
| 620 | IF (LBt.eq.0) THEN
|
|---|
| 621 | Loff=1
|
|---|
| 622 | ELSE
|
|---|
| 623 | Loff=0
|
|---|
| 624 | END IF
|
|---|
| 625 | #ifdef INLINE_2DIO
|
|---|
| 626 | IF (LBk.eq.0) THEN
|
|---|
| 627 | Koff=0
|
|---|
| 628 | ELSE
|
|---|
| 629 | Koff=1
|
|---|
| 630 | END IF
|
|---|
| 631 | #endif
|
|---|
| 632 |
|
|---|
| 633 | !
|
|---|
| 634 | ! Check if the following attributes: "scale_factor", "add_offset", and
|
|---|
| 635 | ! "_FillValue" are present in the input NetCDF variable:
|
|---|
| 636 | !
|
|---|
| 637 | ! If the "scale_value" attribute is present, the data is multiplied by
|
|---|
| 638 | ! this factor after reading.
|
|---|
| 639 | ! If the "add_offset" attribute is present, this value is added to the
|
|---|
| 640 | ! data after reading.
|
|---|
| 641 | ! If both "scale_factor" and "add_offset" attributes are present, the
|
|---|
| 642 | ! data are first scaled before the offset is added.
|
|---|
| 643 | ! If the "_FillValue" attribute is present, the data having this value
|
|---|
| 644 | ! is treated as missing and it is replaced with zero. This feature it
|
|---|
| 645 | ! is usually related with the land/sea masking.
|
|---|
| 646 | !
|
|---|
| 647 | AttName(1)='scale_factor'
|
|---|
| 648 | AttName(2)='add_offset '
|
|---|
| 649 | AttName(3)='_FillValue '
|
|---|
| 650 |
|
|---|
| 651 | CALL netcdf_get_fatt (ng, model, ncname, ncvarid, AttName, &
|
|---|
| 652 | & AttValue, foundit, &
|
|---|
| 653 | & ncid = ncid)
|
|---|
| 654 | IF (FoundError(exit_flag, NoError, __LINE__, &
|
|---|
| 655 | & __FILE__)) THEN
|
|---|
| 656 | nf_fread4d=ioerror
|
|---|
| 657 | RETURN
|
|---|
| 658 | END IF
|
|---|
| 659 |
|
|---|
| 660 | IF (.not.foundit(1)) THEN
|
|---|
| 661 | Afactor=1.0_r8
|
|---|
| 662 | ELSE
|
|---|
| 663 | Afactor=AttValue(1)
|
|---|
| 664 | END IF
|
|---|
| 665 |
|
|---|
| 666 | IF (.not.foundit(2)) THEN
|
|---|
| 667 | Aoffset=0.0_r8
|
|---|
| 668 | ELSE
|
|---|
| 669 | Aoffset=AttValue(2)
|
|---|
| 670 | END IF
|
|---|
| 671 |
|
|---|
| 672 | IF (.not.foundit(3)) THEN
|
|---|
| 673 | Aspval=spval_check
|
|---|
| 674 | ELSE
|
|---|
| 675 | Aspval=AttValue(3)
|
|---|
| 676 | END IF
|
|---|
| 677 |
|
|---|
| 678 | # ifdef DISTRIBUTE
|
|---|
| 679 | !
|
|---|
| 680 | ! Set the number of tile ghost points, Nghost, to scatter in
|
|---|
| 681 | ! distributed-memory applications. If Nghost=0, the ghost points
|
|---|
| 682 | ! are not processed. They will be processed elsewhere by the
|
|---|
| 683 | ! appropriate call to any of the routines in "mp_exchange.F".
|
|---|
| 684 | !
|
|---|
| 685 | # ifdef NO_READ_GHOST
|
|---|
| 686 | Nghost=0
|
|---|
| 687 | # else
|
|---|
| 688 | IF (model.eq.iADM) THEN
|
|---|
| 689 | Nghost=0
|
|---|
| 690 | ELSE
|
|---|
| 691 | Nghost=NghostPoints
|
|---|
| 692 | END IF
|
|---|
| 693 | # endif
|
|---|
| 694 | # endif
|
|---|
| 695 | # if defined READ_WATER && defined MASKING
|
|---|
| 696 | !
|
|---|
| 697 | ! If processing water points only, set number of points and type
|
|---|
| 698 | ! switch.
|
|---|
| 699 | !
|
|---|
| 700 | SELECT CASE (ABS(MyType))
|
|---|
| 701 | CASE (p3dvar)
|
|---|
| 702 | Npts=IOBOUNDS(ng)%xy_psi
|
|---|
| 703 | wtype=p2dvar
|
|---|
| 704 | CASE (r3dvar, w3dvar)
|
|---|
| 705 | Npts=IOBOUNDS(ng)%xy_rho
|
|---|
| 706 | wtype=r2dvar
|
|---|
| 707 | CASE (u3dvar)
|
|---|
| 708 | Npts=IOBOUNDS(ng)%xy_u
|
|---|
| 709 | wtype=u2dvar
|
|---|
| 710 | CASE (v3dvar)
|
|---|
| 711 | Npts=IOBOUNDS(ng)%xy_v
|
|---|
| 712 | wtype=v2dvar
|
|---|
| 713 | CASE DEFAULT
|
|---|
| 714 | Npts=IOBOUNDS(ng)%xy_rho
|
|---|
| 715 | wtype=r2dvar
|
|---|
| 716 | END SELECT
|
|---|
| 717 | NWpts=(Lm(ng)+2)*(Mm(ng)+2)
|
|---|
| 718 | # if !(defined INLINE_2DIO && defined DISTRIBUTE)
|
|---|
| 719 | Npts=Npts*Klen
|
|---|
| 720 | # endif
|
|---|
| 721 | # endif
|
|---|
| 722 | !
|
|---|
| 723 | ! Initialize local array to avoid denormalized numbers. This
|
|---|
| 724 | ! facilitates processing and debugging.
|
|---|
| 725 | !
|
|---|
| 726 | wrk=0.0_r8
|
|---|
| 727 | !
|
|---|
| 728 | !-----------------------------------------------------------------------
|
|---|
| 729 | ! Serial I/O: Read in requested field and scale it.
|
|---|
| 730 | !-----------------------------------------------------------------------
|
|---|
| 731 | !
|
|---|
| 732 | ! Proccess data as 3D slides.
|
|---|
| 733 | !
|
|---|
| 734 | Amin=spval
|
|---|
| 735 | Amax=-spval
|
|---|
| 736 |
|
|---|
| 737 | DO fourth=LBt,UBt
|
|---|
| 738 | IF (MyType.gt.0) THEN
|
|---|
| 739 | start(1)=1
|
|---|
| 740 | total(1)=Ilen
|
|---|
| 741 | start(2)=1
|
|---|
| 742 | total(2)=Jlen
|
|---|
| 743 | start(3)=1
|
|---|
| 744 | total(3)=Klen
|
|---|
| 745 | start(4)=fourth+Loff
|
|---|
| 746 | total(4)=1
|
|---|
| 747 | start(5)=tindex
|
|---|
| 748 | total(5)=1
|
|---|
| 749 | Npts=IJlen
|
|---|
| 750 | # if !(defined INLINE_2DIO && defined DISTRIBUTE)
|
|---|
| 751 | Npts=Npts*Klen
|
|---|
| 752 | # endif
|
|---|
| 753 | # if defined READ_WATER && defined MASKING
|
|---|
| 754 | ELSE
|
|---|
| 755 | start(1)=1+(fourth+Loff-1)*Npts
|
|---|
| 756 | total(1)=Npts
|
|---|
| 757 | start(2)=1
|
|---|
| 758 | total(2)=tindex
|
|---|
| 759 | # endif
|
|---|
| 760 | END IF
|
|---|
| 761 | # if defined INLINE_2DIO && defined DISTRIBUTE
|
|---|
| 762 | !
|
|---|
| 763 | ! If appropriate, process 3D data level by level to reduce memory
|
|---|
| 764 | ! requirements.
|
|---|
| 765 | !
|
|---|
| 766 | DO k=LBk,UBk
|
|---|
| 767 | start(3)=k-Koff+1
|
|---|
| 768 | total(3)=1
|
|---|
| 769 | # endif
|
|---|
| 770 | status=nf90_noerr
|
|---|
| 771 | IF (InpThread) THEN
|
|---|
| 772 | status=nf90_get_var(ncid, ncvarid, wrk, start, total)
|
|---|
| 773 | IF (status.eq.nf90_noerr) THEN
|
|---|
| 774 | DO i=1,Npts
|
|---|
| 775 | IF (ABS(wrk(i)).ge.ABS(Aspval)) THEN
|
|---|
| 776 | wrk(i)=0.0_r8 ! masked with _FillValue
|
|---|
| 777 | ELSE
|
|---|
| 778 | wrk(i)=Ascl*(Afactor*wrk(i)+Aoffset)
|
|---|
| 779 | Amin=MIN(Amin,wrk(i))
|
|---|
| 780 | Amax=MAX(Amax,wrk(i))
|
|---|
| 781 | END IF
|
|---|
| 782 | END DO
|
|---|
| 783 | END IF
|
|---|
| 784 | END IF
|
|---|
| 785 | # ifdef DISTRIBUTE
|
|---|
| 786 | CALL mp_bcasti (ng, model, status)
|
|---|
| 787 | # endif
|
|---|
| 788 | IF (FoundError(status, nf90_noerr, __LINE__, &
|
|---|
| 789 | & __FILE__)) THEN
|
|---|
| 790 | exit_flag=2
|
|---|
| 791 | ioerror=status
|
|---|
| 792 | nf_fread4d=status
|
|---|
| 793 | RETURN
|
|---|
| 794 | END IF
|
|---|
| 795 | !
|
|---|
| 796 | !-----------------------------------------------------------------------
|
|---|
| 797 | ! Serial I/O: Unpack read field.
|
|---|
| 798 | !-----------------------------------------------------------------------
|
|---|
| 799 | !
|
|---|
| 800 | # ifdef DISTRIBUTE
|
|---|
| 801 | # ifdef INLINE_2DIO
|
|---|
| 802 | CALL mp_scatter2d (ng, model, LBi, UBi, LBj, UBj, &
|
|---|
| 803 | & Nghost, MyType, Amin, Amax, &
|
|---|
| 804 | # if defined READ_WATER && defined MASKING
|
|---|
| 805 | & NWpts, SCALARS(ng)%IJwater(:,wtype), &
|
|---|
| 806 | # endif
|
|---|
| 807 | & Npts, wrk, A(:,:,k,fourth))
|
|---|
| 808 | END DO
|
|---|
| 809 | # else
|
|---|
| 810 | CALL mp_scatter3d (ng, model, LBi, UBi, LBj, UBj, LBk, UBk, &
|
|---|
| 811 | & Nghost, MyType, Amin, Amax, &
|
|---|
| 812 | # if defined READ_WATER && defined MASKING
|
|---|
| 813 | & NWpts, SCALARS(ng)%IJwater(:,wtype), &
|
|---|
| 814 | # endif
|
|---|
| 815 | & Npts, wrk, A(:,:,:,fourth))
|
|---|
| 816 | # endif
|
|---|
| 817 | # else
|
|---|
| 818 | IF (MyType.gt.0) THEN
|
|---|
| 819 | ic=0
|
|---|
| 820 | DO k=LBk,UBk
|
|---|
| 821 | DO j=Jmin,Jmax
|
|---|
| 822 | DO i=Imin,Imax
|
|---|
| 823 | ic=ic+1
|
|---|
| 824 | A(i,j,k,fourth)=wrk(ic)
|
|---|
| 825 | END DO
|
|---|
| 826 | END DO
|
|---|
| 827 | END DO
|
|---|
| 828 | # if defined MASKING || defined READ_WATER
|
|---|
| 829 | ELSE
|
|---|
| 830 | ic=0
|
|---|
| 831 | DO k=LBk,UBk
|
|---|
| 832 | DO j=Jmin,Jmax
|
|---|
| 833 | DO i=Imin,Imax
|
|---|
| 834 | IF (Amask(i,j).gt.0.0_r8) THEN
|
|---|
| 835 | ic=ic+1
|
|---|
| 836 | A(i,j,k,fourth)=wrk(ic)
|
|---|
| 837 | ELSE
|
|---|
| 838 | A(i,j,k,fourth)=0.0_r8
|
|---|
| 839 | END IF
|
|---|
| 840 | END DO
|
|---|
| 841 | END DO
|
|---|
| 842 | END DO
|
|---|
| 843 | # endif
|
|---|
| 844 | END IF
|
|---|
| 845 | # endif
|
|---|
| 846 | END DO
|
|---|
| 847 |
|
|---|
| 848 | nf_fread4d=status
|
|---|
| 849 |
|
|---|
| 850 | RETURN
|
|---|
| 851 | END FUNCTION nf_fread4d
|
|---|
| 852 | #endif
|
|---|
| 853 | END MODULE nf_fread4d_mod
|
|---|
