Problem using symbolic sinc?

Hello everybody,

I have tried to get a plot of the absolute value of the Fourier transform for a square pulse using the equation . The code I used is

syms x y w F

A = 1;
tau = 1;
y(x) = piecewise((x<-tau/2), 0, ...
                 (x==tau/2 & x==tau/2), 0.5, ...
                 (x>-tau/2 & x<tau/2), 1, x>0.5, 0);

figure1 = figure;
axes1 = axes('Parent',figure1);
hold(axes1,'on');
fplot(y,[-1.5 1.5],'Parent',axes1,'MarkerSize',6,'LineWidth',3);
ylabel('\Pi(t)','FontSize',26,'FontName','Calibri');
xlabel('t','FontSize',26,'FontName','Calibri');
ylim(axes1,[0 1.5]);
box(axes1,'on');
hold(axes1,'off');
set(axes1,'FontName','Calibri','FontSize',20,'XAxisLocation','origin',...
    'XTick',[-1/2 0 1/2],'XTickLabel',{'-\tau/2','0','\tau/2'}, ...
    'XLimitMethod','tight','YAxisLocation','origin',...
    'YTick',[0 1],'YTickLabel',{'0','A'},'YLimitMethod','tight',...
    'ZLimitMethod','tight');

F = A*tau*sinc(w/(2*pi)*tau);

figure2 = figure;
axes2 = axes('Parent',figure2);
fplot(abs(F),[-12*pi/tau 12*pi/tau],'Parent',axes2,'MarkerSize',6,'LineWidth',3);
ylim(axes2,[-0.5 A*tau*1.2]);
xlim(axes2,[-13*pi/tau 13*pi/tau])
box(axes2,'on');
hold(axes2,'off');
set(axes2,'FontSize',14,'XAxisLocation','origin','XTick', ...
    [-8*pi/tau -6*pi/tau -4*pi/tau -2*pi/tau 0 2*pi/tau 4*pi/tau 6*pi/tau 8*pi/tau], ...
    'XTickLabel',{'-8\pi/\tau','-6\pi/\tau','-4\pi/\tau','-2\pi/\tau','0',...
    '2\pi/\tau','4\pi/\tau','6\pi/\tau','8\pi/\tau'},'XTickLabelRotation',90, ...
    'YGrid','on','YAxisLocation','origin','YTick',[0 A*tau],'YTickLabel',{'0','A\tau'});
ylabel('|F(\omega)|','FontSize',18,'FontName','Calibri');
xlabel('\omega','FontSize',18,'FontName','Calibri');

The plot of seems to be correct but the gray dashed line is a bit surprising.

Than I tried to get the I have tried to get a plot of the absolute value of the Fourier transform for a pulse of sinusoidal function given by the equation using the foolowing code

syms x t

A = 1;
f=8;

g = A*cos(2*pi*f*t);

figure3 = figure;
axes3 = axes('Parent',figure3);
hold(axes3,'on');
% fplot(heaviside(t + 0.5), [-2, 2],'Parent',axes22,'MarkerSize',6,'LineWidth',3)
% fplot(-heaviside(t - 0.5), [-2, 2],'Parent',axes22,'MarkerSize',6,'LineWidth',3)
fplot(t,g*(heaviside(t + 0.5)-heaviside(t - 0.5)),[-2 2],'Parent',axes3,'MarkerSize',6,'LineWidth',2)
ylim(axes3,[-1.5,1.5]);
xlim(axes3,[-1.5,1.5]);
xlabel("t")
ylabel("f(t)")
hold(axes3,'off');
set(axes3,'FontName','Calibri','FontSize',20,'XAxisLocation','origin',...
    'XTick',[-1 -0.5 0 0.5 1],'XTickLabel',{'-\tau','-\tau/2','0','\tau/2','\tau'}, ...
    'XLimitMethod','tight','YAxisLocation','origin',...
    'YTick',[-1 0 1],'YTickLabel',{'-1','0','1'},'YLimitMethod','tight',...
    'ZLimitMethod','tight');

%%

syms w F

f = 8;
w0=2*pi*f;
tau = 1;

F(w) = tau/2*(sinc((w-w0)/(2*pi)*tau)+sinc((w+w0)/(2*pi)*tau));

figure4 = figure;
axes4 = axes('Parent',figure4);
fplot(abs(F),[-5*f*pi/tau 5*f*pi/tau],'Parent',axes4,'MarkerSize',6,'LineWidth',3);
ylim(axes4,[-0.5 tau]);
xlim(axes4,[-6*f*pi/tau 6*f*pi/tau]);
box(axes4,'on');
hold(axes4,'off');
set(axes4,'FontSize',14,'XAxisLocation','origin','XTick', ...
    [-w0/tau 0 w0/tau],'XTickLabel',{'-\omega_0','0','\omega_0'},...
    'YGrid','on','YAxisLocation','origin','YTick',[0 tau/2],...
    'YTickLabel',{'0','\tau/2'});
ylabel('F(\omega)','FontSize',18,'FontName','Calibri');
xlabel('\omega','FontSize',18,'FontName','Calibri');

Again, there are gray dashed lines, this time for the frequencies and .

This meake me think that there is a problem with the sinc function for argument close to or equal 0.

I have done a simple test

>> syms x

>> f = sinc(x)

f =

sin(pi*x)/(x*pi)

>> x = -5:5

x =

-5 -4 -3 -2 -1 0 1 2 3 4 5

>> subs(f)

Error using symengine

Division by zero.

Error in sym/subs>mupadsubs (line 168)

G = mupadmex('symobj::fullsubs',F.s,X2,Y2);

Error in sym/subs (line 153)

G = mupadsubs(F,X,Y);

but

>> sinc(0)

ans =

1

Then I have repeated the calculations for the sinusoidal pule but this time I did not used the symbolic math and I got a bettrer plot.

f = 8;
w0=2*pi*f;
tau = 1;
w = -5*f*pi/tau:0.1*pi:5*f*pi/tau;

% F(w) = tau/2*(sinc((w-w0)/(2*pi)*tau)+sinc((w+w0)/(2*pi)*tau));
F = tau/2*(sinc((w-w0)/(2*pi)*tau)+sinc((w+w0)/(2*pi)*tau));

figure5 = figure;
axes5 = axes('Parent',figure5);
% fplot((abs(F)),[-5*f*pi/tau 5*f*pi/tau],'Parent',axes23,'MarkerSize',6,'LineWidth',3);
plot(w,abs(F),'Parent',axes5,'MarkerSize',6,'LineWidth',3)
ylim(axes5,[0 tau*0.75]);
xlim(axes5,[-6*f*pi/tau 6*f*pi/tau]);
box(axes5,'on');
hold(axes5,'off');
set(axes5,'FontSize',14,'XAxisLocation','origin','XTick', ...
    [-w0/tau 0 w0/tau],'XTickLabel',{'-\omega_0','0','\omega_0'},...
    'YGrid','on','YAxisLocation','origin','YTick',[0 tau/2],...
    'YTickLabel',{'0','\tau/2'});
ylabel('|F(\omega)|','FontSize',18,'FontName','Calibri');
xlabel('\omega','FontSize',18,'FontName','Calibri');

Is there a proble with the sic function for symbolic calculation?

Answer 1

[Paul]

We can get rid of the vertical asymptote with:

syms x y w F
A = 1;
tau = 1;
F = A*tau*sinc(w/(2*pi)*tau);
figure
fplot(abs(F),[-12*pi/tau 12*pi/tau],'ShowPoles','off')
ylim([0 1])

Or

figure
fplot(w,abs(F),[-12*pi/tau 12*pi/tau])
ylim([0 1])

Of course, either approach also removes the aysmptotes for actual poles that might be of interest, should any exist.

Answer 2

[Walter Roberson]

y(x) = piecewise((x<-tau/2), 0, ...
                 (x==tau/2 & x==tau/2), 0.5, ...
                 (x>-tau/2 & x<tau/2), 1, x>0.5, 0);

(x==tau/2 & x==tau/2) is redundant. If you want to test for equality with tau/2 then use a single condition. If you want to test for it being either tau/2 or -tau/2 then use the correct condition, (x==-tau/2 | x == tau/2)

The x>0.5 seems inconsistent with the rest of the tests. It would make more sense if you were testing x>tau/2

Perhaps you want

y(x) = piecewise(abs(x) < tau/2, 1, ...
                 abs(x) == tau/2, 0.5, ...
                 0);

[BLP]

I thought about making the definition of the rectangular pulse shorter, but the my problem at the moment is its the Fourier transform, which I calculate from the equation, not from definiotion - I do not use the sigal. I have included the part with the signal only to ilustrate what is the Fourier transform calculated for. Anyway, thank you very much for the suggestion. I will definitely use your suggestions in my code.

[Walter Roberson]

Instead of subs(f) use

arrayfun(@(X) limit(f, X), x)

This will be slower.

Or you could create a mask

sincs = zeros(size(x));
mask = x ~= 0;
sincs(mask) = subs(f, x(mask));
sincs(~mask) = 1;

[BLP]

Thank you again.

As far as I unederstand your suggwation it will solve the problem with symboli calculation of which by definition should be , but does it make sense to calculate the limit value for using arrayfun? Moreover, does it not mean that I will ned to calculate the values of the function before I will plot it. Is there any other method to plot the function without prior substitution using the fplot function like in the following example

syms w F

F = A*tau*sinc(w/(2*pi)*tau);
fplot(abs(F),[-12*pi/tau 12*pi/tau]);

and get a proper plot?

The solution with mask will work, but seems unnatural to me.

Bay the way, wuld it not be better to use

sincs = ones(size(x));
mask = x ~= 0;
sincs(mask) = subs(f, x(mask));

And again, this is a substitution which I would like to avoid.