The model's own error estimate
The artifact stores enough to say how uncertain it is at a parameter point, not just to produce a waveform. Written 2026-09-02, when it first covered all seven modes over the whole waveform; updated the same day with the owner's decision on the two-route blend (checklist B5), below.
What you get
times, errors = model.get_waveform_error_estimate(mass_ratio, spin1z, spin2z)
delta_A, delta_phi = errors[(2, 2)]
The name and the argument spelling already matched the public surface settled
by checklist A7, so nothing about this call changed there; the times it
returns are the same geometric, peak-aligned samples
get_td_waveform_modes_geometric returns, and
waveform-evaluation.md covers those.
Seven modes, 15775 samples each, over the full inspiral, merger and ringdown.
delta_A is a 1-sigma amplitude uncertainty in the units of the mode content;
delta_phi is delta_A / |c_lm| in radians, and is meaningful ONLY where the
mode stands well above its own envelope.
It is a band, not a number. Three stages, all per-sample:
predict_var (14, 25) variance of each Chebyshev coefficient
14 elements by 25 coefficients
error_envelope (N,) 1 sigma at every sample
mode_errors two (N,) amplitude and phase
The envelope propagates the coefficient variances through the same Chebyshev
evaluation the waveform uses,
env(Phi)^2 = sum_j T_j(x(Phi))^2 Var[C_j] on the element containing each
sample, so the band is wide where the fit is poorly determined and narrow where
it is not.
What it reports, measured
Relative amplitude uncertainty at q = 4, chi = (0.3, -0.2), on the released artifact:
| mode | inspiral | ramp | ringdown |
|---|---|---|---|
| (2,2) | 0.05% | 0.15% | 4.31% |
| (3,3) | 0.23% | 0.41% | 4.03% |
| (4,4) | 0.48% | 0.45% | 14.00% |
It widens by about 86x from inspiral to ringdown on (2,2), which is the expected direction.
Where it comes from, and one trap
The surrogate is a Gaussian process regression, so fitting returns a coefficient AND its variance. Those stored variances are what is reported: the fit's own statement of how well it pinned each coefficient down, which grows away from training simulations.
The merger-ringdown envelope must NOT come from error_envelope. Under the
SWSpH representation the regressed coefficients are the detrended spheroidal
residual, so error_envelope returns the uncertainty of a DIFFERENT quantity,
before re-winding, mixing and the remnant-map term. Wrong, not merely
uncalibrated, and it would look plausible.
spherical_plunge_spheroidal_ringdown_error_envelope is the correct source and
the dispatch is now guarded, so an unrecognised arm raises. The released arm is
SWSpH.
The two envelopes are joined with the SAME blend window the waveform uses,
INSPIRAL_MERGER_RINGDOWN_BLEND, named once for both joins.
remnant_sigma defaults to (1.2e-4, 5.3e-5), sourced exactly to
loo_Mf_rms = 1.1987e-04 and loo_chif_rms = 5.30597e-05 in
paper_data/swsph_remnant_map.npz. Three qualifications: the M_f number
excludes SXS:BBH:1124, which would make it 2.7e-4; they are RMS rather than
maxima; and that file's dphi_tot_lm columns are quoted at t_end = 55.5 M
while the arm reaches +124.5, and that term grows linearly, so those columns
are 2.24x low.
Limits
Calibration, not blindness. The merger-ringdown arm's error is
parameter-space regression essentially in full, e_gpr 3.44e-03 against
e_own 2.55e-05, a 115x ratio, which is exactly the kind of error a GPR
variance is meant to describe. The open question is whether the band is
correctly SIZED. One measurement bears on it: the in-sample full-IMR mismatch
of 7.90e-07 on (2,2) implies about 1.26e-03 relative error while the envelope
reports 4.77e-04, so it is a factor 2.6 OPTIMISTIC, in sample. Out of sample it
will be worse by an unmeasured amount. See [[error-aware-likelihood]], whose
limit 1 carries the retraction of an earlier and wrong framing of this.
The timing term is absent, and it is the largest known omission.
m omega sigma_tau EXCEEDS the reported content phase uncertainty on every
mode, by 1.7x to 28.5x, and by 28.5x on (2,2). Folding it in is blocked on
measurement rather than code: the tau to Phi(t) map is an inversion, so any
folded term is a linearisation with unmeasured second order, the clock and the
content are fitted from the same simulations so they are not independent, and
folding it moves four published numbers by more than an order of magnitude.
Only the diagonal coefficient variances are propagated, so coefficient correlations are discarded.
The reported uncertainty is that of the route carrying the weight, not of
the combination that is served. Owner decision 2026-09-02, checklist B5:
the estimate is NOT propagated through the two-route blend, this is documented
rather than changed, and improvements ride a later release. What it means: on
the action-angle inspiral every mode is served as c = w_D cD + w_R cR, an
inverse-variance combination of the direct route and the PN-ratio route with
one scalar weight per route per mode (aa_api._predict_pieces_in_arm), and
what this method reports on the inspiral is the DIRECT route's envelope over
the direct route's amplitude. Measured at q = 4, chi = (0.3, -0.2):
| mode | (2,2) | (2,1) | (3,3) | (3,2) | (4,4) | (4,3) | (5,5) |
|---|---|---|---|---|---|---|---|
w_D |
0.849 | 0.9976 | 0.9898 | 0.9992 | 0.9988 | 0.9900 | 0.9989 |
So on six modes the served content is the direct route to one percent or better and the band describes what is served; on the dominant mode the ratio route carries 15 percent and the band describes 85 percent of it.
Why it is not simply propagated, because a reader will ask: the two routes are
NOT independent. Both are fitted from the same 260 simulations, the ratio
being the same coorbital coefficient divided by its PN weight, so the textbook
inverse-variance result vD vR / (vD + vR) understates, and it would NARROW a
band already measured a factor 2.6 optimistic in sample (above). Three
further reasons are recorded: propagating moves the four published odd-m
numbers by 0.2 to 5.7 percent; the AdA arm has no route split, so an AA-only
propagation would make the two arms' estimates differ in kind; and the weights
would have to be recomputed outside the serving path, a second copy free to
drift. The fully correlated limit
(w_D sqrt(vD) + w_R sqrt(vR))^2 is a defensible upper bound and is the
interim if one is ever wanted. The proper fix, for a later release, is not to
propagate at all: measure the BLENDED content's own leave-one-out residual and
calibrate one variance for it, exactly as the fourteen existing calibration
factors were produced (commit aad8d7c0), which makes the correlation question
disappear rather than answering it. It costs a leave-one-out campaign.
The envelope is DISCONTINUOUS at the fit's own element boundaries. Found 2026-09-02 while building the paper figure, and new: earlier work established continuity at the ARM JOIN, which is a different thing. Adjacent samples differ by up to a factor 4.11 on (2,2) at t = +2.17 M, beside the merger-ringdown element break at +2.0625 M, and every large jump sits within one grid spacing (0.25 M) of a break, while the served amplitude over the same region changes by at most 0.00977 in log10, which is only the ringdown decay. The cause is structural: the estimate is assembled per element and nothing constrains it to agree across a boundary. Consequences for a caller: differentiating the envelope, or reading a single point near a break, gives an answer that depends on which side of the boundary the sample landed.
The presentation hides the peak. docs/USAGE.md quotes the MEDIAN of these
arrays. For (3,3) the median is 2.362e-05 and the maximum 4.756e-04, twenty
times larger. The envelope-max column carries it, but the median reads as the
answer.
Dynamics uncertainty, which is separate and opt-in
model.get_dynamics_error_estimate(mass_ratio, spin1z, spin2z)
reaches the E, J, k, tau and xpc variances, which are fitted and
stored in the artifact and which nothing could reach before 2026-09-02. It is a
separate method because the product is different: no mode index, different
units, no merger-ringdown counterpart. It costs 4.9 ms against 109 ms for a
waveform, and nothing if you do not call it.
Also true
predict_var warns that a point is outside the training hull at q = 8. The
stored log-q bound is 2.0794264247 against log(8) = 2.0794415417, so the
largest trained mass ratio is 7.999879 and the 8 in NRHJSur3dq8 is 1.51e-05
outside the model's own hull. Owner decision 2026-09-02: leave it and
document it, which this paragraph and
artifacts-and-resolution.md do; the three
options considered are recorded with the model. model.in_hull is the
authority, and q = 7.9998 is inside.
Where the code is
| what | where |
|---|---|
| the three functions | nrhjsurrogate/waveform_error_estimate.py |
| the model-level entry point | nrhjsurrogate/driver/aa_api.py |
| merger-ringdown envelope | nrhjsurrogate/physics/spin_weighted_angular_harmonics.py |
| folded into a likelihood | see [[error-aware-likelihood]] |
| worked example | examples/usage/08_model_error.py |