118 lines
3.4 KiB
Python
118 lines
3.4 KiB
Python
from typing import Callable
|
|
|
|
LAPSE_RATE_DRY = 9.8 / 1000 # degrees C per km
|
|
LAPSE_RATE_MOIST = 4.0 / 1000
|
|
|
|
PRESSURE_MAX = 1050 # millibar
|
|
PRESSURE_MIN = 100
|
|
PRESSURE_STEP = 50
|
|
|
|
def kelvin(c: float) -> float:
|
|
return 273.15 + c
|
|
|
|
def vapor_pressure(dewpoint: float) -> float:
|
|
"""
|
|
Return the pressure, in millibar, of vapor in a parcel of a given
|
|
dewpoint.
|
|
"""
|
|
return 6.11 * 10 ** ((7.5 * dewpoint) / (237.3 + dewpoint))
|
|
|
|
def saturated_vapor_pressure(temp: float) -> float:
|
|
"""
|
|
Return the pressure, in millbar, at which the vapor and condensation of a
|
|
parcel of a given temperature are at equilibrium.
|
|
"""
|
|
return 6.11 * 10 ** ((7.5 * temp) / (237.3 + temp))
|
|
|
|
def mixing_ratio(dewpoint: float, pressure: float) -> float:
|
|
"""
|
|
Return the amount, in kilograms, of water vapor versus dry air in a parcel
|
|
of a given dewpoint and pressure.
|
|
"""
|
|
e = vapor_pressure(dewpoint)
|
|
|
|
return (0.62197 * e) / (pressure - e)
|
|
|
|
def saturated_mixing_ratio(temp: float, pressure: float) -> float:
|
|
"""
|
|
Return the maximum amount, in kilograms, of water vapor a parcel of a
|
|
given temperature and pressure can hold.
|
|
"""
|
|
es = saturated_vapor_pressure(temp)
|
|
|
|
return (0.62197 * es) / (pressure - es)
|
|
|
|
def lcl(temp: float, dewpoint: float) -> float:
|
|
"""
|
|
Return the height, in meters, at which a parcel of the given temperature
|
|
is cooled to the given dewpoint.
|
|
"""
|
|
return (temp - dewpoint) / 0.008
|
|
|
|
def pressure_height(pressure: float) -> float:
|
|
"""
|
|
Return the approximate altitude, in meters, for a given pressure in
|
|
millibar.
|
|
"""
|
|
return (1 - (pressure / 1013.25) ** 0.190284) * 145366.45 * 0.3048
|
|
|
|
def lapse(temp: float, delta: float, rate=LAPSE_RATE_DRY) -> float:
|
|
"""
|
|
Return the temperature of a parcel cooled at the dry lapse rate for a
|
|
given increase in height (in meters).
|
|
"""
|
|
return temp - (rate * delta)
|
|
|
|
def moist_lapse_rate(temp: float, pressure: float) -> float:
|
|
g = 9.8076
|
|
Hv = 2501000
|
|
Rsd = 287
|
|
Rsw = 461.5
|
|
Cpd = 1003.5
|
|
r = saturated_mixing_ratio(temp, pressure)
|
|
T = kelvin(temp)
|
|
|
|
return g * (1 + (Hv * r) / (Rsd * T)) \
|
|
/ (Cpd + (((Hv**2) * r) / (Rsw * (T**2))))
|
|
|
|
def loft_parcel(start_temp: float,
|
|
start_pressure: float,
|
|
lapse_rate: Callable):
|
|
"""
|
|
Loft a parcel of air from a given pressure, at a given temperature,
|
|
yielding a Tuple containing the temperature and pressure of that parcel
|
|
at each increment of elevation. A Callable which, given a temperature,
|
|
dewpoint, pressure value, will be dispatched to obtain the lapse rate
|
|
for each height.
|
|
"""
|
|
start_height = pressure_height(start_pressure)
|
|
sx_last = None
|
|
sy_last = None
|
|
height_last = None
|
|
|
|
temp = start_temp
|
|
pressure = start_pressure
|
|
|
|
while pressure >= PRESSURE_MIN:
|
|
height = pressure_height(pressure)
|
|
|
|
if height_last is None:
|
|
height_last = height
|
|
|
|
try:
|
|
rate = lapse_rate(temp, pressure)
|
|
except OverflowError:
|
|
break
|
|
|
|
temp = lapse(temp, height - height_last, rate)
|
|
|
|
yield temp, pressure
|
|
|
|
height_last = height
|
|
|
|
if pressure == PRESSURE_MIN:
|
|
break
|
|
elif pressure - 10.0 < PRESSURE_MIN:
|
|
pressure = PRESSURE_MIN
|
|
else:
|
|
pressure -= 10.0
|