Skip to content
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

Jump to bottom

Closes #70. Commit changes for PBL height #26

Open
wants to merge 1 commit into
base: develop
Choose a base branch
from
Open

Closes #70. Commit changes for PBL height #26

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
20 changes: 18 additions & 2 deletions GEOSaana_GridComp/GEOSgsi_Coupler/cplr_ensemble.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -397,6 +397,7 @@ subroutine put_geos_ens(this,grd,member,ntindex,pert,iret)
! 2011-10-01 todling - created for testing purposes
! 2012-05-15 el akkraoui/todling - overload for r4/r8; pass time index
! 2015-10-19 todling - add qi/ql/qr/qs and should work for all cases of CV/SV
! 2020-03-09 yzhu - add pblh
!
! input argument list:
! grd - structure variable containing information about grid
Expand Down Expand Up @@ -448,9 +449,9 @@ subroutine put_geos_ens(this,grd,member,ntindex,pert,iret)
type(gsi_bundle) flds

! Declare fields in file (should come from resource file)
integer(i_kind),parameter:: no2d=2
integer(i_kind),parameter:: no2d=3
integer(i_kind),dimension(no2d)::iptr2d
character(len=4), parameter :: ovars2d(no2d) = (/ 'ps ', 'z '/)
character(len=4), parameter :: ovars2d(no2d) = (/ 'ps ', 'z ', 'pblh'/)
integer(i_kind),parameter:: no3d=10
integer(i_kind),dimension(no3d)::iptr3d
character(len=5), parameter :: ovars3d(no3d) = (/ 'u ', 'v ',&
Expand Down Expand Up @@ -511,6 +512,21 @@ subroutine put_geos_ens(this,grd,member,ntindex,pert,iret)
! if(mype==0) write(6,*) myname,': z ens member not in incoming pert'
! call stop2(999)
! endif

! PBLH
call gsi_bundlegetpointer (pert,'pblh',ipnt,istatus)
if(istatus==0) then
call gsi_bundlegetpointer (flds,'pblh',optr2d,istatus)
if(iamsingle) then
optr2d=pert%r2(ipnt)%qr4
else
optr2d=pert%r2(ipnt)%qr8
endif
else
if(mype==0) write(6,*) myname,': pblh ens member not in incoming pert'
call stop2(999)
endif

! SF/VP->U/V
istatus=0
call gsi_bundlegetpointer (pert,'sf',isf,ier);istatus=ier+istatus
Expand Down
2 changes: 2 additions & 0 deletions GEOSaana_GridComp/GEOSgsi_Coupler/g5_pertmod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,6 +65,7 @@ module g5_pertmod
use stepon, only : ak ! GEOS-5 ADM/TLM pressure levels
use stepon, only : bk ! GEOS-5 ADM/TLM pressure levels
use stepon, only : ts
! use stepon, only : pblht
use stepon, only : oro
use stepon, only : job
use stepon, only : nstep
Expand Down Expand Up @@ -1494,6 +1495,7 @@ subroutine gsi2pgcm0_ ( nymd, nhms, xx, which, stat, &
write(fname,'(4a,i3.3)') trim(job), '.', trim(fnxgsi), '_', mycount
endif
call putpert ( job, nymd, nhms, xpert, fvpsasdt, nstep, &
! ak, bk, Ts, pblht, oro, ps, fname, vectype=vectype_, forceflip=.true. )
ak, bk, Ts, oro, ps, fname, vectype=vectype_, forceflip=.true. )
if(ierr/=0)then
stat = 90
Expand Down
60 changes: 55 additions & 5 deletions GEOSaana_GridComp/GEOSgsi_Coupler/geos_StateIO.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -456,19 +456,21 @@ subroutine gsi2pgcm1_ ( xx, xpert, sg, GSI_xGrid, which, stat )
! 15May2010 Todling Update to use GSI_Bundle
! 21Feb2011 Todling Adapt to work with MAPL-SimpleBundle
! 28Oct2013 Todling Rename p3d to prse
! 06Mar2020 Zhu Add pblh
!
!EOP
!-----------------------------------------------------------------------

character(len=*), parameter :: myname_ = myname//'*gsi2pgcm_'

real(r_kind), allocatable, dimension(:,:,:) :: sub_tv,sub_u,sub_v,sub_q,sub_delp
real(r_kind), allocatable, dimension(:,:,:) :: sub_tv,sub_u,sub_v,sub_q,sub_delp,sub_pblh3
real(r_kind), allocatable, dimension(:,:,:) :: sub_oz,sub_cw
real(r_kind), allocatable, dimension(:,:,:) :: sub_qi,sub_ql,sub_qr,sub_qs
real(r_kind), allocatable, dimension(:,:) :: sub_ps
real(r_kind), allocatable, dimension(:,:) :: sub_ps,sub_pblh
real(4), allocatable, dimension(:,:,:) :: wktmp

integer(i_kind) i,j,k,kk,ijk,ij
integer(i_kind) i_u,i_v,i_t,i_q,i_oz,i_cw,i_prse,i_p,i_dp
integer(i_kind) i_u,i_v,i_t,i_q,i_oz,i_cw,i_prse,i_p,i_dp,i_pblh
integer(i_kind) i_qi,i_ql,i_qr,i_qs
integer(i_kind) ierr,istatus,rc,status
character(len=255) :: whatin
Expand Down Expand Up @@ -611,6 +613,26 @@ subroutine gsi2pgcm1_ ( xx, xpert, sg, GSI_xGrid, which, stat )
sub_delp(:,:,k) = sub_ps ! just a trick to use gsi2pert_ w/o interface change
enddo

call gsi_bundlegetpointer(xx,'pblh', i_pblh , istatus)
if (i_pblh>0) then
allocate(sub_pblh(sg%lat2,sg%lon2), stat=ierr )
if ( ierr/=0 ) then
stat = 91
if(mype==ROOT) print*, trim(myname_), ': Alloc(sub_pblh)'
return
end if
do j=1,sg%lon2
do i=1,sg%lat2
sub_pblh(i,j) = xx%r2(i_pblh)%q(i,j)
enddo
enddo
allocate(sub_pblh3(sg%lat2,sg%lon2,nsig), stat=ierr )
sub_pblh3=zero
do k=1,nsig
sub_pblh3(:,:,k) = sub_pblh ! just a trick to use gsi2pert_ w/o interface change
enddo
endif

! Gather from GSI subdomains/Scatter to GCM
! -----------------------------------------
i_p = MAPL_SimpleBundleGetIndex ( xpert, 'ps' , 2, rc=status )
Expand Down Expand Up @@ -648,6 +670,13 @@ subroutine gsi2pgcm1_ ( xx, xpert, sg, GSI_xGrid, which, stat )
i_qs= MAPL_SimpleBundleGetIndex ( xpert, 'qstot',3, rc=status )
call gsi2pert_ ( sub_qs, xpert%r3(i_qs)%qr4, ierr )
endif
if (i_pblh>0) then
i_pblh= MAPL_SimpleBundleGetIndex ( xpert, 'PBLH',2, rc=status )
allocate(wktmp(sg%lat2,sg%lon2,nsig), stat=ierr )
wktmp = zero
call gsi2pert_ ( sub_pblh3, wktmp, ierr )
xpert%r2(i_pblh)%qr4 = wktmp(:,:,1)
endif

! Build surface pressure perturbation - output purposes
! -----------------------------------------------------
Expand All @@ -670,6 +699,9 @@ subroutine gsi2pgcm1_ ( xx, xpert, sg, GSI_xGrid, which, stat )
if (allocated(sub_ql)) deallocate(sub_ql)
if (allocated(sub_qr)) deallocate(sub_qr)
if (allocated(sub_qs)) deallocate(sub_qs)
if (allocated(sub_pblh)) deallocate(sub_pblh)
if (allocated(sub_pblh3)) deallocate(sub_pblh3)
if (allocated(wktmp)) deallocate(wktmp)

CONTAINS

Expand Down Expand Up @@ -1161,6 +1193,7 @@ subroutine gcm2gsi1_ ( xpert, xx, sg, GSI_xGrid, stat )
! 20Feb2011 Todling Adapt to work with MAPL-SimpleBundle
! 30Nov2014 Todling Add some variable exceptions (qi/ql/ph/ts)
! 06May2020 Todling Shuffled pointer check and alloc for increased flexibility
! 06Mar2020 Zhu Add pblh
!
!EOP
!-----------------------------------------------------------------------
Expand All @@ -1171,11 +1204,11 @@ subroutine gcm2gsi1_ ( xpert, xx, sg, GSI_xGrid, stat )
real(r_kind), allocatable, dimension(:,:,:) :: sub_u,sub_v,sub_q,sub_tv
real(r_kind), allocatable, dimension(:,:,:) :: sub_oz,sub_ci,sub_cl
real(r_kind), allocatable, dimension(:,:,:) :: sub_cr,sub_cs
real(r_kind), allocatable, dimension(:,:) :: sub_ps,sub_ph,sub_ts
real(r_kind), allocatable, dimension(:,:) :: sub_ps,sub_ph,sub_ts,sub_pblh
real(r_kind), allocatable, dimension(:,:,:) :: sub_3d

integer(i_kind) i,j,k,ij,ijk,id,ng_d,ng_s,ni,nj
integer(i_kind) i_u,i_v,i_t,i_q,i_oz,i_cw,i_ps,i_ts,i_tsen,i_ph
integer(i_kind) i_u,i_v,i_t,i_q,i_oz,i_cw,i_ps,i_ts,i_tsen,i_ph,i_pblh
integer(i_kind) i_ci,i_cl,i_cr,i_cs
integer(i_kind) j_cr,j_cs
integer(i_kind) ierr,istatus,status
Expand All @@ -1200,6 +1233,7 @@ subroutine gcm2gsi1_ ( xpert, xx, sg, GSI_xGrid, stat )
i_ps = MAPL_SimpleBundleGetIndex ( xpert, 'ps' , 2, rc=status )
i_ts = MAPL_SimpleBundleGetIndex ( xpert, 'ts' , 2, rc=status )
i_ph = MAPL_SimpleBundleGetIndex ( xpert, 'phis' , 2, rc=status )
i_pblh= MAPL_SimpleBundleGetIndex (xpert, 'PBLH' , 2, rc=status )
if (i_u>0.and.i_v>0) then
allocate(sub_u(sg%lat2,sg%lon2,nsig), stat=ierr )
call pert2gsi_ ( xpert%r3(i_u)%qr4 , sub_u , ierr )
Expand Down Expand Up @@ -1238,6 +1272,10 @@ subroutine gcm2gsi1_ ( xpert, xx, sg, GSI_xGrid, stat )
allocate(sub_ts(sg%lat2,sg%lon2))
call pert2gsi2d_ ( xpert%r2(i_ts)%qr4, sub_ts, ierr )
endif
if (i_pblh>0) then
allocate(sub_pblh(sg%lat2,sg%lon2))
call pert2gsi2d_ ( xpert%r2(i_pblh)%qr4, sub_pblh, ierr )
endif
if ( ierr/=0 ) then
stat = 99
if(mype==ROOT) print*, trim(myname_), ': unfinished convertion ...'
Expand Down Expand Up @@ -1428,6 +1466,18 @@ subroutine gcm2gsi1_ ( xpert, xx, sg, GSI_xGrid, stat )
deallocate(sub_ts)
endif

if (allocated(sub_pblh)) then
call gsi_bundlegetpointer(xx,'pblh', i_pblh, ierr)
if (i_pblh>0) then
do j=1,sg%lon2
do i=1,sg%lat2
xx%r2(i_pblh)%q(i,j)= sub_pblh(i,j)
enddo
enddo
endif
deallocate(sub_pblh)
endif

! check on essential vars
if ( ierr/=0 ) then
stat = 99
Expand Down
10 changes: 8 additions & 2 deletions GEOSaana_GridComp/GEOSgsi_Coupler/geos_pertState.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -103,6 +103,7 @@ module geos_pertState
! and it has be verified against a GSI built-in
! adjoint-test in GSI evaljgrad().
! 28Oct2013 Todling Rename p3d to prse
! 06Mar2020 Zhu Add pblh
!
!EOP
!-------------------------------------------------------------------------
Expand Down Expand Up @@ -152,13 +153,14 @@ module geos_pertState
character(len=*),parameter:: oz_GSIpert='oz' ! = O3, in GG=gram/gram
character(len=*),parameter:: cw_GSIpert='cw' ! = QL + QI
character(len=*),parameter:: sst_GSIpert='sst' ! = ?
character(len=*),parameter:: pblh_GSIpert='pblh' ! = ZPBL

character(len=*),dimension(8),parameter:: var3dList_GSIpert = &
(/ u_GSIpert, v_GSIpert,prse_GSIpert, tv_GSIpert, &
tsen_GSIpert, q_GSIpert, oz_GSIpert, cw_GSIpert /)

character(len=*),dimension(2),parameter:: var2dList_GSIpert = &
(/ps_GSIpert, sst_GSIpert /)
character(len=*),dimension(3),parameter:: var2dList_GSIpert = &
(/ps_GSIpert, sst_GSIpert, pblh_GSIpert /)

! These are expected from a PGCM state
character(len=*),parameter:: u_GCMpert='U' ! = u
Expand All @@ -169,11 +171,15 @@ module geos_pertState
character(len=*),parameter:: ql_GCMpert='QL' ! = cw - qi
character(len=*),parameter:: qi_GCMpert='QI' ! = qi (=0)
character(len=*),parameter:: o3_GCMpert='O3' ! = oz, in PPMV
character(len=*),parameter:: pblh_GCMpert='PBLH' ! = zpbl

character(len=*),dimension(8),parameter:: var3dList_GCMpert = &
(/ u_GCMpert, v_GCMpert, tv_GCMpert, dp_GCMpert, &
qv_GCMpert, ql_GCMpert, qi_GCMpert, o3_GCMpert /)

character(len=*),dimension(1),parameter:: var2dList_GCMpert = &
(/ pblh_GCMpert /)

integer(i_kind),parameter:: ROOT =0
real(r_kind), parameter :: kPa_per_Pa = 0.001_r_kind
real(r_kind), parameter :: Pa_per_kPa = 1000._r_kind
Expand Down
Loading