实验二 连续时间信号的频域分析
实验二连续时间信号的频域分析
一、实验目的
1、掌握连续时间周期信号的傅里叶级数的物理意义和分析方法；
2、观察截短傅里叶级数而产生的“Gibbs现象”，了解其特点以及产生的原因；
3、掌握连续时间傅里叶变换的分析方法及其物理意义；
4、掌握各种典型的连续时间非周期信号的频谱特征以及傅里叶变换的主要性质；
5、学习掌握利用Matlab语言编写计算CTFS、CTFT和DTFT的仿真程序，并能利用这些程序对一些典型信号进行频谱分析，验证CTFT、DTFT的若干重要性质。
基本要求：掌握并深刻理傅里叶变换的物理意义，掌握信号的傅里叶变换的计算方法，掌握利用Matlab编程完成相关的傅里叶变换的计算。
二、原理说明
1、连续时间周期信号的傅里叶级数CTFS分析
任何一个周期为T
1
的正弦周期信号，只要满足狄利克利条件，就可以展开成傅里叶级数。
三角傅里叶级数为：
∑∞=+
+
=
1
0 0
)]
sin(
)
cos(
[
)(
k
k
k
t
k
b
t
k
a
a
t
xω
ω 2.1
或： ∑∞
=++=100)cos()(k k k t k c a t x ?ω 2.2 其中1
02T πω=，称为信号的基本频率（Fundamental frequency ），k k b a a ，和，0分别是信号)(t x 的直流分量、余弦分量幅度和正弦分量幅度，k k c ?、为合并同频率项之后各正弦谐波分量的幅度和初相位，它们都是频率0ωk 的函数，绘制出它们与0ωk 之间的图像，称为信号的频谱图（简称“频谱”），k c －0ωk 图像为幅度谱，k ?－0ωk 图像为相位谱。
三角形式傅里叶级数表明，如果一个周期信号x(t)，满足狄里克利条件，就可以被看作是由很多不同频率的互为谐波关系（harmonically related ）的正弦信号所组成，其中每一个不同频率的正弦信号称为正弦谐波分量 (Sinusoid component)，其幅度（amplitude ）为k c 。也可以反过来理解三角傅里叶级数：用无限多个正弦谐波分量可以合成一个任意的非正弦周期信号。
指数形式的傅里叶级数为：
∑∞-∞==
k t jk k e a t x 0)(ω 2.3
其中，k a 为指数形式的傅里叶级数的系数，按如下公式计算： ?
--=2/2/1110)(1T T t jk k dt e t x T a ω 2.4 指数形式的傅里叶级数告诉我们，如果一个周期信号x(t)，满足狄里克利条件，那么，它就可以被看作是由很多不同频率的互为谐波关系（harmonically related ）的周期复指数信号所组成，其中每一个不同频率的周期复指数信号称为基本频率分量，其复幅度
（complex amplitude ）为k a 。这里“复幅度（complex amplitude ）”指的是k a 通常是复数。
上面的傅里叶级数的合成式说明，我们可以用无穷多个不同频率的周期复指数信号来合成任意一个周期信号。然而，用计算机（或任何其它设备）合成一个周期信号，显然不可能做到用无限多个谐波来合成，只能取这些有限个谐波分量来近似合成。
假设谐波项数为N ，则上面的和成式为：
∑-==
N N k t jk k e a t x 0)(ω 2.5
显然，N 越大，所选项数越多，有限项级数合成的结果越逼近原信号x(t)。本实验可以比较直观地了解傅里叶级数的物理意义，并观察到级数中各频率分量对波形的影响包括“Gibbs ”现象：即信号在不连续点附近存在一个幅度大约为9%的过冲，且所选谐波次数越多，过冲点越向不连续点靠近。这一现象在观察周期矩形波信号和周期锯齿波信号时可以看得很清楚。
2、连续时间信号傅里叶变换----CTFT
傅里叶变换在信号分析中具有非常重要的意义，它主要是用来进行信号的频谱分析的。傅里叶变换和其逆变换定义如下：
?
∞
∞--=
dt e t x j X t j ωω)()( 2.6 ?∞∞-=ωωπωd e j X t x t j )(21)( 2.7
连续时间傅里叶变换主要用来描述连续时间非周期信号的频谱。按照教材中的说法，任意非周期信号，如果满足狄里克利条件，那么，它可以被看作是由无穷多个不同频率（这些频率都是非常的接近）的周期复指数信号e j?t 的线性组合构成的，每个频率所对应的周期复指数信号e j?t 称为频率分量（frequency component ），每个频率的幅度为对应频率的|X(j?)|之值，其相位为对应频率的X(j?)的相位。
X(j?)通常为复函数，可以按照复数的极坐标表示方法表示为：
X(j?)=| X(j?)|e j? X(j?)
其中，| X(j?)|称为x(t)的幅度谱，而?X(j?)则称为x(t)的相位谱。
给定一个连续时间非周期信号x(t)，它的频谱也是连续且非周期的。对于连续时间周期信号，也可以用傅里叶变换来表示其频谱，其特点是，连续时间周期信号的傅里叶变换由冲激序列构成的，是离散的——这是连续时间周期信号的傅里叶变换的基本特征。
4.1 傅里叶级数的Matlab 计算
设周期信号x(t)的基本周期为T 1，且满足狄里克利条件，则其傅里叶级数的系数可由
式2.4计算得到。式2.4重写如下：
基本频率为： 102T πω=
对周期信号进行分析时，我们往往只需对其在一个周期内进行分析即可，通常选择主周期（Principle period ）。假定x 1(t)是x(t)中的主周期，则
计算机不能计算无穷多个系数，所以我们假设需要计算的谐波次数为N ，则总的系数个
数为2N+1个。在确定了时间范围和时间变化的步长即T 1和dt 之后，对某一个系数，上述
系数的积分公式可以近似为： ∑?---==n t jk n T T t jk k T dt e t x dt e t x T a 12
/2/11/)()(101
10ωω 对于全部需要的2N+1个系数，上面的计算可以按照矩阵运算实现。Matlab 实现系数计算的程序如下：
dt = 0.01;
T = 2; t = -T/2:dt:T/2; w0 = 2*pi/T;
x1 = input(‘Type in the periodic signal x(t) over one period x1(t)=’);
N = input(‘Type in the number N=’);
k = -N:N; L = 2*N+1;
ak = x1*exp(-j*k*w0*t’)*dt/T;
需要强调的是，时间变量的变化步长dt 的大小对傅里叶级数系数的计算精度的影响非常大，dt 越小，精度越高，但是，计算机计算所花的时间越长。
5、用Matlab 实现CTFT 及其逆变换的计算
5.1 用Matlab 实现CTFT 的计算
Matlab 进行傅里叶变换有两种方法，一种利用符号运算的方法计算，另一种是数值计算，本实验要求采用数值计算的方法来进行傅里叶变换的计算。严格来说，用数值计算的方法计算连续时间信号的傅里叶变换需要有个限定条件，即信号是时限信号（Time limited
signal），也就是当时间|t|大于某个给定时间时其值衰减为零或接近于零，这个条件与前面提到的为什么不能用无限多个谐波分量来合成周期信号的道理是一样的。计算机只能处理有限大小和有限数量的数。
采用数值计算算法的理论依据是：
若信号为时限信号，当时间间隔T取得足够小时，上式可演变为：
上式用Matlab表示为：
X=x*exp(j*t’*w)*T
其中X为信号x(t)的傅里叶变换，w为频率Ω，T为时间步长。
相应的Matlab程序：
T = 0.01; dw = 0.1; %时间和频率变化的步长
t = -10:T:10;
w = -4*pi:dw:4*pi;
X(j?)可以按照下面的矩阵运算来进行：
X=x*exp(-j*t’*?)*T; %傅里叶变换
X1=abs(X); %计算幅度谱
phai=angle(X); %计算相位谱
为了使计算结果能够直观地表现出来，还需要用绘图函数将时间信号x(t)，信号的幅
度谱|X(j?)|和相位谱? X(j?)分别以图形的方式表现出来，并对图形加以适当的标注。
四、实验步骤及内容
6、编写程序，绘制下面的信号的波形图：
其中，?
0= 0.5π，要求将一个图形窗口分割成四个子图，分别绘制cos(?
t)、cos(3?
t)、
cos(5?
t) 和x(t) 的波形图，给图形加title，网格线和x坐标标签。
dt=0.002; %定义时间间隔
t=-2:dt:2; %设置时间起始和截止
w0=0.5*pi; %题目要求，定义wo
x1=cos(w0*t); %输入信号x1
subplot(221); %将信号1的图形放在两行两列中的第一个plot(t,x1),grid on; %输出信号1的图形，并加网格
x2=cos(3*w0*t); %输入信号x2
subplot(222); %将信号2的图形放在两行两列中的第二个plot(t,x2),grid on; %输出信号2的图形，并加网格
x3=cos(5*w0*t); %输入信号3
subplot(223); %将信号3的图形放在两行两列中的第三个
plot(t,x3),grid on; %输出信号3的图形，并加网格
ak=0; %设傅里叶系数的初始值为0
for n=1:9 %循环计算9次傅里叶系数
ak=ak+(-1).^(n+1).*(1/(2*n-1)).*cos((2*n-1).*w0*t);
%计算傅里叶系数
end
subplot(224); %将该图形放在两行两列中的第四个
plot(t,ak),grid on; %输入该图形，以t为横坐标，ak为纵坐标，并加网格7、给定如下两个周期信号：
编写程序，计算x
1(t)和x
2
(t)
共21个系数；
仿照程序Program2_2
12
合成的y
1(t)和y
2
(t)的波形图，以及x
1
(t) 和x
2
(t)的幅度频谱和相位频谱的谱线图。反复
执行该程序，输入不同的N值，观察合成的信号波形中，是否会产生Gibbs现象？为什么？clear,
close all
T=2;dt=0.00001; %定义周期T=2，时间间隔dt为0.00001，
t=-T/2:dt:T/2; %时间轴的起始和截止
x1=u(t+0.2)-u(t-0.2); %输入信号1
x=0; %给x附初值为0
t1=-2*T:dt:2*T; %定义t1时间轴的起始和截止
for m=-2:2
x=x+u(t1+0.2-m*T)-u(t1-0.2-m*T); %定期延长x1(t)形成一个周期性的信号end
subplot(221) %将该图形放在两行两列的第一个
plot(t1,x) %以t1为横坐标，x为纵坐标画图
axis([-2*T,2*T,-0.2,1.2]) %图形的横坐标从-2T到2T，纵坐标从-0.2到1.2
title('x(t)') %命名x（t）
xlabel('Time t') %轴命名Time t
w0=2*pi/T; %定义频率w0
N=10; %谐波组件的数量为10
L=2*N+1;
for k=-N:N;
ak(N+1+k)=x1*exp(-j*k*w0*t')*dt/T; %计算傅里叶系数的值
end
phi=angle(ak); %计算ak的阶段值
y=0;
for q=1:L;
y=y+ak(q)*exp(j*(q-1-N)*2*pi*t1/T); %分析周期性信号y从有限傅里叶级数(t) end;
subp lot(222) %将图形放在两行两列中的第二位
plot(t1,y) %以t1为横坐标，y为纵坐标画图
axis([-2*T,2*T,-0.2,1.2]) %图形的横坐标从-2T到2T，纵坐标从-0.2到1.2
title('有限项级数合成后的信号') %命名‘有限项级数合成后的信号’
xlabel('Time t') %横轴命名Time t
k=-N:N;
subplot(223)
stem(k,abs(ak))
axis([-N,N,0,1]) %限定横纵坐标范围
title('abs(ak)') %命名
xlabel('Frequency index k') %横轴命名
subplot(224)
stem(k,phi) %作图
axis([-N,N,-pi,pi]) %图形的横坐标从-N到N，纵坐标从-PI到PI title('phi') %命名phi
xlabel('Frequency index k') %横轴命名'Frequency index k'
T = 2; dt = 0.00001; t = -2:dt:2; %定义周期，时间间隔，时间轴长度
x1=(1-t).*(u(t)-u(t-1))+(t+1).*(u(t+1)-u(t)); x = 0; %输入信号x1
for m = -1:1
x = x +(1-t+2*m).*(u(t-2*m)-u(t-1-2*m))+(t+1-2*m).*(u(t+1-2*m)-u(t-2*m));
%使信号周期化
end
w0 = 2*pi/T; %定义w0
N = input('Type in the number of the harmonic components N = :'); %
输入N值
L = 2*N+1; %定义L
for k = -N:1:N; %定义区间
ak(N+1+k) = (1/T)*x1*exp(-j*k*w0*t')*dt; %做傅里
叶变换
end
phi = angle(ak); %计算相位谱
y=0;
for q = 1:L;
y = y+ak(q)*exp(j*(-(L-1)/2+q-1)*2*pi*t/T); %周
期化
end;
subp lot(221) %将图形放在两行两列的第一位
plot(t,x), title('The original signal x(t)'), axis([-2,2,-0.2,1.2]), %以t为横坐标x为纵坐标画图，并命名和确定横纵轴的长度
subplot(223),
plot(t,y), title('The synthesis signal y(t)'), axis([-2,2,-0.2,1.2]),
xlabel('Time t'),
subplot(222)
k=-N:N; stem(k,abs(ak),'k.'), title('The amplitude |ak| of x(t)'), axis([-N,N,-0.1,0.6])
subplot(224)
stem(k,phi,'r.'), title('The phase phi(k) of x(t)'), axis([-N,N,-2,2]), xlabel('Index k')
clear,close all
T=2;dt=0.001;t=-T/2:dt:T/2; %定义周期，时间，时间间隔
w=0.5*pi; %定义频率w
N=10; %输入谐波组件的数量
L = 2*N+1; %
x1=cos((pi/2)*t).*(u(t+1)-u(t-1)); %输入信号x1
subplot(241) %将图形放在两行两列的第一位
plot(t,x1) %以t为横坐标，x1为纵坐标画图
for k=-N:1:N %
ak(N+1+k)=(1/T)*x1*exp(-j*k*w*t')*dt; %计算傅里叶级数
end
phi = angle(ak); %计算ak的阶段值
y=0; %
for q = 1:L; %
y = y+ak(q)*exp(j*(-(L-1)/2+q-1)*2*pi*t/T); %分析周期性信号y从有限傅里叶级数(t)
end
subplot(242) %将图形放在两行两列的第二位
plot(t,y) %以t为横坐标，y为纵坐标画图
subplot(243), %将图形放在两行两列的第三位
k=-N:N; stem(k,abs(ak),'k.') %作图
subplot(244) %将图形放在两行两列的第四位
stem(k,phi,'r.') %作图
t=-2:0.01:2; %定义时间轴长度
w=-4*pi:0.1:4*pi; %定义频率w的长度及间隔
X=x1*exp(-j*t'*w)*T; %x1的傅里叶变换
x1=(t+2).*(u(t+2)-u(t+1))+1.*(u(t+1)-u(t-1))+(-t+2).*(u(t-1)-u(t-2)) %输入信号x1
X2=abs(X); %计算幅度谱
phai=angle(X); %计算相位谱
x3=real(X); %X的实部信号赋给信号3
x4=imag(X); % X的虚部信号赋给信号4
subplot(231) %
plot(w,X2); %
title('幅度频谱') %命名为'幅度频谱'
subplot(232) %
plot(w,phai) %
title('相位谱') %命名为'相位谱'
subplot(233) %
plot(t,x1) %
title('x1波形') %命名为(' x1波形')
axis([-2,2,-2,2]); %限定横纵坐标
subplot(234) %
plot(w,X) %
title('傅里叶变换') %命名为‘傅里叶变换’
subplot(235) %
plot(w,x3) %
title('实部') %命名为‘实部’
subplot(236) %
plot(w,x4) %
title('虚部') %命名为‘虚部’
9、从信号分解的角度，谈谈你对傅里叶级数、傅里叶变换及其物理意义和信号频谱概念的理解。
任意连续测量的时序或信号都可以表示为由无限个不同频率的正弦信号叠加而成。傅立叶变换算法则利用直接测量到的原始信号，以累加方式来计算该信号中不同正弦波信号的频率、振幅和相位。傅立叶变换将原来难以处理的时域信号转换成了易于分析的频域信号，可以利用一些工具对这些频域信号进行处理、加工。最后还可以利用傅立叶反变换将这些频域信号转换成时域信号。傅里叶级数有三角函数形式。对于信号频谱是把一个信号分解为正弦信号的集合,得到其正弦信号幅值随角频率变化的分布,称为该信号的频谱。