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

Jump to bottom

Bugfix for where some fields are not in the input increment file #53

Merged
merged 8 commits into from
Sep 17, 2024
Merged

Bugfix for where some fields are not in the input increment file #53

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
590f084
skip possibly missing hydrometeors
CoryMartin-NOAA Sep 10, 2024
723c211
add warning print
CoryMartin-NOAA Sep 10, 2024
9efe71d
Safer code
CoryMartin-NOAA Sep 10, 2024
5c8db46
handle case of missing increment in calc_analysis, add omga to calc_a…
RussTreadon-NOAA Sep 12, 2024
4eda78d
remove omga from calc_analysis
RussTreadon-NOAA Sep 12, 2024
20c6cab
reduce nnciov by one after removing omga
RussTreadon-NOAA Sep 12, 2024
1c512a1
correct erroneous calc_analysis runtime printout
RussTreadon-NOAA Sep 13, 2024
2f38004
add omga to analysis file generated by calc_analysis
RussTreadon-NOAA Sep 17, 2024
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
43 changes: 28 additions & 15 deletions src/netcdf_io/calc_analysis.fd/inc2anl.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@ module inc2anl

integer, parameter :: nincv=13
character(len=7) :: incvars_nemsio(nincv), incvars_netcdf(nincv), incvars_ncio(nincv)
integer, parameter :: nnciov=22
integer, parameter :: nnciov=23
integer, parameter :: naero=14
integer, parameter :: naero_copy=6
character(len=7) :: iovars_netcdf(nnciov), iovars_aero(naero), copyvars_aero(naero_copy)
Expand All @@ -36,7 +36,7 @@ module inc2anl
'delz ', 'dpres ', 'dzdt ', 'grle ', 'hgtsfc ',&
'icmr ', 'o3mr ', 'pressfc', 'rwmr ', 'snmr ',&
'spfh ', 'tmp ', 'ugrd ', 'vgrd ', 'cld_amt',&
'nccice ', 'nconrd '/
'nccice ', 'nconrd ', 'omga '/
data iovars_aero / 'so4 ', 'bc1 ', 'bc2 ', 'oc1 ', 'oc2 ', &
'dust1 ', 'dust2 ', 'dust3 ', 'dust4 ', 'dust5 ',&
'seas1 ', 'seas2 ', 'seas3 ', 'seas4 '/
Expand All @@ -60,11 +60,11 @@ subroutine gen_anl
! increment, add the two together, and write out
! the analysis to a new file
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
use vars_calc_analysis, only: mype, do_aero
use vars_calc_analysis, only: mype, do_aero, jedi
implicit none
! variables local to this subroutine
integer :: i, j, iincvar
logical :: use_increment
logical :: use_increment,readinc


! loop through each variable in the background file
Expand All @@ -85,9 +85,15 @@ subroutine gen_anl
call add_psfc_increment
else
! call generic subroutine for all other fields
if (mype==0) print *, 'Adding Increment to ', iovars_netcdf(i), incvars_netcdf(iincvar)
call add_increment(iovars_netcdf(i), incvars_netcdf(iincvar))
end if
call add_increment(iovars_netcdf(i), incvars_netcdf(iincvar), readinc)
if (mype==0) then
if ((.not.jedi) .or. (jedi.and.readinc)) then
print *, 'Adding Increment to ', iovars_netcdf(i), incvars_netcdf(iincvar)
else
print *, 'Copying from Background ', iovars_netcdf(i)
endif
endif
endif
else
! otherwise just write out what is in the input to the output
if (mype==0) print *, 'Copying from Background ', iovars_netcdf(i)
Expand Down Expand Up @@ -241,7 +247,7 @@ subroutine add_aero_inc(fcstvar, incvar)

end subroutine add_aero_inc

subroutine add_increment(fcstvar, incvar)
subroutine add_increment(fcstvar, incvar, readinc)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! subroutine add_increment
! generic subroutine for adding increment to background
Expand All @@ -250,6 +256,7 @@ subroutine add_increment(fcstvar, incvar)
! fcstvar - input string of netCDF fcst/anal var name
! incvar - input string of netCDF increment var name
! (without _inc suffix added)
! readinc - .true. if read increment, .false otherwise
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
use vars_calc_analysis, only: fcstncfile, anlncfile, incr_file,&
nlat, nlon, nlev, anlfile, use_nemsio_anl, &
Expand All @@ -260,6 +267,7 @@ subroutine add_increment(fcstvar, incvar)
implicit none
! input vars
character(7), intent(in) :: fcstvar, incvar
logical, intent(out):: readinc
! local variables
real, allocatable, dimension(:,:,:) :: work3d_bg
real, allocatable, dimension(:,:) :: work3d_inc_gsi
Expand All @@ -268,6 +276,7 @@ subroutine add_increment(fcstvar, incvar)
integer :: j,jj,k,krev,iret
type(Dataset) :: incncfile

readinc=.true.
if (has_var(fcstncfile, fcstvar)) then
do k=1,nlev
if (mype == levpe(k)) then
Expand All @@ -277,12 +286,16 @@ subroutine add_increment(fcstvar, incvar)
incncfile = open_dataset(incr_file, paropen=.true.)
! JEDI-derived increments have a time dimension, so read to appropriate array
if ( jedi ) then
call read_vardata(incncfile, trim(incvar)//"_inc", work3d_inc_jedi, nslice=k, slicedim=3)
! add increment to background
do j=1,nlat
jj=nlat+1-j ! increment is S->N, history files are N->S
work3d_bg(:,j,1) = work3d_bg(:,j,1) + work3d_inc_jedi(:,jj,1)
end do
if (has_var(incncfile, trim(incvar)//"_inc")) then
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should this check if specific variables are missing, similar to how interp_inc now allows for only the hydrometeors to be missing?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One could limit this check to the specific hydrometeors flagged in the revised driver.F90. A more general approach was taken here to exclude all variables for which no increment is found. What approach shall we take? Thoughts @CatherineThomas-NOAA , @emilyhcliu , @CoryMartin-NOAA

Copy link
Contributor

@emilyhcliu emilyhcliu Sep 12, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The best way is probably to relate the required increments from the state variable table in anavinfo. This requires a bit more coding and related script changes. For now, the proposed change should work and be acceptable.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you @emilyhcliu for your suggestion. I agree with you. We can stick with the easier solution for the time being.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Logical variable readinc added to improve runtime printout. When a JEDI increment file is processed and a specified increment field is not present readinc is set to .false.. The calling routine uses the jedi logical along with readinc to select the appropriate runtime message. If the increment is added to the background the runtime printout is, for example,

nid002740.dogwood.wcoss2.ncep.noaa.gov 0:  Adding Increment to spfh   sphum
nid002740.dogwood.wcoss2.ncep.noaa.gov 0:  Adding Increment to tmp    T
nid002740.dogwood.wcoss2.ncep.noaa.gov 0:  Adding Increment to ugrd   u

Otherwise the printout is

nid002740.dogwood.wcoss2.ncep.noaa.gov 0:  Copying from Background grle
nid002740.dogwood.wcoss2.ncep.noaa.gov 0:  Copying from Background rwmr
nid002740.dogwood.wcoss2.ncep.noaa.gov 0:  Copying from Background snmr

Previously, even though JEDI increments for grle, rwmr, and snmr were not available, the runtime printout was

nid001756.dogwood.wcoss2.ncep.noaa.gov 0:  Adding Increment to grle   grle
nid001756.dogwood.wcoss2.ncep.noaa.gov 0:  Adding Increment to rwmr   rwmr
nid001756.dogwood.wcoss2.ncep.noaa.gov 0:  Adding Increment to snmr   snmr

This change committed at 1c512a1.

call read_vardata(incncfile, trim(incvar)//"_inc", work3d_inc_jedi, nslice=k, slicedim=3)
! add increment to background, otherwise do nothing
do j=1,nlat
jj=nlat+1-j ! increment is S->N, history files are N->S
work3d_bg(:,j,1) = work3d_bg(:,j,1) + work3d_inc_jedi(:,jj,1)
end do
else
readinc = .false.
endif
else
call read_vardata(incncfile, trim(incvar)//"_inc", work3d_inc_gsi, nslice=k, slicedim=3)
RussTreadon-NOAA marked this conversation as resolved.
Show resolved Hide resolved
! add increment to background
Expand All @@ -306,7 +319,7 @@ subroutine add_increment(fcstvar, incvar)
end do
! clean up and close
if ( jedi ) then
deallocate(work3d_inc_jedi)
if (allocated(work3d_inc_jedi)) deallocate(work3d_inc_jedi)
else
deallocate(work3d_inc_gsi)
end if
Expand Down
20 changes: 17 additions & 3 deletions src/netcdf_io/interp_inc.fd/driver.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -81,6 +81,7 @@ program interp_inc
real(8), allocatable :: gi(:,:), gi2(:,:), go(:,:), go2(:,:), go3(:,:)
real(8), allocatable :: send_layer(:), recv_layer(:)

logical :: readvar

! NOTE: u_inc,v_inc must be consecutive
data records /'u_inc', 'v_inc', 'delp_inc', 'delz_inc', 'T_inc', &
Expand Down Expand Up @@ -366,11 +367,24 @@ program interp_inc

if (mype == rec-1) then
print*,'- PROCESS RECORD: ', trim(records(rec))
readvar = .true.

error = nf90_inq_varid(ncid_in, trim(records(rec)), id_var)
call netcdf_err(error, 'inquiring ' // trim(records(rec)) // ' id for file='//trim(infile) )
error = nf90_get_var(ncid_in, id_var, dummy_in)
call netcdf_err(error, 'reading ' // trim(records(rec)) // ' for file='//trim(infile) )
! handle missing hydrometeor increments
if (error .ne. 0) then
if (ANY((/ 'rwmr_inc', 'snmr_inc', 'grle_inc' /) == trim(records(rec)))) then
print *, 'WARNING: ', trim(records(rec)), ' is missing in increment file. Skipping.'
readvar = .false.
else
call netcdf_err(error, 'inquiring ' // trim(records(rec)) // ' id for file='//trim(infile) )
end if
end if
if (readvar) then
error = nf90_get_var(ncid_in, id_var, dummy_in)
call netcdf_err(error, 'reading ' // trim(records(rec)) // ' for file='//trim(infile) )
else
dummy_in(:,:,:) = 0.0
end if

ip = 0 ! bilinear
ipopt = 0
Expand Down
Loading