From d5a02aec196b77202f62d94a28d0a69bb086f65a Mon Sep 17 00:00:00 2001 From: "Christopher H. Jordan" Date: Mon, 6 Jan 2020 10:22:21 +0800 Subject: [PATCH] Add an example for 1D convolution. --- Cargo.toml | 1 + convolution_1d/Cargo.toml | 11 ++++ convolution_1d/README.md | 19 +++++++ convolution_1d/src/lib.rs | 107 +++++++++++++++++++++++++++++++++++++ convolution_1d/src/main.rs | 16 ++++++ 5 files changed, 154 insertions(+) create mode 100644 convolution_1d/Cargo.toml create mode 100644 convolution_1d/README.md create mode 100644 convolution_1d/src/lib.rs create mode 100644 convolution_1d/src/main.rs diff --git a/Cargo.toml b/Cargo.toml index 9335a04..b114f4c 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -3,4 +3,5 @@ members = [ "./", "linear_regression", "k_means", + "convolution_1d", ] \ No newline at end of file diff --git a/convolution_1d/Cargo.toml b/convolution_1d/Cargo.toml new file mode 100644 index 0000000..40bd645 --- /dev/null +++ b/convolution_1d/Cargo.toml @@ -0,0 +1,11 @@ +[package] +name = "convolution_1d" +version = "0.1.0" +authors = ["Christopher H. Jordan "] +edition = "2018" + +# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html + +[dependencies] +float-cmp = "0.6.0" +ndarray = "0.13.0" diff --git a/convolution_1d/README.md b/convolution_1d/README.md new file mode 100644 index 0000000..02ec964 --- /dev/null +++ b/convolution_1d/README.md @@ -0,0 +1,19 @@ +1D Convolution +================== + +An implementation of traditional (i.e. not utilising FFTs) convolution with +ndarray. Designed to match `numpy`'s `convolve` with `mode='same'`. + +The code could be made to work for multi-dimensional arrays, but at present, I +have not investigated how. + +``` +numpy example of convolve with mode='same': +>>> np.convolve([1,2,3],[0,1,0.5], 'same') +array([1. , 2.5, 4. ]) + +ndarray code follows: +data: [1.0, 2.0, 3.0], shape=[3], strides=[1], layout=C | F (0x3), const ndim=1 +window: [0.0, 1.0, 0.5], shape=[3], strides=[1], layout=C | F (0x3), const ndim=1 +convolution: [1.0, 2.5, 4.0], shape=[3], strides=[1], layout=C | F (0x3), const ndim=1 +``` diff --git a/convolution_1d/src/lib.rs b/convolution_1d/src/lib.rs new file mode 100644 index 0000000..8578086 --- /dev/null +++ b/convolution_1d/src/lib.rs @@ -0,0 +1,107 @@ +#[macro_use] +extern crate ndarray; + +use ndarray::prelude::*; + +pub fn convolve(data: ArrayView1, window: ArrayView1) -> Array1 { + let padded = stack![ + Axis(0), + Array1::zeros(window.len() / 2), + data, + Array1::zeros(window.len() / 2) + ]; + let mut w = window.view(); + w.invert_axis(Axis(0)); + + padded + .windows(w.len()) + .into_iter() + .map(|x| (&x * &w).sum()) + .collect() +} + +#[cfg(test)] +mod tests { + use super::*; + use float_cmp::*; + + #[test] + fn convolve_odd_odd() { + let data = array![1., 2., 3.]; + let window = array![0., 1., 0.5]; + let expected = array![1., 2.5, 4.]; + + for (exp, res) in expected.iter().zip(&convolve(data.view(), window.view())) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } + + #[test] + fn convolve_odd_odd2() { + let data = array![1., 2., 3., 4., 5.]; + let window = array![2., 1., 0., 1., 0.5]; + let result = convolve(data.view(), window.view()); + let expected = array![8., 12., 16.5, 9., 5.5]; + + for (exp, res) in expected.iter().zip(&result) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } + + #[test] + fn convolve_even_odd() { + let data = array![1., 2., 3., 4.]; + let window = array![0., 1., 0.5]; + let expected = array![1., 2.5, 4., 5.5]; + + for (exp, res) in expected.iter().zip(&convolve(data.view(), window.view())) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } + + #[test] + fn convolve_even_even() { + let data = array![1., 2., 3., 4.]; + let window = array![1., 0.5]; + let expected = array![1., 2.5, 4., 5.5]; + + for (exp, res) in expected.iter().zip(&convolve(data.view(), window.view())) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } + + #[test] + fn convolve_even_even2() { + let data = array![1., 2., 3., 4.]; + let window = array![1., 0., 1., 0.5]; + let result = convolve(data.view(), window.view()); + let expected = array![2., 4., 6.5, 4.]; + + for (exp, res) in expected.iter().zip(&result) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } + + #[test] + fn convolve_odd_even() { + let data = array![1., 2., 3., 4., 5.]; + let window = array![1., 0.5]; + let expected = array![1., 2.5, 4., 5.5, 7.]; + + for (exp, res) in expected.iter().zip(&convolve(data.view(), window.view())) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } + + #[test] + fn convolve_bigger_window() { + let data = array![1., 2., 3.]; + let window = array![1., 0., 1., 0.5]; + let result = convolve(data.view(), window.view()); + let expected = array![2., 4., 2.5, 4.]; + + for (exp, res) in expected.iter().zip(&result) { + assert!(approx_eq!(f64, *exp, *res, ulps = 2)); + } + } +} diff --git a/convolution_1d/src/main.rs b/convolution_1d/src/main.rs new file mode 100644 index 0000000..05535f7 --- /dev/null +++ b/convolution_1d/src/main.rs @@ -0,0 +1,16 @@ +use ndarray::prelude::*; +use convolution_1d::convolve; + +fn main() { + println!(r#"numpy example of convolve with mode='same': +>>> np.convolve([1,2,3],[0,1,0.5], 'same') +array([1. , 2.5, 4. ]) + +ndarray code follows:"#); + + let data = array![1., 2., 3.]; + let window = array![0., 1., 0.5]; + println!("data: \t\t{:?}", data); + println!("window: \t{:?}", window); + println!("convolution: \t{:?}", convolve(data.view(), window.view())); +}