From 6aa91e7656af1b23eab3cd3106177c5a6ce60d88 Mon Sep 17 00:00:00 2001 From: Graham McIntire Date: Sat, 25 Apr 2026 10:52:42 -0500 Subject: [PATCH] =?UTF-8?q?fix:=20April=202026=20codebase=20review=20?= =?UTF-8?q?=E2=80=94=20address=2013=20bugs=20across=20propagation=20chain?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Each fix is covered by a regression test that fails on `main` and passes on this commit. Round 1 (initial review): * propagation: thread `latitude` into the conditions map so `score_season/4` actually picks up regional multipliers * hrrr_client / fetcher.rs: `nearest_hrrr_hour` rounds DOWN, never at a future cycle that NOAA hasn't published yet * radio: spherical-vector great-circle midpoint replaces the arithmetic mean — anti-meridian paths no longer fold to Greenwich * weather: `reconcile_weather_statuses` scales the longitude band by `1 / cos(lat)` so the bbox stays ~150 km wide at every latitude * radio/maidenhead: clamp 90°/180° below the field-bucket overflow so `from_latlon` never emits invalid characters like 'S' * prop_grid_rs/pipeline: merge HRRR + NEXRAD-derived rain rates and read `best_duct_freq_ghz` into `best_duct_band_ghz` so the Native Duct Boost actually fires * propagation/region (Elixir + Rust): inclusive upper bounds so points exactly at lat_max get the regional multiplier * weather/sounding_params (Elixir + Rust): drop the 10 m gradient floor so HRRR's thin near-surface layers stop hiding sharp ducts * weather/sounding_params: when the profile ends inside a duct, finalize it with the highest sample as the top instead of throwing it away (Rust port already correct) Round 2 (post-fix sweep): * radio + commercial: single canonical haversine in Radio (atan2 form); Commercial delegates instead of carrying a second copy that could disagree at threshold distances * prop_grid_rs/profiles_file: `snap_coords` matches Elixir's step-aware snap (`round(coord/0.125) * 0.125`, then 3-dp round) so Rust-keyed and Elixir-keyed profile maps land on the same cell * weather/grib2/wgrib2: `parse_lon_val_segment` uses `Float.parse` uniformly — wgrib2 dropping the trailing `.0` from a longitude no longer crashes the whole chain step --- bugs.md | 27 +++++++++ lib/microwaveprop/commercial.ex | 19 ++---- lib/microwaveprop/propagation.ex | 3 +- lib/microwaveprop/propagation/region.ex | 2 +- lib/microwaveprop/radio.ex | 38 +++++++++++- lib/microwaveprop/radio/maidenhead.ex | 7 ++- lib/microwaveprop/weather.ex | 17 +++--- lib/microwaveprop/weather/grib2/wgrib2.ex | 22 +++++-- lib/microwaveprop/weather/hrrr_client.ex | 16 ++--- lib/microwaveprop/weather/sounding_params.ex | 20 ++++--- rust/prop_grid_rs/src/fetcher.rs | 27 +++++---- rust/prop_grid_rs/src/pipeline.rs | 59 ++++++++++++++++++- rust/prop_grid_rs/src/profiles_file.rs | 36 ++++++++--- rust/prop_grid_rs/src/region.rs | 15 ++++- rust/prop_grid_rs/src/sounding_params.rs | 5 +- .../microwaveprop/propagation/region_test.exs | 8 +++ test/microwaveprop/propagation_test.exs | 29 +++++++++ test/microwaveprop/radio/maidenhead_test.exs | 12 ++++ test/microwaveprop/radio_test.exs | 25 ++++++++ .../weather/grib2/wgrib2_test.exs | 24 ++++++++ .../weather/grid_snap_property_test.exs | 19 +++--- .../weather/hrrr_client_test.exs | 17 +++--- .../weather/sounding_params_test.exs | 37 ++++++++++++ test/microwaveprop/weather_test.exs | 34 +++++++++++ 24 files changed, 430 insertions(+), 88 deletions(-) create mode 100644 bugs.md diff --git a/bugs.md b/bugs.md new file mode 100644 index 00000000..70b8cb95 --- /dev/null +++ b/bugs.md @@ -0,0 +1,27 @@ +# Identified Bugs and Observations (Post-Fix Review) + +This document reflects the state of the codebase after the fixes for the initial 15 bugs were applied. + +All three items below were verified, covered by failing tests, and fixed on 2026-04-25. + +## 1. Haversine Implementation Inconsistency — FIXED 2026-04-25 +**Files:** `lib/microwaveprop/radio.ex` vs `lib/microwaveprop/commercial.ex` +**Description:** +- `Radio.haversine_km/4` uses the `asin` form: `2 * r * :math.asin(:math.sqrt(a))` +- `Commercial.haversine_km/4` uses the `atan2` form: `2 * r * :math.atan2(:math.sqrt(a), :math.sqrt(1 - a))` +**Impact:** While functionally similar, `atan2` is the industry standard for floating-point stability at large distances. Having two different implementations for the same physical calculation can lead to tiny discrepancies (sub-millimeter) that might cause a contact or commercial link to intermittently cross a distance threshold (e.g., 75 km) depending on which module performs the check. +**Fix:** Switched `Radio.haversine_km/4` to the `atan2` form and made `Commercial`'s private wrapper delegate to it. One canonical implementation, no per-module drift. Numerical-stability regression test added in `test/microwaveprop/radio_test.exs`. + +## 2. Grid Snapping Mathematical Drift — FIXED 2026-04-25 +**Files:** `lib/microwaveprop/propagation/profiles_file.ex` vs `rust/prop_grid_rs/src/profiles_file.rs` +**Description:** +- Elixir `snap/2`: `Float.round(Float.round(lat / step) * step, 3)` +- Rust `snap_coords`: `(lat * 1000.0).round() / 1000.0` +**Impact:** The Elixir version performs a division and multiplication by the grid step (`0.125`) before rounding. The Rust version rounds the raw coordinate to the third decimal place. For coordinates that land exactly on a $0.0625$ boundary (the halfway point between cells), these two methods may snap to different cells, causing the UI to report "No Data" or show the wrong factor breakdown for a clicked point. +**Fix:** Rust `snap_coords` now mirrors the Elixir step-aware snap (`round(coord / 0.125) * 0.125`, then 3-decimal round). Cross-stack parity test in `rust/prop_grid_rs/src/profiles_file.rs` enumerates known-disagreeing inputs. + +## 3. Potential Exception in Wgrib2 Parsing (Theoretical) — FIXED 2026-04-25 +**File:** `lib/microwaveprop/weather/grib2/wgrib2.ex` +**Description:** `parse_lon_val_segment/1` uses `String.to_float(lon_str)` for the longitude, but `Float.parse` for latitude and value. +**Impact:** If `wgrib2` were to ever return an integer string for longitude (e.g. `"235"`) instead of a float (`"235.0"`), `String.to_float` would raise an `ArgumentError`. While `wgrib2` typically outputs floats, using `Float.parse` would be more consistent with the other fields in the same line. +**Fix:** All three numeric fields now go through `Float.parse/1` (in a `with` chain so a malformed segment returns `nil` instead of bringing down the chain step). Function exposed as `@doc false` and covered by direct unit tests for both decimal-formatted and integer-formatted longitude strings. diff --git a/lib/microwaveprop/commercial.ex b/lib/microwaveprop/commercial.ex index d0b4e6b9..4ad1d31d 100644 --- a/lib/microwaveprop/commercial.ex +++ b/lib/microwaveprop/commercial.ex @@ -166,20 +166,11 @@ defmodule Microwaveprop.Commercial do defp link_endpoint(_), do: nil - defp haversine_km(lat1, lon1, lat2, lon2) do - r = 6371.0 - phi1 = lat1 * :math.pi() / 180 - phi2 = lat2 * :math.pi() / 180 - dphi = (lat2 - lat1) * :math.pi() / 180 - dlambda = (lon2 - lon1) * :math.pi() / 180 - - a = - :math.sin(dphi / 2) * :math.sin(dphi / 2) + - :math.cos(phi1) * :math.cos(phi2) * :math.sin(dlambda / 2) * :math.sin(dlambda / 2) - - c = 2 * :math.atan2(:math.sqrt(a), :math.sqrt(1 - a)) - r * c - end + # Single canonical haversine implementation lives in `Microwaveprop.Radio`. + # Two formulae for the same physical distance let identical coordinates + # land on opposite sides of a 75 km link-radius threshold depending on + # which module asked. + defp haversine_km(lat1, lon1, lat2, lon2), do: Microwaveprop.Radio.haversine_km(lat1, lon1, lat2, lon2) defp average([]), do: 0.0 defp average(list), do: Enum.sum(list) / length(list) diff --git a/lib/microwaveprop/propagation.ex b/lib/microwaveprop/propagation.ex index 28fad1a2..6979208b 100644 --- a/lib/microwaveprop/propagation.ex +++ b/lib/microwaveprop/propagation.ex @@ -76,7 +76,7 @@ defmodule Microwaveprop.Propagation do score_with_algorithm(hrrr_profile, valid_time, temp_c, dewpoint_c, derived, latitude, longitude) end - defp score_with_algorithm(hrrr_profile, valid_time, temp_c, dewpoint_c, derived, _latitude, longitude) do + defp score_with_algorithm(hrrr_profile, valid_time, temp_c, dewpoint_c, derived, latitude, longitude) do temp_f = Scorer.c_to_f(temp_c) dewpoint_f = Scorer.c_to_f(dewpoint_c) @@ -89,6 +89,7 @@ defmodule Microwaveprop.Propagation do utc_hour: valid_time.hour, utc_minute: valid_time.minute, month: valid_time.month, + latitude: latitude, longitude: longitude, pressure_mb: hrrr_profile.surface_pressure_mb, prev_pressure_mb: nil, diff --git a/lib/microwaveprop/propagation/region.ex b/lib/microwaveprop/propagation/region.ex index 642d818c..42feeae5 100644 --- a/lib/microwaveprop/propagation/region.ex +++ b/lib/microwaveprop/propagation/region.ex @@ -42,7 +42,7 @@ defmodule Microwaveprop.Propagation.Region do @spec for_point(float, float) :: atom def for_point(lat, lon) do Enum.find_value(@regions, :other, fn {name, {lat_min, lat_max}, {lon_min, lon_max}} -> - if lat >= lat_min and lat < lat_max and lon >= lon_min and lon < lon_max do + if lat >= lat_min and lat <= lat_max and lon >= lon_min and lon <= lon_max do name end end) diff --git a/lib/microwaveprop/radio.ex b/lib/microwaveprop/radio.ex index 5b3e0e22..f49fba75 100644 --- a/lib/microwaveprop/radio.ex +++ b/lib/microwaveprop/radio.ex @@ -368,14 +368,43 @@ defmodule Microwaveprop.Radio do defp extract_latlon(_), do: nil defp build_path_points({lat1, lon1}, {lat2, lon2}) do - mid_lat = (lat1 + lat2) / 2 - mid_lon = (lon1 + lon2) / 2 + {mid_lat, mid_lon} = great_circle_midpoint(lat1, lon1, lat2, lon2) [{lat1, lon1}, {mid_lat, mid_lon}, {lat2, lon2}] end defp build_path_points({lat1, lon1}, _), do: [{lat1, lon1}] defp build_path_points(_, _), do: [] + # Spherical-vector midpoint between two lat/lon points. The arithmetic + # mean of two longitudes folds across the anti-meridian (e.g. + # `(179 + -179) / 2 == 0` back at Greenwich), so paths that cross the + # date line otherwise get a bogus center. + defp great_circle_midpoint(lat1, lon1, lat2, lon2) do + rlat1 = deg_to_rad(lat1) + rlat2 = deg_to_rad(lat2) + rlon1 = deg_to_rad(lon1) + dlon = deg_to_rad(lon2 - lon1) + + bx = :math.cos(rlat2) * :math.cos(dlon) + by = :math.cos(rlat2) * :math.sin(dlon) + + mid_lat = + :math.atan2( + :math.sin(rlat1) + :math.sin(rlat2), + :math.sqrt((:math.cos(rlat1) + bx) ** 2 + by ** 2) + ) + + mid_lon = rlon1 + :math.atan2(by, :math.cos(rlat1) + bx) + + {rad_to_deg(mid_lat), normalize_lon(rad_to_deg(mid_lon))} + end + + defp rad_to_deg(rad), do: rad * 180.0 / :math.pi() + + defp normalize_lon(lon) when lon > 180.0, do: lon - 360.0 + defp normalize_lon(lon) when lon < -180.0, do: lon + 360.0 + defp normalize_lon(lon), do: lon + @earth_radius_km 6371.0 @spec haversine_km(number(), number(), number(), number()) :: float() @@ -389,7 +418,10 @@ defmodule Microwaveprop.Radio do :math.sin(dlat / 2) ** 2 + :math.cos(rlat1) * :math.cos(rlat2) * :math.sin(dlon / 2) ** 2 - 2 * @earth_radius_km * :math.asin(:math.sqrt(a)) + # `atan2(√a, √(1−a))` is the numerically stable form: when `a` + # rounds to slightly above 1 the asin form returns NaN, while + # atan2 stays finite. + 2 * @earth_radius_km * :math.atan2(:math.sqrt(a), :math.sqrt(1 - a)) end @spec backfill_distances([Contact.t()]) :: Postgrex.Result.t() | nil diff --git a/lib/microwaveprop/radio/maidenhead.ex b/lib/microwaveprop/radio/maidenhead.ex index 0064c95c..f20d94da 100644 --- a/lib/microwaveprop/radio/maidenhead.ex +++ b/lib/microwaveprop/radio/maidenhead.ex @@ -53,8 +53,11 @@ defmodule Microwaveprop.Radio.Maidenhead do """ @spec from_latlon(number(), number(), pos_integer()) :: String.t() def from_latlon(lat, lon, precision \\ 6) when is_number(lat) and is_number(lon) and precision in [4, 6, 8] do - lat = lat + 90.0 - lon = lon + 180.0 + # Clamp into the open ranges the field encoder expects. Exactly + # `90.0` / `180.0` would otherwise overflow the trunc into 'S' + # (one past 'R'), producing grids outside the A..R alphabet. + lat = min(max(lat + 90.0, 0.0), 180.0 - 1.0e-9) + lon = min(max(lon + 180.0, 0.0), 360.0 - 1.0e-9) # Field (A-R), 20° lon / 10° lat f1 = trunc(lon / 20.0) diff --git a/lib/microwaveprop/weather.ex b/lib/microwaveprop/weather.ex index 3dd58492..909cbcaa 100644 --- a/lib/microwaveprop/weather.ex +++ b/lib/microwaveprop/weather.ex @@ -366,9 +366,12 @@ defmodule Microwaveprop.Weather do even after their data had landed. This reconciler closes that loop as a single SQL UPDATE, invoked from the hourly enqueuer cron. - Radius encoded as a conservative ±1.5° lat/lon bounding box — the - same rectangle `weather_for_contact/2` uses at `radius_km: 150` - around the equator/mid-latitudes. + Radius encoded as a ±1.5° latitude band; the longitude band is + scaled by `1 / cos(lat)` so the box covers the same physical + east-west distance (~150 km) at every latitude. A fixed 1.5° lon + box collapses to ~110 km at lat 49° and would silently skip + observations the per-contact `weather_for_contact/2` query would + match. """ @spec reconcile_weather_statuses() :: {:ok, non_neg_integer()} def reconcile_weather_statuses do @@ -386,8 +389,8 @@ defmodule Microwaveprop.Weather do AND o.observed_at <= c.qso_timestamp + interval '2 hours' AND s.lat BETWEEN ((c.pos1->>'lat')::float - 1.5) AND ((c.pos1->>'lat')::float + 1.5) - AND s.lon BETWEEN ((c.pos1->>'lon')::float - 1.5) - AND ((c.pos1->>'lon')::float + 1.5) + AND s.lon BETWEEN ((c.pos1->>'lon')::float - 1.5 / GREATEST(cos(radians((c.pos1->>'lat')::float)), 0.01)) + AND ((c.pos1->>'lon')::float + 1.5 / GREATEST(cos(radians((c.pos1->>'lat')::float)), 0.01)) ) OR ( c.pos2 IS NOT NULL @@ -399,8 +402,8 @@ defmodule Microwaveprop.Weather do AND o.observed_at <= c.qso_timestamp + interval '2 hours' AND s.lat BETWEEN ((c.pos2->>'lat')::float - 1.5) AND ((c.pos2->>'lat')::float + 1.5) - AND s.lon BETWEEN ((c.pos2->>'lon')::float - 1.5) - AND ((c.pos2->>'lon')::float + 1.5) + AND s.lon BETWEEN ((c.pos2->>'lon')::float - 1.5 / GREATEST(cos(radians((c.pos2->>'lat')::float)), 0.01)) + AND ((c.pos2->>'lon')::float + 1.5 / GREATEST(cos(radians((c.pos2->>'lat')::float)), 0.01)) ) ) ) diff --git a/lib/microwaveprop/weather/grib2/wgrib2.ex b/lib/microwaveprop/weather/grib2/wgrib2.ex index fc068035..e6650ae0 100644 --- a/lib/microwaveprop/weather/grib2/wgrib2.ex +++ b/lib/microwaveprop/weather/grib2/wgrib2.ex @@ -652,14 +652,24 @@ defmodule Microwaveprop.Weather.Grib2.Wgrib2 do defp extract_var_level(_), do: nil - # Parse "lon=242.958,lat=32.938,val=306.5" - defp parse_lon_val_segment(segment) do + @doc false + # Parse "lon=242.958,lat=32.938,val=306.5". Exposed (with @doc false) + # so the integer-lon edge case (wgrib2 occasionally drops the trailing + # `.0`) can be asserted directly. + def parse_lon_val_segment(segment) do case Regex.run(~r/lon=([\d.]+),lat=([\d.]+),val=([\d.eE+-]+)/, segment) do [_, lon_str, lat_str, val_str] -> - lon = denormalize_lon(String.to_float(lon_str)) - {lat, ""} = Float.parse(lat_str) - {val, ""} = Float.parse(val_str) - {Float.round(lat, 3), Float.round(lon, 3), val} + # `Float.parse/1` accepts both `"243"` and `"243.0"`, where + # `String.to_float/1` only accepts the latter. wgrib2 normally + # emits the trailing `.0` but a build or locale flip that drops + # it would crash the entire chain step. + with {lon_val, ""} <- Float.parse(lon_str), + {lat, ""} <- Float.parse(lat_str), + {val, ""} <- Float.parse(val_str) do + {Float.round(lat, 3), Float.round(denormalize_lon(lon_val), 3), val} + else + _ -> nil + end _ -> nil diff --git a/lib/microwaveprop/weather/hrrr_client.ex b/lib/microwaveprop/weather/hrrr_client.ex index 985578be..fcf7858e 100644 --- a/lib/microwaveprop/weather/hrrr_client.ex +++ b/lib/microwaveprop/weather/hrrr_client.ex @@ -215,16 +215,18 @@ defmodule Microwaveprop.Weather.HrrrClient do end end + @doc """ + Truncates a DateTime down to the nearest published HRRR cycle hour. + + Always rounds DOWN — rounding to the nearest hour used to point at + the next cycle once we passed `:30`, which produced 404/403 errors + when that cycle had not yet been published by NOAA. Callers that need + a later valid_time should reach it via `forecast_hour` instead. + """ @spec nearest_hrrr_hour(DateTime.t()) :: DateTime.t() def nearest_hrrr_hour(dt) do total_seconds = dt.minute * 60 + dt.second - rounded_dt = DateTime.add(dt, -total_seconds, :second) - - if dt.minute >= 30 do - DateTime.add(rounded_dt, 3600, :second) - else - rounded_dt - end + DateTime.add(dt, -total_seconds, :second) end @spec hrrr_url(Date.t(), non_neg_integer(), :surface | :pressure, non_neg_integer()) :: String.t() diff --git a/lib/microwaveprop/weather/sounding_params.ex b/lib/microwaveprop/weather/sounding_params.ex index 976c1930..39586ee0 100644 --- a/lib/microwaveprop/weather/sounding_params.ex +++ b/lib/microwaveprop/weather/sounding_params.ex @@ -140,7 +140,10 @@ defmodule Microwaveprop.Weather.SoundingParams do |> Enum.flat_map(fn [prev, curr] -> dh_km = (curr.h_agl - prev.h_agl) / 1000.0 - if abs(dh_km) < 0.01 do + # Only zero-thickness layers are skipped — the previous 10m guard + # discarded the very thin near-surface levels that HRRR's native + # vertical grid uses to resolve sharp inversions and surface ducts. + if abs(dh_km) < 1.0e-6 do [] else dn = curr.n - prev.n @@ -175,14 +178,17 @@ defmodule Microwaveprop.Weather.SoundingParams do classify_duct_transition(prev, curr, ducts, duct_state) end) - # Handle case where profile ends while still in a duct + # Surface-based ducts that extend past the top of a shallow + # profile still strongly affect ground-to-ground propagation, so + # close the duct using the highest sampled level as the top + # rather than discarding it for "lack of a clear top". ducts = - case in_duct_state do - %{base: _base, base_m: _base_m} -> - # Duct didn't close — discard (no clear top) - ducts + case {in_duct_state, List.last(refract_profile)} do + {%{base: base, base_m: base_m}, %{} = top} -> + {finalized, _} = finalize_duct(ducts, base, base_m, top) + finalized - nil -> + _ -> ducts end diff --git a/rust/prop_grid_rs/src/fetcher.rs b/rust/prop_grid_rs/src/fetcher.rs index 9dc4d28b..c17fead5 100644 --- a/rust/prop_grid_rs/src/fetcher.rs +++ b/rust/prop_grid_rs/src/fetcher.rs @@ -183,14 +183,15 @@ pub fn merge_ranges(mut ranges: Vec) -> Vec { out } +/// Truncates a DateTime down to the nearest published HRRR cycle hour. +/// +/// Always rounds DOWN — rounding up past `:30` used to point at the next +/// cycle, which produced 404/403 errors when that cycle had not yet been +/// published by NOAA. Callers that need a later valid_time should reach +/// it via `forecast_hour` instead. pub fn nearest_hrrr_hour(dt: DateTime) -> DateTime { let total_seconds = dt.minute() as i64 * 60 + dt.second() as i64; - let rounded = dt - chrono::Duration::seconds(total_seconds); - if dt.minute() >= 30 { - rounded + chrono::Duration::hours(1) - } else { - rounded - } + dt - chrono::Duration::seconds(total_seconds) } pub fn pressure_messages_grid() -> Vec<(String, String)> { @@ -722,11 +723,15 @@ mod tests { } #[test] - fn nearest_hour_rounds_properly() { - let base = Utc.with_ymd_and_hms(2026, 4, 19, 15, 29, 0).unwrap(); - assert_eq!(nearest_hrrr_hour(base).hour(), 15); - let up = Utc.with_ymd_and_hms(2026, 4, 19, 15, 30, 0).unwrap(); - assert_eq!(nearest_hrrr_hour(up).hour(), 16); + fn nearest_hour_always_rounds_down() { + // Round-down keeps us from probing a cycle that NOAA has not + // published yet. Callers reach later valid_times via forecast_hour. + let early = Utc.with_ymd_and_hms(2026, 4, 19, 15, 29, 0).unwrap(); + assert_eq!(nearest_hrrr_hour(early).hour(), 15); + let mid = Utc.with_ymd_and_hms(2026, 4, 19, 15, 30, 0).unwrap(); + assert_eq!(nearest_hrrr_hour(mid).hour(), 15); + let late = Utc.with_ymd_and_hms(2026, 4, 19, 15, 59, 0).unwrap(); + assert_eq!(nearest_hrrr_hour(late).hour(), 15); } #[test] diff --git a/rust/prop_grid_rs/src/pipeline.rs b/rust/prop_grid_rs/src/pipeline.rs index c6a5c92b..ae7a3543 100644 --- a/rust/prop_grid_rs/src/pipeline.rs +++ b/rust/prop_grid_rs/src/pipeline.rs @@ -247,12 +247,28 @@ fn cell_to_conditions( .get("PRES:surface") .copied() .map(|pa| pa as f64 / 100.0); - let rain_rate_mmhr = cell + // Mirror Elixir's merged_rain_rate: pick the heavier of HRRR's + // accumulation-derived rate and NEXRAD's reflectivity-derived rate + // so a fast convective cell that hasn't yet shown up in the hourly + // APCP still triggers the rain penalty. + let hrrr_rate = cell .get("APCP:surface") .copied() .map(|mm| mm as f64) - .filter(|v| *v > 0.0); + .unwrap_or(0.0); + let nexrad_rate = scorer::dbz_to_rain_rate_mmhr( + cell.get("nexrad_max_reflectivity_dbz") + .copied() + .map(|v| v as f64), + ); + let merged_rate = hrrr_rate.max(nexrad_rate); + let rain_rate_mmhr = if merged_rate > 0.0 { + Some(merged_rate) + } else { + None + }; let bl_depth_m = cell.get("HPBL:surface").copied().map(|v| v as f64); + let best_duct_band_ghz = cell.get("best_duct_freq_ghz").copied().map(|v| v as f64); let levels: Vec = fetcher::grid_level_keys() .iter() @@ -287,7 +303,7 @@ fn cell_to_conditions( min_refractivity_gradient, bl_depth_m, pwat_mm, - best_duct_band_ghz: None, + best_duct_band_ghz, bulk_richardson: None, }) } @@ -343,6 +359,43 @@ mod tests { assert!(cell_to_conditions(&cell, 32.0, -97.0, &vt).is_none()); } + #[test] + fn nexrad_reflectivity_lifts_rain_rate_above_hrrr_precip() { + // Mirrors the Elixir scorer's rain merge: when HRRR's hourly + // accumulation hasn't caught a fast convective cell but NEXRAD + // sees 45 dBZ overhead, take whichever rate is higher. Without + // this the Rust pipeline silently underestimates rain attenuation. + use chrono::TimeZone; + let mut cell = CellValues::new(); + cell.insert("TMP:2 m above ground".into(), 295.0); + cell.insert("DPT:2 m above ground".into(), 285.0); + cell.insert("APCP:surface".into(), 0.0); + cell.insert("nexrad_max_reflectivity_dbz".into(), 45.0); + + let vt = Utc.with_ymd_and_hms(2026, 6, 15, 18, 0, 0).unwrap(); + let c = cell_to_conditions(&cell, 32.0, -97.0, &vt).unwrap(); + let rate = c.rain_rate_mmhr.expect("nexrad-derived rain rate present"); + assert!(rate > 5.0, "expected >5 mm/hr from 45 dBZ, got {rate}"); + } + + #[test] + fn duct_freq_threads_into_best_duct_band_ghz() { + // After native_duct::merge_duct_grid runs, every cell with a + // detected trapping layer carries `best_duct_freq_ghz`. The + // scorer reads `best_duct_band_ghz` from Conditions for the 15% + // Native Duct Boost. cell_to_conditions must wire one to the + // other or the boost never fires from the Rust pipeline. + use chrono::TimeZone; + let mut cell = CellValues::new(); + cell.insert("TMP:2 m above ground".into(), 295.0); + cell.insert("DPT:2 m above ground".into(), 285.0); + cell.insert("best_duct_freq_ghz".into(), 24.0); + + let vt = Utc.with_ymd_and_hms(2026, 6, 15, 18, 0, 0).unwrap(); + let c = cell_to_conditions(&cell, 32.0, -97.0, &vt).unwrap(); + assert_eq!(c.best_duct_band_ghz, Some(24.0)); + } + #[test] fn rejects_forecast_hour_zero() { use chrono::TimeZone; diff --git a/rust/prop_grid_rs/src/profiles_file.rs b/rust/prop_grid_rs/src/profiles_file.rs index 497aac3e..6f4fae18 100644 --- a/rust/prop_grid_rs/src/profiles_file.rs +++ b/rust/prop_grid_rs/src/profiles_file.rs @@ -136,13 +136,26 @@ pub fn write_atomic( Ok(path) } -/// Round a lat/lon to 3 decimal places — matches Elixir's -/// `Float.round(_, 3)` snap before keying into the profiles map. +/// Snap a lat/lon to the propagation grid step. Mirrors Elixir's +/// `Microwaveprop.Propagation.ProfilesFile.snap/2`: +/// +/// ```elixir +/// Float.round(Float.round(coord / step) * step, 3) +/// ``` +/// +/// Two implementations of "snap" — Rust's old plain 3-decimal round +/// and Elixir's step-aware round — disagreed on every off-grid input, +/// so a coordinate written by Elixir was unreachable from Rust and +/// vice versa. +const GRID_STEP: f64 = 0.125; + pub fn snap_coords(lat: f64, lon: f64) -> (f64, f64) { - ( - (lat * 1000.0).round() / 1000.0, - (lon * 1000.0).round() / 1000.0, - ) + (snap_step(lat), snap_step(lon)) +} + +fn snap_step(coord: f64) -> f64 { + let snapped = (coord / GRID_STEP).round() * GRID_STEP; + (snapped * 1000.0).round() / 1000.0 } /// Convenience helper: build a `rmpv::Value` map from an @@ -239,7 +252,14 @@ mod tests { } #[test] - fn snap_coords_rounds_to_three_decimals() { - assert_eq!(snap_coords(32.8998, -97.04031), (32.900, -97.040)); + fn snap_coords_matches_elixir_step_snap() { + // The Elixir reader writes the profiles map with keys snapped via + // `Float.round(Float.round(coord / 0.125) * 0.125, 3)`. This Rust + // helper has to produce the same keys or `read_point` looks up a + // coordinate that was never written and the UI shows "no data". + assert_eq!(snap_coords(32.8998, -97.04031), (32.875, -97.0)); + assert_eq!(snap_coords(32.07, -97.07), (32.125, -97.125)); + // Round-trip — already-snapped coordinates stay put. + assert_eq!(snap_coords(32.875, -97.0), (32.875, -97.0)); } } diff --git a/rust/prop_grid_rs/src/region.rs b/rust/prop_grid_rs/src/region.rs index 95025ba9..663ca230 100644 --- a/rust/prop_grid_rs/src/region.rs +++ b/rust/prop_grid_rs/src/region.rs @@ -14,7 +14,9 @@ pub enum Region { Other, } -/// Lat min (inclusive), lat max (exclusive), lon min (inclusive), lon max (exclusive). +/// Lat min, lat max, lon min, lon max — all inclusive on both ends so +/// that points sitting exactly on a region's upper boundary still pick +/// up the regional multiplier instead of falling through to `Other`. struct BBox(Region, f64, f64, f64, f64); // Ordered to match Elixir's `Enum.find_value` short-circuit. @@ -32,7 +34,7 @@ const REGIONS: &[BBox] = &[ pub fn for_point(lat: f64, lon: f64) -> Region { for bbox in REGIONS { let BBox(name, lat_min, lat_max, lon_min, lon_max) = *bbox; - if lat >= lat_min && lat < lat_max && lon >= lon_min && lon < lon_max { + if lat >= lat_min && lat <= lat_max && lon >= lon_min && lon <= lon_max { return name; } } @@ -102,4 +104,13 @@ mod tests { // southeast (30..37, -90..-75). Elixir hits gulf_coast first. assert_eq!(for_point(30.5, -90.0), Region::GulfCoast); } + + #[test] + fn includes_upper_boundary_points() { + // CornBelt is 38.0..48.0 lat; PacificNorthwest is 42.0..50.0. + // Half-open bounds dropped points sitting exactly at the upper + // edge into Region::Other and stripped their seasonal multiplier. + assert_ne!(for_point(48.0, -93.0), Region::Other); + assert_ne!(for_point(50.0, -122.0), Region::Other); + } } diff --git a/rust/prop_grid_rs/src/sounding_params.rs b/rust/prop_grid_rs/src/sounding_params.rs index 0e39704c..2eea0eb1 100644 --- a/rust/prop_grid_rs/src/sounding_params.rs +++ b/rust/prop_grid_rs/src/sounding_params.rs @@ -81,7 +81,10 @@ pub fn min_refractivity_gradient(mut levels: Vec) -> Option { let (h0, n0) = pair[0]; let (h1, n1) = pair[1]; let dh_km = (h1 - h0) / 1000.0; - if dh_km.abs() < 0.01 { + // Only zero-thickness layers are skipped. The previous 10 m + // guard discarded the very thin native-level layers HRRR uses + // to resolve sharp surface-based inversions and ducts. + if dh_km.abs() < 1.0e-9 { continue; } let g = (n1 - n0) / dh_km; diff --git a/test/microwaveprop/propagation/region_test.exs b/test/microwaveprop/propagation/region_test.exs index 90fc1e85..239ec0e5 100644 --- a/test/microwaveprop/propagation/region_test.exs +++ b/test/microwaveprop/propagation/region_test.exs @@ -40,6 +40,14 @@ defmodule Microwaveprop.Propagation.RegionTest do test "returns :other for points outside CONUS" do assert Region.for_point(19.0, -155.0) == :other end + + test "includes points sitting exactly on the upper latitude boundary" do + # corn_belt is defined as 38.0..48.0 lat; pacific_northwest as + # 42.0..50.0. The exclusive upper bound dropped points exactly at + # 48.0 / 50.0 into `:other` and stripped their seasonal multiplier. + assert Region.for_point(48.0, -93.0) != :other + assert Region.for_point(50.0, -122.0) != :other + end end describe "seasonal_adjustment/2" do diff --git a/test/microwaveprop/propagation_test.exs b/test/microwaveprop/propagation_test.exs index c1bca60b..ce63106a 100644 --- a/test/microwaveprop/propagation_test.exs +++ b/test/microwaveprop/propagation_test.exs @@ -162,6 +162,35 @@ defmodule Microwaveprop.PropagationTest do assert rain_with_nx == rain_no_nx end + + test "threads latitude into the season factor so regional adjustments fire" do + # August on the Gulf coast (lat ~29) gets the +1.15 regional multiplier + # and Iowa corn belt (lat ~42) gets the 0.80 multiplier on the same + # band/month. If score_grid_point dropped latitude on the floor the two + # season factors would be identical because the underlying base score + # only depends on `month`. + base_profile = %{ + surface_temp_c: 25.0, + surface_dewpoint_c: 18.0, + surface_pressure_mb: 1013.0, + hpbl_m: 500.0, + wind_u: 3.0, + wind_v: 2.0, + cloud_cover_pct: 15.0, + precip_mm: 0.0, + min_refractivity_gradient: -100.0 + } + + valid_time = ~U[2026-08-15 13:00:00Z] + + gulf = Propagation.score_grid_point(base_profile, valid_time, 29.0, -95.0) + iowa = Propagation.score_grid_point(base_profile, valid_time, 42.0, -93.0) + + gulf_season = Enum.find(gulf, &(&1.band_mhz == 10_000)).factors.season + iowa_season = Enum.find(iowa, &(&1.band_mhz == 10_000)).factors.season + + assert gulf_season > iowa_season + end end describe "replace_scores/2" do diff --git a/test/microwaveprop/radio/maidenhead_test.exs b/test/microwaveprop/radio/maidenhead_test.exs index 0e69b530..1f8e33ef 100644 --- a/test/microwaveprop/radio/maidenhead_test.exs +++ b/test/microwaveprop/radio/maidenhead_test.exs @@ -131,5 +131,17 @@ defmodule Microwaveprop.Radio.MaidenheadTest do grid = Maidenhead.from_latlon(32.897, -97.038, 6) assert String.starts_with?(grid, "EM12") end + + test "clamps the upper field boundary at 90 latitude / 180 longitude" do + # Without clamping, lat = 90 + 90 = 180 → trunc(180/10) = 18, which + # the field encoder turns into 'S' — outside the legal A-R range. + # The same overflow happens at lon = 180. Both must round into the + # last valid field instead of producing an invalid grid character. + grid = Maidenhead.from_latlon(90.0, 180.0, 4) + assert Maidenhead.valid?(grid) + <> = grid + assert f1 in ?A..?R + assert f2 in ?A..?R + end end end diff --git a/test/microwaveprop/radio_test.exs b/test/microwaveprop/radio_test.exs index 0b674dc9..e933e180 100644 --- a/test/microwaveprop/radio_test.exs +++ b/test/microwaveprop/radio_test.exs @@ -175,6 +175,15 @@ defmodule Microwaveprop.RadioTest do test "returns 0 for same point" do assert Radio.haversine_km(32.9, -97.0, 32.9, -97.0) == 0.0 end + + test "uses the numerically stable atan2 form" do + # Antipodal points sit on a 6371-km-radius sphere at exactly + # π·r ≈ 20_015 km apart. Both haversine forms agree here, but + # the asin form is fragile when `a` rounds to slightly above 1 + # (sqrt returns >1 and asin returns NaN). atan2 stays finite. + dist = Radio.haversine_km(0.0, 0.0, 0.0, 180.0) + assert_in_delta dist, 20_015, 1.0 + end end describe "backfill_distances/1" do @@ -393,6 +402,22 @@ defmodule Microwaveprop.RadioTest do points = Radio.contact_path_points(contact) assert length(points) == 3 end + + test "computes a great-circle midpoint across the anti-meridian" do + # Two stations on opposite sides of the 180° line. The arithmetic + # mean of -179 and +179 lands at 0° (Greenwich) — clearly wrong; + # the geodesic midpoint sits near the date line at ±180°. + contact = + create_contact(%{ + pos1: %{"lat" => 0.0, "lon" => 179.0}, + pos2: %{"lat" => 0.0, "lon" => -179.0} + }) + + [_, {mid_lat, mid_lon}, _] = Radio.contact_path_points(contact) + + assert_in_delta mid_lat, 0.0, 0.01 + assert abs(mid_lon) > 179.0 + end end describe "change_contact/2" do diff --git a/test/microwaveprop/weather/grib2/wgrib2_test.exs b/test/microwaveprop/weather/grib2/wgrib2_test.exs index f29ae9cd..37440037 100644 --- a/test/microwaveprop/weather/grib2/wgrib2_test.exs +++ b/test/microwaveprop/weather/grib2/wgrib2_test.exs @@ -350,4 +350,28 @@ defmodule Microwaveprop.Weather.Grib2.Wgrib2Test do end end end + + describe "parse_lon_val_segment/1" do + test "parses a typical decimal-formatted segment" do + assert {32.938, lon, 306.5} = + Wgrib2.parse_lon_val_segment("lon=242.958,lat=32.938,val=306.5") + + # `denormalize_lon` flips the 0..360 wgrib2 form into -180..180. + assert_in_delta lon, -117.042, 0.001 + end + + test "tolerates a longitude string with no decimal point" do + # wgrib2 typically formats `lon=242.000` but a build that drops the + # trailing `.0` would crash `String.to_float/1` and bring the entire + # propagation chain step down. The parser should accept either form. + assert {32.938, lon, 306.5} = + Wgrib2.parse_lon_val_segment("lon=243,lat=32.938,val=306.5") + + assert_in_delta lon, -117.0, 0.001 + end + + test "returns nil for unparseable segments" do + assert Wgrib2.parse_lon_val_segment("d=2024092219") == nil + end + end end diff --git a/test/microwaveprop/weather/grid_snap_property_test.exs b/test/microwaveprop/weather/grid_snap_property_test.exs index 44134e7d..bfe4eb16 100644 --- a/test/microwaveprop/weather/grid_snap_property_test.exs +++ b/test/microwaveprop/weather/grid_snap_property_test.exs @@ -70,21 +70,24 @@ defmodule Microwaveprop.Weather.GridSnapPropertyTest do end end - property "output is within 30 minutes of the input" do + property "output is at or before the input by at most one hour" do check all(dt <- datetime_generator()) do result = HrrrClient.nearest_hrrr_hour(dt) - delta_sec = abs(DateTime.diff(result, dt, :second)) - # Input minute < 30 rounds down; >= 30 rounds up. Worst-case - # gap is exactly 30:00. - assert delta_sec <= 30 * 60 + # Round-down only — the result is always the published cycle + # at or before `dt`, never a future cycle that may not exist. + delta_sec = DateTime.diff(dt, result, :second) + assert delta_sec >= 0 + assert delta_sec < 60 * 60 end end - property "rounds minute >= 30 up to next hour" do + property "always rounds down (never points at a future cycle)" do check all(dt <- datetime_generator(minute_range: 30..59)) do result = HrrrClient.nearest_hrrr_hour(dt) - # result hour is strictly greater than input hour (modulo day/month rollover). - assert DateTime.after?(result, dt) + # result is strictly before the input — the input minute > 0 + # ensures a non-zero offset. + assert DateTime.before?(result, dt) + assert result.hour == dt.hour end end end diff --git a/test/microwaveprop/weather/hrrr_client_test.exs b/test/microwaveprop/weather/hrrr_client_test.exs index dc1cbc3e..84940e9a 100644 --- a/test/microwaveprop/weather/hrrr_client_test.exs +++ b/test/microwaveprop/weather/hrrr_client_test.exs @@ -22,7 +22,7 @@ defmodule Microwaveprop.Weather.HrrrClientTest do refute HrrrClient.cycle_available?(~U[2026-04-25 00:00:00Z]) end - test "passes the rounded run_time to the probe" do + test "passes the rounded-down run_time to the probe" do test_pid = self() Application.put_env(:microwaveprop, :hrrr_cycle_available_fn, fn dt -> @@ -31,7 +31,7 @@ defmodule Microwaveprop.Weather.HrrrClientTest do end) HrrrClient.cycle_available?(~U[2026-04-25 00:38:00Z]) - assert_receive {:probed, ~U[2026-04-25 01:00:00Z]} + assert_receive {:probed, ~U[2026-04-25 00:00:00Z]} end end @@ -41,9 +41,12 @@ defmodule Microwaveprop.Weather.HrrrClientTest do ~U[2026-03-28 18:00:00Z] end - test "rounds up when past 30 minutes" do + test "rounds down when past 30 minutes (never points at a future cycle)" do + # Rounding up was the historical behaviour; in production it caused + # 404/403 errors when the next cycle had not yet been published. The + # caller can always reach a later valid_time via `forecast_hour`. assert HrrrClient.nearest_hrrr_hour(~U[2026-03-28 18:45:00Z]) == - ~U[2026-03-28 19:00:00Z] + ~U[2026-03-28 18:00:00Z] end test "stays same when exactly on the hour" do @@ -51,9 +54,9 @@ defmodule Microwaveprop.Weather.HrrrClientTest do ~U[2026-03-28 18:00:00Z] end - test "rounds at exactly 30 minutes" do - result = HrrrClient.nearest_hrrr_hour(~U[2026-03-28 18:30:00Z]) - assert result in [~U[2026-03-28 18:00:00Z], ~U[2026-03-28 19:00:00Z]] + test "rounds down at exactly 30 minutes" do + assert HrrrClient.nearest_hrrr_hour(~U[2026-03-28 18:30:00Z]) == + ~U[2026-03-28 18:00:00Z] end end diff --git a/test/microwaveprop/weather/sounding_params_test.exs b/test/microwaveprop/weather/sounding_params_test.exs index 9bdd366a..fcf36adb 100644 --- a/test/microwaveprop/weather/sounding_params_test.exs +++ b/test/microwaveprop/weather/sounding_params_test.exs @@ -89,6 +89,25 @@ defmodule Microwaveprop.Weather.SoundingParamsTest do assert result.min_refractivity_gradient < 0 end + test "captures sharp sub-10m inversions instead of skipping the layer" do + # HRRR native-level profiles routinely report 0/10/30/50 m thin layers; + # an inversion confined to the lowest 5m is exactly the kind of + # surface-based duct the propagation map cares about, so dropping any + # gradient with `dh < 10 m` quietly hides the strongest signals. + thin_layer_profile = [ + %{"pres" => 1013.0, "hght" => 0.0, "tmpc" => 25.0, "dwpc" => 22.0}, + %{"pres" => 1012.5, "hght" => 5.0, "tmpc" => 24.5, "dwpc" => 5.0}, + %{"pres" => 1010.0, "hght" => 100.0, "tmpc" => 23.0, "dwpc" => 4.0}, + %{"pres" => 950.0, "hght" => 600.0, "tmpc" => 18.0, "dwpc" => 0.0} + ] + + result = SoundingParams.derive(thin_layer_profile) + # The 0→5 m drop in dewpoint creates a strong negative dM/dh and + # therefore a sharp negative refractivity gradient. Without it we + # only see the broad 5→100 m baseline, which is much shallower. + assert result.min_refractivity_gradient < -200 + end + test "detects temperature inversions" do result = SoundingParams.derive(@inversion_profile) assert result.inversions != [] @@ -134,6 +153,24 @@ defmodule Microwaveprop.Weather.SoundingParamsTest do assert result.ducts == [] end + test "keeps a surface duct that has not yet closed by the top of the profile" do + # Three sampled levels with monotonically decreasing modified + # refractivity (M). The profile starts inside a trapping layer + # and the duct never "tops out" before the data runs out. The + # previous code threw the duct away because it never observed + # the layer close — hiding exactly the surface ducts the + # propagation map needs to surface. + open_top_duct = [ + %{"pres" => 1013.0, "hght" => 0.0, "tmpc" => 30.0, "dwpc" => 27.0}, + %{"pres" => 1010.0, "hght" => 20.0, "tmpc" => 33.0, "dwpc" => 5.0}, + %{"pres" => 1005.0, "hght" => 50.0, "tmpc" => 35.0, "dwpc" => -10.0} + ] + + result = SoundingParams.derive(open_top_duct) + assert result.ducting_detected + assert result.ducts != [] + end + test "computes boundary layer depth" do result = SoundingParams.derive(@inversion_profile) # BL depth should be positive and reasonable diff --git a/test/microwaveprop/weather_test.exs b/test/microwaveprop/weather_test.exs index 59146778..fc535151 100644 --- a/test/microwaveprop/weather_test.exs +++ b/test/microwaveprop/weather_test.exs @@ -129,6 +129,40 @@ defmodule Microwaveprop.WeatherTest do assert %{weather_status: :queued} = Repo.get(Contact, contact.id) end + test "scales the longitude band by latitude so high-lat windows still cover ~150 km east/west" do + # At lat 49° (US/Canada border) cos(lat) ≈ 0.656, so a fixed 1.5° + # longitude box collapses to ~110 km — narrower than the 150 km + # window the function documents. A station 1.8° east of the QSO + # is roughly 132 km away, well inside the intended radius, but a + # latitude-blind bounding box would exclude it. + {:ok, station} = + Weather.find_or_create_station(%{ + @station_attrs + | station_code: "KHIGH", + lat: 49.0, + lon: -121.7 + }) + + Weather.upsert_surface_observation(station, %{ + observed_at: ~U[2023-09-17 17:45:00Z], + temp_f: 60.0 + }) + + contact = + create_contact(%{ + station1: "HILAT", + station2: "HILAT2", + pos1: %{"lat" => 49.0, "lon" => -123.5}, + pos2: nil, + grid2: nil + }) + + Radio.set_enrichment_status!([contact.id], :weather_status, :queued) + + {:ok, _n} = Weather.reconcile_weather_statuses() + assert %{weather_status: :complete} = Repo.get(Contact, contact.id) + end + test "leaves :complete and :pending contacts alone" do {:ok, station} = Weather.find_or_create_station(@station_attrs)