- Timestamp:
- 2019-11-26T15:11:43+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/ENHANCE-02_ISF_nemo_TEST_MERGE/src/OCE/SBC/fldread.F90
r10425 r11967 48 48 TYPE, PUBLIC :: FLD_N !: Namelist field informations 49 49 CHARACTER(len = 256) :: clname ! generic name of the NetCDF flux file 50 REAL(wp) :: nfreqh! frequency of each flux file50 REAL(wp) :: freqh ! frequency of each flux file 51 51 CHARACTER(len = 34) :: clvar ! generic name of the variable in the NetCDF flux file 52 52 LOGICAL :: ln_tint ! time interpolation or not (T/F) … … 64 64 CHARACTER(len = 256) :: clrootname ! generic name of the NetCDF file 65 65 CHARACTER(len = 256) :: clname ! current name of the NetCDF file 66 REAL(wp) :: nfreqh! frequency of each flux file66 REAL(wp) :: freqh ! frequency of each flux file 67 67 CHARACTER(len = 34) :: clvar ! generic name of the variable in the NetCDF flux file 68 68 LOGICAL :: ln_tint ! time interpolation or not (T/F) … … 80 80 INTEGER :: nreclast ! last record to be read in the current file 81 81 CHARACTER(len = 256) :: lsmname ! current name of the NetCDF mask file acting as a key 82 INTEGER :: igrd ! grid type for bdy data 83 INTEGER :: ibdy ! bdy set id number 82 ! ! 83 ! ! Variables related to BDY 84 INTEGER :: igrd ! grid type for bdy data 85 INTEGER :: ibdy ! bdy set id number 86 INTEGER, POINTER, DIMENSION(:) :: imap ! Array of integer pointers to 1D arrays 87 LOGICAL :: ltotvel ! total velocity or not (T/F) 88 LOGICAL :: lzint ! T if it requires a vertical interpolation 84 89 END TYPE FLD 85 86 TYPE, PUBLIC :: MAP_POINTER !: Map from input data file to local domain87 INTEGER, POINTER, DIMENSION(:) :: ptr ! Array of integer pointers to 1D arrays88 LOGICAL :: ll_unstruc ! Unstructured (T) or structured (F) boundary data file89 END TYPE MAP_POINTER90 90 91 91 !$AGRIF_DO_NOT_TREAT … … 129 129 CONTAINS 130 130 131 SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset, jpk_bdy, fvl)131 SUBROUTINE fld_read( kt, kn_fsbc, sd, kit, kt_offset ) 132 132 !!--------------------------------------------------------------------- 133 133 !! *** ROUTINE fld_read *** … … 144 144 INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step) 145 145 TYPE(FLD), INTENT(inout), DIMENSION(:) :: sd ! input field related variables 146 TYPE(MAP_POINTER),INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping indices147 146 INTEGER , INTENT(in ), OPTIONAL :: kit ! subcycle timestep for timesplitting option 148 147 INTEGER , INTENT(in ), OPTIONAL :: kt_offset ! provide fields at time other than "now" … … 150 149 ! ! kt_offset = +1 => fields at "after" time level 151 150 ! ! etc. 152 INTEGER , INTENT(in ), OPTIONAL :: jpk_bdy ! number of vertical levels in the BDY data153 LOGICAL , INTENT(in ), OPTIONAL :: fvl ! number of vertical levels in the BDY data154 151 !! 155 152 INTEGER :: itmp ! local variable … … 166 163 REAL(wp) :: ztintb ! ratio applied to before records when doing time interpolation 167 164 CHARACTER(LEN=1000) :: clfmt ! write format 168 TYPE(MAP_POINTER) :: imap ! global-to-local mapping indices169 165 !!--------------------------------------------------------------------- 170 166 ll_firstcall = kt == nit000 … … 175 171 ENDIF 176 172 IF( PRESENT(kt_offset) ) it_offset = kt_offset 177 178 imap%ptr => NULL()179 173 180 174 ! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar … … 188 182 IF( ll_firstcall ) THEN ! initialization 189 183 DO jf = 1, imf 190 IF( PRESENT(map) ) imap = map(jf) 191 IF( PRESENT(jpk_bdy) ) THEN 192 CALL fld_init( kn_fsbc, sd(jf), imap, jpk_bdy, fvl ) ! read each before field (put them in after as they will be swapped) 193 ELSE 194 CALL fld_init( kn_fsbc, sd(jf), imap ) ! read each before field (put them in after as they will be swapped) 195 ENDIF 184 IF( TRIM(sd(jf)%clrootname) == 'NOT USED' ) CYCLE 185 CALL fld_init( kn_fsbc, sd(jf) ) ! read each before field (put them in after as they will be swapped) 196 186 END DO 197 187 IF( lwp ) CALL wgt_print() ! control print … … 202 192 ! 203 193 DO jf = 1, imf ! --- loop over field --- ! 204 194 195 IF( TRIM(sd(jf)%clrootname) == 'NOT USED' ) CYCLE 196 205 197 IF( isecsbc > sd(jf)%nrec_a(2) .OR. ll_firstcall ) THEN ! read/update the after data? 206 207 IF( PRESENT(map) ) imap = map(jf) ! temporary definition of map208 198 209 199 sd(jf)%nrec_b(:) = sd(jf)%nrec_a(:) ! swap before record informations … … 213 203 CALL fld_rec( kn_fsbc, sd(jf), kt_offset = it_offset, kit = kit ) ! update after record informations 214 204 215 ! if kn_fsbc*rdt is larger than nfreqh (which is kind of odd),205 ! if kn_fsbc*rdt is larger than freqh (which is kind of odd), 216 206 ! it is possible that the before value is no more the good one... we have to re-read it 217 207 ! if before is not the last record of the file currently opened and after is the first record to be read … … 222 212 itmp = sd(jf)%nrec_a(1) ! temporary storage 223 213 sd(jf)%nrec_a(1) = sd(jf)%nreclast ! read the last record of the file currently opened 224 CALL fld_get( sd(jf) , imap )! read after data214 CALL fld_get( sd(jf) ) ! read after data 225 215 sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! re-swap before record field 226 216 sd(jf)%nrec_b(1) = sd(jf)%nrec_a(1) ! update before record informations 227 sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - NINT( sd(jf)% nfreqh * 3600) ! assume freq to be in hours in this case217 sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - NINT( sd(jf)%freqh * 3600. ) ! assume freq to be in hours in this case 228 218 sd(jf)%rotn(1) = sd(jf)%rotn(2) ! update before rotate informations 229 219 sd(jf)%nrec_a(1) = itmp ! move back to after record … … 234 224 IF( sd(jf)%ln_tint ) THEN 235 225 236 ! if kn_fsbc*rdt is larger than nfreqh (which is kind of odd),226 ! if kn_fsbc*rdt is larger than freqh (which is kind of odd), 237 227 ! it is possible that the before value is no more the good one... we have to re-read it 238 228 ! if before record is not just just before the after record... … … 240 230 & .AND. sd(jf)%nrec_b(1) /= sd(jf)%nrec_a(1) - 1 ) THEN 241 231 sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) - 1 ! move back to before record 242 CALL fld_get( sd(jf) , imap )! read after data232 CALL fld_get( sd(jf) ) ! read after data 243 233 sd(jf)%fdta(:,:,:,1) = sd(jf)%fdta(:,:,:,2) ! re-swap before record field 244 234 sd(jf)%nrec_b(1) = sd(jf)%nrec_a(1) ! update before record informations 245 sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - NINT( sd(jf)% nfreqh * 3600) ! assume freq to be in hours in this case235 sd(jf)%nrec_b(2) = sd(jf)%nrec_a(2) - NINT( sd(jf)%freqh * 3600. ) ! assume freq to be in hours in this case 246 236 sd(jf)%rotn(1) = sd(jf)%rotn(2) ! update before rotate informations 247 237 sd(jf)%nrec_a(1) = sd(jf)%nrec_a(1) + 1 ! move back to after record … … 268 258 ! year/month/week/day, next year/month/week/day file must exist 269 259 isecend = nsec_year + nsec1jan000 + (nitend - kt) * NINT(rdt) ! second at the end of the run 270 llstop = isecend > sd(jf)%nrec_a(2) 260 llstop = isecend > sd(jf)%nrec_a(2) ! read more than 1 record of next year 271 261 ! we suppose that the date of next file is next day (should be ok even for weekly files...) 272 262 CALL fld_clopn( sd(jf), nyear + COUNT((/llnxtyr /)) , & … … 277 267 CALL ctl_warn('next year/month/week/day file: '//TRIM(sd(jf)%clname)// & 278 268 & ' not present -> back to current year/month/day') 279 CALL fld_clopn( sd(jf) ) ! back to the current year/month/day269 CALL fld_clopn( sd(jf) ) ! back to the current year/month/day 280 270 sd(jf)%nrec_a(1) = sd(jf)%nreclast ! force to read the last record in the current year file 281 271 ENDIF … … 285 275 286 276 ! read after data 287 IF( PRESENT(jpk_bdy) ) THEN 288 CALL fld_get( sd(jf), imap, jpk_bdy, fvl) 289 ELSE 290 CALL fld_get( sd(jf), imap ) 291 ENDIF 277 CALL fld_get( sd(jf) ) 278 292 279 ENDIF ! read new data? 293 280 END DO ! --- end loop over field --- ! … … 296 283 297 284 DO jf = 1, imf ! --- loop over field --- ! 285 ! 286 IF( TRIM(sd(jf)%clrootname) == 'NOT USED' ) CYCLE 298 287 ! 299 288 IF( sd(jf)%ln_tint ) THEN ! temporal interpolation … … 327 316 328 317 329 SUBROUTINE fld_init( kn_fsbc, sdjf , map , jpk_bdy, fvl)318 SUBROUTINE fld_init( kn_fsbc, sdjf ) 330 319 !!--------------------------------------------------------------------- 331 320 !! *** ROUTINE fld_init *** … … 336 325 INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step) 337 326 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 338 TYPE(MAP_POINTER),INTENT(in) :: map ! global-to-local mapping indices339 INTEGER , INTENT(in), OPTIONAL :: jpk_bdy ! number of vertical levels in the BDY data340 LOGICAL , INTENT(in), OPTIONAL :: fvl ! number of vertical levels in the BDY data341 327 !! 342 328 LOGICAL :: llprevyr ! are we reading previous year file? … … 351 337 CHARACTER(LEN=1000) :: clfmt ! write format 352 338 !!--------------------------------------------------------------------- 339 ! 353 340 llprevyr = .FALSE. 354 341 llprevmth = .FALSE. … … 365 352 ! 366 353 IF( sdjf%nrec_a(1) == 0 ) THEN ! we redefine record sdjf%nrec_a(1) with the last record of previous year file 367 IF ( sdjf%nfreqh== -12 ) THEN ! yearly mean354 IF ( NINT(sdjf%freqh) == -12 ) THEN ! yearly mean 368 355 IF( sdjf%cltype == 'yearly' ) THEN ! yearly file 369 356 sdjf%nrec_a(1) = 1 ! force to read the unique record … … 372 359 CALL ctl_stop( "fld_init: yearly mean file must be in a yearly type of file: "//TRIM(sdjf%clrootname) ) 373 360 ENDIF 374 ELSEIF( sdjf%nfreqh== -1 ) THEN ! monthly mean361 ELSEIF( NINT(sdjf%freqh) == -1 ) THEN ! monthly mean 375 362 IF( sdjf%cltype == 'monthly' ) THEN ! monthly file 376 363 sdjf%nrec_a(1) = 1 ! force to read the unique record … … 381 368 llprevyr = .NOT. sdjf%ln_clim ! use previous year file? 382 369 ENDIF 383 ELSE ! higher frequency mean (in hours)370 ELSE ! higher frequency mean (in hours) 384 371 IF ( sdjf%cltype == 'monthly' ) THEN ! monthly file 385 sdjf%nrec_a(1) = NINT( 24 * nmonth_len(nmonth-1) / sdjf%nfreqh )! last record of previous month372 sdjf%nrec_a(1) = NINT( 24. * REAL(nmonth_len(nmonth-1),wp) / sdjf%freqh )! last record of previous month 386 373 llprevmth = .TRUE. ! use previous month file? 387 374 llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file? 388 375 ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ! weekly file 389 376 llprevweek = .TRUE. ! use previous week file? 390 sdjf%nrec_a(1) = NINT( 24 * 7 / sdjf%nfreqh )! last record of previous week377 sdjf%nrec_a(1) = NINT( 24. * 7. / sdjf%freqh ) ! last record of previous week 391 378 isec_week = NINT(rday) * 7 ! add a shift toward previous week 392 379 ELSEIF( sdjf%cltype == 'daily' ) THEN ! daily file 393 sdjf%nrec_a(1) = NINT( 24 / sdjf%nfreqh ) ! last record of previous day380 sdjf%nrec_a(1) = NINT( 24. / sdjf%freqh ) ! last record of previous day 394 381 llprevday = .TRUE. ! use previous day file? 395 382 llprevmth = llprevday .AND. nday == 1 ! use previous month file? 396 383 llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file? 397 384 ELSE ! yearly file 398 sdjf%nrec_a(1) = NINT( 24 * nyear_len(0) / sdjf%nfreqh )! last record of previous year385 sdjf%nrec_a(1) = NINT( 24. * REAL(nyear_len(0),wp) / sdjf%freqh ) ! last record of previous year 399 386 llprevyr = .NOT. sdjf%ln_clim ! use previous year file? 400 387 ENDIF … … 433 420 ! 434 421 ! read before data in after arrays(as we will swap it later) 435 IF( PRESENT(jpk_bdy) ) THEN 436 CALL fld_get( sdjf, map, jpk_bdy, fvl ) 437 ELSE 438 CALL fld_get( sdjf, map ) 439 ENDIF 422 CALL fld_get( sdjf ) 440 423 ! 441 424 clfmt = "(' fld_init : time-interpolation for ', a, ' read previous record = ', i6, ' at time = ', f7.2, ' days')" … … 456 439 !! if sdjf%ln_tint = .FALSE. 457 440 !! nrec_a(1): record number 458 !! nrec_b(2) and nrec_a(2): time of the beginning and end of the record (for print only)441 !! nrec_b(2) and nrec_a(2): time of the beginning and end of the record 459 442 !!---------------------------------------------------------------------- 460 443 INTEGER , INTENT(in ) :: kn_fsbc ! sbc computation period (in time step) … … 484 467 ELSE ; it_offset = 0 485 468 ENDIF 486 IF( PRESENT(kt_offset) ) it_offset = kt_offset469 IF( PRESENT(kt_offset) ) it_offset = kt_offset 487 470 IF( PRESENT(kit) ) THEN ; it_offset = ( kit + it_offset ) * NINT( rdt/REAL(nn_baro,wp) ) 488 471 ELSE ; it_offset = it_offset * NINT( rdt ) 489 472 ENDIF 490 473 ! 491 ! ! =========== !492 IF ( sdjf%nfreqh== -12 ) THEN ! yearly mean493 ! ! =========== !494 ! 495 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record474 ! ! =========== ! 475 IF ( NINT(sdjf%freqh) == -12 ) THEN ! yearly mean 476 ! ! =========== ! 477 ! 478 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record 496 479 ! 497 480 ! INT( ztmp ) … … 505 488 ! forcing record : 1 506 489 ! 507 ztmp = REAL( nsec_year, wp ) / ( REAL( nyear_len(1), wp ) * rday ) + 0.5 &508 &+ REAL( it_offset, wp ) / ( REAL( nyear_len(1), wp ) * rday )490 ztmp = REAL( nsec_year, wp ) / ( REAL( nyear_len(1), wp ) * rday ) + 0.5 & 491 & + REAL( it_offset, wp ) / ( REAL( nyear_len(1), wp ) * rday ) 509 492 sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) 510 493 ! swap at the middle of the year … … 514 497 & INT(ztmp) * INT(rday) * nyear_len(1) + INT(ztmp) * NINT( 0.5 * rday) * nyear_len(2) 515 498 ENDIF 516 ELSE ! no time interpolation499 ELSE ! no time interpolation 517 500 sdjf%nrec_a(1) = 1 518 501 sdjf%nrec_a(2) = NINT(rday) * nyear_len(1) + nsec1jan000 ! swap at the end of the year … … 520 503 ENDIF 521 504 ! 522 ! ! ============ !523 ELSEIF( sdjf%nfreqh== -1 ) THEN ! monthly mean !524 ! ! ============ !525 ! 526 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record505 ! ! ============ ! 506 ELSEIF( NINT(sdjf%freqh) == -1 ) THEN ! monthly mean ! 507 ! ! ============ ! 508 ! 509 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record 527 510 ! 528 511 ! INT( ztmp ) … … 536 519 ! forcing record : nmonth 537 520 ! 538 ztmp = REAL( nsec_month, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday ) + 0.5 &539 & + REAL(it_offset, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday )521 ztmp = REAL( nsec_month, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday ) + 0.5 & 522 & + REAL( it_offset, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday ) 540 523 imth = nmonth + INT( ztmp ) - COUNT((/llbefore/)) 541 524 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/)) … … 551 534 ENDIF 552 535 ! 553 ! ! ================================ !554 ELSE ! higher frequency mean (in hours)555 ! ! ================================ !556 ! 557 ifreq_sec = NINT( sdjf% nfreqh * 3600) ! frequency mean (in seconds)536 ! ! ================================ ! 537 ELSE ! higher frequency mean (in hours) 538 ! ! ================================ ! 539 ! 540 ifreq_sec = NINT( sdjf%freqh * 3600. ) ! frequency mean (in seconds) 558 541 IF( sdjf%cltype(1:4) == 'week' ) isec_week = ksec_week( sdjf%cltype(6:8) ) ! since the first day of the current week 559 542 ! number of second since the beginning of the file … … 565 548 ztmp = ztmp + 0.5 * REAL(kn_fsbc - 1, wp) * rdt + REAL( it_offset, wp ) ! centrered in the middle of sbc time step 566 549 ztmp = ztmp + 0.01 * rdt ! avoid truncation error 567 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record550 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record 568 551 ! 569 552 ! INT( ztmp/ifreq_sec + 0.5 ) … … 579 562 ! 580 563 ztmp= ztmp / REAL(ifreq_sec, wp) + 0.5 581 ELSE ! no time interpolation564 ELSE ! no time interpolation 582 565 ! 583 566 ! INT( ztmp/ifreq_sec ) … … 610 593 ENDIF 611 594 ! 595 IF( .NOT. sdjf%ln_tint ) sdjf%nrec_a(2) = sdjf%nrec_a(2) - 1 ! last second belongs to bext record : *----( 596 ! 612 597 END SUBROUTINE fld_rec 613 598 614 599 615 SUBROUTINE fld_get( sdjf , map, jpk_bdy, fvl)600 SUBROUTINE fld_get( sdjf ) 616 601 !!--------------------------------------------------------------------- 617 602 !! *** ROUTINE fld_get *** … … 620 605 !!---------------------------------------------------------------------- 621 606 TYPE(FLD) , INTENT(inout) :: sdjf ! input field related variables 622 TYPE(MAP_POINTER), INTENT(in ) :: map ! global-to-local mapping indices623 INTEGER , INTENT(in), OPTIONAL :: jpk_bdy ! number of vertical levels in the bdy data624 LOGICAL , INTENT(in), OPTIONAL :: fvl ! number of vertical levels in the bdy data625 607 ! 626 608 INTEGER :: ipk ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) … … 634 616 ipk = SIZE( sdjf%fnow, 3 ) 635 617 ! 636 IF( ASSOCIATED(map%ptr) ) THEN 637 IF( PRESENT(jpk_bdy) ) THEN 638 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), & 639 sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl ) 640 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), & 641 sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl ) 642 ENDIF 643 ELSE 644 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map ) 645 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), map ) 646 ENDIF 647 ENDIF 618 IF( ASSOCIATED(sdjf%imap) ) THEN 619 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), & 620 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint ) 621 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), & 622 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint ) 623 ENDIF 648 624 ELSE IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN 649 625 CALL wgt_list( sdjf, iw ) … … 700 676 END SUBROUTINE fld_get 701 677 702 SUBROUTINE fld_map( num, clvar, dta, nrec, map, igrd, ibdy, jpk_bdy, fvl ) 678 679 SUBROUTINE fld_map( knum, cdvar, pdta, krec, kmap, kgrd, kbdy, ldtotvel, ldzint ) 703 680 !!--------------------------------------------------------------------- 704 681 !! *** ROUTINE fld_map *** … … 707 684 !! using a general mapping (for open boundaries) 708 685 !!---------------------------------------------------------------------- 709 710 USE bdy_oce, ONLY: ln_bdy, idx_bdy, dta_global, dta_global_z, dta_global_dz, dta_global2, dta_global2_z, dta_global2_dz ! workspace to read in global data arrays 711 712 INTEGER , INTENT(in ) :: num ! stream number 713 CHARACTER(LEN=*) , INTENT(in ) :: clvar ! variable name 714 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional) 715 INTEGER , INTENT(in ) :: nrec ! record number to read (ie time slice) 716 TYPE(MAP_POINTER) , INTENT(in ) :: map ! global-to-local mapping indices 717 INTEGER , INTENT(in), OPTIONAL :: igrd, ibdy, jpk_bdy ! grid type, set number and number of vertical levels in the bdy data 718 LOGICAL , INTENT(in), OPTIONAL :: fvl ! grid type, set number and number of vertical levels in the bdy data 719 INTEGER :: jpkm1_bdy! number of vertical levels in the bdy data minus 1 720 !! 721 INTEGER :: ipi ! length of boundary data on local process 722 INTEGER :: ipj ! length of dummy dimension ( = 1 ) 723 INTEGER :: ipk ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk ) 724 INTEGER :: ilendta ! length of data in file 725 INTEGER :: idvar ! variable ID 726 INTEGER :: ib, ik, ji, jj ! loop counters 727 INTEGER :: ierr 728 REAL(wp) :: fv ! fillvalue 729 REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read ! work space for global data 730 REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read_z ! work space for global data 731 REAL(wp), POINTER, DIMENSION(:,:,:) :: dta_read_dz ! work space for global data 732 !!--------------------------------------------------------------------- 733 ! 734 ipi = SIZE( dta, 1 ) 735 ipj = 1 736 ipk = SIZE( dta, 3 ) 737 ! 738 idvar = iom_varid( num, clvar ) 739 ilendta = iom_file(num)%dimsz(1,idvar) 740 741 IF ( ln_bdy ) THEN 742 ipj = iom_file(num)%dimsz(2,idvar) 743 IF( map%ll_unstruc) THEN ! unstructured open boundary data file 744 dta_read => dta_global 745 IF( PRESENT(jpk_bdy) ) THEN 746 IF( jpk_bdy>0 ) THEN 747 dta_read_z => dta_global_z 748 dta_read_dz => dta_global_dz 749 jpkm1_bdy = jpk_bdy-1 750 ENDIF 751 ENDIF 752 ELSE ! structured open boundary file 753 dta_read => dta_global2 754 IF( PRESENT(jpk_bdy) ) THEN 755 IF( jpk_bdy>0 ) THEN 756 dta_read_z => dta_global2_z 757 dta_read_dz => dta_global2_dz 758 jpkm1_bdy = jpk_bdy-1 759 ENDIF 760 ENDIF 761 ENDIF 762 ENDIF 763 764 IF(lwp) WRITE(numout,*) 'Dim size for ', TRIM(clvar),' is ', ilendta 765 IF(lwp) WRITE(numout,*) 'Number of levels for ',TRIM(clvar),' is ', ipk 766 ! 767 SELECT CASE( ipk ) 768 CASE(1) ; 769 CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1 ), nrec ) 770 IF ( map%ll_unstruc) THEN ! unstructured open boundary data file 771 DO ib = 1, ipi 772 DO ik = 1, ipk 773 dta(ib,1,ik) = dta_read(map%ptr(ib),1,ik) 774 END DO 775 END DO 776 ELSE ! we assume that this is a structured open boundary file 777 DO ib = 1, ipi 778 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 779 ji=map%ptr(ib)-(jj-1)*ilendta 780 DO ik = 1, ipk 781 dta(ib,1,ik) = dta_read(ji,jj,ik) 782 END DO 783 END DO 784 ENDIF 686 INTEGER , INTENT(in ) :: knum ! stream number 687 CHARACTER(LEN=*) , INTENT(in ) :: cdvar ! variable name 688 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pdta ! bdy output field on model grid 689 INTEGER , INTENT(in ) :: krec ! record number to read (ie time slice) 690 INTEGER , DIMENSION(:) , INTENT(in ) :: kmap ! global-to-local bdy mapping indices 691 ! optional variables used for vertical interpolation: 692 INTEGER, OPTIONAL , INTENT(in ) :: kgrd ! grid type (t, u, v) 693 INTEGER, OPTIONAL , INTENT(in ) :: kbdy ! bdy number 694 LOGICAL, OPTIONAL , INTENT(in ) :: ldtotvel ! true if total ( = barotrop + barocline) velocity 695 LOGICAL, OPTIONAL , INTENT(in ) :: ldzint ! true if 3D variable requires a vertical interpolation 696 !! 697 INTEGER :: ipi ! length of boundary data on local process 698 INTEGER :: ipj ! length of dummy dimension ( = 1 ) 699 INTEGER :: ipk ! number of vertical levels of pdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 700 INTEGER :: ipkb ! number of vertical levels in boundary data file 701 INTEGER :: idvar ! variable ID 702 INTEGER :: indims ! number of dimensions of the variable 703 INTEGER, DIMENSION(4) :: idimsz ! size of variable dimensions 704 REAL(wp) :: zfv ! fillvalue 705 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zz_read ! work space for global boundary data 706 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read ! work space local data requiring vertical interpolation 707 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read_z ! work space local data requiring vertical interpolation 708 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read_dz ! work space local data requiring vertical interpolation 709 CHARACTER(LEN=1),DIMENSION(3) :: clgrid 710 LOGICAL :: lluld ! is the variable using the unlimited dimension 711 LOGICAL :: llzint ! local value of ldzint 712 !!--------------------------------------------------------------------- 713 ! 714 clgrid = (/'t','u','v'/) 715 ! 716 ipi = SIZE( pdta, 1 ) 717 ipj = SIZE( pdta, 2 ) ! must be equal to 1 718 ipk = SIZE( pdta, 3 ) 719 ! 720 llzint = .FALSE. 721 IF( PRESENT(ldzint) ) llzint = ldzint 722 ! 723 idvar = iom_varid( knum, cdvar, kndims = indims, kdimsz = idimsz, lduld = lluld ) 724 IF( indims == 4 .OR. ( indims == 3 .AND. .NOT. lluld ) ) THEN ; ipkb = idimsz(3) ! xy(zl)t or xy(zl) 725 ELSE ; ipkb = 1 ! xy or xyt 726 ENDIF 727 ! 728 ALLOCATE( zz_read( idimsz(1), idimsz(2), ipkb ) ) ! ++++++++ !!! this can be very big... 729 ! 730 IF( ipk == 1 ) THEN 731 732 IF( ipkb /= 1 ) CALL ctl_stop( 'fld_map : we must have ipkb = 1 to read surface data' ) 733 CALL iom_get ( knum, jpdom_unknown, cdvar, zz_read(:,:,1), krec ) ! call iom_get with a 2D file 734 CALL fld_map_core( zz_read, kmap, pdta ) 785 735 786 736 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 787 737 ! Do we include something here to adjust barotropic velocities ! 788 738 ! in case of a depth difference between bdy files and ! 789 ! bathymetry in the case ln_ full_vel = .false. and jpk_bdy>0?!739 ! bathymetry in the case ln_totvel = .false. and ipkb>0? ! 790 740 ! [as the enveloping and parital cells could change H] ! 791 741 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 792 742 793 CASE DEFAULT ; 794 795 IF( PRESENT(jpk_bdy) .AND. jpk_bdy>0 ) THEN ! boundary data not on model grid: vertical interpolation 796 CALL iom_getatt(num, '_FillValue', fv, cdvar=clvar ) 797 dta_read(:,:,:) = -ABS(fv) 798 dta_read_z(:,:,:) = 0._wp 799 dta_read_dz(:,:,:) = 0._wp 800 CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:jpk_bdy), nrec ) 801 SELECT CASE( igrd ) 802 CASE(1) 803 CALL iom_get ( num, jpdom_unknown, 'gdept', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 804 CALL iom_get ( num, jpdom_unknown, 'e3t', dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 805 CASE(2) 806 CALL iom_get ( num, jpdom_unknown, 'gdepu', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 807 CALL iom_get ( num, jpdom_unknown, 'e3u', dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 808 CASE(3) 809 CALL iom_get ( num, jpdom_unknown, 'gdepv', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 810 CALL iom_get ( num, jpdom_unknown, 'e3v', dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 811 END SELECT 812 813 IF ( ln_bdy ) & 814 CALL fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta) 815 816 ELSE ! boundary data assumed to be on model grid 817 818 CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec ) 819 IF ( map%ll_unstruc) THEN ! unstructured open boundary data file 820 DO ib = 1, ipi 821 DO ik = 1, ipk 822 dta(ib,1,ik) = dta_read(map%ptr(ib),1,ik) 823 END DO 743 ELSE 744 ! 745 CALL iom_get ( knum, jpdom_unknown, cdvar, zz_read(:,:,:), krec ) ! call iom_get with a 3D file 746 ! 747 IF( ipkb /= ipk .OR. llzint ) THEN ! boundary data not on model vertical grid : vertical interpolation 748 ! 749 IF( ipk == jpk .AND. iom_varid(knum,'gdep'//clgrid(kgrd)) /= -1 .AND. iom_varid(knum,'e3'//clgrid(kgrd)) /= -1 ) THEN 750 751 ALLOCATE( zdta_read(ipi,ipj,ipkb), zdta_read_z(ipi,ipj,ipkb), zdta_read_dz(ipi,ipj,ipkb) ) 752 753 CALL fld_map_core( zz_read, kmap, zdta_read ) 754 CALL iom_get ( knum, jpdom_unknown, 'gdep'//clgrid(kgrd), zz_read ) ! read only once? Potential temporal evolution? 755 CALL fld_map_core( zz_read, kmap, zdta_read_z ) 756 CALL iom_get ( knum, jpdom_unknown, 'e3'//clgrid(kgrd), zz_read ) ! read only once? Potential temporal evolution? 757 CALL fld_map_core( zz_read, kmap, zdta_read_dz ) 758 759 CALL iom_getatt(knum, '_FillValue', zfv, cdvar=cdvar ) 760 CALL fld_bdy_interp(zdta_read, zdta_read_z, zdta_read_dz, pdta, kgrd, kbdy, zfv, ldtotvel) 761 DEALLOCATE( zdta_read, zdta_read_z, zdta_read_dz ) 762 763 ELSE 764 IF( ipk /= jpk ) CALL ctl_stop( 'fld_map : this should be an impossible case...' ) 765 WRITE(ctmp1,*) 'fld_map : vertical interpolation for bdy variable '//TRIM(cdvar)//' requires ' 766 IF( iom_varid(knum, 'gdep'//clgrid(kgrd)) == -1 ) CALL ctl_stop( ctmp1//'gdep'//clgrid(kgrd)//' variable' ) 767 IF( iom_varid(knum, 'e3'//clgrid(kgrd)) == -1 ) CALL ctl_stop( ctmp1// 'e3'//clgrid(kgrd)//' variable' ) 768 769 ENDIF 770 ! 771 ELSE ! bdy data assumed to be the same levels as bdy variables 772 ! 773 CALL fld_map_core( zz_read, kmap, pdta ) 774 ! 775 ENDIF ! ipkb /= ipk 776 ENDIF ! ipk == 1 777 778 DEALLOCATE( zz_read ) 779 780 END SUBROUTINE fld_map 781 782 783 SUBROUTINE fld_map_core( pdta_read, kmap, pdta_bdy ) 784 !!--------------------------------------------------------------------- 785 !! *** ROUTINE fld_map_core *** 786 !! 787 !! ** Purpose : inner core of fld_map 788 !!---------------------------------------------------------------------- 789 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pdta_read ! global boundary data 790 INTEGER, DIMENSION(: ), INTENT(in ) :: kmap ! global-to-local bdy mapping indices 791 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pdta_bdy ! bdy output field on model grid 792 !! 793 INTEGER, DIMENSION(3) :: idim_read, idim_bdy ! arrays dimensions 794 INTEGER :: ji, jj, jk, jb ! loop counters 795 INTEGER :: im1 796 !!--------------------------------------------------------------------- 797 ! 798 idim_read = SHAPE( pdta_read ) 799 idim_bdy = SHAPE( pdta_bdy ) 800 ! 801 ! in all cases: idim_bdy(2) == 1 .AND. idim_read(1) * idim_read(2) == idim_bdy(1) 802 ! structured BDY with rimwidth > 1 : idim_read(2) == rimwidth /= 1 803 ! structured BDY with rimwidth == 1 or unstructured BDY: idim_read(2) == 1 804 ! 805 IF( idim_read(2) > 1 ) THEN ! structured BDY with rimwidth > 1 806 DO jk = 1, idim_bdy(3) 807 DO jb = 1, idim_bdy(1) 808 im1 = kmap(jb) - 1 809 jj = im1 / idim_read(1) + 1 810 ji = MOD( im1, idim_read(1) ) + 1 811 pdta_bdy(jb,1,jk) = pdta_read(ji,jj,jk) 824 812 END DO 825 ELSE ! we assume that this is a structured open boundary file 826 DO ib = 1, ipi 827 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 828 ji=map%ptr(ib)-(jj-1)*ilendta 829 DO ik = 1, ipk 830 dta(ib,1,ik) = dta_read(ji,jj,ik) 831 END DO 813 END DO 814 ELSE 815 DO jk = 1, idim_bdy(3) 816 DO jb = 1, idim_bdy(1) ! horizontal remap of bdy data on the local bdy 817 pdta_bdy(jb,1,jk) = pdta_read(kmap(jb),1,jk) 832 818 END DO 833 ENDIF 834 ENDIF ! PRESENT(jpk_bdy) 835 END SELECT 836 837 END SUBROUTINE fld_map 819 END DO 820 ENDIF 821 822 END SUBROUTINE fld_map_core 838 823 839 SUBROUTINE fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta)840 824 825 SUBROUTINE fld_bdy_interp(pdta_read, pdta_read_z, pdta_read_dz, pdta, kgrd, kbdy, pfv, ldtotvel) 841 826 !!--------------------------------------------------------------------- 842 827 !! *** ROUTINE fld_bdy_interp *** … … 847 832 USE bdy_oce, ONLY: idx_bdy ! indexing for map <-> ij transformation 848 833 849 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read ! work space for global data 850 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read_z ! work space for global data 851 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read_dz ! work space for global data 852 REAL(wp) , INTENT(in) :: fv ! fillvalue and alternative -ABS(fv) 853 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional) 854 TYPE(MAP_POINTER) , INTENT(in ) :: map ! global-to-local mapping indices 855 LOGICAL , INTENT(in), OPTIONAL :: fvl ! grid type, set number and number of vertical levels in the bdy data 856 INTEGER , INTENT(in) :: igrd, ibdy, jpk_bdy ! number of levels in bdy data 857 INTEGER , INTENT(in) :: ilendta ! length of data in file 858 !! 859 INTEGER :: ipi ! length of boundary data on local process 860 INTEGER :: ipj ! length of dummy dimension ( = 1 ) 861 INTEGER :: ipk ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk ) 862 INTEGER :: jpkm1_bdy ! number of levels in bdy data minus 1 863 INTEGER :: ib, ik, ikk ! loop counters 864 INTEGER :: ji, jj, zij, zjj ! temporary indices 865 REAL(wp) :: zl, zi, zh ! tmp variable for current depth and interpolation factor 866 REAL(wp) :: fv_alt, ztrans, ztrans_new ! fillvalue and alternative -ABS(fv) 867 CHARACTER (LEN=10) :: ibstr 834 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pdta_read ! data read in bdy file 835 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pdta_read_z ! depth of the data read in bdy file 836 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pdta_read_dz ! thickness of the levels in bdy file 837 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pdta ! output field on model grid (2 dimensional) 838 REAL(wp) , INTENT(in ) :: pfv ! fillvalue of the data read in bdy file 839 LOGICAL , INTENT(in ) :: ldtotvel ! true if toal ( = barotrop + barocline) velocity 840 INTEGER , INTENT(in ) :: kgrd ! grid type (t, u, v) 841 INTEGER , INTENT(in ) :: kbdy ! bdy number 842 !! 843 INTEGER :: ipi ! length of boundary data on local process 844 INTEGER :: ipkb ! number of vertical levels in boundary data file 845 INTEGER :: jb, ji, jj, jk, jkb ! loop counters 846 REAL(wp) :: zcoef 847 REAL(wp) :: zl, zi, zh ! tmp variable for current depth and interpolation factor 848 REAL(wp) :: zfv_alt, ztrans, ztrans_new ! fillvalue and alternative -ABS(pfv) 849 REAL(wp), DIMENSION(jpk) :: zdepth, zdhalf ! level and half-level depth 868 850 !!--------------------------------------------------------------------- 869 851 870 871 ipi = SIZE( dta, 1 ) 872 ipj = SIZE( dta_read, 2 ) 873 ipk = SIZE( dta, 3 ) 874 jpkm1_bdy = jpk_bdy-1 852 ipi = SIZE( pdta, 1 ) 853 ipkb = SIZE( pdta_read, 3 ) 875 854 876 fv_alt = -ABS(fv) ! set _FillValue < 0 as we make use of MAXVAL and MAXLOC later 877 DO ib = 1, ipi 878 zij = idx_bdy(ibdy)%nbi(ib,igrd) 879 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 880 IF(narea==2) WRITE(*,*) 'MAPI', ib, igrd, map%ptr(ib), narea-1, zij, zjj 881 ENDDO 882 ! 883 IF ( map%ll_unstruc ) THEN ! unstructured open boundary data file 884 885 DO ib = 1, ipi 886 DO ik = 1, jpk_bdy 887 IF( ( dta_read(map%ptr(ib),1,ik) == fv ) ) THEN 888 dta_read_z(map%ptr(ib),1,ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data 889 dta_read_dz(map%ptr(ib),1,ik) = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct 855 zfv_alt = -ABS(pfv) ! set _FillValue < 0 as we make use of MAXVAL and MAXLOC later 856 ! 857 WHERE( pdta_read == pfv ) 858 pdta_read_z = zfv_alt ! safety: put fillvalue into external depth field so consistent with data 859 pdta_read_dz = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct 860 ENDWHERE 861 862 DO jb = 1, ipi 863 ji = idx_bdy(kbdy)%nbi(jb,kgrd) 864 jj = idx_bdy(kbdy)%nbj(jb,kgrd) 865 zh = SUM(pdta_read_dz(jb,1,:) ) 866 ! 867 ! Warnings to flag differences in the input and model topgraphy - is this useful/necessary? 868 SELECT CASE( kgrd ) 869 CASE(1) 870 IF( ABS( (zh - ht_n(ji,jj)) / ht_n(ji,jj)) * tmask(ji,jj,1) > 0.01_wp ) THEN 871 WRITE(ctmp1,"(I10.10)") jb 872 CALL ctl_warn('fld_bdy_interp: T depths differ between grids at BDY point '//TRIM(ctmp1)//' by more than 1%') 873 ! IF(lwp) WRITE(numout,*) 'DEPTHT', zh, sum(e3t_n(ji,jj,:), mask=tmask(ji,jj,:)==1), ht_n(ji,jj), jb, jb, ji, jj 874 ENDIF 875 CASE(2) 876 IF( ABS( (zh - hu_n(ji,jj)) * r1_hu_n(ji,jj)) * umask(ji,jj,1) > 0.01_wp ) THEN 877 WRITE(ctmp1,"(I10.10)") jb 878 CALL ctl_warn('fld_bdy_interp: U depths differ between grids at BDY point '//TRIM(ctmp1)//' by more than 1%') 879 ! IF(lwp) WRITE(numout,*) 'DEPTHU', zh, SUM(e3u_n(ji,jj,:), mask=umask(ji,jj,:)==1), SUM(umask(ji,jj,:)), & 880 ! & hu_n(ji,jj), jb, jb, ji, jj, narea-1, pdta_read(jb,1,:) 881 ENDIF 882 CASE(3) 883 IF( ABS( (zh - hv_n(ji,jj)) * r1_hv_n(ji,jj)) * vmask(ji,jj,1) > 0.01_wp ) THEN 884 WRITE(ctmp1,"(I10.10)") jb 885 CALL ctl_warn('fld_bdy_interp: V depths differ between grids at BDY point '//TRIM(ctmp1)//' by more than 1%') 886 ENDIF 887 END SELECT 888 ! 889 SELECT CASE( kgrd ) 890 CASE(1) 891 ! depth of T points: 892 zdepth(:) = gdept_n(ji,jj,:) 893 CASE(2) 894 ! depth of U points: we must not use gdept_n as we don't want to do a communication 895 ! --> copy what is done for gdept_n in domvvl... 896 zdhalf(1) = 0.0_wp 897 zdepth(1) = 0.5_wp * e3uw_n(ji,jj,1) 898 DO jk = 2, jpk ! vertical sum 899 ! zcoef = umask - wumask ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 900 ! ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 901 ! ! 0.5 where jk = mikt 902 !!gm ??????? BUG ? gdept_n as well as gde3w_n does not include the thickness of ISF ?? 903 zcoef = ( umask(ji,jj,jk) - wumask(ji,jj,jk) ) 904 zdhalf(jk) = zdhalf(jk-1) + e3u_n(ji,jj,jk-1) 905 zdepth(jk) = zcoef * ( zdhalf(jk ) + 0.5 * e3uw_n(ji,jj,jk)) & 906 & + (1-zcoef) * ( zdepth(jk-1) + e3uw_n(ji,jj,jk)) 907 END DO 908 CASE(3) 909 ! depth of V points: we must not use gdept_n as we don't want to do a communication 910 ! --> copy what is done for gdept_n in domvvl... 911 zdhalf(1) = 0.0_wp 912 zdepth(1) = 0.5_wp * e3vw_n(ji,jj,1) 913 DO jk = 2, jpk ! vertical sum 914 ! zcoef = vmask - wvmask ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 915 ! ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 916 ! ! 0.5 where jk = mikt 917 !!gm ??????? BUG ? gdept_n as well as gde3w_n does not include the thickness of ISF ?? 918 zcoef = ( vmask(ji,jj,jk) - wvmask(ji,jj,jk) ) 919 zdhalf(jk) = zdhalf(jk-1) + e3v_n(ji,jj,jk-1) 920 zdepth(jk) = zcoef * ( zdhalf(jk ) + 0.5 * e3vw_n(ji,jj,jk)) & 921 & + (1-zcoef) * ( zdepth(jk-1) + e3vw_n(ji,jj,jk)) 922 END DO 923 END SELECT 924 ! 925 DO jk = 1, jpk 926 IF( zdepth(jk) < pdta_read_z(jb,1, 1) ) THEN ! above the first level of external data 927 pdta(jb,1,jk) = pdta_read(jb,1,1) 928 ELSEIF( zdepth(jk) > pdta_read_z(jb,1,ipkb) ) THEN ! below the last level of external data 929 pdta(jb,1,jk) = pdta_read(jb,1,MAXLOC(pdta_read_z(jb,1,:),1)) 930 ELSE ! inbetween: vertical interpolation between jkb & jkb+1 931 DO jkb = 1, ipkb-1 ! when gdept_n(jkb) < zdepth(jk) < gdept_n(jkb+1) 932 IF( ( ( zdepth(jk) - pdta_read_z(jb,1,jkb) ) * ( zdepth(jk) - pdta_read_z(jb,1,jkb+1) ) <= 0._wp ) & 933 & .AND. ( pdta_read_z(jb,1,jkb+1) /= zfv_alt) ) THEN ! linear interpolation between 2 levels 934 zi = ( zdepth(jk) - pdta_read_z(jb,1,jkb) ) / ( pdta_read_z(jb,1,jkb+1) - pdta_read_z(jb,1,jkb) ) 935 pdta(jb,1,jk) = pdta_read(jb,1,jkb) + ( pdta_read (jb,1,jkb+1) - pdta_read (jb,1,jkb) ) * zi 936 ENDIF 937 END DO 938 ENDIF 939 END DO ! jpk 940 ! 941 END DO ! ipi 942 943 IF(kgrd == 2) THEN ! do we need to adjust the transport term? 944 DO jb = 1, ipi 945 ji = idx_bdy(kbdy)%nbi(jb,kgrd) 946 jj = idx_bdy(kbdy)%nbj(jb,kgrd) 947 zh = SUM(pdta_read_dz(jb,1,:) ) 948 ztrans = 0._wp 949 ztrans_new = 0._wp 950 DO jkb = 1, ipkb ! calculate transport on input grid 951 ztrans = ztrans + pdta_read(jb,1,jkb) * pdta_read_dz(jb, 1,jkb) 952 ENDDO 953 DO jk = 1, jpk ! calculate transport on model grid 954 ztrans_new = ztrans_new + pdta(jb,1,jk ) * e3u_n(ji,jj,jk ) * umask(ji,jj,jk) 955 ENDDO 956 DO jk = 1, jpk ! make transport correction 957 IF(ldtotvel) THEN ! bdy data are total velocity so adjust bt transport term to match input data 958 pdta(jb,1,jk) = ( pdta(jb,1,jk) + ( ztrans - ztrans_new ) * r1_hu_n(ji,jj) ) * umask(ji,jj,jk) 959 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 960 pdta(jb,1,jk) = pdta(jb,1,jk) + ( 0._wp - ztrans_new ) * r1_hu_n(ji,jj) * umask(ji,jj,jk) 890 961 ENDIF 891 962 ENDDO 892 ENDDO 893 894 DO ib = 1, ipi 895 zij = idx_bdy(ibdy)%nbi(ib,igrd) 896 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 897 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 898 ! Warnings to flag differences in the input and model topgraphy - is this useful/necessary? 899 SELECT CASE( igrd ) 900 CASE(1) 901 IF( ABS( (zh - ht_n(zij,zjj)) / ht_n(zij,zjj)) * tmask(zij,zjj,1) > 0.01_wp ) THEN 902 WRITE(ibstr,"(I10.10)") map%ptr(ib) 903 CALL ctl_warn('fld_bdy_interp: T depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 904 ! IF(lwp) WRITE(*,*) 'DEPTHT', zh, sum(e3t_n(zij,zjj,:), mask=tmask(zij,zjj,:)==1), ht_n(zij,zjj), map%ptr(ib), ib, zij, zjj 905 ENDIF 906 CASE(2) 907 IF( ABS( (zh - hu_n(zij,zjj)) * r1_hu_n(zij,zjj)) * umask(zij,zjj,1) > 0.01_wp ) THEN 908 WRITE(ibstr,"(I10.10)") map%ptr(ib) 909 CALL ctl_warn('fld_bdy_interp: U depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 910 IF(lwp) WRITE(*,*) 'DEPTHU', zh, sum(e3u_n(zij,zjj,:), mask=umask(zij,zjj,:)==1), sum(umask(zij,zjj,:)), & 911 & hu_n(zij,zjj), map%ptr(ib), ib, zij, zjj, narea-1 , & 912 & dta_read(map%ptr(ib),1,:) 913 ENDIF 914 CASE(3) 915 IF( ABS( (zh - hv_n(zij,zjj)) * r1_hv_n(zij,zjj)) * vmask(zij,zjj,1) > 0.01_wp ) THEN 916 WRITE(ibstr,"(I10.10)") map%ptr(ib) 917 CALL ctl_warn('fld_bdy_interp: V depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 918 ENDIF 919 END SELECT 920 DO ik = 1, ipk 921 SELECT CASE( igrd ) 922 CASE(1) 923 zl = gdept_n(zij,zjj,ik) ! if using in step could use fsdept instead of gdept_n? 924 CASE(2) 925 IF(ln_sco) THEN 926 zl = ( gdept_n(zij,zjj,ik) + gdept_n(zij+1,zjj,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_n? 927 ELSE 928 zl = MIN( gdept_n(zij,zjj,ik), gdept_n(zij+1,zjj,ik) ) 929 ENDIF 930 CASE(3) 931 IF(ln_sco) THEN 932 zl = ( gdept_n(zij,zjj,ik) + gdept_n(zij,zjj+1,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_n? 933 ELSE 934 zl = MIN( gdept_n(zij,zjj,ik), gdept_n(zij,zjj+1,ik) ) 935 ENDIF 936 END SELECT 937 IF( zl < dta_read_z(map%ptr(ib),1,1) ) THEN ! above the first level of external data 938 dta(ib,1,ik) = dta_read(map%ptr(ib),1,1) 939 ELSEIF( zl > MAXVAL(dta_read_z(map%ptr(ib),1,:),1) ) THEN ! below the last level of external data 940 dta(ib,1,ik) = dta_read(map%ptr(ib),1,MAXLOC(dta_read_z(map%ptr(ib),1,:),1)) 941 ELSE ! inbetween : vertical interpolation between ikk & ikk+1 942 DO ikk = 1, jpkm1_bdy ! when gdept_n(ikk) < zl < gdept_n(ikk+1) 943 IF( ( (zl-dta_read_z(map%ptr(ib),1,ikk)) * (zl-dta_read_z(map%ptr(ib),1,ikk+1)) <= 0._wp) & 944 & .AND. (dta_read_z(map%ptr(ib),1,ikk+1) /= fv_alt)) THEN 945 zi = ( zl - dta_read_z(map%ptr(ib),1,ikk) ) / & 946 & ( dta_read_z(map%ptr(ib),1,ikk+1) - dta_read_z(map%ptr(ib),1,ikk) ) 947 dta(ib,1,ik) = dta_read(map%ptr(ib),1,ikk) + & 948 & ( dta_read(map%ptr(ib),1,ikk+1) - dta_read(map%ptr(ib),1,ikk) ) * zi 949 ENDIF 950 END DO 951 ENDIF 952 END DO 953 END DO 954 955 IF(igrd == 2) THEN ! do we need to adjust the transport term? 956 DO ib = 1, ipi 957 zij = idx_bdy(ibdy)%nbi(ib,igrd) 958 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 959 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 960 ztrans = 0._wp 961 ztrans_new = 0._wp 962 DO ik = 1, jpk_bdy ! calculate transport on input grid 963 ztrans = ztrans + dta_read(map%ptr(ib),1,ik) * dta_read_dz(map%ptr(ib),1,ik) 964 ENDDO 965 DO ik = 1, ipk ! calculate transport on model grid 966 ztrans_new = ztrans_new + dta(ib,1,ik) * e3u_n(zij,zjj,ik) * umask(zij,zjj,ik) 967 ENDDO 968 DO ik = 1, ipk ! make transport correction 969 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 970 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hu_n(zij,zjj) ) * umask(zij,zjj,ik) 971 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 972 IF( ABS(ztrans * r1_hu_n(zij,zjj)) > 0.01_wp ) & 973 & CALL ctl_warn('fld_bdy_interp: barotropic component of > 0.01 ms-1 found in baroclinic velocities at') 974 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hu_n(zij,zjj) * umask(zij,zjj,ik) 975 ENDIF 976 ENDDO 963 ENDDO 964 ENDIF 965 966 IF(kgrd == 3) THEN ! do we need to adjust the transport term? 967 DO jb = 1, ipi 968 ji = idx_bdy(kbdy)%nbi(jb,kgrd) 969 jj = idx_bdy(kbdy)%nbj(jb,kgrd) 970 zh = SUM(pdta_read_dz(jb,1,:) ) 971 ztrans = 0._wp 972 ztrans_new = 0._wp 973 DO jkb = 1, ipkb ! calculate transport on input grid 974 ztrans = ztrans + pdta_read(jb,1,jkb) * pdta_read_dz(jb, 1,jkb) 977 975 ENDDO 978 ENDIF 979 980 IF(igrd == 3) THEN ! do we need to adjust the transport term? 981 DO ib = 1, ipi 982 zij = idx_bdy(ibdy)%nbi(ib,igrd) 983 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 984 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 985 ztrans = 0._wp 986 ztrans_new = 0._wp 987 DO ik = 1, jpk_bdy ! calculate transport on input grid 988 ztrans = ztrans + dta_read(map%ptr(ib),1,ik) * dta_read_dz(map%ptr(ib),1,ik) 989 ENDDO 990 DO ik = 1, ipk ! calculate transport on model grid 991 ztrans_new = ztrans_new + dta(ib,1,ik) * e3v_n(zij,zjj,ik) * vmask(zij,zjj,ik) 992 ENDDO 993 DO ik = 1, ipk ! make transport correction 994 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 995 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hv_n(zij,zjj) ) * vmask(zij,zjj,ik) 996 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 997 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hv_n(zij,zjj) * vmask(zij,zjj,ik) 998 ENDIF 999 ENDDO 976 DO jk = 1, jpk ! calculate transport on model grid 977 ztrans_new = ztrans_new + pdta(jb,1,jk ) * e3v_n(ji,jj,jk ) * vmask(ji,jj,jk) 1000 978 ENDDO 1001 ENDIF 1002 1003 ELSE ! structured open boundary file 1004 1005 DO ib = 1, ipi 1006 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 1007 ji=map%ptr(ib)-(jj-1)*ilendta 1008 DO ik = 1, jpk_bdy 1009 IF( ( dta_read(ji,jj,ik) == fv ) ) THEN 1010 dta_read_z(ji,jj,ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data 1011 dta_read_dz(ji,jj,ik) = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct 979 DO jk = 1, jpk ! make transport correction 980 IF(ldtotvel) THEN ! bdy data are total velocity so adjust bt transport term to match input data 981 pdta(jb,1,jk) = ( pdta(jb,1,jk) + ( ztrans - ztrans_new ) * r1_hv_n(ji,jj) ) * vmask(ji,jj,jk) 982 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 983 pdta(jb,1,jk) = pdta(jb,1,jk) + ( 0._wp - ztrans_new ) * r1_hv_n(ji,jj) * vmask(ji,jj,jk) 1012 984 ENDIF 1013 985 ENDDO 1014 ENDDO 1015 1016 1017 DO ib = 1, ipi 1018 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 1019 ji=map%ptr(ib)-(jj-1)*ilendta 1020 zij = idx_bdy(ibdy)%nbi(ib,igrd) 1021 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 1022 zh = SUM(dta_read_dz(ji,jj,:) ) 1023 ! Warnings to flag differences in the input and model topgraphy - is this useful/necessary? 1024 SELECT CASE( igrd ) 1025 CASE(1) 1026 IF( ABS( (zh - ht_n(zij,zjj)) / ht_n(zij,zjj)) * tmask(zij,zjj,1) > 0.01_wp ) THEN 1027 WRITE(ibstr,"(I10.10)") map%ptr(ib) 1028 CALL ctl_warn('fld_bdy_interp: T depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 1029 ! IF(lwp) WRITE(*,*) 'DEPTHT', zh, sum(e3t_n(zij,zjj,:), mask=tmask(zij,zjj,:)==1), ht_n(zij,zjj), map%ptr(ib), ib, zij, zjj 1030 ENDIF 1031 CASE(2) 1032 IF( ABS( (zh - hu_n(zij,zjj)) * r1_hu_n(zij,zjj)) * umask(zij,zjj,1) > 0.01_wp ) THEN 1033 WRITE(ibstr,"(I10.10)") map%ptr(ib) 1034 CALL ctl_warn('fld_bdy_interp: U depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 1035 ENDIF 1036 CASE(3) 1037 IF( ABS( (zh - hv_n(zij,zjj)) * r1_hv_n(zij,zjj)) * vmask(zij,zjj,1) > 0.01_wp ) THEN 1038 WRITE(ibstr,"(I10.10)") map%ptr(ib) 1039 CALL ctl_warn('fld_bdy_interp: V depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 1040 ENDIF 1041 END SELECT 1042 DO ik = 1, ipk 1043 SELECT CASE( igrd ) ! coded for sco - need zco and zps option using min 1044 CASE(1) 1045 zl = gdept_n(zij,zjj,ik) ! if using in step could use fsdept instead of gdept_n? 1046 CASE(2) 1047 IF(ln_sco) THEN 1048 zl = ( gdept_n(zij,zjj,ik) + gdept_n(zij+1,zjj,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_n? 1049 ELSE 1050 zl = MIN( gdept_n(zij,zjj,ik), gdept_n(zij+1,zjj,ik) ) 1051 ENDIF 1052 CASE(3) 1053 IF(ln_sco) THEN 1054 zl = ( gdept_n(zij,zjj,ik) + gdept_n(zij,zjj+1,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_n? 1055 ELSE 1056 zl = MIN( gdept_n(zij,zjj,ik), gdept_n(zij,zjj+1,ik) ) 1057 ENDIF 1058 END SELECT 1059 IF( zl < dta_read_z(ji,jj,1) ) THEN ! above the first level of external data 1060 dta(ib,1,ik) = dta_read(ji,jj,1) 1061 ELSEIF( zl > MAXVAL(dta_read_z(ji,jj,:),1) ) THEN ! below the last level of external data 1062 dta(ib,1,ik) = dta_read(ji,jj,MAXLOC(dta_read_z(ji,jj,:),1)) 1063 ELSE ! inbetween : vertical interpolation between ikk & ikk+1 1064 DO ikk = 1, jpkm1_bdy ! when gdept_n(ikk) < zl < gdept_n(ikk+1) 1065 IF( ( (zl-dta_read_z(ji,jj,ikk)) * (zl-dta_read_z(ji,jj,ikk+1)) <= 0._wp) & 1066 & .AND. (dta_read_z(ji,jj,ikk+1) /= fv_alt)) THEN 1067 zi = ( zl - dta_read_z(ji,jj,ikk) ) / & 1068 & ( dta_read_z(ji,jj,ikk+1) - dta_read_z(ji,jj,ikk) ) 1069 dta(ib,1,ik) = dta_read(ji,jj,ikk) + & 1070 & ( dta_read(ji,jj,ikk+1) - dta_read(ji,jj,ikk) ) * zi 1071 ENDIF 1072 END DO 1073 ENDIF 1074 END DO 1075 END DO 1076 1077 IF(igrd == 2) THEN ! do we need to adjust the transport term? 1078 DO ib = 1, ipi 1079 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 1080 ji=map%ptr(ib)-(jj-1)*ilendta 1081 zij = idx_bdy(ibdy)%nbi(ib,igrd) 1082 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 1083 zh = SUM(dta_read_dz(ji,jj,:) ) 1084 ztrans = 0._wp 1085 ztrans_new = 0._wp 1086 DO ik = 1, jpk_bdy ! calculate transport on input grid 1087 ztrans = ztrans + dta_read(ji,jj,ik) * dta_read_dz(ji,jj,ik) 1088 ENDDO 1089 DO ik = 1, ipk ! calculate transport on model grid 1090 ztrans_new = ztrans_new + dta(ib,1,ik) * e3u_n(zij,zjj,ik) * umask(zij,zjj,ik) 1091 ENDDO 1092 DO ik = 1, ipk ! make transport correction 1093 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 1094 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hu_n(zij,zjj) ) * umask(zij,zjj,ik) 1095 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 1096 dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hu_n(zij,zjj) ) * umask(zij,zjj,ik) 1097 ENDIF 1098 ENDDO 1099 ENDDO 1100 ENDIF 1101 1102 IF(igrd == 3) THEN ! do we need to adjust the transport term? 1103 DO ib = 1, ipi 1104 jj = 1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 1105 ji = map%ptr(ib)-(jj-1)*ilendta 1106 zij = idx_bdy(ibdy)%nbi(ib,igrd) 1107 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 1108 zh = SUM(dta_read_dz(ji,jj,:) ) 1109 ztrans = 0._wp 1110 ztrans_new = 0._wp 1111 DO ik = 1, jpk_bdy ! calculate transport on input grid 1112 ztrans = ztrans + dta_read(ji,jj,ik) * dta_read_dz(ji,jj,ik) 1113 ENDDO 1114 DO ik = 1, ipk ! calculate transport on model grid 1115 ztrans_new = ztrans_new + dta(ib,1,ik) * e3v_n(zij,zjj,ik) * vmask(zij,zjj,ik) 1116 ENDDO 1117 DO ik = 1, ipk ! make transport correction 1118 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 1119 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hv_n(zij,zjj) ) * vmask(zij,zjj,ik) 1120 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 1121 dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hv_n(zij,zjj) ) * vmask(zij,zjj,ik) 1122 ENDIF 1123 ENDDO 1124 ENDDO 1125 ENDIF 1126 1127 ENDIF ! endif unstructured or structured 1128 986 ENDDO 987 ENDIF 988 1129 989 END SUBROUTINE fld_bdy_interp 1130 990 … … 1151 1011 imf = SIZE( sd ) 1152 1012 DO ju = 1, imf 1013 IF( TRIM(sd(ju)%clrootname) == 'NOT USED' ) CYCLE 1153 1014 ill = LEN_TRIM( sd(ju)%vcomp ) 1154 1015 DO jn = 2-COUNT((/sd(ju)%ln_tint/)), 2 … … 1159 1020 iv = -1 1160 1021 DO jv = 1, imf 1022 IF( TRIM(sd(jv)%clrootname) == 'NOT USED' ) CYCLE 1161 1023 IF( TRIM(sd(jv)%vcomp) == TRIM(clcomp) ) iv = jv 1162 1024 END DO … … 1197 1059 LOGICAL, OPTIONAL, INTENT(in ) :: ldstop ! stop if open to read a non-existing file (default = .TRUE.) 1198 1060 ! 1199 LOGICAL :: llprevyr ! are we reading previous year file?1200 LOGICAL :: llprevmth ! are we reading previous month file?1201 INTEGER :: iyear, imonth, iday ! first day of the current file in yyyy mm dd1202 INTEGER :: isec_week ! number of seconds since start of the weekly file1203 INTEGER :: indexyr ! year undex (O/1/2: previous/current/next)1204 INTEGER :: iyear_len, imonth_len ! length (days) of iyear and imonth !1205 CHARACTER(len = 256) :: clname ! temporary file name1061 LOGICAL :: llprevyr ! are we reading previous year file? 1062 LOGICAL :: llprevmth ! are we reading previous month file? 1063 INTEGER :: iyear, imonth, iday ! first day of the current file in yyyy mm dd 1064 INTEGER :: isec_week ! number of seconds since start of the weekly file 1065 INTEGER :: indexyr ! year undex (O/1/2: previous/current/next) 1066 REAL(wp) :: zyear_len, zmonth_len ! length (days) of iyear and imonth ! 1067 CHARACTER(len = 256) :: clname ! temporary file name 1206 1068 !!---------------------------------------------------------------------- 1207 1069 IF( PRESENT(kyear) ) THEN ! use given values … … 1254 1116 ! find the last record to be read -> update sdjf%nreclast 1255 1117 indexyr = iyear - nyear + 1 1256 iyear_len = nyear_len( indexyr)1118 zyear_len = REAL(nyear_len( indexyr ), wp) 1257 1119 SELECT CASE ( indexyr ) 1258 CASE ( 0 ) ; imonth_len = 31! previous year -> imonth = 121259 CASE ( 1 ) ; imonth_len = nmonth_len(imonth)1260 CASE ( 2 ) ; imonth_len = 31! next year -> imonth = 11120 CASE ( 0 ) ; zmonth_len = 31. ! previous year -> imonth = 12 1121 CASE ( 1 ) ; zmonth_len = REAL(nmonth_len(imonth), wp) 1122 CASE ( 2 ) ; zmonth_len = 31. ! next year -> imonth = 1 1261 1123 END SELECT 1262 1124 ! 1263 1125 ! last record to be read in the current file 1264 IF ( sdjf% nfreqh == -12) THEN ; sdjf%nreclast = 1 ! yearly mean1265 ELSEIF( sdjf% nfreqh == -1) THEN ! monthly mean1126 IF ( sdjf%freqh == -12. ) THEN ; sdjf%nreclast = 1 ! yearly mean 1127 ELSEIF( sdjf%freqh == -1. ) THEN ! monthly mean 1266 1128 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = 1 1267 1129 ELSE ; sdjf%nreclast = 12 1268 1130 ENDIF 1269 1131 ELSE ! higher frequency mean (in hours) 1270 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = NINT( 24 * imonth_len / sdjf%nfreqh )1271 ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ; sdjf%nreclast = NINT( 24 * 7 / sdjf%nfreqh )1272 ELSEIF( sdjf%cltype == 'daily' ) THEN ; sdjf%nreclast = NINT( 24 / sdjf%nfreqh )1273 ELSE ; sdjf%nreclast = NINT( 24 * iyear_len / sdjf%nfreqh )1132 IF( sdjf%cltype == 'monthly' ) THEN ; sdjf%nreclast = NINT( 24. * zmonth_len / sdjf%freqh ) 1133 ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ; sdjf%nreclast = NINT( 24. * 7. / sdjf%freqh ) 1134 ELSEIF( sdjf%cltype == 'daily' ) THEN ; sdjf%nreclast = NINT( 24. / sdjf%freqh ) 1135 ELSE ; sdjf%nreclast = NINT( 24. * zyear_len / sdjf%freqh ) 1274 1136 ENDIF 1275 1137 ENDIF … … 1299 1161 ! 1300 1162 DO jf = 1, SIZE(sdf) 1301 sdf(jf)%clrootname = TRIM( cdir )//TRIM( sdf_n(jf)%clname ) 1163 sdf(jf)%clrootname = sdf_n(jf)%clname 1164 IF( TRIM(sdf_n(jf)%clname) /= 'NOT USED' ) sdf(jf)%clrootname = TRIM( cdir )//sdf(jf)%clrootname 1302 1165 sdf(jf)%clname = "not yet defined" 1303 sdf(jf)% nfreqh = sdf_n(jf)%nfreqh1166 sdf(jf)%freqh = sdf_n(jf)%freqh 1304 1167 sdf(jf)%clvar = sdf_n(jf)%clvar 1305 1168 sdf(jf)%ln_tint = sdf_n(jf)%ln_tint … … 1308 1171 sdf(jf)%num = -1 1309 1172 sdf(jf)%wgtname = " " 1310 IF( LEN( TRIM(sdf_n(jf)%wname) ) > 0 ) sdf(jf)%wgtname = TRIM( cdir )// TRIM( sdf_n(jf)%wname )1173 IF( LEN( TRIM(sdf_n(jf)%wname) ) > 0 ) sdf(jf)%wgtname = TRIM( cdir )//sdf_n(jf)%wname 1311 1174 sdf(jf)%lsmname = " " 1312 IF( LEN( TRIM(sdf_n(jf)%lname) ) > 0 ) sdf(jf)%lsmname = TRIM( cdir )// TRIM( sdf_n(jf)%lname )1175 IF( LEN( TRIM(sdf_n(jf)%lname) ) > 0 ) sdf(jf)%lsmname = TRIM( cdir )//sdf_n(jf)%lname 1313 1176 sdf(jf)%vcomp = sdf_n(jf)%vcomp 1314 1177 sdf(jf)%rotn(:) = .TRUE. ! pretend to be rotated -> won't try to rotate data before the first call to fld_get … … 1317 1180 IF( sdf(jf)%cltype(1:4) == 'week' .AND. sdf(jf)%ln_clim ) & 1318 1181 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs ln_clim = .FALSE.') 1319 sdf(jf)%nreclast = -1 ! Set to non zero default value to avoid errors, is updated to meaningful value during fld_clopn 1182 sdf(jf)%nreclast = -1 ! Set to non zero default value to avoid errors, is updated to meaningful value during fld_clopn 1183 sdf(jf)%igrd = 0 1184 sdf(jf)%ibdy = 0 1185 sdf(jf)%imap => NULL() 1186 sdf(jf)%ltotvel = .FALSE. 1187 sdf(jf)%lzint = .FALSE. 1320 1188 END DO 1321 1189 ! … … 1331 1199 DO jf = 1, SIZE(sdf) 1332 1200 WRITE(numout,*) ' root filename: ' , TRIM( sdf(jf)%clrootname ), ' variable name: ', TRIM( sdf(jf)%clvar ) 1333 WRITE(numout,*) ' frequency: ' , sdf(jf)% nfreqh, &1201 WRITE(numout,*) ' frequency: ' , sdf(jf)%freqh , & 1334 1202 & ' time interp: ' , sdf(jf)%ln_tint , & 1335 1203 & ' climatology: ' , sdf(jf)%ln_clim
Note: See TracChangeset
for help on using the changeset viewer.
