q608.m

MATLAB及在电子信息课程中的应用例题索引与程序

%《MATLAB在电子信息课程中的应用》第六章例6.8程序q608
% 非正弦周期信号的频谱分析
% 电子工业出版社出版  陈怀琛 吴大正 高西全合著 2001年10月

clear, format compact
Um=100; T=0.02; w=2*pi*50;
% 方法一：按傅立叶分析定义计算
N=input('取的谐波次数N=  ');	% 取的谐波次数愈高，分段数得愈多
t=linspace(-T/2,T/2,1000);dt=T/999;	% 列出t数组，取100点，有99个间隔dt
u=Um*abs(sin(w*t));		% 在一个周期内生成两个半波的简便方法
for k=0:N						% 循环求系数ak,bk,Ak
   a(k+1)=trapz(u.*cos(k*w*t))*dt/T*2;	% 变量下标不得取0，故将k加1
   b(k+1)=trapz(u.*sin(k*w*t))*dt/T*2;
   A(k+1)=sqrt(a(k+1)^2+b(k+1)^2);
end
[[0:N]',[A(1)/2,A(2:end)]']		% 显示各傅立叶分量，恢复与k的对应关系
stem(0:N,[A(1)/2,A(2:end)])		% 画出各傅立叶分量与k的对应关系
Us11=sqrt(trapz(u.^2)*dt/T)		% 用总功率求出的有效值
Us12=sqrt(A(1)^2/4+sum(A(2:end).^2/2))	% 用各分量功率和求出的有效值
% 方法二：按推导出的半波傅立叶分量公式计算
Us21=Um*sqrt(trapz(sin(w*t).^2)*dt/T)			% 用总功率求出的有效值
% 用前四个分量功率和求出的有效值
Us22=4*Um/pi*sqrt(0.5^2+0.5*sum((1./(4*[1:3].^2-1)).^2))
e=(Us21-Us22)/Us21				% 相对误差