24小时热门版块排行榜    

查看: 1587  |  回复: 5
当前主题已经存档。

zzbei

木虫 (正式写手)

[交流] 自适应matlab程序(LMS和RLS)

%lms算法源程序

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') ;

[ Last edited by 幻影无痕 on 2006-10-26 at 07:47 ]
回复此楼
经历就是一种享受
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

wwfjhx863

0.5

太好了
2楼2006-05-09 07:56:42
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

0.5

3楼2006-05-09 17:21:55
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

wanglong

铁杆木虫 (正式写手)

非常感谢!!!!老大,有啥好东东,多上穿一点
4楼2006-05-11 12:51:08
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

wanglong

铁杆木虫 (正式写手)

0.5

非常感谢,刚才没评价,重来
5楼2006-05-11 12:52:37
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

sanheiliu

木虫 (小有名气)

好用!谢谢!!!
6楼2006-05-16 23:02:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 zzbei 的主题更新
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 时间戳又变了8-15 +8 archvillain 2026-08-15 15/750 2026-08-15 16:23 by 医学老男孩
[考研] 售SCI一区T0P文章,我:8.O.55.1.O.5.4,科目全,可+急 +3 P4JhWrbJGjSe 2026-08-14 5/250 2026-08-15 15:46 by QXtWNz7PI7wZ
[硕博家园] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +3 7lpolszZVXgi 2026-08-14 3/150 2026-08-15 14:53 by QXtWNz7PI7wZ
[基金申请] 有时候,自然基金真的不能太认真 (我的申报经验) +7 majunge000 2026-08-11 9/450 2026-08-15 13:34 by kiz2011
[博后之家] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 k0dTPqJtl0jt 2026-08-14 3/150 2026-08-15 12:41 by KxMI1BYBxWX1
[考研] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +4 7lpolszZVXgi 2026-08-14 8/400 2026-08-15 11:53 by KxMI1BYBxWX1
[考博] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +3 HFw0lei2R37i 2026-08-14 5/250 2026-08-15 11:41 by KxMI1BYBxWX1
[公派出国] 售SCI-T0P文章,我:8O.5.5.1.O.54,科目齐全,可+急 +3 k0dTPqJtl0jt 2026-08-14 5/250 2026-08-15 07:09 by 4wMiSEwB6436
[基金申请] filecode +7 cratir 2026-08-14 11/550 2026-08-15 06:19 by 学员8dgXkO
[公派出国] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +3 7lpolszZVXgi 2026-08-14 5/250 2026-08-15 04:41 by 4wMiSEwB6436
[基金申请] 各位道友,我要去昆明玩几天,回来见。 +7 Tide man 2026-08-14 8/400 2026-08-15 01:11 by arzu_hma
[基金申请] filecode +15 documentary 2026-08-10 17/850 2026-08-14 10:08 by kissu88
[基金申请] 我的国基提前知道中了,可是同事的操作让我实在接受不了,怎么会有这样的人 +10 家与远方 2026-08-10 15/750 2026-08-14 02:08 by 绵羊哥哥
[文学芳草园] 阿姨 +4 汪汪锅 2026-08-09 4/200 2026-08-13 19:43 by arzu_hma
[基金申请] Filecode 又变了,巨变 +3 WH3796 2026-08-12 4/200 2026-08-13 14:13 by 小木虫6752397
[基金申请] 结合人工智能,周易传统文化,filecode打分制来了,3分以上希望很大。 +3 Tide man 2026-08-12 4/200 2026-08-13 08:35 by ZJTJZ
[基金申请] 2019年青年基金涵评意见,大家看看几个A,几个B? +11 Tide man 2026-08-11 11/550 2026-08-13 07:35 by 撸猫猫
[基金申请] 帮忙看看fileCode +7 wwncly 2026-08-10 13/650 2026-08-11 19:36 by 冰心玉壶晴
[基金申请] 据悉今年马上要出结果了 +7 瞬息宇宙 2026-08-10 8/400 2026-08-10 12:42 by Vivilian
[基金申请] 面上项目filecode邪修 +5 西山十月 2026-08-09 7/350 2026-08-10 07:32 by 仁砚薪传
信息提示
请填处理意见