数字信号处理Matlab实现实例有用.docx
《数字信号处理Matlab实现实例有用.docx》由会员分享，可在线阅读，更多相关《数字信号处理Matlab实现实例有用.docx（12页珍藏版）》请在冰豆网上搜索。
数字信号处理Matlab实现实例有用
数字信号处理Matlab_实现实例（有用）
数字信号处理Matlab实现实例离散时间信号与系统
例1-1用MATLAB计算序列-201–13和序列120-1的离
散卷积
解MATLAB程序如下
a[-201-13]
b[120-1]
cconvab
Mlengthc-1
n01M
stemnc
xlabelnylabel幅度
图11给出了卷积结果的图形求得的结果存放在数组c中为-2-413
151-3
例1-2用MATLAB计算差分方程
当输入序列为时的输出结果
MATLAB程序如下
N41
a[08-044036022]
b[107-045-06]
x[1zeros1N-1]
k01N-1
yfilterabx
stemky
xlabelnylabel幅度
图12给出了该差分方程的前41个样点的输出即该系统的单位脉冲响应
1-3用MATLAB计算例1-2差分方程所对应的系统函数的DTFT
1-2差分方程所对应的系统函数为
其DTFT为
用MATLAB计算的程序如下
k256
num[08-044036002]
den[107-045-06]
w0pikpi
hfreqznumdenw
subplot221
plotwpirealhgrid
title实部
xlabel\omega\piylabel幅度
subplot222
plotwpiimaghgrid
title虚部
xlabel\omega\piylabelAmplitude
subplot223
plotwpiabshgrid
title幅度谱
xlabel\omega\piylabel幅值
subplot224
plotwpianglehgrid
title相位谱
xlabel\omega\piylabel弧度
例2-1对连续的单一频率周期信号按采样频率采样截取长度N分别选N
20和N16观察其DFT结果的幅度谱
解此时离散序列即k8用MATLAB计算并作图函数fft用于计算离散傅里
叶变换DFT程序如下
k8
n1[0119]
xa1sin2pin1k
subplot221
plotn1xa1
xlabeltTylabelxn
xk1fftxa1xk1absxk1
subplot222
stemn1xk1
xlabelkylabelXk
n2[0115]
xa2sin2pin2k
subplot223
plotn2xa2
xlabeltTylabelxn
xk2fftxa2xk2absxk2subplot224
stemn2xk2
xlabelkylabelXk计算结果示于图21a和b分别是N20时的截取信号和DFT结果由于截取
了两个半周期频谱出现泄漏c和d分别是N16时的截取信号和DFT结果由
于截取了两个整周期得到单一谱线的频谱上述频谱的误差主要是由于时域中对
信号的非整周期截断产生的频谱泄漏
例2-2用FFT计算两个序列
的互相关函数
解用MATLAB计算程序如下
x[13-112331]
y[21-1120-13]
klengthx
xkfftx2k
ykffty2k
rmrealifftconjxkyk
rm[rmk22krm1k]
m-k1k-1
stemmrm
xlabelmylabel幅度
其计算结果如图22所示
2-3计算两个序列的的互相关函数其中
xn23521–100123530–1–2012ynxn-4enenMATLAB中可以用随机函数rand产生解用MATLAB计算程序如下
x[23521-100123530-1-2012]
y[000023521-100123530-1-2012]
klengthy
erand1k-05
yye
xkfftx2k
ykffty2k
rmrealifftconjxkyk
rm[rmk22krm1k]
m-k1k-1
stemmrm
xlabelmylabel幅度
计算结果如图23a我们看到最大值出现在m4处正好是ynxn的延迟23b是xn自相关函数他和ynyn受到噪声的干扰
第3章无限长单位脉冲响应IIR滤波器的设计方法
例3-1设采样周期T250μs采样频率fs4kHz用脉冲响应不变法和双线性变换法设计一个三阶巴特沃兹滤波器其3dB边界频率为fc1kHz
[BA]butter32pi1000s
[num1den1]impinvarBA4000
[h1w]freqznum1den1
[BA]butter32000025s
[num2den2]bilinearBA4000
[h2w]freqznum2den2
fwpi2000
plotfabsh1-fabsh2-
grid
xlabel频率Hz
ylabel幅值dB
程序中第一个butter的边界频率2π×1000为脉冲响应不变法原型低通滤波器的边界频率第二个butter的边界频率2T2000025为双线性变换法原型低通滤波器的com这两种设计方法所得到的频响虚线为脉冲响应不变法的结果实线为双线性变换法的结果脉冲响应不变法由于混叠效应使得过渡带和阻带的衰减
特性变差并且不存在传输零点同时也看到双线性变换法在z-1即ωπ或f2000Hz处有一个三阶传输零点这个三阶零点正是模拟滤波器在Ω?
处的三阶传输零点通过映射形成的
例3-2设计一数字高通滤波器它的通带为400,500Hz通带内容许有05dB的波动阻带内衰减在小于317Hz的频带内至少为19dB采样频率为1000Hz
wc21000tan2pi40021000
wt21000tan2pi31721000
[Nwn]cheb1ordwcwt0519s
[BA]cheby1N05wnhighs
[numden]bilinearBA1000
[hw]freqznumden
fwpi500
plotf20log10absh
axis[0500-8010]
grid
xlabel
ylabel幅度dB
图32给出了MATLAB计算的结果可以看到模拟滤波器在Ω?
处的三阶零点通过高通变换后出现在ω0z1处这正是高通滤波器所希望得到的
例3-3设计一巴特沃兹带通滤波器其,dB边界频率分别为f2110kHz和f190kHz在阻带f3120kHz处的最小衰减大于,,dB采样频率fs400kHz
w12400tan2pi902400
w22400tan2pi1102400
wr2400tan2pi1202400
[Nwn]buttord[w1w2][0wr]310s
[BA]butterNwns
[numden]bilinearBA400
[hw]freqznumden
fwpi200
plotf20log10absh
axis[40160-3010]
grid
xlabel频率kHz
ylabel幅度dB
图33给出了MATLAB计算的结果可以看出数字滤波器将无穷远点的二阶零点
映射为z?
1的二阶零点数字带通滤波器的极点数是模拟低通滤波器的极点数
的两倍
例3-4一数字滤波器采样频率fs1kHz要求滤除100Hz的干扰其,dB的
边界频率为95Hz和105Hz原型归一化低通滤波器为
w195500
w2105500
[BA]butter1[w1w2]stop
[hw]freqzBA
fwpi500
plotf20log10absh
axis[50150-3010]
grid
xlabel频率Hz
ylabel幅度dB
图34为MATLAB的计算结果
第4章有限长单位脉冲响应FIR滤波器的设计方法2用凯塞窗设计一FIR低通滤波器低通边界频率阻带边界频率阻带衰减
不小于50dB
解首先由过渡带宽和阻带衰减来决定凯塞窗的N和图41给出了以上设计的频率特性a为N30直接截取的频率特性b为凯
塞窗设计的频率特性凯塞窗设计对应的MATLAB程序为
wnkaiser30455
nn[0129]
alfa30-12
hdsin04pinn-alfapinn-alfa
hhdwn
[h1w1]freqzh1
plotw1pi20log10absh1
axis[01-8010]
grid
xlabel归一化频率
ylabel幅度dB
例4-2利用雷米兹交替算法设计一个线性相位低通FIR数字滤波器其指标为通带边界频率fc800Hz阻带边界fr1000Hz通带波动阻带最小衰减At40dB采样频率fs4000Hz
解
在MATLAB中可以用remezord和remez两个函数设计其结果如图42MATLAB程序如下
fedge[8001000]
mval[10]
dev[00559001]
fs4000
[Nfptsmagwt]remezordfedgemvaldevfs
bremezNfptsmagwt
[hw]freqzb1256
plotw2000pi20log10absh
grid
xlabel频率Hz
ylabel幅度dB
函数remezord中的数组fedge为通带和阻带边界频率数组mval是两个边界处的幅值而数组dev是通带和阻带的波动fs是采样频率单位为Hz
第5章数字信号处理系统的实现
例5-1求下列直接型系统函数的零极点并将它转换成二阶节形式
解用MATLAB计算程序如下
num[1-01-03-03-02]
den[101020205]
[zpk]tf2zpnumden
mabsp
dispdispz
disp极点dispp
disp增益系数dispk
soszp2soszpk
disp二阶节disprealsos
zplanenumden
输入到num和den的分别为分子和分母多项式的系数计算求得零极点增益系
数和二阶节的系数
零点
09615
-05730
-0144305850i
-01443-05850i极点
0527606997i
05276-06997i
-0577605635i
-05776-05635i
增益系数
1
二阶节
10000-03885-05509100001155206511
10000028850363010000-1055207679
系统函数的二阶节形式为
极点图见图51
5-2分析五阶椭圆低通滤波器的量化效应其截止频率为04通带纹波为
04dB最小的阻带衰减为50dB对滤波器进行截尾处理时使用函数a2dTm解用以下MATLAB程序分析量化效应
clf
[ba]ellip5045004
[hw]freqzba512
g20log10absh
bqa2dTb5
aqa2dTa5
[hqw]freqzbqaq512
gq20log10abshq
plotwpigbwpigqr
grid
axis[01-805]
xlabel\omega\pi
ylabelGaindB
legend量化前量化后
figure
[z1p1k1]tf2zpba
[z2p2k2]tf2zpbqaq
zplaneplot[z1z2][p1p2]ox
legend量化前的零点量化后的零点量化前的极点量化后的极点
图51a表示系数是无限精度的理想滤波器的频率响应以实线表示以及当滤波器系数截尾到5位时的频率响应以短线表示由图可知系数量化对频带的边缘影响较大经系数量化后增加了通带的波纹幅度减小了过渡带宽并且减小了最小的阻带衰减
图51b给出了系数量化以前和系数量化以后的椭圆低通滤波器的零极点位置由图可知系数的量化会使零极点的位置与它们的理想的标称位置相比发生显著的改变在这个例子中靠近虚轴的零点的位置变动最大并且移向靠它最近的极点的位置只要对程序稍作改变就可以分析舍入量化的影响
为了研究二进制数量化效应对数字滤波器的影响首先需要将十进制表示的滤波器系数转换成二进制数并进行量化二进制数的量化既可以通过截尾法也可以通过舍入法实现我们提供了如下的两个MATLAB程序com这两段程序分别将向量d中的每一个数按二进制数进行截尾或舍入量化量化的精度是小数点以后保留b位量化后返回的向量为beq
functionbeqa2dTdb
beqa2dTdb将十进制数利用截尾法得到b位的二进制数然后将该二进制数再转换为十进制数
m1d1absd
whilefixd10
d1absd2m
mm1
end
beqfixd12b
beqsigndbeq2m-b-1functionbeqa2dRdb
beqa2dRdb将十进制数利用舍入法得到b位的二进制数然后将该二进制数再转换为十进制数
m1d1absd
whilefixd10
d1absd2m
mm1
end
beqfixd12b5
beqsigndbeq2m-b-1第7章多采样率信号处理
例7-1在时域上显示一个信号频率为0042的正弦信号然后以抽取因子3
降采样率并在时域上显示相应的结果比较两者在时域上的特点
MATLAB计算程序如下
M3down-samplingfactor3
fo0042signalfrequency0042
generatetheinputsinusoidalsequence
n0N-1
m0NM-1
xsin2pifom
generatethedown-samplingsquence
yx[1Mlengthx]
subplot211
stemnx1N
title输入序列
xlabel时间n
ylabel幅度
subplot212
stemny
title[输出序列抽取因子为num2strM]
xlabel时间nylabel幅度
图71信号频率为0042
例7-2用汉明窗设计一长度为32的线性相位QMF滤波器组解采用MATLAB设计调用fir2函数设计公共低通滤波器参数缺省即为汉明
窗程序如下
b1fir231[00405055061][11100600]
fork132b2k-1k-1b1k
end
[H1zw]freqzb11256
h1absH1zg120log10h1
[H2zw]freqzb21256
h2absH2zg220log10h2
figure1
plotwpig1-wpig2--
axis[01-10010]
grid
xlabel\omega\piylabel幅度dB
sumh1h1h2h2
d10log10sum
figure2
plotwpidgrid
xlabel\omega\piylabel误差dB
axis[01-0303]
图72a是一个N32的汉明窗设计结果图中实线表示的低通频响虚线
表示它的镜像图72b是基于这种设计方法的分析综合滤波器组的整个频
响从这个图可见重建误差小于?
005dB由于汉明窗设计的频率响应在通带中近乎是平坦的因此最大重建误差发生在这个滤波器的通带
边界和过渡带内
升级会员 