clear all
close all
%channel system order
sysorder = 5
% Number of system points
N=2000
inp = randn(N,1)
n = randn(N,1)
[b,a] = butter(2,0.25)
Gz = tf(b,a,-1)
%This function is submitted to make inverse Z-transform (Matlab central file exchange)
%The first sysorder weight value
%h=ldiv(b,a,sysorder)'
% if you use ldiv this will give h :filter weights to be
h= [0.0976
0.2873
0.3360
0.2210
0.0964]
y = lsim(Gz,inp)
%add some noise
n = n * std(y)/(10*std(n))
d = y + n
totallength=size(d,1)
%Take 60 points for training
N=60
%begin of algorithm
w = zeros ( sysorder , 1 )
for n = sysorder : N
u = inp(n:-1:n-sysorder+1)
y(n)= w' * u
e(n) = d(n) - y(n)
% Start with big mu for speeding the convergence then slow down to reach the correct weights
if n <20
mu=0.32
else
mu=0.15
end
w = w + mu * u * e(n)
end
%check of results
for n = N+1 : totallength
u = inp(n:-1:n-sysorder+1)
y(n) = w' * u
e(n) = d(n) - y(n)
end
hold on
plot(d)
plot(y,'r')
title('System output')
xlabel('Samples')
ylabel('True and estimated output')
figure
semilogy((abs(e)))
title('Error curve')
xlabel('Samples')
ylabel('Error value')
figure
plot(h, 'k+')
hold on
plot(w, 'r*')
legend('Actual weights','Estimated weights')
title('Comparison of the actual weights and the estimated weights')
axis([0 6 0.05 0.35])
% RLS 算法
randn('seed', 0)
rand('seed', 0)
NoOfData = 8000 % Set no of data points used for training
Order = 32% Set the adaptive filter order
Lambda = 0.98% Set the forgetting factor
Delta = 0.001% R initialized to Delta*I
x = randn(NoOfData, 1) % Input assumed to be white
h = rand(Order, 1) % System picked randomly
d = filter(h, 1, x) % Generate output (desired signal)
% Initialize RLS
P = Delta * eye ( Order, Order )
w = zeros ( Order, 1 )
% RLS Adaptation
for n = Order : NoOfData
u = x(n:-1:n-Order+1)
pi_ = u' * P
k = Lambda + pi_ * u
K = pi_'/k
e(n) = d(n) - w' * u
w = w + K * e(n)
PPrime = K * pi_
P = ( P - PPrime ) / Lambda
w_err(n) = norm(h - w)
end
% Plot results
figure
plot(20*log10(abs(e)))
title('Learning Curve')
xlabel('Iteration Number')
ylabel('Output Estimation Error in dB')
figure
semilogy(w_err)
title('Weight Estimation Error')
xlabel('Iteration Number')
ylabel('Weight Error in dB')
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)