Skip to main content

fdars_core/alignment/
set.rs

1//! Set-level alignment operations and elastic decomposition.
2
3use super::pairwise::elastic_align_pair;
4use super::srsf::reparameterize_curve;
5use super::{AlignmentResult, AlignmentSetResult};
6use crate::iter_maybe_parallel;
7use crate::matrix::FdMatrix;
8#[cfg(feature = "parallel")]
9use rayon::iter::ParallelIterator;
10
11/// Result of elastic phase-amplitude decomposition.
12#[derive(Debug, Clone, PartialEq)]
13#[non_exhaustive]
14pub struct DecompositionResult {
15    /// Full alignment result.
16    pub alignment: AlignmentResult,
17    /// Amplitude distance: SRSF distance after alignment.
18    pub d_amplitude: f64,
19    /// Phase distance: geodesic distance of warp from identity.
20    pub d_phase: f64,
21}
22
23/// Align all curves in `data` to a single target curve.
24///
25/// # Arguments
26/// * `data` — Functional data matrix (n × m)
27/// * `target` — Target curve to align to (length m)
28/// * `argvals` — Evaluation points (length m)
29/// * `lambda` — Penalty weight on warp deviation from identity (0.0 = no penalty)
30///
31/// # Returns
32/// [`AlignmentSetResult`] with all warping functions, aligned curves, and distances.
33///
34/// # Examples
35///
36/// ```
37/// use fdars_core::matrix::FdMatrix;
38/// use fdars_core::alignment::align_to_target;
39///
40/// let argvals: Vec<f64> = (0..20).map(|i| i as f64 / 19.0).collect();
41/// let target: Vec<f64> = argvals.iter().map(|&t| (t * 6.0).sin()).collect();
42/// let data = FdMatrix::from_column_major(
43///     (0..60).map(|i| ((i as f64 * 0.05) + 0.1).sin()).collect(),
44///     3, 20,
45/// ).unwrap();
46/// let result = align_to_target(&data, &target, &argvals, 0.0);
47/// assert_eq!(result.aligned_data.shape(), (3, 20));
48/// assert_eq!(result.distances.len(), 3);
49/// ```
50#[must_use = "expensive computation whose result should not be discarded"]
51pub fn align_to_target(
52    data: &FdMatrix,
53    target: &[f64],
54    argvals: &[f64],
55    lambda: f64,
56) -> AlignmentSetResult {
57    let (n, m) = data.shape();
58
59    let results: Vec<AlignmentResult> = iter_maybe_parallel!(0..n)
60        .map(|i| {
61            let fi = data.row(i);
62            elastic_align_pair(target, &fi, argvals, lambda)
63        })
64        .collect();
65
66    let mut gammas = FdMatrix::zeros(n, m);
67    let mut aligned_data = FdMatrix::zeros(n, m);
68    let mut distances = Vec::with_capacity(n);
69
70    for (i, r) in results.into_iter().enumerate() {
71        for j in 0..m {
72            gammas[(i, j)] = r.gamma[j];
73            aligned_data[(i, j)] = r.f_aligned[j];
74        }
75        distances.push(r.distance);
76    }
77
78    AlignmentSetResult {
79        gammas,
80        aligned_data,
81        distances,
82    }
83}
84
85/// Perform elastic phase-amplitude decomposition of two curves.
86///
87/// Returns both the alignment result and the separate amplitude and phase distances.
88///
89/// # Arguments
90/// * `f1` — Target curve (length m)
91/// * `f2` — Curve to decompose against f1 (length m)
92/// * `argvals` — Evaluation points (length m)
93/// * `lambda` — Penalty weight on warp deviation from identity (0.0 = no penalty)
94pub fn elastic_decomposition(
95    f1: &[f64],
96    f2: &[f64],
97    argvals: &[f64],
98    lambda: f64,
99) -> DecompositionResult {
100    let alignment = elastic_align_pair(f1, f2, argvals, lambda);
101    let d_amplitude = alignment.distance;
102    let d_phase = crate::warping::phase_distance(&alignment.gamma, argvals);
103    DecompositionResult {
104        alignment,
105        d_amplitude,
106        d_phase,
107    }
108}
109
110/// Apply stored warps to original curves to produce aligned data.
111pub(super) fn apply_stored_warps(data: &FdMatrix, gammas: &FdMatrix, argvals: &[f64]) -> FdMatrix {
112    let (n, m) = data.shape();
113    let mut aligned = FdMatrix::zeros(n, m);
114    for i in 0..n {
115        let fi = data.row(i);
116        let gamma: Vec<f64> = (0..m).map(|j| gammas[(i, j)]).collect();
117        let f_aligned = reparameterize_curve(&fi, argvals, &gamma);
118        for j in 0..m {
119            aligned[(i, j)] = f_aligned[j];
120        }
121    }
122    aligned
123}