Skip to content

Commit

Permalink
Add superobbing capability for AHI radiance (#1173)
Browse files Browse the repository at this point in the history
TYPE: new feature

KEYWORDS: superobbing, cloud detection, cloud fraction, AHI 

SOURCE: Yali Wu (NCAR/MMM), Jonathan Guerrette (NCAR/MMM).

DESCRIPTION OF CHANGES: This PR adds a new function of superobbing AHI radiances and retrieving cloud fraction based on AHI cloud detection results. The default namelist parameter `ahi_superob_halfwidth` is set to zero, which means no superobbing will be applied.  

ISSUE: no issue

LIST OF MODIFIED FILES: 
M       Registry/registry.var
M       var/da/da_define_structures/da_define_structures.f90
M       var/da/da_monitor/da_rad_diags.f90
M       var/da/da_radiance/da_allocate_rad_iv.inc
M       var/da/da_radiance/da_deallocate_radiance.inc
M       var/da/da_radiance/da_initialize_rad_iv.inc
M       var/da/da_radiance/da_qc_ahi.inc
M       var/da/da_radiance/da_radiance.f90
M       var/da/da_radiance/da_radiance1.f90
M       var/da/da_radiance/da_read_obs_netcdf4ahi_jaxa.inc
M       var/da/da_radiance/da_write_iv_rad_ascii.inc
M       var/da/da_radiance/da_write_oa_rad_ascii.inc

TESTS CONDUCTED:
This capability has been tested in generating one-month AHI obs files for MPAS-JEDI DA purpose. WRFDA regression tests passed.
  • Loading branch information
YaliWu0219 authored May 27, 2020
1 parent 0f61c1f commit ace4756
Show file tree
Hide file tree
Showing 12 changed files with 265 additions and 99 deletions.
1 change: 1 addition & 0 deletions Registry/registry.var
Original file line number Diff line number Diff line change
Expand Up @@ -452,6 +452,7 @@ rconfig integer varbc_nbgerr namelist,wrfvar14 1 5000 - "va
rconfig integer varbc_nobsmin namelist,wrfvar14 1 10 - "varbc_nobsmin" "" ""
rconfig integer use_clddet namelist,wrfvar14 1 2 - "use_clddet" "0: off, 1: mmr, 2: pf, 3: ecmwf" ""
rconfig logical use_clddet_zz namelist,wrfvar14 1 .false. - "use_clddet_zz" "cloud detection scheme from Zhuge X. and Zou X. JAMC, 2016." ""
rconfig integer ahi_superob_halfwidth namelist,wrfvar14 1 0 - "ahi_superob_halfwidth" "" ""
rconfig logical airs_warmest_fov namelist,wrfvar14 1 .false. - "airs_warmest_fov" "" ""
rconfig logical use_satcv namelist,wrfvar14 2 .false. - "use_satcv" "" ""
rconfig logical use_blacklist_rad namelist,wrfvar14 1 .true. - "use_blacklist_rad" "" ""
Expand Down
8 changes: 7 additions & 1 deletion var/da/da_define_structures/da_define_structures.f90
Original file line number Diff line number Diff line change
Expand Up @@ -514,6 +514,10 @@ module da_define_structures
real :: CIRH2O
!real, allocatable :: CIRH2O(:,:,:)
end type clddet_geoir_type
type superob_type
real, allocatable :: tb_obs(:,:)
type(clddet_geoir_type), allocatable :: cld_qc(:)
end type superob_type
type cv_index_type
integer :: ts
integer :: nclouds
Expand Down Expand Up @@ -549,6 +553,7 @@ module da_define_structures
integer, pointer :: cloud_flag(:,:)
integer, pointer :: cloudflag(:)
integer, pointer :: rain_flag(:)
real, allocatable :: cloud_frac(:)
real, pointer :: satzen(:)
real, pointer :: satazi(:)
real, pointer :: solzen(:)
Expand Down Expand Up @@ -629,11 +634,12 @@ module da_define_structures
real, pointer :: ice_coverage(:)
real, pointer :: snow_coverage(:)
integer, pointer :: crtm_climat(:) ! CRTM only
type (clddet_geoir_type), pointer :: cld_qc(:)
integer :: superob_width = 1
type (varbc_info_type) :: varbc_info
type (varbc_type),pointer :: varbc(:)
type (cv_index_type), pointer :: cv_index(:)
type (infa_type) :: info
type (superob_type), allocatable :: superob(:,:)
end type instid_type

type iv_type
Expand Down
28 changes: 22 additions & 6 deletions var/da/da_monitor/da_rad_diags.f90
Original file line number Diff line number Diff line change
Expand Up @@ -58,8 +58,8 @@ program da_rad_diags
integer, dimension(:), allocatable :: landsea_mask, soiltyp, vegtyp
real*4, dimension(:), allocatable :: lat, lon, elv, elev
real*4, dimension(:), allocatable :: ret_clw
real*4, dimension(:), allocatable :: satzen, satazi, t2m, mr2m, u10, v10, ps, ts
real*4, dimension(:), allocatable :: smois, tslb, snowh, vegfra, clwp
real*4, dimension(:), allocatable :: satzen, satazi, solzen, solazi, t2m, mr2m, u10, v10, ps, ts
real*4, dimension(:), allocatable :: smois, tslb, snowh, vegfra, clwp, cloud_frac
integer, dimension(:,:), allocatable :: tb_qc
real*4, dimension(:,:), allocatable :: tb_obs, tb_bak, tb_inv, tb_oma, tb_err, ems, ems_jac
real*4, dimension(:,:), allocatable :: prf_pfull, prf_phalf, prf_t, prf_q, prf_water
Expand Down Expand Up @@ -230,6 +230,8 @@ program da_rad_diags
allocate ( lon(1:total_npixel) )
allocate ( satzen(1:total_npixel) )
allocate ( satazi(1:total_npixel) )
allocate ( solzen(1:total_npixel) )
allocate ( solazi(1:total_npixel) )
allocate ( ret_clw(1:total_npixel) ) !obs retrieved clw
allocate ( t2m(1:total_npixel) )
allocate ( mr2m(1:total_npixel) )
Expand All @@ -246,6 +248,7 @@ program da_rad_diags
allocate ( vegfra(1:total_npixel) )
allocate ( elev(1:total_npixel) )
allocate ( clwp(1:total_npixel) ) !model/guess clwp
allocate ( cloud_frac(1:total_npixel) )
allocate ( tb_obs(1:nchan,1:total_npixel) )
allocate ( tb_bak(1:nchan,1:total_npixel) )
allocate ( tb_inv(1:nchan,1:total_npixel) )
Expand Down Expand Up @@ -341,14 +344,15 @@ program da_rad_diags

npixel_loop: do ipixel = ips, ipe

read(unit=iunit(iproc),fmt='(7x,i7,2x,a19,i6,i3,f6.0,4f8.2,f8.3)',iostat=ios) &
read(unit=iunit(iproc),fmt='(7x,i7,2x,a19,i6,i3,f6.0,6f8.2,f8.3)',iostat=ios) &
n, datestr2(ipixel), scanpos(ipixel), landsea_mask(ipixel), elv(ipixel), &
lat(ipixel), lon(ipixel), satzen(ipixel), satazi(ipixel), ret_clw(ipixel)
lat(ipixel), lon(ipixel), satzen(ipixel), satazi(ipixel), solzen(ipixel), &
solazi(ipixel), ret_clw(ipixel)

read(unit=iunit(iproc),fmt='(14x,9f10.2,3i3,3f10.2)',iostat=ios) &
read(unit=iunit(iproc),fmt='(14x,9f10.2,3i3,3f10.2,f15.5)',iostat=ios) &
t2m(ipixel), mr2m(ipixel), u10(ipixel), v10(ipixel), ps(ipixel), ts(ipixel), &
smois(ipixel), tslb(ipixel), snowh(ipixel), isflg(ipixel), soiltyp(ipixel), &
vegtyp(ipixel), vegfra(ipixel), elev(ipixel), clwp(ipixel)
vegtyp(ipixel), vegfra(ipixel), elev(ipixel), clwp(ipixel), cloud_frac(ipixel)
read(unit=iunit(iproc),fmt='(a)',iostat=ios) buf ! OBS
read(unit=iunit(iproc),fmt='(10f11.2)',iostat=ios) tb_obs(:,ipixel)
read(unit=iunit(iproc),fmt='(a)',iostat=ios) buf ! BAK
Expand Down Expand Up @@ -580,6 +584,8 @@ program da_rad_diags
ios = NF_DEF_VAR(ncid, 'lon', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 'satzen', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 'satazi', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 'solzen', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 'solazi', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 't2m', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 'mr2m', NF_FLOAT, 1, ishape(1), varid)
ios = NF_DEF_VAR(ncid, 'u10', NF_FLOAT, 1, ishape(1), varid)
Expand All @@ -598,6 +604,7 @@ program da_rad_diags
if ( amsr2 ) then
ios = NF_DEF_VAR(ncid, 'ret_clw', NF_FLOAT, 1, ishape(1), varid)
end if
ios = NF_DEF_VAR(ncid, 'cloud_frac', NF_FLOAT, 1, ishape(1), varid)

ios = NF_ENDDEF(ncid)
!
Expand Down Expand Up @@ -768,6 +775,10 @@ program da_rad_diags
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), satzen)
ios = NF_INQ_VARID (ncid, 'satazi', varid)
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), satazi)
ios = NF_INQ_VARID (ncid, 'solzen', varid)
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), solzen)
ios = NF_INQ_VARID (ncid, 'solazi', varid)
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), solazi)
ios = NF_INQ_VARID (ncid, 't2m', varid)
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), t2m)
ios = NF_INQ_VARID (ncid, 'mr2m', varid)
Expand Down Expand Up @@ -802,6 +813,8 @@ program da_rad_diags
ios = NF_INQ_VARID (ncid, 'ret_clw', varid)
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), ret_clw)
end if
ios = NF_INQ_VARID (ncid, 'cloud_frac', varid)
ios = NF_PUT_VARA_REAL(ncid, varid, istart(2), icount(2), cloud_frac)
!
ios = NF_CLOSE(ncid)

Expand All @@ -818,6 +831,8 @@ program da_rad_diags
deallocate ( lon )
deallocate ( satzen )
deallocate ( satazi )
deallocate ( solzen )
deallocate ( solazi )
deallocate ( t2m )
deallocate ( mr2m )
deallocate ( u10 )
Expand All @@ -833,6 +848,7 @@ program da_rad_diags
deallocate ( vegfra )
deallocate ( elev )
deallocate ( clwp )
deallocate ( cloud_frac )
deallocate ( ret_clw )
deallocate ( tb_obs )
deallocate ( tb_bak )
Expand Down
13 changes: 11 additions & 2 deletions var/da/da_radiance/da_allocate_rad_iv.inc
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ subroutine da_allocate_rad_iv (i, nchan, iv)

integer , intent (in) :: i, nchan
type (iv_type) , intent (inout) :: iv
integer :: n
integer :: n, ix, iy
call da_trace_entry("da_allocate_rad_iv")

allocate (iv%instid(i)%info%date_char(iv%instid(i)%num_rad))
Expand Down Expand Up @@ -106,8 +106,17 @@ subroutine da_allocate_rad_iv (i, nchan, iv)
allocate (iv%instid(i)%solazi(iv%instid(i)%num_rad))
allocate (iv%instid(i)%tropt(iv%instid(i)%num_rad))
allocate (iv%instid(i)%gamma_jacobian(nchan,iv%instid(i)%num_rad))
allocate (iv%instid(i)%cloud_frac(iv%instid(i)%num_rad))
if ( use_clddet_zz ) then
allocate (iv%instid(i)%cld_qc(iv%instid(i)%num_rad))
iv%instid(i)%superob_width = 2*ahi_superob_halfwidth+1
allocate (iv%instid(i)%superob(iv%instid(i)%superob_width, &
iv%instid(i)%superob_width))
do iy = 1, iv%instid(i)%superob_width
do ix = 1, iv%instid(i)%superob_width
allocate (iv%instid(i)%superob(ix,iy)%cld_qc(iv%instid(i)%num_rad))
allocate (iv%instid(i)%superob(ix,iy)%tb_obs(nchan,iv%instid(i)%num_rad))
end do
end do
end if
if ( use_rttov_kmatrix .or. use_crtm_kmatrix ) then
allocate(iv%instid(i)%ts_jacobian(nchan,iv%instid(i)%num_rad))
Expand Down
13 changes: 10 additions & 3 deletions var/da/da_radiance/da_deallocate_radiance.inc
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@
type (iv_type), intent(inout) :: iv ! Obs. increment structure.
type (j_type), intent(inout) :: j ! Cost function.

integer :: i,n,ichan
integer :: i,n,ichan,ix,iy

if (trace_use) call da_trace_entry("da_deallocate_radiance")

Expand Down Expand Up @@ -127,10 +127,17 @@
deallocate (iv%instid(i)%satazi)
deallocate (iv%instid(i)%solzen)
deallocate (iv%instid(i)%solazi)
deallocate (iv%instid(i)%tropt)
deallocate (iv%instid(i)%tropt)
deallocate (iv%instid(i)%gamma_jacobian)
deallocate(iv%instid(i)%cloud_frac)
if ( use_clddet_zz ) then
deallocate (iv%instid(i)%cld_qc)
do iy = 1, iv%instid(i)%superob_width
do ix = 1, iv%instid(i)%superob_width
deallocate (iv%instid(i)%superob(ix,iy)%cld_qc)
deallocate (iv%instid(i)%superob(ix,iy)%tb_obs)
end do
end do
deallocate (iv%instid(i)%superob)
end if
if (ANY(use_satcv)) then
if (use_satcv(2)) then
Expand Down
25 changes: 17 additions & 8 deletions var/da/da_radiance/da_initialize_rad_iv.inc
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@ subroutine da_initialize_rad_iv (i, n, iv, p)
integer, intent(in) :: i, n
type(datalink_type), intent(in) :: p
type(iv_type), intent(inout) :: iv
integer :: ix, iy

call da_trace_entry("da_initialize_rad_iv")

Expand Down Expand Up @@ -101,16 +102,24 @@ subroutine da_initialize_rad_iv (i, n, iv, p)
iv%instid(i)%solzen(n) = p%solzen
iv%instid(i)%solazi(n) = p%solazi
iv%instid(i)%tropt(n) = 0.0
iv%instid(i)%cloud_frac(n) = missing_r
! iv%instid(i)%solazi(n) = 0.0
if ( use_clddet_zz ) then
iv%instid(i)%cld_qc(n)%RTCT = p % cld_qc % RTCT
iv%instid(i)%cld_qc(n)%RFMFT = p % cld_qc % RFMFT
iv%instid(i)%cld_qc(n)%TEMPIR = p % cld_qc % TEMPIR
iv%instid(i)%cld_qc(n)%tb_stddev_10 = p % cld_qc % tb_stddev_10
iv%instid(i)%cld_qc(n)%tb_stddev_13 = p % cld_qc % tb_stddev_13
iv%instid(i)%cld_qc(n)%tb_stddev_14 = p % cld_qc % tb_stddev_14
iv%instid(i)%cld_qc(n)%terr_hgt = p % cld_qc % terr_hgt
iv%instid(i)%cld_qc(n)%CIRH2O = p % cld_qc % CIRH2O
if ( allocated ( p % superob ) ) then
do iy = 1, iv%instid(i)%superob_width
do ix = 1, iv%instid(i)%superob_width
iv%instid(i)%superob(ix,iy)%tb_obs(:,n) = p % superob(ix,iy) % tb_obs(:,1)
iv%instid(i)%superob(ix,iy)%cld_qc(n)%RTCT = p % superob(ix,iy) % cld_qc(1) % RTCT
iv%instid(i)%superob(ix,iy)%cld_qc(n)%RFMFT = p % superob(ix,iy) % cld_qc(1) % RFMFT
iv%instid(i)%superob(ix,iy)%cld_qc(n)%TEMPIR = p % superob(ix,iy) % cld_qc(1) % TEMPIR
iv%instid(i)%superob(ix,iy)%cld_qc(n)%tb_stddev_10 = p % superob(ix,iy) % cld_qc(1) % tb_stddev_10
iv%instid(i)%superob(ix,iy)%cld_qc(n)%tb_stddev_13 = p % superob(ix,iy) % cld_qc(1) % tb_stddev_13
iv%instid(i)%superob(ix,iy)%cld_qc(n)%tb_stddev_14 = p % superob(ix,iy) % cld_qc(1) % tb_stddev_14
iv%instid(i)%superob(ix,iy)%cld_qc(n)%terr_hgt = p % superob(ix,iy) % cld_qc(1) % terr_hgt
iv%instid(i)%superob(ix,iy)%cld_qc(n)%CIRH2O = p % superob(ix,iy) % cld_qc(1) % CIRH2O
end do
end do
end if
end if
if ( rtm_option == rtm_option_rttov ) then
iv%instid(i)%surftype(n) = 0
Expand Down
Loading

0 comments on commit ace4756

Please sign in to comment.