///|
pub(all) enum ProfileKind {
Constant
Gaussian
Sinusoid
Ramp
DoubleHump
TopHat
} derive(Eq, Debug, ToJson)
///|
pub(all) struct Profile1D {
kind : ProfileKind
base : Double
amplitude : Double
center : Double
width : Double
wave_number : Double
floor : Double
} derive(Debug, ToJson)
///|
pub fn profile_constant(value : Double) -> Profile1D {
{
kind: Constant,
base: value,
amplitude: 0.0,
center: 0.0,
width: 1.0,
wave_number: 0.0,
floor: value,
}
}
///|
pub fn profile_gaussian(
base : Double,
amplitude : Double,
center : Double,
width : Double,
) -> Profile1D {
{
kind: Gaussian,
base,
amplitude,
center,
width,
wave_number: 0.0,
floor: base,
}
}
///|
pub fn profile_sinusoid(
base : Double,
amplitude : Double,
wave_number : Double,
phase : Double,
) -> Profile1D {
{
kind: Sinusoid,
base,
amplitude,
center: phase,
width: 1.0,
wave_number,
floor: base - amplitude.abs(),
}
}
///|
pub fn profile_ramp(
start : Double,
end : Double,
center : Double,
width : Double,
) -> Profile1D {
{
kind: Ramp,
base: start,
amplitude: end - start,
center,
width,
wave_number: 0.0,
floor: start.min(end),
}
}
///|
pub fn profile_double_hump(
base : Double,
amplitude : Double,
center : Double,
width : Double,
) -> Profile1D {
{
kind: DoubleHump,
base,
amplitude,
center,
width,
wave_number: 0.0,
floor: base,
}
}
///|
pub fn profile_top_hat(
low : Double,
high : Double,
center : Double,
width : Double,
) -> Profile1D {
{
kind: TopHat,
base: high,
amplitude: high - low,
center,
width,
wave_number: 0.0,
floor: low,
}
}
///|
pub fn profile_value(profile : Profile1D, x : Double) -> Double {
match profile.kind {
Constant => profile.base
Gaussian => {
let scale = if profile.width == 0.0 {
0.0
} else {
(x - profile.center) / profile.width
}
profile.base + profile.amplitude * @math.exp(-0.5 * scale * scale)
}
Sinusoid =>
profile.base +
profile.amplitude * @math.cos(profile.wave_number * x + profile.center)
Ramp => {
let t = if profile.width == 0.0 {
0.0
} else {
clamp((x - profile.center) / profile.width + 0.5, 0.0, 1.0)
}
profile.base + profile.amplitude * t
}
DoubleHump => {
let scale = if profile.width == 0.0 {
0.0
} else {
(x - profile.center) / profile.width
}
profile.base +
profile.amplitude * scale * scale * @math.exp(-0.5 * scale * scale)
}
TopHat =>
if (x - profile.center).abs() <= profile.width * 0.5 {
profile.base
} else {
profile.floor
}
}
}
///|
pub fn profile_values(grid : Grid1D, profile : Profile1D) -> Array[Double] {
Array::makei(grid.cells, fn(i) { profile_value(profile, grid.position(i)) })
}
///|
pub fn profile_mean(grid : Grid1D, profile : Profile1D) -> Double {
mean(profile_values(grid, profile))
}
///|
pub fn profile_integral(grid : Grid1D, profile : Profile1D) -> Double {
integrate_trapezoid(grid, profile_values(grid, profile))
}
///|
pub fn normalize_profile(
grid : Grid1D,
values : ArrayView[Double],
target_integral : Double,
) -> Array[Double] {
let current = integrate_trapezoid(grid, values)
if current == 0.0 {
values.to_owned()
} else {
let scale = target_integral / current
values.map(fn(value) { value * scale })
}
}
///|
pub fn add_profile(
a : ArrayView[Double],
b : ArrayView[Double],
) -> Array[Double] {
let output = zeros(a.length())
for i in 0.. Array[Double] {
let output = zeros(a.length())
for i in 0.. Array[Double] {
let output = zeros(grid.cells)
for i in 0.. Double {
min_value(values, 0.0)
}
///|
pub fn profile_maximum(values : ArrayView[Double]) -> Double {
max_value(values, 0.0)
}
///|
pub fn profile_range(values : ArrayView[Double]) -> Double {
profile_maximum(values) - profile_minimum(values)
}
///|
pub fn profile_floor(
values : ArrayView[Double],
floor : Double,
) -> Array[Double] {
values.map(fn(value) { if value < floor { floor } else { value } })
}
///|
pub fn profile_ceiling(
values : ArrayView[Double],
ceiling : Double,
) -> Array[Double] {
values.map(fn(value) { if value > ceiling { ceiling } else { value } })
}