In [1]:

    

import numpy as np, cmath,scipy as sp
import scipy.io
from matplotlib import pyplot as plt

from numpy import pi, sin, cos, exp, sqrt, log, random  #import basic functions from numpy that we'll need

%matplotlib inline


    

In [2]:

    

import seaborn as sns
sns.set_palette('muted')
sns.set_style('darkgrid')


    

In [4]:

    

#2 vectors of random numbers
a = random.randn(10)
b = random.randn(10)

#initialize temporary matrix
pointwise_result = np.zeros(len(a))

for ii in xrange(len(a)):
    pointwise_result[ii] = a[ii] * b[ii]
    
dotproduct = np.sum(pointwise_result)

#above code is useful if unfamiliar with how a dot product works
#following is bit more elegant:
dotproduct = np.sum(a*b)

#better yet
np.dot(a,b)


    

In [17]:

    

#impulse function; all zeros with one 1 in the middle
impfun = np.zeros(100)
impfun[49] = 1

#figure in book uses the following boxcar function
impfun[44:54] = 1

kernel = np.array([1,0.8,.6,.4,.2])

#convolution result using built-in numpy function
numpy_conv_result = np.convolve(impfun,kernel,mode="same")

fig=plt.figure()
fig.add_subplot(share_x=True, share_y = True)

#plot signal
plt.subplot(311)
plt.plot(impfun)

#plot kernel
plt.subplot(312)
plt.plot(kernel,'.-')

#plot result of convolution
plt.subplot(313)
_=plt.plot(numpy_conv_result)


    

In [23]:

    

#zero-pad our signal
zz = np.zeros(len(kernel)-1)
dat4conv = np.concatenate([zz,impfun,zz])

half_of_kernel_size = int(np.ceil((len(kernel) -1)/2))

#initialize convolution output
convolution_result = np.zeros(len(impfun)+len(kernel)-1)

#run convolution (kernel is flipped backwards)
for ti in xrange(len(convolution_result)-half_of_kernel_size):
    convolution_result[ti] = np.sum(dat4conv[ti:ti+len(kernel)]*kernel[::-1])

#cut off edges
convolution_result = convolution_result[half_of_kernel_size:-half_of_kernel_size]

plt.figure()
plt.plot(impfun)
plt.plot(convolution_result,'g')
plt.plot(convolution_result/np.sum(kernel),'r')
plt.plot(numpy_conv_result/np.sum(kernel),'ko')
plt.axis([0,100,-.1,3.1])
_=plt.legend(["raw","unscaled convolution","manual wavelet convolution","numpy convolution"])