Posted by
Gabriel Landini on
Sep 07, 2013; 5:04pm
URL: http://imagej.273.s1.nabble.com/median-filter-that-preserves-peaks-tp5004693p5004716.html
On Saturday 07 Sep 2013 06:19:18 Kenneth Sloan wrote:
> The nice thing about filters like the MedianFilter is that once you have
> committed to that style, you can implement ANY function of the values
> inside the window. This is what separates Median from more
> signal-processing-inspired methods like Mean. Median is a frank "spatial
> domain" operation, with no clear partner in the Fourier domain. Given a
> 5x5 window, say, you have 25 values - and are free to compute "anything you
> want to" from those 25 values.
>
> To completely solve this problem, I propose a small suite of programs, all
> modified in roughly the same way. The idea is that small, local
> fluctuations are "noise" but macroscopic major features are "signal".
>
> The basic idea is: if the central value in a window is the MIN or MAX within
> the window, LEAVE IT ALONE. Otherwise - do the usual thing.
>
> So, XXXFilter becomes:
>
> for (x = 0, width)
> for (y = 0, height)
> {
> gatherWindowAndCompute(x,y,windowW, windowH)
> if (MinOrMax(input(x,y)) output(x,y) = input(x,y)
> else output(x,y) = XXX(window);
> }
>
> I hope that communicates the idea.
>
> Now…to fully implement what the OP wants, I suggest implementing
> MedianFilter and also MeanFilter.
Very interesting. I guess that this would be similar to:
1. compute the local minima and maxima for a kernel size and their grey
values,
2. filter the original image with a median filter
3. replace, in the filtered image, the local minima and maxima of the
original.
This macro implements that:
//------------------8<---------------
// MedianLocalMinMax.txt
// Based on the original idea of Kenneth Sloan
// posted to the IJ mailing list on Sept 7, 2013.
// G. Landini at bham. ac. uk
a=getNumber("Local radius",3);
setBatchMode(true);
b=getTitle();
setOption("BlackBackground", true);
run("Duplicate...", "title=Median");
run("Duplicate...", "title=Maxima");
run("Maximum...", "radius="+a);
run("Image Calculator...", "image1=Maxima operation=Subtract image2=Median
create");
setThreshold(0, 0);
run("Convert to Mask");
selectWindow("Maxima");
run("Close");
selectWindow("Result of Maxima");
rename("LocalMaxima");
selectWindow("Median");
run("Duplicate...", "title=Minima");
run("Minimum...", "radius="+a);
run("Image Calculator...", "image1=Median operation=Subtract image2=Minima
create");
setThreshold(0, 0);
run("Convert to Mask");
selectWindow("Minima");
run("Close");
selectWindow("Result of Median");
rename("LocalMinima");
imageCalculator("Add" , "LocalMaxima","LocalMinima");
selectWindow("LocalMinima");
close();
selectWindow("Median");
run("Median...", "radius="+a);
imageCalculator("Subtract", "Median","LocalMaxima");
imageCalculator("AND", "LocalMaxima", b);
imageCalculator("Add", "Median","LocalMaxima");
selectWindow("LocalMaxima");
close();
selectWindow("Median");
rename("MedianLocalMinMax");
setBatchMode(false);
//------------------8<---------------
Furthermore, it is also possible to replace with the *regional* minima and
maxima instead, as those are regions guaranteed to be surrounded by brighter
and darker pixels respectively and more likely to be peaks and wells in the
image.
Interestingly the regional minima and maxima do not depend on kernel sizes.
The following macro does that (the Domes plugin is part of the morphology
collection, available in my page)
//------------------8<---------------
// MedianRegionalMinMax.txt
// Based on the original idea of Kenneth Sloan
// posted to the IJ mailing list on Sept 7, 2013.
// G. Landini at bham. ac. uk
a=getNumber("Local radius",3);
setBatchMode(true);
b=getTitle();
setOption("BlackBackground", true);
run("Duplicate...", "title=Median");
run("Duplicate...", "title=RMaxima");
run("Domes ", "height=1");
setThreshold(1, 255);
run("Convert to Mask");
selectWindow("Median");
run("Duplicate...", "title=RMinima");
run("Domes ", "height=1 basins");
setThreshold(1, 255);
run("Convert to Mask");
imageCalculator("Add" , "RMaxima","RMinima");
selectWindow("RMinima");
close();
selectWindow("Median");
run("Median...", "radius="+a);
imageCalculator("Subtract", "Median","RMaxima");
imageCalculator("AND", "RMaxima", b);
imageCalculator("Add", "Median","RMaxima");
selectWindow("RMaxima");
close();
selectWindow("Median");
rename("MedianRegionalMinMax");
setBatchMode(false);
//------------------8<---------------
Hope it is useful. Let me know of any bugs.
Cheers
Gabriel
--
ImageJ mailing list:
http://imagej.nih.gov/ij/list.html