ROMS
Loading...
Searching...
No Matches
esmf_wrf_mod::wrf_processimport Interface Reference

Public Member Functions

subroutine wrf_processimport_scalar (grid, model, got, ifield, fieldname, lbi, ubi, lbj, ubj, focn, fdat, rc)
 
subroutine wrf_processimport_vector (grid, model, got, ifield, fieldname, lbi, ubi, lbj, ubj, ubk, focn, fdat, rc)
 

Detailed Description

Definition at line 112 of file esmf_atm_wrf.h.

Member Function/Subroutine Documentation

◆ wrf_processimport_scalar()

subroutine esmf_wrf_mod::wrf_processimport::wrf_processimport_scalar ( type (domain), pointer grid,
type (esmf_gridcomp) model,
logical, dimension(2), intent(in) got,
integer, dimension(2), intent(in) ifield,
character (len=*), dimension(:), intent(in) fieldname,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
real (dp), dimension(lbi:ubi,lbj:ubj), intent(in) focn,
real (dp), dimension(lbi:ubi,lbj:ubj), intent(in) fdat,
integer, intent(out) rc )

Definition at line 3151 of file esmf_atm_wrf.h.

3156!
3157!=======================================================================
3158! !
3159! If both import fields Focn and Fdat are avaliable, it merges !
3160! its values. Otherwise, it loads available data into ouput field, !
3161! Fout. Only sea-water or sea-ice points are processed. !
3162! !
3163! The mergin is needed when the atmosphere and ocean grids are !
3164! incongruent. The DATA model provides values on those grid points !
3165! not covered by the OCEAN model. !
3166! !
3167! On Input: !
3168! !
3169! grid WRF grid state structure (TYPE domain) !
3170! model Gridded component object (TYPE ESMF_GridComp) !
3171! got Switches indicating source and availability of !
3172! imported data (logical vector): !
3173! got(1) OCEAN model switch (T/F) !
3174! got(2) DATA model switch (T/F) !
3175! ifield Imported field index (integer vector) !
3176! ifield(1) OCEAN model field index !
3177! ifield(2) DATA model field index !
3178! FieldName Field shortname (string array) !
3179! ifield Imported field index (integer) !
3180! LBi I-dimension lower bound (integer) !
3181! UBi I-dimension upper bound (integer) !
3182! LBj J-dimension lower bound (integer) !
3183! UBj J-dimension upper bound (integer) !
3184! Focn Imported field from OCEAN model (2D real array) !
3185! Fdat Imported field from DATA model (2D real array) !
3186! !
3187! On Output: !
3188! !
3189! rc Return code (integer) !
3190! !
3191!=======================================================================
3192!
3193 USE module_domain, ONLY : domain
3194 USE strings_mod, ONLY : lowercase
3195!
3196! Imported variable declarations.
3197!
3198 logical, intent(in) :: got(2)
3199!
3200 integer, intent(in) :: ifield(2)
3201 integer, intent(in) :: LBi, UBi, LBj, UBj
3202 integer, intent(out) :: rc
3203!
3204 real (dp), intent(in) :: Focn(LBi:UBi,LBj:UBj)
3205 real (dp), intent(in) :: Fdat(LBi:UBi,LBj:UBj)
3206!
3207 character (len=*), intent(in) :: FieldName(:)
3208!
3209 TYPE (domain), pointer :: grid
3210 TYPE (ESMF_GridComp) :: model
3211!
3212! Local variable declarations.
3213!
3214 logical :: got_dat, got_ocn
3215 logical :: DebugWrite(2) = (/ .false., .false. /)
3216!
3217 integer :: i, ic, is, j, ng
3218 integer :: year, month, day, hour, minutes, seconds, sN, SD
3219 integer :: LakeValue, LandValue
3220 integer :: localDE, localDEcount, localPET, PETcount
3221 integer :: IminP, ImaxP, JminP, JmaxP
3222!
3223 real (dp) :: Fseconds, TimeInDays, Time_Current
3224
3225 real (dp) :: Fval, MyFmax(3), MyFmin(3), Fmin(3), Fmax(3)
3226!
3227 real (dp), pointer :: ptr2d(:,:) => null()
3228!
3229 real (KIND(grid%sst)), pointer :: Fout(:,:) => null()
3230!
3231 character (len=22 ) :: Time_CurrentString
3232
3233 character (len=*), parameter :: MyFile = &
3234 & __FILE__//", WRF_ProcessImport_scalar"
3235!
3236 character (ESMF_MAXSTR) :: cname, fld_string, ofile
3237!
3238 TYPE (ESMF_ArraySpec) :: arraySpec2d
3239 TYPE (ESMF_Clock) :: clock
3240 TYPE (ESMF_Field) :: Fmerge
3241 TYPE (ESMF_StaggerLoc) :: staggerLoc
3242 TYPE (ESMF_Time) :: CurrentTime
3243 TYPE (ESMF_VM) :: vm
3244!
3245!-----------------------------------------------------------------------
3246! Initialize return code flag to success state (no error).
3247!-----------------------------------------------------------------------
3248!
3249 IF (esm_track) THEN
3250 WRITE (trac,'(a,a,i0)') '==> Entering WRF_ProcessImport_scalar',&
3251 & ', PET', petrank
3252 FLUSH (trac)
3253 END IF
3254 rc=esmf_success
3255!
3256!-----------------------------------------------------------------------
3257! Get information about the gridded component.
3258!-----------------------------------------------------------------------
3259!
3260 CALL esmf_gridcompget (model, &
3261 & clock=clock, &
3262 & localpet=localpet, &
3263 & petcount=petcount, &
3264 & vm=vm, &
3265 & name=cname, &
3266 & rc=rc)
3267 IF (esmf_logfounderror(rctocheck=rc, &
3268 & msg=esmf_logerr_passthru, &
3269 & line=__line__, &
3270 & file=myfile)) THEN
3271 RETURN
3272 END IF
3273 ng=grid%grid_id
3274!
3275! Get number of local decomposition elements (DEs). Usually, a single
3276! DE is associated with each Persistent Execution Thread (PETs). Thus,
3277! localDEcount=1.
3278!
3279 CALL esmf_gridget (models(iatmos)%grid(ng), &
3280 & localdecount=localdecount, &
3281 & rc=rc)
3282 IF (esmf_logfounderror(rctocheck=rc, &
3283 & msg=esmf_logerr_passthru, &
3284 & line=__line__, &
3285 & file=myfile)) THEN
3286 RETURN
3287 END IF
3288!
3289! Get current time.
3290!
3291 CALL esmf_clockget (clock, &
3292 & currtime=currenttime, &
3293 & rc=rc)
3294 IF (esmf_logfounderror(rctocheck=rc, &
3295 & msg=esmf_logerr_passthru, &
3296 & line=__line__, &
3297 & file=myfile)) THEN
3298 RETURN
3299 END IF
3300!
3301 CALL esmf_timeget (currenttime, &
3302 & yy=year, &
3303 & mm=month, &
3304 & dd=day, &
3305 & h =hour, &
3306 & m =minutes, &
3307 & s =seconds, &
3308 & sn=sn, &
3309 & sd=sd, &
3310 & rc=rc)
3311 IF (esmf_logfounderror(rctocheck=rc, &
3312 & msg=esmf_logerr_passthru, &
3313 & line=__line__, &
3314 & file=myfile)) THEN
3315 RETURN
3316 END IF
3317!
3318 CALL esmf_timeget (currenttime, &
3319 & s_r8=time_current, &
3320 & timestring=time_currentstring, &
3321 & rc=rc)
3322 IF (esmf_logfounderror(rctocheck=rc, &
3323 & msg=esmf_logerr_passthru, &
3324 & line=__line__, &
3325 & file=myfile)) THEN
3326 RETURN
3327 END IF
3328 fseconds=real(seconds,dp)+real(sn,dp)/real(sd,dp)
3329 timeindays=time_current/86400.0_dp
3330 is=index(time_currentstring, 'T') ! remove 'T' in
3331 IF (is.gt.0) time_currentstring(is:is)=' ' ! ISO 8601 format
3332!
3333!-----------------------------------------------------------------------
3334! Create merged field array.
3335!-----------------------------------------------------------------------
3336!
3337! Set a 2D floating-point array descriptor.
3338!
3339 CALL esmf_arrayspecset (arrayspec2d, &
3340 & typekind=esmf_typekind_r8, &
3341 & rank=2, &
3342 & rc=rc)
3343 IF (esmf_logfounderror(rctocheck=rc, &
3344 & msg=esmf_logerr_passthru, &
3345 & line=__line__, &
3346 & file=myfile)) THEN
3347 RETURN
3348 END IF
3349!
3350! Create 2D merge field from the Grid and arraySpec.
3351!
3352 got_ocn=got(1)
3353 got_dat=got(2)
3354!
3355 IF (.not.got_dat.and.got_ocn) THEN
3356 debugwrite(1)=models(iatmos)%ImportField(ifield(1))%debug_write
3357 fld_string=trim(fieldname(1))
3358 ELSE IF (.not.got_ocn.and.got_dat) THEN
3359 debugwrite(2)=models(iatmos)%ImportField(ifield(2))%debug_write
3360 fld_string=trim(fieldname(2))
3361 ELSE IF (got_ocn.and.got_dat) THEN
3362 debugwrite(1)=models(iatmos)%ImportField(ifield(1))%debug_write
3363 debugwrite(2)=models(iatmos)%ImportField(ifield(2))%debug_write
3364 fld_string=trim(fieldname(1))//'-'//trim(fieldname(2))
3365 END IF
3366 staggerloc=esmf_staggerloc_center
3367!
3368 fmerge=esmf_fieldcreate(models(iatmos)%grid(ng), &
3369 & arrayspec2d, &
3370 & staggerloc=staggerloc, &
3371 & name=trim(fld_string), &
3372 & rc=rc)
3373 IF (esmf_logfounderror(rctocheck=rc, &
3374 & msg=esmf_logerr_passthru, &
3375 & line=__line__, &
3376 & file=myfile)) THEN
3377 RETURN
3378 END IF
3379!
3380! Get pointer for merge array.
3381!
3382 CALL esmf_fieldget (fmerge, &
3383 & farrayptr=ptr2d, &
3384 & rc=rc)
3385 IF (esmf_logfounderror(rctocheck=rc, &
3386 & msg=esmf_logerr_passthru, &
3387 & line=__line__, &
3388 & file=myfile)) THEN
3389 RETURN
3390 END IF
3391 ptr2d=missing_dp
3392!
3393!-----------------------------------------------------------------------
3394! Create pointer to WRF export field target. This logic is awkward but
3395! needed since WRF could have different floating-point representations.
3396! Passing WRF variables as arguments to this routine doesn't work.
3397!-----------------------------------------------------------------------
3398!
3399 SELECT CASE (lowercase(trim(fld_string)))
3400 CASE ('sst', 'dsst', 'sst-dsst', 'dsst-sst')
3401 fout => grid%sst
3402 CASE DEFAULT
3403 IF (localpet.eq.0) THEN
3404 WRITE (cplout,10) trim(fld_string), trim(cinpname)
3405 END IF
3406 rc=esmf_rc_not_found
3407 IF (esmf_logfounderror(rctocheck=rc, &
3408 & msg=esmf_logerr_passthru, &
3409 & line=__line__, &
3410 & file=myfile)) THEN
3411 RETURN
3412 END IF
3413 END SELECT
3414!
3415! Set patch indices range.
3416!
3417 iminp=grid%sp31
3418 imaxp=grid%ep31
3419 jminp=grid%sp33
3420 jmaxp=grid%ep33
3421 IF (grid%ed31.eq.imaxp) THEN
3422 imaxp=imaxp-1
3423 END IF
3424 IF (grid%ed33.eq.jmaxp) THEN
3425 jmaxp=jmaxp-1
3426 END IF
3427!
3428!-----------------------------------------------------------------------
3429! Set WRF mask values:
3430!
3431! lakemask > 0: elsewhere (land, ocean) 1: lake
3432! landmask > 0: elsewhere (ocean, lakes) 1: land
3433!
3434!-----------------------------------------------------------------------
3435!
3436 lakevalue=1
3437 landvalue=1
3438!
3439!-----------------------------------------------------------------------
3440! If only one field is available, load field into output array at
3441! seawater points. Notice that Fout has the same precision as the
3442! WRF variable. It can be single or double precision.
3443!-----------------------------------------------------------------------
3444!
3445 IF (.not.got_dat.and.got_ocn) THEN
3446 myfmin= missing_dp
3447 myfmax=-missing_dp
3448 DO j=jminp,jmaxp
3449 DO i=iminp,imaxp
3450 IF ((int(grid%landmask(i,j)).ne.landvalue).and. &
3451 & (int(grid%lakemask(i,j)).ne.lakevalue)) THEN
3452 fout(i,j)=real(focn(i,j), kind(grid%sst))
3453 END IF
3454 ptr2d(i,j)=real(fout(i,j), dp)
3455 myfmin(1)=min(myfmin(1),fout(i,j))
3456 myfmax(1)=max(myfmax(1),fout(i,j))
3457 END DO
3458 END DO
3459 ELSE IF (.not.got_ocn.and.got_dat) THEN
3460 myfmin= missing_dp
3461 myfmax=-missing_dp
3462 DO j=jminp,jmaxp
3463 DO i=iminp,imaxp
3464 IF ((int(grid%landmask(i,j)).ne.landvalue).and. &
3465 & (int(grid%lakemask(i,j)).ne.lakevalue)) THEN
3466 fout(i,j)=real(fdat(i,j), kind(grid%sst))
3467 END IF
3468 ptr2d(i,j)=real(fout(i,j), dp)
3469 myfmin(1)=min(myfmin(1),fout(i,j))
3470 myfmax(1)=max(myfmax(1),fout(i,j))
3471 END DO
3472 END DO
3473 END IF
3474!
3475!-----------------------------------------------------------------------
3476! Otherwise, merge imported fields.
3477!-----------------------------------------------------------------------
3478!
3479 IF (got_ocn.and.got_dat) THEN
3480!
3481! Merge Focn and Fdat at sea-water and sea-ice points. Notice that
3482! the ESMF regridding will not fill unbounded interpolation points.
3483! Such grid cells still have the pointer initialized value MISSING_dp.
3484! The TOL_dp is used to identify such values. The user has full
3485! control of how the merging is done from the weights coefficients
3486! provided from input NetCDF file specified in "WeightsFile(atmos)".
3487!
3488 myfmin= missing_dp
3489 myfmax=-missing_dp
3490 DO j=jminp,jmaxp
3491 DO i=iminp,imaxp
3492 IF ((int(grid%landmask(i,j)).ne.landvalue).and. &
3493 & (int(grid%lakemask(i,j)).ne.lakevalue)) THEN
3494 IF (abs(fdat(i,j)).lt.tol_dp) THEN
3495 myfmin(2)=min(myfmin(2),fdat(i,j))
3496 myfmax(2)=max(myfmax(2),fdat(i,j))
3497 fval=fdat(i,j) ! initialize with DATA
3498 IF (abs(focn(i,j)).lt.tol_dp) THEN
3499 myfmin(1)=min(myfmin(1),focn(i,j))
3500 myfmax(1)=max(myfmax(1),focn(i,j))
3501 fval=weights(iatmos)%Cdat(i,j)*fval+ &
3502 & weights(iatmos)%Cesm(i,j)*focn(i,j)
3503 END IF
3504 fout(i,j)=real(fval, kind(grid%sst))
3505 ptr2d(i,j)=real(fval, dp)
3506 myfmin(3)=min(myfmin(3),fval)
3507 myfmax(3)=max(myfmax(3),fval)
3508 END IF
3509 END IF
3510 END DO
3511 END DO
3512 END IF
3513!
3514! Get merged field minimun and maximum values.
3515!
3516 IF (got_ocn.and.got_dat) THEN
3517 ic=3
3518 ELSE
3519 ic=1
3520 END IF
3521 CALL esmf_vmallreduce (vm, &
3522 & senddata=myfmin, &
3523 & recvdata=fmin, &
3524 & count=ic, &
3525 & reduceflag=esmf_reduce_min, &
3526 & rc=rc)
3527 IF (esmf_logfounderror(rctocheck=rc, &
3528 & msg=esmf_logerr_passthru, &
3529 & line=__line__, &
3530 & file=myfile)) THEN
3531 RETURN
3532 END IF
3533!
3534 CALL esmf_vmallreduce (vm, &
3535 & senddata=myfmax, &
3536 & recvdata=fmax, &
3537 & count=ic, &
3538 & reduceflag=esmf_reduce_max, &
3539 & rc=rc)
3540 IF (esmf_logfounderror(rctocheck=rc, &
3541 & msg=esmf_logerr_passthru, &
3542 & line=__line__, &
3543 & file=myfile)) THEN
3544 RETURN
3545 END IF
3546!
3547! Report merged field information.
3548!
3549 IF (got_ocn.and.got_dat) THEN
3550 IF ((debuglevel.ge.0).and.(localpet.eq.0)) THEN
3551 WRITE (cplout,20) trim(fld_string), &
3552 & trim(time_currentstring), ng, &
3553 & fmin(1), fmax(1), &
3554 & fmin(2), fmax(2), &
3555 & fmin(3), fmax(3)
3556 END IF
3557 ELSE
3558 IF ((debuglevel.ge.0).and.(localpet.eq.0)) THEN
3559 WRITE (cplout,30) fmin(1), fmax(1)
3560 END IF
3561 END IF
3562!
3563! Debugging: write out merged field into a NetCDF file.
3564!
3565 IF ((debuglevel.ge.3).and.any(debugwrite)) THEN
3566 WRITE (ofile,40) ng, trim(fld_string), &
3567 & year, month, day, hour, minutes, seconds
3568 CALL esmf_fieldwrite (fmerge, &
3569 & trim(ofile), &
3570 & overwrite=.true., &
3571 & rc=rc)
3572 IF (esmf_logfounderror(rctocheck=rc, &
3573 & msg=esmf_logerr_passthru, &
3574 & line=__line__, &
3575 & file=myfile)) THEN
3576 RETURN
3577 END IF
3578 END IF
3579!
3580! Nullify the pointers to ensure that it does not point to a random
3581! part in the memory.
3582!
3583 IF (associated(ptr2d)) nullify (ptr2d)
3584 IF (associated(fout )) nullify (fout)
3585!
3586! Destroy merged field array.
3587!
3588 CALL esmf_fielddestroy (fmerge, &
3589 & nogarbage=.false., &
3590 & rc=rc)
3591 IF (esmf_logfounderror(rctocheck=rc, &
3592 & msg=esmf_logerr_passthru, &
3593 & line=__line__, &
3594 & file=myfile)) THEN
3595 RETURN
3596 END IF
3597!
3598 IF (esm_track) THEN
3599 WRITE (trac,'(a,a,i0)') '<== Exiting WRF_ProcessImport_scalar',&
3600 & ', PET', petrank
3601 FLUSH (trac)
3602 END IF
3603 IF (debuglevel.gt.0) FLUSH (cplout)
3604!
3605 10 FORMAT (/,5x,'WRF_ProcessImport - ', &
3606 & 'unable to find option to import: ',a, &
3607 & /,25x,'check ''Import(atmos)'' in input script: ',a)
3608 20 FORMAT (1x,' WRF_ProcessImport - ESMF: merging field ''',a,'''', &
3609 & t72,a,2x,'Grid ',i2.2, &
3610 & /,19x,'(OcnMin = ', 1p,e15.8,0p, &
3611 & ' OcnMax = ', 1p,e15.8,0p,')', &
3612 & /,19x,'(DatMin = ', 1p,e15.8,0p, &
3613 & ' DatMax = ', 1p,e15.8,0p,')', &
3614 & /,19x,'(OutMin = ', 1p,e15.8,0p, &
3615 & ' OutMax = ', 1p,e15.8,0p,')')
3616 30 FORMAT (19x, '(OutMin = ', 1p,e15.8,0p, &
3617 & ' OutMax = ', 1p,e15.8,0p,') WRF_ProcessImport')
3618 40 FORMAT ('wrf_',i2.2,'_merged_',a,'_',i4.4,2('-',i2.2),'_', &
3619 & i2.2,2('.',i2.2),'.nc')
3620!
3621 RETURN
character(len(sinp)) function, public lowercase(sinp)
Definition strings.F:531

References mod_esmf_esm::cinpname, mod_esmf_esm::cplout, mod_esmf_esm::debuglevel, mod_esmf_esm::esm_track, mod_esmf_esm::iatmos, strings_mod::lowercase(), mod_esmf_esm::missing_dp, mod_esmf_esm::models, mod_esmf_esm::petrank, mod_esmf_esm::tol_dp, mod_esmf_esm::trac, and mod_esmf_esm::weights.

Here is the call graph for this function:

◆ wrf_processimport_vector()

subroutine esmf_wrf_mod::wrf_processimport::wrf_processimport_vector ( type (domain), pointer grid,
type (esmf_gridcomp) model,
logical, dimension(2,ubk), intent(in) got,
integer, dimension(2,ubk), intent(in) ifield,
character (len=*), dimension(2,ubk), intent(in) fieldname,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) ubk,
real (dp), dimension(lbi:ubi,lbj:ubj,ubk), intent(in) focn,
real (dp), dimension(lbi:ubi,lbj:ubj,ubk), intent(in) fdat,
integer, intent(out) rc )

Definition at line 3624 of file esmf_atm_wrf.h.

3629!
3630!=======================================================================
3631! !
3632! If both import vector components Focn and Fdat are avaliable, it !
3633! merges its values. Otherwise, it loads available data into ouput !
3634! component, Uout and Vout. Only sea-water or sea-ice points are !
3635! processed. !
3636! !
3637! The merging is needed when the atmosphere and ocean grids are !
3638! incongruent. The DATA model provides values on those grid points !
3639! not covered by the OCEAN model. !
3640! !
3641! On Input: !
3642! !
3643! grid WRF grid state structure (TYPE domain) !
3644! model Gridded component object (TYPE ESMF_GridComp) !
3645! got Switches indicating source and availability of !
3646! imported data (logical array): !
3647! got(1,1) OCEAN U-component switch (T/F) !
3648! got(1,2) OCEAN V-component switch (T/F) !
3649! got(2,1) DATA U-component switch (T/F) !
3650! got(2,2) DATA V-component switch (T/F) !
3651! ifield Imported vector component indices (integer array) !
3652! ifield(1,1) OCEAN U-component field index !
3653! ifield(1,2) OCEAN V-component field index !
3654! ifield(2,1) DATA U-component field index !
3655! ifield(2,2) DATA V-component field index !
3656! FieldName Vector components shortname (string array) !
3657! ifield Imported vector component indices (integer array) !
3658! LBi I-dimension lower bound (integer) !
3659! UBi I-dimension upper bound (integer) !
3660! LBj J-dimension lower bound (integer) !
3661! UBj J-dimension upper bound (integer) !
3662! UBk K-dimension storing vector component (integer; UBk=2) !
3663! Focn Imported vector from OCEAN model (3D real array) !
3664! Fdat Imported vector from DATA model (3D real array) !
3665! !
3666! On Output: !
3667! !
3668! rc Return code (integer) !
3669! !
3670!=======================================================================
3671!
3672 USE module_domain, ONLY : domain
3673 USE strings_mod, ONLY : lowercase
3674!
3675! Imported variable declarations.
3676!
3677 integer, intent(in) :: LBi, UBi, LBj, UBj, UBk
3678 integer, intent(in) :: ifield(2,UBk)
3679 integer, intent(out) :: rc
3680!
3681 logical, intent(in) :: got(2,UBk)
3682!
3683 real (dp), intent(in) :: Focn(LBi:UBi,LBj:UBj,UBk)
3684 real (dp), intent(in) :: Fdat(LBi:UBi,LBj:UBj,UBk)
3685!
3686 character (len=*), intent(in) :: FieldName(2,UBk)
3687!
3688 TYPE (domain), pointer :: grid
3689 TYPE (ESMF_GridComp) :: model
3690!
3691! Local variable declarations.
3692!
3693 logical :: got_dat, got_ocn
3694 logical :: DebugWrtU(2) = (/ .false., .false. /)
3695 logical :: DebugWrtv(2) = (/ .false., .false. /)
3696!
3697 integer :: i, ic, is, j, ng
3698 integer :: year, month, day, hour, minutes, seconds, sN, SD
3699 integer :: LakeValue, LandValue
3700 integer :: localDE, localDEcount, localPET, PETcount
3701 integer :: IminP, ImaxP, JminP, JmaxP
3702!
3703 real (dp) :: Fseconds, TimeInDays, Time_Current
3704
3705 real (dp) :: MyUmax(3), MyUmin(3), Umin(3), Umax(3), Uval
3706 real (dp) :: MyVmax(3), MyVmin(3), Vmin(3), Vmax(3), Vval
3707!
3708 real (dp), parameter :: MaxOcnVelocity = 10.0_dp ! m/s
3709!
3710 real (dp), pointer :: ptrU2d(:,:) => null()
3711 real (dp), pointer :: ptrV2d(:,:) => null()
3712!
3713 real (KIND(grid%uoce)), pointer :: Uout(:,:) => null()
3714 real (KIND(grid%voce)), pointer :: Vout(:,:) => null()
3715!
3716 character (len=22 ) :: Time_CurrentString
3717
3718 character (len=*), parameter :: MyFile = &
3719 & __FILE__//", WRF_ProcessImport_vector"
3720!
3721 character (ESMF_MAXSTR) :: cname, ofile, U_string, V_string
3722!
3723 TYPE (ESMF_ArraySpec) :: arraySpec2d
3724 TYPE (ESMF_Clock) :: clock
3725 TYPE (ESMF_Field) :: Umerge, Vmerge
3726 TYPE (ESMF_StaggerLoc) :: staggerLoc
3727 TYPE (ESMF_Time) :: CurrentTime
3728 TYPE (ESMF_VM) :: vm
3729!
3730!-----------------------------------------------------------------------
3731! Initialize return code flag to success state (no error).
3732!-----------------------------------------------------------------------
3733!
3734 IF (esm_track) THEN
3735 WRITE (trac,'(a,a,i0)') '==> Entering WRF_ProcessImport_vector',&
3736 & ', PET', petrank
3737 FLUSH (trac)
3738 END IF
3739 rc=esmf_success
3740!
3741!-----------------------------------------------------------------------
3742! Get information about the gridded component.
3743!-----------------------------------------------------------------------
3744!
3745 CALL esmf_gridcompget (model, &
3746 & clock=clock, &
3747 & localpet=localpet, &
3748 & petcount=petcount, &
3749 & vm=vm, &
3750 & name=cname, &
3751 & rc=rc)
3752 IF (esmf_logfounderror(rctocheck=rc, &
3753 & msg=esmf_logerr_passthru, &
3754 & line=__line__, &
3755 & file=myfile)) THEN
3756 RETURN
3757 END IF
3758 ng=grid%grid_id
3759!
3760! Get number of local decomposition elements (DEs). Usually, a single
3761! DE is associated with each Persistent Execution Thread (PETs). Thus,
3762! localDEcount=1.
3763!
3764 CALL esmf_gridget (models(iatmos)%grid(ng), &
3765 & localdecount=localdecount, &
3766 & rc=rc)
3767 IF (esmf_logfounderror(rctocheck=rc, &
3768 & msg=esmf_logerr_passthru, &
3769 & line=__line__, &
3770 & file=myfile)) THEN
3771 RETURN
3772 END IF
3773!
3774! Get current time.
3775!
3776 CALL esmf_clockget (clock, &
3777 & currtime=currenttime, &
3778 & rc=rc)
3779 IF (esmf_logfounderror(rctocheck=rc, &
3780 & msg=esmf_logerr_passthru, &
3781 & line=__line__, &
3782 & file=myfile)) THEN
3783 RETURN
3784 END IF
3785!
3786 CALL esmf_timeget (currenttime, &
3787 & yy=year, &
3788 & mm=month, &
3789 & dd=day, &
3790 & h =hour, &
3791 & m =minutes, &
3792 & s =seconds, &
3793 & sn=sn, &
3794 & sd=sd, &
3795 & rc=rc)
3796 IF (esmf_logfounderror(rctocheck=rc, &
3797 & msg=esmf_logerr_passthru, &
3798 & line=__line__, &
3799 & file=myfile)) THEN
3800 RETURN
3801 END IF
3802!
3803 CALL esmf_timeget (currenttime, &
3804 & s_r8=time_current, &
3805 & timestring=time_currentstring, &
3806 & rc=rc)
3807 IF (esmf_logfounderror(rctocheck=rc, &
3808 & msg=esmf_logerr_passthru, &
3809 & line=__line__, &
3810 & file=myfile)) THEN
3811 RETURN
3812 END IF
3813 fseconds=real(seconds,dp)+real(sn,dp)/real(sd,dp)
3814 timeindays=time_current/86400.0_dp
3815 is=index(time_currentstring, 'T') ! remove 'T' in
3816 IF (is.gt.0) time_currentstring(is:is)=' ' ! ISO 8601 format
3817!
3818!-----------------------------------------------------------------------
3819! Create merged vector components arrays.
3820!-----------------------------------------------------------------------
3821!
3822! Set a 2D floating-point array descriptor.
3823!
3824 CALL esmf_arrayspecset (arrayspec2d, &
3825 & typekind=esmf_typekind_r8, &
3826 & rank=2, &
3827 & rc=rc)
3828 IF (esmf_logfounderror(rctocheck=rc, &
3829 & msg=esmf_logerr_passthru, &
3830 & line=__line__, &
3831 & file=myfile)) THEN
3832 RETURN
3833 END IF
3834!
3835! Create 2D merge vector components from the Grid and arraySpec.
3836!
3837 got_ocn=got(1,1).and.got(1,2)
3838 got_dat=got(2,1).and.got(2,2)
3839!
3840 IF (.not.got_dat.and.got_ocn) THEN
3841 debugwrtu(1)=models(iatmos)%ImportField(ifield(1,1))%debug_write
3842 debugwrtv(1)=models(iatmos)%ImportField(ifield(1,2))%debug_write
3843 u_string=trim(fieldname(1,1))
3844 v_string=trim(fieldname(1,2))
3845 ELSE IF (.not.got_ocn.and.got_dat) THEN
3846 debugwrtu(2)=models(iatmos)%ImportField(ifield(2,1))%debug_write
3847 debugwrtv(2)=models(iatmos)%ImportField(ifield(2,2))%debug_write
3848 u_string=trim(fieldname(2,1))
3849 v_string=trim(fieldname(2,2))
3850 ELSE IF (got_ocn.and.got_dat) THEN
3851 debugwrtu(1)=models(iatmos)%ImportField(ifield(1,1))%debug_write
3852 debugwrtv(1)=models(iatmos)%ImportField(ifield(1,2))%debug_write
3853 debugwrtu(2)=models(iatmos)%ImportField(ifield(2,1))%debug_write
3854 debugwrtv(2)=models(iatmos)%ImportField(ifield(2,2))%debug_write
3855 u_string=trim(fieldname(1,1))//'-'//trim(fieldname(1,2))
3856 v_string=trim(fieldname(2,1))//'-'//trim(fieldname(2,2))
3857 END IF
3858 staggerloc=esmf_staggerloc_center
3859!
3860 umerge=esmf_fieldcreate(models(iatmos)%grid(ng), &
3861 & arrayspec2d, &
3862 & staggerloc=staggerloc, &
3863 & name=trim(u_string), &
3864 & rc=rc)
3865 IF (esmf_logfounderror(rctocheck=rc, &
3866 & msg=esmf_logerr_passthru, &
3867 & line=__line__, &
3868 & file=myfile)) THEN
3869 RETURN
3870 END IF
3871!
3872 vmerge=esmf_fieldcreate(models(iatmos)%grid(ng), &
3873 & arrayspec2d, &
3874 & staggerloc=staggerloc, &
3875 & name=trim(v_string), &
3876 & rc=rc)
3877 IF (esmf_logfounderror(rctocheck=rc, &
3878 & msg=esmf_logerr_passthru, &
3879 & line=__line__, &
3880 & file=myfile)) THEN
3881 RETURN
3882 END IF
3883!
3884! Get pointers for merge arrays.
3885!
3886 CALL esmf_fieldget (umerge, &
3887 & farrayptr=ptru2d, &
3888 & rc=rc)
3889 IF (esmf_logfounderror(rctocheck=rc, &
3890 & msg=esmf_logerr_passthru, &
3891 & line=__line__, &
3892 & file=myfile)) THEN
3893 RETURN
3894 END IF
3895 ptru2d=missing_dp
3896!
3897 CALL esmf_fieldget (vmerge, &
3898 & farrayptr=ptrv2d, &
3899 & rc=rc)
3900 IF (esmf_logfounderror(rctocheck=rc, &
3901 & msg=esmf_logerr_passthru, &
3902 & line=__line__, &
3903 & file=myfile)) THEN
3904 RETURN
3905 END IF
3906 ptrv2d=missing_dp
3907!
3908!-----------------------------------------------------------------------
3909! Create pointers to WRF vector component targets. This logic is
3910! awkward but needed since WRF could have different floating-point
3911! representations. Passing WRF variables as arguments to this routine
3912! doesn't work.
3913!-----------------------------------------------------------------------
3914!
3915 SELECT CASE (lowercase(trim(u_string)))
3916 CASE ('usur', 'dusur', 'usur-dusur', 'dusur-usur')
3917 uout => grid%uoce
3918 CASE DEFAULT
3919 IF (localpet.eq.0) THEN
3920 WRITE (cplout,10) trim(u_string), trim(cinpname)
3921 END IF
3922 rc=esmf_rc_not_found
3923 IF (esmf_logfounderror(rctocheck=rc, &
3924 & msg=esmf_logerr_passthru, &
3925 & line=__line__, &
3926 & file=myfile)) THEN
3927 RETURN
3928 END IF
3929 END SELECT
3930!
3931 SELECT CASE (lowercase(trim(v_string)))
3932 CASE ('vsur', 'dvsur', 'vsur-dvsur', 'dvsur-vsur')
3933 vout => grid%voce
3934 CASE DEFAULT
3935 IF (localpet.eq.0) THEN
3936 WRITE (cplout,10) trim(v_string), trim(cinpname)
3937 END IF
3938 rc=esmf_rc_not_found
3939 IF (esmf_logfounderror(rctocheck=rc, &
3940 & msg=esmf_logerr_passthru, &
3941 & line=__line__, &
3942 & file=myfile)) THEN
3943 RETURN
3944 END IF
3945 END SELECT
3946!
3947! Set patch indices range.
3948!
3949 iminp=grid%sp31
3950 imaxp=grid%ep31
3951 jminp=grid%sp33
3952 jmaxp=grid%ep33
3953 IF (grid%ed31.eq.imaxp) THEN
3954 imaxp=imaxp-1
3955 END IF
3956 IF (grid%ed33.eq.jmaxp) THEN
3957 jmaxp=jmaxp-1
3958 END IF
3959!
3960!-----------------------------------------------------------------------
3961! Set WRF mask values:
3962!
3963! lakemask > 0: elsewhere (land, ocean) 1: lake
3964! landmask > 0: elsewhere (ocean, lakes) 1: land
3965!
3966!-----------------------------------------------------------------------
3967!
3968 lakevalue=1
3969 landvalue=1
3970!
3971!-----------------------------------------------------------------------
3972! If only one vector source is available, load field into output array
3973! at seawater points. Notice that Uout and Vout has the same precision
3974! as the WRF variable. It can be single or double precision.
3975!-----------------------------------------------------------------------
3976!
3977 IF (.not.got_dat.and.got_ocn) THEN
3978 myumin= missing_dp
3979 myumax=-missing_dp
3980 myvmin= missing_dp
3981 myvmax=-missing_dp
3982 DO j=jminp,jmaxp
3983 DO i=iminp,imaxp
3984 IF ((int(grid%landmask(i,j)).ne.landvalue).and. &
3985 & (int(grid%lakemask(i,j)).ne.lakevalue)) THEN
3986 uout(i,j)=real(focn(i,j,1), kind(grid%uoce))
3987 vout(i,j)=real(focn(i,j,2), kind(grid%voce))
3988 END IF
3989 ptru2d(i,j)=real(uout(i,j), dp)
3990 ptrv2d(i,j)=real(vout(i,j), dp)
3991 myumin(1)=min(myumin(1),uout(i,j))
3992 myumax(1)=max(myumax(1),uout(i,j))
3993 myvmin(1)=min(myvmin(1),vout(i,j))
3994 myvmax(1)=max(myvmax(1),vout(i,j))
3995 END DO
3996 END DO
3997 ELSE IF (.not.got_ocn.and.got_dat) THEN
3998 myumin= missing_dp
3999 myumax=-missing_dp
4000 myvmin= missing_dp
4001 myvmax=-missing_dp
4002 DO j=jminp,jmaxp
4003 DO i=iminp,imaxp
4004 IF ((int(grid%landmask(i,j)).ne.landvalue).and. &
4005 & (int(grid%lakemask(i,j)).ne.lakevalue)) THEN
4006 uout(i,j)=real(fdat(i,j,1), kind(grid%uoce))
4007 vout(i,j)=real(fdat(i,j,2), kind(grid%voce))
4008 END IF
4009 ptru2d(i,j)=real(uout(i,j), dp)
4010 ptrv2d(i,j)=real(vout(i,j), dp)
4011 myumin(1)=min(myumin(1),uout(i,j))
4012 myumax(1)=max(myumax(1),uout(i,j))
4013 myvmin(1)=min(myvmin(1),vout(i,j))
4014 myvmax(1)=max(myvmax(1),vout(i,j))
4015 END DO
4016 END DO
4017 END IF
4018!
4019!-----------------------------------------------------------------------
4020! Otherwise, merge imported vector components.
4021!-----------------------------------------------------------------------
4022!
4023 IF (got_ocn.and.got_dat) THEN
4024!
4025! Merge Focn and Fdat at sea-water and sea-ice points. Notice that
4026! the ESMF regridding will not fill unbounded interpolation points.
4027! Such grid cells still have the pointer initialized value MISSING_dp.
4028! The TOL_dp is used to identify such values. The user has full
4029! control of how the merging is done from the weights coefficients
4030! provided from input NetCDF file specified in "WeightsFile(atmos)".
4031!
4032 myumin= missing_dp
4033 myumax=-missing_dp
4034 myvmin= missing_dp
4035 myvmax=-missing_dp
4036 DO j=jminp,jmaxp
4037 DO i=iminp,imaxp
4038 IF ((int(grid%landmask(i,j)).ne.landvalue).and. &
4039 & (int(grid%lakemask(i,j)).ne.lakevalue)) THEN
4040 IF ((abs(fdat(i,j,1)).lt.maxocnvelocity).and. &
4041 & (abs(fdat(i,j,2)).lt.maxocnvelocity)) THEN
4042 myumin(2)=min(myumin(2),fdat(i,j,1))
4043 myumax(2)=max(myumax(2),fdat(i,j,1))
4044 myvmin(2)=min(myvmin(2),fdat(i,j,2))
4045 myvmax(2)=max(myvmax(2),fdat(i,j,2))
4046 uval=fdat(i,j,1) ! initialize with DATA
4047 vval=fdat(i,j,2)
4048 IF ((abs(focn(i,j,1)).lt.maxocnvelocity).and. &
4049 & (abs(focn(i,j,2)).lt.maxocnvelocity)) THEN
4050 myumin(1)=min(myumin(1),focn(i,j,1))
4051 myumax(1)=max(myumax(1),focn(i,j,1))
4052 myvmin(1)=min(myvmin(1),focn(i,j,2))
4053 myvmax(1)=max(myvmax(1),focn(i,j,2))
4054 uval=weights(iatmos)%Cdat(i,j)*uval+ &
4055 & weights(iatmos)%Cesm(i,j)*focn(i,j,1)
4056 vval=weights(iatmos)%Cdat(i,j)*vval+ &
4057 & weights(iatmos)%Cesm(i,j)*focn(i,j,2)
4058 END IF
4059 uout(i,j)=real(uval, kind(grid%uoce))
4060 vout(i,j)=real(vval, kind(grid%voce))
4061 ptru2d(i,j)=real(uval, dp)
4062 ptrv2d(i,j)=real(vval, dp)
4063 myumin(3)=min(myumin(3),uval)
4064 myumax(3)=max(myumax(3),uval)
4065 myvmin(3)=min(myvmin(3),vval)
4066 myvmax(3)=max(myvmax(3),vval)
4067 END IF
4068 END IF
4069 END DO
4070 END DO
4071 END IF
4072!
4073! Get merged vector components minimun and maximum values.
4074!
4075 IF (got_ocn.and.got_dat) THEN
4076 ic=3
4077 ELSE
4078 ic=1
4079 END IF
4080!
4081 CALL esmf_vmallreduce (vm, &
4082 & senddata=myumin, &
4083 & recvdata=umin, &
4084 & count=ic, &
4085 & reduceflag=esmf_reduce_min, &
4086 & rc=rc)
4087 IF (esmf_logfounderror(rctocheck=rc, &
4088 & msg=esmf_logerr_passthru, &
4089 & line=__line__, &
4090 & file=myfile)) THEN
4091 RETURN
4092 END IF
4093!
4094 CALL esmf_vmallreduce (vm, &
4095 & senddata=myumax, &
4096 & recvdata=umax, &
4097 & count=ic, &
4098 & reduceflag=esmf_reduce_max, &
4099 & rc=rc)
4100 IF (esmf_logfounderror(rctocheck=rc, &
4101 & msg=esmf_logerr_passthru, &
4102 & line=__line__, &
4103 & file=myfile)) THEN
4104 RETURN
4105 END IF
4106!
4107 CALL esmf_vmallreduce (vm, &
4108 & senddata=myvmin, &
4109 & recvdata=vmin, &
4110 & count=ic, &
4111 & reduceflag=esmf_reduce_min, &
4112 & rc=rc)
4113 IF (esmf_logfounderror(rctocheck=rc, &
4114 & msg=esmf_logerr_passthru, &
4115 & line=__line__, &
4116 & file=myfile)) THEN
4117 RETURN
4118 END IF
4119!
4120 CALL esmf_vmallreduce (vm, &
4121 & senddata=myvmax, &
4122 & recvdata=vmax, &
4123 & count=ic, &
4124 & reduceflag=esmf_reduce_max, &
4125 & rc=rc)
4126 IF (esmf_logfounderror(rctocheck=rc, &
4127 & msg=esmf_logerr_passthru, &
4128 & line=__line__, &
4129 & file=myfile)) THEN
4130 RETURN
4131 END IF
4132!
4133! Report merged vector components information.
4134!
4135 IF (got_ocn.and.got_dat) THEN
4136 IF ((debuglevel.ge.0).and.(localpet.eq.0)) THEN
4137 WRITE (cplout,20) trim(u_string), &
4138 & trim(time_currentstring), ng, &
4139 & umin(1), umax(1), &
4140 & umin(2), umax(2), &
4141 & umin(3), umax(3)
4142 WRITE (cplout,20) trim(v_string), &
4143 & trim(time_currentstring), ng, &
4144 & vmin(1), vmax(1), &
4145 & vmin(2), vmax(2), &
4146 & vmin(3), vmax(3)
4147 END IF
4148 ELSE
4149 IF ((debuglevel.ge.0).and.(localpet.eq.0)) THEN
4150 WRITE (cplout,30) umin(1), umax(1)
4151 WRITE (cplout,30) vmin(1), vmax(1)
4152 END IF
4153 END IF
4154!
4155! Debugging: write out merged vector components into a NetCDF file.
4156!
4157 IF ((debuglevel.ge.3).and.any(debugwrtu)) THEN
4158 WRITE (ofile,40) ng, trim(u_string), &
4159 & year, month, day, hour, minutes, seconds
4160 CALL esmf_fieldwrite (umerge, &
4161 & trim(ofile), &
4162 & overwrite=.true., &
4163 & rc=rc)
4164 IF (esmf_logfounderror(rctocheck=rc, &
4165 & msg=esmf_logerr_passthru, &
4166 & line=__line__, &
4167 & file=myfile)) THEN
4168 RETURN
4169 END IF
4170 END IF
4171!
4172 IF ((debuglevel.ge.3).and.any(debugwrtv)) THEN
4173 WRITE (ofile,40) ng, trim(v_string), &
4174 & year, month, day, hour, minutes, seconds
4175 CALL esmf_fieldwrite (vmerge, &
4176 & trim(ofile), &
4177 & overwrite=.true., &
4178 & rc=rc)
4179 IF (esmf_logfounderror(rctocheck=rc, &
4180 & msg=esmf_logerr_passthru, &
4181 & line=__line__, &
4182 & file=myfile)) THEN
4183 RETURN
4184 END IF
4185 END IF
4186!
4187! Nullify the pointers to ensure that it does not point to a random
4188! part in the memory.
4189!
4190 IF (associated(ptru2d)) nullify (ptru2d)
4191 IF (associated(ptrv2d)) nullify (ptrv2d)
4192 IF (associated(uout )) nullify (uout)
4193 IF (associated(vout )) nullify (vout)
4194!
4195! Destroy merged vector arrays.
4196!
4197 CALL esmf_fielddestroy (umerge, &
4198 & nogarbage=.false., &
4199 & rc=rc)
4200 IF (esmf_logfounderror(rctocheck=rc, &
4201 & msg=esmf_logerr_passthru, &
4202 & line=__line__, &
4203 & file=myfile)) THEN
4204 RETURN
4205 END IF
4206!
4207 CALL esmf_fielddestroy (vmerge, &
4208 & nogarbage=.false., &
4209 & rc=rc)
4210 IF (esmf_logfounderror(rctocheck=rc, &
4211 & msg=esmf_logerr_passthru, &
4212 & line=__line__, &
4213 & file=myfile)) THEN
4214 RETURN
4215 END IF
4216!
4217 IF (esm_track) THEN
4218 WRITE (trac,'(a,a,i0)') '<== Exiting WRF_ProcessImport_vector',&
4219 & ', PET', petrank
4220 FLUSH (trac)
4221 END IF
4222 IF (debuglevel.gt.0) FLUSH (cplout)
4223!
4224 10 FORMAT (/,5x,'WRF_ProcessImport - ', &
4225 & 'unable to find option to import: ',a, &
4226 & /,25x,'check ''Import(atmos)'' in input script: ',a)
4227 20 FORMAT (1x,' WRF_ProcessImport - ESMF: merging field ''',a,'''', &
4228 & t72,a,2x,'Grid ',i2.2, &
4229 & /,19x,'(OcnMin = ', 1p,e15.8,0p, &
4230 & ' OcnMax = ', 1p,e15.8,0p,')', &
4231 & /,19x,'(DatMin = ', 1p,e15.8,0p, &
4232 & ' DatMax = ', 1p,e15.8,0p,')', &
4233 & /,19x,'(OutMin = ', 1p,e15.8,0p, &
4234 & ' OutMax = ', 1p,e15.8,0p,')')
4235 30 FORMAT (19x, '(OutMin = ', 1p,e15.8,0p, &
4236 & ' OutMax = ', 1p,e15.8,0p,') WRF_ProcessImport')
4237 40 FORMAT ('wrf_',i2.2,'_merged_',a,'_',i4.4,2('-',i2.2),'_', &
4238 & i2.2,2('.',i2.2),'.nc')
4239!
4240 RETURN

References mod_esmf_esm::cinpname, mod_esmf_esm::cplout, mod_esmf_esm::debuglevel, mod_esmf_esm::esm_track, mod_esmf_esm::iatmos, strings_mod::lowercase(), mod_esmf_esm::missing_dp, mod_esmf_esm::models, mod_esmf_esm::petrank, mod_esmf_esm::trac, and mod_esmf_esm::weights.

Here is the call graph for this function:

The documentation for this interface was generated from the following file: