信号与系统实验5连续系统频分析

实验五 连续系统的频域分析及连续信号的采样与重构

一、目的

（1）掌握连续系统频率响应概念

（2）掌握利用MATLAB分析系统频率响应的方法

（3）掌握利用MATLAB实现连续信号采用与重构的方法

二、系统的频率响应

设LTI系统的冲激响应为h(t)，该系统的激励信号为f(t)，则此系统的零状态y(t)响应为

y(t)h(t)*f(t) （5-1）

又设f(t)，h(t)，y(t)的傅立叶变换分别为F(j)，H(j)，Y(j)，根据时域卷积定理，与式（5-1）对应的频域关系为：

Y(j)H(j)F(j) （5-2）

一般地，连续系统的频率响应定义为系统的零状态响应y(t)的傅立叶变换Y(j)与激励信号f(t)的傅立叶变换F(j)之比，即

Y(j)H(j) （5-3）

F(j)通常，H(j)是的复函数，因此，又可将其写为：

H(j)H(j)ej() （5-4）

如果令Y(j)Y(j)ejy()，F(j)F(j)ejf()则应有：

H(j)yY(j) （5-5） F(j)f()()() （5-6）

称H(j)为系统的幅频特性，()为系统的相频特性。

需要注意的是，H(j)是系统的固有属性，它与激励信号f(t)的具体形式无关。求系统的H(j)，当然可以按照式（5-3）的定义求，但在实际工程中往往是给出具体的系统图（如具体电路形式），通过电路分析的方法直接求出H(j)。

通常，H(j)可表示成两个有利多项式B(j)与A(j)的商，即

B(j)b1(j)mb2(j)m1bm1(j)bm （5-7） H(j)nn1A(j)a1(j)a2(j)an1(j)an二、利用MATLAB分析系统频响特性

1、分析方法

MATLAB提供了专门用于连续系统频响H(j)分析的函数freqs()。该函数可以求出系统频响的数值解，并可绘出系统的幅频及相频响应曲线。函数freqs()有如下四种调用格式：

（1）h=freqs(b,a,w)

该调用格式中，b为对应于式（5-7）的向量[b1,b2,,bm]，a为对应于式（5-7）的向量[a1,a2,,an]，w为形如w1:p:w2的冒号运算定义的系统频率响应的频率范围，w1为起始频率，w2为终止频率，p为频率取样间隔。向量h则返回在向量w所定义的频率点上系统频响的样值。

..

.

例如，运行如下命令 a=[1 2 1]; b=[0 1];

h=freqs(b,a,0:0.5:2*pi) %计算0~2频率范围内的频响样值 则运行结果为： h =

Columns 1 through 6

1.0000 0.4800 - 0.6400i 0 - 0.5000i -0.1183 - 0.2840i -0.1200 - 0.1600i -0.0999 - 0.0951i Columns 7 through 12

-0.0800 - 0.0600i -0.0641 - 0.0399i -0.0519 - 0.0277i -0.0426 - 0.0199i -0.0355 - 0.0148i -0.0300 - 0.0113i Column 13 -0.0256 - 0.0088i （2）[h,w]=freqs(b,a)

该调用格式将计算默认频率范围内200个频率点的系统频率响应的样值，并赋值给返回变量h，200个频率点记录在w中。 （3）[h,w]=freqs(b,a,n)

该调用格式将计算默认频率范围内n个频率点的系统频率响应的样值，并赋值给返回变量h，n个频率点记录在w中。 （4）freqs(b,a)

该调用格式并不返回系统频率响应样值，而是以对数坐标的方式绘出系统的幅频响应和相频响应。例如运行如下命令：

a=[1 0.4 1]; b=[1 0 0]; freqs(b,a)

运行结果如图5-1所示。

图5-1 对数坐标的系统幅频及相频响应曲线

..

.

下面通过具体例子说明函数freqs()求解系统频响的方法

例5-1：理想低通滤波器在物理上是不可实现的，但传输特性近似于理想特性的电路却能找到。图5-2是常见的用RLC元件构成的二阶低通滤波器（一般说来，阶数越高，实际滤波器的特性越能接近于理想特性）。设LL0.8H，C0.1F，R2，试用MATLAB的freqs()函数求解该系统频率响应并绘图。

CRV0(t)解：根据原理图，容易写出系统的频率响应e(t)为：

1，将L，C，R的H(j)L图5-2 RLC二阶低通滤波器电路图 12LCjR值代入H(j)的表达式，得

1H(j)H(j)ej() 20.08(j)0.4j1其中：

1 H(j)2410.08

实现求解该系统响应的程序为：

b=[0 0 1]; %生成向量b a=[0.08 0.4

1]; %生成向()arctan0.4210.08量a

[h,w]=freqs(b,a,100); %求系统频响特性 h1=abs(h); %求幅频响应 h2=angle(h); %求相频响应 subplot(211); plot(w,h1); grid

xlabel('角频率(W)'); ylabel('幅度');

title('H(jw)的幅频特性'); subplot(212);

plot(w,h2*180/pi); grid

xlabel('角频率(w)'); ylabel('相位(度)');

title('H(jw)的相频特性'); 运行结果如图5-3所示。

由图5-3可见，当从0开始增大时，该低通滤波器幅度从1降到0，c约为3.5；而()从0°降到-180°，与理论分析结果一致。

..

.

图5-3 RLC二阶低通滤波器的幅频特性及相频特性

2、实验内容

图5-4所示的电路为最平坦幅度型二阶低通滤波器。试用MATLAB程序画出系统

U2(j)响应H(j)的幅度响应及相频响应，并与理论分析的结果进行比较。H(j)的

U1(j)截止频率0? 2H 1Ω

U1(t)1ΩU2(t) 2F 图5-4 二阶低通滤波器电路图

计算传递函数：

H(j)U2(j)U1(j)12j112j112j1 22(j)22j2..

.

幅频特性：H(j)182(222)21444

相频特性：(j)-222arctanarctan2-2ω221-222arctanarctan2-2ω221-2202-2ω2-22022-2ω作图：

w = 0:0.01:10;

h = 1./sqrt(4*w.^4+4); subplot(2,1,1); plot(w, h); hold on;

axis([0,10,0,0.8]);

title('H(jω)的幅频特性'); w = 0:0.01:0.99;

ang = atan(sqrt(2).*w./(w.^2-1))*180/pi; subplot(2,1,2); plot(w, ang); hold on;

axis([0,10,-200,0]);

title('H(jω)的相频特性'); hold on;

w = 1.01:0.01:10;

ang = atan(sqrt(2).*w./(w.^2-1))*180/pi-180; subplot(2,1,2); plot(w, ang); hold off;

Matlab实现：

b=[0 0 1]; %生成向量b a=[2 2*sqrt(2) 2]; %生成向量a [h,w]=freqs(b,a,100); %求系统频响特性 h1=abs(h); %求幅频响应 h2=angle(h); %求相频响应

..

.

subplot(211); plot(w,h1); grid

xlabel('角频率(W)'); ylabel('幅度');

title('H(jw)的幅频特性'); subplot(212);

plot(w,h2*180/pi); grid

xlabel('角频率(w)'); ylabel('相位(度)');

title('H(jw)的相频特性');

三、连续信号的采样与重构

1、信号采样

图5-5给出信号采样原理图 fs(t)f(t) 相乘 Ts(t)

图5-5 信号采样原理图

由图5-5可见，fs(t)f(t)Ts(t)，其中，冲激采样信号Ts(t)的表达式为：

(t)(tnT) （5-8）

Tsns2。设F(j)，Fs(j)分别为f(t)，fs(t)nTs的傅立叶变换，由傅立叶变换的频域卷积定理，可得

11Fs(j)F(j)*s(ns)F[j(ns)] （5-9）

n2Tsn其傅立叶变换为s(ns)，其中s..

.

若设f(t)是带限信号，带宽为m，由式（5-9）可见，f(t)经过采样后的频谱Fs(j)就是将F(j)在频率轴上搬移至0,s,2s,,ns,处（幅度为原频谱的1Ts倍）。因此，当s2m时，频谱不发生混叠；而当s2m时，频谱发生混叠。

应该指出的是，实际信号中，绝大多数都不是严格意义上的带限信号，这时根据实际精度要求来确定信号的带宽m。

例5-2：当采样频率s2m时，称为临界采样，取cm。下列程序实现对信号Sa(t)的采样及由采样信号恢复Sa(t)（见信号恢复小节）。

wm=1; wc=wm; Ts=pi/wm; ws=2*pi/Ts; n=-100:100; nTs=n*Ts; f=sinc(nTs/pi);

Dt=0.005;t=-15:Dt:15;

fa=f*Ts*wc/pi*sinc((wc/pi)*(ones(length(nTs),1)*t-nTs'*ones(1,length(t)))); t1=-15:0.5:15; f1=sinc(t1/pi); subplot(211); stem(t1,f1); xlabel('kTs'); ylabel('f(kTs)');

title('sa(t)=sinc(t/pi)的临界采样信号'); subplot(212); plot(t,fa) xlabel('t'); ylabel('fa(t)');

title('由sa(t)=sinc(t/pi)的临界采样信号重构sa(t)'); grid;

程序运行结果如图5-6所示。 [注意]Sa(t)=sinc(t/pi)

..

.

图5-6 临界采样信号及信号恢复

2、信号重构

设信号f(t)被采样后形成的采样信号为fs(t)，信号的重构是指由fs(t)经过内插处理后，恢复出原来信号f(t)的过程。又称为信号恢复。

若设f(t)是带限信号，带宽为m，经采样后的频谱为Fs(j)。设采样频率s2m，则由式（5-9）知Fs(j)是以s为周期的谱线。现选取一个频率特性

cTsH(j)（其中截止频率c满足mcs）的理想低通滤波器与

2c0Fs(j)相乘，得到的频谱即为原信号的频谱F(j)。

显然，F(j)Fs(j)H(j)，与之对应的时域表达式为

f(t)h(t)*fs(t) （5-10）

而

fs(t)f(t)(tnTs)f(nTs)(tnTs)

nnh(t)F1[H(j)]TsSa(t) cc将h(t)及fs(t)代入式（5-10）得

..

.

TSa(t)f(nT)Sa[(tnT)] （5-11） 式（5-11）即为用f(nTs)求解f(t)的表达式，是利用MATLAB实现信号重构的基本关系式，抽样函数Sa(ct)在此起着内插函数的作用。

f(t)fs(t)*Tscsccnscssint，其F(j)为： t1 F(j)01即f(t)的带宽为m1，为了由f(t)的采样信号fs(t)不失真地重构f(t)，由时域采样定

例5-3：设f(t)Sa(t)理知采样间隔TsSinc(t)。利用MATLAB的抽样函数，取T0.7（过采样）

smn为了比较由采样信号恢复后的信号与原信号的误差，计算两信号的绝对误差。MATLAB实现此过程的程序如下：

wm=1;

wc=1.1*wm; Ts=0.7*pi/wm; ws=2*pi/Ts; n=-100:100; nTs=n*Ts; f=sinc(nTs/pi);

Dt=0.005;t=-15:Dt:15;

fa=f*Ts*wc/pi*sinc((wc/pi)*(ones(length(nTs),1)*t-nTs'*ones(1,length(t)))); error=abs(fa-sinc(t/pi)); t1=-15:0.5:15; f1=sinc(t1/pi); subplot(311); stem(t1,f1); ylabel('f(kTs)');

title('sa(t)=sinc(t/pi)的采样信号'); subplot(312); plot(t,fa) ylabel('fa(t)');

title('由sa(t)=sinc(t/pi)的过采样信号重构sa(t)'); grid;

subplot(313); plot(t,error); ylabel('error(t)');

title('过采样信号与原信号的误差error(t)');

结果如图5-7所示，由图5-7可知，两信号的绝对误差error已在10-6数量级，说明重构信号的精度已经很高。

..

sin(t)来表示Sa(t)，有Sa(t)Sinc(t/)。据此可知： tcTsccf(t)fs(t)*TsSa(ct)f(nT)Sinc[(tnTs)] （5-12） s.

图5-7 过采样信号、重构信号及两信号的绝对误差

将上述程序稍加改动，令m1，cm，Ts1.5（欠采样），其余不变。改动后

程序运行结果如图5-8所示。

由图5-8可见，绝对误差error已大为增加，其原因是因采样信号的频谱混叠，使得在c区域内的频谱相互“干扰”所致。

图5-8 欠采样信号、重构信号及两信号的绝对误差

..

.

3、实验内容

设f(t)e1000t，由于不是严格的带限信号，但其带宽m可根据一定的精度要求做一近似。试根据以下三种情况用MATLAB实现由f(t)采样信号fs(t)重构f(t)并求出两者误差，分析三种情况下的结果。

（1）m5000，cm，Ts/m；

Matlab实现： wm=5000*pi; wc=1*wm; Ts=pi/wm; ws=2*pi/Ts; n=-100:100; nTs=n*Ts;

f=exp(-1000*abs(nTs));

Dt=0.005/10000;t=-15/10000:Dt:15/10000;

fa=f*Ts*wc/pi*sinc((wc/pi)*(ones(length(nTs),1)*t-nTs'*ones(1,length(t)))); error=abs(fa-exp(-1000*abs(t))); t1=-15/10000:0.5/10000:15/10000; f1=exp(-1000*abs(t1)); subplot(311); stem(t1,f1); ylabel('f(kTs)');

title('exp(-1000*abs(t))的采样信号'); subplot(312); plot(t,fa) ylabel('fa(t)');

title('由exp(-1000*abs(t))的采样信号重构'); grid;

subplot(313); plot(t,error); ylabel('error(t)');

title('过采样信号与原信号的误差error(t)'); 作图：

（2）m10000，c1.1m，Ts/m；

Matlab实现： wm=10000*pi;

..

.

wc=1.1*wm; Ts=pi/wm; ws=2*pi/Ts; n=-100:100; nTs=n*Ts;

f=exp(-1000*abs(nTs));

Dt=0.005/10000;t=-15/10000:Dt:15/10000;

fa=f*Ts*wc/pi*sinc((wc/pi)*(ones(length(nTs),1)*t-nTs'*ones(1,length(t)))); error=abs(fa-exp(-1000*abs(t))); t1=-15/10000:0.5/10000:15/10000; f1=exp(-1000*abs(t1)); subplot(311); stem(t1,f1); ylabel('f(kTs)');

title('exp(-1000*abs(t))的采样信号'); subplot(312); plot(t,fa) ylabel('fa(t)');

title('由exp(-1000*abs(t))的采样信号重构'); grid;

subplot(313); plot(t,error); ylabel('error(t)');

title('过采样信号与原信号的误差error(t)') 作图：

（3）m2500，c0.9m，Ts/m；

Matlab实现： wm=2500*pi; wc=0.9*wm; Ts=pi/wm; ws=2*pi/Ts; n=-100:100; nTs=n*Ts;

f=exp(-1000*abs(nTs));

Dt=0.005/10000;t=-15/10000:Dt:15/10000;

fa=f*Ts*wc/pi*sinc((wc/pi)*(ones(length(nTs),1)*t-nTs'*ones(1,length(t)))); error=abs(fa-exp(-1000*abs(t))); t1=-15/10000:0.5/10000:15/10000;

..

.

f1=exp(-1000*abs(t1)); subplot(311); stem(t1,f1); ylabel('f(kTs)');

title('exp(-1000*abs(t))的采样信号'); subplot(312); plot(t,fa) ylabel('fa(t)');

title('由exp(-1000*abs(t))的采样信号重构'); grid;

subplot(313); plot(t,error); ylabel('error(t)');

title('过采样信号与原信号的误差error(t)') 作图：

..