根据傅里叶级数的定义我们知道:
对于任意一个周期为的周期信号
,都可以求出它在三角函数集中各函数中的分量,从而可将
在区间
内表示为三角函数集中各函数的加权和。即:
其中:
借助Matlab编写M文件,以方波信号为例,选取不同的级数项数进行合成,观察与原方波的逼近情况。
源码如下:
%该m文件用于绘制方波信号及其傅立叶级数谐波叠加
clc;clear;close all;
syms x; %定义符号变量x,用于求解an、bn
A = input('请输入方波信号的幅度:');
T = input('请输入方波信号的周期:');
PW = input('请输入方波信号的脉宽:');
c = input('请问您要显示几个周期的方波:');
m = input('请问您需要几次谐波叠加:');
index = 1; %定义方波信号下标
sw = ones(1,c*T*100+1); %预定义方波信号,可提高Matlab运行速率
for N = -c/2:c/2-1 %遍历c个周期
%绘图过程中舍弃重复部分
if N==-c/2
i = 0;
else
i = 0.01;
end
%一个周期内判断何时为幅度为A,何时为0
for t = N*T+i:0.01:(T+N*T)
if((t>=(T-PW)/2 + N*T)&&(t<=((T+PW)/2 + N*T)))
sw(index) = A;
index = index + 1;
else
sw(index) = 0;
index = index + 1;
end
end
end
t = -c*T/2:0.01:c*T/2;
plot(t,sw,'LineWidth',2); %绘制方波信号
xlabel('t');ylabel('幅度');title(['方波信号及其傅立叶级数' num2str(m) '次谐波叠加']);
hold on;
%以下为傅立叶级数展开部分
w = 2*pi/T;
%预定义an、bn
an = ones(1,m);
bn = ones(1,m);
sum = 0;
a0 = 2/T*A*PW;
for n = 1:m
%根据周期的奇偶选择积分区间
if (mod(c,2)==1)
an(n) = 2/T*double(int(A*cos(n*w*x),-PW/2,PW/2));
bn(n) = 2/T*double(int(A*sin(n*w*x),PW/2,PW/2));
else
an(n) = 2/T*double(int(A*cos(n*w*x),(T-PW)/2,(T+PW)/2));
bn(n) = 2/T*double(int(A*sin(n*w*x),(T-PW)/2,(T+PW)/2));
end
ft = a0/2 + an(n)*cos(n*w*t) + bn(n)*sin(n*w*t);
%避免重复添加a0
if n==1
sum = sum + ft;
else
sum = sum + ft - a0/2;
end
end
plot(t,sum,'r','LineWidth',2); %绘制m次谐波叠加后图像
legend('方波信号',[num2str(m) '次谐波叠加']);
在命令窗口输入对应参数:幅度5,周期4,脉宽2,显示5个周期方波,3次谐波叠加。运行该M文件得结果如下:
将叠加次数换为9次,得如下结果:
再将叠加次数换为15次,得如下结果:
对比以上三幅图,可以得到以下结论:
周期信号(本实验中为方波信号)选取有限次傅里叶级数进行叠加。当选取的项数越多,在所合成的波形中出现的峰起越靠近原信号的不连续点。当选取的项数很大时,该峰起值趋于一个常数,大约等于总跳变值的9%(吉布斯效应)。