fdars_core/detrend/
diff.rs1use crate::iter_maybe_parallel;
2use crate::matrix::FdMatrix;
3#[cfg(feature = "parallel")]
4use rayon::iter::ParallelIterator;
5use std::borrow::Cow;
6
7use super::{reassemble_polynomial_results, TrendResult};
8
9pub(super) fn diff_single_curve(
10 curve: &[f64],
11 m: usize,
12 order: usize,
13) -> (Vec<f64>, Vec<f64>, Vec<f64>, f64) {
14 let diff1: Vec<f64> = (0..m - 1).map(|j| curve[j + 1] - curve[j]).collect();
15 let detrended = if order == 2 {
16 (0..diff1.len() - 1)
17 .map(|j| diff1[j + 1] - diff1[j])
18 .collect()
19 } else {
20 diff1.clone()
21 };
22 let initial_values = if order == 2 {
23 vec![curve[0], curve[1]]
24 } else {
25 vec![curve[0]]
26 };
27 let rss: f64 = detrended.iter().map(|&x| x.powi(2)).sum();
28 let new_m = m - order;
29 let mut trend = vec![0.0; m];
30 trend[0] = curve[0];
31 if order == 1 {
32 for j in 1..m {
33 trend[j] = curve[j] - if j <= new_m { detrended[j - 1] } else { 0.0 };
34 }
35 } else {
36 trend = curve.to_vec();
37 }
38 let mut det_full = vec![0.0; m];
39 det_full[..new_m].copy_from_slice(&detrended[..new_m]);
40 (trend, det_full, initial_values, rss)
41}
42
43#[must_use = "returns the detrending result which should not be discarded"]
45pub fn detrend_diff(data: &FdMatrix, order: usize) -> TrendResult {
46 let (n, m) = data.shape();
47 if n == 0 || m <= order || order == 0 || order > 2 {
48 return TrendResult::empty(data, n, m, Cow::Owned(format!("diff{order}")), order);
49 }
50 let results: Vec<(Vec<f64>, Vec<f64>, Vec<f64>, f64)> = iter_maybe_parallel!(0..n)
51 .map(|i| {
52 let curve: Vec<f64> = (0..m).map(|j| data[(i, j)]).collect();
53 diff_single_curve(&curve, m, order)
54 })
55 .collect();
56 let (trend, detrended, coefficients, rss) = reassemble_polynomial_results(results, n, m, order);
57 TrendResult {
58 trend,
59 detrended,
60 method: Cow::Owned(format!("diff{order}")),
61 coefficients: Some(coefficients),
62 rss,
63 n_params: order,
64 }
65}