MATLAB Answers

Calculation of Convolution Integral Over a Positive Domain in MATLAB

131 views (last 30 days)
Ali Pedram
Ali Pedram on 21 Jun 2021
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 results of the two methods don't overlap as it is seen from the figure below.
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 overlap.
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.

Answers (1)

Harikrishnan Balachandran Nair
From my understanding, you have tried convolution using both the inbuilt function in matlab, and your own manually written code, and the results did not seem to overlap for both.
In this example, both the inputs are vectors of length 101, hence the output of convolution 'y' should be of length 201. The third parameter in the 'conv' function determines the shape of the output function. As you have specified it to be 'same', the function will return the central portion of your output, and the length of output will be same as the length of 'x' which is 101 here. Instead , by specifying it to 'full', which is also the default value, you can get the whole output, which will be of length 201. You can refer to the following code .
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),'full')*dx;
plot(x,y(1:101));

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by