smooth.py
Go to the documentation of this file.
00001 """
00002 cookb_signalsmooth.py
00003 
00004 from: http://scipy.org/Cookbook/SignalSmooth
00005 """
00006 
00007 import numpy as np
00008 import matplotlib.pyplot as plt
00009 
00010 def smooth(x, window_len=10, window='hanning'):
00011     """smooth the data using a window with requested size.
00012     
00013     This method is based on the convolution of a scaled window with the signal.
00014     The signal is prepared by introducing reflected copies of the signal 
00015     (with the window size) in both ends so that transient parts are minimized
00016     in the begining and end part of the output signal.
00017     
00018     input:
00019         x: the input signal 
00020         window_len: the dimension of the smoothing window
00021         window: the type of window from 'flat', 'hanning', 'hamming', 'bartlett', 'blackman'
00022             flat window will produce a moving average smoothing.
00023 
00024     output:
00025         the smoothed signal
00026         
00027     example:
00028 
00029     import numpy as np    
00030     t = np.linspace(-2,2,0.1)
00031     x = np.sin(t)+np.random.randn(len(t))*0.1
00032     y = smooth(x)
00033     
00034     see also: 
00035     
00036     numpy.hanning, numpy.hamming, numpy.bartlett, numpy.blackman, numpy.convolve
00037     scipy.signal.lfilter
00038  
00039     TODO: the window parameter could be the window itself if an array instead of a string   
00040     """
00041 
00042     if x.ndim != 1:
00043         raise ValueError, "smooth only accepts 1 dimension arrays."
00044 
00045     if x.size < window_len:
00046         raise ValueError, "Input vector needs to be bigger than window size."
00047 
00048     if window_len < 3:
00049         return x
00050 
00051     if not window in ['flat', 'hanning', 'hamming', 'bartlett', 'blackman']:
00052         raise ValueError, "Window is on of 'flat', 'hanning', 'hamming', 'bartlett', 'blackman'"
00053 
00054     # Hai and Advait: removed things we don't understand. (Feb 26, 2010)
00055     # this code is mirroring the first and last window_len elements.
00056     #s=np.r_[2*x[0]-x[window_len:1:-1], x, 2*x[-1]-x[-1:-window_len:-1]]
00057     # and added in what we do understand
00058     s = x
00059 
00060     if window == 'flat': #moving average
00061         w = np.ones(window_len,'d')
00062     else:
00063         w = getattr(np, window)(window_len)
00064     y = np.convolve(w/w.sum(), s, mode='same')
00065     return y[window_len-1:-window_len+1]
00066 
00067 
00068 #*********** part2: 2d
00069 
00070 from scipy import signal
00071 
00072 def gauss_kern(size, sizey=None):
00073     """ Returns a normalized 2D gauss kernel array for convolutions """
00074     size = int(size)
00075     if not sizey:
00076         sizey = size
00077     else:
00078         sizey = int(sizey)
00079     x, y = np.mgrid[-size:size+1, -sizey:sizey+1]
00080     g = np.exp(-(x**2/float(size) + y**2/float(sizey)))
00081     return g / g.sum()
00082 
00083 def blur_image(im, n, ny=None) :
00084     """ blurs the image by convolving with a gaussian kernel of typical
00085         size n. The optional keyword argument ny allows for a different
00086         size in the y direction.
00087     """
00088     g = gauss_kern(n, sizey=ny)
00089     improc = signal.convolve(im, g, mode='valid')
00090     return(improc)
00091 
00092 
00093 def smooth_demo():
00094     t = np.linspace(-4,4,100)
00095     x = np.sin(t)
00096     xn = x + np.random.randn(len(t)) * 0.1
00097     y = smooth(x)
00098     ws = 31
00099 
00100     plt.subplot(211)
00101     plt.plot(np.ones(ws))
00102 
00103     windows=['flat', 'hanning', 'hamming', 'bartlett', 'blackman']
00104 
00105     plt.hold(True)
00106     for w in windows[1:]:
00107         #eval('plt.plot('+w+'(ws) )')
00108         plt.plot(getattr(np, w)(ws))
00109 
00110     plt.axis([0,30,0,1.1])
00111 
00112     plt.legend(windows)
00113     plt.title("The smoothing windows")
00114     plt.subplot(212)
00115     plt.plot(x)
00116     plt.plot(xn)
00117     for w in windows:
00118         plt.plot(smooth(xn, 10, w))
00119     l = ['original signal', 'signal with noise']
00120     l.extend(windows)
00121     plt.legend(l)
00122     plt.title("Smoothing a noisy signal")
00123     #plt.show()
00124 
00125 
00126 if __name__=='__main__':
00127     
00128     # part 1: 1d
00129     smooth_demo()
00130     
00131     # part 2: 2d
00132     X, Y = np.mgrid[-70:70, -70:70]
00133     Z = np.cos((X**2+Y**2)/200.)+ np.random.normal(size=X.shape)
00134     Z2 = blur_image(Z, 3)
00135     plt.figure()
00136     plt.imshow(Z)
00137     plt.figure()
00138     plt.imshow(Z2)
00139     plt.show()
00140     


2010_biorob_everyday_mechanics
Author(s): Advait Jain, Hai Nguyen, Charles C. Kemp (Healthcare Robotics Lab, Georgia Tech)
autogenerated on Wed Nov 27 2013 11:58:43