[代码]扩展卡尔曼滤波ekf实现

Abstract: 之前写了数学原理,在线标定要用到这个算法,找了些资料实现了一下,这是pitch和yaw角的ekf。


原理可见:[数学]扩展卡尔曼滤波EKF

矩阵的定义

协方差,运动方程,观测方程等

# 状态协方差矩阵 P
self.P = np.diag([0.1, 0.1, 0.5, 0.5])

# 过程噪声协方差 Q
self.Q_base = np.diag([0.001, 0.001, 0.01, 0.01])

# 测量噪声协方差 R
self.R_base = np.diag([0.05**2, 0.03**2])  # sigma_pitch=0.005, sigma_yaw=0.003

# 状态转移矩阵 F
self.F = np.array([[1, 0, dt, 0],
                    [0, 1, 0, dt],
                    [0, 0, 1, 0],
                    [0, 0, 0, 1]])

# 观测矩阵 H
self.H = np.array([[1, 0, 0, 0],
                    [0, 1, 0, 0]])

预测步

与推导一致

def predict(self):
    """预测步"""
    self.x = self.F @ self.x
    self.P = self.F @ self.P @ self.F.T + self.Q_base
    return self.x.copy()

更新

def update(self, z_meas, confidence=1.0):
    """
    更新步
    """
    # 计算新息
    z_pred = self.H @ self.x
    # y是先验观测值减掉通过运动方程转换来的先验状态值
    y = np.array([z_meas[0] - z_pred[0], z_meas[1] - z_pred[1]])
    
    # 角度包裹处理
    y[0] = self.normalize_angle(y[0])
    y[1] = self.normalize_angle(y[1])
    
    # 自适应噪声调整
    R = self.R_base / confidence if self.use_adaptive else self.R_base
    
    # 计算新息协方差
    S = self.H @ self.P @ self.H.T + R
    
    # 卡尔曼增益
    K = self.P @ self.H.T @ np.linalg.inv(S)
    
    # 状态更新
    self.x = self.x + K @ y
    
    # 协方差更新
    I = np.eye(4)
    self.P = (I - K @ self.H) @ self.P @ (I - K @ self.H).T + K @ R @ K.T

模拟

输入十帧的模拟数据:

timestamps = np.array([0.0, 0.033, 0.067, 0.1, 0.133, 0.167, 0.2, 0.233, 0.267, 0.3])

# 观测到的pitch和yaw(单位:弧度),带有噪声
pitch_obs = np.array([0.05, 0.048, 0.045, 0.042, 0.04, 0.038, 0.037, 0.036, 0.035, 0.034])
yaw_obs = np.array([-0.01, -0.008, -0.005, -0.003, 0.0, 0.002, 0.003, 0.004, 0.005, 0.006])

输出:

ekf


Last modified on 2026-03-18