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) and speed of sound. What remains unported from MATLAB’s Atmosphere_and_Refraction (NRLMSISE-00, Jacchia 1971) is accounted for 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. Valid from 0 to ~86 km.
Returns:: state – Atmospheric state containing temperature, pressure, density, and speed of sound.
Return type:: AtmosphereState

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 valid from sea level 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

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 is an approximate inversion of the ISA model, valid primarily in the troposphere.

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).

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. NRLMSISE-00 requires harmonic coefficient tables from NOAA which this library does not distribute; 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.

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; the real model requires harmonic coefficient tables from NOAA that this library does not ship. 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.

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

>>> 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
>>> 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³

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]

Simplified International Reference Ionosphere (IRI) model.

This provides approximate electron density and TEC values based on simplified IRI physics. For accurate predictions, use the full IRI model or external services.

Parameters:

latitude (array_like) – Geodetic latitude in radians.
longitude (array_like) – Geodetic longitude in radians.
altitude (array_like) – Altitude in meters.
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

Notes

This is a simplified empirical model suitable for educational purposes and rough estimates. For operational use, the full IRI-2020 model should be employed.

Examples

>>> 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.