【信号处理】扩展卡尔曼滤波EKF(Matlab代码实现)

💥💥💞💞欢迎来到本博客❤️❤️💥💥

🏆博主优势:🌞🌞🌞博客内容尽量做到思维缜密,逻辑清晰,为了方便读者。

⛳️座右铭:行百里者,半于九十。

📋📋📋本文目录如下:🎁🎁🎁

目录

💥1 概述

📚2 运行结果

🎉3 参考文献

🌈4 Matlab代码实现


💥1 概述

扩展卡尔曼滤波(extended kalman filter,EKF)是在线性最小方差估计的基础上,提出的一种实时递推型的最优估计滤波算法,近年来被逐步应用于工程结构的参数识别研究[1.2.3.4.5]当中。为解决传统EKF算法中其状态向量维度过大导致该算法运行效率变慢、识别精度下降等问题,本文提出了一种改进的EKF算法,引入模态坐标变换,利用模态坐标转换对传统EKF的状态向量进行处理,构建以模态坐标初始值和结构损伤参数的状态向量。由于结构振动响应通常由前几阶的低阶频率成分组成,因此,对复杂结构可以有效缩减状态向量的维数,以保证算法的稳定性及准确性。
 

📚2 运行结果

 

 

 

 

 

部分代码:

%观测站的位置
x1=-10;
y1=0;
x2=10;
y2=0;

r1 = sqrt((Rx-x1)^2+(Ry-y1)^2);
beta1 = atan2((Ry-y1),(Rx-x1));
r2 = sqrt((Rx-x2)^2+(Ry-y2)^2);
beta2 = atan2((Ry-y2),(Rx-x2));

%noise
sigma_u = sqrt(0.0001);     %过程噪声
sigma_R = sqrt(5);        %距离量测噪声
sigma_beta = sqrt(0.0001);    %角度量测噪声

%% Kalman filter CV 2D

%-------Kalman Parameters-------%

A1 = [cos(beta1) -r1*sin(beta1); sin(beta1) r1*cos(beta1)] ;
R1 = A1*[sigma_R^2 0;0 sigma_beta^2]*A1' ;
P1 = [R1(1,1)   R1(1,1)/T     R1(1,2)   R1(1,2)/T
    R1(1,1)/T 2*R1(1,1)/T^2 R1(1,2)/T 2*R1(1,2)/T^2
    R1(1,2)   R1(1,2)/T     R1(2,2)   R1(2,2)/T
    R1(1,2)/T 2*R1(1,2)/T^2 R1(2,2)/T 2*R1(2,2)/T^2 ];
% P1 = 100*eye(4);
A2 = [cos(beta2) -r2*sin(beta2); sin(beta2) r2*cos(beta2)] ;
R2 = A2*[sigma_R^2 0;0 sigma_beta^2]*A2' ;
P2 = [R2(1,1)   R2(1,1)/T     R2(1,2)   R2(1,2)/T
    R2(1,1)/T 2*R2(1,1)/T^2 R2(1,2)/T 2*R2(1,2)/T^2
    R2(1,2)   R2(1,2)/T     R2(2,2)   R2(2,2)/T
    R2(1,2)/T 2*R2(1,2)/T^2 R2(2,2)/T 2*R2(2,2)/T^2 ];
% P2 = 100*eye(4);
%状态转移矩阵
F = [1 T 0 0 
     0 1 0 0
     0 0 1 T
     0 0 0 1];

%过程噪声
B = [T^2/2, 0; T, 0;
     0, T^2/2; 0, T]; %过程噪声分布矩阵
v = sigma_u^2;   %x方向的过程噪声向量//相当于Q
V = B * v * B';
% %观测噪声??
% W = B * noise_x;

%------Data initial-------%
X_real = zeros(4,N);
X = zeros(4,N);

Z1 = zeros(2,N);
Z_polar1 = zeros(2,N);
X_EKF1 = zeros(4,N);
Z2 = zeros(2,N);
Z_polar2 = zeros(2,N);
X_EKF2 = zeros(4,N);

X_CC = zeros(4,N);
X_BC = zeros(4,N);
bias = zeros(8,N,M);

%-------Track Initial-------%
X_real(:,1) = [Rx, vx, Ry, vy]'; %x: km,km/s

X_EKF1(:,1) = [Rx, 0, Ry, 0];
X_EKF2(:,1) = [Rx, 0, Ry, 0];
X_CC(:,1) = [Rx, 0, Ry, 0];
X_BC(:,1) = [Rx, 0, Ry, 0];

%Monto-carlo
for m=1:M
    
    noise_u = randn(2,N).*sigma_u;
    noise_R = randn(1,N).*sigma_R; %过程噪声
    noise_beta = randn(1,N).*sigma_beta; %观测噪声

🎉3 参考文献

部分理论来源于网络,如有侵权请联系删除。

[1]刘志勇,王阿利,王小红.一种基于自构架模糊EKF的目标跟踪方法[J].计算技术与自动化,2022,41(03):14-20.DOI:10.16339/j.cnki.jsjsyzdh.202203003.

[2]胡文强,胡建鹏,吴飞,陆雯霞.基于EKF的初始状态自适应室内融合定位算法[J].传感器与微系统,2022,41(11):147-151.DOI:10.13873/J.1000-9787(2022)11-0147-05.

🌈4 Matlab代码实现

  • 0
    点赞
  • 3
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
以下一个简单的扩展卡滤波(E)的Matlab示例代码: ```matlab% 初始化 x = [0; 0]; % 初始状态向量 [位置; 速度] P = eye(2); % 初始状态协方差矩阵 dt = 0.1; % 时间间隔 % 系统模型 A = [1 dt; 0 1]; % 状态转移矩阵 B = [0.5*dt^2; dt]; % 输入矩阵 H = [1 0]; % 观测矩阵 Q = diag([0.1, 0.1]); % 系统噪声协方差矩阵 R = 1; % 观测噪声方差 % 数据生成 T = 10; % 总时间 N = T/dt; % 数据点个数 u = sin(0:dt:T-dt)'; % 控制输入 x_true = zeros(N, 2); % 实际状态 z = zeros(N, 1); % 观测值 for k = 1:N x_true(k, :) = A*x + B*u(k) + sqrt(Q)*randn(2, 1); % 状态更新 z(k) = H*x_true(k, :) + sqrt(R)*randn; % 观测更新 x = A*x + B*u(k) + sqrt(Q)*randn(2, 1); % 状态预测 end % 扩展卡滤波 x_est = zeros(N, 2); % 估计状态 P_est = zeros(2, 2, N); % 估计状态协方差矩阵 x_est(1, :) = [0; 0]; % 初始估计状态 P_est(:, :, 1) = eye(2); % 初始估计状态协方差矩阵 for k = 2:N % 预测步骤 x_pred = A*x_est(k-1, :) + B*u(k); P_pred = A*P_est(:, :, k-1)*A' + Q; % 更新步骤 K = P_pred*H'/(H*P_pred*H' + R); x_est(k, :) = x_pred' + K*(z(k) - H*x_pred); P_est(:, :, k) = (eye(2) - K*H)*P_pred; end % 可视化结果 t = 0:dt:T-dt; figure; subplot(2, 1, 1); plot(t, x_true(:, 1), 'b', t, x_est(:, 1), 'r'); legend('真实位置', '估计位置'); xlabel('时间'); ylabel('位置'); subplot(2, 1, 2); plot(t, x_true(:, 2), 'b', t, x_est(:, 2), 'r'); legend('真实速度', '估计速度'); xlabel('时间'); ylabel('速度'); ``` 这段代码实现了一个简单的一维运动模型,其中使用了扩展卡滤波来估计位置和速度。代码中的注释会帮助你理解每一步的操作。你可以根据需要进行修改和扩展。希望对你有帮助!

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值