///|
/// Deterministic pseudo-random state for repeatable waveform simulations.
pub struct NoiseGenerator {
  state : Int
  standard_deviation : Double
  mean : Double
}

///|
/// A configurable colored-noise model.
pub enum NoiseColor {
  WhiteNoise
  PinkNoise
  BrownNoise
  VioletNoise
}

///|
/// Return all supported noise colors for simulation configuration.
pub fn supported_noise_colors() -> Array[NoiseColor] {
  [WhiteNoise, PinkNoise, BrownNoise, VioletNoise]
}

///|
/// Parameters for a synthetic event envelope.
pub struct EventEnvelopeConfig {
  onset : Double
  rise_time : Double
  decay_time : Double
  amplitude : Double
  coda_exponent : Double
}

///|
/// Create a deterministic noise generator.
pub fn NoiseGenerator::new(
  seed : Int,
  mean : Double,
  standard_deviation : Double,
) -> NoiseGenerator {
  {
    state: if seed == 0 {
      1
    } else {
      seed.abs()
    },
    mean,
    standard_deviation: standard_deviation.abs(),
  }
}

///|
/// Advance the generator and return a uniform number in `[0, 1)`.
pub fn NoiseGenerator::uniform(
  self : NoiseGenerator,
) -> (NoiseGenerator, Double) {
  let next_state = (self.state * 1664525 + 1013904223) % 2147483647
  let normalized_state = if next_state < 0 {
    next_state + 2147483647
  } else {
    next_state
  }
  (
    {
      state: normalized_state,
      standard_deviation: self.standard_deviation,
      mean: self.mean,
    },
    normalized_state.to_double() / 2147483647.0,
  )
}

///|
/// Return a zero-mean approximately Gaussian sample using Box-Muller.
pub fn NoiseGenerator::normal(
  self : NoiseGenerator,
) -> (NoiseGenerator, Double) {
  let (after_first, first_value) = self.uniform()
  let (after_second, second_value) = after_first.uniform()
  let first = first_value.max(1.0e-12)
  let second = second_value
  let radius = (-2.0 * @math.ln(first)).sqrt()
  let angle = 2.0 * @math.PI * second
  (
    after_second,
    self.mean + self.standard_deviation * radius * @math.cos(angle),
  )
}

///|
/// Generate deterministic Gaussian noise.
pub fn gaussian_noise(
  length : Int,
  seed : Int,
  mean : Double,
  standard_deviation : Double,
) -> Array[Double] {
  let result = Array::new(capacity=length.max(0))
  let mut generator = NoiseGenerator::new(seed, mean, standard_deviation)
  for _ in 0.. Array[Double] {
  let result = Array::new(capacity=length.max(0))
  let mut generator = NoiseGenerator::new(seed, 0.0, 1.0)
  for _ in 0.. Array[Double] {
  let result = Array::make(length.max(0), 0.0)
  if sample_rate <= 0.0 {
    return result
  }
  for index = 0; index < result.length(); index = index + 1 {
    result[index] = amplitude *
      @math.sin(
        2.0 * @math.PI * frequency * index.to_double() / sample_rate + phase,
      )
  }
  result
}

///|
/// Generate a linear chirp with a frequency sweep.
pub fn chirp_signal(
  length : Int,
  sample_rate : Double,
  start_frequency : Double,
  end_frequency : Double,
  amplitude : Double,
) -> Array[Double] {
  let result = Array::make(length.max(0), 0.0)
  if sample_rate <= 0.0 || length <= 0 {
    return result
  }
  let duration = length.to_double() / sample_rate
  let slope = (end_frequency - start_frequency) / duration.max(1.0e-12)
  for index = 0; index < result.length(); index = index + 1 {
    let time = index.to_double() / sample_rate
    let phase = 2.0 *
      @math.PI *
      (start_frequency * time + 0.5 * slope * time * time)
    result[index] = amplitude * @math.sin(phase)
  }
  result
}

///|
/// Generate a Ricker wavelet centered on a sample index.
pub fn ricker_wavelet_noise_model(
  length : Int,
  sample_rate : Double,
  central_frequency : Double,
  center_time : Double,
  amplitude : Double,
) -> Array[Double] {
  let result = Array::make(length.max(0), 0.0)
  if sample_rate <= 0.0 || central_frequency <= 0.0 {
    return result
  }
  for index = 0; index < result.length(); index = index + 1 {
    let time = index.to_double() / sample_rate - center_time
    let x = @math.PI * central_frequency * time
    let square = x * x
    result[index] = amplitude * (1.0 - 2.0 * square) * @math.exp(-square)
  }
  result
}

///|
/// Generate an event envelope with a smooth rise and power-law coda.
pub fn event_envelope(
  length : Int,
  sample_rate : Double,
  config : EventEnvelopeConfig,
) -> Array[Double] {
  let result = Array::make(length.max(0), 0.0)
  if sample_rate <= 0.0 {
    return result
  }
  for index = 0; index < result.length(); index = index + 1 {
    let time = index.to_double() / sample_rate
    if time < config.onset {
      result[index] = 0.0
    } else {
      let elapsed = time - config.onset
      let rise = if config.rise_time <= 0.0 {
        1.0
      } else {
        1.0 - @math.exp(-elapsed / config.rise_time)
      }
      let decay = if config.decay_time <= 0.0 {
        1.0
      } else {
        @math.exp(-elapsed / config.decay_time)
      }
      let coda = @math.pow(1.0 + elapsed, config.coda_exponent)
      result[index] = config.amplitude * rise * decay / coda.max(1.0e-12)
    }
  }
  result
}

///|
/// Add two aligned signals with a configurable noise-free fallback.
pub fn add_noise_scaled(
  signal : Array[Double],
  noise : Array[Double],
  noise_scale : Double,
) -> Array[Double] {
  let result = signal.copy()
  let count = result.length().min(noise.length())
  for index = 0; index < count; index = index + 1 {
    result[index] += noise[index] * noise_scale
  }
  result
}

///|
/// Generate pink noise with a one-pole multi-band approximation.
pub fn colored_noise(
  length : Int,
  seed : Int,
  standard_deviation : Double,
  color : NoiseColor,
) -> Array[Double] {
  let white = gaussian_noise(length, seed, 0.0, standard_deviation)
  let result = Array::make(white.length(), 0.0)
  let mut previous = 0.0
  for index = 0; index < white.length(); index = index + 1 {
    match color {
      WhiteNoise => result[index] = white[index]
      PinkNoise => {
        previous = 0.98 * previous + 0.02 * white[index]
        result[index] = previous + 0.5 * white[index]
      }
      BrownNoise => {
        previous = (previous + white[index] * 0.1)
          .max(-5.0 * standard_deviation)
          .min(5.0 * standard_deviation)
        result[index] = previous
      }
      VioletNoise =>
        result[index] = if index == 0 {
          white[index]
        } else {
          white[index] - white[index - 1]
        }
    }
  }
  let scale = standard_deviation / rms_samples(result).max(1.0e-12)
  result.map(value => value * scale)
}

///|
/// Build a synthetic event by combining a wavelet, envelope, and noise.
pub fn synthetic_event_waveform(
  length : Int,
  sample_rate : Double,
  wavelet_frequency : Double,
  onset : Double,
  amplitude : Double,
  noise_level : Double,
  seed : Int,
) -> Array[Double] {
  let wavelet = ricker_wavelet_noise_model(
    length,
    sample_rate,
    wavelet_frequency,
    onset + 2.0 / wavelet_frequency.max(1.0e-12),
    amplitude,
  )
  let envelope = event_envelope(length, sample_rate, {
    onset,
    rise_time: 0.25,
    decay_time: 8.0,
    amplitude: 1.0,
    coda_exponent: 0.35,
  })
  let noise = colored_noise(length, seed, noise_level, PinkNoise)
  let result = Array::make(length.max(0), 0.0)
  for index = 0; index < result.length(); index = index + 1 {
    result[index] = wavelet[index] * envelope[index] + noise[index]
  }
  result
}

///|
/// Estimate noise level from a quiet prefix and scale a synthetic signal.
pub fn signal_with_target_snr(
  signal : Array[Double],
  noise : Array[Double],
  target_snr_db : Double,
) -> Array[Double] {
  let signal_rms = rms_samples(signal)
  let noise_rms = rms_samples(noise)
  if noise_rms <= 0.0 {
    return signal.copy()
  }
  let desired_noise = signal_rms / @math.pow(10.0, target_snr_db / 20.0)
  add_noise_scaled(signal, noise, desired_noise / noise_rms)
}

///|
/// Generate a deterministic impulse train for calibration tests.
pub fn impulse_train(
  length : Int,
  positions : Array[Int],
  amplitudes : Array[Double],
) -> Array[Double] {
  let result = Array::make(length.max(0), 0.0)
  for index = 0
      index < positions.length() && index < amplitudes.length()
      index = index + 1 {
    let position = positions[index]
    if position >= 0 && position < result.length() {
      result[position] = amplitudes[index]
    }
  }
  result
}

///|
/// Return a synthetic noise model name for reports.
pub fn noise_color_name(color : NoiseColor) -> String {
  match color {
    WhiteNoise => "white"
    PinkNoise => "pink"
    BrownNoise => "brown"
    VioletNoise => "violet"
  }
}

///|
/// Format an event envelope configuration.
pub fn event_envelope_text(config : EventEnvelopeConfig) -> String {
  "onset=" +
  config.onset.to_string() +
  " rise=" +
  config.rise_time.to_string() +
  " decay=" +
  config.decay_time.to_string() +
  " amplitude=" +
  config.amplitude.to_string() +
  " coda_exponent=" +
  config.coda_exponent.to_string()
}