Atmosphere

Standard-atmosphere, thermosphere and ionosphere models for propagation work, plus humidity conversions and the refraction suite: astronomical refraction (Sinclair atmosphere, add/remove), the standard-exponential-model radar refraction conversions (bistatic r-u-v ray tracing, bias approximation, cubature variants), speed of sound, and (since v2.10.0) NRLMSISE-00 itself, shipped as the vendored public-domain reference C implementation compiled into pytcl.atmosphere._nrlmsise00_c with a validated pure-Python fallback. v2.11.0 adds the Jacchia 1971 model (pytcl.atmosphere.jacchia), completing the port of MATLAB’s Atmosphere_and_Refraction models tracked in MATLAB TCL parity inventory.

Atmospheric models for tracking applications.

This module provides standard atmosphere models used for computing temperature, pressure, density, and other properties at various altitudes.

Submodules

models : Standard atmosphere models (US76, ISA) ionosphere : Ionospheric models for GPS/GNSS corrections humidity : Humidity conversions and dew-point calculations refraction : Atmospheric refractivity helpers

Atmospheric Models

US Standard Atmosphere 1976 and ISA density/temperature/pressure models, plus the pressure-altitude, Mach and true-airspeed conversions built on them.

Atmospheric models for tracking applications.

This module provides standard atmosphere models used for computing temperature, pressure, and density at various altitudes.

class pytcl.atmosphere.models.AtmosphereState(temperature, pressure, density, speed_of_sound)[source]

Bases: NamedTuple

Atmospheric state at a given altitude.

Variables:

temperature (float or ndarray) – Temperature in Kelvin.
pressure (float or ndarray) – Pressure in Pascals.
density (float or ndarray) – Density in kg/m³.
speed_of_sound (float or ndarray) – Speed of sound in m/s.

temperature: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 0

pressure: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 1

density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 2

speed_of_sound: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 3

pytcl.atmosphere.models.us_standard_atmosphere_1976(altitude)[source]

Compute atmospheric properties using US Standard Atmosphere 1976.

Parameters:: altitude (array_like) – Geometric altitude in meters. Defined from -5000 m to ~86 km (84852 m geopotential); outside that range the nearest layer’s lapse rate is extrapolated and a warning is raised.
Returns:: state – Atmospheric state containing temperature, pressure, density, and speed of sound.
Return type:: AtmosphereState
Warns:: UserWarning – If altitude is below -5000 m or above the ~86 km ceiling, since the standard’s layer table does not extend past there and the result is extrapolated from the nearest boundary.

Examples

>>> state = us_standard_atmosphere_1976(10000)
>>> round(state.temperature, 3)
223.252
>>> round(state.pressure, 1)
26499.9

Notes

The US Standard Atmosphere 1976 is a model of the Earth’s atmosphere that defines temperature, pressure, and density as functions of altitude. It is defined from -5000 m to approximately 86 km altitude.

References

U.S. Standard Atmosphere, 1976, U.S. Government Printing Office, Washington, D.C., 1976.

pytcl.atmosphere.models.isa_atmosphere(altitude, temperature_offset=0.0)[source]

Compute atmospheric properties using International Standard Atmosphere (ISA).

This is essentially the troposphere portion of US Standard Atmosphere 1976 with an optional temperature offset for non-standard days.

Parameters:

altitude (array_like) – Geometric altitude in meters.
temperature_offset (float, optional) – Temperature offset from ISA conditions in Kelvin (default: 0). Positive values indicate warmer than standard day.

Returns:

state – Atmospheric state.

Return type:

AtmosphereState

Examples

>>> # Standard day at 5000m
>>> state = isa_atmosphere(5000)
>>> # Hot day (+15K) at 5000m
>>> state = isa_atmosphere(5000, temperature_offset=15)

pytcl.atmosphere.models.altitude_from_pressure(pressure)[source]

Compute geometric altitude from pressure (pressure altitude).

Parameters:: pressure (array_like) – Atmospheric pressure in Pascals.
Returns:: altitude – Geometric altitude in meters.
Return type:: ndarray
Warns:: UserWarning – If pressure falls outside the troposphere layer this inversion covers (the pressure at the -5000 m floor down to the pressure at the 11000 m tropopause); the result is clamped to that layer’s boundary rather than extrapolated through the model’s other (differently-sloped) layers.

Examples

>>> # Sea level pressure
>>> bool(abs(altitude_from_pressure(101325)) < 1e-6)
True
>>> # Pressure at approximately 5000m
>>> alt = altitude_from_pressure(54000)
>>> 4800 < alt < 5200
True

Notes

This inverts the closed-form troposphere gradient formula, so it is only exact within the troposphere layer (-5000 m to 11000 m geopotential). Outside that layer’s pressure range the model switches to a different lapse rate that this inversion does not represent.

pytcl.atmosphere.models.mach_number(velocity, altitude)[source]

Compute Mach number from velocity and altitude.

Parameters:

velocity (array_like) – True airspeed in m/s.
altitude (array_like) – Geometric altitude in meters.

Returns:

mach – Mach number.

Return type:

ndarray

Examples

>>> # Aircraft at 300 m/s at sea level
>>> mach_number(300, 0)
0.88...
>>> # Same speed at 10 km altitude (lower speed of sound)
>>> mach_number(300, 10000)
1.00...

pytcl.atmosphere.models.true_airspeed_from_mach(mach, altitude)[source]

Compute true airspeed from Mach number and altitude.

Parameters:

mach (array_like) – Mach number.
altitude (array_like) – Geometric altitude in meters.

Returns:

velocity – True airspeed in m/s.

Return type:

ndarray

Examples

>>> # Mach 0.8 at cruise altitude (10 km)
>>> tas = true_airspeed_from_mach(0.8, 10000)
>>> 230 < tas < 250  # approximately 240 m/s
True
>>> # Supersonic at sea level
>>> true_airspeed_from_mach(1.0, 0)
340.2...

Thermosphere Model

Simplified barometric thermosphere density, temperature and composition model with solar-activity and geomagnetic inputs. Not NRLMSISE-00: usable above ~200 km (within ~2x of published NRLMSISE-00 values), up to 50x wrong below ~86 km where us_standard_atmosphere_1976 should be used. Limits are documented in the module and pinned by validation tests (gh-79). For the real model at any altitude, use pytcl.atmosphere.nrlmsise00.

Barometric thermosphere model.

This is not NRLMSISE-00. The module was named nrlmsise00 and described itself as a high-fidelity NRL model until gh-79. What is implemented here is a set of per-species exponential profiles with temperature-dependent scale heights, driven by F10.7 and Ap, each clamped to a floor.

Deprecated since version 2.11.0: The real NRLMSISE-00 – coefficient tables included, compiled from the reference C – ships as pytcl.atmosphere.nrlmsise00 since v2.10.0. Use it instead; this approximation will be removed in v3.0.0.

It is usable above roughly 200 km, where it agrees with published NRLMSISE-00 densities to within a factor of about two. Below that it is wrong by up to 50x, and pytcl.atmosphere.us_standard_atmosphere_1976() should be used instead. SimplifiedThermosphere carries the measured comparison.

The references below describe the model this one approximates, not the one implemented here. They are retained because the temperature parameterisation and species set follow their structure.

References

Picone, J. M., A. E. Hedin, D. P. Drob, and A. C. Aikin (2002), “NRLMSISE-00 empirical model of the atmosphere: Statistical comparisons and scientific issues,” J. Geophys. Res., 107(A12), 1468, doi:10.1029/2002JA009430
NASA GSFC NRLMSISE-00 Model: https://ccmc.gsfc.nasa.gov/models/nrlmsise00
Drob, D. P., et al. (2008), “An update to the COSPAR International Reference Atmosphere model for the middle atmosphere,” Adv. Space Res., 43(12), 1747-1764

class pytcl.atmosphere.thermosphere.SimplifiedThermosphere(use_meter_altitude=True)[source]

Bases: object

Barometric thermosphere model with solar and geomagnetic coupling.

This is not NRLMSISE-00. It was named NRLMSISE00 and described as “a comprehensive thermosphere model” until gh-79. What this computes is a set of per-species exponential profiles with temperature-dependent scale heights, driven by F10.7 and Ap, with each species density clamped to a floor. The real NRLMSISE-00 ships as pytcl.atmosphere.nrlmsise00 since v2.10.0 – prefer it; this class is deprecated and will be removed in v3.0.0.

Use it above roughly 200 km. Below that it is wrong, in places by more than an order of magnitude:

Altitude	This model	US Standard 1976	Ratio
0 km	0.682 kg/m^3	1.225 kg/m^3	0.56
10 km	0.253	0.414	0.61
30 km	0.0261	0.0184	1.42
50 km	1.044e-3	1.027e-3	1.02
80 km	3.69e-7	1.85e-5	0.02

For altitudes below about 86 km use pytcl.atmosphere.us_standard_atmosphere_1976(), which is validated against the standard. In the thermosphere this model is reasonable: at 400 km with F10.7 = 150 it gives 2.9e-12 kg/m^3, against a published NRLMSISE-00 range of roughly 2-4e-12.

The species floors are visible in the output. h_density returns its floor of 1e8 m^-3 from 0 to 400 km, and he_density its floor of 1e10 m^-3 at 0 and 100 km, because the computed profile falls below the clamp. Those are not measurements of anything.

The model implements: - Temperature profile with solar activity and magnetic coupling - Per-species barometric profiles with floors - Solar flux (F10.7) and magnetic activity (Ap) variations

Parameters:: use_meter_altitude (bool, optional) – If True, expect altitude input in meters. If False, expect km. Default is True (meters).

Examples

>>> import warnings
>>> with warnings.catch_warnings():
...     warnings.simplefilter("ignore", DeprecationWarning)
...     model = SimplifiedThermosphere()
>>> output = model(
...     latitude=np.radians(45),
...     longitude=np.radians(-75),
...     altitude=400_000,  # 400 km
...     year=2024,
...     day_of_year=100,
...     seconds_in_day=43200,
...     f107=150,
...     f107a=150,
...     ap=5
... )
>>> print(f"Density: {output.density:.2e} kg/m³")
Density: 2.85e-12 kg/m³

Notes

This implementation uses empirical correlations for atmospheric properties as a function of geomagnetic and solar activity indices. For highest accuracy, use the original NRLMSISE-00 Fortran code from NASA/NOAA, which includes extensive coefficient tables.

__init__(use_meter_altitude=True)[source]

Initialize the simplified thermosphere model.

__call__(latitude, longitude, altitude, year, day_of_year, seconds_in_day, f107=150.0, f107a=150.0, ap=4.0)[source]

Compute atmospheric density and composition.

Parameters:

latitude (array_like) – Geodetic latitude in radians.
longitude (array_like) – Longitude in radians.
altitude (array_like) – Altitude in meters (or km if use_meter_altitude=False).
year (int) – Year (e.g., 2024).
day_of_year (int) – Day of year (1-366).
seconds_in_day (float) – Seconds since midnight (0-86400).
f107 (float, optional) – 10.7 cm solar flux (daily value, SFU). Default 150.
f107a (float, optional) – 10.7 cm solar flux (81-day average, SFU). Default 150.
ap (float or array_like, optional) – Planetary magnetic index. Can be single value or 8-element array of 3-hour Ap values. Default 4.0.

Returns:

output – Atmospheric properties (density, temperature, composition).

Return type:

ThermosphereState

Notes

The model assumes hydrostatic equilibrium and uses empirical correlations for density and temperature variations.

class pytcl.atmosphere.thermosphere.ThermosphereState(density, temperature, exosphere_temperature, he_density, o_density, n2_density, o2_density, ar_density, h_density, n_density)[source]

Bases: NamedTuple

Output from the simplified thermosphere model (NRLMSISE-00-style interface; NOT NRLMSISE-00 – see module docstring).

Variables:

density (float or ndarray) – Total atmospheric density in kg/m³.
temperature (float or ndarray) – Temperature at altitude (K).
exosphere_temperature (float or ndarray) – Exospheric temperature (K).
he_density (float or ndarray) – Helium density in m⁻³.
o_density (float or ndarray) – Atomic oxygen density in m⁻³.
n2_density (float or ndarray) – N₂ density in m⁻³.
o2_density (float or ndarray) – O₂ density in m⁻³.
ar_density (float or ndarray) – Argon density in m⁻³.
h_density (float or ndarray) – Hydrogen density in m⁻³.
n_density (float or ndarray) – Atomic nitrogen density in m⁻³.

density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 0

temperature: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 1

exosphere_temperature: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 2

he_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 3

o_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 4

n2_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 5

o2_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 6

ar_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 7

h_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 8

n_density: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 9

class pytcl.atmosphere.thermosphere.F107Index(f107, f107a, ap, ap_array=None)[source]

Bases: NamedTuple

Solar activity indices for the simplified thermosphere model (NRLMSISE-00-style interface; NOT NRLMSISE-00).

Variables:

f107 (float) – 10.7 cm solar radio flux (daily, SFU).
f107a (float) – 10.7 cm solar radio flux (81-day average, SFU).
ap (float or ndarray) – Planetary magnetic index (Ap index).
ap_array (ndarray, optional) – Ap values for each 3-hour interval of the day (8 values). Optional 7-element ap history. Nothing derives it from ap – it is a plain default of None, retained for NRLMSISE-00-style input compatibility and unused by SimplifiedThermosphere.

f107: float: Alias for field number 0

f107a: float: Alias for field number 1

ap: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 2

ap_array: ndarray[tuple[Any, ...], dtype[float64]] | None: Alias for field number 3

pytcl.atmosphere.thermosphere.simplified_thermosphere(latitude, longitude, altitude, year, day_of_year, seconds_in_day, f107=150.0, f107a=150.0, ap=4.0)[source]

Compute simplified thermosphere atmospheric properties (NRLMSISE-00-style interface; NOT NRLMSISE-00 – see module docstring).

This is a module-level convenience function wrapping the SimplifiedThermosphere class.

Parameters:

latitude (array_like) – Geodetic latitude in radians.
longitude (array_like) – Longitude in radians.
altitude (array_like) – Altitude in meters.
year (int) – Year (e.g., 2024).
day_of_year (int) – Day of year (1-366).
seconds_in_day (float) – Seconds since midnight (0-86400).
f107 (float, optional) – 10.7 cm solar flux (daily value, SFU). Default 150.
f107a (float, optional) – 10.7 cm solar flux (81-day average, SFU). Default 150.
ap (float or array_like, optional) – Planetary magnetic index. Default 4.0.

Returns:

output – Atmospheric properties.

Return type:

ThermosphereState

Notes

See SimplifiedThermosphere class for more details.

Examples

>>> # ISS altitude (~400 km), magnetic latitude = 40°, quiet geomagnetic activity
>>> import warnings
>>> with warnings.catch_warnings():
...     warnings.simplefilter("ignore", DeprecationWarning)
...     output = simplified_thermosphere(
...     latitude=np.radians(40),
...     longitude=np.radians(-75),
...     altitude=400_000,  # 400 km
...     year=2024,
...     day_of_year=1,
...     seconds_in_day=43200,
...     f107=150,  # Average solar activity
...     f107a=150,
...     ap=5  # Quiet conditions
... )
>>> print(f"Density at ISS: {output.density:.2e} kg/m³")
Density at ISS: 2.89e-12 kg/m³

Jacchia 1971

The Jacchia 1971 reference atmosphere (90 km to 2500 km): exospheric temperature from solar and geomagnetic activity, the empirical temperature profile, and bi-polynomial densities with the semiannual, geomagnetic, and seasonal-latitude corrections.

Jacchia 1971 reference atmosphere.

Port of jacchiaAtmosParam: exospheric temperature from solar flux, geomagnetic activity, and the diurnal bulge; the empirical temperature profile; and density from Montenbruck & Gill’s bi-polynomial fit with the semiannual, geomagnetic, and seasonal-latitude corrections.

Known deviations from the MATLAB original, all deliberate and loud:

The MATLAB function is unrunnable as shipped: it references three Earth-orientation variables (deltaTTUT1, xpyp, dXdY) that are never assigned. The oracle fixtures were captured with a patched copy that takes the solar declination and local hour angle as inputs (see scripts/matlab_capture/capture_jacchia.m); this port computes them with pytcl’s own time-scale chain and a low-precision solar position (accurate to ~0.01 degrees, far below the sensitivity of the diurnal term).
MATLAB computes pressure as rho * R * T / M with the molar mass in g/mol, which yields kPa while documenting Pa. This port returns true pascals (a factor of 1000 larger than MATLAB’s number).
The seasonal-latitude correction’s sin(lat)^3 / |sin(lat)| is 0/0 at the equator; MATLAB silently returns NaN there. This port uses the limit value 0 (the factor is sign(sin lat) * sin^2 lat).
The mean-molecular-mass polynomial behind the pressure output is a fit for the 90-100 km mixing region; at higher altitudes it diverges and MATLAB silently returns a negative pressure. This port emits a RuntimeWarning and returns NaN pressure where the polynomial’s molar mass is not positive.

References

L. G. Jacchia, “Revised Static Models of the Thermosphere and Exosphere with Empirical Temperature Profiles,” SAO Special Report 332, Cambridge, MA, 1971.
O. Montenbruck and E. Gill, Satellite Orbits: Models, Methods and Applications, Springer, 2000 (Tables 3.9 and 3.10).

class pytcl.atmosphere.jacchia.JacchiaState(density, pressure, temperature, exospheric_temperature)[source]

Bases: NamedTuple

Jacchia 1971 atmospheric state.

Variables:

density (float) – Total mass density, kg/m^3.
pressure (float) – Pressure in Pa (ideal-gas law; the MATLAB original returns kPa while documenting Pa – see the module docstring).
temperature (float) – Temperature at the requested altitude, K.
exospheric_temperature (float) – Exospheric temperature, K.

density: float: Alias for field number 0

pressure: float: Alias for field number 1

temperature: float: Alias for field number 2

exospheric_temperature: float: Alias for field number 3

pytcl.atmosphere.jacchia.jacchia_atmos_param(jul1, jul2, point, f107, f107b, kp)[source]

Jacchia 1971 atmospheric density, pressure, and temperature.

Port of jacchiaAtmosParam. Valid for altitudes of 90 km to 2500 km and exospheric temperatures of 500 K to 1900 K (both enforced, as in the original).

Parameters:

jul1 (float) – Two-part Julian date in UTC (days); the full date is the sum.
jul2 (float) – Two-part Julian date in UTC (days); the full date is the sum.
point (array_like) – Geodetic [latitude, longitude, altitude] with latitude and longitude in radians and altitude in meters.
f107 (float) – 10.7 cm solar radio flux averaged over the previous day, in units of 1e-22 W/(m^2 Hz).
f107b (float) – Same flux averaged over three 27-day solar rotations.
kp (float) – Three-hourly planetary geomagnetic index for a time 6.7 hours before the epoch.

Returns:

state – Density (kg/m^3), pressure (Pa), temperature (K), and exospheric temperature (K).

Return type:

JacchiaState

Raises:

ValueError – If the altitude or the resulting exospheric temperature is outside the model’s tabulated domain.

Examples

>>> import numpy as np
>>> state = jacchia_atmos_param(
...     2451545.0, 0.25, [np.radians(40.0), np.radians(-75.0), 400e3],
...     150.0, 150.0, 3.0,
... )
>>> 1e-13 < state.density < 1e-11  # ~ISS-altitude density
True
>>> 600.0 < state.exospheric_temperature < 1400.0
True

NRLMSISE-00

The NRLMSISE-00 empirical atmosphere model: the vendored reference C implementation (compiled at build time) with a machine-precision pure-Python fallback.

The NRLMSISE-00 empirical atmosphere model.

A faithful transcription of the public-domain C reference implementation as vendored (and lightly modified) by the MATLAB TCL (3rd_Party_Libraries/nrlmsise-00-bc9a2fe), which is the validation oracle: every code path here is checked against outputs of that C code compiled directly. The NRL modifications preserved here are the ghp7 altitude return and the MEX-era diagnostic hooks (which become warnings).

The model [1] computes number densities of He, O, N2, O2, Ar, H, N and anomalous O, the total mass density, and the exospheric and local temperatures, for altitudes from the ground into the exosphere.

References

class pytcl.atmosphere.nrlmsise00.NRLMSISEOutput(d, t)[source]

Bases: NamedTuple

Result of the NRLMSISE-00 model.

Variables:

d (ndarray) – The nine densities, indexed as the reference implementation orders them: d[0] He, d[1] O, d[2] N2, d[3] O2, d[4] Ar, d[5] total mass density, d[6] H, d[7] N, d[8] anomalous O. Number densities are per cubic meter and d[5] is kg/m^3 (the SI switch of the reference is always on here, matching the MATLAB wrapper).
t (ndarray) – t[0] the exospheric temperature and t[1] the temperature at altitude, in Kelvin.

d: ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 0

t: ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 1

pytcl.atmosphere.nrlmsise00.nrlmsise00(doy, sec, alt_km, g_lat_deg, g_long_deg, lst, f107a=150.0, f107=150.0, ap=4.0, ap_array=None, effective_density=False)[source]

Evaluate the NRLMSISE-00 model (the gtd7/gtd7d interface).

All 24 model switches are enabled (SI output units), matching the MATLAB TCL wrapper; a 7-element ap_array activates storm mode (switches[9] = -1) exactly as there.

Parameters:

doy (int) – Day of year, counting from 1.
sec (float) – Second of the day (UT).
alt_km (float) – Altitude in kilometers.
g_lat_deg (float) – Geodetic latitude and longitude in degrees.
g_long_deg (float) – Geodetic latitude and longitude in degrees.
lst (float) – Local apparent solar time in hours. The reference documents the relation lst = sec/3600 + g_long/15; it is an independent input here, as in the reference.
f107a (float, optional) – 81-day average 10.7 cm solar flux centered on the day. Default 150.
f107 (float, optional) – Daily 10.7 cm solar flux for the previous day. Default 150.
ap (float, optional) – Daily magnetic index. Default 4.
ap_array (sequence of 7 floats, optional) – Storm-mode Ap history (daily, 3-hourly now/-3h/-6h/-9h, and the 12-33h and 36-57h averages, per the vendored C reference; the MATLAB wrapper header says 36-59h). Activates storm mode.
effective_density (bool, optional) – If True, use gtd7d: include anomalous oxygen in the total mass density d[5] (relevant to drag above ~500 km). Default False (gtd7).

Returns:

result

Return type:

NRLMSISEOutput

Examples

>>> out = nrlmsise00(172, 29000.0, 400.0, 60.0, -70.0, 16.0)
>>> float(round(out.t[0], 2))  # exospheric temperature, K
1250.54
>>> float(round(out.d[5] * 1e12, 4))  # total mass density, kg/m^3
4.0747

Notes

Transcribed from the NRL-vendored reference C implementation and validated against it record-for-record.

pytcl.atmosphere.nrlmsise00.nrlmsise00_alt_for_pressure(doy, sec, press_pa, g_lat_deg, g_long_deg, lst, f107a=150.0, f107=150.0, ap=4.0, ap_array=None)[source]

Altitude at which the atmosphere has the given pressure (ghp7).

Parameters:

doy (int) – As in nrlmsise00().
sec (float) – As in nrlmsise00().
g_lat_deg (float) – As in nrlmsise00().
g_long_deg (float) – As in nrlmsise00().
lst (float) – As in nrlmsise00().
f107a (float) – As in nrlmsise00().
f107 (float) – As in nrlmsise00().
ap (float) – As in nrlmsise00().
ap_array (Sequence[float] | None) – As in nrlmsise00().
press_pa (float) – Pressure in Pascals. Must be positive.

Returns:

alt_km (float) – The altitude in kilometers at which the model pressure matches.
output (NRLMSISEOutput) – The model output at that altitude.

Raises:

ValueError – If the pressure is not positive.

Warns:

UserWarning – Two distinct cases, worded differently so a reader can tell which happened. (1) The Newton iteration behind ghp7 fails outright and returns a non-finite altitude – a genuine non-convergence, the same condition the pure-Python path already raises this warning for. (2) The iteration converges (satisfies its own internal residual test) but to a finite altitude outside the model’s validity range (_GHP7_MIN_VALID_ALT_KM to _GHP7_MAX_VALID_ALT_KM, currently -5.0 to 1000 km) – the compiled backend’s own signal for this is only printf’d to native stdout, invisible to Python, so it is otherwise silent.

Return type:

Tuple[float, NRLMSISEOutput]

Examples

>>> alt, out = nrlmsise00_alt_for_pressure(
...     172, 29000.0, 1000.0, 60.0, -70.0, 16.0)
>>> bool(15.0 < alt < 35.0)  # ~1000 Pa is upper-stratospheric
True

Notes

ghp7 iterates nrlmsise00() with scale-height Newton steps; the pressure is converted to the model’s hPa convention internally. Preserves the NRL modification returning the computed altitude.

pytcl.atmosphere.nrlmsise00.nrlmsise00_gas_temp(day_of_year, second_of_day, lat_lon_alt, ap=4.0, f107=150.0, f107a=150.0, ap_array=None, lst=None)[source]

Atmospheric constituents and temperature at a geodetic location.

The pytcl equivalent of the MATLAB TCL’s NRLMSISE00GasTemp MEX wrapper: geodetic coordinates in radians/meters in, a labeled gas table out.

Parameters:

day_of_year (int) – Day of the year (UTC), counting from 1.
second_of_day (float) – Second of the day, from 0. Values above 86400 are clipped to 86400, as in the MATLAB wrapper.
lat_lon_alt (sequence of 3 floats) – Geodetic latitude and longitude in radians and ellipsoidal height in meters (WGS-84).
ap (float, optional) – Daily geomagnetic index. Default 4.
f107 (float, optional) – Previous-day 10.7 cm solar radio flux. Default 150.
f107a (float, optional) – 81-day average 10.7 cm flux. Default 150.
ap_array (sequence of 7 floats, optional) – Storm-mode Ap history (see nrlmsise00()); activates storm mode.
lst (float, optional) – Local apparent solar time in hours. Default: computed as sec/3600 + lon_deg/15, the approximation the NRLMSISE-00 Users Guide documents and the MATLAB NRLMSISE00Alt4Pres wrapper uses.

Returns:

gas_table (list of (str, float)) – Constituent names (‘He’, ‘O’, ‘N2’, ‘O2’, ‘Ar’, ‘H’, ‘N’, ‘O*’) and their number densities in particles per cubic meter.
t (ndarray, shape (2,)) – t[0] is the exospheric temperature, t[1] the temperature at altitude, in Kelvin.
d (ndarray, shape (9,)) – The raw model densities; d[5] is the total mass density in kg/m^3 (excluding anomalous oxygen), the rest are the number densities of the table.

Return type:

Tuple[List[Tuple[str, float]], ndarray[tuple[Any, …], dtype[float64]], ndarray[tuple[Any, …], dtype[float64]]]

Examples

>>> import math
>>> table, t, d = nrlmsise00_gas_temp(
...     172, 29000.0, [60 * math.pi / 180, -70 * math.pi / 180, 400e3])
>>> table[2][0]
'N2'
>>> bool(t[1] > 500.0)  # thermospheric temperature at 400 km
True

Notes

Deliberate divergences from the MATLAB MEX wrapper, whose upstream defects are documented rather than reproduced: (1) it hardcodes lst = 16 (the formula above sits commented out in its source), ignoring time of day and longitude; here the formula is the default and lst is overridable. (2) Its optional-argument guards test nlhs (the output count, capped at 3) instead of nrhs, so its Ap/F107/F107A inputs are silently ignored and it always runs at defaults; even were they read, the F107A branch assigns into F107, and the storm-mode dimensionality check (mxGetM==7 && mxGetM==1) can never pass. Here all parameters are honored, including storm mode.

pytcl.atmosphere.nrlmsise00.nrlmsise00_pressure_altitude(day_of_year, second_of_day, pressure_pa, lat_lon, ap=4.0, f107=150.0, f107a=150.0, ap_array=None, lst=None)[source]

Altitude of a pressure level at a geodetic location.

The pytcl equivalent of the MATLAB TCL’s NRLMSISE00Alt4Pres MEX wrapper: latitude/longitude in radians and pressure in Pascals in, the altitude in meters plus the labeled gas table there out.

Parameters:

day_of_year (int) – Day of the year (UTC), counting from 1.
second_of_day (float) – Second of the day, from 0; clipped to 86400.
pressure_pa (float) – Pressure in Pascals. Must be positive.
lat_lon (sequence of 2 floats) – Geodetic latitude and longitude in radians (WGS-84).
ap (optional) – As in nrlmsise00_gas_temp().
f107 (optional) – As in nrlmsise00_gas_temp().
f107a (optional) – As in nrlmsise00_gas_temp().
ap_array (optional) – As in nrlmsise00_gas_temp().
lst (optional) – As in nrlmsise00_gas_temp().

Returns:

altitude_m (float) – The altitude in meters at which the model pressure matches.
gas_table (list of (str, float)) – Constituent number densities there (see nrlmsise00_gas_temp()).
t (ndarray, shape (2,)) – Exospheric temperature and temperature at altitude, Kelvin.
d (ndarray, shape (9,)) – The raw model densities at the altitude.

Raises:

ValueError – If pressure_pa is not positive (delegates to nrlmsise00_alt_for_pressure()).

Warns:

UserWarning – As in nrlmsise00_alt_for_pressure(), which this function delegates to – a non-finite or out-of-domain altitude solution warns rather than returning silently.

Return type:

Tuple[float, List[Tuple[str, float]], ndarray[tuple[Any, …], dtype[float64]], ndarray[tuple[Any, …], dtype[float64]]]

Examples

>>> import math
>>> alt_m, table, t, d = nrlmsise00_pressure_altitude(
...     172, 29000.0, 1000.0, [60 * math.pi / 180, -70 * math.pi / 180])
>>> bool(15e3 < alt_m < 35e3)  # ~1000 Pa is upper-stratospheric
True

Notes

The MATLAB MEX wrapper computes lst by the documented formula (its GasTemp sibling hardcodes 16), but shares the nlhs-for-nrhs guard defect: its Ap input takes effect only when exactly four outputs are requested, and its F107/F107A guards (nlhs>4, nlhs>5 against a four-output cap) can never pass, so neither is ever applied. All parameters are honored here.

pytcl.atmosphere.nrlmsise00.uses_compiled_backend()[source]

True when the compiled reference C extension is in use.

The extension is built from the vendored reference implementation (csrc/nrlmsise00); when it cannot be imported this module falls back to the pure-Python transcription, which is validated against the same C code, so results agree to ~1e-14 either way.

Examples

>>> isinstance(uses_compiled_backend(), bool)
True

Humidity and Dew Point

Humidity conversions and dew-point calculations.

Ports of the humidity functions from the MATLAB Tracker Component Library’s Atmosphere_and_Refraction directory: pairwise conversions between absolute, relative and specific humidity, water number density, and the saturation (dew-point) pressure/temperature of water.

Conventions

Temperatures in Kelvin, pressures in Pascals.
Relative humidity is a fraction in [0, 1], not a percent.
Absolute humidity is kilograms of water per cubic meter of air.
Specific humidity is dimensionless; see the definition parameter.

The dew-point algorithms are shared by every function that touches relative humidity:

0 — corrected Clausius-Clapeyron equation (Koutsoyiannis 2012), for use over land or in the upper air.
1 — Magnus-type equation over water (Alduchov & Eskridge 1996), valid for -40 C to +50 C.
2 — Magnus-type equation over ice (Alduchov & Eskridge 1996), valid for -80 C to 0 C.

References

D. Koutsoyiannis, “Clausius-Clapeyron equation and saturation vapour pressure: simple theory reconciled with practice,” European Journal of Physics, vol. 33, no. 2, pp. 295-305, Mar. 2012.
O. A. Alduchov and R. E. Eskridge, “Improved Magnus form approximation of saturation vapor pressure,” Journal of Applied Meteorology, vol. 35, no. 4, pp. 601-609, Apr. 1996.

pytcl.atmosphere.humidity.H2O_MOLAR_MASS: float = 18.015349999999998: Molar mass of water [g/mol], 2*H + O from the 2013 CIAAW standard atomic weights (interval midpoints), matching the MATLAB TCL Constants class.

pytcl.atmosphere.humidity.abs_humid_to_number_density(abs_humid)[source]

Number density of water molecules from absolute humidity.

Port of absHumid2NumberDensH2O.m.

Parameters:: abs_humid (array_like) – Absolute humidity in kilograms of water per cubic meter of air.
Returns:: number_density – Number of water molecules per cubic meter of air.
Return type:: float or ndarray

Examples

>>> nd = abs_humid_to_number_density(0.01)
>>> print(f"{nd:.6e}")
3.342783e+23

pytcl.atmosphere.humidity.abs_humid_to_rel_humid(abs_humid, temperature, algorithm=0)[source]

Convert absolute humidity to relative humidity.

Inverse of rel_humid_to_abs_humid().

Port of absHumid2RelHumid.m.

Parameters:

abs_humid (array_like) – Absolute humidity in kilograms of water per cubic meter of air.
temperature (array_like) – Temperature(s) in Kelvin.
algorithm (int, optional) – Dew-point algorithm; see dew_point_pressure().

Returns:

rel_humid – Relative humidity as a fraction (0 to 1 for physical inputs).

Return type:

float or ndarray

Examples

>>> round(float(abs_humid_to_rel_humid(0.00641652, 288.15)), 6)
0.5

pytcl.atmosphere.humidity.abs_humid_to_spec_humid(abs_humid, dry_air_density, definition=0)[source]

Convert absolute humidity to specific humidity.

Port of absHumid2SpecHumid.m.

Parameters:

abs_humid (array_like) – Absolute humidity in kilograms of water per cubic meter of air.
dry_air_density (array_like) – Mass density of the dry air (not counting the water) in kg/m^3.
definition (int, optional) – 0 (default): specific humidity is the mass density of water over the mass density of dry air (mixing ratio). 1: mass density of water over the total mass density of the air.

Returns:

spec_humid – Specific humidity under the chosen definition (dimensionless).

Return type:

float or ndarray

Examples

>>> round(float(abs_humid_to_spec_humid(0.00641652, 1.225)), 8)
0.00523798
>>> round(float(abs_humid_to_spec_humid(0.00641652, 1.225, definition=1)), 8)
0.00521068

pytcl.atmosphere.humidity.dew_point_pressure(temperature, algorithm=0)[source]

Saturation partial pressure of water for a given temperature.

This is the partial vapor pressure of water that (in equilibrium) cannot be exceeded — the dew-point pressure. Above it, water begins to condense out of gaseous form at this temperature.

Port of dewPointPres4Temp.m.

Parameters:

temperature (array_like) – Temperature(s) in Kelvin.
algorithm (int, optional) – 0 (default) corrected Clausius-Clapeyron equation; 1 Magnus-type over water (-40 C to +50 C); 2 Magnus-type over ice (-80 C to 0 C). See the module docstring.

Returns:

pressure – Saturation pressure(s) of water in Pascals.

Return type:

float or ndarray

Examples

>>> round(float(dew_point_pressure(288.15)), 4)
1706.632
>>> round(float(dew_point_pressure(288.15, algorithm=1)), 4)
1701.9828
>>> round(float(dew_point_pressure(263.15, algorithm=2)), 4)
259.6718

pytcl.atmosphere.humidity.dew_point_temperature(pressure, algorithm=0)[source]

Temperature at which a partial vapor pressure of water saturates.

For a given partial pressure of water, find the temperature at which that pressure is the saturation pressure — the dew-point temperature. Inverse of dew_point_pressure().

Port of dewPointTemp4Pres.m.

Parameters:

pressure (array_like) – Partial vapor pressure(s) of water in Pascals.
algorithm (int, optional) – 0 (default) corrected Clausius-Clapeyron equation (inverted by the fixed-point iteration of Koutsoyiannis Eqs. 44-45, 27 iterations); 1 Magnus-type over water; 2 Magnus-type over ice. See the module docstring.

Returns:

temperature – Dew-point temperature(s) in Kelvin.

Return type:

float or ndarray

Examples

>>> round(float(dew_point_temperature(1706.632)), 4)
288.15
>>> round(float(dew_point_temperature(1701.9828, algorithm=1)), 4)
288.15

pytcl.atmosphere.humidity.number_density_to_abs_humid(number_density)[source]

Absolute humidity from the number density of water molecules.

Port of numberDensH2O2AbsHumid.m.

Parameters:: number_density (array_like) – Number of water molecules per cubic meter of air.
Returns:: abs_humid – Absolute humidity in kilograms of water per cubic meter of air.
Return type:: float or ndarray

Examples

>>> round(number_density_to_abs_humid(3.342783e+23), 8)
0.01

pytcl.atmosphere.humidity.rel_humid_to_abs_humid(rel_humid, temperature, algorithm=0)[source]

Convert relative humidity to absolute humidity.

Assumes the Ideal Gas Law and Dalton’s Law of Partial Pressures: the partial pressure of water is rel_humid times the saturation pressure at temperature, and the corresponding mass density follows from the ideal gas law.

Port of relHumid2AbsHumid.m.

Parameters:

rel_humid (array_like) – Relative humidity as a fraction in [0, 1].
temperature (array_like) – Temperature(s) in Kelvin.
algorithm (int, optional) – Dew-point algorithm; see dew_point_pressure().

Returns:

abs_humid – Absolute humidity in kilograms of water per cubic meter of air.

Return type:

float or ndarray

Examples

>>> round(float(rel_humid_to_abs_humid(0.5, 288.15)), 8)
0.00641652

pytcl.atmosphere.humidity.rel_humid_to_spec_humid(rel_humid, temperature, dry_air_density, definition=0, algorithm=0)[source]

Convert relative humidity to specific humidity.

Port of relHumid2SpecHumid.m.

Parameters:

rel_humid (array_like) – Relative humidity as a fraction in [0, 1].
temperature (array_like) – Temperature(s) in Kelvin.
dry_air_density (array_like) – Mass density of the dry air (not counting the water) in kg/m^3.
definition (int, optional) – Specific-humidity definition; see abs_humid_to_spec_humid().
algorithm (int, optional) – Dew-point algorithm; see dew_point_pressure().

Returns:

spec_humid – Specific humidity under the chosen definition (dimensionless).

Return type:

float or ndarray

Examples

>>> round(float(rel_humid_to_spec_humid(0.5, 288.15, 1.225)), 8)
0.00523798

pytcl.atmosphere.humidity.spec_humid_to_abs_humid(spec_humid, dry_air_density, definition=0)[source]

Convert specific humidity to absolute humidity.

Inverse of abs_humid_to_spec_humid().

Port of specHumid2AbsHumid.m.

Parameters:

spec_humid (array_like) – Specific humidity (dimensionless).
dry_air_density (array_like) – Mass density of the dry air (not counting the water) in kg/m^3.
definition (int, optional) – 0 (default): specific humidity is water density over dry-air density. 1: water density over total air density.

Returns:

abs_humid – Absolute humidity in kilograms of water per cubic meter of air.

Return type:

float or ndarray

Examples

>>> round(float(spec_humid_to_abs_humid(0.00523798, 1.225)), 8)
0.00641653

pytcl.atmosphere.humidity.spec_humid_to_rel_humid(spec_humid, temperature, dry_air_density, definition=0, algorithm=0)[source]

Convert specific humidity to relative humidity.

Inverse of rel_humid_to_spec_humid().

Port of specHumid2RelHumid.m.

Parameters:

spec_humid (array_like) – Specific humidity (dimensionless).
temperature (array_like) – Temperature(s) in Kelvin.
dry_air_density (array_like) – Mass density of the dry air (not counting the water) in kg/m^3.
definition (int, optional) – Specific-humidity definition; see abs_humid_to_spec_humid().
algorithm (int, optional) – Dew-point algorithm; see dew_point_pressure().

Returns:

rel_humid – Relative humidity as a fraction (0 to 1 for physical inputs).

Return type:

float or ndarray

Examples

>>> round(float(spec_humid_to_rel_humid(0.00523798, 288.15, 1.225)), 6)
0.5

Refractivity

Atmospheric refractivity for radar and optical propagation.

Ports from the MATLAB Tracker Component Library’s Atmosphere_and_Refraction directory: refractivity helpers, astronomical refraction (add/remove with the Sinclair atmosphere model), and the standard-exponential-model radar refraction suite (bistatic r-u-v ray tracing, bias approximation, cubature-based Gaussian conversions).

Conventions

Refractivity is N = (n - 1) * 1e6 where n is the index of refraction. The standard exponential atmosphere model takes the refractivity at height h above sea level to be N = Ns * exp(-ce * h) where Ns is the sea-level refractivity and ce is the decay constant returned by atmos_exp_decay_const().

class pytcl.atmosphere.refraction.AstroRefParams(a, b)[source]

Bases: NamedTuple

Constants of the A*tan(z) + B*tan^3(z) astronomical refraction model.

Variables:

a (float) – The tan(z) coefficient in radians.
b (float) – The tan^3(z) coefficient in radians.

a: float: Alias for field number 0

b: float: Alias for field number 1

class pytcl.atmosphere.refraction.AstroRefractionResult(zenith_distance, delta_z)[source]

Bases: NamedTuple

A zenith distance with the refraction correction that produced it.

Variables:

zenith_distance (float or ndarray) – The converted zenith distance(s) in radians: refraction-free for remove_astro_refraction(), refraction-corrupted for add_astro_refraction(). An empty array signals inputs outside the algorithm’s validity region (matching the MATLAB functions, which return empty matrices).
delta_z (float or ndarray) – The refraction correction in radians that was applied, with z_true = z_observed + delta_z.

zenith_distance: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 0

delta_z: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 1

class pytcl.atmosphere.refraction.CubatureConversionResult(mean, covariance)[source]

Bases: NamedTuple

First two moments of a measurement converted by cubature integration.

Variables:

mean (ndarray) – Converted mean(s), shape (3, N).
covariance (ndarray) – Converted covariance matrices, shape (3, 3, N).

mean: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 0

covariance: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 1

class pytcl.atmosphere.refraction.ExpDecayConstResult(ce, delta_n)[source]

Bases: NamedTuple

Exponential-atmosphere decay constant and 1-km refractivity change.

Variables:

ce (float or ndarray) – Decay constant of the refractivity in inverse meters.
delta_n (float or ndarray) – Change in refractivity going 1 km up from sea level (negative).

ce: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 0

delta_n: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 1

class pytcl.atmosphere.refraction.RuvStdRefracResult(z, u_tx, u_tar_rx, u_tar_tx)[source]

Bases: NamedTuple

Refraction-corrupted bistatic r-u-v measurements and ray directions.

Variables:

z (ndarray) – The measurements, shape (3, N) or (4, N) with include_w: bistatic range then the direction cosines of the apparent target direction in the receiver’s local frame.
u_tx (ndarray) – Unit vectors (3, N) in ECEF pointing from the transmitter toward the refraction-corrupted apparent target position.
u_tar_rx (ndarray) – Unit vectors (3, N) of the apparent direction of the receiver as seen by the target.
u_tar_tx (ndarray) – Unit vectors (3, N) of the apparent direction of the transmitter as seen by the target.

z: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 0

u_tx: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 1

u_tar_rx: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 2

u_tar_tx: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 3

class pytcl.atmosphere.refraction.SinclairAtmosResult(n, dndr, temperature, pressure)[source]

Bases: NamedTuple

Atmospheric parameters of the Sinclair model at the queried heights.

Variables:

n (float or ndarray) – Index of refraction.
dndr (float or ndarray) – Derivative of the index of refraction with respect to height, in inverse meters.
temperature (float or ndarray) – Temperature in Kelvin.
pressure (float or ndarray) – Pressure in Pascals.

n: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 0

dndr: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 1

temperature: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 2

pressure: float | ndarray[tuple[Any, ...], dtype[floating]]: Alias for field number 3

class pytcl.atmosphere.refraction.StdRefracBiasResult(delta_r_one_way, delta_theta)[source]

Bases: NamedTuple

Approximate refraction biases of a monostatic radar measurement.

Variables:

delta_r_one_way (float) – Bias in the one-way range in meters (add to the true range to get the measured range).
delta_theta (float) – Bias in the elevation angle in radians.

delta_r_one_way: float: Alias for field number 0

delta_theta: float: Alias for field number 1

pytcl.atmosphere.refraction.add_astro_refraction(algorithm, obs_lat_lon_alt, z_true, rel_humid=0.0, pressure=101325.0, temperature=288.15, wavelength=5.74e-07)[source]

Add atmospheric refraction to a true zenith distance.

The inverse of remove_astro_refraction(): given true (refraction-free) zenith distances of an object outside the atmosphere, compute the refraction-corrupted apparent zenith distances. The inverse problem is solved by a fixed 20 iterations of the forward model, which is generally sufficient for convergence to working precision.

Port of addAstroRefrac.m.

Parameters:

algorithm (int) – Same choices and validity regions as remove_astro_refraction().
obs_lat_lon_alt (array_like) – The observer’s WGS-84 [latitude, longitude, height], radians and meters. Read only by algorithm 0.
z_true (array_like) – True positive zenith distance(s) in radians.
rel_humid (float, optional) – Relative humidity at the observer in [0, 1]. Default 0.
pressure (float, optional) – Pressure at the observer in Pascals. Default 101325 Pa.
temperature (float, optional) – Temperature at the observer in Kelvin. Default 288.15 K.
wavelength (float, optional) – Observation wavelength in meters. Default 0.574e-6 m.

Returns:

result – zenith_distance holds the refraction-corrupted zenith distances, z0 = z_true - delta_z. Both fields are empty arrays when the point falls outside the algorithm’s validity region during iteration.

Return type:

AstroRefractionResult

Examples

>>> obs = [0.61, 0.0, 100.0]
>>> z_t = 1.200706341
>>> z0, dz = add_astro_refraction(0, obs, z_t, 0.5, 101325.0, 288.15)
>>> print(f"{z0:.6f}")
1.200000

pytcl.atmosphere.refraction.approx_refractivity(temperature, pressure, water_vapor_pressure)[source]

Approximate the refractivity of air from temperature and pressure.

Implements Equation 6 in Annex I of ITU-R P.453-11. The approximation does not depend on frequency.

Port of approxRefractivity.m.

References

International Telecommunication Union, “Recommendation ITU-R P.453-11: The radio refractive index: Its formula and refractivity data,” Tech. Rep., Jul. 2015.

Parameters:

temperature (array_like) – Temperature(s) in Kelvin.
pressure (array_like) – Total atmospheric pressure(s) in Pascals (dry pressure plus the partial pressure of water vapor).
water_vapor_pressure (array_like) – Partial pressure(s) of water vapor in Pascals.

Returns:

refractivity – The refractivity N = (n - 1) * 1e6 of the atmosphere.

Return type:

float or ndarray

Examples

>>> round(float(approx_refractivity(288.15, 101325.0, 853.3)), 4)
311.2452

pytcl.atmosphere.refraction.atmos_exp_decay_const(ns)[source]

Decay constant of the exponential refractivity model.

Given the refractivity of air at sea level, obtain the approximate decay constant for refractivity as a function of height per Appendix A of the CRPL Exponential Reference Atmosphere, along with the change in refractivity 1 km above sea level. The refractivity at height h above sea level is then N = ns * exp(-ce * h).

Port of atmosExpDecayConst4Refrac.m.

References

B. R. Bean and G. D. Thayer, CRPL Exponential Reference Atmosphere. Washington, D.C.: U.S. Department of Commerce, National Bureau of Standards, Oct. 1959.

Parameters:: ns (array_like) – Refractivity of air at sea level.
Returns:: result – Named tuple of ce (decay constant, inverse meters) and delta_n (refractivity change 1 km up from sea level), each with the same shape as ns.
Return type:: ExpDecayConstResult

Examples

>>> ce, delta_n = atmos_exp_decay_const(313.0)
>>> print(f"{ce:.6e}")
1.438586e-04
>>> round(delta_n, 4)
-41.9388

pytcl.atmosphere.refraction.cart2ruv_std_refrac(z_c, use_half_range=False, z_tx=None, z_rx=None, m=None, ns=313.0, include_w=False, ce=None, r_e=None, sphere_center=None)[source]

Convert Cartesian points to refraction-corrupted bistatic r-u-v.

Traces rays through the standard exponential atmosphere (N = ns * exp(-ce * h)) over a locally osculating spherical Earth to determine the apparent bistatic range and direction cosines of each target as seen by the receiver. Not suitable for satellite-to-satellite paths grazing the atmosphere, and the ray tracer can fail for paths going too far underground or for targets collocated with the receiver or the transmitter.

Port of Cart2RuvStdRefrac.m.

Parameters:

z_c (array_like) – Cartesian target positions in global ECEF coordinates, shape (3, N) or (3,).
use_half_range (bool, optional) – Whether the bistatic range is halved (one-way range in the monostatic case). Default False.
z_tx (array_like, optional) – Transmitter ECEF position, shape (3,). Default: the origin.
z_rx (array_like, optional) – Receiver ECEF position, shape (3,). Default: the origin.
m (array_like, optional) – 3x3 rotation from global axes to the receiver’s local axes (the receiver boresight is its local z axis). Default: identity.
ns (float, optional) – Refractivity reduced to the reference sphere. Default 313.
include_w (bool, optional) – Include the third direction cosine, making z 4xN. Default False.
ce (float, optional) – Decay constant of the exponential model in inverse meters. Default: derived from ns via the CRPL standard constants (atmos_exp_decay_const()).
r_e (float, optional) – Radius of the spherical-Earth approximation. Default: the osculating sphere at the receiver (pytcl.navigation.geodesy.osculating_sphere()).
sphere_center (array_like, optional) – ECEF offset of the sphere’s center. Defaults to the osculating sphere’s offset when r_e is defaulted, and to zeros when r_e is given.

Returns:

result – The measurements and the apparent ray directions at both ends.

Return type:

RuvStdRefracResult

Examples

>>> z_rx = np.array([6378137.0, 0.0, 0.0])
>>> z_tar = np.array([6428137.0, 100e3, 0.0])
>>> res = cart2ruv_std_refrac(z_tar, True, z_rx, z_rx)
>>> print(f"{res.z[0, 0]:.3f}")
111808.239

pytcl.atmosphere.refraction.cart2ruv_std_refrac_cubature(z_c, sqrt_cov, use_half_range=False, z_tx=None, z_rx=None, m=None, ns=313.0, points=None, weights=None, ce=None, r_e=None, sphere_center=None)[source]

Cubature-based Gaussian conversion of Cartesian states to r-u-v.

Propagates Gaussian state estimates through cart2ruv_std_refrac() by cubature integration, returning the converted means and covariances.

Port of Cart2RuvStdRefracCubature.m.

Parameters:

z_c (array_like) – Cartesian means, shape (3, N) or (3,).
sqrt_cov (array_like) – Lower-triangular square roots of the covariances, shape (3, 3) (shared) or (3, 3, N).
use_half_range (bool) – As in cart2ruv_std_refrac().
z_tx (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
z_rx (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
m (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
ns (float) – As in cart2ruv_std_refrac().
ce (float | None) – As in cart2ruv_std_refrac().
r_e (float | None) – As in cart2ruv_std_refrac().
sphere_center (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
points (array_like, optional) – Cubature points (num_points, 3) and weights for N(0, I). Default: fifth_order_cubature_points().
weights (array_like, optional) – Cubature points (num_points, 3) and weights for N(0, I). Default: fifth_order_cubature_points().

Returns:

result – Converted means (3, N) and covariances (3, 3, N).

Return type:

CubatureConversionResult

Examples

>>> z_rx = np.array([6378137.0, 0.0, 0.0])
>>> z_tar = np.array([6428137.0, 100e3, 0.0])
>>> sr = np.diag([100.0, 100.0, 100.0])
>>> res = cart2ruv_std_refrac_cubature(z_tar, sr, True, z_rx, z_rx)
>>> print(f"{res.mean[0, 0]:.1f}")
111808.3

pytcl.atmosphere.refraction.reduce_std_refrac_to_sphere(n_meas, height, exp_const=0.005577, mult_const=7.32, xatol=0.0001)[source]

Reduce a measured refractivity to the reference sphere.

Given the atmospheric refractivity measured at a height above sea level, determine the equivalent sea-level refractivity under the standard exponential model. The model is scanned over sea-level refractivities in [200, 450]; each sign change brackets a candidate solution refined by bounded scalar minimization, so one or two solutions can be returned.

Port of reduceStdRefrac2Spher.m.

Parameters:

n_meas (float) – The measured refractivity (n - 1) * 1e6.
height (float) – Height of the measurement above mean sea level in meters.
exp_const (float, optional) – Parameters of the decay model deltaN = -mult_const * exp(exp_const * N) per kilometer. Defaults are the CRPL standard values 0.005577 and 7.32.
mult_const (float, optional) – Parameters of the decay model deltaN = -mult_const * exp(exp_const * N) per kilometer. Defaults are the CRPL standard values 0.005577 and 7.32.
xatol (float, optional) – Absolute tolerance of the bounded search (MATLAB fminbnd default 1e-4).

Returns:

ns_values – The candidate sea-level refractivities, shape (num_solutions,). May be empty if the measurement is inconsistent with the model.

Return type:

ndarray

Examples

>>> vals = reduce_std_refrac_to_sphere(300.0, 1000.0)
>>> print(f"{vals[0]:.4f}")
352.1814

pytcl.atmosphere.refraction.remove_astro_refraction(algorithm, obs_lat_lon_alt, z_observed, rel_humid=0.0, pressure=101325.0, temperature=288.15, wavelength=5.74e-07)[source]

Remove atmospheric refraction from an observed zenith distance.

Given refraction-corrupted zenith distances of an object outside the atmosphere seen by a near-surface observer, compute the true (refraction-free) zenith distances using low-precision atmospheric models.

Port of removeAstroRefrac.m.

Parameters:

algorithm (int) – 0 numerical ray integration through the Sinclair atmosphere (the most precise; valid for observed zenith distances up to 100 degrees and observers below the tropopause); 1 the Saastamoinen formula (sea-level observer, zenith distances below 70 degrees); 2 the IAU A*tan(z) + B*tan^3(z) model with Newton-Raphson correction (sea-level observer, any positive zenith distance, degrading near the horizon).
obs_lat_lon_alt (array_like) – The observer’s WGS-84 [latitude, longitude, height], radians and meters. Only algorithm 0 reads it (latitude and height; the longitude does not matter). Ignored by algorithms 1 and 2.
z_observed (array_like) – Refraction-corrupted positive zenith distance(s) in radians, measured down from the local vertical.
rel_humid (float, optional) – Relative humidity at the observer in [0, 1]. Default 0.
pressure (float, optional) – Pressure at the observer in Pascals. Default 101325 Pa.
temperature (float, optional) – Temperature at the observer in Kelvin. Default 288.15 K.
wavelength (float, optional) – Observation wavelength in meters. Default 0.574e-6 m.

Returns:

result – zenith_distance holds the true zenith distances, z_true = z_observed + delta_z. Both fields are empty arrays when any input is outside the algorithm’s validity region (mirroring the MATLAB function’s empty-matrix return).

Return type:

AstroRefractionResult

Raises:

ValueError – If any observed zenith distance is negative, the algorithm is unknown, or (algorithm 0) the observer is above the tropopause.

Examples

>>> obs = [0.61, 0.0, 100.0]
>>> z_true, dz = remove_astro_refraction(0, obs, 1.2, 0.5, 101325.0, 288.15)
>>> print(f"{z_true:.9f} {dz:.4e}")
1.200706341 7.0634e-04
>>> z_true, dz = remove_astro_refraction(2, obs, 1.2, 0.5, 101325.0, 288.15)
>>> print(f"{z_true:.9f}", dz > 0)
1.200705016 True

pytcl.atmosphere.refraction.ruv2cart_std_refrac(z_ruv, use_half_range=False, z_tx=None, z_rx=None, m=None, ns=313.0, ce=None, r_e=None, sphere_center=None, x_max=1000000.0)[source]

Convert refraction-corrupted bistatic r-u-v points to Cartesian.

The inverse of cart2ruv_std_refrac(): shoots a ray from the receiver in the apparent direction through the exponential atmosphere (an initial value problem) and searches along the traced path for the point whose accumulated bistatic range matches the measurement. Fails if the target is collocated with the transmitter or the receiver.

Port of ruv2CartStdRefrac.m. Deviation from MATLAB: the MATLAB function’s near-vertical branch aborts the whole measurement loop (an upstream bug); this port processes remaining measurements.

Parameters:

z_ruv (array_like) – Measurements, shape (3, N) or (4, N) ([r; u; v] or [r; u; v; w]), or a single measurement of shape (3,) or (4,).
use_half_range (bool, optional) – Whether the ranges in z_ruv are halved. Default False.
z_tx (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
z_rx (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
m (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
ns (float) – As in cart2ruv_std_refrac().
ce (float | None) – As in cart2ruv_std_refrac().
r_e (float | None) – As in cart2ruv_std_refrac().
sphere_center (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in cart2ruv_std_refrac().
x_max (float, optional) – Maximum horizontal displacement searched along the ray, in meters. Default 1000 km.

Returns:

z_cart – Cartesian target positions in global ECEF coordinates, (3, N).

Return type:

ndarray

Examples

>>> z_rx = np.array([6378137.0, 0.0, 0.0])
>>> z_tar = np.array([6428137.0, 100e3, 0.0])
>>> ruv = cart2ruv_std_refrac(z_tar, True, z_rx, z_rx).z
>>> back = ruv2cart_std_refrac(ruv, True, z_rx, z_rx)
>>> np.allclose(back[:, 0], z_tar, atol=0.5)
True

pytcl.atmosphere.refraction.ruv2cart_std_refrac_cubature(z_ruv, sqrt_cov, use_half_range=False, z_tx=None, z_rx=None, m=None, ns=313.0, points=None, weights=None, ce=None, r_e=None, sphere_center=None, x_max=1000000.0)[source]

Cubature-based Gaussian conversion of r-u-v measurements to Cartesian.

Propagates Gaussian measurements through ruv2cart_std_refrac() by cubature integration.

Port of ruv2CartStdRefracCubature.m.

Parameters:

z_ruv (array_like) – Measurement means, shape (3, N) or (3,).
sqrt_cov (array_like) – Lower-triangular square roots of the measurement covariances, shape (3, 3) or (3, 3, N).
use_half_range (bool) – As in ruv2cart_std_refrac().
z_tx (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in ruv2cart_std_refrac().
z_rx (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in ruv2cart_std_refrac().
m (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in ruv2cart_std_refrac().
ns (float) – As in ruv2cart_std_refrac().
ce (float | None) – As in ruv2cart_std_refrac().
r_e (float | None) – As in ruv2cart_std_refrac().
sphere_center (_Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None) – As in ruv2cart_std_refrac().
x_max (float) – As in ruv2cart_std_refrac().
points (array_like, optional) – Cubature points and weights for N(0, I). Default: fifth_order_cubature_points().
weights (array_like, optional) – Cubature points and weights for N(0, I). Default: fifth_order_cubature_points().

Returns:

result – Converted means (3, N) and covariances (3, 3, N).

Return type:

CubatureConversionResult

Examples

>>> z_rx = np.array([6378137.0, 0.0, 0.0])
>>> z_tar = z_rx + np.array([1e3, 5e3, 50e3])
>>> ruv = cart2ruv_std_refrac(z_tar, True, z_rx, z_rx).z[:, 0]
>>> sr = np.diag([10.0, 1e-4, 1e-4])
>>> res = ruv2cart_std_refrac_cubature(ruv, sr, True, z_rx, z_rx)
>>> np.allclose(res.mean[:, 0], z_tar, atol=5.0)
True

pytcl.atmosphere.refraction.simple_astro_ref_params(rel_humid=0.0, pressure=101325.0, temperature=288.15, wavelength=5.74e-07)[source]

Constants for the simple A*tan(z) + B*tan^3(z) refraction model.

Port of simpAstroRefParam.m (a MEX wrapper in MATLAB). This function uses computations derived from the IAU SOFA refco routine; it is not itself software provided by or endorsed by SOFA. It differs from the original in taking SI inputs (Pascals, Kelvin, meters) and converting internally to the hPa/Celsius/micrometer units the fit was built for.

Parameters:

rel_humid (float, optional) – Relative humidity at the observer as a fraction in [0, 1]. Default 0.
pressure (float, optional) – Atmospheric pressure at the observer in Pascals. Default 101325 Pa.
temperature (float, optional) – Air temperature at the observer in Kelvin. Default 288.15 K.
wavelength (float, optional) – Observation wavelength in meters; values above 100 micrometers select the radio-frequency fit instead of the optical/IR one. Default 0.574e-6 m (yellow light).

Returns:

params – The coefficients a and b in radians.

Return type:

AstroRefParams

Examples

>>> a, b = simple_astro_ref_params(0.5, 101325.0, 288.15)
>>> print(f"{a:.9e} {b:.9e}")
2.767559220e-04 -3.167276124e-07

pytcl.atmosphere.refraction.sinclair_atmosphere(height, obs_lat_lon_alt, rel_humid=0.0, pressure=101325.0, temperature=288.15, wavelength=5.74e-07, tropopause_height=11000.0)[source]

Atmospheric parameters for the Sinclair refraction model.

A two-layer troposphere/stratosphere model of the index of refraction and its height derivative, used by algorithm 0 of remove_astro_refraction(). Follows Chapter 7.2 of Hohenkerk’s treatment in the Explanatory Supplement and the 1982 HM Nautical Almanac Office technical note.

Port of SinclairAtmos.m.

Parameters:

height (array_like) – Height(s) above the reference ellipsoid in meters at which to evaluate the model.
obs_lat_lon_alt (array_like) – The observer’s WGS-84 [latitude, longitude, height] with angles in radians and height in meters. Only the latitude and height are read.
rel_humid (float, optional) – Relative humidity at the observer in [0, 1]. Default 0.
pressure (float, optional) – Pressure at the observer in Pascals. Default 101325 Pa.
temperature (float, optional) – Temperature at the observer in Kelvin. Default 288.15 K.
wavelength (float, optional) – Observation wavelength in meters. Default 0.574e-6 m.
tropopause_height (float, optional) – Assumed top of the troposphere in meters. Default 11000 m.

Returns:

result – Named tuple of n, dndr, temperature and pressure at each queried height.

Return type:

SinclairAtmosResult

References

C. Y. Hohenkerk and A. T. Sinclair, “The computation of angular atmospheric refraction at large zenith angles,” HM Nautical Almanac Office, Tech. Rep. NAO TN No. 63, Apr. 1985.

Examples

>>> obs = [0.61, 0.0, 100.0]
>>> res = sinclair_atmosphere(1000.0, obs, 0.5, 101325.0, 288.15)
>>> print(f"{res.n:.9f} {res.temperature:.2f}")
1.000254179 282.30

pytcl.atmosphere.refraction.std_refrac_bias_approx(path_length, elevation, radar_height, ns=313.0, ce=None, r_e=None, algorithm=1)[source]

Approximate range and elevation biases due to standard refraction.

For a monostatic radar at a given height observing a target at a given path length and elevation, approximate the offsets that refraction through the standard exponential atmosphere adds to the one-way range and to the elevation angle.

Port of stdRefracBiasApprox.m.

Parameters:

path_length (float) – Length of the refraction-free path to the target in meters.
elevation (float) – Elevation angle of the target above the radar’s local horizontal in radians.
radar_height (float) – Height of the radar above the reference sphere in meters.
ns (float, optional) – Refractivity reduced to the reference sphere. Default 313.
ce (float, optional) – Decay constant in inverse meters. Default: CRPL standard.
r_e (float, optional) – Radius of the reference sphere. Default: the WGS-84 mean radius (2a + b) / 3.
algorithm (int, optional) – 1 (default) numerical ray tracing (the same BVP solve as cart2ruv_std_refrac()); 0 the closed-form Kerce-Blair-Brown approximation, valid for elevations up to 49 degrees.

Returns:

result – The one-way range bias in meters and the elevation bias in radians.

Return type:

StdRefracBiasResult

References

J. C. Kerce, W. D. Blair, and G. C. Brown, “Modeling refraction errors for simulation studies of multisensor target tracking,” Proc. 36th Southeastern Symposium on System Theory, Mar. 2004.

Examples

>>> res = std_refrac_bias_approx(100e3, 0.1, 100.0)
>>> print(f"{res.delta_r_one_way:.4f} {res.delta_theta:.6e}")
15.9561 1.420764e-03
>>> res0 = std_refrac_bias_approx(100e3, 0.1, 100.0, algorithm=0)
>>> print(f"{res0.delta_r_one_way:.4f}")
15.8175

Ionosphere

Ionospheric models for radio propagation and navigation applications.

This module provides ionospheric models used for computing signal delays, electron density profiles, and Total Electron Content (TEC) estimates. These are essential for GPS/GNSS corrections and radio wave propagation.

Models

Klobuchar: GPS broadcast ionospheric model (single-frequency correction)
IRI: International Reference Ionosphere simplified model (simple_iri)
Dual-frequency TEC estimation and TEC-based delay

References

Klobuchar, J.A. (1987). “Ionospheric Time-Delay Algorithm for Single-Frequency GPS Users”. IEEE Transactions on Aerospace and Electronic Systems, AES-23(3), 325-331.

class pytcl.atmosphere.ionosphere.IonosphereState(tec, delay_l1, delay_l2, f_peak, h_peak)[source]

Bases: NamedTuple

Ionospheric state at a given location and time.

Variables:

tec (float or ndarray) – Total Electron Content in TECU (10^16 electrons/m²).
delay_l1 (float or ndarray) – Ionospheric delay at L1 frequency in meters.
delay_l2 (float or ndarray) – Ionospheric delay at L2 frequency in meters.
f_peak (float or ndarray) – Critical frequency of F2 layer in MHz.
h_peak (float or ndarray) – Height of F2 layer peak in km.

tec: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 0

delay_l1: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 1

delay_l2: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 2

f_peak: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 3

h_peak: float | ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 4

class pytcl.atmosphere.ionosphere.KlobucharCoefficients(alpha, beta)[source]

Bases: NamedTuple

Klobuchar ionospheric model coefficients.

These coefficients are broadcast by GPS satellites in the navigation message.

Variables:

alpha (ndarray) – Amplitude coefficients (4 values) in seconds.
beta (ndarray) – Period coefficients (4 values) in seconds.

alpha: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 0

beta: ndarray[tuple[Any, ...], dtype[float64]]: Alias for field number 1

pytcl.atmosphere.ionosphere.klobuchar_delay(latitude, longitude, elevation, azimuth, gps_time, coefficients=None)[source]

Compute ionospheric delay using the Klobuchar model.

The Klobuchar model is the standard GPS broadcast ionospheric correction model. It provides single-frequency ionospheric delay estimates accurate to about 50% RMS.

Parameters:

latitude (array_like) – User geodetic latitude in radians.
longitude (array_like) – User geodetic longitude in radians.
elevation (array_like) – Satellite elevation angle in radians.
azimuth (array_like) – Satellite azimuth angle in radians.
gps_time (array_like) – GPS time of week in seconds.
coefficients (KlobucharCoefficients, optional) – Ionospheric coefficients from GPS navigation message. If None, uses default mid-latitude values.

Returns:

delay – Ionospheric delay in meters (at L1 frequency).

Return type:

ndarray

Examples

>>> # User at 40°N, 105°W, satellite at 45° elevation
>>> delay = klobuchar_delay(
...     np.radians(40), np.radians(-105),
...     np.radians(45), np.radians(180),
...     gps_time=43200  # Noon
... )
>>> delay > 0
True

Notes

The Klobuchar model assumes a thin-shell ionosphere at 350 km altitude and uses a cosine model for diurnal variation. It typically removes about 50% of the ionospheric delay.

References

IS-GPS-200, Interface Specification.

pytcl.atmosphere.ionosphere.dual_frequency_tec(pseudorange_l1, pseudorange_l2)[source]

Compute Total Electron Content from dual-frequency pseudoranges.

This method uses the dispersive nature of the ionosphere to estimate TEC from the difference in L1 and L2 pseudoranges.

Parameters:

pseudorange_l1 (array_like) – L1 pseudorange in meters.
pseudorange_l2 (array_like) – L2 pseudorange in meters.

Returns:

tec – Total Electron Content in TECU (10^16 electrons/m²).

Return type:

ndarray

Examples

>>> # Pseudorange measurements from dual-frequency receiver
>>> p_l1 = 22000000.0  # L1 pseudorange in meters
>>> p_l2 = 22000002.5  # L2 pseudorange (slightly delayed)
>>> tec = dual_frequency_tec(p_l1, p_l2)
>>> tec > 0  # TEC should be positive
True

Notes

The ionospheric delay is proportional to TEC and inversely proportional to frequency squared:

delay = 40.3 * TEC / f²

The difference in delays at L1 and L2 gives:

P2 - P1 = 40.3 * TEC * (1/f1² - 1/f2²)

This is the standard dual-frequency ionospheric correction method.

pytcl.atmosphere.ionosphere.ionospheric_delay_from_tec(tec, frequency=1575420000.0)[source]

Compute ionospheric delay from Total Electron Content.

Parameters:

tec (array_like) – Total Electron Content in TECU (10^16 electrons/m²).
frequency (float, optional) – Signal frequency in Hz. Default is GPS L1.

Returns:

delay – Ionospheric delay in meters.

Return type:

ndarray

Examples

>>> # Typical mid-latitude TEC of 20 TECU
>>> delay = ionospheric_delay_from_tec(20.0)
>>> delay > 0
True
>>> # Delay at L2 is larger than at L1 (lower frequency)
>>> delay_l2 = ionospheric_delay_from_tec(20.0, frequency=F_L2)
>>> delay_l2 > delay
True

Notes

The ionospheric delay for a signal is:: delay = 40.3 * TEC * 10^16 / f²

pytcl.atmosphere.ionosphere.simple_iri(latitude, longitude, altitude, hour, month=6, solar_flux=150.0)[source]

Empirical, latitude-and-time-only approximation of F2-layer TEC.

This is not a reduced IRI model in the sense of evaluating a coarser electron-density profile: it fits the F2 critical frequency (foF2) and total electron content (TEC) from latitude, local hour, month and solar flux alone, using no explicit vertical electron density profile. altitude and longitude are accepted for a call signature that resembles a full IRI evaluation, but neither one affects the result – a warning names both on every call. TEC is identical at 100 km and at 1000 km, and identical at any two longitudes. For an altitude- or longitude-sensitive model, use the full IRI-2020 model or an external service; making this one genuinely altitude/longitude-dependent needs an IRI oracle this repository does not have.

Parameters:

latitude (array_like) – Geodetic latitude in radians.
longitude (array_like) – Geodetic longitude in radians. Accepted but ignored – see Warns.
altitude (array_like) – Altitude in meters. Accepted but ignored – see Warns.
hour (array_like) – Local hour (0-24).
month (int, optional) – Month of year (1-12). Default is 6 (June).
solar_flux (float, optional) – F10.7 solar flux in SFU. Default is 150 (moderate activity).

Returns:

state – Ionospheric state with TEC, delays, and F2 layer parameters.

Return type:

IonosphereState

Warns:

UserWarning – Always: altitude and longitude are accepted but do not affect the returned state.

Notes

This is a simplified empirical model suitable for educational purposes and rough order-of-magnitude estimates only, not for operational use.

Examples

>>> import warnings
>>> with warnings.catch_warnings():
...     warnings.simplefilter("ignore")
...     state = simple_iri(np.radians(40), np.radians(-105), 300e3, 12)
>>> state.tec > 0
True

pytcl.atmosphere.ionosphere.magnetic_latitude(latitude, longitude)[source]

Compute approximate geomagnetic latitude.

Uses a simple dipole approximation with the magnetic pole at approximately 80.5°N, 72.8°W.

Parameters:

latitude (array_like) – Geodetic latitude in radians.
longitude (array_like) – Geodetic longitude in radians.

Returns:

mag_lat – Geomagnetic latitude in radians.

Return type:

ndarray

Examples

>>> import numpy as np
>>> # New York City (40.7°N, 74°W)
>>> mag_lat = magnetic_latitude(np.radians(40.7), np.radians(-74))
>>> np.degrees(mag_lat)
50.19...
>>> # Equator at 0° longitude
>>> mag_lat_eq = magnetic_latitude(0.0, 0.0)
>>> np.abs(np.degrees(mag_lat_eq)) < 15  # Near magnetic equator
True

pytcl.atmosphere.ionosphere.scintillation_index(magnetic_latitude, hour, kp_index=3.0)[source]

Estimate ionospheric scintillation index S4.

Provides a rough estimate of amplitude scintillation based on geomagnetic latitude, local time, and geomagnetic activity.

Parameters:

magnetic_latitude (array_like) – Geomagnetic latitude in radians.
hour (array_like) – Local hour (0-24).
kp_index (float, optional) – Kp geomagnetic activity index (0-9). Default is 3 (moderate).

Returns:

s4 – S4 amplitude scintillation index (0-1).

Return type:

ndarray

Examples

>>> import numpy as np
>>> # Equatorial region at night (high scintillation risk)
>>> s4 = scintillation_index(np.radians(10), 21, kp_index=5.0)
>>> s4 > 0.3  # Moderate to strong scintillation
True
>>> # Mid-latitude during daytime (low scintillation)
>>> s4_low = scintillation_index(np.radians(45), 12, kp_index=1.0)
>>> s4_low < 0.2
True

Notes

S4 > 0.3 indicates moderate scintillation. S4 > 0.6 indicates strong scintillation that may affect receivers.