Resolve "nreal" differences for conventional data

[ilianagenkova]

See email: "data(iamvg,i) out of bounds for some uv obs types" (Russ Treadon, 25 Aug 2022)

R.T.
Part of NCO's IT testing is to build and run apps with "-check all". I've done so for release/gfsda.v16.3.0. gsi.x fails with an out of bounds array error
nid001207.cactus.wcoss2.ncep.noaa.gov 369: forrtl: severe (408): fort: (2): Subscript #1 of the array DATA has value 27 which is greater than the upper bound of 26
Image PC Routine Line Source
gsi.x 00000000058628CF Unknown Unknown Unknown
gsi.x 0000000003DAD6DE w_setupsetupw_mp_ 1817 setupw.f90
gsi.x 0000000003DA57BD w_setup_mp_setupw 1395 setupw.f90
gsi.x 0000000002B4CDED gsi_woper_mp_setu 103 gsi_wOper.F90
gsi.x 0000000001CC6DB7 setuprhsall_ 490 setuprhsall.f90

Line 1817 of setupw.f90 is the imavq line below
! AMVQ Mitigated winds
call nc_diag_metadata("Mitigation_flag_AMVQ", sngl(data(iamvq,i)) )

setupw.f90 sets iamvq=27. Not all uv observation types have 27 elements in the data array. I added prints to the code and reran the above case. The failure occurs when process FL_HDOB uv data. For FL_HDOB uv the data array is declared to be size (26,nobs) ... just like the error message indicates.

Is there an easy way in contents_netcdf_diag_ (in setupw.f90) to identify that we are processing AMV data with a valid value for data(iamvq,i)?

I.G.
I experienced the same issue when I first introduced iamvq=27.
That's because, as you noted, not all uv obs have 27 elements in the data array. I increased the array size in read_satwnd, but I also had to increase the array size in read_prepbufr (winds are read there too). Haixia and I exchanged comments on it in a GitHIb issue (I think) or may be an email.
I was not aware there are other read_* routines that route winds to setupw.f90. Is that the case with FL_HDOB? Then the cdata array has to be increased in the corresponding read routine. Where is FL_HDOB read?

R.T
Yes, the flight level high density data has its own read routine, read_fl_hdob.f90. This data can have a variety of nreal = nele values depending on what data is being processed
if (ltob) then
nreal = 25
if (aircraft_t_bc_pof .or. aircraft_t_bc .or.aircraft_t_bc_ext) nreal=nreal+3
iecol = 2
errmin = half ! set lower bound of ob error for T or Tv
else if (luvob) then
nreal = 26
iecol = 4
errmin = one ! set lower bound of ob error for u,v winds
else if (lspdob) then
nreal = 23
iecol = 4
errmin = one
else if (lqob) then ! set lower bound of ob err for surface wind speed
nreal = 26
iecol = 3
errmin = half ! set lower bound of ob error for moisture (RH)
else if (lpsob) then
nreal = 23
iecol = 5
errmin = one_tenth ! set lower bound of ob error for moisture (RH)
else
write(6,*) ' illegal obs type in read_fl_hdob '
call stop2(94)
end if
if (perturb_obs .and. luvob) nreal = nreal+2
if (perturb_obs .and. (ltob .or. lpsob .or. lqob)) nreal = nreal+1
...
! Allocate array to hold data
allocate(cdata_all(nreal,maxobs))

read_sfcwnd.f90 also processes uv data. It has nreal /= 27. Instead

! nreal keep the dimension of cdata_all
werrmin=one
nsattype=0
nreal=24
...
! Loop over convinfo file entries; operate on matches
allocate(cdata_all(nreal,maxobs),isort(maxobs))

I don't think we want to set nreal=27 for all cases.

I.G.
I see, separate read_routine will cause that.
We can have different nreal for the different obs types (t, q, ps, spd), but all uv data get to setupw. How would you rather handle the differences in nreal for various uv observations? Also worth nothing that different variable may populate the same index in the cdata array, depending on data type, b/c not all data types provide the same info.

R.T.
I doubt we can fix this in a robust way for GFS v16.3.0.
We should not assume a fixed nreal (nele) for a given obs type (e.g., uv) for the reason you state. Different data sources provide different information. Much is common but there are quantities which are unique to a given source. It's not wrong to want to retain this information in the diagnostic file.
Since our diagnostic files are now netcdf I think we can be more flexible in our approach. How best to do this remains an open question. It's my impression the GFS v17 will keep gsi.x. We have time to revisit and refactor our codes to be more robust and efficient. The danger is that this email thread will be forgotten and this issue snags us the next time we pass code to NCO.

[ilianagenkova]

Also see PR #456

[RussTreadon-NOAA]

@ilianagenkova , may we close this issue?

[ilianagenkova]

Yes, thanks for cleaning up old issues!