Possible truncation bug in the binning of observations to baroclinic time steps

[stef]

I'm trying to learn 4dvar and noticed that some observations got rejected, even though the observation locations and time looked fine.

Am I doing the time conversion for observations with the wrong numerical precision, or is there a minor bug in the binning of observations to baroclinic time steps caused by truncation errors?

For an observation made at 2021-01-03 22:11:40, the POSIX timestamp (units seconds since 1970/1/1 00:00) is

to = 1609708300.0

Saving it to the obs file (units of days), it becomes

to/(24*3600) = 18630.88310185185

In output.F at the following lines

Code: Select all

      IF (((time(ng)-0.5_r8*dt(ng)).le.ObsTime(ng)).and.                &
     &    (ObsTime(ng).lt.(time(ng)+0.5_r8*dt(ng)))) THEN
             CALL obs_read (ng, iNLM, .FALSE.)

ObsTime(ng) = 1609708299.9999998, and not 1609708300.0, I guess because of truncation errors:

18630.88310185185 * (24*3600) = 1609708299.9999998

Because of this, for dt=200 the lines in output.F call obs_read() when

time(ng)-0.5_r8*dt(ng) = 1609708100.0000000 (2021-01-03 22:08:20)
time(ng)+0.5_r8*dt(ng) = 1609708300.0000000 (2021-01-03 22:11:40)

but then, obs_write() calculates the window bounds in units of days resulting in

TimeLB 18630.880787037036
TimeUB 18630.883101851850
Tobs(iobs) 18630.883101851850

and the observations get rejected because they are not in the interval

[TimeLB,TimeUB)

In summary, in output.F the following line evaluates to True:

time(ng)-0.5_r8*dt(ng) <= ObsTime(ng)) < time(ng)+0.5_r8*dt(ng)

but in extract_obs2d() (called from obs_write() the following line evaluates to False:

TimeLB <= Tobs(iobs) < TimeUB

I.e. the observations lying exactly *on* the right boundary of the right-open interval get rejected.

Now that I know the problem I can easily fix it by properly centering the observations in the window when I do the binning to superobservations, which is actually what I had intended to do because I read in

Tutorial 10: Building Your Observation Files

that

time - 0.5*dt < ObsTime < time + 0.5*dt

However, I thought that 'time(ng)' is staggered such that it coincides with DSTART.


Which ROMS version did I use? I forked from git commit 58bb0f88fbec7548a3d, a.k.a. 'svn: src927 Wed Feb 8 16:45:30 2023 +0100'. I hope it's not outdated?