Median filtering in constant time ?

[denis]

Keith, folks,

fwiw, Wikipedia Median filter -> http://nomis80.org/ctmf.html ->
a nice paper by Perreault and Hebert, "Median filtering in constant
time" (2007, 6p).
His code for 2d image medians (the 80 % use case ?) with a cython
wrapper
is not "constant time", but increases pretty slowly --

# from: test-png.py
# run: 3 May 2011 15:02 in ~/py/etc/cython/median Denis.local mac
10.4.11 p
read (288, 450, 3) float32 from images/Schnee-Kyoto.png av 0.653
max 1
877 msec ctmedianfilter2 (288, 450, 3) r 1
929 msec ctmedianfilter2 (288, 450, 3) r 8
984 msec ctmedianfilter2 (288, 450, 3) r 16
1048 msec ctmedianfilter2 (288, 450, 3) r 32

I'd put code + wrapper up on https://github.com/kwgoodman/roly
but I'm a complete noob at git ...

Also, Real Code (TM) needs reference test cases; any ideas ?

cheers
-- denis

[Keith Goodman]

On Tue, May 3, 2011 at 7:30 AM, denis <denis...@t-online.de> wrote:
> Keith, folks,
>
>  fwiw, Wikipedia Median filter -> http://nomis80.org/ctmf.html ->
> a nice paper by Perreault and Hebert, "Median filtering in constant
> time" (2007, 6p).
> His code for 2d image medians (the 80 % use case ?) with a cython
> wrapper
> is not "constant time", but increases pretty slowly --

More algorithms! That's exciting.

The three moving medians in roly all work along 1d lines. In fact for
testing they only take 1d input (but can be made to take nd input).
But if you are interested in a nd median we can include that. Roly is
just a testing ground. One which everyone can borrow from.

So far the code in roly is all BSD. The code you link to is GPL. I see
no problem with that. Each file in roly has its own license. I imagine
people will end up taking bits and pieces of roly for their own
projects.

>  # from: test-png.py
>  # run: 3 May 2011 15:02  in ~/py/etc/cython/median  Denis.local  mac
> 10.4.11 p
>  read (288, 450, 3) float32 from images/Schnee-Kyoto.png  av 0.653
> max 1
>  877 msec ctmedianfilter2 (288, 450, 3)  r 1
>  929 msec ctmedianfilter2 (288, 450, 3)  r 8
>  984 msec ctmedianfilter2 (288, 450, 3)  r 16
>  1048 msec ctmedianfilter2 (288, 450, 3)  r 32
>
> I'd put code + wrapper up on https://github.com/kwgoodman/roly
> but I'm a complete noob at git ...

I thought git was difficult to learn coming from svn. But well worth
it. The best documentation I found was http://progit.org/book

> Also, Real Code (TM) needs reference test cases; any ideas ?

I wrote a slow, python "for loop" to calculate the moving median along
1d. You can do that same for the 2d filter.

The 3 moving medians in roly:

roly.linkedlist.move_median
roly.doubleheap.move_median
roly.slow.move_median

Maybe you can name yours (?)

roly.???.move_median2d
roly.slow.move_median2d

[denis]

On May 3, 4:30 pm, denis <denis-bz...@t-online.de> wrote:
> Keith, folks,
>
>  fwiw, Wikipedia Median filter ->http://nomis80.org/ctmf.html->
> a nice paper by Perreault and Hebert, "Median filtering in constant
> time" (2007, 6p).

Note that median2d algorithms operate on 8-bit rgba pixels,
so can use the fantastic trick
"computing the median from a histogram is done in constant time
by summing the values from one end
and stopping when the sum reaches (2r+1)^2 / 2."
Median1d s on the other hand run on floats, doing compares.
Cf. Perreault section IV, O(1) sorts.

Would there be a market for single-pass estimation of quantiles / pdf
aka equi-depth histogram
along the llines of NR3 (boo hiss from some people) p. 435 ?
A cython wrapper would be trivial, but their license may not allow
that, dunno.

cheers
-- denis

[Keith Goodman]

What is a quantile? Does that convert an array of floats to ints where
the smallest n values are assigned 0, next n smallest 1, etc?

[denis]

On May 4, 3:33 pm, Keith Goodman <kwgood...@gmail.com> wrote:

> What is a quantile? Does that convert an array of floats to ints where
> the smallest n values are assigned 0, next n smallest 1, etc?

That's about right. See e.g. Wikipedia Quantile, Digital_audio,
Quantization (image processing).
(Almost all audio is of of course quantized on recording, as on CDs /
MP3,
and many images -- digital cameras, medical images --
are quantized to 8 bits each of R G B and A.
So image filters (scipy.ndimage) do 8bit rgba -> floats -> filter ->
8bits,
usually just uniform steps floor( float * 255 )
but not always.

What's in it for us ? Well, one can quantize floats -> bytes with
arbitrary steps pretty fast
(~ 150 clock cycles in C on my old mac ppc):
a solution looking for problems.

See also Vector_quantization, quantization in 2d 3d 100d, in Wikipedia
and scikits.learn:
tough, interesting, zillions of papers, steps / regions non-uniform.

Don't know if there's a real market for faster scalar quant or
Vector_quant
or just sounds-nice, don't-call-us-we'll-call-you.

cheers
-- denis