使用G2O解决优化问题的简单例子

从零开始一起学习SLAM | 掌握g2o顶点编程套路 点“计算机视觉life”关注,置顶更快接收消息! 小白:师兄,上一次将的g2o框架《从零开始一起学习SLAM | 理解图优化,一步步带你看懂g2o代码》真的很清晰,我现在再去看g2o的那些优化的部分,基本都能看懂了呢! 师兄:那太好啦,以后多练习练习,加深理解 小白:嗯,我开始编程时,发现g2o的顶点和边的定义也非常复杂,光看十四讲里面,就有好几种不同的定义,完全懵圈状态。。。师兄,能否帮我捋... 阅读详情

优化问题描述

假设一个机器人初始起点在0处,然后机器人向前移动,通过编码器测得它向前移动了1m,到达第二个地点x1。接着,又向后返回,编码器测得它向后移动了0.8米。但是,通过闭环检测,发现它回到了原始起点。可以看出,编码器误差导致计算的位姿和观测到有差异,那机器人这几个状态中的位姿到底是怎么样的才最好的满足这些条件呢?

首先构建位姿之间的关系,即图的边:

线性方程组中变量小于方程的个数,要计算出最优的结果,使出杀手锏最小二乘法。先构建残差平方和函数:

为了使残差平方和最小,我们对上面的函数每个变量求偏导,并使得偏导数等于0.

整理得到:

接着矩阵求解线性方程组:

所以调整以后为满足这些边的条件,机器人的位姿为:

可以看到,因为有了x2->x0的闭环检测,几个位置点间才能形成互相约束。这对使用最小二乘解决该优化问题起到了决定性的作用。

该问题描述来源于:https://heyijia.blog.youkuaiyun.com/article/details/47686523

下面利用G2O来解上面的问题,以便理解如何使用G2O。

定义顶点

在该问题中,一个位置点就是图优化中的一个顶点。一个顶点可以包含多个需优化量。比如二维环境下的机器人位置一般是3维的(x,y,theta),即一个顶点有三个需要优化的量。

在此问题中,我们只需优化求解一个一维的距离值。即是,一个顶点只包含一个需优化量。

// 曲线模型的顶点,模板参数:优化变量维度和数据类型
class SimpleVertex : public g2o::BaseVertex<1, double> {
public:
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW
  SimpleVertex(bool fixed = false)
  {
    setToOriginImpl();
    setFixed(fixed);
  }

  SimpleVertex(double e, bool fixed = false)
  {
    _estimate = e;
    setFixed(fixed);
  }

  // 重置
  virtual void setToOriginImpl() override {
    _estimate = 0;
  }

  // 更新
  virtual void oplusImpl(const double *update) override {
    _estimate += *update;
  }

  // 存盘和读盘:留空
  virtual bool read(istream &in) {return true;}

  virtual bool write(ostream &out) const {return true;}
};

上面代码中g2o::BaseVertex<1, double>就表示了优化的量是一维的,数据类型是double。

如果查看TEB中设置的优化量,可以发现它是这样写的:

g2o::BaseVertex<3, PoseSE2 >

TEB中的优化量是三维的,即机器人的位姿(x,y,theta)。此处的PoseSE2便是描述机器人位姿的结构体。

自定义一个顶点通常需要重新实现setToOriginImpl()和oplusImpl(const double *update)这两个函数。至于read(istream &in)和write(ostream &out)是为了保存和读取当前顶点数据的。用不上可留空。

setToOriginImpl()函数用于重置顶点的数据。顶点包含的数据变量是_estimate。该变量的类型即是g2o::BaseVertex<1, double>中设置的double。该函数正是用于重置_estimate和使顶点恢复默认状态。

oplusImpl(const double *update)函数用于叠加优化量的步长。注意有时候这样的叠加操作并不是线性的。比如2D-SLAM中位置步进叠加操作则不是线性的。

在此问题中我们是直接线性叠加一维的距离值。

TEB中的位置叠加也是线性叠加的。放置在下面以作参考。

void plus(const double* pose_as_array)
{
    _position.coeffRef(0) += pose_as_array[0];
    _position.coeffRef(1) += pose_as_array[1];
    _theta = g2o::normalize_theta( _theta + pose_as_array[2] );
}

virtual void oplusImpl(const double* update) override
{
    _estimate.plus(update);
}

另外,顶点是可以设置成固定的。当不需要变动某个顶点时,使用setFixed函数来固定。通常,一个优化问题中,至少需要固定一个顶点,否则所有的顶点都在浮动,优化效果也不会好。

定义边

边即是顶点之间的约束。在该问题中就是两位置间的测量值和观测值之间的差值要趋近于0。这里需要定义的边其实就是下面的等式。

实际上,G2O中边有三种:一元边(g2o::BaseUnaryEdge),二元边(g2o::BaseBinaryEdge)和多元边(g2o::BaseMultiEdge)。

g2o::BaseUnaryEdge<D, E, VertexXi>
g2o::BaseBinaryEdge<D, E, VertexXi, VertexXj>
g2o::BaseMultiEdge<D, E>

不同类型的边有不同数量的模板参数。其中D 是 int 型,表示误差值的维度 (dimension), E 表示测量值的数据类型(即_measurement的类型), VertexXi,VertexXj 分别表示不同顶点的类型。当D为2时,_error的类型变为Eigen::Vector2d,当D为3时,_error的类型变为Eigen::Vector3d。

在上面的约束中,有一个一元边(f1)和三个二元边(f2,f3,f4)。在G2O中可如下定义:

// 误差模型 模板参数:测量值维度,测量值类型,连接顶点类型
class SimpleUnaryEdge : public g2o::BaseUnaryEdge<1, double, SimpleVertex> {
public:
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW

  SimpleUnaryEdge() : BaseUnaryEdge() {}

  // 计算模型误差
  virtual void computeError() override {
    const SimpleVertex *v = static_cast<const SimpleVertex *> (_vertices[0]);
    const double abc = v->estimate();
    _error(0, 0) = _measurement - abc;
  }

  // 计算雅可比矩阵
  virtual void linearizeOplus() override {
    _jacobianOplusXi[0] = -1;//一元边只有一个顶点,所以只有_jacobianOplusXi。_jacobianOplusXi的维度和顶点需优化的维度是一致的
  }

  virtual bool read(istream &in) {}

  virtual bool write(ostream &out) const {}

};

// 误差模型 模板参数:测量值维度,测量值类型,连接顶点类型
class SimpleBinaryEdge : public g2o::BaseBinaryEdge<1, double, SimpleVertex, SimpleVertex> {
public:
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW

  SimpleBinaryEdge() {}

  // 计算模型误差
  virtual void computeError() override {
    const SimpleVertex *v1 = static_cast<const SimpleVertex *> (_vertices[0]);
    const SimpleVertex *v2 = static_cast<const SimpleVertex *> (_vertices[1]);
    const double abc1 = v1->estimate();
    const double abc2 = v2->estimate();
    _error[0] = _measurement - (abc1 - abc2);
  }

  // 计算雅可比矩阵
  virtual void linearizeOplus() override {
    _jacobianOplusXi(0,0) = -1;//二元边有两个顶点,偏差_error对每个顶点的偏导数都需要求解;如果_error是多维的,则每一维的偏导都需要求解。即会出现_jacobianOplusXi(1,0)

    _jacobianOplusXj(0,0) = 1;//因为误差_error只有一维,所以_jacobianOplusXi只有_jacobianOplusXi(0,0)项。此时_jacobianOplusXi[0]与_jacobianOplusXi(0,0)效果等价。
  }

  virtual bool read(istream &in) {}

  virtual bool write(ostream &out) const {}

};

因为我们这里测量值是一维的距离值(double类型)而顶点类型为SimpleVertex,所以边需设置成g2o::BaseUnaryEdge<1, double, SimpleVertex>和g2o::BaseBinaryEdge<1, double, SimpleVertex, SimpleVertex>类型。

自定义边同样需要特别关注两个函数computeError()和linearizeOplus()。

computeError()是用于计算迭代误差的。顶点间的约束正是由误差计算函数构建的。优化时误差项将逐步趋近于0。_error的维度和类型通常由构建的模型决定。比如该问题中误差为距离误差。

TEB的速度边构建的误差项为线速度与最大线速度的差值和角速度与最大角速度的差值。所以这里的误差项是二维的。

  void computeError()
  {
    TEB_ASSERT_MSG(cfg_, "You must call setTebConfig on EdgeVelocity()");
    const VertexPose* conf1 = static_cast<const VertexPose*>(_vertices[0]);
    const VertexPose* conf2 = static_cast<const VertexPose*>(_vertices[1]);
    const VertexTimeDiff* deltaT = static_cast<const VertexTimeDiff*>(_vertices[2]);
    
    const Eigen::Vector2d deltaS = conf2->estimate().position() - conf1->estimate().position();
    
    double dist = deltaS.norm();
    const double angle_diff = g2o::normalize_theta(conf2->theta() - conf1->theta());
    if (cfg_->trajectory.exact_arc_length && angle_diff != 0)
    {
        double radius =  dist/(2*sin(angle_diff/2));
        dist = fabs( angle_diff * radius ); // actual arg length!
    }
    double vel = dist / deltaT->estimate();
    
//     vel *= g2o::sign(deltaS[0]*cos(conf1->theta()) + deltaS[1]*sin(conf1->theta())); // consider direction
    vel *= fast_sigmoid( 100 * (deltaS.x()*cos(conf1->theta()) + deltaS.y()*sin(conf1->theta())) ); // consider direction
    
    const double omega = angle_diff / deltaT->estimate();
  
    _error[0] = penaltyBoundToInterval(vel, -cfg_->robot.max_vel_x_backwards, cfg_->robot.max_vel_x,cfg_->optim.penalty_epsilon);
    _error[1] = penaltyBoundToInterval(omega, cfg_->robot.max_vel_theta,cfg_->optim.penalty_epsilon);

    TEB_ASSERT_MSG(std::isfinite(_error[0]), "EdgeVelocity::computeError() _error[0]=%f _error[1]=%f\n",_error[0],_error[1]);
  }

linearizeOplus()函数里主要是配置雅克比矩阵。当然,G2O是支持自动求导的,该函数可以不实现。优化时由G2O自动处理。但准确的实现可加快优化计算的速度。下面介绍雅克比矩阵该如何计算。

计算雅克比矩阵

雅克比矩阵存储了误差项的每一维相对于顶点各优化成员的偏导数。

一元边

一元边只有一个顶点所以只定义了_jacobianOplusXi,其对应顶点类型VertexXi。注意,_jacobianOplus变量才完整地描述雅克比矩阵,_jacobianOplusXi只是_jacobianOplus中相关部分的引用,以方便配置雅克比矩阵。

// This could be a simple using statement, but in multiple places
// _jacobianOplusXi is used.
template <int D, typename E, typename VertexXi>
class BaseUnaryEdge : public BaseFixedSizedEdge<D, E, VertexXi> {
 public:
  using VertexXiType = VertexXi;
  BaseUnaryEdge() : BaseFixedSizedEdge<D, E, VertexXi>(){};

 protected:
  typename BaseFixedSizedEdge<D, E, VertexXi>::template JacobianType<
      D, VertexXi::Dimension>& _jacobianOplusXi =
      std::get<0>(this->_jacobianOplus);
};

本优化问题的一元边按下面代码的方式计算雅克比矩阵。

// 误差模型 模板参数:误差项维度,测量值类型,连接顶点类型
class SimpleUnaryEdge : public g2o::BaseUnaryEdge<1, double, SimpleVertex> {
public:
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW

  SimpleUnaryEdge() : BaseUnaryEdge() {}

  // 计算曲线模型误差
  virtual void computeError() override {
    const SimpleVertex *v = static_cast<const SimpleVertex *> (_vertices[0]);
    const double abc = v->estimate();
    _error(0, 0) = _measurement - abc;
  }

  // 计算雅可比矩阵
  virtual void linearizeOplus() override {
    _jacobianOplusXi[0] = -1;//一元边只有一个顶点,所以只有_jacobianOplusXi。_jacobianOplusXi的维度和顶点需优化的维度是一致的
  }//此处就是求误差相对于abc变量的导数

  virtual bool read(istream &in) {}

  virtual bool write(ostream &out) const {}

};
二元边

二元边有两个顶点,所以定义了_jacobianOplusXi和_jacobianOplusXj。

// This could be a simple using statement, but in multiple places
// _jacobianOplusXi and _jacobianOplusYi are used.
template <int D, typename E, typename VertexXi, typename VertexXj>
class BaseBinaryEdge : public BaseFixedSizedEdge<D, E, VertexXi, VertexXj> {
 public:
  using VertexXiType = VertexXi;
  using VertexXjType = VertexXj;
  BaseBinaryEdge() : BaseFixedSizedEdge<D, E, VertexXi, VertexXj>(){};

 protected:
  typename BaseFixedSizedEdge<D, E, VertexXi, VertexXj>::template JacobianType<
      D, VertexXi::Dimension>& _jacobianOplusXi =
      std::get<0>(this->_jacobianOplus);
  typename BaseFixedSizedEdge<D, E, VertexXi, VertexXj>::template JacobianType<
      D, VertexXj::Dimension>& _jacobianOplusXj =
      std::get<1>(this->_jacobianOplus);
};

本优化问题的二元边按下面代码的方式计算雅克比矩阵。

// 误差模型 模板参数:误差项维度,测量值类型,连接顶点类型
class SimpleBinaryEdge : public g2o::BaseBinaryEdge<1, double, SimpleVertex, SimpleVertex> {
public:
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW

  SimpleBinaryEdge() {}

  // 计算曲线模型误差
  virtual void computeError() override {
    const SimpleVertex *v1 = static_cast<const SimpleVertex *> (_vertices[0]);
    const SimpleVertex *v2 = static_cast<const SimpleVertex *> (_vertices[1]);
    const double abc1 = v1->estimate();
    const double abc2 = v2->estimate();
    _error[0] = _measurement - (abc1 - abc2);
  }

  // 计算雅可比矩阵
  virtual void linearizeOplus() override {
    _jacobianOplusXi(0,0) = -1;//二元边有两个顶点,偏差_error对每个顶点的偏导数都需要求解;如果_error是多维的,则每一维的偏导都需要求解。即会出现_jacobianOplusXi(1,0)

    _jacobianOplusXj(0,0) = 1;//因为误差_error只有一维,所以_jacobianOplusXi只有_jacobianOplusXi(0,0)项。此时_jacobianOplusXi[0]与_jacobianOplusXi(0,0)效果等价。
  }

  virtual bool read(istream &in) {}

  virtual bool write(ostream &out) const {}

};
多元边

多元边直接配置_jacobianOplus变量。

下面的示例代码是teb中的速度约束,来源于edge_velocity.h文件。

class EdgeVelocity : public BaseTebMultiEdge<2, double>
{
public:
  
  /**
   * @brief Construct edge.
   */	      
  EdgeVelocity()
  {
    this->resize(3); // Since we derive from a g2o::BaseMultiEdge, set the desired number of vertices
  }
  
  /**
   * @brief Actual cost function
   */  
  void computeError()
  {
    TEB_ASSERT_MSG(cfg_, "You must call setTebConfig on EdgeVelocity()");
    const VertexPose* conf1 = static_cast<const VertexPose*>(_vertices[0]);
    const VertexPose* conf2 = static_cast<const VertexPose*>(_vertices[1]);
    const VertexTimeDiff* deltaT = static_cast<const VertexTimeDiff*>(_vertices[2]);
    
    const Eigen::Vector2d deltaS = conf2->estimate().position() - conf1->estimate().position();
    
    double dist = deltaS.norm();
    const double angle_diff = g2o::normalize_theta(conf2->theta() - conf1->theta());
    if (cfg_->trajectory.exact_arc_length && angle_diff != 0)
    {
        double radius =  dist/(2*sin(angle_diff/2));
        dist = fabs( angle_diff * radius ); // actual arg length!
    }
    double vel = dist / deltaT->estimate();
    
//     vel *= g2o::sign(deltaS[0]*cos(conf1->theta()) + deltaS[1]*sin(conf1->theta())); // consider direction
    vel *= fast_sigmoid( 100 * (deltaS.x()*cos(conf1->theta()) + deltaS.y()*sin(conf1->theta())) ); // consider direction
    
    const double omega = angle_diff / deltaT->estimate();
  
    _error[0] = penaltyBoundToInterval(vel, -cfg_->robot.max_vel_x_backwards, cfg_->robot.max_vel_x,cfg_->optim.penalty_epsilon);
    _error[1] = penaltyBoundToInterval(omega, cfg_->robot.max_vel_theta,cfg_->optim.penalty_epsilon);

    TEB_ASSERT_MSG(std::isfinite(_error[0]), "EdgeVelocity::computeError() _error[0]=%f _error[1]=%f\n",_error[0],_error[1]);
  }

#ifdef USE_ANALYTIC_JACOBI
#if 0 //TODO the hardcoded jacobian does not include the changing direction (just the absolute value)
      // Change accordingly...

  /**
   * @brief Jacobi matrix of the cost function specified in computeError().
   */
  void linearizeOplus()
  {
    TEB_ASSERT_MSG(cfg_, "You must call setTebConfig on EdgeVelocity()");
    const VertexPose* conf1 = static_cast<const VertexPose*>(_vertices[0]);
    const VertexPose* conf2 = static_cast<const VertexPose*>(_vertices[1]);
    const VertexTimeDiff* deltaT = static_cast<const VertexTimeDiff*>(_vertices[2]);
    
    Eigen::Vector2d deltaS = conf2->position() - conf1->position();
    double dist = deltaS.norm();
    double aux1 = dist*deltaT->estimate();
    double aux2 = 1/deltaT->estimate();
    
    double vel = dist * aux2;
    double omega = g2o::normalize_theta(conf2->theta() - conf1->theta()) * aux2;
    
    double dev_border_vel = penaltyBoundToIntervalDerivative(vel, -cfg_->robot.max_vel_x_backwards, cfg_->robot.max_vel_x,cfg_->optim.penalty_epsilon);
    double dev_border_omega = penaltyBoundToIntervalDerivative(omega, cfg_->robot.max_vel_theta,cfg_->optim.penalty_epsilon);
    
    _jacobianOplus[0].resize(2,3); // conf1
    _jacobianOplus[1].resize(2,3); // conf2
    _jacobianOplus[2].resize(2,1); // deltaT
    
//  if (aux1==0) aux1=1e-6;
//  if (aux2==0) aux2=1e-6;
    
    if (dev_border_vel!=0)
    {
      double aux3 = dev_border_vel / aux1;
      _jacobianOplus[0](0,0) = -deltaS[0] * aux3; // vel x1
      _jacobianOplus[0](0,1) = -deltaS[1] * aux3; // vel y1	
      _jacobianOplus[1](0,0) = deltaS[0] * aux3; // vel x2
      _jacobianOplus[1](0,1) = deltaS[1] * aux3; // vel y2
      _jacobianOplus[2](0,0) = -vel * aux2 * dev_border_vel; // vel deltaT
    }
    else
    {
      _jacobianOplus[0](0,0) = 0; // vel x1
      _jacobianOplus[0](0,1) = 0; // vel y1	
      _jacobianOplus[1](0,0) = 0; // vel x2
      _jacobianOplus[1](0,1) = 0; // vel y2	
      _jacobianOplus[2](0,0) = 0; // vel deltaT
    }
    
    if (dev_border_omega!=0)
    {
      double aux4 = aux2 * dev_border_omega;
      _jacobianOplus[2](1,0) = -omega * aux4; // omega deltaT
      _jacobianOplus[0](1,2) = -aux4; // omega angle1
      _jacobianOplus[1](1,2) = aux4; // omega angle2
    }
    else
    {
      _jacobianOplus[2](1,0) = 0; // omega deltaT
      _jacobianOplus[0](1,2) = 0; // omega angle1
      _jacobianOplus[1](1,2) = 0; // omega angle2			
    }

    _jacobianOplus[0](1,0) = 0; // omega x1
    _jacobianOplus[0](1,1) = 0; // omega y1
    _jacobianOplus[1](1,0) = 0; // omega x2
    _jacobianOplus[1](1,1) = 0; // omega y2
    _jacobianOplus[0](0,2) = 0; // vel angle1
    _jacobianOplus[1](0,2) = 0; // vel angle2
  }
#endif
#endif
 
  
public:
  
  EIGEN_MAKE_ALIGNED_OPERATOR_NEW

};

上面的代码中_error[0]对应线速度vel误差,_error[1]对应角速度omega误差。所以需要求解vel和omega分别相对于顶点变量x,y,theta的导数。

_jacobianOplus[0].resize(2,3);对应第一个顶点,其中误差项有两维,顶点优化变量有3维,所以雅克比矩阵是一个2x3的矩阵。

_jacobianOplus[1].resize(2,3);对应第二个顶点,其中误差项有两维,顶点优化变量有3维,所以雅克比矩阵是一个2x3的矩阵。

_jacobianOplus[2].resize(2,1);对应第三个顶点,其中误差项有两维,顶点是时间差,只有一维,所以雅克比矩阵是一个2x1的矩阵。

创建优化器

  // 构建图优化,先设定g2o
  typedef g2o::BlockSolver<g2o::BlockSolverTraits<1, 1>> BlockSolverType;  // 每个误差项优化变量维度为1,误差值维度为1
  typedef g2o::LinearSolverDense<BlockSolverType::PoseMatrixType> LinearSolverType; // 线性求解器类型

  // 梯度下降方法,可以从GN, LM, DogLeg 中选
  auto solver = new g2o::OptimizationAlgorithmGaussNewton(
    g2o::make_unique<BlockSolverType>(g2o::make_unique<LinearSolverType>()));
  g2o::SparseOptimizer optimizer;     // 图模型
  optimizer.setAlgorithm(solver);   // 设置求解器
  optimizer.setVerbose(true);       // 打开调试输出

teb中的优化器设置

//! Typedef for the block solver utilized for optimization
typedef g2o::BlockSolverX TEBBlockSolver;

//! Typedef for the linear solver utilized for optimization
typedef g2o::LinearSolverCSparse<TEBBlockSolver::PoseMatrixType> TEBLinearSolver;
//typedef g2o::LinearSolverCholmod<TEBBlockSolver::PoseMatrixType> TEBLinearSolver;

/*
 * initialize g2o optimizer. Set solver settings here.
 * Return: pointer to new SparseOptimizer Object.
 */
std::shared_ptr<g2o::SparseOptimizer> TebOptimalPlanner::initOptimizer()
{
  // Call register_g2o_types once, even for multiple TebOptimalPlanner instances (thread-safe)
  static std::once_flag flag;
  std::call_once(flag, [this](){this->registerG2OTypes();});

  // allocating the optimizer
  std::shared_ptr<g2o::SparseOptimizer> optimizer = std::make_shared<g2o::SparseOptimizer>();
  auto linearSolver = std::make_unique<TEBLinearSolver>(); // see typedef in optimization.h
  linearSolver->setBlockOrdering(true);
  auto blockSolver = std::make_unique<TEBBlockSolver>(std::move(linearSolver));
  g2o::OptimizationAlgorithmLevenberg* solver = new g2o::OptimizationAlgorithmLevenberg(std::move(blockSolver));

  optimizer->setAlgorithm(solver);
  
  optimizer->initMultiThreading(); // required for >Eigen 3.1
  
  return optimizer;
}

从上面的代码对比可以看出,BlockSolver可设成固定的,也可以设成动态的。

下面是本人对于BlockSolver设置的理解,不一定正确。仅为参考。

优化变量的维度通常是明确的,但误差项的维度可能是变化的。每条边的误差项维度可能都不一样。这时应该使用g2o::BlockSolverX,以便能动态适应误差项的维度。

linear solver也是可选的。主要有下面几种可选:

g2o::LinearSolverCSparse
g2o::LinearSolverCholmod
g2o::LinearSolverEigen
g2o::LinearSolverDense
g2o::LinearSolverPCG

不同的linear solver有不同的依赖,需要注意是否已经安装了相应的依赖。

添加顶点

  // 往图中增加顶点
  std::vector<SimpleVertex *> vertexs;

  SimpleVertex *v = new SimpleVertex();
  v->setEstimate(0);
  v->setId(0);
  v->setFixed(true);
  optimizer.addVertex(v);
  vertexs.push_back(v);

  SimpleVertex *v1 = new SimpleVertex();
  v1->setEstimate(1);
  v1->setId(1);
  optimizer.addVertex(v1);
  vertexs.push_back(v1);

  SimpleVertex *v2 = new SimpleVertex();
  v2->setEstimate(0.1);
  v2->setId(2);
  optimizer.addVertex(v2);
  vertexs.push_back(v2);

顶点的id需要确保不能重复。

顶点可使用setEstimate接口设置初始值。

添加边

  // 往图中增加边
  SimpleUnaryEdge *edge = new SimpleUnaryEdge();
  // edge->setId(i);
  edge->setVertex(0, vertexs[0]);                // 设置连接的顶点
  edge->setMeasurement(0);      // 观测数值
  edge->setInformation(Eigen::Matrix<double, 1, 1>::Identity()); // 信息矩阵:协方差矩阵之逆
  optimizer.addEdge(edge);

  SimpleBinaryEdge * edge1 = new SimpleBinaryEdge();
  // edge->setId(i);//id 不设置似乎没有关系,如果设置需要每条边设置成不一样的
  edge1->setVertex(0, vertexs[1]);//这里设置的序号对应的顶点要和边的computeError函数里设定的顶点是一一对应的
  edge1->setVertex(1, vertexs[0]);                // 设置连接的顶点
  edge1->setMeasurement(1);      // 观测数值
  edge1->setInformation(Eigen::Matrix<double, 1, 1>::Identity()*1.0); // 信息矩阵:协方差矩阵之逆
  optimizer.addEdge(edge1);

  SimpleBinaryEdge * edge2 = new SimpleBinaryEdge();
  // edge->setId(i);
  edge2->setVertex(0, vertexs[2]);                // 设置连接的顶点
  edge2->setVertex(1, vertexs[1]);                // 设置连接的顶点
  edge2->setMeasurement(-0.8);      // 观测数值
  edge2->setInformation(Eigen::Matrix<double, 1, 1>::Identity()); // 信息矩阵:协方差矩阵之逆
  optimizer.addEdge(edge2);

  SimpleBinaryEdge * edge3 = new SimpleBinaryEdge();
  // edge->setId(i);
  edge3->setVertex(0, vertexs[2]);                // 设置连接的顶点
  edge3->setVertex(1, vertexs[0]);                // 设置连接的顶点
  edge3->setMeasurement(0);      // 观测数值
  edge3->setInformation(Eigen::Matrix<double, 1, 1>::Identity()); // 信息矩阵:协方差矩阵之逆
  optimizer.addEdge(edge3);

边也有setId接口,但不设置好像也可以。如果设置也需确保不重复。但id号的顺序似乎并没有要求。

使用setVertex接口设置顶点时是有顺序的。这个顺序与边的computeError函数中使用顶点的顺序要对应起来。

setMeasurement接口用于设置内部的_measurement值。

setInformation接口用于设置信息矩阵。信息矩阵是方正矩阵,其行列数由误差项的维度决定。

void TebOptimalPlanner::AddEdgesVelocity()
{
  if (cfg_->robot.max_vel_y == 0) // non-holonomic robot
  {
    if ( cfg_->optim.weight_max_vel_x==0 && cfg_->optim.weight_max_vel_theta==0)
      return; // if weight equals zero skip adding edges!

    int n = teb_.sizePoses();
    Eigen::Matrix<double,2,2> information;
    information(0,0) = cfg_->optim.weight_max_vel_x;
    information(1,1) = cfg_->optim.weight_max_vel_theta;
    information(0,1) = 0.0;
    information(1,0) = 0.0;

    for (int i=0; i < n - 1; ++i)
    {
      EdgeVelocity* velocity_edge = new EdgeVelocity;
      velocity_edge->setVertex(0,teb_.PoseVertex(i));
      velocity_edge->setVertex(1,teb_.PoseVertex(i+1));
      velocity_edge->setVertex(2,teb_.TimeDiffVertex(i));
      velocity_edge->setInformation(information);
      velocity_edge->setTebConfig(*cfg_);
      optimizer_->addEdge(velocity_edge);
    }
  }
  ...
}

最终优化求解的结果为estimated model: 0 0.933333 0.0666667 。基本和标准结果一致。

如果将edge1的信息矩阵如下设置:

edge1->setInformation(Eigen::Matrix<double, 1, 1>::Identity()*10.0);

优化求解的结果为:estimated model: 0 0.990476 0.0952381

可以发现x1更接近1了。说明此时我们更相信编码器测量的从x0到x1的距离值。

完整的测试代码可查看下面的链接:

https://github.com/shoufei403/g2o_learning.git


觉得有用就点赞吧!

我是首飞,一个帮大家填坑的机器人开发攻城狮。

另外在公众号《首飞》内回复“机器人”获取精心推荐的C/C++,Python,Docker,Qt,ROS1/2等机器人行业常用技术资料。

G2O (General Graph Optimization)入门及简单使用 g2o全称是General Graph Optimization,也就是图优化,我们在做SLAM后端或者更加常见的任何优化问题(曲线拟合)都可以使用G2O进行处理。 阅读详情

相关推荐

G2O非线性优化

G2O非线性优化 上一篇文章介绍最小二乘法,本文介绍g2o实现最小二乘法。 g2o即General Graphic Optimization,它是一个基于图优化的库。至于图优化理论,参照半闲居士的博客 接着上回曲线拟合问题,拟合y=Asin(Bx)+Ccos(Dx)y=Asin(Bx)+Ccos(Dx),已知N组数据(xi,yi),i=0,1,⋯N−1(x_i, y_i), i = 0, 1,

STI10hy的博客 8819

TEB算法详解(TebLocalPlannerROS::computeVelocityCommands(3))

第一章主要分析了teb算法的准备条件,包括获取当前位姿与速度、对全局路径的裁减以获取局部路径等以及局部地图的获取等。第二章主要讲述了在获取前置条件后,如何根据前置条件进行位姿优化,teb的路径优化主要是调用了g2o优化算法,以全局路径的点以及理论上点与点之间的运动时间作为g2o的优化顶点,以距离障碍物的距离、机器人最大运动速度、点与点之间的最短路径等约束条件作为g2o的边对局部路径规划的姿态进行了优化。这一章则是主要看一下算法怎么从优化后的点位到运动输出之间的关系。

YiYeZhiNian的博客 3649

g2o图优化简介与基本使用方法

g2o简介 g2o(General Graphic Optimization)是一个基于图优化的库,将非线性优化与图论结合起来的理论,我们可以利用g2o求解任何可以表示为图优化的最小二乘问题。 图优化就是把优化问题表现成图的方式。图由顶点和边组成,其中顶点表示优化变量,边表示误差项,对任意一个非线性?> 最小二乘问题,我们都可以构建与之对应的图。 (注:这里的图是图论意义上的图,可以用概率论里面的定义,贝叶斯图或因子图。) g2o安装 首先安装g2o的依赖 sudo apt install qt5

个人博客 4247

TEB算法3- teb中的g2o问题

针对阿克曼底盘的导航,TEB_local_planner没有直接给出具体的导航策略,而是利用非线性优化问题完成了对该问题的建模。

ゃ﹏淸凨ゞ 3186

视觉SLAM非线性优化之Ceres solver 和 G2o

“ 如果你受苦了,感谢生活,那是它给你的一份感觉; 如果你受苦了,感谢上帝,说明你还活着。” ——伊索 3月就要结束了~ 最近在看视觉SLAM的一些东西...

I-am-Unique 4009

g2o简单测试

开始研究g2o,记录一些过程。 首先使用create_sphere可以生成一个数据文件。然后使用g2o_viewer就可以载入这个文件。 载入以后可以进行优化,优化后是一个球体。 这里g2o_viewer不是万能的可视化窗口,有许多已存在在g2o工程里的数据类型是不知显示的,比如VertexSE3Expmap因此如果希望使用g2o_viewer*来方便调试,有两个办法,第一个办法就是老老实实的

u010566411的专栏 4434

g2o使用总结

需要解决的问题: 如何自定义Vertex和Edge? 如何选择Edge类型?一元还是二元? 如何赋值信息矩阵? 如何设置鲁棒核函数的阈值? 如何选择Vertex设置为固定? 如何边缘化以便稀疏化求解? 如何处理优化结束后outliner? 如何设置根据卡方分布的临界值表对Edge的chi2设置阈值? 1. g2o提供的顶点vertex 1) 李代数位姿 clas...

陋室逢雨 1万+

图优化g2o教程简介

图优化g2o教程简介,包括g2o间接、如何使用,顶点和边的定义,信息矩阵、鲁棒核函数对于优化的结果的影响。

u010196944的博客 2129

【SLAM】G2O优化库超详细解析

通用图优化。为何叫通用呢?G2O的核里带有各种各样的求解器,而它的顶点、边的类型则多种多样。通过自定义顶点和边,事实上,只要一个优化问题能够表达成图,那么就可以用G2O去求解它。常见的,比如bundle adjustment、ICP、数据拟合,都可以用G2O来做。所谓图优化,是把优化问题表现成图的一种方式,这里的图是图论意义上的图。一个图由若干个顶点,以及连着这些顶点的边组成。在这里,我们用顶点表示优化变量,而用边表示残差项。

Yngz_Miao的博客 1万+

理解图优化,一步步带你看懂g2o框架.pdf

理解图优化,一步步带你看懂g2o框架.pdf

如何在win11中用双硬盘或移动硬盘装Ubuntu 20.04 双系统

win11 安装双系统

qq_42817360的博客 5138

G2O优化

G2O优化0.引言0.1.参考0.1.G2O的解决问题0.2.G2O整体框架1. 基本使用1.1.构造 g2o 模型1.2.g2o 类图2.g2o 的顶点(Vertex)2.1.顶点的格式2.2.自定义顶点2.3.将顶点添加到图中3.g2o 中的边3.1.边的格式3.2.自定义边3.3.将边添加到图中4.ORB-SLAM2 中的顶点4.1.相机位姿顶点4.2.地图点坐标顶点4.3.闭环时的 sim3 相机位姿5. ORB-SLAM2 中的边5.1.位姿优化时的重投影误差5.2.BA 优化时的重投影误差5.3

fb_941219的博客 5088

(一)Tushare Pro教程:股票基本数据接口

 一直使用tushare作数据源,这里开个教程大家一起学习。 一、介绍        tushare是一个python股票数据源,现已升级到tushare pro版本,新版本在数据稳定性和获取速率上大幅提升。与wind、东财等收费数据源相比,tushare pro算是金融数据界的良心之作,其为从事金融工程、量化研究等方向的同学们提供了便利。tushare pro还有一点突出的特色就是已经将股票...

虎毛的博客 4万+

G2o,GTSAM,Ceres,Tensorflow优化器的方法比较

首先是一些背景知识: 最速梯度下降法(GD)和牛顿法可以用来优化所有种类的函数。 牛顿高斯和基于牛顿高斯的LM法只能用来优化非线性最小二乘。 SGD是最速梯度下降法的改进,也就是每次迭代随即选取一些样本进行优化。 一般大型优化问题不会用二阶的牛顿法。 Ceres是针对所有种类的函数的优化,所以只能用最速梯度下降法。 G2o是求解的一个标准的最小二乘 GTSAM准确的...

ziliwangmoe的博客 9137

ABtest描述统计功效

我们在面试数据分析时,常常会问到如何确定样本量的问题,会有一个指标:统计功效(power),或者是1-β 统计功效会影响实验的可信度,甚至会影响实验所需的样本量,要了解统计功效,首先要理解两类统计错误 假如我们做了一个AB实验,且两组用户的数据均值为μ1和μ2,我们会有两个初始假设: H0:两个实验组的数据没有差别(μ1=μ2),这个假设也叫“零假设”,也叫原假设 H1:两个实验组之间的数据存在差别(μ1≠μ2),这个假设也叫“非零假设”,也叫备择假设 定义两类错误: 实验本身没有效应,但我们觉

qq_38134556的博客 898

关于RTP时间戳以及播放器对时间戳的处理

首先,了解时间戳几个基本概念: 时间戳单位:时间戳计算的单位不是秒之类的单位,而是由采样频率所代替的单位,这样做的目的就是为了是时间戳单位更为精准。比如说一个音频的采样频率为8000Hz,那么我们可以把时间戳单位设为1 / 8000。 时间戳增量:相邻两帧之间的时间差(以时间戳单位为基准)。 采样频率: 有些地方也叫时钟频率,即 每秒钟抽取样本的次数,例如音频的采...

toshiba689的专栏 3052

一步步带你走入基于图优化的g2o库:General Graphic Optimization

图优化是通过使用图结构来表示和解决优化问题的一种方法。图是由顶点(Vertex)和边(Edge)组成。在SLAM中,搭载传感器的主体的位姿构成顶点,位姿之间的变换关系构成边。一旦图构建完成了,就要主体的位姿去尽量满足这些边构成的约束。对于一个非线性的最小二乘问题,我们可以用顶点来表示优化变量,用边来表示误差项,这样就能构建一个对应的图,我们称这个图为贝叶斯图或者是因子图。可以把SLAM中的图优化问题分解为:构建图、优化图。位于hyper_graph.h中。

2301_76831056的博客 3194

c++通过自然语言处理技术分析语音信号音高

其中,缓冲区大小为2048,最小基频为80Hz,最大基频为1000Hz。通常情况下,音高在85Hz到255Hz之间,因此可以根据实际情况修改最小和最大基频值。C++中可以使用数字信号处理库或语音处理库,比如MATLAB、Praat、YIN算法等等来进行音高分析。对于语音信号的音高分析,可以使用基频提取技术。基频是指一个声音周期的重复率,也就是一个声音波形中最长的周期。通常情况下,人的声音基频范围是85Hz到255Hz。这段代码使用了YIN算法来计算语音信号的基频,得出的结果即为音高。

SYC20110120的博客 1728

G2O使用说明(SLAM中)

G2O使用说明(SLAM) 范例:在使用G2O要明确目标函数是什么,优化变量是几维,误差项是几维。举例说明:现在有很多个点符合y=exp(ax2+bx+c)+wy=exp(ax^{2}+bx+c)+wy=exp(ax2+bx+c)+w,其中w为高斯噪声。那么目标函数就是e(a,b,c)=y−exp(ax2+bx+c)e(a,b,c)=y-exp(ax^{2}+bx+c)e(a,b,c)=y−exp(ax2+bx+c)使得它最小(最小二乘法的原理)。优化变量就是a、b、c是三维的,误差项是e(a,b,c)的值

qq_45369294的博客 2004

基于Simulink使用激光雷达数据跟踪车辆仿真(附源码)

这些模型将长方体的尺寸定义为状态转换期间的常量,并且它们的估计值在滤波器的校正阶段随时间演变。下面的动画显示了在滤波器的预测阶段将具有不同概率的恒定速度和恒定转弯率模型混合在一起的效果。此示例显示了如何使用带有 IMM 过滤器的 JPDA 跟踪器来跟踪使用激光雷达传感器的对象。图像和激光雷达数据读取器从 MAT 文件中读取记录的数据,并分别输出参考图像和点云中点的位置。跟踪器块选中“启用所有轨道输出”和“启用可检测的跟踪 ID 输入”复选框,以从跟踪器输出所有轨道,并将其检测概率作为其状态的函数进行计算。

珞瑜的博客 2566

STM32开发 -- LIS3DH三轴加速度计开发详解

这次要讲一下LIS3DH三轴加速度计的开发使用,之前是有讲过 MPU6050开发 的,但当时是在Linux下做的。额,还是个半成品,哈哈。现在手边有STM32的开发板,可以实实在在的来实现以下卡尔曼滤波了。 一、LIS3DH初识 依照惯例,先看一下手册,大致先了解一下,LIS3DH是什么模块。 下载:LIS3DH datasheet 1、描述 LIS3DH是一款超低功耗高性能三轴线性加速度计,属于...

不积跬步,无以至千里 2万+

g2o的简介、安装、使用说明和曲线拟合例程

g2o(General Graphic Optimization)是基于图优化的库。该工具尤其是在图像优化方面被广泛的应用。本文主要围绕g2o进行展开介绍相关的安装、使用、g2o的结构以及曲线拟合的实例。

非晚非晚的博客 1万+

C-08-实战:搭建一个AI客服机器人(接入微信企微)

客服是一人公司最大的时间杀手。本文教你搭建一个AI客服机器人,能自动回答80%的常见问题,只把复杂问题转给人工。含微信公众号和企业微信接入完整代码。

Johnny2004的专栏 88

Ceres和g2o的配置和使用

上文非线性优化介绍了非线性优化的基本求解方法,并使用C++手动实现了曲线拟合实例。本文介绍ceres和g2o库的配置方法,并通过曲线拟合实例介绍其使用方法。 Ceres安装 Google Ceres 是一个广泛使用的最小二乘问题求解库。在 Ceres 中,我们作为用户,只需按照一定步骤定义待解的优化问题,然后交给求解器计算即可。 Google Ceres库的源码:https://github.com/ceres-solver/ceres-solver 依赖项: ...

陈建驱的博客 2319

最新的sdio 4.0 spec 协议标准 带书签 高清完整PDF版

Part-E1-SDIO-Specification-Ver4.00-Final-120220.pdf:最新的sdio 4.0 spec 协议标准 带书签 高清完整PDF版,包括SDIO 接口的外设,SDIO 接口官方协议等。

迈瑞BS系列全自动生化分析仪LIS协议接口手册(内含HL7和ASTM)

迈瑞BS系列全自动生化分析仪LIS协议接口手册,内含HL7和ASTM)

SLAM——之G2O函数库

1. G2O示例 相较于Ceres而言,G2O函数库相对较为复杂,但是适用面更加广,可以解决较为复杂的重定位问题。Ceres库向通用的最小二乘问题的求解,定义优化问题,设置一些选项,可通过Ceres求解。而图优化,是把优化问题表现成图的一种方式,这里的图是图论意义上的图。一个图由若干个顶点,以及连着这些顶点的边组成。在这里,我们用顶点表示优化变量,而用边表示误差项。 为了使用g2o,首先要将曲线拟合问题抽象成图优化。这个过程中,只要记住节点为优化变量,边为误差项即可。曲线拟合的图优化问题可以画成以下形式:

格物致知 6614

visio安装包

Microsoft Office Visio 是一款专业的制图软件。利用经过更新的形状、协作工具和数据链接图表,创建专业级图表以简化复杂的信息。

【视觉惯性SLAM:六、图优化库(1):g2o的使用指南】

g2o 提供了一些基本顶点类,但也支持用户自定义顶点。自定义顶点需要继承 g2o::BaseVertex 并实现相应的虚函数。示例:自定义一个三维点顶点// 定义一个3维点的顶点public:// 初始化函数,设置顶点初始值// 初始化为 (0, 0, 0)// 更新顶点值// 更新顶点估计值// 可选:保存和加载顶点信息(序列化)

计算机领域博客,与时俱进的思想与技术。 1376

win10 跳过U8环境检测的IE WebControls.rar

win10 安装u8的时候 提示IE WebControls 不通过,这个工具可以欺骗u8的检查。 跳过U8环境检测的IE WebControls。提示导入失败后,重新检测通过。

上一篇: LCM库的简单使用
下一篇: 对pure pursuit算法的理解和改进
首飞爱玩机器人
博客等级 码龄14年 6989粉丝 120原创
评论 3
成就一亿技术人!
拼手气红包6.0元
还能输入1000个字符
 
 条评论被折叠 查看
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付:10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

首飞爱玩机器人

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

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

余额充值