use super::{Layer, Layers};
use crate::{
error::{
AnalysisError::{InvalidInput, MissingProfile},
{AnalysisError, Result},
},
interpolation::{linear_interp, linear_interpolate_sounding},
sounding::Sounding,
};
use itertools::{izip, Itertools};
use metfor::{Celsius, HectoPascal, FREEZING};
use optional::Optioned;
pub fn dendritic_snow_zone(snd: &Sounding) -> Result<Layers> {
temperature_layer(snd, Celsius(-12.0), Celsius(-18.0), HectoPascal(300.0))
}
pub fn hail_growth_zone(snd: &Sounding) -> Result<Layers> {
temperature_layer(snd, Celsius(-10.0), Celsius(-30.0), HectoPascal(1.0))
}
#[inline]
fn temperature_layer(
snd: &Sounding,
warm_side: Celsius,
cold_side: Celsius,
top_pressure: HectoPascal,
) -> Result<Layers> {
let t_profile = snd.temperature_profile();
let p_profile = snd.pressure_profile();
if t_profile.is_empty() || p_profile.is_empty() {
return Err(AnalysisError::MissingProfile);
}
enum Crossing {
Enter(HectoPascal),
Exit(HectoPascal),
Over(HectoPascal, HectoPascal),
In(HectoPascal),
}
fn to_crossing_type(
pnt0: (HectoPascal, Celsius),
pnt1: (HectoPascal, Celsius),
warm_side: Celsius,
cold_side: Celsius,
) -> Option<Crossing> {
let (p0, t0) = pnt0;
let (p1, t1) = pnt1;
debug_assert!(p0 >= p1, "Pressure must decrease with height.");
if t0 < cold_side && t1 >= cold_side && t1 <= warm_side {
let cold_p = linear_interp(cold_side, t0, t1, p0, p1);
Some(Crossing::Enter(cold_p))
} else if t0 > warm_side && t1 <= warm_side && t1 >= cold_side {
let warm_p = linear_interp(warm_side, t0, t1, p0, p1);
Some(Crossing::Enter(warm_p))
} else if (t0 < cold_side && t1 > warm_side) || (t0 > warm_side && t1 < cold_side) {
let warm_p = linear_interp(warm_side, t0, t1, p0, p1);
let cold_p = linear_interp(cold_side, t0, t1, p0, p1);
Some(Crossing::Over(warm_p.max(cold_p), warm_p.min(cold_p)))
} else if t0 > cold_side && t0 < warm_side && t1 > warm_side {
let warm_p = linear_interp(warm_side, t0, t1, p0, p1);
Some(Crossing::Exit(warm_p))
} else if t0 > cold_side && t0 < warm_side && t1 < cold_side {
let cold_p = linear_interp(cold_side, t0, t1, p0, p1);
Some(Crossing::Exit(cold_p))
} else if t0 >= cold_side && t0 <= warm_side && t1 >= cold_side && t1 <= warm_side {
Some(Crossing::In(p0))
} else {
None
}
}
izip!(p_profile, t_profile)
.filter(|(p, t)| p.is_some() && t.is_some())
.map(|(p, t)| (p.unpack(), t.unpack()))
.take_while(move |&(p, _)| p > top_pressure)
.tuple_windows::<(_, _)>()
.filter_map(|(pnt0, pnt1)| to_crossing_type(pnt0, pnt1, warm_side, cold_side))
.scan(None, |bottom_p: &mut Option<_>, crossing_type: Crossing| {
match crossing_type {
Crossing::In(p) => {
if bottom_p.is_none() {
*bottom_p = Some(p);
}
Some(None)
}
Crossing::Enter(p) => {
*bottom_p = Some(p);
Some(None)
}
Crossing::Exit(top_p) => {
debug_assert!(bottom_p.is_some());
Some(Some((bottom_p.take().unwrap(), top_p)))
}
Crossing::Over(p0, p1) => {
debug_assert!(bottom_p.is_none());
Some(Some((p0, p1)))
}
}
})
.filter_map(|opt| opt)
.map(|(bottom_p, top_p)| {
linear_interpolate_sounding(snd, bottom_p).and_then(|bottom| {
linear_interpolate_sounding(snd, top_p).map(|top| Layer { bottom, top })
})
})
.collect()
}
pub fn warm_temperature_layer_aloft(snd: &Sounding) -> Result<Layers> {
warm_layer_aloft(snd, snd.temperature_profile())
}
pub fn warm_wet_bulb_layer_aloft(snd: &Sounding) -> Result<Layers> {
warm_layer_aloft(snd, snd.wet_bulb_profile())
}
#[inline]
fn warm_layer_aloft(snd: &Sounding, t_profile: &[Optioned<Celsius>]) -> Result<Layers> {
let p_profile = snd.pressure_profile();
if t_profile.is_empty() || p_profile.is_empty() {
return Err(AnalysisError::MissingProfile);
}
enum Crossing {
IntoWarmLayer(HectoPascal),
OutOfWarmLayer(HectoPascal),
}
fn crossing_type(
pnt0: (HectoPascal, Celsius),
pnt1: (HectoPascal, Celsius),
) -> Option<Crossing> {
let (p0, t0) = pnt0;
let (p1, t1) = pnt1;
if t0 < FREEZING && t1 >= FREEZING {
let crossing_p = linear_interp(FREEZING, t0, t1, p0, p1);
Some(Crossing::IntoWarmLayer(crossing_p))
} else if t0 >= FREEZING && t1 < FREEZING {
let crossing_p = linear_interp(FREEZING, t0, t1, p0, p1);
Some(Crossing::OutOfWarmLayer(crossing_p))
} else {
None
}
}
izip!(p_profile, t_profile)
.filter(|(p, t)| p.is_some() && t.is_some())
.map(|(p, t)| (p.unpack(), t.unpack()))
.take_while(|&(p, _)| p > HectoPascal(500.0))
.skip_while(|&(_, t)| t > FREEZING)
.tuple_windows::<(_, _)>()
.filter_map(|(pnt0, pnt1)| crossing_type(pnt0, pnt1))
.scan(
None,
|bottom: &mut Option<_>, crossing: Crossing| match crossing {
Crossing::IntoWarmLayer(bottom_p) => {
debug_assert!(bottom.is_none());
*bottom = Some(bottom_p);
Some(None)
}
Crossing::OutOfWarmLayer(top_p) => {
Some(bottom.take().map(|bottom_p| (bottom_p, top_p)))
}
},
)
.filter_map(|opt| opt)
.map(|(bottom_p, top_p)| {
linear_interpolate_sounding(snd, bottom_p).and_then(|bottom| {
linear_interpolate_sounding(snd, top_p).map(|top| Layer { bottom, top })
})
})
.collect()
}
pub fn cold_surface_temperature_layer(snd: &Sounding, warm_layers: &[Layer]) -> Result<Layer> {
cold_surface_layer(snd, snd.temperature_profile(), warm_layers)
}
fn cold_surface_layer(
snd: &Sounding,
t_profile: &[Optioned<Celsius>],
warm_layers: &[Layer],
) -> Result<Layer> {
if warm_layers.is_empty() {
return Err(InvalidInput);
}
let p_profile = snd.pressure_profile();
if t_profile.is_empty() || p_profile.is_empty() {
return Err(MissingProfile);
}
izip!(0usize.., p_profile, t_profile)
.filter(|(_, p, t)| p.is_some() && t.is_some())
.map(|(i, _, t)| (i, t.unpack()))
.take_while(|(_, t)| *t <= FREEZING)
.nth(0)
.and_then(|(i, _)| snd.data_row(i))
.and_then(|bottom| warm_layers.get(0).map(|lyr| (bottom, lyr.bottom)))
.map(|(bottom, top)| Layer { bottom, top })
.ok_or(InvalidInput)
}