Compare commits

..

No commits in common. "c0ed20cce50b02d2b21fdf8e253e6f8f5cf69591" and "796257194b57c66066c45a126dd454f5fae6b7b2" have entirely different histories.

2 changed files with 35 additions and 57 deletions

View file

@ -4,9 +4,12 @@ import cairo
from typing import Callable from typing import Callable
from xmet.sounding import Sounding from xmet.sounding import Sounding
from xmet.thermo import pressure_height, loft_parcel, moist_lapse_rate, \ from xmet.thermo import pressure_height, lapse, moist_lapse_rate, \
LAPSE_RATE_DRY, PRESSURE_MAX, PRESSURE_MIN, \ LAPSE_RATE_DRY
PRESSURE_STEP
PRESSURE_MAX = 1050 # millibar
PRESSURE_MIN = 100
PRESSURE_STEP = 50
PRESSURE_LOG_MAX = math.log(PRESSURE_MAX) PRESSURE_LOG_MAX = math.log(PRESSURE_MAX)
PRESSURE_LOG_MIN = math.log(PRESSURE_MIN) PRESSURE_LOG_MIN = math.log(PRESSURE_MIN)
@ -104,11 +107,26 @@ class SkewTGraph():
start_temp: float, start_temp: float,
start_pressure: float, start_pressure: float,
lapse_rate: Callable): lapse_rate: Callable):
sx_last = None start_height = pressure_height(start_pressure)
sy_last = None sx_last = None
sy_last = None
height_last = None
for level in loft_parcel(start_temp, start_pressure, lapse_rate): temp = start_temp
temp, pressure = level pressure = start_pressure
while pressure >= PRESSURE_MIN:
height = pressure_height(pressure)
if height_last is None:
height_last = height
try:
rate = lapse_rate(temp, temp, pressure)
except OverflowError:
break
temp = lapse(temp, height - height_last, rate)
sx, sy = self.sample_to_screen(temp, pressure) sx, sy = self.sample_to_screen(temp, pressure)
@ -117,8 +135,16 @@ class SkewTGraph():
cr.line_to(x + sx_last, y + sy_last) cr.line_to(x + sx_last, y + sy_last)
cr.stroke() cr.stroke()
sx_last = sx sx_last = sx
sy_last = sy sy_last = sy
height_last = height
if pressure == PRESSURE_MIN:
break
elif pressure - 10.0 < PRESSURE_MIN:
pressure = PRESSURE_MIN
else:
pressure -= 10.0
def draw_adiabats(self, def draw_adiabats(self,
cr: cairo.Context, cr: cairo.Context,

View file

@ -1,12 +1,6 @@
from typing import Callable
LAPSE_RATE_DRY = 9.8 / 1000 # degrees C per km LAPSE_RATE_DRY = 9.8 / 1000 # degrees C per km
LAPSE_RATE_MOIST = 4.0 / 1000 LAPSE_RATE_MOIST = 4.0 / 1000
PRESSURE_MAX = 1050 # millibar
PRESSURE_MIN = 100
PRESSURE_STEP = 50
def kelvin(c: float) -> float: def kelvin(c: float) -> float:
return 273.15 + c return 273.15 + c
@ -73,45 +67,3 @@ def moist_lapse_rate(temp: float, dewpoint: float, pressure: float) -> float:
return g * (1 + (Hv * r) / (Rsd * T)) \ return g * (1 + (Hv * r) / (Rsd * T)) \
/ (Cpd + (((Hv**2) * r) / (Rsw * (T**2)))) / (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, 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