///|
/// An inclusive-exclusive index range clamped to an array length.
pub struct SampleRange {
  start : Int
  end : Int
} derive(Debug, Eq)

///|
/// Clamps an index range without ever producing an invalid array access range.
pub fn clamp_range(length : Int, start : Int, end : Int) -> SampleRange {
  let safe_length = if length < 0 { 0 } else { length }
  let safe_start = if start < 0 {
    0
  } else if start > safe_length {
    safe_length
  } else {
    start
  }
  let safe_end = if end < safe_start {
    safe_start
  } else if end > safe_length {
    safe_length
  } else {
    end
  }
  { start: safe_start, end: safe_end }
}

///|
fn finite_copy(values : Array[Double], range : SampleRange) -> Array[Double] {
  let result : Array[Double] = Array::new(capacity=range.end - range.start)
  for i = range.start; i < range.end; i = i + 1 {
    let value = values[i]
    if !value.is_nan() && !value.is_inf() {
      result.push(value)
    }
  }
  result
}

///|
fn sort_in_place(values : Array[Double]) -> Unit {
  for i = 1; i < values.length(); i = i + 1 {
    let value = values[i]
    let mut j = i
    while j > 0 && values[j - 1] > value {
      values[j] = values[j - 1]
      j = j - 1
    }
    values[j] = value
  }
}

///|
/// Returns the q-th quantile of finite values using linear interpolation.
pub fn quantile(values : Array[Double], q : Double) -> Double {
  let finite = finite_copy(
    values,
    clamp_range(values.length(), 0, values.length()),
  )
  if finite.length() == 0 {
    return 0.0
  }
  sort_in_place(finite)
  let position = q.clamp(min=0.0, max=1.0) * (finite.length() - 1).to_double()
  let lower = position.floor().to_int()
  let upper = position.ceil().to_int()
  if lower == upper {
    finite[lower]
  } else {
    finite[lower] +
    (finite[upper] - finite[lower]) * (position - lower.to_double())
  }
}

///|
/// Returns the median of finite values, or zero when no finite value exists.
pub fn median(values : Array[Double]) -> Double {
  quantile(values, 0.5)
}

///|
/// Returns the median inside an inclusive-exclusive range.
pub fn median_in_range(
  values : Array[Double],
  start : Int,
  end : Int,
) -> Double {
  quantile(copy_range(values, start, end), 0.5)
}

///|
/// Returns the median absolute deviation from the median.
pub fn mad(values : Array[Double]) -> Double {
  let center = median(values)
  let deviations : Array[Double] = Array::new(capacity=values.length())
  for value in values {
    if !value.is_nan() && !value.is_inf() {
      deviations.push((value - center).abs())
    }
  }
  median(deviations)
}

///|
/// Counts values that are neither NaN nor positive/negative infinity.
pub fn count_finite(values : Array[Double]) -> Int {
  let mut count = 0
  for value in values {
    if !value.is_nan() && !value.is_inf() {
      count = count + 1
    }
  }
  count
}

///|
/// Computes the finite mean over a clamped range.
pub fn mean_in_range(values : Array[Double], start : Int, end : Int) -> Double {
  let range = clamp_range(values.length(), start, end)
  let mut sum = 0.0
  let mut count = 0
  for i = range.start; i < range.end; i = i + 1 {
    let value = values[i]
    if !value.is_nan() && !value.is_inf() {
      sum = sum + value
      count = count + 1
    }
  }
  if count == 0 {
    0.0
  } else {
    sum / count.to_double()
  }
}

///|
/// Computes sample variance over a clamped range, ignoring non-finite values.
pub fn variance_in_range(
  values : Array[Double],
  start : Int,
  end : Int,
) -> Double {
  let range = clamp_range(values.length(), start, end)
  let mean = mean_in_range(values, range.start, range.end)
  let mut sum = 0.0
  let mut count = 0
  for i = range.start; i < range.end; i = i + 1 {
    let value = values[i]
    if !value.is_nan() && !value.is_inf() {
      let delta = value - mean
      sum = sum + delta * delta
      count = count + 1
    }
  }
  if count <= 1 {
    0.0
  } else {
    sum / (count - 1).to_double()
  }
}

///|
/// Copies a clamped range. The source array is not changed.
pub fn copy_range(
  values : Array[Double],
  start : Int,
  end : Int,
) -> Array[Double] {
  let range = clamp_range(values.length(), start, end)
  let result : Array[Double] = Array::new(capacity=range.end - range.start)
  for i = range.start; i < range.end; i = i + 1 {
    result.push(values[i])
  }
  result
}