plotting a function using a function

Hello everyone, I like to ask how do i plot a function that integrates a certain function
Cn = integral(function of the sawtooth * exp(-1i*2*pi*n*f*t), 0, T);
so far this is my working code which only plots the sawtooth wave function. I have made multiple changes to it trying to figure out how to plot Cn. I'm still fairly new to matlab. maybe I'm missing something.
close all; figure; t=linspace(0,3,1500); m=linspace(0,3,1500); T=1; f=1/T; amp=1; sawf=amp/2; k=plot(NaN,NaN); 1i; for n=1:1:500 sawf = sawf - (amp/pi)*(1/n)*(sin(2*pi*n*f*t)); set(k, 'XData',t,'YData',sawf); pause(0.01); end
%cn1= integral(sawf*((cos(2*pi*1*f*t))-(-(1i)*sin(2*pi*1*f*t))),0,T); %cn1= sawf*(sin(2*pi*n*f*t)/(2*pi*n*f))+((i)cos(2*pi*n*f*t)/(2*pi*n*f)); %cn = (1./T)*cn1; %set(k, 'XData',m,'YData',cn);

3 commentaires

here's the code in proper format. for some reason cloudflare won't let me post on my pc so i was forced to post using my phone
close all;
figure;
t=linspace(0,3,1500);
m=linspace(0,3,1500);
T=1;
f=1/T;
amp=1;
sawf=amp/2;
k=plot(NaN,NaN);
1i;
for n=1:1:500
sawf = sawf - (amp/pi)*(1/n)*(sin(2*pi*n*f*t));
set(k, 'XData',t,'YData',sawf);
pause(0.01);
end
%cn1= integral(sawf*((cos(2*pi*1*f*t))-(-(1i)*sin(2*pi*1*f*t))),0,T);
%cn1= sawf*(sin(2*pi*n*f*t)/(2*pi*n*f))+((i)cos(2*pi*n*f*t)/(2*pi*n*f));
%cn = (1./T)*cn1;
%set(k, 'XData',m,'YData',cn);
Dyuman Joshi
Dyuman Joshi le 9 Sep 2023
Modifié(e) : Dyuman Joshi le 9 Sep 2023
integral requires a function handle as an integrand. You have not supplied any function, just a numerical value.
Also, I don't understand how this -
%1
cn1= integral(sawf*((cos(2*pi*1*f*t))-(-(1i)*sin(2*pi*1*f*t))),0,T);
results in this -
%2
cn1= sawf*(sin(2*pi*n*f*t)/(2*pi*n*f))+((i)cos(2*pi*n*f*t)/(2*pi*n*f));
or how did you arrive at line 2.
or was it supposed to be n instead of 1 (next to pi inside the trignometric terms) in line 1?
jvfa
jvfa le 9 Sep 2023
Hi, can you elaborate where i where i got wrong on what you said about integral? because it's really a big part of understanding the matlab function of integration that I lack.
those comments don't worry about them they were just my attempts at understanding the function. %1 was originally n then i tried substituting 1 for n to see what the value is sort of doing trial and error. for %2 i manually integrated the thing by hand, given that i have the value for sawf so i can move it out of the integral sign then just integrated exp(-2i*pi*n*f*t)dt instead

Connectez-vous pour commenter.

 Réponse acceptée

Bruno Luong
Bruno Luong le 9 Sep 2023
Modifié(e) : Bruno Luong le 9 Sep 2023
The coefficients Cn = 1i/(2*pi*n).
I'm not very good in symbolic calculation, I never own the tollbox license.
syms t
syms n integer
sawtooth = t;
Cn = simplify(int(sawtooth*exp(-1i*2*pi*n*t), t, 0, 1))
Cn = 
C0 = simplify(int(sawtooth, t, 0, 1)) % funny the above general formula is wrong for n=0
C0 = 
Cnfun = matlabFunction(Cn);
n = -10:10;
Cnval = Cnfun(n);
Cnval(n==0) = subs(C0); % no idea why I need to do this
figure
hold on
plot(n, real(Cnval))
plot(n, imag(Cnval))
legend('real','imag')
xlabel('n')
ylabel('C(n)')
t = linspace(0,3,1000)';
f = @(t) mod(t,1)
f = function_handle with value:
@(t)mod(t,1)
fFourier = sum(Cnval.*exp(1i*2*pi*n.*t),2);
figure
plot(t, f(t), 'r', t, fFourier, 'b')
Warning: Imaginary parts of complex X and/or Y arguments ignored.
xlabel('t')
legend('f', 'Fourier approximation')

5 commentaires

jvfa
jvfa le 9 Sep 2023
I didn't really understand what you said there, also I tried to run your code it's a fourier approximation plot? I'm trying to plot something else i think
Bruno Luong
Bruno Luong le 9 Sep 2023
"I have a homework that asks me to plot Cn with respect to n"
that's the first plot, Cn vs n
Paul
Paul le 9 Sep 2023
Modifié(e) : Paul le 9 Sep 2023
Hi Bruno,
Regarding the concern about the general formula being wrong ...
syms t
syms n integer
sawtooth = t;
Don't use simplify to see what's happening (always safer to use sym(pi) )
Cn = int(sawtooth*exp(-1i*2*sym(pi)*n*t), t, 0, 1)
Cn = 
int comes up with a general solution that doesn't take into account the special case of n = 0. However, note that it's sort of the right answer for n = 0 in the sense that
limit(Cn,n,0)
ans = 
This situation comes up all of the time with Fourier integrals. I'm not sure if it's reasonable to expect that int identify all of the special cases for all parameters in the integrand and then provide a result in terms of piecewise that captures each special case. I've often wondered about this issue.
I don't know the inner workings of simplify, but I'll note that we can do this
[num,den] = numden(Cn)
num = 
den = 
num = simplify(num)
num = 
And then
num/den
ans = 
which is standard for the Symbolic toolbox, i.e., n/n will always be replaced with 1, even if n = 0 is a possibility.
C0 = simplify(int(sawtooth, t, 0, 1)) % funny the above general formula is wrong for n=0
C0 = 
What I do is identify the roots of the denominator of the output of int, before doing any simplification, and then do a piecewise and limit. The roots can be found using solve, for example, but here it's easy to do by hand
Cn = piecewise(n == 0, limit(Cn,n,0), Cn)
Cn = 
Then simplify at the end
Cn(n) = simplify(Cn)
Cn(n) = 
The drawback to this approach is that piecewise is not supported by matlabFunction
try
Cnfun = matlabFunction(Cn);
catch ME
ME.message
end
ans = 'Unable to generate code for piecewise for use in anonymous functions. Use option 'File' to write code for piecewise to a file instead.'
So I just evaluate Cn symbolically and convert to double (if numeric is preferred)
n = -10:10;
%Cnval = Cnfun(n);
Cnval = double(Cn(n));
% subs wasn't needed here
% Cnval(n==0) = subs(C0); % no idea why I need to do this
figure
hold on
stem(n, real(Cnval))
stem(n, imag(Cnval))
legend('real','imag')
xlabel('n')
ylabel('C(n)')
t = linspace(0,3,1000)';
f = @(t) mod(t,1);
fFourier = sum(Cnval.*exp(1i*2*pi*n.*t),2);
figure
plot(t, f(t), 'r', t, fFourier, 'b');
Warning: Imaginary parts of complex X and/or Y arguments ignored.
xlabel('t')
legend('f', 'Fourier approximation')
Bruno Luong
Bruno Luong le 9 Sep 2023
@Paul thanks
jvfa
jvfa le 10 Sep 2023

hi bruno, after studying your comment and code i finally understand. thank you, alsp for paul for the correction

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Programming dans Centre d'aide et File Exchange

Produits

Version

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by