SpectralPlot/src-tauri/src/algorithm/spectraltools.rs

extern crate splines;

use splines::{Spline, Key, Interpolation};
use std::error::Error;

use find_peaks::PeakFinder;

pub fn interpolate_spline<T: Copy + Into<f64>,>(x_t: Vec<T>, y_t: Vec<T>, step: f64) -> Result<Vec<(f64, f64)>, Box<dyn Error>> {
    let x: Vec<f64> = x_t.iter().map(|&x| x.into()).collect();
    let y: Vec<f64> = y_t.iter().map(|&y| y.into()).collect();

    if x.len() != y.len() {
        return Err("x and y must have the same length".into());
    }

    // 创建样条曲线
    let keys: Vec<Key<f64, f64>> = x.iter()
        .zip(y.iter())
        .map(|(&x, &y)| Key::new(x, y, Interpolation::Linear))
        .collect();

    let spline = Spline::from_vec(keys);

    // 计算 x 的最大值和最小值
    let &start = x.iter().min_by(|a, b| a.partial_cmp(b).unwrap()).unwrap();
    let &end = x.iter().max_by(|a, b| a.partial_cmp(b).unwrap()).unwrap();

    // 插值到间隔为 step 的点
    let mut result = Vec::new();
    let mut t = start;

    while t <= end {
        if let Some(value) = spline.clamped_sample(t) {
            result.push((t, value));
        }
        t += step;
    }

    Ok(result)
}

pub fn interpolate_spline_at_points<T: Copy + Into<f64>>(x_t: Vec<T>, y_t: Vec<T>, x_target: Vec<f64>) -> Result<Vec<f64>, Box<dyn Error>> {
    let x: Vec<f64> = x_t.iter().map(|&x| x.into()).collect();
    let y: Vec<f64> = y_t.iter().map(|&y| y.into()).collect();

    if x.len() != y.len() {
        return Err("x and y must have the same length".into());
    }

    // 创建样条曲线
    let keys: Vec<Key<f64, f64>> = x.iter()
        .zip(y.iter())
        .map(|(&x, &y)| Key::new(x, y, Interpolation::Linear))
        .collect();

    let spline = Spline::from_vec(keys);

    // 插值到 x_target 指定的点
    let mut result = Vec::new();

    for &t in x_target.iter() {
        if let Some(value) = spline.clamped_sample(t) {
            result.push(value);
        }
    }

    Ok(result)
}
pub fn find_peek(data:Vec<f64>,minheigh:f64)->Vec<(u32,f64)>{
    let mut fp = PeakFinder::new(&data);
    fp.with_min_prominence(200.);
    fp.with_min_height(minheigh);
    let mut retvec=Vec::new();
    let peaks = fp.find_peaks();
    for p in peaks {
      //  println!("{} {}", p.middle_position(), p.height.unwrap());
        retvec.push((p.middle_position().try_into().unwrap(),p.height.unwrap()));

    }
    retvec
}


#[test]
fn testinterpolate_spline() -> Result<(), Box<dyn Error>> {
    // 示例数据
    let x = vec![0.0,0.5, 0.569, 1.138, 1.707, 2.276, 2.845];
    let y = vec![0.0, 0.4,0.5, 1.0, 0.5, 0.0, -0.5];
    let step = 0.1;

    // 调用插值函数
    let interpolated_values = interpolate_spline(x, y, step)?;

    // 输出结果
    for (xi, yi) in interpolated_values {
        println!("x = {:.3}, y = {:.3}", xi, yi);
    }

    Ok(())
}
#[test]
fn tset_interpolate_spline_at_points() -> Result<(), Box<dyn Error>> {
    let x = vec![0.0,0.5, 0.569, 1.138, 1.707, 2.276, 2.845];
    let y = vec![0.0, 0.4,0.5, 1.0, 0.5, 0.0, -0.5];
    let x_target = vec![0.1, 0.2, 0.3];

    let result = interpolate_spline_at_points(x, y, x_target)?;

    for (yi) in result {
        println!("y = {:.3}",  yi);
    }

    Ok(())
}

use csv::ReaderBuilder;
fn read_csv_to_vec(file_path: &str) -> Result<Vec<f64>, Box<dyn Error>> {
    let mut rdr = ReaderBuilder::new().from_path(file_path)?;
    let mut values = Vec::new();

    for result in rdr.records() {
        let record = result?;
        if let Some(value) = record.get(1) {
            values.push(value.parse::<f64>()?);
        }
    }

    Ok(values)
}

use nalgebra::{DMatrix, DVector};

pub fn compute_weave_coeff(x_data:Vec<f64>,y_data:Vec<f64>)->Vec<f64>{

    assert_eq!(x_data.len(), y_data.len());

    let n = x_data.len();

    // 构建设计矩阵 X 和观测向量 y
    let mut x_matrix = DMatrix::zeros(n, 4); // 三阶多项式有 4 个系数
    let y_vector = DVector::from_vec(y_data.clone());

    for (i, &x) in x_data.iter().enumerate() {
        x_matrix[(i, 0)] = 1.0;       // 常数项
        x_matrix[(i, 1)] = x;         // x
        x_matrix[(i, 2)] = x.powi(2); // x²
        x_matrix[(i, 3)] = x.powi(3); // x³
    }

    // 使用正规方程求解最小二乘问题： (XᵀX)β = Xᵀy
    let xt = x_matrix.transpose();
    let xtx = &xt * &x_matrix;
    let xty = &xt * &y_vector;

    // 求解方程 (XᵀX)β = Xᵀy
    let beta = xtx
        .lu()
        .solve(&xty)
        .expect("无法求解正规方程，可能是矩阵奇异");

    // 输出拟合系数
    println!("拟合的三阶多项式系数：");
    println!("y = {:.4} + {:.4}x + {:.4}x² + {:.4}x³", beta[0], beta[1], beta[2], beta[3]);

    // 示例：使用拟合的多项式进行预测
    let x_test = 6.0;
    let y_pred = beta[0] + beta[1]*x_test + beta[2]*x_test.powi(2) + beta[3]*x_test.powi(3);
    println!("对于 x = {:.2}, 预测的 y = {:.4}", x_test, y_pred);
    let mut retvec=Vec::new();
    for i in 0..4{
        retvec.push(beta[i]);
    }
    retvec
   
}

#[test]
fn test_find_peek(){
    let data = read_csv_to_vec("D:\\06Learn\\rust\\tarui\\myfirst_tauri\\src-tauri\\test0_UP.csv").unwrap();
    let peaks = find_peek(data,10000.0);
    for p in peaks {
        println!("{} {}", p.0, p.1);
    }
}