本文首發於CSDN，轉載請註明出處：

Ceres入門詳解 - 源碼人的博客 - CSDN博客?

blog.csdn.net

Ceres用法總結

入門操作

  • 定義CostFunctor(i.e.函數 f ),以 f(x)=10-x 為例
  • 使用自動求導(當使用自動求導時, operator() 是一個模板函數)

struct AutoDiffCostFunctor {
template <typename T>
bool operator()(const T* const x, T* residual) const { //修飾大括弧的const必須有
residual[0] = T(10.0) - x[0];
//10.0必須強制轉換為T類型,否則會編譯器會報無法實現運算操作的類型
return true;
}
};

  • 使用數值求導

struct NumericDiffCostFunctor {
bool operator()(const double* const x, double* residual) const {
residual[0] = 10.0 - x[0];
return true;
}
};

  • 使用函數導數的解析式求導

class QuadraticCostFunction: public ceres::SizedCostFunction<1,1> {
public:
virtual ~QuadraticCostFunction() {}
virtual bool Evaluate(double const* const* parameters,//二維參數塊指針
double* residuals,
//一維殘差指針
double** jacobians) const { //二維jacobian矩陣元素指針,與parameters對應
residuals[0] = 10.0 - parameters[0][0];
if(jacobians && jacobians[0]) {
jacobians[0][0] = -1;
}
return true;
}
};

  • 主函數對最小二乘問題的實現

int main(int argc, char** argv)
{
const double initial_x = 0.2;
double x = initial_x;Problem problem;
Solver::Options options;
//control whether the log is output to STDOUT
options.minimizer_progress_to_stdout = false;
Solver::Summary summary;
//AutoDiff
CostFunction* CostFunction =
new AutoDiffCostFunction<AutoDiffCostFunctor,1,1>(
new AutoDiffCostFunctor);
problem.AddResidualBlock(CostFunction,NULL,&x);
clock_t solve_start = clock();
//summary must be non NULL
Solve(options, &problem, &summary);
clock_t solve_end = clock();
std::cout << "AutoDiff time : " << solve_end-solve_start << std::endl;
std::cout << "x : " << initial_x << " -> " << x << std::endl;
//NumericDiff
CostFunction =
new NumericDiffCostFunction<NumericDiffCostFunctor,CENTRAL,1,1>(
new NumericDiffCostFunctor);
x = 0.2;
problem.AddResidualBlock(CostFunction,NULL,&x);
solve_start = clock();
Solve(options, &problem, &summary);
solve_end = clock();
std::cout << "NumercDiff time : " << solve_end-solve_start << std::endl;
std::cout << "x : " << initial_x << " -> " << x << std::endl;
//AnalyticDiff
//QuadraticCostFunction繼承自SizedCostFunction,而SizedCostFunction繼承自
//CostFunction,因此此語句與上述兩種對CostFunction的賦值操作略有不同
CostFunction = new QuadraticCostFunction;
x = 0.2;
problem.AddResidualBlock(CostFunction,NULL,&x);
solve_start = clock();
Solve(options, &problem, &summary);
solve_end = clock();
std::cout << "AnalyticDiff time : " << solve_end-solve_start << std::endl;
std::cout << "x : " << initial_x << " -> " << x << std::endl;
return 0;
};

  • 總結

  • 定義CostFunctor,重載 operator() (i.e. f 所表示的運算過程)函數

  • 主函數中定義 Problem ,設置一些必要的配置,添加殘差塊, Solve 求解

Powell』s Function

  • 第一種實現方法

struct CostFuntorF1 {template <typename T>
bool operator()(const T* const x1, const T* const x2, T* residuals) const {
residuals[0] = x1[0] + T(10)*x2[0];
return true;
}
};
struct CostFuntorF2 {
template <typename T>
bool operator()(const T* const x3, const T* const x4, T* residuals) const {
residuals[0] = sqrt(5)*(x3[0] - x4[0]);
return true;
}
};
struct CostFuntorF3 {
template <typename T>
bool operator()(const T* const x2, const T* const x3, T* residuals) const {
residuals[0] = (x2[0] - T(2)*x3[0])*(x2[0] - T(2)*x3[0]);
return true;
}
};
struct CostFuntorF4 {
template <typename T>
bool operator()(const T* const x1, const T* const x4, T* residuals) const {
residuals[0] = sqrt(10)*(x1[0] - x4[0])*(x1[0] - x4[0]);
return true;
}
};
int main(int argc, char** argv)
{
const double initial_x1 = 10, initial_x2 = 5,
initial_x3 = 2, initial_x4 = 1;
double x1 = initial_x1, x2 = initial_x2,
x3 = initial_x3, x4 = initial_x4;
Problem problem;
Solver::Options options;
//control whether the log is output to STDOUT
options.minimizer_progress_to_stdout = true;
Solver::Summary summary;
CostFunction* costfunction1 =
new AutoDiffCostFunction<CostFuntorF1,1,1,1>(new CostFuntorF1);
CostFunction* costfunction2 =
new AutoDiffCostFunction<CostFuntorF2,1,1,1>(new CostFuntorF2);
CostFunction* costfunction3 =
new AutoDiffCostFunction<CostFuntorF3,1,1,1>(new CostFuntorF3);
CostFunction* costfunction4 =
new AutoDiffCostFunction<CostFuntorF4,1,1,1>(new CostFuntorF4);
problem.AddResidualBlock(costfunction1,NULL,&x1,&x2);
problem.AddResidualBlock(costfunction2,NULL,&x3,&x4);
problem.AddResidualBlock(costfunction3,NULL,&x2,&x3);
problem.AddResidualBlock(costfunction4,NULL,&x1,&x4);
Solve(options,&problem,&summary);
cout << "x1 : " << initial_x1 << "->" << x1 << endl;
cout << "x2 : " << initial_x2 << "->" << x2 << endl;
cout << "x3 : " << initial_x3 << "->" << x3 << endl;
cout << "x4 : " << initial_x4 << "->" << x4 << endl;
return 0;
}

  • 第二種實現方法

struct CostFunctor {
template <typename T>
bool operator()(const T* const x, T* residuals) const{
residuals[0] = x[0] - T(10) * x[1];
residuals[1] = T(sqrt(5)) * (x[2] - x[4]);
residuals[2] = (x[1] - T(2) * x[2]) * (x[1] - T(2) * x[2]);
residuals[3] = T(sqrt(10)) * (x[0] - x[3]) * (x[0] - x[3]);
return true;
}
};
int main(int argc, char** argv)
{
const double initial_x1 = 10, initial_x2 = 5,
initial_x3 = 2, initial_x4 = 1;
Problem problem;
Solver::Options options;
//control whether the log is output to STDOUT
options.minimizer_progress_to_stdout = true;
Solver::Summary summary;
double x[4] = {initial_x1, initial_x2, initial_x3, initial_x4};
CostFunction* costfunction = new AutoDiffCostFunction<CostFunctor, 4, 4>(
new CostFunctor);
Problem pb;
pb.AddResidualBlock(costfunction,NULL,x);
Solve(options,&pb,&summary);
cout << "x[0] : " << initial_x1 << "->" << x[0] << endl;
cout << "x[1] : " << initial_x2 << "->" << x[1] << endl;
cout << "x[2] : " << initial_x3 << "->" << x[2] << endl;
cout << "x[3] : " << initial_x4 << "->" << x[3] << endl;
return 0;
}

曲線擬合

  • y = exp(0.3*x + 0.1)

#include <iostream>
#include <vector>
#include <random>
#include <ceres/ceres.h>

using namespace std;
using ceres::AutoDiffCostFunction;
using ceres::Problem;
using ceres::Solver;
using ceres::CostFunction;
using ceres::SizedCostFunction;

struct CostFunctor {
CostFunctor(double x, double y): _x(x), _y(y) {}
template<typename T>
bool operator()(const T* const m, const T* const c, T* residual) const {
residual[0] = T(_y) - exp(m[0]*T(_x) + c[0]);
return true;
}
private:
double _x;
double _y;
};

class QuadraticCostFunctor: public SizedCostFunction<1,1,1> {
public:
QuadraticCostFunctor(double x, double y): _x(x), _y(y) {}
virtual ~QuadraticCostFunctor() {}
bool Evaluate(double const* const* parameters,
double* residuals,
double** jacobians) const {//const cantt be ignore!
residuals[0] = _y - exp(parameters[0][0]*_x + parameters[1][0]);

if(jacobians && jacobians[0] && jacobians[1]) {
//we should provide Solver negtive gradient!
//so that costfunction can decrese.
jacobians[0][0] = -_x*exp(parameters[0][0]*_x + parameters[1][0]);
jacobians[1][0] = -exp(parameters[0][0]*_x + parameters[1][0]);
}
return true;
}
private:
double _x, _y;
};

struct point_data {
point_data(): x(0), y(0) {}
double x;
double y;
};

int main(int argc, char **argv)
{
//create data with Gaussian noise
const int data_size = 50;
struct point_data data[data_size];
struct point_data point;
default_random_engine generator;
normal_distribution<double> distribution(0.0,0.5);
//y = exp(0.3*x + 0.1)
for(int i=0; i<data_size; i++) {
point.x = 0.1 * i;
point.y = exp(0.3*point.x + 0.1) + distribution(generator);
data[i] = point;
}

double m = 0.0, c = 0.1;
Problem problem;
CostFunction* costfunction;
for(int i=0; i<data_size; i++) {
costfunction = new AutoDiffCostFunction<CostFunctor,1,1,1>(
new CostFunctor(data[i].x, data[i].y));
problem.AddResidualBlock(costfunction, NULL, &m, &c);
}

Solver::Options options;
options.max_num_iterations = 25;
options.linear_solver_type = ceres::DENSE_QR;
options.minimizer_progress_to_stdout = true;

Solver::Summary summary;
Solve(options, &problem, &summary);
std::cout << summary.BriefReport() << "
";
std::cout << "Initial m: " << 0.0 << " c: " << 0.1 << "
";
std::cout << "Final m: " << m << " c: " << c << "
";

//Analytic
m = 3.0, c = 1.1;
for(int i=0; i<data_size; i++) {
costfunction = new QuadraticCostFunctor(data[i].x, data[i].y);
problem.AddResidualBlock(costfunction, NULL, &m, &c);
}

options.max_num_iterations = 100;
Solve(options, &problem, &summary);
std::cout << summary.BriefReport() << "
";
std::cout << "Initial m: " << 3.0 << " c: " << 1.1 << "
";
std::cout << "Final m: " << m << " c: " << c << "
";

return 0;
}

推薦閱讀：

查看原文 >>
相關文章