From d5095f935a985fa7a16f34826c21b62394585d6c Mon Sep 17 00:00:00 2001 From: Eun Hee Lee Date: Mon, 24 Jan 2022 22:07:40 -0500 Subject: [PATCH] Updated AMV QC and thinning grids --- .../GSI_GridComp/etc/gmao_global_convinfo.txt | 26 ++++++------ .../GSI_GridComp/etc/gsi.rc.tmpl | 4 +- .../GSI_GridComp/read_satwnd.f90 | 37 +++++++++++++---- GEOSaana_GridComp/GSI_GridComp/setupw.f90 | 41 ++++++++++--------- 4 files changed, 66 insertions(+), 42 deletions(-) diff --git a/GEOSaana_GridComp/GSI_GridComp/etc/gmao_global_convinfo.txt b/GEOSaana_GridComp/GSI_GridComp/etc/gmao_global_convinfo.txt index 8f47fe34..a850eec1 100644 --- a/GEOSaana_GridComp/GSI_GridComp/etc/gmao_global_convinfo.txt +++ b/GEOSaana_GridComp/GSI_GridComp/etc/gmao_global_convinfo.txt @@ -121,12 +121,12 @@ uv 243 57 1 3.0 0 0 0 1.5 15.0 1.4 1.5 0.055000 1 200. 100. 0 0. 0. uv 243 70 1 3.0 0 0 0 1.5 15.0 1.4 1.5 0.055000 1 200. 100. 0 0. 0. uv 244 0 -1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. - uv 244 3 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. - uv 244 4 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. + uv 244 3 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 1 200. 100. 0 0. 0. + uv 244 4 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 1 200. 100. 0 0. 0. uv 244 206 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. uv 244 207 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. - uv 244 209 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. - uv 244 223 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. + uv 244 209 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 1 125. 0. 0 0. 0. + uv 244 223 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 1 125. 0. 0 0. 0. uv 245 0 -1 3.0 0 0 0 2.5 20.0 1.4 1.3 0.005000 1 200. 100. 0 0. 0. uv 245 257 1 3.0 0 0 0 2.5 20.0 1.4 1.3 0.005000 1 200. 100. 0 0. 0. uv 245 259 1 3.0 0 0 0 2.5 20.0 1.4 1.3 0.005000 1 200. 100. 0 0. 0. @@ -140,8 +140,8 @@ uv 247 0 -1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. uv 247 257 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. uv 247 259 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. - uv 247 270 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. - uv 247 271 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 0 0. 0. 0 0. 0. + uv 247 270 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 1 200. 100. 0 0. 0. + uv 247 271 1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.005000 1 200. 100. 0 0. 0. uv 248 0 -1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.000500 0 0. 0. 0 0. 0. uv 249 0 -1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.050500 0 0. 0. 0 0. 0. uv 250 0 -1 3.0 0 0 0 2.5 20.0 1.4 2.5 0.050500 0 0. 0. 0 0. 0. @@ -173,20 +173,20 @@ uv 254 70 1 3.0 0 0 0 1.5 20.0 1.4 1.5 0.050500 1 200. 100. 0 0. 0. uv 256 0 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.000500 0 0. 0. 0 0. 0. uv 257 0 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. - uv 257 783 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. - uv 257 784 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. + uv 257 783 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 1 75. 50. 0 0. 0. + uv 257 784 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 1 100. 50. 0 0. 0. uv 258 0 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. uv 258 783 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. - uv 258 784 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. + uv 258 784 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 1 125. 0. 0 0. 0. uv 259 0 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. uv 259 783 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. - uv 259 784 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 0 0. 0. 0 0. 0. + uv 259 784 1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.005500 1 100. 0. 0 0. 0. uv 260 0 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.000500 0 0. 0. 0 0. 0. uv 260 224 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.000500 0 0. 0. 0 0. 0. uv 260 225 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.000500 0 0. 0. 0 0. 0. uv 270 0 -1 3.0 0 0 0 2.5 20.1 1.4 2.5 0.000500 0 0. 0. 0 0. 0. - uv 280 0 1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 0 0. 0. 0 0. 0. - uv 280 1 1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 0 0. 0. 0 0. 0. + uv 280 0 1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 1 50. 0. 0 0. 0. + uv 280 1 1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 1 50. 0. 0 0. 0. uv 281 0 -1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 0 0. 0. 0 0. 0. uv 282 0 1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 0 0. 0. 0 0. 0. uv 284 0 -1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 0 0. 0. 0 0. 0. @@ -194,7 +194,7 @@ uv 286 0 -1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.005500 0 0. 0. 0 0. 0. uv 287 0 -1 3.0 0 0 0 6.0 6.1 1.4 6.0 0.000500 0 0. 0. 0 0. 0. uv 289 0 1 3.0 0 0 0 5.0 6.1 1.4 5.0 0.000500 0 0. 0. 0 0. 0. - uv 290 4 -1 3.0 0 0 0 5.0 6.1 1.4 5.0 0.000500 1 100. 1200. 0 0. 0. + uv 290 4 -1 3.0 0 0 0 5.0 6.1 1.4 5.0 0.000500 1 150. 1200. 0 0. 0. uv 290 3 1 3.0 0 0 0 5.0 6.1 1.4 5.0 0.000500 1 100. 1200. 0 0. 0. uv 290 5 -1 3.0 0 0 0 5.0 6.1 1.4 5.0 0.000500 1 100. 1200. 0 0. 0. uv 291 0 1 3.0 0 0 0 5.0 6.1 1.4 5.0 0.000500 1 100. 1200. 0 0. 0. diff --git a/GEOSaana_GridComp/GSI_GridComp/etc/gsi.rc.tmpl b/GEOSaana_GridComp/GSI_GridComp/etc/gsi.rc.tmpl index b23162c1..ac4b71f4 100755 --- a/GEOSaana_GridComp/GSI_GridComp/etc/gsi.rc.tmpl +++ b/GEOSaana_GridComp/GSI_GridComp/etc/gsi.rc.tmpl @@ -59,7 +59,7 @@ / &OBSQC dfact=0.75,dfact1=3.0,noiqc=.false.,oberrflg=.true.,c_varqc=0.02,blacklst=.true., - use_poq7=.true.,qc_noirjaco3=.false.,qc_satwnds=.false.,cld_det_dec2bin=.true., + use_poq7=.true.,qc_noirjaco3=.false.,qc_satwnds=.true.,cld_det_dec2bin=.true., half_goesr_err=.false., ! tcp_ermin=0.75,tcp_ermax=0.75, >>>AIRCFT_BIAS<<< @@ -178,7 +178,7 @@ OBS_INPUT:: avcsambufr avhrr metop-b avhrr3_metop-b 0.0 1 0 ncep_avcsam_bufr ! omieffbufr omieff aura omieff_aura 0.0 2 0 test omieffnc omieff aura omieff_aura 0.0 2 0 aura_omieff_nc - ompsnmeffnc ompsnmeff npp ompsnmeff_npp 0.0 2 0 npp_ompsnmeff_nc + ompsnmeffnc ompsnmeff npp ompsnmeff_npp 1.0 2 0 npp_ompsnmeff_nc gsnd1bufr sndrd1 g15 sndrD1_g15 0.0 1 0 ncep_goesfv_bufr gsnd1bufr sndrd2 g15 sndrD2_g15 0.0 1 0 ncep_goesfv_bufr gsnd1bufr sndrd3 g15 sndrD3_g15 0.0 1 0 ncep_goesfv_bufr diff --git a/GEOSaana_GridComp/GSI_GridComp/read_satwnd.f90 b/GEOSaana_GridComp/GSI_GridComp/read_satwnd.f90 index 108f067f..5331d009 100644 --- a/GEOSaana_GridComp/GSI_GridComp/read_satwnd.f90 +++ b/GEOSaana_GridComp/GSI_GridComp/read_satwnd.f90 @@ -70,7 +70,8 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis ! 2018-06-13 Genkova - Goes-16 AMVs use ECMWF QC till new HAM late 2018 ! and OE/2 ! 2021-06-21 Todling - Allow code to bypass hack to half-GOES-R errors -! +! 2021-11-15 Eunhee add QC for MSG IR winds and modified QC for IR winds over land and top air winds removal +! 2022-01-21 Eunhee assign SAZA, SWCM(ich) for monitoring ! ! ! input argument list: @@ -139,6 +140,7 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis real(r_kind),parameter:: r105= 105.0_r_kind real(r_kind),parameter:: r110= 110.0_r_kind real(r_kind),parameter:: r125=125.0_r_kind + real(r_kind),parameter:: r150=150.0_r_kind real(r_kind),parameter:: r200=200.0_r_kind real(r_kind),parameter:: r250=250.0_r_kind real(r_kind),parameter:: r360 = 360.0_r_kind @@ -174,7 +176,7 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis integer(i_kind) nmind,lunin,idate,ilat,ilon,iret,k integer(i_kind) nreal,ithin,iout,ntmp,icount,iiout,ii integer(i_kind) itype,iosub,ixsub,isubsub,iobsub,itypey,ierr - integer(i_kind) qm + integer(i_kind) qm, ich, hamd integer(i_kind) nlevp ! vertical level for thinning integer(i_kind) pflag integer(i_kind) ntest,nvtest @@ -195,7 +197,7 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis real(r_kind) rmesh,ediff,usage,tdiff real(r_kind) u0,v0,uob,vob,dx,dy,dx1,dy1,w00,w10,w01,w11 real(r_kind) dlnpob,ppb,ppb2,qifn,qify,ee,ree,pct1,experr_norm - real(r_kind) woe,dlat,dlon,dlat_earth,dlon_earth + real(r_kind) woe,dlat,dlon,dlat_earth,dlon_earth,saza real(r_kind) dlat_earth_deg,dlon_earth_deg real(r_kind) cdist,disterr,disterrmax,rlon00,rlat00 real(r_kind) vdisterrmax,u00,v00,uob1,vob1 @@ -264,7 +266,7 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis ! Set lower limits for observation errors werrmin=one nsattype=0 - nreal=25 + nreal=28 if(perturb_obs ) nreal=nreal+2 ntread=1 ntmatch=0 @@ -599,6 +601,8 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis ee=r110 qifn=r110 qify=r110 + ich=-1 + hamd=-1 ! Test for BUFR version using lat/lon mnemonics call ufbint(lunin,hdrdat_test,2,1,iret, 'CLAT CLON') @@ -615,7 +619,8 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis hdrdat(3) >100000000.0_r_kind .or. & obsdat(4) > 100000000.0_r_kind) cycle loop_readsb if(ppb >r10000) ppb=ppb/r100 - if (ppb 20.0_r_kind) then call deter_sfc_type(dlat_earth,dlon_earth,t4dv,isflg,tsavg) if(isflg /= 0) cycle loop_readsb + else + call deter_sfc_type(dlat_earth,dlon_earth,t4dv,isflg,tsavg) + if (isflg /= 0 .and. ppb > 700.0_r_kind) cycle loop_readsb ! low over land + !! if (obsdat(1) > 700.0_r_kind) then endif endif - endif + + if(itype ==253) then + if(hdrdat(2) >5.0_r_kind) then + call deter_sfc_type(dlat_earth,dlon_earth,t4dv,isflg,tsavg) + if(isflg /= 0 ) cycle loop_readsb + endif + endif + endif !! convert from wind direction and speed to u,v component uob=-obsdat(4)*sin(obsdat(3)*deg2rad) @@ -1453,10 +1471,13 @@ subroutine read_satwnd(nread,ndata,nodata,infile,obstype,lunout,gstime,twind,sis cdata_all(22,iout)=r_prvstg(1,1) ! provider name cdata_all(23,iout)=r_sprvstg(1,1) ! subprovider name cdata_all(25,iout)=var_jb ! non linear qc parameter + cdata_all(26,iout)=ich ! spectral type(channel) of wind + cdata_all(27,iout)=saza ! sat zenith angle + cdata_all(28,iout)=hamd ! height assignment method if(perturb_obs)then - cdata_all(26,iout)=ran01dom()*perturb_fact ! u perturbation - cdata_all(27,iout)=ran01dom()*perturb_fact ! v perturbation + cdata_all(29,iout)=ran01dom()*perturb_fact ! u perturbation + cdata_all(30,iout)=ran01dom()*perturb_fact ! v perturbation endif enddo loop_readsb diff --git a/GEOSaana_GridComp/GSI_GridComp/setupw.f90 b/GEOSaana_GridComp/GSI_GridComp/setupw.f90 index 8da76b1e..b36addf7 100644 --- a/GEOSaana_GridComp/GSI_GridComp/setupw.f90 +++ b/GEOSaana_GridComp/GSI_GridComp/setupw.f90 @@ -191,6 +191,7 @@ subroutine setupw(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav ! 2018-04-09 pondeca - introduce duplogic to correctly handle the characterization of ! duplicate obs in twodvar_regional applications ! 2020-02-26 todling - reset obsbin from hr to min +! 2021-11-15 Eunhee - Remove the QC for ir winds in the mid atmospheric layer if qc_satwnds=true ! ! ! REMARKS: @@ -263,7 +264,7 @@ subroutine setupw(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav integer(i_kind) idomsfc,isfcr,iskint,iff10 integer(i_kind) msges - integer(i_kind) iswcm,isaza, isccf, qify, qifn + integer(i_kind) iswcm,isaza, ihamd, isccf, qify, qifn real(r_kind) sccf_wavelen real(r_kind),parameter:: rsol=300000000.0_r_kind !speed of light real(r_kind),parameter:: rtomic=1000000.0_r_kind !conv to micron @@ -357,7 +358,8 @@ subroutine setupw(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav iswcm=26 ! spectral type of wind isaza=27 ! saza satellite zen angle - isccf=28 ! spec chan freq (hertz) + ihamd=28 ! height assignment method +! isccf=28 ! spec chan freq (hertz) iptrbu=29 ! index of u perturbation iptrbv=30 ! index of v perturbation @@ -889,20 +891,20 @@ subroutine setupw(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav if(itype ==242 .or. itype ==243 ) then ! visible winds from JMA and EUMETSAT if(presw <700.0_r_kind) error=zero ! no visible winds above 700mb endif - if(itype ==245 ) then - if( presw >399.0_r_kind .and. presw <801.0_r_kind) then !GOES IR winds - error=zero ! no data between 400-800mb - endif - endif - if(itype == 252 .and. presw >499.0_r_kind .and. presw <801.0_r_kind) then ! JMA IR winds - error=zero - endif - if(itype == 253 ) then - if(presw >401.0_r_kind .and. presw <801.0_r_kind) then ! EUMET IR winds - error=zero - endif - endif - if( itype == 246 .or. itype == 250 .or. itype == 254 ) then ! water vapor cloud top +! if(itype ==245 ) then +! if( presw >399.0_r_kind .and. presw <801.0_r_kind) then !GOES IR winds +! error=zero ! no data between 400-800mb +! endif +! endif +! if(itype == 252 .and. presw >499.0_r_kind .and. presw <801.0_r_kind) then ! JMA IR winds +! error=zero +! endif +! if(itype == 253 ) then +! if(presw >401.0_r_kind .and. presw <801.0_r_kind) then ! EUMET IR winds +! error=zero +! endif +! endif + if( itype == 246 .or. itype ==247 .or. itype == 250 .or. itype == 254 ) then ! water vapor cloud top if(presw >399.0_r_kind) error=zero endif if(itype ==257 .and. presw <249.0_r_kind) error=zero @@ -1764,8 +1766,9 @@ subroutine contents_netcdf_diag_(udiag,vdiag) ! Write out in nc diag the extra vars from cdata_all (see read_satwnd.f90) call nc_diag_metadata("SWCM_spec_type", sngl(data(iswcm,i))) call nc_diag_metadata("SAZA_sat_zen_angle", sngl(data(isaza,i))) - sccf_wavelen=(rsol/data(isccf,i))*rtomic !spec chan wavelen(microns) - call nc_diag_metadata("SCCF_chan_wavelen", sngl(sccf_wavelen)) + call nc_diag_metadata("Height_assignment", sngl(data(ihamd,i))) +! sccf_wavelen=(rsol/data(isccf,i))*rtomic !spec chan wavelen(microns) +! call nc_diag_metadata("SCCF_chan_wavelen", sngl(sccf_wavelen)) qify= int(data(ielev,i)/1000.0); qifn= mod(data(ielev,i),1000.0); call nc_diag_metadata("QI_with_FC", sngl(qify)) @@ -1774,7 +1777,7 @@ subroutine contents_netcdf_diag_(udiag,vdiag) ! Write out missing values) call nc_diag_metadata("SWCM_spec_type", sngl(bmiss)) call nc_diag_metadata("SAZA_sat_zen_angle", sngl(bmiss)) - call nc_diag_metadata("SCCF_chan_wavelen", sngl(bmiss)) + call nc_diag_metadata("Height_assignment", sngl(bmiss)) call nc_diag_metadata("QI_with_FC", sngl(bmiss)) call nc_diag_metadata("QI_without_FC", sngl(bmiss)) endif