@@ -60,6 +60,7 @@ static const int NUM_MEAS = NUM_POSES * DoF + NUM_FACTORS * Dim;
6060static const int MAX_ITER = 20 ; // for the solver
6161
6262
63+ // bundle state type to optimize over
6364using BundleT = manif::Bundle<double ,
6465 manif::SE2,
6566 manif::SE2,
@@ -76,35 +77,38 @@ using BundleT = manif::Bundle<double,
7677template <std::size_t XI, std::size_t XJ>
7778void add_pose_factor (
7879 const BundleT & X,
79- const manif::SE2d::Tangent & control,
80+ const SE2Tangentd & control,
8081 const MatrixT & W,
8182 Eigen::Ref<Eigen::Matrix<double , 3 , 1 >> r,
8283 Eigen::Ref<Eigen::Matrix<double , 3 , BundleT::DoF>> J)
8384{
85+ // index start position and length in the DoF of BundleT
8486 constexpr int BegI = manif::internal::intseq_element<XI, BundleT::BegDoF>::value;
8587 constexpr int LenI = manif::internal::intseq_element<XI, BundleT::LenDoF>::value;
8688 constexpr int BegJ = manif::internal::intseq_element<XJ, BundleT::BegDoF>::value;
8789 constexpr int LenJ = manif::internal::intseq_element<XJ, BundleT::LenDoF>::value;
8890
8991 MatrixT J_d_xi, J_d_xj; // Jacobian of motion wrt poses i and j
9092
91- auto d = X.block <XJ>().rminus (X.block <XI>(), J_d_xj, J_d_xi);
93+ auto d = X.element <XJ>().rminus (X.element <XI>(), J_d_xj, J_d_xi);
9294 r = W * (d - control).coeffs ();
9395 J.setZero ();
9496 J.block <3 , LenI>(0 , BegI) = W * J_d_xi;
9597 J.block <3 , LenJ>(0 , BegJ) = W * J_d_xj;
9698}
9799
100+
98101// Insert a landmark measurement factor of landmark LK
99102// from pose XI into the residual-jacobian pair (r, J)
100103template <std::size_t XI, std::size_t LK>
101104void add_beacon_factor (
102105 const BundleT & X,
103- const Eigen::Vector2d & measurement,
104- const Eigen::Matrix2d & S,
106+ const VectorY & measurement,
107+ const MatrixY & S,
105108 Eigen::Ref<Eigen::Matrix<double , 2 , 1 >> r,
106109 Eigen::Ref<Eigen::Matrix<double , 2 , BundleT::DoF>> J)
107110{
111+ // index start position and length in the DoF of BundleT
108112 constexpr int BegX = manif::internal::intseq_element<XI, BundleT::BegDoF>::value;
109113 constexpr int LenX = manif::internal::intseq_element<XI, BundleT::LenDoF>::value;
110114 constexpr int BegLMK = manif::internal::intseq_element<NUM_POSES + LK, BundleT::BegDoF>::value;
@@ -115,7 +119,7 @@ void add_beacon_factor(
115119 MatrixYT J_e_x; // Jacobian of measurement expectation wrt pose
116120 MatrixYB J_e_b; // Jacobian of measurement expectation wrt lmk
117121
118- auto e = X.block <XI>().inverse (J_ix_x).act (X.block <NUM_POSES + LK>().coeffs (), J_e_ix, J_e_b);
122+ auto e = X.element <XI>().inverse (J_ix_x).act (X.element <NUM_POSES + LK>().coeffs (), J_e_ix, J_e_b);
119123 J_e_x = J_e_ix * J_ix_x;
120124 r = S * (e - measurement);
121125 J.setZero ();
@@ -279,15 +283,15 @@ int main()
279283 manif::R2d (landmarks[4 ]));
280284
281285 cout << " prior" << std::showpos << endl;
282- cout << " pose :" << X.block <0 >().translation ().transpose () << " " << X.block <0 >().angle () << endl;
283- cout << " pose :" << X.block <1 >().translation ().transpose () << " " << X.block <1 >().angle () << endl;
284- cout << " pose :" << X.block <2 >().translation ().transpose () << " " << X.block <2 >().angle () << endl;
285-
286- cout << " lmk :" << X.block <3 >() << endl;
287- cout << " lmk :" << X.block <4 >() << endl;
288- cout << " lmk :" << X.block <5 >() << endl;
289- cout << " lmk :" << X.block <6 >() << endl;
290- cout << " lmk :" << X.block <7 >() << endl;
286+ cout << " pose :" << X.element <0 >().translation ().transpose () << " " << X.element <0 >().angle () << endl;
287+ cout << " pose :" << X.element <1 >().translation ().transpose () << " " << X.element <1 >().angle () << endl;
288+ cout << " pose :" << X.element <2 >().translation ().transpose () << " " << X.element <2 >().angle () << endl;
289+
290+ cout << " lmk :" << X.element <3 >() << endl;
291+ cout << " lmk :" << X.element <4 >() << endl;
292+ cout << " lmk :" << X.element <5 >() << endl;
293+ cout << " lmk :" << X.element <6 >() << endl;
294+ cout << " lmk :" << X.element <7 >() << endl;
291295 cout << " -----------------------------------------------" << endl;
292296
293297 // iterate
@@ -301,7 +305,7 @@ int main()
301305 int row = 0 ; // keep track of row in J and r
302306
303307 // first residual: prior
304- r.segment <DoF>(row) = X.block <0 >().lminus (SE2d::Identity (), J.block <DoF, DoF>(row, 0 )).coeffs ();
308+ r.segment <DoF>(row) = X.element <0 >().lminus (SE2d::Identity (), J.block <DoF, DoF>(row, 0 )).coeffs ();
305309 row += DoF;
306310
307311 // motion residuals
@@ -338,7 +342,10 @@ int main()
338342 // compute optimal step
339343 // ATTENTION: This is an expensive step!!
340344 // ATTENTION: Use QR factorization and column reordering for larger problems!!
341- X += BundleT::Tangent (-(J.transpose () * J).inverse () * J.transpose () * r);
345+ const auto dX = (-(J.transpose () * J).inverse () * J.transpose () * r).eval ();
346+
347+ // update estimate
348+ X += BundleT::Tangent (dX);
342349
343350 // DEBUG INFO
344351 cout << " residual norm: " << std::scientific << r.norm () << " , step norm: " << dX.norm () << endl;
@@ -355,15 +362,15 @@ int main()
355362
356363 // solved problem
357364 cout << " posterior" << std::showpos << endl;
358- cout << " pose :" << X.block <0 >().translation ().transpose () << " " << X.block <0 >().angle () << endl;
359- cout << " pose :" << X.block <1 >().translation ().transpose () << " " << X.block <1 >().angle () << endl;
360- cout << " pose :" << X.block <2 >().translation ().transpose () << " " << X.block <2 >().angle () << endl;
361-
362- cout << " lmk :" << X.block <3 >() << endl;
363- cout << " lmk :" << X.block <4 >() << endl;
364- cout << " lmk :" << X.block <5 >() << endl;
365- cout << " lmk :" << X.block <6 >() << endl;
366- cout << " lmk :" << X.block <7 >() << endl;
365+ cout << " pose :" << X.element <0 >().translation ().transpose () << " " << X.element <0 >().angle () << endl;
366+ cout << " pose :" << X.element <1 >().translation ().transpose () << " " << X.element <1 >().angle () << endl;
367+ cout << " pose :" << X.element <2 >().translation ().transpose () << " " << X.element <2 >().angle () << endl;
368+
369+ cout << " lmk :" << X.element <3 >() << endl;
370+ cout << " lmk :" << X.element <4 >() << endl;
371+ cout << " lmk :" << X.element <5 >() << endl;
372+ cout << " lmk :" << X.element <6 >() << endl;
373+ cout << " lmk :" << X.element <7 >() << endl;
367374
368375 cout << " -----------------------------------------------" << endl;
369376
0 commit comments