数字信号实验.docx
《数字信号实验.docx》由会员分享，可在线阅读，更多相关《数字信号实验.docx（19页珍藏版）》请在冰豆网上搜索。
数字信号实验
实验程序：
1产生10点的单位抽样序列δ（n）；
%functionunit_pulse（N）
%unit_pulse.m
N=10；
x=zeros（1,N）;
x
（1）=1;
n=0:
N-1;
figure
（1）;
stem（n,x）;xlabel（'单位抽样序列'）
axis（[-12001.1]）
2产生10点同时移位3位的单位抽样序列δ（n-3）；
%functionshift_unit_pulse（N,k）
%shift_unit_pulse.m
N=10;
k=3;
x=zeros（1,N）;
x（k+1）=1;
n=0:
N-1;
figure
（2）;
stem（n,x）;xlabel（'移位3位的单位抽样序列'）
axis（[-12001.1]）
或function[x,n]=ishift_unit_pulse（n0,ns,nf）
n=[0:
9]；
x=[（n-3）==0]
3产生任意序列
f（n）=8δ（n）+7δ（n-1）+6δ（n-2）+5δ（n-3）+4δ（n-4）+7δ（n-5）；
%functionarbitrary_pulse（N）
%arbitrary_pulse.m
N=10
x=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（[-12009]）
4产生N＝32点的单位阶跃序列
%functionunit_step（N）
%unit_step.m
N=32;
x=ones（1,N）;
n=0:
N-1;
figure（4）;
stem（n,x）;xlabel（'单位阶跃序列'）
axis（[-13201.1]）
5产生斜率为3，n0=4,点数为20点的斜坡序列g（n）=B（n-n0）
%functionslope（N,k,B）
%slope.m
N=20;k=4;B=3;
x=[zeros（1,k）ones（1,N-k）];
fori=1:
N
x（i）=B*x（i）*（i-k）;
end
n=0:
N-1;
figure（5）;
stem（n,x）;xlabel（'斜坡序列'）
axis（[-110090]）
6产生幅度A=3，频率f＝100，初始相位
＝1.2，点数为32点的正弦序列。
%functionsine（N,A,f,fai）
%sine.m
A=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（[-132-3.23.2]）
7产生幅度A=3，角频率ω＝314，点数为32点的复正弦序列。
%functioncomplex_sine（N,A,w）
%complex_sine.m
A=3;w=314;N=32
n=0:
N-1;
x=A*exp（j*w*n）;
figure（7）;
stem（n,x）;xlabel（'复正弦序列'）
axis（[-132-3.23.2]）
8产生幅度A=3，a＝0.7，点数为32点的实指数序列。
%functionreal_exponent（N,A,a）
%real_exponent.m
a=0.7;A=3;
n=0:
N-1;
x=A*a.^n;
figure（8）;
stem（n,x）;xlabel（'实指数序列'）
axis（[-13203.2]）
四：
实验结果和分析
实验程序：
1求序列
的谱分析。
160长度更准确，整周期采样估计更准确。
clearall
closeall
clc
%序列周期为r（2*pi/w）=r（2/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（'nomorlizedfrequency'）
ylabel（'Y（exp（j\omega））'）
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（'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（'nomorlizedfrequency'）
ylabel（'Y（exp（j\omega））'）
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（'Y（k））'）
2求信号
的谱分析，绘出
的幅相特性。
给出采样频率
、采样点数
，谱分辨率
等参数的选择原因及结果。
clearall
closeall
clc
f1=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'）
figure
subplot（2,1,1）%画图
plot（（0:
length（y）-1）*F,abs（y））
xlabel（'frequency/Hz'）
ylabel（'Y（j\Omega）'）
subplot（2,1,2）
plot（（0:
length（y）-1）*F,angle（y））
xlabel（'frequency/Hz'）
ylabel（'\Phi（j\Omega）'）
尝试变换参数，分析信号频谱。
3观察并分析采用不同频率时，对函数
的频谱影响。
（a）：
以
，对其进行采样得到x1（n）。
（b）：
以
，对其进行采样得到x3（n）
a=218.2;
b=50*pi;
fs=1000;%采样频率1000hz
T=1/fs;
Tp=50*0.001;%观察时间50微秒
M=Tp*fs;%采样点数
n=0:
M-1
y=a*exp（-b*n*T）.*sin（b*n*T）%函数表达式
subplot（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点FFT
K=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;%采样频率200hz
T=1/fs;
Tp=50*0.001;
M=Tp*fs;
n=0:
M-1
y=a*exp（-b*n*T）.*sin（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点FFT
K=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）
作用：
计算巴特沃斯数字滤波器的阶数N和3dB截止频率wc,wc为数字频率，单位rad。
说明：
调用参数wp，ws分别为数字滤波器的通带、阻带截止频率的归一化值，要求：
0≤wp≤1，0≤ws≤1。
αp，αs分别为通带最大衰减和组带最小衰减（dB）。
当ws≤wp时，为高通滤波器；当wp和ws为二元矢量时，为带通或带阻滤波器，这时wc也是二元向量。
（2）[N,Ωc]=buttord（Ωp，Ωs，αp，αs，‘s’）
作用：
计算巴特沃斯模拟滤波器的阶数N和3dB截止频率Ωc。
说明：
Ωp，Ωs，Ωc均为实际模拟角频率。
模拟频率f：
每秒经历多少个周期，单位Hz，即1/s，信号的真实频率，可用于模拟信号和数字信号；
模拟角频率Ω：
每秒经历多少弧度，单位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’）
作用：
计算巴特沃斯模拟滤波器系统函数的分子、分母多项式系数向量。
说明：
调用参数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分别为离散系统的系统函数分子、分母多项式的系数向量，返回量H则包含了离散系统频响在0～pi范围内N个频率等分点的值（其中N为正整数），w则包含了范围内N个频率等分点。
调用默认的N时，其值是512。
可以先调用freqz（）函数计算系统的频率响应，然后利用abs（）和angle（）函数及plot（）函数，绘制出系统的频响曲线。
4.滤波器离散化函数：
bilinear（使用双线性变换法把模拟滤波器转换为数字滤波器）
impinvar（使用脉冲响应不变法把模拟滤波器转换为数字滤波器）
四、实验内容及步骤
1.用直接设计法设计BW（巴特沃斯）低通数字滤波器。
采样频率为2000Hz，通带中允许的最大衰减为0.5dB，阻带内的最小衰减为40dB，通带上限临界频率为30Hz，阻带下限临界频率为40Hz。
Ø设计步骤：
（1）确定滤波器的设计指标：
；
（2）运用函数
计算巴特沃斯低通滤波器的阶数N和归一化3db截止频率
；
（3）运用函数
求得低通滤波器的系统函数的分子、分母多项式形式；
（4）作图显示滤波器的幅频特性和相位特性。
2.脉冲响应不变法设计数字滤波器
使用脉冲响应不变法设计数字低通滤波器，其指标为：
通带临界频率0.5
，通带内衰减小于1dB；阻带临界频率0.8
，阻带内衰减大于15dB,采样频率为100Hz。
Ø设计步骤：
（5）确定数字频率指标；
（6）采用脉冲响应不变法求得模拟低通滤波器频率设计指标；
（7）用butterworth设计方法求得模拟低通滤波器的截止频率和阶数；
（8）设计归一化模拟低通滤波器；
（9）利用脉冲响应不变法把模拟滤波器转换为数字滤波器；
（10）画出幅度响应和相位响应图。
3.应用双线性变换方法设计低通数字滤波器
数字低通滤波器的设计指标为：
通带临界频率0.5
，通带内衰减小于1dB；阻带临界频率0.8
，阻带内衰减大于15dB,采样频率为100Hz。
Ø设计步骤：
（11）确定数字频率指标；
（12）采用双线性变换法求得模拟低通滤波器频率设计指标；
（13）用butterworth设计方法求得模拟低通滤波器的截止频率和阶数；
（14）设计归一化模拟低通滤波器；
（15）利用双线性变换法把模拟滤波器转换为数字滤波器；
（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（W*ft/（2*pi）,abs（H））;grid;%绘出频率响应曲线
Øxlabel（'频率/Hz'）;ylabel（'幅值'）;
Øfigure
Øplot（W*ft/（2*pi）,angle（H））;grid;%绘出频率响应曲线
Øxlabel（'频率/Hz'）;ylabel（'相位'）;
Ø
Ø2脉冲响应不变法设计数字滤波器
Øclearall;
Ø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]=butter（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双线性变换法设计数字滤波器
Øclearall;closeall;clc
Øfs=100;
Øwp=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（'双线性变换法'）
五：
实验结果和分析
升级会员 