diff --git a/lib/microwaveprop/weather/skewt_params.ex b/lib/microwaveprop/weather/skewt_params.ex index 8c796bfb..980f1fd2 100644 --- a/lib/microwaveprop/weather/skewt_params.ex +++ b/lib/microwaveprop/weather/skewt_params.ex @@ -357,18 +357,59 @@ defmodule Microwaveprop.Weather.SkewtParams do end) end - defp wet_bulb_c(level) do - # Stipanuk-style iterative approximation: T_wb satisfies - # e_s(T_wb) = e + (p × c_p / (ε × L_v)) × (T - T_wb) - # Solve by bracketing between dewpoint and dry bulb. - t_d = level.dwpc - t = level.tmpc + @doc """ + Wet-bulb temperature (°C) for a level using a Stipanuk-style + iterative bracket between dewpoint and dry bulb. Returns `nil` if + either temperature is missing. - if t_d == nil or t == nil do - nil - else - bisect_wet_bulb(t, t_d, level.pres) - end + Accepts a level map with `tmpc`, `dwpc`, and `pres` keys (the same + normalized shape `derive/1` consumes). + """ + @spec wet_bulb_c(map()) :: float() | nil + def wet_bulb_c(%{tmpc: t, dwpc: t_d, pres: _p}) when t == nil or t_d == nil, do: nil + + def wet_bulb_c(%{tmpc: t, dwpc: t_d, pres: p}) do + bisect_wet_bulb(t, t_d, p) + end + + def wet_bulb_c(_), do: nil + + @doc """ + Mixing ratio (kg/kg) at a given dewpoint and pressure, via the Buck + saturation-vapor formula evaluated at the dewpoint temperature. + Returns `nil` for non-numeric inputs. + """ + @spec mixing_ratio(number() | nil, number() | nil) :: float() | nil + def mixing_ratio(dewpoint_c, pressure_mb) when is_number(dewpoint_c) and is_number(pressure_mb) do + e = saturation_vapor_pressure_mb(dewpoint_c) + @epsilon * e / max(pressure_mb - e, 1.0e-3) + end + + def mixing_ratio(_, _), do: nil + + @doc """ + Virtual temperature (°C) for a parcel/level. Adds the moisture + correction `(1 + 0.61 r)` to the dry-bulb temperature, where `r` is + the mixing ratio in kg/kg. + + Two calling forms: + + * `virtual_temp_c(temp_c, dewpoint_c, pressure_mb)` — environment + virtual temperature (uses the actual mixing ratio from the + dewpoint). + * `virtual_temp_c(temp_c, pressure_mb)` — parcel virtual temperature + assuming saturation (uses the saturation mixing ratio at the + parcel's own temperature). This is what CAPE uses above the LCL. + """ + @spec virtual_temp_c(number(), number(), number()) :: float() + def virtual_temp_c(temp_c, dewpoint_c, pressure_mb) do + r = mixing_ratio(dewpoint_c, pressure_mb) || 0.0 + (temp_c + 273.15) * (1.0 + 0.61 * r) - 273.15 + end + + @spec virtual_temp_c(number(), number()) :: float() + def virtual_temp_c(temp_c, pressure_mb) do + virtual_temp_c(temp_c, temp_c, pressure_mb) end defp bisect_wet_bulb(t, t_d, _p) when t == t_d, do: t diff --git a/lib/microwaveprop_web/live/skewt_live.ex b/lib/microwaveprop_web/live/skewt_live.ex index ac6939ba..124fe6b7 100644 --- a/lib/microwaveprop_web/live/skewt_live.ex +++ b/lib/microwaveprop_web/live/skewt_live.ex @@ -127,7 +127,13 @@ defmodule MicrowavepropWeb.SkewtLive do # 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] || []) + + svg = + SkewtSvg.render(profile, + parcel_trace: spc[:parcel_trace] || [], + critical_levels: critical_levels(spc, profile) + ) + {profile, derived, svg} end end @@ -137,6 +143,52 @@ defmodule MicrowavepropWeb.SkewtLive do cell[:profile] || cell["profile"] || [] end + # Right-edge critical-level markers, mirroring SPC's annotation rail. + # Pressure → label is enough; SkewtSvg uses pressure_y/1 to position. + # 0 °C and wet-bulb-zero come back from SkewtParams as heights, so + # we walk the profile to recover the corresponding pressure. + defp critical_levels(spc, profile) do + Enum.flat_map( + [ + {"LCL", spc[:lcl_pressure_mb]}, + {"LFC", spc[:lfc_pressure_mb]}, + {"EL", spc[:el_pressure_mb]}, + {"0 °C", pressure_at_height(profile, spc[:freezing_level_m])}, + {"WBZ", pressure_at_height(profile, spc[:wbz_m])} + ], + fn + {label, p} when is_number(p) -> [%{label: label, pressure_mb: p}] + _ -> [] + end + ) + end + + defp pressure_at_height(_profile, nil), do: nil + + defp pressure_at_height(profile, target_m) when is_number(target_m) do + profile + |> Enum.map(&extract_pres_hght/1) + |> Enum.filter(&(is_number(&1.pres) and is_number(&1.hght))) + |> Enum.sort_by(& &1.pres, :desc) + |> Enum.chunk_every(2, 1, :discard) + |> Enum.find_value(&interpolate_pressure_at_height(&1, target_m)) + end + + defp extract_pres_hght(lvl) do + %{ + pres: lvl[:pres] || lvl["pres"] || lvl[:pres_mb] || lvl["pres_mb"], + hght: lvl[:hght] || lvl["hght"] || lvl[:hght_m] || lvl["hght_m"] + } + end + + defp interpolate_pressure_at_height([a, b], target_m) do + cond do + a.hght == b.hght -> nil + (a.hght - target_m) * (b.hght - target_m) > 0 -> nil + true -> a.pres + (target_m - a.hght) / (b.hght - a.hght) * (b.pres - a.pres) + end + 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 diff --git a/lib/microwaveprop_web/live/skewt_svg.ex b/lib/microwaveprop_web/live/skewt_svg.ex index 07f335f5..72d44dca 100644 --- a/lib/microwaveprop_web/live/skewt_svg.ex +++ b/lib/microwaveprop_web/live/skewt_svg.ex @@ -10,6 +10,7 @@ defmodule MicrowavepropWeb.SkewtSvg do AGL, not on the pressure-level profile, so a per-level barb column would be misleading. """ + alias Microwaveprop.Weather.SkewtParams # Plot canvas (SVG user units). @width 720 @@ -42,6 +43,7 @@ defmodule MicrowavepropWeb.SkewtSvg do def render(profile, opts \\ []) when is_list(profile) do cleaned = clean_profile(profile) parcel = Keyword.get(opts, :parcel_trace, []) + critical = Keyword.get(opts, :critical_levels, []) plot_w = @right - @left plot_h = @bottom - @top @@ -61,16 +63,108 @@ defmodule MicrowavepropWeb.SkewtSvg do dry_adiabats(), moist_adiabats(), mixing_ratio_lines(), + virtual_temp_trace(cleaned), + wet_bulb_trace(cleaned), profile_traces(cleaned), parcel_trace_line(parcel), + parcel_virtual_trace(parcel), "", frame(), axis_labels(), + critical_level_markers(critical), cursor_layer(), "" ] end + # ── Wet-bulb / virtual-T traces (NOAA SPC supplementary lines) ───── + + defp wet_bulb_trace([]), do: "" + + defp wet_bulb_trace(profile) do + pts = + profile + |> Enum.sort_by(& &1.pres, :desc) + |> Enum.flat_map(fn lvl -> + case SkewtParams.wet_bulb_c(%{ + tmpc: lvl.tmpc, + dwpc: lvl.dwpc, + pres: lvl.pres + }) do + nil -> [] + t_w -> [{temperature_x(t_w, pressure_y(lvl.pres)), pressure_y(lvl.pres)}] + end + end) + |> clip_polyline() + + polyline_or_empty(pts, "wet-bulb-line") + end + + defp virtual_temp_trace([]), do: "" + + defp virtual_temp_trace(profile) do + pts = + profile + |> Enum.sort_by(& &1.pres, :desc) + |> Enum.map(fn lvl -> + t_v = SkewtParams.virtual_temp_c(lvl.tmpc, lvl.dwpc, lvl.pres) + y = pressure_y(lvl.pres) + {temperature_x(t_v, y), y} + end) + |> clip_polyline() + + polyline_or_empty(pts, "virtual-temp-line") + end + + # Parcel ascent uses the saturated mixing ratio at its own + # temperature above the LCL — that's the virtual-T trajectory used + # in the actual CAPE integration. SPC draws this as a thin red + # dashed line slightly to the right of the dry-bulb parcel trace. + defp parcel_virtual_trace([]), do: "" + + defp parcel_virtual_trace(trace) do + pts = + trace + |> Enum.sort_by(& &1.pres, :desc) + |> Enum.map(fn %{pres: p, tmpc: t} -> + t_v = SkewtParams.virtual_temp_c(t, p) + y = pressure_y(p) + {temperature_x(t_v, y), y} + end) + |> clip_polyline() + + polyline_or_empty(pts, "parcel-virtual-line") + end + + # Critical-level annotations on the right edge — each entry is a + # `%{label, pressure_mb}` map; the label appears outside the frame + # at the right so it doesn't obscure the trace, with a short + # in-plot tick. Mirrors SPC's `143 mb`, `158 mb (mix)`, etc. rail. + defp critical_level_markers([]), do: "" + + defp critical_level_markers(levels) do + Enum.flat_map(levels, fn + %{label: label, pressure_mb: p} when is_number(p) -> + y = pressure_y(p) + + if y >= @top and y <= @bottom do + [render_critical_marker(label, p, y)] + else + [] + end + + _ -> + [] + end) + end + + defp render_critical_marker(label, pressure_mb, y) do + [ + ~s||, + ~s|#{label} #{round(pressure_mb)} mb| + ] + end + # Surface-parcel ascent line (dashed grey), matching the SPC viewer's # convention. The trace is a list of `%{pres, tmpc}` from surface up # through every level above the LCL. @@ -147,6 +241,11 @@ defmodule MicrowavepropWeb.SkewtSvg do .temp-line { fill: none; stroke: #ff0000; stroke-width: 2.5; stroke-linejoin: round; stroke-linecap: round; } .dew-line { fill: none; stroke: #00aa00; stroke-width: 2.5; stroke-linejoin: round; stroke-linecap: round; } .parcel-line { fill: none; stroke: #475569; stroke-width: 1.5; stroke-dasharray: 5 4; stroke-linejoin: round; stroke-linecap: round; } + .parcel-virtual-line { fill: none; stroke: #f87171; stroke-width: 1; stroke-dasharray: 2 3; stroke-linejoin: round; stroke-linecap: round; } + .virtual-temp-line { fill: none; stroke: #fb7185; stroke-width: 1; stroke-dasharray: 1 3; stroke-linejoin: round; stroke-linecap: round; } + .wet-bulb-line { fill: none; stroke: #2563eb; stroke-width: 1.4; stroke-linejoin: round; stroke-linecap: round; opacity: 0.85; } + .critical-tick { stroke: currentColor; stroke-width: 1.2; } + .critical-label { font: 10px ui-sans-serif, system-ui, sans-serif; fill: currentColor; fill-opacity: 0.85; } .skewt-cursor-line { stroke: currentColor; stroke-opacity: 0.5; stroke-width: 0.8; stroke-dasharray: 3 3; pointer-events: none; } .skewt-cursor-dot-t { fill: #ff0000; pointer-events: none; } .skewt-cursor-dot-d { fill: #00aa00; pointer-events: none; }