时域GSC原理程序

clear
clc
B=[1,-1,0,0;0,1,-1,0;0,0,0,1];
x1=wavread('noisy0.wav');
x2=wavread('noisy10.wav');
x3=wavread('noisy15.wav');
x4=wavread('noisy20.wav');
d=wavread('clean.wav');%期望信号
M=length(x1);
N=4;%输入信号端数
u=0.008;
sysorder=512;%滤波器阶数
w1=zeros(sysorder,1);
w2=zeros(sysorder,1);
w3=zeros(sysorder,1);
y=zeros(1,sysorder);
W=[w1';w2';w3'];
for n=sysorder+1:M
    X1=x1(n-sysorder+1:1:n);
    X2=x2(n-sysorder+1:1:n);
    X3=x3(n-sysorder+1:1:n);
    X4=x4(n-sysorder+1:1:n);
    X=[X1';X2';X3';X4'];
    Z=B*X;
    Y=W.*Z;
    y(n)=sum(Y(:));
    e(n)=d(n)-y(n);
    W=W+u*e(n)*Z;
end
figure(1)
subplot(2,1,1)
plot(x1);title('输入信号');
subplot(2,1,2)
plot(y);title('输出信号');
wavplay(x1);
wavplay(y);