24小时热门版块排行榜    

查看: 1579  |  回复: 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 的主题更新
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 关于Filecode分析方法 +8 majunge000 2026-08-10 9/450 2026-08-12 22:17 by 虫友是什么虫
[基金申请] 不应该看fileCode +5 且听虎啸 2026-08-12 6/300 2026-08-12 16:26 by Tide man
[论文投稿] 职称评审,求友友推荐见刊最快的期刊 +4 工厂打螺丝 2026-08-08 4/200 2026-08-12 09:18 by hansi2025
[基金申请] 是这周出结果还是下周出结果? +3 yuleib84 2026-08-11 3/150 2026-08-11 21:53 by jnhyjjm
[基金申请] 小木虫上这么多卖论文的,真有人买论文么?感觉没必要啊 +8 Tide man 2026-08-10 9/450 2026-08-11 20:39 by beefly
[基金申请] 有时候,自然基金真的不能太认真 (我的申报经验) +3 majunge000 2026-08-11 4/200 2026-08-11 20:13 by lch2012
[基金申请] 分享一下我之前已中青C的计划书的filecode +3 布布和一二 2026-08-11 4/200 2026-08-11 15:20 by aasahr
[基金申请] 听说今天filecode变了 +26 布布和一二 2026-08-06 49/2450 2026-08-11 13:18 by WH3796
[基金申请] 我的国基提前知道中了,可是同事的操作让我实在接受不了,怎么会有这样的人 +8 家与远方 2026-08-10 13/650 2026-08-11 09:49 by 家与远方
[基金申请] 基金中了 +15 laoda193707 2026-08-06 15/750 2026-08-11 00:11 by jiafei2190
[基金申请] 确定了,国自然21号放榜 +6 布布和一二 2026-08-10 7/350 2026-08-10 19:15 by 2000zf36392
[基金申请] 关于代码变化问题,想知道的进来 +17 且听虎啸 2026-08-07 24/1200 2026-08-10 18:35 by zhangduo2008
[基金申请] 据悉今年马上要出结果了 +7 瞬息宇宙 2026-08-10 8/400 2026-08-10 12:42 by Vivilian
[基金申请] 这样的filecode谁见过 +11 布布和一二 2026-08-08 22/1100 2026-08-10 11:10 by wmfsnow
[基金申请] 关于filecode,很负责任的告诉大家 +6 爱看书的可乐 2026-08-08 7/350 2026-08-08 22:13 by a_niu
[基金申请] 国基金的申报应该改成非等额制,评价高的钱多评价低的钱少,但是增加资助率 +7 a089 2026-08-07 7/350 2026-08-08 18:05 by gltch
[基金申请] 固定端突然变了,今天 +6 archvillain 2026-08-06 10/500 2026-08-07 16:03 by 医学老男孩
[基金申请] filecode变化情况 +6 布布和一二 2026-08-07 22/1100 2026-08-07 14:45 by 且听虎啸
[基金申请] filecode +14 等待解的谜 2026-08-06 19/950 2026-08-07 12:20 by wlwhappy
[基金申请] filecode +8 布布和一二 2026-08-06 11/550 2026-08-06 20:41 by tangpu318
信息提示
请填处理意见