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.

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

Adding I/O for direct analysis of near-surface wind gust for RRFS-based 3DRTMA #730

Merged
Show file tree
Hide file tree
Changes from 8 commits
Commits
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
65 changes: 54 additions & 11 deletions src/gsi/gsi_rfv3io_mod.f90
Original file line number Diff line number Diff line change
Expand Up @@ -61,7 +61,7 @@ module gsi_rfv3io_mod
use rapidrefresh_cldsurf_mod, only: i_use_2mq4b,i_use_2mt4b
use chemmod, only: naero_cmaq_fv3,aeronames_cmaq_fv3,imodes_cmaq_fv3,laeroana_fv3cmaq
use chemmod, only: naero_smoke_fv3,aeronames_smoke_fv3,laeroana_fv3smoke
use rapidrefresh_cldsurf_mod, only: i_howv_3dda
use rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda

implicit none
public type_fv3regfilenameg
Expand Down Expand Up @@ -147,7 +147,7 @@ module gsi_rfv3io_mod
public :: mype_u,mype_v,mype_t,mype_q,mype_p,mype_oz,mype_ql
public :: mype_qi,mype_qr,mype_qs,mype_qg,mype_qnr,mype_w
public :: k_slmsk,k_tsea,k_vfrac,k_vtype,k_stype,k_zorl,k_smc,k_stc
public :: k_snwdph,k_f10m,mype_2d,n2d,k_orog,k_psfc,k_t2m,k_q2m,k_howv
public :: k_snwdph,k_f10m,mype_2d,n2d,k_orog,k_psfc,k_t2m,k_q2m,k_howv,k_gust
public :: ijns,ijns2d,displss,displss2d,ijnz,displsz_g
public :: fv3lam_io_dynmetvars3d_nouv,fv3lam_io_tracermetvars3d_nouv
public :: fv3lam_io_tracerchemvars3d_nouv,fv3lam_io_tracersmokevars3d_nouv
Expand All @@ -158,7 +158,7 @@ module gsi_rfv3io_mod
integer(i_kind) mype_qi,mype_qr,mype_qs,mype_qg,mype_qnr,mype_w

integer(i_kind) k_slmsk,k_tsea,k_vfrac,k_vtype,k_stype,k_zorl,k_smc,k_stc
integer(i_kind) k_snwdph,k_f10m,mype_2d,n2d,k_orog,k_psfc,k_t2m,k_q2m,k_howv
integer(i_kind) k_snwdph,k_f10m,mype_2d,n2d,k_orog,k_psfc,k_t2m,k_q2m,k_howv,k_gust
parameter( &
k_f10m =1, & !fact10
k_stype=2, & !soil_type
Expand All @@ -174,7 +174,8 @@ module gsi_rfv3io_mod
k_q2m =12, & ! 2 m Q
k_orog =13, & !terrain
k_howv =14, & ! significant wave height (aka howv in GSI)
n2d=14 )
k_gust =15, & ! wind gust (aka gust in GSI)
n2d=15 )
logical :: grid_reverse_flag
character(len=max_varname_length),allocatable,dimension(:) :: fv3lam_io_dynmetvars3d_nouv
! copy of cvars3d excluding uv 3-d fields
Expand Down Expand Up @@ -996,6 +997,7 @@ subroutine read_fv3_netcdf_guess(fv3filenamegin)
real(r_kind),dimension(:,:),pointer::ges_t2m=>NULL()
real(r_kind),dimension(:,:),pointer::ges_q2m=>NULL()
real(r_kind),dimension(:,:),pointer::ges_howv=>NULL()
real(r_kind),dimension(:,:),pointer::ges_gust=>NULL()

real(r_kind),dimension(:,:,:),pointer::ges_ql=>NULL()
real(r_kind),dimension(:,:,:),pointer::ges_qi=>NULL()
Expand Down Expand Up @@ -1274,6 +1276,7 @@ subroutine read_fv3_netcdf_guess(fv3filenamegin)
else if(trim(vartem)=='t2m') then
else if(trim(vartem)=='q2m') then
else if(trim(vartem)=='howv') then
else if(trim(vartem)=='gust') then
else
write(6,*)'the metvarname2 ',trim(vartem),' has not been considered yet, stop'
call stop2(333)
Expand All @@ -1294,7 +1297,7 @@ subroutine read_fv3_netcdf_guess(fv3filenamegin)
vartem=trim(name_metvars2d(i))
if(.not.( (trim(vartem)=='ps'.and.fv3sar_bg_opt==0).or.(trim(vartem)=="z") &
.or.(trim(vartem)=="t2m").or.(trim(vartem)=="q2m") &
.or.(trim(vartem)=="howv"))) then ! z is treated separately
.or.(trim(vartem)=="howv").or.(trim(vartem)=="gust"))) then ! z is treated separately
if (ifindstrloc(vardynvars,trim(vartem)) > 0) then
jdynvar=jdynvar+1
fv3lam_io_dynmetvars2d_nouv(jdynvar)=trim(vartem)
Expand Down Expand Up @@ -1557,6 +1560,12 @@ subroutine read_fv3_netcdf_guess(fv3filenamegin)
if (ier/=0) call die(trim(myname),'cannot get pointers for howv, ier=',ier)
endif

!--- wind gust (gust)
if ( i_gust_3dda == 1 ) then
call GSI_BundleGetPointer(GSI_MetGuess_Bundle(it),'gust',ges_gust,istatus ); ier=ier+istatus
if (ier/=0) call die(trim(myname),'cannot get pointers for gust, ier=',ier)
endif

if(mype == 0 ) then
call check(nf90_open(fv3filenamegin(it)%dynvars,nf90_nowrite,loc_id))
call check(nf90_inquire(loc_id,formatNum=ncfmt))
Expand Down Expand Up @@ -1739,7 +1748,8 @@ subroutine read_fv3_netcdf_guess(fv3filenamegin)
endif


call gsi_fv3ncdf2d_read(fv3filenamegin(it),it,ges_z,ges_t2m,ges_q2m,ges_howv)
call gsi_fv3ncdf2d_read(fv3filenamegin(it),it,ges_z,ges_t2m,ges_q2m, &
ges_howv,ges_gust)

if(i_use_2mq4b > 0 .and. i_use_2mt4b > 0 ) then
! Convert 2m guess mixing ratio to specific humidity
Expand Down Expand Up @@ -1975,7 +1985,8 @@ end subroutine gsi_bundlegetpointer_fv3lam_tracerchem_nouv

end subroutine read_fv3_netcdf_guess

subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m, &
ges_howv,ges_gust)
!$$$ subprogram documentation block
! . . . .
! subprogram: gsi_fv3ncdf2d_read
Expand Down Expand Up @@ -2022,6 +2033,7 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
real(r_kind), intent(in),dimension(:,:),pointer::ges_t2m
real(r_kind), intent(in),dimension(:,:),pointer::ges_q2m
real(r_kind), intent(in),dimension(:,:),pointer::ges_howv
real(r_kind), intent(in),dimension(:,:),pointer::ges_gust
type (type_fv3regfilenameg),intent(in) :: fv3filenamegin

character(len=max_varname_length) :: name
Expand All @@ -2036,8 +2048,8 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
integer(i_kind) kk,n,ns,j,ii,jj,mm1
character(len=:),allocatable :: sfcdata !='fv3_sfcdata'
character(len=:),allocatable :: dynvars !='fv3_dynvars'
! for checking the existence of howv in firstguess file
integer(i_kind) id_howv
! for checking the existence of howv/gust in firstguess file
integer(i_kind) id_howv, id_gust
integer(i_kind) iret_bcast

! for io_layout > 1
Expand All @@ -2057,8 +2069,9 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
allocate(work(itotsub*n2d))
allocate( sfcn2d(lat2,lon2,n2d))

!-- initialisation of the array for howv
!-- initialisation of the array for howv/gust
sfcn2d(:,:,k_howv) = zero
sfcn2d(:,:,k_gust) = zero

!-- initialisation of the array for sfc_var_exist
sfc_var_exist = .false.
Expand Down Expand Up @@ -2107,6 +2120,21 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
trim(sfcdata), ', iret, varid = ',iret, id_howv,' (on pe: ', mype,').'
end if
end if
!--- check the existence of wind gust (gust) in 2D FV3-LAM firstguess file
! (similar as done above for howv)
if ( i_gust_3dda == 1 ) then
iret = nf90_inq_varid(gfile_loc,'gust',id_gust)
if ( iret /= nf90_noerr ) then
iret = nf90_inq_varid(gfile_loc,'GUST',id_gust) ! double check with name in uppercase
end if
if ( iret /= nf90_noerr ) then
i_gust_3dda = 0 ! gust does not exist in firstguess, then stop GSI run.
call die('gsi_fv3ncdf2d_read','Warning: CANNOT find gust in firstguess, aborting..., iret = ', iret)
else
write(6,'(1x,A,1x,A,1x,A,1x,I4,1x,I4,1x,A,1x,I4.4,A)') 'gsi_fv3ncdf2d_read:: Found gust in firstguess ', &
trim(sfcdata), ', iret, varid = ',iret, id_gust,' (on pe: ', mype,').'
end if
end if

!!!!!!!!!!!! read in 2d variables !!!!!!!!!!!!!!!!!!!!!!!!!!
do i=ndimensions+1,nvariables
Expand Down Expand Up @@ -2150,6 +2178,9 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
else if( trim(name)=='HOWV'.or.trim(name)=='howv' ) then
k=k_howv
sfc_var_exist(k) = .true.
else if( trim(name)=='GUST'.or.trim(name)=='gust' ) then
k=k_gust
sfc_var_exist(k) = .true.
else
cycle
endif
Expand Down Expand Up @@ -2283,8 +2314,9 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
if(allocated(sfc_fulldomain)) deallocate (sfc_fulldomain)
endif ! mype

!-- broadcast the updated i_howv_3dda to all tasks (!!!!)
!-- broadcast the updated i_howv_3dda, i_gust_3dda to all tasks (!!!!)
call mpi_bcast(i_howv_3dda, 1, mpi_itype, mype_2d, mpi_comm_world, iret_bcast)
call mpi_bcast(i_gust_3dda, 1, mpi_itype, mype_2d, mpi_comm_world, iret_bcast)

!-- broadcast the updated sfc_var_exist to all tasks (!!!!)
call mpi_bcast(sfc_var_exist, n2d, mpi_itype, mype_2d, mpi_comm_world, iret_bcast)
Expand Down Expand Up @@ -2313,6 +2345,9 @@ subroutine gsi_fv3ncdf2d_read(fv3filenamegin,it,ges_z,ges_t2m,ges_q2m,ges_howv)
if ( i_howv_3dda == 1 ) then
if ( sfc_var_exist(k_howv) ) ges_howv(:,:)=sfcn2d(:,:,k_howv)
endif
if ( i_gust_3dda == 1 ) then
if ( sfc_var_exist(k_gust) ) ges_gust(:,:)=sfcn2d(:,:,k_gust)
endif
deallocate (sfcn2d,a)
return
end subroutine gsi_fv3ncdf2d_read
Expand Down Expand Up @@ -3628,6 +3663,7 @@ subroutine wrfv3_netcdf(fv3filenamegin)
real(r_kind),pointer,dimension(:,: ):: ges_t2m =>NULL()
real(r_kind),pointer,dimension(:,: ):: ges_q2m =>NULL()
real(r_kind),pointer,dimension(:,: ):: ges_howv =>NULL()
real(r_kind),pointer,dimension(:,: ):: ges_gust =>NULL()

integer(i_kind) i,k

Expand Down Expand Up @@ -3750,6 +3786,9 @@ subroutine wrfv3_netcdf(fv3filenamegin)
if ( i_howv_3dda == 1 ) then
call GSI_BundleGetPointer (GSI_MetGuess_Bundle(it),'howv',ges_howv,istatus); ier=ier+istatus
endif
if ( i_gust_3dda == 1 ) then
call GSI_BundleGetPointer (GSI_MetGuess_Bundle(it),'gust',ges_gust,istatus); ier=ier+istatus
endif
if (ier/=0) call die('wrfv3_netcdf','cannot get pointers for fv3 met-fields, ier =',ier)

if (laeroana_fv3cmaq) then
Expand Down Expand Up @@ -3964,6 +4003,10 @@ subroutine wrfv3_netcdf(fv3filenamegin)
if ( i_howv_3dda == 1 ) then
call gsi_fv3ncdf_write_sfc(fv3filenamegin,'howv',ges_howv,add_saved)
endif
!-- output analysis of gust
if ( i_gust_3dda == 1 ) then
call gsi_fv3ncdf_write_sfc(fv3filenamegin,'gust',ges_gust,add_saved)
endif

if(allocated(g_prsi)) deallocate(g_prsi)

Expand Down
7 changes: 5 additions & 2 deletions src/gsi/gsimod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -184,7 +184,7 @@ module gsimod
cld_bld_coverage,cld_clr_coverage,&
i_cloud_q_innovation,i_ens_mean,DTsTmax,&
i_T_Q_adjust,l_saturate_bkCloud,l_rtma3d,i_precip_vertical_check, &
corp_howv, hwllp_howv
corp_howv, hwllp_howv, corp_gust, hwllp_gust, oerr_gust
use gsi_metguess_mod, only: gsi_metguess_init,gsi_metguess_final
use gsi_chemguess_mod, only: gsi_chemguess_init,gsi_chemguess_final
use tcv_mod, only: init_tcps_errvals,tcp_refps,tcp_width,tcp_ermin,tcp_ermax
Expand Down Expand Up @@ -1604,6 +1604,9 @@ module gsimod
! = 0.42 meters (default)
! hwllp_howv - real, background error de-correlation length scale of howv
! = 170,000.0 meters (default 170 km)
! corp_gust - real, static background error of gust (stddev error)
! hwllp_gust - real, background error de-correlation length scale of gust
! oerr_gust - real, observation error of gust
!
namelist/rapidrefresh_cldsurf/dfi_radar_latent_heat_time_period, &
metar_impact_radius,metar_impact_radius_lowcloud, &
Expand All @@ -1625,7 +1628,7 @@ module gsimod
cld_bld_coverage,cld_clr_coverage,&
i_cloud_q_innovation,i_ens_mean,DTsTmax, &
i_T_Q_adjust,l_saturate_bkCloud,l_rtma3d,i_precip_vertical_check, &
corp_howv, hwllp_howv
corp_howv, hwllp_howv, corp_gust, hwllp_gust, oerr_gust

! chem(options for gsi chem analysis) :
! berror_chem - .true. when background for chemical species that require
Expand Down
5 changes: 5 additions & 0 deletions src/gsi/m_berror_stats_reg.f90
Original file line number Diff line number Diff line change
Expand Up @@ -313,6 +313,7 @@ subroutine berror_read_wgt_reg(msig,mlat,corz,corp,hwll,hwllp,vz,rlsig,varq,qopt
use mpeu_util,only: getindex
use radiance_mod, only: icloud_cv,n_clouds_fwd,cloud_names_fwd
use chemmod, only: berror_fv3_cmaq_regional,berror_fv3_sd_regional
use rapidrefresh_cldsurf_mod, only: corp_gust, hwllp_gust, l_rtma3d

implicit none

Expand Down Expand Up @@ -830,6 +831,10 @@ subroutine berror_read_wgt_reg(msig,mlat,corz,corp,hwll,hwllp,vz,rlsig,varq,qopt
do i=0,mlat+1
hwllp(i,n)=hwll(i,1,nrf3_q)
end do
if ( l_rtma3d ) then ! tuning stddev and length for analysis of gust in 3drtma test run
ShunLiu-NOAA marked this conversation as resolved.
Show resolved Hide resolved
if ( corp_gust .gt. 0.0_r_kind ) corp(1:mlat, n) = corp_gust
if ( hwllp_gust .gt. 0.0_r_kind ) hwllp(0:mlat+1,n) = hwllp_gust
end if
else if (n==nrf2_vis) then
do i=1,mlat
corp(i,n)=3.0_r_kind
Expand Down
21 changes: 21 additions & 0 deletions src/gsi/rapidrefresh_cldsurf_mod.f90
Original file line number Diff line number Diff line change
Expand Up @@ -197,6 +197,13 @@ module rapidrefresh_cldsurf_mod
! just the reduced static BE of howv. If to make the analysis of howv
! in hyrbid run is as similar as the analysis of howv in pure 3dvar run,
! the static BE of howv used in hybrid run needs to be tuned (inflated actually).
! corp_gust - namelist real, static BE of gust (standard error deviation)
! hwllp_gust - namelist real, static BE de-correlation length scale of gust
! oerr_gust - namelist real, observation error of gust
! i_gust_3dda - integer, control the analysis of gust in 3D analysis (either var or hybrid)
! = 0 (gust-off: default) : no analysis of gust in 3D analysis.
! = 1 (gust-on) : if variable name "gust" is found in anavinfo,
! set it to be 1 to turn on analysis of gust;
!
! attributes:
! language: f90
Expand Down Expand Up @@ -270,6 +277,8 @@ module rapidrefresh_cldsurf_mod
public :: i_precip_vertical_check
public :: corp_howv, hwllp_howv
public :: i_howv_3dda
public :: corp_gust, hwllp_gust, oerr_gust
public :: i_gust_3dda

logical l_hydrometeor_bkio
real(r_kind) dfi_radar_latent_heat_time_period
Expand Down Expand Up @@ -330,6 +339,8 @@ module rapidrefresh_cldsurf_mod
integer(i_kind) i_precip_vertical_check
real(r_kind) :: corp_howv, hwllp_howv
integer(i_kind) :: i_howv_3dda
real(r_kind) :: corp_gust, hwllp_gust, oerr_gust
integer(i_kind) :: i_gust_3dda

contains

Expand Down Expand Up @@ -447,6 +458,10 @@ subroutine init_rapidrefresh_cldsurf
corp_howv = 0.42_r_kind ! 0.42 meters (default)
hwllp_howv = 170000.0_r_kind ! 170,000.0 meters (170km as default for 3DRTMA, 50km is used in 2DRTMA)
i_howv_3dda = 0 ! no analysis of significant wave height (howv) in 3D analysis (default)
corp_gust = -3.00_r_kind ! default to use the value (3.0 m/s) set in m_berror_stats_reg.f90
GangZhao-NOAA marked this conversation as resolved.
Show resolved Hide resolved
hwllp_gust = -170000.0_r_kind ! default to use the value (same as q) set in m_berror_stats_reg.f90
oerr_gust = -1.0_r_kind ! default to use the value (1.0 m/s) set in read_prepbufr.f90
i_gust_3dda = 0 ! no analysis of wind gust (gust) in 3D analysis (default)

!-- searching for specific variable in state variable list (reading from anavinfo)
do i2=1,ns2d
Expand All @@ -456,6 +471,12 @@ subroutine init_rapidrefresh_cldsurf
write(6,'(1x,A,1x,A8,1x,A,1x,I4)')"init_rapidrefresh_cldsurf: anavinfo svars2d (state variable): ",trim(adjustl(svars2d(i2))), " is found in anavinfo, set i_howv_3dda = ", i_howv_3dda
end if
end if
if ( trim(svars2d(i2))=='gust' .or. trim(svars2d(i2))=='GUST' ) then
i_gust_3dda = 1
if ( mype == 0 ) then
write(6,'(1x,A,1x,A8,1x,A,1x,I4)')"init_rapidrefresh_cldsurf: anavinfo svars2d (state variable): ",trim(adjustl(svars2d(i2))), " is found in anavinfo, set i_gust_3dda = ", i_gust_3dda
end if
end if
end do ! i2 : looping over 2-D anasv

return
Expand Down
5 changes: 4 additions & 1 deletion src/gsi/read_prepbufr.f90
Original file line number Diff line number Diff line change
Expand Up @@ -227,7 +227,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
use adjust_cloudobs_mod, only: adjust_convcldobs,adjust_goescldobs
use mpimod, only: npe
use rapidrefresh_cldsurf_mod, only: i_gsdsfc_uselist,i_gsdqc,i_ens_mean
use rapidrefresh_cldsurf_mod, only: l_rtma3d
use rapidrefresh_cldsurf_mod, only: l_rtma3d, oerr_gust
use gsi_io, only: verbose
use phil2, only: denest ! hilbert curve

Expand Down Expand Up @@ -1825,6 +1825,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
gustqm=0
if (kx==188 .or. kx==288 .or. kx==195 .or. kx==295 ) &
call get_gustqm(kx,c_station_id,c_prvstg,c_sprvstg,gustqm)
if ( l_rtma3d ) gustqm = 0 ! skipping get_gustqm for 3drtma run (missing list file)
qm=gustqm
else if(visob) then
visqm=0 ! need to fix this later
Expand Down Expand Up @@ -2556,6 +2557,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
! need to find out gustoe
! gustoe=1.8
gustoe=1.0
if ( l_rtma3d .and. oerr_gust > 0.0_r_kind ) gustoe = oerr_gust
GangZhao-NOAA marked this conversation as resolved.
Show resolved Hide resolved
selev=stnelev
oelev=obsdat(4,k)
if(selev == oelev)oelev=r10+selev
Expand All @@ -2577,6 +2579,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
if ((kx==188).or.(kx==288) .or.(kx==195) .or.(kx==295)) then
! gustoe=2.5
gustoe=1.0
if ( l_rtma3d .and. oerr_gust > 0.0_r_kind ) gustoe = oerr_gust
GangZhao-NOAA marked this conversation as resolved.
Show resolved Hide resolved
windcorr=abs(obsdat(5,k))<1.0 .and. abs(obsdat(6,k))<1.0 .and. obsdat(8,k)>10.0
if (windcorr) gustoe=gustoe*1.5_r_kind

Expand Down
Loading