使用含有纯滞后最小二乘法对含有色噪声的离散模型的参数辨识——Matlab实现

最新推荐文章于 2024-07-27 08:23:39 发布

钩吾之山

最新推荐文章于 2024-07-27 08:23:39 发布

阅读量453

点赞数

文章标签： matlab 最小二乘法算法

本文链接：https://blog.csdn.net/qq_44624573/article/details/129962583

版权

具体原理先顾不上写了，等写完作业再说，先把程序放在这里

%测试
clc
clear
close all
M = 2^10;%数据长度
%%
%产生有色噪声
wnk = idinput(M,'rgs');%白噪声
%plot(wnk)
D=[1 -1.5 0.7 0.1]; C=[1 0.5 0.2]; % D、C多项式的系数
nd=length(D)-1; nc=length(C)-1; % nd、nc为D、C的阶次
wnk0=zeros(nc,1); % 白噪声初值
ek0=zeros(nd,1); % 有色噪声初值
e = zeros(M,1);
for k=1:M
    e(k)=-D(2:nd+1)*ek0+C*[wnk(k);wnk0]; % 产生有色噪声
    %数据更新
    
    for i=nd:-1:2
        ek0(i)=ek0(i-1);
    end
    ek0(1)=e(k);
    for i=nc:-1:2
        wnk0(i)=wnk0(i-1);
    end
    wnk0(1)=wnk(k);
end
plot(e);xlabel('k'); ylabel('噪声幅值'); title('有色噪声序列');
%使用白噪声产生有色噪声的原理以后再补上
%%
u = idinput(M,'prbs');%输入的M序列
u = (u+1)/2;%将变量范围控制在0到1之间
figure
plot((1:M),u)
xlabel('t'); ylabel('输入幅值'); title('M序列');
%A的第二和第三个元素，B的两个元素就是要辨识的对象，一共四个参数
A = [1 -1.5 0.7];%AZ系数
B = [1 0.5];%BZ系数
na = length(A)-1;%AZ的阶次
nb = length(B);%BZ的阶次
lamda = 0.1;%信噪比
zk0=zeros(na,1); % 输入值初值
uk0=zeros(nb,1); % 输出值初值
z = zeros(M,1);%输出序列
for k=1:M

    z(k)=-A(2:na+1)*zk0+B*uk0; % 产生输出


    %数据更新
    zk0(2) = zk0(1);
    zk0(1)=z(k);
    uk0(2) = uk0(1);
    uk0(1)=u(k);
end
z = z+lamda*e;%含噪声的输出序列
figure
plot(z)
title('含噪声的输出序列')

%%
%开始辨识参数
d = 3;%滞后阶次
theta = zeros(4,M+d+na);%辨识参数
theta_T = [-1.5;0.7;1;0.5];%参数标准值
P0 = 10000*eye(na+nb);%Pk初始值
%输入输出序列补0是因为为了解决有色噪声，使用了纯滞后阶次
u_plus = cat(1,zeros(d+na,1),u);
z_plus = cat(1,zeros(na,1),z);
for i = 1:M
    hh = [flip((-1)*u_plus(i:i-1+na));flip(u_plus(i-nb+d+na:i-1+d+na))];%计算h*
    h = [flip((-1)*z_plus(i:i-1+na));flip(u_plus(i+d:i-1+nb+d))];%计算hL
    k = P0*hh*(1/(h'*P0*hh+1));
    theta(:,i+d+na) = theta(:,i+d+na-1)+k*(z_plus(i+na)-h'*theta(:,i+d+na-1));%参数更新
    P0 = (eye(na+nb)-k*h')*P0;
end
figure
plot((1:length(theta)-d-na),theta(:,d+na+1:end),'LineWidth',2)
legend('a1','a2','b2','b2')
grid on
theta_T - theta(:,end)%辨识误差