《数字信号处理》课程研究性学习报告.docx
《《数字信号处理》课程研究性学习报告.docx》由会员分享，可在线阅读，更多相关《《数字信号处理》课程研究性学习报告.docx（16页珍藏版）》请在冰豆网上搜索。
《数字信号处理》课程研究性学习报告
《数字信号处理》课程研究性学习报告
指导教师薛健
时间2014.6
【目的】
（1）掌握IIR和FIR数字滤波器的设计和应用；
（2）掌握多速率信号处理中的基本概念和方法；
（3）学会用Matlab计算小波分解和重建。
（4）了解小波压缩和去噪的基本原理和方法。
【研讨题目】
一、
（1）播放音频信号yourn.wav，确定信号的抽样频率，计算信号的频谱，确定噪声信号的频率范围；
（2）设计IIR数字滤波器，滤除音频信号中的噪声。
通过实验研究
，
的选择对滤波效果及滤波器阶数的影响，给出滤波器指标选择的基本原则，确定你认为最合适的滤波器指标。
（3）设计FIR数字滤波器，滤除音频信号中的噪声。
与
（2）中的IIR数字滤波器，从滤波效果、幅度响应、相位响应、滤波器阶数等方面进行比较。
【设计步骤】
【仿真结果】
【结果分析】
由频谱知噪声频率大于3800Hz。
FIR和IIR都可以实现滤波，但从听觉上讲，人对于听觉不如对图像（视觉）明感，没必要要求线性相位，因此，综合来看选IIR滤波器好一点，因为在同等要求下，IIR滤波器阶数可以做的很低而FIR滤波器阶数太高，自身线性相位的良好特性在此处用处不大。
【自主学习内容】
MATLAB滤波器设计
【阅读文献】
老师课件，教材
【发现问题】（专题研讨或相关知识点学习中发现的问题）：
过渡带的宽度会影响滤波器阶数N
【问题探究】
通过实验，但过渡带越宽时，N越小，滤波器阶数越低，过渡带越窄反之。
这与理论相符合。
【仿真程序】
信号初步处理部分：
[x1,Fs,bits]=wavread（'yourn.wav'）;
sound（x1,Fs）;
y1=fft（x1,1024）;
f=Fs*（0:
511）/1024;
figure
（1）
plot（x1）
title（'原始语音信号时域图谱'）;
xlabel（'timen'）;
ylabel（'magnituden'）;
figure
（2）
freqz（x1）
title（'频率响应图'）
figure（3）
subplot（2,1,1）;
plot（abs（y1（1:
512）））
title（'原始语音信号FFT频谱'）
subplot（2,1,2）;
plot（f,abs（y1（1:
512）））;
title（‘原始语音信号频谱'）
xlabel（'Hz'）;
ylabel（'magnitude'）;
IIR：
fp=2500;fs=3500;
wp=2*pi*fp/FS;
ws=2*pi*fs/FS;
Rp=1;
Rs=15;
Ts=1/Fs;
wp=2*pi*fp/FS;
ws=2*pi*fs/FS;
wp1=2/Ts*tan（wp/2）;
ws1=2/Ts*tan（ws/2）;
t=0:
1/11000:
（size（x1）-1）/11000;
Au=0.03;
d=[Au*cos（2*pi*5000*t）]';
x2=x1+d;
[N,Wn]=buttord（wp1,ws1,Rp,Rs,'s'）;
[Z,P,K]=buttap（N）;
[Bap,Aap]=zp2tf（Z,P,K）;
[b,a]=lp2lp（Bap,Aap,Wn）;
[bz,az]=bilinear（b,a,Fs）;%
[H,W]=freqz（bz,az）;
figure（4）
plot（W*Fs/（2*pi）,abs（H））
grid
xlabel（'频率/Hz'）
ylabel（'频率响应幅度'）
title（'Butterworth'）
f1=filter（bz,az,x2）;
figure（5）
subplot（2,1,1）
plot（t,x2）
title（'滤波前时域波形'）;
subplot（2,1,2）
plot（t,f1）;
title（'滤波后时域波形'）;
sound（f1,FS）;
FIR
[x1,Fs,bits]=wavread（'I:
/dsp_2014_project3/yourn'）;
fp=2500;fs=3500;
wp=2*pi*fp/Fs;
ws=2*pi*fs/Fs;Rs=50;
M=ceil（（Rs-7.95）/（ws-wp）/2.285）;
M=M+mod（M,2）;
beta=0.1102*（Rs-8.7）;
w=kaiser（M+1,beta）;
wc=（wp+ws）/2;
alpha=M/2;
k=0:
M;
hd=（wc/pi）*sinc（（wc/pi）*（k-alpha））;
h=hd.*w';
f1=filter（h,[1],x1）;
[mag,W]=freqz（h,[1]）;
figure
（1）
plot（W*Fs/（2*pi）,abs（mag））;
grid;
xlabel（'频率/Hz'）;
ylabel（'频率响应幅度'）;
title（'Kaiser´窗设计Ⅰ型线性相位FIR低通滤波器'）;
figure
（2）
subplot（2,1,1）
plot（t,x1）
title（'滤波前时域波形'）;
subplot（2,1,2）
plot（t,f1）;
title（'滤波后时域波形'）;
sound（f1,Fs）;
二、
（1）音频信号kdqg24k.wav抽样频率为24kHz，用
y=wavread（'kdqg24k'）;
sound（y,16000）;
播放该信号。
试用频域的方法解释实验中遇到的现象;
（2）设计一数字系统，使得sound（y,16000）可播放出正常的音频信号；讨论滤波器的频率指标、滤波器的的类型（IIR,FIR）对系统的影响。
【仿真结果】
【结果分析】
24K的信号用16K播放，频谱会被拉宽，无法正常播放，通过2倍内插，通过滤波器，然后3倍抽取，得到的信号用16K播放器就能正常播放。
【自主学习内容】
功能：
对时间序列进行重采样。
格式：
1.y=resample（x,p,q）
采用多相滤波器对时间序列进行重采样，得到的序列y的长度为原来的序列x的长度的p/q倍，p和q都为正整数。
此时，默认地采用使用FIR方法设计的抗混叠的低通滤波器。
2.y=resample（x,p,q,n）
采用chebyshevIIR型低通滤波器对时间序列进行重采样，滤波器的长度与n成比例，n缺省值为10.
3.y=resample（x,p,q,n,beta）
beta为设置低通滤波器时使用Kaiser窗的参数，缺省值为5.
4.y=resample（x,p,q,b）
b为重采样过程中滤波器的系数向量。
5.[y,b]=resample（x,p,q）
输出参数b为所使用的滤波器的系数向量。
说明：
x--时间序列
p、q--正整数，指定重采样的长度的倍数。
n--指定所采用的chebyshevIIR型低通滤波器的阶数，滤波器的长度与n成比列。
beta--设计低通滤波器时使用Kaiser窗的参数，缺省值为5.
【阅读文献】
PPt课本
【发现问题】（专题研讨或相关知识点学习中发现的问题）：
采样频率与播放频率之间不是整数倍关系
【问题探究】
此时内插和抽取结合实现正常播放
【仿真程序】
fs=24000;
x1=wavread（'I:
/dsp_2014_project3/kdqg24k.wav'）;
sound（x1,16000）;
y1=fft（x1,1024）;
f1=fs*（0:
511）/1024;
f2=fs/2*3*（0:
511）/1024;
figure
（1）
subplot（2,1,1）;
plot（f1,abs（y1（1:
512）））;
title（'原始语音信号24K正常播放频谱'）
xlabel（'Hz'）;
ylabel（'magnitude'）;
subplot（2,1,2）;
plot（f2,abs（y1（1:
512）））;
title（'原始语音信号16K播放频谱'）
xlabel（'Hz'）;
ylabel（'magnitude'）;
y=resample（x1,2,3）;
sound（y,16000）;
y2=fft（y,1024）;
figure
（2）
subplot（2,1,1）;
plot（f1,abs（y1（1:
512）））;
title（'原始语音信号24K正常播放频谱'）
xlabel（'Hz'）;
ylabel（'magnitude'）;
subplot（2,1,2）;
plot（f2,abs（y2（1:
512）））;
title（'原始语音信号16K经过设计的数字系统后播放频谱'）
xlabel（'Hz'）;
ylabel（'magnitude'）;
三、对连续信号
x（t）=40t2（1-t）4cos（12πt）[0在区间[0，2]均匀抽样1024点得离散信号x[k]，其中n（t）是零均值方差为1的高斯噪声
（1）画出信号x（t）的波形；
（2）计算并画出db7小波的5级小波变换系数；
（3）通过观察小波系数，确定阈值化处理的阈值；
（4）对小波系数进行阈值化处理，画出去噪后的信号波形，求出最大误差和均方误差；
（5）对近似系数和小波系数均进行阈值化处理，画出去噪后的信号波形，求出最大误差和均方误差；
（6）用Haar小波基，重复（3）-（5）；
（7）讨论所得结果。
【仿真结果】
（1）
（2）
（3）T=4.3000t=96.84%
（4）db7小波基
Emabs=0.9166
（5）Haar小波基T=3.100;t=96.18%
Emabs=1.1341
（6）db14小波基T=4.6，t=94.39%
Emabs=0.8965
【结果分析】
选用dbp系列小波基对信号去噪时，p值越大，去噪的效果越好，最大重建误差也小
【自主学习内容】
MATLAB小波变换
【阅读文献】
PPT和课本
【发现问题】（专题研讨或相关知识点学习中发现的问题）：
选用dbp系列小波基对信号去噪时，p值越大，去噪的效果越好，最大重建误差也小
【问题探究】
不同dp小波基不同，会对结果产生不同影响
【仿真程序】
N=4096;
k=linspace（0,2,N）;
nt=randn（size（k））;
x=40*k.^2.*（1-k）.^4.*cos（12*pi*k）.*（0figure;plot（x）;title（'signalwithnoise'）;
dwtmode（'per'）;
[C,L]=wavedec（x,6,'db14'）
figure;plot（k,C）;title（'Waveletcoefficients'）;
M=0;fork=1:
4096;T=4.6;
ifabs（C（1,k））<=T;C（1,k）=0;end
ifC（1,k）~=0;M=M+1;end
end
A1=C.*C;U1=0;A2=x.*x;U2=0;
fork=1:
1024;
U1=U1+A1（1,k）;
U2=U2+A2（1,k）;
end
t=U1/U2
[XD,CXD,LXD]=wden（x,'heursure','s','one',6,'db14'）;
s=waverec（CXD,LXD,'db14'）;
figure;subplot（211）;plot（s）;title（'Ronconstructedsignalwithnoniose'）;
subplot（212）;
plot（x-s）;title（'Error'）;
d=max（x-s）;
升级会员 