///|
/// Time-domain and spectral features for event classification.
pub struct SignalFeatureVector {
  mean : Double
  variance : Double
  rms : Double
  peak : Double
  crest_factor : Double
  skewness : Double
  kurtosis : Double
  zero_crossing_rate : Double
  line_length : Double
  hjorth_activity : Double
  hjorth_mobility : Double
  hjorth_complexity : Double
  spectral_centroid : Double
  spectral_bandwidth : Double
  spectral_entropy : Double
  low_high_ratio : Double
}

///|
/// Compute the third standardized moment.
pub fn signal_skewness(data : Array[Double]) -> Double {
  let deviation = standard_deviation(data)
  if deviation <= 0.0 {
    0.0
  } else {
    let mean = mean_samples(data)
    mean_samples(data.map(value => @math.pow((value - mean) / deviation, 3.0)))
  }
}

///|
/// Compute excess kurtosis.
pub fn signal_kurtosis(data : Array[Double]) -> Double {
  let deviation = standard_deviation(data)
  if deviation <= 0.0 {
    0.0
  } else {
    let mean = mean_samples(data)
    mean_samples(data.map(value => @math.pow((value - mean) / deviation, 4.0))) -
    3.0
  }
}

///|
/// Return waveform line length as a discrete roughness measure.
pub fn signal_line_length(data : Array[Double]) -> Double {
  if data.length() < 2 {
    0.0
  } else {
    let mut total = 0.0
    for index = 1; index < data.length(); index = index + 1 {
      total += (data[index] - data[index - 1]).abs()
    }
    total
  }
}

///|
/// Calculate Hjorth activity, mobility, and complexity.
pub fn hjorth_parameters(data : Array[Double]) -> Array[Double] {
  let activity = variance_samples(data)
  let first = differentiate_samples(data, 1.0)
  let second = differentiate_samples(first, 1.0)
  let first_variance = variance_samples(first)
  let second_variance = variance_samples(second)
  let mobility = if activity <= 0.0 {
    0.0
  } else {
    (first_variance / activity).sqrt()
  }
  let complexity = if first_variance <= 0.0 {
    0.0
  } else {
    (second_variance / first_variance).sqrt() / mobility.max(1.0e-12)
  }
  [activity, mobility, complexity]
}

///|
/// Compute a spectral band ratio using a periodogram.
pub fn feature_band_ratio(
  data : Array[Double],
  sample_rate : Double,
  low_band : SpectralBand,
  high_band : SpectralBand,
) -> Double {
  let spectrum = hann_periodogram(data, sample_rate)
  let low = advanced_band_power(spectrum, low_band)
  let high = advanced_band_power(spectrum, high_band)
  if high <= 0.0 {
    0.0
  } else {
    low / high
  }
}

///|
/// Construct a full event-classification feature vector.
pub fn signal_feature_vector(
  data : Array[Double],
  sample_rate : Double,
  low_band : SpectralBand,
  high_band : SpectralBand,
) -> SignalFeatureVector {
  let spectrum = hann_periodogram(data, sample_rate)
  let hjorth = hjorth_parameters(data)
  let centroid = advanced_spectral_centroid(
    spectrum,
    SpectralBand::new(0.0, sample_rate / 2.0),
  )
  let bandwidth = advanced_spectral_bandwidth(
    spectrum,
    SpectralBand::new(0.0, sample_rate / 2.0),
  )
  let peak = absolute_peak(data)
  let rms = rms_samples(data)
  {
    mean: mean_samples(data),
    variance: variance_samples(data),
    rms,
    peak,
    crest_factor: if rms <= 0.0 {
      0.0
    } else {
      peak / rms
    },
    skewness: signal_skewness(data),
    kurtosis: signal_kurtosis(data),
    zero_crossing_rate: if data.length() < 2 || sample_rate <= 0.0 {
      0.0
    } else {
      zero_crossing_count(data).to_double() /
      (data.length().to_double() / sample_rate)
    },
    line_length: signal_line_length(data),
    hjorth_activity: hjorth[0],
    hjorth_mobility: hjorth[1],
    hjorth_complexity: hjorth[2],
    spectral_centroid: centroid,
    spectral_bandwidth: bandwidth,
    spectral_entropy: advanced_spectral_entropy(
      spectrum,
      SpectralBand::new(0.0, sample_rate / 2.0),
    ),
    low_high_ratio: feature_band_ratio(data, sample_rate, low_band, high_band),
  }
}

///|
/// Calculate a smoothness score from line length and peak amplitude.
pub fn waveform_smoothness(data : Array[Double]) -> Double {
  let peak = absolute_peak(data)
  let length = signal_line_length(data)
  if peak <= 0.0 || data.length() < 2 {
    0.0
  } else {
    (1.0 - length / ((data.length() - 1).to_double() * peak)).max(0.0).min(1.0)
  }
}

///|
/// Return a feature vector as a stable CSV row.
pub fn signal_features_to_csv(features : SignalFeatureVector) -> String {
  "mean,variance,rms,peak,crest,skewness,kurtosis,zcr,line_length,hjorth_activity,hjorth_mobility,hjorth_complexity,spectral_centroid,spectral_bandwidth,spectral_entropy,low_high_ratio\n" +
  features.mean.to_string() +
  "," +
  features.variance.to_string() +
  "," +
  features.rms.to_string() +
  "," +
  features.peak.to_string() +
  "," +
  features.crest_factor.to_string() +
  "," +
  features.skewness.to_string() +
  "," +
  features.kurtosis.to_string() +
  "," +
  features.zero_crossing_rate.to_string() +
  "," +
  features.line_length.to_string() +
  "," +
  features.hjorth_activity.to_string() +
  "," +
  features.hjorth_mobility.to_string() +
  "," +
  features.hjorth_complexity.to_string() +
  "," +
  features.spectral_centroid.to_string() +
  "," +
  features.spectral_bandwidth.to_string() +
  "," +
  features.spectral_entropy.to_string() +
  "," +
  features.low_high_ratio.to_string() +
  "\n"
}

///|
/// Format a compact feature vector summary.
pub fn signal_features_text(features : SignalFeatureVector) -> String {
  "rms=" +
  features.rms.to_string() +
  " peak=" +
  features.peak.to_string() +
  " crest=" +
  features.crest_factor.to_string() +
  " zcr=" +
  features.zero_crossing_rate.to_string() +
  " centroid=" +
  features.spectral_centroid.to_string() +
  " bandwidth=" +
  features.spectral_bandwidth.to_string() +
  " entropy=" +
  features.spectral_entropy.to_string() +
  " low_high=" +
  features.low_high_ratio.to_string()
}