20

I'm throwing this out there in hope that someone will have attempted something this ridiculous before. My goal is to take in an input image, and segment it based upon the standard deviation of a small window around each pixel. Bascially, this should mathematically resemble a gauss or box filter, in that it will be applied to a compile time (or even run-time) user specified window size around each pixel, and the destination array will contain the SD information at each pixel, in an image the same size as the original.

The idea is to do this on an image in HSV space, so that I can easily find regions of homogeneous color (i.e. those with small local SDs in the Hue and Sat planes) and extract them from the image for more in-depth processing.

So the question is, has anyone ever built a custom filter like this before? I don't know how to do the SD in a simple box type filter kernel like the ones used for gauss and blur, so I'm guessing I'll have to use the FilterEngine construct. Also, I forgot to mention I'm doing this in C++.

Your advice and musings are much appreciated.

gankoji
  • 853
  • 1
  • 10
  • 23

2 Answers2

50

Wikipedia has a nice explanation of standard deviation, which you can use to for a standard deviation filter.

Basically, it boils down to blurring the image with a box filter, blurring the square of the image with a box filter, and taking the square root of their difference.

UPDATE: This is probably better shown with the equation from Wikipedia... enter image description here

You can think of the OpenCV blur function as representing the expected value (i.e., E[X] a.k.a. the sample mean) of the neighborhood of interest. The random samples X in this case are represented by image pixels in the local neighborhood. Therefore, by using the above equivalence we have something like sqrt(blur(img^2) - blur(img)^2) in OpenCV. Doing it this way allows you to compute the local means and standard deviations.

Also, just in case you are curious about the mathematical proof. This equivalence is known as the computational formula for variance.

Here is how you can do this in OpenCV:

#include <iostream>
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>

using namespace std;
using namespace cv;

Mat mat2gray(const Mat& src)
{
    Mat dst;
    normalize(src, dst, 0.0, 1.0, NORM_MINMAX);
    return dst;
}

int main()
{
    Mat image = imread("coke-can.jpg", 0);

    Mat image32f;
    image.convertTo(image32f, CV_32F);

    Mat mu;
    blur(image32f, mu, Size(3, 3));

    Mat mu2;
    blur(image32f.mul(image32f), mu2, Size(3, 3));

    Mat sigma;
    cv::sqrt(mu2 - mu.mul(mu), sigma);

    imshow("coke", mat2gray(image32f));
    imshow("mu", mat2gray(mu));
    imshow("sigma",mat2gray(sigma));
    waitKey();
    return 0;
}

This produces the following images:

Original

enter image description here

Mean

enter image description here

Standard Deviation

enter image description here

Hope that helps!

mevatron
  • 13,911
  • 4
  • 55
  • 72
  • 1
    This is phenomenal, thank you! I was having a tough time thinking through the process of summing the squared differences as a single filter kernel, but I see that separating them could work well. One question though, don't you have to take the difference of each element vs. the local mean before you square that difference? The SD formula on the page you mention is (sum((x-mu)^2)^(1/2). I think what you wrote is (sum((x^2)-(mu^2)))^(1/2). Arguably, however, it seems your code has produced rather accurate results. – gankoji Jul 14 '12 at 06:50
  • @gankoji Awesome! Glad you found this sample helpful. Basically, I am using the mathematical equivalence that is shown in the Wiki article. I also updated the answer to hopefully explain this in a more complete fashion. – mevatron Jul 15 '12 at 04:10
  • Ahh, I see the equivalence you mentioned. Sorry, I'm not a statistics buff, I had no idea that was the case! I will try to implement this code on Monday when I return to the office, and let you know how it goes. Thanks again! – gankoji Jul 15 '12 at 16:17
  • @gankoji I updated the answer with a link to the derivation of the equivalence in case you're curious about the math :) – mevatron Jul 16 '12 at 01:04
  • @mevatron hello very nice tutorial....i have to average difference between two edge detected images...very new to opencv i have done in matlab pls help in opencv...with formula – user960439 Mar 25 '13 at 07:12
  • @mevatron : +1 for Theory + Code + Example = "Wonderful!!" :) – G453 Apr 20 '15 at 11:04
  • Sumptuous explanation with example – Jeru Luke Jun 01 '22 at 15:57
1

In case you want to use this in more general way this can produce nan values

Values close to zero can be sometimes "negative".

Mat sigma;
cv::sqrt(mu2 - mu.mul(mu), sigma);

correct way should be

Mat sigma;
cv::sqrt(cv::abs(mu2 - mu.mul(mu)), sigma);