How to plot the continuous convolution result
显示 更早的评论
Hi,
I was trying to plot my fourier transform of a continuous function using fplot, but it doesn't give me the expected result (either a straight line or blank),
can anyone check where I went wrong?
thanks a lot!
there's my code:
clear; clc; clear all;
syms f f1 t
%jammer 1
f1_L = f1 - 2.5*10^6;
f1_H = f1 + 2.5*10^6;
jammer1 = rectangularPulse(f1_L,f1_H,f);
jammer1_time_domain = simplify(ifourier(jammer1, f, t))
% jammer 2
f2 = f1 + 12.5*10^6;
f2_L = f2 - 2.5*10^6;
f2_H = f2 + 2.5*10^6;
jammer2 = rectangularPulse(f2_L,f2_H,f);
jammer2_time_domain = simplify(ifourier(jammer2, f, t))
% signal
f_signal = f1 - 12.5*10^6;
% convolution J1 and J1 on fre domain = multiplicaiton on time domian
fv_1_1_time_domain = jammer1_time_domain * jammer1_time_domain
fv_1_1_freq_domain = simplify(ifourier(jammer1_time_domain * jammer1_time_domain, f, t));
% convolution J1 and J1 and J2
fv_1_1_2 = simplify(fourier(jammer2_time_domain * fv_1_1_time_domain, f, t))
fplot(fv_1_1_freq_domain)
2 个评论
Torsten
2023-2-21
How do you want to plot a function that is 0 everywhere and Inf at t=0 ?
YAMENG SONG
2023-2-21
I think I might do the convolution wrong?
I'd like to convolute jammer_1 * jammer_1*jammer_2 (which are rectangular functions) on freq domain, they why I define:
clear; clc; clear all;
syms f f1 t
f1_L = f1 - 2.5*10^6;
f1_H = f1 + 2.5*10^6;
f2 = f1 + 12.5*10^6;
f2_L = f2 - 2.5*10^6;
f2_H = f2 + 2.5*10^6;
jammer1 = rectangularPulse(f1_L,f1_H,f);
jammer2 = rectangularPulse(f2_L,f2_H,f);
since convolution on freq domian equlals to multiplicaiton on time domian, so I used ifourier to convert them back.
jammer1_time_domain = simplify(ifourier(jammer1, f, t))
jammer2_time_domain = simplify(ifourier(jammer2, f, t))
then I multipliy them together, and converted back to freq domain.
% convolution J1 and J1 on fre domain = multiplicaiton on time domian
fv_1_1_time_domain = jammer1_time_domain * jammer1_time_domain
fv_1_1_2 = simplify(fourier(jammer2_time_domain * fv_1_1_time_domain, f, t))
It shouldn't be a delta function as I got, that's why I'm so confused about.
It should be like what circled out in the picture:
采纳的回答
Hi YAMENG,
Have you considered computing the convolution integrals directly?
Note scaling to get them all to fit on the same plot.
syms f f1 t
f1_L = f1 - 2.5*10^6;
f1_H = f1 + 2.5*10^6;
f2 = f1 + 12.5*10^6;
f2_L = f2 - 2.5*10^6;
f2_H = f2 + 2.5*10^6;
jammer1(f,f1) = rectangularPulse(f1_L,f1_H,f);
jammer2(f,f1) = rectangularPulse(f2_L,f2_H,f);
syms u
j11(f,f1) = int(jammer1(u,f1)*jammer1(f-u,f1),u,-inf,inf);
j211(f,f1) = int(jammer2(u,f1)*j11(f-u,f1),u,-inf,inf);
figure
hold on
fplot(jammer1(f,10e6),[0 6e7])
fplot(jammer2(f,10e6),[0 6e7])
fplot(j11(f,10e6)/1e6,[0 6e7])
fplot(j211(f,10e6)/1e12,[0 6e7])
3 个评论
YAMENG SONG
2023-2-22
YAMENG SONG
2023-2-22
The value of f1 was not specified, so I carried it along as an additional argument. Doing so made it easier to make the plots for a specific value of f1, which I arbitrarily chose as 10e6. There are alternative approaches.
Also, the fourier approach can work as well, if you give Matlab a little bit of help via rewrite. Here, I won't use f1 as an additional argument.
syms f f1 t real
f1_L = f1 - 2.5*10^6;
f1_H = f1 + 2.5*10^6;
f2 = f1 + 12.5*10^6;
f2_L = f2 - 2.5*10^6;
f2_H = f2 + 2.5*10^6;
jammer1(f) = rectangularPulse(f1_L,f1_H,f);
jammer2(f) = rectangularPulse(f2_L,f2_H,f);
J1(t) = rewrite(simplify(ifourier(jammer1(f),f,t)),'exp')
J2(t) = rewrite(simplify(ifourier(jammer2(f),f,t)),'exp')
j211(f) = fourier(J2(t)*J1(t)*J1(t),t,f)
Divide by 1e12 to compare to previous result.
fplot(subs(j211(f),f1,10e6)/1e12,[0 6e7],'Color',[4.940000000000000e-01 1.840000000000000e-01 5.560000000000000e-01])
The reason this result is smaller is that the default paramters for fourier as defined in sympref assume that the frequency is in rad/sec so the ifourier include a 1/(2*pi) out front. Because we used ifourier twice we need to multiply this result by (2*pi)^2, which I think will make the peak of the plot match the peak of the result computed via the convolution integral.
Also, the original code had an error in the last line. Repeating from the Question:
clear all;
syms f f1 t
f1_L = f1 - 2.5*10^6;
f1_H = f1 + 2.5*10^6;
f2 = f1 + 12.5*10^6;
f2_L = f2 - 2.5*10^6;
f2_H = f2 + 2.5*10^6;
jammer1 = rectangularPulse(f1_L,f1_H,f);
jammer2 = rectangularPulse(f2_L,f2_H,f);
jammer1_time_domain = simplify(ifourier(jammer1, f, t))
jammer2_time_domain = simplify(ifourier(jammer2, f, t))
% convolution J1 and J1 on fre domain = multiplicaiton on time domian
fv_1_1_time_domain = jammer1_time_domain * jammer1_time_domain;
This line is incorrect. The f and t arguments are backwards
fv_1_1_2 = simplify(fourier(jammer2_time_domain * fv_1_1_time_domain, f, t))
It should be
fv_1_1_2 = simplify(fourier(jammer2_time_domain * fv_1_1_time_domain, t, f))
But, as can be seen, the result is not useful when jammer*_time_domain are expressed in terms of sin and cos instead of in terms of complex exponentials.
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Calculus 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
