From a5c7dce147bce5043de3b4c067a20993745816a8 Mon Sep 17 00:00:00 2001
From: Graham McIntire <graham@mcintire.me>
Date: Sat, 18 Apr 2026 14:28:43 -0500
Subject: [PATCH] feat(gefs): scaffold extended-horizon forecast ingestion
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

Lays the groundwork for a Day 2-7 propagation outlook driven by NCEP
GEFS ensemble-mean output, picking up where HRRR's 18-hour horizon
leaves off.

- GefsClient: NOMADS URL builder, forecast-hour list (3h to +240,
  6h to +384), ensemble-mean message inventory, and a Magnus
  dewpoint derivation (pgrb2a publishes RH rather than Td)
- gefs_profiles table + GefsProfile schema mirroring hrrr_profiles
  so the scorer can run over rows from either source
- Weather.upsert_gefs_profile/1 and upsert_gefs_profiles_batch/1

No worker or UI yet — those land in follow-up commits.
---
 lib/microwaveprop/weather.ex                  |  41 ++++++
 lib/microwaveprop/weather/gefs_client.ex      | 112 ++++++++++++++++
 lib/microwaveprop/weather/gefs_profile.ex     |  60 +++++++++
 .../20260418191400_create_gefs_profiles.exs   |  33 +++++
 .../weather/gefs_client_test.exs              | 123 ++++++++++++++++++
 .../weather/gefs_profile_test.exs             |  99 ++++++++++++++
 6 files changed, 468 insertions(+)
 create mode 100644 lib/microwaveprop/weather/gefs_client.ex
 create mode 100644 lib/microwaveprop/weather/gefs_profile.ex
 create mode 100644 priv/repo/migrations/20260418191400_create_gefs_profiles.exs
 create mode 100644 test/microwaveprop/weather/gefs_client_test.exs
 create mode 100644 test/microwaveprop/weather/gefs_profile_test.exs

diff --git a/lib/microwaveprop/weather.ex b/lib/microwaveprop/weather.ex
index fb1f4431..cd005add 100644
--- a/lib/microwaveprop/weather.ex
+++ b/lib/microwaveprop/weather.ex
@@ -5,6 +5,7 @@ defmodule Microwaveprop.Weather do
 
   alias Microwaveprop.Propagation.ProfilesFile
   alias Microwaveprop.Repo
+  alias Microwaveprop.Weather.GefsProfile
   alias Microwaveprop.Weather.GridCache
   alias Microwaveprop.Weather.HrrrNativeProfile
   alias Microwaveprop.Weather.HrrrProfile
@@ -727,6 +728,46 @@ defmodule Microwaveprop.Weather do
     |> Map.drop([:profile, :duct_characteristics, :surface_temp_c, :surface_dewpoint_c])
   end
 
+  @spec upsert_gefs_profile(map()) :: {:ok, GefsProfile.t()} | {:error, Ecto.Changeset.t()}
+  def upsert_gefs_profile(attrs) do
+    changeset = GefsProfile.changeset(%GefsProfile{}, attrs)
+
+    if changeset.valid? do
+      Repo.insert(changeset,
+        on_conflict: :nothing,
+        conflict_target: [:lat, :lon, :valid_time]
+      )
+    else
+      {:error, changeset}
+    end
+  end
+
+  @spec upsert_gefs_profiles_batch([map()]) :: {non_neg_integer(), nil}
+  def upsert_gefs_profiles_batch(profiles) do
+    now = DateTime.truncate(DateTime.utc_now(), :second)
+
+    profiles
+    |> Enum.chunk_every(500)
+    |> Enum.reduce({0, nil}, fn chunk, {total_count, _} ->
+      entries =
+        Enum.map(chunk, fn attrs ->
+          Map.merge(attrs, %{
+            id: Ecto.UUID.generate(),
+            inserted_at: now,
+            updated_at: now
+          })
+        end)
+
+      {count, rows} =
+        Repo.insert_all(GefsProfile, entries,
+          on_conflict: :nothing,
+          conflict_target: [:lat, :lon, :valid_time]
+        )
+
+      {total_count + count, rows}
+    end)
+  end
+
   @spec upsert_hrrr_profile(map()) :: {:ok, HrrrProfile.t()} | {:error, Ecto.Changeset.t()}
   def upsert_hrrr_profile(attrs) do
     %HrrrProfile{}
diff --git a/lib/microwaveprop/weather/gefs_client.ex b/lib/microwaveprop/weather/gefs_client.ex
new file mode 100644
index 00000000..542c3676
--- /dev/null
+++ b/lib/microwaveprop/weather/gefs_client.ex
@@ -0,0 +1,112 @@
+defmodule Microwaveprop.Weather.GefsClient do
+  @moduledoc """
+  Client for the NCEP Global Ensemble Forecast System (GEFS) served
+  anonymously from NOMADS at
+  `https://nomads.ncep.noaa.gov/pub/data/nccf/com/gens/prod/`.
+
+  GEFS is the 0.5° / 31-member physics-based global ensemble. It runs
+  four times per day (00/06/12/18Z) with 3-hourly output to +240 h and
+  6-hourly output to +384 h. The `pgrb2ap5/` product family contains
+  the subset of surface and pressure-level variables the propagation
+  scorer needs for extended-horizon (Day 2-7) outlooks, once HRRR's
+  18 h horizon runs out.
+
+  This module exposes the URL, time, and message-inventory helpers.
+  Byte-range idx parsing and HTTP range fetch logic are shared with
+  `Microwaveprop.Weather.HrrrClient` since the GRIB2/idx format is
+  identical; the GEFS-specific pieces are the NOMADS path shape, the
+  3-then-6-hourly forecast-hour list, and the Magnus dewpoint
+  derivation (GEFS pgrb2a publishes RH rather than Td).
+  """
+
+  @base_url "https://nomads.ncep.noaa.gov/pub/data/nccf/com/gens/prod"
+
+  @run_hours [0, 6, 12, 18]
+
+  # Variables the propagation scorer reads. All present in pgrb2a f000-f384
+  # per an April 2026 .idx probe. DPT is absent from GEFS output — derive
+  # it from TMP + RH via Magnus (see `dewpoint_from_rh/2`). HPBL is also
+  # absent; the scorer tolerates a nil HPBL without penalising the score.
+  @surface_messages [
+    %{var: "TMP", level: "2 m above ground"},
+    %{var: "RH", level: "2 m above ground"},
+    %{var: "PRES", level: "surface"},
+    %{var: "PWAT", level: "entire atmosphere (considered as a single layer)"},
+    %{var: "UGRD", level: "10 m above ground"},
+    %{var: "VGRD", level: "10 m above ground"},
+    %{var: "TCDC", level: "entire atmosphere"},
+    %{var: "APCP", level: "surface"}
+  ]
+
+  @doc """
+  Rounds a wall-clock `DateTime` down to the nearest GEFS run cycle
+  (00/06/12/18Z) on the same UTC day.
+  """
+  @spec nearest_run(DateTime.t()) :: DateTime.t()
+  def nearest_run(%DateTime{} = dt) do
+    run_hour = dt.hour |> div(6) |> Kernel.*(6)
+
+    dt
+    |> DateTime.truncate(:second)
+    |> Map.put(:hour, run_hour)
+    |> Map.put(:minute, 0)
+    |> Map.put(:second, 0)
+  end
+
+  @doc """
+  Builds the NOMADS URL for the GEFS ensemble-mean (`geavg`) pgrb2a 0.5°
+  GRIB2 file at the given `run_date`, `run_hour`, and `forecast_hour`.
+  """
+  @spec ensmean_url(Date.t(), non_neg_integer(), non_neg_integer()) :: String.t()
+  def ensmean_url(%Date{} = run_date, run_hour, forecast_hour) when run_hour in @run_hours and forecast_hour >= 0 do
+    date_str = Calendar.strftime(run_date, "%Y%m%d")
+    hh = pad2(run_hour)
+    fff = pad3(forecast_hour)
+
+    "#{@base_url}/gefs.#{date_str}/#{hh}/atmos/pgrb2ap5/geavg.t#{hh}z.pgrb2a.0p50.f#{fff}"
+  end
+
+  @doc "Appends `.idx` to a GRIB2 URL to locate its wgrib2 byte-range index."
+  @spec idx_url(String.t()) :: String.t()
+  def idx_url(grib_url) when is_binary(grib_url), do: grib_url <> ".idx"
+
+  @doc """
+  Returns the list of GEFS forecast hours published per run: 3-hourly
+  from f000 through f240, then 6-hourly from f246 through f384.
+  """
+  @spec forecast_hours() :: [non_neg_integer()]
+  def forecast_hours do
+    short = Enum.to_list(0..240//3)
+    long = Enum.to_list(246..384//6)
+    short ++ long
+  end
+
+  @doc """
+  GRIB2 message descriptors the propagation scorer needs from the
+  GEFS pgrb2a product. Values match the `var` and `level` strings that
+  appear in the file's `.idx` sidecar exactly.
+  """
+  @spec surface_messages() :: [%{var: String.t(), level: String.t()}]
+  def surface_messages, do: @surface_messages
+
+  @doc """
+  Derives dewpoint in °C from air temperature (°C) and relative humidity
+  (%) using the Magnus formula with the Bolton (1980) coefficients.
+
+  GEFS pgrb2a publishes 2 m RH but not 2 m DPT, so downstream code that
+  consumed HRRR's `surface_dewpoint_c` must call this helper to fill in
+  the equivalent value. RH is clamped to a small positive floor to keep
+  the logarithm finite at 0 %.
+  """
+  @spec dewpoint_from_rh(float(), float()) :: float()
+  def dewpoint_from_rh(temp_c, rh_pct) when is_number(temp_c) and is_number(rh_pct) do
+    a = 17.625
+    b = 243.04
+    rh = max(rh_pct, 0.01) / 100.0
+    gamma = :math.log(rh) + a * temp_c / (b + temp_c)
+    b * gamma / (a - gamma)
+  end
+
+  defp pad2(n), do: n |> Integer.to_string() |> String.pad_leading(2, "0")
+  defp pad3(n), do: n |> Integer.to_string() |> String.pad_leading(3, "0")
+end
diff --git a/lib/microwaveprop/weather/gefs_profile.ex b/lib/microwaveprop/weather/gefs_profile.ex
new file mode 100644
index 00000000..393e7901
--- /dev/null
+++ b/lib/microwaveprop/weather/gefs_profile.ex
@@ -0,0 +1,60 @@
+defmodule Microwaveprop.Weather.GefsProfile do
+  @moduledoc """
+  Grid-point atmospheric profile sourced from NCEP GEFS ensemble-mean
+  (`geavg`) pgrb2a 0.5° output. Mirrors the shape of
+  `Microwaveprop.Weather.HrrrProfile` so the existing scorer can run
+  over rows from either source for a given `valid_time`.
+
+  GEFS publishes relative humidity rather than dewpoint, so
+  `surface_dewpoint_c` on these rows is derived by
+  `Microwaveprop.Weather.GefsClient.dewpoint_from_rh/2` at ingest time.
+  HPBL is not in the pgrb2a product; that column intentionally doesn't
+  exist here.
+  """
+
+  use Ecto.Schema
+
+  import Ecto.Changeset
+
+  @primary_key {:id, :binary_id, autogenerate: true}
+  @foreign_key_type :binary_id
+
+  schema "gefs_profiles" do
+    field :valid_time, :utc_datetime
+    field :lat, :float
+    field :lon, :float
+    field :run_time, :utc_datetime
+    field :forecast_hour, :integer
+    field :profile, {:array, :map}
+    field :pwat_mm, :float
+    field :surface_temp_c, :float
+    field :surface_dewpoint_c, :float
+    field :surface_pressure_mb, :float
+    field :surface_refractivity, :float
+    field :min_refractivity_gradient, :float
+    field :ducting_detected, :boolean, default: false
+    field :duct_characteristics, {:array, :map}
+    field :wind_u_mps, :float
+    field :wind_v_mps, :float
+    field :cloud_cover_pct, :float
+    field :precip_mm, :float
+
+    timestamps(type: :utc_datetime)
+  end
+
+  @type t :: %__MODULE__{}
+
+  @required_fields ~w(valid_time lat lon)a
+  @optional_fields ~w(run_time forecast_hour profile pwat_mm surface_temp_c surface_dewpoint_c
+                      surface_pressure_mb surface_refractivity min_refractivity_gradient
+                      ducting_detected duct_characteristics wind_u_mps wind_v_mps
+                      cloud_cover_pct precip_mm)a
+
+  @spec changeset(t() | Ecto.Changeset.t(), map()) :: Ecto.Changeset.t()
+  def changeset(gefs_profile, attrs) do
+    gefs_profile
+    |> cast(attrs, @required_fields ++ @optional_fields)
+    |> validate_required(@required_fields)
+    |> unique_constraint([:lat, :lon, :valid_time])
+  end
+end
diff --git a/priv/repo/migrations/20260418191400_create_gefs_profiles.exs b/priv/repo/migrations/20260418191400_create_gefs_profiles.exs
new file mode 100644
index 00000000..ee54fa95
--- /dev/null
+++ b/priv/repo/migrations/20260418191400_create_gefs_profiles.exs
@@ -0,0 +1,33 @@
+defmodule Microwaveprop.Repo.Migrations.CreateGefsProfiles do
+  use Ecto.Migration
+
+  def change do
+    create table(:gefs_profiles, primary_key: false) do
+      add :id, :binary_id, primary_key: true
+      add :valid_time, :utc_datetime, null: false
+      add :lat, :float, null: false
+      add :lon, :float, null: false
+      add :run_time, :utc_datetime
+      add :forecast_hour, :integer
+      add :profile, {:array, :map}, default: []
+      add :pwat_mm, :float
+      add :surface_temp_c, :float
+      add :surface_dewpoint_c, :float
+      add :surface_pressure_mb, :float
+      add :surface_refractivity, :float
+      add :min_refractivity_gradient, :float
+      add :ducting_detected, :boolean, default: false
+      add :duct_characteristics, {:array, :map}
+      add :wind_u_mps, :float
+      add :wind_v_mps, :float
+      add :cloud_cover_pct, :float
+      add :precip_mm, :float
+
+      timestamps(type: :utc_datetime)
+    end
+
+    create unique_index(:gefs_profiles, [:lat, :lon, :valid_time])
+    create index(:gefs_profiles, [:valid_time])
+    create index(:gefs_profiles, [:run_time])
+  end
+end
diff --git a/test/microwaveprop/weather/gefs_client_test.exs b/test/microwaveprop/weather/gefs_client_test.exs
new file mode 100644
index 00000000..77a66bd4
--- /dev/null
+++ b/test/microwaveprop/weather/gefs_client_test.exs
@@ -0,0 +1,123 @@
+defmodule Microwaveprop.Weather.GefsClientTest do
+  use ExUnit.Case, async: true
+
+  alias Microwaveprop.Weather.GefsClient
+
+  describe "nearest_run/1" do
+    test "rounds down to the same 6-hour cycle within the cycle window" do
+      assert GefsClient.nearest_run(~U[2026-04-18 13:30:00Z]) == ~U[2026-04-18 12:00:00Z]
+    end
+
+    test "snaps to 00Z when given 02:15" do
+      assert GefsClient.nearest_run(~U[2026-04-18 02:15:00Z]) == ~U[2026-04-18 00:00:00Z]
+    end
+
+    test "snaps to 18Z when given 23:59" do
+      assert GefsClient.nearest_run(~U[2026-04-18 23:59:00Z]) == ~U[2026-04-18 18:00:00Z]
+    end
+
+    test "is a no-op when already on a cycle boundary" do
+      assert GefsClient.nearest_run(~U[2026-04-18 06:00:00Z]) == ~U[2026-04-18 06:00:00Z]
+    end
+  end
+
+  describe "ensmean_url/3" do
+    test "builds a pgrb2ap5 ensemble-mean URL for a run + forecast hour" do
+      url = GefsClient.ensmean_url(~D[2026-04-18], 12, 24)
+
+      assert url ==
+               "https://nomads.ncep.noaa.gov/pub/data/nccf/com/gens/prod/gefs.20260418/12/atmos/pgrb2ap5/geavg.t12z.pgrb2a.0p50.f024"
+    end
+
+    test "pads forecast hour to three digits" do
+      url = GefsClient.ensmean_url(~D[2026-04-18], 0, 6)
+      assert String.ends_with?(url, "geavg.t00z.pgrb2a.0p50.f006")
+    end
+
+    test "handles 240-hour forecast horizon" do
+      url = GefsClient.ensmean_url(~D[2026-04-18], 12, 240)
+      assert String.ends_with?(url, "geavg.t12z.pgrb2a.0p50.f240")
+    end
+
+    test "pads single-digit runs" do
+      url = GefsClient.ensmean_url(~D[2026-04-18], 6, 24)
+      assert String.contains?(url, "/06/atmos/")
+      assert String.contains?(url, "geavg.t06z.")
+    end
+  end
+
+  describe "idx_url/1" do
+    test "appends .idx to the GRIB2 URL" do
+      base = GefsClient.ensmean_url(~D[2026-04-18], 12, 24)
+      assert GefsClient.idx_url(base) == base <> ".idx"
+    end
+  end
+
+  describe "forecast_hours/0" do
+    test "returns 3-hourly hours out to 240 then 6-hourly to 384" do
+      hours = GefsClient.forecast_hours()
+      assert Enum.take(hours, 4) == [0, 3, 6, 9]
+      assert 240 in hours
+      assert 246 in hours
+      assert 384 in hours
+      assert List.last(hours) == 384
+      refute 243 in hours
+    end
+  end
+
+  describe "surface_messages/0" do
+    test "lists the variables the scorer needs from the pgrb2a file" do
+      msgs = GefsClient.surface_messages()
+
+      vars = Enum.map(msgs, & &1.var)
+      assert "TMP" in vars
+      assert "RH" in vars
+      assert "PRES" in vars
+      assert "PWAT" in vars
+      assert "UGRD" in vars
+      assert "VGRD" in vars
+      assert "TCDC" in vars
+      assert "APCP" in vars
+    end
+
+    test "2m TMP and RH use the '2 m above ground' level string" do
+      msgs = GefsClient.surface_messages()
+
+      assert %{var: "TMP", level: "2 m above ground"} in msgs
+      assert %{var: "RH", level: "2 m above ground"} in msgs
+    end
+
+    test "10m winds use the '10 m above ground' level string" do
+      msgs = GefsClient.surface_messages()
+
+      assert %{var: "UGRD", level: "10 m above ground"} in msgs
+      assert %{var: "VGRD", level: "10 m above ground"} in msgs
+    end
+  end
+
+  describe "dewpoint_from_rh/2" do
+    test "returns the input temperature when RH is 100%" do
+      # At saturation, Td == T. Magnus returns T exactly.
+      assert_in_delta GefsClient.dewpoint_from_rh(20.0, 100.0), 20.0, 0.05
+    end
+
+    test "returns a cooler dewpoint than temperature when RH is below saturation" do
+      td = GefsClient.dewpoint_from_rh(25.0, 50.0)
+      assert td < 25.0
+      # Tables put Td(25 °C, 50% RH) at ~13.9 °C. Allow ±0.5 °C slack.
+      assert_in_delta td, 13.9, 0.5
+    end
+
+    test "handles very dry air without crashing" do
+      td = GefsClient.dewpoint_from_rh(30.0, 5.0)
+      assert is_float(td)
+      assert td < 0.0
+    end
+
+    test "clamps RH<=0 to a very-dry fallback rather than returning NaN/Inf" do
+      td = GefsClient.dewpoint_from_rh(10.0, 0.0)
+      assert is_float(td)
+      assert td < -30.0
+    end
+  end
+end
diff --git a/test/microwaveprop/weather/gefs_profile_test.exs b/test/microwaveprop/weather/gefs_profile_test.exs
new file mode 100644
index 00000000..be4e85ec
--- /dev/null
+++ b/test/microwaveprop/weather/gefs_profile_test.exs
@@ -0,0 +1,99 @@
+defmodule Microwaveprop.Weather.GefsProfileTest do
+  use Microwaveprop.DataCase, async: true
+
+  alias Microwaveprop.Weather
+  alias Microwaveprop.Weather.GefsProfile
+
+  @valid_time ~U[2026-04-18 18:00:00Z]
+
+  defp base_attrs(overrides \\ %{}) do
+    Map.merge(
+      %{
+        valid_time: @valid_time,
+        lat: 32.9,
+        lon: -97.0,
+        run_time: ~U[2026-04-18 12:00:00Z],
+        forecast_hour: 6,
+        surface_temp_c: 20.0,
+        surface_dewpoint_c: 10.0,
+        surface_pressure_mb: 1013.0,
+        pwat_mm: 20.0,
+        wind_u_mps: 2.0,
+        wind_v_mps: -1.0,
+        cloud_cover_pct: 50.0,
+        precip_mm: 0.0
+      },
+      overrides
+    )
+  end
+
+  describe "changeset/2" do
+    test "is valid with the minimum required fields" do
+      cs = GefsProfile.changeset(%GefsProfile{}, base_attrs())
+      assert cs.valid?
+    end
+
+    test "requires valid_time, lat, and lon" do
+      cs = GefsProfile.changeset(%GefsProfile{}, %{})
+      refute cs.valid?
+      assert %{valid_time: _, lat: _, lon: _} = errors_on(cs)
+    end
+
+    test "accepts an optional pressure-level profile array" do
+      profile = [%{"pres" => 1000.0, "tmpc" => 20.0, "dwpc" => 10.0, "hght" => 100.0}]
+      cs = GefsProfile.changeset(%GefsProfile{}, base_attrs(%{profile: profile}))
+      assert cs.valid?
+      assert Ecto.Changeset.get_field(cs, :profile) == profile
+    end
+  end
+
+  describe "Weather.upsert_gefs_profile/1" do
+    test "inserts a new profile" do
+      assert {:ok, %GefsProfile{} = profile} = Weather.upsert_gefs_profile(base_attrs())
+      assert profile.lat == 32.9
+      assert profile.surface_temp_c == 20.0
+    end
+
+    test "is a no-op on conflict (lat, lon, valid_time)" do
+      {:ok, _} = Weather.upsert_gefs_profile(base_attrs())
+
+      {:ok, second} =
+        Weather.upsert_gefs_profile(base_attrs(%{surface_temp_c: 999.0}))
+
+      # on_conflict: :nothing leaves the existing row alone; the returned
+      # struct just carries the changeset values, but the DB keeps the old
+      # record.
+      persisted = Microwaveprop.Repo.get_by(GefsProfile, lat: 32.9, lon: -97.0, valid_time: @valid_time)
+      assert persisted.surface_temp_c == 20.0
+      assert second.lat == 32.9
+    end
+
+    test "rejects invalid attributes with a changeset error" do
+      assert {:error, %Ecto.Changeset{valid?: false}} = Weather.upsert_gefs_profile(%{})
+    end
+  end
+
+  describe "Weather.upsert_gefs_profiles_batch/1" do
+    test "inserts many rows at once" do
+      attrs =
+        for fh <- [3, 6, 9] do
+          valid = DateTime.add(~U[2026-04-18 12:00:00Z], fh * 3600, :second)
+          base_attrs(%{valid_time: valid, forecast_hour: fh})
+        end
+
+      {count, _} = Weather.upsert_gefs_profiles_batch(attrs)
+      assert count == 3
+
+      rows = Microwaveprop.Repo.all(GefsProfile)
+      assert length(rows) == 3
+    end
+
+    test "skips conflicting rows quietly" do
+      attrs = [base_attrs()]
+      {1, _} = Weather.upsert_gefs_profiles_batch(attrs)
+      {0, _} = Weather.upsert_gefs_profiles_batch(attrs)
+
+      assert Microwaveprop.Repo.aggregate(GefsProfile, :count, :id) == 1
+    end
+  end
+end