36 #include <rtabmap/core/Version.h>
43 #if defined(RTABMAP_G2O) || defined(RTABMAP_ORB_SLAM)
44 #include "g2o/core/sparse_optimizer.h"
45 #include "g2o/core/block_solver.h"
46 #include "g2o/core/factory.h"
47 #include "g2o/core/optimization_algorithm_factory.h"
48 #include "g2o/core/optimization_algorithm_gauss_newton.h"
49 #include "g2o/core/optimization_algorithm_levenberg.h"
50 #include "g2o/core/robust_kernel_impl.h"
57 #include "g2o/types/sba/types_sba.h"
58 #include "g2o/solvers/eigen/linear_solver_eigen.h"
59 #include "g2o/config.h"
60 #include "g2o/types/slam2d/types_slam2d.h"
61 #include "g2o/types/slam3d/types_slam3d.h"
67 #ifdef G2O_HAVE_CSPARSE
68 #include "g2o/solvers/csparse/linear_solver_csparse.h"
70 #include "g2o/solvers/pcg/linear_solver_pcg.h"
71 #ifdef G2O_HAVE_CHOLMOD
72 #include "g2o/solvers/cholmod/linear_solver_cholmod.h"
76 #ifdef RTABMAP_ORB_SLAM
77 #include "g2o/types/types_sba.h"
78 #include "g2o/types/types_six_dof_expmap.h"
79 #include "g2o/solvers/linear_solver_eigen.h"
82 typedef g2o::BlockSolver< g2o::BlockSolverTraits<-1, -1> > SlamBlockSolver;
83 typedef g2o::LinearSolverEigen<SlamBlockSolver::PoseMatrixType> SlamLinearEigenSolver;
85 typedef g2o::LinearSolverPCG<SlamBlockSolver::PoseMatrixType> SlamLinearPCGSolver;
86 #ifdef G2O_HAVE_CSPARSE
87 typedef g2o::LinearSolverCSparse<SlamBlockSolver::PoseMatrixType> SlamLinearCSparseSolver;
89 #ifdef G2O_HAVE_CHOLMOD
90 typedef g2o::LinearSolverCholmod<SlamBlockSolver::PoseMatrixType> SlamLinearCholmodSolver;
100 #if defined(G2O_SRC_DIR) or defined(WIN32)
102 typedef VertexPointXYZ VertexSBAPointXYZ;
106 #if defined(RTABMAP_VERTIGO)
114 #endif // end defined(RTABMAP_G2O) || defined(RTABMAP_ORB_SLAM)
120 #define MULTICAM_OFFSET 10 // 10 means max 10 cameras per pose
126 #if defined(RTABMAP_G2O) || defined(RTABMAP_ORB_SLAM)
135 #ifdef G2O_HAVE_CSPARSE
144 #ifdef G2O_HAVE_CHOLMOD
154 optimizer_(
Parameters::defaultg2oOptimizer()),
155 pixelVariance_(
Parameters::defaultg2oPixelVariance()),
156 robustKernelDelta_(
Parameters::defaultg2oRobustKernelDelta()),
161 if(!g2o::Factory::instance()->knowsTag(
"CACHE_SE3_OFFSET"))
163 #if defined(RTABMAP_G2O_CPP11) && RTABMAP_G2O_CPP11 == 1
164 g2o::Factory::instance()->registerType(
"CACHE_SE3_OFFSET", std::make_unique<g2o::HyperGraphElementCreator<g2o::CacheSE3Offset> >());
166 g2o::Factory::instance()->registerType(
"CACHE_SE3_OFFSET",
new g2o::HyperGraphElementCreator<g2o::CacheSE3Offset>);
185 #ifdef RTABMAP_ORB_SLAM
188 UWARN(
"g2o built with ORB_SLAM has only Eigen solver available, using Eigen=3 instead of %d.",
solver_);
192 #ifndef G2O_HAVE_CHOLMOD
195 UWARN(
"g2o is not built with chmold, so it cannot be used as solver. Using CSparse instead.");
200 #ifndef G2O_HAVE_CSPARSE
203 UWARN(
"g2o is not built with csparse, so it cannot be used as solver. Using PCG instead.");
213 const std::map<int, Transform> & poses,
214 const std::multimap<int, Link> & edgeConstraints,
215 cv::Mat & outputCovariance,
216 std::list<std::map<int, Transform> > * intermediateGraphes,
218 int * iterationsDone)
220 outputCovariance = cv::Mat::eye(6,6,CV_64FC1);
221 std::map<int, Transform> optimizedPoses;
223 UDEBUG(
"Optimizing graph...");
225 #ifndef RTABMAP_VERTIGO
228 UWARN(
"Vertigo robust optimization is not available! Robust optimization is now disabled.");
233 optimizedPoses.clear();
234 if(edgeConstraints.size()>=1 && poses.size()>=2 &&
iterations() > 0 && poses.rbegin()->first > 0)
238 g2o::SparseOptimizer optimizer;
242 g2o::ParameterSE2Offset* odomOffset =
new g2o::ParameterSE2Offset();
244 optimizer.addParameter(odomOffset);
248 g2o::ParameterSE3Offset* odomOffset =
new g2o::ParameterSE3Offset();
250 optimizer.addParameter(odomOffset);
253 #ifdef RTABMAP_G2O_CPP11
255 std::unique_ptr<SlamBlockSolver> blockSolver;
260 auto linearSolver = std::make_unique<SlamLinearEigenSolver>();
261 linearSolver->setBlockOrdering(
false);
262 blockSolver = std::make_unique<SlamBlockSolver>(std::move(linearSolver));
264 #ifdef G2O_HAVE_CHOLMOD
268 auto linearSolver = std::make_unique<SlamLinearCholmodSolver>();
269 linearSolver->setBlockOrdering(
false);
270 blockSolver = std::make_unique<SlamBlockSolver>(std::move(linearSolver));
273 #ifdef G2O_HAVE_CSPARSE
278 auto linearSolver = std::make_unique<SlamLinearCSparseSolver>();
279 linearSolver->setBlockOrdering(
false);
280 blockSolver = std::make_unique<SlamBlockSolver>(std::move(linearSolver));
286 auto linearSolver = std::make_unique<SlamLinearPCGSolver>();
287 blockSolver = std::make_unique<SlamBlockSolver>(std::move(linearSolver));
293 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmGaussNewton(std::move(blockSolver)));
297 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmLevenberg(std::move(blockSolver)));
302 SlamBlockSolver * blockSolver = 0;
307 SlamLinearEigenSolver * linearSolver =
new SlamLinearEigenSolver();
308 linearSolver->setBlockOrdering(
false);
309 blockSolver =
new SlamBlockSolver(linearSolver);
311 #ifdef G2O_HAVE_CHOLMOD
315 SlamLinearCholmodSolver * linearSolver =
new SlamLinearCholmodSolver();
316 linearSolver->setBlockOrdering(
false);
317 blockSolver =
new SlamBlockSolver(linearSolver);
320 #ifdef G2O_HAVE_CSPARSE
324 SlamLinearCSparseSolver* linearSolver =
new SlamLinearCSparseSolver();
325 linearSolver->setBlockOrdering(
false);
326 blockSolver =
new SlamBlockSolver(linearSolver);
332 SlamLinearPCGSolver * linearSolver =
new SlamLinearPCGSolver();
333 blockSolver =
new SlamBlockSolver(linearSolver);
338 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmGaussNewton(blockSolver));
342 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmLevenberg(blockSolver));
346 bool hasGravityConstraints =
false;
349 for(std::multimap<int, Link>::const_iterator
iter=edgeConstraints.begin();
iter!=edgeConstraints.end(); ++
iter)
351 if(
iter->second.from() ==
iter->second.to())
362 hasGravityConstraints =
true;
372 int landmarkVertexOffset = poses.rbegin()->first+1;
373 std::map<int, bool> isLandmarkWithRotation;
375 UDEBUG(
"fill poses to g2o... (rootId=%d hasGravityConstraints=%d isSlam2d=%d)", rootId, hasGravityConstraints?1:0,
isSlam2d()?1:0);
376 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
379 g2o::HyperGraph::Vertex * vertex = 0;
380 int id =
iter->first;
385 g2o::VertexSE2 *
v2 =
new g2o::VertexSE2();
386 v2->setEstimate(g2o::SE2(
iter->second.x(),
iter->second.y(),
iter->second.theta()));
389 UDEBUG(
"Set %d fixed",
id);
397 std::multimap<int, Link>::const_iterator jter=edgeConstraints.find(
id);
398 UASSERT(jter != edgeConstraints.end());
400 if (1 /
static_cast<double>(jter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
402 g2o::VertexPointXY *
v2 =
new g2o::VertexPointXY();
403 v2->setEstimate(Eigen::Vector2d(
iter->second.x(),
iter->second.y()));
405 isLandmarkWithRotation.insert(std::make_pair(
id,
false));
406 id = landmarkVertexOffset -
id;
410 g2o::VertexSE2 *
v2 =
new g2o::VertexSE2();
411 v2->setEstimate(g2o::SE2(
iter->second.x(),
iter->second.y(),
iter->second.theta()));
414 UDEBUG(
"Set %d fixed",
id);
418 isLandmarkWithRotation.insert(std::make_pair(
id,
true));
419 id = landmarkVertexOffset -
id;
431 g2o::VertexSE3 *
v3 =
new g2o::VertexSE3();
437 v3->setEstimate(pose);
438 if(
id == rootId && !hasGravityConstraints)
440 UDEBUG(
"Set %d fixed",
id);
448 std::multimap<int, Link>::const_iterator jter=edgeConstraints.find(
id);
449 UASSERT(jter != edgeConstraints.end());
451 if (1 /
static_cast<double>(jter->second.infMatrix().at<
double>(3,3)) >= 9999.0 ||
452 1 /
static_cast<double>(jter->second.infMatrix().at<
double>(4,4)) >= 9999.0 ||
453 1 /
static_cast<double>(jter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
455 g2o::VertexPointXYZ *
v3 =
new g2o::VertexPointXYZ();
456 v3->setEstimate(Eigen::Vector3d(
iter->second.x(),
iter->second.y(),
iter->second.z()));
458 isLandmarkWithRotation.insert(std::make_pair(
id,
false));
459 id = landmarkVertexOffset -
id;
463 g2o::VertexSE3 *
v3 =
new g2o::VertexSE3();
468 v3->setEstimate(pose);
469 if(
id == rootId && !hasGravityConstraints)
471 UDEBUG(
"Set %d fixed",
id);
475 isLandmarkWithRotation.insert(std::make_pair(
id,
true));
476 id = landmarkVertexOffset -
id;
486 UERROR(
"Could not create vertex for node %d",
id);
496 if(!
isSlam2d() && rootId !=0 && hasGravityConstraints)
498 g2o::VertexSE3*
v1 =
dynamic_cast<g2o::VertexSE3*
>(optimizer.vertex(rootId));
501 g2o::EdgeSE3Prior *
e =
new g2o::EdgeSE3Prior();
507 e->setMeasurement(pose);
511 information(3,3) = information(4,4) = 1;
512 e->setInformation(information);
513 if (!optimizer.addEdge(
e))
516 UERROR(
"Map: Failed adding fixed constraint of rootid %d, set as fixed instead", rootId);
521 UDEBUG(
"Set %d fixed with prior (have gravity constraints)", rootId);
526 UERROR(
"Map: Failed adding fixed constraint of rootid %d (not found in added vertices)", rootId);
530 UDEBUG(
"fill edges to g2o...");
531 #if defined(RTABMAP_VERTIGO)
532 int vertigoVertexId = landmarkVertexOffset - (poses.begin()->first<0?poses.begin()->first-1:0);
534 for(std::multimap<int, Link>::const_iterator
iter=edgeConstraints.begin();
iter!=edgeConstraints.end(); ++
iter)
536 int id1 =
iter->second.from();
537 int id2 =
iter->second.to();
541 g2o::HyperGraph::Edge * edge = 0;
552 id1 = landmarkVertexOffset - id1;
553 id2 = landmarkVertexOffset - id2;
558 if(idTag < 0 && !isLandmarkWithRotation.at(idTag))
561 g2o::VertexPointXY*
v1 = (g2o::VertexPointXY*)optimizer.vertex(id1);
562 priorEdge->setVertex(0,
v1);
567 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
568 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
569 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
570 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
572 priorEdge->setInformation(information);
575 else if (1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
577 g2o::EdgeSE2XYPrior * priorEdge =
new g2o::EdgeSE2XYPrior();
578 g2o::VertexSE2*
v1 = (g2o::VertexSE2*)optimizer.vertex(id1);
579 priorEdge->setVertex(0,
v1);
580 priorEdge->setMeasurement(Eigen::Vector2d(
iter->second.transform().x(),
iter->second.transform().y()));
584 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
585 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
586 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
587 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
589 priorEdge->setInformation(information);
594 g2o::EdgeSE2Prior * priorEdge =
new g2o::EdgeSE2Prior();
595 g2o::VertexSE2*
v1 = (g2o::VertexSE2*)optimizer.vertex(id1);
596 priorEdge->setVertex(0,
v1);
597 priorEdge->setMeasurement(g2o::SE2(
iter->second.transform().x(),
iter->second.transform().y(),
iter->second.transform().theta()));
602 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
603 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
604 information(0,2) =
iter->second.infMatrix().at<
double>(0,5);
605 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
606 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
607 information(1,2) =
iter->second.infMatrix().at<
double>(1,5);
608 information(2,0) =
iter->second.infMatrix().at<
double>(5,0);
609 information(2,1) =
iter->second.infMatrix().at<
double>(5,1);
610 information(2,2) =
iter->second.infMatrix().at<
double>(5,5);
612 priorEdge->setInformation(information);
618 if(idTag < 0 && !isLandmarkWithRotation.at(idTag))
622 g2o::VertexPointXYZ*
v1 = (g2o::VertexPointXYZ*)optimizer.vertex(id1);
623 priorEdge->setVertex(0,
v1);
624 priorEdge->
setMeasurement(Eigen::Vector3d(
iter->second.transform().x(),
iter->second.transform().y(),
iter->second.transform().z()));
629 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
630 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
631 information(0,2) =
iter->second.infMatrix().at<
double>(0,2);
632 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
633 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
634 information(1,2) =
iter->second.infMatrix().at<
double>(1,2);
635 information(2,0) =
iter->second.infMatrix().at<
double>(2,0);
636 information(2,1) =
iter->second.infMatrix().at<
double>(2,1);
637 information(2,2) =
iter->second.infMatrix().at<
double>(2,2);
639 priorEdge->setInformation(information);
642 else if (1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(3,3)) >= 9999.0 ||
643 1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(4,4)) >= 9999.0 ||
644 1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
648 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
649 priorEdge->setVertex(0,
v1);
650 priorEdge->
setMeasurement(Eigen::Vector3d(
iter->second.transform().x(),
iter->second.transform().y(),
iter->second.transform().z()));
655 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
656 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
657 information(0,2) =
iter->second.infMatrix().at<
double>(0,2);
658 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
659 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
660 information(1,2) =
iter->second.infMatrix().at<
double>(1,2);
661 information(2,0) =
iter->second.infMatrix().at<
double>(2,0);
662 information(2,1) =
iter->second.infMatrix().at<
double>(2,1);
663 information(2,2) =
iter->second.infMatrix().at<
double>(2,2);
665 priorEdge->setInformation(information);
671 g2o::EdgeSE3Prior * priorEdge =
new g2o::EdgeSE3Prior();
672 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
673 priorEdge->setVertex(0,
v1);
678 priorEdge->setMeasurement(pose);
683 memcpy(information.
data(),
iter->second.infMatrix().data,
iter->second.infMatrix().total()*
sizeof(
double));
685 priorEdge->setInformation(information);
694 m.head<3>() = Eigen::Vector3d::UnitZ();
702 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
705 priorEdge->setInformation(information);
706 priorEdge->vertices()[0] =
v1;
710 else if(id1<0 || id2 < 0)
715 UASSERT((id1 < 0 && id2 > 0) || (id1 > 0 && id2 < 0));
720 t =
iter->second.transform();
724 t =
iter->second.transform().inverse();
728 id2 = landmarkVertexOffset - id2;
730 #if defined(RTABMAP_VERTIGO)
732 if(this->
isRobust() && isLandmarkWithRotation.at(idTag))
744 v->setId(vertigoVertexId++);
745 UASSERT_MSG(optimizer.addVertex(
v),
uFormat(
"cannot insert switchable vertex %d!?",
v->id()).c_str());
753 prior->setMeasurement(1.0);
757 else if(this->
isRobust() && !isLandmarkWithRotation.at(idTag))
759 UWARN(
"%s cannot be used for landmark constraints without orientation.", Parameters::kOptimizerRobust().
c_str());
765 if(isLandmarkWithRotation.at(idTag))
770 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
771 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
772 information(0,2) =
iter->second.infMatrix().at<
double>(0,5);
773 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
774 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
775 information(1,2) =
iter->second.infMatrix().at<
double>(1,5);
776 information(2,0) =
iter->second.infMatrix().at<
double>(5,0);
777 information(2,1) =
iter->second.infMatrix().at<
double>(5,1);
778 information(2,2) =
iter->second.infMatrix().at<
double>(5,5);
780 #if defined(RTABMAP_VERTIGO)
784 g2o::VertexSE2*
v1 = (g2o::VertexSE2*)optimizer.vertex(id1);
785 g2o::VertexSE2*
v2 = (g2o::VertexSE2*)optimizer.vertex(id2);
791 e->setMeasurement(g2o::SE2(
t.
x(),
t.
y(),
t.theta()));
792 e->setInformation(information);
798 g2o::EdgeSE2 *
e =
new g2o::EdgeSE2();
799 g2o::VertexSE2*
v1 = (g2o::VertexSE2*)optimizer.vertex(id1);
800 g2o::VertexSE2*
v2 = (g2o::VertexSE2*)optimizer.vertex(id2);
805 e->setMeasurement(g2o::SE2(
t.
x(),
t.
y(),
t.theta()));
806 e->setInformation(information);
815 cv::Mat linearCov = cv::Mat(
iter->second.infMatrix(), cv::Range(0,2), cv::Range(0,2)).clone();
816 memcpy(information.
data(), linearCov.data, linearCov.total()*
sizeof(
double));
818 g2o::EdgeSE2PointXY*
e =
new g2o::EdgeSE2PointXY;
819 e->vertices()[0] = optimizer.vertex(id1);
820 e->vertices()[1] = optimizer.vertex(id2);
821 e->setMeasurement(Eigen::Vector2d(
t.
x(),
t.
y()));
822 e->setInformation(information);
829 if(isLandmarkWithRotation.at(idTag))
834 memcpy(information.
data(),
iter->second.infMatrix().data,
iter->second.infMatrix().total()*
sizeof(
double));
839 constraint =
a.linear();
842 #if defined(RTABMAP_VERTIGO)
846 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
847 g2o::VertexSE3*
v2 = (g2o::VertexSE3*)optimizer.vertex(id2);
853 e->setMeasurement(constraint);
854 e->setInformation(information);
860 g2o::EdgeSE3 *
e =
new g2o::EdgeSE3();
861 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
862 g2o::VertexSE3*
v2 = (g2o::VertexSE3*)optimizer.vertex(id2);
867 e->setMeasurement(constraint);
868 e->setInformation(information);
877 cv::Mat linearCov = cv::Mat(
iter->second.infMatrix(), cv::Range(0,3), cv::Range(0,3)).clone();
878 memcpy(information.
data(), linearCov.data, linearCov.total()*
sizeof(
double));
881 g2o::EdgeSE3PointXYZ*
e =
new g2o::EdgeSE3PointXYZ;
882 e->vertices()[0] = optimizer.vertex(id1);
883 e->vertices()[1] = optimizer.vertex(id2);
884 e->setMeasurement(Eigen::Vector3d(
t.
x(),
t.
y(),
t.z()));
885 e->setInformation(information);
894 #if defined(RTABMAP_VERTIGO)
909 v->setId(vertigoVertexId++);
910 UASSERT_MSG(optimizer.addVertex(
v),
uFormat(
"cannot insert switchable vertex %d!?",
v->id()).c_str());
918 prior->setMeasurement(1.0);
929 information(0,0) =
iter->second.infMatrix().at<
double>(0,0);
930 information(0,1) =
iter->second.infMatrix().at<
double>(0,1);
931 information(0,2) =
iter->second.infMatrix().at<
double>(0,5);
932 information(1,0) =
iter->second.infMatrix().at<
double>(1,0);
933 information(1,1) =
iter->second.infMatrix().at<
double>(1,1);
934 information(1,2) =
iter->second.infMatrix().at<
double>(1,5);
935 information(2,0) =
iter->second.infMatrix().at<
double>(5,0);
936 information(2,1) =
iter->second.infMatrix().at<
double>(5,1);
937 information(2,2) =
iter->second.infMatrix().at<
double>(5,5);
940 #if defined(RTABMAP_VERTIGO)
946 g2o::VertexSE2*
v1 = (g2o::VertexSE2*)optimizer.vertex(id1);
947 g2o::VertexSE2*
v2 = (g2o::VertexSE2*)optimizer.vertex(id2);
953 e->setMeasurement(g2o::SE2(
iter->second.transform().x(),
iter->second.transform().y(),
iter->second.transform().theta()));
954 e->setInformation(information);
960 g2o::EdgeSE2 *
e =
new g2o::EdgeSE2();
961 g2o::VertexSE2*
v1 = (g2o::VertexSE2*)optimizer.vertex(id1);
962 g2o::VertexSE2*
v2 = (g2o::VertexSE2*)optimizer.vertex(id2);
967 e->setMeasurement(g2o::SE2(
iter->second.transform().x(),
iter->second.transform().y(),
iter->second.transform().theta()));
968 e->setInformation(information);
977 memcpy(information.
data(),
iter->second.infMatrix().data,
iter->second.infMatrix().total()*
sizeof(
double));
982 constraint =
a.linear();
985 #if defined(RTABMAP_VERTIGO)
991 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
992 g2o::VertexSE3*
v2 = (g2o::VertexSE3*)optimizer.vertex(id2);
998 e->setMeasurement(constraint);
999 e->setInformation(information);
1005 g2o::EdgeSE3 *
e =
new g2o::EdgeSE3();
1006 g2o::VertexSE3*
v1 = (g2o::VertexSE3*)optimizer.vertex(id1);
1007 g2o::VertexSE3*
v2 = (g2o::VertexSE3*)optimizer.vertex(id2);
1010 e->setVertex(0,
v1);
1011 e->setVertex(1,
v2);
1012 e->setMeasurement(constraint);
1013 e->setInformation(information);
1019 if (edge && !optimizer.addEdge(edge))
1022 UERROR(
"Map: Failed adding constraint between %d and %d, skipping", id1, id2);
1028 UDEBUG(
"Initial optimization...");
1029 optimizer.initializeOptimization();
1031 if(!optimizer.verifyInformationMatrices(
true))
1033 UERROR(
"This error can be caused by (1) bad covariance matrix "
1034 "set in odometry messages "
1035 "(see requirements in g2o::OptimizableGraph::verifyInformationMatrices() function) "
1036 "or that (2) PCL and g2o hadn't "
1037 "been built both with or without \"-march=native\" compilation "
1038 "flag (if one library is built with this flag and not the other, "
1039 "this is causing Eigen to not work properly, resulting in segmentation faults).");
1040 return optimizedPoses;
1046 double lastError = 0.0;
1048 if (!optimizer.solver()->init()) {
1049 UERROR(
"g2o: Error while initializing solver");
1050 return optimizedPoses;
1053 if(intermediateGraphes || this->
epsilon() > 0.0)
1057 if(intermediateGraphes)
1061 std::map<int, Transform> tmpPoses;
1064 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
1066 int id =
iter->first;
1069 const g2o::VertexSE2*
v = (
const g2o::VertexSE2*)optimizer.vertex(
id);
1074 Transform t(
v->estimate().translation()[0],
v->estimate().translation()[1],
iter->second.z(),
roll,
pitch,
v->estimate().rotation().angle());
1075 tmpPoses.insert(std::pair<int, Transform>(
id,
t));
1080 UERROR(
"Vertex %d not found!?",
id);
1085 const g2o::OptimizableGraph::Vertex*
v = (
const g2o::OptimizableGraph::Vertex*)optimizer.vertex(landmarkVertexOffset -
id);
1088 if(isLandmarkWithRotation.at(
id))
1090 const g2o::VertexSE2* vSE2 = (
const g2o::VertexSE2*)
v;
1093 Transform t(vSE2->estimate().translation()[0], vSE2->estimate().translation()[1],
iter->second.z(),
roll,
pitch, vSE2->estimate().rotation().angle());
1094 tmpPoses.insert(std::pair<int, Transform>(
id,
t));
1099 const g2o::VertexPointXY* vP = (
const g2o::VertexPointXY*)
v;
1103 tmpPoses.insert(std::pair<int, Transform>(
id,
t));
1109 UERROR(
"Vertex %d not found!?",
id);
1116 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
1118 int id =
iter->first;
1121 const g2o::VertexSE3*
v = (
const g2o::VertexSE3*)optimizer.vertex(
id);
1125 tmpPoses.insert(std::pair<int, Transform>(
id,
t));
1130 UERROR(
"Vertex %d not found!?",
id);
1135 const g2o::OptimizableGraph::Vertex*
v = (
const g2o::OptimizableGraph::Vertex*)optimizer.vertex(landmarkVertexOffset -
id);
1138 if(isLandmarkWithRotation.at(
id))
1140 const g2o::VertexSE3* vSE3 = (
const g2o::VertexSE3*)
v;
1142 tmpPoses.insert(std::pair<int, Transform>(
id,
t));
1147 const g2o::VertexPointXYZ* vP = (
const g2o::VertexPointXYZ*)
v;
1151 tmpPoses.insert(std::pair<int, Transform>(
id,
t));
1157 UERROR(
"Vertex %d not found!?",
id);
1162 intermediateGraphes->push_back(tmpPoses);
1166 g2o::OptimizationAlgorithm::SolverResult
result = optimizer.solver()->solve(
i);
1170 optimizer.computeActiveErrors();
1171 double chi2 = optimizer.activeRobustChi2();
1172 UDEBUG(
"iteration %d: %d nodes, %d edges, chi2: %f",
i, (
int)optimizer.vertices().size(), (
int)optimizer.edges().size(), chi2);
1174 if(
i>0 && optimizer.activeRobustChi2() > 1000000000000.0)
1176 UERROR(
"g2o: Large optimimzation error detected (%f), aborting optimization!");
1177 return optimizedPoses;
1180 if(
result == g2o::OptimizationAlgorithm::Fail)
1182 UERROR(
"g2o: Solver failed, aborting optimization!");
1183 return optimizedPoses;
1186 double errorDelta = lastError - chi2;
1187 if(
i>0 && errorDelta < this->
epsilon())
1191 UDEBUG(
"Negative improvement?! Ignore and continue optimizing... (%f < %f)", errorDelta, this->
epsilon());
1195 UINFO(
"Stop optimizing, not enough improvement (%f < %f)", errorDelta, this->
epsilon());
1199 else if(
i==0 && chi2 < this->
epsilon())
1201 UINFO(
"Stop optimizing, error is already under epsilon (%f < %f)", chi2, this->
epsilon());
1210 optimizer.computeActiveErrors();
1211 UDEBUG(
"%d nodes, %d edges, chi2: %f", (
int)optimizer.vertices().size(), (
int)optimizer.edges().size(), optimizer.activeRobustChi2());
1215 *finalError = lastError;
1219 *iterationsDone = it;
1221 UINFO(
"g2o optimizing end (%d iterations done, error=%f, time = %f s)", it, optimizer.activeRobustChi2(),
timer.ticks());
1223 if(optimizer.activeRobustChi2() > 1000000000000.0)
1225 UERROR(
"g2o: Large optimimzation error detected (%f), aborting optimization!");
1226 return optimizedPoses;
1231 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
1233 int id =
iter->first;
1236 const g2o::VertexSE2*
v = (
const g2o::VertexSE2*)optimizer.vertex(
id);
1241 Transform t(
v->estimate().translation()[0],
v->estimate().translation()[1],
iter->second.z(),
roll,
pitch,
v->estimate().rotation().angle());
1242 optimizedPoses.insert(std::pair<int, Transform>(
id,
t));
1247 UERROR(
"Vertex %d not found!?",
id);
1252 const g2o::OptimizableGraph::Vertex*
v = (
const g2o::OptimizableGraph::Vertex*)optimizer.vertex(landmarkVertexOffset -
id);
1256 if(isLandmarkWithRotation.at(
id))
1258 const g2o::VertexSE2* vSE2 = (
const g2o::VertexSE2*)
v;
1261 Transform t(vSE2->estimate().translation()[0], vSE2->estimate().translation()[1],
iter->second.z(),
roll,
pitch, vSE2->estimate().rotation().angle());
1262 optimizedPoses.insert(std::pair<int, Transform>(
id,
t));
1267 const g2o::VertexPointXY* vP = (
const g2o::VertexPointXY*)
v;
1271 optimizedPoses.insert(std::pair<int, Transform>(
id,
t));
1277 UERROR(
"Vertex %d not found!?",
id);
1282 g2o::VertexSE2*
v = (g2o::VertexSE2*)optimizer.vertex(poses.rbegin()->first);
1286 g2o::SparseBlockMatrix<g2o::MatrixXD> spinv;
1287 optimizer.computeMarginals(spinv,
v);
1288 UINFO(
"Computed marginals = %fs (cols=%d rows=%d, v=%d id=%d)",
t.ticks(), spinv.cols(), spinv.rows(),
v->hessianIndex(), poses.rbegin()->first);
1289 if(
v->hessianIndex() >= 0 &&
v->hessianIndex() < (
int)spinv.blockCols().size())
1291 g2o::SparseBlockMatrix<g2o::MatrixXD>::SparseMatrixBlock *
block = spinv.blockCols()[
v->hessianIndex()].begin()->second;
1293 outputCovariance.at<
double>(0,0) = (*
block)(0,0);
1294 outputCovariance.at<
double>(0,1) = (*
block)(0,1);
1295 outputCovariance.at<
double>(0,5) = (*
block)(0,2);
1296 outputCovariance.at<
double>(1,0) = (*
block)(1,0);
1297 outputCovariance.at<
double>(1,1) = (*
block)(1,1);
1298 outputCovariance.at<
double>(1,5) = (*
block)(1,2);
1299 outputCovariance.at<
double>(5,0) = (*
block)(2,0);
1300 outputCovariance.at<
double>(5,1) = (*
block)(2,1);
1301 outputCovariance.at<
double>(5,5) = (*
block)(2,2);
1303 else if(
v->hessianIndex() < 0)
1305 UWARN(
"Computing marginals: vertex %d has negative hessian index (%d). Cannot compute last pose covariance.", poses.rbegin()->first,
v->hessianIndex());
1309 UWARN(
"Computing marginals: vertex %d has hessian not valid (%d > block size=%d). Cannot compute last pose covariance.", poses.rbegin()->first,
v->hessianIndex(), (
int)spinv.blockCols().size());
1314 UERROR(
"Vertex %d not found!? Cannot compute marginals...", poses.rbegin()->first);
1319 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
1321 int id =
iter->first;
1324 const g2o::VertexSE3*
v = (
const g2o::VertexSE3*)optimizer.vertex(
id);
1328 optimizedPoses.insert(std::pair<int, Transform>(
id,
t));
1333 UERROR(
"Vertex %d not found!?",
id);
1338 const g2o::OptimizableGraph::Vertex*
v = (
const g2o::OptimizableGraph::Vertex*)optimizer.vertex(landmarkVertexOffset -
id);
1342 if(isLandmarkWithRotation.at(
id))
1344 const g2o::VertexSE3* vSE3 = (
const g2o::VertexSE3*)
v;
1346 optimizedPoses.insert(std::pair<int, Transform>(
id,
t));
1351 const g2o::VertexPointXYZ* vP = (
const g2o::VertexPointXYZ*)
v;
1355 optimizedPoses.insert(std::pair<int, Transform>(
id,
t));
1361 UERROR(
"Vertex %d not found!?",
id);
1366 g2o::VertexSE3*
v = (g2o::VertexSE3*)optimizer.vertex(poses.rbegin()->first);
1370 g2o::SparseBlockMatrix<g2o::MatrixXD> spinv;
1371 optimizer.computeMarginals(spinv,
v);
1372 UINFO(
"Computed marginals = %fs (cols=%d rows=%d, v=%d id=%d)",
t.ticks(), spinv.cols(), spinv.rows(),
v->hessianIndex(), poses.rbegin()->first);
1373 if(
v->hessianIndex() >= 0 &&
v->hessianIndex() < (
int)spinv.blockCols().size())
1375 g2o::SparseBlockMatrix<g2o::MatrixXD>::SparseMatrixBlock *
block = spinv.blockCols()[
v->hessianIndex()].begin()->second;
1377 memcpy(outputCovariance.data,
block->data(), outputCovariance.total()*
sizeof(
double));
1379 else if(
v->hessianIndex() < 0)
1381 UWARN(
"Computing marginals: vertex %d has negative hessian index (%d). Cannot compute last pose covariance.", poses.rbegin()->first,
v->hessianIndex());
1383 #ifdef RTABMAP_G2O_CPP11
1386 UWARN(
"Computing marginals: vertex %d has hessian not valid (%d > block size=%d). Cannot compute last pose covariance.", poses.rbegin()->first,
v->hessianIndex(), (
int)spinv.blockCols().size());
1392 UERROR(
"Vertex %d not found!? Cannot compute marginals...", poses.rbegin()->first);
1396 else if(poses.size() == 1 ||
iterations() <= 0)
1398 optimizedPoses = poses;
1402 UWARN(
"This method should be called at least with 1 pose!");
1404 UDEBUG(
"Optimizing graph...end!");
1406 #ifdef RTABMAP_ORB_SLAM
1407 UERROR(
"G2O graph optimization cannot be used with g2o built from ORB_SLAM, only SBA is available.");
1409 UERROR(
"Not built with G2O support!");
1412 return optimizedPoses;
1415 #ifdef RTABMAP_ORB_SLAM
1419 class EdgeSE3Expmap :
public g2o::BaseBinaryEdge<6, g2o::SE3Quat, g2o::VertexSE3Expmap, g2o::VertexSE3Expmap>
1423 EdgeSE3Expmap(): BaseBinaryEdge<6,
g2o::SE3Quat,
g2o::VertexSE3Expmap,
g2o::VertexSE3Expmap>(){}
1424 bool read(std::istream& is)
1429 bool write(std::ostream& os)
const
1436 const g2o::VertexSE3Expmap*
v1 =
dynamic_cast<const g2o::VertexSE3Expmap*
>(_vertices[0]);
1437 const g2o::VertexSE3Expmap*
v2 =
dynamic_cast<const g2o::VertexSE3Expmap*
>(_vertices[1]);
1438 g2o::SE3Quat
delta = _inverseMeasurement * (
v1->estimate().inverse()*
v2->estimate());
1439 _error[0]=
delta.translation().x();
1440 _error[1]=
delta.translation().y();
1441 _error[2]=
delta.translation().z();
1442 _error[3]=
delta.rotation().x();
1443 _error[4]=
delta.rotation().y();
1444 _error[5]=
delta.rotation().z();
1447 virtual void setMeasurement(
const g2o::SE3Quat& meas){
1449 _inverseMeasurement=meas.inverse();
1452 virtual double initialEstimatePossible(
const g2o::OptimizableGraph::VertexSet& , g2o::OptimizableGraph::Vertex* ) {
return 1.;}
1453 virtual void initialEstimate(
const g2o::OptimizableGraph::VertexSet& from_, g2o::OptimizableGraph::Vertex* ){
1454 g2o::VertexSE3Expmap* from =
static_cast<g2o::VertexSE3Expmap*
>(_vertices[0]);
1455 g2o::VertexSE3Expmap* to =
static_cast<g2o::VertexSE3Expmap*
>(_vertices[1]);
1456 if (from_.count(from) > 0)
1457 to->setEstimate((g2o::SE3Quat) from->estimate() * _measurement);
1459 from->setEstimate((g2o::SE3Quat) to->estimate() * _inverseMeasurement);
1462 virtual bool setMeasurementData(
const double* d){
1464 _measurement.fromVector(v);
1465 _inverseMeasurement = _measurement.inverse();
1469 virtual bool getMeasurementData(
double* d)
const{
1471 v = _measurement.toVector();
1475 virtual int measurementDimension()
const {
return 7;}
1477 virtual bool setMeasurementFromState() {
1478 const g2o::VertexSE3Expmap*
v1 =
dynamic_cast<const g2o::VertexSE3Expmap*
>(_vertices[0]);
1479 const g2o::VertexSE3Expmap*
v2 =
dynamic_cast<const g2o::VertexSE3Expmap*
>(_vertices[1]);
1480 _measurement = (
v1->estimate().inverse()*
v2->estimate());
1481 _inverseMeasurement = _measurement.inverse();
1486 g2o::SE3Quat _inverseMeasurement;
1492 const std::map<int, Transform> & poses,
1493 const std::multimap<int, Link> & links,
1494 const std::map<
int, std::vector<CameraModel> > & models,
1495 std::map<int, cv::Point3f> & points3DMap,
1496 const std::map<
int, std::map<int, FeatureBA> > & wordReferences,
1497 std::set<int> * outliers)
1499 std::map<int, Transform> optimizedPoses;
1500 #if defined(RTABMAP_G2O) || defined(RTABMAP_ORB_SLAM)
1501 UDEBUG(
"Optimizing graph...");
1503 optimizedPoses.clear();
1504 if(poses.size()>=2 &&
iterations() > 0 && (models.size() == poses.size() || poses.begin()->first < 0))
1506 g2o::SparseOptimizer optimizer;
1508 #if defined(RTABMAP_G2O_CPP11) && !defined(RTABMAP_ORB_SLAM)
1509 std::unique_ptr<g2o::BlockSolver_6_3::LinearSolverType> linearSolver;
1511 g2o::BlockSolver_6_3::LinearSolverType * linearSolver = 0;
1514 #ifdef RTABMAP_ORB_SLAM
1515 linearSolver =
new g2o::LinearSolverEigen<g2o::BlockSolver_6_3::PoseMatrixType>();
1520 #ifdef RTABMAP_G2O_CPP11
1521 linearSolver = std::make_unique<g2o::LinearSolverEigen<g2o::BlockSolver_6_3::PoseMatrixType> >();
1523 linearSolver =
new g2o::LinearSolverEigen<g2o::BlockSolver_6_3::PoseMatrixType>();
1526 #ifdef G2O_HAVE_CHOLMOD
1530 #ifdef RTABMAP_G2O_CPP11
1531 linearSolver = std::make_unique<g2o::LinearSolverCholmod<g2o::BlockSolver_6_3::PoseMatrixType> >();
1533 linearSolver =
new g2o::LinearSolverCholmod<g2o::BlockSolver_6_3::PoseMatrixType>();
1537 #ifdef G2O_HAVE_CSPARSE
1541 #ifdef RTABMAP_G2O_CPP11
1542 linearSolver = std::make_unique<g2o::LinearSolverCSparse<g2o::BlockSolver_6_3::PoseMatrixType> >();
1544 linearSolver =
new g2o::LinearSolverCSparse<g2o::BlockSolver_6_3::PoseMatrixType>();
1551 #ifdef RTABMAP_G2O_CPP11
1552 linearSolver = std::make_unique<g2o::LinearSolverPCG<g2o::BlockSolver_6_3::PoseMatrixType> >();
1554 linearSolver =
new g2o::LinearSolverPCG<g2o::BlockSolver_6_3::PoseMatrixType>();
1557 #endif // RTABMAP_ORB_SLAM
1559 #ifndef RTABMAP_ORB_SLAM
1562 #ifdef RTABMAP_G2O_CPP11
1563 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmGaussNewton(
1564 std::make_unique<g2o::BlockSolver_6_3>(std::move(linearSolver))));
1566 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmGaussNewton(
new g2o::BlockSolver_6_3(linearSolver)));
1572 #if defined(RTABMAP_G2O_CPP11) && !defined(RTABMAP_ORB_SLAM)
1573 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmLevenberg(
1574 std::make_unique<g2o::BlockSolver_6_3>(std::move(linearSolver))));
1576 optimizer.setAlgorithm(
new g2o::OptimizationAlgorithmLevenberg(
new g2o::BlockSolver_6_3(linearSolver)));
1581 UDEBUG(
"fill poses to g2o...");
1582 for(std::map<int, Transform>::const_iterator
iter=poses.begin();
iter!=poses.end(); ++
iter)
1587 std::map<int, std::vector<CameraModel> >::const_iterator iterModel = models.find(
iter->first);
1588 UASSERT(iterModel != models.end() && !iterModel->second.empty());
1589 for(
size_t i=0;
i<iterModel->second.size(); ++
i)
1591 UASSERT(iterModel->second[
i].isValidForProjection());
1593 Transform camPose =
iter->second * iterModel->second[
i].localTransform();
1597 #ifdef RTABMAP_ORB_SLAM
1598 g2o::VertexSE3Expmap * vCam =
new g2o::VertexSE3Expmap();
1600 g2o::VertexCam * vCam =
new g2o::VertexCam();
1604 #ifdef RTABMAP_ORB_SLAM
1606 vCam->setEstimate(g2o::SE3Quat(
a.linear(),
a.translation()));
1610 iterModel->second[
i].fx(),
1611 iterModel->second[
i].fy(),
1612 iterModel->second[
i].cx(),
1613 iterModel->second[
i].cy(),
1614 iterModel->second[
i].Tx()<0.0?-iterModel->second[
i].Tx()/iterModel->second[
i].fx():
baseline_);
1615 vCam->setEstimate(
cam);
1620 vCam->setFixed((rootId >= 0 &&
iter->first == rootId) || (rootId < 0 && iter->
first != -rootId));
1634 UASSERT_MSG(optimizer.addVertex(vCam),
uFormat(
"cannot insert cam vertex %d (pose=%d)!?", vCam->id(),
iter->first).c_str());
1639 UDEBUG(
"fill edges to g2o...");
1640 for(std::multimap<int, Link>::const_iterator
iter=links.begin();
iter!=links.end(); ++
iter)
1642 if(
iter->second.from() > 0 &&
1643 iter->second.to() > 0 &&
1648 int id1 =
iter->second.from();
1649 int id2 =
iter->second.to();
1653 #ifndef RTABMAP_ORB_SLAM
1654 g2o::HyperGraph::Edge * edge = 0;
1659 m.head<3>() = Eigen::Vector3d::UnitZ();
1669 std::map<int, std::vector<CameraModel> >::const_iterator iterModel = models.find(
iter->first);
1671 UASSERT(iterModel != models.end() && !iterModel->second.empty() && !iterModel->second[0].localTransform().isNull());
1674 priorEdge->setInformation(information);
1675 priorEdge->vertices()[0] =
v1;
1678 if (edge && !optimizer.addEdge(edge))
1681 UERROR(
"Map: Failed adding constraint between %d and %d, skipping", id1, id2);
1682 return optimizedPoses;
1686 else if(id1>0 && id2>0)
1693 memcpy(information.
data(),
iter->second.infMatrix().data,
iter->second.infMatrix().total()*
sizeof(
double));
1697 Transform camLink = models.at(id1)[0].localTransform().
inverse()*
iter->second.transform()*models.at(id2)[0].localTransform();
1704 #ifdef RTABMAP_ORB_SLAM
1705 EdgeSE3Expmap *
e =
new EdgeSE3Expmap();
1706 g2o::VertexSE3Expmap*
v1 = (g2o::VertexSE3Expmap*)optimizer.vertex(id1*
MULTICAM_OFFSET);
1707 g2o::VertexSE3Expmap*
v2 = (g2o::VertexSE3Expmap*)optimizer.vertex(id2*
MULTICAM_OFFSET);
1712 camLink = camPose1 * camPose1 * camLink * camPose2Inv * camPose2Inv;
1714 g2o::EdgeSBACam *
e =
new g2o::EdgeSBACam();
1720 e->setVertex(0,
v1);
1721 e->setVertex(1,
v2);
1723 e->setMeasurement(g2o::SE3Quat(
a.linear(),
a.translation()));
1724 e->setInformation(information);
1726 if (!optimizer.addEdge(
e))
1729 UERROR(
"Map: Failed adding constraint between %d and %d, skipping", id1, id2);
1730 return optimizedPoses;
1736 UDEBUG(
"fill hard edges between camera 0 and other cameras (multicam)...");
1737 for(std::map<
int, std::vector<CameraModel> >::const_iterator
iter=models.begin();
iter!=models.end(); ++
iter)
1739 int id =
iter->first;
1742 for(
size_t i=1;
i<
iter->second.size(); ++
i)
1749 Transform camLink =
iter->second[0].localTransform().inverse()*
iter->second[
i].localTransform();
1750 #ifdef RTABMAP_ORB_SLAM
1751 EdgeSE3Expmap *
e =
new EdgeSE3Expmap();
1752 g2o::VertexSE3Expmap*
v1 = (g2o::VertexSE3Expmap*)optimizer.vertex(
id*
MULTICAM_OFFSET);
1753 g2o::VertexSE3Expmap*
v2 = (g2o::VertexSE3Expmap*)optimizer.vertex(
id*
MULTICAM_OFFSET+
i);
1758 camLink = camPose1 * camPose1 * camLink * camPose2Inv * camPose2Inv;
1760 g2o::EdgeSBACam *
e =
new g2o::EdgeSBACam();
1772 e->setVertex(0,
v1);
1773 e->setVertex(1,
v2);
1775 e->setMeasurement(g2o::SE3Quat(
a.linear(),
a.translation()));
1776 e->setInformation(information);
1778 if (!optimizer.addEdge(
e))
1781 UERROR(
"Map: Failed adding constraint between %d and %d, skipping",
v1->id(),
v2->id());
1782 return optimizedPoses;
1788 UDEBUG(
"fill 3D points to g2o...");
1790 int negVertexOffset = stepVertexId;
1791 if(wordReferences.size() && wordReferences.rbegin()->first>0)
1793 negVertexOffset += wordReferences.rbegin()->first;
1795 UDEBUG(
"stepVertexId=%d, negVertexOffset=%d", stepVertexId, negVertexOffset);
1796 std::list<g2o::OptimizableGraph::Edge*>
edges;
1797 for(std::map<
int, std::map<int, FeatureBA> >::const_iterator
iter = wordReferences.begin();
iter!=wordReferences.end(); ++
iter)
1799 int id =
iter->first;
1800 if(points3DMap.find(
id) != points3DMap.end())
1802 cv::Point3f pt3d = points3DMap.at(
id);
1805 UWARN(
"Ignoring 3D point %d because it has nan value(s)!",
id);
1808 g2o::VertexSBAPointXYZ* vpt3d =
new g2o::VertexSBAPointXYZ();
1810 vpt3d->setEstimate(Eigen::Vector3d(pt3d.x, pt3d.y, pt3d.z));
1813 vpt3d->setId(negVertexOffset +
id*-1);
1817 vpt3d->setId(stepVertexId +
id);
1820 vpt3d->setMarginalized(
true);
1821 optimizer.addVertex(vpt3d);
1826 for(std::map<int, FeatureBA>::const_iterator jter=
iter->second.begin(); jter!=
iter->second.end(); ++jter)
1828 int poseId = jter->first;
1829 int camIndex = jter->second.cameraIndex;
1831 if(poses.find(poseId) != poses.end() && optimizer.vertex(camId) != 0)
1834 double depth = pt.
depth;
1838 g2o::OptimizableGraph::Edge *
e;
1840 #ifdef RTABMAP_ORB_SLAM
1841 g2o::VertexSE3Expmap* vcam =
dynamic_cast<g2o::VertexSE3Expmap*
>(optimizer.vertex(camId));
1842 std::map<int, std::vector<CameraModel> >::const_iterator iterModel = models.find(poseId);
1844 UASSERT(iterModel != models.end() && camIndex<iterModel->
second.size() && iterModel->second[camIndex].isValidForProjection());
1845 baseline = iterModel->second[camIndex].Tx()<0.0?-iterModel->second[camIndex].Tx()/iterModel->second[camIndex].fx():
baseline_;
1847 g2o::VertexCam* vcam =
dynamic_cast<g2o::VertexCam*
>(optimizer.vertex(camId));
1848 baseline = vcam->estimate().baseline;
1854 #ifdef RTABMAP_ORB_SLAM
1855 g2o::EdgeStereoSE3ProjectXYZ*
es =
new g2o::EdgeStereoSE3ProjectXYZ();
1856 float disparity =
baseline * iterModel->second[camIndex].fx() / depth;
1857 Eigen::Vector3d obs( pt.
kpt.pt.x, pt.
kpt.pt.y, pt.
kpt.pt.x-disparity);
1858 es->setMeasurement(obs);
1860 es->setInformation(Eigen::Matrix3d::Identity() / variance);
1861 es->fx = iterModel->second[camIndex].fx();
1862 es->fy = iterModel->second[camIndex].fy();
1863 es->cx = iterModel->second[camIndex].cx();
1864 es->cy = iterModel->second[camIndex].cy();
1868 g2o::EdgeProjectP2SC*
es =
new g2o::EdgeProjectP2SC();
1869 float disparity =
baseline * vcam->estimate().Kcam(0,0) / depth;
1870 Eigen::Vector3d obs( pt.
kpt.pt.x, pt.
kpt.pt.y, pt.
kpt.pt.x-disparity);
1871 es->setMeasurement(obs);
1873 es->setInformation(Eigen::Matrix3d::Identity() / variance);
1881 UDEBUG(
"Stereo camera model detected but current "
1882 "observation (pt=%d to cam=%d, kpt=[%d,%d]) has null depth (%f m), adding "
1883 "mono observation instead.",
1884 vpt3d->id()-stepVertexId, camId, (
int)pt.
kpt.pt.x, (
int)pt.
kpt.pt.y, depth);
1887 #ifdef RTABMAP_ORB_SLAM
1888 g2o::EdgeSE3ProjectXYZ* em =
new g2o::EdgeSE3ProjectXYZ();
1889 Eigen::Vector2d obs( pt.
kpt.pt.x, pt.
kpt.pt.y);
1890 em->setMeasurement(obs);
1891 em->setInformation(Eigen::Matrix2d::Identity() / variance);
1892 em->fx = iterModel->second[camIndex].fx();
1893 em->fy = iterModel->second[camIndex].fy();
1894 em->cx = iterModel->second[camIndex].cx();
1895 em->cy = iterModel->second[camIndex].cy();
1899 g2o::EdgeProjectP2MC* em =
new g2o::EdgeProjectP2MC();
1900 Eigen::Vector2d obs( pt.
kpt.pt.x, pt.
kpt.pt.y);
1901 em->setMeasurement(obs);
1902 em->setInformation(Eigen::Matrix2d::Identity() / variance);
1906 e->setVertex(0, vpt3d);
1907 e->setVertex(1, vcam);
1911 g2o::RobustKernelHuber* kernel =
new g2o::RobustKernelHuber;
1913 e->setRobustKernel(kernel);
1916 optimizer.addEdge(
e);
1923 UDEBUG(
"Initial optimization...");
1924 optimizer.initializeOptimization();
1926 UASSERT(optimizer.verifyInformationMatrices());
1932 int outliersCount = 0;
1933 int outliersCountFar = 0;
1940 optimizer.computeActiveErrors();
1941 double chi2 = optimizer.activeRobustChi2();
1945 UERROR(
"Optimization generated NANs, aborting optimization! Try another g2o's optimizer (current=%d).",
optimizer_);
1946 return optimizedPoses;
1948 UDEBUG(
"iteration %d: %d nodes, %d edges, chi2: %f",
i, (
int)optimizer.vertices().size(), (
int)optimizer.edges().size(), chi2);
1950 if(
i>0 && (optimizer.activeRobustChi2() > 1000000000000.0 || !
uIsFinite(optimizer.activeRobustChi2())))
1952 UWARN(
"g2o: Large optimization error detected (%f), aborting optimization!");
1953 return optimizedPoses;
1960 if((*iter)->level() == 0 && (*iter)->chi2() > (*iter)->robustKernel()->delta())
1962 (*iter)->setLevel(1);
1965 #ifdef RTABMAP_ORB_SLAM
1966 if(
dynamic_cast<g2o::EdgeStereoSE3ProjectXYZ*
>(*
iter) != 0)
1968 d = ((g2o::EdgeStereoSE3ProjectXYZ*)(*
iter))->measurement()[0]-((g2o::EdgeStereoSE3ProjectXYZ*)(*
iter))->measurement()[2];
1972 if(
dynamic_cast<g2o::EdgeProjectP2SC*
>(*
iter) != 0)
1974 d = ((g2o::EdgeProjectP2SC*)(*
iter))->measurement()[0]-((g2o::EdgeProjectP2SC*)(*
iter))->measurement()[2];
1980 if((*iter)->vertex(0)->id() > negVertexOffset)
1982 pt3d = points3DMap.at(negVertexOffset - (*iter)->vertex(0)->id());
1986 pt3d = points3DMap.at((*iter)->vertex(0)->id()-stepVertexId);
1988 ((g2o::VertexSBAPointXYZ*)(*iter)->vertex(0))->setEstimate(Eigen::Vector3d(pt3d.x, pt3d.y, pt3d.z));
1992 outliers->insert((*iter)->vertex(0)->id()-stepVertexId);
2002 optimizer.initializeOptimization(0);
2003 UDEBUG(
"outliers=%d outliersCountFar=%d", outliersCount, outliersCountFar);
2006 UDEBUG(
"g2o optimizing end (%d iterations done, error=%f, outliers=%d/%d (delta=%f) time = %f s)", it, optimizer.activeRobustChi2(), outliersCount, (
int)
edges.size(),
robustKernelDelta_,
timer.ticks());
2008 if(optimizer.activeRobustChi2() > 1000000000000.0)
2010 UWARN(
"g2o: Large optimization error detected (%f), aborting optimization!");
2011 return optimizedPoses;
2015 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
2020 #ifdef RTABMAP_ORB_SLAM
2021 const g2o::VertexSE3Expmap*
v = (
const g2o::VertexSE3Expmap*)optimizer.vertex(camId);
2023 const g2o::VertexCam*
v = (
const g2o::VertexCam*)optimizer.vertex(camId);
2029 #ifdef RTABMAP_ORB_SLAM
2034 t *= models.at(
iter->first)[0].localTransform().inverse();
2039 UERROR(
"Optimized pose %d is null!?!?",
iter->first);
2040 optimizedPoses.clear();
2041 return optimizedPoses;
2048 t =
iter->second.inverse() *
t;
2049 optimizedPoses.insert(std::pair<int, Transform>(
iter->first,
iter->second *
t.to3DoF()));
2053 optimizedPoses.insert(std::pair<int, Transform>(
iter->first,
t));
2058 UERROR(
"Vertex (pose) %d (cam=%d) not found!?",
iter->first, camId);
2065 for(std::map<int, cv::Point3f>::iterator
iter = points3DMap.begin();
iter!=points3DMap.end(); ++
iter)
2067 const g2o::VertexSBAPointXYZ*
v;
2068 int id =
iter->first;
2071 v = (
const g2o::VertexSBAPointXYZ*)optimizer.vertex(negVertexOffset +
id*-1);
2075 v = (
const g2o::VertexSBAPointXYZ*)optimizer.vertex(stepVertexId +
id);
2080 cv::Point3f
p(
v->estimate()[0],
v->estimate()[1],
v->estimate()[2]);
2086 iter->second.x =
iter->second.y =
iter->second.z = std::numeric_limits<float>::quiet_NaN();
2090 else if(poses.size() > 1 && (poses.size() != models.size() && poses.begin()->first > 0))
2092 UERROR(
"This method should be called with size of poses = size camera models!");
2094 else if(poses.size() == 1 ||
iterations() <= 0)
2096 optimizedPoses = poses;
2100 UWARN(
"This method should be called at least with 1 pose!");
2102 UDEBUG(
"Optimizing graph...end!");
2104 UERROR(
"Not built with G2O support!");
2106 return optimizedPoses;
2110 const std::string & fileName,
2111 const std::map<int, Transform> & poses,
2112 const std::multimap<int, Link> & edgeConstraints)
2117 fopen_s(&
file, fileName.c_str(),
"w");
2119 file = fopen(fileName.c_str(),
"w");
2125 setlocale(LC_ALL,
"en_US.UTF-8");
2135 Eigen::Vector3f
v = Eigen::Vector3f::Zero();
2137 fprintf(
file,
"PARAMS_SE3OFFSET %d %f %f %f %f %f %f %f\n",
2149 std::map<int, bool> isLandmarkWithRotation;
2150 for(std::multimap<int, Link>::const_iterator
iter = edgeConstraints.begin();
iter!=edgeConstraints.end(); ++
iter)
2152 int landmarkId =
iter->second.from() < 0?
iter->second.from():
iter->second.to() < 0?
iter->second.to():0;
2153 if(landmarkId != 0 && isLandmarkWithRotation.find(landmarkId) == isLandmarkWithRotation.end())
2157 if (1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
2159 isLandmarkWithRotation.insert(std::make_pair(landmarkId,
false));
2160 UDEBUG(
"Tag %d has no orientation", landmarkId);
2164 isLandmarkWithRotation.insert(std::make_pair(landmarkId,
true));
2165 UDEBUG(
"Tag %d has orientation", landmarkId);
2168 else if (1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(3,3)) >= 9999.0 ||
2169 1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(4,4)) >= 9999.0 ||
2170 1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
2172 isLandmarkWithRotation.insert(std::make_pair(landmarkId,
false));
2173 UDEBUG(
"Tag %d has no orientation", landmarkId);
2177 isLandmarkWithRotation.insert(std::make_pair(landmarkId,
true));
2178 UDEBUG(
"Tag %d has orientation", landmarkId);
2184 int landmarkOffset = poses.size()&&poses.rbegin()->first>0?poses.rbegin()->first:0;
2185 for(std::map<int, Transform>::const_iterator
iter = poses.begin();
iter!=poses.end(); ++
iter)
2192 fprintf(
file,
"VERTEX_SE2 %d %f %f %f\n",
2196 iter->second.theta());
2200 if(
uValue(isLandmarkWithRotation,
iter->first,
false))
2203 fprintf(
file,
"VERTEX_SE2 %d %f %f %f\n",
2204 landmarkOffset-
iter->first,
2207 iter->second.theta());
2212 fprintf(
file,
"VERTEX_XY %d %f %f\n",
2213 landmarkOffset-
iter->first,
2225 fprintf(
file,
"VERTEX_SE3:QUAT %d %f %f %f %f %f %f %f\n",
2237 if(
uValue(isLandmarkWithRotation,
iter->first,
false))
2241 fprintf(
file,
"VERTEX_SE3:QUAT %d %f %f %f %f %f %f %f\n",
2242 landmarkOffset-
iter->first,
2254 fprintf(
file,
"VERTEX_TRACKXYZ %d %f %f %f\n",
2255 landmarkOffset-
iter->first,
2264 int virtualVertexId = landmarkOffset - (poses.size()&&poses.rbegin()->first<0?poses.rbegin()->first:0);
2265 for(std::multimap<int, Link>::const_iterator
iter = edgeConstraints.begin();
iter!=edgeConstraints.end(); ++
iter)
2275 if(
uValue(isLandmarkWithRotation,
iter->first,
false))
2278 fprintf(
file,
"EDGE_SE2 %d %d %f %f %f %f %f %f %f %f %f\n",
2279 iter->second.from()<0?landmarkOffset-
iter->second.from():
iter->second.from(),
2280 iter->second.to()<0?landmarkOffset-
iter->second.to():
iter->second.to(),
2281 iter->second.transform().x(),
2282 iter->second.transform().y(),
2283 iter->second.transform().theta(),
2284 iter->second.infMatrix().at<
double>(0, 0),
2285 iter->second.infMatrix().at<
double>(0, 1),
2286 iter->second.infMatrix().at<
double>(0, 5),
2287 iter->second.infMatrix().at<
double>(1, 1),
2288 iter->second.infMatrix().at<
double>(1, 5),
2289 iter->second.infMatrix().at<
double>(5, 5));
2294 fprintf(
file,
"EDGE_SE2_XY %d %d %f %f %f %f %f\n",
2295 iter->second.from()<0?landmarkOffset-
iter->second.from():
iter->second.from(),
2296 iter->second.to()<0?landmarkOffset-
iter->second.to():
iter->second.to(),
2297 iter->second.transform().x(),
2298 iter->second.transform().y(),
2299 iter->second.infMatrix().at<
double>(0, 0),
2300 iter->second.infMatrix().at<
double>(0, 1),
2301 iter->second.infMatrix().at<
double>(1, 1));
2306 if(
uValue(isLandmarkWithRotation,
iter->first,
false))
2310 fprintf(
file,
"EDGE_SE3 %d %d %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n",
2311 iter->second.from()<0?landmarkOffset-
iter->second.from():
iter->second.from(),
2312 iter->second.to()<0?landmarkOffset-
iter->second.to():
iter->second.to(),
2313 iter->second.transform().x(),
2314 iter->second.transform().y(),
2315 iter->second.transform().z(),
2320 iter->second.infMatrix().at<
double>(0, 0),
2321 iter->second.infMatrix().at<
double>(0, 1),
2322 iter->second.infMatrix().at<
double>(0, 2),
2323 iter->second.infMatrix().at<
double>(0, 3),
2324 iter->second.infMatrix().at<
double>(0, 4),
2325 iter->second.infMatrix().at<
double>(0, 5),
2326 iter->second.infMatrix().at<
double>(1, 1),
2327 iter->second.infMatrix().at<
double>(1, 2),
2328 iter->second.infMatrix().at<
double>(1, 3),
2329 iter->second.infMatrix().at<
double>(1, 4),
2330 iter->second.infMatrix().at<
double>(1, 5),
2331 iter->second.infMatrix().at<
double>(2, 2),
2332 iter->second.infMatrix().at<
double>(2, 3),
2333 iter->second.infMatrix().at<
double>(2, 4),
2334 iter->second.infMatrix().at<
double>(2, 5),
2335 iter->second.infMatrix().at<
double>(3, 3),
2336 iter->second.infMatrix().at<
double>(3, 4),
2337 iter->second.infMatrix().at<
double>(3, 5),
2338 iter->second.infMatrix().at<
double>(4, 4),
2339 iter->second.infMatrix().at<
double>(4, 5),
2340 iter->second.infMatrix().at<
double>(5, 5));
2345 fprintf(
file,
"EDGE_SE3_TRACKXYZ %d %d %d %f %f %f %f %f %f %f %f %f\n",
2346 iter->second.from()<0?landmarkOffset-
iter->second.from():
iter->second.from(),
2347 iter->second.to()<0?landmarkOffset-
iter->second.to():
iter->second.to(),
2349 iter->second.transform().x(),
2350 iter->second.transform().y(),
2351 iter->second.transform().z(),
2352 iter->second.infMatrix().at<
double>(0, 0),
2353 iter->second.infMatrix().at<
double>(0, 1),
2354 iter->second.infMatrix().at<
double>(0, 2),
2355 iter->second.infMatrix().at<
double>(1, 1),
2356 iter->second.infMatrix().at<
double>(1, 2),
2357 iter->second.infMatrix().at<
double>(2, 2));
2363 std::string prefix =
isSlam2d()?
"EDGE_SE2" :
"EDGE_SE3:QUAT";
2364 std::string suffix =
"";
2365 std::string to =
uFormat(
" %d",
iter->second.to());
2382 if (1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
2384 prefix =
"EDGE_PRIOR_SE2_XY";
2389 prefix =
"EDGE_PRIOR_SE2";
2398 if (1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(3,3)) >= 9999.0 ||
2399 1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(4,4)) >= 9999.0 ||
2400 1 /
static_cast<double>(
iter->second.infMatrix().at<
double>(5,5)) >= 9999.0)
2403 prefix =
"EDGE_POINTXYZ_PRIOR";
2409 prefix =
"EDGE_SE3_PRIOR";
2417 fprintf(
file,
"VERTEX_SWITCH %d 1\n", virtualVertexId);
2418 fprintf(
file,
"EDGE_SWITCH_PRIOR %d 1 1.0\n", virtualVertexId);
2419 prefix =
isSlam2d() ?
"EDGE_SE2_SWITCHABLE" :
"EDGE_SE3_SWITCHABLE";
2420 suffix =
uFormat(
" %d", virtualVertexId++);
2429 fprintf(
file,
"%s %d%s%s %f %f %f %f %f %f %f %f %f\n",
2431 iter->second.from(),
2434 iter->second.transform().x(),
2435 iter->second.transform().y(),
2436 iter->second.transform().theta(),
2437 iter->second.infMatrix().at<
double>(0, 0),
2438 iter->second.infMatrix().at<
double>(0, 1),
2439 iter->second.infMatrix().at<
double>(0, 5),
2440 iter->second.infMatrix().at<
double>(1, 1),
2441 iter->second.infMatrix().at<
double>(1, 5),
2442 iter->second.infMatrix().at<
double>(5, 5));
2448 fprintf(
file,
"%s %d%s%s %f %f %f %f %f\n",
2450 iter->second.from(),
2453 iter->second.transform().x(),
2454 iter->second.transform().y(),
2455 iter->second.infMatrix().at<
double>(0, 0),
2456 iter->second.infMatrix().at<
double>(0, 1),
2457 iter->second.infMatrix().at<
double>(1, 1));
2467 fprintf(
file,
"%s %d%s%s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n",
2469 iter->second.from(),
2472 iter->second.transform().x(),
2473 iter->second.transform().y(),
2474 iter->second.transform().z(),
2479 iter->second.infMatrix().at<
double>(0, 0),
2480 iter->second.infMatrix().at<
double>(0, 1),
2481 iter->second.infMatrix().at<
double>(0, 2),
2482 iter->second.infMatrix().at<
double>(0, 3),
2483 iter->second.infMatrix().at<
double>(0, 4),
2484 iter->second.infMatrix().at<
double>(0, 5),
2485 iter->second.infMatrix().at<
double>(1, 1),
2486 iter->second.infMatrix().at<
double>(1, 2),
2487 iter->second.infMatrix().at<
double>(1, 3),
2488 iter->second.infMatrix().at<
double>(1, 4),
2489 iter->second.infMatrix().at<
double>(1, 5),
2490 iter->second.infMatrix().at<
double>(2, 2),
2491 iter->second.infMatrix().at<
double>(2, 3),
2492 iter->second.infMatrix().at<
double>(2, 4),
2493 iter->second.infMatrix().at<
double>(2, 5),
2494 iter->second.infMatrix().at<
double>(3, 3),
2495 iter->second.infMatrix().at<
double>(3, 4),
2496 iter->second.infMatrix().at<
double>(3, 5),
2497 iter->second.infMatrix().at<
double>(4, 4),
2498 iter->second.infMatrix().at<
double>(4, 5),
2499 iter->second.infMatrix().at<
double>(5, 5));
2505 fprintf(
file,
"%s %d%s%s %f %f %f %f %f %f %f %f %f\n",
2507 iter->second.from(),
2510 iter->second.transform().x(),
2511 iter->second.transform().y(),
2512 iter->second.transform().z(),
2513 iter->second.infMatrix().at<
double>(0, 0),
2514 iter->second.infMatrix().at<
double>(0, 1),
2515 iter->second.infMatrix().at<
double>(0, 2),
2516 iter->second.infMatrix().at<
double>(1, 1),
2517 iter->second.infMatrix().at<
double>(1, 2),
2518 iter->second.infMatrix().at<
double>(2, 2));
2522 UINFO(
"Graph saved to %s", fileName.c_str());
2528 UERROR(
"Cannot save to file %s", fileName.c_str());