2014信号与系统实验五(2)

Matlab拓展实验1：周期方波的傅立叶级数分析
一、周期方波信号的产生
1、周期方波函数：square(t,duty)
（1）其中t为时间向量。

（2）duty为占空比，此参数可缺省，缺省时，duty=50。

注：周期方波的占空比为x%，是指一个周期内，高电平持续时间占整个周期时长的x%。

但在square函数中，参数duty是一个0~100的数值。

例如，若需产生占空比为20%的方波，则应设置duty=20。

2、功能描述：square函数产生一个周期为2π、高低电平分别为±1的周期方波，通过适当的编程可以调整为任意周期，任意幅度，任意中心值。

【例1】在时间范围[−4π,4π]产生周期的2π的周期方波。

clear;clf;clc;
step=pi/200; %设置时间步长
t=-4*pi:step:4*pi; %设置时间范围
x1=square(t); %产生占空比为50%的周期方波
duty=20;x2=square(t,duty); %产生占空比为20%的周期方波
subplot(2,1,1);plot(t,x1);
axis([-4*pi,4*pi,-1.5,1.5]);title('占空比为0.5的周期方波');
subplot(2,1,2);plot(t,x2);
axis([-4*pi,4*pi,-1.5,1.5]);title('占空比为0.2的周期方波');
上述程序的运行结果如下：
【思考】如将例1程序中的作图部分修改如下（修改部分以红色标出），这种处理方式称之为横坐标对pi归一化，图形结果会有何不同？这时横坐标的数值标注含义有何变化？
subplot(2,1,1);plot(t/pi,x1);
axis([-4,4,-1.5,1.5]);title('占空比为0.5的周期方波');
subplot(2,1,2);plot(t/pi,x2);
axis([-4,4,-1.5,1.5]);title('占空比为0.2的周期方波');
【实践1】产生周期为1，占空比为30%，高电平为1，低电平为0的周期方波，时间范围为[−2,2]，即能观察到4个完整的周期。

（时间变量的尺度变换）
二、傅立叶级数分解的计算
1、周期信号傅立叶级数展开的理论分析
周期信号的傅立叶级数有若干种形式（如三角形式、指数形式），这里以指数形式的傅立叶级数为例进行说明。

信号x(t)满足，x t=x(t+T)，其中T称为周期（或基波周期，单位s），F=1T称为频率（或基波频率，单位Hz），ω=2πT=2πF称为基波角频率，单位rad s。

（1）指数形式的傅立叶级数：
x t = a k e jkωt ∞
k =−∞
=⋯+a −1e −jωt +a 0+a 1e jωt +a 2e j2ωt +⋯
（2）a k （k =0,±1,±2,⋯）称为傅立叶系数，可通过以下方式计算，
a k =1T x (t )e −jkωt dt T
（3）傅里叶级数展开的物理意义：
周期信号可以表示为直流分量（a 0）与各次谐波（a ±k e ±jkωt ）之和。

其中a ±1e ±jωt 称为基波（或1次谐波），a ±2e ±j2ωt 称为2次谐波，以此类推。

注：虽然表面上看，展开式中的频率有正负之分，但实际上 e jkωt 与 e −jkωt 的频率是完全相同的，这二者合在一起称为k 次谐波。

（4）当x (t )为实信号时，a k 满足共轭对称性，即有a k =a −k ∗。

由于k 次谐波包含了展开式中a ±k e ±jkωt 这两部分，因此对于实信号，k 次谐波的幅度应为2|a k |。

2、周期信号傅立叶级数展开的Matlab 计算
（1）积分计算的Matlab 实现：傅立叶系数a k 的计算必须通过积分运算来完成，Matlab 中关于积分计算的实现有若干种途径，这里推荐一种较为方便的计算方法，即trapz 函数，先来看一个简单的例子。

【例2】计算 sin (t )dt π0
，该积分的正确结果为2，Matlab 可以通过以下代码完成积分运算。

step=pi/200;
%设定时间步长，步长越小，积分的结果越精确 t=0:step:pi;
%设定时间向量，即按步长将积分区间做等间隔划分 x=sin(t); %计算被积信号在时间采样点的数值
result=trapz(t,x);
%计算积分结果，注意时间t 与信号x 的位置前后关系
（2）周期信号的频谱显示：周期信号的频谱是由离散的谱线构成的，谱线的位置（即每根谱线的横坐标）位于基波频率的整数倍上，谱线的高度即为该次谐波的幅度，在指数形式的FS展开式中，为傅立叶系数a k的模的2倍，即2|a k|。

（a）占空比为50%的周期方波（b）周期方波的频谱
【实践2】产生一个周期方波，并计算其傅立叶系数，画出频谱分析图。

此任务可以分成三段程序来完成，分别是：①产生时域信号；②计算傅立叶系数；③绘制频谱分析图。

（1）产生时域信号：产生一个方波信号，周期为1，高电平为1，低电平为0，占空比为50%，时间范围为[0,1]（只含一个周期）。

注：若是绘制信号的时域波形，我们一般会选择一个较长的时间范围，其中包含若干个周期，以便于得出该信号波形的“全貌”。

但如果是用于计算FS系数，则只需要一个周期，这一点非常重要。

参数设置代码如下：
T=1;w=2*pi/T; %设置信号的周期T=1，w为基波角频率
step=T/500; %设置时间步长
t=0:step:T; %设置时间向量
注：这里只给出了周期、角频率、时间步长以及时间向量的设置，方波信号的产生需自行设计。

（2）利用trapz函数计算傅里叶系数a k，其计算公式如下，
a k=1
x(t)e−j kωt dt
T
计算a k（k=1,2,3,⋯,10）的代码如下，
for k=1:10
a(k)=trapz(t,x.*exp(-j*k*w*t))/T;
a(k)=abs(a(k));
end
注：①上述代码计算了k = 1~10次谐波的傅立叶系数。

②方波信号保存在变量x中，x应为一行向量，其长度与时间向量t一致。

③exp为指数函数，j为虚单位，w为基波角频率。

④信号x与指数函数exp产生的结果均为行向量，二者之间是.*运算。

⑤各次谐波的FS系数计算结果，存放于变量a中，a(k)为变量a的第k个值，Matlab规定k必须是从1开始的整数，且必须使用圆括号( )。

⑥由于a(k)的结果一般为复数值，而谐波的幅度是要取模的，其函数为abs()。

⑦循环中没有计算a0，大家自己思考如何计算a0（直流分量），并作相应编程。

提示：编程时，可设置一个变量a0用于保存系数a0的计算结果，但千万不能写成a(0)，为什么？（这里的a0与a是两个完全无关的变量）。

（3）利用stem函数绘制周期方波的频谱分析结果
①stem(n,x) 函数可以绘制离散谱线的图形，其中n为横坐标向量，x为纵坐标向量，这里的横坐标应为频率轴（或角频率）。

②各次谐波的谱线幅度包含于变量a中（这里a是一个长度为10的向量），而各次谐波对应的横坐标为基波角频率的整数倍，思考应如何编程实现。

三、任务与思考
1、完成【实践1】与【实践2】的编程，并书写一个简短的实验报告。

2、对于【实践2】中的周期方波，计算各次谐波的理论幅度值，即求出a k的表达式，进而得出各次谐波相对于基波的幅度比值，对比Matlab的计算结果，看看是否一致？
3、若改变占空比（例如将占空比设置为25%），频谱结果会有何明显的变化？
4、报告内容包括：
①程序代码与运行结果；
②Matlab计算结果与理论结果的分析对比；
③改变占空比引起的频谱变化。