完整word版综合训练.docx

1、完整word版综合训练 数字信号处理课程研究报告 项目: 综合训练课题 班级: 测控142 姓名: 吉宇 学号: 160514205 学年: 20162017 常熟理工学院电气与自动化工程学院课题一：基于MATLAB的语音信号分析与处理1、课题描述录制一段自己的语音信号，并对录制的信号进行采样；画出采样后语音信号的时域波形和频谱波形；根据频率特征，分别设计IIR和FIR低通、高通、带通滤波器，并对声音信号进行滤波处理，回放声音信号，分析比较处理前后声音的变化。二、课题分析语音信号处理主要分成三个部分：语音信号的录制与采样；由于信号会有噪声，因此需要有相滤波器的设计；最后对采样信号进行滤波。三、

2、课题设计 对原始声音信号采样并画出时域波形和频域波形程序设计如下：原始声音信号采集y,fs,nbits=wavread(shengyin.wav); %读取声音sound(y,fs,nbits); %回放声音N=length(y);n=0:N-1;%长度w=2*n*pi/N;y1=fft(y);%对原始信号做傅里叶变化subplot(2,1,1);plot(n,y);title(原始语音信号时域图);xlabel(时间t);ylabel(幅值y);subplot(2,1,2);plot(w/pi,abs(y1);title(原始语音信号频谱);xlabel(频率Hz);ylabel(幅度);设

3、计合适的滤波器前首先需要设计一个原型滤波器，下面分别设计了IIR和FIR低通，带通高通滤波器，设计程序过程如下：原型滤波器设计IIR低通滤波器：Ft=8000;%模拟指标采样周期Fp=1000;%通带截止频率Fs=1200;%阻带截止频率wp=2*pi*Fp/Ft;%数字指标ws=2*pi*Fs/Ft;ft=2*Ft*tan(wp/2);%双线性变化fs=2*Ft*tan(ws/2);n11,wn11=buttord(wp,ws,1,50,s);b11,a11=butter(n11,wn11,s);num11,den11=bilinear(b11,a11,0.5);h,w=freqz(num1

4、1,den11);figure;plot(w*8000*0.5/pi,abs(h);legend(IIR低通滤波器);grid;IIR带通滤波器Fp1=1200;Fp2=3000;Fs1=1000;Fs2=3200;Ft=8000;wp1=tan(pi*Fp1/Ft);wp2=tan(pi*Fp2/Ft);ws1=tan(pi*Fs1/Ft);ws2=tan(pi*Fs2/Ft);w=wp1*wp2/ws2;bw=wp2-wp1;wp=1;ws=(wp1*wp2-w.2)/(bw*w);n12,wn12=buttord(wp,ws,1,50,s);b12,a12=butter(n12,wn12

5、,s);num2,den2=lp2bp(b12,a12,sqrt(wp1*wp2),bw);num12,den12=bilinear(num2,den2,0.5);h,w=freqz(num12,den12);figure;plot(w*8000*0.5/pi,abs(h);axis(0 4500 0 1.5);legend(IIR带通滤波器,Location,NorthWest);grid;IIR高通滤波器：Ft=8000;Fp=4000;Fs=3500;wp1=tan(pi*Fp/Ft);ws1=tan(pi*Fs/Ft);wp=1;ws=wp1*wp/ws1;n13,wn13=cheb1

6、ord(wp,ws,1,50,s);b13,a13=cheby1(n13,1,wn13,s);num,den=lp2hp(b13,a13,wn13);num13,den13=bilinear(num,den,0.5);h,w=freqz(num13,den13);figure;plot(w*21000*0.5/pi,abs(h);legend(IIR高通滤波器,Location,NorthWest);axis(0 11000 0 1.5);grid;FIR低通滤波器：Ft=8000;Fp=1000;Fs=1200;wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.

7、(-rp/20);s=10.(-rs/20);fpts=wp ws;mag=1 0;dev=p s;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);b21=fir1(n21,wn21,kaiser(n21+1,beta);h,w=freqz(b21,1);figure;plot(w*8000*0.5/pi,abs(h);title(FIR低通滤波器,fontweight,bold);grid;FIR带通滤波器Fp1=1200;Fp2=3000;Fs1=1000;Fs2=3200;Ft=8000;wp1=tan(pi*Fp1/Ft);wp2=tan(pi*

8、Fp2/Ft);ws1=tan(pi*Fs1/Ft);ws2=tan(pi*Fs2/Ft);w=wp1*wp2/ws2;bw=wp2-wp1;wp=1;ws=(wp*wp2-w.2)/(bw*w);n22,wn22=buttord(wp,ws,1,50,s);b22,a22=butter(n22,wn22,s);num2,den2=lp2bp(b22,a22,sqrt(wp1*wp2),bw);num22,den22=bilinear(num2,den2,0.5);h,w=freqz(num22,den22);figure;plot(w*8000*0.5/pi,abs(h);axis(0 45

9、00 0 1.5);legend(FIR带通滤波器,Locating,NorthWest);grid;FIR高通滤波器：Ft=8001;Fp=4000;Fs=3500;wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.(-rp/20);s=10.(-rs/20);fpts=ws,wp;mag=0,1;dev=p,s;n23,wn23,beta,ftype=kaiserord(fpts,mag,dev);b23=fir1(n23,wn23,high,kaiser(n23+1,beta);h,w=freqz(b23,1);figure;plot(w*12000*0.5

10、/pi,abs(h);title(fir高通滤波器);axis(2500,5500,0,1.2);grid;对声音信号分别用所设计的滤波器进行滤波，分析滤波前后时域和频域波形的不同，设计程序如下：对声音信号滤波处理后IIR双线性低通滤波：y,fs,nbits=wavread(shengyin.wav);n=length(y);S=fft(y);Ft=8000;Fp=1000;Fs=1200;wp=2*pi*Fp/Ft;ws=2*pi*Fs/Ft;n11,wn11=buttord(wp,ws,1,50,s);b11,a11=butter(n11,wn11,s);num11,den11=bilin

11、ear(b11,a11,0.5);z11=filter(num11,den11,s);sound(z11,fs);m11=fft(z11);figure;subplot(2,2,1);plot(abs(S),g);title(滤波前信号的频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,2);plot(abs(m11),r);title(滤波后信号的频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,3);plot(s);title(滤波前信号的波形,fontwei

12、ght,bold);axis(67000 87000 -0.5 0.5);grid;subplot(2,2,4);plot(z11);title(滤波后信号的波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;IIR双线性带通滤波：y,fs,nbits=wavread(shengyin.wav);n=length(y);S=fft(y);Ft=20000;Fp=700;Fs=1400;wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.(-rp/20);q=10.(-rs/20);fpts=wp ws;mag=1

13、0;dev=p q;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);b21=fir1(n21,wn21,kaiser(n21+1,beta);z21=fftfilt(b21,s);sound(z21,fs);m21=fft(z21);figure(4);subplot(2,2,1);plot(abs(S),g);title(滤波前信号的频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,2);plot(abs(m11),r);title(滤波后信号的频谱,fontweight,bold);a

14、xis(0 15000 0 1000);grid;subplot(2,2,3);plot(s);title(滤波前信号的波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;subplot(2,2,4);plot(z11);title(滤波后信号的波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;IIR双线性高通滤波：y,fs,nbits=wavread(shengyin.wav);n=length(y);S=fft(y);Fp1=1400;Fs1=700;Ft=10000;wp1=tan(

15、pi*Fp1/Ft);ws1=tan(pi*Fs1/Ft);wp=1;ws=wp*wp/ws1;n13,wn13=cheb1ord(wp,ws,1,50,s);b13,a13=cheby1(n13,1,wn13,s);num,den=lp2hp(b13,a13,wn13);num13,den13=bilinear(num,den,0.5);z13=filter(num13,den13,s);sound(z13,fs);m13=fft(z13);figure;subplot(2,2,1);plot(abs(S),g);title(滤波前信号的频谱,fontweight,bold);axis(0

16、15000 0 1000);grid;subplot(2,2,2);plot(abs(m11),r);title(滤波后信号的频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,3);plot(s);title(滤波前信号的波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;subplot(2,2,4);plot(z11);title(滤波后信号的波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;FIR窗函数低通滤波：y,fs

17、,nbits=wavread(shengyin.wav);n=length(y)S=fft(y);Ft=10000;Fp=700;Fs=1400;wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.(-rp/20);q=10.(-rs/20);fpts=wp ws;mag=1 0;dev=p q;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);b21=fir1(n21,wn21,kaiser(n21+1,beta);z21=fftfilt(b21,s);sound(z21,fs);m21=fft(z21);figure(4);

18、subplot(2,2,1);plot(abs(S),g);title(滤波前的信号频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,2);plot(abs(m21),r);title(滤波后的信号频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,3);plot(s);title(滤波前的信号波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;subplot(2,2,4);plot(z21);title(滤波后的

19、信号波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;FIR窗函数带通滤波：y,fs,nbits=wavread(shengyin.wav);n=length(y)S=fft(y);Fp1=1200;Fp2=3000;Fs1=1000;Fs2=3200;Ft=2200;wp1=tan(pi*Fp1/Ft);wp2=tan(pi*Fp2/Ft);ws11=tan(pi*Fs1/Ft);ws2=tan(pi*Fs2/Ft);w=wp1*wp2/ws2;bw=wp2-wp1;wp=1;ws=(wp*wp2-w.2)/(bw*w);n22,wn22

20、=buttord(wp,ws,1,50,s);b22,a22=butter(n22,wn22,s);z22=fftfilt(b22,s);sound(z22,fs);m22=fft(z22);figure;subplot(2,2,1);plot(abs(S),g);title(滤波前的信号频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,2);plot(abs(m22),r);title(滤波后的信号频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,3);plot

21、(s);title(滤波前的信号波形,fontweight,bold);axis(6700 8700 -0.5 0.5);grid;subplot(2,2,4);plot(z22);title(滤波后的信号波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;FIR窗函数高通滤波：y,fs,nbits=wavread(shengyin.wav);n=length(y)S=fft(y);Ft=10000;Fp=1400;Fs=700;wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.(-rp/20);q=10.(-r

22、s/20);fpts=ws wp;mag=1 0;dev=p q;n23,wn23,beta,ftype=kaiserord(fpts,mag,dev);b23=fir1(n23,wn23,high,kaiser(n23+1,beta);z23=fftfilt(b23,s);sound(z23,fs);m23=fft(z23);figure;subplot(2,2,1);plot(abs(S),g);title(滤波前的信号频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,2);plot(abs(m23),r);title(滤波

23、后的信号频谱,fontweight,bold);axis(0 15000 0 1000);grid;subplot(2,2,3);plot(s);title(滤波前的信号波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;subplot(2,2,4);plot(z23);title(滤波后的信号波形,fontweight,bold);axis(67000 87000 -0.5 0.5);grid;课题总结： 本次的数字信号处理综合实验的题目是应用Matlab对语音信号进行频谱分析及滤波，首先通过网络和书籍查找有关本次综合实验的材料，编写相关程

24、序，并通过Matlab软件运行得到相关波形频谱图。 实验中利用双线性变换法设计IIR数字滤波器，利用窗函数设计FIR数字滤波器，可以是低通、高通和带通滤波器的设计，通过设计的滤波器对语音信号进行滤波，再对得出的频谱图进行分析。在实验中得到一些困难，在设计数字滤波器的时候，通带和阻带频率的选取需要满足低通的要求，以及带通允许的最大衰减和阻带应达到最小的衰减。课题二：基于MATLAB的谐波分析FFT本文应用MATLAB来验证定理：方波可用相应频率的基波及其奇次谐波合成。一、公式分析及计算 1.1傅里叶变换的原理 任何具有性质周期为T的波函数()ft都可以表示为三角函数所构成的级数之和，即： （1）

25、 其中：t为时间，为角频率。=2/T（T为周期），第一项1/2 (a0)为直流分量。图1方波所谓周期性函数的傅里叶变换（Fouriertransform）就是将周期性函数张凯成直流分量，基波和所有n次谐波的叠加。 图1所示的方波可以写成函数形式：f(t)= h, (0tT/2) -h,(-T/2t0)在这里，h为常数2。很明显，此方波为奇函数，并且它没有常数项，同时，它是一个周期为T的函数，所以我们可以用傅里叶级数来表示这个函数。我们把它展开，可以得到：1.2傅里叶变换的证明下面，我们要从数学角度来证明为什么公式（3）能成立。由于这是一个奇函数,常数项0a可以用积分函数计算出来：所以其常数项不

26、存在，即a0=0，下面开始计算an与bn:1.3周期信号的分解周期信号是定义在（-,+）区间，每隔一定的时间T，按相同规律重复变化的信号，它可表示为: f(t)=f(t+mT)式中m为任意整数。时间T称为该信号的重复周期，简称为周期。 需要指出的是，只有当周期信号满足狄里赫利条件时，才能展开为傅里叶级数。狄里赫利条件是：1）函数在任意有限区间内连续，或只有有限个第一类间断点（当t从左或右趋于这个间断点时，函数有有限的左极限和右极限）2）在一周期内，函数有有限个极大值或极小值。上式中系数an,bn称为傅里叶系数。为简便式，积分区间（t0,t0+T）取为（ -2T ，2 T ）或（0，T）。考虑到

27、正、余弦函数的正交条件，可得傅里叶系数。以下是公式（2-2），公式（2-3）式中T为函数f(t)的周期，=2/T为角频率，由上述两式，傅里叶系数an和bn 都是n的函数，其中an是n的偶函数，即a(-n)=an；而bn是n的奇函数，既有b(-n)=-bn将式（2-1）中同频率项合并，可写成如下形式:如将式（2-4）的形式化为（2-1）的形式，他们系数之间的关系为:式（2-4）表明，任何满足狄里赫利条件的周期函数可分解为直流和许多余弦（或正 弦）分量。其中第一项 A/2是常数项，它是周期信号中所包涵的直流分量；式中第二项A1cos(1+1）称为基波或一次谐波，它的角频率与原周期信号相同，A1是基

28、波振幅，1 是基波初相角；式中第三项A2cos(2+1）称为二次谐波，它的频率是基波频率的两倍，A2是二次谐波振幅，2 是其初相角。以此类推，还有三次、四次、谐波。一般而言，Ancos(n+1）称为n次谐波，An是n次谐波的振幅，n 是其初相角。式（2-4）表明，周期函数可以分解为各谐波分量.1.4方波的分解 设方波信号f(t)的周期为T ，宽度为2T ，将其展开为傅里叶级数 由式（2-2）和（2-3）可得将它们代入到式（2-1） ，得到信号的傅里叶级数展开式为它只含一、三、五奇次谐波分量。 下图中画出了一个周期的方波组成情况，由图可见，当它包含的谐波分量愈多时，波形就愈接近原来的方波信号f(t)（图中虚线所示），其均方误差愈小，还可以看出，频率较低的谐波，其振幅较大，他们组成方波的主体，而频率较高的高次谐波振幅较小，它们主要影响波形的细节，波形中所包含的高次谐