SLAM系统中,一般的做法是通过P3P或者EPnP等方法进行投影估算出相机位姿,进而构建最小二乘对估计的相机位姿进行优化,从而得到精度比较高的位姿。优化完后可以剔除一些outlier点,这样的话后边再进行投影和优化的时候精度会有所保证。ORB-SLAM系统当中的优化分为下面5个,从优化的范围层次逐步扩大。,

1.位姿优化函数PoseOptimization

一些说明:PoseOptimization优化的是单帧图像的位姿。相机在运动的时候每个图像帧可以对应三维世界中多个地图点,那么此时我们可以知道三维世界中的多个3D地图点,也可以知道图像帧中对应的2D像素点坐标,就可以通过3D-2D的投影关系估计出相机的位姿。PoseOptimization就是对单个图像帧中的多个地图点建立多个一元连接边,构成图进行优化, 优化单个图像帧的SE3位姿。

调用时机:当一个图像帧中的特征点和3D空间中多个地图点匹配,并且有一个“初始位姿”的时候,就可以通过位姿优化函数来对当前帧的位姿进行优化。ORB-SLAM系统当中在Tracking线程中多次调用了位姿优化函数,可以试想一下,毕竟Tracking线程的主要业务就是进行跟踪,一帧一帧图像的位姿描述了相机的轨迹。系统的前两帧图像的位姿都初始化为单位矩阵,在通过2D点匹配建立关联关系后,通过三角化创建了地图点,此时才可以对图像帧的位姿进行优化更新。只要3D-2D匹配关系发生变化,就需要对位姿进行优化。

ORB-SLAM系统中下面几个函数里边都调用了位姿优化函数:

1)Tracking::TrackReferenceKeyFrame()

2)Tracking::TrackWithMotionModel()

3)Tracking::TrackLocalMap()

4)Tracking::Relocalization()

代码:

/*** 这里仅优化帧的位姿,不优化地图点的坐标。(之后会剔除优化后的外点,这个操作不再本函数中做)* (3D-2D)将地图点投影到相机坐标系下的相机平面
*/
int Optimizer::PoseOptimization(Frame *pFrame)
{g2o::SparseOptimizer optimizer;g2o::BlockSolver_6_3::LinearSolverType * linearSolver;linearSolver = new g2o::LinearSolverDense<g2o::BlockSolver_6_3::PoseMatrixType>();g2o::BlockSolver_6_3 * solver_ptr = new g2o::BlockSolver_6_3(linearSolver);g2o::OptimizationAlgorithmLevenberg* solver = new g2o::OptimizationAlgorithmLevenberg(solver_ptr);optimizer.setAlgorithm(solver);int nInitialCorrespondences=0;// Set Frame vertex 设置帧结点//VertexSE3Expmap为李代数位姿g2o::VertexSE3Expmap * vSE3 = new g2o::VertexSE3Expmap();//当前帧的变换矩阵vSE3->setEstimate(Converter::toSE3Quat(pFrame->mTcw));vSE3->setId(0);vSE3->setFixed(false);optimizer.addVertex(vSE3);// Set MapPoint vertices 设置地图点结点const int N = pFrame->N;vector<g2o::EdgeSE3ProjectXYZOnlyPose*> vpEdgesMono;vector<size_t> vnIndexEdgeMono;vpEdgesMono.reserve(N);vnIndexEdgeMono.reserve(N);vector<g2o::EdgeStereoSE3ProjectXYZOnlyPose*> vpEdgesStereo;vector<size_t> vnIndexEdgeStereo;vpEdgesStereo.reserve(N);vnIndexEdgeStereo.reserve(N);const float deltaMono = sqrt(5.991);const float deltaStereo = sqrt(7.815);{unique_lock<mutex> lock(MapPoint::mGlobalMutex);//遍历帧中的MapPointfor(int i=0; i<N; i++){MapPoint* pMP = pFrame->mvpMapPoints[i];if(pMP){// Monocular observation 单目if(pFrame->mvuRight[i]<0){nInitialCorrespondences++;pFrame->mvbOutlier[i] = false;Eigen::Matrix<double,2,1> obs;//获取关键点const cv::KeyPoint &kpUn = pFrame->mvKeysUn[i];//obs中输入关键点的x,y坐标obs << kpUn.pt.x, kpUn.pt.y;/*** EdgeSE3ProjectXYZOnlyPose():PoseEstimation中的重投影误差,将地图点投影到相机坐标系下的相机平面。*///构造pose一元边g2o::EdgeSE3ProjectXYZOnlyPose* e = new g2o::EdgeSE3ProjectXYZOnlyPose();//这里只设置一个顶点//optimizer.vertex(0)返回的就是上边的vSE3结点,也就是当前帧的Tcwe->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(0)));//设置观测数值,obs中为帧中的关键点坐标,也就是三维坐标中物体投影到相机平面的坐标e->setMeasurement(obs);//获取关键点所在层的sigma2值const float invSigma2 = pFrame->mvInvLevelSigma2[kpUn.octave];e->setInformation(Eigen::Matrix2d::Identity()*invSigma2);g2o::RobustKernelHuber* rk = new g2o::RobustKernelHuber;e->setRobustKernel(rk);rk->setDelta(deltaMono);//设置阈值//相机内参e->fx = pFrame->fx;e->fy = pFrame->fy;e->cx = pFrame->cx;e->cy = pFrame->cy;//地图中MapPoint的世界坐标系的坐标值cv::Mat Xw = pMP->GetWorldPos();e->Xw[0] = Xw.at<float>(0);e->Xw[1] = Xw.at<float>(1);e->Xw[2] = Xw.at<float>(2);optimizer.addEdge(e);vpEdgesMono.push_back(e);vnIndexEdgeMono.push_back(i);}else  // Stereo observation 双目{nInitialCorrespondences++;pFrame->mvbOutlier[i] = false;//SET EDGE 设置边Eigen::Matrix<double,3,1> obs;const cv::KeyPoint &kpUn = pFrame->mvKeysUn[i];const float &kp_ur = pFrame->mvuRight[i];//obs中存入的是双目相机中像素坐标下一个关键点的坐标obs << kpUn.pt.x, kpUn.pt.y, kp_ur;//优化变量只有pose,地图点位置固定,是一边元,双目中使用的是EdgeStereoSE3ProjectXYZOnlyPoze()。g2o::EdgeStereoSE3ProjectXYZOnlyPose* e = new g2o::EdgeStereoSE3ProjectXYZOnlyPose();e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(0)));e->setMeasurement(obs);const float invSigma2 = pFrame->mvInvLevelSigma2[kpUn.octave];Eigen::Matrix3d Info = Eigen::Matrix3d::Identity()*invSigma2;e->setInformation(Info);g2o::RobustKernelHuber* rk = new g2o::RobustKernelHuber;e->setRobustKernel(rk);rk->setDelta(deltaStereo);e->fx = pFrame->fx;e->fy = pFrame->fy;e->cx = pFrame->cx;e->cy = pFrame->cy;e->bf = pFrame->mbf;cv::Mat Xw = pMP->GetWorldPos();e->Xw[0] = Xw.at<float>(0);e->Xw[1] = Xw.at<float>(1);e->Xw[2] = Xw.at<float>(2);optimizer.addEdge(e);vpEdgesStereo.push_back(e);vnIndexEdgeStereo.push_back(i);}}}}//至少需要使用三对点进行P3P求解,三对3D-2D匹配点if(nInitialCorrespondences<3)return 0;// We perform 4 optimizations, after each optimization we classify observation as inlier/outlier// At the next optimization, outliers are not included, but at the end they can be classified as inliers again.const float chi2Mono[4]={5.991,5.991,5.991,5.991};const float chi2Stereo[4]={7.815,7.815,7.815, 7.815};const int its[4]={10,10,10,10};    int nBad=0;//优化四次for(size_t it=0; it<4; it++){vSE3->setEstimate(Converter::toSE3Quat(pFrame->mTcw));optimizer.initializeOptimization(0);optimizer.optimize(its[it]);nBad=0;//单目for(size_t i=0, iend=vpEdgesMono.size(); i<iend; i++){g2o::EdgeSE3ProjectXYZOnlyPose* e = vpEdgesMono[i];const size_t idx = vnIndexEdgeMono[i];if(pFrame->mvbOutlier[idx]){e->computeError();}const float chi2 = e->chi2();if(chi2>chi2Mono[it]){                pFrame->mvbOutlier[idx]=true;e->setLevel(1);nBad++;}else{pFrame->mvbOutlier[idx]=false;e->setLevel(0);}if(it==2)e->setRobustKernel(0);}//双目for(size_t i=0, iend=vpEdgesStereo.size(); i<iend; i++){g2o::EdgeStereoSE3ProjectXYZOnlyPose* e = vpEdgesStereo[i];const size_t idx = vnIndexEdgeStereo[i];if(pFrame->mvbOutlier[idx]){e->computeError();}const float chi2 = e->chi2();if(chi2>chi2Stereo[it]){pFrame->mvbOutlier[idx]=true;e->setLevel(1);nBad++;}else{                e->setLevel(0);pFrame->mvbOutlier[idx]=false;}if(it==2)e->setRobustKernel(0);}if(optimizer.edges().size()<10)break;}    // Recover optimized pose and return number of inliersg2o::VertexSE3Expmap* vSE3_recov = static_cast<g2o::VertexSE3Expmap*>(optimizer.vertex(0));g2o::SE3Quat SE3quat_recov = vSE3_recov->estimate();cv::Mat pose = Converter::toCvMat(SE3quat_recov);//设置优化完成后的位姿pFrame->SetPose(pose);return nInitialCorrespondences-nBad;
}

2.优化两帧之间的位姿OptimizeSim3

一些说明:对单目相机来说具有尺度不确定性,当相机运动一段时间后不可避免会产生尺度漂移累积误差增大,当进行闭环检测的时候如果尺度不一致就不能进行很好的对应。相比于SE3来说,相似变换可以提供尺度的表达。那么这时候就可以借助相似变换群sim3来求解当前关键帧和闭环候选帧之间的sim3位姿,同时sim3位姿也可以和SE3进行转换。

要优化的两帧都能观测到多个相同的地图点,并且已经知道有哪些地图点是一致的。此时可以构造一个超定方程组,并对其求解最小化误差优化。优化帧间变化的SIM3位姿与地图的VertexSBAPointXYZ位姿。

调用时机:单目闭环检测中检测出闭环后,对当前关键帧和闭环帧之间计算sim3位姿,当计算出sim3位姿后,就需要调用sim3位姿优化函数。

ORB-SLAM当中在LoopClosing::ComputetSim3()中调用了OptimizeSim3函数。

代码:

int Optimizer::OptimizeSim3(KeyFrame *pKF1, KeyFrame *pKF2, vector<MapPoint *> &vpMatches1, g2o::Sim3 &g2oS12, const float th2, const bool bFixScale)
{g2o::SparseOptimizer optimizer;g2o::BlockSolverX::LinearSolverType * linearSolver;linearSolver = new g2o::LinearSolverDense<g2o::BlockSolverX::PoseMatrixType>();g2o::BlockSolverX * solver_ptr = new g2o::BlockSolverX(linearSolver);g2o::OptimizationAlgorithmLevenberg* solver = new g2o::OptimizationAlgorithmLevenberg(solver_ptr);optimizer.setAlgorithm(solver);// Calibrationconst cv::Mat &K1 = pKF1->mK;const cv::Mat &K2 = pKF2->mK;// Camera posesconst cv::Mat R1w = pKF1->GetRotation();const cv::Mat t1w = pKF1->GetTranslation();const cv::Mat R2w = pKF2->GetRotation();const cv::Mat t2w = pKF2->GetTranslation();// Set Sim3 vertex 设置sim3顶点,将pKF1到pKF2之间的sim3位姿变换加入到图中作为结点0//g2o::VertexSim3Expmap用于创建sim3位姿结点g2o::VertexSim3Expmap * vSim3 = new g2o::VertexSim3Expmap();    vSim3->_fix_scale=bFixScale;vSim3->setEstimate(g2oS12);vSim3->setId(0);vSim3->setFixed(false);//获取相机内参vSim3->_principle_point1[0] = K1.at<float>(0,2);vSim3->_principle_point1[1] = K1.at<float>(1,2);vSim3->_focal_length1[0] = K1.at<float>(0,0);vSim3->_focal_length1[1] = K1.at<float>(1,1);vSim3->_principle_point2[0] = K2.at<float>(0,2);vSim3->_principle_point2[1] = K2.at<float>(1,2);vSim3->_focal_length2[0] = K2.at<float>(0,0);vSim3->_focal_length2[1] = K2.at<float>(1,1);optimizer.addVertex(vSim3);// Set MapPoint verticesconst int N = vpMatches1.size();//pKF1当中所有的地图点const vector<MapPoint*> vpMapPoints1 = pKF1->GetMapPointMatches();//sim3投影边创建vector<g2o::EdgeSim3ProjectXYZ*> vpEdges12;vector<g2o::EdgeInverseSim3ProjectXYZ*> vpEdges21;vector<size_t> vnIndexEdge;vnIndexEdge.reserve(2*N);vpEdges12.reserve(2*N);vpEdges21.reserve(2*N);const float deltaHuber = sqrt(th2);int nCorrespondences = 0;for(int i=0; i<N; i++){if(!vpMatches1[i])continue;//vpMapPoints1中存放的是pKF1当中的地图点MapPoint* pMP1 = vpMapPoints1[i];//vpMatches1当中存放的是pKF1和pKF2两个帧中匹配的地图点MapPoint* pMP2 = vpMatches1[i];const int id1 = 2*i+1;const int id2 = 2*(i+1);const int i2 = pMP2->GetIndexInKeyFrame(pKF2);if(pMP1 && pMP2){if(!pMP1->isBad() && !pMP2->isBad() && i2>=0){//g2o::VertexSBAPointXYZ用于创建特征点空间坐标结点g2o::VertexSBAPointXYZ* vPoint1 = new g2o::VertexSBAPointXYZ();cv::Mat P3D1w = pMP1->GetWorldPos();cv::Mat P3D1c = R1w*P3D1w + t1w;vPoint1->setEstimate(Converter::toVector3d(P3D1c));vPoint1->setId(id1);vPoint1->setFixed(true);optimizer.addVertex(vPoint1);g2o::VertexSBAPointXYZ* vPoint2 = new g2o::VertexSBAPointXYZ();cv::Mat P3D2w = pMP2->GetWorldPos();cv::Mat P3D2c = R2w*P3D2w + t2w;vPoint2->setEstimate(Converter::toVector3d(P3D2c));vPoint2->setId(id2);vPoint2->setFixed(true);optimizer.addVertex(vPoint2);}elsecontinue;}elsecontinue;nCorrespondences++;// Set edge x1 = S12*X2Eigen::Matrix<double,2,1> obs1;const cv::KeyPoint &kpUn1 = pKF1->mvKeysUn[i];obs1 << kpUn1.pt.x, kpUn1.pt.y;//g2o::EdgeSim3ProjectXYZ用于创建重投影误差,优化变量sim3位姿和地图点g2o::EdgeSim3ProjectXYZ* e12 = new g2o::EdgeSim3ProjectXYZ();e12->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(id2)));e12->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(0)));e12->setMeasurement(obs1);const float &invSigmaSquare1 = pKF1->mvInvLevelSigma2[kpUn1.octave];e12->setInformation(Eigen::Matrix2d::Identity()*invSigmaSquare1);g2o::RobustKernelHuber* rk1 = new g2o::RobustKernelHuber;e12->setRobustKernel(rk1);rk1->setDelta(deltaHuber);optimizer.addEdge(e12);// Set edge x2 = S21*X1Eigen::Matrix<double,2,1> obs2;const cv::KeyPoint &kpUn2 = pKF2->mvKeysUn[i2];obs2 << kpUn2.pt.x, kpUn2.pt.y;g2o::EdgeInverseSim3ProjectXYZ* e21 = new g2o::EdgeInverseSim3ProjectXYZ();e21->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(id1)));e21->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(0)));e21->setMeasurement(obs2);float invSigmaSquare2 = pKF2->mvInvLevelSigma2[kpUn2.octave];e21->setInformation(Eigen::Matrix2d::Identity()*invSigmaSquare2);g2o::RobustKernelHuber* rk2 = new g2o::RobustKernelHuber;e21->setRobustKernel(rk2);rk2->setDelta(deltaHuber);optimizer.addEdge(e21);vpEdges12.push_back(e12);vpEdges21.push_back(e21);vnIndexEdge.push_back(i);}// Optimize!optimizer.initializeOptimization();optimizer.optimize(5);// Check inliersint nBad=0;for(size_t i=0; i<vpEdges12.size();i++){g2o::EdgeSim3ProjectXYZ* e12 = vpEdges12[i];g2o::EdgeInverseSim3ProjectXYZ* e21 = vpEdges21[i];if(!e12 || !e21)continue;if(e12->chi2()>th2 || e21->chi2()>th2){size_t idx = vnIndexEdge[i];vpMatches1[idx]=static_cast<MapPoint*>(NULL);optimizer.removeEdge(e12);optimizer.removeEdge(e21);vpEdges12[i]=static_cast<g2o::EdgeSim3ProjectXYZ*>(NULL);vpEdges21[i]=static_cast<g2o::EdgeInverseSim3ProjectXYZ*>(NULL);nBad++;}}int nMoreIterations;if(nBad>0)nMoreIterations=10;elsenMoreIterations=5;if(nCorrespondences-nBad<10)return 0;// Optimize again only with inliersoptimizer.initializeOptimization();optimizer.optimize(nMoreIterations);int nIn = 0;for(size_t i=0; i<vpEdges12.size();i++){g2o::EdgeSim3ProjectXYZ* e12 = vpEdges12[i];g2o::EdgeInverseSim3ProjectXYZ* e21 = vpEdges21[i];if(!e12 || !e21)continue;if(e12->chi2()>th2 || e21->chi2()>th2){size_t idx = vnIndexEdge[i];vpMatches1[idx]=static_cast<MapPoint*>(NULL);}elsenIn++;}// Recover optimized Sim3 恢复优化后的sim3位姿g2o::VertexSim3Expmap* vSim3_recov = static_cast<g2o::VertexSim3Expmap*>(optimizer.vertex(0));g2oS12= vSim3_recov->estimate();return nIn;
}

3.局部BA优化LocalBundleAdjustment

一些说明:局部BA优化,优化的是局部关键帧的位姿和这些局部关键帧可以观测到的地图点的3D坐标。相机在运动过程中,当前时刻的图像帧何其前后的若干个帧是能观测到一些相同的路标点的,也就是说这些帧之间有共视关系。这些共视帧和他们的所有地图点放在一块就体现了“局部”的意思。将局部的所有地图点与可以观测到他们的关键帧放在一起进行优化,将关键帧的SE3位姿和地图点的3D坐标作为待优化量添加到顶点当中。

调用时机:LocalMapping线程当中调用。LocalMapping线程中处理的就是当前关键帧和其共视关键帧、以及这些关键帧能够观测到的地图点之间的连接关系,当然少不了对局部关键帧位姿和地图点3D坐标的更新。

代码:

void Optimizer::LocalBundleAdjustment(KeyFrame *pKF, bool* pbStopFlag, Map* pMap)
{    // Local KeyFrames: First Breath Search from Current Keyframelist<KeyFrame*> lLocalKeyFrames;lLocalKeyFrames.push_back(pKF);pKF->mnBALocalForKF = pKF->mnId;//1.获取局部关键帧列表存放在lLocalKeyFrames中。也就是和当前关键帧有共视关系的关键帧列表。const vector<KeyFrame*> vNeighKFs = pKF->GetVectorCovisibleKeyFrames();for(int i=0, iend=vNeighKFs.size(); i<iend; i++){KeyFrame* pKFi = vNeighKFs[i];pKFi->mnBALocalForKF = pKF->mnId;if(!pKFi->isBad())lLocalKeyFrames.push_back(pKFi);}// Local MapPoints seen in Local KeyFrames//2.创建在局部关键帧中所看到的局部地图点列表,存放在lLocalMapPoints中。list<MapPoint*> lLocalMapPoints;for(list<KeyFrame*>::iterator lit=lLocalKeyFrames.begin() , lend=lLocalKeyFrames.end(); lit!=lend; lit++){vector<MapPoint*> vpMPs = (*lit)->GetMapPointMatches();for(vector<MapPoint*>::iterator vit=vpMPs.begin(), vend=vpMPs.end(); vit!=vend; vit++){MapPoint* pMP = *vit;if(pMP)if(!pMP->isBad())if(pMP->mnBALocalForKF!=pKF->mnId){lLocalMapPoints.push_back(pMP);pMP->mnBALocalForKF=pKF->mnId;}}}// Fixed Keyframes. Keyframes that see Local MapPoints but that are not Local Keyframes//3.创建固定关键帧列表lFixedCameras。这些关键帧可以看到局部地图点,但是这些关键帧并不是局部关键帧。list<KeyFrame*> lFixedCameras;for(list<MapPoint*>::iterator lit=lLocalMapPoints.begin(), lend=lLocalMapPoints.end(); lit!=lend; lit++){map<KeyFrame*,size_t> observations = (*lit)->GetObservations();for(map<KeyFrame*,size_t>::iterator mit=observations.begin(), mend=observations.end(); mit!=mend; mit++){KeyFrame* pKFi = mit->first;if(pKFi->mnBALocalForKF!=pKF->mnId && pKFi->mnBAFixedForKF!=pKF->mnId){//TODO:设置BAFixed标识 pKFi->mnBAFixedForKF=pKF->mnId;if(!pKFi->isBad())lFixedCameras.push_back(pKFi);}}}// Setup optimizerg2o::SparseOptimizer optimizer;g2o::BlockSolver_6_3::LinearSolverType * linearSolver;linearSolver = new g2o::LinearSolverEigen<g2o::BlockSolver_6_3::PoseMatrixType>();g2o::BlockSolver_6_3 * solver_ptr = new g2o::BlockSolver_6_3(linearSolver);g2o::OptimizationAlgorithmLevenberg* solver = new g2o::OptimizationAlgorithmLevenberg(solver_ptr);optimizer.setAlgorithm(solver);if(pbStopFlag)optimizer.setForceStopFlag(pbStopFlag);unsigned long maxKFid = 0;// Set Local KeyFrame vertices// 4. 设置局部关键帧顶点并加入优化器中。for(list<KeyFrame*>::iterator lit=lLocalKeyFrames.begin(), lend=lLocalKeyFrames.end(); lit!=lend; lit++){KeyFrame* pKFi = *lit;g2o::VertexSE3Expmap * vSE3 = new g2o::VertexSE3Expmap();vSE3->setEstimate(Converter::toSE3Quat(pKFi->GetPose()));vSE3->setId(pKFi->mnId);vSE3->setFixed(pKFi->mnId==0);optimizer.addVertex(vSE3);if(pKFi->mnId>maxKFid)maxKFid=pKFi->mnId;}// Set Fixed KeyFrame vertices// 5. 设置固定关键帧顶点,并加入到优化器列表中。for(list<KeyFrame*>::iterator lit=lFixedCameras.begin(), lend=lFixedCameras.end(); lit!=lend; lit++){KeyFrame* pKFi = *lit;g2o::VertexSE3Expmap * vSE3 = new g2o::VertexSE3Expmap();vSE3->setEstimate(Converter::toSE3Quat(pKFi->GetPose()));vSE3->setId(pKFi->mnId);vSE3->setFixed(true);optimizer.addVertex(vSE3);if(pKFi->mnId>maxKFid)maxKFid=pKFi->mnId;}// Set MapPoint vertices// 6.设置MapPoint顶点const int nExpectedSize = (lLocalKeyFrames.size()+lFixedCameras.size())*lLocalMapPoints.size();vector<g2o::EdgeSE3ProjectXYZ*> vpEdgesMono;vpEdgesMono.reserve(nExpectedSize);vector<KeyFrame*> vpEdgeKFMono;vpEdgeKFMono.reserve(nExpectedSize);vector<MapPoint*> vpMapPointEdgeMono;vpMapPointEdgeMono.reserve(nExpectedSize);vector<g2o::EdgeStereoSE3ProjectXYZ*> vpEdgesStereo;vpEdgesStereo.reserve(nExpectedSize);vector<KeyFrame*> vpEdgeKFStereo;vpEdgeKFStereo.reserve(nExpectedSize);vector<MapPoint*> vpMapPointEdgeStereo;vpMapPointEdgeStereo.reserve(nExpectedSize);const float thHuberMono = sqrt(5.991);const float thHuberStereo = sqrt(7.815);// 7.遍历局部地图点列表,设置优化对应的边for(list<MapPoint*>::iterator lit=lLocalMapPoints.begin(), lend=lLocalMapPoints.end(); lit!=lend; lit++){MapPoint* pMP = *lit;g2o::VertexSBAPointXYZ* vPoint = new g2o::VertexSBAPointXYZ();vPoint->setEstimate(Converter::toVector3d(pMP->GetWorldPos()));int id = pMP->mnId+maxKFid+1;vPoint->setId(id);vPoint->setMarginalized(true);optimizer.addVertex(vPoint);const map<KeyFrame*,size_t> observations = pMP->GetObservations();//Set edgesfor(map<KeyFrame*,size_t>::const_iterator mit=observations.begin(), mend=observations.end(); mit!=mend; mit++){KeyFrame* pKFi = mit->first;if(!pKFi->isBad()){                const cv::KeyPoint &kpUn = pKFi->mvKeysUn[mit->second];// Monocular observationif(pKFi->mvuRight[mit->second]<0){Eigen::Matrix<double,2,1> obs;obs << kpUn.pt.x, kpUn.pt.y;g2o::EdgeSE3ProjectXYZ* e = new g2o::EdgeSE3ProjectXYZ();e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(id)));//地图点e->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(pKFi->mnId)));//关键帧e->setMeasurement(obs);const float &invSigma2 = pKFi->mvInvLevelSigma2[kpUn.octave];e->setInformation(Eigen::Matrix2d::Identity()*invSigma2);g2o::RobustKernelHuber* rk = new g2o::RobustKernelHuber;e->setRobustKernel(rk);rk->setDelta(thHuberMono);e->fx = pKFi->fx;e->fy = pKFi->fy;e->cx = pKFi->cx;e->cy = pKFi->cy;optimizer.addEdge(e);vpEdgesMono.push_back(e);vpEdgeKFMono.push_back(pKFi);vpMapPointEdgeMono.push_back(pMP);}else // Stereo observation{Eigen::Matrix<double,3,1> obs;const float kp_ur = pKFi->mvuRight[mit->second];obs << kpUn.pt.x, kpUn.pt.y, kp_ur;g2o::EdgeStereoSE3ProjectXYZ* e = new g2o::EdgeStereoSE3ProjectXYZ();e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(id)));e->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(pKFi->mnId)));e->setMeasurement(obs);const float &invSigma2 = pKFi->mvInvLevelSigma2[kpUn.octave];Eigen::Matrix3d Info = Eigen::Matrix3d::Identity()*invSigma2;e->setInformation(Info);g2o::RobustKernelHuber* rk = new g2o::RobustKernelHuber;e->setRobustKernel(rk);rk->setDelta(thHuberStereo);e->fx = pKFi->fx;e->fy = pKFi->fy;e->cx = pKFi->cx;e->cy = pKFi->cy;e->bf = pKFi->mbf;optimizer.addEdge(e);vpEdgesStereo.push_back(e);vpEdgeKFStereo.push_back(pKFi);vpMapPointEdgeStereo.push_back(pMP);}}}}if(pbStopFlag)if(*pbStopFlag)return;//8.进行优化。这里进行了5次迭代优化。optimizer.initializeOptimization();optimizer.optimize(5);bool bDoMore= true;if(pbStopFlag)if(*pbStopFlag)bDoMore = false;// 9.检查上边5次迭代优化正确点的观测值,并将异常点排除掉。完成后再次进行10次迭代优化。if(bDoMore){// Check inlier observationsfor(size_t i=0, iend=vpEdgesMono.size(); i<iend;i++){g2o::EdgeSE3ProjectXYZ* e = vpEdgesMono[i];MapPoint* pMP = vpMapPointEdgeMono[i];if(pMP->isBad())continue;if(e->chi2()>5.991 || !e->isDepthPositive()){e->setLevel(1);}e->setRobustKernel(0);}for(size_t i=0, iend=vpEdgesStereo.size(); i<iend;i++){g2o::EdgeStereoSE3ProjectXYZ* e = vpEdgesStereo[i];MapPoint* pMP = vpMapPointEdgeStereo[i];if(pMP->isBad())continue;if(e->chi2()>7.815 || !e->isDepthPositive()){e->setLevel(1);}e->setRobustKernel(0);}// Optimize again without the outliers//排除掉不好的点后,再次进行优化optimizer.initializeOptimization(0);optimizer.optimize(10);}vector<pair<KeyFrame*,MapPoint*> > vToErase;vToErase.reserve(vpEdgesMono.size()+vpEdgesStereo.size());// Check inlier observations   //10.优化后进行优化结果检查    for(size_t i=0, iend=vpEdgesMono.size(); i<iend;i++){g2o::EdgeSE3ProjectXYZ* e = vpEdgesMono[i];MapPoint* pMP = vpMapPointEdgeMono[i];if(pMP->isBad())continue;if(e->chi2()>5.991 || !e->isDepthPositive()){KeyFrame* pKFi = vpEdgeKFMono[i];vToErase.push_back(make_pair(pKFi,pMP));}}for(size_t i=0, iend=vpEdgesStereo.size(); i<iend;i++){g2o::EdgeStereoSE3ProjectXYZ* e = vpEdgesStereo[i];MapPoint* pMP = vpMapPointEdgeStereo[i];if(pMP->isBad())continue;if(e->chi2()>7.815 || !e->isDepthPositive()){KeyFrame* pKFi = vpEdgeKFStereo[i];vToErase.push_back(make_pair(pKFi,pMP));}}// Get Map Mutexunique_lock<mutex> lock(pMap->mMutexMapUpdate);if(!vToErase.empty()){for(size_t i=0;i<vToErase.size();i++){KeyFrame* pKFi = vToErase[i].first;MapPoint* pMPi = vToErase[i].second;pKFi->EraseMapPointMatch(pMPi);pMPi->EraseObservation(pKFi);}}// Recover optimized data//11.设置优化后的数据。//Keyframesfor(list<KeyFrame*>::iterator lit=lLocalKeyFrames.begin(), lend=lLocalKeyFrames.end(); lit!=lend; lit++){KeyFrame* pKF = *lit;g2o::VertexSE3Expmap* vSE3 = static_cast<g2o::VertexSE3Expmap*>(optimizer.vertex(pKF->mnId));//优化后的关键帧位姿g2o::SE3Quat SE3quat = vSE3->estimate();pKF->SetPose(Converter::toCvMat(SE3quat));}//Pointsfor(list<MapPoint*>::iterator lit=lLocalMapPoints.begin(), lend=lLocalMapPoints.end(); lit!=lend; lit++){MapPoint* pMP = *lit;g2o::VertexSBAPointXYZ* vPoint = static_cast<g2o::VertexSBAPointXYZ*>(optimizer.vertex(pMP->mnId+maxKFid+1));//优化后的地图点世界坐标pMP->SetWorldPos(Converter::toCvMat(vPoint->estimate()));pMP->UpdateNormalAndDepth();}
}

4.本质图优化OptimizeEssentialGraph

一些说明:函数名称是OptimizeEssentialGraph,那么只有形成了EssentialGraph之后才需要进行优化。EssentialGraph中不但有一般的共视关键帧之间的边,还有闭环边和共视图边。所以其中加入了帧间闭环的考虑,也加入了纠正后的帧的SIM3位姿。所有共视的关键帧、地图点、一般的边和闭环边放在一起进行优化,优化的目标是关键帧的sim3位姿和地图点的坐标。 将优化后的sim3位姿转换为SE3,就得到了相机的位姿,并对地图点的坐标进行投影得到更新后的坐标。

调用时机:LoopClosing线程当中通过sim3调整完关键帧和其闭环帧之间位姿后调用。

代码:

/** 参数解释:* KeyFrame* pLoopKF 为和当前关键帧pCurKF形成闭环的关键帧* KeyFrame* pCurKF  当前关键帧* NonCorrectedSim3 当前关键帧及其共视关键帧最初的位姿* CorrectedSim3  当前关键帧的共视关键帧进行sim3纠正后的位姿* LoopConnections 其中的index为与当前关键帧有共视关系的关键帧*                 值为set<KeyFrame *>表示与当前关键帧的共视关键帧相连接的关键帧
*/
void Optimizer::OptimizeEssentialGraph(Map* pMap, KeyFrame* pLoopKF, KeyFrame* pCurKF,const LoopClosing::KeyFrameAndPose &NonCorrectedSim3,const LoopClosing::KeyFrameAndPose &CorrectedSim3,const map<KeyFrame *, set<KeyFrame *> > &LoopConnections, const bool &bFixScale)
{// Setup optimizerg2o::SparseOptimizer optimizer;optimizer.setVerbose(false);//这里为什么选择7*3的矩阵?因为除了在三个轴方向上的旋转+平移外,还有尺度g2o::BlockSolver_7_3::LinearSolverType * linearSolver =new g2o::LinearSolverEigen<g2o::BlockSolver_7_3::PoseMatrixType>();g2o::BlockSolver_7_3 * solver_ptr= new g2o::BlockSolver_7_3(linearSolver);g2o::OptimizationAlgorithmLevenberg* solver = new g2o::OptimizationAlgorithmLevenberg(solver_ptr);solver->setUserLambdaInit(1e-16);optimizer.setAlgorithm(solver);//获取地图中所有的关键帧和地图点const vector<KeyFrame*> vpKFs = pMap->GetAllKeyFrames();const vector<MapPoint*> vpMPs = pMap->GetAllMapPoints();const unsigned int nMaxKFid = pMap->GetMaxKFid();vector<g2o::Sim3,Eigen::aligned_allocator<g2o::Sim3> > vScw(nMaxKFid+1);vector<g2o::Sim3,Eigen::aligned_allocator<g2o::Sim3> > vCorrectedSwc(nMaxKFid+1);vector<g2o::VertexSim3Expmap*> vpVertices(nMaxKFid+1);const int minFeat = 100;// Set KeyFrame verticesfor(size_t i=0, iend=vpKFs.size(); i<iend;i++){KeyFrame* pKF = vpKFs[i];if(pKF->isBad())continue;g2o::VertexSim3Expmap* VSim3 = new g2o::VertexSim3Expmap();const int nIDi = pKF->mnId;LoopClosing::KeyFrameAndPose::const_iterator it = CorrectedSim3.find(pKF);if(it!=CorrectedSim3.end()){//获取nIDi这个index的关键帧的纠正过的位姿vScw[nIDi] = it->second;VSim3->setEstimate(it->second);}else{Eigen::Matrix<double,3,3> Rcw = Converter::toMatrix3d(pKF->GetRotation());Eigen::Matrix<double,3,1> tcw = Converter::toVector3d(pKF->GetTranslation());g2o::Sim3 Siw(Rcw,tcw,1.0);vScw[nIDi] = Siw;VSim3->setEstimate(Siw);}//pLoopKF为查找到的匹配闭环帧if(pKF==pLoopKF)VSim3->setFixed(true);VSim3->setId(nIDi);VSim3->setMarginalized(false);VSim3->_fix_scale = bFixScale;optimizer.addVertex(VSim3);vpVertices[nIDi]=VSim3;}set<pair<long unsigned int,long unsigned int> > sInsertedEdges;const Eigen::Matrix<double,7,7> matLambda = Eigen::Matrix<double,7,7>::Identity();// Set Loop edgesfor(map<KeyFrame *, set<KeyFrame *> >::const_iterator mit = LoopConnections.begin(), mend=LoopConnections.end(); mit!=mend; mit++){//pKF为和当前关键帧有共视关系的关键帧KeyFrame* pKF = mit->first;const long unsigned int nIDi = pKF->mnId;//mit->second表示和pKF相连接的关键帧const set<KeyFrame*> &spConnections = mit->second;const g2o::Sim3 Siw = vScw[nIDi];const g2o::Sim3 Swi = Siw.inverse();for(set<KeyFrame*>::const_iterator sit=spConnections.begin(), send=spConnections.end(); sit!=send; sit++){const long unsigned int nIDj = (*sit)->mnId;if((nIDi!=pCurKF->mnId || nIDj!=pLoopKF->mnId) && pKF->GetWeight(*sit)<minFeat)continue;const g2o::Sim3 Sjw = vScw[nIDj];const g2o::Sim3 Sji = Sjw * Swi;g2o::EdgeSim3* e = new g2o::EdgeSim3();//两个顶点一个为当前关键帧共视的关键帧,另一个为共视关键帧e->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(nIDj)));e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(nIDi)));e->setMeasurement(Sji);e->information() = matLambda;optimizer.addEdge(e);sInsertedEdges.insert(make_pair(min(nIDi,nIDj),max(nIDi,nIDj)));}}// Set normal edgesfor(size_t i=0, iend=vpKFs.size(); i<iend; i++){KeyFrame* pKF = vpKFs[i];const int nIDi = pKF->mnId;g2o::Sim3 Swi;LoopClosing::KeyFrameAndPose::const_iterator iti = NonCorrectedSim3.find(pKF);if(iti!=NonCorrectedSim3.end())Swi = (iti->second).inverse();elseSwi = vScw[nIDi].inverse();KeyFrame* pParentKF = pKF->GetParent();// Spanning tree edgeif(pParentKF){int nIDj = pParentKF->mnId;g2o::Sim3 Sjw;LoopClosing::KeyFrameAndPose::const_iterator itj = NonCorrectedSim3.find(pParentKF);if(itj!=NonCorrectedSim3.end())Sjw = itj->second;elseSjw = vScw[nIDj];g2o::Sim3 Sji = Sjw * Swi;g2o::EdgeSim3* e = new g2o::EdgeSim3();e->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(nIDj)));e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(nIDi)));e->setMeasurement(Sji);e->information() = matLambda;optimizer.addEdge(e);}// Loop edgesconst set<KeyFrame*> sLoopEdges = pKF->GetLoopEdges();for(set<KeyFrame*>::const_iterator sit=sLoopEdges.begin(), send=sLoopEdges.end(); sit!=send; sit++){KeyFrame* pLKF = *sit;if(pLKF->mnId<pKF->mnId){g2o::Sim3 Slw;LoopClosing::KeyFrameAndPose::const_iterator itl = NonCorrectedSim3.find(pLKF);if(itl!=NonCorrectedSim3.end())Slw = itl->second;elseSlw = vScw[pLKF->mnId];g2o::Sim3 Sli = Slw * Swi;g2o::EdgeSim3* el = new g2o::EdgeSim3();el->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(pLKF->mnId)));el->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(nIDi)));el->setMeasurement(Sli);el->information() = matLambda;optimizer.addEdge(el);}}// Covisibility graph edgesconst vector<KeyFrame*> vpConnectedKFs = pKF->GetCovisiblesByWeight(minFeat);for(vector<KeyFrame*>::const_iterator vit=vpConnectedKFs.begin(); vit!=vpConnectedKFs.end(); vit++){KeyFrame* pKFn = *vit;if(pKFn && pKFn!=pParentKF && !pKF->hasChild(pKFn) && !sLoopEdges.count(pKFn)){if(!pKFn->isBad() && pKFn->mnId<pKF->mnId){if(sInsertedEdges.count(make_pair(min(pKF->mnId,pKFn->mnId),max(pKF->mnId,pKFn->mnId))))continue;g2o::Sim3 Snw;LoopClosing::KeyFrameAndPose::const_iterator itn = NonCorrectedSim3.find(pKFn);if(itn!=NonCorrectedSim3.end())Snw = itn->second;elseSnw = vScw[pKFn->mnId];g2o::Sim3 Sni = Snw * Swi;g2o::EdgeSim3* en = new g2o::EdgeSim3();en->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(pKFn->mnId)));en->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(nIDi)));en->setMeasurement(Sni);en->information() = matLambda;optimizer.addEdge(en);}}}}// Optimize!optimizer.initializeOptimization();optimizer.optimize(20);unique_lock<mutex> lock(pMap->mMutexMapUpdate);// SE3 Pose Recovering. Sim3:[sR t;0 1] -> SE3:[R t/s;0 1]for(size_t i=0;i<vpKFs.size();i++){KeyFrame* pKFi = vpKFs[i];const int nIDi = pKFi->mnId;g2o::VertexSim3Expmap* VSim3 = static_cast<g2o::VertexSim3Expmap*>(optimizer.vertex(nIDi));g2o::Sim3 CorrectedSiw =  VSim3->estimate();vCorrectedSwc[nIDi]=CorrectedSiw.inverse();Eigen::Matrix3d eigR = CorrectedSiw.rotation().toRotationMatrix();Eigen::Vector3d eigt = CorrectedSiw.translation();double s = CorrectedSiw.scale();eigt *=(1./s); //[R t/s;0 1]cv::Mat Tiw = Converter::toCvSE3(eigR,eigt);pKFi->SetPose(Tiw);}// Correct points. Transform to "non-optimized" reference keyframe pose and transform back with optimized posefor(size_t i=0, iend=vpMPs.size(); i<iend; i++){MapPoint* pMP = vpMPs[i];if(pMP->isBad())continue;int nIDr;if(pMP->mnCorrectedByKF==pCurKF->mnId){nIDr = pMP->mnCorrectedReference;}else{KeyFrame* pRefKF = pMP->GetReferenceKeyFrame();nIDr = pRefKF->mnId;}g2o::Sim3 Srw = vScw[nIDr];g2o::Sim3 correctedSwr = vCorrectedSwc[nIDr];cv::Mat P3Dw = pMP->GetWorldPos();Eigen::Matrix<double,3,1> eigP3Dw = Converter::toVector3d(P3Dw);Eigen::Matrix<double,3,1> eigCorrectedP3Dw = correctedSwr.map(Srw.map(eigP3Dw));cv::Mat cvCorrectedP3Dw = Converter::toCvMat(eigCorrectedP3Dw);pMP->SetWorldPos(cvCorrectedP3Dw);pMP->UpdateNormalAndDepth();}
}

5.全局BA优化GlobalBundleAdjustment

一些说明:这个是ORB-SLAM系统当中最大的优化,将全局地图当中所有的关键帧和地图点都放进来一起进行优化。对关键帧的位姿和地图点3D坐标进行优化。

调用时机:LoopClosing线程当中所有工作完成后创建线程进行整体优化。函数调用关系为:RunGlobalBundleAdjustment-->Optimizer::GlobalBundleAdjustemnt-->BundleAdjustment

代码:

RunGlobalBundleAdjustment函数:

void LoopClosing::RunGlobalBundleAdjustment(unsigned long nLoopKF)
{cout << "Starting Global Bundle Adjustment" << endl;int idx =  mnFullBAIdx;Optimizer::GlobalBundleAdjustemnt(mpMap,10,&mbStopGBA,nLoopKF,false);// Update all MapPoints and KeyFrames// Local Mapping was active during BA, that means that there might be new keyframes// not included in the Global BA and they are not consistent with the updated map.// We need to propagate the correction through the spanning tree{unique_lock<mutex> lock(mMutexGBA);if(idx!=mnFullBAIdx)return;if(!mbStopGBA){cout << "Global Bundle Adjustment finished" << endl;cout << "Updating map ..." << endl;mpLocalMapper->RequestStop();// Wait until Local Mapping has effectively stoppedwhile(!mpLocalMapper->isStopped() && !mpLocalMapper->isFinished()){usleep(1000);}// Get Map Mutexunique_lock<mutex> lock(mpMap->mMutexMapUpdate);// Correct keyframes starting at map first keyframelist<KeyFrame*> lpKFtoCheck(mpMap->mvpKeyFrameOrigins.begin(),mpMap->mvpKeyFrameOrigins.end());while(!lpKFtoCheck.empty()){KeyFrame* pKF = lpKFtoCheck.front();const set<KeyFrame*> sChilds = pKF->GetChilds();cv::Mat Twc = pKF->GetPoseInverse();for(set<KeyFrame*>::const_iterator sit=sChilds.begin();sit!=sChilds.end();sit++){KeyFrame* pChild = *sit;if(pChild->mnBAGlobalForKF!=nLoopKF){cv::Mat Tchildc = pChild->GetPose()*Twc;pChild->mTcwGBA = Tchildc*pKF->mTcwGBA;//*Tcorc*pKF->mTcwGBA;pChild->mnBAGlobalForKF=nLoopKF;}lpKFtoCheck.push_back(pChild);}pKF->mTcwBefGBA = pKF->GetPose();pKF->SetPose(pKF->mTcwGBA);lpKFtoCheck.pop_front();}// Correct MapPointsconst vector<MapPoint*> vpMPs = mpMap->GetAllMapPoints();for(size_t i=0; i<vpMPs.size(); i++){MapPoint* pMP = vpMPs[i];if(pMP->isBad())continue;if(pMP->mnBAGlobalForKF==nLoopKF){// If optimized by Global BA, just updatepMP->SetWorldPos(pMP->mPosGBA);}else{// Update according to the correction of its reference keyframeKeyFrame* pRefKF = pMP->GetReferenceKeyFrame();if(pRefKF->mnBAGlobalForKF!=nLoopKF)continue;// Map to non-corrected cameracv::Mat Rcw = pRefKF->mTcwBefGBA.rowRange(0,3).colRange(0,3);cv::Mat tcw = pRefKF->mTcwBefGBA.rowRange(0,3).col(3);cv::Mat Xc = Rcw*pMP->GetWorldPos()+tcw;// Backproject using corrected cameracv::Mat Twc = pRefKF->GetPoseInverse();cv::Mat Rwc = Twc.rowRange(0,3).colRange(0,3);cv::Mat twc = Twc.rowRange(0,3).col(3);pMP->SetWorldPos(Rwc*Xc+twc);}}            mpMap->InformNewBigChange();mpLocalMapper->Release();cout << "Map updated!" << endl;}mbFinishedGBA = true;mbRunningGBA = false;}
}

Optimizer::GlobalBundleAdjustemnt函数:

/*** 优化相机位姿和关键点坐标* 优化中使用的几个数据结构解释如下:* VertexSE3Expmap李代数位姿* VertexSBAPointXYZ 空间点位置* EdgeProjectXYZ2UV 投影方程边*
*/
void Optimizer::GlobalBundleAdjustemnt(Map* pMap, int nIterations, bool* pbStopFlag, const unsigned long nLoopKF, const bool bRobust)
{vector<KeyFrame*> vpKFs = pMap->GetAllKeyFrames();vector<MapPoint*> vpMP = pMap->GetAllMapPoints();BundleAdjustment(vpKFs,vpMP,nIterations,pbStopFlag, nLoopKF, bRobust);
}

Optimizer::BundleAdjustment函数:

/*** vpKFs 地图中所有的关键帧* vpMP  地图中所有的MapPoint* 将3D点投影至图像平面优化投影误差。* 初始化g2o优化器,添加所有的关键帧位姿作为顶点,所有的地图点作为顶点,同时对每一个地图点获取其观测信息,之后建立点-关键帧之间的投影关系作为优化器的边
*/
void Optimizer::BundleAdjustment(const vector<KeyFrame *> &vpKFs, const vector<MapPoint *> &vpMP,int nIterations, bool* pbStopFlag, const unsigned long nLoopKF, const bool bRobust)
{vector<bool> vbNotIncludedMP;vbNotIncludedMP.resize(vpMP.size());g2o::SparseOptimizer optimizer;//第一步:创建一个线性求解器LinearSolver。指定pose维度为6, landmark维度为3g2o::BlockSolver_6_3::LinearSolverType * linearSolver;linearSolver = new g2o::LinearSolverEigen<g2o::BlockSolver_6_3::PoseMatrixType>();//第二步:创建BlockSolver,并用上边定义的线性求解器初始化。这里是个稀疏矩阵块求解器。g2o::BlockSolver_6_3 * solver_ptr = new g2o::BlockSolver_6_3(linearSolver);//第三步:创建总求解器solver,并从GN, LM, DogLeg中选一个梯度下降算法,再用上述块求解器BlockSolver初始化,这里选择ML算法g2o::OptimizationAlgorithmLevenberg* solver = new g2o::OptimizationAlgorithmLevenberg(solver_ptr);//第四步:创建终极大boss稀疏优化器(SparseOptimizer),并进行设置optimizer.setAlgorithm(solver);//设置求解器if(pbStopFlag)optimizer.setForceStopFlag(pbStopFlag);long unsigned int maxKFid = 0;//第五步:定义图的顶点和边,并添加到SparseOptimizer中// Set KeyFrame vertices 设置关键帧顶点并加入优化器for(size_t i=0; i<vpKFs.size(); i++){KeyFrame* pKF = vpKFs[i];if(pKF->isBad())continue;//5.1 VertexSE3Expmap创建相机位姿结点,并加入优化器g2o::VertexSE3Expmap * vSE3 = new g2o::VertexSE3Expmap();vSE3->setEstimate(Converter::toSE3Quat(pKF->GetPose()));//优化的是位姿,也就是变换矩阵vSE3->setId(pKF->mnId);vSE3->setFixed(pKF->mnId==0);optimizer.addVertex(vSE3);//maxKFid用于记录关键帧id的最大值if(pKF->mnId > maxKFid)maxKFid=pKF->mnId;}const float thHuber2D = sqrt(5.99);const float thHuber3D = sqrt(7.815);// Set MapPoint vertices// 设置MapPoint顶点并加入优化器。遍历MapPointsfor(size_t i=0; i<vpMP.size(); i++){MapPoint* pMP = vpMP[i];if(pMP->isBad())continue;//5.2VertexSBAPointXYZ创建特征点空间坐标结点,并加入优化器g2o::VertexSBAPointXYZ* vPoint = new g2o::VertexSBAPointXYZ();//【1】设置待优化空间点3D位置XYZvPoint->setEstimate(Converter::toVector3d(pMP->GetWorldPos()));//id为什么是MapPoint的id加上关键帧的最大id还要加个1?const int id = pMP->mnId + maxKFid+1;//【2】设置Id号vPoint->setId(id);//【3】是否边缘化(以便稀疏化求解)vPoint->setMarginalized(true);optimizer.addVertex(vPoint);//获取能够观察到该MapPoint的关键帧列表const map<KeyFrame*,size_t> observations = pMP->GetObservations();int nEdges = 0;//SET EDGES//设置边并加入优化器。遍历每个MapPoint对应的关键帧,获取该MapPoint在能看到它的关键帧中对应的关键点for(map<KeyFrame*,size_t>::const_iterator mit=observations.begin(); mit!=observations.end(); mit++){KeyFrame* pKF = mit->first;if(pKF->isBad() || pKF->mnId>maxKFid)continue;//计算边数nEdges++;//mit->second为pMP这个MapPoint在pKF这个关键帧中对应的关键点的index。const cv::KeyPoint &kpUn = pKF->mvKeysUn[mit->second];//单目相机时pKF->mvuRight为负数if(pKF->mvuRight[mit->second]<0){Eigen::Matrix<double,2,1> obs;//obs中存储的是关键点坐标obs << kpUn.pt.x, kpUn.pt.y;/*** EdgeSE3ProjectXYZ为(Point-Pose)二元边,既要优化MapPoints的位置,又要优化相机位姿* class  EdgeSE3ProjectXYZ: public  BaseBinaryEdge<2, Vector2d, VertexSBAPointXYZ, VertexSE3Expmap>* 继承于BaseBinaryEdge这个二元边模板类,需要设置的模板参数:* 参数2 :观测值(这里是3D点在像素坐标系下的投影坐标)的维度* 参数Vector :观测值类型,piexl.x,piexl.y* 参数VertexSBAPointXYZ:第一个顶点类型  MapPoint类型* 参数VertexSE3Expmap :第二个顶点类型  关键帧类型* */g2o::EdgeSE3ProjectXYZ* e = new g2o::EdgeSE3ProjectXYZ();//【1】设置第一个顶点,注意该顶点类型与模板参数第一个顶点类型对应e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(id)));//【2】设置第二个顶点e->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(pKF->mnId)));//【3】设置观测值,类型与模板参数对应e->setMeasurement(obs);//kpUn.octave表示该关键点所处的金字塔图的第几层const float &invSigma2 = pKF->mvInvLevelSigma2[kpUn.octave];//【4】设置信息矩阵,也就是协方差矩阵e->setInformation(Eigen::Matrix2d::Identity()*invSigma2);if(bRobust){/*** 【5】设置鲁棒核函数。* 之所以要设置鲁棒核函数是为了平衡误差,不让二范数的误差增加的过快。* 鲁棒核函数里要自己设置delta值,这个delta值是,当误差的绝对值小于等于它的时候,误差函数不变。* */g2o::RobustKernelHuber* rk = new g2o::RobustKernelHuber;e->setRobustKernel(rk);rk->setDelta(thHuber2D);}//【6】设置相机内参e->fx = pKF->fx;e->fy = pKF->fy;e->cx = pKF->cx;e->cy = pKF->cy;optimizer.addEdge(e);}else//双目{Eigen::Matrix<double,3,1> obs;const float kp_ur = pKF->mvuRight[mit->second];obs << kpUn.pt.x, kpUn.pt.y, kp_ur;g2o::EdgeStereoSE3ProjectXYZ* e = new g2o::EdgeStereoSE3ProjectXYZ();e->setVertex(0, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(id)));e->setVertex(1, dynamic_cast<g2o::OptimizableGraph::Vertex*>(optimizer.vertex(pKF->mnId)));e->setMeasurement(obs);const float &invSigma2 = pKF->mvInvLevelSigma2[kpUn.octave];Eigen::Matrix3d Info = Eigen::Matrix3d::Identity()*invSigma2;e->setInformation(Info);if(bRobust){g2o::RobustKernelHuber* rk = new g2o::RobustKernelHuber;e->setRobustKernel(rk);rk->setDelta(thHuber3D);}e->fx = pKF->fx;e->fy = pKF->fy;e->cx = pKF->cx;e->cy = pKF->cy;e->bf = pKF->mbf;optimizer.addEdge(e);}}//如果图中的边数为0,则删除图中的结点if(nEdges==0){optimizer.removeVertex(vPoint);vbNotIncludedMP[i]=true;}else{vbNotIncludedMP[i]=false;}}// Optimize!//第六步:设置优化参数,开始执行优化optimizer.initializeOptimization();//设置迭代次数optimizer.optimize(nIterations);// Recover optimized data//Keyframes//获取优化后keyframe的位姿for(size_t i=0; i<vpKFs.size(); i++){KeyFrame* pKF = vpKFs[i];if(pKF->isBad())continue;g2o::VertexSE3Expmap* vSE3 = static_cast<g2o::VertexSE3Expmap*>(optimizer.vertex(pKF->mnId));g2o::SE3Quat SE3quat = vSE3->estimate();if(nLoopKF==0){pKF->SetPose(Converter::toCvMat(SE3quat));}else{pKF->mTcwGBA.create(4,4,CV_32F);Converter::toCvMat(SE3quat).copyTo(pKF->mTcwGBA);pKF->mnBAGlobalForKF = nLoopKF;}}//Points//获取优化后MapPoint的坐标for(size_t i=0; i<vpMP.size(); i++){if(vbNotIncludedMP[i])continue;MapPoint* pMP = vpMP[i];if(pMP->isBad())continue;g2o::VertexSBAPointXYZ* vPoint = static_cast<g2o::VertexSBAPointXYZ*>(optimizer.vertex(pMP->mnId+maxKFid+1));if(nLoopKF==0){pMP->SetWorldPos(Converter::toCvMat(vPoint->estimate()));pMP->UpdateNormalAndDepth();}else{pMP->mPosGBA.create(3,1,CV_32F);Converter::toCvMat(vPoint->estimate()).copyTo(pMP->mPosGBA);pMP->mnBAGlobalForKF = nLoopKF;}}}

ORB-SLAM当中的优化,在论文的附录中有进行解释。

ORB-SLAM2系统中的优化函数相关推荐

  1. orbslam2可视化_[Ubuntu] ORB SLAM2 编译调试

    ORB SLAM2 是 2015年比较受到关注的一篇文章,它的主要思想是借助 ORB 描述子改进了 Sparse SLAM 的性能,使得其在稳定性和速度上都达到了比较好的程度.从创新性上来讲,它的主要 ...

  2. ORB SLAM2源码解读(三):Frame类

    文章目录 前言 构造函数 双目相机 RGBD相机 单目相机 ExtractORB:提取特征点 ComputeBoW:计算词袋数据 SetPose:设置相机外参 isInFrustum:判断一个MapP ...

  3. Ubuntu下使用单目相机运行ORB SLAM2

    环境:Ubuntu16.04+ROS Kinetic+USB单目摄像头 虽然ORB SLAM2的官方说明中表示没有ROS也可以编译运行,但要实时的跑ORB SLAM2还是需要ROS平台的,所以之前没有 ...

  4. 先进机器人系统中的关键技术

    先进机器人系统中的关键技术 Key technologies coalesce in advanced robotic systems 就在机器人设计进入商业领域,服务于制造业.物流业和服务业之际,概 ...

  5. Linux系统中创建大文件,并作为文件系统使用

    在LInux系统的使用过程中,有时候会遇到诸如某个磁盘分区的大小不够用了,导致其下的文件系统不能正常写入数据.亦或者是系统swap分区太小,不够用或者不满足条件而导致的其他一系列问题.如果我们系统上挂 ...

  6. ip设置 kali 重置_在 Windows 系统中如何重置 TCP/IP 协议堆栈修复网络连接问题

    Internet 在 TCP/IP 协议上工作,如果 TCP/IP 协议堆栈在 Windows 或任何其他操作系统(例如 Linux 或 MacOS)中无法正常工作,则您的 Internet 连接会出 ...

  7. linux ubuntu安装 mono,在Ubuntu 18.04系统中安装Mono及基本使用Mono的方法

    本文介绍在Ubuntu 18.04操作系统中安装Mono及基本使用Mono的方法.Mono是一个基于ECMA/ISO标准开发和运行跨平台应用程序的平台,它是Microsoft .NET框架的免费开源实 ...

  8. ROS系统中实现点云聚类(realsense数据源)

    本文主要介绍ROS系统中如何订阅并解码realsense点云数据,并对点云进行稀疏.去噪.聚类. 环境配置见<ROS系统中从零开始部署YoloV4目标检测算法(3种方式)> 需要安装的第三 ...

  9. win10html5无法播放,win10系统中网页中无法播放视频怎么办

    近日有win10系统用户要通过浏览器来打开网页观看视频的时候,却发现在网页中打开视频的时却无法播放,这是怎么回事呢,经过分析是由于Adobe Flash Player ActiveX插件未安装.版本过 ...

最新文章

  1. linux增量编译不成功,Linux学习笔记-增量编译(Makefile进一步使用)
  2. 74. Search a 2D Matrix (Graph; Divide-and-Conquer)
  3. Django笔记 —— 表单(form)
  4. ios开发--企业帐号发布
  5. 【2013高考作文】重庆作文
  6. 启动dubbo-admin遇到的那些坑
  7. 文件资源管理软件EagleFiler for Mac
  8. 《算法导论》:关于循环不变式
  9. 面试宝典之深度学习面试题(下)
  10. [转]游戏多开的原理
  11. 云桌面终端CT3200,硬件与信号连接
  12. 手写一个简单的mybatis
  13. 为什么用手机拍摄电脑屏幕时会出现波纹?终于明白了!
  14. 2021年6月PMP考试50天备考5A通过经历心得分享
  15. 利用Basemap画世界地图
  16. Nginx下载、安装与使用
  17. 状态空间方程转换传递函数
  18. Expressing Fear(2019/1/4)
  19. python中的结束用语_python结束语句
  20. Kotlin学习教程(一)

热门文章

  1. 怎么看SaaS企业中的收入留存率?
  2. Python读写文件rb,wb,ab模式
  3. YAML实践指南:3:YAML格式检查与转换
  4. latex并排显示多个图片
  5. 用css动态实现圆环百分比分配——初探css3动画
  6. 基于linux的贪吃蛇游戏设计_贪吃蛇还能这么玩?绝对是你从未体验过的全新版本(上)...
  7. Devart Excel Add-ins Crack
  8. 关于程序员日常接单之淘宝运营
  9. ProSci 艾美捷CCR3抗体
  10. 计算机怎样检查视力,电脑视力表同样测视力 测试方法要正确