MATLAB课程设计报告-语音信号的采集与处理.doc

matlab课程设计 学院： 电气信息工程学院 课题： 语音信号的采集与处理 班级： 电 信 1201 组长： 组员： 指导老师： 语音信号的采集与处理一、实验要求 （1） 本课程设计要求掌握 matlab 语言特性、数学运算和图形表示；掌握 matlab 程序设计方法在信息处理方面的一些应用；掌握在 windows 环境下，利用 matlab 进行语音信号采集与处理的基本方法。 （2） 本课程设计分为基础与综合设计题目，其中基础题目为 matlab 基本训练（必做） ；对于综合设计题，需按照任务书要求，对设计过程所进行的有关步骤进行理论分析，并对完成的设计作出评价，总结自己整个设计工作中的经验教训、收获；在课设最后一天提交课程设计报告，报告必须按照统一格式打印，装订成册。二、实验内容1、利用 windows 自带的录音机或者其它录音软件，录制几段语音信号（要有几种不同的声音，要有男声、女声） 。 2、对录制的语音信号进行频谱分析，确定该段语音的主要频率范围，由此频率范围判断该段语音信号的特点（低沉 or 尖锐？） 3、利用采样定理，对该段语音信号进行采样，观察不同采样频率（过采样、欠采样、临界采样） 对信号的影响。 这里涉及到采样下重采样的问题， 请大家思考如何去做。 4、实现语音信号的快放、慢放、倒放、回声、男女变声。 5、对采集到的语音信号进行调制与解调，观测调制与解调前后信号的变化。 6、 对语音信号加噪，然后进行滤波，分析不同的滤波方式对信号的影响。 7、 如果精力和时间允许，鼓励利用 matlab gui 制作语音信号采集与分析演示系统。三、matlab仿真任务一（陈康负责）1. 利用频谱的分析利用windows下的录音机,录制一段男生和女生的语音,存为*.wav的文件。然后在matlab软件平台下,利用函数wavread对语音信号进行采样,记住采样频率和采样点数。首先画出语音信号的时域波形;然后对语音号进行快速傅里叶变换（fft变换）,得到信号的频谱特性，分析语音信号频谱。源程序：fs=44100; x,fs,nbits =wavread(e:gqqqleo) ;derta_fs = fs/length(x);t=0:1/44100:(length(x)-1)/44100;sound(x,fs) ; %对加载的语音信号进行回放figure(1)plot(t,x) %做原始语音信号的时域图形title(原始语音信号);xlabel(time n);ylabel(fuzhi n);n=length(x); %求出语音信号的长度y1=fft(x,n) ; %傅里叶变换w=1/n*0:n-1*fs;figure(2) plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(y1)*2/n); %做原始语音信号的fft频谱图title(原始语音信号fft频谱)频谱分析对比运行结果：女生：频谱分析：从300hz到2000hz都有分布，主要集中在3001200hz，声音尖锐，是女生频谱。男生:频谱分析：从200hz到1000hz都有分布，集中在3001000hz声音低沉，是男生频谱。任务二（詹筱珊负责）利用采样定理，对该段语音信号进行采样，观察不同采样频率（过采样、欠采样、临界采样） 对信号的影响。 这里涉及到采样下重采样的问题， 请大家思考如何去做。 源程序：clear;y,fs,bits=wavread(e:gqqqwoman); %x语音数据fs采样频率bits采样点数 sound(y,fs,bits); %话音回放n = length (y) ; %求出语音信号的长度derta_fs = fs/length(y);y=fft(y,n); %傅里叶变换figure(1);subplot(2,1,1);plot(y);title(原始信号波形); subplot(2,1,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(y)*2/n); title(原始信号频谱);gridy1=decimate(y,2); %采样频率为22050hz，过采样fs1=fs/2;sound(y1,fs1,bits); %话音回放n1 = length (y1) ; %求出语音信号的长度derta_fs = fs/length(y1);y1=fft(y1,n1); %傅里叶变换figure(2);subplot(2,1,1);plot(y1);title(过采样信号波形); subplot(2,1,2); plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(y1)*2/n1);title(过采样信号频谱);gridy2=decimate(y,4); %每4个取一个，采样频率为11025hz接近临界采样fs2=fs/4;sound(y2,fs2,bits); %话音回放n2 = length (y2) ; %求出语音信号的长度derta_fs = fs/length(y2);y2=fft(y2,n2); %傅里叶变换figure(3);subplot(2,1,1);plot(y2);title(临界采样信号波形); subplot(2,1,2); plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(y2)*2/n2);title(临界采样信号频谱);gridy3=decimate(y,8); %抽取函数，频率变为原来1/8，欠采样fs2=fs/8;sound(y2,fs2,bits); %话音回放n3 = length (y3) ; %求出语音信号的长度derta_fs = fs/length(y3);y3=fft(y3,n3); %傅里叶变换figure(4);subplot(2,1,1);plot(y3);title(欠采样信号波形); subplot(2,1,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(y3)*2/n3);title(欠采样信号频谱);grid运行结果：任务三（郭前负责）实现语音信号的快放、慢放、倒放、回声、男女变声源程序：%变速-快放%fs=25600;z,fs,nbits=wavread(e:gqqqleo); %读取声音文件derta_fs = fs/length(z);x=z(:,1); %读入的y矩阵有两列,取第1列y1=voice(x,1.5); %调整voice()第2个参数转换音调,1降调,1为快放，w1降调,1降低;f1升高f=round(f*1000);d=resample(x,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); endend任务四（丁睿智负责）（1）调制：语音信号与高频正弦载波相调制，比较其频谱变化，回放信号，比较是与众语音信号变化。将调制后的信号进行解调，回放信号，比较时域中语音信号变化。所谓调制，就是在传送信号的一方将所要传送的信号附加在高频振荡上，再由天线发射出去。这里高频振荡波就是携带信号的运载工具，也叫载波。振幅调制，就是由调制信号去控制高频载波的振幅，直至随调制信号做线性变化。在线性调制系列中，最先应用的一种幅度调制是全调幅或常规调幅，简称为调幅（am）。在频域中已调波频谱是基带调制信号频谱的线性位移；在时域中，已调波包络与调制信号波形呈线性关系。设正弦载波为：c(t)=acos()式中，a为载波幅度； 为载波角频率； 为载波初始相位(假设=0).调制信号（基带信号）为m(t)。根据调制的定义，振幅调制信号（已调信号）一般可以表示为 ： 设调制信号的频谱为，则已调信号的频谱为：标准调幅波（am）产生原理： 调制信号是只来来自信源的调制信号（基带信号），这些信号可以是模拟的，亦可以是数字的。为首调制的高频振荡信号可称为载波，它可以是正弦波，亦可以是非正弦波(如周期性脉冲序列)。载波由高频信号源直接产生即可，然后经过高频功率放大器进行放大，作为调幅波的载波，调制信号由低频信号源直接产生，二者经过乘法器后即可产生双边带的调幅波。设载波信号的表达式为，调制信号的表达式为，则调幅信号的表达式为 （2） 解调 从高频已调信号中恢复出调制信号的过程称为解调(demodulation )，又称为检波(detection )。对于振幅调制信号，解调(demodulation )就是从它的幅度变化上提取调制信号的过程。解调(demodulation )是调制的逆过程。 可利用乘积型同步检波器实现振幅的解调，让已调信号与本地恢复载波信号相乘并通过低通滤波可获得解调信号。调制与解调：源程序：clear; f1,fs,nbits=wavread(e:gqqqwoman); %读入录音文件dt=1/fs; figure(1); %控制画图窗口1subplot(2,1,1); %将多个图画到一个平面，2行1列，共2个子图，定位第1个图进行操作n=length(f1); derta_fs = fs/length(f1); %设置频率间隔t=0:1/fs:(n-1)/fs; plot(t,f1); % 以t为横坐标，f1为纵坐标绘制连线图title(信息信号的时域波形); %加标题 fy1=fft(f1);w1=0:fs/(n-1):fs;subplot(2,1,2); plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(fy1)*2/n);title(信息信号的频谱); f2=cos(22000*pi*t);fy2 = fft(f2);n2=length(f2);w2=fs/n*0:n-1;figure(2);subplot(2,2,1); plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(fy2)*2/n); title(载波信号的频谱); f1=f1(:,1);f3=f1.*f2;subplot(2,2,2); fy3 = fft(f3);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(fy3)*2/n); title(已调信号的频谱);sound(f3,fs,nbits); %播放语音信号 f4=f3.*f2;subplot(2,1,2);fy4=fft(f4);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(fy4)*2/n);title(解调信号的频谱);sound(f4,fs,nbits); fp1=0; fs1=5000; as1=100;wp1=2*pi*fp1/fs; ws1=2*pi*fs1/fs; bf1=ws1-wp1;wc1=(wp1+ws1)/2;m1=ceil(as1-7.95)/(2.286*bf1)+1;n1=m1+1;beta1=0.1102*(as1-8.7); window=(kaiser(n1,beta1); % 用kaiser窗函数b1=fir1(m1,wc1/pi,window); % fir1是用窗函数法设计线性相位firdf的工具箱函数figure(3);freqz(b1,1,512); title(fir低通滤波器的频率响应); %用于求离散系统频响特性的函数freqzf4_low = filter(b1,1, f4); %一维数字滤波器plot(t,2*f4_low);title(滤波后的解调信号时域波形);sound(f4_low,fs,nbits); f5=fft(f4_low);figure(4);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(f5)*2/n);title(滤波后的解调信号频谱);任务五（吴亚明负责）对语音信号加噪，然后进行滤波，分析不同的滤波方式对信号的影响加噪源程序：y,fs,nbits=wavread (e:gqqqwoman);n = length (y) ; %求出语音信号的长度derta_fs = fs/length(y);noise=0.05*randn(n,2); %随机函数产生噪声si=y+noise; %语音信号加入噪声 sound(si,fs);figure(1);subplot(2,1,1);plot(y);title(原语音信号的时域波形);y=fft(y); %傅里叶变换subplot(2,1,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(y)*2/n); title(原语音信号的频域波形);figure(2);subplot(2,1,1);plot(si);title(加噪语音信号的时域波形);s=fft(si); %傅里叶变换subplot(2,1,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(s)*2/n); title(加噪语音信号的频域波形);低通：y,fs,nbits=wavread (e:gqqqwoman);n = length (y) ; %求出语音信号的长度derta_fs = fs/length(y);noise=0.05*randn(n,2); %随机函数产生噪声s=y+noise; %语音信号加入噪声s=fft(s); %傅里叶变换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); %s域频率响应的参数 num11,den11=bilinear(b11,a11,0.5); %利用双线性变换实现频率响应s域到z域的变换z11=filter(num11,den11,s); %直接滤波器实现sound(z11,fs,nbits);m11=fft(z11); %求滤波后的信号subplot(2,2,1);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(s)*2/n); title(滤波前信号的频谱);grid;subplot(2,2,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(m11)*2/n);title(滤波后信号的频谱);grid;subplot(2,2,3);plot(s); title(滤波前信号的波形);grid;subplot(2,2,4);plot(z11); title(滤波后的信号波形)高通y,fs,nbits=wavread (e:gqqqwoman);n = length (y) ; %求出语音信号的长度derta_fs = fs/length(y);noise=0.001*randn(n,2); %随机函数产生噪声s=y+noise; %语音信号加入噪声s=fft(s); %傅里叶变换fp1=300;fs1=1000;ft=8000;wp1=tan(pi*fp1/ft);ws1=tan(pi*fs1/ft);wp=1;ws=wp1*wp/ws1;n13,wn13=cheb1ord(wp,ws,1,50,s); %模拟的低通滤波器阶数和截止频率b13,a13=cheby1(n13,1,wn13,s); %s域的频率响应的参数num,den=lp2hp(b13,a13,wn13); %s域低通参数转为高通的num13,den13=bilinear(num,den,0.5); %利用双线性变换实现频率响应s域到z域转换z13=filter(num13,den13,s);sound(z13);m13=fft(z13); %求滤波后的信号subplot(2,2,1);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(s)*2/n); title(滤波前信号的频谱);grid;subplot(2,2,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(m13)*2/n);title(滤波后信号的频谱);grid;subplot(2,2,3);plot(s); title(滤波前信号的波形);grid;subplot(2,2,4);plot(z13); title(滤波后的信号波形);带通：y,fs,nbits=wavread (e:gqqqwoman);n = length (y) ; %求出语音信号的长度derta_fs = fs/length(y);noise=0.05*randn(n,2); %随机函数产生噪声s=y+noise; %语音信号加入噪声s=fft(s); %傅里叶变换ft=8000;fp=1000;fs=1200;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;n12,wn12,beta,ftype=kaiserord(fpts,mag,dev);%由kaiserord求滤波器的阶数和截止频率b12=fir1(n12,wn12,kaiser(n12+1,beta); %设计滤波器z12=fftfilt(b12,s);sound(z12,fs,nbits);m12=fft(z12); %求滤波后的信号subplot(2,2,1);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(s)*2/n); title(滤波前信号的频谱);grid;subplot(2,2,2);plot(-fs/2:derta_fs:fs/2-derta_fs,abs(fftshift(m12)*2/n);title(滤波后信号的频谱);grid;subplot(2,2,3);plot(s); title(滤波前信号的波形);grid;subplot(2,2,4);plot(z12); title(滤波后的信号波形);课程设计总结：陈康：为期两周的课程设计让我们对于数字信号处理以及matlab软件有了更深的了解,让我们将课堂上学到的知识融会贯通起来。同时也让我们首先真实地感受到matlab功能的强大，它几乎可以完成任何数学计算；而在整个课程设计过程中，自己也体会到了自主学习的乐趣和意义。遇到问题，自己想都不想就直接问别人，或者干脆放弃是一种解决办法；而自己先思考，实在不会再问别人然后举一反三是另外一种解决方案；而这其中的区别不仅仅是自己的态度，而且决定了你是不是真的通过课程设计学到了一些东西。遇到挫折，这是在所难免的，题目不会、资料找不到、运行之后提示错误等等这些已经是课程设计中的“家常便饭”了，但是只要我们有兴趣，肯去学，就能够搞定这些障碍，因为任何事情都是由不会到会，所以做课程设计也是锻炼自己的耐心的一个过程。 詹筱珊：本次课程设计应用matlab验证时域采样定理，主要侧重于某些理论知识的灵活运用，以及一些关键命令的掌握理解与分析，初步掌握线性系统的设计方法，培养独立的工作能力。加深理解时域采样定理的概念，掌握利用matlab分析系统频率响应的方法，计算在临界采样过采样欠采样三种不同条件下恢复信号的误差，并总结采样频率对信号恢复产生误差的影响，从而验证时域采样定理。一在上述实验当中，首先定义信号时采用了该信号的函数表达式形式。二在matlab中求连续信号的频谱，应用的是傅里叶变换，这样实际运算的仍是对连续信号的采样结果，这里给予了足够高的采样频率，把其作为连续信号来考虑。三实际中对模拟信号进行采样，需要根据最高截止频率fmax,按照采样定理的要求选择采样频率的两倍，即fs2fmax。设计结果中对三种频率时采样分析总结：(1)欠采样即fs2fmax，时域波形恢复过程中已经不能完整表示原信号，有了失真，从频谱上也可看出，同的频谱带互相重叠，已经不能体现原信号频谱的特点了，从而无法得到原来的信号。（2）临界采样即fs=2fmax，时域波形仍然不能恢复完整的原信号，信号只恢复过程中恢复