brk_computer 0.11.0

use super::order_statistics::ExactOrderStats;

/// Sorted sliding window for rolling distribution/median computations.
///
/// Uses sqrt-decomposition for O(sqrt(n)) insert/remove/kth instead of
/// O(n) memmoves with a flat sorted Vec.
pub(crate) struct SlidingWindowSorted {
    sorted: ExactOrderStats,
    prev_start: usize,
}

impl SlidingWindowSorted {
    pub fn with_capacity(cap: usize) -> Self {
        Self {
            sorted: ExactOrderStats::new(cap),
            prev_start: 0,
        }
    }

    /// Reconstruct state from historical data (the elements in [range_start..skip]).
    /// Uses O(n log n) sort + O(n) block construction instead of O(n√n) individual inserts.
    pub fn reconstruct(&mut self, partial_values: &[f64], range_start: usize, skip: usize) {
        self.prev_start = range_start;
        let slice = &partial_values[..skip - range_start];
        if slice.is_empty() {
            return;
        }
        self.sorted = ExactOrderStats::from_unsorted(slice.to_vec());
    }

    /// Add a new value and remove all expired values up to `new_start`.
    pub fn advance(
        &mut self,
        value: f64,
        new_start: usize,
        partial_values: &[f64],
        range_start: usize,
    ) {
        self.sorted.insert(value);

        while self.prev_start < new_start {
            let old = partial_values[self.prev_start - range_start];
            self.sorted.remove(old);
            self.prev_start += 1;
        }
    }

    #[inline]
    pub fn is_empty(&self) -> bool {
        self.sorted.is_empty()
    }

    #[inline]
    pub fn min(&self) -> f64 {
        if self.sorted.is_empty() {
            0.0
        } else {
            self.sorted.first()
        }
    }

    #[inline]
    pub fn max(&self) -> f64 {
        if self.sorted.is_empty() {
            0.0
        } else {
            self.sorted.last()
        }
    }

    /// Extract a percentile (0.0-1.0) using linear interpolation.
    #[inline]
    pub fn percentile(&self, p: f64) -> f64 {
        let len = self.sorted.len();
        if len == 0 {
            return 0.0;
        }
        if len == 1 {
            return self.sorted.kth(0);
        }
        let rank = p * (len - 1) as f64;
        let lo = rank.floor() as usize;
        let hi = rank.ceil() as usize;
        if lo == hi {
            self.sorted.kth(lo)
        } else {
            let frac = rank - lo as f64;
            self.sorted.kth(lo) * (1.0 - frac) + self.sorted.kth(hi) * frac
        }
    }

    /// Extract multiple percentiles in a single pass through the sorted blocks.
    /// Percentiles must be sorted ascending. Returns interpolated values.
    pub fn percentiles<const N: usize>(&self, ps: &[f64; N]) -> [f64; N] {
        let len = self.sorted.len();
        if len == 0 {
            return [0.0; N];
        }
        if len == 1 {
            return [self.sorted.kth(0); N];
        }

        // Collect all unique ranks needed (lo and hi for each percentile)
        let last = (len - 1) as f64;
        let mut rank_set: [usize; 10] = [0; 10];
        let mut rank_count = 0;
        let mut lo_hi: [(usize, usize, f64); N] = [(0, 0, 0.0); N];

        for (i, &p) in ps.iter().enumerate() {
            let rank = p * last;
            let lo = rank.floor() as usize;
            let hi = rank.ceil() as usize;
            lo_hi[i] = (lo, hi, rank - lo as f64);

            rank_set[rank_count] = lo;
            rank_count += 1;
            if hi != lo {
                rank_set[rank_count] = hi;
                rank_count += 1;
            }
        }

        // Sort and deduplicate (interleaved lo/hi values aren't necessarily sorted)
        rank_set[..rank_count].sort_unstable();
        let mut w = 1;
        for r in 1..rank_count {
            if rank_set[r] != rank_set[w - 1] {
                rank_set[w] = rank_set[r];
                w += 1;
            }
        }
        rank_count = w;

        // Single pass through blocks to get all values
        let ranks = &rank_set[..rank_count];
        let mut values = [0.0f64; 10];
        self.sorted.values_at(ranks, &mut values[..rank_count]);

        // Interpolate results
        let mut out = [0.0; N];
        for (i, &(lo, hi, frac)) in lo_hi.iter().enumerate() {
            if lo == hi {
                let ri = ranks.partition_point(|&r| r < lo);
                out[i] = values[ri];
            } else {
                let lo_ri = ranks.partition_point(|&r| r < lo);
                let hi_ri = ranks.partition_point(|&r| r < hi);
                out[i] = values[lo_ri] * (1.0 - frac) + values[hi_ri] * frac;
            }
        }
        out
    }
}