Added driver download and convenience functions

benches/density.rs

@@ -89,5 +89,33 @@ fn entry_points(c: &mut Criterion) {
     group.finish();
 }
 
+/// What the compiled-in driver record costs.
+///
+/// The point of measuring it is that it must not be on the path of a caller who supplies
+/// its own drivers, and must not be large for one who does not.
+#[cfg(feature = "drivers")]
+fn driver_lookup(c: &mut Criterion) {
+    use nrlmsise00::drivers::{self, Drivers};
+
+    let model = Model::new();
+    let record = Drivers::embedded();
+    let when = drivers::mjd(2015, 3, 15, 8, 3, 20.0);
+
+    let mut group = c.benchmark_group("drivers");
+    group.bench_function("one day out of the record", |b| {
+        b.iter(|| record.day(black_box(when)).unwrap())
+    });
+    group.bench_function("a full driver set", |b| {
+        b.iter(|| record.at(black_box(when)).unwrap())
+    });
+    group.bench_function("gtd7d_at, lookup and model together", |b| {
+        b.iter(|| model.gtd7d_at(black_box(when), 400.0, 60.0, -70.0).unwrap())
+    });
+    group.finish();
+}
+
+#[cfg(feature = "drivers")]
+criterion_group!(benches, gtd7d, derivatives, entry_points, driver_lookup);
+#[cfg(not(feature = "drivers"))]
 criterion_group!(benches, gtd7d, derivatives, entry_points);
 criterion_main!(benches);

data/gfz-drivers.txt

@@ -0,0 +1,36 @@
+The solar and geomagnetic record NRLMSISE-00 reads.
+
+Source      https://kp.gfz-potsdam.de/app/files/Kp_ap_Ap_SN_F107_since_1932.txt
+Attribution GFZ German Research Centre for Geosciences and Geomagnetic Observatory Niemegk, CC BY 4.0
+Downloaded  2026-10-02
+Digest      sha256 cfd3c7dcf021f1d4125f8ecc71a8c5edcb13efe228d111fead95b628e3ba5b4e
+            of the file as downloaded, not of the blob next to this note.
+
+Coverage    1957-01-02 to 2026-10-01, 25475 days
+            Modified Julian Day 35840 to 61314
+            The last day a model can be evaluated on is forty days earlier than the last
+            day held, because the 81-day flux average the model wants is centered.
+
+Flux        The observed 10.7 cm radio flux, at the Earth's actual distance from the Sun,
+            which is what the model asks for. The flux adjusted to one astronomical unit
+            is not kept. 56 of 25475 days were not measured and are
+            interpolated from the days either side; they are marked in the blob and
+            reported by Day::flux_measured.
+
+Geomagnetic The daily mean Ap and the eight three-hour ap indices. 153 of
+            25475 days are definitive rather than preliminary, which Day::definitive
+            reports.
+
+Extremes    The flux runs from 53.5 on 1997-03-01 to 938.6 on 2011-03-07.
+            Both ends are in the source and are passed through unaltered. The high end is
+            real: the daily flux includes radio bursts, so a large flare puts it several
+            times above the solar-maximum background, and there are 7 measured
+            days above 400. The low end looks like instrument trouble rather than the Sun —
+            2 measured days fall below 60, each an isolated day between
+            neighbors around 70, and both are in early March. Whether to smooth them is a
+            decision for whoever is using the record, not for the record.
+
+Dropped     The sunspot number, which is the one part of the source under CC BY-NC 4.0,
+            and the adjusted flux.
+
+Regenerate  scripts/fetch-drivers.py

hooks/pre-commit

@@ -1,12 +1,37 @@
 #!/bin/sh
-# Cheap checks, on the whole tree. The expensive ones are in pre-push.
+# The cheap checks, over the whole workspace. The expensive ones are in pre-push.
+#
+# Measured on a warm tree: under a second for all four. A cold target directory costs half
+# a minute once and nothing after.
+#
+# Silent when everything passes. On a failure it names the check and prints what the tool
+# said, which is what a commit from an editor needs in order to be diagnosable.
 set -e
-cargo fmt --check
-cargo clippy --all-targets --all-features -- -D warnings
-cargo test --lib
+
+if ! command -v cargo >/dev/null 2>&1; then
+    echo "pre-commit: cargo is not on PATH" >&2
+    echo "pre-commit: PATH is $PATH" >&2
+    exit 1
+fi
+
+check() {
+    name=$1
+    shift
+    if ! output=$("$@" 2>&1); then
+        echo "pre-commit: $name failed" >&2
+        printf '%s\n' "$output" >&2
+        exit 1
+    fi
+}
+
+check "formatting" cargo fmt --all --check
+check "clippy" cargo clippy --workspace --all-targets --all-features -- -D warnings
+check "the library tests" cargo test --lib
+
 # British spellings keep slipping into generated prose. This file is excluded because it
 # holds the word list, which would otherwise match itself.
-if git diff --cached -U0 -- . ':(exclude)hooks/pre-commit' | grep -nE '^\+.*\b(metre|metres|normalis|stabilis|recognis|licence|modelling|manoeuvre|behaviour|analyse|analysing|centre|colour)' ; then
+if git diff --cached -U0 -- . ':(exclude)hooks/pre-commit' |
+    grep -nE '^\+.*\b(metre|metres|normalis|stabilis|recognis|licence|modelling|manoeuvre|behaviour|analyse|analysing|centre|colour)'; then
     echo "pre-commit: British spelling above; this project is US English" >&2
     exit 1
 fi

hooks/pre-push

@@ -1,6 +1,22 @@
 #!/bin/sh
-# The slow checks: doctests link a binary apiece, and rustdoc is its own pass.
+# The slow checks: every doctest links its own binary, and rustdoc is its own pass.
+#
+# Measured on a warm tree: about a second and a half for all three. One of the tests reads
+# the grid fixture, so this is also what catches a fixture regenerated but not committed.
+#
+# Silent when everything passes.
 set -e
-cargo test
-cargo test --no-default-features
-cargo doc --no-deps
+
+check() {
+    name=$1
+    shift
+    if ! output=$("$@" 2>&1); then
+        echo "pre-push: $name failed" >&2
+        printf '%s\n' "$output" >&2
+        exit 1
+    fi
+}
+
+check "the tests" cargo test --workspace
+check "the tests without default features" cargo test --no-default-features
+check "the docs" cargo doc --no-deps --workspace

readme.md

@@ -154,10 +154,40 @@ own them: `f107`, `f107_average` and the ap history are plain numbers on the inp
 caller with its own record passes them straight in and pays nothing for the convenience it
 is not using.
 
-The `drivers` feature, on by default, is **not implemented yet**. It will add a
-compiled-in record of those measurements and a lookup by date, for callers who would rather
-not source their own. The feature is declared now so that turning it on later is not a
-change to what `default` means.
+The `drivers` feature, on by default, is for callers who would rather not. It compiles in
+GFZ Potsdam's record — 25475 days, 1957 to the present, 535 KB — and looks a date up in it:
+
+```rust
+# #[cfg(feature = "drivers")]
+# fn main() {
+use nrlmsise00::{drivers, Model};
+
+let model = Model::new();
+let when = drivers::mjd(2015, 3, 15, 8, 3, 20.0);   // or hifitime's to_mjd_utc_days()
+let output = model.gtd7d_at(when, 400.0, 60.0, -70.0).expect("inside the record");
+
+assert!((output.total_mass - 2.3e-12).abs() < 0.5e-12);
+# }
+# #[cfg(not(feature = "drivers"))]
+# fn main() {}
+```
+
+That picks the previous day's flux, the 81-day centered average and the seven-element
+magnetic history, and derives the day of the year, the time of day and a consistent local
+solar time. `Drivers::embedded()` exposes the record itself if you want the parts.
+
+Nothing is read or allocated until something asks for a day. Measured: a lookup is 70 ns
+against an evaluation's 1.76 us, the first call is 278 us to accumulate the flux averages,
+and a caller supplying its own drivers touches none of it. Dates are Modified Julian Days
+and the module needs no time library; `drivers::mjd` builds one from a civil date.
+
+The record is measurements, so a date past it is refused rather than extrapolated, and the
+usable range is forty days shorter than the record **at both ends**, because the 81-day
+flux average is centered. Refreshing it is a deliberate act — `scripts/fetch-drivers.py` —
+not something that happens inside an evaluation: a program that fetched its own weather
+would answer differently on different days, and nothing here should acquire a network
+failure mode. `drivers::provenance()` returns the source, license, digest and coverage of
+whatever is compiled in.
 
 ## The model is not this crate's work
 

scripts/fetch-drivers.py

@@ -0,0 +1,257 @@
+#!/usr/bin/env python3
+"""Regenerates data/gfz-drivers.bin, the solar and geomagnetic record the model reads.
+
+Usage: scripts/fetch-drivers.py [path to a local copy of the source file]
+
+The record is committed as a compact binary blob and embedded with include_bytes!, so a
+caller needs no network, no data directory and nothing to parse at startup. The blob is
+numbers only; its provenance is written alongside it in data/gfz-drivers.txt, which the
+crate compiles in as a string so that a program can print where its weather came from.
+
+Source: GFZ Potsdam's combined Kp/ap/Ap/SN/F10.7 record, under CC BY 4.0. Chosen over
+Celestrak's SW-All.csv for two reasons: it states a license where Celestrak states none,
+which matters for a file committed into an MIT crate, and Celestrak's runs decades past
+today on its own predictions, which are not measurements.
+
+The sunspot column is dropped, being the one part of the GFZ file under CC BY-NC 4.0, and
+so is the flux adjusted to one astronomical unit: NRLMSISE-00 asks for the flux at the
+Earth's actual distance from the Sun, which is the observed column.
+
+Unlike the coefficient tables, this record is meant to change — it gains a row a day — so
+there is no digest pinned in the crate. What is recorded instead is the digest of what was
+downloaded, with the date and the coverage, in the provenance file.
+"""
+
+import hashlib
+import struct
+import subprocess
+import sys
+from datetime import date, timedelta
+from pathlib import Path
+
+URL = "https://kp.gfz-potsdam.de/app/files/Kp_ap_Ap_SN_F107_since_1932.txt"
+ATTRIBUTION = (
+    "GFZ German Research Centre for Geosciences and Geomagnetic Observatory Niemegk, "
+    "CC BY 4.0"
+)
+MISSING = -1.0
+
+# Geomagnetic indices reach back to 1932 and the solar flux to February 1947, but the early
+# flux record is sporadic: 105 missing days in 1947 and one hole of 44. From 1957 it is all
+# but complete, which is where this record starts and where the other published driver
+# files start too.
+START = date(1957, 1, 1)
+
+# A run of unmeasured flux days this long is refused rather than interpolated.
+LONGEST_FILL = 7
+
+ROOT = Path(__file__).resolve().parent.parent
+OUT = ROOT / "data" / "gfz-drivers.bin"
+NOTE = ROOT / "data" / "gfz-drivers.txt"
+
+# Modified Julian Day zero. The record is indexed by MJD because a UT day is not always
+# 86400 seconds and MJD still increments by exactly one across a leap second.
+MJD_EPOCH = date(1858, 11, 17)
+
+MAGIC = b"NRLMSISEDRV"
+VERSION = 1
+
+# Field positions in the source, which is whitespace-aligned with 28 columns.
+YEAR, MONTH, DAY = 0, 1, 2
+AP_FIRST, AP_LAST = 15, 22
+AP_DAILY = 23
+F107_OBSERVED = 25
+DEFINITIVE = 27
+
+
+def source_text(argv):
+    if len(argv) > 1:
+        raw = Path(argv[1]).read_bytes()
+        print(f"reading {argv[1]}")
+    else:
+        print(f"downloading {URL}")
+        raw = subprocess.run(
+            ["curl", "-sS", "--fail", "--max-time", "180", "-L", URL],
+            check=True,
+            capture_output=True,
+        ).stdout
+    return raw.decode(), hashlib.sha256(raw).hexdigest()
+
+
+def parse(text):
+    """Daily rows from START onward, with the flux left as None where unmeasured."""
+    rows = []
+    for line in text.splitlines():
+        if line.startswith("#") or not line.strip():
+            continue
+        fields = line.split()
+        if len(fields) != 28:
+            raise SystemExit(f"expected 28 fields, found {len(fields)}: {line[:60]}")
+
+        when = date(int(fields[YEAR]), int(fields[MONTH]), int(fields[DAY]))
+        if when < START:
+            continue
+
+        observed = float(fields[F107_OBSERVED])
+        rows.append(
+            {
+                "date": when,
+                "ap": [int(float(f)) for f in fields[AP_FIRST : AP_LAST + 1]],
+                "ap_daily": int(float(fields[AP_DAILY])),
+                "flux": None if observed <= MISSING else observed,
+                "definitive": int(fields[DEFINITIVE]) == 1,
+                "filled": False,
+            }
+        )
+
+    if not rows:
+        raise SystemExit(f"no rows at or after {START}")
+
+    # The record is addressed by day, so a missing date would shift everything after it.
+    for earlier, later in zip(rows, rows[1:]):
+        if later["date"] - earlier["date"] != timedelta(days=1):
+            raise SystemExit(f"missing day between {earlier['date']} and {later['date']}")
+
+    return rows
+
+
+def trim(rows):
+    """Drops unmeasured days at either end.
+
+    START is a round date rather than a measurement and the flux happens not to have been
+    taken on it. The far end matters more: an unmeasured day there would have to be
+    extrapolated, and refusing to evaluate past the record is precisely so that nothing is.
+    """
+    measured = [index for index, row in enumerate(rows) if row["flux"] is not None]
+    if not measured:
+        raise SystemExit("no measured solar flux at all")
+    return rows[measured[0] : measured[-1] + 1]
+
+
+def fill(rows):
+    """Interpolates the days the flux was not measured, marks them, counts them.
+
+    The gaps are isolated — one or two days, five at the very worst — and the flux moves a
+    few percent over that, so a straight line between the measurements either side is a
+    fair reconstruction. A long gap is not, and is refused.
+    """
+    measured = [index for index, row in enumerate(rows) if row["flux"] is not None]
+    filled = 0
+    for before, after in zip(measured, measured[1:]):
+        span = after - before
+        if span == 1:
+            continue
+        if span - 1 >= LONGEST_FILL:
+            raise SystemExit(
+                f"{span - 1} unmeasured days from {rows[before + 1]['date']}, too many "
+                f"to interpolate"
+            )
+        for offset in range(1, span):
+            weight = offset / span
+            rows[before + offset]["flux"] = (
+                rows[before]["flux"] * (1.0 - weight) + rows[after]["flux"] * weight
+            )
+            rows[before + offset]["filled"] = True
+            filled += 1
+    return filled
+
+
+def mjd(when):
+    return when.toordinal() - MJD_EPOCH.toordinal()
+
+
+def encode(rows, generated):
+    """The blob: a fixed header, then twenty-one bytes a day, little-endian throughout.
+
+    The flux is stored in tenths because that is the precision the source states it to, so
+    nothing is lost and a float is not needed. ap and Ap are whole numbers on the Kp scale,
+    which reaches 400.
+    """
+    blob = bytearray()
+    blob += MAGIC
+    blob += struct.pack("<BiI", VERSION, mjd(rows[0]["date"]), len(rows))
+    blob += struct.pack("<I", int(generated.strftime("%Y%m%d")))
+
+    for row in rows:
+        tenths = round(row["flux"] * 10.0)
+        if not 0 <= tenths <= 0xFFFF:
+            raise SystemExit(f"{row['date']}: flux {row['flux']} does not fit")
+        values = [row["ap_daily"], *row["ap"]]
+        for value in values:
+            if not 0 <= value <= 0xFFFF:
+                raise SystemExit(f"{row['date']}: ap {value} does not fit")
+        flags = (1 if row["definitive"] else 0) | (0 if row["filled"] else 2)
+        blob += struct.pack("<H", tenths)
+        blob += struct.pack("<9H", *values)
+        blob += struct.pack("<B", flags)
+    return bytes(blob)
+
+
+def note(rows, digest, generated, filled):
+    """The provenance, which the crate compiles in as a string."""
+    definitive = sum(1 for row in rows if row["definitive"])
+    highest = max(rows, key=lambda row: row["flux"])
+    lowest = min(rows, key=lambda row: row["flux"])
+    dropouts = [row for row in rows if row["flux"] < 60.0 and not row["filled"]]
+    peaks = [row for row in rows if row["flux"] > 400.0 and not row["filled"]]
+    return f"""The solar and geomagnetic record NRLMSISE-00 reads.
+
+Source      {URL}
+Attribution {ATTRIBUTION}
+Downloaded  {generated.isoformat()}
+Digest      sha256 {digest}
+            of the file as downloaded, not of the blob next to this note.
+
+Coverage    {rows[0]['date'].isoformat()} to {rows[-1]['date'].isoformat()}, \
+{len(rows)} days
+            Modified Julian Day {mjd(rows[0]['date'])} to {mjd(rows[-1]['date'])}
+            The last day a model can be evaluated on is forty days earlier than the last
+            day held, because the 81-day flux average the model wants is centered.
+
+Flux        The observed 10.7 cm radio flux, at the Earth's actual distance from the Sun,
+            which is what the model asks for. The flux adjusted to one astronomical unit
+            is not kept. {filled} of {len(rows)} days were not measured and are
+            interpolated from the days either side; they are marked in the blob and
+            reported by Day::flux_measured.
+
+Geomagnetic The daily mean Ap and the eight three-hour ap indices. {definitive} of
+            {len(rows)} days are definitive rather than preliminary, which Day::definitive
+            reports.
+
+Extremes    The flux runs from {lowest['flux']} on {lowest['date'].isoformat()} to \
+{highest['flux']} on {highest['date'].isoformat()}.
+            Both ends are in the source and are passed through unaltered. The high end is
+            real: the daily flux includes radio bursts, so a large flare puts it several
+            times above the solar-maximum background, and there are {len(peaks)} measured
+            days above 400. The low end looks like instrument trouble rather than the Sun —
+            {len(dropouts)} measured days fall below 60, each an isolated day between
+            neighbors around 70, and both are in early March. Whether to smooth them is a
+            decision for whoever is using the record, not for the record.
+
+Dropped     The sunspot number, which is the one part of the source under CC BY-NC 4.0,
+            and the adjusted flux.
+
+Regenerate  scripts/fetch-drivers.py
+"""
+
+
+def main():
+    text, digest = source_text(sys.argv)
+    rows = trim(parse(text))
+    filled = fill(rows)
+    generated = date.today()
+
+    OUT.parent.mkdir(parents=True, exist_ok=True)
+    OUT.write_bytes(encode(rows, generated))
+    NOTE.write_text(note(rows, digest, generated, filled))
+
+    print(f"digest {digest}")
+    print(f"interpolated {filled} unmeasured days of {len(rows)}")
+    print(
+        f"wrote {OUT.relative_to(ROOT)}: {len(rows)} days, "
+        f"{rows[0]['date']} to {rows[-1]['date']}, {OUT.stat().st_size / 1e3:.0f} KB"
+    )
+
+
+if __name__ == "__main__":
+    main()

src/drivers.rs

@@ -0,0 +1,471 @@
+//! The solar and geomagnetic record the model reads, and lookups into it by date.
+//!
+//! The model's flux and magnetic inputs are measurements of what the Sun did. This crate
+//! does not need to own them — [`crate::Input`] takes them as plain numbers, so a caller
+//! with its own record passes them straight in — but sourcing them is real work, and this
+//! is here for callers who would rather not.
+//!
+//! The record is compiled in, as half a megabyte of [`include_bytes!`], and nothing is
+//! read, parsed or allocated until something asks for a day. A caller who never calls into
+//! this module pays only the bytes, and a build with the `drivers` feature off does not
+//! even have those.
+//!
+//! # What it costs, measured
+//!
+//! | | |
+//! |---|---|
+//! | [`Drivers::embedded`], first call | 278 us, accumulating 25475 flux averages |
+//! | [`Drivers::embedded`], later calls | 29 ns |
+//! | [`Drivers::day`] | 7.2 ns |
+//! | [`Drivers::at`], a full driver set | 70 ns |
+//! | [`crate::Model::gtd7d_at`] against `gtd7d` | 1.79 us against 1.76 us |
+//!
+//! So a lookup is a few percent of an evaluation, the first one is a quarter of a
+//! millisecond, and a caller passing its own drivers to [`crate::Model::gtd7d`] touches
+//! none of it.
+//!
+//! # Dates
+//!
+//! Everything here is addressed by Modified Julian Day, as a float, with no dependency on
+//! a time library. That is deliberate rather than minimal: a UT day is not always 86400
+//! seconds, and MJD still increments by exactly one across a leap second, where a count of
+//! elapsed seconds slips. [`mjd`] builds one from a civil date; `hifitime`'s
+//! `Epoch::to_mjd_utc_days` returns one directly.
+//!
+//! # The record ends, and so does the evaluation
+//!
+//! There are no measurements for tomorrow, so a day past the record is refused rather than
+//! extrapolated. The usable end is earlier than the last day held, because the 81-day flux
+//! average the model wants is *centered*: the last day that has one is forty days before
+//! the last day that has a flux measurement.
+//!
+//! The front of the record is short by the same forty days, for the same reason, so the
+//! usable range is narrower than the record at both ends.
+//!
+//! Refreshing the record is a deliberate act — `scripts/fetch-drivers.py` — and not
+//! something that happens inside an evaluation. Two reasons: a program that fetched its
+//! own weather would give different answers on different days, and nothing here should
+//! acquire a network failure mode.
+//!
+//! # Provenance
+//!
+//! The record is GFZ Potsdam's, used under CC BY 4.0 and attributed to the GFZ German
+//! Research Centre for Geosciences and Geomagnetic Observatory Niemegk. [`provenance`]
+//! returns the full note, including the digest of what was downloaded and the coverage.
+
+use core::fmt;
+use std::sync::OnceLock;
+
+use crate::input::MagneticActivity;
+
+/// The record as generated by `scripts/fetch-drivers.py`.
+static BLOB: &[u8] = include_bytes!("../data/gfz-drivers.bin");
+
+/// What the record is and where it came from.
+static NOTE: &str = include_str!("../data/gfz-drivers.txt");
+
+/// Days either side of a day that its centered flux average spans.
+///
+/// The average is over 81 days: forty each way, and the day itself.
+pub const FLUX_AVERAGE_RADIUS: i64 = 40;
+
+/// Bytes per day in the blob.
+const RECORD_BYTES: usize = 21;
+
+/// Bytes of header before the first day.
+const HEADER: usize = 24;
+
+/// What the record holds for one day.
+#[derive(Clone, Copy, Debug, PartialEq)]
+pub struct Day {
+    /// The observed 10.7 cm solar flux, in solar flux units.
+    ///
+    /// At the Earth's actual distance from the Sun rather than scaled to one astronomical
+    /// unit, which is what the model asks for.
+    pub flux: f64,
+    /// The daily mean geomagnetic index, Ap.
+    pub ap_daily: f64,
+    /// The eight three-hour geomagnetic indices of the day, earliest first.
+    pub ap: [f64; 8],
+    /// Whether the geomagnetic indices are settled rather than preliminary.
+    pub definitive: bool,
+    /// Whether the flux was measured on this day rather than interpolated across a gap.
+    ///
+    /// Fifty-six of the record's twenty-five thousand days are interpolated, none of them
+    /// recent; the gaps are one or two days and the flux moves a few percent over that.
+    pub flux_measured: bool,
+}
+
+/// Everything the model wants for one moment.
+///
+/// What [`Drivers::at`] assembles. The magnetic activity is always the seven-element
+/// history, which is the mode a propagation should use: a storm's effect on density lags
+/// it by hours, and a single daily mean cannot express that.
+#[derive(Clone, Copy, Debug, PartialEq)]
+pub struct DriverSet {
+    /// The previous day's observed flux, which is what the model means by `f107`.
+    pub f107: f64,
+    /// The 81-day average of the observed flux, centered on the day.
+    pub f107_average: f64,
+    /// The geomagnetic history, in the order the model wants it.
+    pub magnetic: MagneticActivity<f64>,
+}
+
+/// A day the record cannot answer for.
+#[derive(Clone, Copy, Debug, PartialEq)]
+pub enum DriverError {
+    /// Before the record can be used.
+    ///
+    /// Which is forty days after it starts, not the day it starts: the centered flux
+    /// average runs off the front of the record just as it runs off the back.
+    BeforeRecord {
+        /// The Modified Julian Day asked for.
+        asked: f64,
+        /// The first Modified Julian Day a model can be evaluated on.
+        first: i64,
+    },
+    /// After the record can be used, which includes every future date.
+    PastRecord {
+        /// The Modified Julian Day asked for.
+        asked: f64,
+        /// The last Modified Julian Day a model can be evaluated on, which is forty days
+        /// before the last day measured.
+        last: i64,
+    },
+}
+
+impl fmt::Display for DriverError {
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        match self {
+            Self::BeforeRecord { asked, first } => {
+                let (year, month, day) = civil(*first);
+                write!(
+                    f,
+                    "MJD {asked} is before the driver record can be used, which is from \
+                     {year:04}-{month:02}-{day:02}; the centered flux average needs forty \
+                     days before the day it is asked about"
+                )
+            }
+            Self::PastRecord { asked, last } => {
+                let (year, month, day) = civil(*last);
+                write!(
+                    f,
+                    "MJD {asked} is past the driver record, which reaches \
+                     {year:04}-{month:02}-{day:02}; drivers are measurements and there are \
+                     none for later, so refresh with scripts/fetch-drivers.py or evaluate \
+                     inside the record"
+                )
+            }
+        }
+    }
+}
+
+impl std::error::Error for DriverError {}
+
+/// The record.
+pub struct Drivers {
+    /// Modified Julian Day of the first day held.
+    first: i64,
+    /// How many days are held.
+    days: usize,
+    /// The 81-day centered flux average per day, `None` where the window runs off an end.
+    centered: Vec<Option<f64>>,
+}
+
+impl fmt::Debug for Drivers {
+    /// The coverage, rather than twenty-five thousand days of it.
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        let (first_year, first_month, first_day) = civil(self.first);
+        let (last_year, last_month, last_day) = civil(self.last_day());
+        write!(
+            f,
+            "Drivers({first_year:04}-{first_month:02}-{first_day:02} to \
+             {last_year:04}-{last_month:02}-{last_day:02}, {} days)",
+            self.days
+        )
+    }
+}
+
+/// What the record is and where it came from, as the fetch script recorded it.
+pub fn provenance() -> &'static str {
+    NOTE
+}
+
+impl Drivers {
+    /// The compiled-in record.
+    ///
+    /// The flux averages are accumulated on the first call and kept; the days themselves
+    /// are read straight out of the compiled-in bytes, so there is nothing to parse and
+    /// nothing that can fail.
+    pub fn embedded() -> &'static Self {
+        static RECORD: OnceLock<Drivers> = OnceLock::new();
+        RECORD.get_or_init(|| {
+            assert_eq!(&BLOB[..11], b"NRLMSISEDRV", "the driver blob is not one");
+            assert_eq!(
+                BLOB[11], 1,
+                "the driver blob is a version this does not read"
+            );
+            let first = i64::from(i32::from_le_bytes(BLOB[12..16].try_into().unwrap()));
+            let days = u32::from_le_bytes(BLOB[16..20].try_into().unwrap()) as usize;
+            assert_eq!(
+                BLOB.len(),
+                HEADER + days * RECORD_BYTES,
+                "the driver blob is the wrong length for {days} days"
+            );
+            let flux: Vec<f64> = (0..days).map(|index| read(index).flux).collect();
+            Drivers {
+                first,
+                days,
+                centered: centered_averages(&flux),
+            }
+        })
+    }
+
+    /// The first Modified Julian Day the record holds.
+    pub fn first_day(&self) -> i64 {
+        self.first
+    }
+
+    /// The last Modified Julian Day the record holds a measurement for.
+    pub fn last_day(&self) -> i64 {
+        self.first + self.days as i64 - 1
+    }
+
+    /// The first Modified Julian Day a model can be evaluated on.
+    ///
+    /// Forty days after [`Drivers::first_day`]. The centered flux average spans forty days
+    /// either side of the day it describes, so it runs off the front of the record as well
+    /// as off the back, and the usable range is shorter than the record at both ends.
+    pub fn first_usable_day(&self) -> i64 {
+        self.first_day() + FLUX_AVERAGE_RADIUS
+    }
+
+    /// The last Modified Julian Day a model can be evaluated on.
+    ///
+    /// Forty days before [`Drivers::last_day`], for the same reason as
+    /// [`Drivers::first_usable_day`].
+    pub fn last_usable_day(&self) -> i64 {
+        self.last_day() - FLUX_AVERAGE_RADIUS
+    }
+
+    /// How many days the record holds.
+    pub fn len(&self) -> usize {
+        self.days
+    }
+
+    /// Whether the record holds no days at all, which the compiled-in one does not.
+    pub fn is_empty(&self) -> bool {
+        self.days == 0
+    }
+
+    /// Whether a model can be evaluated at this Modified Julian Day.
+    pub fn covers(&self, mjd: f64) -> bool {
+        self.index_of(mjd).is_ok()
+    }
+
+    /// Checks that an arc from `start` to `end` is covered throughout.
+    ///
+    /// Only the ends are examined, which is enough: the record is one contiguous run of
+    /// days, so an arc with both ends inside it has nothing missing in between. The order
+    /// of the two does not matter, since a propagation may run backwards.
+    ///
+    /// # Errors
+    ///
+    /// Naming the end that falls outside, and the day the record reaches.
+    pub fn covers_arc(&self, start: f64, end: f64) -> Result<(), DriverError> {
+        self.index_of(start)?;
+        self.index_of(end)?;
+        Ok(())
+    }
+
+    /// What the record holds for the day a Modified Julian Day falls in.
+    ///
+    /// # Errors
+    ///
+    /// If the day is outside what the record can be evaluated on.
+    pub fn day(&self, mjd: f64) -> Result<Day, DriverError> {
+        Ok(read(self.index_of(mjd)?))
+    }
+
+    /// The 81-day centered average of the observed flux around a day.
+    ///
+    /// # Errors
+    ///
+    /// If the day is outside what the record can be evaluated on.
+    pub fn centered_flux(&self, mjd: f64) -> Result<f64, DriverError> {
+        let index = self.index_of(mjd)?;
+        Ok(self.centered[index].expect("inside the usable range"))
+    }
+
+    /// Everything the model wants at a moment.
+    ///
+    /// The flux is the *previous* day's, which is what the model asks for; the average is
+    /// centered on the day itself; and the magnetic history is the seven elements in the
+    /// model's own order, reaching back fifty-seven hours from the three-hour interval the
+    /// moment falls in.
+    ///
+    /// # Errors
+    ///
+    /// If the day is outside what the record can be evaluated on.
+    pub fn at(&self, mjd: f64) -> Result<DriverSet, DriverError> {
+        let index = self.index_of(mjd)?;
+
+        // Three-hour intervals counted from the start of the record, so reaching back
+        // crosses days with no special case.
+        let fraction = mjd - mjd.floor();
+        let interval = index as i64 * 8 + (fraction * 8.0).floor() as i64;
+        let at = |back: i64| -> f64 {
+            let wanted = interval - back;
+            if wanted < 0 {
+                // Before the record. Its first days are 1957 and nothing is propagated
+                // there; holding the earliest value is enough not to produce a NaN.
+                return read(0).ap[0];
+            }
+            let (day, slot) = ((wanted / 8) as usize, (wanted % 8) as usize);
+            read(day.min(self.days - 1)).ap[slot]
+        };
+        let mean = |first_back: i64| -> f64 {
+            (0..8).map(|offset| at(first_back + offset)).sum::<f64>() / 8.0
+        };
+
+        Ok(DriverSet {
+            f107: read(index.saturating_sub(1)).flux,
+            f107_average: self.centered[index].expect("inside the usable range"),
+            magnetic: MagneticActivity::History {
+                ap: [
+                    read(index).ap_daily,
+                    at(0),
+                    at(1),
+                    at(2),
+                    at(3),
+                    // Twelve hours back is four intervals, and the mean runs back from it.
+                    mean(4),
+                    // Thirty-six hours back is twelve.
+                    mean(12),
+                ],
+            },
+        })
+    }
+
+    /// The index of the day a Modified Julian Day falls in.
+    fn index_of(&self, mjd: f64) -> Result<usize, DriverError> {
+        let asked = mjd.floor() as i64;
+        if asked < self.first_usable_day() {
+            return Err(DriverError::BeforeRecord {
+                asked: mjd,
+                first: self.first_usable_day(),
+            });
+        }
+        let index = (asked - self.first) as usize;
+        if index >= self.days || self.centered[index].is_none() {
+            return Err(DriverError::PastRecord {
+                asked: mjd,
+                last: self.last_usable_day(),
+            });
+        }
+        Ok(index)
+    }
+}
+
+/// One day, read out of the compiled-in bytes.
+fn read(index: usize) -> Day {
+    let at = HEADER + index * RECORD_BYTES;
+    let word =
+        |offset: usize| -> u16 { u16::from_le_bytes([BLOB[at + offset], BLOB[at + offset + 1]]) };
+    let mut ap = [0.0; 8];
+    for (slot, value) in ap.iter_mut().enumerate() {
+        *value = f64::from(word(4 + slot * 2));
+    }
+    let flags = BLOB[at + 20];
+    Day {
+        // Stored in tenths, which is the precision the source states the flux to.
+        flux: f64::from(word(0)) / 10.0,
+        ap_daily: f64::from(word(2)),
+        ap,
+        definitive: flags & 1 != 0,
+        flux_measured: flags & 2 != 0,
+    }
+}
+
+/// The 81-day centered average of the flux for every day that can have one.
+fn centered_averages(flux: &[f64]) -> Vec<Option<f64>> {
+    let radius = FLUX_AVERAGE_RADIUS as usize;
+    let window = 2 * radius + 1;
+    if flux.len() < window {
+        return vec![None; flux.len()];
+    }
+
+    // One running sum rather than a sum per day, which at eighty-one days apiece would be
+    // two million additions over the compiled-in record.
+    let mut running: f64 = flux[..window].iter().sum();
+    let mut averages = vec![None; flux.len()];
+    averages[radius] = Some(running / window as f64);
+    for index in radius + 1..flux.len() - radius {
+        running += flux[index + radius] - flux[index - radius - 1];
+        averages[index] = Some(running / window as f64);
+    }
+    averages
+}
+
+/// Modified Julian Day zero as a count of days from the Unix epoch.
+const MJD_UNIX: i64 = 40587;
+
+/// The Modified Julian Day of a civil date and time, UTC.
+///
+/// Days are counted, not seconds, so a leap second does not shift the result.
+///
+/// ```
+/// # use nrlmsise00::drivers::mjd;
+/// // The Unix epoch is MJD 40587, and noon is half a day into it.
+/// assert_eq!(mjd(1970, 1, 1, 0, 0, 0.0), 40587.0);
+/// assert_eq!(mjd(1970, 1, 1, 12, 0, 0.0), 40587.5);
+/// ```
+pub fn mjd(year: i64, month: u32, day: u32, hour: u32, minute: u32, second: f64) -> f64 {
+    let days = days_from_civil(year, month, day) + MJD_UNIX;
+    let fraction = (f64::from(hour) * 3600.0 + f64::from(minute) * 60.0 + second) / 86400.0;
+    days as f64 + fraction
+}
+
+/// The day of the year a Modified Julian Day falls in, 1 to 366.
+pub fn day_of_year(mjd: f64) -> i32 {
+    let day = mjd.floor() as i64;
+    let (year, _, _) = civil(day);
+    (day - (days_from_civil(year, 1, 1) + MJD_UNIX) + 1) as i32
+}
+
+/// The seconds of universal time a Modified Julian Day falls at.
+pub fn seconds_of_day(mjd: f64) -> f64 {
+    (mjd - mjd.floor()) * 86400.0
+}
+
+/// The civil date of a Modified Julian Day, UTC.
+///
+/// Hinnant's algorithm, which is exact over any range an `i64` holds.
+fn civil(mjd: i64) -> (i64, u32, u32) {
+    let z = mjd - MJD_UNIX + 719_468;
+    let era = if z >= 0 { z } else { z - 146_096 } / 146_097;
+    let day_of_era = (z - era * 146_097) as u64;
+    let year_of_era =
+        (day_of_era - day_of_era / 1460 + day_of_era / 36_524 - day_of_era / 146_096) / 365;
+    let year = year_of_era as i64 + era * 400;
+    let day_of_year = day_of_era - (365 * year_of_era + year_of_era / 4 - year_of_era / 100);
+    let shifted_month = (5 * day_of_year + 2) / 153;
+    let day = (day_of_year - (153 * shifted_month + 2) / 5 + 1) as u32;
+    let month = if shifted_month < 10 {
+        shifted_month + 3
+    } else {
+        shifted_month - 9
+    } as u32;
+    (if month <= 2 { year + 1 } else { year }, month, day)
+}
+
+/// Days from 1970-01-01 to a civil date, the inverse of [`civil`].
+fn days_from_civil(year: i64, month: u32, day: u32) -> i64 {
+    let year = if month <= 2 { year - 1 } else { year };
+    let era = if year >= 0 { year } else { year - 399 } / 400;
+    let year_of_era = (year - era * 400) as u64;
+    let shifted_month = if month > 2 { month - 3 } else { month + 9 };
+    let day_of_year = (153 * u64::from(shifted_month) + 2) / 5 + u64::from(day) - 1;
+    let day_of_era = year_of_era * 365 + year_of_era / 4 - year_of_era / 100 + day_of_year;
+    era * 146_097 + day_of_era as i64 - 719_468
+}

src/lib.rs

@@ -3,6 +3,8 @@
 
 pub mod scalar;
 
+#[cfg(feature = "drivers")]
+pub mod drivers;
 mod expansion;
 mod helpers;
 mod input;

src/model.rs

@@ -130,6 +130,69 @@ impl Model {
         self.thermosphere(input, input.altitude_km, &mut work)
     }
 
+    /// Density and temperature at a date, with the drivers looked up for you.
+    ///
+    /// The flux and magnetic activity come from the compiled-in record, the day of the
+    /// year and the time of day come from `mjd`, and the local solar time is derived so
+    /// that it agrees with both. This is the whole of what a caller needs to evaluate the
+    /// model at a moment and a place, and it is [`Model::gtd7d`] underneath, so the total
+    /// mass density is the one a drag calculation wants.
+    ///
+    /// `mjd` is a Modified Julian Day; [`crate::drivers::mjd`] builds one from a civil
+    /// date, and `hifitime`'s `Epoch::to_mjd_utc_days` returns one directly. The epoch is
+    /// `f64` rather than generic because a table lookup has no derivative to carry;
+    /// position stays generic, so a drag partial still comes out.
+    ///
+    /// ```
+    /// use nrlmsise00::{drivers, Model};
+    ///
+    /// let model = Model::new();
+    /// let when = drivers::mjd(2020, 6, 21, 8, 3, 20.0);
+    /// let output = model.gtd7d_at(when, 400.0, 60.0, -70.0).expect("inside the record");
+    /// assert!(output.total_mass > 0.0);
+    /// ```
+    ///
+    /// # Errors
+    ///
+    /// If the date is outside what the record can be evaluated on, which includes every
+    /// date in the future: the drivers are measurements, and there are none for later.
+    #[cfg(feature = "drivers")]
+    pub fn gtd7d_at<T: Scalar>(
+        &self,
+        mjd: f64,
+        altitude_km: T,
+        latitude_deg: T,
+        longitude_deg: T,
+    ) -> Result<Output<T>, crate::drivers::DriverError> {
+        let drivers = crate::drivers::Drivers::embedded().at(mjd)?;
+        let seconds_of_day = T::lift(crate::drivers::seconds_of_day(mjd));
+        Ok(self.gtd7d(&Input {
+            day_of_year: crate::drivers::day_of_year(mjd),
+            seconds_of_day,
+            altitude_km,
+            latitude_deg,
+            longitude_deg,
+            local_solar_time_h: seconds_of_day / 3600.0 + longitude_deg / 15.0,
+            f107_average: T::lift(drivers.f107_average),
+            f107: T::lift(drivers.f107),
+            magnetic: match drivers.magnetic {
+                crate::input::MagneticActivity::History { ap } => {
+                    crate::input::MagneticActivity::History {
+                        ap: ap.map(T::lift),
+                    }
+                }
+                crate::input::MagneticActivity::Daily { ap } => {
+                    crate::input::MagneticActivity::Daily { ap: T::lift(ap) }
+                }
+                crate::input::MagneticActivity::Off { daily_ap } => {
+                    crate::input::MagneticActivity::Off {
+                        daily_ap: T::lift(daily_ap),
+                    }
+                }
+            },
+        }))
+    }
+
     /// `gtd7`: the thermosphere, and below 72.5 km the mesosphere blended onto it.
     fn evaluate<T: Scalar>(&self, input: &Input<T>, work: &mut Work<T>) -> Output<T> {
         let (sw, metric) = (&self.resolved.sw, self.switches.units == Units::Si);

tests/drivers.rs

@@ -0,0 +1,340 @@
+//! The compiled-in driver record, and the lookups into it.
+//!
+//! The record is measurements, so these tests check it against what is known about the
+//! Sun's recent behavior rather than against a second copy of the same numbers. Where a
+//! value is asserted, it is one that can be looked up independently: the solar cycle's
+//! maxima and minima, the flux spike of the 2011 flare, the Gannon storm of May 2024.
+
+#![cfg(feature = "drivers")]
+
+use nrlmsise00::drivers::{self, DriverError, Drivers};
+use nrlmsise00::{MagneticActivity, Model};
+
+#[test]
+fn test_the_record_covers_what_it_says_it_does() {
+    let record = Drivers::embedded();
+
+    // It starts on 1957-01-02: geomagnetic indices reach back to 1932 and the flux to
+    // 1947, but the early flux record is too sporadic to interpolate, and the first
+    // measured day of 1957 is the second.
+    assert_eq!(
+        record.first_day(),
+        drivers::mjd(1957, 1, 2, 0, 0, 0.0) as i64
+    );
+    assert!(record.len() > 25_000, "{} days", record.len());
+    assert!(!record.is_empty());
+
+    // The usable end is forty days before the measured end, because the flux average the
+    // model wants is centered.
+    assert_eq!(
+        record.last_usable_day(),
+        record.last_day() - drivers::FLUX_AVERAGE_RADIUS
+    );
+    // And the front of the record is short by the same forty days, for the same reason.
+    assert_eq!(
+        record.first_usable_day(),
+        record.first_day() + drivers::FLUX_AVERAGE_RADIUS
+    );
+    assert!(record.day(record.first_day() as f64).is_err());
+    assert!(record.day(record.first_usable_day() as f64).is_ok());
+    assert!(record.covers(drivers::mjd(2020, 6, 21, 12, 0, 0.0)));
+    assert!(!record.covers(drivers::mjd(1950, 1, 1, 0, 0, 0.0)));
+
+    // Every day between the ends is present, which is what being indexed by day requires.
+    for day in record.first_usable_day()..=record.last_usable_day() {
+        assert!(record.covers(day as f64), "MJD {day} is missing");
+    }
+}
+
+#[test]
+fn test_a_day_outside_the_record_is_refused_and_says_why() {
+    let record = Drivers::embedded();
+
+    let early = drivers::mjd(1900, 1, 1, 0, 0, 0.0);
+    match record.day(early) {
+        Err(DriverError::BeforeRecord { asked, first }) => {
+            assert_eq!(asked, early);
+            assert_eq!(first, record.first_usable_day());
+        }
+        other => panic!("expected a refusal, got {other:?}"),
+    }
+
+    // The future, which is the case that matters: drivers are measurements.
+    let late = (record.last_day() + 1) as f64;
+    match record.at(late) {
+        Err(error @ DriverError::PastRecord { .. }) => {
+            let message = error.to_string();
+            assert!(message.contains("measurements"), "{message}");
+            assert!(message.contains("fetch-drivers.py"), "{message}");
+        }
+        other => panic!("expected a refusal, got {other:?}"),
+    }
+
+    // And the forty days at the end that have a measurement but no centered average.
+    assert!(record.day((record.last_day() - 1) as f64).is_err());
+    assert!(record.day(record.last_usable_day() as f64).is_ok());
+}
+
+/// The solar cycle is in the record, at the right amplitude and in the right years.
+#[test]
+fn test_the_solar_cycle_is_where_it_should_be() {
+    let record = Drivers::embedded();
+    let flux_in = |year: i64| -> f64 {
+        let (mut sum, mut count) = (0.0, 0);
+        for day in 1..=365 {
+            let when = drivers::mjd(year, 1, 1, 12, 0, 0.0) + f64::from(day - 1);
+            if let Ok(day) = record.day(when) {
+                sum += day.flux;
+                count += 1;
+            }
+        }
+        sum / f64::from(count)
+    };
+
+    // Cycle 23 peaked around 2001-2002 and bottomed out in 2008-2009, which was the
+    // deepest minimum of the space age. Cycle 24 peaked lower, in 2014.
+    let minimum = flux_in(2008);
+    let peak_23 = flux_in(2001);
+    let peak_24 = flux_in(2014);
+    assert!((65.0..75.0).contains(&minimum), "2008 averaged {minimum}");
+    assert!((160.0..200.0).contains(&peak_23), "2001 averaged {peak_23}");
+    assert!((110.0..150.0).contains(&peak_24), "2014 averaged {peak_24}");
+    assert!(peak_23 > peak_24, "cycle 24 was the weaker of the two");
+
+    // The whole record sits inside the range the provenance states, and the two days
+    // that fall below 60 are the two named outliers. Both are isolated dropouts between
+    // neighbors around 70, which is instrument trouble rather than the Sun, and both are
+    // passed through unaltered because the source carries them. If a refresh changes
+    // this, it is worth looking at before the assertion is relaxed.
+    let mut dropouts = Vec::new();
+    for day in record.first_usable_day()..=record.last_usable_day() {
+        let day_record = record.day(day as f64).expect("inside the record");
+        let flux = day_record.flux;
+        assert!((53.0..1000.0).contains(&flux), "MJD {day} has flux {flux}");
+        if flux < 60.0 && day_record.flux_measured {
+            dropouts.push((day, flux));
+        }
+    }
+    assert_eq!(
+        dropouts,
+        vec![
+            (drivers::mjd(1997, 3, 1, 0, 0, 0.0) as i64, 53.5),
+            (drivers::mjd(1998, 3, 2, 0, 0, 0.0) as i64, 56.8),
+        ],
+        "the measured days below 60 sfu"
+    );
+}
+
+/// The flare of 2011-03-07 put the daily flux above 900, which is four times the
+/// solar-maximum background and looks like a parsing error until it is checked.
+#[test]
+fn test_the_2011_flux_spike_is_in_the_record() {
+    let record = Drivers::embedded();
+    let flare = record
+        .day(drivers::mjd(2011, 3, 7, 12, 0, 0.0))
+        .expect("inside the record");
+    assert!(flare.flux > 900.0, "2011-03-07 has flux {}", flare.flux);
+    assert!(flare.flux_measured);
+}
+
+/// The Gannon storm of May 2024 was the largest geomagnetic storm in twenty years.
+#[test]
+fn test_the_gannon_storm_is_in_the_record() {
+    let record = Drivers::embedded();
+    let onset = record
+        .day(drivers::mjd(2024, 5, 10, 12, 0, 0.0))
+        .expect("inside the record");
+    let peak = record
+        .day(drivers::mjd(2024, 5, 11, 12, 0, 0.0))
+        .expect("inside the record");
+
+    // The storm began late on the 10th UTC and the Kp 9 intervals fall on the 11th, which
+    // is why the peak day is the later one. ap reaches 400 at the top of the Kp scale.
+    assert_eq!(onset.ap.iter().cloned().fold(0.0, f64::max), 300.0);
+    assert_eq!(
+        peak.ap.iter().cloned().fold(0.0, f64::max),
+        400.0,
+        "the storm reached the top of the scale on the 11th"
+    );
+    assert!(peak.ap_daily > onset.ap_daily);
+    assert!(peak.ap_daily > 250.0, "daily Ap was {}", peak.ap_daily);
+
+    // A quiet day for contrast.
+    let quiet = record
+        .day(drivers::mjd(2019, 9, 15, 12, 0, 0.0))
+        .expect("inside the record");
+    assert!(quiet.ap_daily < 10.0, "daily Ap was {}", quiet.ap_daily);
+
+    // The storm is denser than the quiet day at the same altitude and place, which is the
+    // whole reason the model reads these.
+    let model = Model::new();
+    let density = |when: f64| {
+        model
+            .gtd7d_at(when, 400.0, 45.0, 0.0)
+            .expect("inside the record")
+            .total_mass
+    };
+    let stormy = density(drivers::mjd(2024, 5, 11, 12, 0, 0.0));
+    let calm = density(drivers::mjd(2019, 9, 15, 12, 0, 0.0));
+    assert!(
+        stormy > 5.0 * calm,
+        "storm {stormy:e} against quiet {calm:e}"
+    );
+}
+
+/// The centered average is the mean of the eighty-one days around the day, and the daily
+/// flux moves around it rather than with it.
+#[test]
+fn test_the_centered_average_is_centered() {
+    let record = Drivers::embedded();
+    let when = drivers::mjd(2015, 3, 15, 0, 0, 0.0);
+
+    let mut sum = 0.0;
+    for offset in -drivers::FLUX_AVERAGE_RADIUS..=drivers::FLUX_AVERAGE_RADIUS {
+        sum += record
+            .day(when + offset as f64)
+            .expect("inside the record")
+            .flux;
+    }
+    let expected = sum / (2.0 * drivers::FLUX_AVERAGE_RADIUS as f64 + 1.0);
+    let average = record.centered_flux(when).expect("inside the record");
+    assert!(
+        (average - expected).abs() < 1e-9,
+        "{average} against {expected}"
+    );
+}
+
+/// The seven magnetic values come out in the order the model wants them.
+#[test]
+fn test_the_magnetic_history_is_in_the_models_order() {
+    let record = Drivers::embedded();
+    // 1500 UT is in the sixth three-hour interval of the day, which runs 15:00 to 18:00.
+    let when = drivers::mjd(2024, 5, 10, 15, 30, 0.0);
+    let set = record.at(when).expect("inside the record");
+    let today = record.day(when).expect("inside the record");
+    let yesterday = record.day(when - 1.0).expect("inside the record");
+
+    let MagneticActivity::History { ap } = set.magnetic else {
+        panic!("a lookup should give the history, which is what a propagation wants");
+    };
+
+    assert_eq!(ap[0], today.ap_daily, "the daily mean comes first");
+    assert_eq!(ap[1], today.ap[5], "then the interval the moment falls in");
+    assert_eq!(ap[2], today.ap[4], "three hours before");
+    assert_eq!(ap[3], today.ap[3], "six hours before");
+    assert_eq!(ap[4], today.ap[2], "nine hours before");
+
+    // Then the mean of the eight intervals from 12 to 33 hours before, which crosses into
+    // the previous day: intervals 1 and 0 of today, then 7 down to 2 of yesterday.
+    let expected: f64 = (today.ap[1]
+        + today.ap[0]
+        + yesterday.ap[7..].iter().sum::<f64>()
+        + yesterday.ap[2..7].iter().sum::<f64>())
+        / 8.0;
+    assert!(
+        (ap[5] - expected).abs() < 1e-12,
+        "{} against {expected}",
+        ap[5]
+    );
+
+    // The flux is the previous day's, as the model asks.
+    assert_eq!(set.f107, yesterday.flux);
+}
+
+/// The dates round-trip, including across a leap year and a century.
+#[test]
+fn test_the_date_arithmetic_round_trips() {
+    for (year, month, day) in [
+        (1957, 1, 2),
+        (1970, 1, 1),
+        (2000, 2, 29),
+        (2024, 2, 29),
+        (2026, 10, 1),
+        (1999, 12, 31),
+    ] {
+        let when = drivers::mjd(year, month, day, 0, 0, 0.0);
+        assert_eq!(when.fract(), 0.0, "{year}-{month}-{day} is not a whole day");
+        // Day of the year, checked against the month lengths.
+        let cumulative: [i32; 12] = [0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334];
+        let leap = (year % 4 == 0 && year % 100 != 0) || year % 400 == 0;
+        let mut expected = cumulative[(month - 1) as usize] + day as i32;
+        if leap && month > 2 {
+            expected += 1;
+        }
+        assert_eq!(drivers::day_of_year(when), expected, "{year}-{month}-{day}");
+    }
+
+    // The time of day comes back out in seconds.
+    let noon = drivers::mjd(2020, 6, 21, 12, 0, 0.0);
+    assert_eq!(drivers::seconds_of_day(noon), 43200.0);
+    let odd = drivers::mjd(2020, 6, 21, 8, 3, 20.0);
+    assert!((drivers::seconds_of_day(odd) - 29000.0).abs() < 1e-6);
+}
+
+/// The convenience that takes a date gives the same answer as assembling the input by hand.
+#[test]
+fn test_the_date_entry_point_matches_assembling_the_input() {
+    use nrlmsise00::Input;
+
+    let model = Model::new();
+    let when = drivers::mjd(2020, 6, 21, 8, 3, 20.0);
+    let set = Drivers::embedded().at(when).expect("inside the record");
+    let seconds = drivers::seconds_of_day(when);
+
+    let by_hand = model.gtd7d(&Input {
+        day_of_year: drivers::day_of_year(when),
+        seconds_of_day: seconds,
+        altitude_km: 400.0,
+        latitude_deg: 60.0,
+        longitude_deg: -70.0,
+        local_solar_time_h: seconds / 3600.0 + -70.0 / 15.0,
+        f107_average: set.f107_average,
+        f107: set.f107,
+        magnetic: set.magnetic,
+    });
+    let by_date = model
+        .gtd7d_at(when, 400.0, 60.0, -70.0)
+        .expect("inside the record");
+    assert_eq!(by_hand.densities(), by_date.densities());
+    assert_eq!(by_hand.temperature, by_date.temperature);
+}
+
+/// A date past the record is refused by the convenience too, rather than extrapolated.
+#[test]
+fn test_the_date_entry_point_refuses_the_future() {
+    let model = Model::new();
+    let far = drivers::mjd(2100, 1, 1, 0, 0, 0.0);
+    assert!(matches!(
+        model.gtd7d_at(far, 400.0, 60.0, -70.0),
+        Err(DriverError::PastRecord { .. })
+    ));
+}
+
+/// Density through the solar cycle, from the record rather than from chosen numbers: the
+/// thermosphere at 400 km is several times denser at solar maximum than at minimum.
+#[test]
+fn test_the_cycle_moves_the_density_the_way_it_should() {
+    let model = Model::new();
+    let at = |year: i64| {
+        model
+            .gtd7d_at(drivers::mjd(year, 6, 21, 12, 0, 0.0), 400.0, 0.0, 0.0)
+            .expect("inside the record")
+            .total_mass
+    };
+    let minimum = at(2008);
+    let maximum = at(2001);
+    assert!(
+        maximum > 4.0 * minimum,
+        "solar maximum {maximum:e} against minimum {minimum:e}"
+    );
+}
+
+/// The provenance says where the record came from, which a program may want to print.
+#[test]
+fn test_the_provenance_is_compiled_in() {
+    let note = drivers::provenance();
+    assert!(note.contains("kp.gfz-potsdam.de"), "no source");
+    assert!(note.contains("CC BY 4.0"), "no license");
+    assert!(note.contains("sha256"), "no digest");
+    assert!(note.contains("fetch-drivers.py"), "no way to regenerate it");
+}