pub fn euclidean_squared_distance_transform(image: &Image<Luma<u8>>) -> Image<Luma<f64>> {
    let (width, height) = image.dimensions();
    let mut result = Image::new(width, height);
    let mut column_envelope = LowerEnvelope::new(height as usize);

    // Compute 1d transforms of each column
    for x in 0..width {
        let source = Column { image, column: x };
        let mut sink = ColumnMut {
            image: &mut result,
            column: x,
        };
        distance_transform_1d_mut(&source, &mut sink, &mut column_envelope);
    }

    let mut row_buffer = vec![0f64; width as usize];
    let mut row_envelope = LowerEnvelope::new(width as usize);

    // Compute 1d transforms of each row
    for y in 0..height {
        for x in 0..width {
            row_buffer[x as usize] = result.get_pixel(x, y)[0];
        }
        let mut sink = Row {
            image: &mut result,
            row: y,
        };
        distance_transform_1d_mut(&row_buffer, &mut sink, &mut row_envelope);
    }

    result
}
Computes the square of the `L2` (Euclidean) distance transform of `image`. Distances are to the
nearest foreground pixel, where a pixel is counted as foreground if it has non-zero value.

Uses the algorithm from [Distance Transforms of Sampled Functions] to achieve time linear
in the size of the image.

[Distance Transforms of Sampled Functions]: https://www.cs.cornell.edu/~dph/papers/dt.pdf
src/distance_transform.rs