Matlab处理音频信号一、问题的提出：数字语音是信号的一种，我们处理数字语音信号，也就是对一种信号的处理，那信号是什么呢？信号是传递信息的函数。
一、问题的提出：数字语音是信号的一种，我们处理数字语音信号，也就是对一种信号的处理，那信号是什么呢？信号是传递信息的函数。
离散时间信号%26mdash;%26mdash;序列%26mdash;%26mdash;可以用图形来表示。
按信号特点的不同，信号可表示成一个或几个独立变量的函数。
例如，图像信号就是空间位置（二元变量）的亮度函数。
一维变量可以是时间，也可以是其他参量，习惯上将其看成时间。
信号有以下几种：（1）连续时间信号：在连续时间范围内定义的信号，但信号的幅值可以是连续数值，也可以是离散数值。
当幅值为连续这一特点情况下又常称为模拟信号。
实际上连续时间信号与模拟信号常常通用，用以说明同一信号。
（2）离时间信号：时间为离散变量的信号，即独立变量时间被量化了。
而幅度仍是连续变化的。
（3）数字信号：时间离散而幅度量化的信号。
语音信号是基于时间轴上的一维数字信号，在这里主要是对语音信号进行频域上的分析。
在信号分析中，频域往往包含了更多的信息。
对于频域来说，大概有8种波形可以让我们分析：矩形方波，锯齿波，梯形波，临界阻尼指数脉冲波形，三角波，余旋波，余旋平方波，高斯波。
对于各种波形，我们都可以用一种方法来分析，就是傅立叶变换：将时域的波形转化到频域来分析。
于是，本课题就从频域的角度对信号进行分析，并通过分析频谱来设计出合适的滤波器。
当然，这些过程的实现都是在MATLAB软件上进行的，MATLAB软件在数字信号处理上发挥了相当大的优势。
二、设计方案：利用MATLAB中的wavread命令来读入（采集）语音信号，将它赋值给某一向量。
再将该向量看作一个普通的信号，对其进行FFT变换实现频谱分析，再依据实际情况对它进行滤波。
对于波形图与频谱图（包括滤波前后的对比图）都可以用MATLAB画出。
我们还可以通过sound命令来对语音信号进行回放，以便在听觉上来感受声音的变化。
选择设计此方案，是对数字信号处理的一次实践。
在数字信号处理的课程学习过程中，我们过多的是理论学习，几乎没有进行实践方面的运用。
这个课题正好是对数字语音处理的一次有利实践，而且语音处理也可以说是信号处理在实际应用中很大众化的一方面。
这个方案用到的软件也是在数字信号处理中非常通用的一个软件%26mdash;%26mdash;MATLAB软件。
所以这个课题的设计过程也是一次数字信号处理在MATLAB中应用的学习过程。
课题用到了较多的MATLAB语句，而由于课题研究范围所限，真正与数字信号有关的命令函数却并不多。
三、主体部分：（一）、语音的录入与打开：[y,fs,bits]=wavread('Blip',[N1 N2]);用于读取语音，采样值放在向量y中，fs表示采样频率(Hz)，bits表示采样位数。
[N1 N2]表示读取从N1点到N2点的值（若只有一个N的点则表示读取前N点的采样值）。
sound(x,fs,bits); 用于对声音的回放。
向量y则就代表了一个信号（也即一个复杂的%26ldquo;函数表达式%26rdquo;）也就是说可以像处理一个信号表达式一样处理这个声音信号。
FFT的MATLAB实现在MATLAB的信号处理工具箱中函数FFT和IFFT用于快速傅立叶变换和逆变换。
下面介绍这些函数。
函数FFT用于序列快速傅立叶变换。
函数的一种调用格式为 y=fft(x)其中，x是序列，y是序列的FFT，x可以为一向量或矩阵，若x为一向量，y是x 的FFT。
且和x相同长度。
若x为一矩阵，则y是对矩阵的每一列向量进行FFT。
如果x长度是2的幂次方，函数fft执行高速基－2FFT算法；否则fft执行一种混合基的离散傅立叶变换算法，计算速度较慢。
函数FFT的另一种调用格式为 y=fft(x,N)式中，x，y意义同前，N为正整数。
函数执行N点的FFT。
若x为向量且长度小于N，则函数将x补零至长度N。
若向量x的长度大于N，则函数截短x使之长度为N。
若x 为矩阵，按相同方法对x 进行处理。
经函数fft求得的序列y一般是复序列，通常要求其幅值和相位。
MATLAB提供求复数的幅值和相位函数：abs，angle，这些函数一般和FFT同时使用。
函数abs(x)用于计算复向量x的幅值，函数angle(x)用于计算复向量的相角，介于和之间，以弧度表示。
函数unwrap(p)用于展开弧度相位角p ,当相位角绝对变化超过时，函数把它扩展至。
用MATLAB工具箱函数fft进行频谱分析时需注意：（1）函数fft返回值y的数据结构对称性若已知序列x=[4,3,2,6,7,8,9,0],求X(k)=DFT[x(n)]。
利用函数fft计算，用MATLAB编程如下：N=8;n=0:N-1;xn=[4 3 2 6 7 8 9 0]';XK=fft(xn)结果为：XK =39.0000-10.7782 + 6.2929i0 - 5.0000i4.7782 - 7.7071i5.00004.7782 + 7.7071i0 + 5.0000i-10.7782 - 6.2929i由程序运行所得结果可见，X(k)和x(n)的维数相同，共有8个元素。
X(k)的第一行元素对应频率值为0，第五行元素对应频率值为Nyquist频率，即标准频率为1.因此第一行至第五行对应的标准频率为0～1。
而第五行至第八行对应的是负频率，其X(k)值是以Nyquist频率为轴对称。
（注：通常表示为Nyquist频率外扩展，标以正值。
）一般而言，对于N点的x(n)序列的FFT是N点的复数序列，其点n=N/2+1对应Nyquist频率，作频谱分析时仅取序列X(k)的前一半，即前N/2点即可。
X(k)的后一半序列和前一半序列时对称的。
（2）频率计算若N点序列x(n)(n=0,1,…,N-1)是在采样频率下获得的。
它的FFT也是N点序列，即X(k)(k=0,1,2,…,N-1)，则第k点所对应实际频率值为f=k*f /N.（3）作FFT分析时，幅值大小与FFT选择点数有关，但不影响分析结果。
2、设计内容：（1）下面的一段程序是语音信号在MATLAB中的最简单表现，它实现了语音的读入打开，以及绘出了语音信号的波形频谱图。
[x,fs,bits]=wavread('ding.wav',[1024 5120]);sound(x,fs,bits);X=fft(x,4096);magX=abs(X);angX=angle(X);subplot(221);plot(x);title('原始信号波形');subplot(222);plot(X); title('原始信号频谱');subplot(223);plot(magX);title('原始信号幅值');subplot(224);plot(angX);title('原始信号相位');程序运行可以听到声音，得到的图形为：（2）定点分析：已知一个语音信号，数据采样频率为100Hz，试分别绘制N＝128点DFT的幅频图和N＝1024点DFT幅频图。
编程如下：x=wavread('ding.wav');sound(x);fs=100;N=128;y=fft(x,N);magy=abs(y);f=(0:length(y)-1)'*fs/length(y); subplot(221);plot(f,magy);xlabel('频率（Hz）');ylabel('幅值'); title('N=128(a)');gridsubplot(222);plot(f(1:N/2),magy(1:N/2)); xlabel('频率（Hz）');ylabel('幅值'); title('N=128(b)');gridfs=100;N=1024;y=fft(x,N);magy=abs(y);f=(0:length(y)-1)'*fs/length(y); subplot(223);plot(f,magy);xlabel('频率（Hz）');ylabel('幅值');title('N=1024(c)');gridsubplot(224);plot(f(1:N/2),magy(1:N/2));xlabel('频率（Hz）');ylabel('幅值');title('N=1024(d)');grid运行结果如图：上图(a)、(b)为N=128点幅频谱图，(c)、(d)为N=1024点幅频谱图。
由于采样频率f =100Hz，故Nyquist频率为 50Hz。
(a)、(c)是0～100Hz频谱图，(b)、(d)是0～50Hz 频谱图。
由(a)或(c)可见，整个频谱图是以Nyquist频率为轴对称的。
因此利用fft对信号作频谱分析，只要考察0～Nyquist频率（采样频率一半）范围的幅频特性。
比较(a)和(c)或(b)和(d)可见，幅值大小与fft选用点数N有关，但只要点数N足够不影响研究结果。
从上图幅频谱可见，信号中包括15Hz和40Hz的正弦分量。
（3）若信号长度T=25.6s，即抽样后x(n)点数为T/Ts=256，所得频率分辨率为Hz，以此观察数据长度N的变化对DTFT分辨率的影响：编程如下：[x,fs,bits]=wavread('ding.wav');N=256;f=0:fs/N:fs/2-1/N;X=fft(x);X=abs(X);subplot(211)plot(f(45:60),X(45:60));gridxlabel('Hz'),ylabel('|H(ejw)|')%数据长度N扩大4倍后观察信号频谱N=N*4;f=0:fs/N:fs/2-1/N;X=fft(x);X=abs(X);subplot(212)plot(f(45*4:4*60),X(4*45:4*60));gridxlabel('Hz'),ylabel('|H(ejw)|')结果如图：（三）、滤波器设计：1、相关原理：设计数字滤波器的任务就是寻求一个因果稳定的线性时不变系统，并使系统函数H(z)具有指定的频率特性。
数字滤波器从实现的网络结构或者从单位冲激响应分类，可以分成无限长单位冲激响应（IIR）数字滤波器和有限长单位冲激响应（FIR）数字滤波器。