From f632ea83bdd128d228eba1d17b7c2e0739e02774 Mon Sep 17 00:00:00 2001
From: Graham McIntire <graham@mcintire.me>
Date: Wed, 15 Apr 2026 14:43:47 -0500
Subject: [PATCH] Sporadic-E physics module + Ionosphere.nearest_foes
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

Adds the missing link between the GIRO ionosonde data we just started
ingesting and the VHF/UHF band scoring:

* `Microwaveprop.Propagation.SporadicE` — ITU-R P.534-6 / Davies 1990
  thin-layer Es MUF (sec(i) · foEs with h=110 km), plus an `es_score/3`
  that maps (foEs, target band, hop distance) into a 0-100 single-hop
  propagation likelihood. Calibrated against literature: 50 MHz Es at
  2000 km needs foEs ≳ 5.5 MHz (routine summer); 144 MHz needs
  foEs ≳ 15.8 MHz (rare Jun/Jul peaks only); 440 MHz Es does not occur
  at physical foEs values. Multi-hop Es and Es-scatter are separate
  factors and explicitly out of scope.

* `Ionosphere.nearest_foes/3` — given a lat/lon, returns the latest
  observation from the nearest polled GIRO station (Millstone Hill or
  Alpena for now), with a configurable staleness cutoff (default 2h).
  Returns `{:error, :stale}` or `{:error, :no_data}` so callers can
  choose whether to apply the Es factor at all.

Neither is wired into the grid scorer yet — that's a separate commit
so the integration can be reviewed on its own.
---
 lib/microwaveprop/ionosphere.ex               | 64 ++++++++++++++
 lib/microwaveprop/propagation/sporadic_e.ex   | 84 ++++++++++++++++++
 .../ionosphere/ionosphere_test.exs            | 48 +++++++++++
 .../propagation/sporadic_e_test.exs           | 86 +++++++++++++++++++
 4 files changed, 282 insertions(+)
 create mode 100644 lib/microwaveprop/propagation/sporadic_e.ex
 create mode 100644 test/microwaveprop/propagation/sporadic_e_test.exs

diff --git a/lib/microwaveprop/ionosphere.ex b/lib/microwaveprop/ionosphere.ex
index 80b8e6e8..fb51bdb9 100644
--- a/lib/microwaveprop/ionosphere.ex
+++ b/lib/microwaveprop/ionosphere.ex
@@ -49,4 +49,68 @@ defmodule Microwaveprop.Ionosphere do
         limit: 1
     )
   end
+
+  # Station coordinates for the nearest-lookup. Must stay in sync with
+  # `IonosphereFetchWorker.stations/0` — these are the stations we're
+  # actively polling. Kept here (rather than imported) to avoid a
+  # context → worker dependency.
+  @stations [
+    %{code: "MHJ45", lat: 42.6, lon: -71.5},
+    %{code: "AL945", lat: 45.1, lon: -83.6}
+  ]
+
+  @default_max_age_seconds 2 * 3600
+
+  @doc """
+  Look up the nearest polled ionosonde station's latest observation for
+  the given lat/lon.
+
+  Options:
+  - `:max_age_seconds` — reject observations older than this (default 2h).
+    Returns `{:error, :stale}` if the nearest station has no fresh data.
+
+  Returns `{:ok, observation}`, `{:error, :stale}`, or `{:error, :no_data}`.
+  """
+  @spec nearest_foes(number(), number(), keyword()) ::
+          {:ok, Observation.t()} | {:error, :stale | :no_data}
+  def nearest_foes(lat, lon, opts \\ []) do
+    max_age = Keyword.get(opts, :max_age_seconds, @default_max_age_seconds)
+    nearest_station = nearest_station_to(lat, lon)
+
+    case latest_observation(nearest_station.code) do
+      nil ->
+        {:error, :no_data}
+
+      %Observation{valid_time: vt} = obs ->
+        age = DateTime.diff(DateTime.utc_now(), vt, :second)
+
+        if age <= max_age do
+          {:ok, obs}
+        else
+          {:error, :stale}
+        end
+    end
+  end
+
+  defp nearest_station_to(lat, lon) do
+    Enum.min_by(@stations, fn s -> haversine_km(lat, lon, s.lat, s.lon) end)
+  end
+
+  # Standard haversine distance in km. Good enough for station-pick
+  # ranking — we're not computing hop geometry here.
+  defp haversine_km(lat1, lon1, lat2, lon2) do
+    r = 6371.0
+    dlat = :math.pi() * (lat2 - lat1) / 180
+    dlon = :math.pi() * (lon2 - lon1) / 180
+    lat1_rad = :math.pi() * lat1 / 180
+    lat2_rad = :math.pi() * lat2 / 180
+
+    a =
+      :math.sin(dlat / 2) * :math.sin(dlat / 2) +
+        :math.cos(lat1_rad) * :math.cos(lat2_rad) *
+          :math.sin(dlon / 2) * :math.sin(dlon / 2)
+
+    c = 2 * :math.atan2(:math.sqrt(a), :math.sqrt(1 - a))
+    r * c
+  end
 end
diff --git a/lib/microwaveprop/propagation/sporadic_e.ex b/lib/microwaveprop/propagation/sporadic_e.ex
new file mode 100644
index 00000000..4375220d
--- /dev/null
+++ b/lib/microwaveprop/propagation/sporadic_e.ex
@@ -0,0 +1,84 @@
+defmodule Microwaveprop.Propagation.SporadicE do
+  @moduledoc """
+  Sporadic-E propagation physics for VHF scoring.
+
+  Given foEs (the E-layer critical frequency measured by an ionosonde,
+  in MHz) and a ground-range distance in km, computes the single-hop
+  Es MUF via the thin-layer obliquity factor from ITU-R P.534-6 /
+  Davies 1990 *Ionospheric Radio*, Eq. 6.26:
+
+      MUF_Es = foEs · sec(i)
+
+  where `sec(i) = √(1 + (D / (2h))²)` with the reflection height
+  `h = 110 km` (standard thin-E approximation).
+
+  The `es_score/3` function maps a `(foEs, target band, distance)`
+  triple into a 0–100 propagation likelihood for **single-hop Es
+  only**. Multi-hop Es (which extends the range but needs consecutive
+  intense Es patches) and Es-scatter are not modelled here — those
+  are separate factors.
+
+  Calibration notes:
+  - 50 MHz Es at 2000 km needs foEs ≳ 5.5 MHz (routine summer Es).
+  - 144 MHz Es at 2000 km needs foEs ≳ 15.8 MHz (rare, observed during
+    major Jun/Jul peaks).
+  - 440 MHz Es does not occur at physically observed foEs values.
+  - Single-hop Es has a minimum useful range of ~500 km (steeper rays
+    escape the layer rather than reflect) and a maximum of ~2500 km
+    (the geometric horizon limit for a 110 km reflection height).
+  """
+
+  # Standard E-layer reflection height for the thin-layer approximation.
+  @e_layer_height_km 110.0
+
+  @min_hop_km 500
+  @max_hop_km 2500
+
+  @doc """
+  Single-hop Es MUF in MHz. Returns 0.0 for nil/non-positive inputs.
+  """
+  @spec single_hop_muf(number() | nil, number() | nil) :: float()
+  def single_hop_muf(nil, _), do: 0.0
+  def single_hop_muf(_, nil), do: 0.0
+  def single_hop_muf(fo_es_mhz, _) when fo_es_mhz <= 0, do: 0.0
+  def single_hop_muf(_, distance_km) when distance_km <= 0, do: 0.0
+
+  def single_hop_muf(fo_es_mhz, distance_km) do
+    sec_i = :math.sqrt(1 + :math.pow(distance_km / (2 * @e_layer_height_km), 2))
+    fo_es_mhz * sec_i
+  end
+
+  @doc """
+  Returns a 0–100 score for single-hop Es propagation of `band_mhz`
+  over a ground-range `distance_km`, given the measured `fo_es_mhz`.
+
+  - Short paths (< 500 km) and long paths (> 2500 km) are outside the
+    single-hop window → 0.
+  - `nil` or 0 foEs → 0.
+  - Otherwise the score ramps with `MUF_Es / band_mhz`:
+      ≥ 1.20 → 100 (comfortably above the critical threshold)
+      ≥ 1.00 → 80  (just above, marginal Es)
+      ≥ 0.90 → 50  (close, occasional fringe)
+      ≥ 0.70 → 20  (weak fringe, unreliable)
+      below  → 0
+  """
+  @spec es_score(number() | nil, number(), number() | nil) :: integer()
+  def es_score(nil, _band, _distance), do: 0
+  def es_score(_fo_es, _band, nil), do: 0
+  def es_score(fo_es_mhz, _band, _distance) when fo_es_mhz <= 0, do: 0
+
+  def es_score(_fo_es, _band, distance_km) when distance_km < @min_hop_km or distance_km > @max_hop_km, do: 0
+
+  def es_score(fo_es_mhz, band_mhz, distance_km) do
+    muf = single_hop_muf(fo_es_mhz, distance_km)
+    ratio = muf / band_mhz
+
+    cond do
+      ratio >= 1.20 -> 100
+      ratio >= 1.00 -> 80
+      ratio >= 0.90 -> 50
+      ratio >= 0.70 -> 20
+      true -> 0
+    end
+  end
+end
diff --git a/test/microwaveprop/ionosphere/ionosphere_test.exs b/test/microwaveprop/ionosphere/ionosphere_test.exs
index 8392b3a1..599d3e3f 100644
--- a/test/microwaveprop/ionosphere/ionosphere_test.exs
+++ b/test/microwaveprop/ionosphere/ionosphere_test.exs
@@ -80,4 +80,52 @@ defmodule Microwaveprop.IonosphereTest do
       assert Ionosphere.latest_observation("MHJ45") == nil
     end
   end
+
+  describe "nearest_foes/3" do
+    setup do
+      recent = DateTime.truncate(DateTime.utc_now(), :second)
+
+      {:ok, _} =
+        Ionosphere.upsert_observations("MHJ45", [
+          %{valid_time: recent, fo_es_mhz: 7.5, fo_f2_mhz: 8.0, mufd_mhz: 25.0}
+        ])
+
+      {:ok, _} =
+        Ionosphere.upsert_observations("AL945", [
+          %{valid_time: recent, fo_es_mhz: 12.3, fo_f2_mhz: 7.5, mufd_mhz: 24.0}
+        ])
+
+      :ok
+    end
+
+    test "returns the nearest station's latest foEs for a grid point in New England" do
+      # 42°N, -72°W is essentially on Millstone Hill (42.6, -71.5).
+      assert {:ok, obs} = Ionosphere.nearest_foes(42.0, -72.0)
+      assert obs.station_code == "MHJ45"
+      assert obs.fo_es_mhz == 7.5
+    end
+
+    test "returns the nearest station's latest foEs for a grid point in Michigan" do
+      # 45°N, -84°W is next to Alpena (45.1, -83.6).
+      assert {:ok, obs} = Ionosphere.nearest_foes(45.0, -84.0)
+      assert obs.station_code == "AL945"
+      assert obs.fo_es_mhz == 12.3
+    end
+
+    test "returns {:error, :stale} when the nearest station's latest value is older than max_age" do
+      stale = ~U[2020-01-01 00:00:00Z]
+
+      Microwaveprop.Repo.update_all(
+        Observation,
+        set: [valid_time: stale]
+      )
+
+      assert {:error, :stale} = Ionosphere.nearest_foes(42.0, -72.0, max_age_seconds: 3600)
+    end
+
+    test "returns {:error, :no_data} when no observations have been ingested" do
+      Microwaveprop.Repo.delete_all(Observation)
+      assert {:error, :no_data} = Ionosphere.nearest_foes(42.0, -72.0)
+    end
+  end
 end
diff --git a/test/microwaveprop/propagation/sporadic_e_test.exs b/test/microwaveprop/propagation/sporadic_e_test.exs
new file mode 100644
index 00000000..6bfb873e
--- /dev/null
+++ b/test/microwaveprop/propagation/sporadic_e_test.exs
@@ -0,0 +1,86 @@
+defmodule Microwaveprop.Propagation.SporadicETest do
+  use ExUnit.Case, async: true
+
+  alias Microwaveprop.Propagation.SporadicE
+
+  describe "single_hop_muf/2" do
+    test "returns 0.0 for nil or non-positive inputs" do
+      assert SporadicE.single_hop_muf(nil, 2000) == 0.0
+      assert SporadicE.single_hop_muf(5.0, nil) == 0.0
+      assert SporadicE.single_hop_muf(0.0, 2000) == 0.0
+      assert SporadicE.single_hop_muf(5.0, 0) == 0.0
+      assert SporadicE.single_hop_muf(-1.0, 2000) == 0.0
+    end
+
+    test "matches the thin-layer secant-of-incidence formula for known distances" do
+      # h = 110 km. sec(i) = sqrt(1 + (D/(2h))^2). Values cross-checked
+      # against literature (Davies 1990, Ionospheric Radio Eq. 6.26).
+      # At D=2000 km, sec(i) ≈ 9.146, so foEs × 9.146 MHz.
+      assert_in_delta SporadicE.single_hop_muf(10.0, 2000), 10.0 * 9.146, 0.05
+      assert_in_delta SporadicE.single_hop_muf(5.0, 1000), 5.0 * 4.654, 0.05
+      assert_in_delta SporadicE.single_hop_muf(8.0, 1500), 8.0 * 6.891, 0.05
+    end
+
+    test "degenerate short range approaches foEs (vertical incidence)" do
+      # As D → 0 the ray is nearly vertical and MUF → foEs.
+      assert_in_delta SporadicE.single_hop_muf(5.0, 1), 5.0, 0.001
+    end
+  end
+
+  describe "es_score/3" do
+    test "is 0 when foEs is nil or zero" do
+      assert SporadicE.es_score(nil, 144, 2000) == 0
+      assert SporadicE.es_score(0.0, 144, 2000) == 0
+    end
+
+    test "is 0 for paths shorter than the single-hop minimum" do
+      # Es geometry doesn't support paths much under ~500 km (the
+      # elevation angle is too steep for the layer to reflect a VHF
+      # signal regardless of foEs).
+      assert SporadicE.es_score(20.0, 144, 100) == 0
+      assert SporadicE.es_score(20.0, 144, 400) == 0
+    end
+
+    test "is 0 for paths beyond the single-hop maximum (~2500 km)" do
+      # Beyond ~2500 km the geometry exits the single-hop window. We
+      # don't model multi-hop Es here — that's a different factor.
+      assert SporadicE.es_score(20.0, 144, 3000) == 0
+    end
+
+    test "scores 50 MHz as ROUTINE when foEs is typical summer Es (6-8 MHz) at 2000 km" do
+      # foEs 6 MHz × sec(i) at 2000 km ≈ 54.9 MHz, just above 50 MHz.
+      # 50 MHz Es at 2000 km with foEs=6 is a classic summer opening —
+      # should score high.
+      score = SporadicE.es_score(6.0, 50, 2000)
+      assert score >= 80
+    end
+
+    test "scores 144 MHz as NONE when foEs is typical summer Es" do
+      # foEs=6 at 2000 km ≈ 54.9 MHz MUF — well below 144 MHz.
+      # Typical summer Es doesn't propagate 2m.
+      assert SporadicE.es_score(6.0, 144, 2000) == 0
+    end
+
+    test "scores 144 MHz as VIABLE when foEs reaches intense-Es levels (16+ MHz) at 2000 km" do
+      # foEs=16 × 9.15 ≈ 146 MHz MUF → 2m Es is feasible. These are
+      # rare but documented events (June/July peaks).
+      score = SporadicE.es_score(16.0, 144, 2000)
+      assert score >= 50
+    end
+
+    test "scores 440 MHz as NONE even at extreme foEs (sporadic-E does not reach 70cm)" do
+      # foEs would need to be ~48 MHz to get 440 MHz via single-hop Es
+      # — way beyond any physically observed value.
+      assert SporadicE.es_score(25.0, 440, 2000) == 0
+    end
+
+    test "higher foEs strictly increases the score at a given band/distance (monotonic)" do
+      weak = SporadicE.es_score(10.0, 50, 2000)
+      mid = SporadicE.es_score(14.0, 50, 2000)
+      strong = SporadicE.es_score(20.0, 50, 2000)
+      assert weak <= mid
+      assert mid <= strong
+      assert strong == 100
+    end
+  end
+end