///|
/// Noise and integrity measures for one station channel.
pub struct StationQualityReport {
  station : String
  channel : String
  sample_count : Int
  sample_rate : Double
  duration_seconds : Double
  completeness : Double
  clipping_fraction : Double
  flatline_fraction : Double
  gap_count : Int
  longest_gap : Int
  noise_rms : Double
  signal_rms : Double
  signal_to_noise_db : Double
  score : Double
  status : String
}

///|
/// A time interval used in station availability accounting.
pub struct AvailabilityInterval {
  start_time : Double
  end_time : Double
  duration : Double
  complete : Bool
}

///|
/// Noise metrics computed from a baseline segment.
pub struct NoiseMetrics {
  mean : Double
  rms : Double
  standard_deviation : Double
  peak : Double
  zero_crossings_per_second : Double
  dominant_frequency : Double
  spectral_flatness : Double
}

///|
/// Aggregate quality over multiple station channels.
pub struct NetworkQualitySummary {
  channels : Int
  healthy_channels : Int
  mean_completeness : Double
  mean_score : Double
  worst_station : String
  worst_score : Double
}

///|
/// Count exact consecutive repeats, a useful flatline proxy.
pub fn flatline_count(data : Array[Double], tolerance : Double) -> Int {
  if data.length() < 2 {
    0
  } else {
    let limit = tolerance.abs()
    let mut count = 0
    for index = 1; index < data.length(); index = index + 1 {
      if (data[index] - data[index - 1]).abs() <= limit {
        count += 1
      }
    }
    count
  }
}

///|
/// Return the longest consecutive run of near-identical values.
pub fn longest_flatline(data : Array[Double], tolerance : Double) -> Int {
  let limit = tolerance.abs()
  let mut longest = if data.is_empty() { 0 } else { 1 }
  let mut current = if data.is_empty() { 0 } else { 1 }
  for index = 1; index < data.length(); index = index + 1 {
    if (data[index] - data[index - 1]).abs() <= limit {
      current += 1
    } else {
      current = 1
    }
    if current > longest {
      longest = current
    }
  }
  longest
}

///|
/// Estimate clipping fraction at a symmetric full-scale limit.
pub fn clipping_fraction(data : Array[Double], full_scale : Double) -> Double {
  if data.is_empty() || full_scale <= 0.0 {
    0.0
  } else {
    data.filter(value => value.abs() >= full_scale).length().to_double() /
    data.length().to_double()
  }
}

///|
/// Estimate observed completeness relative to an expected sample count.
pub fn completeness_fraction(observed : Int, expected : Int) -> Double {
  if expected <= 0 {
    0.0
  } else {
    observed.to_double() / expected.to_double().min(1.0)
  }
}

///|
/// Compute baseline noise metrics.
pub fn baseline_noise_metrics(
  data : Array[Double],
  sample_rate : Double,
) -> NoiseMetrics {
  let spectrum = hann_periodogram(data, sample_rate)
  let band = SpectralBand::new(0.0, sample_rate.max(0.0) / 2.0)
  {
    mean: mean_samples(data),
    rms: rms_samples(data),
    standard_deviation: standard_deviation(data),
    peak: absolute_peak(data),
    zero_crossings_per_second: if sample_rate <= 0.0 {
      0.0
    } else {
      zero_crossing_count(data).to_double() /
      (data.length().to_double() / sample_rate).max(1.0e-12)
    },
    dominant_frequency: dominant_frequency(spectrum),
    spectral_flatness: advanced_spectral_flatness(spectrum, band),
  }
}

///|
/// Compute a signal-to-noise ratio between two arrays.
pub fn arrays_snr_db(noise : Array[Double], signal : Array[Double]) -> Double {
  let noise_rms = rms_samples(noise)
  let signal_rms = rms_samples(signal)
  if noise_rms <= 0.0 {
    0.0
  } else {
    20.0 * @math.log10(signal_rms / noise_rms)
  }
}

///|
/// Build availability intervals for a fixed sample period and gap mask.
pub fn availability_intervals(
  start_time : Double,
  sample_rate : Double,
  sample_count : Int,
  gaps : Array[GapInterval],
) -> Array[AvailabilityInterval] {
  let result = Array::new()
  let duration = if sample_rate <= 0.0 {
    0.0
  } else {
    sample_count.to_double() / sample_rate
  }
  if sample_count <= 0 || sample_rate <= 0.0 {
    return result
  }
  let mut cursor = start_time
  for gap in gaps {
    let gap_start = gap.start_time.max(start_time)
    let gap_end = gap.end_time.min(start_time + duration)
    if gap_start > cursor {
      result.push({
        start_time: cursor,
        end_time: gap_start,
        duration: gap_start - cursor,
        complete: true,
      })
    }
    if gap_end > gap_start {
      result.push({
        start_time: gap_start,
        end_time: gap_end,
        duration: gap_end - gap_start,
        complete: false,
      })
      cursor = gap_end
    }
  }
  if cursor < start_time + duration {
    result.push({
      start_time: cursor,
      end_time: start_time + duration,
      duration: start_time + duration - cursor,
      complete: true,
    })
  }
  result
}

///|
/// Score a station channel from completeness, noise, and integrity signals.
pub fn station_quality_score(
  completeness : Double,
  clipping : Double,
  flatline : Double,
  signal_to_noise : Double,
) -> Double {
  let completeness_score = completeness.max(0.0).min(1.0)
  let clipping_score = (1.0 - clipping.max(0.0).min(1.0) * 4.0).max(0.0)
  let flatline_score = (1.0 - flatline.max(0.0).min(1.0) * 2.0).max(0.0)
  let snr_score = ((signal_to_noise + 10.0) / 40.0).max(0.0).min(1.0)
  0.4 * completeness_score +
  0.25 * clipping_score +
  0.2 * flatline_score +
  0.15 * snr_score
}

///|
/// Convert a score to an operator-facing health status.
pub fn quality_status(score : Double) -> String {
  if score >= 0.85 {
    "healthy"
  } else if score >= 0.65 {
    "degraded"
  } else if score >= 0.4 {
    "review"
  } else {
    "unhealthy"
  }
}

///|
/// Build a quality report from a trace and a baseline window.
pub fn station_quality_report(
  trace : Trace,
  expected_samples : Int,
  full_scale : Double,
  flatline_tolerance : Double,
  baseline : Array[Double],
) -> StationQualityReport {
  let gaps = find_gaps(trace.data, 0.0, trace.sampling_rate, trace.start_time)
  let mut longest_gap = 0
  for gap in gaps {
    if gap.length > longest_gap {
      longest_gap = gap.length
    }
  }
  let completeness = if expected_samples <= 0 {
    1.0
  } else {
    (trace.npts().to_double() / expected_samples.to_double()).min(1.0)
  }
  let clipped = clipping_fraction(trace.data, full_scale)
  let flatline = if trace.npts() == 0 {
    0.0
  } else {
    flatline_count(trace.data, flatline_tolerance).to_double() /
    trace.npts().to_double()
  }
  let noise = baseline_noise_metrics(baseline, trace.sampling_rate)
  let signal_rms = rms_samples(trace.data)
  let snr = arrays_snr_db(baseline, trace.data)
  let score = station_quality_score(completeness, clipped, flatline, snr)
  {
    station: trace.station,
    channel: trace.channel,
    sample_count: trace.npts(),
    sample_rate: trace.sampling_rate,
    duration_seconds: trace.duration(),
    completeness,
    clipping_fraction: clipped,
    flatline_fraction: flatline,
    gap_count: gaps.length(),
    longest_gap,
    noise_rms: noise.rms,
    signal_rms,
    signal_to_noise_db: snr,
    score,
    status: quality_status(score),
  }
}

///|
/// Compare two channels by normalized waveform agreement.
pub fn waveform_agreement(
  left : Array[Double],
  right : Array[Double],
) -> Double {
  let count = left.length().min(right.length())
  if count == 0 {
    return 0.0
  }
  let left_mean = mean_samples(slice_samples(left, 0, count))
  let right_mean = mean_samples(slice_samples(right, 0, count))
  let mut covariance = 0.0
  let mut left_energy = 0.0
  let mut right_energy = 0.0
  for index = 0; index < count; index = index + 1 {
    let a = left[index] - left_mean
    let b = right[index] - right_mean
    covariance += a * b
    left_energy += a * a
    right_energy += b * b
  }
  let denominator = (left_energy * right_energy).sqrt()
  if denominator <= 0.0 {
    0.0
  } else {
    covariance / denominator
  }
}

///|
/// Detect a possible packet loss by looking for a zero run.
pub fn longest_gap_samples(data : Array[Double], sentinel : Double) -> Int {
  let gaps = find_sentinel_gaps(data, sentinel)
  let mut longest = 0
  for gap in gaps {
    if gap.1 - gap.0 > longest {
      longest = gap.1 - gap.0
    }
  }
  longest
}

///|
/// Aggregate quality reports and identify the weakest channel.
pub fn summarize_network_quality(
  reports : Array[StationQualityReport],
) -> NetworkQualitySummary {
  let mut completeness = 0.0
  let mut score = 0.0
  let mut healthy = 0
  let mut worst_station = ""
  let mut worst_score = 1.0e300
  for report in reports {
    completeness += report.completeness
    score += report.score
    if report.status == "healthy" {
      healthy += 1
    }
    if report.score < worst_score {
      worst_score = report.score
      worst_station = report.station + "." + report.channel
    }
  }
  let divisor = reports.length().to_double().max(1.0)
  {
    channels: reports.length(),
    healthy_channels: healthy,
    mean_completeness: completeness / divisor,
    mean_score: score / divisor,
    worst_station,
    worst_score: if reports.is_empty() {
      0.0
    } else {
      worst_score
    },
  }
}

///|
/// Export quality reports for dashboards.
pub fn station_quality_to_csv(reports : Array[StationQualityReport]) -> String {
  let mut output = "station,channel,samples,rate,duration,completeness,clipping,flatline,gaps,longest_gap,noise_rms,signal_rms,snr_db,score,status\n"
  for report in reports {
    output = output +
      report.station +
      "," +
      report.channel +
      "," +
      report.sample_count.to_string() +
      "," +
      report.sample_rate.to_string() +
      "," +
      report.duration_seconds.to_string() +
      "," +
      report.completeness.to_string() +
      "," +
      report.clipping_fraction.to_string() +
      "," +
      report.flatline_fraction.to_string() +
      "," +
      report.gap_count.to_string() +
      "," +
      report.longest_gap.to_string() +
      "," +
      report.noise_rms.to_string() +
      "," +
      report.signal_rms.to_string() +
      "," +
      report.signal_to_noise_db.to_string() +
      "," +
      report.score.to_string() +
      "," +
      report.status +
      "\n"
  }
  output
}

///|
/// Format a network quality summary for the command line.
pub fn network_quality_summary_text(summary : NetworkQualitySummary) -> String {
  "channels=" +
  summary.channels.to_string() +
  " healthy=" +
  summary.healthy_channels.to_string() +
  " completeness=" +
  summary.mean_completeness.to_string() +
  " score=" +
  summary.mean_score.to_string() +
  " worst=" +
  summary.worst_station +
  " worst_score=" +
  summary.worst_score.to_string()
}