数字信号处理MATLAB仿真.docx
《数字信号处理MATLAB仿真.docx》由会员分享，可在线阅读，更多相关《数字信号处理MATLAB仿真.docx（14页珍藏版）》请在冰豆网上搜索。
数字信号处理MATLAB仿真
实验一数字信号处理的Matlab仿真
一、实验目的
1、掌握连续信号及其MATLAB实现方法；
2、掌握离散信号及其MATLAB实现方法
3、掌握离散信号的基本运算方法，以及MATLAB实现
4、了解离散傅里叶变换的MATLAB实现
5、了解IIR数字滤波器设计
6、了解FIR数字滤波器设计1
二、实验设备
计算机，Matlab软件
三、实验内容
（一）、连续信号及其MATLAB实现
1、单位冲击信号
例1.1：
t=1/A=50时，单位脉冲序列的MATLAB实现程序如下：
clearall;
t1=-0.5:
0.001:
0;
A=50;
A1=1/A;
n1=length（t1）;
u1=zeros（1,n1）;
t2=0:
0.001:
A1;
t0=0;
u2=A*stepfun（t2,t0）;
t3=A1:
0.001:
1;
n3=length（t3）;
u3=zeros（1,n3）;
t=[t1t2t3];
u=[u1u2u3];
plot（t,u）
axis（[-0.510A+2]）
2、任意函数
例1.2：
用MATLAB画出如下表达式的脉冲序列
clearall;
t=-2:
1:
3;
N=length（t）;
x=zeros（1,N）;
x
（1）=0.4;
x
（2）=0.8
x（3）=1.2;
x（4）=1.5;
x（5）=1.0;
x（6）=0.7;
stem（t,x）;
axis（[-2.23.201.7]）
3、单位阶跃函数
例1.3：
用MATLAB实现单位阶跃函数
clearall;
t=-0.5:
0.001:
1;
t0=0;
u=stepfun（t,t0）;
plot（t,u）
axis（[-0.51-0.21.2]）
4、斜坡函数
例1.4：
用MATLAB实现g（t）=3（t-1）
clearall;
t=0:
0.01:
3;
B=3;
t0=1;
u=stepfun（t,t0）;
n=length（t）;
fori=1:
n
u（i）=B*u（i）*（t（i）-t0）;
end
plot（t,u）
axis（[-0.23.1-0.26.2]）
5、实指数函数
例1.5：
用MATLAB实现
clearall;
t=0:
0.001:
3;
A=3;
a=0.5;
u=A*exp（a*t）;
plot（t,u）
axis（[-0.23.1-0.214]）
6、正弦函数
例1.6：
用MATLAB实现正弦函数f（t）=3cos（10πt+1）
clearall;
t=-0.5:
0.001:
1;
A=3;
f=5;
fai=1;
u=A*sin（2*pi*f*t+fai）;
plot（t,u）
axis（[-0.51-3.23.2]）
（二）、离散信号及其MATLAB实现
1、单位冲激序列
例2.1：
用MATLAB产生64点的单位冲激序列
clearall;
N=64;
x=zeros（1,N）;
x
（1）=1;
xn=0:
N-1;
stem（xn,x）
axis（[-16501.1]）
2、任意序列
例2.2：
用MATLAB画出如下表达式的脉冲序列
clearall;
N=8;
x=zeros（1,N）;
x
（1）=8.0;
x
（2）=3.4
x（3）=1.8;
x（4）=5.6;
x（5）=2.9;
x（6）=0.7;
xn=0:
N-1;
stem（xn,x）
axis（[-1808.2]）
3、单位阶跃序列
例2.3：
用MATLAB实现单位阶跃函数
clearall;
N=32;
x=ones（1,N）;
xn=0:
N-1;
stem（xn,x）
axis（[-13201.1]）
4、斜坡序列
例2.4：
用MATLAB实现g（n）=3（n-4）点数为32的斜坡序列
clearall;
N=32;
k=4
B=3;
t0=1;
x=[zeros（1,k）ones（1,N-k）];
fori=1:
N
x（i）=B*x（i）*（i-k）;
end
xn=0:
N-1;
stem（xn,x）
axis（[-132090]）
5、正弦序列
例2.5：
用MATLAB实现幅度A=3，频率f=100，初始相位Φ=1.2，点数为32的正弦信号
clearall;
N=32;
A=3;
f=100;
fai=1.2;
xn=0:
N-1;
x=A*sin（2*pi*f*（xn/N）+fai）;
stem（xn,x）
axis（[-132-3.23.2]）
6、实指数序列
例2.6：
用MATLAB实现
，点数为32的实指数序列
clearall;
N=32;
A=3;
a=0.7;
xn=0:
N-1;
x=A*a.^xn;
stem（xn,x）
7、复指数序列
例2.7：
用MATLAB实现幅度A=3，a=0.7，角频率ω=314，点数为32的实指数序列
clearall;
N=32;
A=3;
a=0.7;
w=314;
xn=0:
N-1;
x=A*exp（（a+j*w）*xn）;
stem（xn,x）
8、随机序列
利用MATLAB产生两种随机信号：
rand（1,N）在区间上产生N点均匀分布的随机序列
randn（1,N）产生均值为0，方差为1的高斯随机序列，即白噪声序列
例2.8：
用MATLAB产生点数为32的均匀分布的随机序列与高斯随机序列
clearall;
N=32;
x_rand=rand（1,N）;
x_randn=randn（1,N）;
xn=0:
N-1;
figure
（1）;
stem（xn,x_rand）
figure
（2）;
stem（xn,x_randn）
（三）、离散信号的基本运算
1、信号的延迟
给定离散信号x（n），若信号y（n）定义为：
y（n）=x（n-k），那么y（n）是信号x（n）在时间轴上右移k个抽样周期得到的新序列。
例3.1：
正弦序列y（n）=sin（100n）右移3个抽样周期后所得的序列，MATLAB程序如下：
clearall;
N=32;
w=100;
k=3;
x1=zeros（1,k）;
xn=0:
N-1;
x2=sin（100*xn）;
figure
（1）
stem（xn,x2）
x=[x1x2];
axis（[-1N-1.11.1]）
N=N+k;
xn=0:
N-1;
figure
（2）
stem（xn,x）
axis（[-1N-1.11.1]）
2、信号相加
若信号
，值得注意的是当序列
和
的长度不相等或者位置不对应时，首先应该使两者的位置对齐，然后通过zeros函数左右补零使其长度相等后再相加
例3.2：
用MATLAB实现两序列相加
clearall;
n1=0:
3
x1=[20.50.91];
figure
（1）
stem（n1,x1）
axis（[-1802.1]）
n2=0:
7
x2=[00.10.20.30.40.50.60.7];
figure
（2）
stem（n2,x2）
axis（[-1800.8]）
n=0:
7;
x1=[x1zeros（1,8-length（n1））];
x2=[zeros（1,8-length（n2））,x2];
x=x1+x2;
figure（3）
stem（n,x）
axis（[-1802.1]）
3、信号相乘
信号序列
和
相乘所得信号
的表达式为：
这是样本与样本之间的点乘运算，在MATLAB中可采用“.*”来实现，但是在信号序列相乘之前，应对其做与相加运算一样的操作。
例3.3：
用MATLAB实现上例中两序列相乘
clearall;
n1=0:
3
x1=[20.50.91];
figure
（1）
stem（n1,x1）
axis（[-1802.1]）
n2=0:
7
x2=[00.10.20.30.40.50.60.7];
figure
（2）
stem（n2,x2）
axis（[-1800.8]）
n=0:
7;
x1=[x1zeros（1,8-length（n1））];
x2=[zeros（1,8-length（n2））,x2];
x=x1.*x2;
figure（3）
stem（n,x）
axis（[-1800.35]）
4、信号翻转
信号翻转的表达式为：
y（n）=x（-n），在MATLAB中可以用fliplr函数实现此操作
例3.4：
用MATLAB实现“信号相加”中的
序列翻转
clearall;
n=0:
3
x1=[20.50.91];
x=fliplr（x1）;
stem（n,x）
axis（[-1402.1]）
5、信号和
对于N点信号
，其和的定义为：
例3.5：
用MATLAB实现“信号相加”中的
序列和
clearall;
n=0:
3
x1=[20.50.91];
x=sum（x1）
6、信号积
对于N点信号
，其积的定义为：
例3.5：
用MATLAB实现“信号相加”中的
序列积
clearall;
n=0:
3
x1=[20.50.91];
x=prod（x1）
（四）、离散傅里叶变换的MATLAB实现
例４：
若
是一个N=32的有限序列，利用MATLAB计算它的DFT并画出图形。
N=32;
n=0:
N-1;
xn=cos（pi*n/6）;
k=0:
N-1;
WN=exp（-j*2*pi/N）;
nk=n’*k;
WNnk=WN.^nk;
Xk=xn*WNnk;
figure
（1）
stem（n,xn）
figure
（2）
stem（k,abs（Xk））
在MATLAB中，可以直接利用内部函数fft来实现FFT算法，该函数是机器语言，而不是MATLAB指令写成的，执行速度很快。
常用格式为：
y=fft（x）
y=fft（x,N）
（五）、IIR数字滤波器设计
1、基于巴特沃斯法直接设计IIR数字滤波器
例５.1：
设计一个10阶的带通巴特沃斯数字滤波器，带通频率为100Hz到200Hz，采样频率为1000Hz，绘出该滤波器的幅频于相频特性，以及其冲击响应图
clearall;
N=10;
Wn=[100200]/500;
[b,a]=butter（N,Wn,’bandpass’）;
freqz（b,a,128,1000）
figure
（2）
[y,t]=impz（b,a,101）;
stem（t,y）
2、基于切比雪夫法直接设计IIR数字滤波器
例5.2：
设计一个切比雪夫Ⅰ型数字低通滤波器，要求：
Ws=200Hz,Wp=100Hz,Rp=3dB,Rs=30dB,Fs=1000Hz
clearall;
Wp=100;
Rp=3;
Ws=200;
Rs=30;
Fs=1000;
[N,Wn]=cheb1ord（Wp/（Fs/2）,Ws/（Fs/2）,Rp,Rs）;
[b,a]=cheby1（N,Rp,Wn）;
freqz（b,a,512,1000）;
例5.3：
设计一个切比雪夫Ⅱ型数字带通滤波器，要求带通范围100-250Hz，带阻上限为300Hz，下限为50Hz，通带内纹波小于3dB，阻带纹波为30dB，抽样频率为1000Hz，并利用最小的阶次实现。
clearall;
Wpl=100;
Wph=250;
Wp=[Wpl,Wph];
Rp=3;
Wsl=50;
Wsh=300;
Ws=[Wsl,Wsh];
Rs=30;
Fs=1000;
[N,Wn]=cheb2ord（Wp/（Fs/2）,Ws/（Fs/2）,Rp,Rs）;
[b,a]=cheby2（N,Rp,Wn）;
freqz（b,a,512,1000）;
（六）、FIR数字滤波器设计
1、、在MATLAB中产生窗函数十分简单：
（1）矩形窗（RectangleWindow）
调用格式：
w=boxcar（n），根据长度n产生一个矩形窗w。
（2）三角窗（TriangularWindow）
调用格式：
w=triang（n），根据长度n产生一个三角窗w。
（3）汉宁窗（HanningWindow）
调用格式：
w=hanning（n），根据长度n产生一个汉宁窗w。
（4）海明窗（HammingWindow）
调用格式：
w=hamming（n），根据长度n产生一个海明窗w。
（5）布拉克曼窗（BlackmanWindow）
调用格式：
w=blackman（n），根据长度n产生一个布拉克曼窗w。
（6）恺撒窗（KaiserWindow）
调用格式：
w=kaiser（n,beta），根据长度n和影响窗函数旁瓣的β参数产生一个恺撒窗w。
2、基于窗函数的FIR滤波器设计
利用MATLAB提供的函数firl来实现
调用格式：
firl（n,Wn,’ftype’,Window），n为阶数、Wn是截止频率（如果输入是形如[W1W2]
的矢量时，本函数将设计带通滤波器，其通带为W1<ω通-省略该参数、高通-ftype=high、带阻-ftype=stop）、Window是窗函数。
例6.1：
设计一个长度为8的线性相位FIR滤波器。
其理想幅频特性满足
Window=boxcar（8）;
b=fir1（7,0.4,Window）;
freqz（b,1）
例6.2：
设计线性相位带通滤波器，其长度N=15，上下边带截止频率分别为W1=0.3π，w2=0.5
π
Window=blackman（16）;
b=fir1（15,[0.30.5],Window）;
freqz（b,1）
例6.3：
MATLAB中的chirp.mat文件中存储信号ｙ的数据，该信号的大部分号能量集中在Fs/4（或二分之一奈奎斯特）以上，试设计一个34阶的FIR高通滤波器，滤除频率低于Fs/4的信号成分，其中滤波器的截止频率为0.48，阻带衰减为30dB，滤波器窗采用切比雪夫窗
clearall;
loadchirp
window=chebwin（35,30）;
b=fir1（34,0.48,’high’,window）;
yfit=filter（b,1,y）;
[Py,fy]=pburg（y,10,512,Fs）;
[Pyfit,fyfit]=pburg（yfit,10,512,Fs）;
plot（fy,10*log10（Py）,’.’,fyfit,10*log10（Pyfit））;
gridon
ylabel（‘幅度（dB）’）
xlabel（‘频率（Hz’）
legend（‘滤波前的线性调频信号’,‘滤波后的线性调频信号’）
升级会员 