vector_ta/indicators/moving_averages/

sama.rs

#[cfg(feature = "python")]
use numpy::PyUntypedArrayMethods;
#[cfg(feature = "python")]
use numpy::{IntoPyArray, PyArray1, PyArrayMethods, PyReadonlyArray1, PyReadonlyArray2};
#[cfg(feature = "python")]
use pyo3::exceptions::PyValueError;
#[cfg(feature = "python")]
use pyo3::prelude::*;
#[cfg(feature = "python")]
use pyo3::types::PyDict;

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
use serde::{Deserialize, Serialize};
#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
use wasm_bindgen::prelude::*;

#[cfg(all(feature = "python", feature = "cuda"))]
use crate::cuda::moving_averages::sama_wrapper::DeviceArrayF32Sama;
#[cfg(all(feature = "python", feature = "cuda"))]
use crate::cuda::{cuda_available, moving_averages::CudaSama};
use crate::utilities::data_loader::{source_type, Candles};
#[cfg(all(feature = "python", feature = "cuda"))]
use crate::utilities::dlpack_cuda::export_f32_cuda_dlpack_2d;
use crate::utilities::enums::Kernel;
use crate::utilities::helpers::{
    alloc_with_nan_prefix, detect_best_batch_kernel, detect_best_kernel, init_matrix_prefixes,
    make_uninit_matrix,
};
#[cfg(feature = "python")]
use crate::utilities::kernel_validation::validate_kernel;
#[cfg(all(feature = "python", feature = "cuda"))]
use cust::memory::DeviceBuffer;

#[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
use core::arch::x86_64::*;

#[cfg(not(target_arch = "wasm32"))]
use rayon::prelude::*;

use std::convert::AsRef;
use std::error::Error;
use std::mem::MaybeUninit;
use thiserror::Error;

impl<'a> AsRef<[f64]> for SamaInput<'a> {
    #[inline(always)]
    fn as_ref(&self) -> &[f64] {
        match &self.data {
            SamaData::Slice(slice) => slice,
            SamaData::Candles { candles, source } => source_type(candles, source),
        }
    }
}

#[derive(Debug, Clone)]
pub enum SamaData<'a> {
    Candles {
        candles: &'a Candles,
        source: &'a str,
    },
    Slice(&'a [f64]),
}

#[derive(Debug, Clone)]
pub struct SamaOutput {
    pub values: Vec<f64>,
}

#[derive(Debug, Clone)]
#[cfg_attr(
    all(target_arch = "wasm32", feature = "wasm"),
    derive(Serialize, Deserialize)
)]
pub struct SamaParams {
    pub length: Option<usize>,
    pub maj_length: Option<usize>,
    pub min_length: Option<usize>,
}

impl Default for SamaParams {
    fn default() -> Self {
        Self {
            length: Some(200),
            maj_length: Some(14),
            min_length: Some(6),
        }
    }
}

#[derive(Debug, Clone)]
pub struct SamaInput<'a> {
    pub data: SamaData<'a>,
    pub params: SamaParams,
}

impl<'a> SamaInput<'a> {
    #[inline]
    pub fn from_candles(c: &'a Candles, s: &'a str, p: SamaParams) -> Self {
        Self {
            data: SamaData::Candles {
                candles: c,
                source: s,
            },
            params: p,
        }
    }

    #[inline]
    pub fn from_slice(sl: &'a [f64], p: SamaParams) -> Self {
        Self {
            data: SamaData::Slice(sl),
            params: p,
        }
    }

    #[inline]
    pub fn with_default_candles(c: &'a Candles) -> Self {
        Self::from_candles(c, "close", SamaParams::default())
    }

    #[inline]
    pub fn get_length(&self) -> usize {
        self.params.length.unwrap_or(200)
    }

    #[inline]
    pub fn get_maj_length(&self) -> usize {
        self.params.maj_length.unwrap_or(14)
    }

    #[inline]
    pub fn get_min_length(&self) -> usize {
        self.params.min_length.unwrap_or(6)
    }
}

#[derive(Clone, Debug)]
pub struct SamaBuilder {
    length: Option<usize>,
    maj_length: Option<usize>,
    min_length: Option<usize>,
    kernel: Kernel,
}

impl Default for SamaBuilder {
    fn default() -> Self {
        Self {
            length: None,
            maj_length: None,
            min_length: None,
            kernel: Kernel::Auto,
        }
    }
}

impl SamaBuilder {
    #[inline(always)]
    pub fn new() -> Self {
        Self::default()
    }

    #[inline(always)]
    pub fn length(mut self, val: usize) -> Self {
        self.length = Some(val);
        self
    }

    #[inline(always)]
    pub fn maj_length(mut self, val: usize) -> Self {
        self.maj_length = Some(val);
        self
    }

    #[inline(always)]
    pub fn min_length(mut self, val: usize) -> Self {
        self.min_length = Some(val);
        self
    }

    #[inline(always)]
    pub fn kernel(mut self, k: Kernel) -> Self {
        self.kernel = k;
        self
    }

    #[inline(always)]
    pub fn apply(self, c: &Candles) -> Result<SamaOutput, SamaError> {
        let p = SamaParams {
            length: self.length,
            maj_length: self.maj_length,
            min_length: self.min_length,
        };
        let i = SamaInput::from_candles(c, "close", p);
        sama_with_kernel(&i, self.kernel)
    }

    #[inline(always)]
    pub fn apply_slice(self, d: &[f64]) -> Result<SamaOutput, SamaError> {
        let p = SamaParams {
            length: self.length,
            maj_length: self.maj_length,
            min_length: self.min_length,
        };
        let i = SamaInput::from_slice(d, p);
        sama_with_kernel(&i, self.kernel)
    }

    #[inline(always)]
    pub fn into_stream(self) -> Result<SamaStream, SamaError> {
        let p = SamaParams {
            length: self.length,
            maj_length: self.maj_length,
            min_length: self.min_length,
        };
        SamaStream::try_new(p)
    }
}

#[derive(Debug, Error)]
pub enum SamaError {
    #[error("sama: Input data slice is empty.")]
    EmptyInputData,

    #[error("sama: All values are NaN.")]
    AllValuesNaN,

    #[error("sama: Invalid period: period = {period}, data length = {data_len}")]
    InvalidPeriod { period: usize, data_len: usize },

    #[error("sama: Not enough valid data: needed = {needed}, valid = {valid}")]
    NotEnoughValidData { needed: usize, valid: usize },

    #[error("sama: Output length mismatch: expected = {expected}, got = {got}")]
    OutputLengthMismatch { expected: usize, got: usize },

    #[error("sama: Invalid range expansion: start = {start}, end = {end}, step = {step}")]
    InvalidRange {
        start: usize,
        end: usize,
        step: usize,
    },

    #[error("sama: Invalid kernel for batch path: {0:?}")]
    InvalidKernelForBatch(Kernel),
}

#[inline(always)]
pub fn sama(input: &SamaInput) -> Result<SamaOutput, SamaError> {
    sama_with_kernel(input, Kernel::Auto)
}

#[inline(always)]
pub fn sama_with_kernel(input: &SamaInput, kernel: Kernel) -> Result<SamaOutput, SamaError> {
    let (data, length, maj_length, min_length, first, chosen) = sama_prepare(input, kernel)?;

    let mut out = alloc_with_nan_prefix(data.len(), first);
    sama_compute_into(
        data, length, maj_length, min_length, first, chosen, &mut out,
    );
    Ok(SamaOutput { values: out })
}

#[cfg(not(all(target_arch = "wasm32", feature = "wasm")))]
pub fn sama_into(input: &SamaInput, out: &mut [f64]) -> Result<(), SamaError> {
    let (data, length, maj_length, min_length, first, chosen) = sama_prepare(input, Kernel::Auto)?;

    if out.len() != data.len() {
        return Err(SamaError::OutputLengthMismatch {
            expected: data.len(),
            got: out.len(),
        });
    }

    let warm = first.min(out.len());
    for i in 0..warm {
        out[i] = f64::from_bits(0x7ff8_0000_0000_0000);
    }

    sama_compute_into(data, length, maj_length, min_length, first, chosen, out);

    Ok(())
}

#[inline(always)]
pub fn sama_into_slice(dst: &mut [f64], input: &SamaInput, kern: Kernel) -> Result<(), SamaError> {
    let (data, length, maj_length, min_length, first, chosen) = sama_prepare(input, kern)?;

    if dst.len() != data.len() {
        return Err(SamaError::OutputLengthMismatch {
            expected: data.len(),
            got: dst.len(),
        });
    }

    sama_compute_into(data, length, maj_length, min_length, first, chosen, dst);

    for v in &mut dst[..first] {
        *v = f64::NAN;
    }
    Ok(())
}

#[inline(always)]
fn sama_prepare<'a>(
    input: &'a SamaInput,
    kernel: Kernel,
) -> Result<(&'a [f64], usize, usize, usize, usize, Kernel), SamaError> {
    let data: &[f64] = input.as_ref();
    let n = data.len();

    if n == 0 {
        return Err(SamaError::EmptyInputData);
    }

    let first = data
        .iter()
        .position(|x| !x.is_nan())
        .ok_or(SamaError::AllValuesNaN)?;

    let length = input.get_length();
    let maj_length = input.get_maj_length();
    let min_length = input.get_min_length();

    if length + 1 > n || length == 0 {
        return Err(SamaError::InvalidPeriod {
            period: length,
            data_len: n,
        });
    }

    if maj_length == 0 || min_length == 0 {
        return Err(SamaError::InvalidPeriod {
            period: 0,
            data_len: n,
        });
    }

    let valid = n - first;

    if valid < 1 {
        return Err(SamaError::NotEnoughValidData { needed: 1, valid });
    }

    let chosen = match kernel {
        Kernel::Auto => detect_best_kernel(),
        k => k,
    };

    Ok((data, length, maj_length, min_length, first, chosen))
}

fn sama_compute_into(
    data: &[f64],
    length: usize,
    maj_length: usize,
    min_length: usize,
    first: usize,
    kernel: Kernel,
    out: &mut [f64],
) {
    unsafe {
        #[cfg(all(target_arch = "wasm32", target_feature = "simd128"))]
        {
            if matches!(kernel, Kernel::Scalar | Kernel::ScalarBatch) {
                sama_simd128(data, length, maj_length, min_length, first, out);
                return;
            }
        }

        match kernel {
            Kernel::Scalar | Kernel::ScalarBatch => {
                sama_scalar(data, length, maj_length, min_length, first, out)
            }
            #[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
            Kernel::Avx2 | Kernel::Avx2Batch => {
                sama_avx2(data, length, maj_length, min_length, first, out)
            }
            #[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
            Kernel::Avx512 | Kernel::Avx512Batch => {
                sama_avx512(data, length, maj_length, min_length, first, out)
            }
            #[cfg(not(all(feature = "nightly-avx", target_arch = "x86_64")))]
            Kernel::Avx2 | Kernel::Avx2Batch | Kernel::Avx512 | Kernel::Avx512Batch => {
                sama_scalar(data, length, maj_length, min_length, first, out)
            }
            _ => unreachable!(),
        }
    }
}

#[inline]
pub fn sama_scalar(
    data: &[f64],
    length: usize,
    maj_length: usize,
    min_length: usize,
    first: usize,
    out: &mut [f64],
) {
    let n = data.len();
    if n == 0 {
        return;
    }

    let maj_alpha = 2.0 / (maj_length as f64 + 1.0);
    let min_alpha = 2.0 / (min_length as f64 + 1.0);
    let delta = min_alpha - maj_alpha;

    let cap = length + 1;
    let mut max_idx = vec![0usize; cap];
    let mut min_idx = vec![0usize; cap];
    let mut max_head = 0usize;
    let mut min_head = 0usize;
    let mut max_len = 0usize;
    let mut min_len = 0usize;

    let mut sama_val = f64::NAN;

    for i in first..n {
        let p = data[i];
        if p.is_nan() {
            out[i] = f64::NAN;
            continue;
        }

        let wstart = i.saturating_sub(length);

        while max_len > 0 {
            let idx = max_idx[max_head];
            if idx >= wstart {
                break;
            }
            max_head += 1;
            if max_head == cap {
                max_head = 0;
            }
            max_len -= 1;
        }

        while min_len > 0 {
            let idx = min_idx[min_head];
            if idx >= wstart {
                break;
            }
            min_head += 1;
            if min_head == cap {
                min_head = 0;
            }
            min_len -= 1;
        }

        while max_len > 0 {
            let mut last_pos = max_head + max_len - 1;
            if last_pos >= cap {
                last_pos -= cap;
            }
            let last_idx = max_idx[last_pos];
            let last_val = data[last_idx];
            if last_val <= p {
                max_len -= 1;
            } else {
                break;
            }
        }
        let mut ins_pos_max = max_head + max_len;
        if ins_pos_max >= cap {
            ins_pos_max -= cap;
        }
        max_idx[ins_pos_max] = i;
        max_len += 1;

        while min_len > 0 {
            let mut last_pos = min_head + min_len - 1;
            if last_pos >= cap {
                last_pos -= cap;
            }
            let last_idx = min_idx[last_pos];
            let last_val = data[last_idx];
            if last_val >= p {
                min_len -= 1;
            } else {
                break;
            }
        }
        let mut ins_pos_min = min_head + min_len;
        if ins_pos_min >= cap {
            ins_pos_min -= cap;
        }
        min_idx[ins_pos_min] = i;
        min_len += 1;

        let hh = data[max_idx[max_head]];
        let ll = data[min_idx[min_head]];

        let denom = hh - ll;
        let c = (p + p) - (hh + ll);
        let mult = if denom > 0.0 { c.abs() / denom } else { 0.0 };

        let a = mult.mul_add(delta, maj_alpha);
        let alpha = a * a;

        if sama_val.is_nan() {
            sama_val = p;
        } else {
            let diff = p - sama_val;
            sama_val = diff.mul_add(alpha, sama_val);
        }

        out[i] = sama_val;
    }
}

#[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
#[inline]
pub fn sama_avx2(
    data: &[f64],
    length: usize,
    maj_length: usize,
    min_length: usize,
    first: usize,
    out: &mut [f64],
) {
    sama_scalar(data, length, maj_length, min_length, first, out);
}

#[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
#[inline]
pub fn sama_avx512(
    data: &[f64],
    length: usize,
    maj_length: usize,
    min_length: usize,
    first: usize,
    out: &mut [f64],
) {
    sama_scalar(data, length, maj_length, min_length, first, out);
}

#[cfg(all(target_arch = "wasm32", target_feature = "simd128"))]
#[inline]
pub fn sama_simd128(
    data: &[f64],
    length: usize,
    maj_length: usize,
    min_length: usize,
    first: usize,
    out: &mut [f64],
) {
    sama_scalar(data, length, maj_length, min_length, first, out);
}

#[derive(Debug, Clone)]
pub struct SamaBatchRange {
    pub length: (usize, usize, usize),
    pub maj_length: (usize, usize, usize),
    pub min_length: (usize, usize, usize),
}

#[derive(Debug, Clone)]
pub struct SamaBatchOutput {
    pub values: Vec<f64>,
    pub combos: Vec<SamaParams>,
    pub rows: usize,
    pub cols: usize,
}

impl SamaBatchOutput {
    #[inline]
    pub fn row_for_params(&self, p: &SamaParams) -> Option<usize> {
        self.combos.iter().position(|c| {
            c.length.unwrap_or(200) == p.length.unwrap_or(200)
                && c.maj_length.unwrap_or(14) == p.maj_length.unwrap_or(14)
                && c.min_length.unwrap_or(6) == p.min_length.unwrap_or(6)
        })
    }

    #[inline]
    pub fn values_for(&self, p: &SamaParams) -> Option<&[f64]> {
        self.row_for_params(p).map(|row| {
            let start = row * self.cols;
            &self.values[start..start + self.cols]
        })
    }
}

#[derive(Clone, Debug, Default)]
pub struct SamaBatchBuilder {
    range: SamaBatchRange,
    kernel: Kernel,
}

impl SamaBatchBuilder {
    pub fn new() -> Self {
        Self::default()
    }

    pub fn kernel(mut self, k: Kernel) -> Self {
        self.kernel = k;
        self
    }

    #[inline]
    pub fn length_range(mut self, start: usize, end: usize, step: usize) -> Self {
        self.range.length = (start, end, step);
        self
    }

    #[inline]
    pub fn maj_length_range(mut self, start: usize, end: usize, step: usize) -> Self {
        self.range.maj_length = (start, end, step);
        self
    }

    #[inline]
    pub fn min_length_range(mut self, start: usize, end: usize, step: usize) -> Self {
        self.range.min_length = (start, end, step);
        self
    }

    #[inline]
    pub fn length_static(mut self, v: usize) -> Self {
        self.range.length = (v, v, 0);
        self
    }

    #[inline]
    pub fn maj_length_static(mut self, v: usize) -> Self {
        self.range.maj_length = (v, v, 0);
        self
    }

    #[inline]
    pub fn min_length_static(mut self, v: usize) -> Self {
        self.range.min_length = (v, v, 0);
        self
    }

    pub fn apply_slice(self, data: &[f64]) -> Result<SamaBatchOutput, SamaError> {
        sama_batch_with_kernel(data, &self.range, self.kernel)
    }

    pub fn apply_candles(self, c: &Candles, src: &str) -> Result<SamaBatchOutput, SamaError> {
        self.apply_slice(source_type(c, src))
    }

    pub fn with_default_slice(data: &[f64], k: Kernel) -> Result<SamaBatchOutput, SamaError> {
        SamaBatchBuilder::new().kernel(k).apply_slice(data)
    }

    pub fn with_default_candles(c: &Candles) -> Result<SamaBatchOutput, SamaError> {
        SamaBatchBuilder::new()
            .kernel(Kernel::Auto)
            .length_static(200)
            .maj_length_static(14)
            .min_length_static(6)
            .apply_candles(c, "close")
    }
}

impl Default for SamaBatchRange {
    fn default() -> Self {
        Self {
            length: (200, 449, 1),
            maj_length: (14, 14, 0),
            min_length: (6, 6, 0),
        }
    }
}

#[inline(always)]
fn axis_usize((start, end, step): (usize, usize, usize)) -> Result<Vec<usize>, SamaError> {
    if start == end || step == 0 {
        return Ok(vec![start]);
    }
    if start < end {
        if step == 0 {
            return Ok(vec![start]);
        }
        let mut v = Vec::new();
        let mut x = start;
        while x <= end {
            v.push(x);
            match x.checked_add(step) {
                Some(nx) => x = nx,
                None => break,
            }
        }
        if v.is_empty() {
            return Err(SamaError::InvalidRange { start, end, step });
        }
        Ok(v)
    } else {
        if step == 0 {
            return Ok(vec![start]);
        }
        let mut v = Vec::new();
        let mut x = start;
        loop {
            v.push(x);
            if x <= end {
                break;
            }
            x = x.saturating_sub(step);
            if x == 0 && end > 0 && x < end {
                break;
            }
            if v.len() > 1 && *v.last().unwrap() == x {
                break;
            }
            if x == 0 && end == 0 {
                break;
            }
            if x < end {
                break;
            }
        }
        if v.is_empty() {
            return Err(SamaError::InvalidRange { start, end, step });
        }
        Ok(v)
    }
}

#[inline(always)]
fn expand_grid_sama(r: &SamaBatchRange) -> Result<Vec<SamaParams>, SamaError> {
    let lens = axis_usize(r.length)?;
    let maj = axis_usize(r.maj_length)?;
    let min = axis_usize(r.min_length)?;
    if lens.is_empty() || maj.is_empty() || min.is_empty() {
        return Err(SamaError::InvalidRange {
            start: 0,
            end: 0,
            step: 0,
        });
    }
    let mut out = Vec::with_capacity(
        lens.len()
            .saturating_mul(maj.len())
            .saturating_mul(min.len()),
    );
    for &l in &lens {
        for &j in &maj {
            for &m in &min {
                out.push(SamaParams {
                    length: Some(l),
                    maj_length: Some(j),
                    min_length: Some(m),
                });
            }
        }
    }
    Ok(out)
}

pub fn sama_batch_with_kernel(
    data: &[f64],
    sweep: &SamaBatchRange,
    k: Kernel,
) -> Result<SamaBatchOutput, SamaError> {
    let kernel = match k {
        Kernel::Auto => detect_best_batch_kernel(),
        other if other.is_batch() => other,
        other => return Err(SamaError::InvalidKernelForBatch(other)),
    };

    let simd = match kernel {
        Kernel::Avx512Batch => Kernel::Avx512,
        Kernel::Avx2Batch => Kernel::Avx2,
        Kernel::ScalarBatch => Kernel::Scalar,
        _ => unreachable!(),
    };

    sama_batch_par_slice(data, sweep, simd)
}

#[inline(always)]
pub fn sama_batch_par_slice(
    data: &[f64],
    sweep: &SamaBatchRange,
    kern: Kernel,
) -> Result<SamaBatchOutput, SamaError> {
    sama_batch_inner(data, sweep, kern, true)
}

#[inline(always)]
pub fn sama_batch_slice(
    data: &[f64],
    sweep: &SamaBatchRange,
    kern: Kernel,
) -> Result<SamaBatchOutput, SamaError> {
    sama_batch_inner(data, sweep, kern, false)
}

#[inline(always)]
fn sama_batch_inner(
    data: &[f64],
    sweep: &SamaBatchRange,
    kern: Kernel,
    parallel: bool,
) -> Result<SamaBatchOutput, SamaError> {
    let combos = expand_grid_sama(sweep)?;
    let cols = data.len();
    let rows = combos.len();
    if cols == 0 {
        return Err(SamaError::EmptyInputData);
    }

    let total = rows.checked_mul(cols).ok_or(SamaError::InvalidRange {
        start: rows,
        end: cols,
        step: 0,
    })?;

    let mut buf_mu = make_uninit_matrix(rows, cols);

    let first = data
        .iter()
        .position(|x| !x.is_nan())
        .ok_or(SamaError::AllValuesNaN)?;
    let warm: Vec<usize> = combos.iter().map(|_| first).collect();
    init_matrix_prefixes(&mut buf_mu, cols, &warm);

    let mut guard = core::mem::ManuallyDrop::new(buf_mu);
    let out: &mut [f64] =
        unsafe { core::slice::from_raw_parts_mut(guard.as_mut_ptr() as *mut f64, guard.len()) };

    sama_batch_inner_into(data, &combos, first, kern, parallel, out)?;

    let values = unsafe {
        Vec::from_raw_parts(
            guard.as_mut_ptr() as *mut f64,
            guard.len(),
            guard.capacity(),
        )
    };

    Ok(SamaBatchOutput {
        values,
        combos,
        rows,
        cols,
    })
}

#[inline(always)]
fn sama_batch_inner_into(
    data: &[f64],
    combos: &[SamaParams],
    first: usize,
    kern: Kernel,
    parallel: bool,
    out_flat: &mut [f64],
) -> Result<(), SamaError> {
    if combos.is_empty() {
        return Ok(());
    }

    let rows = combos.len();
    let cols = data.len();

    if cols - first < 1 {
        return Err(SamaError::NotEnoughValidData {
            needed: 1,
            valid: cols - first,
        });
    }

    use std::collections::HashMap;
    let mut uniq_lengths: Vec<usize> = Vec::new();
    uniq_lengths.reserve(combos.len());
    for prm in combos.iter() {
        let l = prm.length.unwrap_or(200);
        if !uniq_lengths.contains(&l) {
            uniq_lengths.push(l);
        }
    }

    fn build_rolling_extrema(data: &[f64], length: usize) -> (Vec<f64>, Vec<f64>) {
        let n = data.len();
        let cap = length + 1;
        let mut max_idx = vec![0usize; cap];
        let mut min_idx = vec![0usize; cap];
        let mut max_head = 0usize;
        let mut min_head = 0usize;
        let mut max_len = 0usize;
        let mut min_len = 0usize;
        let mut hh = vec![f64::NAN; n];
        let mut ll = vec![f64::NAN; n];

        for i in 0..n {
            let p = data[i];
            let wstart = i.saturating_sub(length);

            while max_len > 0 {
                let idx = max_idx[max_head];
                if idx >= wstart {
                    break;
                }
                max_head += 1;
                if max_head == cap {
                    max_head = 0;
                }
                max_len -= 1;
            }
            while min_len > 0 {
                let idx = min_idx[min_head];
                if idx >= wstart {
                    break;
                }
                min_head += 1;
                if min_head == cap {
                    min_head = 0;
                }
                min_len -= 1;
            }

            if !p.is_nan() {
                while max_len > 0 {
                    let last_pos = (max_head + max_len - 1) % cap;
                    let last_idx = max_idx[last_pos];
                    if data[last_idx] <= p {
                        max_len -= 1;
                    } else {
                        break;
                    }
                }
                let ins_pos_max = (max_head + max_len) % cap;
                max_idx[ins_pos_max] = i;
                max_len += 1;

                while min_len > 0 {
                    let last_pos = (min_head + min_len - 1) % cap;
                    let last_idx = min_idx[last_pos];
                    if data[last_idx] >= p {
                        min_len -= 1;
                    } else {
                        break;
                    }
                }
                let ins_pos_min = (min_head + min_len) % cap;
                min_idx[ins_pos_min] = i;
                min_len += 1;
            }

            if max_len > 0 && min_len > 0 {
                hh[i] = data[max_idx[max_head]];
                ll[i] = data[min_idx[min_head]];
            }
        }
        (hh, ll)
    }

    let mut hh_map: HashMap<usize, Vec<f64>> = HashMap::with_capacity(uniq_lengths.len());
    let mut ll_map: HashMap<usize, Vec<f64>> = HashMap::with_capacity(uniq_lengths.len());
    for &l in &uniq_lengths {
        let (hh, ll) = build_rolling_extrema(data, l);
        hh_map.insert(l, hh);
        ll_map.insert(l, ll);
    }

    let do_row = |row: usize, row_dst: &mut [f64]| {
        let prm = &combos[row];
        let length = prm.length.unwrap_or(200);
        let maj_length = prm.maj_length.unwrap_or(14);
        let min_length = prm.min_length.unwrap_or(6);

        let maj_alpha = 2.0 / (maj_length as f64 + 1.0);
        let min_alpha = 2.0 / (min_length as f64 + 1.0);
        let delta = min_alpha - maj_alpha;

        let hh = &hh_map[&length];
        let ll = &ll_map[&length];

        let mut sama_val = f64::NAN;
        for i in first..cols {
            let p = data[i];
            if p.is_nan() {
                row_dst[i] = f64::NAN;
                continue;
            }

            let h = hh[i];
            let l = ll[i];
            let denom = h - l;
            let c = (p + p) - (h + l);
            let mult = if denom > 0.0 { c.abs() / denom } else { 0.0 };
            let a = mult.mul_add(delta, maj_alpha);
            let alpha = a * a;

            if sama_val.is_nan() {
                sama_val = p;
            } else {
                let diff = p - sama_val;
                sama_val = diff.mul_add(alpha, sama_val);
            }
            row_dst[i] = sama_val;
        }
    };

    if parallel {
        #[cfg(not(target_arch = "wasm32"))]
        {
            use rayon::prelude::*;
            out_flat
                .par_chunks_mut(cols)
                .enumerate()
                .for_each(|(row, dst)| do_row(row, dst));
        }
        #[cfg(target_arch = "wasm32")]
        {
            for (row, dst) in out_flat.chunks_mut(cols).enumerate() {
                do_row(row, dst);
            }
        }
    } else {
        for (row, dst) in out_flat.chunks_mut(cols).enumerate() {
            do_row(row, dst);
        }
    }

    Ok(())
}

#[derive(Debug, Clone)]
pub struct SamaStream {
    length: usize,
    maj_length: usize,
    min_length: usize,

    maj_alpha: f64,
    min_alpha: f64,
    delta: f64,

    cap: usize,
    buf: Vec<f64>,
    head: usize,
    tick: usize,

    max_idx: Vec<usize>,
    min_idx: Vec<usize>,
    max_head: usize,
    min_head: usize,
    max_len: usize,
    min_len: usize,

    sama_val: f64,
}

impl SamaStream {
    pub fn try_new(params: SamaParams) -> Result<Self, SamaError> {
        let length = params.length.unwrap_or(200);
        let maj = params.maj_length.unwrap_or(14);
        let min = params.min_length.unwrap_or(6);
        if length == 0 || maj == 0 || min == 0 {
            return Err(SamaError::InvalidPeriod {
                period: 0,
                data_len: 0,
            });
        }

        let cap = length + 1;
        let maj_alpha = 2.0 / (maj as f64 + 1.0);
        let min_alpha = 2.0 / (min as f64 + 1.0);

        Ok(Self {
            length,
            maj_length: maj,
            min_length: min,

            maj_alpha,
            min_alpha,
            delta: min_alpha - maj_alpha,

            cap,
            buf: vec![f64::NAN; cap],
            head: 0,
            tick: 0,

            max_idx: vec![0usize; cap],
            min_idx: vec![0usize; cap],
            max_head: 0,
            min_head: 0,
            max_len: 0,
            min_len: 0,

            sama_val: f64::NAN,
        })
    }

    #[inline]
    pub fn update(&mut self, value: f64) -> Option<f64> {
        let t = self.tick;
        let cap = self.cap;

        let oldest = t.saturating_sub(self.length);

        while self.max_len > 0 {
            let idx = self.max_idx[self.max_head];
            if idx >= oldest {
                break;
            }
            self.max_head += 1;
            if self.max_head == cap {
                self.max_head = 0;
            }
            self.max_len -= 1;
        }
        while self.min_len > 0 {
            let idx = self.min_idx[self.min_head];
            if idx >= oldest {
                break;
            }
            self.min_head += 1;
            if self.min_head == cap {
                self.min_head = 0;
            }
            self.min_len -= 1;
        }

        self.buf[self.head] = value;

        if value.is_nan() {
            self.head += 1;
            if self.head == cap {
                self.head = 0;
            }
            self.tick = t.wrapping_add(1);
            return Some(f64::NAN);
        }

        while self.max_len > 0 {
            let back_pos = (self.max_head + self.max_len - 1) % cap;
            let back_idx = self.max_idx[back_pos];
            let back_val = self.buf[back_idx % cap];
            if back_val <= value {
                self.max_len -= 1;
            } else {
                break;
            }
        }
        let ins_pos_max = (self.max_head + self.max_len) % cap;
        self.max_idx[ins_pos_max] = t;
        self.max_len += 1;

        while self.min_len > 0 {
            let back_pos = (self.min_head + self.min_len - 1) % cap;
            let back_idx = self.min_idx[back_pos];
            let back_val = self.buf[back_idx % cap];
            if back_val >= value {
                self.min_len -= 1;
            } else {
                break;
            }
        }
        let ins_pos_min = (self.min_head + self.min_len) % cap;
        self.min_idx[ins_pos_min] = t;
        self.min_len += 1;

        let hh = self.buf[self.max_idx[self.max_head] % cap];
        let ll = self.buf[self.min_idx[self.min_head] % cap];

        let denom = hh - ll;
        let c = (value + value) - (hh + ll);
        let mult = if denom > 0.0 { c.abs() / denom } else { 0.0 };

        let a = mult.mul_add(self.delta, self.maj_alpha);
        let alpha = a * a;

        if self.sama_val.is_nan() {
            self.sama_val = value;
        } else {
            let diff = value - self.sama_val;
            self.sama_val = diff.mul_add(alpha, self.sama_val);
        }

        self.head += 1;
        if self.head == cap {
            self.head = 0;
        }
        self.tick = t.wrapping_add(1);

        Some(self.sama_val)
    }

    #[inline]
    pub fn next(&mut self, value: f64) -> f64 {
        self.update(value).unwrap_or(f64::NAN)
    }

    #[inline]
    pub fn reset(&mut self) {
        self.buf.fill(f64::NAN);
        self.head = 0;
        self.tick = 0;
        self.max_len = 0;
        self.min_len = 0;
        self.max_head = 0;
        self.min_head = 0;
        self.sama_val = f64::NAN;
    }
}

#[cfg(feature = "python")]
#[pyfunction(name = "sama")]
#[pyo3(signature = (data, length, maj_length, min_length, kernel=None))]
pub fn sama_py<'py>(
    py: Python<'py>,
    data: PyReadonlyArray1<'py, f64>,
    length: usize,
    maj_length: usize,
    min_length: usize,
    kernel: Option<&str>,
) -> PyResult<Bound<'py, PyArray1<f64>>> {
    let slice_in = data.as_slice()?;
    let kern = validate_kernel(kernel, false)?;
    let params = SamaParams {
        length: Some(length),
        maj_length: Some(maj_length),
        min_length: Some(min_length),
    };
    let input = SamaInput::from_slice(slice_in, params);
    let result_vec: Vec<f64> = py
        .allow_threads(|| sama_with_kernel(&input, kern).map(|o| o.values))
        .map_err(|e| PyValueError::new_err(e.to_string()))?;
    Ok(result_vec.into_pyarray(py))
}

#[cfg(feature = "python")]
#[pyfunction(name = "sama_batch")]
#[pyo3(signature = (data, length_range, maj_length_range, min_length_range, kernel=None))]
pub fn sama_batch_py<'py>(
    py: Python<'py>,
    data: PyReadonlyArray1<'py, f64>,
    length_range: (usize, usize, usize),
    maj_length_range: (usize, usize, usize),
    min_length_range: (usize, usize, usize),
    kernel: Option<&str>,
) -> PyResult<Bound<'py, PyDict>> {
    use numpy::{IntoPyArray, PyArray1, PyArrayMethods};
    use pyo3::types::PyDict;

    let slice_in = data.as_slice()?;
    let sweep = SamaBatchRange {
        length: length_range,
        maj_length: maj_length_range,
        min_length: min_length_range,
    };

    let combos = expand_grid_sama(&sweep).map_err(|e| PyValueError::new_err(e.to_string()))?;
    let rows = combos.len();
    let cols = slice_in.len();

    let out_arr = unsafe { PyArray1::<f64>::new(py, [rows * cols], false) };
    let slice_out = unsafe { out_arr.as_slice_mut()? };

    let kern = validate_kernel(kernel, true)?;

    let combos = py
        .allow_threads(|| {
            let mapped = match kern {
                Kernel::Auto => detect_best_batch_kernel(),
                k => k,
            };
            let simd = match mapped {
                Kernel::Avx512Batch => Kernel::Avx512,
                Kernel::Avx2Batch => Kernel::Avx2,
                Kernel::ScalarBatch => Kernel::Scalar,
                _ => Kernel::Scalar,
            };

            let first = slice_in
                .iter()
                .position(|x| !x.is_nan())
                .ok_or(SamaError::AllValuesNaN)?;
            sama_batch_inner_into(slice_in, &combos, first, simd, true, slice_out).map(|_| combos)
        })
        .map_err(|e| PyValueError::new_err(e.to_string()))?;

    let dict = PyDict::new(py);
    dict.set_item("values", out_arr.reshape((rows, cols))?)?;
    dict.set_item(
        "lengths",
        combos
            .iter()
            .map(|p| p.length.unwrap_or(200) as u64)
            .collect::<Vec<_>>()
            .into_pyarray(py),
    )?;
    dict.set_item(
        "maj_lengths",
        combos
            .iter()
            .map(|p| p.maj_length.unwrap_or(14) as u64)
            .collect::<Vec<_>>()
            .into_pyarray(py),
    )?;
    dict.set_item(
        "min_lengths",
        combos
            .iter()
            .map(|p| p.min_length.unwrap_or(6) as u64)
            .collect::<Vec<_>>()
            .into_pyarray(py),
    )?;
    Ok(dict.into())
}

#[cfg(all(feature = "python", feature = "cuda"))]
#[pyfunction(name = "sama_cuda_batch_dev")]
#[pyo3(signature = (data_f32, length_range=(200, 200, 0), maj_length_range=(14, 14, 0), min_length_range=(6, 6, 0), device_id=0))]
pub fn sama_cuda_batch_dev_py(
    py: Python<'_>,
    data_f32: PyReadonlyArray1<'_, f32>,
    length_range: (usize, usize, usize),
    maj_length_range: (usize, usize, usize),
    min_length_range: (usize, usize, usize),
    device_id: usize,
) -> PyResult<DeviceArrayF32SamaPy> {
    if !cuda_available() {
        return Err(PyValueError::new_err("CUDA not available"));
    }

    let slice_in = data_f32.as_slice()?;
    let sweep = SamaBatchRange {
        length: length_range,
        maj_length: maj_length_range,
        min_length: min_length_range,
    };

    let inner = py.allow_threads(|| {
        let cuda = CudaSama::new(device_id).map_err(|e| PyValueError::new_err(e.to_string()))?;
        cuda.sama_batch_dev(slice_in, &sweep)
            .map_err(|e| PyValueError::new_err(e.to_string()))
    })?;

    Ok(DeviceArrayF32SamaPy { inner })
}

#[cfg(all(feature = "python", feature = "cuda"))]
#[pyfunction(name = "sama_cuda_many_series_one_param_dev")]
#[pyo3(signature = (data_tm_f32, length, maj_length, min_length, device_id=0))]
pub fn sama_cuda_many_series_one_param_dev_py(
    py: Python<'_>,
    data_tm_f32: PyReadonlyArray2<'_, f32>,
    length: usize,
    maj_length: usize,
    min_length: usize,
    device_id: usize,
) -> PyResult<DeviceArrayF32SamaPy> {
    if !cuda_available() {
        return Err(PyValueError::new_err("CUDA not available"));
    }
    if length == 0 || maj_length == 0 || min_length == 0 {
        return Err(PyValueError::new_err(
            "length, maj_length, and min_length must be positive",
        ));
    }

    let flat = data_tm_f32.as_slice()?;
    let shape = data_tm_f32.shape();
    let series_len = shape[0];
    let num_series = shape[1];
    let params = SamaParams {
        length: Some(length),
        maj_length: Some(maj_length),
        min_length: Some(min_length),
    };

    let inner = py.allow_threads(|| {
        let cuda = CudaSama::new(device_id).map_err(|e| PyValueError::new_err(e.to_string()))?;
        cuda.sama_many_series_one_param_time_major_dev(flat, num_series, series_len, &params)
            .map_err(|e| PyValueError::new_err(e.to_string()))
    })?;

    Ok(DeviceArrayF32SamaPy { inner })
}

#[cfg(feature = "python")]
#[pyclass(name = "SamaStream")]
pub struct SamaStreamPy {
    stream: SamaStream,
}

#[cfg(feature = "python")]
#[pymethods]
impl SamaStreamPy {
    #[new]
    fn new(length: usize, maj_length: usize, min_length: usize) -> PyResult<Self> {
        let params = SamaParams {
            length: Some(length),
            maj_length: Some(maj_length),
            min_length: Some(min_length),
        };
        let stream =
            SamaStream::try_new(params).map_err(|e| PyValueError::new_err(e.to_string()))?;
        Ok(SamaStreamPy { stream })
    }

    fn update(&mut self, value: f64) -> Option<f64> {
        self.stream.update(value)
    }
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[wasm_bindgen]
pub fn sama_js(
    data: &[f64],
    length: usize,
    maj_length: usize,
    min_length: usize,
) -> Result<Vec<f64>, JsValue> {
    let params = SamaParams {
        length: Some(length),
        maj_length: Some(maj_length),
        min_length: Some(min_length),
    };
    let input = SamaInput::from_slice(data, params);
    let mut out = vec![0.0; data.len()];
    sama_into_slice(&mut out, &input, detect_best_kernel())
        .map_err(|e| JsValue::from_str(&e.to_string()))?;
    Ok(out)
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[derive(Serialize, Deserialize)]
pub struct SamaBatchConfig {
    pub length_range: (usize, usize, usize),
    pub maj_length_range: (usize, usize, usize),
    pub min_length_range: (usize, usize, usize),
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[derive(Serialize, Deserialize)]
pub struct SamaBatchJsOutput {
    pub values: Vec<f64>,
    pub combos: Vec<SamaParams>,
    pub rows: usize,
    pub cols: usize,
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[wasm_bindgen(js_name = sama_batch)]
pub fn sama_batch_unified_js(data: &[f64], config: JsValue) -> Result<JsValue, JsValue> {
    let cfg: SamaBatchConfig = serde_wasm_bindgen::from_value(config)
        .map_err(|e| JsValue::from_str(&format!("Invalid config: {}", e)))?;
    let sweep = SamaBatchRange {
        length: cfg.length_range,
        maj_length: cfg.maj_length_range,
        min_length: cfg.min_length_range,
    };
    let out = sama_batch_inner(data, &sweep, detect_best_kernel(), false)
        .map_err(|e| JsValue::from_str(&e.to_string()))?;
    serde_wasm_bindgen::to_value(&SamaBatchJsOutput {
        values: out.values,
        combos: out.combos,
        rows: out.rows,
        cols: out.cols,
    })
    .map_err(|e| JsValue::from_str(&format!("Serialization error: {}", e)))
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[wasm_bindgen]
pub fn sama_alloc(len: usize) -> *mut f64 {
    let mut v = Vec::<f64>::with_capacity(len);
    let p = v.as_mut_ptr();
    std::mem::forget(v);
    p
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[wasm_bindgen]
pub fn sama_free(ptr: *mut f64, len: usize) {
    unsafe {
        let _ = Vec::from_raw_parts(ptr, len, len);
    }
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[wasm_bindgen]
pub fn sama_into(
    in_ptr: *const f64,
    out_ptr: *mut f64,
    len: usize,
    length: usize,
    maj_length: usize,
    min_length: usize,
) -> Result<(), JsValue> {
    if in_ptr.is_null() || out_ptr.is_null() {
        return Err(JsValue::from_str("null pointer"));
    }

    unsafe {
        let data = std::slice::from_raw_parts(in_ptr, len);
        let params = SamaParams {
            length: Some(length),
            maj_length: Some(maj_length),
            min_length: Some(min_length),
        };
        let input = SamaInput::from_slice(data, params);

        if in_ptr == out_ptr {
            let mut tmp = vec![0.0; len];
            sama_into_slice(&mut tmp, &input, detect_best_kernel())
                .map_err(|e| JsValue::from_str(&e.to_string()))?;
            let out = std::slice::from_raw_parts_mut(out_ptr, len);
            out.copy_from_slice(&tmp);
        } else {
            let out = std::slice::from_raw_parts_mut(out_ptr, len);

            sama_into_slice(out, &input, detect_best_kernel())
                .map_err(|e| JsValue::from_str(&e.to_string()))?;
        }
    }
    Ok(())
}

#[cfg(all(target_arch = "wasm32", feature = "wasm"))]
#[wasm_bindgen]
pub fn sama_batch_into(
    in_ptr: *const f64,
    out_ptr: *mut f64,
    len: usize,
    length_start: usize,
    length_end: usize,
    length_step: usize,
    maj_start: usize,
    maj_end: usize,
    maj_step: usize,
    min_start: usize,
    min_end: usize,
    min_step: usize,
) -> Result<usize, JsValue> {
    if in_ptr.is_null() || out_ptr.is_null() {
        return Err(JsValue::from_str("null pointer"));
    }
    unsafe {
        let data = std::slice::from_raw_parts(in_ptr, len);
        let sweep = SamaBatchRange {
            length: (length_start, length_end, length_step),
            maj_length: (maj_start, maj_end, maj_step),
            min_length: (min_start, min_end, min_step),
        };
        let combos = expand_grid_sama(&sweep).map_err(|_| JsValue::from_str("Invalid range"))?;
        let rows = combos.len();
        let cols = len;
        let out = std::slice::from_raw_parts_mut(out_ptr, rows * cols);
        let first = data
            .iter()
            .position(|x| !x.is_nan())
            .ok_or_else(|| JsValue::from_str("All NaN"))?;
        sama_batch_inner_into(data, &combos, first, detect_best_kernel(), false, out)
            .map_err(|e| JsValue::from_str(&e.to_string()))?;
        Ok(rows)
    }
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::skip_if_unsupported;
    use crate::utilities::data_loader::read_candles_from_csv;
    use std::error::Error;

    fn check_sama_accuracy(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let input = SamaInput::from_candles(&candles, "close", SamaParams::default());
        let result = sama_with_kernel(&input, kernel)?;

        assert_eq!(result.values.len(), candles.close.len());

        let valid_count = result.values.iter().filter(|&&v| !v.is_nan()).count();
        assert!(
            valid_count > 0,
            "[{}] SAMA should produce valid values",
            test_name
        );

        Ok(())
    }

    #[cfg(not(all(target_arch = "wasm32", feature = "wasm")))]
    #[test]
    fn test_sama_into_matches_api() -> Result<(), Box<dyn Error>> {
        let mut data = Vec::with_capacity(256);
        for _ in 0..5 {
            data.push(f64::NAN);
        }
        for i in 0..251 {
            let x = (i as f64).sin() * 0.12345 + 50.0 + ((i % 11) as f64) * 0.001;
            data.push(x);
        }

        let input = SamaInput::from_slice(&data, SamaParams::default());
        let baseline = sama(&input)?.values;

        let mut out = vec![0.0; data.len()];
        sama_into(&input, &mut out)?;

        assert_eq!(baseline.len(), out.len());
        fn eq_or_both_nan(a: f64, b: f64) -> bool {
            (a.is_nan() && b.is_nan()) || (a == b) || (a - b).abs() <= 1e-12
        }
        for (i, (&a, &b)) in baseline.iter().zip(out.iter()).enumerate() {
            assert!(
                eq_or_both_nan(a, b),
                "mismatch at index {}: api={} into={}",
                i,
                a,
                b
            );
        }

        Ok(())
    }

    fn check_sama_partial_params(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let default_params = SamaParams {
            length: None,
            maj_length: None,
            min_length: None,
        };
        let input = SamaInput::from_candles(&candles, "close", default_params);
        let output = sama_with_kernel(&input, kernel)?;
        assert_eq!(output.values.len(), candles.close.len());

        Ok(())
    }

    fn check_sama_default_candles(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let input = SamaInput::with_default_candles(&candles);
        match input.data {
            SamaData::Candles { source, .. } => assert_eq!(source, "close"),
            _ => panic!("Expected SamaData::Candles"),
        }
        let output = sama_with_kernel(&input, kernel)?;
        assert_eq!(output.values.len(), candles.close.len());

        Ok(())
    }

    fn check_sama_zero_period(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let input_data = [10.0, 20.0, 30.0];
        let params = SamaParams {
            length: Some(0),
            maj_length: None,
            min_length: None,
        };
        let input = SamaInput::from_slice(&input_data, params);
        let res = sama_with_kernel(&input, kernel);
        assert!(
            res.is_err(),
            "[{}] SAMA should fail with zero period",
            test_name
        );
        Ok(())
    }

    fn check_sama_period_exceeds_length(
        test_name: &str,
        kernel: Kernel,
    ) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let data_small = [10.0, 20.0, 30.0];
        let params = SamaParams {
            length: Some(10),
            maj_length: None,
            min_length: None,
        };
        let input = SamaInput::from_slice(&data_small, params);
        let res = sama_with_kernel(&input, kernel);
        assert!(
            res.is_err(),
            "[{}] SAMA should fail with period exceeding length",
            test_name
        );
        Ok(())
    }

    fn check_sama_very_small_dataset(
        test_name: &str,
        kernel: Kernel,
    ) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let single_point = [42.0];
        let params = SamaParams::default();
        let input = SamaInput::from_slice(&single_point, params);
        let res = sama_with_kernel(&input, kernel);
        assert!(
            res.is_err(),
            "[{}] SAMA should fail with insufficient data",
            test_name
        );
        Ok(())
    }

    fn check_sama_empty_input(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let empty: [f64; 0] = [];
        let params = SamaParams::default();
        let input = SamaInput::from_slice(&empty, params);
        let res = sama_with_kernel(&input, kernel);
        assert!(
            matches!(res, Err(SamaError::EmptyInputData)),
            "[{}] SAMA should fail with empty input",
            test_name
        );
        Ok(())
    }

    fn check_sama_all_nan(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let nan_data = [f64::NAN, f64::NAN, f64::NAN];
        let params = SamaParams::default();
        let input = SamaInput::from_slice(&nan_data, params);
        let res = sama_with_kernel(&input, kernel);
        assert!(
            matches!(res, Err(SamaError::AllValuesNaN)),
            "[{}] SAMA should fail with all NaN values",
            test_name
        );
        Ok(())
    }

    fn check_sama_reinput(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let first_params = SamaParams {
            length: Some(50),
            maj_length: Some(14),
            min_length: Some(6),
        };
        let first_input = SamaInput::from_candles(&candles, "close", first_params);
        let first_result = sama_with_kernel(&first_input, kernel)?;

        let second_params = SamaParams {
            length: Some(50),
            maj_length: Some(14),
            min_length: Some(6),
        };
        let second_input = SamaInput::from_slice(&first_result.values, second_params);
        let second_result = sama_with_kernel(&second_input, kernel)?;

        assert_eq!(second_result.values.len(), first_result.values.len());
        Ok(())
    }

    fn check_sama_nan_handling(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let params = SamaParams::default();
        let input = SamaInput::from_candles(&candles, "close", params);
        let output = sama_with_kernel(&input, kernel)?;

        let first_valid = candles.close.iter().position(|x| !x.is_nan()).unwrap_or(0);
        let warmup = first_valid + input.get_length();

        for (i, &val) in output
            .values
            .iter()
            .enumerate()
            .skip(warmup.min(output.values.len()))
        {
            assert!(
                !val.is_nan(),
                "[{}] Unexpected NaN at index {}",
                test_name,
                i
            );
        }
        Ok(())
    }

    fn check_sama_streaming(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let params = SamaParams::default();

        let batch_input = SamaInput::from_candles(&candles, "close", params.clone());
        let batch_result = sama_with_kernel(&batch_input, kernel)?;

        let mut stream = SamaStream::try_new(params)?;
        let mut stream_results = Vec::with_capacity(candles.close.len());
        for &price in &candles.close {
            stream_results.push(stream.next(price));
        }

        assert_eq!(batch_result.values.len(), stream_results.len());

        for (i, (&batch_val, &stream_val)) in batch_result
            .values
            .iter()
            .zip(stream_results.iter())
            .enumerate()
        {
            if batch_val.is_nan() && stream_val.is_nan() {
                continue;
            }
            if !batch_val.is_nan() && !stream_val.is_nan() {
                let diff = (batch_val - stream_val).abs();
                assert!(
                    diff < 1e-9,
                    "[{}] Stream mismatch at index {}: batch={}, stream={}, diff={}",
                    test_name,
                    i,
                    batch_val,
                    stream_val,
                    diff
                );
            }
        }
        Ok(())
    }

    fn check_batch_sweep(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let range = SamaBatchRange {
            length: (190, 210, 1),
            maj_length: (12, 16, 1),
            min_length: (4, 8, 1),
        };

        let output = sama_batch_with_kernel(&candles.close, &range, kernel)?;

        let expected_count = 21 * 5 * 5;
        assert_eq!(
            output.rows, expected_count,
            "[{}] Expected {} batch results",
            test_name, expected_count
        );
        assert_eq!(
            output.values.len(),
            expected_count * candles.close.len(),
            "[{}] Expected flattened array size",
            test_name
        );

        Ok(())
    }

    fn check_batch_default_row(test_name: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test_name);
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let range = SamaBatchRange {
            length: (200, 200, 0),
            maj_length: (14, 14, 0),
            min_length: (6, 6, 0),
        };

        let output = sama_batch_with_kernel(&candles.close, &range, kernel)?;

        assert_eq!(
            output.rows, 1,
            "[{}] Should have 1 result for default params",
            test_name
        );
        assert_eq!(output.cols, candles.close.len());

        Ok(())
    }

    macro_rules! generate_all_sama_tests {
        ($($test_fn:ident),*) => {
            paste::paste! {
                $(
                    #[test]
                    fn [<$test_fn _scalar>]() -> Result<(), Box<dyn Error>> {
                        $test_fn(stringify!([<$test_fn _scalar>]), Kernel::Scalar)
                    }
                )*

                #[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
                $(
                    #[test]
                    fn [<$test_fn _avx2>]() -> Result<(), Box<dyn Error>> {
                        $test_fn(stringify!([<$test_fn _avx2>]), Kernel::Avx2)
                    }

                    #[test]
                    fn [<$test_fn _avx512>]() -> Result<(), Box<dyn Error>> {
                        $test_fn(stringify!([<$test_fn _avx512>]), Kernel::Avx512)
                    }
                )*

                #[cfg(all(target_arch = "wasm32", target_feature = "simd128"))]
                $(
                    #[test]
                    fn [<$test_fn _simd128>]() -> Result<(), Box<dyn Error>> {
                        $test_fn(stringify!([<$test_fn _simd128>]), Kernel::Scalar)
                    }
                )*
            }
        }
    }

    macro_rules! gen_batch_tests {
        ($fn_name:ident) => {
            paste::paste! {
                #[test]
                fn [<$fn_name _scalar>]() -> Result<(), Box<dyn Error>> {
                    $fn_name(stringify!([<$fn_name _scalar>]), Kernel::ScalarBatch)
                }

                #[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
                #[test]
                fn [<$fn_name _avx2>]() -> Result<(), Box<dyn Error>> {
                    $fn_name(stringify!([<$fn_name _avx2>]), Kernel::Avx2Batch)
                }

                #[cfg(all(feature = "nightly-avx", target_arch = "x86_64"))]
                #[test]
                fn [<$fn_name _avx512>]() -> Result<(), Box<dyn Error>> {
                    $fn_name(stringify!([<$fn_name _avx512>]), Kernel::Avx512Batch)
                }

                #[test]
                fn [<$fn_name _auto_detect>]() -> Result<(), Box<dyn Error>> {
                    $fn_name(stringify!([<$fn_name _auto_detect>]), Kernel::Auto)
                }
            }
        };
    }

    generate_all_sama_tests!(
        check_sama_accuracy,
        check_sama_partial_params,
        check_sama_default_candles,
        check_sama_zero_period,
        check_sama_period_exceeds_length,
        check_sama_very_small_dataset,
        check_sama_empty_input,
        check_sama_all_nan,
        check_sama_reinput,
        check_sama_nan_handling,
        check_sama_streaming
    );

    gen_batch_tests!(check_batch_sweep);
    gen_batch_tests!(check_batch_default_row);

    #[cfg(debug_assertions)]
    fn check_batch_no_poison(test: &str, kernel: Kernel) -> Result<(), Box<dyn Error>> {
        skip_if_unsupported!(kernel, test);

        let file = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let c = read_candles_from_csv(file)?;

        let cfgs = vec![
            (190, 200, 1, 12, 14, 1, 4, 8, 1),
            (200, 200, 0, 14, 14, 0, 6, 6, 0),
            (195, 205, 5, 10, 16, 2, 3, 9, 2),
        ];

        for (ls, le, lstep, js, je, jstep, ms, me, mstep) in cfgs {
            let out = SamaBatchRange {
                length: (ls, le, lstep),
                maj_length: (js, je, jstep),
                min_length: (ms, me, mstep),
            };
            let res = sama_batch_with_kernel(&c.close, &out, kernel)?;
            for &v in &res.values {
                if v.is_nan() {
                    continue;
                }
                let b = v.to_bits();
                assert_ne!(b, 0x1111_1111_1111_1111);
                assert_ne!(b, 0x2222_2222_2222_2222);
                assert_ne!(b, 0x3333_3333_3333_3333);
            }
        }
        Ok(())
    }

    #[cfg(debug_assertions)]
    gen_batch_tests!(check_batch_no_poison);

    #[cfg(not(debug_assertions))]
    fn check_batch_no_poison_scalar() -> Result<(), Box<dyn Error>> {
        Ok(())
    }

    #[cfg(all(target_arch = "wasm32", target_feature = "simd128"))]
    #[test]
    fn test_sama_simd128_correctness() {
        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0];
        let params = SamaParams::default();
        let input = SamaInput::from_slice(&data, params);
        let scalar = sama_with_kernel(&input, Kernel::Scalar).unwrap();
        let simd = sama_with_kernel(&input, Kernel::Scalar).unwrap();
        assert_eq!(scalar.values.len(), simd.values.len());
        for (a, b) in scalar.values.iter().zip(simd.values.iter()) {
            assert!((a - b).abs() < 1e-10);
        }
    }

    #[cfg(debug_assertions)]
    #[test]
    fn test_sama_no_poison_values() -> Result<(), Box<dyn Error>> {
        let file_path = "src/data/2018-09-01-2024-Bitfinex_Spot-4h.csv";
        let candles = read_candles_from_csv(file_path)?;

        let input = SamaInput::from_candles(&candles, "close", SamaParams::default());
        let output = sama(&input)?;

        for &v in &output.values {
            if v.is_nan() {
                continue;
            }
            let b = v.to_bits();

            assert_ne!(
                b, 0x11111111_11111111,
                "Found poison value 0x11111111_11111111"
            );
            assert_ne!(
                b, 0x22222222_22222222,
                "Found poison value 0x22222222_22222222"
            );
            assert_ne!(
                b, 0x33333333_33333333,
                "Found poison value 0x33333333_33333333"
            );
            assert_ne!(
                b, 0xDEADBEEF_DEADBEEF,
                "Found poison value 0xDEADBEEF_DEADBEEF"
            );
            assert_ne!(
                b, 0xFEEEFEEE_FEEEFEEE,
                "Found poison value 0xFEEEFEEE_FEEEFEEE"
            );
        }
        Ok(())
    }

    #[test]
    fn test_sama_stream_incremental() -> Result<(), Box<dyn Error>> {
        let params = SamaParams {
            length: Some(10),
            maj_length: Some(5),
            min_length: Some(3),
        };

        let mut stream = SamaStream::try_new(params)?;
        let data = vec![
            10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 110.0, 120.0,
        ];

        let mut results = Vec::new();
        for &val in &data {
            let result = stream.next(val);
            if !result.is_nan() {
                results.push(result);
            }
        }

        assert!(
            !results.is_empty(),
            "Stream should produce results immediately"
        );

        Ok(())
    }

    #[test]
    fn sama_pine_parity_head_start() -> Result<(), Box<dyn Error>> {
        let mut long = vec![0.0; 5000];
        for i in 0..long.len() {
            long[i] = (i as f64).sin() + (i as f64 * 0.01).cos();
        }

        let pine_params = SamaParams::default();
        let pine_out = SamaInput::from_slice(&long, pine_params.clone());
        let a = sama_with_kernel(&pine_out, Kernel::Scalar)?.values;

        let tail = &long[2000..];
        let pine_like_tail = SamaInput::from_slice(tail, pine_params);
        let b = sama_with_kernel(&pine_like_tail, Kernel::Scalar)?.values;

        let tol = 0.1;
        for (i, (&x, &y)) in a[2000..].iter().zip(b.iter()).enumerate().skip(100) {
            if x.is_finite() && y.is_finite() {
                assert!((x - y).abs() < tol, "i={}, |Δ|={}", i, (x - y).abs());
            }
        }

        Ok(())
    }
}

#[cfg(all(feature = "proptest", not(target_arch = "wasm32")))]
mod prop_tests {
    use super::*;
    use proptest::prelude::*;

    proptest! {
        #[test]
        fn sama_properties(data in prop::collection::vec(-1e6f64..1e6, 5..400),
                           len in 2usize..64,
                           maj in 2usize..64,
                           min in 2usize..64) {


            if data.len() <= len {
                return Ok(());
            }

            let params = SamaParams {
                length: Some(len),
                maj_length: Some(maj),
                min_length: Some(min),
            };
            let input = SamaInput::from_slice(&data, params);
            let SamaOutput { values: out } = sama_with_kernel(&input, Kernel::Scalar).unwrap();


            let first = data.iter().position(|x| !x.is_nan()).unwrap_or(0);
            let warm = first + len;
            for i in warm..data.len() {
                let wstart = i - len;
                let window = &data[wstart..=i];
                if window.iter().all(|v| v.is_finite()) {
                    let y = out[i];



                    prop_assert!(
                        y.is_finite(),
                        "Output {} at index {} is not finite",
                        y, i
                    );
                }
            }
        }
    }
}

#[cfg(all(feature = "python", feature = "cuda"))]
#[pyclass(module = "ta_indicators.cuda", unsendable)]
pub struct DeviceArrayF32SamaPy {
    pub(crate) inner: DeviceArrayF32Sama,
}

#[cfg(all(feature = "python", feature = "cuda"))]
#[pymethods]
impl DeviceArrayF32SamaPy {
    #[getter]
    fn __cuda_array_interface__<'py>(&self, py: Python<'py>) -> PyResult<Bound<'py, PyDict>> {
        let d = PyDict::new(py);
        d.set_item("shape", (self.inner.rows, self.inner.cols))?;
        d.set_item("typestr", "<f4")?;

        d.set_item(
            "strides",
            (
                self.inner.cols * std::mem::size_of::<f32>(),
                std::mem::size_of::<f32>(),
            ),
        )?;

        let len = self.inner.rows.checked_mul(self.inner.cols).unwrap_or(0);
        let ptr_usize = if len == 0 {
            0usize
        } else {
            self.inner.device_ptr() as usize
        };
        d.set_item("data", (ptr_usize, false))?;

        d.set_item("version", 3)?;
        Ok(d)
    }

    fn __dlpack_device__(&self) -> (i32, i32) {
        (2, self.inner.device_id as i32)
    }

    #[pyo3(signature = (stream=None, max_version=None, dl_device=None, copy=None))]
    fn __dlpack__<'py>(
        &mut self,
        py: Python<'py>,
        stream: Option<pyo3::PyObject>,
        max_version: Option<pyo3::PyObject>,
        dl_device: Option<pyo3::PyObject>,
        copy: Option<pyo3::PyObject>,
    ) -> PyResult<PyObject> {
        let (kdl, alloc_dev) = self.__dlpack_device__();
        if let Some(dev_obj) = dl_device.as_ref() {
            if let Ok((dev_ty, dev_id)) = dev_obj.extract::<(i32, i32)>(py) {
                if dev_ty != kdl || dev_id != alloc_dev {
                    let wants_copy = copy
                        .as_ref()
                        .and_then(|c| c.extract::<bool>(py).ok())
                        .unwrap_or(false);
                    if wants_copy {
                        return Err(PyValueError::new_err(
                            "device copy not implemented for __dlpack__",
                        ));
                    } else {
                        return Err(PyValueError::new_err("dl_device mismatch for __dlpack__"));
                    }
                }
            }
        }

        let _ = stream;

        let dummy =
            DeviceBuffer::from_slice(&[]).map_err(|e| PyValueError::new_err(e.to_string()))?;
        let device_id = self.inner.device_id;
        let ctx = self.inner.ctx.clone();
        let inner = std::mem::replace(
            &mut self.inner,
            DeviceArrayF32Sama {
                buf: dummy,
                rows: 0,
                cols: 0,
                ctx,
                device_id,
            },
        );

        let rows = inner.rows;
        let cols = inner.cols;
        let buf = inner.buf;

        let max_version_bound = max_version.map(|obj| obj.into_bound(py));

        export_f32_cuda_dlpack_2d(py, buf, rows, cols, alloc_dev, max_version_bound)
    }
}