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

Add assimilation of GLM flash-extent density (FED) observations to EnKF #654

Merged
merged 2 commits into from
Dec 7, 2023
Merged

Add assimilation of GLM flash-extent density (FED) observations to EnKF #654

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
92126e3
Add assimilation of GLM flash-extent density (FED) observations, with…
daviddowellNOAA Nov 3, 2023
431acca
Merge branch 'NOAA-EMC:develop' into GSI_EnKF_FED
daviddowellNOAA Dec 5, 2023
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
2 changes: 1 addition & 1 deletion src/enkf/enkf.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -51,7 +51,7 @@ module enkf
! NH, tropics and SH, and in the horizontal, vertical and time dimensions,
! using the namelist parameters corrlengthnh, corrlengthtr, corrlengthsh,
! lnsigcutoffnh, lnsigcutofftr, lnsigcutoffsh (lnsigcutoffsatnh,
! lnsigcutoffsattr, lnsigcutoffsatsh for satellite obs, similar for ps obs)
! lnsigcutoffsattr, lnsigcutoffsatsh for satellite obs, similar for ps and fed obs)
! obtimelnh, obtimeltr, obtimelsh. The length scales should be given in km for the
! horizontal, hours for time, and 'scale heights' (units of -log(p/pref)) in the
! vertical. The function used for localization (function taper)
Expand Down
3 changes: 3 additions & 0 deletions src/enkf/enkf_obs_sensitivity.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@ module enkf_obs_sensitivity
use params, only: efsoi_flag,latbound,nlevs,nanals,datestring, &
lnsigcutoffsatnh,lnsigcutoffsattr,lnsigcutoffsatsh, &
lnsigcutoffpsnh,lnsigcutoffpstr,lnsigcutoffpssh, &
lnsigcutofffednh,lnsigcutofffedtr,lnsigcutofffedsh, &
lnsigcutoffnh,lnsigcutofftr,lnsigcutoffsh, &
corrlengthnh,corrlengthtr,corrlengthsh, &
obtimelnh,obtimeltr,obtimelsh,letkf_flag, &
Expand Down Expand Up @@ -292,6 +293,8 @@ subroutine read_ob_sens
lnsigl(nob) = latval(deglat,lnsigcutoffsatnh,lnsigcutoffsattr,lnsigcutoffsatsh)
else if (obtype(nob)(1:3) == ' ps') then
lnsigl(nob) = latval(deglat,lnsigcutoffpsnh,lnsigcutoffpstr,lnsigcutoffpssh)
else if (obtype(nob)(1:3) == 'fed') then
lnsigl(nob) = latval(deglat,lnsigcutofffednh,lnsigcutofffedtr,lnsigcutofffedsh)
else
lnsigl(nob)=latval(deglat,lnsigcutoffnh,lnsigcutofftr,lnsigcutoffsh)
end if
Expand Down
7 changes: 7 additions & 0 deletions src/enkf/enkf_obsmod.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,8 @@ module enkf_obsmod
lnsigcutoffnh, lnsigcutoffsh, lnsigcutofftr, corrlengthnh,&
corrlengthtr, corrlengthsh, obtimelnh, obtimeltr, obtimelsh,&
lnsigcutoffsatnh, lnsigcutoffsatsh, lnsigcutoffsattr,&
lnsigcutofffednh, lnsigcutofffedsh, lnsigcutofffedtr,&
corrlengthfednh, corrlengthfedtr, corrlengthfedsh, &
varqc, huber, zhuberleft, zhuberright, modelspace_vloc, &
lnsigcutoffpsnh, lnsigcutoffpssh, lnsigcutoffpstr, neigv, &
lnsigcutoffrdrnh, lnsigcutoffrdrsh, lnsigcutoffrdrtr,&
Expand Down Expand Up @@ -276,6 +278,8 @@ subroutine readobs()
lnsigl(nob) = latval(deglat,lnsigcutoffsatnh,lnsigcutoffsattr,lnsigcutoffsatsh)
else if (obtype(nob)(1:3) == ' ps') then
lnsigl(nob) = latval(deglat,lnsigcutoffpsnh,lnsigcutoffpstr,lnsigcutoffpssh)
else if (obtype(nob)(1:3) == 'fed') then
lnsigl(nob) = latval(deglat,lnsigcutofffednh,lnsigcutofffedtr,lnsigcutofffedsh)
else if ( (obtype(nob)(1:3) == 'dbz' .or. obtype(nob)(1:3) == ' rw') .and. l_use_enkf_directZDA ) then
lnsigl(nob) = latval(deglat,lnsigcutoffrdrnh,lnsigcutoffrdrtr,lnsigcutoffrdrsh)
else
Expand All @@ -293,6 +297,9 @@ subroutine readobs()
if ( (obtype(nob)(1:3) == 'dbz' .or. obtype(nob)(1:3) == ' rw') .and. l_use_enkf_directZDA ) then
corrlengthsq(nob)=latval(deglat,corrlengthrdrnh,corrlengthrdrtr,corrlengthrdrsh)**2
end if
if (obtype(nob)(1:3) == 'fed') then
corrlengthsq(nob)=latval(deglat,corrlengthfednh,corrlengthfedtr,corrlengthfedsh)**2
end if
obtimel(nob)=latval(deglat,obtimelnh,obtimeltr,obtimelsh)
end do

Expand Down
16 changes: 16 additions & 0 deletions src/enkf/innovstats.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,7 @@ subroutine print_innovstats(obfit,obsprd)
nobsspd_nh,nobsspd_sh,nobsspd_tr,&
nobsgps_nh,nobsgps_sh,nobsgps_tr,&
nobsdbz_nh,nobsdbz_sh,nobsdbz_tr,&
nobsfed_nh,nobsfed_sh,nobsfed_tr,&
nobsrw_nh,nobsrw_sh,nobsrw_tr,&
nobsq_nh,nobsq_sh,nobsq_tr,nobswnd_nh,nobswnd_sh,nobswnd_tr,&
nobsoz_nh,nobsoz_sh,nobsoz_tr,nobsps_sh,nobsps_nh,nobsps_tr,nob
Expand All @@ -67,6 +68,9 @@ subroutine print_innovstats(obfit,obsprd)
sumdbz_nh,biasdbz_nh,sumdbz_spread_nh,sumdbz_oberr_nh,&
sumdbz_sh,biasdbz_sh,sumdbz_spread_sh,sumdbz_oberr_sh,&
sumdbz_tr,biasdbz_tr,sumdbz_spread_tr,sumdbz_oberr_tr,&
sumfed_nh,biasfed_nh,sumfed_spread_nh,sumfed_oberr_nh,&
sumfed_sh,biasfed_sh,sumfed_spread_sh,sumfed_oberr_sh,&
sumfed_tr,biasfed_tr,sumfed_spread_tr,sumfed_oberr_tr,&
sumrw_nh,biasrw_nh,sumrw_spread_nh,sumrw_oberr_nh,&
sumrw_sh,biasrw_sh,sumrw_spread_sh,sumrw_oberr_sh,&
sumrw_tr,biasrw_tr,sumrw_spread_tr,sumrw_oberr_tr,&
Expand Down Expand Up @@ -112,6 +116,9 @@ subroutine print_innovstats(obfit,obsprd)
nobsdbz_nh = 0
nobsdbz_sh = 0
nobsdbz_tr = 0
nobsfed_nh = 0
nobsfed_sh = 0
nobsfed_tr = 0
nobsrw_nh = 0
nobsrw_sh = 0
nobsrw_tr = 0
Expand Down Expand Up @@ -168,6 +175,12 @@ subroutine print_innovstats(obfit,obsprd)
sumdbz_nh,biasdbz_nh,sumdbz_spread_nh,sumdbz_oberr_nh,nobsdbz_nh,&
sumdbz_sh,biasdbz_sh,sumdbz_spread_sh,sumdbz_oberr_sh,nobsdbz_sh,&
sumdbz_tr,biasdbz_tr,sumdbz_spread_tr,sumdbz_oberr_tr,nobsdbz_tr)
else if (obtype(nob)(1:3) == 'fed') then
call obstats(obfit(nob),oberrvar_orig(nob),&
obsprd(nob),obloclat(nob),&
sumfed_nh,biasfed_nh,sumfed_spread_nh,sumfed_oberr_nh,nobsfed_nh,&
sumfed_sh,biasfed_sh,sumfed_spread_sh,sumfed_oberr_sh,nobsfed_sh,&
sumfed_tr,biasfed_tr,sumfed_spread_tr,sumfed_oberr_tr,nobsfed_tr)
else if (obtype(nob)(1:3) == ' rw') then
call obstats(obfit(nob),oberrvar_orig(nob),&
obsprd(nob),obloclat(nob),&
Expand Down Expand Up @@ -216,6 +229,9 @@ subroutine print_innovstats(obfit,obsprd)
call printstats(' all dbz',sumdbz_nh,biasdbz_nh,sumdbz_spread_nh,sumdbz_oberr_nh,nobsdbz_nh,&
sumdbz_sh,biasdbz_sh,sumdbz_spread_sh,sumdbz_oberr_sh,nobsdbz_sh,&
sumdbz_tr,biasdbz_tr,sumdbz_spread_tr,sumdbz_oberr_tr,nobsdbz_tr)
call printstats(' all fed',sumfed_nh,biasfed_nh,sumfed_spread_nh,sumfed_oberr_nh,nobsfed_nh,&
sumfed_sh,biasfed_sh,sumfed_spread_sh,sumfed_oberr_sh,nobsfed_sh,&
sumfed_tr,biasfed_tr,sumfed_spread_tr,sumfed_oberr_tr,nobsfed_tr)
call printstats(' all rw',sumrw_nh,biasq_nh,sumrw_spread_nh,sumrw_oberr_nh,nobsrw_nh,&
sumrw_sh,biasrw_sh,sumrw_spread_sh,sumrw_oberr_sh,nobsrw_sh,&
sumrw_tr,biasrw_tr,sumrw_spread_tr,sumrw_oberr_tr,nobsrw_tr)
Expand Down
15 changes: 15 additions & 0 deletions src/enkf/params.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -124,6 +124,8 @@ module params
real(r_single),public :: lnsigcutoffnh,lnsigcutofftr,lnsigcutoffsh,&
lnsigcutoffsatnh,lnsigcutoffsattr,lnsigcutoffsatsh,&
lnsigcutoffpsnh,lnsigcutoffpstr,lnsigcutoffpssh
real(r_single),public :: corrlengthfednh,corrlengthfedtr,corrlengthfedsh, &
lnsigcutofffednh,lnsigcutofffedtr,lnsigcutofffedsh
real(r_single),public :: corrlengthrdrnh,corrlengthrdrtr,corrlengthrdrsh, &
lnsigcutoffrdrnh,lnsigcutoffrdrtr,lnsigcutoffrdrsh
real(r_single),public :: analpertwtnh,analpertwtsh,analpertwttr,sprd_tol,saterrfact
Expand Down Expand Up @@ -261,6 +263,8 @@ module params
lnsigcutoffnh,lnsigcutofftr,lnsigcutoffsh,&
lnsigcutoffsatnh,lnsigcutoffsattr,lnsigcutoffsatsh,&
lnsigcutoffpsnh,lnsigcutoffpstr,lnsigcutoffpssh,&
corrlengthfednh,corrlengthfedsh,corrlengthfedtr,&
lnsigcutofffednh,lnsigcutofffedsh,lnsigcutofffedtr,&
fgfileprefixes,fgsfcfileprefixes,anlfileprefixes, &
anlsfcfileprefixes,incfileprefixes,incsfcfileprefixes,&
statefileprefixes,statesfcfileprefixes, &
Expand Down Expand Up @@ -317,6 +321,10 @@ subroutine read_namelist()
corrlengthrdrnh = 10
corrlengthrdrtr = 10
corrlengthrdrsh = 10
! corrlength (km) for GLM flash extent density
corrlengthfednh = 30_r_single
corrlengthfedtr = 30_r_single
corrlengthfedsh = 30_r_single
! read in localization length scales from an external file.
readin_localization = .false.
! min and max inflation.
Expand All @@ -341,6 +349,9 @@ subroutine read_namelist()
lnsigcutoffrdrnh = 0.2_r_single ! value for radar
lnsigcutoffrdrtr = 0.2_r_single ! value for radar
lnsigcutoffrdrsh = 0.2_r_single ! value for radar
lnsigcutofffednh = 2._r_single ! value for GLM flash extent density
lnsigcutofffedtr = 2._r_single ! value for GLM flash extent density
lnsigcutofffedsh = 2._r_single ! value for GLM flash extent density
! ob time localization
obtimelnh = 1.e10_r_single
obtimeltr = 1.e10_r_single
Expand Down Expand Up @@ -813,6 +824,10 @@ subroutine read_namelist()
corrlengthrdrnh = corrlengthrdrnh * 1.e3_r_single/rearth
corrlengthrdrtr = corrlengthrdrtr * 1.e3_r_single/rearth
corrlengthrdrsh = corrlengthrdrsh * 1.e3_r_single/rearth
! rescale covariance localization length for GLM FED
corrlengthfednh = corrlengthfednh * 1.e3_r_single/rearth
corrlengthfedtr = corrlengthfedtr * 1.e3_r_single/rearth
corrlengthfedsh = corrlengthfedsh * 1.e3_r_single/rearth

! convert targe area boundary into radians
tar_minlat = tar_minlat * deg2rad
Expand Down
12 changes: 9 additions & 3 deletions src/enkf/readconvobs.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -42,9 +42,9 @@ module readconvobs


!> observation types to read from netcdf files
integer(i_kind), parameter :: nobtype = 11
integer(i_kind), parameter :: nobtype = 12
character(len=3), dimension(nobtype), parameter :: obtypes = (/' t', ' q', ' ps', ' uv', 'tcp', &
'gps', 'spd', ' pw', ' dw', ' rw', 'dbz' /)
'gps', 'spd', ' pw', ' dw', ' rw', 'dbz', 'fed' /)

contains

Expand Down Expand Up @@ -79,7 +79,7 @@ subroutine get_num_convobs_bin(obspath,datestring,num_obs_tot,num_obs_totdiag,id
integer(i_kind) :: iunit, nchar, nreal, ii, mype, ios, idate, i, ipe, ioff0
integer(i_kind),dimension(2) :: nn,nobst, nobsps, nobsq, nobsuv, nobsgps, &
nobstcp,nobstcx,nobstcy,nobstcz,nobssst, nobsspd, nobsdw, nobsrw, nobspw, &
nobsdbz
nobsdbz, nobsfed
character(8),allocatable,dimension(:):: cdiagbuf
real(r_single),allocatable,dimension(:,:)::rdiagbuf
real(r_kind) :: errorlimit,errorlimit2,error,pres,obmax
Expand All @@ -104,6 +104,7 @@ subroutine get_num_convobs_bin(obspath,datestring,num_obs_tot,num_obs_totdiag,id
nobspw = 0
nobsgps = 0
nobsdbz = 0
nobsfed = 0
nobstcp = 0; nobstcx = 0; nobstcy = 0; nobstcz = 0
init_pass = .true.
peloop: do ipe=0,npefiles
Expand Down Expand Up @@ -187,6 +188,9 @@ subroutine get_num_convobs_bin(obspath,datestring,num_obs_tot,num_obs_totdiag,id
else if (obtype == 'dbz') then
nobsdbz = nobsdbz + nn
num_obs_tot = num_obs_tot + nn(2)
else if (obtype == 'fed') then
nobsfed = nobsfed + nn
num_obs_tot = num_obs_tot + nn(2)
else if (obtype == 'gps') then
nobsgps = nobsgps + nn
num_obs_tot = num_obs_tot + nn(2)
Expand Down Expand Up @@ -231,6 +235,7 @@ subroutine get_num_convobs_bin(obspath,datestring,num_obs_tot,num_obs_totdiag,id
write(6,100) 'dw',nobsdw(1),nobsdw(2)
write(6,100) 'rw',nobsrw(1),nobsrw(2)
write(6,100) 'dbz',nobsdbz(1),nobsdbz(2)
write(6,100) 'fed',nobsfed(1),nobsfed(2)
write(6,100) 'tcp',nobstcp(1),nobstcp(2)
if (nobstcx(2) .gt. 0) then
write(6,100) 'tcx',nobstcx(1),nobstcx(2)
Expand Down Expand Up @@ -1075,6 +1080,7 @@ subroutine get_convobs_data_bin(obspath, datestring, nobs_max, nobs_maxdiag, &
if (obtype == ' t' .or. obtype == ' uv' .or. obtype == ' ps' .or. &
obtype == 'tcp' .or. obtype == ' q' .or. obtype == 'spd' .or. &
obtype == 'sst' .or. obtype == ' rw' .or. obtype == 'dbz' .or. &
obtype == 'fed' .or. &
obtype == 'gps' .or. obtype == ' dw' .or. obtype == ' pw') then

! direct reflectivitiy DA has a different routine for dbz obs.
Expand Down