prop/docs/plans/2026-04-15-merra2-historical-backfill.md

NARR Historical Backfill Implementation Plan

For Claude: REQUIRED SUB-SKILL: Use superpowers:executing-plans to implement this plan task-by-task.

History: This file originally specified a NASA MERRA-2 path. While executing Phase 1 it turned out the MERRA-2 AWS bucket (gesdisc-cumulus-prod-protected) requires Earthdata credentials, contradicting the initial research. While investigating alternatives I discovered NCEP NARR is anonymously available at NCEI with native GRIB1 + wgrib1-format .inv files. Plan was rewritten around NARR. A second correction landed during the spike: wgrib2 cannot decode GRIB1 messages (rd_grib2_seq_msg, grib1 message ignored, use wgrib), so we do need a separate decoder. cdo handles GRIB1 + Lambert-conformal projection

• nearest-neighbor remap natively, so the plan now uses cdo for decoding. Filename kept for git continuity.

Goal: Replace the broken ERA5/CDS historical backfill with NCEP NARR fetched anonymously from NCEI, parsed via the existing wgrib2 pipeline, and inserted into the existing era5_profiles table.

Architecture: New Microwaveprop.Weather.NarrClient mirrors the existing HrrrClient byte-range fetch pattern: download the small .inv file, parse it into a map of {var, level} => {byte_offset, length}, issue HTTP Range: bytes=START-END GETs for the ~10 records we need, write each record to a per-record temp file, then run cdo -merge to combine them into a single GRIB1 and cdo -outputtab,name,lev,value -remapnn,lon=X_lat=Y to extract values at the QSO point. New NarrFetchWorker is a single-stage Oban worker (no submit/poll split — there's no upstream queue) that takes a contact_id, picks the nearest 3-hourly analysis hour to the QSO timestamp, calls NarrClient, and inserts an era5_profiles row. ContactWeatherEnqueueWorker dispatches NarrFetchWorker for pre-2014 contacts in place of Era5FetchWorker. The era5_profiles table is reused 1:1 — schema rename can be a follow-up if we feel like it.

Tech Stack: Elixir, Ecto, Oban, Req for HTTPS byte-range fetches, cdo (Climate Data Operators) CLI for GRIB1 record handling and Lambert-conformal nearest-neighbor extraction. The existing wgrib2 CLI does not decode GRIB1 messages so we add cdo to the Dockerfile (one new system dep, ~30–60 MB image growth via libnetcdf/libhdf5/libudunits2 transitive deps).

cdo parameter-number → semantic-name lookup (verified empirically against NARR 2010-06-15 12Z):

cdo name	NCEP GRIB1 param	Meaning	Where in NARR
var1	1	PRES (Pressure)	sfc / 10 m / hybrid
var7	7	HGT (Geopotential height)	pressure levels
var11	11	TMP (Temperature)	sfc / 2 m / 10 m / pressure levels
var17	17	DPT (Dew point)	2 m only
var33	33	UGRD (u wind)	10 m / pressure levels
var34	34	VGRD (v wind)	10 m / pressure levels
var51	51	SPFH (Specific humidity)	2 m / pressure levels
var54	54	PWAT (Precipitable water)	atmos col
var221	221	HPBL (PBL height)	sfc

cdo's lev column for pressure-level records is in Pa (85000 = 850 mb). For 2 m above ground the lev is 2, for 10 m the lev is 10, for surface/atmos col the lev is 0. The parser must use both the var ID and the lev value to disambiguate (e.g. var11 lev=2 is 2 m TMP, var11 lev=85000 is 850 mb TMP).

Why NARR over MERRA-2 / ERA5:

• Native GRIB1 with wgrib1-format .inv (already verified — narr-a_221_20100615_1200_000.inv parses cleanly)
• Accept-Ranges: bytes confirmed on https://www.ncei.noaa.gov/data/north-american-regional-reanalysis/access/3-hourly/... — same byte-range fetch pattern HrrrClient already uses
• Anonymous HTTPS — no auth, no quota, no job queue
• 32 km Lambert conformal grid covering CONUS (Grid 221), 3-hourly, 1979-01-01 → 2014-10-02 — exactly the pre-HRRR window
• DPT is a direct surface dewpoint record, no specific-humidity → Td derivation
• 56 MB per analysis hour total, but a per-point fetch only pulls ~10 records ≈ 1.5 MB max
• Same atmospheric variables we need on the same surface + pressure-level layout as the existing era5_profiles schema

Spike findings (2026-04-15, verified by hand):

• For 2010-06-15 12Z at DFW (32.9°N, -97.0°W): TMP_2m=296.94 K (23.8°C), DPT_2m=294.49 K (21.3°C), HPBL=703.5 m, PRES_sfc=99383.6 Pa (993.8 mb), PWAT=45.1 mm, TMP_850mb=291.52 K (18.4°C), HGT_850mb=1547.8 m, SPFH_850mb=10.4 g/kg — all physically sane Texas summer values
• cdo -merge is the right operator for combining GRIB1 records of different variables. cdo -cat is for time concatenation only and silently drops records with mismatched levels. Shell cat produces an invalid composite (GRIB message error). Always use cdo -merge.

Verified facts (gathered live during plan rewrite):

• File URL pattern: https://www.ncei.noaa.gov/data/north-american-regional-reanalysis/access/3-hourly/{YYYYMM}/{YYYYMMDD}/narr-a_221_{YYYYMMDD}_{HHMM}_000.grb
• Inventory file: same path with .inv extension
• HHMM ∈ {0000, 0300, 0600, 0900, 1200, 1500, 1800, 2100} (3-hourly analyses)
• Inventory format: record:offset:D=YYYYMMDDHH:VAR:level:kpds=...:anl:...:"description [units]
• Records we need (verified present in real inventory for 2010-06-15 12Z):
  • 15: HPBL:sfc
  • 265: PRES:sfc
  • 288: TMP:2 m above gnd
  • 290: DPT:2 m above gnd
  • 293: UGRD:10 m above gnd
  • 294: VGRD:10 m above gnd
  • 317: PWAT:atmos col
  • Pressure levels 100, 125, 150, 175, 200, 225, 250, 275, 300, 350, 400, 450, 500, 550, 600, 650, 700, 725, 750, 775, 800, 825, 850, 875, 900, 925, 950, 975, 1000 mb each have TMP, HGT, SPFH records

---

Phase 1: Spike — DONE 2026-04-15

Verified by hand. Findings recorded in the spike-findings block above. Key deltas from the original assumption:

• wgrib2 cannot read GRIB1 messages — it errors out with rd_grib2_seq_msg, grib1 message ignored, use wgrib.
• cdo handles GRIB1 + Lambert-conformal projection + nearest-neighbor remap natively. cdo -merge for combining records, cdo -outputtab,name,lev,value -remapnn,lon=X_lat=Y for point extraction.
• cdo emits variable names as varNNN (NCEP GRIB1 parameter numbers). Mapping table in the architecture section.
• All 5 surface vars + 3 pressure-level vars at 850 mb extract to sane physical values for DFW June 12Z.
• Adding cdo to the Dockerfile is now a Phase 2 prerequisite (see Task 2.5).

---

Phase 2: Build NarrClient (TDD)

Task 2.5 (PREREQUISITE): Add cdo to Dockerfile runtime stage

Files:

• Modify: Dockerfile:137-141 (final runtime apt-get install)

Add cdo to the apt-get install line. cdo on Debian trixie pulls in libnetcdf, libhdf5, libudunits2 transitively. Verify with docker build --target final then docker run --rm <image> cdo --version.

RUN apt-get update \
  && apt-get install -y --no-install-recommends \
     libstdc++6 openssl libncurses6 locales ca-certificates snmp \
     libgfortran5 libaec0 zlib1g libpng16-16t64 libopenjp2-7 \
     cdo \
  && rm -rf /var/lib/apt/lists/*

This is a prerequisite because every other Phase 2/3 task depends on cdo being available at runtime. Local dev already has cdo from brew install cdo during the spike — no local change needed there. Commit this Dockerfile change before the NarrClient tasks so prod CI builds with the new image.

Task 2: NarrClient.url_for/2 + url_for_inventory/2

Files:

• Create: lib/microwaveprop/weather/narr_client.ex
• Create: test/microwaveprop/weather/narr_client_test.exs

Step 1: Failing test

defmodule Microwaveprop.Weather.NarrClientTest do
  use ExUnit.Case, async: true
  alias Microwaveprop.Weather.NarrClient

  describe "url_for/1" do
    test "builds the NCEI NARR HTTPS URL for a given analysis time" do
      vt = ~U[2010-06-15 12:00:00Z]
      assert NarrClient.url_for(vt) ==
        "https://www.ncei.noaa.gov/data/north-american-regional-reanalysis/access/3-hourly/201006/20100615/narr-a_221_20100615_1200_000.grb"
    end

    test "raises on non-3-hourly analysis times (NARR is 00/03/06/.../21Z only)" do
      vt = ~U[2010-06-15 13:00:00Z]
      assert_raise ArgumentError, fn -> NarrClient.url_for(vt) end
    end
  end

  describe "url_for_inventory/1" do
    test "swaps .grb for .inv" do
      vt = ~U[2010-06-15 12:00:00Z]
      assert NarrClient.url_for_inventory(vt) ==
        "https://www.ncei.noaa.gov/data/north-american-regional-reanalysis/access/3-hourly/201006/20100615/narr-a_221_20100615_1200_000.inv"
    end
  end

  describe "snap_to_analysis_hour/1" do
    test "rounds an arbitrary timestamp down to the nearest 3-hourly NARR slot" do
      assert NarrClient.snap_to_analysis_hour(~U[2010-06-15 13:42:00Z]) == ~U[2010-06-15 12:00:00Z]
      assert NarrClient.snap_to_analysis_hour(~U[2010-06-15 14:30:00Z]) == ~U[2010-06-15 12:00:00Z]
      assert NarrClient.snap_to_analysis_hour(~U[2010-06-15 15:00:00Z]) == ~U[2010-06-15 15:00:00Z]
      assert NarrClient.snap_to_analysis_hour(~U[2010-06-15 23:59:59Z]) == ~U[2010-06-15 21:00:00Z]
    end
  end
end

Step 2: Run, see fail.

Step 3: Implement minimally.

Step 4: Run, green.

Step 5: Commit

Task 3: NarrClient.fetch_inventory/1

Files: modify lib/microwaveprop/weather/narr_client.ex + test file.

Step 1: Failing test using a fixture inventory file.

Capture a real NARR .inv file from Phase 1 and commit it as a fixture:

cp /tmp/narr.inv test/fixtures/narr/narr-a_221_20100615_1200_000.inv

Then test the parser:

describe "parse_inventory/2" do
  test "parses the inv text into a map of {var, level} => {offset, length}" do
    inv = File.read!("test/fixtures/narr/narr-a_221_20100615_1200_000.inv")
    file_size = 56_221_430

    {:ok, index} = NarrClient.parse_inventory(inv, file_size)

    # Spot-check a surface record
    assert {37_836_396, _len_1} = index[{"TMP", "2 m above gnd"}]
    # Spot-check a pressure-level record
    assert {_off, _len} = index[{"TMP", "850 mb"}]
    # File-end record's length is computed against file_size
    assert {_off, len_last} = index[{"PWAT", "atmos col"}]
    assert len_last > 0
  end
end

Step 2: Run, see fail.

Step 3: Implement.

@spec parse_inventory(String.t(), pos_integer()) :: {:ok, %{{String.t(), String.t()} => {non_neg_integer(), pos_integer()}}}
def parse_inventory(inv_text, file_size) do
  lines =
    inv_text
    |> String.split("\n", trim: true)
    |> Enum.map(&parse_line/1)
    |> Enum.reject(&is_nil/1)

  index =
    lines
    |> Enum.zip(Enum.drop(lines, 1) ++ [{nil, file_size, nil, nil}])
    |> Map.new(fn {{_n, off, var, level}, {_, next_off, _, _}} ->
      {{var, level}, {off, next_off - off}}
    end)

  {:ok, index}
end

defp parse_line(line) do
  case String.split(line, ":") do
    [num, off, _date, var, level | _] ->
      {String.to_integer(num), String.to_integer(off), var, level}
    _ ->
      nil
  end
end

Step 4: Run, green.

Step 5: Add fetch_inventory/1 that wraps parse_inventory with an Req.get of the .inv URL — this part needs a Req.Test.stub test to avoid hitting NCEI in the test suite.

Step 6: Commit.

Task 4: NarrClient.fetch_profile_at/2

Files: modify narr_client.ex + test.

This is the big one. Given a valid_time and a {lat, lon} point, fetch all the records we need via byte-range, write each record to its own per-record temp file, run cdo -merge to combine them, run cdo -outputtab,name,lev,value -remapnn,lon=X_lat=Y to extract values at the QSO point, and return a %{surface_temp_c: ..., profile: [...], ...} map shaped exactly like what Era5BatchClient.build_profile_attrs/5 produces.

Step 1: Failing test using a real captured composite GRIB1 fixture.

Capture a tiny composite GRIB1 (~1 MB total) containing all needed records for one analysis hour from the spike, commit as test/fixtures/narr/narr_dfw_2010-06-15_12z.grb. The test reads the fixture, calls NarrClient.extract_profile_from_file/3 (a public-but-internal helper exposed for testability), and asserts the returned map has all expected keys with sane values.

test "extract_profile_from_file/3 builds a profile from a captured fixture" do
  profile = NarrClient.extract_profile_from_file(
    "test/fixtures/narr/narr_dfw_2010-06-15_12z.grb",
    32.9,
    -97.0
  )

  # Surface — values come from the spike (DFW 2010-06-15 12Z)
  assert_in_delta profile.surface_temp_c, 23.8, 1.0
  assert_in_delta profile.surface_dewpoint_c, 21.3, 1.0
  assert profile.surface_pressure_mb > 950 and profile.surface_pressure_mb < 1050
  assert profile.hpbl_m >= 0 and profile.hpbl_m < 5000
  assert profile.pwat_mm > 30 and profile.pwat_mm < 60

  # Profile array
  assert length(profile.profile) >= 5
  Enum.each(profile.profile, fn level ->
    assert is_float(level["pres"])
    assert is_float(level["tmpc"])
    assert is_float(level["hght"])
  end)
end

Step 2: Run, see fail.

Step 3: Implement.

Two layers:

• extract_profile_from_file/3 — testable, takes a composite GRIB1 path + lat + lon, runs cdo -outputtab,name,lev,value -remapnn,lon=#{lon}_lat=#{lat} <path>, parses the text output, maps varNNN+lev to semantic field names via the lookup table in the Architecture section, builds the era5_profiles-shaped map. Reuses ThetaE.dewpoint_from_spfh/2 for pressure-level dewpoint derivation (NARR has DPT only at 2 m).
• fetch_profile_at/2 — non-testable (does network), calls fetch_inventory/1, picks the records to fetch, byte-range-fetches each into its own /tmp/narr_<tag>_<i>.grb file, calls cdo -merge to combine them into /tmp/narr_<tag>_merged.grb, calls extract_profile_from_file/3, then deletes all temp files.

cdo invocation pattern (verified in spike):

cdo -merge /tmp/r1.grb /tmp/r2.grb /tmp/r3.grb /tmp/merged.grb
cdo -outputtab,name,lev,value -remapnn,lon=-97.0_lat=32.9 /tmp/merged.grb

Output format (one line per record):

#    name    lev    value 
   var11      2  296.939 
   var17      2 294.4932 
  var221      0  703.519 
    var1      0 99383.64 
   var54      0 45.10547 
   var11  85000 291.5186 
    var7  85000 1547.793 
   var51  85000 0.01040225

Parse: skip the header line (starts with #), split each data line on whitespace, ignore empty fields, take [name, lev, value], look up name in the parameter table, build the row.

Step 4: Run, green.

Step 5: Commit.

---

Phase 3: NarrFetchWorker (TDD)

Task 5: Build NarrFetchWorker

Files:

• Create: lib/microwaveprop/workers/narr_fetch_worker.ex
• Create: test/microwaveprop/workers/narr_fetch_worker_test.exs

Per-contact Oban worker. Args: %{contact_id: id}. Looks up the contact, picks the nearest 3-hourly analysis hour to qso_timestamp, calls NarrClient.fetch_profile_at(snapped_vt, {contact.lat, contact.lon}), inserts an era5_profiles row via Era5BatchClient.bulk_insert_profiles/1, marks the contact as enriched.

Step 1: Failing test with Req.Test.stub for the network calls. Asserts the era5_profiles row is created and the contact's status is updated.

Step 2: Implement.

Step 3: Commit.

Task 6: Add :narr queue to dev.exs and runtime.exs

Files:

• Modify: config/dev.exs (queue list)
• Modify: config/runtime.exs (queue list)

narr: 6 in both. Network-bound, so high concurrency is fine.

Task 7: Wire ContactWeatherEnqueueWorker to NarrFetchWorker

Files:

• Modify: lib/microwaveprop/workers/contact_weather_enqueue_worker.ex (Era5FetchWorker.new(%{...}) call site, currently around line 196)
• Modify: matching test file

Step 1: Update the failing test to assert NarrFetchWorker is enqueued instead of Era5FetchWorker.

Step 2: Implement the swap. The pre-2014 cutoff stays the same — anything with qso_timestamp < 2014-10-02 00:00 UTC gets NarrFetchWorker, after that gets HrrrFetchWorker.

Step 3: Run all tests.

Step 4: Commit.

Task 8: End-to-end smoke test in dev

Files: none — manual verification.

Step 1: Pick one pre-2014 contact in dev with hrrr_status == :unavailable, manually enqueue NarrFetchWorker for it, and verify a row lands in era5_profiles and the propagation map shows weather data for that historical timestamp.

---

Phase 4: Decommission ERA5

Task 9: Permanently document Era5* deprecation in runtime.exs

Files:

• Modify: config/runtime.exs (era5* paused-queue comment block)

Update the existing "TEMPORARILY PAUSED while we diagnose" comment to "DEPRECATED in favor of NARR via NarrFetchWorker — see plan 2026-04-15. Slated for deletion in a follow-up PR."

Task 10 (follow-up PR, do not execute in this pass): Drop Era5* code paths

This is a separate PR after NARR has been running cleanly for a few days. Scope:

• Drop the era5_cds_jobs table (Ecto migration)
• Delete Era5SubmitWorker, Era5PollWorker, Era5FetchWorker, Era5MonthBatchWorker, Era5BatchClient, Era5Client, Era5CdsJob schema, and the era5_* queue config blocks
• Maybe rename era5_profiles → historical_profiles (or leave it — the misnomer is a footgun but harmless)

---

Open Questions / Risks

1. wgrib2 -lola on the Lambert conformal grid. Resolved by spike: wgrib2 doesn't read GRIB1 at all, so we use cdo's -remapnn,lon=X_lat=Y instead. cdo handles the Lambert conformal projection natively.

2. NARR records sometimes have different metadata strings across years. The 1979 records may have slightly different inventory strings than 2014 records (kpds tables changed once or twice). The Phase 1 spike used 2010 — also test 1995 and 2013 in Phase 3 Task 8 to confirm the parser handles the full date range.

3. NARR is 3-hourly, ERA5 was hourly. That's a 3× temporal coarsening, but the propagation algorithm doesn't care — it picks the nearest sample to a QSO timestamp anyway. Acceptable.

4. No tile-month batching. ERA5 had it because CDS jobs were month-tile shaped. NARR fetches are per-point, much smaller, much faster — no batching needed. The flow is simpler: enqueue worker → worker calls client → client byte-range-fetches → done.

5. NARR file 56 MB total per analysis hour, but byte-range fetches only pull ~1.5 MB max for our 10 records. NCEI's HTTP server supports ranges (verified). If for some reason a range request fails, fall back to the full file download — but expect this to be rare.

6. narr-a vs narr-b. narr-a is the analysis (what we want). narr-b is a 3-hour forecast interpolation (filling in the off-hour times). Use narr-a exclusively — analyses are higher quality than forecast interpolations.