Skip to content

Fix TimeOperLinInterp: interpolate between the states bracketing each obs - #81

Open
pera-noaa wants to merge 2 commits into
JCSDA:developfrom
pera-noaa:bugfix/timeoper-linear-interp
Open

pera-noaa wants to merge 2 commits into
JCSDA:developfrom
pera-noaa:bugfix/timeoper-linear-interp

Conversation

@pera-noaa

Copy link
Copy Markdown

Description

TimeOperLinInterp does not interpolate in time. Since #2688 (6577fd0) it computes w*gv + (1-w)*gv on a single path per obs, so it returns nearest-state H(x). The existing timeoper tests miss this because their GeoVaLs are the same at every time.

Fix:

  • locations() requests two paths per obs, at the bracketing state times (windowSub apart). A path on an excluded window bound is moved 1 s inside.
  • NL/TL store w*before + (1-w)*after in the reduced GeoVaLs. The AD is the exact adjoint. computeReducedVars is overridden, so filters see interpolated profiles.
  • GeoVaLs with one path per obs pass through unchanged, so the existing tests are untouched.
  • The constructor aborts if windowSub does not divide the window.

The new test, opr_timeoper_two_states, needs no data files. It fills GeoVaLs as nearest-state GetValues would, from three states that differ in time. It then checks H(x), the TL and the adjoint against hand-computed values.

Use with default (nearest) GetValues. GetValues' own TL/AD are nearest-only (GetValues.h:569, 764), so this is the only NL/TL/AD-consistent linear time interpolation.

Issue(s) addressed

None (issues are disabled on JCSDA/ufo).

Dependencies

None.

Impact

None known: public JCSDA and NOAA-EMC code has no TimeOperLinInterp configs outside ufo's tests. Any existing user now gets interpolated H(x).

Manual Testing Instructions

Orion, spack-stack 2.1.1 ue-oneapi-2025.3.1, Release; ufo-data develop 55298a5.

Commit Test Result
42c5ed5 opr_timeoper_two_states FAIL: H(x) 280 vs 281
611442f opr_timeoper_two_states PASS (adjoint 6705.33 = 6705.33)
611442f opr_timeoper, opr_timeoper_pe4, linopr_timeoper, coding_norms PASS
ctest -R 'opr_timeoper|linopr_timeoper' --output-on-failure

Checklist

  • I have performed a self-review of my own code
  • I have made corresponding changes to the documentation (none exists yet; can add to jedi-docs)
  • I have run the unit tests before creating the PR

🤖 Generated with Claude Code

The existing timeoper test reads single-time GeoVaLs, so it passes whether or
not the operator interpolates in time. The new test fills GeoVaLs the way
oops::GetValues does by default (nearest state in time) from three states
3 h apart whose values differ, and checks H(x), the TL and the adjoint against
hand-computed linear interpolation. It fails on the current operator, which
returns the nearest-state value.
… obs

Since #2688 (6577fd0) removed GeoVaLs::split, simulateObs, setTrajectory and
simulateObsTL have computed w0*gv + w1*gv from two copies of the same GeoVaLs,
and locations() has returned one path per obs at the obs time. With w0 + w1 = 1
the operator returned the nearest-state H(x) and did no time interpolation.

locations() now requests two interpolation paths per obs, at the model states
before and after it (windowSub apart; a path on an excluded window bound is
moved one second inside). The forward, trajectory and TL operators store
w0*before + w1*after in the reduced GeoVaLs format and run the wrapped operator
on it; the AD is the exact adjoint. computeReducedVars is overridden so filters
and bias predictors see the interpolated profiles. GeoVaLs with one path per
obs (e.g. read from single-time test files) are passed through unchanged.
@pera-noaa

Copy link
Copy Markdown
Author

@wsmigaj @fmahebert — Bo Huang suggested you as reviewers. Wojciech, this touches the multi-path GeoVaLs path from #2688; the new test fails on develop and passes with the fix. Happy to have it carried into JCSDA-internal however suits you.

@MarekWlasak

MarekWlasak commented Oct 9, 2026 •

Copy link
Copy Markdown
Contributor

Note that there is now separate code in OOPS for time-interpolation of state GeoVaLs

@fmahebert

Copy link
Copy Markdown
Contributor

The OOPS-side capability @MarekWlasak mentioned is enabled via time interpolation: linear in the yaml for the observer. See the QG yaml hofx_tinterp.yaml for an example.

I'm not familiar with the details of the TimeOperLinInterp operator, but I wonder if the OOPS capability would work for you, and whether the operator is perhaps redundant.

@pera-noaa

Copy link
Copy Markdown
Author

@fmahebert @MarekWlasak Thanks. For H(x) the OOPS option does the same job, but its TL and adjoint are still nearest-state: GetValues widens the window only in the nonlinear path (L395), while processTL and processAD keep halfWidth_ = tstep/2 (L628, L805). So a 4D-Var with time interpolation: linear pairs an interpolated H(x) with a nearest-state TL and adjoint; TimeOperLinInterp is consistent in all three. Once the OOPS TL/AD interpolate, I'd be glad to retire it.

@MarekWlasak

Copy link
Copy Markdown
Contributor

@pera-noaa You are right. I would be happy for this code to retire too if we can get oops to do optionally the TL and AD part. The problem there, I think, (not 100% sure) is that in oops there were problems with having the TL and AD parts of this with 4DEnVar.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants