-
Notifications
You must be signed in to change notification settings - Fork 1
/
gauss.py
62 lines (51 loc) · 1.96 KB
/
gauss.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
# -*- coding: utf-8 -*-
"""
Created on Tue Aug 28 01:05:13 2018
@author: zhang
"""
#!/usr/bin/env python
"""Module providing functionality surrounding gaussian function.
"""
SVN_REVISION = '$LastChangedRevision: 16541 $'
import sys
import numpy
def gaussian2(size, sigma):
"""Returns a normalized circularly symmetric 2D gauss kernel array
f(x,y) = A.e^{-(x^2/2*sigma^2 + y^2/2*sigma^2)} where
A = 1/(2*pi*sigma^2)
as define by Wolfram Mathworld
http://mathworld.wolfram.com/GaussianFunction.html
"""
A = 1/(2.0*numpy.pi*sigma**2)
x, y = numpy.mgrid[-size//2 + 1:size//2 + 1, -size//2 + 1:size//2 + 1]
g = A*numpy.exp(-((x**2/(2.0*sigma**2))+(y**2/(2.0*sigma**2))))
return g
def fspecial_gauss(size, sigma):
"""Function to mimic the 'fspecial' gaussian MATLAB function
"""
x, y = numpy.mgrid[-size//2 + 1:size//2 + 1, -size//2 + 1:size//2 + 1]
g = numpy.exp(-((x**2 + y**2)/(2.0*sigma**2)))
return g/g.sum()
def main():
"""Show simple use cases for functionality provided by this module."""
from mpl_toolkits.mplot3d.axes3d import Axes3D
import pylab
argv = sys.argv
if len(argv) != 3:
print >>sys.stderr, 'usage: python -m pim.sp.gauss size sigma'
sys.exit(2)
size = int(argv[1])
sigma = float(argv[2])
x, y = numpy.mgrid[-size//2 + 1:size//2 + 1, -size//2 + 1:size//2 + 1]
fig = pylab.figure()
fig.suptitle('Some 2-D Gauss Functions')
ax = fig.add_subplot(2, 1, 1, projection='3d')
ax.plot_surface(x, y, fspecial_gauss(size, sigma), rstride=1, cstride=1,
linewidth=0, antialiased=False, cmap=pylab.jet())
ax = fig.add_subplot(2, 1, 2, projection='3d')
ax.plot_surface(x, y, gaussian2(size, sigma), rstride=1, cstride=1,
linewidth=0, antialiased=False, cmap=pylab.jet())
pylab.show()
return 0
if __name__ == '__main__':
sys.exit(main())