CODEWITHSUNDEEP
functionmatlabmathsignal-processingconvolution
Ali Pedram
Ali Pedram

Reputation: 155

Calculation of Convolution Integral Over a Positive Domain in MATLAB

Suppose we want to compute the convolution of sine and a Gaussian functions from 0 to 10. Below you can see a code that does this using two methods. In the first method the MATLAB function conv() is used and in the second method it is calculated directly.

clear all
clc

x = linspace(0, 10, 101);
dx = x(2) - x(1);

a = @(x) sin(x);
b = @(x) -exp(-x.^2);

y = conv(a(x), b(x),'same') * dx;

t = linspace(-10, 10, 100); 
z = zeros(size(x));
for i = 1:length(x)
    uu = a(t).*b(x(i)-t);
    z(i) = trapz(t,uu);
end

figure(1)
hold on
plot(x, y, 'DisplayName','conv()')
plot(x, z, 'DisplayName','direct')

The two methods give different results as it is seen from the figure below.

enter image description here

However if we use conv() from -10 to 10 and after doing the calculation only keep the part of the result that corresponds to the interval 0 to 10, the two methods give the same results.

clear all
clc

x = linspace(0, 10, 101);
xf = linspace(-10,10,201);
dx = x(2) - x(1);
x0 = find(xf==0);

a = @(x) sin(x);
b = @(x) -exp(-x.^2);

y = conv(a(xf), b(xf),'same') * dx;
y = y(x0:end);

t = linspace(-10, 10, 100); 
z = zeros(size(x));
for i = 1:length(x)
    uu = a(t).*b(x(i)-t);
    z(i) = trapz(t,uu);
end

figure(1)
hold on
plot(x, y, 'DisplayName','conv()')
plot(x, z, 'DisplayName','direct')

As it is seen from the figure below the results are the same.

enter image description here

I have two questions:

1- Why conv() doesn't work for the positive domain and only when the negative numbers are included we get the right result?

2- Is there a way to use only the function conv() and the interval 0 to 10 and get the same result? I want to use convolution integrals in a very complicated code and the direct method is very slow so I want to use conv(). However the structure of the code is such that I can't use the negative domain as I have used in this example.

Upvotes: 1

Views: 955

Answers (1)

Cris Luengo
Cris Luengo

Reputation: 60780

Plot the two sampled signals you’re convolving. You’ll notice that you have only half a Gaussian. When using the (-10,10) interval, you sample the full Gaussian, and hence get a correct result.

To get the right result, you need to shift the Gaussian so that it is fully sampled:

y = conv(a(x), b(x-5),'same') * dx;

If you shift the Gaussian such that its origin is in the middle sample, the output function will not seem shifted, conv assumes the origin is in the middle of the arrays.

Upvotes: 2

Related Questions