Skip to content

Commit

Permalink
Merge pull request #12 from imeka/morphological
Browse files Browse the repository at this point in the history
Ball and generic kernels
  • Loading branch information
nilgoyette committed Nov 8, 2022
2 parents cbcc026 + 76d32c9 commit 0a2623f
Show file tree
Hide file tree
Showing 3 changed files with 160 additions and 55 deletions.
59 changes: 55 additions & 4 deletions src/lib.rs
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@
//! The `ndarray-image` crate provides multidimensional image processing for `ArrayBase`,
//! the *n*-dimensional array data structure provided by [`ndarray`].

use ndarray::{Array, Array3, ArrayBase, Data, Dimension, Ix3, ShapeBuilder};
use ndarray::{arr3, Array, Array3, ArrayBase, ArrayView3, Data, Dimension, Ix3, ShapeBuilder};

mod filters;
mod interpolation;
Expand All @@ -30,12 +30,63 @@ pub use pad::{pad, pad_to, PadMode};
pub type Mask = Array3<bool>;

/// 3D common kernels. Also called Structuring Element.
#[derive(Clone, Copy, Debug, PartialEq)]
pub enum Kernel3d {
#[derive(Clone, Copy, PartialEq)]
pub enum Kernel3d<'a> {
/// Diamond/star kernel (center and sides).
///
/// Equivalent to `generate_binary_structure(3, 1)`.
Star,
/// 3x3 cube.
/// Ball kernel (center and sides).
///
/// Equivalent to `generate_binary_structure(3, 2)`.
Ball,
/// 3x3x3 cube.
///
/// Equivalent to `generate_binary_structure(3, 3)`.
Full,
/// Generic kernel of any 3D size.
///
/// The generic kernels are incredibly slower on all morphological operations.
Generic(ArrayView3<'a, bool>),
}

impl<'a> std::fmt::Debug for Kernel3d<'a> {
fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
match *self {
Kernel3d::Star => write!(f, "Star {:?}", self.dim()),
Kernel3d::Ball => write!(f, "Ball {:?}", self.dim()),
Kernel3d::Full => write!(f, "Full {:?}", self.dim()),
Kernel3d::Generic(k) => write!(f, "Generic {:?}", k.dim()),
}
}
}

impl<'a> Kernel3d<'a> {
/// Return the kernel dimension.
pub fn dim(&self) -> (usize, usize, usize) {
match *self {
Kernel3d::Star | Kernel3d::Ball | Kernel3d::Full => (3, 3, 3),
Kernel3d::Generic(k) => k.dim(),
}
}

/// Return the actual 3D array.
pub fn array(&self) -> Array3<bool> {
match *self {
Kernel3d::Star => arr3(&[
[[false, false, false], [false, true, false], [false, false, false]],
[[false, true, false], [true, true, true], [false, true, false]],
[[false, false, false], [false, true, false], [false, false, false]],
]),
Kernel3d::Ball => arr3(&[
[[false, true, false], [true, true, true], [false, true, false]],
[[true, true, true], [true, true, true], [true, true, true]],
[[false, true, false], [true, true, true], [false, true, false]],
]),
Kernel3d::Full => Array3::from_elem((3, 3, 3), true),
Kernel3d::Generic(k) => k.to_owned(),
}
}
}

/// Utilitary function that returns a new *n*-dimensional array of dimension `shape` with the same
Expand Down
125 changes: 75 additions & 50 deletions src/morphology.rs
Original file line number Diff line number Diff line change
Expand Up @@ -2,42 +2,54 @@ use ndarray::{s, ArrayBase, ArrayView3, ArrayViewMut3, Data, Ix3, Zip};

use crate::{array_like, dim_minus, Kernel3d, Mask};

impl Kernel3d {
impl<'a> Kernel3d<'a> {
/// Erode the mask when at least one of the values doesn't respect the kernel.
///
/// An erosion is defined either as `all(!(!w & k))` or `!any(!w & k)`. Thus, an empty kernel
/// will always produce a full mask.
#[rustfmt::skip]
fn erode(&self, from: ArrayView3<bool>, into: ArrayViewMut3<bool>) {
let windows = from.windows(self.dim());
match *self {
Kernel3d::Full => {
// TODO This is incredibly ugly but this is much faster (3x) than the Zip
//Zip::from(from.windows((3, 3, 3)))
// .map_assign_into(into, |w| Zip::from(w).all(|&m| m));
Zip::from(from.windows((3, 3, 3))).map_assign_into(into, |w| {
w[(0, 0, 0)] && w[(0, 0, 1)] && w[(0, 0, 2)]
&& w[(0, 1, 0)] && w[(0, 1, 1)] && w[(0, 1, 2)]
&& w[(0, 2, 0)] && w[(0, 2, 1)] && w[(0, 2, 2)]
&& w[(1, 0, 0)] && w[(1, 0, 1)] && w[(1, 0, 2)]
&& w[(1, 1, 0)] && w[(1, 1, 1)] && w[(1, 1, 2)]
&& w[(1, 2, 0)] && w[(1, 2, 1)] && w[(1, 2, 2)]
&& w[(2, 0, 0)] && w[(2, 0, 1)] && w[(2, 0, 2)]
&& w[(2, 1, 0)] && w[(2, 1, 1)] && w[(2, 1, 2)]
&& w[(2, 2, 0)] && w[(2, 2, 1)] && w[(2, 2, 2)]
})
}
Kernel3d::Star => Zip::from(from.windows((3, 3, 3))).map_assign_into(into, |w| {
Kernel3d::Full => Zip::from(windows).map_assign_into(into, |w| {
// This is incredibly ugly but this is much faster (3x) than the Zip
// |w| Zip::from(w).all(|&m| m)
w[(0, 0, 0)] && w[(0, 0, 1)] && w[(0, 0, 2)]
&& w[(0, 1, 0)] && w[(0, 1, 1)] && w[(0, 1, 2)]
&& w[(0, 2, 0)] && w[(0, 2, 1)] && w[(0, 2, 2)]
&& w[(1, 0, 0)] && w[(1, 0, 1)] && w[(1, 0, 2)]
&& w[(1, 1, 0)] && w[(1, 1, 1)] && w[(1, 1, 2)]
&& w[(1, 2, 0)] && w[(1, 2, 1)] && w[(1, 2, 2)]
&& w[(2, 0, 0)] && w[(2, 0, 1)] && w[(2, 0, 2)]
&& w[(2, 1, 0)] && w[(2, 1, 1)] && w[(2, 1, 2)]
&& w[(2, 2, 0)] && w[(2, 2, 1)] && w[(2, 2, 2)]
}),
Kernel3d::Ball => Zip::from(windows).map_assign_into(into, |w| {
w[(0, 0, 1)]
&& w[(0, 1, 0)] && w[(0, 1, 1)] && w[(0, 1, 2)]
&& w[(0, 2, 1)]
&& w[(1, 0, 0)] && w[(1, 0, 1)] && w[(1, 0, 2)]
&& w[(1, 1, 0)] && w[(1, 1, 1)] && w[(1, 1, 2)]
&& w[(1, 2, 0)] && w[(1, 2, 1)] && w[(1, 2, 2)]
&& w[(2, 0, 1)]
&& w[(2, 1, 0)] && w[(2, 1, 1)] && w[(2, 1, 2)]
&& w[(2, 2, 1)]
}),
Kernel3d::Star => Zip::from(windows).map_assign_into(into, |w| {
// This ugly condition is equivalent to
// *mask = !w.iter().zip(&kernel).any(|(w, k)| !w & k)
// but it's around 5x faster because there's no branch misprediction
w[(0, 1, 1)]
&& w[(1, 0, 1)]
&& w[(1, 1, 0)]
&& w[(1, 1, 1)]
&& w[(1, 1, 2)]
&& w[(1, 2, 1)]
&& w[(2, 1, 1)]
// but it's extremely faster because there's no branch misprediction
w[(0, 1, 1)]
&& w[(1, 0, 1)]
&& w[(1, 1, 0)] && w[(1, 1, 1)] && w[(1, 1, 2)]
&& w[(1, 2, 1)]
&& w[(2, 1, 1)]
}),
Kernel3d::Generic(kernel) => Zip::from(windows).map_assign_into(into, |w| {
// TODO Maybe use Zip::any when available
// !Zip::from(w).and(kernel).any(|(&w, &k)| !w & k)
Zip::from(w).and(kernel).all(|&w, &k| !(!w & k))
})
}
}

Expand All @@ -47,34 +59,47 @@ impl Kernel3d {
/// mask.
#[rustfmt::skip]
fn dilate(&self, from: ArrayView3<bool>, into: ArrayViewMut3<bool>) {
let windows = from.windows(self.dim());
match *self {
Kernel3d::Full => {
// TODO This is incredibly ugly but this is much faster (6x) than the Zip
//Zip::from(from.windows((3, 3, 3))).map_assign_into(into, |w| w.iter().any(|&w| w))
Zip::from(from.windows((3, 3, 3))).map_assign_into(into, |w| {
w[(0, 0, 0)] || w[(0, 0, 1)] || w[(0, 0, 2)]
|| w[(0, 1, 0)] || w[(0, 1, 1)] || w[(0, 1, 2)]
|| w[(0, 2, 0)] || w[(0, 2, 1)] || w[(0, 2, 2)]
|| w[(1, 0, 0)] || w[(1, 0, 1)] || w[(1, 0, 2)]
|| w[(1, 1, 0)] || w[(1, 1, 1)] || w[(1, 1, 2)]
|| w[(1, 2, 0)] || w[(1, 2, 1)] || w[(1, 2, 2)]
|| w[(2, 0, 0)] || w[(2, 0, 1)] || w[(2, 0, 2)]
|| w[(2, 1, 0)] || w[(2, 1, 1)] || w[(2, 1, 2)]
|| w[(2, 2, 0)] || w[(2, 2, 1)] || w[(2, 2, 2)]
})
}
Kernel3d::Star => Zip::from(from.windows((3, 3, 3))).map_assign_into(into, |w| {
Kernel3d::Full => Zip::from(windows).map_assign_into(into, |w| {
// This is incredibly ugly but this is much faster (6x) than the Zip
// |w| w.iter().any(|&w| w))
w[(0, 0, 0)] || w[(0, 0, 1)] || w[(0, 0, 2)]
|| w[(0, 1, 0)] || w[(0, 1, 1)] || w[(0, 1, 2)]
|| w[(0, 2, 0)] || w[(0, 2, 1)] || w[(0, 2, 2)]
|| w[(1, 0, 0)] || w[(1, 0, 1)] || w[(1, 0, 2)]
|| w[(1, 1, 0)] || w[(1, 1, 1)] || w[(1, 1, 2)]
|| w[(1, 2, 0)] || w[(1, 2, 1)] || w[(1, 2, 2)]
|| w[(2, 0, 0)] || w[(2, 0, 1)] || w[(2, 0, 2)]
|| w[(2, 1, 0)] || w[(2, 1, 1)] || w[(2, 1, 2)]
|| w[(2, 2, 0)] || w[(2, 2, 1)] || w[(2, 2, 2)]
}),
Kernel3d::Ball => Zip::from(windows).map_assign_into(into, |w| {
w[(0, 0, 1)]
|| w[(0, 1, 0)] || w[(0, 1, 1)] || w[(0, 1, 2)]
|| w[(0, 2, 1)]
|| w[(1, 0, 0)] || w[(1, 0, 1)] || w[(1, 0, 2)]
|| w[(1, 1, 0)] || w[(1, 1, 1)] || w[(1, 1, 2)]
|| w[(1, 2, 0)] || w[(1, 2, 1)] || w[(1, 2, 2)]
|| w[(2, 0, 1)]
|| w[(2, 1, 0)] || w[(2, 1, 1)] || w[(2, 1, 2)]
|| w[(2, 2, 1)]
}),
Kernel3d::Star => Zip::from(windows).map_assign_into(into, |w| {
// This ugly condition is equivalent to
// *mask = w.iter().zip(&kernel).any(|(w, k)| w & k)
// |(w, k)| w & k
// but it's around 5x faster because there's no branch misprediction
w[(0, 1, 1)]
|| w[(1, 0, 1)]
|| w[(1, 1, 0)]
|| w[(1, 1, 1)]
|| w[(1, 1, 2)]
|| w[(1, 2, 1)]
|| w[(2, 1, 1)]
w[(0, 1, 1)]
|| w[(1, 0, 1)]
|| w[(1, 1, 0)] || w[(1, 1, 1)] || w[(1, 1, 2)]
|| w[(1, 2, 1)]
|| w[(2, 1, 1)]
}),
Kernel3d::Generic(kernel) => Zip::from(windows).map_assign_into(into, |w| {
// TODO Use Zip::any when available
// Zip::from(w).and(kernel).any(|idx, &w, &k| w & k)
w.iter().zip(&kernel).any(|(w, k)| w & k)
})
}
}
}
Expand Down
31 changes: 30 additions & 1 deletion tests/morphology.rs
Original file line number Diff line number Diff line change
Expand Up @@ -49,10 +49,27 @@ fn test_binary_erosion_hole() {
gt.slice_mut(s![5, 5, 4..7]).fill(false);

assert_eq!(gt, binary_erosion(&mask, Kernel3d::Star, 1));
assert_eq!(gt, binary_erosion(&mask, Kernel3d::Generic(Kernel3d::Star.array().view()), 1));
}

#[test] // Results verified with the `binary_erosion` function from SciPy. (v1.7.0)
fn test_binary_erosion_cube_kernel() {
fn test_binary_erosion_ball_kernel() {
let mut mask = Mask::from_elem((11, 11, 11), true);
mask[(5, 5, 5)] = false;

let mut gt = Mask::from_elem((11, 11, 11), false);
let (width, height, depth) = dim_minus(&mask, 1);
gt.slice_mut(s![1..width, 1..height, 1..depth]).fill(true);
// Remove the ball shape in the image center.
gt.slice_mut(s![4..7, 4..7, 4..7]).fill(false);
gt.slice_mut(s![4..7; 2, 4..7; 2, 4..7; 2]).fill(true);

assert_eq!(gt, binary_erosion(&mask, Kernel3d::Ball, 1));
assert_eq!(gt, binary_erosion(&mask, Kernel3d::Generic(Kernel3d::Ball.array().view()), 1));
}

#[test] // Results verified with the `binary_erosion` function from SciPy. (v1.7.0)
fn test_binary_erosion_full_kernel() {
let mut mask = Mask::from_elem((11, 11, 11), true);
mask[(5, 5, 5)] = false;

Expand All @@ -63,6 +80,7 @@ fn test_binary_erosion_cube_kernel() {
gt.slice_mut(s![4..7, 4..7, 4..7]).fill(false);

assert_eq!(gt, binary_erosion(&mask, Kernel3d::Full, 1));
assert_eq!(gt, binary_erosion(&mask, Kernel3d::Generic(Kernel3d::Full.array().view()), 1));
}

#[test] // Results verified with the `binary_dilation` function from SciPy. (v1.7.0)
Expand Down Expand Up @@ -98,18 +116,29 @@ fn test_binary_dilation_plain() {
gt.slice_mut(s![2..w - 1, h - 1, 2..d - 1]).fill(true);

assert_eq!(gt, binary_dilation(&mask.view(), Kernel3d::Star, 1));
assert_eq!(gt, binary_dilation(&mask, Kernel3d::Generic(Kernel3d::Star.array().view()), 1));

let mut mask = Mask::from_elem((w, h, d), false);
mask.slice_mut(s![4, 4, 4..]).fill(true);
let mut gt = Mask::from_elem((w, h, d), false);
gt.slice_mut(s![3..6, 3..6, 3..]).fill(true);
gt.slice_mut(s![3..6; 2, 3..6; 2, 3]).fill(false);
assert_eq!(gt, binary_dilation(&mask.view(), Kernel3d::Ball, 1));
assert_eq!(gt, binary_dilation(&mask, Kernel3d::Generic(Kernel3d::Ball.array().view()), 1));

let mut mask = Mask::from_elem((w, h, d), false);
mask[(4, 4, 4)] = true;
let mut gt = Mask::from_elem((w, h, d), false);
gt.slice_mut(s![2.., 2.., 2..]).fill(true);
assert_eq!(gt, binary_dilation(&mask.view(), Kernel3d::Full, 2));
assert_eq!(gt, binary_dilation(&mask, Kernel3d::Generic(Kernel3d::Full.array().view()), 2));

let mut mask = Mask::from_elem((w, h, d), false);
mask[(4, 5, 5)] = true;
let mut gt = Mask::from_elem((w, h, d), false);
gt.slice_mut(s![1.., 2.., 2..]).fill(true);
assert_eq!(gt, binary_dilation(&mask.view(), Kernel3d::Full, 3));
assert_eq!(gt, binary_dilation(&mask, Kernel3d::Generic(Kernel3d::Full.array().view()), 3));
}

#[test] // Results verified with the `binary_dilation` function from SciPy. (v1.7.0)
Expand Down

0 comments on commit 0a2623f

Please sign in to comment.