1、数字信号实验实验程序: 1 产生10点的单位抽样序列(n);%function unit_pulse(N)% unit_pulse.mN=10;x=zeros(1,N);x(1)=1;n=0:N-1;figure(1);stem(n,x);xlabel(单位抽样序列)axis(-1 20 0 1.1)2 产生10点同时移位3位的单位抽样序列(n-3);%function shift_unit_pulse(N,k)% shift_unit_pulse.mN=10;k=3;x=zeros(1,N);x(k+1)=1;n=0:N-1;figure(2);stem(n,x);xlabel(移位3位的单
2、位抽样序列)axis(-1 20 0 1.1)或function x, n=i shift_unit_pulse (n0,ns,nf)n=0:9;x=(n-3)=03 产生任意序列f(n)=8(n)+7(n-1)+6(n-2) +5(n-3)+ 4(n-4)+7(n-5);%function arbitrary_pulse(N)% arbitrary_pulse.mN=10x=zeros(1,N);x(1)=8;x(2)=7;x(3)=6;x(4)=5;x(5)=4;x(6)=7;n=0:N-1;figure(3);stem(n,x);xlabel(任意序列f(n)axis(-1 20 0 9
3、)4 产生N32点的单位阶跃序列%function unit_step(N)% unit_step.mN=32;x=ones(1,N);n=0:N-1;figure(4);stem(n,x);xlabel(单位阶跃序列)axis(-1 32 0 1.1)5 产生斜率为3,n0=4,点数为20点的斜坡序列g(n)=B(n-n0)%function slope(N,k,B)% slope.mN=20;k=4;B=3;x=zeros(1,k) ones(1,N-k);for i=1:N x(i)=B*x(i)*(i-k);endn=0:N-1;figure(5);stem(n,x);xlabel(斜
4、坡序列)axis(-1 10 0 90)6 产生幅度A=3,频率f100,初始相位1.2,点数为32点的正弦序列。%function sine(N,A,f,fai)% sine.mA=3;f=100;fai=1.2;N=32;n=0:N-1;x=A*sin(2*pi*f*(n/N)+fai);figure(6);stem(n,x);xlabel(正弦序列)axis(-1 32 -3.2 3.2)7 产生幅度A=3,角频率314,点数为32点的复正弦序列。%function complex_sine(N,A,w)% complex_sine.mA=3; w=314;N=32n=0:N-1;x=A
5、*exp(j*w*n);figure(7);stem(n,x);xlabel(复正弦序列)axis(-1 32 -3.2 3.2)8 产生幅度A=3,a0.7,点数为32点的实指数序列。%function real_exponent(N,A,a)% real_exponent.ma=0.7;A=3;n=0:N-1;x=A*a.n;figure(8);stem(n,x);xlabel(实指数序列)axis(-1 32 0 3.2)四:实验结果和分析实验程序:1 求序列的谱分析。160长度更准确,整周期采样估计更准确。clear allclose allclc%序列周期为r(2*pi/w)=r(2
6、/0.225)=80/9%选定序列的长度为80点可准确估计N=160;n=0:N-1;x=cos(0.225*pi*n);%生成信号y=fft(x); %求信号的频谱subplot(2,1,1) %画图plot(n,x)xlabel(n)ylabel(x(n)subplot(2,1,2)plot(n/N*2,abs(y)xlabel(nomorlized frequency)ylabel(Y(exp(jomega)xlim(0,1)figure;subplot(2,1,1)stem(n,x)xlabel(n)ylabel(x(n)subplot(2,1,2)stem(n,abs(y)xlabe
7、l(k)ylabel(Y(k)N=256; %点数不同,估计准确度不同n=0:N-1;x=cos(0.225*pi*n);y=fft(x);figure;subplot(2,1,1)plot(n,x)xlabel(n)ylabel(x(n)subplot(2,1,2)plot(n/N*2,abs(y)xlabel(nomorlized frequency)ylabel(Y(exp(jomega)xlim(0,1)figure;subplot(2,1,1)stem(n,x)xlabel(n)ylabel(x(n)subplot(2,1,2)stem(n,abs(y)xlabel(k)ylabel
8、(Y(k)2求信号的谱分析,绘出的幅相特性。给出采样频率、采样点数,谱分辨率等参数的选择原因及结果。clear allclose allclcf1=0.5; f2=1.25;fs=12.5;ts=1/fs;F=0.25;Tp=1/F;t=0:ts:Tp;xt=cos(2*pi*f1*t)+cos(2*pi*f2*t);y=fft(xt); %求信号的频谱plot(t,xt)xlabel(time/s)ylabel(x(t)title(signal)figuresubplot(2,1,1) %画图plot(0:length(y)-1)*F,abs(y)xlabel(frequency/Hz)yl
9、abel(Y(jOmega)subplot(2,1,2)plot(0:length(y)-1)*F,angle(y)xlabel(frequency/Hz)ylabel(Phi(jOmega)尝试变换参数,分析信号频谱。3 观察并分析采用不同频率时,对函数的频谱影响。(a):以,对其进行采样得到x1(n)。 (b):以,对其进行采样得到x3(n) a=218.2;b=50*pi;fs=1000; %采样频率1000hzT=1/fs;Tp=50*0.001; %观察时间50微秒M=Tp*fs; %采样点数n=0:M-1y=a*exp(-b*n*T).*sin(b*n*T ) %函数表达式subp
10、lot(3,2,1)stem(n,y)xlabel(n);ylabel(xa(nT);title(fs=1000HZ);axis(0,M,-2,1.2*max(abs(y)yk=T*fft(y,M) %M点FFTK=0:M-1;fk=K/Tp;subplot(3,2,2)plot(fk,abs(yk)xlabel(f(Hz);ylabel(幅度);title(T*FFT,fs=1000HZ);axis(0,fs,0,1.2*max(abs(yk) fs=200; %采样频率200hzT=1/fs;Tp=50*0.001; M=Tp*fs; n=0:M-1y=a*exp(-b*n*T).*sin
11、(b*n*T )subplot(3,2,5)stem(n,y)xlabel(n);ylabel(xa(nT);title(fs=200HZ);axis(0,M,-2,1.2*max(abs(y);yk=T*fft(y,M) %M点FFTK=0:M-1;fk=K/Tp;subplot(3,2,6)plot(fk,abs(yk)xlabel(f(Hz);ylabel(幅度);title(T*FFT,fs=200HZ);axis(0,fs,0,1.2*max(abs(yk)四:实验结果和分析 三、实验所采用的功能函数1. 巴特沃斯滤波器阶数选择函数(1)N,wc=buttord(wp,ws,p,s)
12、作用: 计算巴特沃斯数字滤波器的阶数N和3dB截止频率wc, wc为数字频率,单位rad。说明: 调用参数wp,ws分别为数字滤波器的通带、阻带截止频率的归一化值,要求:0wp1,0ws1。p,s分别为通带最大衰减和组带最小衰减(dB)。当wswp时,为高通滤波器;当wp和ws为二元矢量时,为带通或带阻滤波器,这时wc也是二元向量。(2)N,c=buttord(p,s,p,s,s)作用: 计算巴特沃斯模拟滤波器的阶数N和3dB截止频率c。说明:p,s,c均为实际模拟角频率。模拟频率f:每秒经历多少个周期,单位Hz,即1/s,信号的真实频率,可用于模拟信号和数字信号;模拟角频率:每秒经历多少弧度
13、,单位rad/s,通常只于模拟信号;数字频率w:每个采样点间隔之间的弧度,单位rad,通常只用于数字信号。关系:=2pi*f;w = *T=2pi*f/F。(F=1/Ts为采样频率,Ts为采样间隔)2. 完整巴特沃斯滤波器设计函数(1)格式: b,a=butter(N,wc,ftype) 作用: 计算N阶巴特沃斯数字滤波器系统函数分子、分母多项式的系数向量b、a。说明: 调用参数N和wc分别为巴特沃斯数字滤波器的阶数和3dB截止频率的归一化值,一般是调用buttord格式(1)计算N和wc。系数b、a是按照z-1的升幂排列。(2)格式:B,A=butter(N,c,ftype,s) 作用:计算
14、巴特沃斯模拟滤波器系统函数的分子、分母多项式系数向量。 说明:调用参数N和c分别为巴特沃斯模拟滤波器的阶数和3dB截止频率(实际角频率),可调用buttord(2)格式计算N和c。系数B、A按s的正降幂排列。 tfype为滤波器的类型: ftype=high时,高通;c只有1个值。 ftype=stop时,带阻;c=cl,cu,分别为带阻滤波器的通带3dB下截止频率和上截止频率。 ftype缺省时:若c只有1个值,则默认为低通;若c有2个值,则默认为带通;其通带频率区间cl cu。3. 求离散系统频响特性的函数freqz()格式:H,w=freqz(b,a,N)说明:b和a分别为离散系统的系统
15、函数分子、分母多项式的系数向量,返回量H则包含了离散系统频响在 0pi范围内N个频率等分点的值(其中N为正整数),w则包含了范围内N个频率等分点。调用默认的N时,其值是512。可以先调用freqz()函数计算系统的频率响应,然后利用abs()和angle()函数及plot()函数,绘制出系统的频响曲线。4. 滤波器离散化函数:bilinear(使用双线性变换法把模拟滤波器转换为数字滤波器)impinvar(使用脉冲响应不变法把模拟滤波器转换为数字滤波器)四、实验内容及步骤1. 用直接设计法设计BW(巴特沃斯)低通数字滤波器。 采样频率为2000Hz,通带中允许的最大衰减为0.5dB,阻带内的最
16、小衰减为40dB,通带上限临界频率为30Hz,阻带下限临界频率为40Hz。 设计步骤:(1) 确定滤波器的设计指标:;(2) 运用函数计算巴特沃斯低通滤波器的阶数N和归一化3db截止频率;(3) 运用函数求得低通滤波器的系统函数的分子、分母多项式形式;(4) 作图显示滤波器的幅频特性和相位特性。2. 脉冲响应不变法设计数字滤波器使用脉冲响应不变法设计数字低通滤波器,其指标为:通带临界频率0.5,通带内衰减小于1dB;阻带临界频率0.8,阻带内衰减大于15dB,采样频率为100Hz。 设计步骤:(5) 确定数字频率指标;(6) 采用脉冲响应不变法求得模拟低通滤波器频率设计指标;(7) 用butt
17、erworth设计方法求得模拟低通滤波器的截止频率和阶数;(8) 设计归一化模拟低通滤波器;(9) 利用脉冲响应不变法把模拟滤波器转换为数字滤波器;(10) 画出幅度响应和相位响应图。3. 应用双线性变换方法设计低通数字滤波器数字低通滤波器的设计指标为:通带临界频率0.5,通带内衰减小于1dB;阻带临界频率0.8,阻带内衰减大于15dB,采样频率为100Hz。 设计步骤:(11) 确定数字频率指标;(12) 采用双线性变换法求得模拟低通滤波器频率设计指标;(13) 用butterworth设计方法求得模拟低通滤波器的截止频率和阶数;(14) 设计归一化模拟低通滤波器;(15) 利用双线性变换法
18、把模拟滤波器转换为数字滤波器;(16) 画出幅度响应和相位响应图。 实验程序: 1 直接设计法设计BW(巴特沃斯)低通数字滤波器 fp=40; %带通截止频率 fs=30; %阻通截止频率 ft=200; %采样频率 rp=0.5; rs=40; wp=fp/(ft/2); %利用Nyquist频率进行归一化 ws=fs/(ft/2); n,wc=buttord(wp,ws,rp,rs); %求数字滤波器的最小阶数和截止频率 b,a=butter(n,wc); %设计低通数字滤波器系数b,a H,W=freqz(b,a); %求系统频响特性,W为数字角频率,单位rad figure; plot
19、(W*ft/(2*pi),abs(H);grid; %绘出频率响应曲线 xlabel(频率/Hz);ylabel(幅值); figure plot(W*ft/(2*pi),angle(H);grid; %绘出频率响应曲线 xlabel(频率/Hz);ylabel(相位); 2 脉冲响应不变法设计数字滤波器 clear all; fs=100 wp=0.5*pi; ws=0.8*pi; rp=1; rs=15; Wp=wp*fs; %由数字角频率转换为模拟角频率(脉冲响应不变法) Ws=ws*fs; n,wc=buttord(Wp,Ws,rp,rs,s); %选择滤波器的最小阶数 b,a=but
20、ter(n,wc,s); bz,az=impinvar(b,a,fs); %脉冲相应不变法变换为数字滤波器 H,W= freqz(bz,az); %求解数字滤波器的频率响应 plot(W*fs/(2*pi),abs(H); grid; xlabel(频率/hz);ylabel(幅值/dB); title(脉冲响应不变变换法) figure; plot(W/pi,20*log10(abs(H); grid; xlabel(归一化频率);ylabel(幅值/dB); title(脉冲响应不变变换法) 3 双线性变换法设计数字滤波器 clear all;close all;clc fs=100; w
21、p=0.5*pi; ws=0.8*pi; rp=1; rs=15; Wp=2*fs*tan(wp/2); %由数字角频率转换为模拟角频率(脉冲响应不变法) Ws=2*fs*tan(ws/2); n,wc=buttord(Wp,Ws,rp,rs,s); %选择滤波器的最小阶数 b,a=butter(n,wc,s); bz,az=bilinear(b,a,fs); %脉冲相应不变法变换为数字滤波器 H,W= freqz(bz,az); %求解数字滤波器的频率响应 plot(W*fs/(2*pi),20*log10(abs(H); grid; xlabel(频率/hz);ylabel(幅值/dB); title(双线性变换法) figure; plot(W/pi,20*log10(abs(H); grid; xlabel(归一化频率);ylabel(幅值/dB); title(双线性变换法)五:实验结果和分析
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1