I am trying to learn a bit of signal processing , specifically using Python. Here's a sample code I wrote.
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import deconvolve
a = np.linspace(-1,1,50)
b = np.linspace(-1,1,50)**2
c = np.convolve(a,b,mode='same')
quotient,remainder = deconvolve(c,b);
plt.plot(a/max(a),"g")
plt.plot(b/max(b),"r")
plt.plot(c/max(c),"b")
plt.plot(remainder/max(remainder),"k")
#plt.plot(quotient/max(quotient),"k")
plt.legend(['a_original','b_original','convolution_a_b','deconvolution_a_b'])
In my understanding, the deconvolution of the convoluted array should return exactly the same array 'a' since I am using 'b' as the filter. This is clearly not the case as seen from the plots below.
I am not really sure if my mathematical understanding of deconvolution is wrong or if there is something wrong with the code. Any help is greatly appreciated!
You are using mode='same'
, and this seems not to be compatible with scipy deconvolve. Try with mode='full'
, it should work much better.
Here the corrected code:
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import deconvolve
a = np.linspace(-1,1,50)
b = np.linspace(-1,1,50)**2
c = np.convolve(a,b,mode='full')
quotient,remainder = deconvolve(c,b)
plt.plot(a,"g")
plt.plot(b,"r")
plt.plot(c,"b")
plt.plot(quotient,"k")
plt.xlim(0,50)
plt.ylim(-6,2)
plt.legend(['a_original','b_original','convolution_a_b','deconvolution_c_b'])