- Timestamp:
- 2016-06-24T09:50:27+02:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/NERC/dev_r3874_FASTNEt/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90
r3851 r6736 7 7 !! ! 05-2008 (S. Alderson) Modified for Interpolation in memory 8 8 !! ! from input grid to model grid 9 !! ! 04-2013 (J. Harle) Addition to interpolate bdy data onto 10 !! ! model grid 9 11 !!---------------------------------------------------------------------- 10 12 … … 27 29 28 30 PUBLIC fld_map ! routine called by tides_init 29 PUBLIC fld_read, fld_fill ! called by sbc... modules30 31 31 32 TYPE, PUBLIC :: FLD_N !: Namelist field informations … … 58 59 ! ! into the WGTLIST structure 59 60 CHARACTER(len = 34) :: vcomp ! symbolic name for a vector component that needs rotation 60 LOGICAL, DIMENSION(2) :: rotn ! flag to indicate whether before/after field has been rotated 61 INTEGER :: nreclast ! last record to be read in the current file 61 LOGICAL :: rotn ! flag to indicate whether field has been rotated 62 62 END TYPE FLD 63 63 … … 98 98 !$AGRIF_END_DO_NOT_TREAT 99 99 100 PUBLIC fld_read, fld_fill ! called by sbc... modules 101 100 102 !!---------------------------------------------------------------------- 101 103 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 105 107 CONTAINS 106 108 107 SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset)109 SUBROUTINE fld_read( kt, kn_fsbc, sd, map, jit, time_offset, jpk_1) 108 110 !!--------------------------------------------------------------------- 109 111 !! *** ROUTINE fld_read *** … … 120 122 INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step) 121 123 TYPE(FLD), INTENT(inout), DIMENSION(:) :: sd ! input field related variables 122 TYPE(MAP_POINTER),INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping ind ices123 INTEGER , INTENT(in ), OPTIONAL :: kit ! subcycle timestep for timesplitting option124 INTEGER , INTENT(in ), OPTIONAL :: kt_offset ! provide fields at time other than "now"125 ! kt_offset = -1 => fields at "before" time level126 ! kt_offset = +1 => fields at "after" time level127 !etc.128 !!129 INTEGER :: itmp ! temporary variable124 TYPE(MAP_POINTER),INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping index 125 INTEGER , INTENT(in ), OPTIONAL :: jit ! subcycle timestep for timesplitting option 126 INTEGER , INTENT(in ), OPTIONAL :: time_offset ! provide fields at time other than "now" 127 ! time_offset = -1 => fields at "before" time level 128 ! time_offset = +1 => fields at "after" time levels 129 ! etc. 130 INTEGER , INTENT(in ), OPTIONAL :: jpk_1 ! 131 !! 130 132 INTEGER :: imf ! size of the structure sd 131 133 INTEGER :: jf ! dummy indices 134 INTEGER :: ireclast ! last record to be read in the current year file 132 135 INTEGER :: isecend ! number of second since Jan. 1st 00h of nit000 year at nitend 133 136 INTEGER :: isecsbc ! number of seconds between Jan. 1st 00h of nit000 year and the middle of sbc time step 134 INTEGER :: it _offset! local time offset variable137 INTEGER :: itime_add ! local time offset variable 135 138 LOGICAL :: llnxtyr ! open next year file? 136 139 LOGICAL :: llnxtmth ! open next month file? … … 140 143 REAL(wp) :: ztintb ! ratio applied to before records when doing time interpolation 141 144 CHARACTER(LEN=1000) :: clfmt ! write format 142 TYPE(MAP_POINTER) :: imap ! global-to-local mapping indices 143 !!--------------------------------------------------------------------- 144 ll_firstcall = kt == nit000 145 IF( PRESENT(kit) ) ll_firstcall = ll_firstcall .and. kit == 1 146 147 it_offset = 0 148 IF( PRESENT(kt_offset) ) it_offset = kt_offset 149 150 imap%ptr => NULL() 151 145 !!--------------------------------------------------------------------- 146 ll_firstcall = .false. 147 IF( PRESENT(jit) ) THEN 148 IF(kt == nit000 .and. jit == 1) ll_firstcall = .true. 149 ELSE 150 IF(kt == nit000) ll_firstcall = .true. 151 ENDIF 152 153 itime_add = 0 154 IF( PRESENT(time_offset) ) itime_add = time_offset 155 152 156 ! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar 153 IF( present(kit) ) THEN ! ignore kn_fsbc in this case 154 isecsbc = nsec_year + nsec1jan000 + (kit+it_offset)*NINT( rdt/REAL(nn_baro,wp) ) 155 ELSE ! middle of sbc time step 156 isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) + it_offset * NINT(rdttra(1)) 157 IF( present(jit) ) THEN 158 ! ignore kn_fsbc in this case 159 isecsbc = nsec_year + nsec1jan000 + (jit+itime_add)*rdt/REAL(nn_baro,wp) 160 ELSE 161 isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) + itime_add * rdttra(1) ! middle of sbc time step 157 162 ENDIF 158 163 imf = SIZE( sd ) 159 164 ! 160 165 IF( ll_firstcall ) THEN ! initialization 161 DO jf = 1, imf 162 IF( PRESENT(map) ) imap = map(jf) 163 CALL fld_init( kn_fsbc, sd(jf), imap ) ! read each before field (put them in after as they will be swapped) 164 END DO 166 IF( PRESENT(map) ) THEN 167 DO jf = 1, imf 168 IF( PRESENT(jpk_1) ) THEN 169 CALL fld_init( kn_fsbc, sd(jf), map(jf)%ptr, jpk_1 ) ! read each before field (put them in after as they will be swapped) 170 ELSE 171 CALL fld_init( kn_fsbc, sd(jf), map(jf)%ptr ) ! read each before field (put them in after as they will be swapped) 172 ENDIF 173 END DO 174 ELSE 175 DO jf = 1, imf 176 CALL fld_init( kn_fsbc, sd(jf) ) ! read each before field (put them in after as they will be swapped) 177 END DO 178 ENDIF 165 179 IF( lwp ) CALL wgt_print() ! control print 180 CALL fld_rot( kt, sd ) ! rotate vector fiels if needed 166 181 ENDIF 167 182 ! ! ====================================== ! … … 171 186 DO jf = 1, imf ! --- loop over field --- ! 172 187 173 IF( isecsbc > sd(jf)%nrec_a(2) .OR. ll_firstcall ) THEN ! read/update the after data? 174 175 IF( PRESENT(map) ) imap = map(jf) ! temporary definition of map 176 177 sd(jf)%nrec_b(:) = sd(jf)%nrec_a(:) ! swap before record informations 178 sd(jf)%rotn(1) = sd(jf)%rotn(2) ! swap before rotate informations 179 IF( sd(jf)%ln_tint ) sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! swap before record field 180 181 CALL fld_rec( kn_fsbc, sd(jf), kt_offset = it_offset, kit = kit ) ! update after record informations 182 183 ! if kn_fsbc*rdttra is larger than nfreqh (which is kind of odd), 184 ! it is possible that the before value is no more the good one... we have to re-read it 185 ! if before is not the last record of the file currently opened and after is the first record to be read 186 ! in a new file which means after = 1 (the file to be opened corresponds to the current time) 187 ! or after = nreclast + 1 (the file to be opened corresponds to a future time step) 188 IF( .NOT. ll_firstcall .AND. sd(jf)%ln_tint .AND. sd(jf)%nrec_b(1) /= sd(jf)%nreclast & 189 & .AND. MOD( sd(jf)%nrec_a(1), sd(jf)%nreclast ) == 1 ) THEN 190 itmp = sd(jf)%nrec_a(1) ! temporary storage 191 sd(jf)%nrec_a(1) = sd(jf)%nreclast ! read the last record of the file currently opened 192 CALL fld_get( sd(jf), imap ) ! read after data 193 sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! re-swap before record field 194 sd(jf)%nrec_b(1) = sd(jf)%nrec_a(1) ! update before record informations 195 sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - sd(jf)%nfreqh * 3600 ! assume freq to be in hours in this case 196 sd(jf)%rotn(1) = sd(jf)%rotn(2) ! update before rotate informations 197 sd(jf)%nrec_a(1) = itmp ! move back to after record 188 IF( isecsbc > sd(jf)%nrec_a(2) .OR. ll_firstcall ) THEN ! read/update the after data? 189 190 IF( sd(jf)%ln_tint ) THEN ! swap before record field and informations 191 sd(jf)%nrec_b(:) = sd(jf)%nrec_a(:) 192 !CDIR COLLAPSE 193 sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) 198 194 ENDIF 199 195 200 CALL fld_clopn( sd(jf) ) ! Do we need to open a new year/month/week/day file? 201 196 IF( PRESENT(jit) ) THEN 197 CALL fld_rec( kn_fsbc, sd(jf), time_offset=itime_add, jit=jit ) ! update record informations 198 ELSE 199 CALL fld_rec( kn_fsbc, sd(jf), time_offset=itime_add ) ! update record informations 200 ENDIF 201 202 ! do we have to change the year/month/week/day of the forcing field?? 202 203 IF( sd(jf)%ln_tint ) THEN 203 204 ! if kn_fsbc*rdttra is larger than nfreqh (which is kind of odd),205 ! it is possible that the before value is no more the good one... we have to re-read it206 ! if before record is not just just before the after record...207 IF( .NOT. ll_firstcall .AND. MOD( sd(jf)%nrec_a(1), sd(jf)%nreclast ) /= 1 &208 & .AND. sd(jf)%nrec_b(1) /= sd(jf)%nrec_a(1) - 1 ) THEN209 sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) - 1 ! move back to before record210 CALL fld_get( sd(jf), imap ) ! read after data211 sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! re-swap before record field212 sd(jf)%nrec_b(1) = sd(jf)%nrec_a(1) ! update before record informations213 sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - sd(jf)%nfreqh * 3600 ! assume freq to be in hours in this case214 sd(jf)%rotn(1) = sd(jf)%rotn(2) ! update before rotate informations215 sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) + 1 ! move back to after record216 ENDIF217 218 ! do we have to change the year/month/week/day of the forcing field??219 204 ! if we do time interpolation we will need to open next year/month/week/day file before the end of the current 220 205 ! one. If so, we are still before the end of the year/month/week/day when calling fld_rec so sd(jf)%nrec_a(1) 221 206 ! will be larger than the record number that should be read for current year/month/week/day 207 208 ! last record to be read in the current file 209 IF ( sd(jf)%nfreqh == -12 ) THEN ; ireclast = 1 ! yearly mean 210 ELSEIF( sd(jf)%nfreqh == -1 ) THEN ! monthly mean 211 IF( sd(jf)%cltype == 'monthly' ) THEN ; ireclast = 1 212 ELSE ; ireclast = 12 213 ENDIF 214 ELSE ! higher frequency mean (in hours) 215 IF( sd(jf)%cltype == 'monthly' ) THEN ; ireclast = 24 * nmonth_len(nmonth) / sd(jf)%nfreqh 216 ELSEIF( sd(jf)%cltype(1:4) == 'week' ) THEN ; ireclast = 24 * 7 / sd(jf)%nfreqh 217 ELSEIF( sd(jf)%cltype == 'daily' ) THEN ; ireclast = 24 / sd(jf)%nfreqh 218 ELSE ; ireclast = 24 * nyear_len( 1 ) / sd(jf)%nfreqh 219 ENDIF 220 ENDIF 221 222 222 ! do we need next file data? 223 IF( sd(jf)%nrec_a(1) > sd(jf)%nreclast ) THEN224 225 sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) - sd(jf)%nreclast !226 227 IF( .NOT. ( sd(jf)%ln_clim .AND. sd(jf)%cltype == 'yearly' ) ) THEN ! close/open the current/new file228 223 IF( sd(jf)%nrec_a(1) > ireclast ) THEN 224 225 sd(jf)%nrec_a(1) = 1 ! force to read the first record of the next file 226 227 IF( .NOT. sd(jf)%ln_clim ) THEN ! close the current file and open a new one. 228 229 229 llnxtmth = sd(jf)%cltype == 'monthly' .OR. nday == nmonth_len(nmonth) ! open next month file? 230 230 llnxtyr = sd(jf)%cltype == 'yearly' .OR. (nmonth == 12 .AND. llnxtmth) ! open next year file? … … 235 235 isecend = nsec_year + nsec1jan000 + (nitend - kt) * NINT(rdttra(1)) ! second at the end of the run 236 236 llstop = isecend > sd(jf)%nrec_a(2) ! read more than 1 record of next year 237 ! we suppose that the date of next file is next day (should be ok even for weekly files...) 237 238 238 CALL fld_clopn( sd(jf), nyear + COUNT((/llnxtyr /)) , & 239 239 & nmonth + COUNT((/llnxtmth/)) - 12 * COUNT((/llnxtyr /)), & … … 243 243 CALL ctl_warn('next year/month/week/day file: '//TRIM(sd(jf)%clname)// & 244 244 & ' not present -> back to current year/month/day') 245 CALL fld_clopn( sd(jf) ) ! back to the current year/month/day246 sd(jf)%nrec_a(1) = sd(jf)%nreclast ! force to read the last record in the current year file245 CALL fld_clopn( sd(jf), nyear, nmonth, nday ) ! back to the current year/month/day 246 sd(jf)%nrec_a(1) = ireclast ! force to read the last record to be read in the current year file 247 247 ENDIF 248 248 249 249 ENDIF 250 ENDIF ! open need next file? 251 252 ENDIF ! temporal interpolation? 250 ENDIF 251 252 ELSE 253 ! if we are not doing time interpolation, we must change the year/month/week/day of the file just after 254 ! switching to the NEW year/month/week/day. If it is the case, we are at the beginning of the 255 ! year/month/week/day when calling fld_rec so sd(jf)%nrec_a(1) = 1 256 IF( sd(jf)%nrec_a(1) == 1 .AND. .NOT. ( sd(jf)%ln_clim .AND. sd(jf)%cltype == 'yearly' ) ) & 257 & CALL fld_clopn( sd(jf), nyear, nmonth, nday ) 258 ENDIF 253 259 254 260 ! read after data 255 CALL fld_get( sd(jf), imap ) 256 257 ENDIF ! read new data? 261 IF( PRESENT(map) ) THEN 262 IF( PRESENT(jpk_1) ) THEN 263 CALL fld_get( sd(jf), map(jf)%ptr, jpk_1) 264 ELSE 265 CALL fld_get( sd(jf), map(jf)%ptr) 266 ENDIF 267 ELSE 268 CALL fld_get( sd(jf) ) 269 ENDIF 270 271 ENDIF 258 272 END DO ! --- end loop over field --- ! 259 273 260 CALL fld_rot( kt, sd ) ! rotate vector before/now/after fields if needed274 CALL fld_rot( kt, sd ) ! rotate vector fiels if needed 261 275 262 276 DO jf = 1, imf ! --- loop over field --- ! … … 268 282 WRITE(numout, clfmt) TRIM( sd(jf)%clvar ), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 269 283 & sd(jf)%nrec_b(1), sd(jf)%nrec_a(1), REAL(sd(jf)%nrec_b(2),wp)/rday, REAL(sd(jf)%nrec_a(2),wp)/rday 270 WRITE(numout, *) 'it _offset is : ',it_offset284 WRITE(numout, *) 'itime_add is : ',itime_add 271 285 ENDIF 272 286 ! temporal interpolation weights … … 295 309 296 310 297 SUBROUTINE fld_init( kn_fsbc, sdjf, map )311 SUBROUTINE fld_init( kn_fsbc, sdjf, map , jpk_1 ) 298 312 !!--------------------------------------------------------------------- 299 313 !! *** ROUTINE fld_init *** … … 304 318 INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step) 305 319 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 306 TYPE(MAP_POINTER),INTENT(in) :: map ! global-to-local mapping indices 320 INTEGER , INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping indices 321 INTEGER , INTENT(in), OPTIONAL :: jpk_1 ! global-to-local mapping indices 307 322 !! 308 323 LOGICAL :: llprevyr ! are we reading previous year file? … … 317 332 CHARACTER(LEN=1000) :: clfmt ! write format 318 333 !!--------------------------------------------------------------------- 334 335 ! some default definitions... 336 sdjf%num = 0 ! default definition for non-opened file 337 IF( sdjf%ln_clim ) sdjf%clname = TRIM( sdjf%clrootname ) ! file name defaut definition, never change in this case 319 338 llprevyr = .FALSE. 320 339 llprevmth = .FALSE. … … 323 342 isec_week = 0 324 343 344 IF( sdjf%cltype(1:4) == 'week' .AND. nn_leapy == 0 ) & 345 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdjf%clrootname)//') needs nn_leapy = 1') 346 IF( sdjf%cltype(1:4) == 'week' .AND. sdjf%ln_clim ) & 347 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdjf%clrootname)//') needs ln_clim = .FALSE.') 348 325 349 ! define record informations 326 350 CALL fld_rec( kn_fsbc, sdjf, ldbefore = .TRUE. ) ! return before values in sdjf%nrec_a (as we will swap it later) … … 336 360 llprevyr = .NOT. sdjf%ln_clim ! use previous year file? 337 361 ELSE 338 CALL ctl_stop( "fld_init: yearly mean file must be in a yearly type of file: "//TRIM(sdjf%cl rootname) )362 CALL ctl_stop( "fld_init: yearly mean file must be in a yearly type of file: "//TRIM(sdjf%clname) ) 339 363 ENDIF 340 364 ELSEIF( sdjf%nfreqh == -1 ) THEN ! monthly mean … … 343 367 llprevmth = .TRUE. ! use previous month file? 344 368 llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file? 369 ! IF (lwp) write(numout,*) sdjf%clvar,'AFTER', sdjf%nrec_a(1), sdjf%nrec_a(2), sdjf%clname 345 370 ELSE ! yearly file 346 371 sdjf%nrec_a(1) = 12 ! force to read december mean … … 367 392 ENDIF 368 393 ENDIF 369 !370 394 IF ( sdjf%cltype(1:4) == 'week' ) THEN 371 395 isec_week = isec_week + ksec_week( sdjf%cltype(6:8) ) ! second since the beginning of the week … … 383 407 ! if previous year/month/day file does not exist, we switch to the current year/month/day 384 408 IF( llprev .AND. sdjf%num <= 0 ) THEN 385 CALL ctl_warn( 'previous year/month/week/day file: '//TRIM(sdjf%cl rootname)// &409 CALL ctl_warn( 'previous year/month/week/day file: '//TRIM(sdjf%clname)// & 386 410 & ' not present -> back to current year/month/week/day' ) 387 411 ! we force to read the first record of the current year/month/day instead of last record of previous year/month/day 388 412 llprev = .FALSE. 389 413 sdjf%nrec_a(1) = 1 390 CALL fld_clopn( sdjf )414 CALL fld_clopn( sdjf, nyear, nmonth, nday ) 391 415 ENDIF 392 416 393 IF( llprev ) THEN ! check if the record sdjf%nrec_a(1) exists in the file417 IF( llprev ) THEN ! check if the last record sdjf%nrec_n(1) exists in the file 394 418 idvar = iom_varid( sdjf%num, sdjf%clvar ) ! id of the variable sdjf%clvar 395 419 IF( idvar <= 0 ) RETURN … … 398 422 ENDIF 399 423 400 ! read before data in after arrays(as we will swap it later) 401 CALL fld_get( sdjf, map ) 424 ! read before data 425 IF( PRESENT(map) ) THEN 426 IF( PRESENT(jpk_1) ) THEN 427 CALL fld_get( sdjf, map , jpk_1) ! read before values in after arrays(as we will swap it later) 428 ELSE 429 CALL fld_get( sdjf, map ) ! read before values in after arrays(as we will swap it later) 430 ENDIF 431 ELSE 432 CALL fld_get( sdjf ) ! read before values in after arrays(as we will swap it later) 433 ENDIF 402 434 403 435 clfmt = "('fld_init : time-interpolation for ', a, ' read previous record = ', i4, ' at time = ', f7.2, ' days')" 404 436 IF(lwp) WRITE(numout, clfmt) TRIM(sdjf%clvar), sdjf%nrec_a(1), REAL(sdjf%nrec_a(2),wp)/rday 405 437 438 IF( llprev ) CALL iom_close( sdjf%num ) ! force to close previous year file (-> redefine sdjf%num to 0) 439 406 440 ENDIF 441 442 ! make sure current year/month/day file is opened 443 IF( sdjf%num <= 0 ) THEN 444 ! 445 IF ( sdjf%cltype(1:4) == 'week' ) THEN 446 isec_week = ksec_week( sdjf%cltype(6:8) ) ! second since the beginning of the week 447 llprevmth = isec_week > nsec_month ! longer time since beginning of the week than the month 448 llprevyr = llprevmth .AND. nmonth == 1 449 ELSE 450 isec_week = 0 451 llprevmth = .FALSE. 452 llprevyr = .FALSE. 453 ENDIF 454 ! 455 iyear = nyear - COUNT((/llprevyr /)) 456 imonth = nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /)) 457 iday = nday + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week / NINT(rday) 458 ! 459 CALL fld_clopn( sdjf, iyear, imonth, iday ) 460 ENDIF 407 461 ! 408 462 END SUBROUTINE fld_init 409 463 410 464 411 SUBROUTINE fld_rec( kn_fsbc, sdjf, ldbefore, kit, kt_offset )465 SUBROUTINE fld_rec( kn_fsbc, sdjf, ldbefore, jit, time_offset ) 412 466 !!--------------------------------------------------------------------- 413 467 !! *** ROUTINE fld_rec *** … … 423 477 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 424 478 LOGICAL , INTENT(in ), OPTIONAL :: ldbefore ! sent back before record values (default = .FALSE.) 425 INTEGER , INTENT(in ), OPTIONAL :: kit ! index of barotropic subcycle479 INTEGER , INTENT(in ), OPTIONAL :: jit ! index of barotropic subcycle 426 480 ! used only if sdjf%ln_tint = .TRUE. 427 INTEGER , INTENT(in ), OPTIONAL :: kt_offset! Offset of required time level compared to "now"428 !time level in units of time steps.481 INTEGER , INTENT(in ), OPTIONAL :: time_offset ! Offset of required time level compared to "now" 482 ! time level in units of time steps. 429 483 !! 430 484 LOGICAL :: llbefore ! local definition of ldbefore … … 433 487 INTEGER :: ifreq_sec ! frequency mean (in seconds) 434 488 INTEGER :: isec_week ! number of seconds since the start of the weekly file 435 INTEGER :: it _offset! local time offset variable489 INTEGER :: itime_add ! local time offset variable 436 490 REAL(wp) :: ztmp ! temporary variable 437 491 !!---------------------------------------------------------------------- … … 443 497 ENDIF 444 498 ! 445 it_offset = 0 446 IF( PRESENT(kt_offset) ) it_offset = kt_offset 447 IF( PRESENT(kit) ) THEN ; it_offset = ( kit + it_offset ) * NINT( rdt/REAL(nn_baro,wp) ) 448 ELSE ; it_offset = it_offset * NINT( rdttra(1) ) 449 ENDIF 499 itime_add = 0 500 IF( PRESENT(time_offset) ) itime_add = time_offset 450 501 ! 451 502 ! ! =========== ! … … 465 516 ! forcing record : 1 466 517 ! 467 ztmp = REAL( nday, wp ) / REAL( nyear_len(1), wp ) + 0.5 + REAL( it_offset, wp ) 518 ztmp = REAL( nday, wp ) / REAL( nyear_len(1), wp ) + 0.5 519 IF( PRESENT(jit) ) THEN 520 ztmp = ztmp + (jit+itime_add)*rdt/REAL(nn_baro,wp) 521 ELSE 522 ztmp = ztmp + itime_add*rdttra(1) 523 ENDIF 468 524 sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) 469 525 ! swap at the middle of the year … … 493 549 ! forcing record : nmonth 494 550 ! 495 ztmp = REAL( nday, wp ) / REAL( nmonth_len(nmonth), wp ) + 0.5 + REAL( it_offset, wp ) 551 ztmp = REAL( nday, wp ) / REAL( nmonth_len(nmonth), wp ) + 0.5 552 IF( PRESENT(jit) ) THEN 553 ztmp = ztmp + (jit+itime_add)*rdt/REAL(nn_baro,wp) / ( REAL( nmonth_len(nmonth), wp )* 86400. ) 554 ELSE 555 ztmp = ztmp + itime_add*rdttra(1) / ( REAL( nmonth_len(nmonth), wp ) * 86400. ) 556 ENDIF 496 557 imth = nmonth + INT( ztmp ) - COUNT((/llbefore/)) 497 558 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) … … 499 560 ENDIF 500 561 sdjf%nrec_a(2) = nmonth_half( imth ) + nsec1jan000 ! swap at the middle of the month 562 ! IF (lwp) write(numout,*) sdjf%clvar, sdjf%nrec_a(1), sdjf%nrec_a(2), nday, nmonth, itime_add, & 563 ! rdttra(1), COUNT((/llbefore/)), ztmp, nmonth_half( imth ), & 564 ! nsec1jan000, REAL( nmonth_len(nmonth), wp ) 501 565 ELSE ! no time interpolation 502 566 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nrec_a(1) = 1 … … 519 583 ELSE ; ztmp = REAL(nsec_year ,wp) ! since 00h on Jan 1 of the current year 520 584 ENDIF 521 ztmp = ztmp + 0.5 * REAL(kn_fsbc - 1, wp) * rdttra(1) + REAL( it_offset, wp ) ! centrered in the middle of sbc time step 522 ztmp = ztmp + 0.01 * rdttra(1) ! avoid truncation error 585 ztmp = ztmp + 0.5 * REAL(kn_fsbc - 1, wp) * rdttra(1) ! shift time to be centrered in the middle of sbc time step 586 ztmp = ztmp + 0.01 * rdttra(1) ! add 0.01 time step to avoid truncation error 587 IF( PRESENT(jit) ) THEN 588 ztmp = ztmp + (jit+itime_add)*rdt/REAL(nn_baro,wp) 589 ELSE 590 ztmp = ztmp + itime_add*rdttra(1) 591 ENDIF 523 592 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record 524 593 ! 525 ! INT( ztmp/ifreq_sec + 0.5)594 ! INT( ztmp ) 526 595 ! /|\ 527 596 ! 2 | *-----( … … 529 598 ! 0 |--( 530 599 ! |--+--|--+--|--+--|--> time 531 ! 0 /|\ 1 /|\ 2 /|\ 3 (ztmp/ifreq_sec)600 ! 0 /|\ 1 /|\ 2 /|\ 3 (nsec_year/ifreq_sec) or (nsec_month/ifreq_sec) 532 601 ! | | | 533 602 ! | | | … … 537 606 ELSE ! no time interpolation 538 607 ! 539 ! INT( ztmp/ifreq_sec)608 ! INT( ztmp ) 540 609 ! /|\ 541 610 ! 2 | *-----( … … 543 612 ! 0 |-----( 544 613 ! |--+--|--+--|--+--|--> time 545 ! 0 /|\ 1 /|\ 2 /|\ 3 (ztmp/ifreq_sec)614 ! 0 /|\ 1 /|\ 2 /|\ 3 (nsec_year/ifreq_sec) or (nsec_month/ifreq_sec) 546 615 ! | | | 547 616 ! | | | … … 550 619 ztmp= ztmp / REAL(ifreq_sec, wp) 551 620 ENDIF 552 sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) ! record n umber to be read621 sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) ! record nomber to be read 553 622 554 623 iendrec = ifreq_sec * sdjf%nrec_a(1) + nsec1jan000 ! end of this record (in second) … … 569 638 570 639 571 SUBROUTINE fld_get( sdjf, map )640 SUBROUTINE fld_get( sdjf, map, jpk_1 ) 572 641 !!--------------------------------------------------------------------- 573 642 !! *** ROUTINE fld_get *** … … 576 645 !!---------------------------------------------------------------------- 577 646 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 578 TYPE(MAP_POINTER),INTENT(in) :: map ! global-to-local mapping indices 647 INTEGER , INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping indices 648 INTEGER , INTENT(in), OPTIONAL :: jpk_1 ! number of levels in bdy data 579 649 !! 580 650 INTEGER :: ipk ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 581 651 INTEGER :: iw ! index into wgts array 582 INTEGER :: ipdom ! index of the domain 583 !!--------------------------------------------------------------------- 584 ! 652 !!--------------------------------------------------------------------- 653 585 654 ipk = SIZE( sdjf%fnow, 3 ) 586 ! 587 IF( ASSOCIATED(map%ptr) ) THEN 588 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map%ptr ) 589 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), map%ptr ) 655 656 IF( PRESENT(map) ) THEN 657 IF( PRESENT(jpk_1) ) THEN 658 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map, jpk_1 ) 659 IF(lwp) WRITE(numout,*) 'in get 2' 660 CALL flush(numout) 661 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), map, jpk_1 ) 662 IF(lwp) WRITE(numout,*) 'in get 1' 663 CALL flush(numout) 664 ENDIF 665 ELSE 666 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map ) 667 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), map ) 668 ENDIF 590 669 ENDIF 591 670 ELSE IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN … … 595 674 ENDIF 596 675 ELSE 597 IF( SIZE(sdjf%fnow, 1) == jpi ) THEN ; ipdom = jpdom_data598 ELSE ; ipdom = jpdom_unknown599 ENDIF600 676 SELECT CASE( ipk ) 601 CASE(1) 602 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_a(1) )603 ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,1 ), sdjf%nrec_a(1) )677 CASE(1) 678 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_a(1) ) 679 ELSE ; CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fnow(:,:,1 ), sdjf%nrec_a(1) ) 604 680 ENDIF 605 681 CASE DEFAULT 606 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1) )607 ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1) )682 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1) ) 683 ELSE ; CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1) ) 608 684 ENDIF 609 685 END SELECT 610 686 ENDIF 611 687 ! 612 sdjf%rotn (2)= .false. ! vector not yet rotated688 sdjf%rotn = .false. ! vector not yet rotated 613 689 614 690 END SUBROUTINE fld_get 615 691 616 SUBROUTINE fld_map( num, clvar, dta, nrec, map )617 !!--------------------------------------------------------------------- 618 !! *** ROUTINE fld_ map***692 SUBROUTINE fld_map( num, clvar, dta, nrec, map, jpk_1 ) 693 !!--------------------------------------------------------------------- 694 !! *** ROUTINE fld_get *** 619 695 !! 620 696 !! ** Purpose : read global data from file and map onto local data 621 697 !! using a general mapping (for open boundaries) 698 !! 699 !! 12-04-13 updated to include interpolation of boundary 700 !! data from non-native vertical grid 622 701 !!---------------------------------------------------------------------- 623 702 #if defined key_bdy 624 USE bdy_oce, ONLY: dta_global, dta_global 2! workspace to read in global data arrays703 USE bdy_oce, ONLY: dta_global, dta_global_1, dta_global_2, idx_bdy ! workspace to read in global data arrays 625 704 #endif 626 INTEGER , INTENT(in ) :: num ! stream number 627 CHARACTER(LEN=*) , INTENT(in ) :: clvar ! variable name 628 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional) 629 INTEGER , INTENT(in ) :: nrec ! record number to read (ie time slice) 630 INTEGER, DIMENSION(:) , INTENT(in ) :: map ! global-to-local mapping indices 631 !! 632 INTEGER :: ipi ! length of boundary data on local process 633 INTEGER :: ipj ! length of dummy dimension ( = 1 ) 634 INTEGER :: ipk ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk ) 635 INTEGER :: ilendta ! length of data in file 636 INTEGER :: idvar ! variable ID 637 INTEGER :: ib, ik, ji, jj ! loop counters 705 706 INTEGER , INTENT(in ) :: num ! stream number 707 CHARACTER(LEN=*) , INTENT(in ) :: clvar ! variable name 708 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional) 709 INTEGER , INTENT(in ) :: nrec ! record number to read (ie time slice) 710 INTEGER, DIMENSION(:) , INTENT(in ) :: map ! global-to-local mapping indices 711 INTEGER , INTENT(in), OPTIONAL :: jpk_1 ! number of levels in bdy data 712 INTEGER :: jpkm1_1 ! number of levels in bdy data minus 1 713 !! 714 INTEGER :: ipi ! length of boundary data on local process 715 INTEGER :: ipj ! length of dummy dimension ( = 1 ) 716 INTEGER :: ipk, ipkm1 ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk ) 717 INTEGER :: ilendta ! length of data in file 718 INTEGER :: idvar ! variable ID 719 INTEGER :: ib, ik, ikk! loop counters 638 720 INTEGER :: ierr 639 REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read ! work space for global data 721 INTEGER :: igrd, ib_bdy 722 REAL(wp) :: zl, zi ! tmp variable for current depth and interpolation factor 723 REAL(wp) :: fv, fv_alt ! fillvalue and alternative -ABS(fv) 724 REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read ! work space for global data 725 REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read_1 ! work space for BDY data from file 726 REAL(wp), POINTER, DIMENSION(:,:) :: dta_read_2 ! work space for BDY depth data from file 640 727 !!--------------------------------------------------------------------- 641 728 729 #if defined key_bdy 730 dta_read => dta_global 731 IF( PRESENT(jpk_1) ) THEN 732 IF( jpk_1>0 ) THEN 733 dta_read_1 => dta_global_1 734 dta_read_2 => dta_global_2 735 jpkm1_1 = jpk_1 - 1 736 ENDIF 737 ENDIF 738 igrd = 1 ! T/S only so far 739 ib_bdy = 1 ! and only one bdy file 740 #endif 741 642 742 ipi = SIZE( dta, 1 ) 643 743 ipj = 1 644 744 ipk = SIZE( dta, 3 ) 745 ipkm1 = ipk - 1 645 746 646 747 idvar = iom_varid( num, clvar ) 647 748 ilendta = iom_file(num)%dimsz(1,idvar) 648 649 #if defined key_bdy650 ipj = iom_file(num)%dimsz(2,idvar)651 IF (ipj == 1) THEN ! we assume that this is a structured open boundary file652 dta_read => dta_global653 ELSE654 dta_read => dta_global2655 ENDIF656 #endif657 658 749 IF(lwp) WRITE(numout,*) 'Dim size for ',TRIM(clvar),' is ', ilendta 659 750 IF(lwp) WRITE(numout,*) 'Number of levels for ',TRIM(clvar),' is ', ipk 751 CALL flush(numout) 660 752 661 753 SELECT CASE( ipk ) 662 CASE(1) ; CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1 ), nrec ) 663 CASE DEFAULT ; CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec ) 754 CASE(1) 755 CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1 ), nrec ) 756 CASE DEFAULT 757 #if defined key_bdy 758 IF( PRESENT(jpk_1) ) THEN ! boundary data not on model grid: veritcal interpolation 759 IF( jpk_1>0 ) THEN 760 IF( lwp )THEN 761 WRITE(numout,*) 'BDY: interpolate T & S data onto new vertical mesh' 762 ENDIF 763 ! 764 ! gather data from file along with depth and _FillValue info 765 ! 766 CALL iom_get ( num, jpdom_unknown, clvar, dta_read_1(1:ilendta,1:ipj,1:jpk_1), nrec ) 767 CALL iom_get ( num, jpdom_unknown, 'deptht', dta_read_2(1:ilendta,1:jpk_1) ) 768 CALL iom_getatt(num, '_FillValue', fv, cdvar=clvar ) 769 ! 770 fv_alt = -ABS(fv) ! set _FillValue < 0 as we make use of MAXVAL and MAXLOC later 771 ! 772 DO ib = 1, ipi 773 DO ik = 1, ipk 774 IF( ( dta_read_1(map(ib),1,ik) == fv ) ) THEN 775 dta_read_2(map(ib),ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data 776 ENDIF 777 dta_read(map(ib),1,ik) = fv_alt ! put fillvalue into new field as if all goes well all wet points will be replaced 778 ENDDO 779 ENDDO ! had to use map in this loop ?? tried looping over ib but failed !! investigate TODO 780 ! 781 DO ib = 1, ipi 782 DO ik = 1, ipk 783 zl = gdept_1(idx_bdy(ib_bdy)%nbi(ib,igrd),idx_bdy(ib_bdy)%nbj(ib,igrd),ik) ! if using in step could use fsdept instead of gdept_1? 784 IF( zl < dta_read_2(map(ib),1) ) THEN ! above the first level of external data 785 dta_read(map(ib),1,ik) = dta_read_1(map(ib),1,1) 786 ELSEIF( zl > MAXVAL(dta_read_2(map(ib),:),1) ) THEN ! below the last level of external data 787 dta_read(map(ib),1,ik) = dta_read_1(map(ib),1,MAXLOC(dta_read_2(map(ib),:),1)) 788 ELSE ! inbetween : vertical interpolation between ikk & ikk+1 789 DO ikk = 1, ipkm1 ! when gdept_1(ikk) < zl < gdept_1(ikk+1) 790 IF( ( (zl-dta_read_2(map(ib),ikk)) * (zl-dta_read_2(map(ib),ikk+1)) <= 0._wp) & 791 & .AND. (dta_read_2(map(ib),ikk+1) /= fv_alt)) THEN 792 zi = ( zl - dta_read_2(map(ib),ikk) ) / (dta_read_2(map(ib),ikk+1)-dta_read_2(map(ib),ikk)) 793 dta_read(map(ib),1,ik) = dta_read_1(map(ib),1,ikk) + & 794 & ( dta_read_1(map(ib),1,ikk+1) - dta_read_1(map(ib),1,ikk) ) * zi 795 ENDIF 796 END DO 797 ENDIF 798 END DO 799 END DO 800 ! 801 IF(lwp) WRITE(numout,*) 'BDY: finished interpolating T & S data onto new vertical mesh' 802 ! 803 ENDIF ! is jpk_1 > 0 804 ELSE ! must be on model grid already 805 CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec ) 806 ENDIF ! end PRESENT jpk_1 807 #else 808 CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec ) 809 #endif 664 810 END SELECT 665 811 ! 666 IF (ipj==1) THEN 667 DO ib = 1, ipi 668 DO ik = 1, ipk 669 dta(ib,1,ik) = dta_read(map(ib),1,ik) 670 END DO 812 DO ib = 1, ipi 813 DO ik = 1, ipk 814 dta(ib,1,ik) = dta_read(map(ib),1,ik) 671 815 END DO 672 ELSE ! we assume that this is a structured open boundary file 673 DO ib = 1, ipi 674 jj=1+floor(REAL(map(ib)-1)/REAL(ilendta)) 675 ji=map(ib)-(jj-1)*ilendta 676 DO ik = 1, ipk 677 dta(ib,1,ik) = dta_read(ji,jj,ik) 678 END DO 679 END DO 680 ENDIF 816 END DO 681 817 682 818 END SUBROUTINE fld_map … … 692 828 TYPE(FLD), INTENT(inout), DIMENSION(:) :: sd ! input field related variables 693 829 !! 694 INTEGER :: ju, jv,jk,jn! loop indices830 INTEGER :: ju, jv, jk ! loop indices 695 831 INTEGER :: imf ! size of the structure sd 696 832 INTEGER :: ill ! character length … … 707 843 DO ju = 1, imf 708 844 ill = LEN_TRIM( sd(ju)%vcomp ) 709 DO jn = 2-COUNT((/sd(ju)%ln_tint/)), 2 710 IF( ill > 0 .AND. .NOT. sd(ju)%rotn(jn) ) THEN ! find vector rotations required 711 IF( sd(ju)%vcomp(1:1) == 'U' ) THEN ! east-west component has symbolic name starting with 'U' 712 ! look for the north-south component which has same symbolic name but with 'U' replaced with 'V' 713 clcomp = 'V' // sd(ju)%vcomp(2:ill) ! works even if ill == 1 714 iv = -1 715 DO jv = 1, imf 716 IF( TRIM(sd(jv)%vcomp) == TRIM(clcomp) ) iv = jv 717 END DO 718 IF( iv > 0 ) THEN ! fields ju and iv are two components which need to be rotated together 719 DO jk = 1, SIZE( sd(ju)%fnow, 3 ) 720 IF( sd(ju)%ln_tint )THEN 721 CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->i', utmp(:,:) ) 722 CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->j', vtmp(:,:) ) 723 sd(ju)%fdta(:,:,jk,jn) = utmp(:,:) ; sd(iv)%fdta(:,:,jk,jn) = vtmp(:,:) 724 ELSE 725 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->i', utmp(:,:) ) 726 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->j', vtmp(:,:) ) 727 sd(ju)%fnow(:,:,jk ) = utmp(:,:) ; sd(iv)%fnow(:,:,jk ) = vtmp(:,:) 728 ENDIF 729 END DO 730 sd(ju)%rotn(jn) = .TRUE. ! vector was rotated 731 IF( lwp .AND. kt == nit000 ) WRITE(numout,*) & 732 & 'fld_read: vector pair ('//TRIM(sd(ju)%clvar)//', '//TRIM(sd(iv)%clvar)//') rotated on to model grid' 733 ENDIF 734 ENDIF 735 ENDIF 736 END DO 845 IF( ill > 0 .AND. .NOT. sd(ju)%rotn ) THEN ! find vector rotations required 846 IF( sd(ju)%vcomp(1:1) == 'U' ) THEN ! east-west component has symbolic name starting with 'U' 847 ! look for the north-south component which has same symbolic name but with 'U' replaced with 'V' 848 clcomp = 'V' // sd(ju)%vcomp(2:ill) ! works even if ill == 1 849 iv = -1 850 DO jv = 1, imf 851 IF( TRIM(sd(jv)%vcomp) == TRIM(clcomp) ) iv = jv 852 END DO 853 IF( iv > 0 ) THEN ! fields ju and iv are two components which need to be rotated together 854 DO jk = 1, SIZE( sd(ju)%fnow, 3 ) 855 IF( sd(ju)%ln_tint )THEN 856 CALL rot_rep( sd(ju)%fdta(:,:,jk,2), sd(iv)%fdta(:,:,jk,2), 'T', 'en->i', utmp(:,:) ) 857 CALL rot_rep( sd(ju)%fdta(:,:,jk,2), sd(iv)%fdta(:,:,jk,2), 'T', 'en->j', vtmp(:,:) ) 858 sd(ju)%fdta(:,:,jk,2) = utmp(:,:) ; sd(iv)%fdta(:,:,jk,2) = vtmp(:,:) 859 ELSE 860 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->i', utmp(:,:) ) 861 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->j', vtmp(:,:) ) 862 sd(ju)%fnow(:,:,jk ) = utmp(:,:) ; sd(iv)%fnow(:,:,jk ) = vtmp(:,:) 863 ENDIF 864 END DO 865 sd(ju)%rotn = .TRUE. ! vector was rotated 866 IF( lwp .AND. kt == nit000 ) WRITE(numout,*) & 867 & 'fld_read: vector pair ('//TRIM(sd(ju)%clvar)//', '//TRIM(sd(iv)%clvar)//') rotated on to model grid' 868 ENDIF 869 ENDIF 870 ENDIF 737 871 END DO 738 872 ! … … 749 883 !!---------------------------------------------------------------------- 750 884 TYPE(FLD) , INTENT(inout) :: sdjf ! input field related variables 751 INTEGER , OPTIONAL, INTENT(in ) :: kyear ! year value752 INTEGER , OPTIONAL, INTENT(in ) :: kmonth ! month value753 INTEGER , OPTIONAL, INTENT(in ) :: kday ! day value885 INTEGER , INTENT(in ) :: kyear ! year value 886 INTEGER , INTENT(in ) :: kmonth ! month value 887 INTEGER , INTENT(in ) :: kday ! day value 754 888 LOGICAL, OPTIONAL, INTENT(in ) :: ldstop ! stop if open to read a non-existing file (default = .TRUE.) 755 !! 756 LOGICAL :: llprevyr ! are we reading previous year file? 757 LOGICAL :: llprevmth ! are we reading previous month file? 758 INTEGER :: iyear, imonth, iday ! first day of the current file in yyyy mm dd 759 INTEGER :: isec_week ! number of seconds since start of the weekly file 760 INTEGER :: indexyr ! year undex (O/1/2: previous/current/next) 761 INTEGER :: iyear_len, imonth_len ! length (days) of iyear and imonth ! 762 CHARACTER(len = 256):: clname ! temporary file name 763 !!---------------------------------------------------------------------- 764 IF( PRESENT(kyear) ) THEN ! use given values 765 iyear = kyear 766 imonth = kmonth 767 iday = kday 768 ELSE ! use current day values 769 IF ( sdjf%cltype(1:4) == 'week' ) THEN ! find the day of the beginning of the week 770 isec_week = ksec_week( sdjf%cltype(6:8) ) ! second since the beginning of the week 771 llprevmth = isec_week > nsec_month ! longer time since beginning of the week than the month 772 llprevyr = llprevmth .AND. nmonth == 1 773 ELSE 774 isec_week = 0 775 llprevmth = .FALSE. 776 llprevyr = .FALSE. 777 ENDIF 778 iyear = nyear - COUNT((/llprevyr /)) 779 imonth = nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /)) 780 iday = nday + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week / NINT(rday) 781 ENDIF 782 889 !!---------------------------------------------------------------------- 890 891 IF( sdjf%num /= 0 ) CALL iom_close( sdjf%num ) ! close file if already open 783 892 ! build the new filename if not climatological data 784 clname=TRIM(sdjf%clrootname)785 ! 786 ! note that sdjf%ln_clim is is only acting on the presence of the year in the file name893 sdjf%clname=TRIM(sdjf%clrootname) 894 ! 895 ! note that sdjf%ln_clim is is only acting on presence of the year in the file 787 896 IF( .NOT. sdjf%ln_clim ) THEN 788 WRITE( clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), iyear ! add year789 IF( sdjf%cltype /= 'yearly' ) WRITE( clname, '(a,"m" ,i2.2)' ) TRIM( clname ), imonth ! add month897 WRITE(sdjf%clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), kyear ! add year 898 IF( sdjf%cltype /= 'yearly' ) WRITE(sdjf%clname, '(a,"m" ,i2.2)' ) TRIM( sdjf%clname ), kmonth ! add month 790 899 ELSE 791 900 ! build the new filename if climatological data 792 IF( sdjf%cltype /= 'yearly' ) WRITE( clname, '(a,"_m",i2.2)' ) TRIM( sdjf%clrootname ), imonth ! add month901 IF( sdjf%cltype /= 'yearly' ) WRITE(sdjf%clname, '(a,"_m",i2.2)' ) TRIM( sdjf%clrootname ), kmonth ! add month 793 902 ENDIF 794 903 IF( sdjf%cltype == 'daily' .OR. sdjf%cltype(1:4) == 'week' ) & 795 & WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname ), iday ! add day 796 ! 797 IF( TRIM(clname) /= TRIM(sdjf%clname) .OR. sdjf%num == 0 ) THEN ! new file to be open 798 799 sdjf%clname = TRIM(clname) 800 IF( sdjf%num /= 0 ) CALL iom_close( sdjf%num ) ! close file if already open 801 CALL iom_open( sdjf%clname, sdjf%num, ldstop = ldstop, ldiof = LEN(TRIM(sdjf%wgtname)) > 0 ) 802 803 ! find the last record to be read -> update sdjf%nreclast 804 indexyr = iyear - nyear + 1 805 iyear_len = nyear_len( indexyr ) 806 SELECT CASE ( indexyr ) 807 CASE ( 0 ) ; imonth_len = 31 ! previous year -> imonth = 12 808 CASE ( 1 ) ; imonth_len = nmonth_len(imonth) 809 CASE ( 2 ) ; imonth_len = 31 ! next year -> imonth = 1 810 END SELECT 811 812 ! last record to be read in the current file 813 IF ( sdjf%nfreqh == -12 ) THEN ; sdjf%nreclast = 1 ! yearly mean 814 ELSEIF( sdjf%nfreqh == -1 ) THEN ! monthly mean 815 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = 1 816 ELSE ; sdjf%nreclast = 12 817 ENDIF 818 ELSE ! higher frequency mean (in hours) 819 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = 24 * imonth_len / sdjf%nfreqh 820 ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ; sdjf%nreclast = 24 * 7 / sdjf%nfreqh 821 ELSEIF( sdjf%cltype == 'daily' ) THEN ; sdjf%nreclast = 24 / sdjf%nfreqh 822 ELSE ; sdjf%nreclast = 24 * iyear_len / sdjf%nfreqh 823 ENDIF 824 ENDIF 825 826 ENDIF 827 ! 904 & WRITE(sdjf%clname, '(a,"d" ,i2.2)' ) TRIM( sdjf%clname ), kday ! add day 905 ! 906 CALL iom_open( sdjf%clname, sdjf%num, ldstop = ldstop, ldiof = LEN(TRIM(sdjf%wgtname)) > 0 ) 907 ! 828 908 END SUBROUTINE fld_clopn 829 909 … … 847 927 DO jf = 1, SIZE(sdf) 848 928 sdf(jf)%clrootname = TRIM( cdir )//TRIM( sdf_n(jf)%clname ) 849 sdf(jf)%clname = "not yet defined"850 929 sdf(jf)%nfreqh = sdf_n(jf)%nfreqh 851 930 sdf(jf)%clvar = sdf_n(jf)%clvar … … 853 932 sdf(jf)%ln_clim = sdf_n(jf)%ln_clim 854 933 sdf(jf)%cltype = sdf_n(jf)%cltype 855 sdf(jf)%num = -1 856 sdf(jf)%wgtname = " " 934 sdf(jf)%wgtname = " " 857 935 IF( LEN( TRIM(sdf_n(jf)%wname) ) > 0 ) sdf(jf)%wgtname = TRIM( cdir )//TRIM( sdf_n(jf)%wname ) 858 sdf(jf)%vcomp = sdf_n(jf)%vcomp 859 sdf(jf)%rotn(:) = .TRUE. ! pretend to be rotated -> won't try to rotate data before the first call to fld_get 860 IF( sdf(jf)%cltype(1:4) == 'week' .AND. nn_leapy == 0 ) & 861 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs nn_leapy = 1') 862 IF( sdf(jf)%cltype(1:4) == 'week' .AND. sdf(jf)%ln_clim ) & 863 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs ln_clim = .FALSE.') 936 sdf(jf)%vcomp = sdf_n(jf)%vcomp 937 sdf(jf)%rotn = .TRUE. 864 938 END DO 865 939
Note: See TracChangeset
for help on using the changeset viewer.