forked from wrf-model/WRF
-
Notifications
You must be signed in to change notification settings - Fork 1
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
New Capability of GPM GMI Radiance Data Assimilation (wrf-model#1585)
KEYWORDS: GPM-GMI, Radiance SOURCE: Dongmei Xu (NUIST/NCAR) DESCRIPTION OF CHANGES: Add new capability to assimilate GPM-GMI Level-1 data in HDF5 format, which can be downloaded from https://gpm.nasa.gov/data/ LIST OF MODIFIED FILES: M Registry/registry.var M var/build/depend.txt 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_rad.inc M var/da/da_radiance/da_radiance.f90 M var/da/da_radiance/da_radiance1.f90 M var/da/da_radiance/da_radiance_init.inc M var/da/da_radiance/da_setup_radiance_structures.inc M var/da/da_radiance/module_radiance.f90 M var/da/da_setup_structures/da_setup_obs_structures.inc M var/da/da_setup_structures/da_setup_structures.f90 M var/run/VARBC.in A var/da/da_radiance/da_qc_gmi.inc A var/da/da_radiance/da_read_obs_hdf5gmi.inc A var/run/radiance_info/gpm-1-gmi.info TESTS CONDUCTED: WRFDA regression tests passed. RELEASE NOTE: A new capability to assimilate GPM-GMI radiance data Shen, et al.,2021: Assimilation of GPM Microwave Imager Radiance data with the WRF Hybrid 3DEnVar System for the Prediction of Typhoon Chan-hom (2015), Atmospheric Research. 251, 105422.
- Loading branch information
1 parent
1b7735a
commit 38183b5
Showing
17 changed files
with
1,071 additions
and
15 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,239 @@ | ||
subroutine da_qc_gmi (it, i, nchan, ob, iv) | ||
|
||
!--------------------------------------------------------------------------- | ||
! Purpose: perform quality control for gmi data. | ||
!--------------------------------------------------------------------------- | ||
|
||
implicit none | ||
|
||
integer, intent(in) :: it ! outer loop count | ||
integer, intent(in) :: i ! sensor index. | ||
integer, intent(in) :: nchan ! number of channel | ||
type (y_type), intent(in) :: ob ! Observation structure. | ||
type (iv_type), intent(inout) :: iv ! O-B structure. | ||
|
||
real, parameter :: clwcutofx(1:14) = & | ||
(/0.350, 0.350, 0.350, 0.350, 0.350, 0.350, 0.300, & | ||
0.300, 0.250, 0.250, 0.100, 0.100, 0.020, 0.020 /) | ||
|
||
! local variables | ||
integer :: n,k,isflg,ios,fgat_rad_unit | ||
real :: si, coscon, sincon, bearaz, sun_zenith, sun_azimuth, sun_glint | ||
integer :: ngood(nchan),nrej(nchan),nrej_omb_abs(nchan), & | ||
nrej_omb_std(nchan), & | ||
nrej_clw(nchan), num_proc_domain, & | ||
nrej_land, nrej_glint | ||
|
||
character(len=30) :: filename | ||
real :: c37_mean | ||
|
||
if (trace_use) call da_trace_entry("da_qc_gmi") | ||
|
||
ngood(:) = 0 | ||
nrej(:) = 0 | ||
nrej_omb_abs(:) = 0 | ||
nrej_omb_std(:) = 0 | ||
nrej_clw(:) = 0 | ||
nrej_land = 0 | ||
nrej_glint = 0 | ||
num_proc_domain = 0 | ||
|
||
coscon = cos( (90.0 - 55.0)*deg_to_rad) | ||
sincon = sin( (90.0 - 55.0)*deg_to_rad) | ||
|
||
do n= iv%instid(i)%info%n1,iv%instid(i)%info%n2 | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
num_proc_domain = num_proc_domain + 1 | ||
|
||
if ( crtm_cloud ) then | ||
! calculate c37_mean | ||
c37_mean = 1.0-(ob%instid(i)%tb(6,n)-ob%instid(i)%tb(7,n)+ & | ||
iv%instid(i)%tb_xb(6,n)-iv%instid(i)%tb_xb(7,n))/ & | ||
(2.0*(iv%instid(i)%tb_xb_clr(6,n)-iv%instid(i)%tb_xb_clr(7,n))) | ||
end if | ||
|
||
! 0.0 initialise QC by flags assuming good obs | ||
!----------------------------------------------------------------- | ||
iv%instid(i)%tb_qc(:,n) = qc_good | ||
|
||
! 1.0 reject all channels over land | ||
!----------------------------------------------------------------- | ||
isflg = iv%instid(i)%isflg(n) !model surface type at ob location | ||
if (only_sea_rad) then | ||
if (isflg > 0 .or. iv%instid(i)%landsea_mask(n) == 0) then | ||
iv%instid(i)%tb_qc(:,n) = qc_bad | ||
nrej_land = nrej_land + 1 | ||
end if | ||
end if | ||
|
||
|
||
|
||
! 2.0 check sun_glint angle for 10.6 GHz | ||
! It should be >= 25 degrees | ||
!----------------------------------------------------------------- | ||
bearaz = 270.0 - iv%instid(i)%satazi(n) | ||
sun_zenith = iv%instid(i)%solzen(n) | ||
sun_azimuth = 90.0 - iv%instid(i)%solazi(n) | ||
|
||
sun_glint = acos(coscon * cos( bearaz ) * cos( sun_zenith ) * cos( sun_azimuth ) + & | ||
coscon * sin( bearaz ) * cos( sun_zenith ) * sin( sun_azimuth ) + & | ||
sincon * sin( sun_zenith )) * rad_to_deg | ||
|
||
if (sun_glint < 25.0) then | ||
! apply only to 10.6 GHz for both V & H polarizations | ||
iv%instid(i)%tb_qc(1:2,n) = qc_bad | ||
nrej_glint = nrej_glint + 1 | ||
end if | ||
|
||
! 3.0 check iuse | ||
!----------------------------------------------------------------- | ||
do k = 1, nchan | ||
if (satinfo(i)%iuse(k) .eq. -1) & | ||
iv%instid(i)%tb_qc(k,n) = qc_bad | ||
end do | ||
|
||
! 4.0 check cloud | ||
!----------------------------------------------------------------- | ||
if (.not. crtm_cloud ) then | ||
do k = 1, nchan | ||
! clw check | ||
! channel dependent criteria | ||
if( iv%instid(i)%clw(n) < 0 .or. & !bad obs retrieved clw | ||
iv%instid(i)%clw(n) > clwcutofx(k) .or. & !obs retrieved clw | ||
iv%instid(i)%clwp(n) > clwcutofx(k) ) then !model/guess clwp | ||
iv%instid(i)%tb_qc(k,n) = qc_bad | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
nrej_clw(k) = nrej_clw(k) + 1 | ||
end if | ||
end do | ||
end if | ||
|
||
! assigning obs errors | ||
if (.not. crtm_cloud ) then | ||
do k = 1, nchan | ||
if (use_error_factor_rad) then | ||
iv%instid(i)%tb_error(k,n) = & | ||
satinfo(i)%error_std(k)*satinfo(i)%error_factor(k) | ||
else | ||
iv%instid(i)%tb_error(k,n) = satinfo(i)%error_std(k) | ||
end if | ||
end do ! nchan | ||
|
||
else !crtm_cloud | ||
! do k = 1, nchan | ||
! ! clw check--remove clw <0 bad quality pixel | ||
! ! channel dependent criteria | ||
! if( iv%instid(i)%clw(n) < 0 ) then !bad obs retrieved clw | ||
! iv%instid(i)%tb_qc(k,n) = qc_bad | ||
! if (iv%instid(i)%info%proc_domain(1,n)) & | ||
! nrej_clw(k) = nrej_clw(k) + 1 | ||
! end if | ||
! end do !nchan | ||
|
||
! symmetric error model, Geer and Bauer (2011) | ||
do k = 1, nchan | ||
if (c37_mean.lt.0.05) then | ||
iv%instid(i)%tb_error(k,n)= satinfo(i)%error_std(k) | ||
else if (c37_mean.ge.0.05.and.c37_mean.lt.0.5) then | ||
iv%instid(i)%tb_error(k,n)= satinfo(i)%error_std(k)+ & | ||
(c37_mean-0.05)*(satinfo(i)%error_cld(k)-satinfo(i)%error_std(k))/(0.5-0.05) | ||
else | ||
iv%instid(i)%tb_error(k,n)= satinfo(i)%error_cld(k) | ||
end if | ||
end do ! nchan | ||
|
||
end if | ||
|
||
! 5.0 check innovation | ||
!----------------------------------------------------------------- | ||
if (.not. crtm_cloud ) then | ||
! absolute departure check | ||
do k = 1, nchan | ||
if ( k <= 2 .or. k == 6 .or. k == 7) then | ||
if (abs(iv%instid(i)%tb_inv(k,n)) > 6.0) then | ||
iv%instid(i)%tb_qc(k,n) = qc_bad | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
nrej_omb_abs(k) = nrej_omb_abs(k) + 1 | ||
end if | ||
else if ( k == 4 .or. k == 5 ) then | ||
if (abs(iv%instid(i)%tb_inv(k,n)) > 8.0) then | ||
iv%instid(i)%tb_qc(k,n) = qc_bad | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
nrej_omb_abs(k) = nrej_omb_abs(k) + 1 | ||
end if | ||
else | ||
if (abs(iv%instid(i)%tb_inv(k,n)) > 10.0) then | ||
iv%instid(i)%tb_qc(k,n) = qc_bad | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
nrej_omb_abs(k) = nrej_omb_abs(k) + 1 | ||
end if | ||
end if | ||
end do ! nchan | ||
end if | ||
|
||
do k = 1, nchan | ||
! relative departure check | ||
if (abs(iv%instid(i)%tb_inv(k,n)) > 3.0*iv%instid(i)%tb_error(k,n)) then | ||
iv%instid(i)%tb_qc(k,n) = qc_bad | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
nrej_omb_std(k) = nrej_omb_std(k) + 1 | ||
end if | ||
|
||
! final QC decsion | ||
if (iv%instid(i)%tb_qc(k,n) == qc_bad) then | ||
iv%instid(i)%tb_error(k,n) = 500.0 | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
nrej(k) = nrej(k) + 1 | ||
else | ||
if (iv%instid(i)%info%proc_domain(1,n)) & | ||
ngood(k) = ngood(k) + 1 | ||
end if | ||
end do ! nchan | ||
|
||
end do ! end loop pixel | ||
|
||
! Do inter-processor communication to gather statistics. | ||
call da_proc_sum_int (num_proc_domain) | ||
call da_proc_sum_ints (nrej_clw) | ||
call da_proc_sum_int (nrej_land) | ||
call da_proc_sum_int (nrej_glint) | ||
call da_proc_sum_ints (nrej_omb_abs) | ||
call da_proc_sum_ints (nrej_omb_std) | ||
call da_proc_sum_ints (nrej) | ||
call da_proc_sum_ints (ngood) | ||
|
||
if (rootproc) then | ||
if (num_fgat_time > 1) then | ||
write(filename,'(i2.2,a,i2.2)') it,'_qcstat_'//trim(iv%instid(i)%rttovid_string)//'_',iv%time | ||
else | ||
write(filename,'(i2.2,a)') it,'_qcstat_'//trim(iv%instid(i)%rttovid_string) | ||
end if | ||
|
||
call da_get_unit(fgat_rad_unit) | ||
open(fgat_rad_unit,file=trim(filename),form='formatted',iostat=ios) | ||
if (ios /= 0) then | ||
write(unit=message(1),fmt='(A,A)') 'error opening the output file ', filename | ||
call da_error(__FILE__,__LINE__,message(1:1)) | ||
end if | ||
|
||
write(fgat_rad_unit, fmt='(/a/)') ' Quality Control Statistics for '//iv%instid(i)%rttovid_string | ||
if(num_proc_domain > 0) write(fgat_rad_unit,'(a20,i7)') ' num_proc_domain = ', num_proc_domain | ||
if(nrej_land > 0)write(fgat_rad_unit,'(a20,i7)') ' nrej_land = ', nrej_land | ||
if(nrej_glint> 0)write(fgat_rad_unit,'(a20,i7)') ' nrej_glint = ', nrej_glint | ||
write(fgat_rad_unit,'(a20)') ' nrej_omb_abs(:) = ' | ||
write(fgat_rad_unit,'(10i7)') nrej_omb_abs(:) | ||
write(fgat_rad_unit,'(a20)') ' nrej_omb_std(:) = ' | ||
write(fgat_rad_unit,'(10i7)') nrej_omb_std(:) | ||
write(fgat_rad_unit,'(a20)') ' nrej_clw(:) = ' | ||
write(fgat_rad_unit,'(10i7)') nrej_clw(:) | ||
write(fgat_rad_unit,'(a20)') ' nrej(:) = ' | ||
write(fgat_rad_unit,'(10i7)') nrej(:) | ||
write(fgat_rad_unit,'(a20)') ' ngood(:) = ' | ||
write(fgat_rad_unit,'(10i7)') ngood(:) | ||
|
||
close(fgat_rad_unit) | ||
call da_free_unit(fgat_rad_unit) | ||
end if | ||
if (trace_use) call da_trace_exit("da_qc_gmi") | ||
|
||
end subroutine da_qc_gmi |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Oops, something went wrong.