当前位置: 首页 > 工具软件 > Ceres Solver > 使用案例 >

Ceres使用

曾永新
2023-12-01

定义:ceres是一款非线性优化问题的数值求解器。

手册地址:
http://www.ceres-solver.org/nnls_modeling.html#_CPPv2N5ceres29AutoDiffLocalParameterizationE

1. Problem

Ceres的求解过程包括构建最小二乘和求解最小二乘问题两部分,其中构建最小二乘问题的相关方法均包含在Ceres::Problem类中,涉及的成员函数主要包括Problem::AddResidualBlock()和Problem::AddParameterBlock()。

1.1 Problem::AddResidualBlock( )

AddResidualBlock()顾名思义主要用于向Problem类传递残差模块的信息,函数原型如下,传递的参数主要包括代价函数模块、损失函数模块和参数模块。

ResidualBlockId Problem::AddResidualBlock(CostFunction *cost_function, 
										  LossFunction *loss_function,
										  double *x0, double *x1, ...)

注:x0、x1为估计参数

1.1.1 CostFunction

代价函数:包含了参数模块的维度信息,内部使用仿函数定义误差函数的计算方式。AddResidualBlock( )函数会检测传入的参数模块是否和代价函数模块中定义的维数一致,维度不一致时程序会强制退出。

ceres提供了许多种CostFunction模板,较为常用的包括以下三种:

  • 自动导数(AutoDiffCostFunction):由ceres自行决定导数的计算方式,最常用的求导方式。
  • 数值导数(NumericDiffCostFunction):由用户手动编写导数的数值求解形式,通常在残差函数的计算使用无法直接调用的库函数,导致调用AutoDiffCostFunction类构建时使用;但手动编写的精度和计算效率不如模板类,因此不到不得已,官方并不建议使用该方法。
  • 解析导数(Analytic Derivatives):当导数存在闭合解析形式时使用,用于可基于CostFunciton基类自行编写;但由于需要自行管理残差和雅克比矩阵,除非闭合解具有具有明显的精度和效率优势,否则同样不建议使用。
    AutoDiffCostFunction为模板类,构造函数如下:
ceres::AutoDiffCostFunction<CostFunctor, int residualDim, int paramDim>(CostFunctor* functor);

模板参数依次为仿函数(functor)类型CostFunctor,残差维数residualDim和待优化变量维数paramDim,接受参数类型为仿函数指针CostFunctor*。

ceres::CostFunction* cost_function=new ceres::AutoDiffCostFunction<
											ReprojectionError3D, 2, 4, 3, 3>(new ReprojectionError3D(observed_x,observed_y));
problem.AddResidualBlock(cost_function, NULL, Q, T, Position); 	

仿函数CostFunctor

仿函数的本质为结构体struct或者类class,由于重载了()运算符,使得其能够具有和函数一样的调用行为,因此被称为仿函数。ceres中采用仿函数来表示残差的计算过程。

  • 重载操作符()(必有)
    操作符()是一个模板方法,返回值为bool型,接受参数依次为待优化变量和残差变量,且待优化变量的传入方式应和Probelm::AddResidualBlock()一致。
struct ReprojectionError3D
{
	ReprojectionError3D(double _observed_u, double _observed_v)
		:observed_u(_observed_u), observed_v(_observed_v)
		{}
	template <typename T>
	bool operator()(const T* const camera_R, const T* const camera_T, const T* point, T* residuals) const  
	{
		T p[3];
		ceres::QuaternionRotatePoint(camera_R, point, p);
		p[0] += camera_T[0]; p[1] += camera_T[1]; p[2] += camera_T[2];
		T xp = p[0] / p[2];
    	T yp = p[1] / p[2];
    	residuals[0] = xp - T(observed_u);
    	residuals[1] = yp - T(observed_v);
    	return true;
	}
	static ceres::CostFunction* Create(const double observed_x,const double observed_y) 
	{
	  return (new ceres::AutoDiffCostFunction<
	          ReprojectionError3D, 2, 4, 3, 3>(
	          	new ReprojectionError3D(observed_x,observed_y)));
	}
	double observed_u;
	double observed_v;
};

1.1.2 LossFunction

损失函数:用于处理参数中含有野值的情况,避免错误量测对估计的影响,常用参数包括HuberLoss、CauchyLoss等;该参数可以取NULL或nullptr,此时损失函数为单位函数。

1.2 Problem::AddParameterBlock( )

用户在调用AddResidualBlock( )时其实已经隐式地向Problem传递了参数模块,但在一些情况下,需要用户显示地向Problem传入参数模块(通常出现在需要对优化参数进行重新参数化的情况)。Ceres提供了Problem::AddParameterBlock( )函数用于用户显式传递参数模块:

void Problem::AddParameterBlock(double *values, int size)

void Problem::AddParameterBlock(double *values, int size, LocalParameterization *local_parameterization)

注:values表示优化变量,size表示优化变量的维度。
其中,第一种函数原型除了会增加一些额外的参数检查之外,功能上和隐式传递参数并没有太大区别。第二种函数原型则会额外传入LocalParameterization参数,用于重构优化参数的维数,这里我们着重讲解一下LocalParameterization类。

1.2.1 LocalParameterization

LocalParameterization是在优化Manifold(流形)上的变量时需要考虑的,Manifold上变量是过参数的,即Manifold上变量的维度大于其自由度。这会导致Manifold上变量各个量之间存在约束,如果直接对这些量求导、优化,那么这就是一个有约束的优化,实现困难。为了解决这个问题,在数学上对Manifold在当前变量值处形成的切空间求导,在切空间上优化,最后投影回Manifold。

对于SLAM问题,广泛遇到的Manifold是旋转,旋转仅需要3个量,但实际运用中涉及到万向锁问题,在更高维空间表达旋转,四元数就是在维度4表达3个自由度的三维空间的旋转。

bool ComputeJacobian()计算得到一个4*3的矩阵(global_to_local),含义是Manifold上变量对Tangent Space上变量的导数,在ceres::CostFunction处提供residuals对Manifold上变量的倒数,乘以这个矩阵,之后变就变成了对Tangent Space上变量的导数。

problem.AddParameterBlock(quaternion, 4);// 直接传递4维参数

ceres::LocalParameterization* local_param = new ceres::QuaternionParameterization();
problem.AddParameterBlock(quaternion, 4, local_param)//重构参数,优化时实际使用的是3维的等效旋转矢量

除了上面提到的QuaternionParameterization外,ceres还提供下述预定义LocalParameterization子类,具体可查手册。

1.2.2 自定义LocalParameterization

LocalParaneterization本身是一个虚基类,详细定义如下。用户可以自行定义自己需要使用的子类,或使用Ceres预先定义好的子类。

1.3 其他成员函数

Probelm还提供了其他关于ResidualBlock和ParameterBlock的函数,例如获取模块维数、判断是否存在模块、存在的模块数目等,这里只列出几个比较重要的函数,完整的列表参见ceres API:

// 设定对应的参数模块在优化过程中保持不变
void Problem::SetParameterBlockConstant(double *values)
// 设定对应的参数模块在优化过程中可变
void Problem::SetParameterBlockVariable(double *values)
// 设定优化下界
void Problem::SetParameterLowerBound(double *values, int index, double lower_bound)
// 设定优化上界
void Problem::SetParameterUpperBound(double *values, int index, double upper_bound)
// 该函数紧跟在参数赋值后,在给定的参数位置求解Problem,给出当前位置处的cost、梯度以及Jacobian矩阵;
bool Problem::Evaluate(const Problem::EvaluateOptions &options, 
					   double *cost, vector<double>* residuals, 
					   vector<double> *gradient, CRSMatrix *jacobian)

2. Solver

选取合适的求解器

2.1 Solver::Options

Ceres的参数主要有三类,一类通用参数,比如迭代次数什么的;第二类是和优化算法的参数;第三类是和线性求解器(在信任域算法中被使用)有关的参数。
常用通用参数如下:

  • options.max_solver_time_in_seconds
    默认值:1e6
    最长运行时间,单位为秒。
  • options.linear_solver_type
    线性求解器的类型,用于计算Levenberg-Marquardt算法每次迭代中线性最小二乘问题的解。 如果编译Ceres时加入了SuiteSparse或CXSparse或Eigen的稀疏Cholesky分解选项,则默认为SPARSE_NORMAL_CHOLESKY,否则为DENSE_QR。

2.2 Solver::Summary

输出优化过程及结果

  • termination_type
    求解器终止的原因:
enum TerminationType {
  CONVERGENCE,
  NO_CONVERGENCE,
  FAILURE,
  USER_SUCCESS,
  USER_FAILURE
};
  • final_cost
    优化后的最终代价

2.3 使用

ceres::Solve(options, &problem, &summary);
 类似资料: