Understanding convolution, 'same' might give wrong result

112 views (last 30 days)
Silke on 7 Nov 2017
Commented: Silke on 10 Nov 2017
I am having a problem with the convolution function. I am trying to convolute two functions, one is bi-exponential decay and the other one a pulse shape. I did it like that:
t = 0:250e-12:1e-6;
a1 = 0.6;
a2 = 0.4;
D1 = 25e-9;
D2 = 350e-9;
filename1 = 'pulse410.txt';
pathname1 = 'D:\My Documents\PulseShapes\';
pulse410 = importfile1([pathname1, filename1]);
[P,I] = max(pulse410(:,1));
pulse410s=spline(pulse410(:,1), pulse410(:,2), t); % pulse410s contains the same amount of data points as t.
t_in = find(t==P);
for ii = t_in:numel(t)
pulse410s(ii)=0;
end
f1 = conv(a1 * exp(-t/D1) + a2 * exp(-t/D2), pulse410s./max(pulse410s), 'same') ;
In the picture below, data 1 is the bi-exponential decay, data 2 is pulse410s, data3 is the convoluted signal. Well, I would expect the convolution of the two signals to look different. Especially, the drop at 0.5e-6 looks weird to me. Is this weird result due to me using the 'same'-option in the convolution? How can I obtain a good convolution by keeping the amount of data points the same before and after the convolution? Thanks for your help!

David Goodmanson on 7 Nov 2017
Hi Silke,
Yes, this is a direct result of using the 'same' option. If you take a look at the convolution without that option, all 'same' is doing is taking the middle section, between 1/4 and 3/4 of the x axis and throwing away the rest. But some of the rest has needed information.
When you do a convolution of arrays with n1 and n2 points, the result has n1+n2+1 points. Right now you are zerofilling the pulse to give it the same number of points as the exponentials, 4001 points. There is no need to do that. It appears that the pulse has a nonzero extent of about 65 points after interpolation. If you keep an array of that size and convolve it with the exponentials, you will have an array of about 4067 points or so with the same information.
Starting from somewhere in negative time, the exponentials have a sharp rise at t=0 which is real, and then a sharp drop at the end because your time record runs out. The second drop is an artifact. When you do the convolution to get 4067 or whatever points, the part at the beginning is real and the part at the end is not. So if you want the convolution to have the same length as the signal, the best approximation is to take the first 4001 points in the convolution.
Silke on 8 Nov 2017
Yes, this is actually what I was looking for. But I am not sure if it is really required for my goal. I am quite puzzled myself at the moment. Well, I have put the convolution problem now in a bigger context (fitting a function). If you have time/interest, it would be great if you could have a look at it. https://nl.mathworks.com/matlabcentral/answers/365865-nlinfit-gives-no-result Thank you!

Eric on 7 Nov 2017
Your convolution looks correct, at least for what you have done, and yes, it is partially because of the 'same' option. The other part of the reason you are seeing the behavior you are is due to your extending the pulse waveform by appending zeros to it.
Let us assume I have a pulse that is 10ms long and my decay is simply e^-t. Here are several different ways to use MATLAB's convolution function:
t = 0:0.001:1;
f0 = conv(exp(-t),[ones(10,1)]); %The default, which uses 'full'
f1 = conv(exp(-t),[ones(10,1)],'same'); %Just the pulse, 'same'
f2 = conv([ones(10,1)],exp(-t),'same'); %Just the pulse, 'same' but inputs reversed
f3 = conv(exp(-t),[ones(10,1); zeros(length(t)-10,1)],'same'); %Pulse appended with zeros, 'same'
f4 = conv(exp(-t),[ones(10,1); zeros(length(t)-10,1)]); %Pulse appended with zeros, 'full'
f5 = conv(exp(-t),[ones(10,1)],'valid'); %Just the pulse, 'valid'
What you will find when comparing f1 and f2 to f0, or f3 to f4, is that the argument 'same' basically takes the center n points from the full convolution, where n is the length of the first vector. This is expected and how it is defined in the documentation. Now, when comparing the plots of f0 and f4, you can see that the 'full' result of f4 is also extended by a number of zeros relative to f0. As such, when you call for 'same' using the f4 method, your center points will be different than if you had not appended the zeros (f0). MATLAB will do any zero padding for you internally, so there is no need to do it yourself. P.S. There is also the option of 'valid', (f5) which may interest you, which will return only those results where it did not need to use zero-padding. Basically, It will remove the sharp drops at the edges of the convolution. If you are concerned about it, the output vector's length can be obtained based on the lengths of the inputs (so adding zeros will change this). See conv:shape for more info.
Silke on 10 Nov 2017
Hi Eric and David,
thanks for your reply. Well, I am actually interested in getting an answer with the correct units, as I am trying to fit measurement data. So I really believe that I need to multiply with dx. Thanks a lot for your help!