Skip to content

Saving the original reported observation information if they are changed by terrain-match adjustment in subroutine setupt #868

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

Merged
ShunLiu-NOAA merged 3 commits into from
Apr 28, 2025
Merged

Saving the original reported observation information if they are changed by terrain-match adjustment in subroutine setupt #868

Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
20935a5
Minor modification in setupt.f90:
GangZhao-NOAA Apr 24, 2025
667600e
Merge remote-tracking branch 'upstream/develop' into feature/saving_o…
GangZhao-NOAA Apr 24, 2025
6153cb9
Saving the marker index in the binary format obs-diag file of T obs t…
GangZhao-NOAA Apr 25, 2025
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
74 changes: 74 additions & 0 deletions src/gsi/setupt.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -286,6 +286,7 @@ subroutine setupt(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav
real(r_kind),dimension(34) :: ptablt
real(r_single),allocatable,dimension(:,:)::rdiagbuf
real(r_single),allocatable,dimension(:,:)::rdiagbufp
real(r_kind), allocatable,dimension(:,:)::data_kept


real(r_kind),dimension(nsig):: prsltmp2
Expand All @@ -298,6 +299,7 @@ subroutine setupt(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav
integer(i_kind) regime
integer(i_kind) idomsfc,iskint,iff10,isfcr
integer(i_kind) ibb,ikk,idddd
integer(i_kind) i_dat

integer(i_kind),dimension(nobs):: buddyuse

Expand Down Expand Up @@ -375,7 +377,26 @@ subroutine setupt(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav

! call GSD terrain match for surface temperature observation
if(l_gsd_terrain_match_surftobs) then
if ( l_rtma3d ) then !keeping original reported info changed in gsd_terrain_match_surfTobs
if ( allocated(data_kept) ) deallocate(data_kept)
allocate(data_kept(5,nobs))
data_kept(:,:) = -9999.0_r_kind
data_kept(1,:) = data( 5,:) !tob : temperature obsevation (tv or tsen) (itob=5)
data_kept(2,:) = data( 4,:) !dlnpob : log of reported pressure (cb) (ipres=4)
data_kept(3,:) = data(19,:) !stnelev : station elevation (istnelv=19)
data_kept(4,:) = data( 1,:) !toe : obs error (ier=1)
data_kept(5,:) = 0.0_r_kind !marker : =0, if obs info is not changed in gsd_terrain_match_surfTobs
end if
call gsd_terrain_match_surfTobs(mype,nele,nobs,data)
if ( l_rtma3d ) then !checking if original reported info is changed in gsd_terrain_match_surfTobs
do i_dat=1,nobs
if ( abs(data_kept(1,i_dat)-data( 5,i_dat)) > tiny_r_kind .or. &
abs(data_kept(2,i_dat)-data( 4,i_dat)) > tiny_r_kind .or. &
abs(data_kept(3,i_dat)-data(19,i_dat)) > tiny_r_kind ) then
data_kept(5,:) = 1.0_r_kind !marker : =1, if obs info is changed in gsd_terrain_match_surfTobs
end if
end do
end if
endif

! index information for data array (see reading routine)
Expand Down Expand Up @@ -502,6 +523,9 @@ subroutine setupt(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav
allocate(cprvstg(nobs),csprvstg(nobs)) ! provider/subprovider info
if(l_pbl_pseudo_surfobst) allocate(cprvstgp(nobs*3),csprvstgp(nobs*3)) ! provider of pseudo obs
endif
if ( l_rtma3d .and. l_gsd_terrain_match_surftobs ) then
nreal = nreal + 4 ! to save the reported obs info (tob, lnpob, stnelev and change-marker)
end if
if (save_jacobian) then
nnz = 2 ! number of non-zero elements in dH(x)/dx profile
nind = 1
Expand Down Expand Up @@ -1382,6 +1406,8 @@ subroutine setupt(obsLL,odiagLL,lunin,mype,bwork,awork,nele,nobs,is,conv_diagsav
if(l_pbl_pseudo_surfobst) deallocate(cdiagbufp,rdiagbufp)
end if

if ( allocated(data_kept) ) deallocate(data_kept)


! End of routine

Expand Down Expand Up @@ -1593,6 +1619,8 @@ end subroutine init_netcdf_diag_
subroutine contents_binary_diag_(odiag)
type(obs_diag),pointer,intent(in):: odiag

real(r_kind) :: prest_kept

cdiagbuf(ii) = station_id ! station id

rdiagbuf(1,ii) = ictype(ikx) ! observation type
Expand Down Expand Up @@ -1668,6 +1696,19 @@ subroutine contents_binary_diag_(odiag)
csprvstg(ii) = c_sprvstg ! subprovider name
endif

! saving the original reported obs info changed in gsd_terrain_match_surfTobs
if ( l_rtma3d .and. l_gsd_terrain_match_surftobs ) then
idia = idia + 1
rdiagbuf(idia,ii) = data_kept(1,i) ! original reported tob
idia = idia + 1
prest_kept = r10*exp(data_kept(2,i)) ! cb --> mb
rdiagbuf(idia,ii) = prest_kept ! original reported log of pressure (cb)
idia = idia + 1
rdiagbuf(idia,ii) = data_kept(3,i) ! original reported station elevation
idia = idia + 1
rdiagbuf(idia,ii) = data_kept(5,i) ! marker if obs info is changed (0:no; 1: changed)
end if

if (save_jacobian) then
call writearray(dhx_dx, rdiagbuf(idia+1:nreal,ii))
idia = idia + size(dhx_dx)
Expand Down Expand Up @@ -1746,6 +1787,20 @@ subroutine contents_binary_diagp_(odiag)
! r_sprvstg = data(isprvd,i)
csprvstgp(iip) = '88888888' !c_sprvstg ! subprovider name
endif

! saving the original reported obs info changed in gsd_terrain_match_surfTobs
! (pseudo obs does not call gsd_terrain_match_surfTobs, but because same "nreal"
! is used for both true obs and pseudo obs, save the "fake" info for pseudo obs)
if ( l_rtma3d .and. l_gsd_terrain_match_surftobs ) then
idia = idia + 1
rdiagbufp(idia,iip) = -9999._r_single ! original reported tob
idia = idia + 1
rdiagbufp(idia,iip) = -9999._r_single ! original reported log of pressure (cb)
idia = idia + 1
rdiagbufp(idia,iip) = -9999._r_single ! original reported station elevation
idia = idia + 1
rdiagbufp(idia,iip) = -9999._r_single ! marker whether obs info is changed or not
end if
!----

if (save_jacobian) then
Expand All @@ -1762,6 +1817,7 @@ subroutine contents_netcdf_diag_(odiag)
real(r_single),parameter:: missing = -9.99e9_r_single

real(r_kind),dimension(miter) :: obsdiag_iuse
real(r_kind) :: prest_kept

call nc_diag_metadata("Station_ID", station_id )
call nc_diag_metadata("Observation_Class", obsclass )
Expand Down Expand Up @@ -1841,6 +1897,15 @@ subroutine contents_netcdf_diag_(odiag)
call nc_diag_metadata("Subprovider_Name", c_sprvstg )
endif

! saving the original reported obs info changed in gsd_terrain_match_surfTobs
if ( l_rtma3d .and. l_gsd_terrain_match_surftobs ) then
call nc_diag_metadata_to_single("Observation_reported", data_kept(1,i) )
prest_kept = r10*exp(data_kept(2,i)) ! cb --> mb
call nc_diag_metadata_to_single("Pressure_reported", prest_kept )
call nc_diag_metadata_to_single("Station_Elevation_reported", data_kept(3,i) )
call nc_diag_metadata_to_single("ObsInfo_UnChanged_marker", data_kept(5,i) )
end if

if (save_jacobian) then
call nc_diag_data2d("Observation_Operator_Jacobian_stind", dhx_dx%st_ind)
call nc_diag_data2d("Observation_Operator_Jacobian_endind", dhx_dx%end_ind)
Expand Down Expand Up @@ -1909,6 +1974,15 @@ subroutine contents_netcdf_diagp_(odiag)
call nc_diag_metadata("Subprovider_Name", "88888888" )
endif

! saving the original reported obs info changed in gsd_terrain_match_surfTobs
! (pseudo obs does not call gsd_terrain_match_surfTobs, but because same "nreal"
! is used for both true obs and pseudo obs, save the "fake" info for pseudo obs)
if ( l_rtma3d .and. l_gsd_terrain_match_surftobs ) then
call nc_diag_metadata("Observation_reported", missing )
call nc_diag_metadata("Pressure_reported", missing )
call nc_diag_metadata("Station_Elevation_reported", missing )
call nc_diag_metadata("ObsInfo_UnChanged_marker", missing )
end if
!----
if (save_jacobian) then
call nc_diag_data2d("Observation_Operator_Jacobian_stind", dhx_dx%st_ind)
Expand Down