数字信号处理备用总程序.docx

数字信号处理备用总程序.docx

《数字信号处理备用总程序.docx》由会员分享，可在线阅读，更多相关《数字信号处理备用总程序.docx（44页珍藏版）》请在冰豆网上搜索。

数字信号处理备用总程序

实验一

1.在MATLAB中实现

序列，显示范围

（函数命名为impseq（n0,n1,n2））。

并利用该函数实现序列：

；

程序如下：

%函数impseq（n0,n1,n2）

functiony=impseq（n0,n1,n2）

n=n1:

1:

n2%横坐标

y=[（n-n0）==0]%生成离散信号y（n）

%脚本文件

%调用以上函数

n=-3:

10;

y=2*impseq（3,-3,10）+impseq（6,-3,10）;

stem（n,y,'k','filled'）

title（'单位脉冲序列'）

xlabel（'时间（n）'）;

ylabel（'幅度y（n）'）;

仿真结果：

2.在MATLAB中实现

序列，显示范围

（函数命名为stepseq（n0,n1,n2））。

并利用该函数实现序列：

程序如下：

%函数stepseq（n0,n1,n2）

functiony=stepseq（n0,n1,n2）

n=n1:

1:

n2;

y=[（n-n0）>=0];

%脚本文件

%调用以上函数

n=-5:

20;

y=stepseq（2,-5,20）+stepseq（-2,-5,20）;

stem（n,y,'k',’filled’）

title（'单位阶跃序列'）

xlabel（'时间（n）'）;

ylabel（'幅度y（n）'）;

仿真结果：

3.在MATLAB中利用数组运算符“.^”来实现一个实指数序列。

如：

程序如下：

%函数zhishu（a,n1,n2）

functiony=zhishu（a,n1,n2）

n=n1:

n2;

y=（a）.^n;

%脚本文件

%调用以上函数

y=zhishu（0.3,0,50）;

n=0:

50;

stem（n,y,'k','filled'）

title（'实指数序列'）

xlabel（'时间（n）'）;

ylabel（'幅度y（n）'）;

仿真结果：

4.在MATLAB中用函数sin或cos产生正余弦序列，如：

程序如下：

%脚本文件

n=0:

20;

x=11*sin（0.3*pi*n+0.2*pi）+5*cos（0.3*pi*n）;

stem（n,x,'k','filled'）

title（'正余弦序列'）

xlabel（'时间（n）'）;

ylabel（'幅度y（n）'）;

仿真结果：

5.已知

，试显示

在

区间的波形。

程序如下：

%函数yuxian（a,n1,n2）

functiony=yuxian（a,n1,n2）

n=n1:

n2;

y=3*cos（2*pi/10*（n-a））;

%脚本文件

%调用以上函数

n=0:

20;

y1=yuxian（0,0,20）

y2=yuxian（3,0,20）

y3=yuxian（-3,0,20）

subplot（3,1,1）

stem（n,y1,'r'）

title（'x（n）'）

subplot（3,1,2）

stem（n,y2,'g'）

title（'x（n-3）'）

subplot（3,1,3）

stem（n,y3,'y'）

title（'x（n+3）'）

仿真结果：

6.参加运算的两个序列维数不同，已知

，

，求

。

程序如下：

%函数u（n0,n1,n2）

functiony=u（n0,n1,n2）

n=n1:

1:

n2;

y=[（n-n0）>=0];

%脚本文件

%调用以上函数

n1=-4:

6;

n2=-5:

8;

x1=u（-2,-4,6）;

x2=u（4,-5,8）;

y1=[0x100];

y2=x2;

y=y1+y2;

stem（n2,y,'k','filled'）

title（'x（n）=x1（n）+x2（n）'）

xlabel（'时间（n）'）;

ylabel（'幅度y（n）'）;

仿真结果：

实验2

1.利用interp1函数重构采样信号

已知一个连续时间信号

，

，取最高有限带宽频率

，对x（t）进行等间隔采样，采样频率为fs，要求：

（1）画出原连续时间信号x（t）的波形；

程序如下

f0=1;%设置正弦波信号的频率

T0=1/f0;

t=0:

0.01:

3*T0;

xt=sin（2*pi*f0*t）+1/3*sin（6*pi*f0*t）;%产生正弦波信号

fs=15*f0;%设置采样频率

Ts=1/fs;

n=0:

1:

3*T0/Ts;%采样点数

xn=1/3*sin（2*pi*f0*n*Ts）;%产生采样信号

subplot（2,1,1）

plot（t,xt）;%绘制原连续信号

title（'正弦信号的采样'）

subplot（2,1,2）

stem（n,xn）;%绘制采样信号

仿真波形：

（2）设

三种情况下对连续信号分别进行采样，画出采样信号波形，并利用内插公式重建原信号，对结果进行分析。

程序1：

f0=1;%设置正弦波信号的频率

T0=1/f0;

t=0:

0.01:

3*T0;

xt=sin（2*pi*f0*t）+1/3*sin（6*pi*f0*t）;%产生正弦波信号

subplot（3,1,1）

plot（t,xt）;%绘制原连续信号

fs=15*f0;%设置采样频率

Ts=1/fs;

n=0:

1:

3*T0/Ts;%采样点数

xn=sin（2*pi*f0*n*Ts）+1/3*sin（6*pi*f0*n*Ts）;%产生采样信号

tn=0:

Ts:

3*T0;%采样时间

yt=interp1（tn,xn,t,'spline'）;%内插恢复信号

subplot（3,1,2）

stem（n,xn）;%绘制采样信号

title（'采样信号的恢复fs=3fh'）

subplot（3,1,3）

plot（t,yt）;%绘制内插恢复信号

仿真波形1：

程序2：

f0=1;%设置正弦波信号的频率

T0=1/f0;

t=0:

0.01:

3*T0;

xt=sin（2*pi*f0*t）+1/3*sin（6*pi*f0*t）;%产生正弦波信号

subplot（3,1,1）

plot（t,xt）;%绘制原连续信号

fs=10*f0;%设置采样频率

Ts=1/fs;

n=0:

1:

3*T0/Ts;%采样点数

xn=sin（2*pi*f0*n*Ts）+1/3*sin（6*pi*f0*n*Ts）;%产生采样信号

tn=0:

Ts:

3*T0;%采样时间

yt=interp1（tn,xn,t,'spline'）;%内插恢复信号

subplot（3,1,2）

stem（n,xn）;%绘制采样信号

title（'采样信号的恢复fs=3fh'）

subplot（3,1,3）

plot（t,yt）;%绘制内插恢复信号

仿真波形2：

程序3：

f0=1;T0=1/f0;

t=0:

0.01:

3*T0;

xt=sin（2*pi*f0*t）+1/3*sin（6*pi*f0*t）;subplot（3,1,1）

plot（t,xt）

title（'原正弦信号'）

fs=5*f0;

Ts=1/fs;

n=0:

1:

3*T0/Ts;

xn=sin（2*pi*f0*n*Ts）+1/3*sin（6*pi*f0*n*Ts）;

yt=interp1（tn,xn,t,'spline'）;subplot（3,1,2）

stem（n,xn）;

title（'采样信号的恢复fs=fh'）

subplot（3,1,3）

plot（t,yt）;

仿真波形3：

分析：

当fs>=2fh进行采样时，在进行恢复时可以原样恢复；当fs<2fh时，采样恢复时不能够恢复原图像。

2.利用conv函数计算两个有限长序列的卷积

已知序列x（n）和h（n），试求其卷积结果y（n）。

程序：

xn=[5,9,3,6,-8];

hn=[18,7,5,20,11,14,9];

n=-4:

6;

yn=conv（hn,xn）

stem（n,yn）

title（'卷积运算'）

仿真波形：

3.利用impz和dstep子函数求解离散系统的单位脉冲响应和单位阶跃响应

已知线性常系数差分方程为

，求系统的单位脉冲响应和单位阶跃响应，并画出图形（

）。

程序：

pulse=[1,zeros（1,50）];

b=[1];

a=[1,-0.6,0.8];

h1=impz（b,a,50）;

subplot（2,1,1）,stem（h1）,title（'impzfunction'）;

n=0:

50;

pulse2=n>=0;

h2=dstep（b,a,50）;

subplot（2,1,2）,stem（h2）,title（'dstepfunction'）;

仿真波形：

4.在MATLAB中利用filtic和filter子函数求解离散系统的单位脉冲响应

已知线性常系数差分方程为

（1）若

，输入为

，计算系统的响应，并画出图形（

）；

程序：

pulse=[1,zeros（1,25）];

b=[1];

a=[1,-0.9,0.5];

n=0:

25;

h1=filter（b,a,pulse）;

stem（n,h1）

title（'y（n）=0;n<0'）;

仿真波形：

（2）若

，重做

（1）；

程序：

pulse=[1,zeros（1,25）];

b=[1];

a=[1,-0.9,0.5];

n=0:

25;

h1=filter（b,a,pulse）;

subplot（2,1,1）

stem（n,h1）

title（'y（n）=0;n<0'）;

Y0=[1]

zi0=filtic（b,a,Y0）

h2=filter（b,a,pulse,zi0）

subplot（2,1,2）

stem（n,h2）

title（'y（-1）=1'）;

仿真波形：

（3）若

，

，重做

（1）；

程序：

pulse=[1,zeros（1,25）];

b=[1];

a=[1,-0.9,0.5];

n=0:

25;

h1=filter（b,a,pulse）;

subplot（2,1,1）

stem（n,h1）

title（'y（n）=0;n<0'）;

Y1=[2,1]

zi=filtic（b,a,Y1）

h2=filter（b,a,pulse,zi）

subplot（2,1,2）

stem（n,h2）

title（'y（-2）=1,y（-1）=2'）;

仿真波形：

（4）分析该系统的稳定性。

答：

对于任意有界的输入其零状态响应是有界的，在本实验中，n趋于无穷时是收敛的，所以该系统是稳定的。

实验3

1.已知序列

，画出由离散时间傅里叶变换（DTFT）求得的幅度谱

和相位谱

图形。

程序：

x=[0,1,2,3,4,5,6,7];%定义序列x

N=length（x）;%求x的列长

n=0:

N-1;

w=linspace（-2*pi,2*pi,500）;%给出频率ω及其范围

X=x*exp（-j*n'*w）;%求DTFT

subplot（2,1,1）

plot（w/pi,abs（X））;title（'幅度谱'）;axis（[-2,2,min（abs（X））,max（abs（X））]）;%画幅度谱

subplot（2,1,2）

plot（w/pi,angle（X））;title（'相位谱'）;axis（[-2,2,min（angle（X））,max（angle（X））]）;%画相位谱

仿真结果：

2.已知周期序列的主值

，求

周期重复次数为6次时的DFS。

要求画出原主值序列、周期序列及其傅里叶级数变换对应的

图形。

程序：

x=[0,1,2,3,4,5,6,7];

n1=0:

7;

N=length（x）;n=0:

6*N-1;k=0:

6*N-1;

x1=x（mod（n,N）+1）;%序列周期重复6次

Xk=x1*exp（-j*2*pi/N）.^（n'*k）;%求DFS

subplot（2,2,1）,stem（n1,x）;title（'原主值序列'）;

subplot（2,2,2）,stem（n,x1）;title（'主值序列周期重复6次'）;axis（[0,6*N,min（x1）,max（x1）]）;

subplot（2,2,3）,stem（k,abs（Xk））;

title（'序列周期重复6次的DFS幅频特性'）;axis（[0,6*N,min（abs（Xk））,max（abs（Xk））]）;

subplot（2,2,4）,stem（k,angle（Xk））;

title（'序列周期重复6次的DFS相频特性'）;axis（[0,6*N,min（angle（Xk））,max（angle（Xk））]）;

仿真结果：

3.已知

，求

的DFT和

的IDFT。

要求：

画出序列离散傅里叶变换对应的

图形、原信号

与

的傅里叶逆变换

的图形，并进行比较。

（1）程序：

x=[01234567];

N=length（x）;n=0:

N-1;k=0:

N-1;

Xk=x*exp（-j*2*pi/N）.^（n'*k）;%求DFT

x1=（Xk*exp（j*2*pi/N）.^（n'*k））/N;%求IDFT

subplot（2,2,1）,stem（k,abs（Xk））;title（'|X（k）|'）;

subplot（2,2,2）,stem（k,angle（Xk））;title（'arg|X（k）|'）;

subplot（2,2,3）,stem（n,x）;title（'x（n）'）;

subplot（2,2,4）,stem（n,x1）;title（'IDFT[X（k）]'）;

仿真结果：

比较分析：

经IDFT进行反变换所得的x（n）与原序列x（n）一致

4.利用FFT计算序列的频谱

已知序列

，

，

（1）求序列

、

、

和

的离散傅里叶变换

值；

程序：

x=[8,4,2,1];

N1=length（x）;n1=0:

N1-1;k1=0:

N1-1;

Xk=x*exp（-j*2*pi/N1）.^（n1'*k）%求X（k）

g=[8,4,2,1,0,0,0,0];

N2=length（g）;n2=0:

N2-1;k2=0:

N2-1;

Gk=g*exp（-j*2*pi/N2）.^（n2'*k2）%求G（k）

y=[8,0,4,0,2,0,1,0];

N3=length（y）;n3=0:

N3-1;k3=0:

N3-1;

Yk=y*exp（-j*2*pi/N3）.^（n3'*k3）%求Y（k）

h=[8,4,2,1,8,4,2,1];

N4=length（h）;n4=0:

N4-1;k4=0:

N4-1;

Hk=h*exp（-j*2*pi/N4）.^（n4'*k4）%求H（k）

subplot（2,2,1）,stem（k,Xk）;title（'X（k）'）;

subplot（2,2,2）,stem（k,Gk）;title（'G（k）'）;

subplot（2,2,3）,stem（k,Yk）;title（'Y（k）'）;

subplot（2,2,4）,stem（k,Hk）;title（'H（k）'）;

仿真结果：

值如下：

Xk=

Columns1through4

15.00006.0000-3.0000i5.0000-0.0000i6.0000+3.0000i

Columns5through8

15.0000+0.0000i6.0000-3.0000i5.0000-0.0000i6.0000+3.0000i

Gk=

Columns1through4

15.000010.1213-5.5355i6.0000-3.0000i5.8787-1.5355i

Columns5through8

5.0000-0.0000i5.8787+1.5355i6.0000+3.0000i10.1213+5.5355i

Yk=

Columns1through4

15.00006.0000-3.0000i5.0000-0.0000i6.0000+3.0000i

Columns5through8

15.0000+0.0000i6.0000-3.0000i5.0000-0.0000i6.0000+3.0000i

Hk=

Columns1through4

30.0000-0.0000-0.0000i12.0000-6.0000i-0.0000-0.0000i

Columns5through8

10.0000+0.0000i0.0000-0.0000i12.0000+6.0000i0.0000-0.0000i

>>

（2）画出

、

、

以及

图形，并进行比较。

程序：

x=[8,4,2,1];

N1=length（x）;n1=0:

N1-1;k1=0:

N1-1;

Xk=x*exp（-j*2*pi/N1）.^（n1'*k）;%求X（k）

subplot（4,2,1）,stem（k,abs（Xk））;title（'|X（k）|'）;

subplot（4,2,2）,stem（k,angle（Xk））;title（'arg|X（k）|'）;

g=[8,4,2,1,0,0,0,0];

N2=length（g）;n2=0:

N2-1;k2=0:

N2-1;

Gk=g*exp（-j*2*pi/N2）.^（n2'*k2）;%求G（k）

subplot（4,2,3）,stem（k2,abs（Gk））;title（'|G（k）|'）;

subplot（4,2,4）,stem（k2,angle（Gk））;title（'arg|G（k）|'）;

y=[8,0,4,0,2,0,1,0];

N3=length（y）;n3=0:

N3-1;k3=0:

N3-1;

Yk=y*exp（-j*2*pi/N3）.^（n3'*k3）;%求Y（k）

subplot（4,2,5）,stem（k3,abs（Yk））;title（'|Y（k）|'）;

subplot（4,2,6）,stem（k3,angle（Yk））;title（'arg|Y（k）|'）;

h=[8,4,2,1,8,4,2,1];

N4=length（h）;n4=0:

N4-1;k4=0:

N4-1;

Hk=h*exp（-j*2*pi/N4）.^（n4'*k4）;%求H（k）

subplot（4,2,7）,stem（k4,abs（Hk））;title（'|H（k）|'）;

subplot（4,2,8）,stem（k4,angle（Hk））;title（'arg|H（k）|'）;

仿真结果：

比较分析：

序列的DFT是序列频谱的等间隔采样。

G（k）与X（k）的频谱是相对应的，G（k）比X（k）频谱间隔小，谱线密；Y（k）是X（k）的重复，y（n）与x（n）相比，相当于改变了对信号的取样频率，H（k）与X（k）相比较，改变了谱间隔及谱线的比例。

实验4

4.1已知离散时间系统的系统函数为

求系统的零极点，画出零极点分布图，分析系统的因果稳定性。

程序如下：

B=[0.20.10.30.10.2];

A=[1-1.11.5-0.70.3];

r1=roots（B）%求分子多项式的根，即系统的零点

r2=roots（A）%求分母多项式的根，即系统的极点

figure

（1）

zplane（B,A）;%调用zplane函数画零极点图

仿真结果：

系统的零极点：

零点：

r1=

-0.5000+0.8660i

-0.5000-0.8660i

0.2500+0.9682i

0.2500-0.9682i

极点：

r2=

0.2367+0.8915i

0.2367-0.8915i

0.3133+0.5045i

0.3133-0.5045i

系统的因果稳定性分析：

因果稳定系统的充要条件是系统函数的极点都位于Z平面单位圆内部，不包括单位圆。

由系统零极图可知该系统的全部极点都在单位圆内，所以该系统因果稳定。

4.2已知离散时间系统的系统函数为

画出系统在

频率范围内的幅频响应

和相频响应

图形。

程序如下：

B=[0.20.10.30.10.2];

A=[1-1.11.5-0.70.3];

[H,w]=freqz（B,A）;%求离散系统频响特性

Hf=abs（H）;%取幅度

Hx=angle（H）;%取相角

figure

（1）

subplot（2,1,1）

plot（w,Hf）;%画幅度谱

title（'幅频特性曲线'）;

subplot（2,1,2）

plot（w,Hx）;%画相位谱

title（'相频特性曲线'）;

仿真结果：

4.3已知系统的单位脉冲响应

，利用freqz函数画出系统在

频率范围内的幅频响应

和相频响应

图形。

程序如下：

B=[8,4,2,1];

A=[1,0,0,0];

N=1024;

[H,w]=freqz（B,A,N,'whole'）;%求hn的离散时间傅里叶变换Hk

subplot（2,1,1）,plot（w/pi,abs（H））;

title（'离散时间傅里叶变换后的幅度谱'）;

subplot（2,1,2）,plot（w/pi,angle（H））;

title（'离散时间傅里叶变换后的相位谱'）;

仿真结果：

4.4已知离散时间系统的系统函数为

画出系统在

频率范围内的幅频响应和相频响应图形，并与课本上例2-9进行对比。

程序如下：

B=conv（[1,-exp（pi/2*1i）],[1,-exp（-pi/2*1i）]）;

A=conv（[1,-0.8*exp（pi/4*1i）],[1,-0.8*exp（-pi/4*1i）]）;

[H,w]=freqz（B,A）;

subplot（2,1,1）

plot（w,abs（H））;

title（'幅频特性'）;

subplot（2,1,2）

plot（w,angle（H））

title（'相频特性'）;

仿真结果：

4.5已知系统单位脉冲响应为

，如果输入为

，利用圆周卷积定理求系统输出

。

（1）用FFT实现：

n1=0:

19;

hn=cos（0.5*n