语音信号处理与滤波(数字信号处理课程设计)

《数字信号处理》

课程设计报告

题目：语音信号处理与滤波

专业：通信工程

学号：__________

姓等_______________

指导教师：周老师

成绩：________________________

二。一二年一月三日

一、课程设计目的

综合运用数字信号处理的理论知识进行频谱分析和滤波器设计，通过理论推

导得出相应结论，再利用MATLAB作为编程工具进行计算机实现，从而加深对

所学知识的理解。

二、课程设计基本要求及分组

1、对知识点的掌握要求：

(1)熟悉离散信号和系统的时域特性。

(2)掌握序列快速傅里叶变换FFT方法。

(3)学会MATLAB的使用，掌握MATLAB的程序设计方法。

(4)利用MATLAB对语音信号进行频谱分析。

(5)掌握MATLAB设计FIR和IIR数字滤波器的方法。

三、详细设计过程

1、问题描述：录制一段自己的语音信号，取不同的数据点对语音信号进行频谱

分析；对所有数据进行插值和抽取处理，改变抽样率再对信号进行频谱分析;

设计FIR和HR数字滤波器，并对被抽样后的语音信号进行滤波，分析滤波

后信号的时域和频域特征，回放语音信号。

2、详细操作步骤和部分运行结果

(1)采集语音信号：

利用Windows下的录音机，录制一段自己的话音，时间控制在4s左右；然

后在MATLAB软件平台下，利用函数wavread对语音信号进行采样，记住

采样频率和采样点数。通过使用wavread函数，理解采样频率、采样位数等

概念。

wavread函数调用格式如下：

y二wavread(file),读取file所规定的wav文件，返回采样值放在向量y中。

[y,fs,iibilb]=waviead(file),采样值放在向量y中，fs表示采样频率(Hz),

bits表示采样位数。

y二wavread(file,N),读取前N点的采样值放在向量y中。

y=wavread(file,[Nl,N2]),读取从N1点到N2点的采样值放在向量y中。

录制并存储语音的程序代码：(函数名’luyin'被主函数"nain，调用)

fs=8000;

t=4;

a=msgbox(,请说话..J);

x=wavrecord(t*fs,fsz1z*Double');名录制语音信号

close(a);

wavwrite(xzfs,，D:\MATLAB\R2010b\work\数字信号课程设计杜\

课程设计.wav，)"存储录音信号

原始语音时域信号波形图：

程序代码：

[xl,fszbits]=wavread(，D:\MATLAB\R2010b\work\数字信号课程

设计'课程设计.wav，)；

sound(xl,fs,bits);名播放语音

figure(1);

subplot(lz1,1);plot(xl);%做原始语音信号的时域图形

title(1原始语音信号作品,);

xlabel('时间n*);ylabel(*n*);

运行结果：

(2)对语音信号进行频谱分析:

在MATLAB中，利用函数fft对信号进行快速傅里叶变换，得到信号的频

谱特性。首先画出语音信号的时域波形，然后对语音信号进行频谱分析。

程序代码：

[xl,fs,bits]=wavread(，D:\MATLAB\R2010b\work\数字信号课程设计

忠、课程设计.wav，)；

yl=fft(xlz1024);

f=fs*(0:511)/1024;

figure(1)

yl=fft(xl);%®length(xl)^I^FFT

yl=fftshift(yl);金平移，是频率中心为0

dertafs=fs/length(xl);的设置频谱的间隔，分辨率

subplot(2,1,1);plot([-fs/2:derta_fs:fs/2-derta_fs],abs(yl));

生画出原始语音信号的频谱图

title(，原始语音信号的频谱杜作品，)；

subplot(2,l,2);plot(abs(yl(l:1024)));title「原始语音信号FFT频谱

杜作品D;

xlabel(，点数N，);

ylabel(,幅值

原始语音信号频谱和FFT转换后的频谱图

原始语音信号的频谱杜艺忠作品

600------------1-------------1-------------1-------------1------------1-------------1-

400-

点数N

(3)设计数字滤波器和画出频率响应：

根据语音信号的特点给出有关滤波器的性能指标：

1低通滤波器性能指标:fp=1000Hz,fc=1200Hz,As=100dB,Ap=ldB；

2高通滤波器性能指标：fc=4800Hz,fp=5000Hz,As=100dB,Ap=ldB；

3带通滤波器性能指标：fpl=1200Hz,fp2=3000Hz,fcl=1000Hz,

fc2=3200Hz,As=100dB,Ap=ldBo

首先用窗函数法(矩形窗(Rectangularwindow)、三角窗(Triangularwindow)、

汉宇窗(Hanningwindow)、海明窗(Hammingwindow)、布拉克曼窗(Blackman

window)、切比雪夫窗(Chebyshevwindow)>巴特里特窗(Bartlettwindow)及凯

塞窗(Kaiserwindow)。)设计上面要求的三种滤波器。在MATLAB中，利用

函数firl设计FIR滤波器；然后再用双线性变换法设计上面要求的三种滤波

器，利用函数butte、chcbyl和cilip设计HR滤波器；最后，利用MATLAB

中的函数freqz画出各种滤波器的频率响应。

(4)HR滤波器的设计过程

HR低通滤波器：

程序代码：()

fs=8000;

T=l/fs;

duration=3;

[y,fs,bits]=wavread(*D:\MATLAB\R2010b\work\数字信号课程设计

杜'课程设计.wavD；

soundview(*D:\MATLAB\R2010b\work\数字信号课程设计杜、课程设

计.wav');

N=length(y);

t=(1:N)/fs;

df=fs/N；

nl=l:N/2;

f=(nl-1)*df;

subplot(2,1,1);

plot(y);

title。原始信号时域的图形，)

xlabel。时间/1)

ylabel('幅值')

subplot(2,1,2);

q-fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，原始信号频谱图形，)

xlabel。频率/Hz，)

ylabel(,幅值/dbD

fs=8000;

fp=1000;

fst=1200;

ap=l;

ast=100;

wp=2*fp/fs;

wst=2*fst/fs;

[N,wc]=cheblord(wp,wst,ap,ast);

[num,den]=chebyl(N,ap,wc,'low');

[h,w]=freqz(num,den);

subplot(2,1,1);

plot(w/pi,abs(h));grid;

xlabel('\omega/\pi');

ylabel。振幅(幅值)，)；

title('契比雪夫I型低通滤波器的幅频响应，);

subplot(2,1,2);

plot(w/pi,20*logl0(abs(h)));grid;

xlabel('\omega/\pi');

ylabel「振幅i分贝)');

title。契比雪夫I型低通滤波器的幅频响应D;

yO=filter(num,den,y);

fftwav=fft(yO);

subplot(2,2,1);

plot(y);

title(,原始信号的波形，)

xlabel。时间/tJ)

ylabel(，幅值，)

subplot(2Z2,2);

plot(yO);

title(，通过契比雪夫I型低通滤波器后的信号图形，)

xlabel(，时间/t')

ylabel(1幅值,)

subplot(2,2,3);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title「滤波前的信号频谱图，)

xlabel(,频率/Hz，)

ylabel「幅值/db，)

subplot(2,2,4);

plot(f,20*logl0(abs(fftwav(nl))));

title。滤波后信号的频谱图，)

xlabel(,频率/HzD

ylabel「幅值/db，)

wavwrite(yO,fs,'D:\MATLAB\R2010b\work\数字信号课程设计\

课程设计工工R_LP.wav，);

[yO,fs,bits]=wavread('D:\MATLAB\R2010b\work\数字信号课程设计

杜\课程设计工工R_LP.wav，);卷语音回放

soundview(yO,fsj课程设计工工R_LP•wav

运行结果如下图:

原始信号的波形通过契比雪夫I型低通滤波器后的信号图先

12

0.51

I0胆

和0

log

-0.5-1

-1-2-------------------------------------

0IM123401234

4

时间a时间外,in

滤波前的信号频谱囱川滤波后信号的频谱囱1口

100I------------------------------------------100

IIR高通滤波器设计

程序代码：

fS=8000;

T=l/fs;

duration=3;

[y,fs,bits]=wavread(*D:\MATLAB\R2010b\work\^(^:,fa^iSW'®i+

忠、课程设计.waW);

soundview(，D:\MATLAB\R2010b\work\数字信号课程设计杜艺忠、课程设

计.wav*);

N=length(y);

t=(l:N)/fs;

df=fs/N;

nl=l:N/2;

f=(nl-1)*df;

subplot(2,1,1);

plot(y);

title「原始信号时域的图形，)

xlabel(，时间/t/)

ylabel。幅值')

subplot(2,1,2);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，原始信号频谱图形，)

xlabel(‘频率/Hz')

ylabel('幅值/db')

Fs=16000;

fp=4000;

fst=3800;

ap=l;

ast=100;

wp=2*fp/Fs;

wst=2*fst/Fs;

(Nzwc]=cheblord(wp,wst,ap,ast);

[num,den]=chebyl(N,ap,wc,'high');

(h,w]=freqz(num,den);

subplot(2,1,1);

plot(w/pi,abs(h));grid;

xlabel('\omega/\pi');

ylabel「振幅(幅值)，);

title「契比雪夫I型高通滤波器的幅频响应1);

subplot(2,1,2);

plot(w/pi,20*logl0(abs(h)));grid;

xlabel('\omega/\pi');

ylabel「振幅(分乂)，)；

title「契比雪夫I型高通滤波器的幅频响应，);

yO=filter(num,den,y);

fftwav=fft(yO);

subplot(2,2,1);

plot(y);

title(，原始信号的波形，)

xlabel(，时间/L)

ylabel(1')

subplot(2,2,2);

plot(yO);

title(，通过契比雪夫I型高通滤波器后的信号图形，)

xlabel「时间/丁)

ylabel。幅值')

subplot(2,2,3);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，滤波前的信号频谱图，)

xlabel「频率/Hz1

ylabel「幅值/db)

subplot(2,2,4);

plot(f,20*logl0(abs(fftwav(nl))));

title(，滤波后信号的频谱图，)

xlabel「频率/Hz，)

ylabel('iB{fi/db')

wavwrite(yOzfs,'D:\MATLAB\R2010b\work\数字信号课程设计'课程

iSitlIR_HP.wav');

[yO,fs,bits]=wavread(，D:\MATLAB\R2010b\work\数字信号课程设iI，

忠、课程设计工工R_HP.wav，)；名语音回放

soundview(yO,fs,1课程.设计工工R_HP.wav');

程序运行结果如下图:

原始信号的波形通过契比雪夫I型高通滤波器后的信号图先

胆胆

log悭

时间/t104

滤波前的信号频谱囱

10050

IIR带通滤波器的设计

程序代码：

fs=8000;

T=l/fs;

duration=3;

[y,fs,bits]=wavread(1D:\MATLAB\R2010b\work\数字信号课程设计

忠、课程设计.wav1);

soundview「D:\MATLAB\R2010b\work\数’力：号课程设计杜'课程设

计.wav');

N=length(y);

t=(1:N)/fs;

df=fs/N;

nl=l:N/2;

f=(nl-1)*df;

subplot(2,1,1);

plot(y);

title(,原始信号时域的图形,)

xlabel「时间/1，)

ylabelL幅值,)

subplot(2,1,2);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，原始信号频谱图形，)

xlabel「频率/Hz，)

ylabel(，幅值/db，)

fs=8000;

ast-100;

ap=l;

wp=[12003000]*2/fs;

wst=[10003200)*2/fs;

[Nzwn]=cheblord(wp,wst,apzast);

(num,den]=chebyl(N,ap,wn);

[h,w]=freqz(num,den);

subplot(2,1,1);

plot(w/pi,abs(h));grid;

xlabel('\omega/\pi');

ylabel(，振幅(幅值)，);

title(，契比雪夫I型带通滤波器的幅频响应，);

subplot(2,1,2);

plot(w/pi,20*logl0(abs(h)));grid;

xlabel('\omega/\pi');

ylabel「振幅(分贝),);

title(，契比雪夫I型带通滤波器的幅频响应，)；

yl=filter(num,den,y);

fftwav=fft(yl);

subplot(2,2,1);

plot(y);

title(，原始信号的波形，)

xlabel「时间/L)

ylabel。幅值,)

subplot(2,2,2);

plot(yl);

title(，通过工工R带通滤波器后的信号图形，)

xlabell时间/1)

ylabel「幅值，)

subplot(2,2,3);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，滤波前的信号频谱图，)

xlabel「频率/Hz，)

ylabelC幅值/db1)

subplot(2,2,4)；

plot(f,20*logl0(abs(fftwav(nl))));

title(，带通工工R滤波器滤波后信号的频谱图，)

xlabel「频率/Hz，)

ylabel('幅值/db')

wavwrite(yl,fsJD:\MATLAB\R2010b\work\数字信号课程设计杜、课程

设计工IR_BP.wav，);

[yl,fs,bits]=wavread「D:\MATLAB\R2010b\work\数字信号课程设计

忠、课程设计工工R_BP.wav，)；传语音回放)

soundview(yl,fs,'课程设计工工R.BP.wav');

程序运行结果：

原始信号的波形通过炽带通滤波器后的信号图形

1

0.5

胆胆

悭ipy0

-0.5

-1

01234

带通IIR滤波器滤卷9信号的威普图

50

(5)FIR滤波器的设计

FIR低通滤波器：

程序代码：

fs=8000;

T=l/fs;

duration=3;

f

[y,fszsbits]=wavread(D:\MATLAB\R2010b\work\数字信号课程设计

杜艺'课程设计.wav，)；

soundview('D:\MATLAB\R2010b\work\数字信号课程设计杜\课程设

计.wav*);

N=length(y);

t=(1:N)/fs;

df=fs/N;

nl=l:N/2;

f=(nl-l)*df;

subplot(2,1,1);

plot(y);

title「原始信号时域的图形,)

xlabel「时恒/丁)

ylabel「幅值')

subplot(2,1,2);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，原始信号频谱图形，)

xlabel。频率/Hz，)

ylabel「福值/db，)

fs=8000;

fp=1000;

fb=1200;

ap=l;

ast=100;

wp=2*pi*fp/fs;

wst=2*pi*fb/fs;

deltaW=wst-wp;

wc=(wst+wp)/2;

N=ceil(2*6.4*pi/deltaW);

hdWindow=ideallp(wc,N);%---ideally()函数(II：系统自有函数)在系

统安装目录的WORK子目录中，ideallp.m

wdWindow=Kaiser(N,10.056);

hr=wdWindow.*hdWindow';

n=0:N-1;

subplot(2,2,1)；

stem(n,wdWindow)；

xlabel(,时间i)

ylabel('振幅作品，)

title「凯泽窗，)

subplot(2,2,2);

stem(n,hr);

xlabel(,时间i)

ylabel(，振幅,)

title(，凯泽窗冲激响应，)

[H,W]=freqz(hr,1);

subplot(2,2,3);

plot(W/pizabs(H));

xlabel('\omega/\pi')

ylabel(，绝对振幅，)

title。低通滤波器幅频特性(1)，)

subplot(2,2,4);

plot(W/pi,20*logl0(abs(H)));

xlabel('\omega/\pi')

ylabel「相对振幅，)

title。低通滤波器幅频特性(2)，)

yO=fftfilt(hr,y);

fftwav=fft(yO);

subplot(2,2,1);

plot(y);

title。原始信号的波形，)

xlabel('IbJfHj/t')

ylabel。幅值，)

subplot(2,2,2);

plot(yO);

title(，通过凯泽低通滤波器后的信号图形，)

xlabel('时间/t')

ylabel「幅值，)

subplot(2,2,3);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title「滤波前的信号频谱图，)

xlabel('频率/Hz')

ylabel(,i|g{fi；db,)

subplot(2,2,4);

plot(f,20*logl0(abs(fftwav(nl))));

title(,滤波后信号的频谱图,)

xlabel('频率/Hz

ylabel「幅值/db，)

wavwrite（yO,fsJD:\MATLAB\R2010b\work\数字信号课程设计艺'

课程设计F工1\_1^.讨@口）;

[yO,fs,bits]=wavread（*D:\MATLAB\R2010b\work\数字号课程设计

杜、课程设计F工R_LP.wav，）；密语音同放

soundview（yO,fsJ课程设计F工R_LP.wav，）;

程序运行结果：

原始信号的波形通过凯泽低通滤波器后的信号图形

2

1

胆胆0

।亚log-

-1

时同/t.4

滤波后信号的频谱囱1口n

100

FIR高通滤波器：

程序代码：

fs=8000;

T=l/fs;

duration=3;

[y,fs,sbits]=wavread(，D:\MATLAB\R2010b\work\数字信号课程设计总'课

程设计.wav，)；

soundview「D:\MATLAB\R2010b\work\数字信号课程设计杜、课程设

计.wav')；

N=length(y);

t=(1:N)/fs;

df=fs/N;

nl=l:N/2;

f=(nl-l)*df;

subplot(2,1,1);

plot(y);

titled原始信号时域的图形,)

xlabel('WlHj/t')

ylabel「幅值')

subplot(2,1,2);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，原始信号频谱图形，)

xlabol(,频j率/Hz')

ylabel「幅值/db，)

Fs=16000;

fp=4000;

fst=3800;

ast=100;

ap=l;

wp=2*pi*fp/Fs;

wst-2*pi*fst/F3;

wc=(wst+wp)/2;

deltaW=wp-wst;

N0=ceil(2*6.4*pi/deltaW);%——查表7-3(P342)凯泽窗

N=N0+mod(N0+l,2);---确保N为奇数

N=ceil(2*6.4*pi/deltaW);

hdWindow=ideallp(pi,N)-ideallp(wc,N);%---idealip()函数(非系统自

有函数)在系统安装目录的WORK子目录中，ideallp.ru

wdWindow=Kaiser(N,10.056);

hr=wdWindow.*hdWindow,;

n=0:N-l;

subplot(2,2,1);

stem(nzwdWindow);

xlabel(,时间，)

ylabel「振幅，)

title(，凯泽窗，)

subplot(2,2,2);

stem(n,hr);

xlabel。时间D

ylabel「振幅，)

title(，凯泽窗冲激响应，)

[H,W]=freqz(hr,1);

subplot(2,2,3);

plot(W/pizabs(H));

xlabel(*\omega/\pi')

ylabel「绝对振幅,)

title(，凯泽高通滤波器幅频特性(1)1)

subplot(2,2,4);

plot(W/pi,20*logl0(abs(H)))；

xlabel(1\omega/\pi')

ylabel「相对振幅，)

title(，凯泽高通滤波器幅频特性(2)，)

yO=fftfilt(hr,y);

fftwav=fft(yO);

subplot(2,2,1);

plot(y);

title(，原始信号的波形，)

xlabel('时间/t')

ylabel。幅值’)

subplot(2,2,2);

plot(yO);

title('通过凯泽高通滤波器后的信号图形()

xlabel。时间/L)

ylabel「巾鬲值，)

subplot(2,2,3);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，滤波前的信号频谱图，)

xlabel「频率/Hz')

ylabel(，幅值/db，)

subplot(2,2,4);

plot(f,20*logl0(abs(fftwav(nl))));

title(•凯泽而通滤波器流波后信号的频谱图')

xlabel(,频率/Hz')

ylabel。幅值/db，)

wavwrite(yO,fsJD:\MATLAB\R2010b\work\数字信号课程设计杜'课程

SitFIR_HP.wav');

[yO,fs,bits]=wavread「D:\MATLAB\R2010Q\work\数字信号课程设计

忠〜课程设计三工1弁^3/)外语音回放

soundview(yO,fs,,课程设计FIR_HP.wav，);

程序运行结果：

原始信号的波形通过凯泽图通滤波器后的信力图形

1

0.5

胆0

logU

-0.5

1234

时间/t104

滤波前的信号频谱囹

100

5o

pq

通

o

lols

频率/Hz

FIR带通滤波器:

程序代码:

fs=8000;

T=l/fs;

duration=3;

ly,ts,sbitsJ=wavread('D:\MATLAB\R2UlUQ\work\数字信号课程设计'课

程设计.wav];

soundview(*D:\MATLAB\R2010b\work\数字信号课程设计杜、课程设

计.wav');

N=length(y);

t=(1:N)/fs;

df=fs/N;

nl=l:N/2;

f=(nl-l)*df;

subplot(2,1,1);

plot(y);

title(，原始信号时域的图形，)

xlabel('时间/t*)

ylabel「幅值，)

subplot(2,1,2);

q=fft(y);

plot(f,20*logl0(abs(q(nl))));

title(，原始信号频谱图形，)

xlabel「频率/Hz，)

ylabel(,幅值/db，)

fs=8000;

fpl=1200;

fp2=3000;

fstl=1000;

fst2=3200;

ast=l00;

ap=l;

wpl=2*pi*fpl/fs;

wp2=2*pi*fp2/fs;

wstl=2*pi*fstl/fs;

wst2=2*pi*fst2/fs;

wcl=(wstl+wpl)/2;

wc2=(wst2+wp2)/2;

deltaW=min((wpl-wstl),(wst2-wp2));

NO-ceil(2*6.4*pi/deltaW);%——瓷表7-3(P342)0U单窗

N=N0+mod(N0+l,2);%——确保N为奇数

N=ceil(2*6.4*pi/deltaW);

hdWindow=ideallp(wc2,N)-ideallp(wcl,N);%--ideallp()函数(系统

自有函数)在系统安装目录的WORK子目录中，ideallp.m

wdWindow=Kaiser(N,10.056);

hr=wdWindow.*hdWindow1;

n=0:N-1;

subplot(2,2,1);

stem(n,wdwindow);

xlabel('时间1)

ylabel「振幅D

title(,凯泽窗,)

subplot(2,2,2);

stem(n,hr);

xlabel。时间D

ylabel「振幅,)

title(，凯泽窗冲激响应，)

[H,W]=freqz(hr,1);

subplot(2,2,3);

plot(W/pi,abs(H));

xlabel('\omega/\pi*)

ylabel(，绝对振幅i)

title(，带通滤波器幅频特性(1)')

subplot(2,2,4);

plot(W/pi,20*logl0(abs(H)));

xlabel(*\omega/\pi*)

ylabel「相对振幅，）

title「带通滤波器幅频特性（2），）

yO=fftfilt(hr,y);

fftwav=fft(yO);

subplot(2,2,1);

plot(y);

titlet原始信号的波形，)

xlabel(1时间/t1)

ylabel(1幅值')

subplot(2,2,2);

plot(yO);

title(，通过带通滤波器后的信号图形，)

xlabel(，时间/ti*)

ylabel(，幅值，)

subplot(2,2,3);

q=fft(y)；

plot(f,20*logl0(aba(q(nl))))；

title(,滤波前的信号频谱图1)

xlabel「频率/Hz')

ylabel(，幅值/db，)

subplot(2,2,4);

plot(f,20*logl0(abs(fftwav(nl))));

title(，带通滤波器滤波后信号的频谱图，)

xlabel频率/Hz')

ylabel('幅值/db')

wavwrite(yUztsJD:\MATLAB\R2UlUb\work\数字信号课程设计'课程

设ilFIR_BP.wav');

[yO,fs,bits]=wavread('D:\MATLAB\R2010b\work\数字信弓课程设计杜艺

忠、课程设计FIR_BP.wav];%语音回放

soundview(yO,fs,,课程设计F工R_BP•wav1);

程序运行结果：

原始信号的波形通过带通滤波器后的信号图形

11

0.50.5

|0胆

log0

-0.5-0.5

-1-10

4234

时间£fl/时间九翱图

滤波前的信号频谱囹仞带通滤波器滤波后信号的

1001-------------5--0----------------------------

(6)设计系统界面