sn.m

用matlab仿真维纳滤波器

%维纳滤波器的计算机实现

clear
L=input('请输入信号样本个数L=');
N=input('请输入滤波器的阶数N=');

%产生一个W1(n)均值零方差1-a^2的均匀分布的白噪声
a=0.95;
p=1-a^2;
u=rand(1,L);
u=u-mean(u);
d=var(u);
power_u=sqrt(1/d*p);
W1=u*power_u;
figure;
plot(W1(L-100:L));
grid on;
ylabel('W(n)');
xlabel('n');
title('W(n)为均值零 方差为1-a^2的均匀分布的白噪声')
W_power=var(W1)
W_e=mean(W1)
pause;


%产生一个s(n)
n=2;
a=0.95;
s=zeros(1,L);
s(1)=0;

while n<=L
       s(n)=a*s(n-1)+W1(n);
        n=n+1;
    end
    figure;
plot(s(L-100:L));grid on;
ylabel('S(n)');
xlabel('n');
title('信号S(n)');
pause;


%产生一个V1(n)均值零方差1的均匀分布的白噪声
p=1;
u=rand(1,L);
u=u-mean(u);
d=var(u);
power_u=sqrt(1/d*p);
V1=u*power_u;
figure;
plot(V1(L-100:L));
grid on;
V_power=var(V1)
V_e=mean(V1)
title('V(n)均值零 方差为1的均匀分布的白噪声')
ylabel('V(n)');
xlabel('n');
pause;

%在同一坐标绘出最后个s(n)和X(n)
n=1;
while n<=L
      X(n)=s(n)+V1(n);
        n=n+1;
    end
figure;
plot(X(L-100:L),'--');
grid on;
hold on;
xlabel('n');
ylabel('X(n)=S(n)+V(n)');
pause;
plot(s(L-100:L),'r');
xlabel('n');
ylabel('S(n)红色');
title('S(n)与X(n)的比较');
hold off;
pause;

%利用L个s(n)和X(n)估计Rxx和rxs
Mlag=N;
rxx=xcorr(X,Mlag,'biased');% 2*N-1 length
rxs=xcorr(X,s,Mlag,'biased');
for i=1:N
   for j=1:N
     Rxx(i,j)=rxx(N+1+abs(i-j));
 end
end
sum=zeros(1,N);
for i=0:(N-1)
    for j=1:L-i
      sum(i+1)=sum(i+1)+X(j)*s(j+i)';
   end  
  Rxs(i+1)=sum(i+1)/(L-i);
end


%估计的h(n)
h=inv(Rxx)*Rxs';


%理论H(n)
H=zeros(1,N);
u1=ones(1,N);
for n=1:N
  H(n)=0.238*(0.724)^n*u1(n);
end
figure;
plot(h,'--');
hold on;
grid on;
title('估计的h(n)');
pause;
plot(H,'r');
hold off;
title('理论的H(n)和估计的h(n)的比较')
xlabel('n');
ylabel('h(n)和H(n)红色');
pause;


%利用公式（1-15），进行理想的维娜滤波得到L个Sl(n),绘图
Sl=zeros(1,L);
for n=2:L
     Sl(n)=0.724*Sl(n-1)+0.238*X(n);
end
figure;
plot(s(L-100:L),'--');
hold on;
grid on;
title('原信号S(n)');
pause;
plot(Sl(L-100:L),'r');
hold off;
title('理想的维娜滤波得到Sl(n)红色,和S(n)比较');
pause;


%由h(n)根据公式（1-6）对x(n)进行过滤得L个Sr(n)，将最后100个s(n)和Sr(n)绘于同一坐标内
Sr=zeros(1,L);
for n=N+1:L
    for m=1:N
        %if m>=n
           % Sr(n)=Sr(n-1);
           %else
    Sr(n)=Sr(n)+h(m)*X(n-m);
    %end
end
end
figure;
plot(s(L-100:L),'--');
hold on;
grid on;
title('原信号S(n)');
pause;
plot(Sr(L-100:L),'r');
hold off;
title('对x(n)进行过滤得L个Sr(n)红色,和S(n)比较');
pause;


%由L个x (n)，s ( n) ， 和Sl(n)和Sr(n)，根据公式统计ex,el,er
ex1=zeros(1,L);
for i=2:L
ex1(i)=ex1(i-1)+(X(i)-s(i))^2;
end
ex=ex1(L)/L

ei1=zeros(1,L);
for i=2:L
ei1(i)=ei1(i-1)+(Sl(i)-s(i))^2;
end
ei=ei1(L)/L

er1=zeros(1,L);
for i=2:L
er1(i)=er1(i-1)+(Sr(i)-s(i))^2;
end
er=er1(L)/L