5 #ifndef __pinocchio_multibody_liegroup_special_euclidean_operation_hpp__ 6 #define __pinocchio_multibody_liegroup_special_euclidean_operation_hpp__ 10 #include "pinocchio/macros.hpp" 11 #include "pinocchio/math/quaternion.hpp" 12 #include "pinocchio/spatial/fwd.hpp" 13 #include "pinocchio/utils/static-if.hpp" 14 #include "pinocchio/spatial/se3.hpp" 15 #include "pinocchio/multibody/liegroup/liegroup-base.hpp" 17 #include "pinocchio/multibody/liegroup/vector-space.hpp" 18 #include "pinocchio/multibody/liegroup/cartesian-product.hpp" 19 #include "pinocchio/multibody/liegroup/special-orthogonal.hpp" 23 template<
int Dim,
typename Scalar,
int Options = 0>
27 template<
int Dim,
typename Scalar,
int Options>
31 template<
typename _Scalar,
int _Options>
44 template<
typename _Scalar,
int _Options>
46 :
public LieGroupBase <SpecialEuclideanOperationTpl<2,_Scalar,_Options> >
54 typedef Eigen::Matrix<Scalar,2,2,Options>
Matrix2;
55 typedef Eigen::Matrix<Scalar,2,1,Options>
Vector2;
57 template<
typename TangentVector,
typename Matrix2Like,
typename Vector2Like>
58 static void exp(
const Eigen::MatrixBase<TangentVector> &
v,
59 const Eigen::MatrixBase<Matrix2Like> & R,
60 const Eigen::MatrixBase<Vector2Like> &
t)
63 EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(Matrix2Like, 2, 2);
64 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector2Like, 2);
67 const Scalar omega =
v(2);
68 Scalar cv,sv;
SINCOS(omega, &sv, &cv);
74 vcross -= -
v(1)*R.col(0) +
v(0)*R.col(1);
76 Scalar omega_abs = math::fabs(omega);
88 template<
typename Matrix2Like,
typename Vector2Like,
typename Matrix3Like>
90 const Eigen::MatrixBase<Vector2Like> &
t,
91 const Eigen::MatrixBase<Matrix3Like> & M,
94 EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(Matrix2Like, 2, 2);
95 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector2Like, 2);
105 Mout.template topLeftCorner<2,2>() = R.transpose();
106 Mout.template topRightCorner<2,1>() = tinv;
107 Mout.template bottomLeftCorner<1,2>().
setZero();
111 Mout.template topLeftCorner<2,2>() += R.transpose();
112 Mout.template topRightCorner<2,1>() += tinv;
116 Mout.template topLeftCorner<2,2>() -= R.transpose();
117 Mout.template topRightCorner<2,1>() -= tinv;
121 assert(
false &&
"Wrong Op requesed value");
129 template<
typename Matrix2Like,
typename Vector2Like,
typename TangentVector>
130 static void log(
const Eigen::MatrixBase<Matrix2Like> & R,
131 const Eigen::MatrixBase<Vector2Like> &
p,
132 const Eigen::MatrixBase<TangentVector> &
v)
134 EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(Matrix2Like, 2, 2);
135 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector2Like, 2);
143 const Scalar1 tabs = math::fabs(t);
144 const Scalar1 t2 = t*
t;
145 Scalar1 st,ct;
SINCOS(tabs, &st, &ct);
148 1 - t2/12 - t2*t2/720,
151 vout.template head<2>().noalias() = alpha *
p;
152 vout(0) += t/2 *
p(1);
153 vout(1) += -t/2 *
p(0);
157 template<
typename Matrix2Like,
typename Vector2Like,
typename JacobianOutLike>
158 static void Jlog(
const Eigen::MatrixBase<Matrix2Like> & R,
159 const Eigen::MatrixBase<Vector2Like> &
p,
160 const Eigen::MatrixBase<JacobianOutLike> &
J)
162 EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(Matrix2Like, 2, 2);
163 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector2Like, 2);
171 const Scalar1 tabs = math::fabs(t);
172 Scalar1 alpha, alpha_dot;
174 Scalar1 st,ct;
SINCOS(t, &st, &ct);
175 Scalar1 inv_2_1_ct = 0.5 / (1-ct);
181 - t / 6 - t2*t / 180,
182 (st-t) * inv_2_1_ct);
185 V(0,0) =
V(1,1) = alpha;
189 Jout.template topLeftCorner <2,2>().noalias() =
V * R;
190 Jout.template topRightCorner<2,1>() << alpha_dot*p[0] + p[1]/2, -
p(0)/2 + alpha_dot*
p(1);
191 Jout.template bottomLeftCorner<1,2>().
setZero();
217 return std::string(
"SE(2)");
220 template <
class ConfigL_t,
class ConfigR_t,
class Tangent_t>
222 const Eigen::MatrixBase<ConfigR_t> &
q1,
223 const Eigen::MatrixBase<Tangent_t> &
d)
225 Matrix2 R0, R1; Vector2 t0, t1;
228 Matrix2
R (R0.transpose() * R1);
229 Vector2
t (R0.transpose() * (t1 - t0));
234 template <ArgumentPosition arg,
class ConfigL_t,
class ConfigR_t,
class JacobianOut_t>
236 const Eigen::MatrixBase<ConfigR_t> &
q1,
237 const Eigen::MatrixBase<JacobianOut_t>&
J)
const 239 Matrix2 R0, R1; Vector2 t0, t1;
242 Matrix2
R (R0.transpose() * R1);
243 Vector2
t (R0.transpose() * (t1 - t0));
250 Vector2 pcross (
q1(1) -
q0(1),
q0(0) -
q1(0));
253 J0.template topLeftCorner <2,2> ().noalias() = -
R.transpose();
254 J0.template topRightCorner<2,1> ().noalias() = R1.transpose() * pcross;
255 J0.template bottomLeftCorner <1,2> ().
setZero();
257 J0.applyOnTheLeft(J1);
258 }
else if (arg ==
ARG1) {
263 template <
class ConfigIn_t,
class Velocity_t,
class ConfigOut_t>
265 const Eigen::MatrixBase<Velocity_t> &
v,
266 const Eigen::MatrixBase<ConfigOut_t> & qout)
275 out.template head<2>().noalias() = R0 * t + t0;
276 out.template tail<2>().noalias() = R0 * R.col(0);
279 template <
class Config_t,
class Jacobian_t>
281 const Eigen::MatrixBase<Jacobian_t> &
J)
283 assert(J.rows() ==
nq() && J.cols() ==
nv() &&
"J is not of the right dimension");
291 Jout.template topLeftCorner<2,2>() << c_theta, -s_theta, s_theta, c_theta;
292 Jout.template bottomRightCorner<2,1>() << -s_theta, c_theta;
295 template <
class Config_t,
class Tangent_t,
class JacobianOut_t>
297 const Eigen::MatrixBase<Tangent_t> &
v,
298 const Eigen::MatrixBase<JacobianOut_t>&
J,
307 toInverseActionMatrix(R, t, Jout, op);
310 template <
class Config_t,
class Tangent_t,
class JacobianOut_t>
312 const Eigen::MatrixBase<Tangent_t> &
v,
313 const Eigen::MatrixBase<JacobianOut_t> &
J,
319 Eigen::Matrix<Scalar,6,6> Jtmp6;
325 Jout << Jtmp6.template topLeftCorner<2,2>(), Jtmp6.template topRightCorner<2,1>(),
326 Jtmp6.template bottomLeftCorner<1,2>(), Jtmp6.template bottomRightCorner<1,1>();
329 Jout.template topLeftCorner<2,2>() += Jtmp6.template topLeftCorner<2,2>();
330 Jout.template topRightCorner<2,1>() += Jtmp6.template topRightCorner<2,1>();
331 Jout.template bottomLeftCorner<1,2>() += Jtmp6.template bottomLeftCorner<1,2>();
332 Jout.template bottomRightCorner<1,1>() += Jtmp6.template bottomRightCorner<1,1>();
335 Jout.template topLeftCorner<2,2>() -= Jtmp6.template topLeftCorner<2,2>();
336 Jout.template topRightCorner<2,1>() -= Jtmp6.template topRightCorner<2,1>();
337 Jout.template bottomLeftCorner<1,2>() -= Jtmp6.template bottomLeftCorner<1,2>();
338 Jout.template bottomRightCorner<1,1>() -= Jtmp6.template bottomRightCorner<1,1>();
341 assert(
false &&
"Wrong Op requesed value");
346 template <
class Config_t,
class Tangent_t,
class JacobianIn_t,
class JacobianOut_t>
348 const Eigen::MatrixBase<Tangent_t> &
v,
349 const Eigen::MatrixBase<JacobianIn_t> & Jin,
350 const Eigen::MatrixBase<JacobianOut_t> & J_out)
const 357 Vector2 tinv = (R.transpose() *
t).reverse();
360 Jout.template topRows<2>().noalias() = R.transpose() * Jin.template topRows<2>();
361 Jout.template topRows<2>().noalias() += tinv * Jin.template bottomRows<1>();
362 Jout.template bottomRows<1>() = Jin.template bottomRows<1>();
365 template <
class Config_t,
class Tangent_t,
class JacobianIn_t,
class JacobianOut_t>
367 const Eigen::MatrixBase<Tangent_t> &
v,
368 const Eigen::MatrixBase<JacobianIn_t> & Jin,
369 const Eigen::MatrixBase<JacobianOut_t> & J_out)
const 374 Eigen::Matrix<Scalar,6,6> Jtmp6;
377 Jout.template topRows<2>().noalias()
378 = Jtmp6.template topLeftCorner<2,2>() * Jin.template topRows<2>();
379 Jout.template topRows<2>().noalias()
380 += Jtmp6.template topRightCorner<2,1>() * Jin.template bottomRows<1>();
381 Jout.template bottomRows<1>().noalias()
382 = Jtmp6.template bottomLeftCorner<1,2>()* Jin.template topRows<2>();
383 Jout.template bottomRows<1>().noalias()
384 += Jtmp6.template bottomRightCorner<1,1>() * Jin.template bottomRows<1>();
388 template <
class Config_t,
class Tangent_t,
class Jacobian_t>
390 const Eigen::MatrixBase<Tangent_t> &
v,
391 const Eigen::MatrixBase<Jacobian_t> &
J)
const 398 Vector2 tinv = (R.transpose() *
t).reverse();
401 Jout.template topRows<2>() = R.transpose() * Jout.template topRows<2>();
403 Jout.template topRows<2>().noalias() += tinv * Jout.template bottomRows<1>();
406 template <
class Config_t,
class Tangent_t,
class Jacobian_t>
408 const Eigen::MatrixBase<Tangent_t> &
v,
409 const Eigen::MatrixBase<Jacobian_t> &
J)
const 414 Eigen::Matrix<Scalar,6,6> Jtmp6;
417 Jout.template topRows<2>()
418 = Jtmp6.template topLeftCorner<2,2>() * Jout.template topRows<2>();
419 Jout.template topRows<2>().noalias()
420 += Jtmp6.template topRightCorner<2,1>() * Jout.template bottomRows<1>();
421 Jout.template bottomRows<1>()
422 = Jtmp6.template bottomRightCorner<1,1>() * Jout.template bottomRows<1>();
423 Jout.template bottomRows<1>().noalias()
424 += Jtmp6.template bottomLeftCorner<1,2>() * Jout.template topRows<2>();
427 template <
class Config_t>
433 template <
class Config_t>
437 const Scalar norm =
Scalar(qin.template tail<2>().norm());
439 return abs(norm -
Scalar(1.0)) < prec;
442 template <
class Config_t>
445 R2crossSO2_t().random(qout);
448 template <
class ConfigL_t,
class ConfigR_t,
class ConfigOut_t>
450 const Eigen::MatrixBase<ConfigR_t> & upper,
451 const Eigen::MatrixBase<ConfigOut_t> & qout)
454 R2crossSO2_t ().randomConfiguration(lower, upper, qout);
457 template <
class ConfigL_t,
class ConfigR_t>
459 const Eigen::MatrixBase<ConfigR_t> &
q1,
462 return R2crossSO2_t().isSameConfiguration(q0, q1, prec);
467 template<
typename Matrix2Like,
typename Vector2Like,
typename Vector4Like>
469 const Eigen::MatrixBase<Vector2Like> &
t,
470 const Eigen::MatrixBase<Vector4Like> &
q)
472 EIGEN_STATIC_ASSERT_SAME_MATRIX_SIZE(Matrix2Like, Matrix2);
473 EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(Vector2Like, Vector2);
484 template<
typename _Scalar,
int _Options>
497 template<
typename _Scalar,
int _Options>
499 :
public LieGroupBase <SpecialEuclideanOperationTpl<3,_Scalar,_Options> >
534 return std::string(
"SE(3)");
537 template <
class ConfigL_t,
class ConfigR_t,
class Tangent_t>
539 const Eigen::MatrixBase<ConfigR_t> &
q1,
540 const Eigen::MatrixBase<Tangent_t> &
d)
542 ConstQuaternionMap_t quat0 (q0.derived().template tail<4>().
data());
544 ConstQuaternionMap_t quat1 (q1.derived().template tail<4>().
data());
548 =
log6(
SE3(quat0.matrix(), q0.derived().template head<3>()).
inverse()
549 *
SE3(quat1.matrix(), q1.derived().template head<3>())).toVector();
553 template <ArgumentPosition arg,
class ConfigL_t,
class ConfigR_t,
class JacobianOut_t>
555 const Eigen::MatrixBase<ConfigR_t> &
q1,
556 const Eigen::MatrixBase<JacobianOut_t> &
J)
const 561 ConstQuaternionMap_t quat0 (q0.derived().template tail<4>().
data());
563 ConstQuaternionMap_t quat1 (q1.derived().template tail<4>().
data());
566 Matrix3 R0(quat0.matrix()), R1 (quat1.matrix());
569 const SE3
M (
SE3(R0, q0.template head<3>()).
inverse()
570 *
SE3(R1, q1.template head<3>()));
576 const Vector3 p1_p0 = R1.transpose()*(q1.template head<3>() - q0.template head<3>());
579 J0.template bottomRightCorner<3,3> ().noalias() = J0.template topLeftCorner <3,3> ().noalias() = -
M.rotation().transpose();
580 J0.template topRightCorner<3,3> ().noalias() =
skew(p1_p0) *
M.rotation().transpose();
581 J0.template bottomLeftCorner<3,3> ().
setZero();
582 J0.applyOnTheLeft(J1);
584 else if (arg ==
ARG1) {
589 template <
class ConfigIn_t,
class Velocity_t,
class ConfigOut_t>
591 const Eigen::MatrixBase<Velocity_t> &
v,
592 const Eigen::MatrixBase<ConfigOut_t> & qout)
595 Quaternion_t
const quat(q.derived().template tail<4>());
597 QuaternionMap_t res_quat (out.template tail<4>().data());
601 SE3 M0 (quat.matrix(), q.derived().template head<3>());
603 SE3 M1 (M0 *
exp6(mref_v));
605 out.template head<3>() = M1.translation();
607 const Scalar dot_product = res_quat.dot(quat);
608 for(Eigen::DenseIndex k = 0; k < 4; ++k)
611 -res_quat.coeffs().coeff(k),
612 res_quat.coeffs().coeff(k));
621 template <
class Config_t,
class Jacobian_t>
623 const Eigen::MatrixBase<Jacobian_t> &
J)
625 assert(J.rows() ==
nq() && J.cols() ==
nv() &&
"J is not of the right dimension");
634 ConstQuaternionMap_t quat_map(q.derived().template tail<4>().
data());
636 Jout.template topLeftCorner<3,3>() = quat_map.toRotationMatrix();
639 typedef Eigen::Matrix<Scalar,4,3,JacobianPlainType::Options|Eigen::RowMajor> Jacobian43;
641 Jacobian43 Jexp3QuatCoeffWise;
653 if(quat_map.coeffs()[3] >=
Scalar(0))
659 template <
class Config_t,
class Tangent_t,
class JacobianOut_t>
661 const Eigen::MatrixBase<Tangent_t> &
v,
662 const Eigen::MatrixBase<JacobianOut_t>&
J,
679 assert(
false &&
"Wrong Op requesed value");
684 template <
class Config_t,
class Tangent_t,
class JacobianOut_t>
686 const Eigen::MatrixBase<Tangent_t> &
v,
687 const Eigen::MatrixBase<JacobianOut_t>&
J,
702 assert(
false &&
"Wrong Op requesed value");
707 template <
class Config_t,
class Tangent_t,
class JacobianIn_t,
class JacobianOut_t>
709 const Eigen::MatrixBase<Tangent_t> &
v,
710 const Eigen::MatrixBase<JacobianIn_t> & Jin,
711 const Eigen::MatrixBase<JacobianOut_t> & J_out)
const 714 Eigen::Matrix<Scalar,6,6> Jtmp6;
717 Jout.template topRows<3>().noalias()
718 = Jtmp6.template topLeftCorner<3,3>() * Jin.template topRows<3>();
719 Jout.template topRows<3>().noalias()
720 += Jtmp6.template topRightCorner<3,3>() * Jin.template bottomRows<3>();
721 Jout.template bottomRows<3>().noalias()
722 = Jtmp6.template bottomRightCorner<3,3>() * Jin.template bottomRows<3>();
725 template <
class Config_t,
class Tangent_t,
class JacobianIn_t,
class JacobianOut_t>
727 const Eigen::MatrixBase<Tangent_t> &
v,
728 const Eigen::MatrixBase<JacobianIn_t> & Jin,
729 const Eigen::MatrixBase<JacobianOut_t> & J_out)
const 732 Eigen::Matrix<Scalar,6,6> Jtmp6;
735 Jout.template topRows<3>().noalias()
736 = Jtmp6.template topLeftCorner<3,3>() * Jin.template topRows<3>();
737 Jout.template topRows<3>().noalias()
738 += Jtmp6.template topRightCorner<3,3>() * Jin.template bottomRows<3>();
739 Jout.template bottomRows<3>().noalias()
740 = Jtmp6.template bottomRightCorner<3,3>() * Jin.template bottomRows<3>();
744 template <
class Config_t,
class Tangent_t,
class Jacobian_t>
746 const Eigen::MatrixBase<Tangent_t> &
v,
747 const Eigen::MatrixBase<Jacobian_t> & J_out)
const 750 Eigen::Matrix<Scalar,6,6> Jtmp6;
754 Jout.template topRows<3>()
755 = Jtmp6.template topLeftCorner<3,3>() * Jout.template topRows<3>();
756 Jout.template topRows<3>().noalias()
757 += Jtmp6.template topRightCorner<3,3>() * Jout.template bottomRows<3>();
758 Jout.template bottomRows<3>()
759 = Jtmp6.template bottomRightCorner<3,3>() * Jout.template bottomRows<3>();
762 template <
class Config_t,
class Tangent_t,
class Jacobian_t>
764 const Eigen::MatrixBase<Tangent_t> &
v,
765 const Eigen::MatrixBase<Jacobian_t> & J_out)
const 768 Eigen::Matrix<Scalar,6,6> Jtmp6;
771 Jout.template topRows<3>()
772 = Jtmp6.template topLeftCorner<3,3>() * Jout.template topRows<3>();
773 Jout.template topRows<3>().noalias()
774 += Jtmp6.template topRightCorner<3,3>() * Jout.template bottomRows<3>();
775 Jout.template bottomRows<3>()
776 = Jtmp6.template bottomRightCorner<3,3>() * Jout.template bottomRows<3>();
779 template <
class ConfigL_t,
class ConfigR_t>
781 const Eigen::MatrixBase<ConfigR_t> &
q1)
784 difference_impl(q0, q1, t);
785 return t.squaredNorm();
788 template <
class Config_t>
795 template <
class Config_t>
801 return abs(norm -
Scalar(1.0)) < prec;
804 template <
class Config_t>
807 R3crossSO3_t().random(qout);
810 template <
class ConfigL_t,
class ConfigR_t,
class ConfigOut_t>
812 const Eigen::MatrixBase<ConfigR_t> & upper,
813 const Eigen::MatrixBase<ConfigOut_t> & qout)
816 R3crossSO3_t ().randomConfiguration(lower, upper, qout);
819 template <
class ConfigL_t,
class ConfigR_t>
821 const Eigen::MatrixBase<ConfigR_t> &
q1,
824 return R3crossSO3_t().isSameConfiguration(q0, q1, prec);
830 #endif // ifndef __pinocchio_multibody_liegroup_special_euclidean_operation_hpp__ VectorSpaceOperationTpl< 2, Scalar, Options > R2_t
JointCollectionTpl const Eigen::MatrixBase< ConfigVectorType > const Eigen::MatrixBase< TangentVectorType > & v
void dIntegrateTransport_dv_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianIn_t > &Jin, const Eigen::MatrixBase< JacobianOut_t > &J_out) const
static void exp(const Eigen::MatrixBase< TangentVector > &v, const Eigen::MatrixBase< Matrix2Like > &R, const Eigen::MatrixBase< Vector2Like > &t)
int nv(const JointModelTpl< Scalar, Options, JointCollectionTpl > &jmodel)
Visit a JointModelTpl through JointNvVisitor to get the dimension of the joint tangent space...
static void integrate_impl(const Eigen::MatrixBase< ConfigIn_t > &q, const Eigen::MatrixBase< Velocity_t > &v, const Eigen::MatrixBase< ConfigOut_t > &qout)
static bool isSameConfiguration_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1, const Scalar &prec)
void dDifference_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1, const Eigen::MatrixBase< JacobianOut_t > &J) const
Eigen::Matrix< Scalar, NQ, 1, Options > ConfigVector_t
static void integrateCoeffWiseJacobian_impl(const Eigen::MatrixBase< Config_t > &q, const Eigen::MatrixBase< Jacobian_t > &J)
#define PINOCCHIO_ASSERT_MATRIX_SPECIFIC_SIZE(type, M, nrows, ncols)
Ensure that a matrix (or vector) is of correct size (compile-time and run-time assertion) ...
static void toInverseActionMatrix(const Eigen::MatrixBase< Matrix2Like > &R, const Eigen::MatrixBase< Vector2Like > &t, const Eigen::MatrixBase< Matrix3Like > &M, const AssignmentOperatorType op)
static void integrateCoeffWiseJacobian_impl(const Eigen::MatrixBase< Config_t > &q, const Eigen::MatrixBase< Jacobian_t > &J)
static bool isSameConfiguration_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1, const Scalar &prec)
Eigen::Matrix< Scalar, NV, NV, Options > JacobianMatrix_t
int nq(const JointModelTpl< Scalar, Options, JointCollectionTpl > &jmodel)
Visit a JointModelTpl through JointNqVisitor to get the dimension of the joint configuration space...
static void dIntegrate_dv_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianOut_t > &J, const AssignmentOperatorType op=SETTO)
void Jlog3(const Scalar &theta, const Eigen::MatrixBase< Vector3Like > &log, const Eigen::MatrixBase< Matrix3Like > &Jlog)
Derivative of log3.
void dIntegrateTransport_dv_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< Jacobian_t > &J_out) const
Eigen::Matrix< Scalar, 2, 1, Options > Vector2
Eigen::Matrix< Scalar, NV, 1, Options > TangentVector_t
Eigen::Map< const Quaternion_t > ConstQuaternionMap_t
static std::string name()
void firstOrderNormalize(const Eigen::QuaternionBase< D > &q)
static void normalize_impl(const Eigen::MatrixBase< Config_t > &qout)
static void dIntegrate_dv_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianOut_t > &J, const AssignmentOperatorType op=SETTO)
void randomConfiguration_impl(const Eigen::MatrixBase< ConfigL_t > &lower, const Eigen::MatrixBase< ConfigR_t > &upper, const Eigen::MatrixBase< ConfigOut_t > &qout) const
void Jlog6(const SE3Tpl< Scalar, Options > &M, const Eigen::MatrixBase< Matrix6Like > &Jlog)
Derivative of log6.
MotionTpl< Scalar, Options > log6(const SE3Tpl< Scalar, Options > &M)
Log: SE3 -> se3.
#define PINOCCHIO_LIE_GROUP_TPL_PUBLIC_INTERFACE(Derived)
Eigen::Matrix< Scalar, 2, 2, Options > Matrix2
void Jexp6(const MotionDense< MotionDerived > &nu, const Eigen::MatrixBase< Matrix6Like > &Jexp)
Derivative of exp6 Computed as the inverse of Jlog6.
traits< LieGroupDerived >::Scalar Scalar
void dIntegrateTransport_dv_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< Jacobian_t > &J) const
JointCollectionTpl const Eigen::MatrixBase< ConfigVectorType > & q
#define PINOCCHIO_EIGEN_CONST_CAST(TYPE, OBJ)
Macro for an automatic const_cast.
static void normalize_impl(const Eigen::MatrixBase< Config_t > &qout)
bool isUnitary(const Eigen::MatrixBase< MatrixLike > &mat, const typename MatrixLike::RealScalar &prec=Eigen::NumTraits< typename MatrixLike::Scalar >::dummy_precision())
Check whether the input matrix is Unitary within the given precision.
void Jexp3CoeffWise(const Eigen::MatrixBase< Vector3Like > &v, const Eigen::MatrixBase< Matrix43Like > &Jexp)
Derivative of where is a small perturbation of at identity.
static ConfigVector_t neutral()
void dIntegrateTransport_dq_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< Jacobian_t > &J_out) const
static Index nv()
Get dimension of Lie Group tangent space.
SpecialOrthogonalOperationTpl< 2, Scalar, Options > SO2_t
void forwardKinematics(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, DataTpl< Scalar, Options, JointCollectionTpl > &data, const Eigen::MatrixBase< ConfigVectorType > &q)
Update the joint placements according to the current joint configuration.
void dIntegrateTransport_dq_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianIn_t > &Jin, const Eigen::MatrixBase< JacobianOut_t > &J_out) const
Eigen::Matrix< typename QuaternionLike::Scalar, 3, 1, PINOCCHIO_EIGEN_PLAIN_TYPE(typename QuaternionLike::Vector3)::Options > log3(const Eigen::QuaternionBase< QuaternionLike > &quat, typename QuaternionLike::Scalar &theta)
Same as log3 but with a unit quaternion as input.
void random_impl(const Eigen::MatrixBase< Config_t > &qout) const
SE3Tpl inverse() const
aXb = bXa.inverse()
Eigen::NumTraits< Scalar >::Real RealScalar
static Index nv()
Get dimension of Lie Group tangent space.
JointCollectionTpl const Eigen::MatrixBase< ConfigVectorIn1 > const Eigen::MatrixBase< ConfigVectorIn2 > & q1
void SINCOS(const S1 &a, S2 *sa, S3 *ca)
Computes sin/cos values of a given input scalar.
static bool isNormalized_impl(const Eigen::MatrixBase< Config_t > &qin, const Scalar &prec)
void dDifference_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1, const Eigen::MatrixBase< JacobianOut_t > &J) const
static Scalar squaredDistance_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1)
traits< SE3Tpl >::Vector3 Vector3
void dIntegrateTransport_dq_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianIn_t > &Jin, const Eigen::MatrixBase< JacobianOut_t > &J_out) const
SE3Tpl< Scalar, Options > Transformation_t
Main pinocchio namespace.
void dIntegrateTransport_dv_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianIn_t > &Jin, const Eigen::MatrixBase< JacobianOut_t > &J_out) const
ToVectorConstReturnType toVector() const
void random_impl(const Eigen::MatrixBase< Config_t > &qout) const
JointCollectionTpl const Eigen::MatrixBase< ConfigVectorIn1 > & q0
void assignQuaternion(Eigen::QuaternionBase< D > &quat, const Eigen::MatrixBase< Matrix3 > &R)
static ConfigVector_t neutral()
PINOCCHIO_EIGEN_PLAIN_TYPE(ConfigVectorType) integrate(const ModelTpl< Scalar
Integrate a configuration vector for the specified model for a tangent vector during one unit time...
void randomConfiguration_impl(const Eigen::MatrixBase< ConfigL_t > &lower, const Eigen::MatrixBase< ConfigR_t > &upper, const Eigen::MatrixBase< ConfigOut_t > &qout) const
Eigen::Map< Quaternion_t > QuaternionMap_t
bool isNormalized(const Eigen::QuaternionBase< Quaternion > &quat, const typename Quaternion::Coefficients::RealScalar &prec=Eigen::NumTraits< typename Quaternion::Coefficients::RealScalar >::dummy_precision())
Check whether the input quaternion is Normalized within the given precision.
#define PINOCCHIO_DEFAULT_QUATERNION_NORM_TOLERANCE_VALUE
static void dIntegrate_dq_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianOut_t > &J, const AssignmentOperatorType op=SETTO)
CartesianProductOperation< R2_t, SO2_t > R2crossSO2_t
static void Jlog(const Eigen::MatrixBase< Matrix2Like > &R, const Eigen::MatrixBase< Vector2Like > &p, const Eigen::MatrixBase< JacobianOutLike > &J)
static void difference_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1, const Eigen::MatrixBase< Tangent_t > &d)
static void difference_impl(const Eigen::MatrixBase< ConfigL_t > &q0, const Eigen::MatrixBase< ConfigR_t > &q1, const Eigen::MatrixBase< Tangent_t > &d)
static void forwardKinematics(const Eigen::MatrixBase< Matrix2Like > &R, const Eigen::MatrixBase< Vector2Like > &t, const Eigen::MatrixBase< Vector4Like > &q)
void setZero(std::vector< MatType, Eigen::aligned_allocator< MatType > > &Ms)
Common traits structure to fully define base classes for CRTP.
if_then_else_impl< LhsType, RhsType, ThenType, ElseType >::ReturnType if_then_else(const ComparisonOperators op, const LhsType &lhs_value, const RhsType &rhs_value, const ThenType &then_value, const ElseType &else_value)
static std::string name()
void skew(const Eigen::MatrixBase< Vector3 > &v, const Eigen::MatrixBase< Matrix3 > &M)
Computes the skew representation of a given 3d vector, i.e. the antisymmetric matrix representation o...
void normalize(const ModelTpl< Scalar, Options, JointCollectionTpl > &model, const Eigen::MatrixBase< ConfigVectorType > &qout)
Normalize a configuration vector.
SE3Tpl< Scalar, Options > SE3
CartesianProductOperation< VectorSpaceOperationTpl< 3, Scalar, Options >, SpecialOrthogonalOperationTpl< 3, Scalar, Options > > R3crossSO3_t
static bool isNormalized_impl(const Eigen::MatrixBase< Config_t > &qin, const Scalar &prec)
SE3Tpl< typename MotionDerived::Scalar, PINOCCHIO_EIGEN_PLAIN_TYPE(typename MotionDerived::Vector3)::Options > exp6(const MotionDense< MotionDerived > &nu)
Exp: se3 -> SE3.
traits< SE3Tpl >::Matrix3 Matrix3
static void integrate_impl(const Eigen::MatrixBase< ConfigIn_t > &q, const Eigen::MatrixBase< Velocity_t > &v, const Eigen::MatrixBase< ConfigOut_t > &qout)
void dIntegrateTransport_dq_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< Jacobian_t > &J) const
static void log(const Eigen::MatrixBase< Matrix2Like > &R, const Eigen::MatrixBase< Vector2Like > &p, const Eigen::MatrixBase< TangentVector > &v)
Eigen::Quaternion< Scalar, Options > Quaternion_t
static void dIntegrate_dq_impl(const Eigen::MatrixBase< Config_t > &, const Eigen::MatrixBase< Tangent_t > &v, const Eigen::MatrixBase< JacobianOut_t > &J, const AssignmentOperatorType op=SETTO)