Height Adjustment of Background 10-m Wind Gust When Analyzing Near-surface Wind Gust observations in GSI for HRRR-3DRTMA

src/gsi/cplr_read_wrf_mass_guess.f90

@@ -1394,7 +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 rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda, i_sfcrough_fgs
     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
@@ -1447,7 +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) i_howv, i_gust, i_rough
     integer(i_kind) ier, istatus
     integer(i_kind) n_actual_clouds
     integer(i_kind) iv,n_gocart_var
@@ -1566,6 +1566,7 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
        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( i_sfcrough_fgs > 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'
@@ -1758,6 +1759,12 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
          write(identity(i),'("record ",i3,"--gust")')i
          jsig_skip(i)=0 ; igtype(i)=1
        end if
+  !  for surface roughness (sfcrough is after tsk/q2/soilt1/th2, and before cloud hydrometers) 
+       if ( i_sfcrough_fgs >0 ) then
+         i=i+1 ; i_rough=i                                               ! roughtness
+         write(identity(i),'("record ",i3,"--rough")')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
@@ -2283,7 +2290,12 @@ subroutine read_wrf_mass_netcdf_guess_wrf(this,mype)
                 psfc_this=(psfc_this_dry-pt_ll)*q_integral(j,i)+pt_ll+q_integralc4h(j,i)
                 ges_ps_it(j,i)=one_tenth*psfc_this   ! convert from mb to cb
                 sno(j,i,it)=real(all_loc(j,i,i_0+i_sno),r_kind)
-                sfc_rough(j,i,it)=rough_default
+  ! surface roughness (ZNT)
+                if ( i_sfcrough_fgs >0 ) then
+                    sfc_rough(j,i,it) = real(all_loc(j,i,i_0+i_rough),r_kind)
+                else
+                    sfc_rough(j,i,it) = rough_default
+                end if
                 if(i_use_2mt4b > 0 ) then
                    ges_t2m_it(j,i)=real(all_loc(j,i,i_0+i_th2),r_kind)
   ! convert from potential to sensible temperature

src/gsi/cplr_regional_io.f90

@@ -89,7 +89,7 @@ subroutine convert_regional_guess_wrf(this,mype,ctph0,stph0,tlm0)
       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
+      use rapidrefresh_cldsurf_mod, only: i_howv_mask, i_sfcrough_fgs
 
       implicit none
 
@@ -150,11 +150,12 @@ subroutine convert_regional_guess_wrf(this,mype,ctph0,stph0,tlm0)
          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)
+         call mpi_bcast(i_sfcrough_fgs, 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,').'
+            write(6,'(1x,A,4(2x,I4),2x,A6,I6.6,A2)')                                   &
+              ' in convert_regional_guess, i_howv_3dda i_gust_3dda i_howv_mask i_sfcrough_fgs=',     &
+              i_howv_3dda,i_gust_3dda,i_howv_mask,i_sfcrough_fgs,'  (pe=',mype,').'
          end if
 
       elseif (cmaq_regional) then

src/gsi/cplr_wrf_netcdf_interface.f90

@@ -67,7 +67,7 @@ subroutine convert_netcdf_mass_wrf(this)
     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 rapidrefresh_cldsurf_mod, only: i_howv_mask, i_sfcrough_fgs
     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
@@ -1141,6 +1141,35 @@ subroutine convert_netcdf_mass_wrf(this)
           end if
        endif  ! i_gust_3dda (reading 2D 10-m wind gust from netcdf-format background) 
 
+!      Reading Surface Roughness (ZNT) in firstguess (it is dumped in sigf file after gust and before cloud variables.)
+       if ( i_sfcrough_fgs >= 1 ) then
+          rmse_var='ZNT'
+          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),' 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
+          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 ZNT(sfc_rough)=',maxval(field2),minval(field2)
+             write(iunit)field2   !ZNT  (Surface Roughness Length)
+             i_sfcrough_fgs = 1
+          else
+             i_sfcrough_fgs = 0
+             write(6,'(1x,A,1x,I4)') &
+               'convert_netcdf_mass_wrf:: ZNT (sfc_roughness) is NOT in firstguess. Re-set i_sfcrough_fgs=',i_sfcrough_fgs
+          end if
+       end if ! i_sfcrough_fgs (reading surface roughness 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, &
@@ -2556,7 +2585,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 rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda, i_sfcrough_fgs
     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
@@ -3174,8 +3203,14 @@ subroutine update_netcdf_mass_wrf(this)
             start_index,end_index1,               & !mem
             start_index,end_index1,               & !pat
             ierr                                 )
-    endif  ! i_howv_3dda (wave height)
+    endif  ! i_gust_3dda (wind gust)
 
+!   Reading surface roughtness from binary siganl (only when i_sfcrough_fgs = 1),
+!    and no need to write it to netcdf analysis file, since it is not updated.
+    if ( i_sfcrough_fgs == 1 ) then
+       read(iunit)   field2   !ZNT (surface roughness)
+       if(print_verbose)write(6,*)'update_netcdf_mass_wrf:: max,min ZNT(sfc_rough)=',maxval(field2),minval(field2)
+    end if ! i_sfcrough_fgs (surface roughness)
 
     if (l_hydrometeor_bkio .and. n_actual_clouds>0) then
       do k=1,nsig_regional

src/gsi/cplr_wrwrfmassa.f90

@@ -1877,7 +1877,7 @@ subroutine wrwrfmassa_netcdf_wrf(this,mype)
     use rapidrefresh_cldsurf_mod, only: l_hydrometeor_bkio,l_gsd_soilTQ_nudge,&
          i_use_2mq4b,i_use_2mt4b
     use rapidrefresh_cldsurf_mod, only: i_howv_3dda, i_gust_3dda
-    use rapidrefresh_cldsurf_mod, only: i_howv_mask
+    use rapidrefresh_cldsurf_mod, only: i_howv_mask, i_sfcrough_fgs
 
     use chemmod, only: laeroana_gocart,wrf_pm2_5
     use gsi_bundlemod, only: GSI_BundleGetPointer
@@ -2971,6 +2971,15 @@ subroutine wrwrfmassa_netcdf_wrf(this,mype)
 
     end if ! i_gust_3dda > 0
 
+  ! ZNT: surface roughness (no updating to it)
+    if ( i_sfcrough_fgs > 0 ) then
+       if(mype == 0) then
+  !         This is surface roughness
+          read(lendian_in)temp1
+          write(lendian_out)temp1
+       end if
+    end if
+
   !
   ! for saving cloud analysis results
     if(l_hydrometeor_bkio .and. n_actual_clouds>0) then

src/gsi/gsimod.F90

@@ -184,7 +184,8 @@ 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_gust, hwllp_gust, oerr_gust, i_howv_mask
+                            corp_howv, hwllp_howv, corp_gust, hwllp_gust, oerr_gust, i_howv_mask, &
+                            i_sfcrough_fgs
   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
@@ -1619,6 +1620,10 @@ module gsimod
 !      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
+!      i_sfcrough_fgs - integer, option to control the read-in of surface roughness from firstguess
+!                           = 0 : do not read surface roughness from firstguess,
+!                                 and use the default value instead (default)
+!                           = 1 : read surface roughness from firstguess and use it in analysis
 !
   namelist/rapidrefresh_cldsurf/dfi_radar_latent_heat_time_period, &
                                 metar_impact_radius,metar_impact_radius_lowcloud, &
@@ -1640,7 +1645,8 @@ 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_gust, hwllp_gust, oerr_gust, i_howv_mask
+                                corp_howv, hwllp_howv, corp_gust, hwllp_gust, oerr_gust, i_howv_mask, &
+                                i_sfcrough_fgs
 
 ! chem(options for gsi chem analysis) :
 !     berror_chem       - .true. when background  for chemical species that require

src/gsi/rapidrefresh_cldsurf_mod.f90

@@ -221,6 +221,10 @@ module rapidrefresh_cldsurf_mod
 !                          = 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;
+!      i_sfcrough_fgs - integer, namelist option to control the read-in and usage of surface roughness in firstguess
+!                          = 0 : do not read surface roughness from firstguess,
+!                                and use the default value instead (default)
+!                          = 1 : read surface roughness from firstguess and use it in analysis
 !
 ! attributes:
 !   language: f90
@@ -297,6 +301,7 @@ module rapidrefresh_cldsurf_mod
   public :: i_howv_mask
   public :: corp_gust, hwllp_gust, oerr_gust
   public :: i_gust_3dda
+  public :: i_sfcrough_fgs
 
   logical l_hydrometeor_bkio
   real(r_kind)  dfi_radar_latent_heat_time_period
@@ -360,6 +365,7 @@ module rapidrefresh_cldsurf_mod
   integer(i_kind)   :: i_howv_mask
   real(r_kind)      :: corp_gust, hwllp_gust, oerr_gust
   integer(i_kind)   :: i_gust_3dda
+  integer(i_kind)   :: i_sfcrough_fgs
 
 contains
 
@@ -495,6 +501,8 @@ subroutine init_rapidrefresh_cldsurf
 
     i_gust_3dda         = 0                           ! no analysis of wind gust (gust) in 3D analysis (default)
 
+    i_sfcrough_fgs      = 0                           ! do not read surface roughness from firstguess (default)
+
 !-- searching for specific variable in state variable list (reading from anavinfo)
     do i2=1,ns2d
       if ( trim(svars2d(i2))=='howv' .or. trim(svars2d(i2))=='HOWV'   ) then

src/gsi/read_prepbufr.f90

@@ -452,7 +452,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
 
   integer:: icase,klev,ikkk,tkk
   real:: diffhgt,diffuu,diffvv
-  integer,dimension(3)::kcount
+  integer,dimension(5)::kcount
 
   real(r_double),dimension(3,1500):: fcstdat
   logical print_verbose
@@ -2066,6 +2066,12 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
                        call find_wind_height(c_prvstg,c_sprvstg,windsensht,kcount)
                     endif
                  endif
+                 !retrieve wind sensor height for mesonet gustob only when running 3DRTMA
+                 if (l_rtma3d)  then
+                    if ( gustob .and. (kx==188.or.kx==195) )  then
+                       call find_wind_height(c_prvstg,c_sprvstg,windsensht,kcount)
+                    endif
+                 endif
                  if(i_gsdqc==2) then  ! filter bad 2-m dew point and  0 mesonet wind obs
                     if (kx==288.or.kx==295) then ! for mesonet wind
                        if(abs(obsdat(5,k))<0.01_r_kind .and. abs(obsdat(6,k))<0.01_r_kind) usage=115._r_kind
@@ -2572,7 +2578,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
                    if ((kx==280).or.(kx==180)) oelev=r20+selev
                    if ((kx==299).or.(kx==199)) oelev=r20+selev
                    if ((kx==282).or.(kx==182)) oelev=r20+selev
-                   if (((kx==295).or.(kx==288).or.(kx==195).or.(kx==188)).and.twodvar_regional) then
+                   if (((kx==295).or.(kx==288).or.(kx==195).or.(kx==188)).and.(twodvar_regional.or.l_rtma3d)) then
                       !account for mesonet wind sensor height
                       oelev=windsensht+selev
                    end if
@@ -3309,7 +3315,7 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
   call destroy_aircraft_rjlists
   if(i_gsdsfc_uselist==1) call destroy_gsd_sfcuselist
   if (lhilbert) call destroy_hilbertcurve
-  if (twodvar_regional) then
+  if (twodvar_regional .or. l_rtma3d) then
      call destroy_ndfdgrid
      call destroy_windht_lists
   endif
@@ -3323,6 +3329,10 @@ subroutine read_prepbufr(nread,ndata,nodata,infile,obstype,lunout,twindin,sis,&
   if (twodvar_regional .and. (uvob .or. gustob .or. spdob))  then
     write(6,*) 'kcount values from find wind height = ',kcount
   end if
+  if ( l_rtma3d  .and. gustob )  then
+     write(6,'(1x,A,1x,A,A,5(1x,I9))') 'read_prepbufr:: for obstype=',            &
+       trim(adjustl(obstype)), ': kcount values from find wind height = ',kcount
+  end if