Skip to content

Restored BLEND_QV_AT_TP option for Replay #1148

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
sdrabenh merged 3 commits into from
Sep 23, 2025
Merged

Restored BLEND_QV_AT_TP option for Replay #1148

Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
085de2e
Restored BLEND_QV_AT_TP option
mmanyin Aug 28, 2025
026f382
Merge branch 'develop' into feature/mmanyin/restore_blend_option
sdrabenh Sep 23, 2025
de7ede2
Merge branch 'develop' into feature/mmanyin/restore_blend_option
sdrabenh Sep 23, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
108 changes: 102 additions & 6 deletions GEOSmkiau_GridComp/GEOS_mkiauGridComp.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -96,6 +96,7 @@ subroutine SetServices ( GC, RC )
type (ESMF_Config) :: CF

logical :: BLEND_AT_PBL
logical :: BLEND_QV_AT_TP
#ifdef PYMKIAU_INTEGRATION
! IEEE trapping see below
logical :: halting_mode(5)
Expand Down Expand Up @@ -125,6 +126,11 @@ subroutine SetServices ( GC, RC )
call MAPL_GetResource(MAPL, BLEND_AT_PBL, LABEL="REPLAY_BLEND_AT_PBL:", default=.FALSE., RC=status)
VERIFY_(STATUS)

call MAPL_GetResource(MAPL, BLEND_QV_AT_TP, LABEL="REPLAY_BLEND_QV_AT_TP:", default=.FALSE., RC=status)
VERIFY_(STATUS)

if ( BLEND_AT_PBL ) BLEND_QV_AT_TP = .FALSE.

! Set the Run entry points (phase 1 for regular IAU and phase 2 for clearing
! --------------------------------------------------------------------------

Expand Down Expand Up @@ -251,6 +257,17 @@ subroutine SetServices ( GC, RC )
VERIFY_(STATUS)
endif

if( BLEND_QV_AT_TP ) then
call MAPL_AddImportSpec(GC, &
SHORT_NAME = 'TROPP_BLENDED', &
LONG_NAME = 'tropopause_pressure_based_on_blended_estimate', &
UNITS = 'Pa', &
DIMS = MAPL_DimsHorzOnly, &
VLOCATION = MAPL_VLocationNone, &
RC=STATUS )
VERIFY_(STATUS)
endif


! !EXPORT STATE:
! --------------
Expand Down Expand Up @@ -597,6 +614,7 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC )
real, pointer, dimension(:,:,:) :: pdum1 => null()
real, pointer, dimension(:,:,:) :: pdum2 => null()
real, pointer, dimension(:,:) :: blnpp => null()
real, pointer, dimension(:,:) :: tropp => null()

real, allocatable, dimension(:,:,:) :: du_fix
real, allocatable, dimension(:,:,:) :: dv_fix
Expand Down Expand Up @@ -662,6 +680,7 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC )
real :: FACP1, FACP0, FACM1, FACM2
real :: DAMPBEG, DAMPEND
logical :: BLEND_AT_PBL
logical :: BLEND_QV_AT_TP
integer :: i,j,L,n
integer :: nt,nvars,natts
integer :: nymd, nhms
Expand Down Expand Up @@ -894,6 +913,9 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC )
call MAPL_GetResource(MAPL, BLEND_AT_PBL, LABEL="REPLAY_BLEND_AT_PBL:", default=.FALSE., RC=status)
VERIFY_(STATUS)

call MAPL_GetResource(MAPL, BLEND_QV_AT_TP, LABEL="REPLAY_BLEND_QV_AT_TP:", default=.FALSE., RC=status)
VERIFY_(STATUS)

CREMAP = ESMF_UtilStringUpperCase(CREMAP)
FIXWIND = ESMF_UtilStringUpperCase(FIXWIND)
DOWINDFIX = trim(FIXWIND)=="YES"
Expand Down Expand Up @@ -2145,7 +2167,7 @@ subroutine handleANA_
! **** with option to blend QV specially, starting at tropopause. ****
! **********************************************************************

if( DAMPBEG.ne.DAMPEND .or. BLEND_AT_PBL ) then
if( DAMPBEG.ne.DAMPEND .or. BLEND_AT_PBL .or. BLEND_QV_AT_TP ) then

if(first .and. MAPL_AM_I_ROOT()) then
if(DAMPBEG.ne.DAMPEND) then
Expand All @@ -2158,9 +2180,21 @@ subroutine handleANA_
else
print *, 'No blending based on PBL'
endif
if(BLEND_QV_AT_TP) then
if(BLEND_AT_PBL) then
print *, 'Blending at PBL supercedes QV blending at TROPP'
else
print *, 'Blending ANA and BKG QV based on TROPP'
endif
else
print *, 'No blending of QV based on TROPP'
endif
print *
endif

! Enforce priority now that warning has been logged
if ( BLEND_AT_PBL ) BLEND_QV_AT_TP = .FALSE.

if( BLEND_AT_PBL ) then
allocate ( pdum1(IMbkg,JMbkg,1) )
allocate ( pdum2(IM, JM, 1) )
Expand All @@ -2178,13 +2212,30 @@ subroutine handleANA_
endif
blnpp => pdum2(:,:,1)
endif
if( BLEND_QV_AT_TP ) then
allocate ( pdum1(IMbkg,JMbkg,1) )
allocate ( pdum2(IM, JM, 1) )
pdum1=0.0

call MAPL_GetPointer(import, ptr2d, 'TROPP_BLENDED', RC=STATUS)
VERIFY_(STATUS)
pdum1(:,:,1) = ptr2d

if (trim(GRIDINC)=="ANA" .and. do_transforms) then
call mkiau_internal_state%bkg2ana_regridder%regrid(pdum1, pdum2, RC=STATUS)
VERIFY_(STATUS)
else
pdum2=pdum1
endif
tropp => pdum2(:,:,1)
endif

call blend ( ple_ana,u_ana,v_ana,t_ana,q_ana,o3_ana, &
ple_bkg,u_bkg,v_bkg,t_bkg,q_bkg,o3_bkg, &
im,jm,LMbkg, DAMPBEG,DAMPEND, BLEND_AT_PBL, &
blnpp=blnpp )
BLEND_QV_AT_TP, blnpp=blnpp, tropp=tropp )

if( BLEND_AT_PBL ) then
if( BLEND_AT_PBL .or. BLEND_QV_AT_TP ) then
deallocate ( pdum1 )
deallocate ( pdum2 )
endif
Expand Down Expand Up @@ -2628,24 +2679,30 @@ end subroutine CLEAR
subroutine blend ( plea,ua,va,ta,qa,oa, &
pleb,ub,vb,tb,qb,ob, &
im,jm,lm, pabove,pbelow, &
BLEND_AT_PBL, blnpp )
BLEND_AT_PBL, &
BLEND_QV_AT_TP, &
blnpp, tropp )

! Blends Anaylsis and Background values.
! This routine is called if pabove /= pbelow or BLEND_AT_PBL
! This routine is called if pabove /= pbelow or BLEND_AT_PBL or BLEND_QV_AT_TP
! ***************************************************************************

implicit none
integer, intent(IN) :: im,jm,lm
real, intent(IN) :: pabove,pbelow
logical, intent(IN) :: BLEND_AT_PBL
logical, intent(IN) :: BLEND_QV_AT_TP

! Background values
real, intent(IN) :: pleb(im,jm,lm+1)
real, intent(IN) :: ub(im,jm,lm)
real, intent(IN) :: vb(im,jm,lm)
real, intent(IN) :: tb(im,jm,lm)
real, intent(IN) :: qb(im,jm,lm)
real, intent(IN) :: ob(im,jm,lm)

! IN: Anaylsis values
! OUT: Blended values
real, intent(INOUT) :: plea(im,jm,lm+1)
real, intent(INOUT) :: ua(im,jm,lm)
real, intent(INOUT) :: va(im,jm,lm)
Expand All @@ -2654,6 +2711,7 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
real, intent(INOUT) :: oa(im,jm,lm)

real, intent(IN), optional, pointer :: blnpp(:,:) ! blending pressure when BLEND_AT_PBL is TRUE
real, intent(IN), optional, pointer :: tropp(:,:) ! Tropopause Pressure used when BLEND_QV_AT_TP is TRUE

! Locals
! ------
Expand All @@ -2670,6 +2728,9 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
real pabove_BL,pbelow_BL
real bl_press

real pabove_QV,pbelow_QV ! compute from tropp
real tp_press

real alf,eps,p
integer i,j,L

Expand Down Expand Up @@ -2709,7 +2770,7 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
ua(i,j,L) = ub(i,j,L) + alf*( ua(i,j,L)- ub(i,j,L) )
va(i,j,L) = vb(i,j,L) + alf*( va(i,j,L)- vb(i,j,L) )
oa(i,j,L) = ob(i,j,L) + alf*( oa(i,j,L)- ob(i,j,L) )
qa(i,j,L) = qb(i,j,L) + alf*( qa(i,j,L)- qb(i,j,L) )
IF (.NOT. BLEND_QV_AT_TP) qa(i,j,L) = qb(i,j,L) + alf*( qa(i,j,L)- qb(i,j,L) )
enddo
enddo
enddo
Expand Down Expand Up @@ -2779,6 +2840,41 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
enddo
endif

! Use Analysis values of water vapor in the troposphere
! Relax to Background values in the stratosphere
! -----------------------------------------------------
if ( BLEND_QV_AT_TP ) then
do j=1,jm
do i=1,im

IF ( tropp(i,j) == MAPL_UNDEF ) THEN
tp_press = 100.0 * 100.0 ! 100 hPa
ELSE
tp_press = tropp(i,j)
ENDIF

pabove_QV = tp_press * 0.5
pbelow_QV = tp_press * 1.0

do L=1,lm
p = 0.5*( pleb(i,j,L)+pleb(i,j,L+1) )
if( p.le.pabove_QV ) then
alf = 0.0 ! use the background value
else if( p.gt.pabove_QV .and. p.le.pbelow_QV ) then
alf = (LOG(p) -LOG(pabove_QV))/ &
(LOG(pbelow_QV)-LOG(pabove_QV))
else
alf = 1.0 ! use the analysis value
endif

qa(i,j,L) = qb(i,j,L) + alf*( qa(i,j,L)- qb(i,j,L) )

enddo

enddo
enddo
endif

return
end subroutine blend

Expand Down