Nephele Weather-Sensing Mathematical Model Sheet

Version: 2026-07-14 Scope: LEO/GNSS-to-ground geometry, atmospheric propagation, complex RF reception, multilayer weather attribution, and 3-D/4-D reconstruction.

This is the governing-math reference for the Nephele weather-sensing goal. It connects the current interactive digital twin to the phase-accurate scientific model required for calibrated field measurements.

0. Scientific interpretation and status labels

Nephele presently combines public catalog orbits, public numerical-weather-model fields, explicit RF assumptions, and synthetic tomographic observations. It is a physically informed simulator, not yet a live atmospheric retrieval instrument.

• [Implemented] Equation is evaluated in the current application.
• [Approximation] Useful engineering or visualization model, but not a calibrated observable.
• [Scientific target] Required for phase-accurate field measurements or defensible retrievals.
• [Measured input required] Cannot be inferred from public orbit/weather catalogs alone.

The end-to-end chain is

where the recovered atmospheric state may include water vapor, liquid water, ice, rain, and wind, together with posterior uncertainty.

Quick equation index

1. Symbols, units, and sign convention

Nephele should use SI units internally except where an ITU recommendation explicitly uses GHz, km, hPa, mm/h, or g/m$^3$.

Phasor convention

Use

with propagation factor $e^{-jk_0\tilde n L}$. Then positive excess path produces negative phase and positive $\kappa$ produces attenuation:

For a loss $A$ stated in power dB,

Do not mix field-amplitude and power ratios.

2. Orbit, Earth, station, and line-of-sight geometry

2.1 Catalog-orbit propagation

[Implemented] For a catalog OMM/TLE record,

Nephele uses satellite.js, converts the propagated state to Earth-fixed coordinates, and interpolates between orbit anchor samples spaced by 5 s. The globe can animate every display frame, but this does not mean that a satellite is newly observed every frame. Catalog elements are a propagated estimate and are not phase-grade ephemerides.

[Scientific target] Use precise orbit, clock, attitude, antenna phase-center, and Earth-orientation products. Solve light time rather than using same-epoch Euclidean positions:

Iterate until $t_{tx}$ and $\tau$ converge. Earth rotation/Sagnac and relativistic terms must use one consistent reference-frame convention.

2.2 WGS-84 geodetic coordinates to ECEF

For geodetic latitude $\varphi$, longitude $\ell$, and ellipsoidal height $h$,

A nominally fixed ground station has $\mathbf v_g=0$ in ECEF, but in an inertial frame

before loading, tides, and tectonic velocity are added.

2.3 Range, range rate, elevation, and azimuth

In a common frame,

[Implemented] The current twin estimates $\dot R$ by a centered finite difference of Earth-fixed range.

The ECEF-to-local-ENU rotation at the station is

The geometric line of sight is above the ideal horizon when $E>0$. A practical station also needs a terrain/building horizon mask and a link-budget threshold.

2.4 Satellite and ground antenna orientation

Let $\hat{\mathbf b}_t$ and $\hat{\mathbf b}_r$ be transmitter and receiver boresights. For the propagation directions $\hat{\mathbf u}_t$ and $\hat{\mathbf u}_r$,

[Approximation, implemented] Nephele uses the smooth reference pattern

This is 3 dB down at half the stated full beamwidth and has a 30 dB floor. Replace it with measured, frequency-dependent, polarized patterns when available.

For polarization vectors $\mathbf e_t$ and $\mathbf e_r$,

The simple linear-polarization case is $\eta_{pol}=\cos^2\psi$.

2.5 Ray/shell and ray/voxel intersections

For a straight ray

the crossing of a spherical height shell $r_h=R_E+h$ satisfies

The physical root inside $[0,R]$ gives the shell pierce point. Intersections with the lower and upper shells give each layer's slant length $L_k$.

For a locally horizontal slab and elevation $E$,

[Implemented] The selected-link gas, rain, and cloud models use the exact straight-ray distance between spherical height shells. The 3-D array operator clips each straight local-ENU ray against every crossed voxel and stores exact intersection lengths in sparse CSR form. The LOCAL_ENU approximation is separate from the intersection algorithm: production inversion quantifies corner curvature and rejects domains above its declared 500 m tangent-plane tolerance; the coarse browser OSSE allows at most 1,100 m and displays the actual sag.

2.6 Refracted ray

[Scientific target] In a spatially varying refractive index, the geometric-optics ray obeys

At low elevation, use a bent-ray shooting solution through the 3-D refractivity field. Straight rays become increasingly biased near strong boundary-layer gradients, ducting, terrain, and multipath.

3. Unified complex RF channel

3.1 Direct one-way channel

For one direct path, the complex baseband voltage transfer is

Therefore

or in dB,

The free-space path loss is

equivalently, for $f$ in GHz and $R$ in km,

3.2 Multipath and wideband channel

For resolved paths $p$,

The channel impulse response is

A single scalar phase is insufficient for a wideband OFDM signal near dispersive gas lines or in the ionosphere. Evaluate $H(f,t)$ on the actual pilot/subcarrier frequencies.

3.3 Atmospheric amplitude and phase along one ray

Write the total power attenuation and excess phase path as

The unwrapped propagation phase is

For phase/group index $n_\phi$ and group index $n_g$,

Phase delay and group delay are not interchangeable in dispersive media.

4. Neutral atmosphere: humidity, refractivity, gas, and bending

4.1 Relative humidity to vapor pressure

Let $t_C$ be temperature in degrees Celsius and $RH$ in percent. Nephele implements the ITU-R P.453 phase-specific saturation-vapor-pressure form

Pressure and $e_s$ are in hPa.

4.2 Neutral radio refractivity

[Implemented] ITU-R P.453-14 gives

The neutral excess path is

Nephele requests up to 44 public NOAA weather-model pressure surfaces from 1000 to 10 hPa and derives $N$ only where the returned thermodynamic inputs are valid. The usable count and ceiling depend on model and location. These are model samples, not direct soundings or independent observed degrees of freedom.

4.3 Gaseous absorption and dispersion

[Implemented in the offline research estimator] The ITU-R P.676-13 line-by-line oxygen and water-vapor attenuation model gives

satrf/src/satrf/spectroscopy.py implements the complete P.676-13 oxygen/water-vapor line tables, line strengths, widths, continua, and attenuation equations. Its numerical output was independently checked against ITU-Rpy 0.4.0. The nonlinear joint estimator combines this attenuation with P.453 neutral delay and P.840-9 liquid absorption.

[Implemented in the additive high-fidelity Python path] The piecewise WGS-84 ray and frozen-ray array inverse use the real dispersive part of P.676 under one phase/group convention, with P.453/P.840, curved-ray bending, attenuation, antenna/hardware outputs, and explicit rank/DOFS abstention. It accepts either calibrated medium-only quantities or raw resolved carrier phase, group ToF, and RSSI after full-covariance whitening and exact cross-link clock/gain nullspace projection. The latter is differential-only: wrapped phase, unresolved integer ambiguity, static-state Doppler assimilation, absolute nuisance recovery, and atmospheric emission/brightness temperature remain scientific targets. See HIGH_FIDELITY_REFRACTION.md.

[Implemented in the interactive path] lib/atmospheric-propagation.ts evaluates the complete P.676-13 oxygen and water-vapour line tables and continua for every supplied pressure-level state. It converts total pressure, temperature, and relative humidity with P.453-14, then integrates gas loss and neutral excess path over caller-supplied spherical-shell segment lengths. Repository float64 golden vectors cover scalar kernels and complete column integrals. The public NWP column requests 44 pressure surfaces but retains only the usable fields returned by the explicit HRRR or GFS model, so the UI labels its actual ceiling and never treats the requested layer count as independent observational resolution.

5. Liquid cloud, multiple cloud layers, rain, ice, and sky noise

5.1 ITU liquid-cloud attenuation

For small liquid droplets under the Rayleigh approximation, ITU-R P.840-9 gives

Here $M$ is LWC in g/m$^3$ and $K_l$ has units $(\mathrm{dB/km})/(\mathrm{g/m^3})$. The small-droplet model is specified for liquid cloud/fog over 1--200 GHz; it is not an ice-cloud model.

5.2 Double-Debye water coefficient used by Nephele

[Implemented] With $f$ in GHz, $233\le T\le303$ K, and

For layer $k$,

5.3 Cloud carrier-phase model

P.840 standardizes attenuation; it does not alone provide a complete coherent cloud phase observable. A consistent dilute-droplet extension uses complex Maxwell-Garnett mixing.

Let the liquid volume fraction be

and

Then

The complex cloud-only field ratio is

Thus layer $k$ contributes

[Approximation, implemented] The current display uses the real, first-order dilute form

It wraps displayed phase into $[-180^\circ,180^\circ)$, although inversion must retain unwrapped phase and cycle-slip flags.

5.4 Current public-weather cloud-layer construction

[Approximation, implemented] The public pressure-level feed supplies cloud cover, temperature, height, and wind but not LWC on those levels. Nephele therefore creates no more than four scenario slabs from levels whose cover $C_k\ge12\%$.

Layer boundaries are midpoints between adjacent pressure-level heights, with thickness clipped to 0.25--3 km. Scenario LWC is

This creates separately inspectable layers but does not make their water contents observed. A field system should ingest cloud liquid/ice mixing ratios, radar/radiometer constraints, or retrieve them from calibrated RF observables.

5.5 Rain attenuation

ITU-R P.838-3 gives

For horizontal/vertical coefficients $k_H,k_V,\alpha_H,\alpha_V$, path elevation $E$, and polarization tilt $\tau$,

[Implemented] Nephele evaluates the P.838-3 coefficient fits. For an instantaneous grid rain rate, it forms a declared scenario column from the instantaneous pressure-level 0 °C isotherm plus 0.36 km and integrates $A_R=\gamma_RL_R$ through that curved-Earth column. That instantaneous height is not labeled P.839 climatology. P.839 rain-height and P.618 effective-path attribution are available only when the caller explicitly supplies both a climatological rain height and the $R_{0.01}$ rate. Missing rain rate, unbracketed freezing level, unresolved station MSL height, or below-horizon geometry returns unavailable rather than dry weather.

P.838 does not define a complete coherent rain phase, particle-size distribution, depolarization, or Doppler spectrum. [Scientific target] Use a complex forward-scattering model driven by drop-size distribution, temperature, shape, canting, and polarization; its extinction must be checked against P.838.

5.6 Ice, snow, hail, and mixed phase

[Scientific target] Ice response depends on ice-water content, particle size distribution, habit, density, orientation, aggregation, melting, and frequency. Use documented scattering databases or T-matrix/discrete-dipole models. Do not apply the liquid P.840 coefficient to ice. Until constrained, ice/mixed-phase layers must carry an explicit model-error term.

Polarimetric differential phase is

under the common weather-radar convention for one-way specific differential phase $K_{DP}$. The exact factor and sign must match the adopted PRO/radar convention. Differential phase can suppress common geometry/clock terms and is substantially more diagnostic of oriented nonspherical hydrometeors than a single scalar carrier phase.

5.7 Atmospheric emission and sky noise

Attenuating atmosphere is also emissive. A useful equivalent-sky approximation is

The receiver should propagate this through system noise temperature/G/T rather than showing attenuation alone.

6. Phase and Doppler decomposition

6.1 Total received phase

The tracked carrier phase can be written

$N_{amb}$ is the integer-cycle ambiguity. Absolute cloud phase is not recoverable until clocks, geometry, hardware, ambiguity, and non-cloud media are estimated or differenced.

6.2 One-way geometric Doppler

The instantaneous residual frequency is

The first-order geometric term is

[Implemented] Nephele evaluates this term for every selected satellite/station pair.

6.3 Neutral/cloud excess Doppler

For any phase excess path,

For a moving ray parameterized by $\mathbf x(\xi,t)$,

This separates medium evolution/advection from the rapidly moving LEO ray footprint.

[Approximation, implemented] The displayed cloud-motion scenario uses

This is an exploratory frozen-flow phase-screen sensitivity, not a retrieved wind product.

6.4 Moving-scatterer Doppler

Let $\hat{\mathbf k}_{in}$ point from transmitter to scatterer and $\hat{\mathbf k}_{out}$ from scatterer to receiver. A moving hydrometeor/scatterer contributes approximately

Its magnitude approaches $2|v_r|/\lambda$ for monostatic backscatter and approaches zero for direct forward propagation. Therefore direct transparent-cloud Doppler is not equivalent to weather-radar Doppler.

[Implemented simulation] The array mode now places a synthetic rain-scattering cell inside the dominant modeled cloud layer and generates a delay-gated bistatic Doppler observation for every usable satellite--station pair. A weighted least-squares inverse estimates the local east, north, and upward hydrometeor-velocity components from the scattering vectors only. The truth generator combines scenario horizontal wind, vertical air motion, and terminal fall speed; those truth values are withheld from the inverse and used solely to report vector error. Rank-deficient geometries and fewer than four informative observations are rejected. This is a proof-of-concept observing-system simulation, not a claim that an unmodified communications receiver exposes range-resolved precipitation echoes.

6.5 Ionospheric phase and Doppler

For slant TEC,

the first-order group and carrier-phase paths are

Thus

Dual/multi-frequency observations separate the leading $f^{-2}$ ionosphere from approximately nondispersive neutral delay away from absorption lines.

6.6 Error propagation from orbit knowledge

Range and range-rate errors directly contaminate weather phase and Doppler:

At 12 GHz, a 1 m/s range-rate error is about 40 Hz. A 1 Hz geometric-Doppler residual requires roughly 2.5 cm/s range-rate knowledge before clock, oscillator, and hardware errors are considered.

7. RF-to-digital receiver emulation

7.1 RF, downconversion, and complex samples

For transmitted complex waveform $x(t)$,

Mixing with local oscillator $f_{LO}$ and low-pass filtering gives

where $\Delta f=f_c+f_D-f_{LO}$ and $q[n]$ is quantization error.

The receiver should first wipe off predicted geometric Doppler, then estimate residual carrier/phase on known pilots:

7.2 Noise, SNR, $C/N_0$, and quantization

At approximately 290 K,

For an ideal $b$-bit ADC driven near full scale,

[Implemented] Nephele combines thermal and quantization SNR as

For bit rate $R_b$,

7.3 QPSK decode model

[Approximation, implemented] Nephele assumes Gray-coded coherent QPSK with $R_b=2R_s$ and

The demonstration decoder declares lock if line of sight exists,

Its output bytes are synthetic. A real system needs the actual waveform, pilots, coding, AGC state, sample clock, phase-noise model, beam/handover metadata, and raw calibrated I/Q.

8. What one moving satellite-to-ground link can retrieve

At observation time $t_i$, a direct link generally measures a path integral:

For horizontally uniform layers with thicknesses $\Delta h_j$,

so every matrix row is a scalar multiple of the same thickness vector:

Therefore $\operatorname{rank}(\mathbf A)=1$ in the ideal plane-parallel case. A single direct channel can robustly constrain a total column or total attenuation, but it cannot uniquely separate multiple vertical layers merely by taking many closely related samples.

The moving LEO nevertheless provides useful samples

which can show changing integrated weather along the pass. A 2-D along-track curtain becomes possible if the cloud field is horizontally structured and a frozen-flow/advection prior is credible. Unique 3-D vertical attribution still requires crossing rays, frequency/polarization diversity, external priors, or occultation geometry.

The sample interval should satisfy both receiver dynamics and atmospheric evolution. If the highest residual fluctuation of interest is $f_{max}$,

while a tomography window $T_w$ should keep evolution error below the measurement error:

More display samples do not create more independent atmospheric information if the rows of the observation Jacobian remain collinear.

9. Multiple satellites, stations, layers, and 3-D/4-D tomography

9.1 Linear ray-to-voxel model

For ray $m$ and voxel $v$, let $A_{mv}=\ell_{mv}$ be the in-voxel path length. Then

For slant water vapor (SWV), $x_v$ is water-vapor density. For cloud attenuation at one frequency,

For several physical mechanisms,

This joint multi-frequency/polarization problem is nonlinear when rain microphysics, ice, clocks, and wind are unknown.

9.2 Multilayer contribution accounting

For known layer geometry, every ray's modeled contributions can be kept separate:

This answers “how much did each modeled layer contribute?” It does not prove that the measured total uniquely determines each term. Measured layer separation requires an adequately ranked Jacobian

with independent elevation, azimuth, time, frequency, polarization, and/or external constraints.

9.3 Current synthetic cloud truth

[Implemented, synthetic] Nephele constructs an advected multilayer LWC field for testing. For voxel $v$ and layer $k$,

Each Gaussian center is translated using the layer wind. The synthetic observation is

with a deterministic grid-relative perturbation. The browser OSSE offsets paths by as much as 20% of the smallest horizontal cell and 15% of a vertical cell, plus a declared 2% gain bias. The inversion never receives $\mathbf x_{truth}$ directly.

9.4 ART

For ray row $\mathbf a_i^T$, non-negative algebraic reconstruction is

[Implemented] Nephele uses $\lambda=0.72$, then applies either neighbor smoothing or a smoothed total-variation step.

9.5 SIRT

Define

[Implemented] The update is

9.6 MART

[Implemented] For row prediction $\hat y_i=\mathbf a_i^T\mathbf x$ and $s_i=\sum_jA_{ij}$,

9.7 MAP/regularized scientific inversion

[Implemented offline research estimator; reduced form shown] Estimate the atmospheric and nuisance state by

$\boldsymbol\beta$ contains satellite/station clocks, EIRP and gain offsets, beam state, wet-antenna loss, multipath, and other nuisance parameters. The implemented OSSE instantiates clock/path, receiver-gain, and beam/handover states. Wet-antenna and sparse multipath terms remain model-discrepancy stresses, not estimated state blocks.

9.8 4-D advection/evolution

An atmospheric scalar $x$ may evolve as

The simplest frozen-flow model is

[Implemented in both paths, with different roles] The browser uses horizontal semi-Lagrangian back-tracing and bilinear interpolation as a prior to each independent display solve. The offline research estimator places the same operator inside one weak-constraint six-epoch optimization. Formation, dissipation, vertical motion, turbulence, phase changes, and precipitation are model error, not captured by pure translation.

9.9 Implemented nonlinear joint thermodynamic, multi-epoch inverse

The original joint_4d.py wet-refractivity/ liquid-water inverse remains a precursor and regression benchmark. The active research estimator is thermodynamic_4d.py. At epoch $t$ it estimates $4\times4\times4$ trilinear control fields for temperature $T$, relative humidity $H$, and cloud liquid $M$, mapped to a $5\times5\times4$ reporting grid. A deterministic ISA-like dry-pressure profile $p_d(z)$ is a declared structural model, not a measured or historical weather input.

P.453 closes RH, vapor density, and refractivity through

The first equation is solved by a fixed-point iteration because the P.453 enhancement factor uses total pressure while P.676 separates dry and vapor pressure. Each frequency has

using exact P.676-13 and P.840-9 coefficients. With ray-length vector $\mathbf a_{it}$ in km and line-of-sight unit vector $\hat{\mathbf l}_{it}$,

The signed nuisance state therefore contains receiver/satellite path-equivalent clocks, three-axis satellite orbit error, receiver gain, and beam/handover gain. Bounds enforce $0\le H\le100$, $M\ge0$, and temperature within 15 K of the structural profile. A physically scaled sparse Gauss--Newton solve minimizes the nonlinear RF likelihood plus absolute-profile, spatial, semi-Lagrangian dynamic, and nuisance random-walk penalties. It receives no prior weather measurement.

Six epochs, 768 paths, and one delay plus six attenuation channels at 12/22.235/30/40/50/54 GHz produce 5,376 observations. On untouched seeds 904--905, mean RMSE is 0.620 deg C, 2.829 RH percentage points, and 0.406 mm PWV; liquid RMSE is 0.0062 g/m3. A no-RF standard profile gives 1.307 deg C, 5.814 points, and 3.499 mm. The likelihood Jacobian is full column rank on seed 904 (1,152/1,152) but has condition number $2.04\times10^7$: bounds and regularization are essential, and output voxels are not independently resolved.

A staged post-freeze audit adds unchanged base-model seeds 917--924 and 925--928. Across all 14 base trials (including 904--905), temperature/RH/PWV RMSE is $0.679\pm0.103$ deg C, $3.203\pm0.432$ points, and $0.433\pm0.056$ mm; 13/14 trials pass all three gates. Seed 921 is the single failure at 4.039 RH points. A subsequent physics-only E-optimal screen selected 175 GHz under a 40 dB/10-km loss constraint, but paired untouched seeds 925--928 reject it: PWV improves from 0.423 to 0.212 mm while temperature/RH degrade from 0.636/3.003 to 0.883/4.642, all-gate passes fall from 4/4 to 1/4, and mean solve time rises from 44 to 130 s. The active model therefore retains the six-band design. Local Jacobian rank and singular values are necessary diagnostics, not sufficient model-selection criteria for the global nonlinear weak-constraint inverse.

The seed-906 ablation gives 2.613 deg C/15.385 points/2.438 mm without nuisance estimation, 0.671/3.533/0.368 for joint static, 0.629/3.300/0.388 for joint 4D, and 0.503/3.250/0.388 for an oracle-calibrated dynamic solve. Temporal coupling improves temperature, RH, and liquid over joint static, but not PWV on this one ablation seed.

[Implemented robust cycle] For a same-link group $G_l$ containing delay and all six attenuation bands, the standardized residual magnitude and group-Huber weight are

All rows in $G_l$ receive $\sqrt{w_l}$. The frozen cycle uses $\delta=1.5$ and $w_{min}=0.05$; robust serving requires weighted reduced $\chi^2\le3$, at most 10% flagged links, and at least 75% effective weight. Under 5% coherent bursts, the Gaussian solve fails at 1.364 deg C/6.602 points/0.700 mm, while robust weighting recovers 0.591/3.224/0.373, with 94.7% flag precision and 100% recall. Thirty-percent link loss, doubled nuisance, wrong wind, and the combined robust stress also remain below the three requested gates. See the frozen thermodynamic_4d_robustness.json.

9.10 Posterior uncertainty, resolution, and observability

For a locally linear model with prior covariance $\mathbf B$ and observation covariance $\mathbf R$,

The averaging/resolution matrix is

and degrees of freedom for signal are

Use singular values, condition number, vertical/horizontal averaging kernels, posterior correlation length, and $DFS$ to report actual resolving power. Voxel size alone is not resolution.

A network-design objective can maximize incremental information,