% 本实验中利用近似方法，即最佳FIR维纳滤波方法，在计算机上实现随机信号的维纳滤波。
% w(n)是零均值,方差为(1-a^2)的均匀分布白噪声
% s(n)为真实信号：s(n)=a*s(n-1)+w(n)
% v(n)是与s(n)互不相关的均匀分布白噪声，其均值为零，方差为1
% x(n)为接收到的添加了白噪声的信号：x = s + v
% si(n) 为使用理想维纳滤波器滤波得到结果
% y(n) 为使用近似FIR维纳滤波器滤波得到的结果
% If you have any problem, please email me:
qijiaxing@live.cn
% This program is written by xndpd/p>
% 2009.11.20
clc
clear all
N = 20; % length of the FIR
filter
a = 0.95;
% white noise with mean of 0 and var of (1-a^2)
w = sqrt( 12 * (1 - a^2)) * ( rand(1,L) -
0.5 );
% true signal: s(n)=a*s(n-1)+w(n)
s = zeros(1,L); s(1) = w(1);
for ii = 2:欢喜的铃铛/p>
s(ii) =
a * s(ii-1) + w(ii);
end
% white noise with mean of 0 and var of 1
v = sqrt( 12 ) * ( rand(1,L) - 0.5 );
% received signal: x = s + v
x = s + v;
% r_xx is the autocorrelation of x
r_xx = xcorr( x );
% R_xx is the N-dimentional autocorrelation matrix of x
R_xx = zeros( N );
for ii = 1 : N
R_xx( :
, ii ) = r_xx( L+1-ii : L+1-ii+N-1 )';
end
% r_xs is the cross-correlation of x and
s
r_xs = xcorr( x , s );
r_xs = r_xs( L : L+N-1 )';
% according to R_xx * h_FIR = R_xs
% we can calculate the h_FIR
h_FIR = inv(R_xx) * r_xs;
h_FIR = h_FIR';
% y is the reslut signal filtered by h_FIR
y = cconv( h_FIR , x, L );
% si is the result signal filtered by ideal Wiener filter
h
n = 0:L-1;
h = 0.238 * ( 0.724 .^ n );
si = cconv( h , x, L );
%% plot the true signal s(n) and the received signal x(n) of
the last 100
% points
t = L - 99 : L;
figure(1);
plot( t , x(t), '--b' , t , s(t) , 'r' );
legend( 'x(n)' , 's(n)' , 0 );
title('The true signal s(n) and the received signal
x(n)');
xlabel('n');ylabel('Signal Amplitude');
%% plot the ideal impulze response h(n) and the FIR impulze
response
% h_FIR(n) of N points
n = 1:N;
figure(2);
plot( n,h(n),'--b',n,h_FIR(n),'r' );
legend('The Ideal Wiener Filter h(n)','Approximate FIR Filter
h_FIR(n)',0);
title('The Impulze Response of Ideal Wiener Filter h(n) and
the Impulze Response h_FIR(n)');
xlabel('N points');ylabel('Amplitude');
%% plot the true signal s(n) and the result signal y(n)
filtered by
% h_FIR
figure(3);
plot( t, s(t), '--b', t, y(t), 'r' );
legend('True Signal s(n)','Result Signal y(n)',0);
title('The True Signal s(n) and The Result Signal y(n) Filter
by h_FIR');
xlabel('n');ylabel('Signal Amplitude');
%% plot the true signal s(n) and the result signal si(n)
filtered by h
figure(4);
plot( t, s(t), '--b', t, si(t), 'r' );
legend('True Signal s(n)','Result Signal si(n)',0);
title('The True Signal s(n) and The Result Signal si(n)
Filtered by h');
xlabel('n');ylabel('Signal Amplitude');