回声信号的产生与消除.doc

数字信号处理课程设计回声信号的产生与消除姓 名 张针海 学 号 10300123 专业 电子信息工程 指导教师 樊玲 年 级 10级电信2班 日期 2013 .5 . 25 【摘要】 本课程是利用Windows下的录音机，录制一段自己不小于10s的语音，然后在Matlab软件平台下，利用函数wavread对语音信号进行采样，并记录采样频率和采样点数。在抽样信号的基础上，通过采样后的的信号与原信号实现一次及多次延迟、叠加产生回波信号，再使用Matlab绘出有回声及无回声语音信号的时域波形和频谱图。再分别用频率抽样法设计的FIR滤波器和冲激相应不变法设计设计的IIR滤波器消除回声，并记录滤波器的频域响应，再绘制滤波后信号的时域波形和频谱，并对前后信号进行对比，分析信号的变化。关键词 录音 matlab 采样 滤波 抽样Abstract this course is to use a tape recorder to record voice under Windows, a section of their own not less than 10s, then in Matlab software platform, sampling of the speech signal using the function wavread, and record the sampling frequency and sampling points. Based on the sampling signal, through its implementation of single and multiple superposition delay, echo, and use Matlab to draw the echo and echo free speech signal time-domain waveform and spectrum. FIR filter respectively by frequency sampling design method and impulse corresponding invariant IIR filter design to eliminate echo, and record the response of the filter in frequency domain, and then draw the time-domain waveform and spectrum of the filtered signal, and compared before and after the signal, analysis of signal changes目录1 设计目的及要求31.1设计回音目的及要求31.2设计滤波器目的及要求31.2.1 FIR滤波器31. 2. 2 巴特沃兹滤波器31. 2. 3 距离估计要求42 设计原理43设计内容431语音采集.43. 2信号分析433制作回音53.4设计滤波器及滤波83. 4. 1 设计FIR滤波器及滤波83.4.1.1单回声的滤波83. 4.2设计巴特沃兹滤波器及滤波123.4.2.1设计巴特沃斯数字低通滤波器1235估算距离143.5.1通过理论计算法143.5.2程序返回测量法144总结155、参考文献151 设计目的及要求1.1设计回音目的及要求现代通信中回波是影响通信质量的噪声，本课程设计是在matble库元件中搜索一段不小于10s的录音，再利用函数wavread对语音信号进行采样，并自身实现一次及多次延迟、叠加产生回波信号，再使用Matlab绘出有回声及无回声语音信号的时域波形和频谱图。在此过程中必须灵活运用matlab中 wavread函数对信号进行采样，同时加深对声频信号中噪声的认识。1.2设计滤波器目的及要求1.2.1 FIR滤波器 FIR滤波器是有限长冲激响应滤波器，它是通信，语音与图象处理模式识别及频谱分析等应用中一种基本的处理部件。他可以满足滤波器幅度和相位特性的严格要求避免模拟滤波器无法克服的电压漂移，温度漂移和噪声等问题。数字滤波器的作用是滤除信号中某一部分频率分量。信号经过滤波处理。就是相当于信号频谱与滤波响应相乘的结果。时域上看就是信号与滤波器的冲激响应卷积的结果。有限长冲激响应滤波器可以保证任何幅频特性的同时具有严格的相频特性，同时其单位冲激响是有限的。没有输入到输出的反馈，是稳定的系统。1. 2. 2 巴特沃兹滤波器 欲消除回声，则需设计滤波器，绘出滤波器的频域响应，绘出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化，理想低通滤波器的一个非因果系统，如果允许低通滤波器的通带和阻带之间有一定的过渡带，且通带和阻带允许有一定的衰减，我们就可以用物理可实现的系统去逼近理想低通滤波器的频率特性，从而获得较好的滤波效果；巴特沃兹滤波器就是工程中常用的频率响应逼近理想低通滤波器的物理可实现系统，还原出原有音频信号。1. 2. 3 距离估计要求从信号y中估计反射物的距离，即回声延迟的时间，可理解为估算原有信号中的N值。也就是，估测y(n)中的原始声音信号x(n)与其延时衰减分量kx(n-N)的相关联的程度。2 设计原理 语音采集： 通过搜索在pc机的中matble库函数中chimes.wav，绘制其时域波形，对此音频信号用FFT作谱分析 。 制作回音： 利用MATLAB软件中的wavread函数进行采样分析得到的源噪声，通过叠 加到原始语音信号中，模仿语音信号被污染，并对其进行频谱分析； 消除回音： 设计FIR和IIR数字滤波器的滤波器，对带有回声的声音信号进行滤波，恢复出原始信号。绘制所设计滤波器的幅频和相频特性，及滤波后的信号的时域波形和频谱图，分析滤波后信号的时域和频域特征，回放语音信号。 估算距离： 计算信号y中估计反射物的距离，可理解为估计延时序列中的N值。即估计y(n)中的原始声音信号x(n)与其延时衰减分量kx(n-N)的相关联的程度；相关是指两个确定信号或两个随机信号之间的相互关系，对于随机信号，信号一般是不确定的，但是通过对它的规律进行统计，它们的相关函数往往是确定的，故在随机信号处理中，可以用相关函数来描述一个平稳随机信号的统计特性，因而可通过相关分析法估算反射物的距离。3设计内容31语音采集 读取本地matble库中的自带的文件chimes.Wav。3. 2信号分析在matbal中绘制其时域波形，利用wavread()函数将其提取出来绘制其时域波形对此音频信号用FFT作谱分析，用plot()函数绘制其图形。%声音信号的提取fs=22050;x,fs,Nbits=wavread(D:chimes.wav); %把语音信号进行加载入Matlab仿真软件平台中wavplay(x,fs);% 回放语音信号。或者sound(x,fs)figure(1);N=length(x);%求语音信号的长度subplot(3,1,1);plot(x(1:N);title(原始信号波形);y=fft(x,N);%傅立叶变换subplot(3,1,2);plot(0:N-1)/N*fs,abs(y);title(原始信号幅值);subplot(3,1,3);plot(angle(y);title(原始信号相位);2，产生的原信号的波形，以及其幅度、相位谱如图1所示：图1 原信号波形图33制作回音1 参数的设置： 因为人耳能分辨出的声音延迟至少是0.1s，因此，最小延迟量不能小于0.1s。在此先先延迟时间为0.2s,即最小延迟量N=0.2*fs=0.2*16000=3200。 在已有声音信号x的基础上产生带回声的声音信号，可以表达为在于娜信号的基础上叠加其延时的分量。假设只有一个回声的情况下，可简化其模型为 y(n)=x(n)+kx(n-m) k为反射系数；m为延迟时间。 这里设m=3200 ； k=0.412 声音延迟： ,利用矩阵置零产生x的延迟，以及得到y信号：为了保证图像的完整性，对读取的信号先延长4000个采样点，将原始信号延长3200个采样点，然后再在后面补上800个点，得到如下代码： 1 %一次回波的产生x,fs=wavread(D:chimes.wav);%把语音信号进行加载入Matlab仿真软件平台中。N=length(x);%语音信号的长度。x1=x(1:N);x2=x(1:N);x1=x1,zeros(1,4000);%zeros(1,4000)产生1行3000列全零矩阵加到x1后面。x2=zeros(1,3200),0.4*x2,zeros(1,800);%N+3200+x=N+4000，得x=800y=x1+x2;%加入回音的信号。figure(2);subplot(3,1,1);plot(y(1:3200.+N);title(含回声信号波形);y1=fft(y);subplot(3,1,2);plot(abs(y1);title(含回声信号幅值);subplot(3,1,3);plot(angle(y1);title(含回声信号相位);sound(y,fs);加入回声以后的波形，如图2所示：图2 含一次回声波形图2%多次回声的产生x,fs=wavread(D:chimes.wav);%把语音信号进行加载入Matlab仿真软件平台中。N=length(x);%语音信号的长度。x1=x(1:N);x2=x(1:N);x3=x(1:N);x1=x1,zeros(1,6400);%zeros(1,4000)产生1行3000列全零矩阵加到x1后面。x2=zeros(1,3200),0.4*x2,zeros(1,3200);%N+3200+x=N+4000，得x=800x3=zeros(1,3200),0.7*x3,zeros(1,3200);y=x1+x2+x3;%加入回音的信号。figure(2);subplot(3,1,1);plot(y(1:3200.+N);title(含回声信号波形);h=fft(y,N);%傅立叶变换subplot(3,1,2);plot(0:N-1)/N*fs,abs(h);title(含回声信号幅值);subplot(3,1,3);plot(angle(y);title(含回声信号相位);sound(y,fs);程序设计如图3所示：图3 含多次回声波形图3.4设计滤波器及滤波3. 4.1 设计FIR滤波器及滤波3.4.1.1单回声的滤波%单个回声的滤波x,fs=wavread(D:chimes.wav);%把语音信号进行加载入Matlab仿真软件平台中。N=length(x);%语音信号的长度。x1=x(1:N);x2=x(1:N);x3=x(1:N);subplot(6,1,1);plot(x);title(原信号波形);h=fft(x,N);%傅立叶变换subplot(6,1,2);plot(0:N-1)/N*fs,abs(h);title(原始信号幅值);x1=x1,zeros(1,4000);%zeros(1,4000)产生1行4000列全零矩阵加到x1后面。x2=zeros(1,3200),0.8*x2,zeros(1,800);%N+3200+x=N+4000，得x=800y=x1+x2;%加入回音的信号。subplot(6,1,3);plot(y(1:3200.+N);title(含回声信号波形);y1=fft(y);subplot(6,1,4);plot(abs(y1);title(含回声信号幅值);b=1;a=zeros(1,N);a(1)=1;a(3201)=0.4;z=filter(b,a,y);subplot(6,1,5);plot(z);title(滤波后信号波形);z1=fft(z);subplot(6,1,6); %3行 1列排列第1个图plot(abs(z1); %画出abs（z1）title(滤波后信号的幅值 ); %绘图区标题sound(z1,fs);程序设计后运行结果如图4所示：图4 滤除单个回波后的波形图 3.4.1.2多次回声的滤波 %多个回声的滤除x,fs=wavread(D:chimes.wav);%把语音信号进行加载入Matlab仿真软件平台中。N=length(x);%语音信号的长度。x1=x(1:N);x2=x(1:N);x3=x(1:N);subplot(6,1,1);plot(x);title(原信号波形);h=fft(x,N);%傅立叶变换subplot(6,1,2);plot(0:N-1)/N*fs,abs(h);title(原始信号幅值);x1=x1,zeros(1,2400);%zeros(1,4000)产生1行3000列全零矩阵加到x1后面。x2=zeros(1,3200),0.4*x2,zeros(1,3200);%N+3200+x=N+6400，得x=3200x3=zeros(1,5200),0.7*x3,zeros(1,5200);y=x1+x2+x3;%加入回音的信号。subplot(6,1,3);plot(y(1:3200.+N);title(含回声信号波形);y1=fft(y);subplot(6,1,4);plot(abs(y1);%plot(0:N-1)/N*fs,abs(y1);title(含回声信号幅值);%subplot(6,1,5);%plot(angle(y1);%title(含回声信号相位);b=1;a(1)=1;a(3201)=0.4;a(5201)=0.7;z=filter(b,a,y);subplot(6,1,5);plot(z);title(滤波后信号波形);z1=fft(z);subplot(6,1,6); %3行 1列排列第1个图%plot(0:N-1)/N*fs,abs(z1);plot(abs(z1); %画出abs（z2title(滤波后信号的幅值 ); %绘图区标题sound(z1,fs);程序设计后运行结果如图5所示：:图5 滤除多个回波后的波形3.4. 2设计巴特沃兹滤波器及滤波3.4.2.1设计巴特沃斯数字低通滤波器fp=2100; fs=8000; Fs=20000; Rp=0.5;Rs=30; T=1/Fs; %设计指标 W1p=fp/Fs*2;W1s=fs/Fs*2;%求归一化频率 N,Wn=buttord(W1p,W1s,Rp,Rs,s); %确定butterworth的最小介数N和频率参数Wnz,p,k=buttap(N); %设计模拟低通原型的零极点增益参数 bp,ap=zp2tf(z,p,k); %将零极点增益转换成分子分母参数 bs,as=lp2lp(bp,ap,Wn*pi*Fs);%将低通原型转换为模拟低通bz,az=impinvar(bs,as,Fs); %用脉冲响应不变法进行模数变换 sys=tf(bz,az,T);%给出传输函数H(Z) H,W=freqz(bz,az,512,Fs); %生成频率响应参数 subplot(2,1,1);plot(W,20*log10(abs(H); %绘制幅频响应 grid on; %加坐标网格 xlabel(频率/Hz);ylabel(振幅/dB);subplot(2,1,2); plot(W,abs(H);grid on; xlabel(频率/Hz);ylabel(振幅/H)； 运行后的结果所得波形如图6所示： 图6 巴特沃兹滤波器图运行结果如下： N =4 bz = 0.0000 0.0999 0.1914 0.0252 az= 1.0000 -1.4336 1.0984 -0.4115 0.0627可以得出：只需编程，结果非常直观35估算距离3.5.1通过理论计算法 （1） %计算单个回声方法 m=3200，而采样频率为fs=25050Hz，可得延时t=m/fs=0.145s 设反射物距离为s，声速为v，则有2*s=v*t 得：s=v*t/2；声速大概在340m/s左右，把t=0.13s,v=340m/s代入求得s=24.65m，即反射物的距离为24.65米。3.5.2程序返回测量法（2）采用相关分析法从带有多个回声的声音信号中估计反射物的距离。 %相关函数法测障碍物距离运行后的结果所得波形如图7所示：h=xcorr(y);figure(4); %建立图形窗口subplot(4,1,1) %4行1列第1个图plot(abs(h); %画出abs（h）title(最值 ); %绘图区标题h1=h(17800:18000); %h1矩阵r1,t1=max(h1); %最大值t1=t1+17800-1;subplot(4,1,2); %4行1列第2个图plot(h1) %画出h1title(点1); %绘图区标题h2=h(21000:22000); %h2矩阵r2,t2=max(h2); %最大值t2=t2+21000-t1;subplot(4,1,3); %4行1列第3个图plot(h2); %画出h3title(点2); %绘图区标题h3=h(23000:23200); %h3矩