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

Adding I/O interfaces in GSI for Analysis of Significant Wave Height (HOWV) and near-surface Wind Gust (GUST) for WRF-ARW based 3DRTMA #835

Open
wants to merge 6 commits into
base: develop
Choose a base branch
from
+510 −17
Open

Adding I/O interfaces in GSI for Analysis of Significant Wave Height (HOWV) and near-surface Wind Gust (GUST) for WRF-ARW based 3DRTMA #835

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
9e2a3ff
Adding code to analyze wave height and wind gust for 3DRTMA when using
GangZhao-NOAA Jan 30, 2025
796b173
if howv (or gust) is set in anavinfo, but its data is not found in
GangZhao-NOAA Feb 5, 2025
d87be10
Merge remote-tracking branch 'upstream/develop' into feature/rtma3d_w…
GangZhao-NOAA Feb 5, 2025
47b3115
A trivial fix to a typo in comment.
GangZhao-NOAA Feb 5, 2025
5ce20c4
Merge remote-tracking branch 'upstream/develop' into feature/rtma3d_w…
GangZhao-NOAA Feb 14, 2025
94bed49
Clear off the new changelog entries introduced by PR#835, and a few
GangZhao-NOAA Feb 14, 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
34 changes: 34 additions & 0 deletions src/gsi/cplr_read_wrf_mass_guess.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1394,6 +1394,7 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
aerotot_guess,init_aerotot_guess,wrf_pm2_5,aero_ratios
use rapidrefresh_cldsurf_mod, only: l_hydrometeor_bkio,l_gsd_soiltq_nudge
use rapidrefresh_cldsurf_mod, only: i_use_2mq4b,i_use_2mt4b
use rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda
use wrf_mass_guess_mod, only: soil_temp_cld,isli_cld,ges_xlon,ges_xlat,ges_tten,create_cld_grids
use gsi_bundlemod, only: GSI_BundleGetPointer
use gsi_metguess_mod, only: gsi_metguess_get,GSI_MetGuess_Bundle
Expand Down Expand Up @@ -1446,6 +1447,7 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
integer(i_kind) i_qc,i_qi,i_qr,i_qs,i_qg,i_qnr,i_qni,i_qnc,i_w,i_dbz
integer(i_kind) kqc,kqi,kqr,kqs,kqg,kqnr,kqni,kqnc,i_xlon,i_xlat,i_tt,ktt
integer(i_kind) i_th2,i_q2,i_soilt1,ksmois,ktslb
integer(i_kind) i_howv, i_gust
integer(i_kind) ier, istatus
integer(i_kind) n_actual_clouds
integer(i_kind) iv,n_gocart_var
Expand All @@ -1469,6 +1471,8 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
real(r_kind), pointer :: ges_tv_it (:,:,:)=>NULL()
real(r_kind), pointer :: ges_q_it (:,:,:)=>NULL()
real(r_kind), pointer :: ges_w_it (:,:,:)=>NULL()
real(r_kind), pointer :: ges_howv_it (:,:)=>NULL()
real(r_kind), pointer :: ges_gust_it (:,:)=>NULL()

real(r_kind), pointer :: ges_qc (:,:,:)=>NULL()
real(r_kind), pointer :: ges_qi (:,:,:)=>NULL()
Expand Down Expand Up @@ -1560,6 +1564,8 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
if(l_gsd_soilTQ_nudge) num_mass_fields=num_mass_fields+2*(nsig_soil-1)+1
if(i_use_2mt4b > 0 ) num_mass_fields=num_mass_fields + 2
if(i_use_2mq4b > 0 .and. i_use_2mt4b <=0 ) num_mass_fields=num_mass_fields + 1
if( i_howv_3dda > 0 ) num_mass_fields = num_mass_fields + 1
if( i_gust_3dda > 0 ) num_mass_fields = num_mass_fields + 1

if (laeroana_gocart .and. wrf_pm2_5 ) then
if(mype==0) write(6,*)'laeroana_gocart canoot be both true'
Expand Down Expand Up @@ -1740,6 +1746,18 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
write(identity(i),'("record ",i3,"--th2(",i2,")")')i,k
jsig_skip(i)=0 ; igtype(i)=1
endif
! for wave height (howv is after tsk/q2/soilt1/th2, and before cloud hydrometers)
if ( i_howv_3dda >0 ) then
i=i+1 ; i_howv=i ! howv
write(identity(i),'("record ",i3,"--howv")')i
jsig_skip(i)=0 ; igtype(i)=1
end if
! for wind gust (gust is after tsk/q2/soilt1/th2, and before cloud hydrometers)
if ( i_gust_3dda >0 ) then
i=i+1 ; i_gust=i ! gust
write(identity(i),'("record ",i3,"--gust")')i
jsig_skip(i)=0 ; igtype(i)=1
end if
! for cloud array
if(l_hydrometeor_bkio .and. n_actual_clouds>0) then
i_qc=i+1
Expand Down Expand Up @@ -1944,6 +1962,14 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
call GSI_BundleGetPointer ( GSI_MetGuess_Bundle(it), 't2m',ges_t2m_it, istatus );ier=ier+istatus
if (ier/=0) call die(trim(myname),'cannot get pointers for t2m,ier =',ier)
endif
if ( i_howv_3dda >0 ) then
call GSI_BundleGetPointer ( GSI_MetGuess_Bundle(it), 'howv',ges_howv_it,istatus );ier=ier+istatus
if (ier/=0) call die(trim(myname),'cannot get pointers for met-field: howv, ier =',ier)
end if
if ( i_gust_3dda >0 ) then
call GSI_BundleGetPointer ( GSI_MetGuess_Bundle(it), 'gust',ges_gust_it,istatus );ier=ier+istatus
if (ier/=0) call die(trim(myname),'cannot get pointers for met-field: gust, ier =',ier)
end if
if (l_gsd_soilTQ_nudge) then
call GSI_BundleGetPointer ( GSI_MetGuess_Bundle(it), 'tskn',ges_tsk_it, istatus );ier=ier+istatus
call GSI_BundleGetPointer ( GSI_MetGuess_Bundle(it), 'tsoil',ges_soilt1_it,istatus);ier=ier+istatus
Expand Down Expand Up @@ -2273,6 +2299,14 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
ges_q2_it(j,i)=real(all_loc(j,i,i_0+i_q2),r_kind)
ges_q2_it(j,i)=ges_q2_it(j,i)/(one+ges_q2_it(j,i))
endif
! wave height (howv)
if ( i_howv_3dda >0 ) then
ges_howv_it(j,i) = real(all_loc(j,i,i_0+i_howv),r_kind)
end if
! wind gust (gust)
if ( i_gust_3dda >0 ) then
ges_gust_it(j,i) = real(all_loc(j,i,i_0+i_gust),r_kind)
end if
! for cloud analysis
if(l_hydrometeor_bkio .and. n_actual_clouds>0) then
soil_temp_cld(j,i,it)=soil_temp(j,i,it)
Expand Down
12 changes: 11 additions & 1 deletion src/gsi/cplr_regional_io.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,7 @@ subroutine convert_regional_guess_wrf(this,mype,ctph0,stph0,tlm0)
!$$$ end documentation block

use kinds, only: i_kind,r_kind
use mpimod, only: mpi_integer4,mpi_rtype
use mpimod, only: mpi_integer4,mpi_rtype,mpi_itype
use gridmod, only: wrf_mass_regional,wrf_nmm_regional,&
nems_nmmb_regional,cmaq_regional,&
twodvar_regional,netcdf
Expand All @@ -88,6 +88,8 @@ subroutine convert_regional_guess_wrf(this,mype,ctph0,stph0,tlm0)
use mpimod, only: mpi_comm_world,ierror
use wrf_params_mod, only: update_pint,cold_start
use gsi_io, only: verbose
use rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda
use rapidrefresh_cldsurf_mod, only: i_howv_mask

implicit none

Expand Down Expand Up @@ -145,7 +147,15 @@ subroutine convert_regional_guess_wrf(this,mype,ctph0,stph0,tlm0)
end if
call mpi_barrier(mpi_comm_world,ierror)
call mpi_bcast(byte_swap,1,mpi_integer4,0,mpi_comm_world,ierror)
call mpi_bcast(i_howv_3dda, 1, mpi_itype, 0, mpi_comm_world, ierror)
call mpi_bcast(i_gust_3dda, 1, mpi_itype, 0, mpi_comm_world, ierror)
call mpi_bcast(i_howv_mask, 1, mpi_itype, 0, mpi_comm_world, ierror)
if(print_verbose)write(6,*)' in convert_regional_guess, for wrf arw binary input, byte_swap=',byte_swap
if (mype <= 1 .and. print_verbose) then
write(6,'(1x,A,3(2x,I4),2x,A6,I6.6,A2)') &
' in convert_regional_guess, i_howv_3dda i_gust_3dda i_howv_mask =', &
i_howv_3dda,i_gust_3dda,i_howv_mask,' (pe=',mype,').'
end if

elseif (cmaq_regional) then
if (mype==0) then
Expand Down
175 changes: 175 additions & 0 deletions src/gsi/cplr_wrf_netcdf_interface.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -66,11 +66,14 @@ subroutine convert_netcdf_mass_wrf(this)
use gsi_4dvar, only: nhr_assimilation
use rapidrefresh_cldsurf_mod, only: l_hydrometeor_bkio,l_gsd_soilTQ_nudge
use rapidrefresh_cldsurf_mod, only: i_use_2mt4b,i_use_2mq4b
use rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda
use rapidrefresh_cldsurf_mod, only: i_howv_mask
use gsi_metguess_mod, only: gsi_metguess_get
use chemmod, only: laeroana_gocart, ppmv_conv,wrf_pm2_5
use gsi_chemguess_mod, only: gsi_chemguess_get
use gridmod, only: wrf_mass_hybridcord
use netcdf_mod, only: nc_check
use mpeu_util, only: die, perr

use wrf_vars_mod, only : w_exist, dbz_exist
use constants, only: zero
Expand All @@ -91,6 +94,12 @@ subroutine convert_netcdf_mass_wrf(this)
integer(i_kind) :: dh1

integer(i_kind) :: iunit

character(len=120) :: filename_lakemask ! lakemask.bin
integer(i_kind) :: iunit_lakemask
logical :: l_iunit_opened
character(len=200) :: derr_msg ! error message used in die
character(len=24),parameter :: myname_ = 'convert_netcdf_mass_wrf'

integer(i_kind) :: i,j,k
integer(i_kind) :: ndim1
Expand Down Expand Up @@ -1022,6 +1031,116 @@ subroutine convert_netcdf_mass_wrf(this)
write(iunit)field2 !TH2
endif

! Reading Significant Wave Height (HOWV) in firstguess
if(i_howv_3dda >= 1) then ! if howv is found in anavinfo
rmse_var='HOWV'
call ext_ncd_get_var_info (dh1,trim(rmse_var),ndim1,ordering,staggering, &
start_index,end_index, WrfType, ierr )
if(print_verbose)then
write(6,*)' rmse_var=',trim(rmse_var)
write(6,*)' ordering=',ordering
write(6,*)' WrfType,ierr=',WrfType,ierr
write(6,*)' ndim1=',ndim1
write(6,*)' staggering=',staggering
write(6,*)' start_index=',start_index
write(6,*)' end_index=',end_index
end if
if(ierr == 0) then
call ext_ncd_read_field(dh1,DateStr1,TRIM(rmse_var), &
field2,WRF_REAL,0,0,0,ordering, &
staggering, dimnames , &
start_index,end_index, & !dom
start_index,end_index, & !mem
start_index,end_index, & !pat
ierr )
if(print_verbose)write(6,*)' max,min ',trim(adjustl(rmse_var)),'=',maxval(field2),minval(field2)
write(iunit)field2 !HOWV (2D wave height)
i_howv_3dda = 2 ! howv was found both in anavinfo and firstguess
else
i_howv_3dda = 0 ! skipping analysis of this variable
derr_msg='Warning: Error when get info in firstguess for var : '//trim(rmse_var)// &
' ==> Stop GSI analysis ... ierr='
! write(6,'(1x,A,A2,A,1x,I6)')trim(adjustl(myname_)),'::',trim(adjustl(derr_msg)),ierr
! call stop2(ierr)
call die(trim(adjustl(myname_)),trim(adjustl(derr_msg)),ierr)
end if

! reading lakemask for screening off wave height over lake (and/or land) area
if(i_howv_mask .eq. 2) then
rmse_var='LAKEMASK' ! see the code for rmse_var='XLAND'
call ext_ncd_get_var_info (dh1,trim(rmse_var),ndim1,ordering,staggering, &
start_index,end_index, WrfType, ierr )
if(ierr .eq. 0) then
if(print_verbose)then
write(6,*)' rmse_var = ',trim(rmse_var),' ndim1=',ndim1
write(6,*)' WrfType = ',WrfType,' WRF_REAL=',WRF_REAL,'ierr = ',ierr !DEDE
write(6,*)' ordering = ',trim(ordering),' staggering = ',trim(staggering)
write(6,*)' start_index = ',start_index,' end_index = ',end_index
end if
call ext_ncd_read_field(dh1,DateStr1,TRIM(rmse_var), &
field2,WRF_REAL,0,0,0,ordering, &
staggering, dimnames , &
start_index,end_index, & !dom
start_index,end_index, & !mem
start_index,end_index, & !pat
ierr )
if(print_verbose)then
write(6,*)' max,min lakemask=',maxval(field2),minval(field2)
write(6,*)' lakemask(1,1),lakemask(nlon,1)=',field2(1,1),field2(nlon_regional,1)
write(6,*)' lakemask(1,nlat),lakemask(nlon,nlat)=', &
field2(1,nlat_regional),field2(nlon_regional,nlat_regional)
end if
! dumping out lakemask into an temporary binary file
! (which would be used for howv mask when updating the analysis of howv)
filename_lakemask="lakemask.bin"; iunit_lakemask=115; l_iunit_opened=.true.
inquire(unit=iunit_lakemask, opened=l_iunit_opened)
if(l_iunit_opened) iunit_lakemask=iunit_lakemask+5000
if(print_verbose) write(6,*)' I/O unit for writing lakemask: ', iunit_lakemask
open(iunit_lakemask,file=trim(adjustl(filename_lakemask)),form='unformatted')
write(iunit_lakemask)field2 !LAKEMASK (1=lake, 0=water)
close(iunit_lakemask)
else
write(6,'(1x,A)')' No LAKEMASK data in firstguess data. No lakemask is used for HOWV mask. reset i_howv_mask=1'
i_howv_mask=1 ! i_howv_mask is changed on pe 0, so it must be broadcast to all PEs.
endif ! if error in reading lakemask
endif ! i_howv_mask = 2
endif ! i_howv_3dda (reading 2D wave height from netcdf-format background)

! Reading 10-m Wind Gust (GUST) in firstguess
if(i_gust_3dda >= 1) then ! if gust is found in anavinfo
rmse_var='GUST'
call ext_ncd_get_var_info (dh1,trim(rmse_var),ndim1,ordering,staggering, &
start_index,end_index, WrfType, ierr )
if(print_verbose)then
write(6,*)' rmse_var=',trim(rmse_var)
write(6,*)' ordering=',ordering
write(6,*)' WrfType,ierr=',WrfType,ierr
write(6,*)' ndim1=',ndim1
write(6,*)' staggering=',staggering
write(6,*)' start_index=',start_index
write(6,*)' end_index=',end_index
end if
if(ierr == 0) then
call ext_ncd_read_field(dh1,DateStr1,TRIM(rmse_var), &
field2,WRF_REAL,0,0,0,ordering, &
staggering, dimnames , &
start_index,end_index, & !dom
start_index,end_index, & !mem
start_index,end_index, & !pat
ierr )
if(print_verbose)write(6,*)'convert_netcdf_mass_wrf:: max,min ',trim(adjustl(rmse_var)),'=',maxval(field2),minval(field2)
write(iunit)field2 !GUST (2D Wind Gust)
i_gust_3dda = 2 ! gust was found both in anavinfo and firstguess
else
i_gust_3dda = 0 ! skipping analysis of this variable
derr_msg='Warning: Error when get info in firstguess for var : '//trim(rmse_var)// &
' ==> Stop GSI analysis ... ierr='
! write(6,'(1x,A,A2,A,1x,I6)')trim(adjustl(myname_)),'::',trim(adjustl(derr_msg)),ierr
! call stop2(ierr)
call die(trim(adjustl(myname_)),trim(adjustl(derr_msg)),ierr)
end if
endif ! i_gust_3dda (reading 2D 10-m wind gust from netcdf-format background)

if(l_hydrometeor_bkio .and. n_actual_clouds>0) then
rmse_var='QCLOUD'
call ext_ncd_get_var_info (dh1,trim(rmse_var),ndim1,ordering,staggering, &
Expand Down Expand Up @@ -2437,6 +2556,7 @@ subroutine update_netcdf_mass_wrf(this)
use constants, only: h300,tiny_single
use rapidrefresh_cldsurf_mod, only: l_hydrometeor_bkio,l_gsd_soilTQ_nudge
use rapidrefresh_cldsurf_mod, only: i_gsdcldanal_type
use rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda
use gsi_metguess_mod, only: gsi_metguess_get,GSI_MetGuess_Bundle
use rapidrefresh_cldsurf_mod, only: i_use_2mt4b,i_use_2mq4b
use gsi_bundlemod, only: GSI_BundleGetPointer
Expand Down Expand Up @@ -3002,6 +3122,61 @@ subroutine update_netcdf_mass_wrf(this)
ierr )
endif

! writing analysis of significantg wave height (howv) into firstguess (analysis) file
if(i_howv_3dda==2) then
read(iunit) field2 !HOWV (wave height)
rmse_var='HOWV'
if(print_verbose)write(6,'(1x,4A,2(1x,F15.7))') &
trim(adjustl(myname_)),'::max,min ',trim(rmse_var),'=', &
maxval(field2),minval(field2)
call ext_ncd_get_var_info (dh1,trim(rmse_var),ndim1,ordering,staggering, &
start_index,end_index1, WrfType, ierr )
if(print_verbose)then
write(6,*)' rmse_var=',trim(rmse_var)
write(6,*)' ordering=',ordering
write(6,*)' WrfType,WRF_REAL=',WrfType,WRF_REAL
write(6,*)' ndim1=',ndim1
write(6,*)' staggering=',staggering
write(6,*)' start_index=',start_index
write(6,*)' end_index1=',end_index1
end if
call ext_ncd_write_field(dh1,DateStr1,TRIM(rmse_var), &
field2,WRF_REAL,0,0,0,ordering, &
staggering, dimnames , &
start_index,end_index1, & !dom
start_index,end_index1, & !mem
start_index,end_index1, & !pat
ierr )
endif ! i_howv_3dda (wave height)

! writing analysis of 10-m wind gust (gust) into firstguess (analysis) file
if(i_gust_3dda==2) then
read(iunit) field2 !GUST (wind gust)
rmse_var='GUST'
if(print_verbose)write(6,'(1x,4A,2(1x,F15.7))') &
trim(adjustl(myname_)),'::max,min ',trim(rmse_var),'=', &
maxval(field2),minval(field2)
call ext_ncd_get_var_info (dh1,trim(rmse_var),ndim1,ordering,staggering, &
start_index,end_index1, WrfType, ierr )
if(print_verbose)then
write(6,*)' rmse_var=',trim(rmse_var)
write(6,*)' ordering=',ordering
write(6,*)' WrfType,WRF_REAL=',WrfType,WRF_REAL
write(6,*)' ndim1=',ndim1
write(6,*)' staggering=',staggering
write(6,*)' start_index=',start_index
write(6,*)' end_index1=',end_index1
end if
call ext_ncd_write_field(dh1,DateStr1,TRIM(rmse_var), &
field2,WRF_REAL,0,0,0,ordering, &
staggering, dimnames , &
start_index,end_index1, & !dom
start_index,end_index1, & !mem
start_index,end_index1, & !pat
ierr )
endif ! i_howv_3dda (wave height)


if (l_hydrometeor_bkio .and. n_actual_clouds>0) then
do k=1,nsig_regional
read(iunit)((field3(i,j,k),i=1,nlon_regional),j=1,nlat_regional) ! Qc
Expand Down
Loading