simplecv/filter/
mod.rs

1//! 2D filters for image processing.
2
3use ndarray::prelude::*;
4use ndarray::{Data, DataMut};
5
6/// Representing the border type for filters.
7///
8/// Following border types are supported:
9/// * Constant(v): a constant, i.e., vvvv|abcdefgh|vvvv
10/// * Reflect: reflect the image, i.e., edcb|abcdefgh|gfed
11/// * Replicate: copy the value at the border, i.e., aaaa|abcdefgh|hhhh
12#[derive(Copy, Clone)]
13pub enum BorderType {
14    Constant(f64),
15    Reflect,
16    Replicate
17}
18
19/// Compute the source location of the outside point.
20///
21/// This function is used by `filter`. 
22/// For example, when the border type is `Reflect`, 
23/// ```
24/// use simplecv::filter::*;
25/// let nx = border_interpolate(-2, 10, BorderType::Reflect).unwrap();
26/// assert_eq!(nx, 2);
27/// ```
28/// 
29/// The function return None when border type is Constant.
30pub fn border_interpolate(p:i32, len:usize, border: BorderType) -> Option<usize> {
31    fn abs(a: i32) -> i32{
32        if a < 0 {
33            -a
34        } else {
35            a
36        }
37    }
38    match border{
39        BorderType::Constant(_) => None,
40        BorderType::Reflect => Some(abs(p % (len as i32)) as usize),
41        BorderType::Replicate => {
42            if p < 0 {
43                Some(0usize)
44            } else {
45                Some(len - 1)
46            }
47        }
48    }
49}
50
51/// Get the value of an image at a location which may be outside the image.
52///
53/// This function is used by `filter`.
54fn access_img_border<S>(src: &ArrayBase<S, Ix2>, x:i32, y:i32, border: BorderType) -> f64 
55    where S:Data<Elem=f64> 
56{
57    if x >= 0 && y >= 0 && x < src.shape()[0] as i32 && y < src.shape()[1] as i32 {
58        src[[x as usize, y as usize]]
59    }
60    else {
61        match border {
62            BorderType::Constant(v) => v,
63            BorderType::Reflect|BorderType::Replicate => {
64                let nx = border_interpolate(x, src.shape()[0], border).unwrap();
65                let ny = border_interpolate(y, src.shape()[1], border).unwrap();
66                src[[nx, ny]]
67            }
68        }
69    }
70}
71
72/// Apply a linear filter to the source image.
73///
74/// `out` must have the exactly same shape of `src`. Both the `src` and `kernel` 
75/// should be 2D array. For more channels you may need to write a wrapper by yourself.
76///
77/// The method of dealing with border situation is selected by `border`. By setting
78/// `border=Reflect`, you will get the default result of OpenCV.
79///
80pub fn filter_<S, T, K>(src: &ArrayBase<S, Ix2>, kernel: &ArrayBase<K, Ix2>, border: BorderType, 
81               out:&mut ArrayBase<T, Ix2>)
82    where S: Data<Elem=f64>, T: DataMut<Elem=f64>, K: Data<Elem=f64>
83{
84    let kh = kernel.shape()[0];
85    let kw = kernel.shape()[1];
86    let kcx = (kh / 2) as i32; // kernel center x
87    let kcy = (kw / 2) as i32; // kernel center x
88    let height = src.shape()[0];
89    let width = src.shape()[1];
90    for i in 0..height {
91        for j in 0..width {
92            let mut val = 0.0f64;
93            for ki in 0..kh {
94                for kj in 0..kw {
95                    let sx = i as i32 + ki as i32 - kcx;
96                    let sy = j as i32 + kj as i32 - kcy;
97                    let sval = access_img_border(src, sx, sy, border);
98                    val = val + sval * kernel[[ki, kj]];
99                }
100            }
101            out[[i, j]] = val;
102        }
103    }
104}
105
106/// Apply a linear filter to the source image.
107///
108/// The method of dealing with border situation is selected by `border`. By setting
109/// `border=Reflect`, you will get the default result of OpenCV.
110/// 
111/// # Example
112/// ```
113///  use simplecv::filter::*;
114///  use ndarray::arr2;
115///  let A = arr2(&[[0.0, 0.0, 0.0, 0.0, 0.0],
116///                 [0.0, 1.0, 1.0, 1.0, 0.0],
117///                 [0.0, 1.0, 1.0, 1.0, 0.0],
118///                 [0.0, 1.0, 1.0, 1.0, 0.0],
119///                 [0.0, 0.0, 0.0, 0.0, 0.0]]);
120///  let kernel = arr2(&[[1.0, 2.0, 1.0],
121///                      [2.0, 4.0, 2.0],
122///                      [1.0, 2.0, 1.0]]);
123///  let target = arr2(&[[1.0, 3.0, 4.0, 3.0, 1.0],
124///                      [3.0, 9.0, 12.0, 9.0, 3.0],
125///                      [4.0, 12.0, 16.0, 12.0, 4.0],
126///                      [3.0, 9.0, 12.0, 9.0, 3.0],
127///                      [1.0, 3.0, 4.0, 3.0, 1.0]]);
128///  let output = filter(&A, &kernel, BorderType::Constant(0.0));
129///  assert_eq!(target, output);
130/// ```
131///
132pub fn filter<S, K>(src: &ArrayBase<S, Ix2>, kernel: &ArrayBase<K, Ix2>, border: BorderType) -> Array<f64, Ix2>
133    where S: Data<Elem=f64>, K: Data<Elem=f64>
134{
135    let shape = src.shape();
136    let height = shape[0] as usize;
137    let width = shape[1] as usize;
138    let mut buffer = Array::zeros((height, width));
139    filter_(src, kernel, border, &mut buffer);
140    buffer
141}
142
143/// Generate a Gaussian kernel with the simplest method.
144pub fn gaussian_kernel_generator(ksize: usize) -> Array<f64, Ix2>{
145    fn sqr_dis(dx:i32, dy:i32) -> i32{
146        dx * dx + dy * dy
147    }
148    let cx = (ksize / 2) as i32;
149    let cy = (ksize / 2) as i32;
150    let mut kernel = Array::zeros((ksize, ksize));
151    for x in 0..ksize {
152        for y in 0..ksize{
153            let dist = sqr_dis(cx - x as i32, cy - y as i32);
154            kernel[[x, y]] = -dist as f64 / 2.0;
155        }
156    }
157    kernel.map_inplace(|x| *x = x.exp());
158    let normalized_constant = kernel.fold(0.0, |acc, x| acc + x);
159    kernel /= normalized_constant;
160    kernel
161}
162
163/// Smooth the image with a gaussian kernel.
164///
165/// The output buffer should be allocated by users.
166/// * `ksize`: is the kernel size. 
167/// * `border`: how to deal with the border.
168pub fn gaussian_smooth_<S, T>(src: &ArrayBase<S, Ix2>, ksize:usize, border: BorderType, out:&mut ArrayBase<T, Ix2>) 
169    where S: Data<Elem=f64>, T:DataMut<Elem=f64>
170{
171    let kernel = gaussian_kernel_generator(ksize);
172    filter_(src, &kernel, border, out);
173}
174
175/// Smooth the image with a gaussian kernel.
176///
177/// * `ksize`: is the kernel size. 
178/// * `border`: how to deal with the border.
179pub fn gaussian_smooth<S>(src: &ArrayBase<S, Ix2>, ksize:usize, border: BorderType) -> Array<f64, Ix2>
180    where S: Data<Elem=f64>
181{
182    let kernel = gaussian_kernel_generator(ksize);
183    filter(src, &kernel, border)
184}
185
186/// Smooth the image with a mean kernel.
187///
188/// The output buffer should be allocated by users.
189/// * `ksize`: is the kernel size. 
190/// * `border`: how to deal with the border.
191pub fn mean_smooth_<S, T>(src: &ArrayBase<S, Ix2>, ksize:usize, border:BorderType, out: &mut ArrayBase<T, Ix2>)
192    where S: Data<Elem=f64>, T:DataMut<Elem=f64>
193{
194    let kernel = Array::ones((ksize, ksize)) / ((ksize * ksize) as f64);
195    filter_(src, &kernel, border, out);
196}
197
198/// Smooth the image with a mean kernel.
199///
200/// * `ksize`: is the kernel size. 
201/// * `border`: how to deal with the border.
202///
203/// # Example
204/// ```
205/// let a = ndarray::arr2(&[[0.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 0.0]]);
206/// let smoothed = simplecv::filter::mean_smooth(&a, 3,
207///                     simplecv::filter::BorderType::Constant(0.0));
208/// let diff = simplecv::utils::max_diff(&smoothed, &(ndarray::Array::ones((3, 3)) / 9.0));
209/// assert!(diff < 1e-4);
210/// ```
211///
212pub fn mean_smooth<S>(src: &ArrayBase<S, Ix2>, ksize: usize, border:BorderType) -> Array<f64, Ix2>
213    where S: Data<Elem=f64>
214{
215    let kernel = Array::ones((ksize, ksize)) / ((ksize * ksize) as f64);
216    filter(src, &kernel, border)
217}
218
219/// Sobel operator implementation.
220///
221/// The output buffer should be allocated by users.
222///
223/// When kernel size is 3, the classical Sobel filter is applied. Read 
224/// [OpenCV Sobel()](https://docs.opencv.org/3.4/d4/d86/group__imgproc__filter.html#gacea54f142e81b6758cb6f375ce782c8d)
225/// for more details. Only `dx=1, dy=0` and `dx=0, dy=1` are supported now.
226///
227/// * `ksize`: the kernel size. Currently, only `ksize=3` is supported.
228/// * `dx`: order of the derivative x.
229/// * `dy`: order of the derivative y.
230/// * `border`: border type.
231pub fn sobel_<S, T>(src: &ArrayBase<S, Ix2>, ksize: usize, dx: u32, dy: u32, border: BorderType, out: &mut ArrayBase<T, Ix2>)
232    where S: Data<Elem=f64>, T: DataMut<Elem=f64>
233{
234    assert!(ksize==3, "Only ksize=3 is supported in sobel_() now.");
235    assert!(dx + dy == 1, "Only first order gradient of one direction is supported in sobel_() now.");
236    let kernel = arr2(&[[-1.0, 0.0, 1.0], [-2.0, 0.0, 2.0], [-1.0, 0.0, 1.0]]);
237    if dx == 1 {
238        filter_(src, &kernel, border, out);
239    } else{
240        filter_(src, &kernel.t(), border, out);
241    }
242}
243
244/// Sobel operator implementation.
245///
246/// When kernel size is 3, the classical Sobel filter is applied. Read 
247/// [OpenCV Sobel()](https://docs.opencv.org/3.4/d4/d86/group__imgproc__filter.html#gacea54f142e81b6758cb6f375ce782c8d)
248/// for more details. Only `dx=1, dy=0` and `dx=0, dy=1` are supported now.
249///
250/// * `ksize`: the kernel size. Currently, only `ksize=3` is supported.
251/// * `dx`: order of the derivative x.
252/// * `dy`: order of the derivative y.
253/// * `border`: border type.
254pub fn sobel<S>(src: &ArrayBase<S, Ix2>, ksize: usize, dx: u32, dy: u32, border: BorderType) -> Array<f64, Ix2>
255    where S: Data<Elem=f64>
256{
257    let mut buffer = Array::zeros((src.shape()[0], src.shape()[1]));
258    sobel_(src, ksize, dx, dy, border, &mut buffer);
259    buffer
260}
261
262/// Get the norm of image processed by a Sobel operation.
263///
264/// Currently `norm=-1, 1, 2` are supported, where -1 means the infinty norm (max of absolute value). 
265/// This function can be used to obtain the edge of original image.
266///
267/// This function is implemented by `sobel()`. First order derivative of x and y direction are 
268/// used for computing the gradient.
269/// * `ksize`: the kernel size.
270/// * `norm`: the norm used for computation.
271/// * `border`: border type.
272
273pub fn sobel_norm<S>(src: &ArrayBase<S, Ix2>, ksize: usize, norm: i32, border: BorderType) -> Array<f64, Ix2>
274    where S: Data<Elem=f64>
275{
276    let gx = sobel(src, ksize, 1, 0, border);
277    let gy = sobel(src, ksize, 0, 1, border);
278    let gnorm = match norm {
279        2 => {
280            let t = gx.mapv(|x| x.powi(2)) + gy.mapv(|x| x.powi(2));
281            t.mapv(f64::sqrt)
282        }
283        1 => { gx.mapv(|x| x.abs()) + gy.mapv(|x| x.abs()) }
284        -1 => {
285            let mut buffer = Array::zeros((src.shape()[0], src.shape()[1]));
286            for x in 0..src.shape()[0] {
287                for y in 0..src.shape()[1] {
288                    buffer[[x, y]] = super::utils::max(gx[[x, y]].abs(), gy[[x, y]].abs());
289                }
290            }
291            buffer
292        }
293        _ => {
294            panic!(format!("Norm = {} is not supported by sobel_norm()!", norm));
295        }
296    };
297    gnorm
298}
299
300pub mod canny;
301pub use canny::canny_edge;
302pub use canny::canny_edge_;