Troubles with median

[mileslucas]

I was checking this package out and seeing if it would be able to replace some simple median-filtering code I've got in a package. Previously, I was using ImageFiltering.jl's mapwindow(median!, ...) to accomplish this. My solution using StaticKernels is

using StaticKernels
k = Kernel{(-5:5,-5:5)}(median!)
map(k, extend(img, StaticKernels.ExtensionReplicate()))

This fails due to the lack of Base.iterate for the window view passed to median! (similarly fails for median. If I try Kernel{...}(w -> median(Tuple(w))) the runtime explodes (like, I have to sigint the REPL).

Am I missing something here?

[stev47]

A couple of notes here:

• An 11x11 kernel is something I would consider "big" and StaticKernels might not be the right choice for this, since we just loop over all windows in a dumb way. There are better algorithms for your specific usecase which are less dependent of the window size (see e.g. the Herk-Gil-Werman algorithm which scales logarithmically for median)
• Making the window view iterable is definitely a near-term goal from which I've held off up until now since it always introduced allocations before Julia v1.5.
• Kernel{...}(w -> median(Tuple(w))) works, but for your kernel size you are running into the julia exponential compile-time problem for nested for-loops (there is a warning in the readme). This should eventually be fixed within Julia, but I'm inclined to just make the boundary windows non-statically sized as a workaround (you then get slower runtime though).
• In contrast to Images.jl we are not copying the image window into a buffer, so mutating the window will modify the underying image and potentially mess up your results (since windows overlap).
  So you ideally want a non-allocating non-mutating median algorithm, which doesn't seem to be implemented in Statistics currently.

[mileslucas]

Wow, what a thorough, well-put, and intelligent response. Thanks, this answers a lot of questions, even ones I didn't know I had! I suppose I'll be keeping an eye out on this package and let it mature and seep into the various ecosystems.

[RoyiAvital]

Is there a non allocating non mutating median without pre allocated buffer?
Maybe an idea is that map() will allow sending kwargs into the kernel function?
So in that case a buffer can be used.

[stev47]

There is one in this old pull request.

I would want to avoid changing map syntax too much, but you can definitely use a buffer in your window function already by creating a closure.

[RoyiAvital]

What's blocking the pull request?

I am not sure I get what you mean, maybe an example?
Will closure be optimal performance wise?

[stev47]

It says "Review required".

I'm saying that you can reference a buffer from within your window function like this:

const buf = zeros(3)
k = Kernel{(-1:1,)}(@inline w -> (copyto!(buf, Tuple(w)); median!(buf))

It should be as fast as you get using a heap-allocated buffer, but please benchmark for yourself.

[RoyiAvital]

We now have JuliaArrays/StaticArrays.jl#973 merged.