爽一把手写Bundle Adjustment

本文介绍了在视觉SLAM中的重要技术Bundle Adjustment(BA),它是一个大规模最小二乘优化问题,用于优化相机位姿和路标点的估计。通过Levenberg-Marquardt优化算法,详细讲解了如何构建和求解增量方程,最后提供了使用Ceres优化库的示例和完整的BA代码链接。
摘要由CSDN通过智能技术生成

今天是传说中的5.20,同时也是本人最近特别钟爱的神剧《权力的游戏》大结局的日子,作为一名伪技术宅,这么重大的日子当然要写一篇文来纪念一下。所以就有了这篇手写Bundle Adjustment的记录文。

Bundle Adjustment(后面简写成BA)在视觉SLAM里面的作用不用多说,自然是非常重要的。BA主要求解的是一个关于相机位姿以及路标点的空间坐标的大规模最小二乘优化问题,每一个相机位姿都与多个路标点关联,而每一个路标点也同时和多个相机位姿相关联,由此形成复杂的网状的相互约束,通过调整相机位姿和路标点的位置估计,使得所有路标点在所有相机中的重投影误差最小。下面我们就来亲自动手尝试一下求解Bundle Adjustment。相关的理论可以参考《视觉SLAM十四讲》第六章和第十章,这里仅简单略过。

1  Levenberg-Marquardt优化

考虑简单的最小二乘优化问题:

  \min_{\textbf{x}}{\frac{1}{2}\left \| f(\textbf{x}) \right \|_2^2}

f(\textbf{x})进行一阶泰勒展开:

  f(\textbf{x})\approx f(\textbf{x})+J(\textbf{x})\Delta \textbf{x}

同时通过拉格朗日乘子添加对增量的约束,得到以下形式:

  \min_{\textbf{x}}{\frac{1}{2}\left \| f(\textbf{x})+J(\textbf{x})\Delta \textbf{x} \right \|_2^2+\frac{\lambda}{2}\left \| D\Delta \textbf{x} \right \|^2}

\Delta \textbf{x}求导并令导数等于0,可得以下方程

(J(\textbf{x})^TJ(\textbf{x})+\lambda D^TD)\Delta \textbf{x}=-J(\textbf{x})^Tf(\textbf{x})

求解以上增量方程,可得\Delta \textbf{x}的值,进而令\textbf{x}=\textbf{x}+\Delta \textbf{x}更新对优化变量\textbf{x}的估计。

在实际应用中通常取D=I或者D^TD=diag(J(\textbf{x})^TJ(\textbf{x})),同时还涉及到对\lambda大小的调整。

在应用Levenberg-Marquardt优化算法时,最主要的步骤就是求函数值f(\textbf{x})和雅可比J(\textbf{x}),以及求解增量方程得出\Delta \textbf{x}

下面尝试一下Ceres优化库的一个曲线拟合的实例y=e^{mx+c}。拟合用的数据点通过对真实曲线y=e^{0.3x+0.1}进行采样并附加标准差\sigma=0.2的高斯噪声生成。

#include<random>
#include<vector>

//用来生成高斯噪声
std::default_random_engine e;
std::normal_distribution<double> n(0, 0.2);
//用来存储数据点
std::vector<std::pair<double, double>> data;
//生成数据点,假设x的范围是0-5
for(double x = 0; x < 5; x += 0.05)
{
	double y = exp(0.3*x + 0.1)+n(e);
	data.push_back(std::make_pair(x, y));
}

根据以上数据点拟合曲线的最小二乘优化问题为

  \min_{m,c}{\frac{1}{2}\sum_k (e^{mx_k+c}-y_k)^2}=\min_{m,c}{\frac{1}{2}\left \| \begin{bmatrix}e^{mx_1+c}-y_1\\e^{mx_2+c}-y_2\\\dots\\e^{mx_N+c}-y_N \end{bmatrix} \right \|^2_2}=\min_{m,c}{\frac{1}{2}\left \| F(m,c) \right \|^2_2}

优化的目标是求得最佳的mc的值,使得\frac{1}{2}\left \| F(m,c) \right \|^2_2值最小。

Bundle Adjustment 是一种用于同时估计多个相机位置,姿态和三维点云的优化算法。下面是一个用Python编的简单的Bundle Adjustment玩具程序: ```python import numpy as np import scipy.sparse as sp import scipy.sparse.linalg as splinalg # 生成测试数据 num_cameras = 3 num_points = 5 num_observations = 10 noise_scale = 0.1 # 相机的位置和姿态 camera_positions = np.random.randn(num_cameras, 3) camera_orientations = np.random.randn(num_cameras, 3) # 三维点云 points_3d = np.random.randn(num_points, 3) # 观测到的2D点 observed_points = np.zeros((num_observations, 2)) camera_indices = np.zeros(num_observations, dtype=int) point_indices = np.zeros(num_observations, dtype=int) for i in range(num_observations): camera_idx = np.random.randint(num_cameras) point_idx = np.random.randint(num_points) observed_points[i] = np.random.randn(2) + \ np.dot(camera_orientations[camera_idx], points_3d[point_idx]) + \ camera_positions[camera_idx] observed_points[i] += np.random.randn(2) * noise_scale camera_indices[i] = camera_idx point_indices[i] = point_idx # 构建稀疏矩阵 A = sp.lil_matrix((num_observations * 2, num_cameras * 6 + num_points * 3)) b = np.zeros(num_observations * 2) for i in range(num_observations): camera_idx = camera_indices[i] point_idx = point_indices[i] point_3d = points_3d[point_idx] camera_pos = camera_positions[camera_idx] camera_orient = camera_orientations[camera_idx] # 计算重投影误差 projected_point = np.dot(camera_orient, point_3d) + camera_pos projected_point /= projected_point[2] projected_point = projected_point[:2] error = observed_points[i] - projected_point # 构建雅可比矩阵 J_camera = np.hstack((point_3d, np.zeros(3), np.eye(3))) J_point = np.hstack((np.zeros(2), camera_orient[:2].reshape(2, 1) * point_3d.reshape(1, 3))) J = np.vstack((J_camera, J_point)) A[i * 2: (i+1) * 2, camera_idx * 6: (camera_idx+1) * 6] = J[:, :6] A[i * 2: (i+1) * 2, num_cameras * 6 + point_idx * 3: num_cameras * 6 + (point_idx+1) * 3] = J[:, 6:] b[i * 2: (i+1) * 2] = error # 优化 x0 = np.hstack((camera_positions.ravel(), camera_orientations.ravel(), points_3d.ravel())) res = splinalg.lsmr(A, b, x0=x0) optimized_positions = res[0][:num_cameras * 3].reshape(num_cameras, 3) optimized_orientations = res[0][num_cameras * 3: num_cameras * 6].reshape(num_cameras, 3) optimized_points_3d = res[0][num_cameras * 6:].reshape(num_points, 3) print("Original Camera Positions:\n", camera_positions) print("Original Camera Orientations:\n", camera_orientations) print("Original 3D Points:\n", points_3d) print("Optimized Camera Positions:\n", optimized_positions) print("Optimized Camera Orientations:\n", optimized_orientations) print("Optimized 3D Points:\n", optimized_points_3d) ``` 这个程序生成了一些随机的相机位置、姿态和三维点云,然后随机生成了一些观测到的2D点,并添加了一些高斯噪声。程序使用Bundle Adjustment来估计相机和3D点云的位置,并输出优化后的结果。
评论 5
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值