Suppose that $G(i,j)$ is a Gaussian decay function on the distance between points $i$ and $j$ of an image. In addition, $D(i,j)$ is the difference between the VALUES of the image at those points. Now, at every point $i$ of an image I need to compute the summation of the neighboring $G(i,j) D(i,j) $ , this is, the product of both functions. This is a convolution-like operation: for every point I add the contributions of a function at every neighboring point. However, I do not know if that can be expressed as a convolution because G and D are functions on different variables: the distance (a raw gaussian filter) and the image intensity difference. I need help at expressing the summation as a combination of convolution operations in order to easily code it in a program without the need for loops. In addition, it would allow me to possibly do it in the frequency domain. Any help will be appreciated!
-
$\begingroup$ Welcome to DSP.SE! Do you have links to, or more specific descriptions of, what $G$ and $D$ actually are mathematically? $\endgroup$– PhononCommented Jan 2, 2014 at 6:12
-
$\begingroup$ I think that you are trying to implement bilateral filtering. Am I correct? For speedups, try the integral histogram approach. $\endgroup$– visoftCommented Jan 2, 2014 at 8:57
-
$\begingroup$ Hello guys, thank you. It is not bilateral filtering. What I need is to apply a very complex filter (I could simply code it with a couple of loops but I want to simplify it). One of the difficulties I find is that I bump into a Gaussian decay G, i.e. it convolves penalizing distance from the center of the kernel as usual, but it is multiplied by a gaussian decay on the intensity difference between the points, that's D. It is not a convolution I think because I'm comparing the intensities at both sites and subtracting. $\endgroup$– user116773Commented Jan 2, 2014 at 18:05
-
$\begingroup$ The summation is outside the product. I mean, for every point I take G times H of every neighbor and then add them. $\endgroup$– user116773Commented Jan 2, 2014 at 18:06
-
$\begingroup$ Well, maybe what I am asking is equivalent to a bilateral filter, right? $\endgroup$– user116773Commented Jan 2, 2014 at 18:09
2 Answers
It is equivalent to a bilateral filter. You can move the kernel to each local part of image to do the dot product, then sum them up to get the result in a loop. There is also a fast implementation with O(1) in matlab programming for your reference. And the relevant publications are:
K.N. Chaudhury, D. Sage, and M. Unser, "Fast O(1) bilateral filtering using trigonometric range kernels," IEEE Transactions on Image Processing, vol. 20, no. 11, 2011.
K.N. Chaudhury, "Acceleration of the shiftable O(1) algorithm for bilateral filtering and non-local means," arXiv:1203.5128v1.
I assume the Gaussian kernel is position shift-invariant. So normalization means the sum of elements inside the kernel equals to 1. For a 3*3 Gaussian kernel, it is
1/16 2/16 1/16
2/16 4/16 2/16
1/16 2/16 1/16
The bilateral filter can also be normalized with the K factor of integral of the product of Gaussian kernel and intensity kernel. You can try both to observe the effect of normalization on your image data.
-
$\begingroup$ Thanks, that's great. And what if I do not need to normalize the kernels? What I need is the "weighted" convolution as explained above but I do not need to normalize the kernel (which in fact is not a constant normalization factor in bilateral filtering). $\endgroup$ Commented Jan 2, 2014 at 18:36
-
$\begingroup$ Please see my update. Would that work for you? $\endgroup$ Commented Jan 2, 2014 at 18:43
-
$\begingroup$ That is okay for the Gaussian kernel. But what about the other term, the one that is a Gaussian decay on the intensity difference? The bilateral filter, as far as I know (maybe I'm not correct) normalizes the whole kernel which is a product of two terms. I need, instead, to compute the added contributions of neighbors without normalization (even if the value is biased). Maybe I can simply avoid a normalization factor, but I do not know if that fast method allows to do that in its formulation. $\endgroup$ Commented Jan 2, 2014 at 19:00
-
$\begingroup$ Yes bilateral filter can also be scaled with both kernels. See my update for the link. I would suggest you try both and see which one works better. Thanks $\endgroup$ Commented Jan 3, 2014 at 13:51
Let * denote convolution operation, $g$ denote Gaussian kernel, &h& some other kernel and I be your image.
Convolution is a linear operation: $g*h*I = g*(h*I)$
Gradient operation is also a linear one.
$\nabla (g*I) = (\nabla g) * I$
For this reason if you first differentiate your kernel then apply it you will end up with smooth gradients. Actually, this is what Sobel has achieved:
http://en.wikipedia.org/wiki/Sobel_operator
If you further smooth the image prior to applying Sobel, this actually would have the same effect of increasing kernel size.
For further info on difference and Gaussian tricks, I would also suggest you to check http://en.wikipedia.org/wiki/Difference_of_Gaussians
