Modified observation errors for radar-reflectivity observations in pr…

…ecipitation (#650)

This PR addresses issue 649: Reflectivity high bias resulting from EnVar
radar-reflectivity data assimilation.

Observation errors for radar-reflectivity observations are changed as
follows:
(1) Increase default observation error standard deviation from 5.0 dBZ
to 10.0 dBZ for reflectivity observations in precipitation (i.e.,
observations >= 5 dBZ). For non-precipitation observations (< 5 dBZ),
keep the existing 5.0 dBZ error standard deviation.
(2) For reflectivity observations in precipitation, further inflate the
error for observations that fail the gross error check by a factor of
1.0-2.0. This change will be combined with a stricter gross error check,
implemented through a separate PR to the regional workflow. Also, the
gross error check won't be applied to non-precipitation reflectivity
observations.
The changes described above were initially discussed by David Dowell,
Jacob Carley, and Sho Yokota in emails on 11 August 2023.

The proposed changes were tested in a prototype CONUS RRFSv1 for a
summer 2022 retrospective period.

src/gsi/gsimod.F90

@@ -26,7 +26,7 @@ module gsimod
   use obsmod, only: doradaroneob,dofedoneob,oneoblat,oneoblon,oneobheight,oneobvalue,oneobddiff,oneobradid,&
      radar_no_thinning,ens_hx_dbz_cut,static_gsi_nopcp_dbz,rmesh_dbz,&
      rmesh_vr,zmesh_dbz,zmesh_vr,if_vterminal, if_model_dbz,if_model_fed,innov_use_model_fed,if_vrobs_raw,if_use_w_vr,&
-     minobrangedbz,maxobrangedbz,maxobrangevr,maxtiltvr,missing_to_nopcp,&
+     minobrangedbz,maxobrangedbz,maxobrangevr,maxtiltvr,inflate_dbz_obserr,missing_to_nopcp,&
      ntilt_radarfiles,whichradar,&
      minobrangevr,maxtiltdbz,mintiltvr,mintiltdbz,l2rwthin,hurricane_radar,&
      r_hgt_fed
@@ -744,6 +744,10 @@ module gsimod
 !     optconv - downweighting option for iasi and cris for moisture channels to
 !     improve convergence.  default 0.0 (no change).  Larger number improves
 !     convergence.
+!     inflate_dbz_obserr - logical that controls inflation of reflectivity ob error
+!                          for obs that exceed gross error magnitude
+!                          if true, inflate ob error
+!                          if false, reject ob
 !
 !     NOTE:  for now, if in regional mode, then iguess=-1 is forced internally.
 !            add use of guess file later for regional mode.
@@ -779,7 +783,7 @@ module gsimod
        use_gfs_stratosphere,pblend0,pblend1,step_start,diag_precon,lrun_subdirs,&
        use_sp_eqspace,lnested_loops,lsingleradob,thin4d,use_readin_anl_sfcmask,&
        luse_obsdiag,id_drifter,id_ship,verbose,print_obs_para,lsingleradar,singleradar,lnobalance, &
-       missing_to_nopcp,minobrangedbz,minobrangedbz,maxobrangedbz,&
+       inflate_dbz_obserr,missing_to_nopcp,minobrangedbz,minobrangedbz,maxobrangedbz,&
        maxobrangevr,maxtiltvr,whichradar,doradaroneob,dofedoneob,oneoblat,&
        oneoblon,oneobheight,oneobvalue,oneobddiff,oneobradid,&
        rmesh_vr,zmesh_dbz,zmesh_vr, ntilt_radarfiles, whichradar,&

src/gsi/obsmod.F90

@@ -480,7 +480,8 @@ module obsmod
   ! ==== DBZ DA ===
   public :: ntilt_radarfiles
   public :: whichradar
-  public :: vr_dealisingopt, if_vterminal, if_model_dbz, inflate_obserr, if_vrobs_raw, if_use_w_vr, l2rwthin
+  public :: vr_dealisingopt, if_vterminal, if_model_dbz, if_vrobs_raw, if_use_w_vr, l2rwthin
+  public :: inflate_dbz_obserr
 
   public :: doradaroneob,oneoblat,oneoblon
   public :: oneobddiff,oneobvalue,oneobheight,oneobradid
@@ -634,7 +635,8 @@ module obsmod
 
   logical ::  ta2tb
   logical ::  doradaroneob,dofedoneob
-  logical :: vr_dealisingopt, if_vterminal, if_model_dbz,if_model_fed, innov_use_model_fed,inflate_obserr, if_vrobs_raw, if_use_w_vr, l2rwthin
+  logical :: vr_dealisingopt, if_vterminal, if_model_dbz,if_model_fed, innov_use_model_fed, if_vrobs_raw, if_use_w_vr, l2rwthin
+  logical :: inflate_dbz_obserr
   character(4) :: whichradar,oneobradid
   real(r_kind) :: oneoblat,oneoblon,oneobddiff,oneobvalue,oneobheight
   logical :: radar_no_thinning
@@ -769,7 +771,7 @@ subroutine init_obsmod_dflts
     if_model_dbz=.false.
     if_model_fed=.false.
     innov_use_model_fed=.false.
-    inflate_obserr=.false.
+    inflate_dbz_obserr=.false.
     whichradar="KKKK"
 
     oneobradid="KKKK"

src/gsi/read_dbz_nc.f90

@@ -412,9 +412,15 @@ subroutine read_dbz_nc(nread,ndata,nodata,infile,lunout,obstype,sis,hgtl_full,no
 ! changed to hard-coded value for now; dbznoise used for two different purposes in this subroutine:
 !                   (1) threshold for lowest reflectivity value considered to be an observation and 
 !                   (2) ob error
-        thiserr = 5.0_r_kind
-
-
+
+!       Specify a larger error standard deviation for reflectivity observations in precipitation
+!       than for reflectivity observations that indicate a lack of preciptation.
+        if( dbzQC(i,j,k) < 5.0_r_kind ) then
+          thiserr = 5.0_r_kind
+        else
+          thiserr = 10.0_r_kind
+        end if
+
         nread = nread + 1
 
  !####################       Data thinning       ###################

src/gsi/setupdbz.f90

@@ -150,7 +150,7 @@ subroutine setupdbz(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,radardbz_d
   use convinfo, only: nconvtype,cermin,cermax,cgross,cvar_b,cvar_pg,ictype
   use convinfo, only: icsubtype
   use m_dtime, only: dtime_setup, dtime_check
-  use obsmod, only   : if_model_dbz, inflate_obserr
+  use obsmod, only   : if_model_dbz, inflate_dbz_obserr
   use setupdbz_lib, only:hx_dart,jqr_dart,jqs_dart,jqg_dart 
   use gridmod, only: wrf_mass_regional,nems_nmmb_regional, fv3_regional
   use sparsearr, only: sparr2, new, size, writearray, fullarray
@@ -1260,8 +1260,11 @@ subroutine setupdbz(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,radardbz_d
          end if
 
      else
-        if (ratio > cgross(ikx) .or. ratio_errors < tiny_r_kind) then
-           if ( inflate_obserr .and. (ratio-cgross(ikx)) <= cgross(ikx) .and. ratio_errors >= tiny_r_kind) then 
+
+!       Apply gross error check only to reflectivity observations in precipitation (>= 5 dBZ).
+        if ( ( (data(idbzob,i) >= 5_r_kind) .and. (ratio > cgross(ikx)) ) .or. (ratio_errors < tiny_r_kind) ) then
+
+           if ( inflate_dbz_obserr .and. (ratio-cgross(ikx)) <= cgross(ikx) .and. ratio_errors >= tiny_r_kind) then 
            ! Since radar reflectivity can be very different from the model background
            ! good observations may be rejected during this QC step.  However, if these observations
            ! are allowed through, they can yield problems with convergence.  Therefore the error