pub fn suppress_non_maximum_reference<I, C>(image: &I, radius: u32) -> Image<Luma<C>>
    where
        I: GenericImage<Pixel = Luma<C>>,
        C: Primitive + Ord,
    {
        let (width, height) = image.dimensions();
        let mut out = Image::new(width, height);
        out.copy_from(image, 0, 0).unwrap();

        let iradius = radius as i32;
        let iheight = height as i32;
        let iwidth = width as i32;

        // We update zero values from out as we go, so to check intensities
        // we need to read values from the input image.
        for y in 0..height {
            for x in 0..width {
                let intensity = image.get_pixel(x, y)[0];
                let mut is_max = true;

                let y_lower = cmp::max(0, y as i32 - iradius);
                let y_upper = cmp::min(y as i32 + iradius + 1, iheight);
                let x_lower = cmp::max(0, x as i32 - iradius);
                let x_upper = cmp::min(x as i32 + iradius + 1, iwidth);

                for py in y_lower..y_upper {
                    for px in x_lower..x_upper {
                        let v = image.get_pixel(px as u32, py as u32)[0];
                        // Handle intensity tiebreaks lexicographically
                        let candidate_is_lexically_earlier = (px as u32, py as u32) < (x, y);
                        if v > intensity || (v == intensity && candidate_is_lexically_earlier) {
                            is_max = false;
                            break;
                        }
                    }
                }

                if !is_max {
                    out.put_pixel(x, y, Luma([C::zero()]));
                }
            }
        }

        out
    }
Reference implementation of suppress_non_maximum. Used to validate
the (presumably faster) actual implementation.
src/suppress.rs