- Timestamp:
- 2010-10-12T20:49:32+02:00 (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/fldread.F90
r2004 r2236 15 15 USE oce ! ocean dynamics and tracers 16 16 USE dom_oce ! ocean space and time domain 17 USE ioipsl, ONLY : ymds2ju, ju2ymds ! for calendar 17 18 USE phycst ! ??? 18 19 USE in_out_manager ! I/O manager … … 29 30 LOGICAL :: ln_tint ! time interpolation or not (T/F) 30 31 LOGICAL :: ln_clim ! climatology or not (T/F) 31 CHARACTER(len = 7) :: cltype ! type of data file 'daily', 'monthly' or yearly'32 CHARACTER(len = 8) :: cltype ! type of data file 'daily', 'monthly' or yearly' 32 33 CHARACTER(len = 34) :: wname ! generic name of a NetCDF weights file to be used, blank if not 33 34 CHARACTER(len = 34) :: vcomp ! symbolic component name if a vector that needs rotation … … 43 44 LOGICAL :: ln_tint ! time interpolation or not (T/F) 44 45 LOGICAL :: ln_clim ! climatology or not (T/F) 45 CHARACTER(len = 7) :: cltype ! type of data file 'daily', 'monthly' or yearly'46 CHARACTER(len = 8) :: cltype ! type of data file 'daily', 'monthly' or yearly' 46 47 INTEGER :: num ! iom id of the jpfld files to be read 47 48 INTEGER :: nswap_sec ! swapping time in second since Jan. 1st 00h of nit000 year … … 78 79 INTEGER, DIMENSION(:,:,:), POINTER :: data_jpj ! array of source integers 79 80 REAL(wp), DIMENSION(:,:,:), POINTER :: data_wgt ! array of weights on model grid 80 REAL(wp), DIMENSION(:,: ), POINTER :: fly_dta ! array of values on input grid81 REAL(wp), DIMENSION(:,: ), POINTER :: col2 ! temporary array for reading in columns81 REAL(wp), DIMENSION(:,:,:), POINTER :: fly_dta ! array of values on input grid 82 REAL(wp), DIMENSION(:,:,:), POINTER :: col2 ! temporary array for reading in columns 82 83 END TYPE WGT 83 84 … … 93 94 !! OPA 9.0 , LOCEAN-IPSL (2006) 94 95 !! $Id$ 95 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)96 !! Software governed by the CeCILL licence (NEMOGCM/License_CeCILL.txt) 96 97 !!---------------------------------------------------------------------- 97 98 … … 159 160 160 161 ! last record to be read in the current file 161 IF( sd(jf)%nfreqh == -1 ) THEN ; ireclast = 12 162 IF( sd(jf)%nfreqh == -1 ) THEN 163 IF( sd(jf)%cltype == 'monthly' ) THEN ; ireclast = 1 164 ELSE ; ireclast = 12 165 ENDIF 162 166 ELSE 163 IF( sd(jf)%cltype == 'monthly' ) THEN ; ireclast = 24 * nmonth_len(nmonth) / sd(jf)%nfreqh 164 ELSEIF( sd(jf)%cltype == 'daily' ) THEN ; ireclast = 24 / sd(jf)%nfreqh 165 ELSE ; ireclast = 24 * nyear_len( 1 ) / sd(jf)%nfreqh 167 IF( sd(jf)%cltype == 'monthly' ) THEN ; ireclast = 24 * nmonth_len(nmonth) / sd(jf)%nfreqh 168 ELSEIF( sd(jf)%cltype(1:4) == 'week' ) THEN ; ireclast = 24.* 7 / sd(jf)%nfreqh 169 ELSEIF( sd(jf)%cltype == 'daily' ) THEN ; ireclast = 24 / sd(jf)%nfreqh 170 ELSE ; ireclast = 24 * nyear_len( 1 ) / sd(jf)%nfreqh 166 171 ENDIF 167 172 ENDIF … … 207 212 IF( LEN(TRIM(sd(jf)%wgtname)) > 0 ) THEN 208 213 CALL wgt_list( sd(jf), kw ) 209 DO jk = 1, ipk 210 CALL fld_interp( sd(jf)%num, sd(jf)%clvar, kw, sd(jf)%fdta(:,:,jk,2), sd(jf)%nrec_a(1) ) 211 END DO 214 ipk = SIZE(sd(jf)%fnow,3) 215 IF( sd(jf)%ln_tint ) THEN 216 CALL fld_interp( sd(jf)%num, sd(jf)%clvar , kw , ipk, sd(jf)%fdta(:,:,:,2) , sd(jf)%nrec_a(1) ) 217 ELSE 218 CALL fld_interp( sd(jf)%num, sd(jf)%clvar , kw , ipk, sd(jf)%fnow(:,:,:) , sd(jf)%nrec_a(1) ) 219 ENDIF 212 220 ELSE 213 IF( ipk == 1 ) THEN 214 CALL iom_get( sd(jf)%num, jpdom_data, sd(jf)%clvar, sd(jf)%fdta(:,:,1,2), sd(jf)%nrec_a(1) ) 215 ELSE 216 CALL iom_get( sd(jf)%num, jpdom_data, sd(jf)%clvar, sd(jf)%fdta(:,:,:,2), sd(jf)%nrec_a(1) ) 217 ENDIF 221 SELECT CASE( SIZE(sd(jf)%fnow,3) ) 222 CASE(1) 223 IF( sd(jf)%ln_tint ) THEN 224 CALL iom_get( sd(jf)%num, jpdom_data, sd(jf)%clvar, sd(jf)%fdta(:,:,1,2), sd(jf)%nrec_a(1) ) 225 ELSE 226 CALL iom_get( sd(jf)%num, jpdom_data, sd(jf)%clvar, sd(jf)%fnow(:,:,1) , sd(jf)%nrec_a(1) ) 227 ENDIF 228 CASE(jpk) 229 IF( sd(jf)%ln_tint ) THEN 230 CALL iom_get( sd(jf)%num, jpdom_data, sd(jf)%clvar, sd(jf)%fdta(:,:,:,2), sd(jf)%nrec_a(1) ) 231 ELSE 232 CALL iom_get( sd(jf)%num, jpdom_data, sd(jf)%clvar, sd(jf)%fnow(:,:,:) , sd(jf)%nrec_a(1) ) 233 ENDIF 234 END SELECT 218 235 ENDIF 219 236 sd(jf)%rotn(2) = .FALSE. … … 249 266 IF( kf > 0 ) THEN 250 267 !! fields jf,kf are two components which need to be rotated together 251 DO nf = 1,2 268 IF( sd(jf)%ln_tint )THEN 269 DO nf = 1,2 270 !! check each time level of this pair 271 IF( .NOT. sd(jf)%rotn(nf) .AND. .NOT. sd(kf)%rotn(nf) ) THEN 272 utmp(:,:) = 0.0 273 vtmp(:,:) = 0.0 274 ! 275 ipk = SIZE( sd(kf)%fnow(:,:,:) ,3 ) 276 DO jk = 1,ipk 277 CALL rot_rep( sd(jf)%fdta(:,:,jk,nf),sd(kf)%fdta(:,:,jk,nf),'T', 'en->i', utmp(:,:) ) 278 CALL rot_rep( sd(jf)%fdta(:,:,jk,nf),sd(kf)%fdta(:,:,jk,nf),'T', 'en->j', vtmp(:,:) ) 279 sd(jf)%fdta(:,:,jk,nf) = utmp(:,:) 280 sd(kf)%fdta(:,:,jk,nf) = vtmp(:,:) 281 ENDDO 282 ! 283 sd(jf)%rotn(nf) = .TRUE. 284 sd(kf)%rotn(nf) = .TRUE. 285 IF( lwp .AND. kt == nit000 ) & 286 WRITE(numout,*) 'fld_read: vector pair (', & 287 TRIM(sd(jf)%clvar),',',TRIM(sd(kf)%clvar), & 288 ') rotated on to model grid' 289 ENDIF 290 END DO 291 ELSE 252 292 !! check each time level of this pair 253 293 IF( .NOT. sd(jf)%rotn(nf) .AND. .NOT. sd(kf)%rotn(nf) ) THEN … … 255 295 vtmp(:,:) = 0.0 256 296 ! 257 DO jk = 1, SIZE( sd(kf)%fdta, 3 ) 258 CALL rot_rep( sd(jf)%fdta(:,:,jk,nf),sd(kf)%fdta(:,:,jk,nf),'T', 'en->i', utmp(:,:) ) 259 CALL rot_rep( sd(jf)%fdta(:,:,jk,nf),sd(kf)%fdta(:,:,jk,nf),'T', 'en->j', vtmp(:,:) ) 260 sd(jf)%fdta(:,:,jk,nf) = utmp(:,:) 261 sd(kf)%fdta(:,:,jk,nf) = vtmp(:,:) 262 END DO 297 ipk = SIZE( sd(kf)%fnow(:,:,:) ,3 ) 298 DO jk = 1,ipk 299 CALL rot_rep( sd(jf)%fnow(:,:,jk),sd(kf)%fnow(:,:,jk),'T', 'en->i', utmp(:,:) ) 300 CALL rot_rep( sd(jf)%fnow(:,:,jk),sd(kf)%fnow(:,:,jk),'T', 'en->j', vtmp(:,:) ) 301 sd(jf)%fnow(:,:,jk) = utmp(:,:) 302 sd(kf)%fnow(:,:,jk) = vtmp(:,:) 303 ENDDO 263 304 ! 264 305 sd(jf)%rotn(nf) = .TRUE. … … 269 310 ') rotated on to model grid' 270 311 ENDIF 271 END DO312 ENDIF 272 313 ENDIF 273 314 ENDIF … … 301 342 ENDIF 302 343 !CDIR COLLAPSE 303 sd(jf)%fnow(:,:,:) = sd(jf)%fdta(:,:,:,2) ! piecewise constant field304 305 344 ENDIF 306 345 ! … … 326 365 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 327 366 !! 328 LOGICAL :: llprevyr ! are we reading previous year file? 329 LOGICAL :: llprevmth ! are we reading previous month file? 330 LOGICAL :: llprevday ! are we reading previous day file? 331 LOGICAL :: llprev ! llprevyr .OR. llprevmth .OR. llprevday 332 INTEGER :: idvar ! variable id 333 INTEGER :: inrec ! number of record existing for this variable 367 LOGICAL :: llprevyr ! are we reading previous year file? 368 LOGICAL :: llprevmth ! are we reading previous month file? 369 LOGICAL :: llprevweek ! are we reading previous week file? 370 LOGICAL :: llprevday ! are we reading previous day file? 371 LOGICAL :: llprev ! llprevyr .OR. llprevmth .OR. llprevday 372 INTEGER :: idvar ! variable id 373 INTEGER :: inrec ! number of record existing for this variable 334 374 INTEGER :: kwgt 335 INTEGER :: jk ! vertical loop variable 336 INTEGER :: ipk ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 375 INTEGER :: jk !vertical loop variable 376 INTEGER :: ipk !number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 377 INTEGER :: iyear, imonth, iday ! first day of the current file in yyyy mm dd 378 INTEGER :: isec_week ! number of seconds since start of the weekly file 337 379 CHARACTER(LEN=1000) :: clfmt ! write format 338 380 !!--------------------------------------------------------------------- 339 381 340 382 ! some default definitions... 341 383 sdjf%num = 0 ! default definition for non-opened file 342 384 IF( sdjf%ln_clim ) sdjf%clname = TRIM( sdjf%clrootname ) ! file name defaut definition, never change in this case 343 llprevyr = .FALSE. 344 llprevmth = .FALSE. 345 llprevday = .FALSE. 385 llprevyr = .FALSE. 386 llprevmth = .FALSE. 387 llprevweek = .FALSE. 388 llprevday = .FALSE. 389 isec_week = 0 346 390 347 391 ! define record informations … … 365 409 llprevmth = .NOT. sdjf%ln_clim ! use previous month file? 366 410 llprevyr = llprevmth .AND. nmonth == 1 ! use previous year file? 411 ELSE IF ( sdjf%cltype(1:4) == 'week' ) THEN !weekly file 412 isec_week = 86400 * 7 413 sdjf%nrec_b(1) = 24. / sdjf%nfreqh * 7 ! last record of previous weekly file 367 414 ELSEIF( sdjf%cltype == 'daily' ) THEN ! daily file 368 415 sdjf%nrec_b(1) = 24 / sdjf%nfreqh ! last record of previous day … … 376 423 ENDIF 377 424 ENDIF 378 llprev = llprevyr .OR. llprevmth .OR. llprev day425 llprev = llprevyr .OR. llprevmth .OR. llprevweek .OR. llprevday 379 426 380 427 CALL fld_clopn( sdjf, nyear - COUNT((/llprevyr /)) , & 381 428 & nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /)), & 382 429 & nday - COUNT((/llprevday/)) + nmonth_len(nmonth-1) * COUNT((/llprevmth/)), .NOT. llprev ) 383 430 431 IF ( sdjf%cltype(1:4) == 'week' ) THEN 432 isec_week = ksec_week( sdjf%cltype(6:8) ) 433 if(lwp)write(numout,*)'cbr test2 isec_week = ',isec_week 434 llprevmth = ( isec_week .GT. nsec_month ) 435 llprevyr = llprevmth .AND. nmonth==1 436 ENDIF 437 ! 438 iyear = nyear - COUNT((/llprevyr /)) 439 imonth = nmonth - COUNT((/llprevmth/)) + 12* COUNT((/llprevyr /)) 440 iday = nday - COUNT((/llprevday/)) + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - INT( isec_week )/86400 441 ! 442 CALL fld_clopn( sdjf , iyear , imonth , iday , .NOT. llprev ) 443 384 444 ! if previous year/month/day file does not exist, we switch to the current year/month/day 385 445 IF( llprev .AND. sdjf%num <= 0 ) THEN … … 399 459 400 460 ! read before data into sdjf%fdta(:,:,2) because we will swap data in the following part of fld_read 401 ipk = SIZE( sdjf%fdta, 3 )402 461 IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN 403 462 CALL wgt_list( sdjf, kwgt ) 404 DO jk = 1, ipk 405 CALL fld_interp( sdjf%num, sdjf%clvar, kwgt, sdjf%fdta(:,:,jk,2), sdjf%nrec_b(1) ) 406 END DO 463 ipk = SIZE(sdjf%fnow,3) 464 IF( sdjf%ln_tint ) THEN 465 CALL fld_interp( sdjf%num, sdjf%clvar, kwgt, ipk, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1) ) 466 ELSE 467 CALL fld_interp( sdjf%num, sdjf%clvar, kwgt, ipk, sdjf%fnow(:,:,:) , sdjf%nrec_a(1) ) 468 ENDIF 407 469 ELSE 408 IF( ipk == 1 ) THEN 409 CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_b(1) ) 410 ELSE 411 CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_b(1) ) 412 ENDIF 470 SELECT CASE( SIZE(sdjf%fnow,3) ) 471 CASE(1) 472 IF( sdjf%ln_tint ) THEN 473 CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fdta(:,:,1,2), sdjf%nrec_b(1) ) 474 ELSE 475 CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fnow(:,:,1) , sdjf%nrec_b(1) ) 476 ENDIF 477 CASE(jpk) 478 IF( sdjf%ln_tint ) THEN 479 CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_b(1) ) 480 ELSE 481 CALL iom_get( sdjf%num, jpdom_data, sdjf%clvar, sdjf%fnow(:,:,:) , sdjf%nrec_b(1) ) 482 ENDIF 483 END SELECT 413 484 ENDIF 414 485 sdjf%rotn(2) = .FALSE. … … 422 493 423 494 IF( sdjf%num <= 0 ) CALL fld_clopn( sdjf, nyear, nmonth, nday ) ! make sure current year/month/day file is opened 495 ! make sure current year/month/day file is opened 496 IF( sdjf%num == 0 ) THEN 497 isec_week = 0 498 llprevyr = .FALSE. 499 llprevmth = .FALSE. 500 llprevweek = .FALSE. 501 ! 502 IF ( sdjf%cltype(1:4) == 'week' ) THEN 503 isec_week = ksec_week( sdjf%cltype(6:8) ) 504 llprevmth = ( isec_week .GT. nsec_month ) 505 llprevyr = llprevmth .AND. nmonth==1 506 ENDIF 507 ! 508 iyear = nyear - COUNT((/llprevyr /)) 509 imonth = nmonth - COUNT((/llprevmth/)) + 12* COUNT((/llprevyr /)) 510 iday = nday + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week/86400 511 ! 512 CALL fld_clopn( sdjf, iyear, imonth, iday ) 513 ENDIF 424 514 425 515 sdjf%nswap_sec = nsec_year + nsec1jan000 - 1 ! force read/update the after data in the following part of fld_read 426 516 517 427 518 END SUBROUTINE fld_init 428 519 … … 442 533 REAL(wp) :: ztmp ! temporary variable 443 534 INTEGER :: ifreq_sec ! frequency mean (in seconds) 535 INTEGER :: isec_week ! number of seconds since the start of the weekly file 444 536 !!---------------------------------------------------------------------- 445 537 ! … … 458 550 ! forcing record : nmonth 459 551 ! 552 ztmp = 0.e0 460 553 ztmp = REAL( nday, wp ) / REAL( nmonth_len(nmonth), wp ) + 0.5 461 554 ELSE … … 468 561 ENDIF 469 562 470 sdjf%nrec_a(:) = (/ irec, nmonth_half(irec) + nsec1jan000 /) ! define after record number and time 471 irec = irec - 1 ! move back to previous record 472 sdjf%nrec_b(:) = (/ irec, nmonth_half(irec) + nsec1jan000 /) ! define before record number and time 563 IF( sdjf%cltype == 'monthly' ) THEN 564 565 sdjf%nrec_b(:) = (/ 0, nmonth_half(irec - 1 ) + nsec1jan000 /) 566 sdjf%nrec_a(:) = (/ 1, nmonth_half(irec ) + nsec1jan000 /) 567 568 IF( ztmp == 1. ) THEN 569 sdjf%nrec_b(1) = 1 570 sdjf%nrec_a(1) = 2 571 ENDIF 572 573 ELSE 574 575 sdjf%nrec_a(:) = (/ irec, nmonth_half(irec) + nsec1jan000 /) ! define after record number and time 576 irec = irec - 1 ! move back to previous record 577 sdjf%nrec_b(:) = (/ irec, nmonth_half(irec) + nsec1jan000 /) ! define before record number and time 578 579 ENDIF 473 580 ! 474 581 ELSE ! higher frequency mean (in hours) 475 582 ! 476 583 ifreq_sec = sdjf%nfreqh * 3600 ! frequency mean (in seconds) 584 IF( sdjf%cltype(1:4) == 'week' ) isec_week = ksec_week( sdjf%cltype(6:8)) !since the first day of the current week 477 585 ! number of second since the beginning of the file 478 IF( sdjf%cltype == 'monthly' ) THEN ; ztmp = REAL(nsec_month,wp) ! since 00h on the 1st day of the current month 479 ELSEIF( sdjf%cltype == 'daily' ) THEN ; ztmp = REAL(nsec_day ,wp) ! since 00h of the current day 480 ELSE ; ztmp = REAL(nsec_year ,wp) ! since 00h on Jan 1 of the current year 586 IF( sdjf%cltype == 'monthly' ) THEN ; ztmp = REAL(nsec_month ,wp) ! since 00h on the 1st day of the current month 587 ELSEIF( sdjf%cltype(1:4) == 'week' ) THEN ; ztmp = REAL(isec_week ,wp) ! since the first day of the current week 588 ELSEIF( sdjf%cltype == 'daily' ) THEN ; ztmp = REAL(nsec_day ,wp) ! since 00h of the current day 589 ELSE ; ztmp = REAL(nsec_year ,wp) ! since 00h on Jan 1 of the current year 481 590 ENDIF 482 591 IF( sdjf%ln_tint ) THEN ! time interpolation, shift by 1/2 record … … 514 623 ! after record index and second since Jan. 1st 00h of nit000 year 515 624 sdjf%nrec_a(:) = (/ irec, ifreq_sec * irec - ifreq_sec / 2 + nsec1jan000 /) 516 IF( sdjf%cltype == 'monthly' ) & ! add the number of seconds between 00h Jan 1 and the end of previous month625 IF( sdjf%cltype == 'monthly' ) & ! add the number of seconds between 00h Jan 1 and the end of previous month 517 626 sdjf%nrec_a(2) = sdjf%nrec_a(2) + isecd * SUM(nmonth_len(1:nmonth -1)) ! ok if nmonth=1 518 IF( sdjf%cltype == 'daily' ) & ! add the number of seconds between 00h Jan 1 and the end of previous day 627 IF( sdjf%cltype(1:4) == 'week' ) & ! add the number of seconds between 00h Jan 1 and the end of previous week 628 sdjf%nrec_a(2) = sdjf%nrec_a(2) + ( nsec_year - isec_week ) 629 IF( sdjf%cltype == 'daily' ) & ! add the number of seconds between 00h Jan 1 and the end of previous day 519 630 sdjf%nrec_a(2) = sdjf%nrec_a(2) + isecd * ( nday_year - 1 ) 520 631 … … 522 633 irec = irec - 1. ! move back to previous record 523 634 sdjf%nrec_b(:) = (/ irec, ifreq_sec * irec - ifreq_sec / 2 + nsec1jan000 /) 524 IF( sdjf%cltype == 'monthly' ) & ! add the number of seconds between 00h Jan 1 and the end of previous month635 IF( sdjf%cltype == 'monthly' ) & ! add the number of seconds between 00h Jan 1 and the end of previous month 525 636 sdjf%nrec_b(2) = sdjf%nrec_b(2) + isecd * SUM(nmonth_len(1:nmonth -1)) ! ok if nmonth=1 526 IF( sdjf%cltype == 'daily' ) & ! add the number of seconds between 00h Jan 1 and the end of previous day 637 IF( sdjf%cltype(1:4) == 'week' ) & ! add the number of seconds between 00h Jan 1 and the end of previous week 638 sdjf%nrec_b(2) = sdjf%nrec_b(2) + ( nsec_year - isec_week ) 639 IF( sdjf%cltype == 'daily' ) & ! add the number of seconds between 00h Jan 1 and the end of previous day 527 640 sdjf%nrec_b(2) = sdjf%nrec_b(2) + isecd * ( nday_year - 1 ) 528 641 … … 545 658 !! ** Method : 546 659 !!---------------------------------------------------------------------- 547 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 548 INTEGER , INTENT(in ) :: kyear ! year value 549 INTEGER , INTENT(in ) :: kmonth ! month value 550 INTEGER , INTENT(in ) :: kday ! day value 551 LOGICAL , INTENT(in ), OPTIONAL :: ldstop ! stop if open to read a non-existing file (default = .TRUE.) 660 TYPE(FLD), INTENT(inout) :: sdjf ! input field related variables 661 INTEGER , INTENT(in ) :: kyear ! year value 662 INTEGER , INTENT(in ) :: kmonth ! month value 663 INTEGER , INTENT(in ) :: kday ! day value 664 LOGICAL , INTENT(in ), OPTIONAL :: ldstop ! stop if open to read a non-existing file (default = .TRUE.) 665 INTEGER :: iyear, imonth, iday ! firt day of the current week in yyyy mm dd 666 REAL(wp) :: zsec, zjul !temp variable 552 667 553 668 IF( sdjf%num /= 0 ) CALL iom_close( sdjf%num ) ! close file if already open 554 669 ! build the new filename if not climatological data 555 IF( .NOT. sdjf%ln_clim ) THEN ; WRITE(sdjf%clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), kyear ! add year 556 IF( sdjf%cltype /= 'yearly' ) WRITE(sdjf%clname, '(a,"m" ,i2.2)' ) TRIM( sdjf%clname ), kmonth ! add month 557 IF( sdjf%cltype == 'daily' ) WRITE(sdjf%clname, '(a,"d" ,i2.2)' ) TRIM( sdjf%clname ), kday ! add day 670 sdjf%clname=TRIM(sdjf%clrootname) 671 ! 672 IF( sdjf%cltype(1:4) == 'week' .AND. nn_leapy==0 )CALL ctl_stop( 'fld_clopn: weekly file and nn_leapy=0 are not compatible' ) 673 ! 674 IF( .NOT. sdjf%ln_clim ) THEN 675 WRITE(sdjf%clname, '(a,"_y",i4.4)' ) TRIM( sdjf%clrootname ), kyear ! add year 676 IF( sdjf%cltype /= 'yearly' ) & 677 & WRITE(sdjf%clname, '(a,"m" ,i2.2)' ) TRIM( sdjf%clname ), kmonth ! add month 678 IF( sdjf%cltype == 'daily' .OR. sdjf%cltype(1:4) == 'week' ) & 679 & WRITE(sdjf%clname, '(a,"d" ,i2.2)' ) TRIM( sdjf%clname ), kday ! add day 680 ELSE 681 ! build the new filename if climatological data 682 IF( sdjf%cltype == 'monthly' ) WRITE(sdjf%clname, '(a,"_m" ,i2.2)' ) TRIM( sdjf%clrootname ), kmonth ! add month 558 683 ELSE 559 684 ! build the new filename if climatological data … … 610 735 & ' pairing : ' , TRIM( sdf(jf)%vcomp ), & 611 736 & ' data type: ' , sdf(jf)%cltype 737 call flush(numout) 612 738 END DO 613 739 ENDIF … … 707 833 INTEGER :: inum ! temporary logical unit 708 834 INTEGER :: id ! temporary variable id 835 INTEGER :: ipk ! temporary vertical dimension 709 836 CHARACTER (len=5) :: aname 710 837 INTEGER , DIMENSION(3) :: ddims … … 869 996 ! SA: +3 stencil is a patch to avoid out-of-bound computation in some configuration. 870 997 ! a more robust solution will be given in next release 871 ALLOCATE( ref_wgts(nxt_wgt)%fly_dta(ref_wgts(nxt_wgt)%jpiwgt+3, ref_wgts(nxt_wgt)%jpjwgt+3) ) 872 IF( ref_wgts(nxt_wgt)%cyclic ) ALLOCATE( ref_wgts(nxt_wgt)%col2(2,ref_wgts(nxt_wgt)%jpjwgt+3) ) 998 ipk = SIZE(sd%fnow,3) 999 ALLOCATE( ref_wgts(nxt_wgt)%fly_dta(ref_wgts(nxt_wgt)%jpiwgt+3, ref_wgts(nxt_wgt)%jpjwgt+3 ,ipk) ) 1000 IF( ref_wgts(nxt_wgt)%cyclic ) ALLOCATE( ref_wgts(nxt_wgt)%col2(2,ref_wgts(nxt_wgt)%jpjwgt+3,ipk) ) 873 1001 874 1002 nxt_wgt = nxt_wgt + 1 … … 880 1008 END SUBROUTINE fld_weight 881 1009 882 SUBROUTINE fld_interp(num, clvar, kw, dta, nrec)1010 SUBROUTINE fld_interp(num, clvar, kw, kk, dta, nrec) 883 1011 !!--------------------------------------------------------------------- 884 1012 !! *** ROUTINE fld_interp *** … … 892 1020 CHARACTER(LEN=*), INTENT(in) :: clvar ! variable name 893 1021 INTEGER, INTENT(in) :: kw ! weights number 894 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj) :: dta ! output field on model grid 1022 INTEGER, INTENT(in) :: kk ! vertical dimension of kk 1023 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,kk) :: dta ! output field on model grid 895 1024 INTEGER, INTENT(in) :: nrec ! record number to read (ie time slice) 896 1025 !! 897 INTEGER, DIMENSION( 2) :: rec1,recn ! temporary arrays for start and length1026 INTEGER, DIMENSION(3) :: rec1,recn ! temporary arrays for start and length 898 1027 INTEGER :: jk, jn, jm ! loop counters 899 1028 INTEGER :: ni, nj ! lengths … … 918 1047 rec1(1) = MAX( jpimin-1, 1 ) 919 1048 rec1(2) = MAX( jpjmin-1, 1 ) 1049 rec1(3) = 1 920 1050 recn(1) = MIN( jpiwid+2, ref_wgts(kw)%ddims(1)-rec1(1)+1 ) 921 1051 recn(2) = MIN( jpjwid+2, ref_wgts(kw)%ddims(2)-rec1(2)+1 ) 1052 recn(3) = kk 922 1053 923 1054 !! where we need to read it to … … 927 1058 jpj2 = jpj1 + recn(2) - 1 928 1059 929 ref_wgts(kw)%fly_dta(:,:) = 0.0 930 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2), nrec, rec1, recn) 1060 ref_wgts(kw)%fly_dta(:,:,:) = 0.0 1061 SELECT CASE( SIZE(ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:),3) ) 1062 CASE(1) 1063 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,1), nrec, rec1, recn) 1064 CASE(jpk) 1065 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%fly_dta(jpi1:jpi2,jpj1:jpj2,:), nrec, rec1, recn) 1066 END SELECT 931 1067 932 1068 !! first four weights common to both bilinear and bicubic 933 1069 !! note that we have to offset by 1 into fly_dta array because of halo 934 dta(:,: ) = 0.01070 dta(:,:,:) = 0.0 935 1071 DO jk = 1,4 936 DO jn = 1, jpj937 DO jm = 1, jpi1072 DO jn = 1, nlcj 1073 DO jm = 1,nlci 938 1074 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 939 1075 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 940 dta(jm,jn ) = dta(jm,jn) + ref_wgts(kw)%data_wgt(jm,jn,jk) * ref_wgts(kw)%fly_dta(ni+1,nj+1)1076 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk) * ref_wgts(kw)%fly_dta(ni+1,nj+1,jk) 941 1077 END DO 942 1078 END DO … … 947 1083 !! fix up halo points that we couldnt read from file 948 1084 IF( jpi1 == 2 ) THEN 949 ref_wgts(kw)%fly_dta(jpi1-1,: ) = ref_wgts(kw)%fly_dta(jpi1,:)1085 ref_wgts(kw)%fly_dta(jpi1-1,:,:) = ref_wgts(kw)%fly_dta(jpi1,:,:) 950 1086 ENDIF 951 1087 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 952 ref_wgts(kw)%fly_dta(jpi2+1,: ) = ref_wgts(kw)%fly_dta(jpi2,:)1088 ref_wgts(kw)%fly_dta(jpi2+1,:,:) = ref_wgts(kw)%fly_dta(jpi2,:,:) 953 1089 ENDIF 954 1090 IF( jpj1 == 2 ) THEN 955 ref_wgts(kw)%fly_dta(:,jpj1-1 ) = ref_wgts(kw)%fly_dta(:,jpj1)1091 ref_wgts(kw)%fly_dta(:,jpj1-1,:) = ref_wgts(kw)%fly_dta(:,jpj1,:) 956 1092 ENDIF 957 1093 IF( jpj2 + jpjmin - 1 == ref_wgts(kw)%ddims(2)+1 .AND. jpj2 .lt. jpjwid+2 ) THEN 958 ref_wgts(kw)%fly_dta(:,jpj2+1 ) = 2.0*ref_wgts(kw)%fly_dta(:,jpj2) - ref_wgts(kw)%fly_dta(:,jpj2-1)1094 ref_wgts(kw)%fly_dta(:,jpj2+1,:) = 2.0*ref_wgts(kw)%fly_dta(:,jpj2,:) - ref_wgts(kw)%fly_dta(:,jpj2-1,:) 959 1095 ENDIF 960 1096 … … 969 1105 IF( jpi1 == 2 ) THEN 970 1106 rec1(1) = ref_wgts(kw)%ddims(1) - 1 971 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col2(:,jpj1:jpj2), nrec, rec1, recn) 972 ref_wgts(kw)%fly_dta(jpi1-1,jpj1:jpj2) = ref_wgts(kw)%col2(ref_wgts(kw)%offset+1,jpj1:jpj2) 1107 SELECT CASE( SIZE( ref_wgts(kw)%col2(:,jpj1:jpj2,:),3) ) 1108 CASE(1) 1109 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col2(:,jpj1:jpj2,1), nrec, rec1, recn) 1110 CASE(jpk) 1111 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col2(:,jpj1:jpj2,:), nrec, rec1, recn) 1112 END SELECT 1113 ref_wgts(kw)%fly_dta(jpi1-1,jpj1:jpj2,:) = ref_wgts(kw)%col2(ref_wgts(kw)%offset+1,jpj1:jpj2,:) 973 1114 ENDIF 974 1115 IF( jpi2 + jpimin - 1 == ref_wgts(kw)%ddims(1)+1 ) THEN 975 1116 rec1(1) = 1 976 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col2(:,jpj1:jpj2), nrec, rec1, recn) 977 ref_wgts(kw)%fly_dta(jpi2+1,jpj1:jpj2) = ref_wgts(kw)%col2(2-ref_wgts(kw)%offset,jpj1:jpj2) 1117 SELECT CASE( SIZE( ref_wgts(kw)%col2(:,jpj1:jpj2,:),3) ) 1118 CASE(1) 1119 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col2(:,jpj1:jpj2,1), nrec, rec1, recn) 1120 CASE(jpk) 1121 CALL iom_get( num, jpdom_unknown, clvar, ref_wgts(kw)%col2(:,jpj1:jpj2,:), nrec, rec1, recn) 1122 END SELECT 1123 ref_wgts(kw)%fly_dta(jpi2+1,jpj1:jpj2,:) = ref_wgts(kw)%col2(2-ref_wgts(kw)%offset,jpj1:jpj2,:) 978 1124 ENDIF 979 1125 ENDIF … … 981 1127 ! gradient in the i direction 982 1128 DO jk = 1,4 983 DO jn = 1, jpj984 DO jm = 1, jpi1129 DO jn = 1, nlcj 1130 DO jm = 1,nlci 985 1131 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 986 1132 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 987 dta(jm,jn ) = dta(jm,jn) + ref_wgts(kw)%data_wgt(jm,jn,jk+4) * 0.5 * &988 (ref_wgts(kw)%fly_dta(ni+2,nj+1 ) - ref_wgts(kw)%fly_dta(ni,nj+1))1133 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+4) * 0.5 * & 1134 (ref_wgts(kw)%fly_dta(ni+2,nj+1,:) - ref_wgts(kw)%fly_dta(ni,nj+1,:)) 989 1135 END DO 990 1136 END DO … … 993 1139 ! gradient in the j direction 994 1140 DO jk = 1,4 995 DO jn = 1, jpj996 DO jm = 1, jpi1141 DO jn = 1, nlcj 1142 DO jm = 1,nlci 997 1143 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 998 1144 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 999 dta(jm,jn ) = dta(jm,jn) + ref_wgts(kw)%data_wgt(jm,jn,jk+8) * 0.5 * &1000 (ref_wgts(kw)%fly_dta(ni+1,nj+2 ) - ref_wgts(kw)%fly_dta(ni+1,nj))1145 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+8) * 0.5 * & 1146 (ref_wgts(kw)%fly_dta(ni+1,nj+2,:) - ref_wgts(kw)%fly_dta(ni+1,nj,:)) 1001 1147 END DO 1002 1148 END DO … … 1009 1155 ni = ref_wgts(kw)%data_jpi(jm,jn,jk) 1010 1156 nj = ref_wgts(kw)%data_jpj(jm,jn,jk) 1011 dta(jm,jn ) = dta(jm,jn) + ref_wgts(kw)%data_wgt(jm,jn,jk+12) * 0.25 * ( &1012 (ref_wgts(kw)%fly_dta(ni+2,nj+2 ) - ref_wgts(kw)%fly_dta(ni ,nj+2)) - &1013 (ref_wgts(kw)%fly_dta(ni+2,nj ) - ref_wgts(kw)%fly_dta(ni ,nj)))1157 dta(jm,jn,:) = dta(jm,jn,:) + ref_wgts(kw)%data_wgt(jm,jn,jk+12) * 0.25 * ( & 1158 (ref_wgts(kw)%fly_dta(ni+2,nj+2,:) - ref_wgts(kw)%fly_dta(ni ,nj+2,:)) - & 1159 (ref_wgts(kw)%fly_dta(ni+2,nj ,:) - ref_wgts(kw)%fly_dta(ni ,nj ,:))) 1014 1160 END DO 1015 1161 END DO … … 1019 1165 1020 1166 END SUBROUTINE fld_interp 1021 1167 1168 FUNCTION ksec_week( cdday ) 1169 !!--------------------------------------------------------------------- 1170 !! *** FUNCTION kshift_week *** 1171 !! 1172 !! ** Purpose : 1173 !! 1174 !! ** Method : 1175 !!--------------------------------------------------------------------- 1176 CHARACTER(len=*), INTENT(in) :: cdday !3 first letters of the first day of the weekly file 1177 !! 1178 INTEGER :: ksec_week ! output variable 1179 INTEGER :: ijul !temp variable 1180 INTEGER :: ishift !temp variable 1181 CHARACTER(len=3),DIMENSION(7) :: cl_week 1182 !!---------------------------------------------------------------------- 1183 cl_week = (/"sun","sat","fri","thu","wed","tue","mon"/) 1184 DO ijul=1,7 1185 IF( cl_week(ijul)==TRIM(cdday) ) EXIT 1186 ENDDO 1187 IF( ijul .GT. 7 ) CALL ctl_stop( 'ksec_week: wrong day for sdjf%cltype(6:8): ',TRIM(cdday) ) 1188 ! 1189 ishift = ( ijul ) * 86400 1190 ! 1191 ksec_week = nsec_week + ishift 1192 ksec_week = MOD( ksec_week , 86400*7 ) 1193 if(lwp)write(numout,*)'cbr ijul ksec_week ',ijul,ksec_week 1194 ! 1195 END FUNCTION ksec_week 1196 1022 1197 END MODULE fldread
Note: See TracChangeset
for help on using the changeset viewer.
