MATLAB与信息系统仿真课程设计报告.doc

matlab与信息系统仿真课程设计报告课题：语音信号的采集与处理课程：信号与系统、数字信号处理学院：xx学院班级：组员：指导老师:完成日期：2015年1月16日1目录1、 引言.22、 正文.22.1 设计要求与设备.22.2 语音信号的频谱分析.22.2.1 设计任务.22.2.2 设计原理.22.2.3 设计结果.32.3 语音信号的采样.62.3.1 设计任务.62.3.2 设计原理.72.3.3 设计结果.72.4 语音信号的转变.82.4.1 设计任务.82.4.2 设计原理.82.4.3 设计结果.82.5 语音信号的调制与解调.102.5.1 设计任务.102.5.2 设计原理.102.5.3 设计结果.102.6 加噪语音信号的滤波.142.6.1 设计任务.142.6.2 设计原理.142.6.3 设计结果.143、 心得体会.184、 参考文献.20附录：程序代码.21程序代码:program1.21程序代码:program2.24程序代码:program3.25程序代码:program4.27程序代码:program5.281、 引言语音信号采集与处理是用数字信号处理技术对语音信号进行处理的方式。它的目的是想通过处理得到一些反应语音信号重要特征的语音参数以辨别个体语音的区别，再者是通过使用matlab实验软件来对语音信号进行各种变化以便更好地掌握快速傅立叶变换（fft）、采样定理（奈奎斯特采样定律）、调制与解调、各种滤波器的原理等。因此，在语音信号处理之前我们要先了解语音信号的特征，以便进一步的实验验证，还要掌握数字信号处理技术的应用。对语音信号进行时域与频域的分析，正确地应用快速傅立叶变换的方法，但我们平常所接触到的大多数是一些结构单一、周期型的函数信号，例如正余弦函数信号、冲激函数和阶跃函数，而对于复杂多变的语音信号，我们该如何驾驭?这便是我们这次的课程设计需要思考的内容。2、 正文2.1 设计要求与设备综合设计：语音信号的采集与处理1. 利用 windows 自带的录音机或者其它录音软件，录制几段语音信号（要有几种不同的声音，要有男声、女声）。 2. 对录制的语音信号进行频谱分析，确定该段语音的主要频率范围，由此频率范围判断该段语音信号的特点（低沉 or 尖锐？） 3. 利用采样定理，对该段语音信号进行采样，观察不同采样频率（过采样、欠采样、临界采样）对信号的影响。 这里涉及到采样下重采样的问题，请大家思考如何去做。 4. 实现语音信号的快放、慢放、倒放、回声、男女变声。 5. 对采集到的语音信号进行调制与解调，观测调制与解调前后信号的变化。 6. 对语音信号加噪，然后进行滤波，分析不同的滤波方式对信号的影响。 设计设备清单：1 pc机一台2 matlab软件2.2 语音信号的频谱分析2.2.1 设计任务对录制的语音信号进行频谱分析，确定该段语音的主要频率范围，由此频率范围判断该段语音信号的特点。（低沉 or 尖锐？）2.2.2 设计原理使用matlab提供的函数来实现fft最为简捷。matlab提供内部函数fft来快速计算序列的离散傅立叶变换，其格式如下：式中，x是时域序列；n是点数；x是变换所得的频域序列。如果x的长度小于n，则在其后补零，使之成为长度为n的序列。也可以省略n，这时采用的就是x的长度。如果x是一个矩阵，则计算x中每列的n点dft。函数fft是matlab的内建函数，是用机器语言而不是用matlab的指令写成的。因此不能在命令行输入type命令来查看程序代码。函数fft给出的频域序列以坐标原点为始点。有时候，希望以坐标原点作为频域序列的中点。为此，可用函数fftshift来计算离散傅立叶变换。通过采用快速傅立叶变换（fft）函数来得到语音信号的频谱可直接看出该段语音的主要频率范围，通过男女声频谱可直接做出比较，得到男女语音信号的特点：在低频部分，在男声语音信号相对于女声语音信号所占的比例大；在高频部分，在女声语音信号相对于男声语音信号所占的比例大；为了能够让阅读者更加直观的发现男女语音信号在高低频率的区别，我们计算了男女语音信号各自在高低频段所占比率，以图示的形式展示。此频谱分析中包含语音信号的时域分析、频域分析、降采样和男女信号各自在高低频段所占比例的对比。2.2.3 设计结果(见附录：程序代码program1)1) 男女语音信号的时域分析图：2) 男女语音信号的频域分析图：3) 降采样：4) 男女信号各自在高低频段所占比例的对比：2.3 语音信号的采样2.3.1 设计任务利用采样定理，对该段语音信号进行采样，观察不同采样频率（过采样、欠采样、临界采样）对信号的影响。 （这里涉及到采样下重采样的问题，请大家思考如何去做。 ）2.3.2 设计原理数字信号处理的对象，是在采样时钟的控制下，通过ad转换器以一定的采样率对模拟信号进行采样而得到的。根据采样定理，采样率必须大于模拟信号最高频率的2倍（奈奎斯特采样定律）。然而在许多情况下，要求对信号以不同的频率进行采样。这时需要对采样数据进行处理。在matlab中可以利用抽取或内插方法改变采样率。其中，resample函数是以p/q倍速度对向量x进行重采样，该函数调用的格式如下：y=resample(x,p,q)；或者采用downsample函数，它的功能是对时间序列重采样，在原时间序列中等间隔地取出一些项，得到新的序列，格式为y=downsample(x,n)，指从第一项开始，等间隔n对信号x采样，得到新序列y。这里改变间隔n的大小来观察新序列的频谱是否出现失真，若频谱失真，则处于欠采样（采样频率低于信号带宽的2倍），在即将失真前的那一个间隔i为临界采样（采样频率等于信号带宽的2倍），小于i时为过采样（采样频率高于信号带宽的2倍），但语音信号本身在录入电脑时已处于过采样阶段，此处可省略过采样。2.3.3 设计结果(见附录：程序代码program2)1) 原始语音信号和采样后的语音信号的时域对比2) 原始语音信号和采样后的语音信号的频域对比2.4 语音信号的转变2.4.1 设计任务 实现语音信号的快放、慢放、倒放、回声、男女变声。 2.4.2 设计原理变速：改变语音的播放速度也就是改变采样间隔（即改变了采样频率），但是这个频实现率依然要在2f(nyquist rate)之上，否则就会产生失真。倒放：使用flipud()倒置语音矩阵，逆序输出音频即可。回声：回声在时域上幅值减小了，频域上的特征不变，只需要把原信号添加一个延时（delay）和对时域的幅度添加一个参数.然后和原信号叠加即可获得回响效果。男女变声：我们使用另外写的voice(x,f)函数实现音调转换，x为需要转换的声音，通过抽取插值更改采样率来改变基频，当f1时音调降低;f1音调升高。然后再进行时长整合使语音文件恢复原来的时长。时长整合使用重叠叠加算法来实现。经过我们试听，转换效果还是很好的。2.4.3 设计结果(见附录：程序代码program3)1) 加入回声的波形2) 转换前与转换后语音信号波形的对比3) 转换前与转换后语音信号频谱的对比2.5 语音信号的调制与解调2.5.1 设计任务采集到的语音信号进行调制与解调，观测调制与解调前后信号的变化。2.5.2 设计原理调制与解调是互为相反功能的信号频谱搬移过程。在发送端，通过调制将传输信号频谱搬移到指定传输信道的频段上，以便于传输、信道复用以及干扰抑制；在接收端，再以相反的过程（解调）将传输信号恢复出来。 通带模拟调制解调的模型：该设计中我们采用的是线性调制：双边带调制（dsb），设输入被调信号为m(t)，则调制输出y(t)为，其中为载波频率（hz）；是载波初始相位（rad），此处，我们将取0。双边带调制解调的优点是功率高效，节省了载波功率，调制效率100%，缺点是频带利用并不高效，多传了一倍的频带信息，且不能用包络检波，需用相干检波，较复杂。 dsb调制解调的框图： 2.5.3 设计结果(见附录：程序代码program4)1) 语音信号的时域波形2) 语音信号的频谱3) 载波信号的频谱4) 已调信号的频谱5) 解调信号的频谱6) 低通滤波器的特性曲线7) 滤波后的语音信号的时域图和频谱2.6 加噪语音信号的滤波2.6.1 设计任务对语音信号加噪，然后进行滤波，分析不同的滤波方式对信号的影响。2.6.2 设计原理该设计加入的噪声是noise=0.1*sin(2*pi*fn*t)（fn=3000）可供更改噪声频率。此处我们简化了噪声的复杂度，没有将噪声信号noise=randn(size(t)随机化,是降低了滤波的难度，但基本原理是一致的。这里我们设计了3种滤波器对语音进行处理，分别为椭圆数字型低通滤波、椭圆数字型高通滤波和巴特沃斯型数字低通滤波器，绘制出相应的幅度、相位谱图，滤波后的波形、频谱图。椭圆数字型低通滤波：调用格式:n,wn=ellipord(wp,ws,rp,rs)b,a=ellip(n,rp,rs,wp);椭圆数字型高通滤波：调用格式:n,wn=ellipord(wp,ws,rp,rs);b,a=ellip(n,rp,rs,wp,high);巴特沃斯型数字低通滤波器：调用格式:n,wn=buttord(wp,ws,rp,rs);b,a=butter(n,wn,low);2.6.3 设计结果(见附录：程序代码program5)1) 原始语音信号与加噪信号时域图的对比2) 原始语音信号与加噪信号频谱的对比3) 椭圆低通滤波器的特性曲线4) 椭圆低通滤波后的语音信号时域图和频谱5) 椭圆高通滤波器的特性曲线6) 椭圆高通滤波后的语音信号时域图和频谱7) 巴特沃斯型数字低通滤波器的特性曲线8) 巴特沃斯型数字低通滤波器后的语音信号时域图和频谱3、 心得体会(自我感言)4、 参考文献1 数字信号处理：理论与实践/王大伦，王志新，王康编著.北京：清华大学出版社，2010.22 通信原理/樊昌信，曹丽娜编著.7版.北京：国防工业出版社，2012.113 matlab信号处理/刘波，文忠，曾涯编著.北京：电子工业出版社，2006.14 matlab/simulink通信系统建模与仿真实例分析/邵玉斌编著.北京：清华大学出版社，2008.65 数字信号处理/张小虹主编.2版.北京：机械工业出版社，2008.86 信号与系统分析/和卫星，许波主编.西安：西安电子科技大学出版社，2007.3附录：程序代码（各个程序代码独立）1 程序代码：program1%频谱分析clear;%fs1=fs2=44100hzx1,fs1,bits =wavread(d:yllvoice女声.wav) ;x2,fs2,bits =wavread(d:yllvoice男声.wav) ;%sound(x1,fs1,bits); %对语音信号进行回放%sound(x2,fs2,bits); %对语音信号进行回放m1=length(x1);m2=length(x2); %语音信号的长度n1=round(log2(m1);n2=round(log2(m2);n1=2n1;n2=2n2;t1=0:1/fs1:(m1-1)/fs1;t2=0:1/fs2:(m2-1)/fs2;%时域频谱图 figure(1)subplot(211);plot(t1,x1);xlabel(时间(t1);ylabel(幅值); %原始语音信号的时域图形title(女声信号x1n);subplot(212);plot(t2,x2);xlabel(时间(t2);ylabel(幅值); %原始语音信号的时域图形title(男声信号x2n);figure(2) xf1=fft(x1,n1);f1=0:n1-1*fs1/n1;amp1=abs(xf1)/(n1/2);amp1(1,:)=amp1(1,:)/2;subplot(211);plot(f1,amp1);xlabel(频率（f1）);ylabel(|x1(f)|);%原始语音信号的fft频谱图 title(女声信号fft频谱);%axis(0 5000 0 0.003);%可加可不加，方便从频谱中看出男女声的差别xf2=fft(x2,n2);f2=0:n2-1*fs2/n2;amp2=abs(xf2)/(n2/2);amp2(1,:)=amp2(1,:)/2;subplot(212);plot(f2,amp2);xlabel(频率（f2）);ylabel(|x2(f)|); %原始语音信号的fft频谱图 title(男声信号fft频谱);%axis(0 5000 0 0.003);%降采样figure(3)%采样频率为原来的1/2时的频谱图x11=x1(1:2:length(x1);x21=x2(1:2:length(x2); %降采样m11=length(x11);m21=length(x21); %语音信号的长度n11=round(log2(m11);n21=round(log2(m21);n11=2n11;n21=2n21;xf11=fft(x11,n11);f11=0:n11-1*fs1/n11;amp11=abs(xf11)/(n11/2);amp11(1,:)=amp11(1,:)/2;subplot(211);plot(f11,amp11);xlabel(频率（f1）);ylabel(|x11(f)|); %原始语音信号的fft频谱图 title(1/2时的女声信号fft频谱);xf21=fft(x21,n21);f21=0:n21-1*fs2/n21;amp21=abs(xf21)/(n21/2);amp21(1,:)=amp21(1,:)/2;subplot(212);plot(f21,amp21);xlabel(频率（f2）);ylabel(|x21(f)|); %原始语音信号的fft频谱图 title(1/2时的男声信号fft频谱);figure(4)%采样频率为原来的1/4时的频谱图x12=x1(1:4:length(x1);x22=x2(1:4:length(x2); %降采样m12=length(x12);m22=length(x22); %语音信号的长度n12=round(log2(m12);n22=round(log2(m22);n12=2n12;n22=2n22;xf12=fft(x12,n12);f12=0:n12-1*fs1/n12;amp12=abs(xf12)/(n12/2);amp12(1,:)=amp12(1,:)/2;subplot(211);plot(f12,amp12);xlabel(频率（f1）);ylabel(|x12(f)|); %原始语音信号的fft频谱图 title(1/4时的女声信号fft频谱);xf22=fft(x22,n22);f22=0:n22-1*fs2/n22;amp22=abs(xf22)/(n22/2);amp22(1,:)=amp22(1,:)/2;subplot(212);plot(f22,amp22);xlabel(频率（f2）);ylabel(|x22(f)|); %原始语音信号的fft频谱图 title(1/4时的男声信号fft频谱);figure(5)%采样频率为原来的1/8时的频谱图x13=x1(1:8:length(x1);x23=x2(1:8:length(x2); %降采样m13=length(x13);m23=length(x23); %语音信号的长度n13=round(log2(m13);n23=round(log2(m23);n13=2n13;n23=2n23;xf13=fft(x13,n13);f13=0:n13-1*fs1/n13;amp13=abs(xf13)/(n13/2);amp13(1,:)=amp13(1,:)/2;subplot(211);plot(f13,amp13);xlabel(频率（f1）);ylabel(|x13(f)|); %原始语音信号的fft频谱图 title(1/8时的女声信号fft频谱);xf23=fft(x23,n23);f23=0:n23-1*fs2/n23;amp23=abs(xf23)/(n23/2);amp23(1,:)=amp23(1,:)/2;subplot(212);plot(f23,amp23);xlabel(频率（f2）);ylabel(|x23(f)|); %原始语音信号的fft频谱图 title(1/8时的男声信号fft频谱);figure(6)%采样频率为原来的1/16时的频谱图x14=x1(1:16:length(x1);x24=x2(1:16:length(x2); %降采样m14=length(x14);m24=length(x24); %语音信号的长度n14=round(log2(m14);n24=round(log2(m24);n14=2n14;n24=2n24;xf14=fft(x14,n14);f14=0:n14-1*fs1/n14;amp14=abs(xf14)/(n14/2);amp14(1,:)=amp14(1,:)/2;subplot(211);plot(f14,amp14);xlabel(频率（f1）);ylabel(|x14(f)|); %原始语音信号的fft频谱图 title(1/16时的女声信号fft频谱);xf24=fft(x24,n24);f24=0:n24-1*fs2/n24;amp24=abs(xf24)/(n24/2);amp24(1,:)=amp24(1,:)/2;subplot(212);plot(f24,amp24);xlabel(频率（f2）);ylabel(|x24(f)|); %原始语音信号的fft频谱图 title(1/16时的男声信号fft频谱);%计算男声和女声信号在高频、低频所占的比例%低频部分的比较xf1=fft(x1,n1);xf2=fft(x2,n2);%低频信号量sum11=0;sum21=0;for i=50:500 sum11=sum11+abs(xf1(i); sum21=sum21+abs(xf2(i);end%总信号量sum1=0;sum2=0;for i1=1:length(xf1) sum1=sum1+abs(xf1(i1); endfor i2=1:length(xf2) sum2=sum2+abs(xf2(i2); end%所占比例n11=sum11/sum1;n21=sum21/sum2;%高频部分的比较%高频信号量sum12=0;sum22=0;for i3=15000:20000 sum12=sum12+abs(xf1(i3); sum22=sum22+abs(xf2(i3);end%所占比例n12=sum12/sum1;n22=sum22/sum2;%显示结果figure(7)axis(-2 2 -2 2);text(-1.5,1.6,女声的总信号量sum1=,num2str(sum1);text(-1.5,1.4,男声的总信号量sum2=,num2str(sum2);text(-1.5,1.2,女声的低频信号量sum11=,num2str(sum11);text(-1.5,1.0,男声的低频信号量sum21=,num2str(sum21);text(-1.5,0.8,女声的高频信号量sum12=,num2str(sum12);text(-1.5,0.6,男声的高频信号量sum22=,num2str(sum22);text(-1.5,0.4,在低频部分女声所占比例n11=,num2str(n11);text(-1.5,0.2,在低频部分男声所占比例n21=,num2str(n21);text(-1.5,0.0,在高频部分女声所占比例n12=,num2str(n12);text(-1.5,-0.2,在高频部分男声所占比例n22=,num2str(n22);text(-1.5,-0.8,在低频部分女声n11在高频部分男声n22);drawnow;sum1%女声的总信号量sum2%男声的总信号量sum11%女声的低频信号量sum21%男声的低频信号量sum12%女声的高频信号量sum22%男声的高频信号量n11%在低频部分女声所占比例n21%在低频部分男声所占比例n12%在高频部分女声所占比例n22%在高频部分男声所占比例2 程序代码：program2%对该段语音信号进行采样，观察不同采样频率（欠采样、临界采样）clear;x,fs,bits =wavread(d:yllvoice女声.wav) ;%sound(x,fs,bits); %对语音信号进行回放m=length(x); %语音信号的长度n=round(log2(m);n=2n;t=0:1/fs:(m-1)/fs;%临界采样i=10;%i10,属于欠采样，从频谱上看出；i1降调,1降低;f1升高f=round(f*1000);d=resample(x0,f,1000);%时长整合使语音文件恢复原来时长 w=400; wov=w/2; kmax=w*2; wsim=wov; xdecim=8; kdecim=2; x=d; f=f/1000; ss =w-wov; xpts = size(x,2);ypts = round(xpts / f); y = zeros(1, ypts); xfwin = (1:wov)/(wov+1); ovix = (1-wov):0; newix = 1:(w-wov);simix = (1:xdecim:wsim) - wsim; padx = zeros(1, wsim), x, zeros(1,kmax+w-wov);y(1:wsim) = x(1:wsim); lastxpos = 0; km = 0; for ypos = wsim:ss:(ypts-w) xpos = round(f * ypos); kmpred = km + (xpos - lastxpos); lastxpos = xpos; if (kmpred = kmax) km = kmpred; else ysim = y(ypos + simix); rxy = zeros(1, kmax+1); rxx = zeros(1, kmax+1); kmin = 0; for k = kmin:kdecim:kmax xsim = padx(wsim + xpos + k + simix); rxx(k+1) = norm(xsim); rxy(k+1) = (ysim * xsim); end rxy = (rxx = 0).*rxy./(rxx+(rxx=0); km = min(find(rxy = max(rxy)-1); end xabs = xpos+km; y(ypos+ovix) = (1-xfwin).*y(ypos+ovix) + (xfwin.*padx(wsim+xabs+ovix); y(ypos+newix) = padx(wsim+xabs+newix); endend4 程序代码：program4%调制解调双边带调制dsbclear; x,fs,bits=wavread(d:yllvoice女声.wav); n=length(x); %语音信号的长度t=0:1/fs:(n-1)/fs;f=0:n-1*fs/n; figure(1) plot(t,x);xlabel(时间(t);ylabel(幅值); title(语音信号的时域波形);figure(2)fy1=fft(x,n); %求语音信号的快速傅立叶变换amp1=abs(fy1)/(n/2);amp1(1,:)=amp1(1,:)/2; plot(f,amp1);xlabel(频率（f）);ylabel(|x(f)|);title(语音信号的频谱); fn=11000;yc=cos(2*fn*pi*t); %载波figure(3)fyc=fft(yc,n); %求载波的快速傅立叶变换amp2=abs(fyc)/(n/2);amp2(1,:)=amp2(1,:)/2; plot(f,amp2);xlabel(频率（f）);ylabel(|yc(f)|);title(载波信号的频谱);x=x(:,1); %取矩阵的单边，即单边道y3=x.*yc; %语音信号赋予载波之上figure(4)fy3 = fft(y3);amp3=abs(abs(fy3)/(n/2);amp3(1,:)=amp3(1,:)/2;plot(f,amp3);xlabel(频率（f）);ylabel(|y3(f)|); title(已调信号的频谱);%abs(abs(fy3)：取模平方%sound(y3,fs,nbits);y4=y3.*yc; %相干解调figure(5)fy4=fft(y4);amp4=abs(abs(fy4)/(n/2);amp4(1,:)=amp4(1,:)/2;plot(f,amp4);xlabel(频率（f）);ylabel(|y4(f)|);title(解调信号的频谱);%sound(y4,fs,bits);%巴特沃斯型数字低通滤波器fp5=800;fs5=1200;rp=3;rs=25;ws5=fs5/(fs/2);wp5=fp5/(fs/2);n5,wn5=buttord(wp5,ws5,rp,rs);b5,a5=butter(n5,wn5,low);y5=filter(b5,a5,y4);figure(6)freqz(b5,a5);title(低通滤波器的特性曲线);figure(7)subplot(211);plot(t,y5);xlabel(时间(t);ylabel(幅值);title(巴特沃斯型数字低通滤波器后的语音信号时域图);axis(0,t(end),min(y5),1.2*max(y5); %坐标范围 amp5=abs(fft(y5)/(n/2);amp5(1,:)=amp5(1,:)/2;subplot(212);plot(amp5);xlabel(频率（f）);ylabel(|y5(f)|);title(巴特沃斯型数字低通滤波器后的信号频谱图);%sound(y5,fs,bits); %对滤波过后的语音信号进行回放5 程序代码：program5%加噪滤波clear;%fs=44100hzx,fs,bits =wavread(d:yllvoice女声.wav) ;%sound(x,fs,bits); %对语音信号进行回放m=length(x); %语音信号的长度n=round(log2(m);n=2n;t=0:1/fs:(m-1)/fs;x=x(:,1);x=x; %取矩阵的单边f=0:n-1/n*fs;fn=3000;y=x+0.1*sin(2*pi*fn*t); %加噪信号，噪声为随机白噪声 %sound(y,fs,bits); %对加噪语音信号进行回放%加噪信号与原始信号的时域图和频域图对比figure(1) subplot(211);plot(t,x);xlabel(时间(t);ylabel(幅值); %原始语音信号的时域图形title(原始语音信号xn);subplot(212);plot(t,y);xlabel(时间(t);ylabel(幅值); %加噪语音信号的时域图形title(加噪语音信号yn); figure(2)xf=fft(x,n);amp0=abs(xf)/(n/2);amp0(1,:)=amp0(1,:)/2;subplot(211);plot(f,amp0);xlabel(频率（f）);ylabel(|x(f)|); %原始语音信号的fft频谱图 title(原始语音信号fft频谱);yf=fft(y,n);amp01=abs(yf)/(n/2);amp01(1,:)=amp01(1,:)/2;subplot(212);plot(f,amp01);xlabel(频率（f）);ylabel(|y(f)|); %加噪语音信号的fft频谱图 title(加噪语音信号fft频谱);%滤波%椭圆低通滤波器fp1=800;fs1=1500;rp=1;rs=40; %低通滤波器通带截止频率800hz和阻带截止频率1500hzwp1=2*fp1/fs; ws1=2*fs1/fs;n1,wn1=ellipord(wp1,ws1,rp,rs); %计算椭圆低通模拟滤波器的阶数和通带边界频率b1,a1=ellip(n1,rp,rs,wp1); %计算低通滤波器模拟滤波器系统函数系数y1=filter(b1,a1,y); %滤波器软件实现%低通滤波器绘图部分figure(3);freqz(b1,a1);t