1use crate::simd;
2
3#[must_use]
4pub fn filter_valid(x: &[f64], taps: &[f64]) -> Vec<f64> {
5 if taps.is_empty() || x.len() < taps.len() {
6 return Vec::new();
7 }
8 let n = x.len() - taps.len() + 1;
9 let mut out = vec![0.0; n];
10 for (k, &t) in taps.iter().enumerate() {
11 simd::axpy(t, &x[k..k + n], &mut out);
12 }
13 out
14}
15
16#[must_use]
17pub fn filter_same(x: &[f64], taps: &[f64]) -> Vec<f64> {
18 assert!(taps.len() % 2 == 1, "centred filtering needs an odd tap count");
19 let r = taps.len() / 2;
20 if x.is_empty() {
21 return Vec::new();
22 }
23 let mut padded = Vec::with_capacity(x.len() + 2 * r);
24 padded.extend(std::iter::repeat_n(x[0], r));
25 padded.extend_from_slice(x);
26 padded.extend(std::iter::repeat_n(x[x.len() - 1], r));
27 filter_valid(&padded, taps)
28}
29
30#[must_use]
31pub fn filter_2d(img: &[f64], w: usize, h: usize, htaps: &[f64], vtaps: &[f64]) -> Vec<f64> {
32 let mut tmp = vec![0.0; w * h];
33 for r in 0..h {
34 tmp[r * w..(r + 1) * w].copy_from_slice(&filter_same(&img[r * w..(r + 1) * w], htaps));
35 }
36 let rv = vtaps.len() / 2;
37 let mut out = vec![0.0; w * h];
38 for r in 0..h {
39 let orow = &mut out[r * w..(r + 1) * w];
40 for (k, &t) in vtaps.iter().enumerate() {
41 let src = (r + k).saturating_sub(rv).min(h - 1);
42 simd::axpy(t, &tmp[src * w..(src + 1) * w], orow);
43 }
44 }
45 out
46}
47
48#[cfg(test)]
49mod tests {
50 use super::*;
51
52 #[test]
53 fn filters_match_naive() {
54 let x: Vec<f64> = (0..37).map(|i| (i as f64 * 0.5).sin()).collect();
55 let taps = [0.25, 0.5, 0.25];
56 let got = filter_valid(&x, &taps);
57 for (i, g) in got.iter().enumerate() {
58 let want: f64 = taps.iter().enumerate().map(|(k, t)| t * x[i + k]).sum();
59 assert!((g - want).abs() < 1e-12);
60 }
61 let same = filter_same(&x, &taps);
62 assert_eq!(same.len(), x.len());
63 assert!((same[0] - (0.75 * x[0] + 0.25 * x[1])).abs() < 1e-12);
64 let img = vec![1.0; 6 * 5];
65 let f = filter_2d(&img, 6, 5, &taps, &taps);
66 assert!(f.iter().all(|v| (v - 1.0).abs() < 1e-12));
67 }
68}