Example: Runway Crosswind & Wind-Shift Operability

An airport decides whether a runway is usable over the coming day. Two things close it: the crosswind component of the wind (its speed times the sine of the angle between the wind and the runway heading) exceeding an aircraft limit, and abrupt wind shifts between hours. Both are driven by a wind that evolves over the day — and wind direction is a circular quantity that wanders across North and often switches between two prevailing regimes. The full runnable source is examples/runway_crosswind/model.py.

This example is the process side of prophys end-to-end: a coupled stochastic wind, a circular direction with a bimodal step, and operability read off the marginal of the day’s worst crosswind — its expectation, its tail, and the exceedance probability that decides downtime.

Step 1 — The wind as a coupled process

Wind speed is a mean-reverting Ornstein–Uhlenbeck process. Wind direction is a circular_autoregressive() process reverting to the prevailing heading, with a two-component Gaussian mixture innovation for the backing/veering regimes — a step distribution that is bimodal, not Gaussian. Sampling such a process is exactly what a mixture’s derived quantile makes possible (a bare mixture had no quantile before):

import prophys as prp

speed = prp.ornstein_uhlenbeck(init=8.0, mean=8.0, theta=0.5, sigma=2.5, n_steps=24)

direction = prp.circular_autoregressive(
    init=250.0, mean=250.0, phi=0.85,
    innovation=prp.Mixture([prp.Gaussian(-6.0, 4.0), prp.Gaussian(6.0, 4.0)],
                          logits=[0.0, 0.0]),
    n_steps=24, period=360.0,
)

The circular geometry matters: reverting toward the prevailing 250° from a direction near North follows the short arc across the 0/360 seam, and every sampled hour stays a valid heading in [0, 360) — visible as the smooth wander in the left panel above.

The two are coupled into one JointProcess so gusty hours and shifting hours co-occur, then exposed as a single joint random variable whose components share the underlying draw:

wind_process = prp.JointProcess([speed, direction], corr_raw=[0.35])
wind = prp.JointRandomVariable("wind", wind_process, components=["speed", "direction"])

speed_path = prp.clip(wind["speed"], 0.0, None)   # (..., 24) m/s
direction_path = wind["direction"]                 # (..., 24) deg

Step 2 — Operability read off the wind

The instantaneous crosswind on a runway of a given heading is \(V\,\lvert\sin(\theta - \text{heading})\rvert\), and the operability question is about the worst crosswind of the day — a maximum along the time axis, observed through small anemometer noise so it is a proper UncertainAttribute:

def crosswind_component(heading_deg):
    return prp.abs_(speed_path * prp.sin(prp.deg2rad(direction_path - heading_deg)))

max_crosswind = prp.UncertainAttribute(
    "max_crosswind_07_25",
    prp.Gaussian(mean=time_max(crosswind_component(250.0)), sigma=anemometer_sigma),
    unit="m/s",
)

The largest hour-to-hour wind shift is a circular first difference, so a swing across North is measured as the short arc rather than a spurious near-360° jump:

from prophys.processes import Circular
hourly_shift = prp.time_diff(direction_path, geometry=Circular(360.0))

Step 3 — Expectation, tail, and downtime probability

Compiling marginalizes the wind process by Monte Carlo. For each runway the model reports the expected worst crosswind, its Conditional Value-at-Risk (the mean of the worst 5 % of days — the number that actually governs operations), and the downtime probability that the day’s worst crosswind exceeds the limit. The exceedance is computed two ways that agree: directly by Monte Carlo, and from a Gaussian fitted to the sampled maxima via its first-class cdf (\(P(X > \ell) = 1 - F(\ell)\)):

compiled = model.compile(mode="eval", n_samples=4000)

e     = compiled.expectation("max_crosswind_07_25")
cvar  = compiled.cvar("max_crosswind_07_25", alpha=0.95)
p_mc  = 1.0 - compiled.prob("max_crosswind_07_25", CROSSWIND_LIMIT)

samples = compiled.sample("max_crosswind_07_25", key, (4000,))
fitted  = prp.Gaussian(float(samples.mean()), float(samples.std()))
p_cdf   = 1.0 - float(fitted.cdf(CROSSWIND_LIMIT).evaluate({}))

For the modelled airport the runway aligned with the prevailing wind sees a low expected crosswind and essentially no downtime, while the near-perpendicular runway straddles the limit with ~50 % downtime (the two distributions in the right panel above) — so the model picks the aligned runway.

Primitives exercised

Category	Used for
ornstein_uhlenbeck(), circular_autoregressive()	Wind speed (mean-reverting) and direction (circular, seam-correct)
Mixture	Bimodal backing/veering hourly step — usable as an innovation via its derived quantile
JointProcess, RandomVariable	Coupling gusty and shifting hours; exposing a shared joint draw
time_diff() with Circular	Seam-safe hour-to-hour wind-shift magnitude
expectation / cvar / prob and Gaussian cdf/quantile	Expected & tail crosswind; downtime probability by Monte Carlo and by CDF

Running it

python examples/runway_crosswind/model.py