您的位置：360文档中心› 《数字信号处理》课程研究性学习报告解读

《数字信号处理》课程研究性学习报告解读

《数字信号处理》课程研究性学习报告解读
《数字信号处理》课程研究性学习报告解读

《数字信号处理》课程研究性学习报告

指导教师薛健

时间2014.6

【目的】

(1) 掌握IIR 和FIR 数字滤波器的设计和应用； (2) 掌握多速率信号处理中的基本概念和方法 ； (3) 学会用Matlab 计算小波分解和重建。 (4)了解小波压缩和去噪的基本原理和方法。

【研讨题目】 一、

(1)播放音频信号 yourn.wav ，确定信号的抽样频率，计算信号的频谱，确定噪声信号的频率范围； (2)设计IIR 数字滤波器，滤除音频信号中的噪声。通过实验研究s P ,ΩΩ，s P ,A A 的选择对滤波效果及滤波器阶数的影响，给出滤波器指标选择的基本原则，确定你认为最合适的滤波器指标。

（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('time n');

ylabel('magnitude n');

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）

Wavelet coefficients

（

3）T=4.3000 t=96.84%

(4)db7小波基

Ronconstructed signal with no niose

Error

Emabs=0.9166

(5)Haar 小波基 T=3.100; t=96.18%

Wavelet coefficients

Ronconstructed signal with no niose

Error

Emabs=1.1341

(6)db14小波基T=4.6，t=94.39%

Wavelet coefficients

Ronconstructed signal with no niose

Error

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).*(0

figure;plot(x);title('signal with noise');

dwtmode('per');

[C,L] = wavedec(x,6,'db14')

figure;plot(k,C);title('Wavelet coefficients');

M=0;for k=1:4096;T=4.6;

if abs(C(1,k))<=T;C(1,k)=0;end

if C(1,k)~=0;M=M+1; end

end

A1=C.*C;U1=0;A2=x.*x;U2=0;

for k=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('Ronconstructed signal with no niose'); subplot(212);

plot(x-s);title('Error');

d=max(x-s);

相关主题
相关文档
最新文档