Hello!
I wrote a program in Matlab, which draws a spectrum for a periodic signal and a graph of the convergence of this signal. But for some reason, for some signals the correct graph of convergence is obtained, and for some it is not, although the algorithm is the same (see the photo). I would be very grateful if anyone would advise how to fix it.
syms tk T=input('Vvedi period T sek. ='); n=input('skolko harmonik nado poschitatj? = '); s1=input('Opishi s1(t) funkciju = '); s2=input('Opishi s2(t) funkciju(esli net-pishi 0) = '); s3=input('Opishi s3(t) funkciju(esli net-pishi 0) = '); %======================================================================== %=================== пределы интегрирования ================== a=input('s1(t) funkcija nachinajetsa(moment vremeni): = '); b=input('s1(t) funkcija konchajetsa(moment vremeni): = '); c=input('s2(t) funkcija nachinajetsa(moment vremeni): = '); d=input('s2(t) funkcija konchaetsa(moment vremeni): = '); % e=input('s3(t) funkcija nachinajetsa(moment vremeni): = '); f=input('s3(t) funkcija konchaetsa(moment vremeni): = '); %======================================================================== Lidzkomponente=(1/T)*int(s1,t,a,b)+(1/T)*int(s2,t,c,d)+(1/T)*int(s3,t,e,f) %Lidzkomponente - это постоянная составляющая %======================================================================== %======================= A koef ===================== for i=1:n, A(i)=(2/T)*int(s1*cos(((i*2*pi)/T)*t),t,a,b)+... (2/T)*int(s2*cos(((i*2*pi)/T)*t),t,c,d)+... (2/T)*int(s3*cos(((i*2*pi)/T)*t),t,e,f); end; A=double(A)% Rezultātu izvade %======================================================================== %======================= B koef ===================== for i=1:n, B(i)=(2/T)*int(s1*sin(((i*2*pi)/T)*t),t,a,b)+... (2/T)*int(s2*sin(((i*2*pi)/T)*t),t,c,d)+... (2/T)*int(s3*sin(((i*2*pi)/T)*t),t,e,f); end; B=double(B) %======================================================================== %==================== считаем амплитуды гармоник ================= for i=1:n, C(i)=sqrt(A(i)^2+B(i)^2); end; C=double(C) %======================================================================== %=========================== спектр ========================== FigHandle = figure('Position', [170,120, 1350,550]); subplot(1,2,1); x=[0,1/T:1/T:i/T]; CC=[Lidzkomponente,C]; stem(x,CC,'LineWidth',2) set(gca,'xtick', [1/T:1/T:i/T]); title('Frequency spectrum');xlabel('f,Hz');ylabel('Amplitude'); grid on %======================================================================== %=========================== считаем фазы ========================= for i=1:n, Phase_in_degrees(i)=-atand(B(i)/A(i)); end; Phase_in_degrees(isnan(Phase_in_degrees))=0 %========================================================================= %=========================== фазовый спектр =============================== subplot(1,2,2) x=1/T:1/T:i/T; y=Phase_in_degrees; stem(x,y,'LineWidth',2); title('Phase spectrum');xlabel('f,Hz');ylabel('\Psi,°') set(gca,'xtick', [1/T:1/T:i/T]); grid on %========================================================================= % сходимость функции solis=T/1000; % шаг tt=0:solis:3*T;% временная ось Phs_in_pi=deg2rad(Phase_in_degrees); for i=1:n, harm=C(i)*cos(((i*2*pi/T).*tt)+Phs_in_pi(i)); H{i}=harm; end, catH=cat(3,H{:}); Fr=sum(catH,3); figure(3) Frr=[Lidzkomponente+Fr]; plot(tt,Frr);grid on;xlabel('t,sec') ylabel('S(t)') 
