wiener滤波程序.docx
yyhxcorr=xcorr(x(l:4096);size(yyhxcorr);=yyhxcorr (4096:4595);yyhdcorr=xcorr(d(l:4096), x(l:4096);size(yyhdcorr);B=yyhdcorr(4096:4595);M=500;yyhresult=wienerf i 1 ter (x, A, B, M) ;%用维纳滤波进行去噪yyhresult=yyhresult(300:8192+299);fq=fft(yyhresult);subplot (3, 1, 3);f=Fs*(0:4095)/8192;plot(f, abs(fq(l :4096);titleC经过维纳滤波后语音信号的频域图形')；xlabel (,频率 f,);ylabel (, FFT');wienerfilter.mfunction y=wienerfilter(x, Rxx, Rxd, N)%进行维纳滤波%是输入信号，R是输入信号的自相关向量%R是输入信号和理想信号的的互相关向量，N是维纳滤波器的长度%输出y是输入信号通过维纳滤波器进行维纳滤波后的输出h=yulewalker (Rxx, Rxd, N) ;%求解维纳滤波器系数t=conv (x, h) ;%进行滤波Lh=length(h) ;%得到滤波器的长度Lx=length(x) ;%得到输入信号的长度y=t (double(uintl6(Lh2) :Lx+double(uint 16(Lh/2)-1) ;%输出序列 y 的长度和输入序列x的长度相同%以下是维纳滤波器系数的求解function h=yulewalker(A, B, M)%求解Yule-Walker方程%A是接收信号的自相关向量为Rxx(0),Rxx(l), Rxx(M-l)%B是接收信号和没有噪声干扰信号的互相关向量为Rxd(0),Rxd(l),Rxd(M-l)%M是滤波器的长度%h保存滤波器的系数%例如 A=6, 5, 4, 3, 2, 1 ;B=100, 90, 120, 50, 80, 200 ;M=6;%求解出 h=26. 4286 -20. 000050. 0000 -50. 0000 -45. 000081. 4286,Tl=zeros(l,M) ;%T1存放中间方程的解向量T2=zeros(l,M) ;%T2存放中间方程的解向量T1(1)=B(1)/A(1);T2(l)=A(2)A(l)jX=zeros(1, M);for i=2:M-ltempl=0;temp2=0;for j=l:i-ltempl=templ+A(i-j+l)*T1 (j);temp2=temp2+A(i-j+l) *T2(j);endX(i) = (B(i)-tempi)/(A(l)-temp2);for j=l:i-lX(j)=Tl(j)-X(i)*T2(j)jendfor j=l:iTl(j)=X(j);endtemp1=0;temp2=0;for j=l:i-ltempl=templ+A(j+l)*T2(j);temp2=temp2÷A(j÷l)*T2(i-j);endX(l) = (A(i + l)- temp 1)/(A(1)-1emp2);for j=2:iX(j)=T2(j-l)-X(l)*T2(i-j+l)jendfor j=l:iT2(j)=X(j);endendtempl=0;temp2=0;for j=l:M-1templ=templ+A(M-j÷l)*T1(j);temp2=temp2+A(M-j+l)*T2(j);endX (M) = (B (M) -tempi) / (A (1) -temp2);for j=l:M-1X(j)=Tl(j)-X(M)*T2(j);endh=X;
