Convolutions in 1d in rust

CONTRIBUTORS.md

@@ -64,4 +64,6 @@ This file lists everyone, who contributed to this repo and wanted to show up her
 - K. Shudipto Amin
 - Peanutbutter_Warrior
 - Thijs Raymakers 
+- eengamer2007
 - Michael Ciccotosto-Camp
+

contents/convolutions/1d/1d.md

@@ -58,6 +58,8 @@ With this in mind, we can almost directly transcribe the discrete equation into
 [import:63-84, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import:20-31, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import:69-83, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 The easiest way to reason about this code is to read it as you might read a textbook.
@@ -193,6 +195,8 @@ Here it is again for clarity:
 [import:63-84, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import:20-31, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import:69-83, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 Here, the main difference between the bounded and unbounded versions is that the output array size is smaller in the bounded case.
@@ -205,6 +209,8 @@ For an unbounded convolution, the function would be called with a the output arr
 [import:96-97, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import:41-42, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import:12-12, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 On the other hand, the bounded call would set the output array size to simply be the length of the signal
@@ -216,6 +222,8 @@ On the other hand, the bounded call would set the output array size to simply be
 [import:98-99, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import:44-45, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import:13-13, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 Finally, as we mentioned before, it is possible to center bounded convolutions by changing the location where we calculate the each point along the filter.
@@ -228,6 +236,8 @@ This can be done by modifying the following line:
 [import:71-71, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import:25-25, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import:74-74, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 Here, `j` counts from `i-length(filter)` to `i`.
@@ -264,6 +274,8 @@ In code, this typically amounts to using some form of modulus operation, as show
 [import:38-61, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import:5-17, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import:85-99, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 This is essentially the same as before, except for the modulus operations, which allow us to work on a periodic domain.
@@ -283,6 +295,8 @@ For the code associated with this chapter, we have used the convolution to gener
 [import, lang:"csharp"](code/csharp/1DConvolution.cs)
 {% sample lang="py" %}
 [import, lang:"python"](code/python/1d_convolution.py)
+{% sample lang="rs" %}
+[import, lang:"rust"](code/rust/1d_convolution.rs)
 {% endmethod %}
 
 At a test case, we have chosen to use two sawtooth functions, which should produce the following images:

contents/convolutions/1d/code/rust/1d_convolution.rs

@@ -0,0 +1,99 @@
+use std::fs::File;
+use std::io::Write;
+use std::cmp::max;
+
+fn main() {
+    let x = normalize(create_sawtooth(200));
+    let y = normalize(create_sawtooth(200));
+
+    let len_x = x.len();
+    let len_y = y.len();
+
+    let full_linear_output = convolve_linear(&x, &y, len_x + len_y - 1);
+    let simple_linear_output = convolve_linear(&x, &y, len_x);
+    let cyclic_output = convolve_cyclic(&x, &y);
+
+    // Save the convolutions to plot them.
+    // The way I do it is a little weird but it is to store the data the same way as the other programs
+    let mut full_file = File::create("full_linear.dat").unwrap();
+    for i in full_linear_output {
+        writeln!(full_file, "{i}").unwrap();
+    }
+    full_file.sync_all().unwrap();
+
+    let mut simple_file = File::create("simple_linear.dat").unwrap();
+    for i in simple_linear_output {
+        writeln!(simple_file, "{i}").unwrap();
+    }
+    simple_file.sync_all().unwrap();
+
+    let mut cyclic_file = File::create("cyclic.dat").unwrap();
+    for i in cyclic_output {
+        writeln!(cyclic_file, "{i}").unwrap();
+    }
+    cyclic_file.sync_all().unwrap();
+}
+
+// Generates a sawtooth function with a given length.
+fn create_sawtooth(length: usize) -> Vec<f64> {
+    let mut array: Vec<f64> = Vec::with_capacity(length);
+    for i in 0..length {
+        array.push((i+1) as f64 / length as f64); // divide by length for normalization
+    }
+    array
+}
+
+// Normalizes the given array.
+fn normalize(array: Vec<f64>) -> Vec<f64> {
+    let norm = norm(&array);
+    let mut output: Vec<f64> = Vec::with_capacity(array.len());
+
+    for value in array {
+        output.push(value / norm);
+    }
+
+    output
+}
+
+// Calculates the norm of an array.
+fn norm(array: &[f64]) -> f64 {
+    let sum = array.iter().map(|i| i * i).sum::<f64>();
+    sum.sqrt()
+}
+
+
+// Modulus function that handles negative values correctly.
+// Assumes that y >= 0.
+fn modulus(x: isize, y: usize) -> usize {((x%y as isize) as usize + y) % y}
+
+fn convolve_linear(signal: &Vec<f64>, filter: &Vec<f64>, output_size: usize) -> Vec<f64> {
+    let mut output = Vec::with_capacity(output_size);
+
+    for i in 0..(output_size as isize) {
+        let mut sum: f64 = 0.;
+        for j in max(0, i - filter.len() as isize)..=i {
+            if j < signal.len() as isize && (i - j) < filter.len() as isize {
+                sum += signal[j as usize] * filter[(i-j) as usize];
+            }
+        }
+        output.push(sum);
+    }
+
+    output
+}
+
+fn convolve_cyclic(signal: &Vec<f64>, filter: &Vec<f64>) -> Vec<f64> {
+    let output_size = max(signal.len(), filter.len());
+
+    let mut output: Vec<f64> = Vec::with_capacity(output_size);
+    for i in 0..output_size {
+        let mut sum: f64 = 0.;
+        for j in 0..output_size {
+            if modulus((i - j) as isize, output_size) < filter.len() {
+                sum += signal[modulus((j - 1) as isize, output_size)] * filter[modulus((i - j) as isize, output_size)];
+            }
+        }
+        output.push(sum);
+    }
+    output
+}

contents/convolutions/1d/code/rust/Cargo.toml

@@ -0,0 +1,12 @@
+[package]
+name = "rust"
+version = "0.1.0"
+edition = "2021"
+
+# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html
+
+[dependencies]
+
+[[bin]]
+path = "./1d_convolution.rs"
+name = "main"