Graham McIntire
e40ade3b19
Single search bar accepts an address, Maidenhead grid square, or callsign; resolves it to lat/lon via Geocoder / Maidenhead / CallsignLocation (cheapest classification first), snaps to the HRRR grid via the existing ProfilesFile reader, and renders an SVG Skew-T-Log-P with isobars, skewed isotherms, dry/moist adiabats, and mixing-ratio lines. A button row picks between every available forecast hour (now → f18); default is the most recent analysis. Right pane lists derived stability and refractivity stats from SoundingParams.derive (PWAT, BL depth, K-index, lifted index, min dN/dh, ducting status + per-duct base/top/ΔM). Renderer is server-side SVG so the page works without JS and serializes into the LiveView payload as a single tag. Wind barbs are deliberately omitted — HRRR's persisted profile carries wind only at 10 m AGL.
365 lines
11 KiB
Elixir
365 lines
11 KiB
Elixir
defmodule MicrowavepropWeb.SkewtSvg do
|
|
@moduledoc """
|
|
Render an HRRR vertical profile as a Skew-T-Log-P diagram in SVG.
|
|
|
|
HRRR's persisted profile stores `pres`, `hght`, `tmpc`, `dwpc` per
|
|
pressure level — same canonical shape `Microwaveprop.Weather.SoundingParams`
|
|
consumes. This renderer ignores anything else.
|
|
|
|
Wind barbs are deliberately omitted: HRRR persists wind only at 10 m
|
|
AGL, not on the pressure-level profile, so a per-level barb column
|
|
would be misleading.
|
|
"""
|
|
|
|
# Plot canvas (SVG user units).
|
|
@width 720
|
|
@height 640
|
|
@left 60
|
|
@right 660
|
|
@top 30
|
|
@bottom 580
|
|
|
|
@p_bottom 1050.0
|
|
@p_top 100.0
|
|
@t_min -40.0
|
|
@t_max 40.0
|
|
|
|
# Skew factor — pixels of horizontal offset added per pixel of
|
|
# height. 1.0 puts isotherms at 45°.
|
|
@skew 1.0
|
|
|
|
@doc """
|
|
Render a Skew-T diagram for `profile` (a list of `%{"pres", "hght",
|
|
"tmpc", "dwpc"}` maps) and return an inline SVG string suitable for
|
|
`Phoenix.HTML.raw/1`.
|
|
"""
|
|
@spec render([map()]) :: iodata()
|
|
def render(profile) when is_list(profile) do
|
|
cleaned = clean_profile(profile)
|
|
|
|
[
|
|
~s|<svg viewBox="0 0 #{@width} #{@height}" xmlns="http://www.w3.org/2000/svg" class="w-full h-auto" role="img" aria-label="Skew-T-Log-P diagram">|,
|
|
style_block(),
|
|
isobars(),
|
|
isotherms(),
|
|
dry_adiabats(),
|
|
moist_adiabats(),
|
|
mixing_ratio_lines(),
|
|
frame(),
|
|
axis_labels(),
|
|
profile_traces(cleaned),
|
|
"</svg>"
|
|
]
|
|
end
|
|
|
|
defp style_block do
|
|
"""
|
|
<defs>
|
|
<style>
|
|
.frame { fill: none; stroke: currentColor; stroke-width: 1.2; }
|
|
.grid-major { fill: none; stroke: currentColor; stroke-opacity: 0.35; stroke-width: 0.6; }
|
|
.grid-minor { fill: none; stroke: currentColor; stroke-opacity: 0.18; stroke-width: 0.4; }
|
|
.grid-iso { fill: none; stroke: #b91c1c; stroke-opacity: 0.30; stroke-width: 0.5; }
|
|
.grid-dry { fill: none; stroke: #b45309; stroke-opacity: 0.30; stroke-width: 0.4; stroke-dasharray: 4 3; }
|
|
.grid-moist { fill: none; stroke: #047857; stroke-opacity: 0.30; stroke-width: 0.4; stroke-dasharray: 1 4; }
|
|
.grid-mix { fill: none; stroke: #6d28d9; stroke-opacity: 0.30; stroke-width: 0.4; stroke-dasharray: 2 2; }
|
|
.axis { font: 10px ui-sans-serif, system-ui, sans-serif; fill: currentColor; }
|
|
.iso-label { font: 9px ui-sans-serif, system-ui, sans-serif; fill: #b91c1c; fill-opacity: 0.7; }
|
|
.pres-label { font: 9px ui-sans-serif, system-ui, sans-serif; fill: currentColor; fill-opacity: 0.7; }
|
|
.temp-line { fill: none; stroke: #dc2626; stroke-width: 2; stroke-linejoin: round; stroke-linecap: round; }
|
|
.dew-line { fill: none; stroke: #16a34a; stroke-width: 2; stroke-linejoin: round; stroke-linecap: round; }
|
|
</style>
|
|
</defs>
|
|
"""
|
|
end
|
|
|
|
# ── Coordinate transforms ───────────────────────────────────────────
|
|
|
|
@doc false
|
|
def pressure_y(p) when is_number(p) do
|
|
# High pressure → bottom of plot, low pressure → top.
|
|
@bottom -
|
|
(@bottom - @top) * (:math.log(@p_bottom) - :math.log(p)) /
|
|
(:math.log(@p_bottom) - :math.log(@p_top))
|
|
end
|
|
|
|
@doc false
|
|
def temperature_x(t, y) when is_number(t) and is_number(y) do
|
|
plot_w = @right - @left
|
|
base_x = @left + (t - @t_min) / (@t_max - @t_min) * plot_w
|
|
base_x + @skew * (@bottom - y)
|
|
end
|
|
|
|
defp clamp_x(x), do: x |> max(@left) |> min(@right)
|
|
|
|
defp clip_segment({x1, y1}, {x2, y2}) do
|
|
cond do
|
|
x1 < @left and x2 < @left -> nil
|
|
x1 > @right and x2 > @right -> nil
|
|
true -> {{clamp_x(x1), y1}, {clamp_x(x2), y2}}
|
|
end
|
|
end
|
|
|
|
# ── Frame & axes ────────────────────────────────────────────────────
|
|
|
|
defp frame do
|
|
~s|<rect class="frame" x="#{@left}" y="#{@top}" width="#{@right - @left}" height="#{@bottom - @top}" />|
|
|
end
|
|
|
|
defp axis_labels do
|
|
pressures = [1000, 850, 700, 500, 400, 300, 200, 100]
|
|
|
|
pres_labels =
|
|
for p <- pressures do
|
|
y = pressure_y(p)
|
|
~s|<text class="pres-label" x="#{@left - 6}" y="#{y + 3}" text-anchor="end">#{p}</text>|
|
|
end
|
|
|
|
temps = -40..40//10
|
|
|
|
temp_labels =
|
|
for t <- temps do
|
|
x = temperature_x(t, @bottom)
|
|
|
|
if x >= @left and x <= @right do
|
|
~s|<text class="axis" x="#{x}" y="#{@bottom + 14}" text-anchor="middle">#{t}</text>|
|
|
else
|
|
""
|
|
end
|
|
end
|
|
|
|
[
|
|
pres_labels,
|
|
temp_labels,
|
|
~s|<text class="axis" x="#{@left - 38}" y="#{(@top + @bottom) / 2}" text-anchor="middle" transform="rotate(-90 #{@left - 38} #{(@top + @bottom) / 2})">Pressure (mb)</text>|,
|
|
~s|<text class="axis" x="#{(@left + @right) / 2}" y="#{@bottom + 30}" text-anchor="middle">Temperature (°C)</text>|
|
|
]
|
|
end
|
|
|
|
# ── Background grid ─────────────────────────────────────────────────
|
|
|
|
defp isobars do
|
|
pressures = [1000, 925, 850, 700, 500, 400, 300, 250, 200, 150, 100]
|
|
|
|
for p <- pressures do
|
|
y = pressure_y(p)
|
|
cls = if p in [1000, 850, 700, 500, 300, 200, 100], do: "grid-major", else: "grid-minor"
|
|
~s|<line class="#{cls}" x1="#{@left}" y1="#{y}" x2="#{@right}" y2="#{y}" />|
|
|
end
|
|
end
|
|
|
|
defp isotherms do
|
|
# Skewed isotherms every 10°C across the plot. Each line goes from
|
|
# (T at bottom) to (T at top), with skew already baked into
|
|
# temperature_x. Clip to the plot box.
|
|
for t <- -100..40//10 do
|
|
render_isotherm(t)
|
|
end
|
|
end
|
|
|
|
defp render_isotherm(t) do
|
|
p1 = {temperature_x(t, @bottom), @bottom}
|
|
p2 = {temperature_x(t, @top), @top}
|
|
|
|
case clip_segment(p1, p2) do
|
|
nil -> ""
|
|
clipped -> isotherm_with_label(t, clipped)
|
|
end
|
|
end
|
|
|
|
defp isotherm_with_label(t, {{x1, y1}, {x2, y2}}) do
|
|
line = ~s|<line class="grid-iso" x1="#{r(x1)}" y1="#{r(y1)}" x2="#{r(x2)}" y2="#{r(y2)}" />|
|
|
[line, isotherm_label(t, x1)]
|
|
end
|
|
|
|
defp isotherm_label(t, x1) when t in [-40, -20, 0, 20, 40] do
|
|
if x1 > @left + 2 and x1 < @right - 12 do
|
|
~s|<text class="iso-label" x="#{r(x1) + 2}" y="#{@bottom - 2}">#{t}°</text>|
|
|
else
|
|
""
|
|
end
|
|
end
|
|
|
|
defp isotherm_label(_t, _x1), do: ""
|
|
|
|
defp dry_adiabats do
|
|
# Potential-temperature curves: T_K(p) = θ * (p / 1000)^(R/cp).
|
|
# Sample every theta = -30..150°C every 10° at 50-mb pressure steps.
|
|
for theta_c <- -30..150//10 do
|
|
theta_k = theta_c + 273.15
|
|
|
|
pts =
|
|
sample_pressures()
|
|
|> Enum.map(fn p ->
|
|
t_k = theta_k * :math.pow(p / 1000.0, 0.2854)
|
|
t_c = t_k - 273.15
|
|
y = pressure_y(p)
|
|
x = temperature_x(t_c, y)
|
|
{x, y}
|
|
end)
|
|
|> clip_polyline()
|
|
|
|
polyline_or_empty(pts, "grid-dry")
|
|
end
|
|
end
|
|
|
|
defp moist_adiabats do
|
|
# Saturated pseudo-adiabats — quick approximation by integrating
|
|
# dT/dp via the moist-adiabatic lapse rate. Good enough to put
|
|
# green dashes in the right place; not used for parcel ascent.
|
|
for theta_e_c <- 0..40//4 do
|
|
pts =
|
|
(theta_e_c + 273.15)
|
|
|> moist_adiabat_curve()
|
|
|> Enum.map(fn {p, t_c} ->
|
|
y = pressure_y(p)
|
|
x = temperature_x(t_c, y)
|
|
{x, y}
|
|
end)
|
|
|> clip_polyline()
|
|
|
|
polyline_or_empty(pts, "grid-moist")
|
|
end
|
|
end
|
|
|
|
defp moist_adiabat_curve(theta_e_k) do
|
|
# Integrate from 1000 mb upward using dT/dz = Γ_m and dp/dz from
|
|
# hydrostatic. Cheap approximation good enough for plotting.
|
|
starting_t = theta_e_k - 273.15
|
|
|
|
sample_pressures()
|
|
|> Enum.reduce({[], starting_t, 1000.0}, fn p, {acc, t, prev_p} ->
|
|
dp = p - prev_p
|
|
# Δz from hydrostatic with mean T_K ≈ t + 273
|
|
mean_t_k = max(t + 273.15, 200.0)
|
|
dz = -dp * 287.0 * mean_t_k / (9.81 * p)
|
|
gamma_m = moist_lapse_c_per_m(t, p)
|
|
new_t = t + dz * gamma_m
|
|
{[{p, new_t} | acc], new_t, p}
|
|
end)
|
|
|> elem(0)
|
|
|> Enum.reverse()
|
|
end
|
|
|
|
defp moist_lapse_c_per_m(t_c, p_mb) do
|
|
t_k = t_c + 273.15
|
|
# Saturation mixing ratio (g/g) via Buck's eq.
|
|
es = 6.1121 * :math.exp((18.678 - t_c / 234.5) * (t_c / (257.14 + t_c)))
|
|
qs = 0.622 * es / (p_mb - es)
|
|
lv = 2.5e6
|
|
cp = 1004.0
|
|
rd = 287.0
|
|
rv = 461.0
|
|
g = 9.81
|
|
num = g * (1.0 + lv * qs / (rd * t_k))
|
|
den = cp + lv * lv * qs / (rv * t_k * t_k)
|
|
-num / den
|
|
end
|
|
|
|
defp mixing_ratio_lines do
|
|
# Constant-mixing-ratio lines (g/kg). Each line is the locus of
|
|
# (T_d, p) for a given saturation mixing ratio.
|
|
for w <- [0.4, 1, 2, 4, 7, 10, 16, 24, 32] do
|
|
pts =
|
|
sample_pressures()
|
|
|> Enum.filter(&(&1 >= 400))
|
|
|> Enum.map(fn p ->
|
|
# Solve for T from w = 0.622 * e_s(T) / (p - e_s).
|
|
e = w / 1000.0 * p / (0.622 + w / 1000.0)
|
|
t_d = dewpoint_from_vapor(e)
|
|
y = pressure_y(p)
|
|
x = temperature_x(t_d, y)
|
|
{x, y}
|
|
end)
|
|
|> clip_polyline()
|
|
|
|
polyline_or_empty(pts, "grid-mix")
|
|
end
|
|
end
|
|
|
|
# Inverse Buck — find T given vapor pressure e (hPa).
|
|
defp dewpoint_from_vapor(e) when e <= 0, do: -80.0
|
|
|
|
defp dewpoint_from_vapor(e) do
|
|
# Magnus form: T = 243.5 * ln(e/6.112) / (17.67 - ln(e/6.112))
|
|
ratio = :math.log(e / 6.112)
|
|
243.5 * ratio / (17.67 - ratio)
|
|
end
|
|
|
|
defp sample_pressures do
|
|
# 1050 down to 100 in 25-mb steps.
|
|
100..1050//25 |> Enum.to_list() |> Enum.reverse()
|
|
end
|
|
|
|
defp clip_polyline(points) do
|
|
Enum.filter(points, fn {x, y} -> x >= @left - 200 and x <= @right + 200 and y >= @top and y <= @bottom end)
|
|
end
|
|
|
|
defp polyline_or_empty([], _cls), do: ""
|
|
|
|
defp polyline_or_empty(points, cls) do
|
|
pts =
|
|
Enum.map_join(points, " ", fn {x, y} -> "#{r(x)},#{r(y)}" end)
|
|
|
|
~s|<polyline class="#{cls}" points="#{pts}" />|
|
|
end
|
|
|
|
# ── Profile traces ──────────────────────────────────────────────────
|
|
|
|
defp profile_traces([]), do: ""
|
|
|
|
defp profile_traces(profile) do
|
|
sorted = Enum.sort_by(profile, & &1.pres, :desc)
|
|
|
|
temp_pts =
|
|
sorted
|
|
|> Enum.map(fn p ->
|
|
y = pressure_y(p.pres)
|
|
x = temperature_x(p.tmpc, y)
|
|
{x, y}
|
|
end)
|
|
|> Enum.filter(fn {_x, y} -> y >= @top and y <= @bottom end)
|
|
|
|
dew_pts =
|
|
sorted
|
|
|> Enum.filter(& &1.dwpc)
|
|
|> Enum.map(fn p ->
|
|
y = pressure_y(p.pres)
|
|
x = temperature_x(p.dwpc, y)
|
|
{x, y}
|
|
end)
|
|
|> Enum.filter(fn {_x, y} -> y >= @top and y <= @bottom end)
|
|
|
|
[
|
|
polyline_or_empty(dew_pts, "dew-line"),
|
|
polyline_or_empty(temp_pts, "temp-line")
|
|
]
|
|
end
|
|
|
|
defp clean_profile(profile) do
|
|
profile
|
|
|> Enum.map(&normalize_level/1)
|
|
|> Enum.filter(&(&1.pres && &1.tmpc))
|
|
end
|
|
|
|
defp normalize_level(%{} = level) do
|
|
%{
|
|
pres: pick_num(level, ["pres", "pres_mb", :pres, :pres_mb]),
|
|
hght: pick_num(level, ["hght", "hght_m", :hght, :hght_m]),
|
|
tmpc: pick_num(level, ["tmpc", :tmpc]),
|
|
dwpc: pick_num(level, ["dwpc", :dwpc])
|
|
}
|
|
end
|
|
|
|
defp pick_num(map, keys) do
|
|
keys
|
|
|> Enum.find_value(fn k -> Map.get(map, k) end)
|
|
|> coerce_num()
|
|
end
|
|
|
|
defp coerce_num(n) when is_number(n), do: n * 1.0
|
|
defp coerce_num(_), do: nil
|
|
|
|
defp r(n) when is_float(n), do: Float.round(n, 1)
|
|
defp r(n), do: n
|
|
end
|