《数字信号处理》课程设计

设计题目一: 语音信号的处理与滤波（难度系数：0.8）

———— 广西科技大学

目录

摘要 ............................................. 3 一、课程设计目的： ............................... 3 二、课程设计内容： ............................... 3 三、基本步骤： ................................... 3 四、滤波流程图 ................................... 5 五、具体内容与实现步骤： ......................... 6

1．录制、存储、回放声音信号 ........................................................ 6 2．分别取8000个和16000个数据进行频谱分析 .......................... 7 3．对声音进行滤波 ............................................................................ 9

六、学习心得 .................................... 15

摘要：

本次课程设计主要是利用MATLAB工具对歌曲语音进行录制，然后设计IIR滤波器，并通过所设计的滤波器进行语音信号滤波以及分析等，对以前所学的知识《信号与线性系统》、《MATLAB7.0》进行一个系统的归纳和进一步理解，并且在实际应用操作当中体会《数字信号处理》这门课程的相关知识。

一、课程设计目的：

综合运用本课程的理论知识进行频谱分析以及滤波器设计，通过理论推导得出相应结论，并利用MATLAB作为工具进行实现，从而复习巩固课堂所学的理论知识，提高对所学知识的综合应用能力，并从实践上初步实现对数字信号的处理。

二、课程设计内容：

1．熟悉并掌握MATLAB中有关声音（wave）录制、播放、存储和读取的函数。

2．在MATLAB环境中，使用声音相关函数录制2秒左右自己的声音，抽样率是8000Hz/s。（考虑如何解决一个实际问题：录制刚开始时，常会出现实际发出声音落后录制动作半拍的现象，如何排除对这些无效点的采样？）

3．分别取8000个和16000个数据进行频谱分析，得到幅度和相位谱，比较二者异同并分析原因。

4．针对电话信道（最高3500Hz），设计一个FIR或IIR滤波器进行滤波，

把抽样率转变为7000Hz/s，并进行频谱分析，得到幅度和相位谱。 5．把处理后的所有数据储存为声音文件，与原始声音进行比较。

三、基本步骤：

1．语音信号的采集

熟悉并掌握MATLAB中有关声音（wave）录制、播放、存储和读取的函数，在MATLAB环境中，使用声音相关函数录制3秒左右自己的声音。（考虑如何解决一个实际问题：录制刚开始时，常会出现实际发出声音落后录制动作半拍的现象，如何排除对这些无效点的采样？）

2．语音信号的频谱分析

要求首先画出语音信号的时域波形；然后对语音号进行快速傅里叶变换，得到信号的频谱特性，从而加深对频谱特性的理解。

3．设计数字滤波器并画出其频率响应 给出各滤波器的性能指标：

(1)低通滤波器性能指标：fb＝1000 Hz，fc＝1200 Hz，As＝100dB，Ap＝1dB。 (2)高通滤波器性能指标：fc＝4800 Hz，fb＝5000 Hz As＝100dB，Ap＝1dB。 (3)带通滤波器性能指标：fp1＝1200 Hz，fp2＝3000 Hz，fsc1＝1000 Hz，fsc2

＝3200Hz，As＝100dB，Ap＝1dB。

要求用窗函数法和双线性变换法设计上面要求的3种滤波器。在Matlab中，使用函数butter，cheby1等函数设计IIR滤波器；利用Matlab中的函数freqz画出各滤波器的频率响应。

4．用滤波器对信号进行滤波

要求用自己设计的各种滤波器分别对采集的信号进行滤波，在Matlab中，FIR滤波器利用函数fftfilt对信号进行滤波，IIR滤波器利用函数filter对信号进行滤波。

5．比较滤波前后语音信号的波形及频谱

要求在一个窗口同时画出滤波前后的波形及频谱。 6．回放语音信号

在Matlab中，使用相关声音函数对声音进行回放。辨别滤波前后的声音有何变化。

四、滤波（低通、带通、高通）流程图如下：

观察滤波后 的信号频谱图（与滤波 前频谱图对 比，播放滤波后语音） 画出信号频谱图并对信号滤波 由信号频率及所要取得的信号频率确定滤波器技术指标 由技术指标确定滤波器类型 设计IIR滤波器

五、具体内容与实现步骤：

（1）录制、存储、回放声音信号（韩红---《青藏高原》）。 代码如下： fs=8000;%采样频率 t=5; %录音时间

fprintf('5秒钟录音中……'); x=wavrecord(t*fs, fs,'double'); %wavrecord为声音采集函数，开始录制声音

wavwrite(x,fs,'F:\\shengyin\\yin.wav');%存放所录声音••••••••• fprintf('\\n录音结束\\n');

y=wavread('F:\\shengyin\\yin.wav');%读取此前声音信号路径 sound(y);%播放声音信号

y1=wavread('F:\\shengyin\\yin.wav',[1 40000]); %读取声音文件 subplot(2,1,1);

plot(y1); title('原始声音信号时域图形'); xlabel('时间'); ylabel('幅值'); grid;

y2=wavread('F:\\shengyin\\yin.wav',[5000 30000]); %读取声音文件 subplot(2,1,2); plot(y2);

title('截取后的有效声音信号时域图形'); xlabel('时间'); ylabel('幅值'); grid;

录制刚开始时，出现一些幅值比较小的点，后面也出现一些。出现这种现象的原因主要是录音开始时，未播放声音所致，出现了一些无效点，而后面出现的无效的点，主要是已经没有声音了，但如果我们一直播放声音再开始录音就不会出现这类情况。我们首先要想办法要去掉这些无效点

（2）分别取8000个和16000个数据进行频谱分析，得到幅度和相位谱，比较二者异同并分析原因。 代码如下：

y=wavread('F:\\shengyin\\yin.wav'); %读取声音文件 d=fft(y,8000); %快速傅里叶变换

%对语音信号做8000点的FFT变换，FFT是离散傅立叶变换的快速算法，可以将一个信号变换到频域。有些信号在时域上是很难看出什么特征的，但是如果变换到频域之后，就很容易看出特征了。 k=0:7999;wk=2*k/8000;%关于π的归一化 figure();%窗口图

subplot(2,2,1); plot(wk,abs(d)); %绘制幅频特性图 axis([0 1 0 100]);%关于x\\y轴的取值 xlabel('ω/π'); %关于x轴的标注

ylabel('幅度谱'); %关于y轴的标注 title('取8000个数据进行频谱分析'); %备注标题 grid; %打网格 subplot(2,2,2);

plot(wk,angle(d)); %绘制相频特性图 xlabel('ω/π'); %关于x轴的标注 ylabel('相位谱'); %关于y轴的标注 axis([0 1 -4 4]); ;%关于x\\y轴的取值

title('取8000个数据进行频谱分析'); %备注标题 grid; %打网格

y=wavread('F:\\shengyin\\yin.wav ');%读取声音文件 d=fft(y,16000);

k=0:15999;wk=2*k/16000;%关于π的归一化 subplot(2,2,3); plot(wk,abs(d)); %绘制幅频特性图

axis([0 1 0 200]); ;%关于x\\y轴的取值 xlabel('ω/π'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('取16000个数据进行频谱分析'); %备注标题 grid; %打网格 subplot(2,2,4);

plot(wk,angle(d)); %绘制相频特性图 xlabel('ω/π'); %关于x轴的标注 ylabel('相位谱');%关于y轴的标注

axis([0 1 -4 4]); ;%关于x\\y轴的取值 title('取16000个数据进行频谱分析'); %备注标题 grid %打网格

取8000点进行频谱分析和取16000点进行频谱分析所得的幅度谱和相位

谱的图形基本上是相同的；不同之处是8000点图形线比16000点的图形线梳。

（3）IIR滤波器进行滤波，把抽样率为fs=7000Hz/s，并进行频谱分析，得到幅度和相位谱。

低通滤波器： 代码如下： fs=7000;%采样频率

y=wavread('F:\\shengyin\\yin.wav '); %读取声音文件 wp=2*1200/fs; % Wp为归一化通带截止频率（0到1之间的整数）ws=2*2000/fs; %Ws为归一化阻带截止频率 Rp=0.5; %Rp为单位为分贝的通带波纹系数 Rs=60; %Rs为最小阻带衰减

[N,Wn]=cheb1ord(wp,ws,Rp,Rs);%求阶数N和滤波器截止频率Wn [num,den]=cheby1(N,Rp,Wn); % 传输分子和分母的系数

[h,w]=freqz(num,den);%求频率响应 figure();%窗口图

f=fft(y,1024); %快速傅里叶变换 k=0:1023;wk=2*k/1024;%关于π的归一化 subplot(2,2,1)

plot(wk,abs(f));%绘制幅频特性图 axis([0 1 0 15]); %关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('原语音幅度谱'); subplot(2,2,2)

plot(wk,angle(f)); %绘制相频特性图 xlabel('ω/π'); %关于x轴的标注 ylabel('相位谱'); %关于y轴的标注 title('原语音相位谱');

axis([0 1 -4 4]); %关于x\\y轴的取值 w=filter(num,den,y); % 滤波•

wavwrite(w,fs,'F:\\shengyin\\yin3.wav');%••写入滤波后的语音文件••••••••••

y1=wavread('F:\\shengyin\\yin3.wav'); %•读取滤波后的语音文件• f1=fft(w,1024); %快速傅里叶变换 k=0:1023;wk=2*k/1024;%关于π的归一化 subplot(2,2,3) plot(wk,abs(f1));

axis([0 1 0 15]); %关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('低通滤波后语音幅度谱'); subplot(2,2,4) plot(wk,angle(f1));

axis([0 1 -4 4]);%关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注

title('低通滤波后语音相位谱'); grid;%打网格

sound(y1);%播放滤波后声音

带通滤波器： fs=7000;%采样频率

Wp=[0.3 0.5];% Wp为归一化通带截止频率（0到1之间的整数） Ws=[0.2 0.6];%Ws为归一化阻带截止频率 Rp=0.5; %Rp为单位为分贝的通带波纹系数 Rs=60; %Rs为最小阻带衰减

[N,Wn]=cheb1ord(Wp,Ws,Rp,Rs);%求阶数N和滤波器截止频率Wn [num,den]=cheby1(N,Rp,Wn); %传输分子和分母的系数 [h,w]=freqz(num,den); %求频率响应

y3=wavread('F:\\shengyin\\yin.wav'); %•••••••••••••••• figure();%窗口图

f3=fft(y3,1024); %快速傅里叶变换 k=0:1023;wk=2*k/1024;%关于π的归一化 subplot(2,2,1)

plot(wk,abs(f3));

axis([0 1 0 15]); %关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('原语音幅度谱'); subplot(2,2,2) plot(wk,angle(f3));

axis([0 1 -4 4]);%关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('原语音相位谱'); grid;%打网格

wavwrite(w,fs,'F:\\shengyin\\yin5.wav'); w=filter(num,den,y3); %•••¨

wavwrite(w,fs,'F:\\shengyin\\yin5.wav');%•••••••••••• w=wavread('F:\\shengyin\\yin5.wav'); %•读取滤波后的语音文件• f1=fft(w,1024); %快速傅里叶变换 k=0:1023;wk=2*k/1024;%关于π的归一化 subplot(2,2,3) plot(wk,abs(f1));

axis([0 1 0 3]); %关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('带通滤波后语音幅度谱'); subplot(2,2,4) plot(wk,angle(f1));

axis([0 1 -4 4]);%关于x\\y轴的取值 xlabel('ω/π•'); %关于x轴的标注 ylabel('幅度谱'); %关于y轴的标注 title('带通滤波后语音相位谱'); grid;%打网格

sound(w);%播放滤波后声音

高通滤波器：

fs=7000;%采样频率

Wp=0.4; % Wp为归一化通带截止频率（0到1之间的整数） Ws=0.5; %Ws为归一化阻带截止频率 Rp=0.5;%Rp为单位为分贝的通带波纹系数 Rs=100;%Rs为最小阻带衰减

[N,Wn]=cheb1ord(wp,ws,Rp,Rs);%求阶数N和滤波器截止频率Wn [num,den]=cheby1(N,Rp,Wn,'high'); % 传输分子和分母的系数[h,w]=freqz(num,den);%求频率响应

y3=wavread('F:\\shengyin\\yin.wav '); %读取声音文件 figure()%窗口图••••••••••• f=fft(y,1024); %快速傅里叶变换 k=0:1023;wk=2*k/1024;%关于π的归一化 subplot(2,2,1)

plot(wk,abs(f));xlabel('ω/π'); %绘制幅频特性图 axis([0 1 0 15]);%关于x\\y轴的取值

ylabel('振幅(幅值)'); title('原语音幅度谱'); %原语音幅度谱 subplot(2,2,2)

plot(wk,angle(f)); %绘制相频特性图 axis([0 1 -4 4]);%关于x\\y轴的取值 xlabel('ω/π'); ylabel('振幅(幅值)');

title('原语音相位谱'); %原语音幅度谱 w=filter(num,den,y); %滤波

wavwrite(w,fs,'F:\\shengyin\\yin3.wav');%存储录音信号 y1=wavread('F:\\shengyin\\yin3.wav'); %读取滤波后语音信号 f1=fft(w,1024); %快速傅里叶变换 subplot(2,2,3)

plot(wk,abs(f1)); %绘制幅频特性图 axis([0 1 0 2]);%关于x\\y轴的取值 xlabel('ω/π');

ylabel('振幅(幅值)'); % 低通滤波后语音幅度谱 title('高通滤波后的幅度谱'); subplot(2,2,4)

plot(wk,angle(f1)); %绘制相频特性图 axis([0 1 -4 4]);%关于x\\y轴的取值 xlabel('ω/π'); ylabel('振幅(幅值)');

title('高通滤波后语音相位谱'); %原语音幅度谱 grid;%打网格

sound(y1);%播放滤波后声音

设计界面如下：

六、学习心得

这门《数字信号处理教程》课是我们通信工程的专业课，它与《信号与线

性系统分析》有紧密的联系，打好《信号与线性系统分析》基础对学起《数字信号处理教程》来是事半功倍的。因此，在做这次课程设计时我也时常翻起《信号与线性系统分析》结合相关知识来复习的。这门学科主要以数学为基础，较为抽象，因此我们在学习这门课的时候，学要细心和耐心。

通过这次课程设计，学到了很多东西。

首先，在做这次课程设计前我又把此前上过的课程《MATLAB7.0从入门与精通》重新复习了一遍，进一步深入地了解了相关函数的用法。

其次，通过本次课程设计，再结合课堂上老师所教授的具体分析方法，我加深了对“信号”的理解。对信号的采集，处理，传输，显示，存储和分析等有了一个系统的掌握和理解。在做课程设计的过程中会遇到很多的困难，在困难面前不要放弃，只要有细心和耐心，坚持下去会达到想要的设计结果。

最后，做课程设计贵在持之以恒，遇到困难总是在所难免，不能因一时的困难而放弃，做课程设计重要的是一种过程，在过程中学到了许多学习的方法，独立学习能力和与同学交流精神。