数字信号处理陈后金版作业
数字信号处理作业：
基于MATLAB 的有噪声的语音信号分析与处理设计
一、题目要求：
1） 选择一个语音信号作为分析对象，或录制一段语音信号；
2） 对语音信号进行采样，画出采样后语音信号的时域波形和频谱图； 3） 利用MATLAB 中的正弦函数产生噪声加入到语音信号中，使语音信号被污染，然后进行频谱分析；
4） 设计FIR 和IIR 数字滤波器，并对被噪声污染的语音信号进行滤波，画出滤波前后信号的时域波形和频谱，并对滤波前后的信号进行比较，分析信号的变化；
5） 进行信号频谱的搬移，进行变声实验 6） 回放语音信号
二、题目的实现方法
1) 语音信号的采集
利用PC 机上的声卡和WINDOWS 操作系统可以进行数字信号的采集。将话筒输入计算机的语音输入插口上, 启动录音机。按下录音按钮, 接着对话筒说话“语音信号处理”,说完后停止录音, 屏幕左侧将显示所录声音的长度。点击放音按钮, 可以实现所录音的重现。以文件名“sound ”保存入g :\the sound文件夹中。可以看到, 文件存储器的后缀默认为. wav ,这是WINDOWS 操作系统规定的声音文件存的标准。
2) 语音信号的时频分析
Matlab 软件平台下，利用wavread 函数对语音信号进行采样，记住采样频率和采样点数
Wavread 函数调用格式
y=wavread(file)%读取file 所规定的wav 文件，返回采样值放在向量y 中。 [y,Ft,nbits]=wavread(file) %采样值放在向量y 中，fs 表示采样频率（Hz ），nbits 表示采样位数。
y=wavread(file，N)%读取钱N 点的采样值放在向量y 中。
y=wavread(file，[N1,N2])%读取从N1到N2点的采样值放在向量y 中。 对语音信号shengyin.wav 进行采样其程序如下：
[y,Ft,nbits]=wavered (shengyin); %把语音信号进行加载入Matlab 仿真软件平台中
Ft = 22050 nbits = 16
首先画出语音信号的时域波形，然后对语音信号进行频谱分析。在matlab 中利用fft 对信号进行快速傅里叶变换，得到信号的频谱特性。
其程序如下：
[x,Ft,bits]=wavread('shengyin');%读入语音信号, 返回采样频率Fs x_length=length(x);%计算语音信号的长度
sound(x,Ft,bits);%以Fs 的频率播放语音信号,Fs 的改变可以改变声调，以16位的精
度重新播放语音信号播放
subplot(311);
plot(x);%x是采集得到的离散声音信号 xlabel('抽样时刻t'); title('原始信号波形'); grid;
fz=fftshift(fft(x,x_length));%计算语音离散信号x 的x_length点的DFT ，并重新排列,
把DFT 移位到对称位置；
W=linspace(-pi,pi,x_length); f=W*Fs/(2*pi); subplot(312);
plot(f,abs(fz));%画出幅频图 xlabel('频率');
title('原始信号的频谱'); grid;
subplot(313); plot(f,angle(fz)); xlabel('频率');
title('原始信号的相位'); grid;
程序运行结果如下图1：
原始信号波形
10
-1
[1**********]00
40005000抽样时刻t 原始信号的频谱
[1**********]0
200100
0-4000
-3000
-2000
-1000
01000频率
原始信号的相位
2000
3000
4000
50
-5-4000
-3000-2000-1000
0频率
[**************]0
图1 原始信号的波形
3) 语音信号加噪与频谱分析
用正弦函数sin(wt)对原始语音信号进行加噪声处理。
程序如下：
clf;%清除原先曲线波形 t=0:length(x)-1;
Fn=1500;%干扰正弦信号的频率 omiga=2*pi*Fn/Ft;
noise1=0.1*sin(omiga*t);%加入正弦干扰信号，频率1.5khz x1=x+noise1'; subplot(311); plot(t,x1); grid;
xlabel('时间（t ）');
title('加噪声的信号波形'); fz1=fftshift(fft(x1,length(x))); subplot(312);
plot(f,abs(fz1));%画出加入高频干扰之后的语音信号幅频图 xlabel('频率');
title('加噪声的信号频谱'); grid;
subplot(313); plot(f,angle(fz)); xlabel('频率');
title('加噪声的信号相位'); grid on;
sound(x1,Ft);%播放加入干扰信号后的声音信号
程序运行结果如下图2： 说明：由程序运行结果可以看到，加入噪声后信号的时域信号波形里面已经有噪声成分，可以看到频谱图中两个幅值很大的尖峰，其频率是1500Hz ，说明程序加噪声已经成功，播放加噪声的语音信号时听到“滴滴”的噪声信号。
根据语音信号的特点给出有关滤波器的新能指标： 4) 设计FIR 数字滤波器
根据语音信号的特点, 设计FIR 数字低通滤波器，则滤波器的的性能指标
fp=1000Hz，fc=1200Hz，As=50db ,Ap=1dB
在Matlab 中，可以利用函数fir1设计FIR 滤波器，利用函数butter,cheby1和ellip 设计IIR 滤波器，利用Matlab 中的函数freqz 画出各步步器的频率响应。
分析如下：函数fir1默认的设计滤波器的方法为窗函数法，其中可选的窗函数有Rectangular Barlrtt Hamming Hann Blackman窗，其相应的都有实现函数。 函数butter,cheby1和ellip 设计IIR 滤波器时都是默认的双线性变换法，所以在设计滤波器时只需要代入相应的实现函数即可。
加噪声的信号波形
10
-1
[1**********]00
40005000时间（t ）
加噪声的信号频谱
[1**********]0
500
-4000
-3000-2000-1000
01000频率
加噪声的信号相位
[1**********]0
50
-5-4000
-3000-2000-1000
0频率
[**************]0
图2 加入噪声后的信号的波形
FIR 低通滤波器设计程序
用窗函数设计低通滤波器的程序如下：
Ft=22050;%信号的采样频率 Fp=1000; Fs=1200; wp=2*Fp/Ft; ws=2*Fs/Ft; rp=1; rs=50;
p=1-10.^(-rp/20); %通带阻带波纹 s=10.^(-rs/20); fpts=[wp ws]; mag=[1 0]; dev=[p s];
[n21,wn21,beta,ftype]=kaiserord(fpts,mag,dev);%由kaiserord 求滤波器的阶数和截止频率 b21=fir1(n21,wn21,Kaiser(n21+1,beta)); %由fir1设计滤波器 [h,w]=freqz(b21,1); %得到频率响应 plot(w/pi,abs(h));
title('FIR低通滤波器'); grid;
程序运行结果如下图3：
FIR 低通滤波器
1.4
1.2
1
0.8
0.6
0.4
0.2
00.10.20.30.40.50.60.70.80.91
图3 FIR低通滤波器的频率响应曲线图
5) 用滤波器对加噪语音信号进行滤波
用自己设计的各滤波器分别对加噪的语音信号进行滤波，在Matlab 中，FIR 滤波器利用函数fftfilt 对信号进行滤波，IIR 滤波器利用函数filter 对信号进行滤波。
窗函数法：
低通滤波器
z21=fftfilt(b21,s);
滤波程序如下 z21=fftfilt(b21,x1); sound(z21);
m21=fft(z21); %求滤波后的信号 subplot(2,2,1); plot(abs(x1),'g');
title('滤波前信号的频谱'); grid;
subplot(2,2,2); plot(abs(m21),'r');
title('滤波后信号的频谱'); grid;
subplot(2,2,3); plot(x1);
title('滤波前信号的波形');
grid;
subplot(2,2,4); plot(z21);
title('滤波后的信号波形'); grid;
程序运行结果如下图4：
滤波前信号的频谱
500400
300
100
200100
-4000-2000
2000
4000
500-4000200150
滤波后信号的频谱
-[1**********]00
滤波前信号的波形10.5
0-0.5-10
2000
4000
6000
8000
10.5
0-0.5-10
滤波后的信号波形
[**************]0
图4 滤波前后信号波形对比
分析：加入噪声后回放的声音与原始的语音信号有明显的不同, 其伴随较尖锐的干扰啸叫声。从含噪语音信号的频谱图中可以看出是含噪的语音信号的频谱, 在整个频域范围内分布均匀, 这正是干扰所造成的。通过滤波前后的对比，低通滤波后效果很好。由此可见，语音信号主要分布在低频，而噪声主要分布在高频。 6） 语音信号的变声实验
变声程序如下： clc;
deltaf=150;
deltaw=2*pi*deltaf/Ft;%频移大小 n=0:x_length-1;
xk1=2*cos(deltaw*n).*x';
yk1=fftshift(fft(xk1,x_length));%计算语音离散信号x 的yk1_length点的DFT ，并重新排列, 把DFT 移位到对称位置； subplot(411); plot(f,abs(yk1));
xlabel('频率'); ylabel('幅度'); grid on;
subplot(412); plot(angle(yk1)); xlabel('频率'); ylabel('相位');
axis([0 8000 -5 5]); subplot(413); plot(x);
xlabel('原始语音信号'); subplot(414); plot(xk1);
xlabel('频移后语音信号'); sound(xk1,fs); grid;
程序运行结果如下图5：
频移后的频谱
400200
0-400010
-1
-3000-2000-1000
0频率
[**************]0
010002000
[1**********]0
原始语音信号
[1**********]0
20
-2
010002000
[1**********]0频移后语音信号
[1**********]0
图5 频移后的频谱图和语音图
分析：程序中利用余弦函数对原语音信号进行了频谱的搬移，由图可以看到频移大小为150Hz ，程序运行的结果可以明显听到声音发生了改变，但是没有实现男声变女声的结果。
三、心得体会
我感受到只有在了解课本知识的前提下，才能更好的应用这个工具；并且熟练的应用MATLAB 也可以很好的加深我对课程的理解，为我带来方便。这次设计使我了解了MATLAB 的使用方法，学会分析滤波器的优劣和性能，提高了分析和动手实践能力。同时我相信，进一步加强对MATLAB 的学习与研究对我今后的学习将会起到很大的帮助。
通过这次作业使我懂得了理论与实际相结合是很重要的，只有理论知识是远远不够的，只有把所学的理论知识与实践相结合起来，从理论中得出结论，才能真正为社会服务，从而提高自己的实际动手能力和独立思考的能力。在设计的过程中遇到问题，可以说得是困难重重，这毕竟第一次做的，难免会遇到过各种各样的问题，同时在设计的过程中发现了自己的不足之处，对以前所学过的知识理解得不够深刻，掌握得不够牢固。
总的来说，在写作业的过程中遇到了很多问题，在同学和书本的帮助下，终于一一解决了，有点小小的成就感，终于觉得平时所学的知识有了实用的价值，达到了理论与实际相结合的目的，不仅学到了不少知识，而且锻炼了自己的能力，使自己对以后的路有了更加清楚的认识，同时，对未来有了更多的信心。