- Timestamp:
- 2020-09-29T12:41:06+02:00 (4 years ago)
- Location:
- NEMO/branches/2020/r12377_ticket2386
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/SBC/fldread.F90 (modified) (50 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/r12377_ticket2386
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev @HEADext/AGRIF5 ^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 8 9 9 # SETTE 10 ^/utils/CI/sette@ HEADsette10 ^/utils/CI/sette@13507 sette
-
- Property svn:externals
-
NEMO/branches/2020/r12377_ticket2386/src/OCE/SBC/fldread.F90
r12511 r13540 53 53 LOGICAL :: ln_tint ! time interpolation or not (T/F) 54 54 LOGICAL :: ln_clim ! climatology or not (T/F) 55 CHARACTER(len = 8) :: cl type! type of data file 'daily', 'monthly' or yearly'55 CHARACTER(len = 8) :: clftyp ! type of data file 'daily', 'monthly' or yearly' 56 56 CHARACTER(len = 256) :: wname ! generic name of a NetCDF weights file to be used, blank if not 57 57 CHARACTER(len = 34) :: vcomp ! symbolic component name if a vector that needs rotation … … 69 69 LOGICAL :: ln_tint ! time interpolation or not (T/F) 70 70 LOGICAL :: ln_clim ! climatology or not (T/F) 71 CHARACTER(len = 8) :: cltype ! type of data file 'daily', 'monthly' or yearly' 71 CHARACTER(len = 8) :: clftyp ! type of data file 'daily', 'monthly' or yearly' 72 CHARACTER(len = 1) :: cltype ! nature of grid-points: T, U, V... 73 REAL(wp) :: zsgn ! -1. the sign change across the north fold, = 1. otherwise 72 74 INTEGER :: num ! iom id of the jpfld files to be read 73 INTEGER , DIMENSION(2) :: nrec_b ! before record (1: index, 2: second since Jan. 1st 00h of nit000 year) 74 INTEGER , DIMENSION(2) :: nrec_a ! after record (1: index, 2: second since Jan. 1st 00h of nit000 year) 75 INTEGER , ALLOCATABLE, DIMENSION(: ) :: nrecsec ! 76 REAL(wp), ALLOCATABLE, DIMENSION(:,:,: ) :: fnow ! input fields interpolated to now time step 77 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: fdta ! 2 consecutive record of input fields 75 INTEGER , DIMENSION(2,2) :: nrec ! before/after record (1: index, 2: second since Jan. 1st 00h of yr nit000) 76 INTEGER :: nbb ! index of before values 77 INTEGER :: naa ! index of after values 78 INTEGER , ALLOCATABLE, DIMENSION(:) :: nrecsec ! 79 REAL(wp), POINTER, DIMENSION(:,:,: ) :: fnow ! input fields interpolated to now time step 80 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: fdta ! 2 consecutive record of input fields 78 81 CHARACTER(len = 256) :: wgtname ! current name of the NetCDF weight file acting as a key 79 82 ! ! into the WGTLIST structure … … 127 130 !! * Substitutions 128 131 # include "do_loop_substitute.h90" 132 # include "domzgr_substitute.h90" 129 133 !!---------------------------------------------------------------------- 130 134 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 156 160 INTEGER :: jf ! dummy indices 157 161 INTEGER :: isecsbc ! number of seconds between Jan. 1st 00h of nit000 year and the middle of sbc time step 162 INTEGER :: ibb, iaa ! shorter name for sd(jf)%nbb and sd(jf)%naa 158 163 LOGICAL :: ll_firstcall ! true if this is the first call to fld_read for this set of fields 159 164 REAL(wp) :: zt_offset ! local time offset variable … … 203 208 IF( TRIM(sd(jf)%clrootname) == 'NOT USED' ) CYCLE 204 209 ! 210 ibb = sd(jf)%nbb ; iaa = sd(jf)%naa 211 ! 205 212 IF( sd(jf)%ln_tint ) THEN ! temporal interpolation 206 213 IF(lwp .AND. kt - nit000 <= 100 ) THEN … … 208 215 & "', records b/a: ', i6.4, '/', i6.4, ' (days ', f9.4,'/', f9.4, ')')" 209 216 WRITE(numout, clfmt) TRIM( sd(jf)%clvar ), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 210 & 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)/rday217 & sd(jf)%nrec(1,ibb), sd(jf)%nrec(1,iaa), REAL(sd(jf)%nrec(2,ibb),wp)/rday, REAL(sd(jf)%nrec(2,iaa),wp)/rday 211 218 WRITE(numout, *) ' zt_offset is : ',zt_offset 212 219 ENDIF 213 220 ! temporal interpolation weights 214 ztinta = REAL( isecsbc - sd(jf)%nrec _b(2), wp ) / REAL( sd(jf)%nrec_a(2) - sd(jf)%nrec_b(2), wp )221 ztinta = REAL( isecsbc - sd(jf)%nrec(2,ibb), wp ) / REAL( sd(jf)%nrec(2,iaa) - sd(jf)%nrec(2,ibb), wp ) 215 222 ztintb = 1. - ztinta 216 sd(jf)%fnow(:,:,:) = ztintb * sd(jf)%fdta(:,:,:, 1) + ztinta * sd(jf)%fdta(:,:,:,2)223 sd(jf)%fnow(:,:,:) = ztintb * sd(jf)%fdta(:,:,:,ibb) + ztinta * sd(jf)%fdta(:,:,:,iaa) 217 224 ELSE ! nothing to do... 218 225 IF(lwp .AND. kt - nit000 <= 100 ) THEN … … 220 227 & "', record: ', i6.4, ' (days ', f9.4, ' <-> ', f9.4, ')')" 221 228 WRITE(numout, clfmt) TRIM(sd(jf)%clvar), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 222 & sd(jf)%nrec _a(1), REAL(sd(jf)%nrec_b(2),wp)/rday, REAL(sd(jf)%nrec_a(2),wp)/rday229 & sd(jf)%nrec(1,iaa), REAL(sd(jf)%nrec(2,ibb),wp)/rday, REAL(sd(jf)%nrec(2,iaa),wp)/rday 223 230 ENDIF 224 231 ENDIF … … 250 257 ! 251 258 CALL fld_clopn( sdjf ) 252 sdjf%nrec _a(:) = (/ 1, nflag /) ! default definition to force flp_update to read the file.259 sdjf%nrec(:,sdjf%naa) = (/ 1, nflag /) ! default definition to force flp_update to read the file. 253 260 ! 254 261 END SUBROUTINE fld_init … … 261 268 !! ** Purpose : Compute 262 269 !! if sdjf%ln_tint = .TRUE. 263 !! nrec _a: record number and its time (nrec_b is obtained from nrec_awhen swapping)270 !! nrec(:,iaa): record number and its time (nrec(:,ibb) is obtained from nrec(:,iaa) when swapping) 264 271 !! if sdjf%ln_tint = .FALSE. 265 !! nrec _a(1): record number266 !! nrec _b(2) and nrec_a(2): time of the beginning and end of the record272 !! nrec(1,iaa): record number 273 !! nrec(2,ibb) and nrec(2,iaa): time of the beginning and end of the record 267 274 !!---------------------------------------------------------------------- 268 275 INTEGER , INTENT(in ) :: ksecsbc ! … … 270 277 INTEGER , OPTIONAL, INTENT(in ) :: Kmm ! ocean time level index 271 278 ! 272 INTEGER :: ja ! end of this record (in seconds) 273 !!---------------------------------------------------------------------- 274 ! 275 IF( ksecsbc > sdjf%nrec_a(2) ) THEN ! --> we need to update after data 279 INTEGER :: ja ! end of this record (in seconds) 280 INTEGER :: ibb, iaa ! shorter name for sdjf%nbb and sdjf%naa 281 !!---------------------------------------------------------------------- 282 ibb = sdjf%nbb ; iaa = sdjf%naa 283 ! 284 IF( ksecsbc > sdjf%nrec(2,iaa) ) THEN ! --> we need to update after data 276 285 277 ! find where is the new after record... (it is not necessary sdjf%nrec _a(1)+1 )278 ja = sdjf%nrec _a(1)286 ! find where is the new after record... (it is not necessary sdjf%nrec(1,iaa)+1 ) 287 ja = sdjf%nrec(1,iaa) 279 288 DO WHILE ( ksecsbc >= sdjf%nrecsec(ja) .AND. ja < sdjf%nreclast ) ! Warning: make sure ja <= sdjf%nreclast in this test 280 289 ja = ja + 1 … … 283 292 284 293 ! if ln_tint and if the new after is not ja+1, we need also to update after data before the swap 285 ! so, after the swap, sdjf%nrec _b(2) will still be the closest value located just before ksecsbc286 IF( sdjf%ln_tint .AND. ( ja > sdjf%nrec _a(1) + 1 .OR. sdjf%nrec_a(2) == nflag ) ) THEN287 sdjf%nrec _a(:) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! update nrec_awith before information288 CALL fld_get( sdjf, Kmm ) ! read after data that will be used as before data294 ! so, after the swap, sdjf%nrec(2,ibb) will still be the closest value located just before ksecsbc 295 IF( sdjf%ln_tint .AND. ( ja > sdjf%nrec(1,iaa) + 1 .OR. sdjf%nrec(2,iaa) == nflag ) ) THEN 296 sdjf%nrec(:,iaa) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! update nrec(:,iaa) with before information 297 CALL fld_get( sdjf, Kmm ) ! read after data that will be used as before data 289 298 ENDIF 290 299 … … 309 318 ! if ln_tint and if after is not the first record, we must (potentially again) update after data before the swap 310 319 IF( sdjf%ln_tint .AND. ja > 1 ) THEN 311 IF( sdjf%nrecsec(0) /= nflag ) THEN ! no trick used: after file is not the current file312 sdjf%nrec _a(:) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! update nrec_awith before information313 CALL fld_get( sdjf, Kmm ) ! read after data that will be used as before data320 IF( sdjf%nrecsec(0) /= nflag ) THEN ! no trick used: after file is not the current file 321 sdjf%nrec(:,iaa) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! update nrec(:,iaa) with before information 322 CALL fld_get( sdjf, Kmm ) ! read after data that will be used as before data 314 323 ENDIF 315 324 ENDIF … … 317 326 ENDIF 318 327 319 IF( sdjf%ln_tint ) THEN 320 ! Swap data 321 sdjf%nrec_b(:) = sdjf%nrec_a(:) ! swap before record informations 322 sdjf%rotn(1) = sdjf%rotn(2) ! swap before rotate informations 323 sdjf%fdta(:,:,:,1) = sdjf%fdta(:,:,:,2) ! swap before record field 324 ELSE 325 sdjf%nrec_b(:) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! only for print 328 IF( sdjf%ln_tint ) THEN ! Swap data 329 sdjf%nbb = sdjf%naa ! swap indices 330 sdjf%naa = 3 - sdjf%naa ! = 2(1) if naa == 1(2) 331 ELSE ! No swap 332 sdjf%nrec(:,ibb) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! only for print 326 333 ENDIF 327 334 328 335 ! read new after data 329 sdjf%nrec _a(:) = (/ ja, sdjf%nrecsec(ja) /) ! update nrec_aas it is used by fld_get330 CALL fld_get( sdjf, Kmm ) ! read after data (with nrec_ainformations)336 sdjf%nrec(:,sdjf%naa) = (/ ja, sdjf%nrecsec(ja) /) ! update nrec(:,naa) as it is used by fld_get 337 CALL fld_get( sdjf, Kmm ) ! read after data (with nrec(:,naa) informations) 331 338 332 339 ENDIF … … 345 352 ! 346 353 INTEGER :: ipk ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 354 INTEGER :: iaa ! shorter name for sdjf%naa 347 355 INTEGER :: iw ! index into wgts array 348 INTEGER :: ipdom ! index of the domain349 356 INTEGER :: idvar ! variable ID 350 357 INTEGER :: idmspc ! number of spatial dimensions 351 358 LOGICAL :: lmoor ! C1D case: point data 352 !!--------------------------------------------------------------------- 353 ! 354 ipk = SIZE( sdjf%fnow, 3 ) 355 ! 356 IF( ASSOCIATED(sdjf%imap) ) THEN 357 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), & 358 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint, Kmm ) 359 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), & 360 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint, Kmm ) 361 ENDIF 362 ELSE IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN 359 REAL(wp), DIMENSION(:,:,:), POINTER :: dta_alias ! short cut 360 !!--------------------------------------------------------------------- 361 iaa = sdjf%naa 362 ! 363 IF( sdjf%ln_tint ) THEN ; dta_alias => sdjf%fdta(:,:,:,iaa) 364 ELSE ; dta_alias => sdjf%fnow(:,:,: ) 365 ENDIF 366 ipk = SIZE( dta_alias, 3 ) 367 ! 368 IF( ASSOCIATED(sdjf%imap) ) THEN ! BDY case 369 CALL fld_map( sdjf%num, sdjf%clvar, dta_alias(:,:,:), sdjf%nrec(1,iaa), & 370 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint, Kmm ) 371 ELSE IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN ! On-the-fly interpolation 363 372 CALL wgt_list( sdjf, iw ) 364 IF( sdjf%ln_tint ) THEN ; CALL fld_interp( sdjf%num, sdjf%clvar, iw, ipk, sdjf%fdta(:,:,:,2), & 365 & sdjf%nrec_a(1), sdjf%lsmname ) 366 ELSE ; CALL fld_interp( sdjf%num, sdjf%clvar, iw, ipk, sdjf%fnow(:,:,: ), & 367 & sdjf%nrec_a(1), sdjf%lsmname ) 368 ENDIF 369 ELSE 370 IF( SIZE(sdjf%fnow, 1) == jpi ) THEN ; ipdom = jpdom_data 371 ELSE ; ipdom = jpdom_unknown 372 ENDIF 373 CALL fld_interp( sdjf%num, sdjf%clvar, iw, ipk, dta_alias(:,:,:), sdjf%nrec(1,iaa), sdjf%lsmname ) 374 CALL lbc_lnk( 'fldread', dta_alias(:,:,:), sdjf%cltype, sdjf%zsgn, kfillmode = jpfillcopy ) 375 ELSE ! default case 373 376 ! C1D case: If product of spatial dimensions == ipk, then x,y are of 374 377 ! size 1 (point/mooring data): this must be read onto the central grid point 375 378 idvar = iom_varid( sdjf%num, sdjf%clvar ) 376 379 idmspc = iom_file ( sdjf%num )%ndims( idvar ) 377 IF( iom_file( sdjf%num )%luld( idvar ) ) idmspc = idmspc - 1 378 lmoor = ( idmspc == 0 .OR. PRODUCT( iom_file( sdjf%num )%dimsz( 1:MAX(idmspc,1) ,idvar ) ) == ipk ) 379 ! 380 SELECT CASE( ipk ) 381 CASE(1) 382 IF( lk_c1d .AND. lmoor ) THEN 383 IF( sdjf%ln_tint ) THEN 384 CALL iom_get( sdjf%num, sdjf%clvar, sdjf%fdta(2,2,1,2), sdjf%nrec_a(1) ) 385 CALL lbc_lnk( 'fldread', sdjf%fdta(:,:,1,2),'Z',1. ) 386 ELSE 387 CALL iom_get( sdjf%num, sdjf%clvar, sdjf%fnow(2,2,1 ), sdjf%nrec_a(1) ) 388 CALL lbc_lnk( 'fldread', sdjf%fnow(:,:,1 ),'Z',1. ) 389 ENDIF 390 ELSE 391 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_a(1) ) 392 ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,1 ), sdjf%nrec_a(1) ) 393 ENDIF 394 ENDIF 395 CASE DEFAULT 396 IF(lk_c1d .AND. lmoor ) THEN 397 IF( sdjf%ln_tint ) THEN 398 CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, sdjf%fdta(2,2,:,2), sdjf%nrec_a(1) ) 399 CALL lbc_lnk( 'fldread', sdjf%fdta(:,:,:,2),'Z',1. ) 400 ELSE 401 CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, sdjf%fnow(2,2,: ), sdjf%nrec_a(1) ) 402 CALL lbc_lnk( 'fldread', sdjf%fnow(:,:,: ),'Z',1. ) 403 ENDIF 404 ELSE 405 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1) ) 406 ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1) ) 407 ENDIF 408 ENDIF 409 END SELECT 410 ENDIF 411 ! 412 sdjf%rotn(2) = .false. ! vector not yet rotated 380 IF( iom_file( sdjf%num )%luld( idvar ) ) idmspc = idmspc - 1 ! id of the last spatial dimension 381 lmoor = ( idmspc == 0 .OR. PRODUCT( iom_file( sdjf%num )%dimsz( 1:MAX(idmspc,1) ,idvar ) ) == ipk ) 382 ! 383 IF( lk_c1d .AND. lmoor ) THEN 384 CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, dta_alias(2,2,:), sdjf%nrec(1,iaa) ) ! jpdom_unknown -> no lbc_lnk 385 CALL lbc_lnk( 'fldread', dta_alias(:,:,:), 'T', 1., kfillmode = jpfillcopy ) 386 ELSE 387 CALL iom_get( sdjf%num, jpdom_global, sdjf%clvar, dta_alias(:,:,:), sdjf%nrec(1,iaa), & 388 & sdjf%cltype, sdjf%zsgn, kfill = jpfillcopy ) 389 ENDIF 390 ENDIF 391 ! 392 sdjf%rotn(iaa) = .false. ! vector not yet rotated 413 393 ! 414 394 END SUBROUTINE fld_get … … 446 426 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read_z ! work space local data requiring vertical interpolation 447 427 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read_dz ! work space local data requiring vertical interpolation 448 CHARACTER(LEN=1),DIMENSION(3) :: cl grid428 CHARACTER(LEN=1),DIMENSION(3) :: cltype 449 429 LOGICAL :: lluld ! is the variable using the unlimited dimension 450 430 LOGICAL :: llzint ! local value of ldzint 451 431 !!--------------------------------------------------------------------- 452 432 ! 453 cl grid= (/'t','u','v'/)433 cltype = (/'t','u','v'/) 454 434 ! 455 435 ipi = SIZE( pdta, 1 ) … … 486 466 IF( ipkb /= ipk .OR. llzint ) THEN ! boundary data not on model vertical grid : vertical interpolation 487 467 ! 488 IF( ipk == jpk .AND. iom_varid(knum,'gdep'//cl grid(kgrd)) /= -1 .AND. iom_varid(knum,'e3'//clgrid(kgrd)) /= -1 ) THEN468 IF( ipk == jpk .AND. iom_varid(knum,'gdep'//cltype(kgrd)) /= -1 .AND. iom_varid(knum,'e3'//cltype(kgrd)) /= -1 ) THEN 489 469 490 470 ALLOCATE( zdta_read(ipi,ipj,ipkb), zdta_read_z(ipi,ipj,ipkb), zdta_read_dz(ipi,ipj,ipkb) ) 491 471 492 472 CALL fld_map_core( zz_read, kmap, zdta_read ) 493 CALL iom_get ( knum, jpdom_unknown, 'gdep'//cl grid(kgrd), zz_read ) ! read only once? Potential temporal evolution?473 CALL iom_get ( knum, jpdom_unknown, 'gdep'//cltype(kgrd), zz_read ) ! read only once? Potential temporal evolution? 494 474 CALL fld_map_core( zz_read, kmap, zdta_read_z ) 495 CALL iom_get ( knum, jpdom_unknown, 'e3'//cl grid(kgrd), zz_read ) ! read only once? Potential temporal evolution?475 CALL iom_get ( knum, jpdom_unknown, 'e3'//cltype(kgrd), zz_read ) ! read only once? Potential temporal evolution? 496 476 CALL fld_map_core( zz_read, kmap, zdta_read_dz ) 497 477 … … 503 483 IF( ipk /= jpk ) CALL ctl_stop( 'fld_map : this should be an impossible case...' ) 504 484 WRITE(ctmp1,*) 'fld_map : vertical interpolation for bdy variable '//TRIM(cdvar)//' requires ' 505 IF( iom_varid(knum, 'gdep'//cl grid(kgrd)) == -1 ) CALL ctl_stop( ctmp1//'gdep'//clgrid(kgrd)//' variable' )506 IF( iom_varid(knum, 'e3'//cl grid(kgrd)) == -1 ) CALL ctl_stop( ctmp1// 'e3'//clgrid(kgrd)//' variable' )485 IF( iom_varid(knum, 'gdep'//cltype(kgrd)) == -1 ) CALL ctl_stop( ctmp1//'gdep'//cltype(kgrd)//' variable' ) 486 IF( iom_varid(knum, 'e3'//cltype(kgrd)) == -1 ) CALL ctl_stop( ctmp1// 'e3'//cltype(kgrd)//' variable' ) 507 487 508 488 ENDIF … … 632 612 zdhalf(jk) = zdhalf(jk-1) + e3v(ji,jj,jk-1,Kmm) 633 613 zdepth(jk) = zcoef * ( zdhalf(jk ) + 0.5_wp * e3vw(ji,jj,jk,Kmm)) & 634 &+ (1._wp-zcoef) * ( zdepth(jk-1) + e3vw(ji,jj,jk,Kmm))614 + (1._wp-zcoef) * ( zdepth(jk-1) + e3vw(ji,jj,jk,Kmm)) 635 615 END DO 636 616 END SELECT … … 727 707 CHARACTER (LEN=100) :: clcomp ! dummy weight name 728 708 REAL(wp), DIMENSION(jpi,jpj) :: utmp, vtmp ! temporary arrays for vector rotation 709 REAL(wp), DIMENSION(:,:,:), POINTER :: dta_u, dta_v ! short cut 729 710 !!--------------------------------------------------------------------- 730 711 ! … … 746 727 END DO 747 728 IF( iv > 0 ) THEN ! fields ju and iv are two components which need to be rotated together 729 IF( sd(ju)%ln_tint ) THEN ; dta_u => sd(ju)%fdta(:,:,:,jn) ; dta_v => sd(iv)%fdta(:,:,:,jn) 730 ELSE ; dta_u => sd(ju)%fnow(:,:,: ) ; dta_v => sd(iv)%fnow(:,:,: ) 731 ENDIF 748 732 DO jk = 1, SIZE( sd(ju)%fnow, 3 ) 749 IF( sd(ju)%ln_tint )THEN 750 CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->i', utmp(:,:) ) 751 CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->j', vtmp(:,:) ) 752 sd(ju)%fdta(:,:,jk,jn) = utmp(:,:) ; sd(iv)%fdta(:,:,jk,jn) = vtmp(:,:) 753 ELSE 754 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->i', utmp(:,:) ) 755 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->j', vtmp(:,:) ) 756 sd(ju)%fnow(:,:,jk ) = utmp(:,:) ; sd(iv)%fnow(:,:,jk ) = vtmp(:,:) 757 ENDIF 733 CALL rot_rep( dta_u(:,:,jk), dta_v(:,:,jk), 'T', 'en->i', utmp(:,:) ) 734 CALL rot_rep( dta_u(:,:,jk), dta_v(:,:,jk), 'T', 'en->j', vtmp(:,:) ) 735 dta_u(:,:,jk) = utmp(:,:) ; dta_v(:,:,jk) = vtmp(:,:) 758 736 END DO 759 737 sd(ju)%rotn(jn) = .TRUE. ! vector was rotated … … 801 779 802 780 ! current file parameters 803 IF( sdjf%cl type(1:4) == 'week' ) THEN! find the day of the beginning of the current week804 isecwk = ksec_week( sdjf%cl type(6:8) ) ! seconds between the beginning of the week and half of current time step805 llprevmt = isecwk > nsec_month 781 IF( sdjf%clftyp(1:4) == 'week' ) THEN ! find the day of the beginning of the current week 782 isecwk = ksec_week( sdjf%clftyp(6:8) ) ! seconds between the beginning of the week and half of current time step 783 llprevmt = isecwk > nsec_month ! longer time since beginning of the current week than the current month 806 784 llprevyr = llprevmt .AND. nmonth == 1 807 785 iyr = nyear - COUNT((/llprevyr/)) 808 786 imt = nmonth - COUNT((/llprevmt/)) + 12 * COUNT((/llprevyr/)) 809 787 idy = nday + nmonth_len(nmonth-1) * COUNT((/llprevmt/)) - isecwk / idaysec 810 isecwk = nsec_year - isecwk ! seconds between 00h jan 1st of current year and current week beginning788 isecwk = nsec_year - isecwk ! seconds between 00h jan 1st of current year and current week beginning 811 789 ELSE 812 790 iyr = nyear … … 818 796 ! previous file parameters 819 797 IF( llprev ) THEN 820 IF( sdjf%cl type(1:4) == 'week' ) THEN! find the day of the beginning of previous week821 isecwk = isecwk + 7 * idaysec ! seconds between the beginning of previous week and half of the time step822 llprevmt = isecwk > nsec_month 798 IF( sdjf%clftyp(1:4) == 'week' ) THEN ! find the day of the beginning of previous week 799 isecwk = isecwk + 7 * idaysec ! seconds between the beginning of previous week and half of the time step 800 llprevmt = isecwk > nsec_month ! longer time since beginning of the previous week than the current month 823 801 llprevyr = llprevmt .AND. nmonth == 1 824 802 iyr = nyear - COUNT((/llprevyr/)) 825 803 imt = nmonth - COUNT((/llprevmt/)) + 12 * COUNT((/llprevyr/)) 826 804 idy = nday + nmonth_len(nmonth-1) * COUNT((/llprevmt/)) - isecwk / idaysec 827 isecwk = nsec_year - isecwk ! seconds between 00h jan 1st of current year and previous week beginning805 isecwk = nsec_year - isecwk ! seconds between 00h jan 1st of current year and previous week beginning 828 806 ELSE 829 idy = nday - COUNT((/ sdjf%cl type== 'daily' /))830 imt = nmonth - COUNT((/ sdjf%cl type== 'monthly' .OR. idy == 0 /))831 iyr = nyear - COUNT((/ sdjf%cl type== 'yearly' .OR. imt == 0 /))807 idy = nday - COUNT((/ sdjf%clftyp == 'daily' /)) 808 imt = nmonth - COUNT((/ sdjf%clftyp == 'monthly' .OR. idy == 0 /)) 809 iyr = nyear - COUNT((/ sdjf%clftyp == 'yearly' .OR. imt == 0 /)) 832 810 IF( idy == 0 ) idy = nmonth_len(imt) 833 811 IF( imt == 0 ) imt = 12 … … 838 816 ! next file parameters 839 817 IF( llnext ) THEN 840 IF( sdjf%cl type(1:4) == 'week' ) THEN! find the day of the beginning of next week841 isecwk = 7 * idaysec - isecwk ! seconds between half of the time step and the beginning of next week818 IF( sdjf%clftyp(1:4) == 'week' ) THEN ! find the day of the beginning of next week 819 isecwk = 7 * idaysec - isecwk ! seconds between half of the time step and the beginning of next week 842 820 llnextmt = isecwk > ( nmonth_len(nmonth)*idaysec - nsec_month ) ! larger than the seconds to the end of the month 843 821 llnextyr = llnextmt .AND. nmonth == 12 … … 845 823 imt = nmonth + COUNT((/llnextmt/)) - 12 * COUNT((/llnextyr/)) 846 824 idy = nday - nmonth_len(nmonth) * COUNT((/llnextmt/)) + isecwk / idaysec + 1 847 isecwk = nsec_year + isecwk ! seconds between 00h jan 1st of current year and next week beginning825 isecwk = nsec_year + isecwk ! seconds between 00h jan 1st of current year and next week beginning 848 826 ELSE 849 idy = nday + COUNT((/ sdjf%cl type== 'daily' /))850 imt = nmonth + COUNT((/ sdjf%cl type== 'monthly' .OR. idy > nmonth_len(nmonth) /))851 iyr = nyear + COUNT((/ sdjf%cl type== 'yearly' .OR. imt == 13 /))827 idy = nday + COUNT((/ sdjf%clftyp == 'daily' /)) 828 imt = nmonth + COUNT((/ sdjf%clftyp == 'monthly' .OR. idy > nmonth_len(nmonth) /)) 829 iyr = nyear + COUNT((/ sdjf%clftyp == 'yearly' .OR. imt == 13 /)) 852 830 IF( idy > nmonth_len(nmonth) ) idy = 1 853 831 IF( imt == 13 ) imt = 1 … … 866 844 IF ( NINT(sdjf%freqh) == -12 ) THEN ; ireclast = 1 ! yearly mean: consider only 1 record 867 845 ELSEIF( NINT(sdjf%freqh) == -1 ) THEN ! monthly mean: 868 IF( sdjf%cl type== 'monthly' ) THEN ; ireclast = 1 ! consider that the file has 1 record846 IF( sdjf%clftyp == 'monthly' ) THEN ; ireclast = 1 ! consider that the file has 1 record 869 847 ELSE ; ireclast = 12 ! consider that the file has 12 record 870 848 ENDIF 871 849 ELSE ! higher frequency mean (in hours) 872 IF( sdjf%cl type== 'monthly' ) THEN ; ireclast = NINT( 24. * REAL(nmonth_len(indexmt), wp) / sdjf%freqh )873 ELSEIF( sdjf%cl type(1:4) == 'week' ) THEN ; ireclast = NINT( 24. * 7. / sdjf%freqh )874 ELSEIF( sdjf%cl type== 'daily' ) THEN ; ireclast = NINT( 24. / sdjf%freqh )850 IF( sdjf%clftyp == 'monthly' ) THEN ; ireclast = NINT( 24. * REAL(nmonth_len(indexmt), wp) / sdjf%freqh ) 851 ELSEIF( sdjf%clftyp(1:4) == 'week' ) THEN ; ireclast = NINT( 24. * 7. / sdjf%freqh ) 852 ELSEIF( sdjf%clftyp == 'daily' ) THEN ; ireclast = NINT( 24. / sdjf%freqh ) 875 853 ELSE ; ireclast = NINT( 24. * REAL( nyear_len(indexyr), wp) / sdjf%freqh ) 876 854 ENDIF … … 890 868 sdjf%nrecsec(1) = sdjf%nrecsec(0) + nyear_len( indexyr ) * idaysec 891 869 ELSEIF( NINT(sdjf%freqh) == -1 ) THEN ! monthly mean: 892 IF( sdjf%cl type== 'monthly' ) THEN ! monthly file870 IF( sdjf%clftyp == 'monthly' ) THEN ! monthly file 893 871 sdjf%nrecsec(0 ) = nsec1jan000 + nmonth_beg(indexmt ) 894 872 sdjf%nrecsec(1 ) = nsec1jan000 + nmonth_beg(indexmt+1) … … 898 876 ENDIF 899 877 ELSE ! higher frequency mean (in hours) 900 IF( sdjf%cl type== 'monthly' ) THEN ; istart = nsec1jan000 + nmonth_beg(indexmt)901 ELSEIF( sdjf%cl type(1:4) == 'week' ) THEN ; istart = nsec1jan000 + isecwk902 ELSEIF( sdjf%cl type== 'daily' ) THEN ; istart = nsec1jan000 + nmonth_beg(indexmt) + ( idy - 1 ) * idaysec878 IF( sdjf%clftyp == 'monthly' ) THEN ; istart = nsec1jan000 + nmonth_beg(indexmt) 879 ELSEIF( sdjf%clftyp(1:4) == 'week' ) THEN ; istart = nsec1jan000 + isecwk 880 ELSEIF( sdjf%clftyp == 'daily' ) THEN ; istart = nsec1jan000 + nmonth_beg(indexmt) + ( idy - 1 ) * idaysec 903 881 ELSEIF( indexyr == 0 ) THEN ; istart = nsec1jan000 - nyear_len( 0 ) * idaysec 904 882 ELSEIF( indexyr == 2 ) THEN ; istart = nsec1jan000 + nyear_len( 1 ) * idaysec … … 941 919 IF( sdjf%num <= 0 .OR. .NOT. sdjf%ln_clim ) THEN 942 920 IF( sdjf%num > 0 ) CALL iom_close( sdjf%num ) ! close file if already open 943 CALL iom_open( sdjf%clname, sdjf%num, ldstop = llstop, ldiof = LEN (TRIM(sdjf%wgtname)) > 0 )921 CALL iom_open( sdjf%clname, sdjf%num, ldstop = llstop, ldiof = LEN_TRIM(sdjf%wgtname) > 0 ) 944 922 ENDIF 945 923 ! … … 963 941 ENDIF 964 942 ! 965 CALL iom_open( sdjf%clname, sdjf%num, ldiof = LEN (TRIM(sdjf%wgtname)) > 0 )943 CALL iom_open( sdjf%clname, sdjf%num, ldiof = LEN_TRIM(sdjf%wgtname) > 0 ) 966 944 ! 967 945 ENDIF … … 996 974 sdf(jf)%ln_tint = sdf_n(jf)%ln_tint 997 975 sdf(jf)%ln_clim = sdf_n(jf)%ln_clim 998 sdf(jf)%cltype = sdf_n(jf)%cltype 976 sdf(jf)%clftyp = sdf_n(jf)%clftyp 977 sdf(jf)%cltype = 'T' ! by default don't do any call to lbc_lnk in iom_get 978 sdf(jf)%zsgn = 1. ! by default don't do change signe across the north fold 999 979 sdf(jf)%num = -1 980 sdf(jf)%nbb = 1 ! start with before data in 1 981 sdf(jf)%naa = 2 ! start with after data in 2 1000 982 sdf(jf)%wgtname = " " 1001 983 IF( LEN( TRIM(sdf_n(jf)%wname) ) > 0 ) sdf(jf)%wgtname = TRIM( cdir )//sdf_n(jf)%wname … … 1004 986 sdf(jf)%vcomp = sdf_n(jf)%vcomp 1005 987 sdf(jf)%rotn(:) = .TRUE. ! pretend to be rotated -> won't try to rotate data before the first call to fld_get 1006 IF( sdf(jf)%cl type(1:4) == 'week' .AND. nn_leapy == 0 ) &988 IF( sdf(jf)%clftyp(1:4) == 'week' .AND. nn_leapy == 0 ) & 1007 989 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs nn_leapy = 1') 1008 IF( sdf(jf)%cl type(1:4) == 'week' .AND. sdf(jf)%ln_clim ) &990 IF( sdf(jf)%clftyp(1:4) == 'week' .AND. sdf(jf)%ln_clim ) & 1009 991 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs ln_clim = .FALSE.') 1010 992 sdf(jf)%nreclast = -1 ! Set to non zero default value to avoid errors, is updated to meaningful value during fld_clopn … … 1032 1014 WRITE(numout,*) ' weights: ' , TRIM( sdf(jf)%wgtname ), & 1033 1015 & ' pairing: ' , TRIM( sdf(jf)%vcomp ), & 1034 & ' data type: ' , sdf(jf)%cl type, &1016 & ' data type: ' , sdf(jf)%clftyp , & 1035 1017 & ' land/sea mask:' , TRIM( sdf(jf)%lsmname ) 1036 1018 call flush(numout) … … 1050 1032 !!---------------------------------------------------------------------- 1051 1033 TYPE( FLD ), INTENT(in ) :: sd ! field with name of weights file 1052 INTEGER , INTENT( inout) :: kwgt ! index of weights1034 INTEGER , INTENT( out) :: kwgt ! index of weights 1053 1035 ! 1054 1036 INTEGER :: kw, nestid ! local integer 1055 LOGICAL :: found ! local logical1056 1037 !!---------------------------------------------------------------------- 1057 1038 ! 1058 1039 !! search down linked list 1059 1040 !! weights filename is either present or we hit the end of the list 1060 found = .FALSE.1061 1041 ! 1062 1042 !! because agrif nest part of filenames are now added in iom_open … … 1068 1048 #endif 1069 1049 DO kw = 1, nxt_wgt-1 1070 IF( TRIM(ref_wgts(kw)%wgtname) == TRIM(sd%wgtname).AND. &1071 ref_wgts(kw)%nestid == nestid) THEN1050 IF( ref_wgts(kw)%wgtname == sd%wgtname .AND. & 1051 ref_wgts(kw)%nestid == nestid) THEN 1072 1052 kwgt = kw 1073 found = .TRUE. 1074 EXIT 1053 RETURN 1075 1054 ENDIF 1076 1055 END DO 1077 IF( .NOT.found ) THEN 1078 kwgt = nxt_wgt 1079 CALL fld_weight( sd ) 1080 ENDIF 1056 kwgt = nxt_wgt 1057 CALL fld_weight( sd ) 1081 1058 ! 1082 1059 END SUBROUTINE wgt_list … … 1121 1098 TYPE( FLD ), INTENT(in) :: sd ! field with name of weights file 1122 1099 !! 1123 INTEGER :: j n! dummy loop indices1100 INTEGER :: ji,jj,jn ! dummy loop indices 1124 1101 INTEGER :: inum ! local logical unit 1125 1102 INTEGER :: id ! local variable id … … 1127 1104 INTEGER :: zwrap ! local integer 1128 1105 LOGICAL :: cyclical ! 1129 CHARACTER (len=5) :: aname !1130 INTEGER , DIMENSION( :), ALLOCATABLE:: ddims1131 INTEGER , DIMENSION(jpi,jpj) :: data_src1106 CHARACTER (len=5) :: clname ! 1107 INTEGER , DIMENSION(4) :: ddims 1108 INTEGER :: isrc 1132 1109 REAL(wp), DIMENSION(jpi,jpj) :: data_tmp 1133 1110 !!---------------------------------------------------------------------- … … 1142 1119 !! current weights file 1143 1120 1144 !! open input data file (non-model grid) 1145 CALL iom_open( sd%clname, inum, ldiof = LEN(TRIM(sd%wgtname)) > 0 ) 1146 1147 !! get dimensions: we consider 2D data as 3D data with vertical dim size = 1 1148 IF( SIZE(sd%fnow, 3) > 0 ) THEN 1149 ALLOCATE( ddims(4) ) 1150 ELSE 1151 ALLOCATE( ddims(3) ) 1152 ENDIF 1153 id = iom_varid( inum, sd%clvar, ddims ) 1154 1155 !! close it 1156 CALL iom_close( inum ) 1121 !! get data grid dimensions 1122 id = iom_varid( sd%num, sd%clvar, ddims ) 1157 1123 1158 1124 !! now open the weights file 1159 1160 1125 CALL iom_open ( sd%wgtname, inum ) ! interpolation weights 1161 1126 IF( inum > 0 ) THEN … … 1193 1158 !! two possible cases: bilinear (4 weights) or bicubic (16 weights) 1194 1159 id = iom_varid(inum, 'src05', ldstop=.FALSE.) 1195 IF( id <= 0) THEN 1196 ref_wgts(nxt_wgt)%numwgt = 4 1197 ELSE 1198 ref_wgts(nxt_wgt)%numwgt = 16 1199 ENDIF 1200 1201 ALLOCATE( ref_wgts(nxt_wgt)%data_jpi(jpi,jpj,4) ) 1202 ALLOCATE( ref_wgts(nxt_wgt)%data_jpj(jpi,jpj,4) ) 1203 ALLOCATE( ref_wgts(nxt_wgt)%data_wgt(jpi,jpj,ref_wgts(nxt_wgt)%numwgt) ) 1160 IF( id <= 0 ) THEN ; ref_wgts(nxt_wgt)%numwgt = 4 1161 ELSE ; ref_wgts(nxt_wgt)%numwgt = 16 1162 ENDIF 1163 1164 ALLOCATE( ref_wgts(nxt_wgt)%data_jpi(Nis0:Nie0,Njs0:Nje0,4) ) 1165 ALLOCATE( ref_wgts(nxt_wgt)%data_jpj(Nis0:Nie0,Njs0:Nje0,4) ) 1166 ALLOCATE( ref_wgts(nxt_wgt)%data_wgt(Nis0:Nie0,Njs0:Nje0,ref_wgts(nxt_wgt)%numwgt) ) 1204 1167 1205 1168 DO jn = 1,4 1206 aname = ' '1207 WRITE(aname,'(a3,i2.2)') 'src',jn1208 data_tmp(:,:) = 01209 CALL iom_get ( inum, jpdom_data, aname, data_tmp(:,:) )1210 data_src(:,:) = INT(data_tmp(:,:))1211 ref_wgts(nxt_wgt)%data_jpj(:,:,jn) = 1 + (data_src(:,:)-1)/ ref_wgts(nxt_wgt)%ddims(1)1212 ref_wgts(nxt_wgt)%data_jpi(:,:,jn) = data_src(:,:) - ref_wgts(nxt_wgt)%ddims(1)*(ref_wgts(nxt_wgt)%data_jpj(:,:,jn)-1)1169 WRITE(clname,'(a3,i2.2)') 'src',jn 1170 CALL iom_get ( inum, jpdom_global, clname, data_tmp(:,:), cd_type = 'Z' ) ! no call to lbc_lnk 1171 DO_2D( 0, 0, 0, 0 ) 1172 isrc = NINT(data_tmp(ji,jj)) - 1 1173 ref_wgts(nxt_wgt)%data_jpi(ji,jj,jn) = 1 + MOD(isrc, ref_wgts(nxt_wgt)%ddims(1)) 1174 ref_wgts(nxt_wgt)%data_jpj(ji,jj,jn) = 1 + isrc / ref_wgts(nxt_wgt)%ddims(1) 1175 END_2D 1213 1176 END DO 1214 1177 1215 1178 DO jn = 1, ref_wgts(nxt_wgt)%numwgt 1216 aname = ' ' 1217 WRITE(aname,'(a3,i2.2)') 'wgt',jn 1218 ref_wgts(nxt_wgt)%data_wgt(:,:,jn) = 0.0 1219 CALL iom_get ( inum, jpdom_data, aname, ref_wgts(nxt_wgt)%data_wgt(:,:,jn) ) 1179 WRITE(clname,'(a3,i2.2)') 'wgt',jn 1180 CALL iom_get ( inum, jpdom_global, clname, data_tmp(:,:), cd_type = 'Z' ) ! no call to lbc_lnk 1181 DO_2D( 0, 0, 0, 0 ) 1182 ref_wgts(nxt_wgt)%data_wgt(ji,jj,jn) = data_tmp(ji,jj) 1183 END_2D 1220 1184 END DO 1221 1185 CALL iom_close (inum) 1222 1186 1223 1187 ! find min and max indices in grid 1224 ref_wgts(nxt_wgt)%botleft( 1) = MINVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:))1225 ref_wgts(nxt_wgt)%botleft( 2) = MINVAL(ref_wgts(nxt_wgt)%data_jpj(:,:,:))1188 ref_wgts(nxt_wgt)%botleft( 1) = MINVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:)) 1189 ref_wgts(nxt_wgt)%botleft( 2) = MINVAL(ref_wgts(nxt_wgt)%data_jpj(:,:,:)) 1226 1190 ref_wgts(nxt_wgt)%topright(1) = MAXVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:)) 1227 1191 ref_wgts(nxt_wgt)%topright(2) = MAXVAL(ref_wgts(nxt_wgt)%data_jpj(:,:,:)) … … 1247 1211 CALL ctl_stop( ' fld_weight : unable to read the file ' ) 1248 1212 ENDIF 1249 1250 DEALLOCATE (ddims )1251 1213 ! 1252 1214 END SUBROUTINE fld_weight … … 1281 1243 SELECT CASE( SIZE(zfieldo(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:),3) ) 1282 1244 CASE(1) 1283 CALL iom_get( inum, jpdom_unknown, 'LSM', zslmec1(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,1), 1, rec1_lsm, recn_lsm) 1245 CALL iom_get( inum, jpdom_unknown, 'LSM', zslmec1(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,1), & 1246 & 1, kstart = rec1_lsm, kcount = recn_lsm) 1284 1247 CASE DEFAULT 1285 CALL iom_get( inum, jpdom_unknown, 'LSM', zslmec1(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), 1, rec1_lsm, recn_lsm) 1248 CALL iom_get( inum, jpdom_unknown, 'LSM', zslmec1(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), & 1249 & 1, kstart = rec1_lsm, kcount = recn_lsm) 1286 1250 END SELECT 1287 1251 CALL iom_close( inum ) … … 1326 1290 !! D. Delrosso INGV 1327 1291 !!---------------------------------------------------------------------- 1328 INTEGER , INTENT(in ) :: ileni,ilenj ! lengths1329 REAL , DIMENSION (ileni,ilenj), INTENT(in ) :: zfieldn ! array of forcing field with undeff for land points1330 REAL , DIMENSION (ileni,ilenj), INTENT( out) :: zfield ! array of forcing field1331 ! 1332 REAL 1333 REAL 1334 REAL 1335 REAL 1336 LOGICAL , DIMENSION (ileni,ilenj,8) :: ll_msknan3d ! logical mask for undeff detection1337 LOGICAL , DIMENSION (ileni,ilenj) :: ll_msknan2d ! logical mask for undeff detection1292 INTEGER , INTENT(in ) :: ileni,ilenj ! lengths 1293 REAL(wp), DIMENSION (ileni,ilenj), INTENT(in ) :: zfieldn ! array of forcing field with undeff for land points 1294 REAL(wp), DIMENSION (ileni,ilenj), INTENT( out) :: zfield ! array of forcing field 1295 ! 1296 REAL(wp) , DIMENSION (ileni,ilenj) :: zmat1, zmat2, zmat3, zmat4 ! local arrays 1297 REAL(wp) , DIMENSION (ileni,ilenj) :: zmat5, zmat6, zmat7, zmat8 ! - - 1298 REAL(wp) , DIMENSION (ileni,ilenj) :: zlsm2d ! - - 1299 REAL(wp) , DIMENSION (ileni,ilenj,8) :: zlsm3d ! - - 1300 LOGICAL , DIMENSION (ileni,ilenj,8) :: ll_msknan3d ! logical mask for undeff detection 1301 LOGICAL , DIMENSION (ileni,ilenj) :: ll_msknan2d ! logical mask for undeff detection 1338 1302 !!---------------------------------------------------------------------- 1339 1303 zmat8 = eoshift( zfieldn , SHIFT=-1 , BOUNDARY = (/zfieldn(:,1)/) , DIM=2 ) … … 1356 1320 1357 1321 1358 SUBROUTINE fld_interp( num, clvar, kw, kk, dta, & 1359 & nrec, lsmfile) 1322 SUBROUTINE fld_interp( num, clvar, kw, kk, dta, nrec, lsmfile) 1360 1323 !!--------------------------------------------------------------------- 1361 1324 !! *** ROUTINE fld_interp *** … … 1375 1338 INTEGER, DIMENSION(3) :: rec1_lsm, recn_lsm ! temporary arrays for start and length in case of seaoverland 1376 1339 INTEGER :: ii_lsm1,ii_lsm2,ij_lsm1,ij_lsm2 ! temporary indices 1377 INTEGER :: jk, jn, jm, jir, jjr ! loop counters 1340 INTEGER :: ji, jj, jk, jn, jir, jjr ! loop counters 1341 INTEGER :: ipk 1378 1342 INTEGER :: ni, nj ! lengths 1379 1343 INTEGER :: jpimin,jpiwid ! temporary indices … … 1386 1350 REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: ztmp_fly_dta ! local array of values on input grid 1387 1351 !!---------------------------------------------------------------------- 1352 ipk = SIZE(dta, 3) 1388 1353 ! 1389 1354 !! for weighted interpolation we have weights at four corners of a box surrounding … … 1415 1380 1416 1381 1417 IF( LEN ( TRIM(lsmfile)) > 0 ) THEN1382 IF( LEN_TRIM(lsmfile) > 0 ) THEN 1418 1383 !! indeces for ztmp_fly_dta 1419 1384 ! -------------------------- … … 1445 1410 CASE(1) 1446 1411 CALL iom_get( num, jpdom_unknown, clvar, ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,1), & 1447 & nrec, rec1_lsm,recn_lsm)1412 & nrec, kstart = rec1_lsm, kcount = recn_lsm) 1448 1413 CASE DEFAULT 1449 1414 CALL iom_get( num, jpdom_unknown, clvar, ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), & 1450 & nrec, rec1_lsm,recn_lsm)1415 & nrec, kstart = rec1_lsm, kcount = recn_lsm) 1451 1416 END SELECT 1452 1417 CALL apply_seaoverland(lsmfile,ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), & … … 1468 1433 1469 1434 ref_wgts(kw)%fly_dta(:,:,:) = 0.0 1470 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:), nrec, rec1,recn)1435 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:), nrec, kstart = rec1, kcount = recn) 1471 1436 ENDIF 1472 1437 … … 1474 1439 !! first four weights common to both bilinear and bicubic 1475 1440 !! data_jpi, data_jpj have already been shifted to (1,1) corresponding to botleft 1476 !! note that we have to offset by 1 into fly_dta array because of halo 1477 dta(:,:,:) = 0.0 1478 DO jk = 1,4 1479 DO jn = 1, jpj 1480 DO jm = 1,jpi 1481 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1482 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1483 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk) * ref_wgts(kw)%fly_dta(ni+1,nj+1,:) 1484 END DO 1485 END DO 1441 !! note that we have to offset by 1 into fly_dta array because of halo added to fly_dta (rec1 definition) 1442 dta(:,:,:) = 0._wp 1443 DO jn = 1,4 1444 DO_3D( 0, 0, 0, 0, 1,ipk ) 1445 ni = ref_wgts(kw)%data_jpi(ji,jj,jn) + 1 1446 nj = ref_wgts(kw)%data_jpj(ji,jj,jn) + 1 1447 dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn) * ref_wgts(kw)%fly_dta(ni,nj,jk) 1448 END_3D 1486 1449 END DO 1487 1450 1488 1451 IF(ref_wgts(kw)%numwgt .EQ. 16) THEN 1489 1452 1490 !! fix up halo points that we couldnt read from file 1491 IF( jpi1 == 2 ) THEN 1492 ref_wgts(kw)%fly_dta(jpi1-1,:,:) = ref_wgts(kw)%fly_dta(jpi1,:,:) 1493 ENDIF 1494 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 1495 ref_wgts(kw)%fly_dta(jpi2+1,:,:) = ref_wgts(kw)%fly_dta(jpi2,:,:) 1496 ENDIF 1497 IF( jpj1 == 2 ) THEN 1498 ref_wgts(kw)%fly_dta(:,jpj1-1,:) = ref_wgts(kw)%fly_dta(:,jpj1,:) 1499 ENDIF 1500 IF( jpj2 + jpjmin - 1 == ref_wgts(kw)%ddims(2)+1 .AND. jpj2 .lt. jpjwid+2 ) THEN 1501 ref_wgts(kw)%fly_dta(:,jpj2+1,:) = 2.0*ref_wgts(kw)%fly_dta(:,jpj2,:) - ref_wgts(kw)%fly_dta(:,jpj2-1,:) 1502 ENDIF 1503 1504 !! if data grid is cyclic we can do better on east-west edges 1505 !! but have to allow for whether first and last columns are coincident 1506 IF( ref_wgts(kw)%cyclic ) THEN 1507 rec1(2) = MAX( jpjmin-1, 1 ) 1508 recn(1) = 1 1509 recn(2) = MIN( jpjwid+2, ref_wgts(kw)%ddims(2)-rec1(2)+1 ) 1510 jpj1 = 2 + rec1(2) - jpjmin 1511 jpj2 = jpj1 + recn(2) - 1 1512 IF( jpi1 == 2 ) THEN 1513 rec1(1) = ref_wgts(kw)%ddims(1) - ref_wgts(kw)%overlap 1514 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, rec1, recn) 1515 ref_wgts(kw)%fly_dta(jpi1-1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:) 1516 ENDIF 1517 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 1518 rec1(1) = 1 + ref_wgts(kw)%overlap 1519 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, rec1, recn) 1520 ref_wgts(kw)%fly_dta(jpi2+1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:) 1521 ENDIF 1522 ENDIF 1523 1524 ! gradient in the i direction 1525 DO jk = 1,4 1526 DO jn = 1, jpj 1527 DO jm = 1,jpi 1528 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1529 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1530 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+4) * 0.5 * & 1531 (ref_wgts(kw)%fly_dta(ni+2,nj+1,:) - ref_wgts(kw)%fly_dta(ni,nj+1,:)) 1532 END DO 1533 END DO 1534 END DO 1535 1536 ! gradient in the j direction 1537 DO jk = 1,4 1538 DO jn = 1, jpj 1539 DO jm = 1,jpi 1540 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1541 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1542 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+8) * 0.5 * & 1543 (ref_wgts(kw)%fly_dta(ni+1,nj+2,:) - ref_wgts(kw)%fly_dta(ni+1,nj,:)) 1544 END DO 1545 END DO 1546 END DO 1547 1548 ! gradient in the ij direction 1549 DO jk = 1,4 1550 DO jn = 1, jpj 1551 DO jm = 1,jpi 1552 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1553 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1554 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+12) * 0.25 * ( & 1555 (ref_wgts(kw)%fly_dta(ni+2,nj+2,:) - ref_wgts(kw)%fly_dta(ni ,nj+2,:)) - & 1556 (ref_wgts(kw)%fly_dta(ni+2,nj ,:) - ref_wgts(kw)%fly_dta(ni ,nj ,:))) 1557 END DO 1558 END DO 1453 !! fix up halo points that we couldnt read from file 1454 IF( jpi1 == 2 ) THEN 1455 ref_wgts(kw)%fly_dta(jpi1-1,:,:) = ref_wgts(kw)%fly_dta(jpi1,:,:) 1456 ENDIF 1457 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 1458 ref_wgts(kw)%fly_dta(jpi2+1,:,:) = ref_wgts(kw)%fly_dta(jpi2,:,:) 1459 ENDIF 1460 IF( jpj1 == 2 ) THEN 1461 ref_wgts(kw)%fly_dta(:,jpj1-1,:) = ref_wgts(kw)%fly_dta(:,jpj1,:) 1462 ENDIF 1463 IF( jpj2 + jpjmin - 1 == ref_wgts(kw)%ddims(2)+1 .AND. jpj2 .LT. jpjwid+2 ) THEN 1464 ref_wgts(kw)%fly_dta(:,jpj2+1,:) = 2.0*ref_wgts(kw)%fly_dta(:,jpj2,:) - ref_wgts(kw)%fly_dta(:,jpj2-1,:) 1465 ENDIF 1466 1467 !! if data grid is cyclic we can do better on east-west edges 1468 !! but have to allow for whether first and last columns are coincident 1469 IF( ref_wgts(kw)%cyclic ) THEN 1470 rec1(2) = MAX( jpjmin-1, 1 ) 1471 recn(1) = 1 1472 recn(2) = MIN( jpjwid+2, ref_wgts(kw)%ddims(2)-rec1(2)+1 ) 1473 jpj1 = 2 + rec1(2) - jpjmin 1474 jpj2 = jpj1 + recn(2) - 1 1475 IF( jpi1 == 2 ) THEN 1476 rec1(1) = ref_wgts(kw)%ddims(1) - ref_wgts(kw)%overlap 1477 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, kstart = rec1, kcount = recn) 1478 ref_wgts(kw)%fly_dta(jpi1-1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:) 1479 ENDIF 1480 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 1481 rec1(1) = 1 + ref_wgts(kw)%overlap 1482 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, kstart = rec1, kcount = recn) 1483 ref_wgts(kw)%fly_dta(jpi2+1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:) 1484 ENDIF 1485 ENDIF 1486 ! 1487 !!$ DO jn = 1,4 1488 !!$ DO_3D( 0, 0, 0, 0, 1,ipk ) 1489 !!$ ni = ref_wgts(kw)%data_jpi(ji,jj,jn) + 1 1490 !!$ nj = ref_wgts(kw)%data_jpj(ji,jj,jn) + 1 1491 !!$ dta(ji,jj,jk) = dta(ji,jj,jk) & 1492 !!$ ! gradient in the i direction 1493 !!$ & + ref_wgts(kw)%data_wgt(ji,jj,jn+4) * 0.5_wp * & 1494 !!$ & (ref_wgts(kw)%fly_dta(ni+1,nj ,jk) - ref_wgts(kw)%fly_dta(ni-1,nj ,jk)) & 1495 !!$ ! gradient in the j direction 1496 !!$ & + ref_wgts(kw)%data_wgt(ji,jj,jn+8) * 0.5_wp * & 1497 !!$ & (ref_wgts(kw)%fly_dta(ni ,nj+1,jk) - ref_wgts(kw)%fly_dta(ni ,nj-1,jk)) & 1498 !!$ ! gradient in the ij direction 1499 !!$ & + ref_wgts(kw)%data_wgt(ji,jj,jn+12) * 0.25_wp * & 1500 !!$ & ((ref_wgts(kw)%fly_dta(ni+1,nj+1,jk) - ref_wgts(kw)%fly_dta(ni-1,nj+1,jk)) - & 1501 !!$ & (ref_wgts(kw)%fly_dta(ni+1,nj-1,jk) - ref_wgts(kw)%fly_dta(ni-1,nj-1,jk))) 1502 !!$ END_3D 1503 !!$ END DO 1504 ! 1505 DO jn = 1,4 1506 DO_3D( 0, 0, 0, 0, 1,ipk ) 1507 ni = ref_wgts(kw)%data_jpi(ji,jj,jn) 1508 nj = ref_wgts(kw)%data_jpj(ji,jj,jn) 1509 ! gradient in the i direction 1510 dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+4) * 0.5_wp * & 1511 & (ref_wgts(kw)%fly_dta(ni+2,nj+1,jk) - ref_wgts(kw)%fly_dta(ni ,nj+1,jk)) 1512 END_3D 1513 END DO 1514 DO jn = 1,4 1515 DO_3D( 0, 0, 0, 0, 1,ipk ) 1516 ni = ref_wgts(kw)%data_jpi(ji,jj,jn) 1517 nj = ref_wgts(kw)%data_jpj(ji,jj,jn) 1518 ! gradient in the j direction 1519 dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+8) * 0.5_wp * & 1520 & (ref_wgts(kw)%fly_dta(ni+1,nj+2,jk) - ref_wgts(kw)%fly_dta(ni+1,nj ,jk)) 1521 END_3D 1522 END DO 1523 DO jn = 1,4 1524 DO_3D( 0, 0, 0, 0, 1,ipk ) 1525 ni = ref_wgts(kw)%data_jpi(ji,jj,jn) 1526 nj = ref_wgts(kw)%data_jpj(ji,jj,jn) 1527 ! gradient in the ij direction 1528 dta(ji,jj,jk) = dta(ji,jj,jk) + ref_wgts(kw)%data_wgt(ji,jj,jn+12) * 0.25_wp * ( & 1529 & (ref_wgts(kw)%fly_dta(ni+2,nj+2,jk) - ref_wgts(kw)%fly_dta(ni ,nj+2,jk)) - & 1530 & (ref_wgts(kw)%fly_dta(ni+2,nj ,jk) - ref_wgts(kw)%fly_dta(ni ,nj ,jk))) 1531 END_3D 1559 1532 END DO 1560 1533 ! … … 1583 1556 IF( .NOT. sdjf%ln_clim ) THEN 1584 1557 WRITE(clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), kyear ! add year 1585 IF( sdjf%cl type/= 'yearly' ) WRITE(clname, '(a, "m",i2.2)' ) TRIM( clname ), kmonth ! add month1558 IF( sdjf%clftyp /= 'yearly' ) WRITE(clname, '(a, "m",i2.2)' ) TRIM( clname ), kmonth ! add month 1586 1559 ELSE 1587 1560 ! build the new filename if climatological data 1588 IF( sdjf%cl type/= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( sdjf%clrootname ), kmonth ! add month1589 ENDIF 1590 IF( sdjf%cl type == 'daily' .OR. sdjf%cltype(1:4) == 'week' ) &1561 IF( sdjf%clftyp /= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( sdjf%clrootname ), kmonth ! add month 1562 ENDIF 1563 IF( sdjf%clftyp == 'daily' .OR. sdjf%clftyp(1:4) == 'week' ) & 1591 1564 & WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname ), kday ! add day 1592 1565 … … 1612 1585 IF( cl_week(ijul) == TRIM(cdday) ) EXIT 1613 1586 END DO 1614 IF( ijul .GT. 7 ) CALL ctl_stop( 'ksec_week: wrong day for sdjf%cl type(6:8): '//TRIM(cdday) )1587 IF( ijul .GT. 7 ) CALL ctl_stop( 'ksec_week: wrong day for sdjf%clftyp(6:8): '//TRIM(cdday) ) 1615 1588 ! 1616 1589 ishift = ijul * NINT(rday)
Note: See TracChangeset
for help on using the changeset viewer.
