实验三 周期信号的傅里叶级数分析
是由于高频部分信号的损失。这就导致在构成有跃变的连续时间周期函数时，在跃变点的附近存在一个幅度大约为9%的过冲，且所选谐波次数越多，过冲点越向不连续点靠近。这种现象称为Gibbs现象。
若在计算机上编程对周期函数f(t)进行傅里叶级数展开，必须对函数f(t)作等间距抽样。若抽样周期为Ts，且令T1?NTs，则??k?1?k2?，(3.3)式离散化为 NTs?a0??2???2???f?nTs??f?t?|t?nTs????akcos?nk??bksin?nk?? (3.12)
2k?1??N??N??时间抽样后，(3.11)式离散化为
Ka0??2???2???f?n?????akcos?nk??bksin?nk?? (3.13)
2k?1?NN????? 将上式与(3.6)式进行比较，实际傅里叶级数展开式f?n?与f?nTs?之间的误差为
??n,K????akcos?k?K??2?2??nk?bksinnk? (3.14) NN?n?t?nTs?和截断频率K??c?K?1?的函数。
上式表明，实际展开后的误差是时间
下图给出了一个方波信号展开成有限长傅里叶级数后，在跃变点的附近产生的Gibbs现象，
而且不连续的跃变点也扩展成了有一定上升空间的连续函数。
图3.1 方波展开成有限长傅里叶级数后，在跃变点的附近产生Gibbs现象
3. 吉布斯现象的消除
为了消除这种频域截断形成的Gibbs现象，通常不采用矩形窗作截断处理，而是采用汉宁窗、海明窗或三角窗等进行加权计算。
1）以0点为中心的汉宁(Hanning)窗（也称为升余弦窗）定义为
?1?2?k???1?cos?k?K/2w?k???2? (3.15) K?otherwise?0?17
实验三 周期信号的傅里叶级数分析
2) 以0点为中心的海明(Hamming)窗定义为
2?kk?K/2??0.54?0.46cos (3.16) w?k???Kotherwise?0?3) 以0点为中心的三角(Bartlett)窗定义为
?2kk?K/2?1? (3.17) w?k???Kotherwise??0下图列出了四种窗的图像，可进行比较：
图3.2 几种加权窗函数的比较
图3.3 用Hanning窗加权后方波傅里叶级数展开的跃变点附近的Gibbs现象的消除
例如图3.1中的方波信号展开式用Hanning窗加权截断后，图像如图3.3所示，显然Gibbs现象已经基本消除。
采用频域Hanning窗加权或Hamming窗加权的方法进行截断，与矩形截断相比，可以减弱或消除Gibbs现象，但不会减小由于频域截断产生的误差，反而因加权导致截取区域内频谱发生变化，增大了误差。
三、实验内容
1. 绘制周期T1=1、幅度E=1的对称方波的前十项傅里叶级数的系数(三角函数形式)，并
用前五项恢复原信号。
解 本题中函数为偶对称，因而正弦分量全部为0，由式(3.3)可写出前K项(未计入直流项)
18
实验三 周期信号的傅里叶级数分析
级数的合成公式
f?t0?cos??0t0??cos?K?0t0????1??a0??f?t??t???1??a?????????cos?t??t?cosK?t??t00000???????1? ???????????????????????????ft?T??t1cos?t?T??t?cosK?t?T??t001001?01????aK?对应程序：
clc
clear all
E=1; %定义方波幅度 T1=1; %定义方波周期 omg0=2*pi/T1; %定义基频
N=1000; %定义时域抽样点数
t=linspace(-T1/2,T1/2-T1/N,N)'; %生成时域抽样点 f=0*t; %初始化时域信号 f(:)=-E/2;
f(t>-T1/4&t k1=-10; %确定系数的起止下标 k2=10; k=[k1:k2]'; F=1/N*exp(-j*kron(k*omg0,t.'))*f; %指数形式的傅里叶级数系数 a0=F(11); %直流分量 本来是a0/2，这里为了方便，直接表示成a0 ak=F(12:21)+F(10:-1:1); fs=cos(kron(t,[0:5]*omg0))*[a0;ak(1:5)]; %合成原函数 subplot(2,1,1) plot(t,f,t,fs,'--') subplot(2,1,2) stem(k(11:21),[a0,ak']) %三角形式傅里叶系数图只有正半轴 2. 画出如下图对称方波（取E=1、T=1），并采用有限项傅里叶级数对原函数进行逼近， 画出对称方波的1、3、5、7、9、11次谐波的傅里叶级数合成波形，观察吉布斯现象。 19 实验三 周期信号的傅里叶级数分析 f(t)-T/4T/4T 对应的程序以及结果图分别显示如下： sum=0; t=-3:0.01:3; E=1;T=1;ta=T/2;w=2*3.14159/T; for n=1:1 fn=(2*E*ta/T)*sin(w*ta*n/2)/(w*ta*n/2); f=(E*ta/T)+cos(n*w*t)*fn-E/2; %由于没有直流分量，将原公式中减去直流分量的差值得 f?t???ET??2?EESa?n?0?/2?cos?n?0t?? 2n?1T?sum=sum+f; end plot(t,sum) 20 实验三 周期信号的傅里叶级数分析 sum=0; t=-3:0.01:3; E=1;T=1;ta=T/2;w=2*3.14159/T; for n=1:3 fn=(2*E*ta/T)*sin(w*ta*n/2)/(w*ta*n/2); f=(E*ta/T)+cos(n*w*t)*fn-E/2; sum =sum+f; end plot(t,sum) 21