How I proved a from-scratch harmonic tide engine is accurate

You wrote a tide engine. It computes a number: 1.42 m at 3:10 PM. How do you know
that number is right? You can't measure the future ocean. You can't diff against
"the truth" because the truth is a prediction too. This is the whole problem with a
from-scratch harmonic tide predictor — the output looks authoritative the day you
first run it, and authoritative-looking is not the same as correct.

TL;DR — Validate a harmonic tide engine on two axes. Algorithm: golden
vectors against a trusted oracle (@neaps/tide-predictor) — every layer to
floating-point agreement. Reality: head-to-head against the tide authority's
own published predictions — NOAA CO-OPS at Friday Harbor, max 7.9 min / 3.5 cm
across 12 highs and lows, inside the ±15 min / ±0.15 m you'd accept from a printed
table. And before you trust any of it: if your predictions read higher than the
printed tide tables, that's usually not a bug — it's the chart-datum offset.
Jump to the numbers.

The trap that makes a correct engine look wrong

Here's the first thing that will convince you your engine is broken when it isn't.

Compute the harmonic tide height for a station, line it up against the printed tide
tables for the same minute, and your number comes out higher — consistently,
by a fixed amount. The instinct is to hunt for the bug. There isn't one.

A harmonic sum is referenced to the datum its constituents were derived against —
typically mean sea level (MSL), the average of all water levels. Printed tide
tables and nautical charts are referenced to chart datum: MLLW (mean lower
low water) in the US, LAT (lowest astronomical tide) in Canada. Chart datum sits
below MSL by design — it's a low-water reference so that charted depths are
conservative. So a raw harmonic prediction reads higher than the table by exactly the
distance between the two datums.

The fix is a single per-station offset. In this engine it's one optional field on the
station:

let station = Station(
    constituents: [HarmonicConstituent(name: "M2", amplitude: 0.96, phase: 128) /* … */],
    offset: 1.387  // datum offset (e.g. MSL → MLLW), optional
)

And that offset isn't a fudge factor — it's published metadata. NOAA ships the datum
values per station; the offset is just their difference:

const fh = stations.find((s) => s.id === 'noaa/9449880');
const offset = fh.datums.MSL - fh.datums.MLLW; // shift MSL-relative harmonics to chart datum MLLW

The apps that "match the printed number" without applying this are the ones actually
carrying an error — they've hidden a datum mismatch by accident, or they never summed
the harmonics in the first place. Once you know the offset exists, it's a knob, not a
mystery. (This is a subtlety of the harmonic-prediction path — computing heights
from constituents versus reading a printed table. Hydrographic services like NOAA and
CHS serve their own predictions at their own stated datum; they're not wrong, they're
just answering a different question than a from-scratch sum does.)

That's the trap. Now the actual validation.

Method: golden vectors against an oracle

You can't test a predictor against the ocean. You can test it against another
predictor that's already trusted — and if the two agree to floating-point precision,
your port of the algorithm is faithful by construction.

The oracle here is @neaps/tide-predictor
@0.10.0, a well-exercised JavaScript harmonic engine. The Swift engine is a port of
its algorithm, so the bar isn't "close" — it's the same number. A Node tool
generates golden fixtures from Neaps; swift test replays them and asserts agreement
layer by layer:

swift test                    # golden validation against the Neaps reference
node tools/gen-golden.mjs     # regenerate golden fixtures from @neaps/tide-predictor
node tools/gen-realworld.mjs  # refresh the NOAA real-world fixture

The point of going layer by layer is that a tide prediction is a stack of
computations, and an error in any one of them hides inside a plausible-looking final
curve. Test the whole pipeline only end-to-end and a sign error in the node
corrections can cancel against a phase error somewhere else for the 48 hours you
happened to check. So each layer gets pinned independently:

Layer	Check	Tolerance	Result
Astronomy	mean longitudes + node angles, 8 times across the 18.6-yr nodal cycle	1e-6°	pass
Node corrections	IHO f/u, 17 base constituents × 3 times	1e-6	pass
Constituents	V₀ + compound f/u, ~39 constituents × 2 times	1e-6	pass
Prediction	48 h height series, mixed-tide set	< 1e-6 m	pass
Extremes	hi/lo count, kind, time, height	60 s / 0.02 m	pass

A few of these deserve a note on why they're the layers that break:

• Astronomy across the 18.6-year nodal cycle. The moon's orbital plane precesses over 18.6 years, and the mean longitudes and node angle are sampled across that full cycle — not just at one date — because a bug in the node angle is invisible near one epoch and large near another.
• Node corrections (the f/u factors). These are the per-constituent amplitude (f) and phase (u) adjustments that account for that same nodal modulation. They're the most common place a harmonic engine quietly drifts, which is why they get their own fixture instead of only being checked through the final height.
• Extremes, not just the height curve. Sampling heights on a grid and calling the peaks "high tide" gives you the wrong minute. The extremes layer solves for the actual turning points and is checked on count, kind (high vs low), time, and height separately — the tolerances there (60 s, 0.02 m) are the tight ones.

Every layer lands at floating-point agreement (1e-6). The algorithm is faithful to
the reference. But "faithful to the reference" only proves the port is correct — it
says nothing about whether the reference matches the real world. That's the second
axis.

The datum offset is a step in the harness, not a footnote — it's the only thing standing between "agrees to a few centimetres" and "off by the whole MSL-to-MLLW gap."

The receipt: vs NOAA's own predictions

The stronger test is against a completely independent authority: NOAA's own published
predictions, straight from the CO-OPS Data Retrieval
API. This is a real head-to-head —
feed the engine the published harmonic constituents for a station, then compare its
highs and lows to what NOAA's own internal engine says for the same station and days.

Station: Friday Harbor, WA (NOAA 9449880). Constituents come from
@neaps/tide-database (sourced from
NOAA, public domain, bundled offline). The comparison target is fetched live from
CO-OPS at MLLW, GMT:

const url = `https://api.tidesandcurrents.noaa.gov/api/prod/datagetter`
  + `?begin_date=20260715&end_date=20260717`
  + `&station=9449880&product=predictions&datum=MLLW&interval=hilo`
  + `&units=metric&time_zone=gmt&format=json`;

Twelve highs and lows over 2026-07-15…17:

Metric	Engine vs NOAA	Navigational tolerance
Max time error	7.9 min	±15 min
Max height error	3.5 cm	±0.15 m

The engine reproduces NOAA's published tide tables to a few minutes and a few
centimetres — inside half the tolerance you'd accept from a paper table, on both axes,
computed entirely offline from the bundled constants. That's the claim that matters:
not "agrees with the code I ported from," but "agrees with the tide authority."

The honest residual — why ~8 minutes

I'm not going to tell you the error is zero, because it isn't, and the residual is
worth understanding rather than papering over.

The ~8 min / 3.5 cm gap against CO-OPS is expected. NOAA's operational engine
isn't the same engine — it uses a different node-correction epoch and a somewhat
different constituent set than the Neaps algorithm this port follows. Two correct
harmonic engines, fed the same station's constants but differing in those internals,
will land a few minutes apart. That's not error accumulating toward a wrong answer;
it's the known spread between two legitimate methods, and it sits an order of magnitude
inside navigational tolerance.

The tell that it's method-spread and not a bug: against the Neaps oracle — same
algorithm — the agreement is 1e-6, essentially exact. Against NOAA — different
algorithm — it's minutes. If the port were broken, the first number wouldn't be
floating-point tight. The two results are consistent with exactly one story: a faithful
port of one correct method, differing from another correct method by the amount two
correct methods differ.

A note on data sourcing

The constituents themselves are public. NOAA harmonic constants are public domain and
bundled with the engine (~3400 stations via @neaps/tide-database), which is what makes
offline prediction possible at all — no network, no key, no live dependency, predictions
years ahead from data that ships with the app. For Canadian waters, CHS publishes
predictions (not just constants) through its IWLS API, fetched and cached online.

One gotcha worth flagging for anyone bundling multi-source constants: constituent
naming isn't standardized. NOAA writes NU2, MM, RHO; CHS and other sources use
their own variants for the same physical constituents. Feed a station's published
constants straight in and the ones whose names don't match your catalog silently drop
out of the sum — a prediction that's subtly, plausibly wrong. The engine's catalog
resolves 83 aliases to canonical names so published constants from any source
predict correctly.

Close

A tide number you computed yourself is only worth anything if you can show your work.
Two axes do it: golden vectors prove the algorithm is a faithful port, and a
head-to-head against the authority's own predictions proves the port matches reality —
with the datum offset applied so you're comparing to the right zero. The engine is
MIT-licensed and the validation is reproducible from the fixtures:
slackwater-engine. If it's
off in your home waters, the fixtures are right there — check it and send a number back.

Not for navigation. Predictions are astronomical estimates and don't account for
weather, surge, or local effects. Carry official tables and charts.