- Timestamp:
- 2020-11-27T17:26:33+01:00 (4 years ago)
- Location:
- NEMO/branches/2020/tickets_icb_1900
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/SBC/fldread.F90 (modified) (47 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/tickets_icb_1900
- Property svn:externals
-
NEMO/branches/2020/tickets_icb_1900/src/OCE/SBC/fldread.F90
r13237 r13899 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 … … 157 160 INTEGER :: jf ! dummy indices 158 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 159 163 LOGICAL :: ll_firstcall ! true if this is the first call to fld_read for this set of fields 160 164 REAL(wp) :: zt_offset ! local time offset variable … … 204 208 IF( TRIM(sd(jf)%clrootname) == 'NOT USED' ) CYCLE 205 209 ! 210 ibb = sd(jf)%nbb ; iaa = sd(jf)%naa 211 ! 206 212 IF( sd(jf)%ln_tint ) THEN ! temporal interpolation 207 213 IF(lwp .AND. kt - nit000 <= 100 ) THEN … … 209 215 & "', records b/a: ', i6.4, '/', i6.4, ' (days ', f9.4,'/', f9.4, ')')" 210 216 WRITE(numout, clfmt) TRIM( sd(jf)%clvar ), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 211 & 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)/rday212 WRITE(numout, *) ' zt_offset is : ',zt_offset217 & 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 218 IF( zt_offset /= 0._wp ) WRITE(numout, *) ' zt_offset is : ', zt_offset 213 219 ENDIF 214 220 ! temporal interpolation weights 215 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 ) 216 222 ztintb = 1. - ztinta 217 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) 218 224 ELSE ! nothing to do... 219 225 IF(lwp .AND. kt - nit000 <= 100 ) THEN … … 221 227 & "', record: ', i6.4, ' (days ', f9.4, ' <-> ', f9.4, ')')" 222 228 WRITE(numout, clfmt) TRIM(sd(jf)%clvar), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 223 & 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 224 230 ENDIF 225 231 ENDIF … … 251 257 ! 252 258 CALL fld_clopn( sdjf ) 253 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. 254 260 ! 255 261 END SUBROUTINE fld_init … … 262 268 !! ** Purpose : Compute 263 269 !! if sdjf%ln_tint = .TRUE. 264 !! 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) 265 271 !! if sdjf%ln_tint = .FALSE. 266 !! nrec _a(1): record number267 !! 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 268 274 !!---------------------------------------------------------------------- 269 275 INTEGER , INTENT(in ) :: ksecsbc ! … … 271 277 INTEGER , OPTIONAL, INTENT(in ) :: Kmm ! ocean time level index 272 278 ! 273 INTEGER :: ja ! end of this record (in seconds) 274 !!---------------------------------------------------------------------- 275 ! 276 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 277 285 278 ! find where is the new after record... (it is not necessary sdjf%nrec _a(1)+1 )279 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) 280 288 DO WHILE ( ksecsbc >= sdjf%nrecsec(ja) .AND. ja < sdjf%nreclast ) ! Warning: make sure ja <= sdjf%nreclast in this test 281 289 ja = ja + 1 … … 284 292 285 293 ! if ln_tint and if the new after is not ja+1, we need also to update after data before the swap 286 ! so, after the swap, sdjf%nrec _b(2) will still be the closest value located just before ksecsbc287 IF( sdjf%ln_tint .AND. ( ja > sdjf%nrec _a(1) + 1 .OR. sdjf%nrec_a(2) == nflag ) ) THEN288 sdjf%nrec _a(:) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! update nrec_awith before information289 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 290 298 ENDIF 291 299 … … 310 318 ! if ln_tint and if after is not the first record, we must (potentially again) update after data before the swap 311 319 IF( sdjf%ln_tint .AND. ja > 1 ) THEN 312 IF( sdjf%nrecsec(0) /= nflag ) THEN ! no trick used: after file is not the current file313 sdjf%nrec _a(:) = (/ ja-1, sdjf%nrecsec(ja-1) /) ! update nrec_awith before information314 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 315 323 ENDIF 316 324 ENDIF … … 318 326 ENDIF 319 327 320 IF( sdjf%ln_tint ) THEN 321 ! Swap data 322 sdjf%nrec_b(:) = sdjf%nrec_a(:) ! swap before record informations 323 sdjf%rotn(1) = sdjf%rotn(2) ! swap before rotate informations 324 sdjf%fdta(:,:,:,1) = sdjf%fdta(:,:,:,2) ! swap before record field 325 ELSE 326 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 327 333 ENDIF 328 334 329 335 ! read new after data 330 sdjf%nrec _a(:) = (/ ja, sdjf%nrecsec(ja) /) ! update nrec_aas it is used by fld_get331 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) 332 338 333 339 ENDIF … … 346 352 ! 347 353 INTEGER :: ipk ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 354 INTEGER :: iaa ! shorter name for sdjf%naa 348 355 INTEGER :: iw ! index into wgts array 349 INTEGER :: ipdom ! index of the domain350 356 INTEGER :: idvar ! variable ID 351 357 INTEGER :: idmspc ! number of spatial dimensions 352 358 LOGICAL :: lmoor ! C1D case: point data 353 !!--------------------------------------------------------------------- 354 ! 355 ipk = SIZE( sdjf%fnow, 3 ) 356 ! 357 IF( ASSOCIATED(sdjf%imap) ) THEN 358 IF( sdjf%ln_tint ) THEN ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), & 359 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint, Kmm ) 360 ELSE ; CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1), & 361 & sdjf%imap, sdjf%igrd, sdjf%ibdy, sdjf%ltotvel, sdjf%lzint, Kmm ) 362 ENDIF 363 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 364 372 CALL wgt_list( sdjf, iw ) 365 IF( sdjf%ln_tint ) THEN ; CALL fld_interp( sdjf%num, sdjf%clvar, iw, ipk, sdjf%fdta(:,:,:,2), & 366 & sdjf%nrec_a(1), sdjf%lsmname ) 367 ELSE ; CALL fld_interp( sdjf%num, sdjf%clvar, iw, ipk, sdjf%fnow(:,:,: ), & 368 & sdjf%nrec_a(1), sdjf%lsmname ) 369 ENDIF 370 ELSE 371 IF( SIZE(sdjf%fnow, 1) == jpi ) THEN ; ipdom = jpdom_data 372 ELSE ; ipdom = jpdom_unknown 373 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 374 376 ! C1D case: If product of spatial dimensions == ipk, then x,y are of 375 377 ! size 1 (point/mooring data): this must be read onto the central grid point 376 378 idvar = iom_varid( sdjf%num, sdjf%clvar ) 377 379 idmspc = iom_file ( sdjf%num )%ndims( idvar ) 378 IF( iom_file( sdjf%num )%luld( idvar ) ) idmspc = idmspc - 1 379 lmoor = ( idmspc == 0 .OR. PRODUCT( iom_file( sdjf%num )%dimsz( 1:MAX(idmspc,1) ,idvar ) ) == ipk ) 380 ! 381 SELECT CASE( ipk ) 382 CASE(1) 383 IF( lk_c1d .AND. lmoor ) THEN 384 IF( sdjf%ln_tint ) THEN 385 CALL iom_get( sdjf%num, sdjf%clvar, sdjf%fdta(2,2,1,2), sdjf%nrec_a(1) ) 386 CALL lbc_lnk( 'fldread', sdjf%fdta(:,:,1,2),'Z',1.0_wp ) 387 ELSE 388 CALL iom_get( sdjf%num, sdjf%clvar, sdjf%fnow(2,2,1 ), sdjf%nrec_a(1) ) 389 CALL lbc_lnk( 'fldread', sdjf%fnow(:,:,1 ),'Z',1.0_wp ) 390 ENDIF 391 ELSE 392 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_a(1) ) 393 ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,1 ), sdjf%nrec_a(1) ) 394 ENDIF 395 ENDIF 396 CASE DEFAULT 397 IF(lk_c1d .AND. lmoor ) THEN 398 IF( sdjf%ln_tint ) THEN 399 CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, sdjf%fdta(2,2,:,2), sdjf%nrec_a(1) ) 400 CALL lbc_lnk( 'fldread', sdjf%fdta(:,:,:,2),'Z',1.0_wp ) 401 ELSE 402 CALL iom_get( sdjf%num, jpdom_unknown, sdjf%clvar, sdjf%fnow(2,2,: ), sdjf%nrec_a(1) ) 403 CALL lbc_lnk( 'fldread', sdjf%fnow(:,:,: ),'Z',1.0_wp ) 404 ENDIF 405 ELSE 406 IF( sdjf%ln_tint ) THEN ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1) ) 407 ELSE ; CALL iom_get( sdjf%num, ipdom, sdjf%clvar, sdjf%fnow(:,:,: ), sdjf%nrec_a(1) ) 408 ENDIF 409 ENDIF 410 END SELECT 411 ENDIF 412 ! 413 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 414 393 ! 415 394 END SUBROUTINE fld_get … … 447 426 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read_z ! work space local data requiring vertical interpolation 448 427 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zdta_read_dz ! work space local data requiring vertical interpolation 449 CHARACTER(LEN=1),DIMENSION(3) :: cl grid428 CHARACTER(LEN=1),DIMENSION(3) :: cltype 450 429 LOGICAL :: lluld ! is the variable using the unlimited dimension 451 430 LOGICAL :: llzint ! local value of ldzint 452 431 !!--------------------------------------------------------------------- 453 432 ! 454 cl grid= (/'t','u','v'/)433 cltype = (/'t','u','v'/) 455 434 ! 456 435 ipi = SIZE( pdta, 1 ) … … 487 466 IF( ipkb /= ipk .OR. llzint ) THEN ! boundary data not on model vertical grid : vertical interpolation 488 467 ! 489 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 490 469 491 470 ALLOCATE( zdta_read(ipi,ipj,ipkb), zdta_read_z(ipi,ipj,ipkb), zdta_read_dz(ipi,ipj,ipkb) ) 492 471 493 472 CALL fld_map_core( zz_read, kmap, zdta_read ) 494 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? 495 474 CALL fld_map_core( zz_read, kmap, zdta_read_z ) 496 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? 497 476 CALL fld_map_core( zz_read, kmap, zdta_read_dz ) 498 477 … … 504 483 IF( ipk /= jpk ) CALL ctl_stop( 'fld_map : this should be an impossible case...' ) 505 484 WRITE(ctmp1,*) 'fld_map : vertical interpolation for bdy variable '//TRIM(cdvar)//' requires ' 506 IF( iom_varid(knum, 'gdep'//cl grid(kgrd)) == -1 ) CALL ctl_stop( ctmp1//'gdep'//clgrid(kgrd)//' variable' )507 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' ) 508 487 509 488 ENDIF … … 728 707 CHARACTER (LEN=100) :: clcomp ! dummy weight name 729 708 REAL(wp), DIMENSION(jpi,jpj) :: utmp, vtmp ! temporary arrays for vector rotation 709 REAL(wp), DIMENSION(:,:,:), POINTER :: dta_u, dta_v ! short cut 730 710 !!--------------------------------------------------------------------- 731 711 ! … … 747 727 END DO 748 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 749 732 DO jk = 1, SIZE( sd(ju)%fnow, 3 ) 750 IF( sd(ju)%ln_tint )THEN 751 CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->i', utmp(:,:) ) 752 CALL rot_rep( sd(ju)%fdta(:,:,jk,jn), sd(iv)%fdta(:,:,jk,jn), 'T', 'en->j', vtmp(:,:) ) 753 sd(ju)%fdta(:,:,jk,jn) = utmp(:,:) ; sd(iv)%fdta(:,:,jk,jn) = vtmp(:,:) 754 ELSE 755 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->i', utmp(:,:) ) 756 CALL rot_rep( sd(ju)%fnow(:,:,jk ), sd(iv)%fnow(:,:,jk ), 'T', 'en->j', vtmp(:,:) ) 757 sd(ju)%fnow(:,:,jk ) = utmp(:,:) ; sd(iv)%fnow(:,:,jk ) = vtmp(:,:) 758 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(:,:) 759 736 END DO 760 737 sd(ju)%rotn(jn) = .TRUE. ! vector was rotated … … 802 779 803 780 ! current file parameters 804 IF( sdjf%cl type(1:4) == 'week' ) THEN! find the day of the beginning of the current week805 isecwk = ksec_week( sdjf%cl type(6:8) ) ! seconds between the beginning of the week and half of current time step806 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 807 784 llprevyr = llprevmt .AND. nmonth == 1 808 785 iyr = nyear - COUNT((/llprevyr/)) 809 786 imt = nmonth - COUNT((/llprevmt/)) + 12 * COUNT((/llprevyr/)) 810 787 idy = nday + nmonth_len(nmonth-1) * COUNT((/llprevmt/)) - isecwk / idaysec 811 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 812 789 ELSE 813 790 iyr = nyear … … 819 796 ! previous file parameters 820 797 IF( llprev ) THEN 821 IF( sdjf%cl type(1:4) == 'week' ) THEN! find the day of the beginning of previous week822 isecwk = isecwk + 7 * idaysec ! seconds between the beginning of previous week and half of the time step823 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 824 801 llprevyr = llprevmt .AND. nmonth == 1 825 802 iyr = nyear - COUNT((/llprevyr/)) 826 803 imt = nmonth - COUNT((/llprevmt/)) + 12 * COUNT((/llprevyr/)) 827 804 idy = nday + nmonth_len(nmonth-1) * COUNT((/llprevmt/)) - isecwk / idaysec 828 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 829 806 ELSE 830 idy = nday - COUNT((/ sdjf%cl type== 'daily' /))831 imt = nmonth - COUNT((/ sdjf%cl type== 'monthly' .OR. idy == 0 /))832 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 /)) 833 810 IF( idy == 0 ) idy = nmonth_len(imt) 834 811 IF( imt == 0 ) imt = 12 … … 839 816 ! next file parameters 840 817 IF( llnext ) THEN 841 IF( sdjf%cl type(1:4) == 'week' ) THEN! find the day of the beginning of next week842 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 843 820 llnextmt = isecwk > ( nmonth_len(nmonth)*idaysec - nsec_month ) ! larger than the seconds to the end of the month 844 821 llnextyr = llnextmt .AND. nmonth == 12 … … 846 823 imt = nmonth + COUNT((/llnextmt/)) - 12 * COUNT((/llnextyr/)) 847 824 idy = nday - nmonth_len(nmonth) * COUNT((/llnextmt/)) + isecwk / idaysec + 1 848 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 849 826 ELSE 850 idy = nday + COUNT((/ sdjf%cl type== 'daily' /))851 imt = nmonth + COUNT((/ sdjf%cl type== 'monthly' .OR. idy > nmonth_len(nmonth) /))852 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 /)) 853 830 IF( idy > nmonth_len(nmonth) ) idy = 1 854 831 IF( imt == 13 ) imt = 1 … … 867 844 IF ( NINT(sdjf%freqh) == -12 ) THEN ; ireclast = 1 ! yearly mean: consider only 1 record 868 845 ELSEIF( NINT(sdjf%freqh) == -1 ) THEN ! monthly mean: 869 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 870 847 ELSE ; ireclast = 12 ! consider that the file has 12 record 871 848 ENDIF 872 849 ELSE ! higher frequency mean (in hours) 873 IF( sdjf%cl type== 'monthly' ) THEN ; ireclast = NINT( 24. * REAL(nmonth_len(indexmt), wp) / sdjf%freqh )874 ELSEIF( sdjf%cl type(1:4) == 'week' ) THEN ; ireclast = NINT( 24. * 7. / sdjf%freqh )875 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 ) 876 853 ELSE ; ireclast = NINT( 24. * REAL( nyear_len(indexyr), wp) / sdjf%freqh ) 877 854 ENDIF … … 891 868 sdjf%nrecsec(1) = sdjf%nrecsec(0) + nyear_len( indexyr ) * idaysec 892 869 ELSEIF( NINT(sdjf%freqh) == -1 ) THEN ! monthly mean: 893 IF( sdjf%cl type== 'monthly' ) THEN ! monthly file870 IF( sdjf%clftyp == 'monthly' ) THEN ! monthly file 894 871 sdjf%nrecsec(0 ) = nsec1jan000 + nmonth_beg(indexmt ) 895 872 sdjf%nrecsec(1 ) = nsec1jan000 + nmonth_beg(indexmt+1) … … 899 876 ENDIF 900 877 ELSE ! higher frequency mean (in hours) 901 IF( sdjf%cl type== 'monthly' ) THEN ; istart = nsec1jan000 + nmonth_beg(indexmt)902 ELSEIF( sdjf%cl type(1:4) == 'week' ) THEN ; istart = nsec1jan000 + isecwk903 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 904 881 ELSEIF( indexyr == 0 ) THEN ; istart = nsec1jan000 - nyear_len( 0 ) * idaysec 905 882 ELSEIF( indexyr == 2 ) THEN ; istart = nsec1jan000 + nyear_len( 1 ) * idaysec … … 942 919 IF( sdjf%num <= 0 .OR. .NOT. sdjf%ln_clim ) THEN 943 920 IF( sdjf%num > 0 ) CALL iom_close( sdjf%num ) ! close file if already open 944 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 ) 945 922 ENDIF 946 923 ! … … 964 941 ENDIF 965 942 ! 966 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 ) 967 944 ! 968 945 ENDIF … … 997 974 sdf(jf)%ln_tint = sdf_n(jf)%ln_tint 998 975 sdf(jf)%ln_clim = sdf_n(jf)%ln_clim 999 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 1000 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 1001 982 sdf(jf)%wgtname = " " 1002 983 IF( LEN( TRIM(sdf_n(jf)%wname) ) > 0 ) sdf(jf)%wgtname = TRIM( cdir )//sdf_n(jf)%wname … … 1005 986 sdf(jf)%vcomp = sdf_n(jf)%vcomp 1006 987 sdf(jf)%rotn(:) = .TRUE. ! pretend to be rotated -> won't try to rotate data before the first call to fld_get 1007 IF( sdf(jf)%cl type(1:4) == 'week' .AND. nn_leapy == 0 ) &988 IF( sdf(jf)%clftyp(1:4) == 'week' .AND. nn_leapy == 0 ) & 1008 989 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs nn_leapy = 1') 1009 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 ) & 1010 991 & CALL ctl_stop('fld_clopn: weekly file ('//TRIM(sdf(jf)%clrootname)//') needs ln_clim = .FALSE.') 1011 992 sdf(jf)%nreclast = -1 ! Set to non zero default value to avoid errors, is updated to meaningful value during fld_clopn … … 1033 1014 WRITE(numout,*) ' weights: ' , TRIM( sdf(jf)%wgtname ), & 1034 1015 & ' pairing: ' , TRIM( sdf(jf)%vcomp ), & 1035 & ' data type: ' , sdf(jf)%cl type, &1016 & ' data type: ' , sdf(jf)%clftyp , & 1036 1017 & ' land/sea mask:' , TRIM( sdf(jf)%lsmname ) 1037 1018 call flush(numout) … … 1051 1032 !!---------------------------------------------------------------------- 1052 1033 TYPE( FLD ), INTENT(in ) :: sd ! field with name of weights file 1053 INTEGER , INTENT( inout) :: kwgt ! index of weights1034 INTEGER , INTENT( out) :: kwgt ! index of weights 1054 1035 ! 1055 1036 INTEGER :: kw, nestid ! local integer 1056 LOGICAL :: found ! local logical1057 1037 !!---------------------------------------------------------------------- 1058 1038 ! 1059 1039 !! search down linked list 1060 1040 !! weights filename is either present or we hit the end of the list 1061 found = .FALSE.1062 1041 ! 1063 1042 !! because agrif nest part of filenames are now added in iom_open … … 1069 1048 #endif 1070 1049 DO kw = 1, nxt_wgt-1 1071 IF( TRIM(ref_wgts(kw)%wgtname) == TRIM(sd%wgtname).AND. &1072 ref_wgts(kw)%nestid == nestid) THEN1050 IF( ref_wgts(kw)%wgtname == sd%wgtname .AND. & 1051 ref_wgts(kw)%nestid == nestid) THEN 1073 1052 kwgt = kw 1074 found = .TRUE. 1075 EXIT 1053 RETURN 1076 1054 ENDIF 1077 1055 END DO 1078 IF( .NOT.found ) THEN 1079 kwgt = nxt_wgt 1080 CALL fld_weight( sd ) 1081 ENDIF 1056 kwgt = nxt_wgt 1057 CALL fld_weight( sd ) 1082 1058 ! 1083 1059 END SUBROUTINE wgt_list … … 1122 1098 TYPE( FLD ), INTENT(in) :: sd ! field with name of weights file 1123 1099 !! 1124 INTEGER :: j n! dummy loop indices1100 INTEGER :: ji,jj,jn ! dummy loop indices 1125 1101 INTEGER :: inum ! local logical unit 1126 1102 INTEGER :: id ! local variable id … … 1128 1104 INTEGER :: zwrap ! local integer 1129 1105 LOGICAL :: cyclical ! 1130 CHARACTER (len=5) :: aname !1131 INTEGER , DIMENSION( :), ALLOCATABLE:: ddims1132 INTEGER , DIMENSION(jpi,jpj) :: data_src1106 CHARACTER (len=5) :: clname ! 1107 INTEGER , DIMENSION(4) :: ddims 1108 INTEGER :: isrc 1133 1109 REAL(wp), DIMENSION(jpi,jpj) :: data_tmp 1134 1110 !!---------------------------------------------------------------------- … … 1143 1119 !! current weights file 1144 1120 1145 !! open input data file (non-model grid) 1146 CALL iom_open( sd%clname, inum, ldiof = LEN(TRIM(sd%wgtname)) > 0 ) 1147 1148 !! get dimensions: we consider 2D data as 3D data with vertical dim size = 1 1149 IF( SIZE(sd%fnow, 3) > 0 ) THEN 1150 ALLOCATE( ddims(4) ) 1151 ELSE 1152 ALLOCATE( ddims(3) ) 1153 ENDIF 1154 id = iom_varid( inum, sd%clvar, ddims ) 1155 1156 !! close it 1157 CALL iom_close( inum ) 1121 !! get data grid dimensions 1122 id = iom_varid( sd%num, sd%clvar, ddims ) 1158 1123 1159 1124 !! now open the weights file 1160 1161 1125 CALL iom_open ( sd%wgtname, inum ) ! interpolation weights 1162 1126 IF( inum > 0 ) THEN … … 1194 1158 !! two possible cases: bilinear (4 weights) or bicubic (16 weights) 1195 1159 id = iom_varid(inum, 'src05', ldstop=.FALSE.) 1196 IF( id <= 0) THEN 1197 ref_wgts(nxt_wgt)%numwgt = 4 1198 ELSE 1199 ref_wgts(nxt_wgt)%numwgt = 16 1200 ENDIF 1201 1202 ALLOCATE( ref_wgts(nxt_wgt)%data_jpi(jpi,jpj,4) ) 1203 ALLOCATE( ref_wgts(nxt_wgt)%data_jpj(jpi,jpj,4) ) 1204 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) ) 1205 1167 1206 1168 DO jn = 1,4 1207 aname = ' '1208 WRITE(aname,'(a3,i2.2)') 'src',jn1209 data_tmp(:,:) = 01210 CALL iom_get ( inum, jpdom_data, aname, data_tmp(:,:) )1211 data_src(:,:) = INT(data_tmp(:,:))1212 ref_wgts(nxt_wgt)%data_jpj(:,:,jn) = 1 + (data_src(:,:)-1)/ ref_wgts(nxt_wgt)%ddims(1)1213 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 1214 1176 END DO 1215 1177 1216 1178 DO jn = 1, ref_wgts(nxt_wgt)%numwgt 1217 aname = ' ' 1218 WRITE(aname,'(a3,i2.2)') 'wgt',jn 1219 ref_wgts(nxt_wgt)%data_wgt(:,:,jn) = 0.0 1220 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 1221 1184 END DO 1222 1185 CALL iom_close (inum) 1223 1186 1224 1187 ! find min and max indices in grid 1225 ref_wgts(nxt_wgt)%botleft( 1) = MINVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:))1226 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(:,:,:)) 1227 1190 ref_wgts(nxt_wgt)%topright(1) = MAXVAL(ref_wgts(nxt_wgt)%data_jpi(:,:,:)) 1228 1191 ref_wgts(nxt_wgt)%topright(2) = MAXVAL(ref_wgts(nxt_wgt)%data_jpj(:,:,:)) … … 1248 1211 CALL ctl_stop( ' fld_weight : unable to read the file ' ) 1249 1212 ENDIF 1250 1251 DEALLOCATE (ddims )1252 1213 ! 1253 1214 END SUBROUTINE fld_weight … … 1282 1243 SELECT CASE( SIZE(zfieldo(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:),3) ) 1283 1244 CASE(1) 1284 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) 1285 1247 CASE DEFAULT 1286 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) 1287 1250 END SELECT 1288 1251 CALL iom_close( inum ) … … 1357 1320 1358 1321 1359 SUBROUTINE fld_interp( num, clvar, kw, kk, dta, & 1360 & nrec, lsmfile) 1322 SUBROUTINE fld_interp( num, clvar, kw, kk, dta, nrec, lsmfile) 1361 1323 !!--------------------------------------------------------------------- 1362 1324 !! *** ROUTINE fld_interp *** … … 1376 1338 INTEGER, DIMENSION(3) :: rec1_lsm, recn_lsm ! temporary arrays for start and length in case of seaoverland 1377 1339 INTEGER :: ii_lsm1,ii_lsm2,ij_lsm1,ij_lsm2 ! temporary indices 1378 INTEGER :: jk, jn, jm, jir, jjr ! loop counters 1340 INTEGER :: ji, jj, jk, jn, jir, jjr ! loop counters 1341 INTEGER :: ipk 1379 1342 INTEGER :: ni, nj ! lengths 1380 1343 INTEGER :: jpimin,jpiwid ! temporary indices … … 1387 1350 REAL(wp),DIMENSION(:,:,:), ALLOCATABLE :: ztmp_fly_dta ! local array of values on input grid 1388 1351 !!---------------------------------------------------------------------- 1352 ipk = SIZE(dta, 3) 1389 1353 ! 1390 1354 !! for weighted interpolation we have weights at four corners of a box surrounding … … 1416 1380 1417 1381 1418 IF( LEN ( TRIM(lsmfile)) > 0 ) THEN1382 IF( LEN_TRIM(lsmfile) > 0 ) THEN 1419 1383 !! indeces for ztmp_fly_dta 1420 1384 ! -------------------------- … … 1446 1410 CASE(1) 1447 1411 CALL iom_get( num, jpdom_unknown, clvar, ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,1), & 1448 & nrec, rec1_lsm,recn_lsm)1412 & nrec, kstart = rec1_lsm, kcount = recn_lsm) 1449 1413 CASE DEFAULT 1450 1414 CALL iom_get( num, jpdom_unknown, clvar, ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), & 1451 & nrec, rec1_lsm,recn_lsm)1415 & nrec, kstart = rec1_lsm, kcount = recn_lsm) 1452 1416 END SELECT 1453 1417 CALL apply_seaoverland(lsmfile,ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:), & … … 1469 1433 1470 1434 ref_wgts(kw)%fly_dta(:,:,:) = 0.0 1471 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) 1472 1436 ENDIF 1473 1437 … … 1475 1439 !! first four weights common to both bilinear and bicubic 1476 1440 !! data_jpi, data_jpj have already been shifted to (1,1) corresponding to botleft 1477 !! note that we have to offset by 1 into fly_dta array because of halo 1478 dta(:,:,:) = 0.0 1479 DO jk = 1,4 1480 DO jn = 1, jpj 1481 DO jm = 1,jpi 1482 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1483 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1484 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk) * ref_wgts(kw)%fly_dta(ni+1,nj+1,:) 1485 END DO 1486 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 1487 1449 END DO 1488 1450 1489 1451 IF(ref_wgts(kw)%numwgt .EQ. 16) THEN 1490 1452 1491 !! fix up halo points that we couldnt read from file 1492 IF( jpi1 == 2 ) THEN 1493 ref_wgts(kw)%fly_dta(jpi1-1,:,:) = ref_wgts(kw)%fly_dta(jpi1,:,:) 1494 ENDIF 1495 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 1496 ref_wgts(kw)%fly_dta(jpi2+1,:,:) = ref_wgts(kw)%fly_dta(jpi2,:,:) 1497 ENDIF 1498 IF( jpj1 == 2 ) THEN 1499 ref_wgts(kw)%fly_dta(:,jpj1-1,:) = ref_wgts(kw)%fly_dta(:,jpj1,:) 1500 ENDIF 1501 IF( jpj2 + jpjmin - 1 == ref_wgts(kw)%ddims(2)+1 .AND. jpj2 .lt. jpjwid+2 ) THEN 1502 ref_wgts(kw)%fly_dta(:,jpj2+1,:) = 2.0*ref_wgts(kw)%fly_dta(:,jpj2,:) - ref_wgts(kw)%fly_dta(:,jpj2-1,:) 1503 ENDIF 1504 1505 !! if data grid is cyclic we can do better on east-west edges 1506 !! but have to allow for whether first and last columns are coincident 1507 IF( ref_wgts(kw)%cyclic ) THEN 1508 rec1(2) = MAX( jpjmin-1, 1 ) 1509 recn(1) = 1 1510 recn(2) = MIN( jpjwid+2, ref_wgts(kw)%ddims(2)-rec1(2)+1 ) 1511 jpj1 = 2 + rec1(2) - jpjmin 1512 jpj2 = jpj1 + recn(2) - 1 1513 IF( jpi1 == 2 ) THEN 1514 rec1(1) = ref_wgts(kw)%ddims(1) - ref_wgts(kw)%overlap 1515 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, rec1, recn) 1516 ref_wgts(kw)%fly_dta(jpi1-1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:) 1517 ENDIF 1518 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 1519 rec1(1) = 1 + ref_wgts(kw)%overlap 1520 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col(:,jpj1:jpj2,:), nrec, rec1, recn) 1521 ref_wgts(kw)%fly_dta(jpi2+1,jpj1:jpj2,:) = ref_wgts(kw)%col(1,jpj1:jpj2,:) 1522 ENDIF 1523 ENDIF 1524 1525 ! gradient in the i direction 1526 DO jk = 1,4 1527 DO jn = 1, jpj 1528 DO jm = 1,jpi 1529 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1530 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1531 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+4) * 0.5 * & 1532 (ref_wgts(kw)%fly_dta(ni+2,nj+1,:) - ref_wgts(kw)%fly_dta(ni,nj+1,:)) 1533 END DO 1534 END DO 1535 END DO 1536 1537 ! gradient in the j direction 1538 DO jk = 1,4 1539 DO jn = 1, jpj 1540 DO jm = 1,jpi 1541 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1542 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1543 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+8) * 0.5 * & 1544 (ref_wgts(kw)%fly_dta(ni+1,nj+2,:) - ref_wgts(kw)%fly_dta(ni+1,nj,:)) 1545 END DO 1546 END DO 1547 END DO 1548 1549 ! gradient in the ij direction 1550 DO jk = 1,4 1551 DO jn = 1, jpj 1552 DO jm = 1,jpi 1553 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1554 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1555 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+12) * 0.25 * ( & 1556 (ref_wgts(kw)%fly_dta(ni+2,nj+2,:) - ref_wgts(kw)%fly_dta(ni ,nj+2,:)) - & 1557 (ref_wgts(kw)%fly_dta(ni+2,nj ,:) - ref_wgts(kw)%fly_dta(ni ,nj ,:))) 1558 END DO 1559 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 1560 1532 END DO 1561 1533 ! … … 1584 1556 IF( .NOT. sdjf%ln_clim ) THEN 1585 1557 WRITE(clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), kyear ! add year 1586 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 1587 1559 ELSE 1588 1560 ! build the new filename if climatological data 1589 IF( sdjf%cl type/= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( sdjf%clrootname ), kmonth ! add month1590 ENDIF 1591 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' ) & 1592 1564 & WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname ), kday ! add day 1593 1565 … … 1613 1585 IF( cl_week(ijul) == TRIM(cdday) ) EXIT 1614 1586 END DO 1615 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) ) 1616 1588 ! 1617 1589 ishift = ijul * NINT(rday)
Note: See TracChangeset
for help on using the changeset viewer.
