数字信号处理实验讲义4修改版.docx
- 文档编号:18086155
- 上传时间:2023-08-13
- 格式:DOCX
- 页数:28
- 大小:554.09KB
数字信号处理实验讲义4修改版.docx
《数字信号处理实验讲义4修改版.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验讲义4修改版.docx(28页珍藏版)》请在冰点文库上搜索。
数字信号处理实验讲义4修改版
实验一、离散信号的matlab实现
一、实验目的
1、熟悉matlab软件,学会matlab语言的编写
2、使用matlab软件产生一些常见的离散信号
3、掌握用matlab软件作信号的相关分析
二、实验环境
计算机操作系统、matlab软件
三、实验内容
1、用matlab程序产生下列离散信号或连续信号,并画出其波形。
a单位抽样序列的产生
参考程序:
N=100;
x=zeros(1,N);产生一个1行N列值全为0的矩阵,如看成数组x
(1)-x(100)都为0
x
(1)=1;
n=0:
N-1;
stem(n,x);
xlabel(‘n’);
ylabel(‘x(n)’);
title(‘单位抽样序列’)
产生序列
参考程序:
N=100;
x=zeros(1,N);
k=20;
x(k+1)=1;
xn=0:
N-1;
stem(xn,x);
b单位阶跃序列的产生
参考程序:
N=32;
x=ones(1,N);产生一个1行N列值全为1的矩阵
n=0:
N-1;
stem(n,x);
产生序列
参考程序:
N=32;
k=20;
x1=zeros(1,k);
x2=ones(1,N-k);
x=[x1,x2];
xn=0:
N-1;
stem(xn,x);
c模拟信号
,以t=0.01n(n=0:
N-1)进行采样后的离散信号。
参考程序:
N=128;n=[0:
N-1];
t=0.01*n;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
figure
(1);
subplot(211);
stem(t,x);
subplot(212);
stem(n,x);
d产生一个sinc(t)=sint/t抽样函数
参考程序:
n=200;
step=4*pi/n;
t=-2*pi:
step:
2*pi;
y=sinc(t);
plot(t,y,t,zeros(size(t)));%同时画出y(t)和横轴
gridon;
plot(t,y,t,zeros(size(t)),zeros(size(y)),y);%同时画出y(t)和横轴、纵轴
e方波信号square(t)square(t,duty)产生周期是2pi,幅度为正负1的方波,duty占空比,高电平跟整个周期的比值
参考程序:
t=0:
0.01:
2*pi;
y=square(t,50);
plot(t,y);
试产生一个周期为2pi,高低电平分别为半个周期的方波信号
2、相关分析去除噪声
x(n)=sin(2*pi*n)+u(n)噪声为高斯分布白噪声,使用相关分析去除噪声,噪声1功率为1,噪声2功率为0.1
%rxy=xcorr(x,y);
%rx=xcorr(x,Mlag,'flag')Mlag表示rx的单边长度,总长度为2Mlag+1,flag---'biased'rx(m)/N--unbiasedrx(m)/(N-abs(m))
参考程序:
N=500;
p1=1;
p2=0.1;
f=1/8;
Mlag=60;
u=randn(1,N);
u2=u*sqrt(p2);
n=[0:
N-1];
s=sin(2*pi*f*n);
x1=u(1:
N)+s;
rx1=xcorr(x1,Mlag,'biased');
subplot(311);
plot(n,x1);%叠加噪声后的正弦信号
subplot(312);
plot(-Mlag:
Mlag,rx1);
x2=u2(1:
N)+s;
rx2=xcorr(x2,Mlag,'biased');
subplot(313);
plot(-Mlag:
Mlag,rx2);
实验二、离散信号的傅里叶变换
一、实验目的
1、进一步熟悉matlab软件的使用,熟悉matlab的编程语言
2、用matlab语言编写程序进行离散信号的傅里叶分析
二、实验原理
设离散序列
,长度为N,其DTFT定义为:
在实际计算中无法取到无限长序列,通常通过无限长序列加窗作有限长序列的DTFT。
三、实验内容
试求序列
,其中n=0,1,……,N-1,N=12的DTFT。
并画出其幅频特性曲线和相频特性曲线。
若N=24,36,120,其幅频特性曲线和相频特性曲线如何变化,为什么?
参考程序:
N=12;
n=[0:
1:
N-1];
xn=cos(n*pi/6);
w=0:
0.01:
2*pi;
M=length(w);
xejw=zeros(1,M);
fori=0:
1:
N-1
xejw=xejw+cos(i*pi/6)*exp(-j*w*i);
end;
xr=abs(xejw);
xphase=angle(xejw);
subplot(211);
plot(w,xr);
xlabel('w');
ylabel('幅度');
title('幅频特性曲线');
subplot(212);
plot(w,xphase);
xlabel('w');
ylabel('相位');
title('相频特性曲线');
N=12:
N=24:
N=36:
N=120:
实验三、FFT频谱分析及应用
一.实验目的
1、通过实验,加深对FFT的理解,熟悉FFT子程序
2、熟悉应用FFT对典型信号进行频谱分析的方法
二、实验原理与方法
在各种信号序列中,有限长序列占有重要地位,对有限长序列,可以利用离散傅里叶变换(DFT)进行分析。
DFT不但可以很好地反映序列的频谱特性,而且易于用快速算法(FFT)在计算机上实现。
设离散序列
,长度为N,其DFT定义为:
k=0,1,2……,N-1
有限长序列的DFT是其Z变换在单位圆上的等距采样,因此可以用于序列的谱分析。
FFT是DFT的一种快速算法,它是对变换式进行一次次分解,使其成为若干小点数的组合,从而减少运算量。
常用的FFT是以2为基数的,其长度
。
它的效率高,程序简单,使用方便。
当要变换的序列长度不等于2的整数次方时,为了使用以2为基数的FFT,可以用末尾补零的方法,使其长度延长至2的整数次方。
在MATLAB中,可以用函数FFT来实现。
格式:
y=FFT(x)
y=FFT(x,n)
说明:
y=FFT(x)为利用FFT算法计算矢量的离散傅里叶变换,当x为矩阵时,y为矩阵x每一列的FFT。
当x的长度为2的整数次方时,则FFT函数采用基2的FFT算法,否则采用稍慢的混合基算法。
y=FFT(x,n)采用n点FFT。
当x的长度小于n时,FFT函数在x的尾部补零,以构成n点数据;当x的长度大于n时,FFT会截断序列x。
三、实验内容
1、模拟信号
,以
进行采样。
求
(1)N=40点FFT的幅度频谱,从图中,能否观察出信号的2个频率分量?
(2)提高采样点数,如N=128,再求该信号的幅度频谱,此时幅度频谱发生了什么变化?
信号的2个模拟频率和数字频率各为多少?
FFT频谱分析结果与理论上是否一致?
2、研究高密度频谱与高分辨率频谱
设有连续信号
以采样频率
对信号
采样,分析下列三种情况的幅频特性。
(1)采样数据长度N=16点,做N=16的FFT,并画出幅频特性
(2)采集数据长度N=16点,补零到256点,做256点的FFT,并画出幅频特性
(3)采样数据长度N=256点,做N=256的FFT,并画出幅频特性。
参考程序一:
N=40;
n=[0:
N-1];
t=0.01*n;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
k=[0:
1:
N-1];
w=2*pi/N*k;
X=fft(x,N);
magx=abs(X);
subplot(211);stem(n,x);title('signalx(n)');
subplot(212);plot(w/pi,magx);title('fftN=40');
xlabel('频率(单位:
pi)');
ylabel('|X|');
grid;
取N=128时,程序参考以上程序.
参考程序二:
(1)N=16;
n=[0:
N-1];
fs=32000;
T=1/32000;
t=T*n;
x=cos(2*pi*6500*t)+cos(2*pi*7000*t)+cos(2*pi*9000*t);
k=[0:
1:
N-1];
w=2*pi/N*k;
X=fft(x,N);
magx=abs(X);
subplot(211);stem(n,x);title('signalx(n)');
subplot(212);plot(w/pi,magx);title('fftN=16');
xlabel('频率(单位:
pi)');
ylabel('|X|');
grid;
(2)N=16;
n=[0:
N-1];
fs=32000;
T=1/32000;
t=T*n;
x=cos(2*pi*6500*t)+cos(2*pi*7000*t)+cos(2*pi*9000*t);
X=fft(x,256);
k=[0:
1:
255];
w=2*pi/256*k;
magx=abs(X);
subplot(211);stem(n,x);title('signalx(n)');
subplot(212);plot(w/pi,magx);title('fftN=16');
xlabel('频率(单位:
pi)');
ylabel('|X|');
title(‘幅频特性’);
grid;
(3)
N=256;
n=[0:
N-1];
fs=32000;
T=1/32000;
t=T*n;
x=cos(2*pi*6500*t)+cos(2*pi*7000*t)+cos(2*pi*9000*t);
k=[0:
1:
N-1];
w=2*pi/N*k;
X=fft(x,N);
magx=abs(X);
subplot(211);stem(n,x);title('signalx(n)');
subplot(212);plot(w/pi,magx);title('fftN=16');
xlabel('频率(单位:
pi)');
ylabel('|X|');
grid;
实验四、IIR数字滤波器的设计
一实验目的
1、掌握脉冲响应不变法和双线性变换法设计IIR数字滤波器的具体方法和原理,熟悉双线性变换法和脉冲响应不变法设计低通,带通IIR数字滤波器的计算机编程.
2、观察双线性变换法和脉冲响应不变法设计的数字滤波器的频域特性,了解双线性变换法和脉冲响应不变法的特点和区别.
3、熟悉Butterworth滤波器和Chebyshev滤波器的频率特性.
二实验原理与方法
本实验利用模拟滤波器设计IIR数字滤波器,这是IIR数字滤波器设计最常用的方法.利用模拟滤波器设计,需要将模拟域的
转换为数字域的
最常用的转换方法为脉冲响应不变法和双线性变换法.
1、脉冲响应不变法
用数字滤波器的单位脉冲响应序列
模仿模拟滤波器的冲激响应
让
正好等于
的采样值,即:
其中T为采样间隔,如果以
及
分别表示
的拉氏变换及
的z变换,则:
在MATLAB中,可用函数impinvar实现从模拟滤波器到数字滤波器的脉冲响应不变映射,调用格式为:
[b,a]=impinvar(c,d,fs)
[b,a]=impinvar(c,d)
其中,c、d分别为模拟滤波器的分子和分母多项式系数向量,fs为采样频率(Hz),缺省值fs=1Hz,b、a分别为数字滤波器分子和分母多项式系数向量。
例:
已知
,T=0.1,利用脉冲响应不变法求H(z)
MATLAB程序:
c=[21];d=[123];T=0.1;
[b,a]=impinvar(c,d,1/T);
执行结果:
b=0.2-0.1881
a=1-1.79160.8187
数字滤波器为:
2、双线性变换法
S平面与z平面之间满足下列映射关系:
S平面的虚轴映射在z平面的单位圆上,s平面的左半平面完全映射到z平面的单位圆内。
双线性变换不存在频率混叠问题。
在MATLAB中,提供了一个叫做bilinear的函数实现这种映射,调用格式为:
[b,a]=bilinear(c,d,fs)
双线性变换是一种非线性变换,即
,这种非线性引起的幅频特性畸变可通过预畸得到校正。
3、一般设计步骤
(1)给定技术指标转换为模拟低通原型设计性能指标。
(2)估计满足性能指标的模拟低通原型阶数和截止频率。
利用MATLAB中buttord、cheb1ord、cheb2ord、ellipord等函数,调用格式如:
[n,wn]=buttord(wp,ws,rp,rs,’s’)
(3)设计模拟低通原型
利用MATLAB中buttap、cheb1ap、cheb2ap、ellipap等函数,调用格式如:
[z,p,k]=buttap(n)
采用上述函数所得到原型滤波器的传递函数为零点、极点、增益形式,需要和函数
[c,d]=zp2tf(z,p,k)配合使用,以转化为多项式形式。
(4)由模拟低通原型经频率变换获得模拟低通、高通、带通或带阻滤波器。
利用MATLAB中lp2lp、lp2hp、lp2bp、lp2bs等函数,调用格式如:
[c1,d1]=lp2lp(c,d,wn)
(5)利用脉冲响应不变法或者双线性不变法,实现模拟滤波器到数字滤波器的映射。
说明:
MATLAB信号处理工具箱还提供了模拟滤波器设计的完全工具函数:
butter、cheby1、cheby2、ellip、besself。
用户只需一次调用就可自动完成以上设计步骤中的3-4步,调用格式如:
[c,d]=butter(n,wn,’ftype’,’s’)
三、实验内容
1、已知fp=0.3kHz,rp=1.2db,fs=0.2kHz,rs=20db,T=1ms,利用双线性变换法设计一个chebyshevI型数字高通滤波器,观察通带损耗和阻带衰减是否满足要求。
按照要求写出matlab程序,并附上所设计的滤波器传递函数H(z)及相应的幅频特性曲线。
2、要求fp=0.2kHz,rp=1db,fs=0.3kHz,rs=25db,T=1ms,分别用脉冲响应不变法和双线性变换法设计一个butterworth数字低通滤波器,观察所设计数字滤波器的幅频特性曲线,记录带宽和衰减量,检查是否满足要求。
比较这两种方法的优缺点。
程序1:
rp=1.5;rs=40;T=0.001;fp=0.3;fs=0.2;wp=2*pi*fp;ws=2*pi*fs;%wp,ws为数字频率
wp1=(2/T)*tan(wp/2);ws1=(2/T)*tan(ws/2);%预畸为模拟频率
[n,wn]=cheb1ord(wp1,ws1,rp,rs,'s');%求阶数n,3db截止频率wn
[z,p,k]=cheb1ap(n,rp);%模拟低通原型
[c,d]=zp2tf(z,p,k);%转化为多项式形式
[c1,d1]=lp2hp(c,d,wn);%模拟低通到模拟高通
[b,a]=bilinear(c1,d1,1/T);%双线性变换法设计出数字滤波器
[db,mag,pha,grd,w]=freqz_m(b,a);%作数字滤波器的幅频响应和相频响应
subplot(211);
plot(w/pi,mag);
title('幅频特性');
xlabel('w(/pi)');
ylabel('|H(jw)|');
subplot(212);
plot(w/pi,pha/pi);
title('相频特性');
xlabel('w(/pi)');
ylabel('pha(/pi)');
freqz_m函数如下:
function[db,mag,pha,grd,w]=freqz_m(b,a)
%滤波器的幅值响应(相对、绝对)
%Usage:
[db,mag,pha,grd,w]=freqz_m(b,a)
%w采样频率;b系统函数H(z)的分子项(对FIR,b=h)
%a系统函数H(z)的分母项(对FIR,a=1)
[H,w]=freqz(b,a); %500点的复频响应
mag=abs(H); %绝对幅值响应
db=20*log10(mag/max(mag)); %相对幅值响应
pha=angle(H); %相位响应
grd=grpdelay(b,a,w); %群延迟响应
程序2:
fp=0.2;fs=0.3;
wp=2*pi*fp;ws=2*pi*fs;%数字频率
rp=1;rs=25;T=0.001;
wp1=tan(wp/2);ws1=tan(ws/2);%预变形为模拟频率
[n,wn]=buttord(wp1,ws1,rp,rs,'s');%设计模拟滤波器阶数
[c,d]=butter(n,wn,'s');
[bz,az]=bilinear(c,d,1/2);%双线性变换法由s到z
[h,w]=freqz(bz,az,1024,1000);
plot(w,abs(h));gridon;
程序3:
fp=200;fs=300;
rp=1;rs=25;T=0.001;
wp=2*pi*fp;ws=2*pi*fs;%模拟频率
[n,wn]=buttord(wp,ws,rp,rs,'s');%设计模拟滤波器的阶数
[z,p,k]=buttap(n);
[bp,ap]=zp2tf(z,p,k);
[bs,as]=lp2lp(bp,ap,wp);
[bz,az]=impinvar(bs,as,1000);
%[bz,az]=butter(n,wn/pi);%设计数字滤波器
[db,mag,pha,grd,w]=freqz_m(b,a);%作出数字滤波器的频率特性图
subplot(211);
plot(w/pi,mag);
title('幅频特性');
xlabel('w(/pi)');
ylabel('|H(jw)|');
subplot(212);
plot(w/pi,pha/pi);
title('相频特性');
xlabel('w(/pi)');
ylabel('pha(/pi)');
实验五、FIR数字滤波器的设计
一、实验目的:
1、掌握用窗函数法设计FIR滤波器的原理和方法,熟悉相应的计算机编程。
2、了解用频率采样法设计FIR滤波器的原理和实现。
3、熟悉线性相位FIR滤
波器的幅频特性和相频特性。
4、了解不同窗函数对滤波器性能的影响。
二、实验原理与方法
1、窗函数法:
窗函数法设计线性相位FIR滤波器步骤:
(1)确定数字滤波器的性能要求:
临界频率{wk},滤波器单位脉冲响应长度N。
(2)根据性能要求,合理选择h(n)的奇、偶对称性,从而确定理想频率响应
的幅频特性和相位特性。
(3)求理想单位脉冲响应
。
在实际计算中,可对
按M等间隔采样,并对其求IDFT得
,用
代替
。
(4)选择合适的窗函数w(n),由
求所设计的FIR滤波器单位脉冲响应。
(5)求
,分析其幅频特性,若不满足要求,可适当改变窗函数形式或长度N,重复上述设计过程,以得到满意的结果。
窗函数的傅氏变换
的主瓣决定了
的过渡带宽。
的旁瓣大小和多少决定了
在通带和阻带范围内的波动幅度。
2.频率采样法
频率采样法是从频域出发,将给定的
加以等间隔采样:
然后以此
作为实际数字滤波器的频率特性的采样值
,即令
由H(k)通过IDFT可得有限长序列h(n)
代入到z变换中可得:
用频率采样法设计线性相位FIR滤波器的一般步骤为:
(1)由设计要求选择滤波器的种类
(2)根据线性相位的约束条件,确定
和
,进而得到
(3)将
代入
内插公式得到所设计滤波器的频率响应。
关于第3步,在MATLAB中可由函数
h=real(ifft(H,N))和[db,mag,pha,grd,w]=freqz_m(h,l)实现。
三、实验内容:
1、用窗函数法设计一线性相位FIR低通滤波器,设计指标为:
(1)选择一个合适的窗函数,取N=15,确定脉冲响应,并给出所设计的滤波器的频率响应图,分析是否满足设计要求。
(2)若取N=45,重复这一设计,观察幅频和相位特性变化,分析长度N变化的影响。
(3)保持N=45不变,改变窗函数,(如由hamming窗变为blackman窗),观察并记录窗函数对滤波器幅频特性的影响,比较两种窗的特点。
2、用kaiser窗设计一个数字带通滤波器,设计指标为:
低阻带:
低通带:
高通带:
高阻带:
3、用频率采样法设计一个低通滤波器:
问:
(1)采样点数N=33,过渡带设置1个采样点,H(k)=0.5,最小阻带衰减为多少,是否满足设计要求?
(2)采样点数N=34,过渡带设置2个采样点,
最小阻带衰减为多少,是否满足设计要求?
参考程序1:
N=15;M=128;
b1=fir1(N,0.3,boxcar(N+1));
b2=fir1(N,0.3,hamming(N+1));
h1=freqz(b1,1,M);h2=freqz(b2,1,M);
t=0:
N;
subplot(221);stem(t,b2,'.');grid;
f=0:
0.5/M:
0.5-0.5/M;
subplot(222);
plot(f*2,abs(h1),'b-',f*2,abs(h2),'g-');grid;
db1=20*log10(abs(h1)/1);db2=20*log10(abs(h2)/1);
subplot(224);
subplot(224);plot(f*2,db1,'b-',f*2,db2,'g-');
参考程序1(方法二):
N=15;M=128;
b1=fir1(N,0.3,boxcar(N+1));
b2=fir1(N,0.3,hamming(N+1));
[h1,w]=freqz(b1,1,M);h2=freqz(b2,1,M);
t=0:
N;
subplot(221);stem(t,b2,'.');grid;
f=0:
0.5/M:
0.5-0.5/M;
subplot(222);
plot(w/pi,abs(h1),'b-',w/pi,abs(h2),'g-');grid;
db1=20*log10(abs(h1)/1);db2=20*log10(abs(h2)/1);
subplot(224);
subplot(224);plot(w/pi,db1,'b-',w/pi,db2,'g-');
参考程序2:
N=15;M=128;
b2=fir1(N,[0.35,0.65],kaiser(N+1));
h2=freqz(b2,1,M);
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号 处理 实验 讲义 修改
![提示](https://static.bingdoc.com/images/bang_tan.gif)