Simpliest way to generate a 1D gaussian kernel
Asked Answered
T

2

11

I'm wondering what would be the easiest way to generate a 1D gaussian kernel in python given the filter length. I think that the idea is to evaluate the normal distribution for the values of the vector [-filter-length,...,filter_length], is it correct?

So far, I've done this, but I don't know why it is not correct:

result = np.zeros( filter_length )

mid = filter_length/2
result=[(1/(sigma*np.sqrt(2*np.pi)))*(1/(numpy.exp((i**2)/(2*sigma**2)))) for i in range(-mid,mid+1)]  

return result

where sigma is the standard deviation, which is a parameter. filter-length is also a parameter.

It's incorrect because I get, for example, for length=3 and sigma=math.sqrt(1.0/2/math.log(2))

[0.23485931967491286, 0.46971863934982572, 0.23485931967491286]

And it should be:

[0.25, 0.5, 0.25]

So, is there any problem of rounding? I don't know what is going on...

Edit I think that I should truncate somehow

Problem Solved The problem was that I wasn't normalizing. I had to divide the vector by the sum of all its components.

Trinitarianism answered 16/2, 2013 at 23:0 Comment(5)
what value is sigma for the results you give?Bystreet
for sigma=1 i think the values should be [0.24,0.40,0.24]. i am unsure why you expect [0.25,0.5,0.25].Bystreet
no, it's for sigma=math.sqrt(1.0/2/math.log(2))Trinitarianism
I have added what I think is the correct sigma to my answer, at least with the Gaussian formula you use. Why did you choose sigma=math.sqrt(1.0/2/math.log(2))?Soler
[0.25, 0.5, 0.25] is a triangular kernel, not a Gaussian. You cannot make a Gaussian kernel this short, it will never share the good properties of the Gaussian kernel. You need more samples. See here for details: crisluengo.net/archives/695Monetary
S
5

I am not very firm with numpy syntax, but if you convolve a kernel with a dirac impulse, you get the same kernel as output.

So you could simply use the inbuild scipy.ndimage.filters.gaussian_filter1d function, and use this array as input: [ 0, 0, 0, ... 0, 1, 0, ...0, 0, 0]

The output should be a gaussian kernel, with a value of 1 at its peak. (replace 1 with the maximum you want in your desired kernel)

So in essence, you will get the Gaussian kernel that gaussian_filter1d function uses internally as the output. This should be the simplest and least error-prone way to generate a Gaussian kernel, and you can use the same approach to generate a 2d kernel, with the respective scipy 2d function. Of course if the goal is to do it from scratch, then this approach is only good as a reference

In regards to your equation:
to get [..., 0.5, ...] as the output with your formula, you need to solve
(1/(sigma*np.sqrt(2*np.pi)) = 0.5
so the correct sigma should be
sigma = math.sqrt(2*1/np.pi)

Soler answered 16/2, 2013 at 23:9 Comment(3)
reading those docs, you'd want mode='constant' and also it's not clear to me whether the amplitude will be 1, or whether the filter has an integrated value of 1. so you might need to divide by max(...).Bystreet
what do you mean that you get the same kernel?Trinitarianism
Yes, you get the same kernel as output that the gaussian_filter1d function uses internally. I am pretty sure that this is the simplest way to generate a 1D Gaussian kernel. Of course, if you want to generate the kernel from scratch as an exercise, you will need a different approach. But you could at least use this method as a reference to compare to your outputSoler
R
0

For those that like to have a fully working copy/paste code example:

import numpy as np

filter_length = 3
sigma=math.sqrt(1.0/2/math.log(2))
result = np.zeros( filter_length )

mid = int(filter_length/2)
result=[(1/(sigma*np.sqrt(2*np.pi)))*(1/(np.exp((i**2)/(2*sigma**2)))) for i in range(-mid,mid+1)]  
sumresult = np.sum(result)
print(result/sumresult)

[0.25 0.5 0.25]

Raynaraynah answered 19/3 at 8:44 Comment(1)
If you’re going to normalize, you can leave the 1/(sigma*np.sqrt(2*np.pi)) computation out. It’s no longer doing anything.Monetary

© 2022 - 2024 — McMap. All rights reserved.