A fast, iterative, correct approach to Stackblur, resulting in a very smooth and high-quality output, with no edge bleeding

So I'm using your crate in my app and hitting an issue while trying to blur a gdk_pixbuf::Pixbufs, that is for certain images the blurred result is grayed out. I believe that it is related to the size of the image, but I'm not entirely too sure. This is an example that exhibits this issue

use gdk_pixbuf::{Colorspace, Pixbuf};
use stackblur_iter::{blur_argb, imgref::ImgRefMut};

fn main() {
    let pixbuf = Pixbuf::from_file("unknown.png").unwrap();
    assert!(pixbuf.n_channels() == 3); // This assertion exists because the Pixbuf I'm blurring in my program should always lack an alpha channel

    let pixels = pixbuf.pixel_bytes().unwrap().to_vec();
    let width = pixbuf.width() as usize;
    let height = pixbuf.height() as usize;

    let mut pixels = match pixbuf.n_channels() {
        3 => blur_rgb(pixels, width, height, 5),
        _ => unreachable!(),
    };

    Pixbuf::from_mut_slice(
        &mut pixels,
        Colorspace::Rgb,
        false,
        8,
        pixbuf.width(),
        pixbuf.height(),
        pixbuf.rowstride(),
    )
    .savev("blurred.png", "png", &[])
    .unwrap();
}

fn blur_rgb(pixels: Vec<u8>, width: usize, height: usize, radius: usize) -> Vec<u8> {
    assert!(
        pixels.len() % 3 == 0,
        "The pixel buffer's length should be a multiple of 3, but it was '{}'.",
        pixels.len()
    );

    let mut pixels = pixels
        .chunks_exact(3)
        .map(|chunk| {
            let red = chunk[0];
            let green = chunk[1];
            let blue = chunk[2];

            u32::from_be_bytes([0xff, red, green, blue])
        })
        .collect::<Vec<_>>();

    let mut img = ImgRefMut::new(&mut pixels, width, height);

    blur_argb(&mut img, radius);

    pixels
        .into_iter()
        .flat_map(|pixel| {
            let [_, r, g, b] = pixel.to_be_bytes();

            [r, g, b]
        })
        .collect()
}

[package]
name = "stackblurtest"
version = "0.1.0"
edition = "2021"

# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html

[dependencies]
stackblur-iter = { git = "https://github.com/LoganDark/stackblur-iter/", rev = "53be4ee685f2bc1fbbb79b9ab02ee6679de85694" }
gdk-pixbuf = "0.15.11"

The Pixbuf mentioned above stores its data in an [rgb][rgb][rgb]... format and I'm converting it to ARGB u32s as that's pretty convenient because it allows me to use the blur_argb method provided by your crate and not have to write my own RGB type.

For the record I don't believe this to be an issue with me using 0xff for every alpha byte in the intermediary buffer as I replaced it with a randomly generated byte using the rand crate and the output didn't look like it changed.

This snippet fails due to undefined behavior when run under miri:

use stackblur_iter::imgref::ImgRefMut;

fn main() {
    let mut buf = vec![
        0, 0, 0, 0,   0, 0, 0, 0,   0, 0, 0, 0,
        0, 0, 0, 0,   0, 0, 0, 0,   0, 0, 0, 0,
        0, 0, 0, 0,   0, 0, 0, 0,   0, 0, 0, 0,
    ];

    let mut img_ref = ImgRefMut::new(&mut buf, 3, 3);
    stackblur_iter::blur_argb(&mut img_ref, 1);
}

stackblur-iter = "0.1.5" rustc --version: rustc 1.64.0-nightly (f8588549c 2022-07-18)

Error:

$ cargo miri run
    Finished dev [unoptimized + debuginfo] target(s) in 0.13s
     Running `/Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/bin/cargo-miri /Users/marko/.target-cargo/miri/x86_64-apple-darwin/debug/test-stackblur`
error: Undefined Behavior: trying to reborrow from <3782> for SharedReadOnly permission at alloc1660[0x0], but that tag does not exist in the borrow stack for this location
   --> /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/imgref-iter-0.3.3/src/iter/windows/ptr.rs:135:13
    |
135 |         let len = (*self.0).len();
    |                   ^^^^^^^^^^^^^^^
    |                   |
    |                   trying to reborrow from <3782> for SharedReadOnly permission at alloc1660[0x0], but that tag does not exist in the borrow stack for this location
    |                   this error occurs as part of a reborrow at alloc1660[0x0..0xc]
    |
    = help: this indicates a potential bug in the program: it performed an invalid operation, but the Stacked Borrows rules it violated are still experimental
    = help: see https://github.com/rust-lang/unsafe-code-guidelines/blob/master/wip/stacked-borrows.md for further information
help: <3782> was created by a retag at offsets [0x0..0x90]
   --> src/main.rs:11:5
    |
11  |     stackblur_iter::blur_argb(&mut img_ref, 1);
    |     ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
help: <3782> was later invalidated at offsets [0x8..0xc]
   --> src/main.rs:11:5
    |
11  |     stackblur_iter::blur_argb(&mut img_ref, 1);
    |     ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
    = note: backtrace:
    = note: inside `imgref_iter::iter::windows::ptr::IterWindowsPtr::<u32>::window` at /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/imgref-iter-0.3.3/src/iter/windows/ptr.rs:135:13
    = note: inside closure at /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/imgref-iter-0.3.3/src/iter/windows/ptr.rs:145:51
    = note: inside `std::option::Option::<usize>::map::<imgref_iter::iter::generic::ptr::IterPtr<u32>, [closure@<imgref_iter::iter::windows::ptr::IterWindowsPtr<u32> as std::iter::Iterator>::next::{closure#0}]>` at /Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/lib/rustlib/src/rust/library/core/src/option.rs:929:29
    = note: inside `<imgref_iter::iter::windows::ptr::IterWindowsPtr<u32> as std::iter::Iterator>::next` at /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/imgref-iter-0.3.3/src/iter/windows/ptr.rs:145:3
    = note: inside `<imgref_iter::iter::windows::IterWindows<u32> as std::iter::Iterator>::next` at /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/imgref-iter-0.3.3/src/iter/windows/mod.rs:43:3
    = note: inside `<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::adapters::zip::ZipImpl<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>>::next` at /Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/lib/rustlib/src/rust/library/core/src/iter/adapters/zip.rs:155:21
    = note: inside `<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::next` at /Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/lib/rustlib/src/rust/library/core/src/iter/adapters/zip.rs:84:9
    = note: inside `<for<'r> fn(&'r mut std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>) -> std::option::Option<<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::Item> {<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::next} as std::ops::FnOnce<(&mut std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>,)>>::call_once - shim(for<'r> fn(&'r mut std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>) -> std::option::Option<<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::Item> {<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::next})` at /Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/lib/rustlib/src/rust/library/core/src/ops/function.rs:248:5
    = note: inside `std::iter::adapters::chain::and_then_or_clear::<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>, (imgref_iter::iter::generic::ptr::IterPtrMut<u32>, imgref_iter::iter::generic::Iter<u32>), for<'r> fn(&'r mut std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>) -> std::option::Option<<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::Item> {<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>> as std::iter::Iterator>::next}>` at /Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/lib/rustlib/src/rust/library/core/src/iter/adapters/chain.rs:287:13
    = note: inside `<std::iter::Chain<std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>, std::iter::Zip<imgref_iter::iter::windows::ptr::IterWindowsPtrMut<u32>, imgref_iter::iter::windows::IterWindows<u32>>> as std::iter::Iterator>::next` at /Users/marko/.rustup/toolchains/nightly-x86_64-apple-darwin/lib/rustlib/src/rust/library/core/src/iter/adapters/chain.rs:50:9
    = note: inside `stackblur_iter::blur::<u32, stackblur_iter::color::Argb<stackblur_iter::color::serial::StackBlurrableU32>, [closure@stackblur_iter::blur_argb::{closure#0}], fn(stackblur_iter::color::Argb<stackblur_iter::color::serial::StackBlurrableU32>) -> u32 {stackblur_iter::color::Argb::<stackblur_iter::color::serial::StackBlurrableU32>::to_u32}>` at /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/stackblur-iter-0.1.5/src/lib.rs:88:23
    = note: inside `stackblur_iter::blur_argb` at /Users/marko/.cargo/registry/src/github.com-1ecc6299db9ec823/stackblur-iter-0.1.5/src/lib.rs:233:2
note: inside `main` at src/main.rs:11:5
   --> src/main.rs:11:5
    |
11  |     stackblur_iter::blur_argb(&mut img_ref, 1);
    |     ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

note: some details are omitted, run with `MIRIFLAGS=-Zmiri-backtrace=full` for a verbose backtrace

error: aborting due to previous error

Combines #4 and #5

Closes #3

Using both at the same time now actually outperforms rayon on its own.

This is surprising.

Currently a lane count of 8 is the fastest on my machine, but this is tunable.

This is all thanks to imgref-iter's experimental SIMD iterators.

It's great to see a correct stackblur implementation in Rust! I switched out the algorithm I use in i3lockr with yours and it takes about ~100 ms whereas the multithreaded C version of stackblur I was using takes ~30 ms. However, the C version is naive and clobbers the gamma. The 100 ms is including sRGB to linear round trip. I expect with multithreading yours may be very close to 30 ms, if not faster.

~~This makes the library slower in single-threaded. I doubt it will be any different in multi-threaded. I implemented it on top of the rayon branch but it's not going to actually support rayon until SIMD is a speedup, which it is apparently not.~~

~~Would appreciate help getting it to be faster.~~

New performance analysis TBD

I made changes in imgref-iter to make this possible.

Behold:

stackblur_1024       ... bench:  18,514,600 ns/iter (+/- 3,095,887)
stackblur_128        ... bench:  11,743,780 ns/iter (+/- 2,533,038)
stackblur_16         ... bench:  11,109,190 ns/iter (+/- 1,912,941)
blur_srgb_1024       ... bench:  17,162,590 ns/iter (+/- 3,337,943)
blur_srgb_128        ... bench:  14,260,590 ns/iter (+/- 2,214,481)
blur_srgb_16         ... bench:  13,107,140 ns/iter (+/- 2,408,900)
par_blur_srgb_1024   ... bench:   4,309,255 ns/iter (+/- 352,918)
par_blur_srgb_128    ... bench:   3,979,880 ns/iter (+/- 305,431)
par_blur_srgb_16     ... bench:   4,063,585 ns/iter (+/- 310,803)

Fixes #3

So there's this thing called libdivide that I don't fully understand, but basically it makes it possible to perform integer divides without actually using division instructions, which is a lot faster.

The fastdivide crate has a really good explanation:

Division is a very costly operation for your CPU. You may have noticed that when the divisor is known at compile time, your compiler transforms the operations into a cryptic combination of a multiplication and bitshift.

The key idea is that, rather than computing

N / D

It is faster to compute (with k sufficiently large)

N * ( 2^k / D ) / (2^k)

If D is known in advance, (2^k / D) can be precomputed by the compiler.

Unfortunately if the divisor is unknown at compile time, the compiler cannot use this trick.

The point of fastdivide is to apply the same trick by letting you precompute a DivideU64 object.

Dividing takes up a lot of the CPU time spent blurring. If divisions can be eliminated, the blur could become a lot faster.