prop/lib/microwaveprop_web/live/skewt_live.ex
Graham McIntire be1174914f
feat(skewt): NOAA SPC line colors + sounding-parameter parity
Two changes that bring /skewt visually and analytically into parity
with the NOAA SPC sounding viewer
(https://www.spc.noaa.gov/exper/soundings/):

NOAA-style line set
-------------------

Updated the SkewtSvg style block to match the public SPC chart:

  * Temperature trace ........ pure red (#ff0000) at 2.5 px stroke
  * Dewpoint trace ........... pure green (#00aa00) at 2.5 px stroke
  * Isotherms ................ dashed teal (#14b8a6)
  * Dry adiabats ............. solid orange (#f59e0b)
  * Moist (saturated) adiabats dashed sky-blue (#38bdf8)
  * Mixing-ratio lines ....... dotted darker green (#16a34a), with
                               labels at 700 mb showing the g/kg
                               value, drawn only below 600 mb (matches
                               the SPC chart's lower-band number row)
  * Parcel ascent ............ new dashed grey (#475569) line, drawn
                               whenever SkewtParams.derive returns a
                               non-empty parcel_trace

The mixing-ratio set is the SPC standard (0.4 / 0.7 / 1 / 2 / 3 / 5 /
8 / 12 / 16 / 20 g/kg) instead of the prior arbitrary set.

SkewtParams: SPC parcel-theory + indices module
-----------------------------------------------

New `Microwaveprop.Weather.SkewtParams` derives the parameter set that
ships in the SPC `*.txt` companion file:

  * LCL pressure / temperature / height (Bolton 1980 eq 22)
  * LFC pressure / height
  * EL pressure / height
  * SBCAPE, SBCIN (J/kg, virtual-temperature integration with
    surface parcel)
  * 3 km CAPE
  * Lifted Index (parcel ↔ environment T at 500 mb)
  * Freezing level (T = 0 °C height) and Wet-Bulb Zero
  * DCAPE (downdraft CAPE — picks the min-θₑ source level in
    400–700 mb, descends moist-adiabatically to surface)
  * Layer lapse rates: sfc–3 km, 700–500 mb, 850–500 mb
  * `parcel_trace` — pressure/temperature trajectory used by
    SkewtSvg to draw the dashed parcel-ascent overlay

Wind-derived indices (SRH, BWD, Bunkers motion, SCP, STP, SHIP) are
deliberately *not* produced and the LiveView surfaces a one-line note
explaining why: HRRR persists wind only at 10 m AGL, so per-pressure-
level shear can't be computed from this data source.

LiveView panel restructured into grouped sections (Surface, Convective
parcel, Levels, Lapse rates, Moisture/stability, Refractivity) so the
layout matches the SPC viewer's order. The Levels rows show pressure
*and* height ("872 mb · 1,140 m") for quick cross-reference between
the diagram axis and the readout.

Tests: 8/8 new SkewtParamsTest cases (LCL Bolton check, saturated-air
edge case, full field-set presence, summer profile produces SBCAPE > 0
and LI < 0, capped profile produces SBCAPE = 0 and LI > 0, freezing
level interpolated correctly, sfc–3 km lapse rate, empty profile no
crash). 25/25 across all skewt-related test files. mix test green
(2,902 / 2,902 + 221 properties), mix credo --strict clean, mix
assets.build clean.
2026-04-25 13:26:41 -05:00

424 lines
14 KiB
Elixir
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

defmodule MicrowavepropWeb.SkewtLive do
@moduledoc """
`/skewt` interactive Skew-T-Log-P diagram.
The user types an address, Maidenhead grid square, or callsign into
the search bar; the page resolves it to lat/lon, snaps to the HRRR
grid, and renders the latest analysis profile as a Skew-T plus the
derived stability / refractivity stats. A time selector switches
between the analysis hour and each available forecast hour out to
HRRR's f18 horizon.
"""
use MicrowavepropWeb, :live_view
alias Microwaveprop.Propagation.ProfilesFile
alias Microwaveprop.Weather.HrrrProfileLookup
alias Microwaveprop.Weather.SkewtParams
alias Microwaveprop.Weather.SoundingParams
alias MicrowavepropWeb.SkewtLocationResolver
alias MicrowavepropWeb.SkewtSvg
@impl true
def mount(_params, _session, socket) do
{:ok,
assign(socket,
page_title: "Skew-T",
query: "",
location: nil,
error: nil,
valid_times: [],
selected_valid_time: nil,
profile: nil,
derived: nil,
svg: nil
)}
end
@impl true
def handle_params(params, _uri, socket) do
query = String.trim(params["q"] || "")
requested_time = parse_time(params["valid_time"])
{:noreply, refresh(socket, query, requested_time)}
end
@impl true
def handle_event("search", %{"q" => q}, socket) do
{:noreply, push_patch(socket, to: ~p"/skewt?#{[q: q]}")}
end
def handle_event("select_time", %{"valid_time" => iso}, socket) do
params = [q: socket.assigns.query, valid_time: iso]
{:noreply, push_patch(socket, to: ~p"/skewt?#{params}")}
end
# ── Refresh logic ──────────────────────────────────────────────────
defp refresh(socket, "", _time) do
assign(socket,
query: "",
location: nil,
error: nil,
valid_times: [],
selected_valid_time: nil,
profile: nil,
derived: nil,
svg: nil
)
end
defp refresh(socket, query, requested_time) do
case SkewtLocationResolver.resolve(query) do
{:ok, %{lat: lat, lon: lon} = location} ->
valid_times = available_valid_times(lat, lon)
chosen = pick_valid_time(valid_times, requested_time)
{profile, derived, svg} =
case chosen do
nil ->
{nil, nil, nil}
time ->
load_profile(time, lat, lon)
end
assign(socket,
query: query,
location: location,
error: nil,
valid_times: valid_times,
selected_valid_time: chosen,
profile: profile,
derived: derived,
svg: svg
)
{:error, reason} ->
assign(socket,
query: query,
location: nil,
error: reason,
valid_times: [],
selected_valid_time: nil,
profile: nil,
derived: nil,
svg: nil
)
end
end
defp load_profile(valid_time, lat, lon) do
cell = ProfilesFile.read_point(valid_time, lat, lon) || HrrrProfileLookup.read_point_near(valid_time, lat, lon)
case cell do
nil ->
{nil, nil, nil}
%{} = cell ->
profile = extract_profile(cell)
if profile == [] do
{nil, nil, nil}
else
sounding = SoundingParams.derive(profile)
spc = SkewtParams.derive(profile)
# Merge so the LiveView template only has to read one map.
# SoundingParams keys (atoms) and SkewtParams keys (atoms)
# are disjoint by construction.
derived = Map.merge(sounding || %{}, spc)
svg = SkewtSvg.render(profile, parcel_trace: spc[:parcel_trace] || [])
{profile, derived, svg}
end
end
end
defp extract_profile(cell) do
cell[:profile] || cell["profile"] || []
end
defp available_valid_times(lat, lon) do
# The chain produces analysis + 18 forecast hours every wall-clock
# hour. By the time the next chain finishes the previous chain's
# analysis is up to 70 min old, and during a missed cycle it can be
# ~2 h old. A 1-hour past cutoff would leave the page empty while
# the chain catches up — keep 3 h so a stale analysis still shows.
now = DateTime.utc_now()
past_cutoff = DateTime.add(now, -3 * 3600, :second)
future_cutoff = DateTime.add(now, 18 * 3600, :second)
on_disk =
Enum.filter(ProfilesFile.list_valid_times(), fn t ->
DateTime.compare(t, past_cutoff) != :lt and DateTime.compare(t, future_cutoff) != :gt
end)
case on_disk do
[] ->
# On-disk store is empty for this window — common in dev (no NFS
# mount, no local Rust pipeline) and during missed-cycle windows
# in prod. Fall back to whatever HRRR profiles the DB has near
# the location, regardless of how old. Better to render historic
# data than a "no profiles stored" stub.
HrrrProfileLookup.list_valid_times_near(lat, lon)
times ->
Enum.sort(times, DateTime)
end
end
defp pick_valid_time([], _requested), do: nil
defp pick_valid_time(times, nil) do
# Default = the most recent valid_time at or before "now". If
# everything is in the future (cold start), fall back to the
# earliest available.
now = DateTime.utc_now()
times
|> Enum.filter(fn t -> DateTime.compare(t, now) != :gt end)
|> case do
[] -> List.first(times)
past -> Enum.max(past, DateTime)
end
end
defp pick_valid_time(times, %DateTime{} = requested) do
case Enum.find(times, &(DateTime.compare(&1, requested) == :eq)) do
nil -> pick_valid_time(times, nil)
match -> match
end
end
defp parse_time(nil), do: nil
defp parse_time(""), do: nil
defp parse_time(iso) when is_binary(iso) do
case DateTime.from_iso8601(iso) do
{:ok, dt, _} -> DateTime.truncate(dt, :second)
_ -> nil
end
end
# ── Render ─────────────────────────────────────────────────────────
@impl true
def render(assigns) do
~H"""
<Layouts.app flash={@flash} current_scope={@current_scope} max_width="max-w-6xl">
<.header>
Skew-T-Log-P
<:subtitle>
Vertical atmospheric profile from the latest HRRR analysis (or any
available forecast hour) for an address, grid square, or callsign.
</:subtitle>
</.header>
<form phx-submit="search" class="mt-4 flex gap-2">
<input
type="text"
name="q"
value={@query}
placeholder="EM12kp, W5ISP, or '123 Main St, Plano TX'"
class="input input-bordered flex-1"
autocomplete="off"
autofocus
/>
<button type="submit" class="btn btn-primary">Plot</button>
</form>
<%= if @error do %>
<div role="alert" class="alert alert-error mt-4">
<p>{@error}</p>
</div>
<% end %>
<%= if @location do %>
<div class="mt-4 flex flex-wrap items-baseline gap-4 text-sm">
<span class="font-semibold">{@location.label}</span>
<span class="opacity-70">
{fmt_coord(@location.lat)}, {fmt_coord(@location.lon)}
</span>
</div>
<% end %>
<%= if @valid_times != [] do %>
<div class="mt-4 flex flex-wrap gap-2">
<%= for t <- @valid_times do %>
<button
type="button"
phx-click="select_time"
phx-value-valid_time={DateTime.to_iso8601(t)}
class={[
"btn btn-xs",
if(@selected_valid_time && DateTime.compare(t, @selected_valid_time) == :eq,
do: "btn-primary",
else: "btn-outline"
)
]}
>
{time_button_label(t, @valid_times)}
</button>
<% end %>
</div>
<% end %>
<%= if @svg do %>
<div class="mt-6 grid grid-cols-1 lg:grid-cols-3 gap-6">
<div
id="skewt-container"
phx-hook="Skewt"
phx-update="ignore"
data-profile={Jason.encode!(@profile)}
data-geometry={Jason.encode!(SkewtSvg.geometry())}
class="lg:col-span-2 bg-base-100 rounded-box p-3 border border-base-300"
>
{Phoenix.HTML.raw(@svg)}
</div>
<div class="bg-base-200 rounded-box p-4 text-sm space-y-4">
<h3 class="font-semibold text-base">Sounding parameters</h3>
<%= if @derived do %>
<%= for {section_title, rows} <- derived_sections(@derived) do %>
<div class="space-y-1">
<h4 class="font-semibold text-xs uppercase tracking-wide opacity-70">
{section_title}
</h4>
<dl class="grid grid-cols-2 gap-x-3 gap-y-1">
<%= for {label, value} <- rows do %>
<dt class="opacity-70">{label}</dt>
<dd class="font-mono text-right">{value}</dd>
<% end %>
</dl>
</div>
<% end %>
<%= if @derived[:ducting_detected] && @derived[:ducts] != [] do %>
<div class="pt-2 border-t border-base-300 space-y-1">
<h4 class="font-semibold text-xs uppercase tracking-wide opacity-70">
Ducts
</h4>
<%= for duct <- @derived.ducts do %>
<div class="font-mono text-xs">
base {round(duct["base"])} m → top {round(duct["top"])} m, ΔM {duct["strength"]}
</div>
<% end %>
</div>
<% end %>
<p class="pt-2 border-t border-base-300 text-xs opacity-60 leading-relaxed">
Wind-derived indices (SRH, BWD, Bunkers motion, SCP/STP/SHIP)
aren't shown — HRRR persists wind only at 10 m AGL, so per-
level shear can't be computed from this data source.
</p>
<% end %>
</div>
</div>
<% else %>
<%= if @location && @selected_valid_time do %>
<div role="alert" class="alert alert-warning mt-6">
<p>
No HRRR profile is stored for this point at {fmt_time(@selected_valid_time)}.
Profiles are written by the hourly chain and only kept for the active
forecast window.
</p>
</div>
<% else %>
<%= if @location && @valid_times == [] do %>
<div role="alert" class="alert alert-info mt-6">
<p>
No HRRR profiles are currently stored. Wait for the next hourly
chain to publish, or check the propagation pipeline status.
</p>
</div>
<% end %>
<% end %>
<% end %>
</Layouts.app>
"""
end
defp fmt_coord(c) when is_number(c), do: :erlang.float_to_binary(c * 1.0, decimals: 3)
defp fmt_coord(_), do: ""
defp fmt_time(%DateTime{} = t) do
Calendar.strftime(t, "%Y-%m-%d %H:%MZ")
end
defp time_button_label(%DateTime{} = t, all_times) do
base = Enum.min(all_times, DateTime)
delta_h = div(DateTime.diff(t, base, :second), 3600)
IO.iodata_to_binary("f#{:io_lib.format("~2..0B", [delta_h])} · #{Calendar.strftime(t, "%H:%MZ")}")
end
# Grouped derived-parameter sections rendered on the right side of
# the page. Mirrors the SPC viewer's "parcel / level / lapse rate /
# moisture / refractivity" layout so users moving between the two
# tools don't have to context-switch.
defp derived_sections(d) do
[
{"Surface",
[
{"T", fmt_temp(d[:surface_temp_c])},
{"Td", fmt_temp(d[:surface_dewpoint_c])},
{"P", fmt_num(d[:surface_pressure_mb], "mb")}
]},
{"Convective parcel",
[
{"SBCAPE", fmt_num(d[:sbcape], "J/kg")},
{"SBCIN", fmt_num(d[:sbcin], "J/kg")},
{"3 km CAPE", fmt_num(d[:cape_3km], "J/kg")},
{"Lifted index", fmt_num(d[:lifted_index], "")},
{"DCAPE", fmt_num(d[:dcape], "J/kg")}
]},
{"Levels",
[
{"LCL", fmt_pressure_height(d[:lcl_pressure_mb], d[:lcl_height_m])},
{"LFC", fmt_pressure_height(d[:lfc_pressure_mb], d[:lfc_height_m])},
{"EL", fmt_pressure_height(d[:el_pressure_mb], d[:el_height_m])},
{"0 °C", fmt_height(d[:freezing_level_m])},
{"Wet-bulb 0", fmt_height(d[:wbz_m])}
]},
{"Lapse rates",
[
{"Sfc3 km", fmt_num(d[:lapse_rate_sfc_3km_c_per_km], "°C/km")},
{"700500 mb", fmt_num(d[:lapse_rate_700_500_c_per_km], "°C/km")},
{"850500 mb", fmt_num(d[:lapse_rate_850_500_c_per_km], "°C/km")}
]},
{"Moisture / stability",
[
{"PWAT", fmt_num(d[:precipitable_water_mm], "mm")},
{"BL depth", fmt_num(d[:boundary_layer_depth_m], "m")},
{"K-index", fmt_num(d[:k_index], "")}
]},
{"Refractivity",
[
{"Surface N", fmt_num(d[:surface_refractivity], "")},
{"Min dN/dh", fmt_num(d[:min_refractivity_gradient], "/km")},
{"Ducting", if(d[:ducting_detected], do: "yes", else: "no")}
]}
]
end
defp fmt_temp(nil), do: ""
defp fmt_temp(v) when is_number(v), do: "#{Float.round(v * 1.0, 1)} °C"
defp fmt_num(nil, _suffix), do: ""
defp fmt_num(v, ""), do: format_value(v)
defp fmt_num(v, suffix), do: "#{format_value(v)} #{suffix}"
defp fmt_height(nil), do: ""
defp fmt_height(v) when is_number(v) do
"#{round(v)} m / #{round(v * 3.28084)} ft"
end
defp fmt_pressure_height(nil, _), do: ""
defp fmt_pressure_height(p, nil), do: "#{round(p)} mb"
defp fmt_pressure_height(p, h) when is_number(p) and is_number(h) do
"#{round(p)} mb · #{round(h)} m"
end
defp format_value(v) when is_float(v), do: v |> Float.round(1) |> :erlang.float_to_binary(decimals: 1)
defp format_value(v) when is_integer(v), do: Integer.to_string(v)
defp format_value(v), do: to_string(v)
end