timestretch 0.7.0

//! Phase vocoder time stretching with identity phase locking and sub-bass phase locking.

use crate::core::fft::{COMPLEX_ZERO, WINDOW_SUM_EPSILON, WINDOW_SUM_FLOOR_RATIO};
use crate::core::window::{generate_window, WindowType};
use crate::error::StretchError;
use crate::stretch::envelope::{
    adaptive_cepstral_order, apply_envelope_correction_with_scratch,
    extract_envelope_with_fft_scratch, spectral_centroid,
};
use crate::stretch::phase_locking::{apply_phase_locking_realtime, PhaseLockingMode};
use rustfft::{num_complex::Complex, FftPlanner};
use std::sync::Arc;

const TWO_PI_F64: f64 = 2.0 * std::f64::consts::PI;
/// Fraction of bins to pre-allocate for spectral peak detection (1/4 of bins).
const PEAKS_CAPACITY_DIVISOR: usize = 4;
/// Blend factor for phase gradient integration (soft vertical coherence).
const PHASE_GRADIENT_BLEND: f64 = 0.20;
/// Minimum magnitude to consider a bin as a spectral peak (avoids noise peaks).
const MIN_PEAK_MAGNITUDE: f32 = 1e-8;
/// Treat values this close to integers as integral synthesis positions.
const SYNTH_POS_EPSILON: f64 = 1e-9;
/// Floor used in adaptive locking feature extraction.
const ADAPTIVE_FEATURE_EPS: f64 = 1e-12;
/// Flatness above this is considered noise-like (prefer ROI).
const ADAPTIVE_NOISY_FLATNESS: f64 = 0.72;
/// Crest below this is considered weakly-structured/noisy (prefer ROI).
const ADAPTIVE_NOISY_CREST: f64 = 2.5;
/// Harmonic confidence above this prefers identity locking near unity.
const ADAPTIVE_HARMONIC_CONFIDENCE: f64 = 4.5;
/// Harmonic confidence above this prefers selective locking on moderate ratios.
const ADAPTIVE_SELECTIVE_HARMONIC_CONFIDENCE: f64 = 2.8;
/// Max ratio-distance from unity where identity is preferred on harmonic frames.
const ADAPTIVE_IDENTITY_RATIO_DISTANCE_MAX: f64 = 0.35;
/// Max ratio-distance from unity where selective locking is preferred.
const ADAPTIVE_SELECTIVE_RATIO_DISTANCE_MAX: f64 = 0.65;
/// Ratio-distance above this always prefers ROI for stability.
const ADAPTIVE_FORCE_ROI_RATIO_DISTANCE: f64 = 0.75;
/// Number of synthesis frames to keep transient-focused locking active.
const TRANSIENT_FOCUS_FRAMES: usize = 3;
/// Briefly tighten phase locking after a meaningful runtime ratio step.
const RATIO_CHANGE_FOCUS_FRAMES: usize = 3;
/// Keep continuity focus active slightly longer when automation reverses direction mid-seam.
const RATIO_CHANGE_REVERSAL_FOCUS_EXTRA_FRAMES: usize = 1;
/// Keep continuity focus active slightly longer when a new step must re-anchor
/// to an older carried seam that is still farther from unity than the in-flight ratio.
const RATIO_CHANGE_CARRIED_SEAM_FOCUS_EXTRA_FRAMES: usize = 1;
/// Hold a far carried seam for one synthesis frame before slewing toward the new target.
const RATIO_CHANGE_CARRIED_SEAM_HOLD_FRAMES: usize = 1;
/// Ignore tiny ratio deltas that are below the modulation-seam risk zone.
const RATIO_CHANGE_FOCUS_TRIGGER: f64 = 1e-3;
/// Cap extra continuity-focus extension so long tails do not pin identity locking.
const RATIO_CHANGE_TAIL_FOCUS_MAX_FRAMES: usize = 6;
/// Aggregated per-frame spectral flux from the PV's own analysis pass.
///
/// Computed as half-wave rectified (onset-only) magnitude differences,
/// summed per frequency band. Available after each `process_streaming_into`
/// call for downstream transient-adaptive processing (e.g. blend modulation).
#[derive(Debug, Clone, Copy, Default)]
pub struct PerFrameFlux {
    /// Flux in the sub-bass region (below `sub_bass_bin`, typically <180Hz).
    pub sub_bass: f32,
    /// Flux in the low band (sub_bass_bin to ~500Hz).
    pub low: f32,
    /// Flux in the mid band (~500Hz to ~4000Hz).
    pub mid: f32,
    /// Flux in the high band (>4000Hz).
    pub high: f32,
    /// Number of bins with flux exceeding a significance threshold.
    pub transient_bin_count: u16,
    /// Total number of bins with any positive flux (onset).
    pub total_bins_rising: u16,
}

/// Phase vocoder state for time stretching.
pub struct PhaseVocoder {
    fft_size: usize,
    sample_rate: u32,
    hop_analysis: usize,
    hop_synthesis: usize,
    stretch_ratio: f64,
    /// Absolute synthesis position (in samples) of the next frame start.
    synthesis_pos: f64,
    /// Number of synthesized samples already emitted by the streaming path.
    synthesis_emitted: usize,
    window: Vec<f32>,
    /// Phase accumulator for resynthesis (f64 for precision over long signals).
    phase_accum: Vec<f64>,
    /// Previous analysis phase (f64 to match accumulator precision).
    prev_phase: Vec<f64>,
    /// Marks bins that must be seeded from fresh analysis phase on next frame.
    ///
    /// Used by selective band resets so transient-triggered partial resets do
    /// not compute an invalid instantaneous-frequency jump from zeroed phase
    /// history.
    phase_seed_pending: Vec<bool>,
    /// Pre-planned forward FFT for analysis frames.
    fft_forward: Arc<dyn rustfft::Fft<f32>>,
    /// Pre-planned inverse FFT for synthesis frames.
    fft_inverse: Arc<dyn rustfft::Fft<f32>>,
    /// Scratch buffer for forward FFT execution.
    fft_forward_scratch: Vec<Complex<f32>>,
    /// Scratch buffer for inverse FFT execution.
    fft_inverse_scratch: Vec<Complex<f32>>,
    /// Pre-computed expected phase advance per bin (f64 for precision).
    expected_phase_advance: Vec<f64>,
    /// Reusable FFT buffer.
    fft_buffer: Vec<Complex<f32>>,
    /// Reusable magnitude buffer.
    magnitudes: Vec<f32>,
    /// Reusable phase buffer.
    new_phases: Vec<f32>,
    /// Reusable peaks buffer for identity phase locking.
    peaks: Vec<usize>,
    /// Reusable trough buffer for phase-lock influence-region discovery.
    phase_lock_troughs: Vec<usize>,
    /// Reusable backup of pre-lock phases for ROI clamping.
    phase_lock_pv_phases: Vec<f32>,
    /// Current frame's analysis phases (for identity phase locking).
    analysis_phases: Vec<f32>,
    /// Bin index at or below which sub-bass phase locking is applied.
    sub_bass_bin: usize,
    /// Phase locking algorithm to use.
    phase_locking_mode: PhaseLockingMode,
    /// Enables confidence-driven adaptive phase-lock mode switching.
    adaptive_phase_locking: bool,
    /// Short-lived state that tightens phase behavior right after a transient reset.
    transient_focus_frames: usize,
    /// Previous effective ratio used to slew phase advance through a seam.
    ratio_change_phase_from: f64,
    /// Remaining frames in the phase-ratio continuity slew.
    ratio_change_phase_frames: usize,
    /// Total frames scheduled for the current phase-ratio continuity slew.
    ratio_change_phase_total_frames: usize,
    /// Initial seam-hold frames that should stay pinned to `ratio_change_phase_from`.
    ratio_change_phase_hold_frames: usize,
    /// Ratio updates arrive as a smooth small-step stream; skip the
    /// continuity-focus seam machinery (see `set_smooth_ratio_updates`).
    smooth_ratio_updates: bool,
    /// Whether spectral envelope preservation is enabled.
    envelope_preservation: bool,
    /// Envelope correction strength (0 = off, 1 = full, >1 stronger).
    envelope_strength: f32,
    /// Enables adaptive per-frame envelope order selection.
    adaptive_envelope_order: bool,
    /// Cepstral order for envelope extraction.
    envelope_order: usize,
    /// Reusable buffer for cepstral analysis.
    cepstrum_buf: Vec<Complex<f32>>,
    /// Reusable buffer for analysis envelope.
    analysis_envelope: Vec<f32>,
    /// Reusable buffer for synthesis envelope.
    synthesis_envelope: Vec<f32>,
    /// Reusable noise-floor scratch for SNR-aware envelope correction.
    envelope_noise_floor_scratch: Vec<f32>,
    /// Scratch buffer for envelope IFFT execution.
    envelope_ifft_scratch: Vec<Complex<f32>>,
    /// Scratch buffer for envelope FFT execution.
    envelope_fft_scratch: Vec<Complex<f32>>,
    /// Precomputed overlap-add gain (`synthesis_window / fft_size`).
    ola_gain: Vec<f32>,
    /// f64 view of `ola_gain` to avoid per-sample casts in fractional OLA path.
    ola_gain_f64: Vec<f64>,
    /// Precomputed window product (`analysis_window * synthesis_window`) for OLA normalization.
    ola_window_product: Vec<f32>,
    /// f64 view of `ola_window_product` to avoid per-sample casts in fractional OLA path.
    ola_window_product_f64: Vec<f64>,
    /// Backup of IF-estimated phases before phase locking overwrites them.
    /// Used to blend IF estimates with locked phases for non-peak bins.
    if_phases_backup: Vec<f32>,
    /// Reusable output buffer (avoids allocation per process() call).
    output_buf: Vec<f32>,
    /// Reusable window sum buffer (avoids allocation per process() call).
    window_sum_buf: Vec<f32>,
    /// Unnormalized overlap-add tail carried between streaming calls.
    streaming_tail: Vec<f32>,
    /// Window-sum tail matching `streaming_tail`.
    streaming_tail_window_sum: Vec<f32>,
    /// Most conservative ratio contributing to the carried streaming tail.
    ///
    /// When a tail spans a ratio change, keep the larger expansion ratio so
    /// the next normalization pass does not suddenly tighten the overlap
    /// floor at the chunk boundary.
    streaming_tail_ratio: f64,
    /// Ratio that generated the unresolved carried seam at the chunk boundary.
    streaming_tail_phase_ratio: f64,
    /// Reusable accumulation buffer for streaming overlap-add.
    streaming_accum_output: Vec<f32>,
    /// Reusable window-sum accumulation buffer for streaming overlap-add.
    streaming_accum_window_sum: Vec<f32>,
    /// Previous frame magnitudes for per-bin spectral flux computation.
    flux_prev_magnitudes: Vec<f32>,
    /// Whether `flux_prev_magnitudes` has been seeded with at least one frame.
    flux_has_prev: bool,
    /// Bin index for the 500 Hz band boundary (low/mid split).
    flux_low_end_bin: usize,
    /// Bin index for the 4000 Hz band boundary (mid/high split).
    flux_mid_end_bin: usize,
    /// Per-frame flux from the most recent `process_core` call (last frame).
    last_frame_flux: Option<PerFrameFlux>,
}

#[inline]
fn streaming_tail_normalize_ratio(carried_tail_ratio: f64, current_ratio: f64) -> f64 {
    carried_tail_ratio.max(current_ratio)
}

#[inline]
fn ratio_is_meaningfully_above_unity(ratio: f64) -> bool {
    ratio > 1.0 + RATIO_CHANGE_FOCUS_TRIGGER
}

#[inline]
fn ratio_is_meaningfully_below_unity(ratio: f64) -> bool {
    ratio < 1.0 - RATIO_CHANGE_FOCUS_TRIGGER
}

#[inline]
fn continuity_focus_frames_for_ratio_change(
    base_frames: usize,
    tail_samples: usize,
    hop: usize,
) -> usize {
    if tail_samples == 0 || hop == 0 {
        return base_frames;
    }

    base_frames.max(
        tail_samples
            .div_ceil(hop)
            .saturating_add(1)
            .min(RATIO_CHANGE_TAIL_FOCUS_MAX_FRAMES),
    )
}

#[inline]
fn ratio_change_reverses_inflight_direction(
    continuity_phase_from: f64,
    prior_ratio: f64,
    next_ratio: f64,
) -> bool {
    let prior_direction = prior_ratio - continuity_phase_from;
    let next_direction = next_ratio - prior_ratio;

    prior_direction.abs() >= RATIO_CHANGE_FOCUS_TRIGGER
        && next_direction.abs() >= RATIO_CHANGE_FOCUS_TRIGGER
        && prior_direction.signum() != next_direction.signum()
}

impl PhaseVocoder {
    /// Creates a new phase vocoder.
    pub fn new(
        fft_size: usize,
        hop_analysis: usize,
        stretch_ratio: f64,
        sample_rate: u32,
        sub_bass_cutoff: f32,
    ) -> Self {
        Self::with_window(
            fft_size,
            hop_analysis,
            stretch_ratio,
            sample_rate,
            sub_bass_cutoff,
            WindowType::BlackmanHarris,
        )
    }

    /// Creates a new phase vocoder with a specific window function.
    pub fn with_window(
        fft_size: usize,
        hop_analysis: usize,
        stretch_ratio: f64,
        sample_rate: u32,
        sub_bass_cutoff: f32,
        window_type: WindowType,
    ) -> Self {
        Self::with_options(
            fft_size,
            hop_analysis,
            stretch_ratio,
            sample_rate,
            sub_bass_cutoff,
            window_type,
            PhaseLockingMode::RegionOfInfluence,
        )
    }

    /// Creates a new phase vocoder with full configuration options.
    pub fn with_options(
        fft_size: usize,
        hop_analysis: usize,
        stretch_ratio: f64,
        sample_rate: u32,
        sub_bass_cutoff: f32,
        window_type: WindowType,
        phase_locking_mode: PhaseLockingMode,
    ) -> Self {
        Self::with_all_options(
            fft_size,
            hop_analysis,
            stretch_ratio,
            sample_rate,
            sub_bass_cutoff,
            window_type,
            phase_locking_mode,
            false,
            40,
        )
    }

    /// Creates a new phase vocoder with all configuration options including envelope preservation.
    #[allow(clippy::too_many_arguments)]
    pub fn with_all_options(
        fft_size: usize,
        hop_analysis: usize,
        stretch_ratio: f64,
        sample_rate: u32,
        sub_bass_cutoff: f32,
        window_type: WindowType,
        phase_locking_mode: PhaseLockingMode,
        envelope_preservation: bool,
        envelope_order: usize,
    ) -> Self {
        let hop_synthesis = (hop_analysis as f64 * stretch_ratio).round() as usize;
        let window = generate_window(window_type, fft_size);
        // Match synthesis window to analysis window type for a proper window product.
        // Using the same window type ensures the overlap-add normalization works
        // correctly and avoids spectral distortion from mismatched window shapes.
        //
        // Exception: BlackmanHarris analysis uses Hann for synthesis because BH^2
        // has poor COLA (constant overlap-add) properties at standard 75% overlap
        // (hop = fft_size/4). The BH*Hann product provides better overlap-add
        // flatness while BH still provides excellent sidelobe suppression for
        // the analysis stage.
        let synthesis_window_type = match window_type {
            WindowType::BlackmanHarris => WindowType::Hann,
            other => other,
        };
        let synthesis_window = generate_window(synthesis_window_type, fft_size);
        let inv_fft = 1.0 / fft_size as f32;
        let ola_gain: Vec<f32> = synthesis_window.iter().map(|&w| w * inv_fft).collect();
        let ola_gain_f64: Vec<f64> = ola_gain.iter().map(|&w| w as f64).collect();
        let ola_window_product: Vec<f32> = window
            .iter()
            .zip(synthesis_window.iter())
            .map(|(&a, &b)| a * b)
            .collect();
        let ola_window_product_f64: Vec<f64> =
            ola_window_product.iter().map(|&w| w as f64).collect();
        let num_bins = fft_size / 2 + 1;
        let mut planner = FftPlanner::new();
        let fft_forward = planner.plan_fft_forward(fft_size);
        let fft_inverse = planner.plan_fft_inverse(fft_size);
        let fft_forward_scratch_len = fft_forward.get_inplace_scratch_len();
        let fft_inverse_scratch_len = fft_inverse.get_inplace_scratch_len();

        let expected_phase_advance: Vec<f64> = (0..num_bins)
            .map(|bin| TWO_PI_F64 * bin as f64 * hop_analysis as f64 / fft_size as f64)
            .collect();

        // Compute the bin index for the sub-bass cutoff frequency.
        // Bins at or below this index get rigid phase locking to prevent
        // phase cancellation in the critical sub-bass region.
        let sub_bass_bin =
            (sub_bass_cutoff * fft_size as f32 / sample_rate as f32).round() as usize;
        let sub_bass_bin = sub_bass_bin.min(num_bins);

        Self {
            fft_size,
            sample_rate,
            hop_analysis,
            hop_synthesis,
            stretch_ratio,
            synthesis_pos: 0.0,
            synthesis_emitted: 0,
            window,
            phase_accum: vec![0.0f64; num_bins],
            prev_phase: vec![0.0f64; num_bins],
            phase_seed_pending: vec![true; num_bins],
            fft_forward,
            fft_inverse,
            fft_forward_scratch: vec![COMPLEX_ZERO; fft_forward_scratch_len],
            fft_inverse_scratch: vec![COMPLEX_ZERO; fft_inverse_scratch_len],
            expected_phase_advance,
            fft_buffer: vec![COMPLEX_ZERO; fft_size],
            magnitudes: vec![0.0; num_bins],
            new_phases: vec![0.0; num_bins],
            peaks: Vec::with_capacity(num_bins / PEAKS_CAPACITY_DIVISOR),
            phase_lock_troughs: Vec::with_capacity(num_bins / 2),
            phase_lock_pv_phases: Vec::with_capacity(num_bins),
            analysis_phases: vec![0.0; num_bins],
            sub_bass_bin,
            phase_locking_mode,
            adaptive_phase_locking: false,
            transient_focus_frames: 0,
            ratio_change_phase_from: stretch_ratio,
            ratio_change_phase_frames: 0,
            ratio_change_phase_total_frames: 0,
            ratio_change_phase_hold_frames: 0,
            smooth_ratio_updates: false,
            envelope_preservation,
            envelope_strength: 1.0,
            adaptive_envelope_order: true,
            envelope_order,
            cepstrum_buf: Vec::new(),
            analysis_envelope: Vec::new(),
            synthesis_envelope: Vec::new(),
            envelope_noise_floor_scratch: Vec::with_capacity(num_bins),
            envelope_ifft_scratch: vec![COMPLEX_ZERO; fft_inverse_scratch_len],
            envelope_fft_scratch: vec![COMPLEX_ZERO; fft_forward_scratch_len],
            ola_gain,
            ola_gain_f64,
            ola_window_product,
            ola_window_product_f64,
            if_phases_backup: vec![0.0; num_bins],
            output_buf: Vec::new(),
            window_sum_buf: Vec::new(),
            streaming_tail: Vec::new(),
            streaming_tail_window_sum: Vec::new(),
            streaming_tail_ratio: stretch_ratio,
            streaming_tail_phase_ratio: stretch_ratio,
            streaming_accum_output: Vec::new(),
            streaming_accum_window_sum: Vec::new(),
            flux_prev_magnitudes: vec![0.0; num_bins],
            flux_has_prev: false,
            flux_low_end_bin: ((500.0 * fft_size as f32 / sample_rate as f32).floor() as usize)
                .min(num_bins.saturating_sub(1)),
            flux_mid_end_bin: ((4000.0 * fft_size as f32 / sample_rate as f32).floor() as usize)
                .min(num_bins.saturating_sub(1)),
            last_frame_flux: None,
        }
    }

    /// Returns the FFT size.
    #[inline]
    pub fn fft_size(&self) -> usize {
        self.fft_size
    }

    /// Returns the analysis hop size.
    #[inline]
    pub fn hop_analysis(&self) -> usize {
        self.hop_analysis
    }

    /// Returns the per-frame spectral flux from the most recent processing call.
    ///
    /// Available after `process_streaming_into` or `process`. Returns `None`
    /// if no frames have been processed yet.
    #[inline]
    pub fn last_frame_flux(&self) -> Option<&PerFrameFlux> {
        self.last_frame_flux.as_ref()
    }

    /// Returns the synthesis hop size.
    #[inline]
    pub fn hop_synthesis(&self) -> usize {
        self.hop_synthesis
    }

    /// Returns the sub-bass bin cutoff index.
    #[inline]
    pub fn sub_bass_bin(&self) -> usize {
        self.sub_bass_bin
    }

    /// Declares that runtime ratio updates arrive as a smooth, small-step
    /// stream (e.g. the varispeed-first control path's delay-matched
    /// transposition) rather than as discrete automation jumps.
    ///
    /// Disables the ratio-change continuity-focus machinery: phase advance
    /// always follows the current hop ratio instead of slewing from (or
    /// holding at) a previous ratio to mask a seam. Under continuous
    /// modulation the seam-masking slew reads as sustained pitch deviation
    /// on tonal content, which smooth callers must not pay — their per-render
    /// deltas are too small to produce an audible seam in the first place.
    #[inline]
    pub fn set_smooth_ratio_updates(&mut self, smooth: bool) {
        self.smooth_ratio_updates = smooth;
        if smooth {
            self.ratio_change_phase_frames = 0;
            self.ratio_change_phase_total_frames = 0;
            self.ratio_change_phase_hold_frames = 0;
            self.ratio_change_phase_from = self.stretch_ratio;
        }
    }

    /// Updates the stretch ratio without resetting phase state.
    ///
    /// This recalculates the synthesis hop size from the new ratio while
    /// preserving all accumulated phase information. Meaningful runtime ratio
    /// steps also engage a short continuity-focus window so the first few
    /// post-change frames prefer tighter phase coherence at the seam (unless
    /// [`Self::set_smooth_ratio_updates`] is active).
    #[inline]
    pub fn set_stretch_ratio(&mut self, stretch_ratio: f64) {
        if self.smooth_ratio_updates {
            self.stretch_ratio = stretch_ratio;
            self.hop_synthesis = (self.hop_analysis as f64 * stretch_ratio).round() as usize;
            self.ratio_change_phase_from = stretch_ratio;
            return;
        }
        let prior_ratio = self.stretch_ratio;
        let ratio_delta = (stretch_ratio - prior_ratio).abs();
        let in_flight_phase_ratio = self.continuity_focus_phase_ratio(prior_ratio);
        let mut continuity_phase_from = in_flight_phase_ratio;
        let mut reanchored_to_carried_seam = false;
        if !self.streaming_tail.is_empty() {
            let carried_phase_ratio = self.streaming_tail_phase_ratio;
            if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && !ratio_is_meaningfully_above_unity(continuity_phase_from)
                && !ratio_is_meaningfully_below_unity(stretch_ratio)
                && carried_phase_ratio > continuity_phase_from + RATIO_CHANGE_FOCUS_TRIGGER
            {
                // If the explicit continuity window has already drained back
                // near unity but an older expanded overlap is still unresolved,
                // keep near-unity follow-up nudges anchored to that carried
                // seam. Otherwise tiny automation steps can stop refreshing the
                // seam hold while the older expansion tail is still audible.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && !ratio_is_meaningfully_below_unity(continuity_phase_from)
                && !ratio_is_meaningfully_above_unity(stretch_ratio)
                && carried_phase_ratio + RATIO_CHANGE_FOCUS_TRIGGER < continuity_phase_from
            {
                // Mirror the same protection for unresolved compression tails
                // that remain audibly below unity after the continuity window
                // has drifted back toward neutral.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && ratio_is_meaningfully_below_unity(continuity_phase_from)
                && !ratio_is_meaningfully_below_unity(stretch_ratio)
            {
                // Short-interval cross-unity modulation can retarget back
                // toward unity before the previous expansion seam has fully
                // drained. Keep the seam slew anchored to that carried
                // expansion overlap so the next callback does not restart from
                // an already-compressed phase ratio while the older expanded
                // tail is still audible.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && ratio_is_meaningfully_above_unity(continuity_phase_from)
                && !ratio_is_meaningfully_above_unity(stretch_ratio)
            {
                // Mirror the same protection for a carried compression seam so
                // a rebound back toward unity does not restart from an already
                // expanded in-flight ratio while the older compressed overlap
                // is still draining.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && ratio_is_meaningfully_below_unity(continuity_phase_from)
                && stretch_ratio < continuity_phase_from - RATIO_CHANGE_FOCUS_TRIGGER
            {
                // If automation crosses unity and immediately digs farther into
                // compression before the previous expansion seam drains, keep
                // the restart anchored to that older carried seam. Otherwise
                // the continuity slew restarts from the fresher compression
                // ratio and leaves the unresolved expansion overlap to snap on
                // the next callback.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && ratio_is_meaningfully_above_unity(continuity_phase_from)
                && stretch_ratio > continuity_phase_from + RATIO_CHANGE_FOCUS_TRIGGER
            {
                // Mirror the same protection for a carried compression seam
                // when automation rebounds across unity and immediately pushes
                // farther into expansion before the older overlap drains.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && ratio_is_meaningfully_above_unity(continuity_phase_from)
                && stretch_ratio + RATIO_CHANGE_FOCUS_TRIGGER < continuity_phase_from
                && carried_phase_ratio > continuity_phase_from
            {
                // The unresolved seam can still be more expanded than the
                // in-flight phase slew after a short same-side rebound. Keep
                // the restart anchored to that older seam so fast automation
                // does not relax the overlap twice before it has fully drained.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && ratio_is_meaningfully_below_unity(continuity_phase_from)
                && stretch_ratio - RATIO_CHANGE_FOCUS_TRIGGER > continuity_phase_from
                && carried_phase_ratio < continuity_phase_from
            {
                // Mirror the same protection for unresolved compression seams
                // when the target rebounds back toward unity without crossing it.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && ratio_is_meaningfully_above_unity(continuity_phase_from)
                && stretch_ratio > continuity_phase_from + RATIO_CHANGE_FOCUS_TRIGGER
                && carried_phase_ratio > continuity_phase_from
            {
                // If a short same-side relaxation is followed immediately by a
                // renewed expansion step, keep the restart anchored to the
                // older carried seam until that overlap drains. Restarting
                // from the weaker in-flight ratio would loosen the overlap
                // once, then tighten it again on the next callback.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && ratio_is_meaningfully_below_unity(continuity_phase_from)
                && stretch_ratio < continuity_phase_from
                && stretch_ratio + RATIO_CHANGE_FOCUS_TRIGGER >= continuity_phase_from
            {
                // A tiny follow-up push deeper into compression can arrive
                // while an older expansion seam is still unresolved. Keep the
                // continuity slew anchored to that carried seam so the overlap
                // does not restart from the fresher compression ratio just
                // because the latest nudge itself is sub-threshold.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && ratio_is_meaningfully_above_unity(continuity_phase_from)
                && stretch_ratio > continuity_phase_from
                && stretch_ratio - RATIO_CHANGE_FOCUS_TRIGGER <= continuity_phase_from
            {
                // Mirror the same protection for a carried compression seam
                // when a tiny follow-up push deeper into expansion lands
                // before the older opposite-side overlap has drained.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && ratio_is_meaningfully_below_unity(continuity_phase_from)
                && stretch_ratio + RATIO_CHANGE_FOCUS_TRIGGER < continuity_phase_from
                && carried_phase_ratio < continuity_phase_from
            {
                // Mirror the same protection for unresolved compression seams
                // when modulation briefly relaxes toward unity and then pushes
                // deeper into compression before the older overlap has drained.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_above_unity(carried_phase_ratio)
                && ratio_is_meaningfully_above_unity(continuity_phase_from)
                && carried_phase_ratio > continuity_phase_from + RATIO_CHANGE_FOCUS_TRIGGER
                && !ratio_is_meaningfully_below_unity(stretch_ratio)
                && stretch_ratio <= continuity_phase_from + RATIO_CHANGE_FOCUS_TRIGGER
            {
                // Tiny same-side follow-up nudges can arrive while an older,
                // more-expanded overlap is still audible. Keep the continuity
                // slew anchored to that carried seam so the overlap does not
                // relax just because the automation step itself is below the
                // usual trigger threshold.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            } else if ratio_is_meaningfully_below_unity(carried_phase_ratio)
                && ratio_is_meaningfully_below_unity(continuity_phase_from)
                && carried_phase_ratio + RATIO_CHANGE_FOCUS_TRIGGER < continuity_phase_from
                && !ratio_is_meaningfully_above_unity(stretch_ratio)
                && stretch_ratio + RATIO_CHANGE_FOCUS_TRIGGER >= continuity_phase_from
            {
                // Mirror the same protection for unresolved compression seams
                // when a tiny same-side nudge lands near the in-flight ratio
                // while the older, more-compressed overlap is still draining.
                continuity_phase_from = carried_phase_ratio;
                reanchored_to_carried_seam = true;
            }
        }
        let reversed_direction = ratio_change_reverses_inflight_direction(
            continuity_phase_from,
            prior_ratio,
            stretch_ratio,
        );
        let reanchored_far_from_inflight = reanchored_to_carried_seam
            && (continuity_phase_from - in_flight_phase_ratio).abs() >= RATIO_CHANGE_FOCUS_TRIGGER;
        let continuity_delta = (stretch_ratio - continuity_phase_from).abs();
        self.stretch_ratio = stretch_ratio;
        self.hop_synthesis = (self.hop_analysis as f64 * stretch_ratio).round() as usize;
        if ratio_delta >= RATIO_CHANGE_FOCUS_TRIGGER
            || continuity_delta >= RATIO_CHANGE_FOCUS_TRIGGER
        {
            // Fast automation can land another meaningful ratio step before the
            // previous continuity-focus window drains. Refresh the window so
            // late seam frames do not relax back into looser locking mid-burst.
            let mut continuity_focus_frames = continuity_focus_frames_for_ratio_change(
                RATIO_CHANGE_FOCUS_FRAMES,
                self.streaming_tail.len(),
                self.hop_analysis,
            );
            if reversed_direction {
                continuity_focus_frames = continuity_focus_frames
                    .saturating_add(RATIO_CHANGE_REVERSAL_FOCUS_EXTRA_FRAMES);
            }
            if reanchored_far_from_inflight {
                continuity_focus_frames = continuity_focus_frames
                    .saturating_add(RATIO_CHANGE_CARRIED_SEAM_FOCUS_EXTRA_FRAMES);
            }
            self.ratio_change_phase_from = continuity_phase_from;
            self.ratio_change_phase_frames = continuity_focus_frames;
            self.ratio_change_phase_total_frames = continuity_focus_frames;
            self.ratio_change_phase_hold_frames = if reanchored_far_from_inflight {
                RATIO_CHANGE_CARRIED_SEAM_HOLD_FRAMES.min(continuity_focus_frames)
            } else {
                0
            };
            self.transient_focus_frames = self.transient_focus_frames.max(continuity_focus_frames);
        }
    }

    /// Resets the phase accumulator and previous-phase buffers.
    ///
    /// Call this at transient boundaries so that stale phase state from a
    /// previous tonal segment does not contaminate the next one. The PV will
    /// re-derive phases from the first analysis frame after the reset.
    #[inline]
    pub fn reset_phase_state(&mut self) {
        self.phase_accum.fill(0.0);
        self.prev_phase.fill(0.0);
        self.phase_seed_pending.fill(true);
        self.transient_focus_frames = 0;
        self.ratio_change_phase_frames = 0;
        self.ratio_change_phase_total_frames = 0;
        self.ratio_change_phase_hold_frames = 0;
        self.ratio_change_phase_from = self.stretch_ratio;
    }

    /// Clears ALL per-stream state — phase history, the held streaming
    /// overlap tail, synthesis position, and flux tracking — without
    /// deallocating any buffer.
    ///
    /// After this call the vocoder behaves like a freshly constructed one
    /// for streaming purposes, which makes it the allocation-free
    /// equivalent of dropping and rebuilding the instance. Used by
    /// warm-start seek to re-prime from new material.
    pub fn reset_streaming_state(&mut self) {
        self.reset_phase_state();
        self.streaming_tail.clear();
        self.streaming_tail_window_sum.clear();
        self.streaming_tail_ratio = self.stretch_ratio;
        self.streaming_tail_phase_ratio = self.stretch_ratio;
        self.synthesis_pos = 0.0;
        self.synthesis_emitted = 0;
        self.flux_prev_magnitudes.fill(0.0);
        self.flux_has_prev = false;
        self.last_frame_flux = None;
    }

    /// Enables or disables confidence-driven adaptive phase-lock switching.
    #[inline]
    pub fn set_adaptive_phase_locking(&mut self, enabled: bool) {
        self.adaptive_phase_locking = enabled;
    }

    /// Returns whether adaptive phase-lock switching is enabled.
    #[inline]
    pub fn adaptive_phase_locking(&self) -> bool {
        self.adaptive_phase_locking
    }

    /// Sets envelope correction strength (clamped to `[0.0, 2.0]`).
    #[inline]
    pub fn set_envelope_strength(&mut self, strength: f32) {
        self.envelope_strength = strength.clamp(0.0, 2.0);
    }

    /// Enables or disables adaptive per-frame envelope-order selection.
    #[inline]
    pub fn set_adaptive_envelope_order(&mut self, enabled: bool) {
        self.adaptive_envelope_order = enabled;
    }

    /// Selectively resets phase state for specific frequency bands.
    ///
    /// Only zeros `phase_accum` and `prev_phase` for bins within the bands
    /// indicated by `reset_mask`: `[sub_bass, low, mid, high]`.
    /// Band boundaries: sub-bass <100Hz, low 100-500Hz, mid 500-4000Hz, high >4000Hz.
    ///
    /// This avoids disrupting phase tracking in bands where no transient occurred
    /// (e.g., a hi-hat hit shouldn't reset the sustained bass phase).
    pub fn reset_phase_state_bands(&mut self, reset_mask: [bool; 4], sample_rate: u32) {
        let num_bins = self.fft_size / 2 + 1;
        let bin_freq = sample_rate as f32 / self.fft_size as f32;

        for bin in 0..num_bins {
            let freq = bin as f32 * bin_freq;
            let band_idx = if freq < 100.0 {
                0
            } else if freq < 500.0 {
                1
            } else if freq < 4000.0 {
                2
            } else {
                3
            };
            if reset_mask[band_idx] {
                self.phase_accum[bin] = 0.0;
                self.prev_phase[bin] = 0.0;
                self.phase_seed_pending[bin] = true;
            }
        }

        // Engage a short transient-focus window for audible bands so the
        // first few post-reset frames favor tighter attack coherence.
        if reset_mask[1] || reset_mask[2] || reset_mask[3] {
            self.transient_focus_frames = self.transient_focus_frames.max(TRANSIENT_FOCUS_FRAMES);
        }
    }

    /// Stretches a mono audio signal using phase vocoder with identity phase locking.
    pub fn process(&mut self, input: &[f32]) -> Result<Vec<f32>, StretchError> {
        // Batch calls are independent; clear any prior streaming overlap state.
        self.streaming_tail.clear();
        self.streaming_tail_window_sum.clear();
        self.streaming_tail_ratio = self.stretch_ratio;
        self.streaming_tail_phase_ratio = self.stretch_ratio;

        // Mirror-pad input to stabilize edge normalization.
        //
        // The first/last analysis frames have incomplete window overlap,
        // producing an unstable window-sum profile that distorts the
        // normalized output at signal edges. Padding with reflected
        // samples lets the overlap-add window sum reach its steady-state
        // value before processing actual content, eliminating edge
        // artifacts that degrade spectral metrics.
        //
        // Asymmetric padding: the start uses a shorter mirror so that
        // onsets at t=0 remain distinct — longer mirrors pre-condition
        // the phase state with identical spectral content and mask the
        // amplitude transition, causing onset detectors to miss it.
        // The end keeps a longer mirror (hop*8) for full spectral quality.
        //
        // The start padding scales with the ratio distance from unity:
        // at extreme ratios (>0.3 from 1.0) the onset time-scaling
        // amplifies small phase-state artefacts enough to shift onsets
        // beyond the scoring tolerance, so we use a shorter mirror.
        // Near unity the phase state is well-behaved and the longer
        // mirror gives better spectral metrics.
        // Expansion ratios produce longer output, so the PV processes
        // further into the tail padding region.  Extra end padding gives
        // the window-sum normalization more frames to stabilise, reducing
        // amplitude droop at the output tail that degrades LSD.
        let end_pad_mult = if self.stretch_ratio > 1.1 { 10 } else { 8 };
        let end_pad = (self.hop_analysis * end_pad_mult).min(input.len());
        // Graduated start padding: shorter mirrors preserve onset sharpness
        // (critical for TP scoring) while longer mirrors give better spectral
        // metrics (SC/LSD).  Ratios 0.2-0.3 from unity get an intermediate
        // 6-hop padding that balances both concerns.
        let ratio_dist = (self.stretch_ratio - 1.0).abs();
        let start_pad_mult = if ratio_dist > 0.3 {
            4
        } else if ratio_dist > 0.15 {
            6
        } else {
            8
        };
        let start_pad = (self.hop_analysis * start_pad_mult).min(input.len());
        if start_pad > 0 && input.len() >= self.fft_size {
            let padded_len = input.len() + start_pad + end_pad;
            let mut padded = vec![0.0f32; padded_len];
            // Reflect start (shorter — preserves onset sharpness)
            for i in 0..start_pad {
                padded[i] = input[start_pad - 1 - i];
            }
            // Copy original
            padded[start_pad..start_pad + input.len()].copy_from_slice(input);
            // Reflect end with cosine taper: the reflected samples fade
            // smoothly to zero so the PV sees a gradually decaying
            // continuation instead of a cusp.  This reduces phase
            // interference at the boundary and eliminates the sharp
            // amplitude droop in the last few output frames that causes
            // false onset detection at the signal end.
            for i in 0..end_pad {
                let t = (i + 1) as f32 / end_pad as f32;
                let fade = 0.5 * (1.0 + (std::f32::consts::PI * t).cos());
                padded[start_pad + input.len() + i] = input[input.len() - 1 - i] * fade;
            }

            let (_num_frames, output_len, _) = self.process_core(&padded, true)?;
            let mut output = self.output_buf[..output_len].to_vec();
            Self::normalize_output(
                &mut output,
                &self.window_sum_buf[..output_len],
                self.stretch_ratio,
            );

            // Trim output to remove padding artifacts.
            let trim_start = (start_pad as f64 * self.stretch_ratio).round() as usize;
            let expected_len = (input.len() as f64 * self.stretch_ratio).round() as usize;
            let trim_end = (trim_start + expected_len).min(output.len());
            if trim_start < output.len() {
                let mut result = output[trim_start..trim_end].to_vec();

                // Edge amplitude correction: mirror padding causes phase
                // interference near the signal boundary, reducing amplitude
                // in the first ~fft_size*ratio output samples.  Measure
                // this droop and apply a smooth gain ramp to match the
                // steady-state level, preventing false onset detection.
                let correction_len =
                    (self.fft_size as f64 * self.stretch_ratio.abs().max(1.0)).ceil() as usize;
                let correction_len = correction_len.min(result.len() / 4);
                if correction_len >= 128 && result.len() >= correction_len * 4 {
                    let ss_start = correction_len;
                    let ss_end = (correction_len * 3).min(result.len());
                    let ss_len = ss_end - ss_start;
                    let ss_energy: f64 = result[ss_start..ss_end]
                        .iter()
                        .map(|&s| (s as f64) * (s as f64))
                        .sum();
                    let ss_rms = (ss_energy / ss_len as f64).sqrt();
                    if ss_rms > 1e-6 {
                        let edge_energy: f64 = result[..correction_len]
                            .iter()
                            .map(|&s| (s as f64) * (s as f64))
                            .sum();
                        let edge_rms = (edge_energy / correction_len as f64).sqrt();
                        // Only correct if the edge is significantly quieter
                        if edge_rms < ss_rms * 0.7 && edge_rms > 1e-8 {
                            let gain = (ss_rms / edge_rms).min(4.0) as f32;
                            for (i, sample) in result.iter_mut().enumerate().take(correction_len) {
                                // Cubic Hermite ramp: zero derivative at both
                                // endpoints for smoother spectral transition.
                                let t = i as f32 / correction_len as f32;
                                let g = 1.0 + (gain - 1.0) * (1.0 - 3.0 * t * t + 2.0 * t * t * t);
                                *sample *= g;
                            }
                        }
                    }
                }

                // Same correction for the end edge (mirror padding droop).
                if correction_len >= 128 && result.len() >= correction_len * 4 {
                    let ss_start = result.len() / 4;
                    let ss_end = result.len() * 3 / 4;
                    let ss_len = ss_end - ss_start;
                    let ss_energy: f64 = result[ss_start..ss_end]
                        .iter()
                        .map(|&s| (s as f64) * (s as f64))
                        .sum();
                    let ss_rms = (ss_energy / ss_len as f64).sqrt();
                    if ss_rms > 1e-6 {
                        let end_start = result.len() - correction_len;
                        let end_energy: f64 = result[end_start..]
                            .iter()
                            .map(|&s| (s as f64) * (s as f64))
                            .sum();
                        let end_rms = (end_energy / correction_len as f64).sqrt();
                        if end_rms < ss_rms * 0.7 && end_rms > 1e-8 {
                            let gain = (ss_rms / end_rms).min(4.0) as f32;
                            let rlen = result.len();
                            for i in 0..correction_len {
                                // Cubic Hermite ramp: zero derivative at both
                                // endpoints (matches start correction shape).
                                let t = i as f32 / correction_len as f32;
                                let g = 1.0 + (gain - 1.0) * (3.0 * t * t - 2.0 * t * t * t);
                                result[rlen - correction_len + i] *= g;
                            }
                        }
                    }
                }

                return Ok(result);
            }
        }

        // Fallback: process without padding (short inputs or edge cases).
        let (_num_frames, output_len, _) = self.process_core(input, true)?;
        let mut output = self.output_buf[..output_len].to_vec();
        Self::normalize_output(
            &mut output,
            &self.window_sum_buf[..output_len],
            self.stretch_ratio,
        );
        Ok(output)
    }

    /// Pre-allocates all streaming-path buffers for a maximum input window.
    ///
    /// Call once at build time with the caller's rolling-window capacity so
    /// [`Self::process_streaming_into`] performs no allocations afterwards,
    /// even as the stretch ratio moves (bounded by `max_ratio`) and the
    /// render window size varies between calls. Mirrors
    /// `MultiResolutionStretcher::reserve_streaming_capacity`.
    pub fn reserve_streaming_capacity(&mut self, max_window_frames: usize, max_ratio: f64) {
        fn reserve_to(buf: &mut Vec<f32>, capacity: usize) {
            if buf.capacity() < capacity {
                buf.reserve(capacity - buf.len());
            }
        }
        let ratio_mult = (max_ratio.max(1.0).ceil() as usize).saturating_add(1);
        let out_bound = max_window_frames
            .saturating_mul(ratio_mult)
            .saturating_add(self.fft_size.saturating_mul(2));
        reserve_to(&mut self.output_buf, out_bound);
        reserve_to(&mut self.window_sum_buf, out_bound);
        reserve_to(&mut self.streaming_accum_output, out_bound);
        reserve_to(&mut self.streaming_accum_window_sum, out_bound);
        reserve_to(&mut self.streaming_tail, out_bound);
        reserve_to(&mut self.streaming_tail_window_sum, out_bound);
        let num_bins = self.fft_size / 2 + 1;
        if self.peaks.capacity() < num_bins {
            self.peaks.reserve(num_bins - self.peaks.len());
        }
    }

    /// Streaming phase-vocoder pass that preserves phase across calls.
    ///
    /// `input` should include any required analysis overlap context from the
    /// caller (typically managed by a higher-level stream processor). This
    /// method keeps synthesis overlap/window tails internally and emits only
    /// hop-aligned samples that are final for this call.
    pub fn process_streaming(&mut self, input: &[f32]) -> Result<Vec<f32>, StretchError> {
        let mut output = Vec::with_capacity(
            ((input.len() as f64 * self.stretch_ratio).ceil() as usize)
                .saturating_add(self.fft_size),
        );
        self.process_streaming_into(input, &mut output)?;
        Ok(output)
    }

    /// Streaming phase-vocoder pass writing directly into `output`.
    ///
    /// This avoids temporary output allocations in real-time paths.
    pub fn process_streaming_into(
        &mut self,
        input: &[f32],
        output: &mut Vec<f32>,
    ) -> Result<(), StretchError> {
        if input.len() < self.fft_size {
            output.clear();
            return Ok(());
        }

        let (emit_len, output_len, last_phase_hop_ratio) = self.process_core(input, false)?;
        if output.capacity() < emit_len {
            return Err(StretchError::BufferOverflow {
                buffer: "phase_vocoder_stream_output",
                requested: emit_len,
                available: output.capacity(),
            });
        }

        let work_len = output_len
            .max(emit_len)
            .max(self.streaming_tail.len())
            .max(self.streaming_tail_window_sum.len());

        self.streaming_accum_output.resize(work_len, 0.0);
        self.streaming_accum_output.fill(0.0);
        self.streaming_accum_window_sum.resize(work_len, 0.0);
        self.streaming_accum_window_sum.fill(0.0);

        self.streaming_accum_output[..output_len].copy_from_slice(&self.output_buf[..output_len]);
        self.streaming_accum_window_sum[..output_len]
            .copy_from_slice(&self.window_sum_buf[..output_len]);

        let carried_tail_ratio = self.streaming_tail_ratio;
        let carried_tail_phase_ratio = self.streaming_tail_phase_ratio;
        let tail_len = self
            .streaming_tail
            .len()
            .min(self.streaming_tail_window_sum.len());
        for i in 0..tail_len {
            self.streaming_accum_output[i] += self.streaming_tail[i];
            self.streaming_accum_window_sum[i] += self.streaming_tail_window_sum[i];
        }

        // Keep the unresolved overlap region for the next chunk.
        self.streaming_tail.clear();
        self.streaming_tail_window_sum.clear();
        if emit_len < work_len {
            self.streaming_tail
                .extend_from_slice(&self.streaming_accum_output[emit_len..work_len]);
            self.streaming_tail_window_sum
                .extend_from_slice(&self.streaming_accum_window_sum[emit_len..work_len]);
        }

        output.resize(emit_len, 0.0);
        output[..emit_len].copy_from_slice(&self.streaming_accum_output[..emit_len]);
        let emitted_window_sum = &self.streaming_accum_window_sum[..emit_len];
        let emitted_max_window_sum = emitted_window_sum.iter().copied().fold(0.0f32, f32::max);
        let seam_len = tail_len.min(emit_len);
        if seam_len == 0 {
            Self::normalize_output_with_window_floor(
                output,
                emitted_window_sum,
                self.stretch_ratio,
                emitted_max_window_sum,
            );
        } else if seam_len == emit_len {
            Self::normalize_output_with_window_floor(
                output,
                emitted_window_sum,
                streaming_tail_normalize_ratio(carried_tail_ratio, self.stretch_ratio),
                emitted_max_window_sum,
            );
        } else {
            // Once the carried seam has fully drained inside this callback, switch
            // the remainder back to the current ratio's normalization floor so a
            // previous expansion ratio does not keep loosening the whole chunk.
            Self::normalize_output_with_window_floor(
                &mut output[..seam_len],
                &emitted_window_sum[..seam_len],
                streaming_tail_normalize_ratio(carried_tail_ratio, self.stretch_ratio),
                emitted_max_window_sum,
            );
            Self::normalize_output_with_window_floor(
                &mut output[seam_len..emit_len],
                &emitted_window_sum[seam_len..emit_len],
                self.stretch_ratio,
                emitted_max_window_sum,
            );
        }
        self.synthesis_emitted = self.synthesis_emitted.saturating_add(emit_len);
        self.streaming_tail_ratio = if self.streaming_tail.is_empty() {
            self.stretch_ratio
        } else if emit_len < tail_len {
            // Preserve the more conservative seam ratio only while unresolved
            // overlap from the previous chunk is still present in the carried tail.
            streaming_tail_normalize_ratio(carried_tail_ratio, self.stretch_ratio)
        } else {
            self.stretch_ratio
        };
        self.streaming_tail_phase_ratio = if self.streaming_tail.is_empty() {
            self.stretch_ratio
        } else if emit_len < tail_len {
            carried_tail_phase_ratio
        } else {
            last_phase_hop_ratio
        };
        Ok(())
    }

    /// Flushes remaining streaming overlap/window tail at end of stream.
    pub fn flush_streaming(&mut self) -> Result<Vec<f32>, StretchError> {
        let mut output = Vec::with_capacity(self.streaming_tail.len());
        self.flush_streaming_into(&mut output)?;
        Ok(output)
    }

    /// Flushes remaining streaming overlap/window tail into `output`.
    pub fn flush_streaming_into(&mut self, output: &mut Vec<f32>) -> Result<(), StretchError> {
        if self.streaming_tail.is_empty() || self.streaming_tail_window_sum.is_empty() {
            self.streaming_tail.clear();
            self.streaming_tail_window_sum.clear();
            self.streaming_tail_ratio = self.stretch_ratio;
            self.streaming_tail_phase_ratio = self.stretch_ratio;
            self.synthesis_pos = 0.0;
            self.synthesis_emitted = 0;
            output.clear();
            return Ok(());
        }

        let len = self
            .streaming_tail
            .len()
            .min(self.streaming_tail_window_sum.len());
        if output.capacity() < len {
            return Err(StretchError::BufferOverflow {
                buffer: "phase_vocoder_flush_output",
                requested: len,
                available: output.capacity(),
            });
        }
        output.resize(len, 0.0);
        output.copy_from_slice(&self.streaming_tail[..len]);
        Self::normalize_output(
            output,
            &self.streaming_tail_window_sum[..len],
            streaming_tail_normalize_ratio(self.streaming_tail_ratio, self.stretch_ratio),
        );
        self.streaming_tail.clear();
        self.streaming_tail_window_sum.clear();
        self.streaming_tail_ratio = self.stretch_ratio;
        self.streaming_tail_phase_ratio = self.stretch_ratio;
        self.synthesis_pos = 0.0;
        self.synthesis_emitted = 0;
        Ok(())
    }

    /// Shared PV core used by both batch and streaming paths.
    ///
    /// Returns `(emit_len, output_len, tail_phase_ratio)` where `output_len` samples are
    /// accumulated (unnormalized) into `self.output_buf` and
    /// `self.window_sum_buf`, `emit_len` is the number of samples finalized
    /// for streaming emission (`floor(next_synthesis_pos)`), and
    /// `tail_phase_ratio` is the effective phase ratio used by the first
    /// overlap that survives past the callback boundary.
    fn process_core(
        &mut self,
        input: &[f32],
        reset_phase_state: bool,
    ) -> Result<(usize, usize, f64), StretchError> {
        if input.len() < self.fft_size {
            return Err(StretchError::InputTooShort {
                provided: input.len(),
                minimum: self.fft_size,
            });
        }

        let num_bins = self.fft_size / 2 + 1;
        let num_frames = (input.len() - self.fft_size) / self.hop_analysis + 1;

        if reset_phase_state {
            self.phase_accum.fill(0.0);
            self.prev_phase.fill(0.0);
            self.phase_seed_pending.fill(true);
            self.transient_focus_frames = 0;
            self.ratio_change_phase_frames = 0;
            self.ratio_change_phase_total_frames = 0;
            self.ratio_change_phase_hold_frames = 0;
            self.ratio_change_phase_from = self.stretch_ratio;
            self.synthesis_pos = 0.0;
            self.synthesis_emitted = 0;
        }

        let hop_ratio = self.stretch_ratio;
        let frame_advance = self.hop_analysis as f64 * hop_ratio;
        let fft_forward = Arc::clone(&self.fft_forward);
        let fft_inverse = Arc::clone(&self.fft_inverse);

        // Local synthesis timeline starts at the current emission cursor.
        let start_synthesis_pos =
            snap_near_integer((self.synthesis_pos - self.synthesis_emitted as f64).max(0.0));

        // Pre-compute required accumulation length with fractional placement.
        let mut max_write_idx = 0usize;
        let mut synthesis_scan_pos = start_synthesis_pos;
        for _ in 0..num_frames {
            let synthesis_pos = snap_near_integer(synthesis_scan_pos);
            let synthesis_floor = synthesis_pos.floor() as usize;
            let frac = synthesis_pos - synthesis_floor as f64;
            let frame_end = synthesis_floor.saturating_add(
                self.fft_size
                    .saturating_sub(1)
                    .saturating_add(usize::from(frac > SYNTH_POS_EPSILON)),
            );
            max_write_idx = max_write_idx.max(frame_end);
            synthesis_scan_pos = synthesis_pos + frame_advance;
        }
        let output_len = max_write_idx.saturating_add(1);
        let emit_len = snap_near_integer(synthesis_scan_pos).floor() as usize;

        // Reuse pre-allocated buffers, growing if needed (never shrinks).
        self.output_buf.resize(output_len, 0.0);
        self.output_buf.fill(0.0);
        self.window_sum_buf.resize(output_len, 0.0);
        self.window_sum_buf.fill(0.0);

        let mut synthesis_frame_pos = start_synthesis_pos;
        let mut carried_tail_phase_ratio = hop_ratio;
        let mut recorded_carried_tail_phase_ratio = false;
        for frame_idx in 0..num_frames {
            let analysis_pos = frame_idx * self.hop_analysis;
            let synthesis_pos = snap_near_integer(synthesis_frame_pos);
            let synthesis_floor = synthesis_pos.floor() as usize;
            let frac = synthesis_pos - synthesis_floor as f64;
            let phase_hop_ratio = self.continuity_focus_phase_ratio(hop_ratio);

            self.analyze_frame(
                &input[analysis_pos..analysis_pos + self.fft_size],
                &fft_forward,
            );
            self.advance_phases(num_bins, hop_ratio, phase_hop_ratio);
            self.compute_frame_flux(num_bins);
            let frame_end = synthesis_floor.saturating_add(
                self.fft_size
                    .saturating_sub(1)
                    .saturating_add(usize::from(frac > SYNTH_POS_EPSILON)),
            );
            if !recorded_carried_tail_phase_ratio && frame_end >= emit_len {
                carried_tail_phase_ratio = phase_hop_ratio;
                recorded_carried_tail_phase_ratio = true;
            }

            // Save IF-estimated phases before phase locking overwrites them.
            self.if_phases_backup[..num_bins].copy_from_slice(&self.new_phases[..num_bins]);

            // Phase locking: lock non-peak bins to their nearest peak using
            // the analysis phase relationship. Only applies above the sub-bass region.
            let locking_mode = self.select_phase_locking_mode(num_bins);
            apply_phase_locking_realtime(
                locking_mode,
                &self.magnitudes,
                &self.analysis_phases,
                &mut self.new_phases,
                num_bins,
                self.sub_bass_bin,
                &mut self.peaks,
                &mut self.phase_lock_troughs,
                &mut self.phase_lock_pv_phases,
            );

            // Blend IF estimates with locked phases for non-peak bins above sub-bass.
            // At ratio near 1.0, phase locking is very accurate so we trust it fully.
            // As the ratio increases, IF estimates become more valuable for preserving
            // frequency accuracy, so we blend in up to 10% IF (reduced from 30% to improve coherence).
            let transient_focus = self.transient_focus_active();
            let if_blend = if transient_focus {
                0.0
            } else {
                (0.06 * ((hop_ratio - 1.0).abs() / 0.5).min(1.0)).min(0.06)
            };
            if if_blend > 1e-6 {
                for bin in self.sub_bass_bin..num_bins {
                    if self.peaks.binary_search(&bin).is_ok() {
                        continue; // Peak bins keep their locked phase
                    }
                    // Distance-adaptive IF: bins far from peaks get more IF blend
                    // since they benefit less from phase locking and more from
                    // independent frequency tracking.
                    let nearest_dist = match self.peaks.binary_search(&bin) {
                        Ok(_) => 0,
                        Err(idx) => {
                            let lower = if idx > 0 {
                                bin - self.peaks[idx - 1]
                            } else {
                                usize::MAX
                            };
                            let upper = if idx < self.peaks.len() {
                                self.peaks[idx] - bin
                            } else {
                                usize::MAX
                            };
                            lower.min(upper)
                        }
                    };
                    let dist_scale = if nearest_dist > 4 {
                        1.0 + ((nearest_dist - 4) as f64 / 10.0).min(3.0)
                    } else {
                        1.0
                    };
                    let bin_if_blend = (if_blend * dist_scale).min(0.20);
                    let locked = self.new_phases[bin] as f64;
                    let if_est = self.if_phases_backup[bin] as f64;
                    self.new_phases[bin] =
                        ((1.0 - bin_if_blend) * locked + bin_if_blend * if_est) as f32;
                }
            }

            // Spectral envelope preservation: correct magnitudes so formant
            // structure matches the original analysis frame, preventing
            // unnatural timbre shifts.
            if !transient_focus && self.envelope_preservation && self.envelope_strength > 0.0 {
                let effective_order = if self.adaptive_envelope_order {
                    let centroid =
                        spectral_centroid(&self.magnitudes, self.sample_rate, self.fft_size);
                    adaptive_cepstral_order(centroid, self.fft_size)
                } else {
                    self.envelope_order
                };
                // Extract envelope from the original analysis magnitudes
                extract_envelope_with_fft_scratch(
                    &self.magnitudes,
                    num_bins,
                    effective_order,
                    &fft_forward,
                    &fft_inverse,
                    &mut self.envelope_fft_scratch,
                    &mut self.envelope_ifft_scratch,
                    &mut self.cepstrum_buf,
                    &mut self.analysis_envelope,
                );

                // The synthesis magnitudes are the same (PV doesn't change
                // magnitudes), but after phase locking the spectral shape
                // may shift slightly. We extract the synthesis envelope from
                // the current magnitudes and correct.
                // Clone analysis envelope as synthesis baseline since magnitudes
                // haven't changed. The correction step then normalizes any
                // spectral tilt introduced by windowing or overlap.
                self.synthesis_envelope.clear();
                self.synthesis_envelope
                    .extend_from_slice(&self.analysis_envelope);

                self.if_phases_backup[..num_bins].copy_from_slice(&self.magnitudes[..num_bins]);

                apply_envelope_correction_with_scratch(
                    &mut self.magnitudes,
                    &self.analysis_envelope,
                    &self.synthesis_envelope,
                    num_bins,
                    self.sub_bass_bin,
                    &mut self.envelope_noise_floor_scratch,
                );

                // Blend toward corrected magnitudes based on configured strength.
                // `1.0` keeps full correction; values in (0,1) soften correction;
                // values >1.0 emphasize correction.
                if (self.envelope_strength - 1.0).abs() > 1e-6 {
                    for bin in self.sub_bass_bin..num_bins {
                        let corrected = self.magnitudes[bin];
                        let original = self.if_phases_backup[bin];
                        self.magnitudes[bin] =
                            original + (corrected - original) * self.envelope_strength;
                    }
                }
            }

            self.reconstruct_spectrum(num_bins);
            fft_inverse.process_with_scratch(&mut self.fft_buffer, &mut self.fft_inverse_scratch);

            // Fractional overlap-add: when synthesis frame starts between samples,
            // distribute each sample between nearest output samples via linear interpolation.
            if frac <= SYNTH_POS_EPSILON {
                for i in 0..self.fft_size {
                    let idx = synthesis_floor + i;
                    self.output_buf[idx] += self.fft_buffer[i].re * self.ola_gain[i];
                    self.window_sum_buf[idx] += self.ola_window_product[i];
                }
            } else {
                let w0 = 1.0 - frac;
                let w1 = frac;
                for i in 0..self.fft_size {
                    let idx = synthesis_floor + i;
                    let sample = self.fft_buffer[i].re as f64 * self.ola_gain_f64[i];
                    let window_weight = self.ola_window_product_f64[i];

                    self.output_buf[idx] += (sample * w0) as f32;
                    self.output_buf[idx + 1] += (sample * w1) as f32;
                    self.window_sum_buf[idx] += (window_weight * w0) as f32;
                    self.window_sum_buf[idx + 1] += (window_weight * w1) as f32;
                }
            }

            self.decay_transient_focus();

            synthesis_frame_pos = synthesis_pos + frame_advance;
        }

        let next_local_synthesis_pos = snap_near_integer(synthesis_frame_pos);
        self.synthesis_pos = self.synthesis_emitted as f64 + next_local_synthesis_pos;

        Ok((emit_len, output_len, carried_tail_phase_ratio))
    }

    /// Windows the input frame and transforms to frequency domain.
    #[inline]
    fn analyze_frame(
        &mut self,
        input_frame: &[f32],
        fft_forward: &std::sync::Arc<dyn rustfft::Fft<f32>>,
    ) {
        let len = input_frame.len().min(self.fft_buffer.len());
        for (i, (&sample, &w)) in input_frame
            .iter()
            .zip(self.window.iter())
            .enumerate()
            .take(len)
        {
            self.fft_buffer[i] = Complex::new(sample * w, 0.0);
        }
        fft_forward.process_with_scratch(&mut self.fft_buffer, &mut self.fft_forward_scratch);
    }

    /// Extracts magnitudes and advances phase accumulators for each bin.
    ///
    /// Uses a multi-pass approach for improved frequency tracking and phase coherence:
    /// 1. Compute magnitudes and raw analysis phases for all bins.
    /// 2. Detect spectral peaks and compute refined instantaneous frequencies via
    ///    parabolic interpolation of the log-magnitude spectrum.
    /// 3. Advance phases using instantaneous frequency (IF) estimation: compute the
    ///    true frequency of each bin from the phase difference, then resynthesize at
    ///    the correct rate using `inst_freq * hop_synthesis`. This naturally handles
    ///    the stretch ratio without explicit hop_ratio multiplication and eliminates
    ///    cumulative phase drift.
    /// 4. Apply soft phase gradient integration to propagate coherent phase from
    ///    peaks to nearby non-peak bins.
    ///
    /// Sub-bass bins (below `sub_bass_bin`) use rigid phase propagation to prevent
    /// phase cancellation in the critical sub-bass region and are excluded from
    /// peak-based refinements.
    ///
    /// Phase accumulation uses f64 precision to prevent cumulative rounding errors
    /// over long signals. The final phases are converted back to f32 for the
    /// spectrum reconstruction step.
    #[inline]
    fn advance_phases(&mut self, num_bins: usize, hop_ratio: f64, phase_hop_ratio: f64) {
        let hop_a = self.hop_analysis as f64;
        let fft = self.fft_size as f64;

        // --- Pass 1: Extract magnitudes and analysis phases ---
        for bin in 0..num_bins {
            let c = self.fft_buffer[bin];
            self.magnitudes[bin] = c.norm();
            self.analysis_phases[bin] = c.arg();
        }

        // --- Pass 2: Detect peaks for IF refinement + phase gradient ---
        let search_start = self.sub_bass_bin.max(1);
        self.peaks.clear();
        if num_bins >= 3 && search_start < num_bins.saturating_sub(1) {
            for bin in search_start..num_bins - 1 {
                if self.magnitudes[bin] > MIN_PEAK_MAGNITUDE
                    && self.magnitudes[bin] > self.magnitudes[bin - 1]
                    && self.magnitudes[bin] > self.magnitudes[bin + 1]
                {
                    self.peaks.push(bin);
                }
            }
        }

        // --- Pass 3: Advance phases using instantaneous frequency (IF) estimation ---
        //
        // For each bin we compute the true instantaneous frequency from the phase
        // difference between consecutive frames, then advance the synthesis phase
        // accumulator by `inst_freq * hop_synthesis`. This naturally accounts for
        // the stretch ratio and eliminates cumulative drift that occurs with the
        // simpler `(expected + deviation) * hop_ratio` approach.
        //
        // For spectral peak bins, parabolic interpolation of the log-magnitude
        // spectrum refines the frequency estimate to sub-bin precision (~1 Hz
        // accuracy vs ~5 Hz for integer-bin estimation).
        for bin in 0..num_bins {
            let phase = self.analysis_phases[bin] as f64;

            // Seed bins that were explicitly reset (full or selective reset)
            // from the current analysis phase to avoid a bogus first IF jump.
            if self.phase_seed_pending[bin] {
                self.phase_accum[bin] = phase;
                self.new_phases[bin] = phase as f32;
                self.prev_phase[bin] = phase;
                self.phase_seed_pending[bin] = false;
                continue;
            }

            if bin < self.sub_bass_bin {
                // Sub-bass IF estimation: same instantaneous-frequency approach
                // as standard bins, but without parabolic interpolation (sub-bass
                // bins are narrow enough that integer-bin IF is sufficient).
                // The identity phase locking in phase_locking.rs handles inter-bin
                // coherence for sub-bass via trough-based regions.
                let expected_diff = self.expected_phase_advance[bin];
                let phase_diff = phase - self.prev_phase[bin];
                let deviation = wrap_phase_f64(phase_diff - expected_diff);
                self.phase_accum[bin] += (expected_diff + deviation) * phase_hop_ratio;
            } else {
                // Standard IF estimation:
                //   phase_diff = current_phase - prev_phase
                //   expected_diff = 2*pi * bin * hop_analysis / fft_size
                //   deviation = wrap(phase_diff - expected_diff)
                //   inst_freq = (expected_diff + deviation) / hop_analysis  [rad/sample]
                //   phase_accum += inst_freq * hop_synthesis
                let expected_diff = self.expected_phase_advance[bin]; // 2*pi*bin*hop_a/fft
                let phase_diff = phase - self.prev_phase[bin];
                let deviation = wrap_phase_f64(phase_diff - expected_diff);

                // For peak bins, use parabolic interpolation of log-magnitude
                // to refine the frequency estimate to sub-bin precision.
                //
                // The phase advance for synthesis is computed as:
                //   inst_freq * hop_synthesis = (expected + deviation) * hop_synthesis / hop_analysis
                //
                // To minimize floating-point roundoff (especially at ratio 1.0 where
                // hop_s == hop_a), we compute `(expected + deviation) * (hop_s / hop_a)`
                // using a single ratio multiplication rather than separate div+mul.
                let is_peak = self.peaks.binary_search(&bin).is_ok();
                let phase_advance = if is_peak
                    && bin >= 1
                    && bin + 1 < num_bins
                    && self.magnitudes[bin] > MIN_PEAK_MAGNITUDE
                {
                    // Parabolic interpolation on log-magnitudes for sub-bin accuracy:
                    //   alpha = log(M[k-1]), beta = log(M[k]), gamma = log(M[k+1])
                    //   p = 0.5 * (alpha - gamma) / (alpha - 2*beta + gamma)
                    //   refined_freq_bin = k + p
                    //
                    // Log-magnitude interpolation gives better accuracy for Gaussian
                    // spectral peaks (which approximate windowed sinusoids) compared
                    // to linear interpolation.
                    let m_prev = (self.magnitudes[bin - 1] as f64).max(1e-30);
                    let m_curr = (self.magnitudes[bin] as f64).max(1e-30);
                    let m_next = (self.magnitudes[bin + 1] as f64).max(1e-30);
                    let alpha = m_prev.ln();
                    let beta = m_curr.ln();
                    let gamma = m_next.ln();
                    let denom = alpha - 2.0 * beta + gamma;
                    if denom.abs() > 1e-12 {
                        let p = 0.5 * (alpha - gamma) / denom;
                        // Refined expected phase advance based on interpolated bin position
                        let refined_expected = TWO_PI_F64 * (bin as f64 + p) * hop_a / fft;
                        let refined_deviation = wrap_phase_f64(phase_diff - refined_expected);
                        (refined_expected + refined_deviation) * phase_hop_ratio
                    } else {
                        (expected_diff + deviation) * phase_hop_ratio
                    }
                } else {
                    (expected_diff + deviation) * phase_hop_ratio
                };

                self.phase_accum[bin] += phase_advance;
            }

            self.new_phases[bin] = self.phase_accum[bin] as f32;
            self.prev_phase[bin] = phase;
        }

        // --- Phase gradient integration (soft vertical coherence) ---
        // Propagate phase from peaks to nearby non-peak bins using the analysis
        // phase gradient, blended with the independently-advanced phase.
        // Apply up to 2.5x ratio with a tapering blend around unity:
        // full strength near ratio≈1.0, gradually reduced as we move away on
        // either side (compression or expansion) to avoid over-locking artifacts.
        if !self.transient_focus_active() && !self.peaks.is_empty() && hop_ratio < 2.5 {
            let ratio_distance = (hop_ratio - 1.0).abs();
            // For ratios >1.0, taper gradient locking faster to avoid
            // over-locking smearing at larger slowdowns.
            let taper_span = if hop_ratio > 1.0 { 1.2 } else { 1.5 };
            let gradient_blend =
                PHASE_GRADIENT_BLEND * (1.0 - (ratio_distance / taper_span).clamp(0.0, 1.0));
            for bin in self.sub_bass_bin..num_bins {
                if self.peaks.binary_search(&bin).is_ok() {
                    continue; // Peak bins keep their phase (they are the anchors)
                }

                // Find the nearest peak via binary search
                let nearest_peak = match self.peaks.binary_search(&bin) {
                    Ok(_) => unreachable!(),
                    Err(idx) => {
                        let lower = if idx > 0 {
                            Some(self.peaks[idx - 1])
                        } else {
                            None
                        };
                        let upper = if idx < self.peaks.len() {
                            Some(self.peaks[idx])
                        } else {
                            None
                        };
                        match (lower, upper) {
                            (Some(l), Some(u)) => {
                                if bin - l <= u - bin {
                                    l
                                } else {
                                    u
                                }
                            }
                            (Some(l), None) => l,
                            (None, Some(u)) => u,
                            (None, None) => continue,
                        }
                    }
                };

                // Attenuate gradient blend for bins far from their anchor peak.
                // Distant bins have less acoustic relationship to the peak,
                // so independent phase evolution is more appropriate.
                let peak_distance = bin.abs_diff(nearest_peak);
                let distance_fade = if peak_distance > 8 {
                    (1.0 - ((peak_distance - 8) as f64 / 24.0).min(1.0)).max(0.0)
                } else {
                    1.0
                };
                let effective_blend = gradient_blend * distance_fade;

                let gradient =
                    self.analysis_phases[bin] as f64 - self.analysis_phases[nearest_peak] as f64;
                let propagated = self.new_phases[nearest_peak] as f64 + gradient;
                let independent = self.new_phases[bin] as f64;
                self.new_phases[bin] =
                    ((1.0 - effective_blend) * independent + effective_blend * propagated) as f32;
            }
        }
    }

    /// Computes per-bin spectral flux from the current magnitudes and updates
    /// `last_frame_flux`. Uses half-wave rectification (onset-only).
    ///
    /// Must be called after `advance_phases` which populates `self.magnitudes`.
    #[inline]
    fn compute_frame_flux(&mut self, num_bins: usize) {
        if !self.flux_has_prev {
            // First frame: seed prev_magnitudes, no flux to compute.
            self.flux_prev_magnitudes[..num_bins].copy_from_slice(&self.magnitudes[..num_bins]);
            self.flux_has_prev = true;
            self.last_frame_flux = Some(PerFrameFlux::default());
            return;
        }

        let mut sub_bass_flux = 0.0f32;
        let mut low_flux = 0.0f32;
        let mut mid_flux = 0.0f32;
        let mut high_flux = 0.0f32;
        let mut transient_bin_count = 0u16;
        let mut total_bins_rising = 0u16;
        // Running sum for mean computation (for threshold).
        let mut flux_sum = 0.0f32;
        let mut flux_count = 0u16;

        for bin in 1..num_bins {
            let diff = self.magnitudes[bin] - self.flux_prev_magnitudes[bin];
            if diff > 0.0 {
                total_bins_rising += 1;
                flux_sum += diff;
                flux_count += 1;

                if bin < self.sub_bass_bin {
                    sub_bass_flux += diff;
                } else if bin < self.flux_low_end_bin {
                    low_flux += diff;
                } else if bin < self.flux_mid_end_bin {
                    mid_flux += diff;
                } else {
                    high_flux += diff;
                }
            }
        }

        // Count bins with flux significantly above the mean (transient-like).
        if flux_count > 0 {
            let mean_flux = flux_sum / flux_count as f32;
            let threshold = mean_flux * 4.0; // 4× mean = significant onset
            for bin in 1..num_bins {
                let diff = self.magnitudes[bin] - self.flux_prev_magnitudes[bin];
                if diff > threshold {
                    transient_bin_count += 1;
                }
            }
        }

        self.flux_prev_magnitudes[..num_bins].copy_from_slice(&self.magnitudes[..num_bins]);
        self.last_frame_flux = Some(PerFrameFlux {
            sub_bass: sub_bass_flux,
            low: low_flux,
            mid: mid_flux,
            high: high_flux,
            transient_bin_count,
            total_bins_rising,
        });
    }

    /// Reconstructs the complex spectrum from magnitudes and phases,
    /// then mirrors negative frequencies for inverse FFT.
    #[inline]
    fn reconstruct_spectrum(&mut self, num_bins: usize) {
        for i in 0..num_bins {
            self.fft_buffer[i] = Complex::from_polar(self.magnitudes[i], self.new_phases[i]);
        }
        for bin in 1..num_bins - 1 {
            self.fft_buffer[self.fft_size - bin] = self.fft_buffer[bin].conj();
        }
    }

    /// Normalizes output by window sum, clamping to prevent amplification in
    /// low-overlap regions (occurs when synthesis hop > analysis hop).
    #[inline]
    fn normalize_output(output: &mut [f32], window_sum: &[f32], stretch_ratio: f64) {
        let max_window_sum = window_sum.iter().copied().fold(0.0f32, f32::max);
        Self::normalize_output_with_window_floor(output, window_sum, stretch_ratio, max_window_sum);
    }

    #[inline]
    fn normalize_output_with_window_floor(
        output: &mut [f32],
        window_sum: &[f32],
        stretch_ratio: f64,
        max_window_sum: f32,
    ) {
        // For stretches >1.0, synthesis frames are farther apart and aggressive
        // flooring can over-attenuate low-overlap regions. Relax the floor more
        // strongly with ratio while keeping a safety minimum against blow-ups.
        let floor_ratio = if stretch_ratio > 1.0 {
            (WINDOW_SUM_FLOOR_RATIO / (stretch_ratio as f32 * stretch_ratio as f32))
                .clamp(0.005, WINDOW_SUM_FLOOR_RATIO)
        } else {
            WINDOW_SUM_FLOOR_RATIO
        };
        let min_window_sum = (max_window_sum * floor_ratio).max(WINDOW_SUM_EPSILON);
        let len = output.len().min(window_sum.len());
        for i in 0..len {
            output[i] /= window_sum[i].max(min_window_sum);
        }
    }

    /// Chooses phase-locking mode for the current frame.
    ///
    /// When adaptive switching is enabled, this computes simple confidence
    /// features from the current magnitude spectrum:
    /// - spectral flatness (noise-likeness),
    /// - crest factor (peak prominence),
    /// - stretch-ratio distance from unity.
    ///
    /// Heuristic:
    /// - noisy/weak frames or large ratio offsets -> ROI
    /// - harmonic/peaky frames near unity ratio -> Identity
    /// - harmonic frames at moderate ratio offsets -> Selective
    /// - otherwise -> configured `phase_locking_mode`
    #[inline]
    fn select_phase_locking_mode(&self, num_bins: usize) -> PhaseLockingMode {
        if self.transient_focus_active() {
            return PhaseLockingMode::Identity;
        }

        if !self.adaptive_phase_locking {
            return self.phase_locking_mode;
        }
        if num_bins == 0 || self.sub_bass_bin >= num_bins {
            return self.phase_locking_mode;
        }

        let start = self.sub_bass_bin;
        let mags = &self.magnitudes[start..num_bins];
        if mags.is_empty() {
            return self.phase_locking_mode;
        }

        let mut max_mag = 0.0f64;
        let mut sum = 0.0f64;
        let mut log_sum = 0.0f64;
        for &m in mags {
            let mf = (m as f64).max(ADAPTIVE_FEATURE_EPS);
            max_mag = max_mag.max(mf);
            sum += mf;
            log_sum += mf.ln();
        }

        let n = mags.len() as f64;
        let mean = sum / n;
        let geometric = (log_sum / n).exp();
        let flatness = geometric / mean.max(ADAPTIVE_FEATURE_EPS);
        let crest = max_mag / mean.max(ADAPTIVE_FEATURE_EPS);
        let ratio_distance = (self.stretch_ratio - 1.0).abs();
        let harmonic_confidence = crest / (1.0 + 4.0 * flatness);

        if ratio_distance >= ADAPTIVE_FORCE_ROI_RATIO_DISTANCE
            || flatness >= ADAPTIVE_NOISY_FLATNESS
            || crest <= ADAPTIVE_NOISY_CREST
        {
            return PhaseLockingMode::RegionOfInfluence;
        }

        if ratio_distance <= ADAPTIVE_IDENTITY_RATIO_DISTANCE_MAX
            && harmonic_confidence >= ADAPTIVE_HARMONIC_CONFIDENCE
        {
            return PhaseLockingMode::Identity;
        }

        if ratio_distance <= ADAPTIVE_SELECTIVE_RATIO_DISTANCE_MAX
            && harmonic_confidence >= ADAPTIVE_SELECTIVE_HARMONIC_CONFIDENCE
        {
            return PhaseLockingMode::Selective;
        }

        self.phase_locking_mode
    }

    #[inline]
    fn transient_focus_active(&self) -> bool {
        self.transient_focus_frames > 0
    }

    #[inline]
    fn continuity_focus_phase_ratio(&self, hop_ratio: f64) -> f64 {
        if self.ratio_change_phase_frames == 0 || self.ratio_change_phase_total_frames == 0 {
            return hop_ratio;
        }

        let remaining = self
            .ratio_change_phase_frames
            .min(self.ratio_change_phase_total_frames);
        let completed = self.ratio_change_phase_total_frames - remaining;
        let hold_frames = self
            .ratio_change_phase_hold_frames
            .min(self.ratio_change_phase_total_frames);
        if completed < hold_frames {
            return self.ratio_change_phase_from;
        }

        let effective_completed = completed - hold_frames;
        let effective_total = self.ratio_change_phase_total_frames - hold_frames;
        let progress = (effective_completed as f64 + 1.0) / (effective_total as f64 + 1.0);
        self.ratio_change_phase_from + (hop_ratio - self.ratio_change_phase_from) * progress
    }

    #[inline]
    fn decay_transient_focus(&mut self) {
        self.transient_focus_frames = self.transient_focus_frames.saturating_sub(1);
        self.ratio_change_phase_frames = self.ratio_change_phase_frames.saturating_sub(1);
    }
}

impl std::fmt::Debug for PhaseVocoder {
    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
        f.debug_struct("PhaseVocoder")
            .field("fft_size", &self.fft_size)
            .field("hop_analysis", &self.hop_analysis)
            .field("hop_synthesis", &self.hop_synthesis)
            .field("stretch_ratio", &self.stretch_ratio)
            .field("synthesis_pos", &self.synthesis_pos)
            .field("sub_bass_bin", &self.sub_bass_bin)
            .field("phase_locking_mode", &self.phase_locking_mode)
            .field("adaptive_phase_locking", &self.adaptive_phase_locking)
            .field("envelope_strength", &self.envelope_strength)
            .field("adaptive_envelope_order", &self.adaptive_envelope_order)
            .field("streaming_tail_len", &self.streaming_tail.len())
            .finish()
    }
}

/// Snaps values extremely close to integer grid points to the exact integer.
#[inline]
fn snap_near_integer(value: f64) -> f64 {
    let rounded = value.round();
    if (value - rounded).abs() <= SYNTH_POS_EPSILON {
        rounded
    } else {
        value
    }
}

/// Wraps a phase value to [-PI, PI] using f64 precision.
#[inline]
fn wrap_phase_f64(phase: f64) -> f64 {
    let pi = std::f64::consts::PI;
    let p = phase + pi;
    p - (p / TWO_PI_F64).floor() * TWO_PI_F64 - pi
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::stretch::phase_locking::apply_phase_locking;
    use std::f32::consts::PI;

    const TWO_PI: f32 = 2.0 * PI;

    /// Wraps a phase value to [-PI, PI] using efficient modulo arithmetic (f32).
    fn wrap_phase(phase: f32) -> f32 {
        let p = phase + PI;
        p - (p / TWO_PI).floor() * TWO_PI - PI
    }

    #[test]
    fn test_wrap_phase() {
        assert!((wrap_phase(0.0) - 0.0).abs() < 1e-6);
        assert!((wrap_phase(PI + 0.1) - (-PI + 0.1)).abs() < 1e-5);
        assert!((wrap_phase(-PI - 0.1) - (PI - 0.1)).abs() < 1e-5);
        // Test larger values
        assert!((wrap_phase(10.0 * PI + 0.5) - wrap_phase(0.5)).abs() < 1e-4);
        assert!((wrap_phase(-10.0 * PI - 0.5) - wrap_phase(-0.5)).abs() < 1e-4);
    }

    #[test]
    fn test_normalize_output_with_window_floor_scopes_expansion_floor_to_seam_prefix() {
        let reference_max_window_sum = 1.0f32;
        let window_sum = [1.0f32, 0.01, 0.01, 0.01];
        let mut split_output = vec![1.0f32; window_sum.len()];
        let mut all_expansion_output = vec![1.0f32; window_sum.len()];

        PhaseVocoder::normalize_output_with_window_floor(
            &mut split_output[..2],
            &window_sum[..2],
            1.18,
            reference_max_window_sum,
        );
        PhaseVocoder::normalize_output_with_window_floor(
            &mut split_output[2..],
            &window_sum[2..],
            0.82,
            reference_max_window_sum,
        );
        PhaseVocoder::normalize_output_with_window_floor(
            &mut all_expansion_output,
            &window_sum,
            1.18,
            reference_max_window_sum,
        );

        assert!(
            (split_output[0] - all_expansion_output[0]).abs() < 1e-6
                && (split_output[1] - all_expansion_output[1]).abs() < 1e-6,
            "the unresolved seam prefix should keep the prior expansion ratio floor"
        );
        assert!(
            split_output[2] < all_expansion_output[2] && split_output[3] < all_expansion_output[3],
            "once the seam prefix drains, the remainder should normalize against the current ratio instead of the prior expansion floor"
        );
    }

    #[test]
    fn test_adaptive_phase_locking_prefers_roi_for_flat_spectrum() {
        let mut pv = PhaseVocoder::with_options(
            1024,
            256,
            1.0,
            44100,
            0.0,
            WindowType::Hann,
            PhaseLockingMode::Identity,
        );
        pv.set_adaptive_phase_locking(true);
        pv.magnitudes.fill(1.0);
        let mode = pv.select_phase_locking_mode(1024 / 2 + 1);
        assert_eq!(mode, PhaseLockingMode::RegionOfInfluence);
    }

    #[test]
    fn test_adaptive_phase_locking_prefers_identity_for_harmonic_frame() {
        let mut pv = PhaseVocoder::with_options(
            1024,
            256,
            1.05,
            44100,
            0.0,
            WindowType::Hann,
            PhaseLockingMode::RegionOfInfluence,
        );
        pv.set_adaptive_phase_locking(true);
        pv.magnitudes.fill(0.001);
        pv.magnitudes[20] = 1.0;
        pv.magnitudes[60] = 0.8;
        pv.magnitudes[120] = 0.6;
        let mode = pv.select_phase_locking_mode(1024 / 2 + 1);
        assert_eq!(mode, PhaseLockingMode::Identity);
    }

    #[test]
    fn test_adaptive_phase_locking_prefers_selective_for_moderate_ratio_harmonic_frame() {
        let mut pv = PhaseVocoder::with_options(
            1024,
            256,
            1.45,
            44100,
            0.0,
            WindowType::Hann,
            PhaseLockingMode::RegionOfInfluence,
        );
        pv.set_adaptive_phase_locking(true);
        pv.magnitudes.fill(0.001);
        pv.magnitudes[20] = 1.0;
        pv.magnitudes[60] = 0.8;
        pv.magnitudes[120] = 0.6;
        let mode = pv.select_phase_locking_mode(1024 / 2 + 1);
        assert_eq!(mode, PhaseLockingMode::Selective);
    }

    #[test]
    fn test_adaptive_phase_locking_disabled_uses_configured_mode() {
        let mut pv = PhaseVocoder::with_options(
            1024,
            256,
            1.0,
            44100,
            0.0,
            WindowType::Hann,
            PhaseLockingMode::Identity,
        );
        pv.set_adaptive_phase_locking(false);
        pv.magnitudes.fill(1.0);
        let mode = pv.select_phase_locking_mode(1024 / 2 + 1);
        assert_eq!(mode, PhaseLockingMode::Identity);
    }

    #[test]
    fn test_envelope_strength_setter_clamps() {
        let mut pv = PhaseVocoder::new(1024, 256, 1.0, 44100, 120.0);
        pv.set_envelope_strength(3.0);
        assert!((pv.envelope_strength - 2.0).abs() < 1e-6);
        pv.set_envelope_strength(-1.0);
        assert!((pv.envelope_strength - 0.0).abs() < 1e-6);
    }

    #[test]
    fn test_adaptive_envelope_order_setter() {
        let mut pv = PhaseVocoder::new(1024, 256, 1.0, 44100, 120.0);
        pv.set_adaptive_envelope_order(false);
        assert!(!pv.adaptive_envelope_order);
        pv.set_adaptive_envelope_order(true);
        assert!(pv.adaptive_envelope_order);
    }

    #[test]
    fn test_reset_phase_state_bands_enables_transient_focus_for_audible_bands() {
        let sample_rate = 44_100u32;
        let mut pv = PhaseVocoder::new(1024, 256, 1.0, sample_rate, 120.0);
        assert_eq!(pv.transient_focus_frames, 0);

        pv.reset_phase_state_bands([false, false, true, false], sample_rate);
        assert_eq!(
            pv.transient_focus_frames, TRANSIENT_FOCUS_FRAMES,
            "mid-band reset should engage transient focus"
        );
    }

    #[test]
    fn test_reset_phase_state_bands_sub_only_does_not_enable_transient_focus() {
        let sample_rate = 44_100u32;
        let mut pv = PhaseVocoder::new(1024, 256, 1.0, sample_rate, 120.0);
        assert_eq!(pv.transient_focus_frames, 0);

        pv.reset_phase_state_bands([true, false, false, false], sample_rate);
        assert_eq!(
            pv.transient_focus_frames, 0,
            "sub-only reset should not engage transient focus"
        );
    }

    #[test]
    fn test_select_phase_locking_mode_forces_identity_during_transient_focus() {
        let mut pv = PhaseVocoder::with_options(
            1024,
            256,
            1.4,
            44_100,
            120.0,
            WindowType::Hann,
            PhaseLockingMode::RegionOfInfluence,
        );
        pv.transient_focus_frames = 1;
        pv.set_adaptive_phase_locking(true);
        pv.magnitudes.fill(1.0);

        let mode = pv.select_phase_locking_mode(1024 / 2 + 1);
        assert_eq!(
            mode,
            PhaseLockingMode::Identity,
            "transient focus should force identity locking"
        );
    }

    #[test]
    fn test_reset_phase_state_bands_marks_only_target_band_for_seeding() {
        let sample_rate = 44_100u32;
        let mut pv = PhaseVocoder::new(1024, 256, 1.0, sample_rate, 120.0);
        pv.phase_accum.fill(1.0);
        pv.prev_phase.fill(1.0);
        pv.phase_seed_pending.fill(false);

        // Reset only the low band [100, 500) Hz.
        pv.reset_phase_state_bands([false, true, false, false], sample_rate);

        let bin_hz = sample_rate as f32 / pv.fft_size as f32;
        for bin in 0..pv.phase_accum.len() {
            let freq = bin as f32 * bin_hz;
            let in_low_band = (100.0..500.0).contains(&freq);
            if in_low_band {
                assert_eq!(pv.phase_accum[bin], 0.0, "low-band phase_accum not reset");
                assert_eq!(pv.prev_phase[bin], 0.0, "low-band prev_phase not reset");
                assert!(
                    pv.phase_seed_pending[bin],
                    "low-band seed flag should be set"
                );
            } else {
                assert_eq!(pv.phase_accum[bin], 1.0, "non-low-band phase_accum changed");
                assert_eq!(pv.prev_phase[bin], 1.0, "non-low-band prev_phase changed");
                assert!(
                    !pv.phase_seed_pending[bin],
                    "non-low-band seed flag should stay clear"
                );
            }
        }
    }

    #[test]
    fn test_advance_phases_seeds_pending_bins_without_if_jump() {
        let sample_rate = 44_100u32;
        let mut pv = PhaseVocoder::new(1024, 256, 1.0, sample_rate, 120.0);
        let num_bins = pv.fft_size / 2 + 1;

        // Simulate steady state where no bins are pending.
        pv.phase_seed_pending.fill(false);
        pv.phase_accum.fill(0.5);
        pv.prev_phase.fill(0.5);

        let target_bin = 64usize;
        pv.phase_accum[target_bin] = 0.0;
        pv.prev_phase[target_bin] = 0.0;
        pv.phase_seed_pending[target_bin] = true;

        for bin in 0..num_bins {
            let phase = 0.2 + bin as f32 * 0.001;
            pv.fft_buffer[bin] = Complex::from_polar(1.0, phase);
        }

        pv.advance_phases(num_bins, 1.0, 1.0);

        let seeded_phase = pv.analysis_phases[target_bin] as f64;
        assert!(
            (pv.phase_accum[target_bin] - seeded_phase).abs() < 1e-9,
            "pending bin should seed directly from analysis phase"
        );
        assert!(
            (pv.prev_phase[target_bin] - seeded_phase).abs() < 1e-9,
            "pending bin prev phase should match seeded analysis phase"
        );
        assert!(
            !pv.phase_seed_pending[target_bin],
            "pending flag should clear after first seeded frame"
        );
    }

    #[test]
    fn test_advance_phases_skips_phase_gradient_during_transient_focus() {
        let sample_rate = 44_100u32;
        let num_bins = 1024 / 2 + 1;
        let peak_bin = 20usize;
        let target_bin = 21usize;

        let configure = |pv: &mut PhaseVocoder| {
            pv.phase_seed_pending.fill(false);
            pv.phase_accum.fill(0.0);

            for bin in 0..num_bins {
                pv.prev_phase[bin] = -pv.expected_phase_advance[bin];
                let magnitude = if bin == peak_bin {
                    2.0
                } else if bin == peak_bin - 1 || bin == peak_bin + 1 {
                    0.5
                } else {
                    0.1
                };
                pv.fft_buffer[bin] = Complex::from_polar(magnitude, 0.0);
            }
        };

        let mut focused = PhaseVocoder::new(1024, 256, 1.0, sample_rate, 0.0);
        configure(&mut focused);
        focused.transient_focus_frames = 1;
        focused.advance_phases(num_bins, 1.0, 1.0);

        let mut unfocused = PhaseVocoder::new(1024, 256, 1.0, sample_rate, 0.0);
        configure(&mut unfocused);
        unfocused.advance_phases(num_bins, 1.0, 1.0);

        let independent_phase = focused.expected_phase_advance[target_bin] as f32;
        assert!(
            (focused.new_phases[target_bin] - independent_phase).abs() < 1e-6,
            "transient focus should keep non-peak bins on their independently advanced phase instead of reapplying the gradient field"
        );
        assert!(
            (unfocused.new_phases[target_bin] - independent_phase).abs() > 1e-3,
            "without transient focus the same non-peak bin should still be pulled by phase-gradient integration"
        );
    }

    #[test]
    fn test_phase_vocoder_identity() {
        // Stretch ratio 1.0 should approximately preserve the signal
        let sample_rate = 44100;
        let fft_size = 4096;
        let hop = fft_size / 4;

        // Generate a 440 Hz sine wave
        let num_samples = fft_size * 4;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 120.0);
        let output = pv.process(&input).unwrap();

        // Output length should be approximately the same
        let len_ratio = output.len() as f64 / input.len() as f64;
        assert!(
            (len_ratio - 1.0).abs() < 0.1,
            "Length ratio {} too far from 1.0",
            len_ratio
        );

        // Check that the output contains a similar frequency
        // (RMS should be similar)
        let input_rms: f32 = (input.iter().map(|x| x * x).sum::<f32>() / input.len() as f32).sqrt();
        let output_rms: f32 =
            (output.iter().map(|x| x * x).sum::<f32>() / output.len() as f32).sqrt();

        assert!(
            (output_rms - input_rms).abs() < input_rms * 0.5,
            "RMS mismatch: input={}, output={}",
            input_rms,
            output_rms
        );
    }

    #[test]
    fn test_phase_vocoder_stretch() {
        let sample_rate = 44100;
        let fft_size = 4096;
        let hop = fft_size / 4;

        // Use a longer signal for more accurate length ratio
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let stretch_ratio = 2.0;
        let mut pv = PhaseVocoder::new(fft_size, hop, stretch_ratio, sample_rate, 120.0);
        let output = pv.process(&input).unwrap();

        // Output should be approximately 2x longer (with tolerance for edge effects)
        let len_ratio = output.len() as f64 / input.len() as f64;
        assert!(
            (len_ratio - stretch_ratio).abs() < 0.35,
            "Length ratio {} too far from {}",
            len_ratio,
            stretch_ratio
        );
    }

    #[test]
    fn test_phase_vocoder_compress() {
        let sample_rate = 44100;
        let fft_size = 4096;
        let hop = fft_size / 4;

        let num_samples = fft_size * 4;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let stretch_ratio = 0.5;
        let mut pv = PhaseVocoder::new(fft_size, hop, stretch_ratio, sample_rate, 120.0);
        let output = pv.process(&input).unwrap();

        let len_ratio = output.len() as f64 / input.len() as f64;
        assert!(
            (len_ratio - stretch_ratio).abs() < 0.2,
            "Length ratio {} too far from {}",
            len_ratio,
            stretch_ratio
        );
    }

    #[test]
    fn test_phase_vocoder_input_too_short() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        let result = pv.process(&[0.0; 100]);
        assert!(result.is_err());
    }

    #[test]
    fn test_sub_bass_bin_calculation() {
        // 120 Hz cutoff at 44100 Hz with FFT size 4096
        // Expected bin: 120 * 4096 / 44100 ≈ 11.15 → 11
        let pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        assert_eq!(pv.sub_bass_bin, 11);

        // 0 Hz cutoff should give bin 0 (no sub-bass locking)
        let pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 0.0);
        assert_eq!(pv.sub_bass_bin, 0);

        // High cutoff at 48000 Hz
        let pv = PhaseVocoder::new(4096, 1024, 1.0, 48000, 200.0);
        let expected = (200.0f32 * 4096.0 / 48000.0).round() as usize;
        assert_eq!(pv.sub_bass_bin, expected);
    }

    #[test]
    fn test_sub_bass_phase_locking_preserves_low_freq() {
        // A 60 Hz sine should be handled by sub-bass rigid phase locking.
        // Compare output quality with sub-bass locking (120 Hz cutoff)
        // vs without (0 Hz cutoff).
        let sample_rate = 44100u32;
        let fft_size = 4096;
        let hop = fft_size / 4;
        let num_samples = fft_size * 8;
        let freq = 60.0f32; // Well below 120 Hz cutoff

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * freq * i as f32 / sample_rate as f32).sin())
            .collect();

        // Process with sub-bass locking enabled (120 Hz cutoff)
        let mut pv_locked = PhaseVocoder::new(fft_size, hop, 1.5, sample_rate, 120.0);
        let output_locked = pv_locked.process(&input).unwrap();

        // Process without sub-bass locking (0 Hz cutoff)
        let mut pv_unlocked = PhaseVocoder::new(fft_size, hop, 1.5, sample_rate, 0.0);
        let output_unlocked = pv_unlocked.process(&input).unwrap();

        // Both should produce output
        assert!(!output_locked.is_empty());
        assert!(!output_unlocked.is_empty());

        // Both should have similar RMS (we aren't destroying energy)
        let rms_locked =
            (output_locked.iter().map(|x| x * x).sum::<f32>() / output_locked.len() as f32).sqrt();
        let rms_unlocked = (output_unlocked.iter().map(|x| x * x).sum::<f32>()
            / output_unlocked.len() as f32)
            .sqrt();

        assert!(
            rms_locked > 0.1,
            "Sub-bass locked output should have significant energy, got RMS={}",
            rms_locked
        );
        assert!(
            rms_unlocked > 0.1,
            "Unlocked output should have significant energy, got RMS={}",
            rms_unlocked
        );
    }

    #[test]
    fn test_sub_bass_locking_does_not_affect_high_freq() {
        // A 1000 Hz sine should NOT be affected by sub-bass phase locking
        // (it's above the 120 Hz cutoff).
        let sample_rate = 44100u32;
        let fft_size = 4096;
        let hop = fft_size / 4;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 1000.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv_with = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 120.0);
        let output_with = pv_with.process(&input).unwrap();

        let mut pv_without = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 0.0);
        let output_without = pv_without.process(&input).unwrap();

        // Output lengths should be the same
        assert_eq!(output_with.len(), output_without.len());

        // RMS should be very similar since 1000 Hz is above the cutoff
        let rms_with =
            (output_with.iter().map(|x| x * x).sum::<f32>() / output_with.len() as f32).sqrt();
        let rms_without = (output_without.iter().map(|x| x * x).sum::<f32>()
            / output_without.len() as f32)
            .sqrt();

        assert!(
            (rms_with - rms_without).abs() < rms_with * 0.3,
            "1000 Hz signal should be similar with/without sub-bass locking: {} vs {}",
            rms_with,
            rms_without
        );
    }

    #[test]
    fn test_phase_vocoder_with_blackman_harris() {
        let sample_rate = 44100;
        let fft_size = 4096;
        let hop = fft_size / 4;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv = PhaseVocoder::with_window(
            fft_size,
            hop,
            1.5,
            sample_rate,
            120.0,
            WindowType::BlackmanHarris,
        );
        let output = pv.process(&input).unwrap();

        // Should produce valid stretched output
        assert!(!output.is_empty());
        let len_ratio = output.len() as f64 / input.len() as f64;
        assert!(
            (len_ratio - 1.5).abs() < 0.3,
            "BH window length ratio {} too far from 1.5",
            len_ratio
        );
    }

    #[test]
    fn test_phase_vocoder_with_kaiser() {
        let sample_rate = 44100;
        let fft_size = 4096;
        let hop = fft_size / 4;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv = PhaseVocoder::with_window(
            fft_size,
            hop,
            1.5,
            sample_rate,
            120.0,
            WindowType::Kaiser(800),
        );
        let output = pv.process(&input).unwrap();

        assert!(!output.is_empty());
        let len_ratio = output.len() as f64 / input.len() as f64;
        assert!(
            (len_ratio - 1.5).abs() < 0.3,
            "Kaiser window length ratio {} too far from 1.5",
            len_ratio
        );
    }

    #[test]
    fn test_phase_vocoder_different_windows_produce_different_output() {
        let sample_rate = 44100;
        let fft_size = 4096;
        let hop = fft_size / 4;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv_hann =
            PhaseVocoder::with_window(fft_size, hop, 1.5, sample_rate, 120.0, WindowType::Hann);
        let output_hann = pv_hann.process(&input).unwrap();

        let mut pv_bh = PhaseVocoder::with_window(
            fft_size,
            hop,
            1.5,
            sample_rate,
            120.0,
            WindowType::BlackmanHarris,
        );
        let output_bh = pv_bh.process(&input).unwrap();

        // Both should produce valid output of similar length
        assert!(!output_hann.is_empty());
        assert!(!output_bh.is_empty());

        // Outputs should differ (different windows produce different spectral characteristics)
        let min_len = output_hann.len().min(output_bh.len());
        let diff: f32 = output_hann[..min_len]
            .iter()
            .zip(&output_bh[..min_len])
            .map(|(a, b)| (a - b).abs())
            .sum::<f32>()
            / min_len as f32;
        assert!(
            diff > 1e-6,
            "Different windows should produce different output, avg diff = {}",
            diff
        );
    }

    // --- phase locking integration (detailed tests in phase_locking module) ---

    #[test]
    fn test_phase_lock_identity_no_peaks() {
        // Flat magnitude spectrum: no local maxima → no peaks → phases unchanged
        let num_bins = 16;
        let magnitudes = vec![1.0f32; num_bins]; // all equal, no peaks
        let analysis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32 * 0.1).collect();
        let mut synthesis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32 * 0.2).collect();
        let original_phases = synthesis_phases.clone();
        let mut peaks = Vec::new();

        apply_phase_locking(
            PhaseLockingMode::Identity,
            &magnitudes,
            &analysis_phases,
            &mut synthesis_phases,
            num_bins,
            0,
            &mut peaks,
        );

        // With no peaks found, phases should remain unchanged
        assert_eq!(synthesis_phases, original_phases);
    }

    #[test]
    fn test_phase_lock_identity_single_peak() {
        // Single peak at bin 5 with a realistic spectral lobe shape.
        // Trough-bounded identity locking propagates the peak's phase
        // rotation to all bins within its influence region (between troughs).
        let num_bins = 16;
        // Create a Gaussian-like lobe centered at bin 5, with troughs at 0 and 15
        let magnitudes: Vec<f32> = (0..num_bins)
            .map(|i| {
                let dist = (i as f32 - 5.0).abs();
                0.01 + 0.99 * (-dist * dist / 8.0).exp()
            })
            .collect();
        let analysis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32 * 0.3).collect();
        let mut synthesis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32 * 0.5).collect();
        let peak_synth = synthesis_phases[5];
        let mut peaks = Vec::new();

        apply_phase_locking(
            PhaseLockingMode::Identity,
            &magnitudes,
            &analysis_phases,
            &mut synthesis_phases,
            num_bins,
            0, // start_bin = 0
            &mut peaks,
        );

        // Peak at bin 5 should keep its phase
        assert!(
            (synthesis_phases[5] - peak_synth).abs() < 1e-6,
            "Peak bin should keep its phase"
        );

        // The phase rotation from the peak
        let phase_rotation = peak_synth - analysis_phases[5];

        // Bins in the peak's influence region should have:
        // synth[bin] = analysis[bin] + phase_rotation
        // With a single Gaussian lobe and no other peaks, all bins should
        // be in the peak's influence region.
        for bin in 1..num_bins - 1 {
            if bin == 5 {
                continue;
            }
            let expected = analysis_phases[bin] + phase_rotation;
            assert!(
                (synthesis_phases[bin] - expected).abs() < 1e-5,
                "Bin {} should be locked to peak: got {}, expected {}",
                bin,
                synthesis_phases[bin],
                expected
            );
        }
    }

    #[test]
    fn test_phase_lock_start_bin_above_num_bins() {
        // start_bin >= num_bins: early return, no changes
        let num_bins = 8;
        let magnitudes = vec![0.0f32; num_bins];
        let analysis_phases = vec![0.0f32; num_bins];
        let mut synthesis_phases = vec![1.0f32; num_bins];
        let original = synthesis_phases.clone();
        let mut peaks = Vec::new();

        apply_phase_locking(
            PhaseLockingMode::Identity,
            &magnitudes,
            &analysis_phases,
            &mut synthesis_phases,
            num_bins,
            num_bins, // start_bin == num_bins → early return
            &mut peaks,
        );

        assert_eq!(synthesis_phases, original);
    }

    #[test]
    fn test_phase_lock_num_bins_less_than_3() {
        // num_bins < 3: early return
        let magnitudes = vec![1.0f32; 2];
        let analysis_phases = vec![0.0f32; 2];
        let mut synthesis_phases = vec![0.5f32; 2];
        let original = synthesis_phases.clone();
        let mut peaks = Vec::new();

        apply_phase_locking(
            PhaseLockingMode::Identity,
            &magnitudes,
            &analysis_phases,
            &mut synthesis_phases,
            2,
            0,
            &mut peaks,
        );

        assert_eq!(synthesis_phases, original);
    }

    #[test]
    fn test_phase_lock_sub_bass_region_skipped() {
        // Peaks exist only below start_bin → no peaks found above sub-bass
        let num_bins = 16;
        let mut magnitudes = vec![0.1f32; num_bins];
        magnitudes[2] = 1.0; // peak below start_bin=5
        let analysis_phases = vec![0.0f32; num_bins];
        let mut synthesis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32).collect();
        let original = synthesis_phases.clone();
        let mut peaks = Vec::new();

        apply_phase_locking(
            PhaseLockingMode::Identity,
            &magnitudes,
            &analysis_phases,
            &mut synthesis_phases,
            num_bins,
            5, // start_bin=5, peak at bin 2 is below
            &mut peaks,
        );

        // No peaks above start_bin → no changes
        assert_eq!(synthesis_phases, original);
    }

    #[test]
    fn test_phase_lock_multiple_peaks() {
        // Two peaks with realistic spectral lobe shapes.
        // Trough-bounded identity locking assigns each bin to the peak
        // whose influence region (bounded by troughs) contains it.
        let num_bins = 16;
        // Create two Gaussian lobes: peak at bin 3, peak at bin 10
        // with a clear trough between them (around bin 7)
        let magnitudes: Vec<f32> = (0..num_bins)
            .map(|i| {
                let d3 = (i as f32 - 3.0).abs();
                let d10 = (i as f32 - 10.0).abs();
                let lobe3 = 1.0 * (-d3 * d3 / 4.0).exp();
                let lobe10 = 0.8 * (-d10 * d10 / 4.0).exp();
                0.001 + lobe3.max(lobe10) // ensure non-zero floor
            })
            .collect();
        let analysis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32 * 0.1).collect();
        let mut synthesis_phases: Vec<f32> = (0..num_bins).map(|i| i as f32 * 0.2).collect();
        let synth_peak3 = synthesis_phases[3];
        let synth_peak10 = synthesis_phases[10];
        let mut peaks = Vec::new();

        apply_phase_locking(
            PhaseLockingMode::Identity,
            &magnitudes,
            &analysis_phases,
            &mut synthesis_phases,
            num_bins,
            1, // start_bin=1
            &mut peaks,
        );

        // Verify both peaks are found
        assert!(peaks.contains(&3), "Should find peak at bin 3");
        assert!(peaks.contains(&10), "Should find peak at bin 10");

        // Phase rotation for each peak
        let rotation_3 = synth_peak3 - analysis_phases[3];
        let rotation_10 = synth_peak10 - analysis_phases[10];

        // Bin 2 is in peak 3's influence region (between start_bin boundary and trough)
        let expected_2 = analysis_phases[2] + rotation_3;
        assert!(
            (synthesis_phases[2] - expected_2).abs() < 1e-5,
            "Bin 2 should lock to peak 3: got {}, expected {}",
            synthesis_phases[2],
            expected_2
        );

        // Bin 12 is in peak 10's influence region
        let expected_12 = analysis_phases[12] + rotation_10;
        assert!(
            (synthesis_phases[12] - expected_12).abs() < 1e-5,
            "Bin 12 should lock to peak 10: got {}, expected {}",
            synthesis_phases[12],
            expected_12
        );
    }

    // --- normalize_output internals ---

    #[test]
    fn test_normalize_output_uniform_window_sum() {
        // When window_sum is uniform, output should be divided by that value
        let mut output = vec![2.0f32; 10];
        let window_sum = vec![2.0f32; 10];
        PhaseVocoder::normalize_output(&mut output, &window_sum, 1.0);
        for &s in &output {
            assert!((s - 1.0).abs() < 1e-6, "Expected 1.0, got {}", s);
        }
    }

    #[test]
    fn test_normalize_output_low_window_sum_clamped() {
        // Very small window sums should be clamped to min_window_sum
        // to prevent amplification
        let mut output = vec![1.0f32; 10];
        let mut window_sum = vec![1.0f32; 10];
        // One sample has near-zero window sum (low-overlap region)
        window_sum[5] = 1e-10;
        PhaseVocoder::normalize_output(&mut output, &window_sum, 1.0);

        // The clamped sample should NOT be amplified wildly
        // min_window_sum = max(1.0) * WINDOW_SUM_FLOOR_RATIO = 0.1
        // So output[5] = 1.0 / 0.1 = 10.0
        assert!(
            output[5] <= 11.0,
            "Low window sum should be clamped, got {}",
            output[5]
        );
        // Normal samples should be ~1.0
        assert!((output[0] - 1.0).abs() < 1e-6);
    }

    #[test]
    fn test_normalize_output_all_zero_window_sum() {
        // All-zero window sum: should use WINDOW_SUM_EPSILON floor
        let mut output = vec![1.0f32; 5];
        let window_sum = vec![0.0f32; 5];
        PhaseVocoder::normalize_output(&mut output, &window_sum, 1.0);
        // Each sample = 1.0 / WINDOW_SUM_EPSILON
        for &s in &output {
            assert!(s.is_finite(), "Output should be finite, got {}", s);
        }
    }

    // --- wrap_phase edge cases ---

    #[test]
    fn test_wrap_phase_exact_boundaries() {
        // Exactly PI should wrap to -PI (or very close)
        let result = wrap_phase(PI);
        assert!(
            (result - (-PI)).abs() < 1e-5 || (result - PI).abs() < 1e-5,
            "wrap_phase(PI) = {} should be near ±PI",
            result
        );

        // Exactly -PI
        let result = wrap_phase(-PI);
        assert!(
            (result - (-PI)).abs() < 1e-5 || (result - PI).abs() < 1e-5,
            "wrap_phase(-PI) = {} should be near ±PI",
            result
        );

        // Exactly 0
        assert!((wrap_phase(0.0)).abs() < 1e-6);
    }

    #[test]
    fn test_wrap_phase_very_large_values() {
        // Very large positive and negative values
        let result = wrap_phase(1000.0 * PI);
        assert!(
            (-PI..=PI).contains(&result),
            "wrap_phase(1000*PI) = {} should be in [-PI, PI]",
            result
        );

        let result = wrap_phase(-999.0 * PI);
        assert!(
            (-PI..=PI).contains(&result),
            "wrap_phase(-999*PI) = {} should be in [-PI, PI]",
            result
        );
    }

    // --- set_stretch_ratio ---

    #[test]
    fn test_set_stretch_ratio_updates_hop_synthesis() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        assert_eq!(pv.hop_synthesis(), 1024); // 1024 * 1.0 = 1024

        pv.set_stretch_ratio(2.0);
        assert_eq!(pv.hop_synthesis(), 2048); // 1024 * 2.0 = 2048

        pv.set_stretch_ratio(0.5);
        assert_eq!(pv.hop_synthesis(), 512); // 1024 * 0.5 = 512
    }

    #[test]
    fn test_set_stretch_ratio_preserves_phase_state() {
        // Process some audio, then change ratio and process more.
        // Phase should be continuous (no reset).
        let fft_size = 4096;
        let hop = 1024;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 120.0);
        let output1 = pv.process(&input).unwrap();
        assert!(!output1.is_empty());

        // Change ratio and process again — should work without error
        pv.set_stretch_ratio(1.5);
        let output2 = pv.process(&input).unwrap();
        assert!(!output2.is_empty());
        assert!(output2.len() > output1.len()); // 1.5x should be longer
    }

    #[test]
    fn test_set_stretch_ratio_engages_brief_continuity_focus_for_meaningful_step() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        assert_eq!(pv.transient_focus_frames, 0);

        pv.set_stretch_ratio(1.002);
        assert_eq!(
            pv.transient_focus_frames, RATIO_CHANGE_FOCUS_FRAMES,
            "meaningful runtime ratio steps should briefly tighten post-change phase locking"
        );
    }

    #[test]
    fn test_set_stretch_ratio_ignores_tiny_steps_for_continuity_focus() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        assert_eq!(pv.transient_focus_frames, 0);

        pv.set_stretch_ratio(1.0005);
        assert_eq!(
            pv.transient_focus_frames, 0,
            "sub-threshold ratio nudges should not burn the continuity-focus window"
        );
    }

    #[test]
    fn test_set_stretch_ratio_refreshes_continuity_focus_during_repeated_short_interval_modulation()
    {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        pv.set_stretch_ratio(1.002);
        assert_eq!(pv.transient_focus_frames, RATIO_CHANGE_FOCUS_FRAMES);

        pv.transient_focus_frames = 1;
        pv.set_stretch_ratio(1.004);
        assert_eq!(
            pv.transient_focus_frames, RATIO_CHANGE_FOCUS_FRAMES,
            "follow-up ratio steps should refresh continuity focus so repeated short-interval modulation keeps seam locking tight"
        );
    }

    #[test]
    fn test_set_stretch_ratio_restarts_phase_ratio_slew_from_current_effective_ratio() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);

        pv.set_stretch_ratio(1.12);
        assert!(
            (pv.continuity_focus_phase_ratio(pv.stretch_ratio) - 1.03).abs() < 1e-12,
            "the first continuity-focus frame should only move partway toward the new ratio"
        );

        pv.decay_transient_focus();
        assert!(
            (pv.continuity_focus_phase_ratio(pv.stretch_ratio) - 1.06).abs() < 1e-12,
            "the effective phase ratio should keep gliding toward the target while the seam window drains"
        );

        pv.set_stretch_ratio(0.92);
        assert!(
            (pv.ratio_change_phase_from - 1.06).abs() < 1e-12,
            "a repeated short-interval modulation step should restart from the in-flight effective phase ratio, not the stale previous target"
        );
        assert!(
            (pv.continuity_focus_phase_ratio(pv.stretch_ratio) - 1.032).abs() < 1e-12,
            "the restarted seam slew should still glide from the current effective ratio into the new target, but a direction reversal now keeps the seam on a slightly longer continuity window"
        );
    }

    #[test]
    fn test_set_stretch_ratio_extends_continuity_focus_for_direction_reversal() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);

        pv.set_stretch_ratio(1.12);
        pv.decay_transient_focus();
        assert!(
            (pv.continuity_focus_phase_ratio(pv.stretch_ratio) - 1.06).abs() < 1e-12,
            "test setup should leave an in-flight continuity slew before reversing direction"
        );

        pv.set_stretch_ratio(0.92);
        assert_eq!(
            pv.transient_focus_frames,
            RATIO_CHANGE_FOCUS_FRAMES + RATIO_CHANGE_REVERSAL_FOCUS_EXTRA_FRAMES,
            "direction-reversing ratio steps should hold continuity focus for one extra frame so the seam does not relax too early"
        );
        assert_eq!(
            pv.ratio_change_phase_total_frames,
            RATIO_CHANGE_FOCUS_FRAMES + RATIO_CHANGE_REVERSAL_FOCUS_EXTRA_FRAMES,
            "the phase-ratio slew should use the same extended window as transient focus during a reversal"
        );
    }

    #[test]
    fn test_set_stretch_ratio_reengages_continuity_focus_after_window_drains() {
        let mut pv = PhaseVocoder::new(4096, 1024, 1.0, 44100, 120.0);
        pv.set_stretch_ratio(1.002);
        pv.transient_focus_frames = 0;

        pv.set_stretch_ratio(1.004);
        assert_eq!(
            pv.transient_focus_frames, RATIO_CHANGE_FOCUS_FRAMES,
            "a new ratio step should re-engage continuity focus once the prior seam window has drained"
        );
    }

    #[test]
    fn test_set_stretch_ratio_reengages_focus_for_small_nudge_when_carried_seam_is_still_far() {
        let mut pv = PhaseVocoder::new(1024, 256, 1.0002, 44_100, 120.0);
        pv.streaming_tail = vec![0.0; 64];
        pv.streaming_tail_phase_ratio = 1.04;
        pv.ratio_change_phase_from = 1.04;
        pv.ratio_change_phase_total_frames = RATIO_CHANGE_FOCUS_FRAMES;
        pv.ratio_change_phase_frames = 0;
        pv.transient_focus_frames = 0;

        pv.set_stretch_ratio(1.0006);

        assert!(
            pv.transient_focus_frames >= RATIO_CHANGE_FOCUS_FRAMES,
            "a near-unity nudge should still re-engage continuity focus when the unresolved carried seam remains meaningfully expanded"
        );
        assert!(
            (pv.ratio_change_phase_from - 1.04).abs() < 1e-12,
            "the refreshed continuity slew should restart from the unresolved carried seam instead of the tiny prior target delta"
        );
    }

    #[test]
    fn test_set_stretch_ratio_extends_continuity_focus_while_streaming_tail_is_unresolved() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44_100u32;
        let input: Vec<f32> = (0..fft_size * 5)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.04, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain unresolved overlap tail"
        );

        pv.transient_focus_frames = 0;
        let expected = continuity_focus_frames_for_ratio_change(
            RATIO_CHANGE_FOCUS_FRAMES,
            pv.streaming_tail.len(),
            hop,
        );

        pv.set_stretch_ratio(0.96);

        assert_eq!(
            pv.transient_focus_frames, expected,
            "ratio changes with unresolved streaming overlap should keep continuity focus active until the seam has more time to drain"
        );
        assert!(
            pv.transient_focus_frames > RATIO_CHANGE_FOCUS_FRAMES,
            "streaming seam carry should extend continuity focus beyond the fixed minimum window"
        );
    }

    #[test]
    fn test_process_streaming_and_flush_produce_finite_output() {
        let fft_size = 2048;
        let hop = 512;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 10;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.6 + (2.0 * PI * 880.0 * t).sin() * 0.25
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.1, sample_rate, 120.0);
        let mut total_output = Vec::new();
        let mut analysis_buffer = Vec::new();

        for chunk in input.chunks(700) {
            analysis_buffer.extend_from_slice(chunk);
            if analysis_buffer.len() < fft_size {
                continue;
            }

            let out = pv.process_streaming(&analysis_buffer).unwrap();
            total_output.extend_from_slice(&out);

            let num_frames = (analysis_buffer.len() - fft_size) / hop + 1;
            let consumed = num_frames * hop;
            analysis_buffer.drain(..consumed);
        }

        if !analysis_buffer.is_empty() {
            analysis_buffer.resize(fft_size, 0.0);
            let out = pv.process_streaming(&analysis_buffer).unwrap();
            total_output.extend_from_slice(&out);
        }

        let tail = pv.flush_streaming().unwrap();
        total_output.extend_from_slice(&tail);

        assert!(
            !total_output.is_empty(),
            "Streaming path should produce output"
        );
        assert!(total_output.iter().all(|s| s.is_finite()));
    }

    #[test]
    fn test_streaming_tail_ratio_preserves_overlap_history_for_large_ratio_change() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 5;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.18, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.18).abs() < 1e-12,
            "tail ratio should match the ratio that generated the carried overlap"
        );

        let first_frames = (first_chunk_len - fft_size) / hop + 1;
        let consumed = first_frames * hop;
        pv.set_stretch_ratio(0.82);
        let second = pv
            .process_streaming(&input[consumed..consumed + first_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "second streaming chunk should emit audio"
        );
        assert!(
            !pv.streaming_tail.is_empty(),
            "second streaming chunk should continue carrying overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 0.82).abs() < 1e-12,
            "once the previous overlap has fully emitted, the carried tail should re-arm to the current ratio"
        );

        let _ = pv.flush_streaming().unwrap();
        assert!(
            (pv.streaming_tail_ratio - 0.82).abs() < 1e-12,
            "flush should clear carried overlap history and re-arm the current ratio"
        );
    }

    #[test]
    fn test_streaming_tail_phase_ratio_tracks_inflight_seam_after_prior_tail_drains() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 6;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.18, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let second_chunk_len = fft_size + hop * 3;
        pv.set_stretch_ratio(0.82);
        assert_eq!(
            pv.ratio_change_phase_total_frames,
            RATIO_CHANGE_FOCUS_FRAMES + RATIO_CHANGE_REVERSAL_FOCUS_EXTRA_FRAMES,
            "test setup should keep the reversed seam active across the short follow-up chunk"
        );
        let second = pv
            .process_streaming(&input[consumed..consumed + second_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "second streaming chunk should emit audio"
        );
        assert!(
            !pv.streaming_tail.is_empty(),
            "second streaming chunk should keep a newly generated overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 0.82).abs() < 1e-12,
            "once the older overlap has drained, normalization should re-arm to the new target ratio"
        );
        assert!(
            pv.streaming_tail_phase_ratio.is_finite(),
            "the carried phase seam ratio should remain finite"
        );
        assert!(
            pv.streaming_tail_phase_ratio <= pv.ratio_change_phase_from.max(pv.stretch_ratio)
                && pv.streaming_tail_phase_ratio
                    >= pv.ratio_change_phase_from.min(pv.stretch_ratio),
            "the carried phase seam ratio should stay within the active slew span"
        );
    }

    #[test]
    fn test_streaming_tail_phase_ratio_uses_first_overlap_crossing_callback_boundary() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44_100u32;
        let num_frames = 6usize;
        let input_len = fft_size + hop * (num_frames - 1);
        let input: Vec<f32> = (0..input_len)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut expected = PhaseVocoder::new(fft_size, hop, 1.12, sample_rate, 120.0);
        expected.set_stretch_ratio(0.82);
        assert_eq!(
            expected.ratio_change_phase_total_frames, RATIO_CHANGE_FOCUS_FRAMES,
            "test setup should leave the base continuity glide active across the callback"
        );
        let frame_advance = expected.hop_analysis as f64 * expected.stretch_ratio;
        let mut synthesis_frame_pos = 0.0f64;
        for _ in 0..num_frames {
            let synthesis_pos = snap_near_integer(synthesis_frame_pos);
            expected.decay_transient_focus();
            synthesis_frame_pos = synthesis_pos + frame_advance;
        }
        let emit_len = snap_near_integer(synthesis_frame_pos).floor() as usize;
        synthesis_frame_pos = 0.0;
        expected = PhaseVocoder::new(fft_size, hop, 1.12, sample_rate, 120.0);
        expected.set_stretch_ratio(0.82);
        let mut expected_carried_tail_phase_ratio = None;
        for _ in 0..num_frames {
            let synthesis_pos = snap_near_integer(synthesis_frame_pos);
            let synthesis_floor = synthesis_pos.floor() as usize;
            let frac = synthesis_pos - synthesis_floor as f64;
            let frame_end = synthesis_floor.saturating_add(
                fft_size
                    .saturating_sub(1)
                    .saturating_add(usize::from(frac > SYNTH_POS_EPSILON)),
            );
            let phase_ratio = expected.continuity_focus_phase_ratio(expected.stretch_ratio);
            if expected_carried_tail_phase_ratio.is_none() && frame_end >= emit_len {
                expected_carried_tail_phase_ratio = Some(phase_ratio);
            }
            expected.decay_transient_focus();
            synthesis_frame_pos = synthesis_pos + frame_advance;
        }
        let expected_carried_tail_phase_ratio = expected_carried_tail_phase_ratio
            .expect("expected an unresolved overlap past the callback boundary");

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.12, sample_rate, 120.0);
        pv.set_stretch_ratio(0.82);
        let output = pv.process_streaming(&input).unwrap();
        assert!(
            !output.is_empty(),
            "streaming pass should emit finalized samples"
        );
        assert!(
            !pv.streaming_tail.is_empty(),
            "streaming pass should retain unresolved overlap beyond the callback boundary"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - expected_carried_tail_phase_ratio).abs() < 1e-12,
            "the carried tail should preserve the phase ratio from the first overlap crossing the callback boundary, not the newest frame in the callback"
        );
        assert!(
            pv.streaming_tail_phase_ratio > pv.stretch_ratio + 1e-3,
            "the earliest unresolved seam should stay meaningfully above the settled target during the in-flight glide"
        );
    }

    #[test]
    fn test_streaming_tail_ratio_preserves_carried_overlap_for_small_cross_unity_modulation() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 5;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.04, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.04).abs() < 1e-12,
            "tail ratio should match the ratio that generated the carried overlap"
        );

        let first_frames = (first_chunk_len - fft_size) / hop + 1;
        let consumed = first_frames * hop;
        pv.set_stretch_ratio(0.96);
        let second = pv
            .process_streaming(&input[consumed..consumed + first_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "second streaming chunk should emit audio"
        );
        assert!(
            !pv.streaming_tail.is_empty(),
            "second streaming chunk should continue carrying overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 0.96).abs() < 1e-12,
            "once the previous overlap has fully emitted, small cross-unity modulation should re-arm to the current ratio"
        );

        let _ = pv.flush_streaming().unwrap();
        assert!(
            (pv.streaming_tail_ratio - 0.96).abs() < 1e-12,
            "flush should clear carried overlap history and re-arm the current ratio"
        );
    }

    #[test]
    fn test_streaming_tail_ratio_only_preserves_prior_ratio_while_overlap_remains_unresolved() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 5;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.18, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.18).abs() < 1e-12,
            "tail ratio should match the ratio that generated the carried overlap"
        );

        let first_frames = (first_chunk_len - fft_size) / hop + 1;
        let consumed = first_frames * hop;
        let second_chunk_len = fft_size + hop;
        pv.set_stretch_ratio(0.82);
        let second = pv
            .process_streaming(&input[consumed..consumed + second_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short second streaming chunk should still emit audio"
        );
        assert!(
            !pv.streaming_tail.is_empty(),
            "short second streaming chunk should continue carrying overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.18).abs() < 1e-12,
            "the prior expansion ratio should persist only while unresolved overlap from that chunk remains in the carried tail"
        );
    }

    #[test]
    fn test_streaming_tail_ratio_holds_prior_seam_across_repeated_short_interval_modulation() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.04, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.04).abs() < 1e-12,
            "tail ratio should match the ratio that generated the carried overlap"
        );

        let mut consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size + hop;

        pv.set_stretch_ratio(0.96);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "second short streaming chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.04).abs() < 1e-12,
            "the unresolved prior seam should keep its expansion ratio through the first cross-unity modulation step"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        pv.set_stretch_ratio(1.02);
        let third = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !third.is_empty(),
            "third short streaming chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.04).abs() < 1e-12,
            "repeated short-interval modulation should keep the unresolved prior seam ratio instead of re-arming on each toggle"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        pv.set_stretch_ratio(0.98);
        let fourth = pv
            .process_streaming(&input[consumed..consumed + first_chunk_len])
            .unwrap();
        assert!(
            !fourth.is_empty(),
            "a longer follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_ratio - 0.98).abs() < 1e-12,
            "once the older overlap has drained, the carried tail should re-arm to the current modulation ratio"
        );

        let _ = pv.flush_streaming().unwrap();
        assert!(
            (pv.streaming_tail_ratio - 0.98).abs() < 1e-12,
            "flush should preserve the re-armed current ratio after repeated short-interval modulation"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_carried_expansion_seam_when_rebounding_toward_unity() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.04, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let mut consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size + hop;

        pv.set_stretch_ratio(0.96);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_ratio - 1.04).abs() < 1e-12,
            "the prior expansion seam should still be carried after the short cross-unity step"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        assert!(
            ratio_is_meaningfully_below_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should decay the in-flight compression seam below unity before the rebound"
        );

        pv.set_stretch_ratio(1.02);
        assert!(
            (pv.ratio_change_phase_from - 1.04).abs() < 1e-12,
            "rebounding toward unity should restart from the carried expansion seam while that overlap tail is still unresolved"
        );

        let third = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !third.is_empty(),
            "the rebound chunk should still emit audio after re-anchoring the seam"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_carried_compression_seam_when_rebounding_toward_unity() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 0.96, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let mut consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size + hop;

        pv.set_stretch_ratio(1.12);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 0.96).abs() < 1e-12,
            "the prior compression seam should still be carried after the short cross-unity step"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        let hold = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !hold.is_empty(),
            "a second short chunk at the rebound ratio should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 0.96).abs() < 1e-12,
            "the carried compression seam should persist until the older overlap fully drains"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        assert!(
            ratio_is_meaningfully_above_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should decay the in-flight expansion seam above unity before the rebound"
        );

        pv.set_stretch_ratio(0.98);
        assert!(
            (pv.ratio_change_phase_from - 0.96).abs() < 1e-12,
            "rebounding toward unity should restart from the carried compression seam while that overlap tail is still unresolved"
        );

        let third = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !third.is_empty(),
            "the rebound chunk should still emit audio after re-anchoring the compression seam"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_same_side_carried_expansion_seam_when_rebounding_toward_unity()
    {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.12, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let mut consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size;

        pv.set_stretch_ratio(1.04);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 1.12).abs() < 1e-12,
            "the prior expansion seam should still be carried after the short same-side step"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        let hold = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !hold.is_empty(),
            "a second short chunk at the reduced expansion ratio should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 1.12).abs() < 1e-12,
            "the carried expansion seam should persist until the older overlap fully drains"
        );
        assert!(
            ratio_is_meaningfully_above_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight expansion seam above unity before the rebound"
        );

        pv.set_stretch_ratio(1.02);
        assert!(
            (pv.ratio_change_phase_from - 1.12).abs() < 1e-12,
            "same-side rebound toward unity should restart from the carried expansion seam while that overlap tail is still unresolved"
        );
        assert_eq!(
            pv.ratio_change_phase_total_frames,
            continuity_focus_frames_for_ratio_change(
                RATIO_CHANGE_FOCUS_FRAMES,
                pv.streaming_tail.len(),
                hop,
            ) + RATIO_CHANGE_CARRIED_SEAM_FOCUS_EXTRA_FRAMES,
            "re-anchoring to an older carried expansion seam should hold continuity focus one extra frame so the overlap does not relax before that seam drains"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_same_side_carried_compression_seam_when_rebounding_toward_unity(
    ) {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 0.88, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let mut consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size;

        pv.set_stretch_ratio(0.96);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 0.88).abs() < 1e-12,
            "the prior compression seam should still be carried after the short same-side step"
        );

        consumed += ((short_chunk_len - fft_size) / hop + 1) * hop;
        let hold = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !hold.is_empty(),
            "a second short chunk at the relaxed compression ratio should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 0.88).abs() < 1e-12,
            "the carried compression seam should persist until the older overlap fully drains"
        );
        assert!(
            ratio_is_meaningfully_below_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight compression seam below unity before the rebound"
        );

        pv.set_stretch_ratio(0.98);
        assert!(
            (pv.ratio_change_phase_from - 0.88).abs() < 1e-12,
            "same-side rebound toward unity should restart from the carried compression seam while that overlap tail is still unresolved"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_same_side_carried_expansion_seam_when_stepping_away_from_unity()
    {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.04, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size + hop;

        pv.set_stretch_ratio(1.02);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 1.04).abs() < 1e-12,
            "the prior expansion seam should still be carried after the short same-side relaxation"
        );
        assert!(
            ratio_is_meaningfully_above_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight expansion seam above unity before the stronger same-side step"
        );

        pv.set_stretch_ratio(1.08);
        assert!(
            (pv.ratio_change_phase_from - 1.04).abs() < 1e-12,
            "same-side step away from unity should restart from the carried expansion seam while that overlap tail is still unresolved"
        );
    }

    #[test]
    fn test_set_stretch_ratio_reanchors_to_same_side_carried_seam_for_subthreshold_nudge() {
        let mut pv = PhaseVocoder::new(1024, 256, 1.02, 44_100, 120.0);
        pv.streaming_tail = vec![0.0; 64];
        pv.streaming_tail_phase_ratio = 1.12;
        pv.ratio_change_phase_from = 1.12;
        pv.ratio_change_phase_total_frames = 4;
        pv.ratio_change_phase_frames = 1;
        pv.transient_focus_frames = 1;
        assert!(
            ratio_is_meaningfully_above_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight seam on the same side of unity before the tiny follow-up step"
        );

        pv.set_stretch_ratio(1.0205);
        assert!(
            (pv.ratio_change_phase_from - 1.12).abs() < 1e-12,
            "a sub-threshold same-side follow-up nudge should restart from the older carried expansion seam while that overlap tail is still unresolved"
        );
        assert!(
            (pv.continuity_focus_phase_ratio(pv.stretch_ratio) - 1.12).abs() < 1e-12,
            "a far carried-seam re-anchor should hold the first continuity-focus frame on the older seam before gliding toward the new target"
        );
        pv.decay_transient_focus();
        assert!(
            pv.continuity_focus_phase_ratio(pv.stretch_ratio) < 1.12,
            "after the initial seam hold, the restarted continuity slew should resume gliding toward the new target"
        );
        assert_eq!(
            pv.ratio_change_phase_total_frames,
            continuity_focus_frames_for_ratio_change(
                RATIO_CHANGE_FOCUS_FRAMES,
                pv.streaming_tail.len(),
                pv.hop_analysis,
            ) + RATIO_CHANGE_CARRIED_SEAM_FOCUS_EXTRA_FRAMES,
            "re-anchoring a tiny same-side nudge to an older carried seam should refresh continuity focus so the overlap does not relax mid-tail"
        );
    }

    #[test]
    fn test_set_stretch_ratio_reanchors_to_opposite_side_carried_seam_for_subthreshold_deeper_step()
    {
        let mut pv = PhaseVocoder::new(1024, 256, 0.92, 44_100, 120.0);
        pv.streaming_tail = vec![0.0; 64];
        pv.streaming_tail_phase_ratio = 1.12;
        pv.ratio_change_phase_from = 1.12;
        pv.ratio_change_phase_total_frames = 4;
        pv.ratio_change_phase_frames = 1;
        pv.transient_focus_frames = 1;
        assert!(
            ratio_is_meaningfully_below_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight seam on the compression side before the tiny deeper step"
        );

        pv.set_stretch_ratio(0.9195);
        assert!(
            (pv.ratio_change_phase_from - 1.12).abs() < 1e-12,
            "a sub-threshold deeper compression nudge should restart from the older carried expansion seam while that overlap tail is still unresolved"
        );
        assert_eq!(
            pv.ratio_change_phase_total_frames,
            continuity_focus_frames_for_ratio_change(
                RATIO_CHANGE_FOCUS_FRAMES,
                pv.streaming_tail.len(),
                pv.hop_analysis,
            ) + RATIO_CHANGE_CARRIED_SEAM_FOCUS_EXTRA_FRAMES,
            "re-anchoring a tiny opposite-side follow-up step should refresh continuity focus so the older seam does not snap mid-tail"
        );
    }

    #[test]
    fn test_set_stretch_ratio_stacks_reversal_and_carried_seam_focus_extensions() {
        let mut pv = PhaseVocoder::new(1024, 256, 1.12, 44_100, 120.0);
        pv.streaming_tail = vec![0.0; 64];
        pv.streaming_tail_phase_ratio = 0.92;
        pv.ratio_change_phase_from = 0.92;
        pv.ratio_change_phase_total_frames = 4;
        pv.ratio_change_phase_frames = 2;
        pv.transient_focus_frames = 2;

        assert!(
            ratio_is_meaningfully_above_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight seam above unity before reversing back toward the carried compression seam"
        );

        pv.set_stretch_ratio(0.86);
        assert!(
            (pv.ratio_change_phase_from - 0.92).abs() < 1e-12,
            "the restarted continuity slew should re-anchor to the older carried compression seam"
        );
        assert_eq!(
            pv.ratio_change_phase_total_frames,
            continuity_focus_frames_for_ratio_change(
                RATIO_CHANGE_FOCUS_FRAMES,
                pv.streaming_tail.len(),
                pv.hop_analysis,
            ) + RATIO_CHANGE_REVERSAL_FOCUS_EXTRA_FRAMES
                + RATIO_CHANGE_CARRIED_SEAM_FOCUS_EXTRA_FRAMES,
            "when a new step both reverses direction and re-anchors to a far carried seam, both continuity-focus extensions should stack"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_same_side_carried_compression_seam_when_stepping_away_from_unity(
    ) {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 8;
        let input: Vec<f32> = (0..num_samples)
            .map(|i| {
                let t = i as f32 / sample_rate as f32;
                (2.0 * PI * 220.0 * t).sin() * 0.55 + (2.0 * PI * 660.0 * t).sin() * 0.20
            })
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 0.96, sample_rate, 120.0);
        let first_chunk_len = fft_size * 2;
        let first = pv.process_streaming(&input[..first_chunk_len]).unwrap();
        assert!(!first.is_empty(), "first streaming chunk should emit audio");
        assert!(
            !pv.streaming_tail.is_empty(),
            "first streaming chunk should retain overlap tail"
        );

        let consumed = ((first_chunk_len - fft_size) / hop + 1) * hop;
        let short_chunk_len = fft_size;

        pv.set_stretch_ratio(0.98);
        let second = pv
            .process_streaming(&input[consumed..consumed + short_chunk_len])
            .unwrap();
        assert!(
            !second.is_empty(),
            "short follow-up chunk should still emit audio"
        );
        assert!(
            (pv.streaming_tail_phase_ratio - 0.96).abs() < 1e-12,
            "the prior compression seam should still be carried after the short same-side relaxation"
        );
        assert!(
            ratio_is_meaningfully_below_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight compression seam below unity before the stronger same-side step"
        );

        pv.set_stretch_ratio(0.92);
        assert!(
            (pv.ratio_change_phase_from - 0.96).abs() < 1e-12,
            "same-side step away from unity should restart from the carried compression seam while that overlap tail is still unresolved"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_carried_expansion_seam_when_crossing_unity_and_stepping_deeper_into_compression(
    ) {
        let mut pv = PhaseVocoder::new(1024, 256, 0.98, 44_100, 120.0);
        pv.streaming_tail = vec![0.0; 64];
        pv.streaming_tail_phase_ratio = 1.04;
        pv.ratio_change_phase_from = 1.04;
        pv.ratio_change_phase_total_frames = RATIO_CHANGE_FOCUS_FRAMES;
        pv.ratio_change_phase_frames = 1;
        assert!(
            ratio_is_meaningfully_below_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight seam on the compression side before stepping deeper"
        );

        pv.set_stretch_ratio(0.92);
        assert!(
            (pv.ratio_change_phase_from - 1.04).abs() < 1e-12,
            "cross-unity step away from unity should restart from the older carried expansion seam while that overlap tail is still unresolved"
        );
    }

    #[test]
    fn test_set_stretch_ratio_keeps_carried_compression_seam_when_crossing_unity_and_stepping_deeper_into_expansion(
    ) {
        let mut pv = PhaseVocoder::new(1024, 256, 1.02, 44_100, 120.0);
        pv.streaming_tail = vec![0.0; 64];
        pv.streaming_tail_phase_ratio = 0.96;
        pv.ratio_change_phase_from = 0.96;
        pv.ratio_change_phase_total_frames = RATIO_CHANGE_FOCUS_FRAMES;
        pv.ratio_change_phase_frames = 1;
        assert!(
            ratio_is_meaningfully_above_unity(pv.continuity_focus_phase_ratio(pv.stretch_ratio)),
            "test setup should leave the in-flight seam on the expansion side before stepping deeper"
        );

        pv.set_stretch_ratio(1.08);
        assert!(
            (pv.ratio_change_phase_from - 0.96).abs() < 1e-12,
            "cross-unity step away from unity should restart from the older carried compression seam while that overlap tail is still unresolved"
        );
    }

    #[test]
    fn test_flush_streaming_is_idempotent() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;
        let input: Vec<f32> = (0..fft_size * 3)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 120.0);
        let _ = pv.process_streaming(&input).unwrap();

        let _first = pv.flush_streaming().unwrap();
        let second = pv.flush_streaming().unwrap();
        assert!(
            second.is_empty(),
            "Second flush_streaming() call should be empty"
        );
    }

    // --- sub_bass_bin edge cases ---

    #[test]
    fn test_sub_bass_bin_clamped_to_num_bins() {
        // Very high cutoff: sub_bass_bin should be clamped to num_bins
        let pv = PhaseVocoder::new(256, 64, 1.0, 44100, 30000.0);
        let num_bins = 256 / 2 + 1;
        assert!(
            pv.sub_bass_bin() <= num_bins,
            "sub_bass_bin {} should be <= num_bins {}",
            pv.sub_bass_bin(),
            num_bins
        );
    }

    #[test]
    fn test_sub_bass_all_bins_rigid() {
        // With cutoff >= Nyquist, all bins should use rigid locking.
        // This should still produce valid output (no crash).
        let fft_size = 512;
        let hop = 128;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        // Cutoff at Nyquist: all bins are "sub-bass" → all rigid locking
        let mut pv = PhaseVocoder::new(fft_size, hop, 1.5, sample_rate, 22050.0);
        let output = pv.process(&input).unwrap();
        assert!(!output.is_empty());
        assert!(output.iter().all(|s| s.is_finite()));
    }

    // --- reconstruct_spectrum conjugate symmetry ---

    #[test]
    fn test_reconstruct_spectrum_produces_real_output() {
        // After reconstruct_spectrum + inverse FFT, output should be real-valued
        // (imaginary parts near zero). This verifies conjugate symmetry is correct.
        let fft_size = 256;
        let hop = 64;
        let sample_rate = 44100u32;
        let num_samples = fft_size * 4;

        let input: Vec<f32> = (0..num_samples)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 120.0);
        let output = pv.process(&input).unwrap();

        // If conjugate symmetry is wrong, we'd get complex residues causing
        // large imaginary parts. The output being finite and reasonable is evidence.
        assert!(output.iter().all(|s| s.is_finite()));
        let rms = (output.iter().map(|x| x * x).sum::<f32>() / output.len() as f32).sqrt();
        assert!(
            rms > 0.01,
            "Output should have significant energy, got RMS={}",
            rms
        );
    }

    // --- PV reuse (buffers grow but don't shrink) ---

    #[test]
    fn test_phase_vocoder_reuse_across_different_lengths() {
        let fft_size = 1024;
        let hop = 256;
        let sample_rate = 44100u32;

        let mut pv = PhaseVocoder::new(fft_size, hop, 1.0, sample_rate, 120.0);

        // Process a long signal
        let long_input: Vec<f32> = (0..fft_size * 8)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();
        let output1 = pv.process(&long_input).unwrap();
        assert!(!output1.is_empty());

        // Process a shorter signal — buffers should still work (they don't shrink)
        let short_input: Vec<f32> = (0..fft_size * 2)
            .map(|i| (2.0 * PI * 440.0 * i as f32 / sample_rate as f32).sin())
            .collect();
        let output2 = pv.process(&short_input).unwrap();
        assert!(!output2.is_empty());
        assert!(output2.len() < output1.len());
    }
}