diff --git a/include/hpp/fcl/BV/AABB.h b/include/hpp/fcl/BV/AABB.h index f71d1a0ca6e08b201435bd8cf9e1d96ecf1e2dcb..522c47f42b5ea96a4c01f25b41f6cc487034c33b 100644 --- a/include/hpp/fcl/BV/AABB.h +++ b/include/hpp/fcl/BV/AABB.h @@ -62,8 +62,8 @@ public: } /// @brief Creating an AABB with two endpoints a and b - AABB(const Vec3f& a, const Vec3f&b) : min_(min(a, b)), - max_(max(a, b)) + AABB(const Vec3f& a, const Vec3f&b) : min_(a.cwiseMin(b)), + max_(a.cwiseMax(b)) { } @@ -74,8 +74,8 @@ public: } /// @brief Creating an AABB contains three points - AABB(const Vec3f& a, const Vec3f& b, const Vec3f& c) : min_(min(min(a, b), c)), - max_(max(max(a, b), c)) + AABB(const Vec3f& a, const Vec3f& b, const Vec3f& c) : min_(a.cwiseMin(b).cwiseMin(c)), + max_(a.cwiseMax(b).cwiseMax(c)) { } @@ -124,8 +124,8 @@ public: return false; } - overlap_part.min_ = max(min_, other.min_); - overlap_part.max_ = min(max_, other.max_); + overlap_part.min_ = min_.cwiseMax(other.min_); + overlap_part.max_ = max_.cwiseMin(other.max_); return true; } @@ -214,7 +214,7 @@ public: /// @brief whether two AABB are equal inline bool equal(const AABB& other) const { - return min_.equal(other.min_) && max_.equal(other.max_); + return isEqual(min_, other.min_) && isEqual(max_, other.max_); } /// @brief expand the half size of the AABB by delta, and keep the center unchanged. diff --git a/include/hpp/fcl/BV/BV_node.h b/include/hpp/fcl/BV/BV_node.h index 9d2f5cff04f4cca126103826ae744895c4d922d8..52a40baa4567c98e55354c3dc29731965ddff256 100644 --- a/include/hpp/fcl/BV/BV_node.h +++ b/include/hpp/fcl/BV/BV_node.h @@ -105,31 +105,31 @@ struct BVNode : public BVNodeBase Vec3f getCenter() const { return bv.center(); } /// @brief Access the orientation of the BV - Matrix3f getOrientation() const { return Matrix3f::getIdentity(); } + Matrix3f getOrientation() const { return Matrix3f::Identity(); } }; template<> inline Matrix3f BVNode<OBB>::getOrientation() const { - return Matrix3f(bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], + return (Matrix3f() << bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], - bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]); + bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]).finished(); } template<> inline Matrix3f BVNode<RSS>::getOrientation() const { - return Matrix3f(bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], + return (Matrix3f() << bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], - bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]); + bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]).finished(); } template<> inline Matrix3f BVNode<OBBRSS>::getOrientation() const { - return Matrix3f(bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], + return (Matrix3f() << bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], bv.obb.axis[0][1], bv.obb.axis[1][1], bv.obb.axis[2][1], - bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]); + bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]).finished(); } diff --git a/include/hpp/fcl/BVH/BVH_model.h b/include/hpp/fcl/BVH/BVH_model.h index cb85a70aee63466200c38eb55c2e28a8eeefc8d6..5794850c94a6c89af55b49aa7a2cfbf3c879359d 100644 --- a/include/hpp/fcl/BVH/BVH_model.h +++ b/include/hpp/fcl/BVH/BVH_model.h @@ -222,13 +222,12 @@ public: Matrix3f computeMomentofInertia() const { - Matrix3f C(0, 0, 0, - 0, 0, 0, - 0, 0, 0); + Matrix3f C = Matrix3f::Zero(); - Matrix3f C_canonical(1/60.0, 1/120.0, 1/120.0, - 1/120.0, 1/60.0, 1/120.0, - 1/120.0, 1/120.0, 1/60.0); + Matrix3f C_canonical; + C_canonical << 1/60.0, 1/120.0, 1/120.0, + 1/120.0, 1/60.0, 1/120.0, + 1/120.0, 1/120.0, 1/60.0; for(int i = 0; i < num_tris; ++i) { @@ -236,15 +235,11 @@ public: const Vec3f& v1 = vertices[tri[0]]; const Vec3f& v2 = vertices[tri[1]]; const Vec3f& v3 = vertices[tri[2]]; - Matrix3f A(v1, v2, v3); + Matrix3f A; A << v1.transpose(), v2.transpose(), v3.transpose(); C += A.derived().transpose() * C_canonical * A * (v1.cross(v2)).dot(v3); } - FCL_REAL trace_C = C(0, 0) + C(1, 1) + C(2, 2); - - return Matrix3f(trace_C - C(0, 0), -C(0, 1), -C(0, 2), - -C(1, 0), trace_C - C(1, 1), -C(1, 2), - -C(2, 0), -C(2, 1), trace_C - C(2, 2)); + return C.trace() * Matrix3f::Identity() - C; } public: diff --git a/include/hpp/fcl/BVH/BV_splitter.h b/include/hpp/fcl/BVH/BV_splitter.h index c10670fdbd0d45f74e826d321f4f72a8f72ba3b6..c8178c885f6ac3ab1999b3971f9ebad862c749c2 100644 --- a/include/hpp/fcl/BVH/BV_splitter.h +++ b/include/hpp/fcl/BVH/BV_splitter.h @@ -76,7 +76,7 @@ class BVSplitter : public BVSplitterBase<BV> { public: - BVSplitter(SplitMethodType method) : split_method(method) + BVSplitter(SplitMethodType method) : split_vector(0,0,0), split_method(method) { } diff --git a/include/hpp/fcl/broadphase/broadphase_SaP.h b/include/hpp/fcl/broadphase/broadphase_SaP.h index 21b07bb6618a70cb2f7377e829af6507a86ecdf4..b0b230766bab0a03a81a0c742c3ddf41f4f7028a 100644 --- a/include/hpp/fcl/broadphase/broadphase_SaP.h +++ b/include/hpp/fcl/broadphase/broadphase_SaP.h @@ -164,13 +164,13 @@ protected: else return aabb->cached.min_; } - inline Vec3f::U getVal(size_t i) const + inline Vec3f::Scalar getVal(size_t i) const { if(minmax) return aabb->cached.max_[i]; else return aabb->cached.min_[i]; } - inline Vec3f::U& getVal(size_t i) + inline Vec3f::Scalar& getVal(size_t i) { if(minmax) return aabb->cached.max_[i]; else return aabb->cached.min_[i]; diff --git a/include/hpp/fcl/ccd/motion.h b/include/hpp/fcl/ccd/motion.h index c9ddc3881632d59b926f818950a6ba7965c69a3c..6277abb6dcb018a6652ebfb02f460e820f8313a0 100644 --- a/include/hpp/fcl/ccd/motion.h +++ b/include/hpp/fcl/ccd/motion.h @@ -166,7 +166,6 @@ public: } // set tm - Matrix3f I(1, 0, 0, 0, 1, 0, 0, 0, 1); // R(t) = R(t0) + R'(t0) (t-t0) + 1/2 R''(t0)(t-t0)^2 + 1 / 6 R'''(t0) (t-t0)^3 + 1 / 24 R''''(l)(t-t0)^4; t0 = 0.5 /// 1. compute M(1/2) Vec3f Rt0 = (Rd[0] + Rd[1] * 23 + Rd[2] * 23 + Rd[3]) * (1 / 48.0); @@ -182,7 +181,7 @@ public: Matrix3f hatWt0; hat(hatWt0, Wt0); Matrix3f hatWt0_sqr = hatWt0 * hatWt0; - Matrix3f Mt0 = I + hatWt0 * sintheta0 + hatWt0_sqr * (1 - costheta0); + Matrix3f Mt0 = Matrix3f::Identity() + hatWt0 * sintheta0 + hatWt0_sqr * (1 - costheta0); /// 2. compute M'(1/2) Vec3f dRt0 = (-Rd[0] - Rd[1] * 5 + Rd[2] * 5 + Rd[3]) * (1 / 8.0); @@ -340,7 +339,7 @@ public: TaylorModel sin_model(getTimeInterval()); generateTaylorModelForSinFunc(sin_model, angular_vel, 0); - TMatrix3 delta_R = hat_axis * sin_model - hat_axis * hat_axis * (cos_model - 1) + Matrix3f(1, 0, 0, 0, 1, 0, 0, 0, 1); + TMatrix3 delta_R = hat_axis * sin_model - hat_axis * hat_axis * (cos_model - 1) + Matrix3f::Identity(); TaylorModel a(getTimeInterval()), b(getTimeInterval()), c(getTimeInterval()); generateTaylorModelForLinearFunc(a, 0, linear_vel * axis[0]); @@ -494,7 +493,7 @@ public: TaylorModel sin_model(getTimeInterval()); generateTaylorModelForSinFunc(sin_model, angular_vel, 0); - TMatrix3 delta_R = hat_angular_axis * sin_model - hat_angular_axis * hat_angular_axis * (cos_model - 1) + Matrix3f(1, 0, 0, 0, 1, 0, 0, 0, 1); + TMatrix3 delta_R = hat_angular_axis * sin_model - hat_angular_axis * hat_angular_axis * (cos_model - 1) + Matrix3f::Identity(); TaylorModel a(getTimeInterval()), b(getTimeInterval()), c(getTimeInterval()); generateTaylorModelForLinearFunc(a, 0, linear_vel[0]); diff --git a/include/hpp/fcl/collision_object.h b/include/hpp/fcl/collision_object.h index fb140008088f4aae81f26cb2760311478dc4ce8e..ec9ecdce999484cd39ff19f0c539d5b13dc407f1 100644 --- a/include/hpp/fcl/collision_object.h +++ b/include/hpp/fcl/collision_object.h @@ -134,15 +134,15 @@ public: Vec3f com = computeCOM(); FCL_REAL V = computeVolume(); - return Matrix3f(C(0, 0) - V * (com[1] * com[1] + com[2] * com[2]), - C(0, 1) + V * com[0] * com[1], - C(0, 2) + V * com[0] * com[2], - C(1, 0) + V * com[1] * com[0], - C(1, 1) - V * (com[0] * com[0] + com[2] * com[2]), - C(1, 2) + V * com[1] * com[2], - C(2, 0) + V * com[2] * com[0], - C(2, 1) + V * com[2] * com[1], - C(2, 2) - V * (com[0] * com[0] + com[1] * com[1])); + return (Matrix3f() << C(0, 0) - V * (com[1] * com[1] + com[2] * com[2]), + C(0, 1) + V * com[0] * com[1], + C(0, 2) + V * com[0] * com[2], + C(1, 0) + V * com[1] * com[0], + C(1, 1) - V * (com[0] * com[0] + com[2] * com[2]), + C(1, 2) + V * com[1] * com[2], + C(2, 0) + V * com[2] * com[0], + C(2, 1) + V * com[2] * com[1], + C(2, 2) - V * (com[0] * com[0] + com[1] * com[1])).finished(); } }; @@ -207,7 +207,7 @@ public: else { Vec3f center = t.transform(cgeom->aabb_center); - Vec3f delta(cgeom->aabb_radius); + Vec3f delta(Vec3f::Constant(cgeom->aabb_radius)); aabb.min_ = center - delta; aabb.max_ = center + delta; } diff --git a/include/hpp/fcl/eigen/plugins/ccd/interval-vector.h b/include/hpp/fcl/eigen/plugins/ccd/interval-vector.h index 1675c87653042995fa45de61594cebdbda68dfb4..9ec6217fef914991214ae3157237b0d8da40e228 100644 --- a/include/hpp/fcl/eigen/plugins/ccd/interval-vector.h +++ b/include/hpp/fcl/eigen/plugins/ccd/interval-vector.h @@ -1,11 +1,4 @@ template <typename Derived> -IVector3& operator=(const FclType<Derived>& other) -{ - const Vec3f& tmp = other.fcl(); - setValue (tmp); - return *this; -} -template <typename Derived> IVector3& operator=(const Eigen::MatrixBase<Derived>& other) { const Vec3f& tmp (other); diff --git a/include/hpp/fcl/eigen/taylor_operator.h b/include/hpp/fcl/eigen/taylor_operator.h index fa02f2b7aea2b7454b255e0c0cb93e35cf673a2c..1de50b3173aa4531a95677a6dc1a6b36a6db3cc0 100644 --- a/include/hpp/fcl/eigen/taylor_operator.h +++ b/include/hpp/fcl/eigen/taylor_operator.h @@ -14,9 +14,9 @@ template<> struct TaylorReturnType<1> { typedef TVector3 type; typedef Vec3f template<> struct TaylorReturnType<3> { typedef TMatrix3 type; typedef Matrix3f eigen_type; }; template<typename Derived> -typename TaylorReturnType<Derived::ColsAtCompileTime>::type operator * (const FclType<Derived>& v, const TaylorModel& a) +typename TaylorReturnType<Derived::ColsAtCompileTime>::type operator * (const Eigen::MatrixBase<Derived>& v, const TaylorModel& a) { - const typename TaylorReturnType<Derived::ColsAtCompileTime>::eigen_type b = v.fcl(); + const typename TaylorReturnType<Derived::ColsAtCompileTime>::eigen_type b = v.derived(); return b * a; } diff --git a/include/hpp/fcl/math/matrix_3f.h b/include/hpp/fcl/math/matrix_3f.h index b04d816b96a9360eeac190d9b2b2f295d181b6c9..5226b53774cf4ecc85e8b1fa95a67fb890b53092 100644 --- a/include/hpp/fcl/math/matrix_3f.h +++ b/include/hpp/fcl/math/matrix_3f.h @@ -48,13 +48,14 @@ namespace fcl { #if FCL_HAVE_EIGEN - typedef Eigen::FclMatrix<FCL_REAL, 3> Matrix3f; + typedef Eigen::Matrix<FCL_REAL, 3, 3> Matrix3f; #elif FCL_HAVE_SSE typedef Matrix3fX<details::sse_meta_f12> Matrix3f; #else typedef Matrix3fX<details::Matrix3Data<FCL_REAL> > Matrix3f; #endif +#if ! FCL_HAVE_EIGEN static inline std::ostream& operator << (std::ostream& o, const Matrix3f& m) { o << "[(" << m(0, 0) << " " << m(0, 1) << " " << m(0, 2) << ")(" @@ -62,6 +63,7 @@ static inline std::ostream& operator << (std::ostream& o, const Matrix3f& m) << m(2, 0) << " " << m(2, 1) << " " << m(2, 2) << ")]"; return o; } +#endif @@ -73,7 +75,7 @@ public: Matrix3f Sigma; /// @brief Variations along the eign axes - Matrix3f::U sigma[3]; + Matrix3f::Scalar sigma[3]; /// @brief Eigen axes of the variation matrix Vec3f axis[3]; diff --git a/include/hpp/fcl/math/tools.h b/include/hpp/fcl/math/tools.h new file mode 100644 index 0000000000000000000000000000000000000000..77cbfe622b37cf87605304579ef0fe8d20be0fbc --- /dev/null +++ b/include/hpp/fcl/math/tools.h @@ -0,0 +1,268 @@ +/* + * Software License Agreement (BSD License) + * + * Copyright (c) 2011-2014, Willow Garage, Inc. + * Copyright (c) 2014-2015, Open Source Robotics Foundation + * All rights reserved. + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * * Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * * Redistributions in binary form must reproduce the above + * copyright notice, this list of conditions and the following + * disclaimer in the documentation and/or other materials provided + * with the distribution. + * * Neither the name of Open Source Robotics Foundation nor the names of its + * contributors may be used to endorse or promote products derived + * from this software without specific prior written permission. + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS + * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE + * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, + * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, + * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; + * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER + * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT + * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN + * ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE + * POSSIBILITY OF SUCH DAMAGE. + */ + +/** \author Joseph Mirabel */ + +#ifndef FCL_MATH_TOOLS_H +#define FCL_MATH_TOOLS_H + +#include <hpp/fcl/deprecated.hh> +#include <hpp/fcl/config.hh> +#include <hpp/fcl/config-fcl.hh> + +#include <Eigen/Dense> +#include <Eigen/Geometry> + +#include <cmath> +#include <iostream> +#include <limits> + +#define FCL_CCD_INTERVALVECTOR_PLUGIN <hpp/fcl/eigen/plugins/ccd/interval-vector.h> +#define FCL_CCD_MATRIXVECTOR_PLUGIN <hpp/fcl/eigen/plugins/ccd/interval-matrix.h> + +namespace fcl +{ + +template<typename Derived> +static inline typename Derived::Scalar triple( + const Eigen::MatrixBase<Derived>& x, + const Eigen::MatrixBase<Derived>& y, + const Eigen::MatrixBase<Derived>& z) +{ + return x.derived().dot(y.derived().cross(z.derived())); +} + + +/* +template<typename Derived, typename OtherDerived> +static inline const typename Eigen::BinaryReturnType<const Derived, const OtherDerived>::Min + min(const Eigen::MatrixBase<Derived>& x, const Eigen::MatrixBase<OtherDerived>& y) +{ + return typename Eigen::BinaryReturnType<const Derived, const OtherDerived>::Min (x.derived(), y.derived()); +} + +template<typename Derived, typename OtherDerived> +static inline const typename Eigen::BinaryReturnType<const Derived, const OtherDerived>::Max + max(const Eigen::MatrixBase<Derived>& x, const Eigen::MatrixBase<OtherDerived>& y) +{ + return typename Eigen::BinaryReturnType<const Derived, const OtherDerived>::Max (x.derived(), y.derived()); +} + +template<typename Derived> +static inline const typename Eigen::UnaryReturnType<const Derived>::Abs +abs(const Eigen::MatrixBase<Derived>& x) +{ + return typename Eigen::UnaryReturnType<const Derived>::Abs (x.derived()); +} */ + +template<typename Derived1, typename Derived2, typename Derived3> +void generateCoordinateSystem( + const Eigen::MatrixBase<Derived1>& _w, + const Eigen::MatrixBase<Derived2>& _u, + const Eigen::MatrixBase<Derived3>& _v) +{ + typedef typename Derived1::Scalar T; + + Eigen::MatrixBase<Derived1>& w = const_cast < Eigen::MatrixBase<Derived1>& > (_w); + Eigen::MatrixBase<Derived2>& u = const_cast < Eigen::MatrixBase<Derived2>& > (_u); + Eigen::MatrixBase<Derived3>& v = const_cast < Eigen::MatrixBase<Derived3>& > (_v); + + T inv_length; + if(std::abs(w[0]) >= std::abs(w[1])) + { + inv_length = (T)1.0 / sqrt(w[0] * w[0] + w[2] * w[2]); + u[0] = -w[2] * inv_length; + u[1] = (T)0; + u[2] = w[0] * inv_length; + v[0] = w[1] * u[2]; + v[1] = w[2] * u[0] - w[0] * u[2]; + v[2] = -w[1] * u[0]; + } + else + { + inv_length = (T)1.0 / sqrt(w[1] * w[1] + w[2] * w[2]); + u[0] = (T)0; + u[1] = w[2] * inv_length; + u[2] = -w[1] * inv_length; + v[0] = w[1] * u[2] - w[2] * u[1]; + v[1] = -w[0] * u[2]; + v[2] = w[0] * u[1]; + } +} + +/* ----- Start Matrices ------ */ +template<typename Derived, typename OtherDerived> +void hat(const Eigen::MatrixBase<Derived>& mat, const Eigen::MatrixBase<OtherDerived>& vec) HPP_FCL_DEPRECATED; + +template<typename Derived, typename OtherDerived> +void hat(const Eigen::MatrixBase<Derived>& mat, const Eigen::MatrixBase<OtherDerived>& vec) +{ + EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(Derived, 3, 3); + EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(OtherDerived, 3, 1); + const_cast< Eigen::MatrixBase<Derived>& >(mat) << 0, -vec[2], vec[1], vec[2], 0, -vec[0], -vec[1], vec[0], 0; +} + +template<typename Derived, typename OtherDerived> +void relativeTransform(const Eigen::MatrixBase<Derived>& R1, const Eigen::MatrixBase<OtherDerived>& t1, + const Eigen::MatrixBase<Derived>& R2, const Eigen::MatrixBase<OtherDerived>& t2, + const Eigen::MatrixBase<Derived>& R , const Eigen::MatrixBase<OtherDerived>& t) +{ + const_cast< Eigen::MatrixBase<Derived>& >(R) = R1.transpose() * R2; + const_cast< Eigen::MatrixBase<OtherDerived>& >(t) = R1.transpose() * (t2 - t1); +} + +/// @brief compute the eigen vector and eigen vector of a matrix. dout is the eigen values, vout is the eigen vectors +template<typename Derived, typename Vector> +void eigen(const Eigen::MatrixBase<Derived>& m, typename Derived::Scalar dout[3], Vector* vout) +{ + typedef typename Derived::Scalar Scalar; + Derived R(m.derived()); + int n = 3; + int j, iq, ip, i; + Scalar tresh, theta, tau, t, sm, s, h, g, c; + int nrot; + Scalar b[3]; + Scalar z[3]; + Scalar v[3][3] = {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}; + Scalar d[3]; + + for(ip = 0; ip < n; ++ip) + { + b[ip] = d[ip] = R(ip, ip); + z[ip] = 0; + } + + nrot = 0; + + for(i = 0; i < 50; ++i) + { + sm = 0; + for(ip = 0; ip < n; ++ip) + for(iq = ip + 1; iq < n; ++iq) + sm += std::abs(R(ip, iq)); + if(sm == 0.0) + { + vout[0] << v[0][0], v[0][1], v[0][2]; + vout[1] << v[1][0], v[1][1], v[1][2]; + vout[2] << v[2][0], v[2][1], v[2][2]; + dout[0] = d[0]; dout[1] = d[1]; dout[2] = d[2]; + return; + } + + if(i < 3) tresh = 0.2 * sm / (n * n); + else tresh = 0.0; + + for(ip = 0; ip < n; ++ip) + { + for(iq= ip + 1; iq < n; ++iq) + { + g = 100.0 * std::abs(R(ip, iq)); + if(i > 3 && + std::abs(d[ip]) + g == std::abs(d[ip]) && + std::abs(d[iq]) + g == std::abs(d[iq])) + R(ip, iq) = 0.0; + else if(std::abs(R(ip, iq)) > tresh) + { + h = d[iq] - d[ip]; + if(std::abs(h) + g == std::abs(h)) t = (R(ip, iq)) / h; + else + { + theta = 0.5 * h / (R(ip, iq)); + t = 1.0 /(std::abs(theta) + std::sqrt(1.0 + theta * theta)); + if(theta < 0.0) t = -t; + } + c = 1.0 / std::sqrt(1 + t * t); + s = t * c; + tau = s / (1.0 + c); + h = t * R(ip, iq); + z[ip] -= h; + z[iq] += h; + d[ip] -= h; + d[iq] += h; + R(ip, iq) = 0.0; + for(j = 0; j < ip; ++j) { g = R(j, ip); h = R(j, iq); R(j, ip) = g - s * (h + g * tau); R(j, iq) = h + s * (g - h * tau); } + for(j = ip + 1; j < iq; ++j) { g = R(ip, j); h = R(j, iq); R(ip, j) = g - s * (h + g * tau); R(j, iq) = h + s * (g - h * tau); } + for(j = iq + 1; j < n; ++j) { g = R(ip, j); h = R(iq, j); R(ip, j) = g - s * (h + g * tau); R(iq, j) = h + s * (g - h * tau); } + for(j = 0; j < n; ++j) { g = v[j][ip]; h = v[j][iq]; v[j][ip] = g - s * (h + g * tau); v[j][iq] = h + s * (g - h * tau); } + nrot++; + } + } + } + for(ip = 0; ip < n; ++ip) + { + b[ip] += z[ip]; + d[ip] = b[ip]; + z[ip] = 0.0; + } + } + + std::cerr << "eigen: too many iterations in Jacobi transform." << std::endl; + + return; +} + +template<typename Derived, typename OtherDerived> +typename Derived::Scalar quadraticForm(const Eigen::MatrixBase<Derived>& R, const Eigen::MatrixBase<OtherDerived>& v) +{ + return v.dot(R * v); +} + +template<typename Derived, typename OtherDerived> +bool isEqual(const Eigen::MatrixBase<Derived>& lhs, const Eigen::MatrixBase<OtherDerived>& rhs, const FCL_REAL tol = std::numeric_limits<FCL_REAL>::epsilon()*100) +{ + return ((lhs - rhs).array().abs() < tol).all(); +} + +template <typename Derived> +inline Derived& setEulerZYX(const Eigen::MatrixBase<Derived>& R, FCL_REAL eulerX, FCL_REAL eulerY, FCL_REAL eulerZ) +{ + const_cast< Eigen::MatrixBase<Derived>& >(R).noalias() = ( + Eigen::AngleAxisd (eulerZ, Eigen::Vector3d::UnitZ()) * + Eigen::AngleAxisd (eulerY, Eigen::Vector3d::UnitY()) * + Eigen::AngleAxisd (eulerX, Eigen::Vector3d::UnitX()) + ).toRotationMatrix(); + return const_cast< Eigen::MatrixBase<Derived>& >(R).derived(); +} + +template <typename Derived> +inline Derived& setEulerYPR(const Eigen::MatrixBase<Derived>& R, FCL_REAL yaw, FCL_REAL pitch, FCL_REAL roll) +{ + return setEulerZYX(R, roll, pitch, yaw); +} + +} + + +#endif diff --git a/include/hpp/fcl/math/vec_3f.h b/include/hpp/fcl/math/vec_3f.h index f81e6da128d6ed07db88736b87c7e5d817e71d60..1ecc46b184afbb92e26697190ad0d4ff714d8b79 100644 --- a/include/hpp/fcl/math/vec_3f.h +++ b/include/hpp/fcl/math/vec_3f.h @@ -42,7 +42,9 @@ #include <hpp/fcl/data_types.h> #if FCL_HAVE_EIGEN -# include <hpp/fcl/eigen/math_details.h> +#include <Eigen/Core> +#include <Eigen/Geometry> +#include <hpp/fcl/math/tools.h> #elif FCL_HAVE_SSE # include <hpp/fcl/simd/math_simd_details.h> #else @@ -50,7 +52,6 @@ #endif #if FCL_HAVE_EIGEN -# include <hpp/fcl/eigen/vec_3fx.h> #else # include <hpp/fcl/math/vec_3fx.h> #endif @@ -64,7 +65,7 @@ namespace fcl { #if FCL_HAVE_EIGEN - typedef Eigen::FclMatrix<FCL_REAL, 1> Vec3f; + typedef Eigen::Matrix<FCL_REAL, 3, 1> Vec3f; #elif FCL_HAVE_SSE typedef Vec3fX<details::sse_meta_f4> Vec3f; #else diff --git a/include/hpp/fcl/narrowphase/gjk.h b/include/hpp/fcl/narrowphase/gjk.h index bd5924085d4306360bcb3eb21979c15571eac64f..8ad919aeda95db8175191483e223db64e48af871 100644 --- a/include/hpp/fcl/narrowphase/gjk.h +++ b/include/hpp/fcl/narrowphase/gjk.h @@ -130,6 +130,8 @@ struct GJK Vec3f d; /// @brieg support vector (i.e., the furthest point on the shape along the support direction) Vec3f w; + + SimplexV () : d(Vec3f::Zero()), w(Vec3f::Zero()) {} }; struct Simplex @@ -220,6 +222,8 @@ private: SimplexF* l[2]; // the pre and post faces in the list size_t e[3]; size_t pass; + + SimplexF () : n(Vec3f::Zero()) {}; }; struct SimplexList diff --git a/include/hpp/fcl/narrowphase/narrowphase.h b/include/hpp/fcl/narrowphase/narrowphase.h index 482afc76098e186d5107dff1e948aaf18302ecf5..2929408d2853279c9d38eff1f4839423c23757d4 100644 --- a/include/hpp/fcl/narrowphase/narrowphase.h +++ b/include/hpp/fcl/narrowphase/narrowphase.h @@ -453,7 +453,7 @@ struct GJKSolver_indep details::MinkowskiDiff shape; shape.shapes[0] = &s1; shape.shapes[1] = &s2; - shape.toshape1 = tf2.getRotation().transposeTimes(tf1.getRotation()); + shape.toshape1 = tf2.getRotation().transpose() * tf1.getRotation(); shape.toshape0 = tf1.inverseTimes(tf2); details::GJK gjk(gjk_max_iterations, gjk_tolerance); @@ -468,7 +468,7 @@ struct GJKSolver_indep details::EPA::Status epa_status = epa.evaluate(gjk, -guess); if(epa_status != details::EPA::Failed) { - Vec3f w0; + Vec3f w0 (Vec3f::Zero()); for(size_t i = 0; i < epa.result.rank; ++i) { w0 += shape.support(epa.result.c[i]->d, 0) * epa.result.p[i]; @@ -517,7 +517,7 @@ struct GJKSolver_indep details::EPA::Status epa_status = epa.evaluate(gjk, -guess); if(epa_status != details::EPA::Failed) { - Vec3f w0; + Vec3f w0 (Vec3f::Zero()); for(size_t i = 0; i < epa.result.rank; ++i) { w0 += shape.support(epa.result.c[i]->d, 0) * epa.result.p[i]; @@ -551,7 +551,7 @@ struct GJKSolver_indep details::MinkowskiDiff shape; shape.shapes[0] = &s; shape.shapes[1] = &tri; - shape.toshape1 = tf2.getRotation().transposeTimes(tf1.getRotation()); + shape.toshape1 = tf2.getRotation().transpose() * tf1.getRotation(); shape.toshape0 = tf1.inverseTimes(tf2); details::GJK gjk(gjk_max_iterations, gjk_tolerance); @@ -566,7 +566,7 @@ struct GJKSolver_indep details::EPA::Status epa_status = epa.evaluate(gjk, -guess); if(epa_status != details::EPA::Failed) { - Vec3f w0; + Vec3f w0 (Vec3f::Zero()); for(size_t i = 0; i < epa.result.rank; ++i) { w0 += shape.support(epa.result.c[i]->d, 0) * epa.result.p[i]; @@ -598,7 +598,7 @@ struct GJKSolver_indep details::MinkowskiDiff shape; shape.shapes[0] = &s1; shape.shapes[1] = &s2; - shape.toshape1 = tf2.getRotation().transposeTimes(tf1.getRotation()); + shape.toshape1 = tf2.getRotation().transpose() * tf1.getRotation(); shape.toshape0 = tf1.inverseTimes(tf2); details::GJK gjk(gjk_max_iterations, gjk_tolerance); @@ -607,7 +607,7 @@ struct GJKSolver_indep if(gjk_status == details::GJK::Valid) { - Vec3f w0, w1; + Vec3f w0 (Vec3f::Zero()), w1 (Vec3f::Zero()); for(size_t i = 0; i < gjk.getSimplex()->rank; ++i) { FCL_REAL p = gjk.getSimplex()->p[i]; @@ -659,7 +659,7 @@ struct GJKSolver_indep if(gjk_status == details::GJK::Valid) { - Vec3f w0, w1; + Vec3f w0 (Vec3f::Zero()), w1 (Vec3f::Zero()); for(size_t i = 0; i < gjk.getSimplex()->rank; ++i) { FCL_REAL p = gjk.getSimplex()->p[i]; @@ -700,7 +700,7 @@ struct GJKSolver_indep details::MinkowskiDiff shape; shape.shapes[0] = &s; shape.shapes[1] = &tri; - shape.toshape1 = tf2.getRotation().transposeTimes(tf1.getRotation()); + shape.toshape1 = tf2.getRotation().transpose() * tf1.getRotation(); shape.toshape0 = tf1.inverseTimes(tf2); details::GJK gjk(gjk_max_iterations, gjk_tolerance); @@ -709,7 +709,7 @@ struct GJKSolver_indep if(gjk_status == details::GJK::Valid) { - Vec3f w0, w1; + Vec3f w0 (Vec3f::Zero()), w1 (Vec3f::Zero()); for(size_t i = 0; i < gjk.getSimplex()->rank; ++i) { FCL_REAL p = gjk.getSimplex()->p[i]; diff --git a/include/hpp/fcl/shape/geometric_shapes.h b/include/hpp/fcl/shape/geometric_shapes.h index 9c76ec11b5ad22e7450b2de916ed34bfb2a7553e..7703561aa6908b4624b5cfb03fa4e84318e446d1 100644 --- a/include/hpp/fcl/shape/geometric_shapes.h +++ b/include/hpp/fcl/shape/geometric_shapes.h @@ -106,9 +106,9 @@ public: FCL_REAL a2 = side[0] * side[0] * V; FCL_REAL b2 = side[1] * side[1] * V; FCL_REAL c2 = side[2] * side[2] * V; - return Matrix3f((b2 + c2) / 12, 0, 0, - 0, (a2 + c2) / 12, 0, - 0, 0, (a2 + b2) / 12); + return (Matrix3f() << (b2 + c2) / 12, 0, 0, + 0, (a2 + c2) / 12, 0, + 0, 0, (a2 + b2) / 12).finished(); } }; @@ -132,9 +132,7 @@ public: Matrix3f computeMomentofInertia() const { FCL_REAL I = 0.4 * radius * radius * computeVolume(); - return Matrix3f(I, 0, 0, - 0, I, 0, - 0, 0, I); + return I * Matrix3f::Identity(); } FCL_REAL computeVolume() const @@ -176,9 +174,9 @@ public: FCL_REAL ix = v_cyl * lz * lz / 12.0 + 0.25 * v_cyl * radius + 0.4 * v_sph * radius * radius + 0.25 * v_sph * lz * lz; FCL_REAL iz = (0.5 * v_cyl + 0.4 * v_sph) * radius * radius; - return Matrix3f(ix, 0, 0, - 0, ix, 0, - 0, 0, iz); + return (Matrix3f() << ix, 0, 0, + 0, ix, 0, + 0, 0, iz).finished(); } }; @@ -214,9 +212,9 @@ public: FCL_REAL ix = V * (0.1 * lz * lz + 3 * radius * radius / 20); FCL_REAL iz = 0.3 * V * radius * radius; - return Matrix3f(ix, 0, 0, - 0, ix, 0, - 0, 0, iz); + return (Matrix3f() << ix, 0, 0, + 0, ix, 0, + 0, 0, iz).finished(); } Vec3f computeCOM() const @@ -256,9 +254,9 @@ public: FCL_REAL V = computeVolume(); FCL_REAL ix = V * (3 * radius * radius + lz * lz) / 12; FCL_REAL iz = V * radius * radius / 2; - return Matrix3f(ix, 0, 0, - 0, ix, 0, - 0, 0, iz); + return (Matrix3f() << ix, 0, 0, + 0, ix, 0, + 0, 0, iz).finished(); } }; @@ -342,13 +340,12 @@ public: Matrix3f computeMomentofInertia() const { - Matrix3f C(0, 0, 0, - 0, 0, 0, - 0, 0, 0); + Matrix3f C = Matrix3f::Zero(); - Matrix3f C_canonical(1/60.0, 1/120.0, 1/120.0, - 1/120.0, 1/60.0, 1/120.0, - 1/120.0, 1/120.0, 1/60.0); + Matrix3f C_canonical; + C_canonical << 1/60.0, 1/120.0, 1/120.0, + 1/120.0, 1/60.0, 1/120.0, + 1/120.0, 1/120.0, 1/60.0; int* points_in_poly = polygons; int* index = polygons + 1; @@ -369,7 +366,7 @@ public: int e_second = index[(j+1)%*points_in_poly]; const Vec3f& v1 = points[e_first]; const Vec3f& v2 = points[e_second]; - Matrix3f A(v1, v2, v3); // this is A' in the original document + Matrix3f A; A << v1.transpose(), v2.transpose(), v3.transpose(); // this is A' in the original document C += A.derived().transpose() * C_canonical * A * (v1.cross(v2)).dot(v3); } @@ -377,12 +374,7 @@ public: index = points_in_poly + 1; } - FCL_REAL trace_C = C(0, 0) + C(1, 1) + C(2, 2); - - return Matrix3f(trace_C - C(0, 0), -C(0, 1), -C(0, 2), - -C(1, 0), trace_C - C(1, 1), -C(1, 2), - -C(2, 0), -C(2, 1), trace_C - C(2, 2)); - + return C.trace() * Matrix3f::Identity() - C; } Vec3f computeCOM() const diff --git a/src/BV/AABB.cpp b/src/BV/AABB.cpp index 0f66b45915cf78239eabf3a7ab8cefa52560912c..82f71598a5416e6b42a8dbf9d07eee5b3ff49ce0 100644 --- a/src/BV/AABB.cpp +++ b/src/BV/AABB.cpp @@ -43,8 +43,8 @@ namespace fcl { -AABB::AABB() : min_(std::numeric_limits<FCL_REAL>::max()), - max_(-std::numeric_limits<FCL_REAL>::max()) +AABB::AABB() : min_(Vec3f::Constant(std::numeric_limits<FCL_REAL>::max())), + max_(Vec3f::Constant(-std::numeric_limits<FCL_REAL>::max())) { } diff --git a/src/BV/OBB.cpp b/src/BV/OBB.cpp index d69b4cf6147f092c6d266b30313e9aa4b2eccc79..bf045f803d1ca176e3bf95c6a923373401f7f96f 100644 --- a/src/BV/OBB.cpp +++ b/src/BV/OBB.cpp @@ -75,7 +75,7 @@ inline OBB merge_largedist(const OBB& b1, const OBB& b2) computeVertices(b2, vertex + 8); Matrix3f M; Vec3f E[3]; - Matrix3f::U s[3] = {0, 0, 0}; + Matrix3f::Scalar s[3] = {0, 0, 0}; // obb axes Vec3f& R0 = b.axis[0]; @@ -178,48 +178,57 @@ bool obbDisjoint(const Matrix3f& B, const Vec3f& T, const Vec3f& a, const Vec3f& register FCL_REAL t, s; const FCL_REAL reps = 1e-6; - Matrix3f Bf = abs(B); - Bf += reps; + Matrix3f Bf (B.array().abs() + reps); + // Bf += reps; // if any of these tests are one-sided, then the polyhedra are disjoint // A1 x A2 = A0 t = ((T[0] < 0.0) ? -T[0] : T[0]); - if(t > (a[0] + Bf.dotX(b))) + // if(t > (a[0] + Bf.dotX(b))) + if(t > (a[0] + Bf.row(0).dot(b))) return true; // B1 x B2 = B0 - s = B.transposeDotX(T); + // s = B.transposeDotX(T); + s = B.col(0).dot(T); t = ((s < 0.0) ? -s : s); - if(t > (b[0] + Bf.transposeDotX(a))) + // if(t > (b[0] + Bf.transposeDotX(a))) + if(t > (b[0] + Bf.col(0).dot(a))) return true; // A2 x A0 = A1 t = ((T[1] < 0.0) ? -T[1] : T[1]); - if(t > (a[1] + Bf.dotY(b))) + // if(t > (a[1] + Bf.dotY(b))) + if(t > (a[1] + Bf.row(1).dot(b))) return true; // A0 x A1 = A2 t =((T[2] < 0.0) ? -T[2] : T[2]); - if(t > (a[2] + Bf.dotZ(b))) + // if(t > (a[2] + Bf.dotZ(b))) + if(t > (a[2] + Bf.row(2).dot(b))) return true; // B2 x B0 = B1 - s = B.transposeDotY(T); + // s = B.transposeDotY(T); + s = B.col(1).dot(T); t = ((s < 0.0) ? -s : s); - if(t > (b[1] + Bf.transposeDotY(a))) + // if(t > (b[1] + Bf.transposeDotY(a))) + if(t > (b[1] + Bf.col(1).dot(a))) return true; // B0 x B1 = B2 - s = B.transposeDotZ(T); + // s = B.transposeDotZ(T); + s = B.col(2).dot(T); t = ((s < 0.0) ? -s : s); - if(t > (b[2] + Bf.transposeDotZ(a))) + // if(t > (b[2] + Bf.transposeDotZ(a))) + if(t > (b[2] + Bf.col(2).dot(a))) return true; // A0 x B0 @@ -312,8 +321,9 @@ bool obbDisjointAndLowerBoundDistance (const Matrix3f& B, const Vec3f& T, b [0] + b [1] + b [2]); FCL_REAL breakDistance2 = breakDistance * breakDistance; - Matrix3f Bf = abs(B); - Bf += reps; + // Matrix3f Bf = abs(B); + // Bf += reps; + Matrix3f Bf (B.array().abs() + reps); squaredLowerBoundDistance = 0; // if any of these tests are one-sided, then the polyhedra are disjoint @@ -321,7 +331,7 @@ bool obbDisjointAndLowerBoundDistance (const Matrix3f& B, const Vec3f& T, // A1 x A2 = A0 t = ((T[0] < 0.0) ? -T[0] : T[0]); - diff = t - (a[0] + Bf.dotX(b)); + diff = t - (a[0] + Bf.row(0).dot(b)); if (diff > 0) { squaredLowerBoundDistance += diff*diff; } @@ -329,7 +339,7 @@ bool obbDisjointAndLowerBoundDistance (const Matrix3f& B, const Vec3f& T, // A2 x A0 = A1 t = ((T[1] < 0.0) ? -T[1] : T[1]); - diff = t - (a[1] + Bf.dotY(b)); + diff = t - (a[1] + Bf.row(1).dot(b)); if (diff > 0) { squaredLowerBoundDistance += diff*diff; } @@ -337,7 +347,7 @@ bool obbDisjointAndLowerBoundDistance (const Matrix3f& B, const Vec3f& T, // A0 x A1 = A2 t =((T[2] < 0.0) ? -T[2] : T[2]); - diff = t - (a[2] + Bf.dotZ(b)); + diff = t - (a[2] + Bf.row(2).dot(b)); if (diff > 0) { squaredLowerBoundDistance += diff*diff; } @@ -346,28 +356,28 @@ bool obbDisjointAndLowerBoundDistance (const Matrix3f& B, const Vec3f& T, return true; // B1 x B2 = B0 - s = B.transposeDotX(T); + s = B.col(0).dot(T); t = ((s < 0.0) ? -s : s); - diff = t - (b[0] + Bf.transposeDotX(a)); + diff = t - (b[0] + Bf.col(0).dot(a)); if (diff > 0) { squaredLowerBoundDistance += diff*diff; } // B2 x B0 = B1 - s = B.transposeDotY(T); + s = B.col(1).dot(T); t = ((s < 0.0) ? -s : s); - diff = t - (b[1] + Bf.transposeDotY(a)); + diff = t - (b[1] + Bf.col(1).dot(a)); if (diff > 0) { squaredLowerBoundDistance += diff*diff; } // B0 x B1 = B2 - s = B.transposeDotZ(T); + s = B.col(2).dot(T); t = ((s < 0.0) ? -s : s); - diff = t - (b[2] + Bf.transposeDotZ(a)); + diff = t - (b[2] + Bf.col(2).dot(a)); if (diff > 0) { squaredLowerBoundDistance += diff*diff; } @@ -537,9 +547,10 @@ bool OBB::overlap(const OBB& other) const /// First compute the rotation part, then translation part Vec3f t = other.To - To; // T2 - T1 Vec3f T(t.dot(axis[0]), t.dot(axis[1]), t.dot(axis[2])); // R1'(T2-T1) - Matrix3f R(axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), axis[0].dot(other.axis[2]), - axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), axis[1].dot(other.axis[2]), - axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), axis[2].dot(other.axis[2])); + Matrix3f R; + R << axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), axis[0].dot(other.axis[2]), + axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), axis[1].dot(other.axis[2]), + axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), axis[2].dot(other.axis[2]); return !obbDisjoint(R, T, extent, other.extent); } @@ -551,12 +562,13 @@ bool OBB::overlap(const OBB& other) const /// First compute the rotation part, then translation part Vec3f t = other.To - To; // T2 - T1 Vec3f T(t.dot(axis[0]), t.dot(axis[1]), t.dot(axis[2])); // R1'(T2-T1) - Matrix3f R(axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), - axis[0].dot(other.axis[2]), - axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), - axis[1].dot(other.axis[2]), - axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), - axis[2].dot(other.axis[2])); + Matrix3f R; + R << axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), + axis[0].dot(other.axis[2]), + axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), + axis[1].dot(other.axis[2]), + axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), + axis[2].dot(other.axis[2]); return !obbDisjointAndLowerBoundDistance (R, T, extent, other.extent, sqrDistLowerBound); @@ -619,13 +631,15 @@ FCL_REAL OBB::distance(const OBB& /*other*/, Vec3f* /*P*/, Vec3f* /*Q*/) const bool overlap(const Matrix3f& R0, const Vec3f& T0, const OBB& b1, const OBB& b2) { - Matrix3f R0b2(R0.dotX(b2.axis[0]), R0.dotX(b2.axis[1]), R0.dotX(b2.axis[2]), - R0.dotY(b2.axis[0]), R0.dotY(b2.axis[1]), R0.dotY(b2.axis[2]), - R0.dotZ(b2.axis[0]), R0.dotZ(b2.axis[1]), R0.dotZ(b2.axis[2])); + Matrix3f R0b2; + R0b2 << R0.row(0).dot(b2.axis[0]), R0.row(0).dot(b2.axis[1]), R0.row(0).dot(b2.axis[2]), + R0.row(1).dot(b2.axis[0]), R0.row(1).dot(b2.axis[1]), R0.row(1).dot(b2.axis[2]), + R0.row(2).dot(b2.axis[0]), R0.row(2).dot(b2.axis[1]), R0.row(2).dot(b2.axis[2]); - Matrix3f R(R0b2.transposeDotX(b1.axis[0]), R0b2.transposeDotY(b1.axis[0]), R0b2.transposeDotZ(b1.axis[0]), - R0b2.transposeDotX(b1.axis[1]), R0b2.transposeDotY(b1.axis[1]), R0b2.transposeDotZ(b1.axis[1]), - R0b2.transposeDotX(b1.axis[2]), R0b2.transposeDotY(b1.axis[2]), R0b2.transposeDotZ(b1.axis[2])); + Matrix3f R; + R << R0b2.col(0).dot(b1.axis[0]), R0b2.col(1).dot(b1.axis[0]), R0b2.col(2).dot(b1.axis[0]), + R0b2.col(0).dot(b1.axis[1]), R0b2.col(1).dot(b1.axis[1]), R0b2.col(2).dot(b1.axis[1]), + R0b2.col(0).dot(b1.axis[2]), R0b2.col(1).dot(b1.axis[2]), R0b2.col(2).dot(b1.axis[2]); Vec3f Ttemp = R0 * b2.To + T0 - b1.To; Vec3f T(Ttemp.dot(b1.axis[0]), Ttemp.dot(b1.axis[1]), Ttemp.dot(b1.axis[2])); @@ -636,13 +650,15 @@ bool overlap(const Matrix3f& R0, const Vec3f& T0, const OBB& b1, const OBB& b2) bool overlap(const Matrix3f& R0, const Vec3f& T0, const OBB& b1, const OBB& b2, FCL_REAL& sqrDistLowerBound) { - Matrix3f R0b2(R0.dotX(b2.axis[0]), R0.dotX(b2.axis[1]), R0.dotX(b2.axis[2]), - R0.dotY(b2.axis[0]), R0.dotY(b2.axis[1]), R0.dotY(b2.axis[2]), - R0.dotZ(b2.axis[0]), R0.dotZ(b2.axis[1]), R0.dotZ(b2.axis[2])); - - Matrix3f R(R0b2.transposeDotX(b1.axis[0]), R0b2.transposeDotY(b1.axis[0]), R0b2.transposeDotZ(b1.axis[0]), - R0b2.transposeDotX(b1.axis[1]), R0b2.transposeDotY(b1.axis[1]), R0b2.transposeDotZ(b1.axis[1]), - R0b2.transposeDotX(b1.axis[2]), R0b2.transposeDotY(b1.axis[2]), R0b2.transposeDotZ(b1.axis[2])); + Matrix3f R0b2; + R0b2 << R0.row(0).dot(b2.axis[0]), R0.row(0).dot(b2.axis[1]), R0.row(0).dot(b2.axis[2]), + R0.row(1).dot(b2.axis[0]), R0.row(1).dot(b2.axis[1]), R0.row(1).dot(b2.axis[2]), + R0.row(2).dot(b2.axis[0]), R0.row(2).dot(b2.axis[1]), R0.row(2).dot(b2.axis[2]); + + Matrix3f R; + R << R0b2.col(0).dot(b1.axis[0]), R0b2.col(1).dot(b1.axis[0]), R0b2.col(2).dot(b1.axis[0]), + R0b2.col(0).dot(b1.axis[1]), R0b2.col(1).dot(b1.axis[1]), R0b2.col(2).dot(b1.axis[1]), + R0b2.col(0).dot(b1.axis[2]), R0b2.col(1).dot(b1.axis[2]), R0b2.col(2).dot(b1.axis[2]); Vec3f Ttemp = R0 * b2.To + T0 - b1.To; Vec3f T(Ttemp.dot(b1.axis[0]), Ttemp.dot(b1.axis[1]), Ttemp.dot(b1.axis[2])); diff --git a/src/BV/RSS.cpp b/src/BV/RSS.cpp index 607d5cbc83aad11b1159df26930c83eb21bbeff1..435a1b2afe91da051ce8defdf079aeab7c3717f7 100644 --- a/src/BV/RSS.cpp +++ b/src/BV/RSS.cpp @@ -135,7 +135,7 @@ FCL_REAL rectDistance(const Matrix3f& Rab, Vec3f const& Tab, const FCL_REAL a[2] bA0_dot_B1 = b[1] * A0_dot_B1; bA1_dot_B1 = b[1] * A1_dot_B1; - Vec3f Tba = Rab.transposeTimes(Tab); + Vec3f Tba (Rab.transpose() * Tab); Vec3f S; FCL_REAL t, u; @@ -842,9 +842,10 @@ bool RSS::overlap(const RSS& other) const Vec3f T(t.dot(axis[0]), t.dot(axis[1]), t.dot(axis[2])); /// Now compute R1'R2 - Matrix3f R(axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), axis[0].dot(other.axis[2]), + Matrix3f R; + R << axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), axis[0].dot(other.axis[2]), axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), axis[1].dot(other.axis[2]), - axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), axis[2].dot(other.axis[2])); + axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), axis[2].dot(other.axis[2]); FCL_REAL dist = rectDistance(R, T, l, other.l); return (dist <= (r + other.r)); @@ -854,15 +855,17 @@ bool overlap(const Matrix3f& R0, const Vec3f& T0, const RSS& b1, const RSS& b2) { // ROb2 = R0 . b2 // where b2 = [ b2.axis [0] | b2.axis [1] | b2.axis [2] ] - Matrix3f R0b2(R0.dotX(b2.axis[0]), R0.dotX(b2.axis[1]), R0.dotX(b2.axis[2]), - R0.dotY(b2.axis[0]), R0.dotY(b2.axis[1]), R0.dotY(b2.axis[2]), - R0.dotZ(b2.axis[0]), R0.dotZ(b2.axis[1]), R0.dotZ(b2.axis[2])); + Matrix3f R0b2; + R0b2 << R0.row(0).dot(b2.axis[0]), R0.row(0).dot(b2.axis[1]), R0.row(0).dot(b2.axis[2]), + R0.row(1).dot(b2.axis[0]), R0.row(1).dot(b2.axis[1]), R0.row(1).dot(b2.axis[2]), + R0.row(2).dot(b2.axis[0]), R0.row(2).dot(b2.axis[1]), R0.row(2).dot(b2.axis[2]); // (1 0 0)^T R0b2^T axis [0] = (1 0 0)^T b2^T R0^T axis [0] // R = b2^T RO^T b1 - Matrix3f R(R0b2.transposeDotX(b1.axis[0]), R0b2.transposeDotY(b1.axis[0]), R0b2.transposeDotZ(b1.axis[0]), - R0b2.transposeDotX(b1.axis[1]), R0b2.transposeDotY(b1.axis[1]), R0b2.transposeDotZ(b1.axis[1]), - R0b2.transposeDotX(b1.axis[2]), R0b2.transposeDotY(b1.axis[2]), R0b2.transposeDotZ(b1.axis[2])); + Matrix3f R; + R << R0b2.col(0).dot(b1.axis[0]), R0b2.col(1).dot(b1.axis[0]), R0b2.col(2).dot(b1.axis[0]), + R0b2.col(0).dot(b1.axis[1]), R0b2.col(1).dot(b1.axis[1]), R0b2.col(2).dot(b1.axis[1]), + R0b2.col(0).dot(b1.axis[2]), R0b2.col(1).dot(b1.axis[2]), R0b2.col(2).dot(b1.axis[2]); Vec3f Ttemp = R0 * b2.Tr + T0 - b1.Tr; Vec3f T(Ttemp.dot(b1.axis[0]), Ttemp.dot(b1.axis[1]), Ttemp.dot(b1.axis[2])); @@ -1081,7 +1084,7 @@ RSS RSS::operator + (const RSS& other) const Matrix3f M; // row first matrix Vec3f E[3]; // row first eigen-vectors - Matrix3f::U s[3] = {0, 0, 0}; + Matrix3f::Scalar s[3] = {0, 0, 0}; getCovariance(v, NULL, NULL, NULL, 16, M); eigen(M, s, E); @@ -1113,9 +1116,10 @@ FCL_REAL RSS::distance(const RSS& other, Vec3f* P, Vec3f* Q) const // First compute the rotation part, then translation part Vec3f t = other.Tr - Tr; // T2 - T1 Vec3f T(t.dot(axis[0]), t.dot(axis[1]), t.dot(axis[2])); // R1'(T2-T1) - Matrix3f R(axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), axis[0].dot(other.axis[2]), - axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), axis[1].dot(other.axis[2]), - axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), axis[2].dot(other.axis[2])); + Matrix3f R; + R << axis[0].dot(other.axis[0]), axis[0].dot(other.axis[1]), axis[0].dot(other.axis[2]), + axis[1].dot(other.axis[0]), axis[1].dot(other.axis[1]), axis[1].dot(other.axis[2]), + axis[2].dot(other.axis[0]), axis[2].dot(other.axis[1]), axis[2].dot(other.axis[2]); FCL_REAL dist = rectDistance(R, T, l, other.l, P, Q); dist -= (r + other.r); @@ -1124,13 +1128,15 @@ FCL_REAL RSS::distance(const RSS& other, Vec3f* P, Vec3f* Q) const FCL_REAL distance(const Matrix3f& R0, const Vec3f& T0, const RSS& b1, const RSS& b2, Vec3f* P, Vec3f* Q) { - Matrix3f R0b2(R0.dotX(b2.axis[0]), R0.dotX(b2.axis[1]), R0.dotX(b2.axis[2]), - R0.dotY(b2.axis[0]), R0.dotY(b2.axis[1]), R0.dotY(b2.axis[2]), - R0.dotZ(b2.axis[0]), R0.dotZ(b2.axis[1]), R0.dotZ(b2.axis[2])); - - Matrix3f R(R0b2.transposeDotX(b1.axis[0]), R0b2.transposeDotY(b1.axis[0]), R0b2.transposeDotZ(b1.axis[0]), - R0b2.transposeDotX(b1.axis[1]), R0b2.transposeDotY(b1.axis[1]), R0b2.transposeDotZ(b1.axis[1]), - R0b2.transposeDotX(b1.axis[2]), R0b2.transposeDotY(b1.axis[2]), R0b2.transposeDotZ(b1.axis[2])); + Matrix3f R0b2; + R0b2 << R0.row(0).dot(b2.axis[0]), R0.row(0).dot(b2.axis[1]), R0.row(0).dot(b2.axis[2]), + R0.row(1).dot(b2.axis[0]), R0.row(1).dot(b2.axis[1]), R0.row(1).dot(b2.axis[2]), + R0.row(2).dot(b2.axis[0]), R0.row(2).dot(b2.axis[1]), R0.row(2).dot(b2.axis[2]); + + Matrix3f R; + R << R0b2.col(0).dot(b1.axis[0]), R0b2.col(1).dot(b1.axis[0]), R0b2.col(2).dot(b1.axis[0]), + R0b2.col(0).dot(b1.axis[1]), R0b2.col(1).dot(b1.axis[1]), R0b2.col(2).dot(b1.axis[1]), + R0b2.col(0).dot(b1.axis[2]), R0b2.col(1).dot(b1.axis[2]), R0b2.col(2).dot(b1.axis[2]); Vec3f Ttemp = R0 * b2.Tr + T0 - b1.Tr; diff --git a/src/BVH/BVH_utility.cpp b/src/BVH/BVH_utility.cpp index c748d022047a76e403f7cfd6616b22c999606f44..e06b1a6e68f97debcdb1b22698554ea96afa6908 100644 --- a/src/BVH/BVH_utility.cpp +++ b/src/BVH/BVH_utility.cpp @@ -106,8 +106,8 @@ void BVHExpand(BVHModel<RSS>& model, const Variance3f* ucs, FCL_REAL r = 1.0) void getCovariance(Vec3f* ps, Vec3f* ps2, Triangle* ts, unsigned int* indices, int n, Matrix3f& M) { - Vec3f S1; - Vec3f S2[3]; + Vec3f S1 (Vec3f::Zero()); + Vec3f S2[3] = { Vec3f::Zero(), Vec3f::Zero(), Vec3f::Zero() }; if(ts) { diff --git a/src/BVH/BV_fitter.cpp b/src/BVH/BV_fitter.cpp index 6eda96e41d5452195142976343d87d9e73eabb77..548bc4f4b0b4053238038a6736b46b5bec9cc0e8 100644 --- a/src/BVH/BV_fitter.cpp +++ b/src/BVH/BV_fitter.cpp @@ -49,7 +49,7 @@ static const double invCosA = 2.0 / sqrt(3.0); static const double sinA = 0.5; static const double cosA = sqrt(3.0) / 2.0; -static inline void axisFromEigen(Vec3f eigenV[3], Matrix3f::U eigenS[3], Vec3f axis[3]) +static inline void axisFromEigen(Vec3f eigenV[3], Matrix3f::Scalar eigenS[3], Vec3f axis[3]) { int min, mid, max; if(eigenS[0] > eigenS[1]) { max = 0; min = 1; } @@ -138,7 +138,7 @@ void fitn(Vec3f* ps, int n, OBB& bv) { Matrix3f M; Vec3f E[3]; - Matrix3f::U s[3] = {0, 0, 0}; // three eigen values + Matrix3f::Scalar s[3] = {0, 0, 0}; // three eigen values getCovariance(ps, NULL, NULL, NULL, n, M); eigen(M, s, E); @@ -224,7 +224,7 @@ void fitn(Vec3f* ps, int n, RSS& bv) { Matrix3f M; // row first matrix Vec3f E[3]; // row first eigen-vectors - Matrix3f::U s[3] = {0, 0, 0}; + Matrix3f::Scalar s[3] = {0, 0, 0}; getCovariance(ps, NULL, NULL, NULL, n, M); eigen(M, s, E); @@ -341,7 +341,7 @@ void fitn(Vec3f* ps, int n, kIOS& bv) { Matrix3f M; Vec3f E[3]; - Matrix3f::U s[3] = {0, 0, 0}; // three eigen values; + Matrix3f::Scalar s[3] = {0, 0, 0}; // three eigen values; getCovariance(ps, NULL, NULL, NULL, n, M); eigen(M, s, E); @@ -523,7 +523,7 @@ OBB BVFitter<OBB>::fit(unsigned int* primitive_indices, int num_primitives) Matrix3f M; // row first matrix Vec3f E[3]; // row first eigen-vectors - Matrix3f::U s[3]; // three eigen values + Matrix3f::Scalar s[3]; // three eigen values getCovariance(vertices, prev_vertices, tri_indices, primitive_indices, num_primitives, M); eigen(M, s, E); @@ -541,7 +541,7 @@ OBBRSS BVFitter<OBBRSS>::fit(unsigned int* primitive_indices, int num_primitives OBBRSS bv; Matrix3f M; Vec3f E[3]; - Matrix3f::U s[3]; + Matrix3f::Scalar s[3]; getCovariance(vertices, prev_vertices, tri_indices, primitive_indices, num_primitives, M); eigen(M, s, E); @@ -572,7 +572,7 @@ RSS BVFitter<RSS>::fit(unsigned int* primitive_indices, int num_primitives) Matrix3f M; // row first matrix Vec3f E[3]; // row first eigen-vectors - Matrix3f::U s[3]; // three eigen values + Matrix3f::Scalar s[3]; // three eigen values getCovariance(vertices, prev_vertices, tri_indices, primitive_indices, num_primitives, M); eigen(M, s, E); axisFromEigen(E, s, bv.axis); @@ -600,7 +600,7 @@ kIOS BVFitter<kIOS>::fit(unsigned int* primitive_indices, int num_primitives) Matrix3f M; // row first matrix Vec3f E[3]; // row first eigen-vectors - Matrix3f::U s[3]; + Matrix3f::Scalar s[3]; getCovariance(vertices, prev_vertices, tri_indices, primitive_indices, num_primitives, M); eigen(M, s, E); diff --git a/src/broadphase/broadphase_SSaP.cpp b/src/broadphase/broadphase_SSaP.cpp index a08bf5a92daf707463ae55f2eef3a4ebc4aa5d7b..4a7cee1497818907623386f466f5f6578e59dba9 100644 --- a/src/broadphase/broadphase_SSaP.cpp +++ b/src/broadphase/broadphase_SSaP.cpp @@ -353,7 +353,7 @@ bool SSaPCollisionManager::distance_(CollisionObject* obj, void* cdata, Distance } else // need more loop { - if(dummy_vector.equal(obj->getAABB().max_)) + if(isEqual(dummy_vector, obj->getAABB().max_)) dummy_vector = dummy_vector + delta; else dummy_vector = dummy_vector * 2 - obj->getAABB().max_; diff --git a/src/broadphase/broadphase_SaP.cpp b/src/broadphase/broadphase_SaP.cpp index e60c0a352fdfb1f084e68da9535c6ad3c71605b9..0fb886ea5fd8b7a1ad38770d6cab134a9e892f9b 100644 --- a/src/broadphase/broadphase_SaP.cpp +++ b/src/broadphase/broadphase_SaP.cpp @@ -120,9 +120,9 @@ void SaPCollisionManager::registerObjects(const std::vector<CollisionObject*>& o for(size_t coord = 0; coord < 3; ++coord) { std::sort(endpoints.begin(), endpoints.end(), - boost::bind(std::less<Vec3f::U>(), - boost::bind(static_cast<Vec3f::U (EndPoint::*)(size_t) const >(&EndPoint::getVal), _1, coord), - boost::bind(static_cast<Vec3f::U (EndPoint::*)(size_t) const >(&EndPoint::getVal), _2, coord))); + boost::bind(std::less<Vec3f::Scalar>(), + boost::bind(static_cast<Vec3f::Scalar (EndPoint::*)(size_t) const >(&EndPoint::getVal), _1, coord), + boost::bind(static_cast<Vec3f::Scalar (EndPoint::*)(size_t) const >(&EndPoint::getVal), _2, coord))); endpoints[0]->prev[coord] = NULL; endpoints[0]->next[coord] = endpoints[1]; @@ -513,9 +513,9 @@ bool SaPCollisionManager::collide_(CollisionObject* obj, void* cdata, CollisionC // compute stop_pos by binary search, this is cheaper than check it in while iteration linearly std::vector<EndPoint*>::const_iterator res_it = std::upper_bound(velist[axis].begin(), velist[axis].end(), &dummy, - boost::bind(std::less<Vec3f::U>(), - boost::bind(static_cast<Vec3f::U (EndPoint::*)(size_t) const>(&EndPoint::getVal), _1, axis), - boost::bind(static_cast<Vec3f::U (EndPoint::*)(size_t) const>(&EndPoint::getVal), _2, axis))); + boost::bind(std::less<Vec3f::Scalar>(), + boost::bind(static_cast<Vec3f::Scalar (EndPoint::*)(size_t) const>(&EndPoint::getVal), _1, axis), + boost::bind(static_cast<Vec3f::Scalar (EndPoint::*)(size_t) const>(&EndPoint::getVal), _2, axis))); EndPoint* end_pos = NULL; if(res_it != velist[axis].end()) @@ -579,9 +579,9 @@ bool SaPCollisionManager::distance_(CollisionObject* obj, void* cdata, DistanceC std::vector<EndPoint*>::const_iterator res_it = std::upper_bound(velist[axis].begin(), velist[axis].end(), &dummy, - boost::bind(std::less<Vec3f::U>(), - boost::bind(static_cast<Vec3f::U (EndPoint::*)(size_t) const>(&EndPoint::getVal), _1, axis), - boost::bind(static_cast<Vec3f::U (EndPoint::*)(size_t) const>(&EndPoint::getVal), _2, axis))); + boost::bind(std::less<Vec3f::Scalar>(), + boost::bind(static_cast<Vec3f::Scalar (EndPoint::*)(size_t) const>(&EndPoint::getVal), _1, axis), + boost::bind(static_cast<Vec3f::Scalar (EndPoint::*)(size_t) const>(&EndPoint::getVal), _2, axis))); EndPoint* end_pos = NULL; if(res_it != velist[axis].end()) diff --git a/src/broadphase/hierarchy_tree.cpp b/src/broadphase/hierarchy_tree.cpp index e1d48d989a0eb0789b7d525a9bb1002a700bf27a..c738be797e74cf21d12020dd068262cbee5c1d8a 100644 --- a/src/broadphase/hierarchy_tree.cpp +++ b/src/broadphase/hierarchy_tree.cpp @@ -73,7 +73,7 @@ bool HierarchyTree<AABB>::update(NodeBase<AABB>* leaf, const AABB& bv_, const Ve { AABB bv(bv_); if(leaf->bv.contain(bv)) return false; - Vec3f marginv(margin); + Vec3f marginv(Vec3f::Constant(margin)); bv.min_ -= marginv; bv.max_ += marginv; if(vel[0] > 0) bv.max_[0] += vel[0]; diff --git a/src/ccd/interval_matrix.cpp b/src/ccd/interval_matrix.cpp index 927afec22b158e84f9eed56dc2b1a57fcbd4b6cd..a2ad5811d7e94654b41c25ffdf7b82f929803987 100644 --- a/src/ccd/interval_matrix.cpp +++ b/src/ccd/interval_matrix.cpp @@ -124,16 +124,16 @@ Vec3f IMatrix3::getRowHigh(size_t i) const Matrix3f IMatrix3::getLow() const { - return Matrix3f(v_[0][0][0], v_[0][1][0], v_[0][2][0], - v_[1][0][0], v_[1][1][0], v_[1][2][0], - v_[2][0][0], v_[2][1][0], v_[2][2][0]); + return (Matrix3f() << v_[0][0][0], v_[0][1][0], v_[0][2][0], + v_[1][0][0], v_[1][1][0], v_[1][2][0], + v_[2][0][0], v_[2][1][0], v_[2][2][0]).finished(); } Matrix3f IMatrix3::getHigh() const { - return Matrix3f(v_[0][0][1], v_[0][1][1], v_[0][2][1], - v_[1][0][1], v_[1][1][1], v_[1][2][1], - v_[2][0][1], v_[2][1][1], v_[2][2][1]); + return (Matrix3f() << v_[0][0][1], v_[0][1][1], v_[0][2][1], + v_[1][0][1], v_[1][1][1], v_[1][2][1], + v_[2][0][1], v_[2][1][1], v_[2][2][1]).finished(); } IMatrix3 IMatrix3::operator * (const Matrix3f& m) const diff --git a/src/collision_node.cpp b/src/collision_node.cpp index 02be38e001ea37e404b6d848e16b7f490b76c107..5d1bddb032985ecc9c93df0450ace5a8064ea5f3 100644 --- a/src/collision_node.cpp +++ b/src/collision_node.cpp @@ -67,10 +67,10 @@ void collide2(MeshCollisionTraversalNodeOBB* node, BVHFrontList* front_list) Matrix3f Rtemp, R; Vec3f Ttemp, T; Rtemp = node->R * node->model2->getBV(0).getOrientation(); - R = node->model1->getBV(0).getOrientation().transposeTimes(Rtemp); + R = node->model1->getBV(0).getOrientation().transpose() * Rtemp; Ttemp = node->R * node->model2->getBV(0).getCenter() + node->T; Ttemp -= node->model1->getBV(0).getCenter(); - T = node->model1->getBV(0).getOrientation().transposeTimes(Ttemp); + T = node->model1->getBV(0).getOrientation().transpose() * Ttemp; collisionRecurse(node, 0, 0, R, T, front_list); } diff --git a/src/intersect.cpp b/src/intersect.cpp index bcf3d10db0a818a854ccf4123a34d45914827447..7991ae4a0d29d60cb55c2200a3b1a88f7cd187dd 100644 --- a/src/intersect.cpp +++ b/src/intersect.cpp @@ -1112,7 +1112,7 @@ bool Intersect::buildTrianglePlane(const Vec3f& v1, const Vec3f& v2, const Vec3f { Vec3f n_ = (v2 - v1).cross(v3 - v1); FCL_REAL norm2 = n_.squaredNorm(); - if (n > 0) + if (norm2 > 0) { *n = n_ / sqrt(norm2); *t = n_.dot(v1); @@ -1126,7 +1126,7 @@ bool Intersect::buildEdgePlane(const Vec3f& v1, const Vec3f& v2, const Vec3f& tn { Vec3f n_ = (v2 - v1).cross(tn); FCL_REAL norm2 = n_.squaredNorm(); - if (n > 0) + if (norm2 > 0) { *n = n_ / sqrt(norm2); *t = n_.dot(v1); diff --git a/src/math/transform.cpp b/src/math/transform.cpp index 5fb31596972f611ade8f790c2f9f1f058f1b65ab..c5f406d2a8a01229e62e5750814bdb02ffd2c205 100644 --- a/src/math/transform.cpp +++ b/src/math/transform.cpp @@ -338,7 +338,8 @@ Quaternion3f inverse(const Quaternion3f& q) void Quaternion3f::fromEuler(FCL_REAL a, FCL_REAL b, FCL_REAL c) { Matrix3f R; - R.setEulerYPR(a, b, c); + // R.setEulerYPR(a, b, c); + setEulerYPR(R, a, b, c); fromRotation(R); } diff --git a/src/narrowphase/gjk.cpp b/src/narrowphase/gjk.cpp index 6e01e2021ba875e23a48d857ef1ae7c6add6d57e..fc65947a3c8dc9bd961a584dc9d613d80c63b602 100644 --- a/src/narrowphase/gjk.cpp +++ b/src/narrowphase/gjk.cpp @@ -142,7 +142,7 @@ Vec3f getSupport(const ShapeBase* shape, const Vec3f& dir) const Convex* convex = static_cast<const Convex*>(shape); FCL_REAL maxdot = - std::numeric_limits<FCL_REAL>::max(); Vec3f* curp = convex->points; - Vec3f bestv; + Vec3f bestv(0,0,0); for(int i = 0; i < convex->num_points; ++i, curp+=1) { FCL_REAL dot = dir.dot(*curp); @@ -166,7 +166,7 @@ Vec3f getSupport(const ShapeBase* shape, const Vec3f& dir) void GJK::initialize() { - ray = Vec3f(); + ray = Vec3f::Zero(); nfree = 0; status = Failed; current = 0; @@ -201,7 +201,9 @@ GJK::Status GJK::evaluate(const MinkowskiDiff& shape_, const Vec3f& guess) simplices[0].rank = 0; ray = guess; - appendVertex(simplices[0], (ray.squaredNorm() > 0) ? -ray : Vec3f(1, 0, 0)); + // appendVertex(simplices[0], (ray.squaredNorm() > 0) ? -ray : Vec3f(1, 0, 0)); + if (ray.squaredNorm() > 0) appendVertex(simplices[0], -ray); + else appendVertex(simplices[0], Vec3f(1, 0, 0)); simplices[0].p[0] = 1; ray = simplices[0].c[0]->w; lastw[0] = lastw[1] = lastw[2] = lastw[3] = ray; // cache previous support points, the new support point will compare with it to avoid too close support points @@ -266,7 +268,7 @@ GJK::Status GJK::evaluate(const MinkowskiDiff& shape_, const Vec3f& guess) if(project_res.sqr_distance >= 0) { next_simplex.rank = 0; - ray = Vec3f(); + ray = Vec3f(0,0,0); current = next; for(size_t i = 0; i < curr_simplex.rank; ++i) { @@ -303,13 +305,13 @@ GJK::Status GJK::evaluate(const MinkowskiDiff& shape_, const Vec3f& guess) void GJK::getSupport(const Vec3f& d, SimplexV& sv) const { - sv.d = normalize(d); + sv.d = d.normalized(); sv.w = shape.support(sv.d); } void GJK::getSupport(const Vec3f& d, const Vec3f& v, SimplexV& sv) const { - sv.d = normalize(d); + sv.d = d.normalized(); sv.w = shape.support(sv.d, v); } @@ -333,7 +335,7 @@ bool GJK::encloseOrigin() { for(size_t i = 0; i < 3; ++i) { - Vec3f axis; + Vec3f axis(Vec3f::Zero()); axis[i] = 1; appendVertex(*simplex, axis); if(encloseOrigin()) return true; @@ -349,7 +351,7 @@ bool GJK::encloseOrigin() Vec3f d = simplex->c[1]->w - simplex->c[0]->w; for(size_t i = 0; i < 3; ++i) { - Vec3f axis; + Vec3f axis(0,0,0); axis[i] = 1; Vec3f p = d.cross(axis); if(p.squaredNorm() > 0) diff --git a/src/narrowphase/narrowphase.cpp b/src/narrowphase/narrowphase.cpp index 3e1601ec3a88a08039c4e2d50f1cda8afb68d852..befa60902cfa77b950d5662dd8b0488a20b8f9a6 100644 --- a/src/narrowphase/narrowphase.cpp +++ b/src/narrowphase/narrowphase.cpp @@ -89,7 +89,7 @@ bool sphereCapsuleIntersect(const Sphere& s1, const Transform3f& tf1, if (distance > 0) return false; - Vec3f normal = diff.normalize() * - FCL_REAL(1); + Vec3f normal = diff.normalized() * - FCL_REAL(1); if (distance < 0 && penetration_depth) *penetration_depth = -distance; @@ -304,7 +304,7 @@ bool sphereTriangleIntersect(const Sphere& s, const Transform3f& tf, if(distance_sqr > 0) { FCL_REAL distance = std::sqrt(distance_sqr); - if(normal_) *normal_ = normalize(contact_to_center); + if(normal_) *normal_ = contact_to_center.normalized(); if(contact_points) *contact_points = contact_point; if(penetration_depth) *penetration_depth = -(radius - distance); } @@ -795,16 +795,16 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, FCL_REAL s, s2, l; int invert_normal, code; - Vec3f p = T2 - T1; // get vector from centers of box 1 to box 2, relative to box 1 - Vec3f pp = R1.transposeTimes(p); // get pp = p relative to body 1 + Vec3f p (T2 - T1); // get vector from centers of box 1 to box 2, relative to box 1 + Vec3f pp (R1.transpose() * p); // get pp = p relative to body 1 // get side lengths / 2 Vec3f A = side1 * 0.5; Vec3f B = side2 * 0.5; // Rij is R1'*R2, i.e. the relative rotation between R1 and R2 - Matrix3f R = R1.transposeTimes(R2); - Matrix3f Q = abs(R); + Matrix3f R (R1.transpose() * R2); + Matrix3f Q (R.cwiseAbs()); // for all 15 possible separating axes: @@ -827,7 +827,7 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, // separating axis = u1, u2, u3 tmp = pp[0]; - s2 = std::abs(tmp) - (Q.dotX(B) + A[0]); + s2 = std::abs(tmp) - (Q.row(0).dot(B) + A[0]); if(s2 > 0) { *return_code = 0; return 0; } if(s2 > s) { @@ -839,7 +839,7 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, } tmp = pp[1]; - s2 = std::abs(tmp) - (Q.dotY(B) + A[1]); + s2 = std::abs(tmp) - (Q.row(1).dot(B) + A[1]); if(s2 > 0) { *return_code = 0; return 0; } if(s2 > s) { @@ -851,7 +851,7 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, } tmp = pp[2]; - s2 = std::abs(tmp) - (Q.dotZ(B) + A[2]); + s2 = std::abs(tmp) - (Q.row(2).dot(B) + A[2]); if(s2 > 0) { *return_code = 0; return 0; } if(s2 > s) { @@ -863,8 +863,8 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, } // separating axis = v1, v2, v3 - tmp = R2.transposeDotX(p); - s2 = std::abs(tmp) - (Q.transposeDotX(A) + B[0]); + tmp = R2.col(0).dot(p); + s2 = std::abs(tmp) - (Q.col(0).dot(A) + B[0]); if(s2 > 0) { *return_code = 0; return 0; } if(s2 > s) { @@ -875,8 +875,8 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, code = 4; } - tmp = R2.transposeDotY(p); - s2 = std::abs(tmp) - (Q.transposeDotY(A) + B[1]); + tmp = R2.col(1).dot(p); + s2 = std::abs(tmp) - (Q.col(1).dot(A) + B[1]); if(s2 > 0) { *return_code = 0; return 0; } if(s2 > s) { @@ -887,8 +887,8 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, code = 5; } - tmp = R2.transposeDotZ(p); - s2 = std::abs(tmp) - (Q.transposeDotZ(A) + B[2]); + tmp = R2.col(2).dot(p); + s2 = std::abs(tmp) - (Q.col(2).dot(A) + B[2]); if(s2 > 0) { *return_code = 0; return 0; } if(s2 > s) { @@ -901,7 +901,7 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, FCL_REAL fudge2(1.0e-6); - Q += fudge2; + Q.array() += fudge2; Vec3f n; FCL_REAL eps = std::numeric_limits<FCL_REAL>::epsilon(); @@ -1098,7 +1098,7 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, for(int j = 0; j < 3; ++j) { - sign = (R1.transposeDot(j, normal) > 0) ? 1 : -1; + sign = (R1.col(j).dot(normal) > 0) ? 1 : -1; pa += R1.col(j) * (A[j] * sign); } @@ -1107,7 +1107,7 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, for(int j = 0; j < 3; ++j) { - sign = (R2.transposeDot(j, normal) > 0) ? -1 : 1; + sign = (R2.col(j).dot(normal) > 0) ? -1 : 1; pb += R2.col(j) * (B[j] * sign); } @@ -1163,8 +1163,8 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, else normal2 = -normal; - nr = Rb->transposeTimes(normal2); - anr = abs(nr); + nr = Rb->transpose() * normal2; + anr = nr.cwiseAbs(); // find the largest compontent of anr: this corresponds to the normal // for the indident face. the other axis numbers of the indicent face @@ -1233,17 +1233,17 @@ int boxBox2(const Vec3f& side1, const Matrix3f& R1, const Vec3f& T1, // find the four corners of the incident face, in reference-face coordinates FCL_REAL quad[8]; // 2D coordinate of incident face (x,y pairs) FCL_REAL c1, c2, m11, m12, m21, m22; - c1 = Ra->transposeDot(code1, center); - c2 = Ra->transposeDot(code2, center); + c1 = Ra->col(code1).dot(center); + c2 = Ra->col(code2).dot(center); // optimize this? - we have already computed this data above, but it is not // stored in an easy-to-index format. for now it's quicker just to recompute // the four dot products. Vec3f tempRac = Ra->col(code1); - m11 = Rb->transposeDot(a1, tempRac); - m12 = Rb->transposeDot(a2, tempRac); + m11 = Rb->col(a1).dot(tempRac); + m12 = Rb->col(a2).dot(tempRac); tempRac = Ra->col(code2); - m21 = Rb->transposeDot(a1, tempRac); - m22 = Rb->transposeDot(a2, tempRac); + m21 = Rb->col(a1).dot(tempRac); + m22 = Rb->col(a2).dot(tempRac); FCL_REAL k1 = m11 * (*Sb)[a1]; FCL_REAL k2 = m21 * (*Sb)[a1]; @@ -1433,9 +1433,9 @@ bool boxHalfspaceIntersect(const Box& s1, const Transform3f& tf1, const Matrix3f& R = tf1.getRotation(); const Vec3f& T = tf1.getTranslation(); - Vec3f Q = R.transposeTimes(new_s2.n); + Vec3f Q = R.transpose() * new_s2.n; Vec3f A(Q[0] * s1.side[0], Q[1] * s1.side[1], Q[2] * s1.side[2]); - Vec3f B = abs(A); + Vec3f B = A.cwiseAbs(); FCL_REAL depth = 0.5 * (B[0] + B[1] + B[2]) - new_s2.signedDistance(T); return (depth >= 0); @@ -1455,9 +1455,9 @@ bool boxHalfspaceIntersect(const Box& s1, const Transform3f& tf1, const Matrix3f& R = tf1.getRotation(); const Vec3f& T = tf1.getTranslation(); - Vec3f Q = R.transposeTimes(new_s2.n); + Vec3f Q = R.transpose() * new_s2.n; Vec3f A(Q[0] * s1.side[0], Q[1] * s1.side[1], Q[2] * s1.side[2]); - Vec3f B = abs(A); + Vec3f B = A.cwiseAbs(); FCL_REAL depth = 0.5 * (B[0] + B[1] + B[2]) - new_s2.signedDistance(T); if(depth < 0) return false; @@ -1908,9 +1908,9 @@ bool boxPlaneIntersect(const Box& s1, const Transform3f& tf1, const Matrix3f& R = tf1.getRotation(); const Vec3f& T = tf1.getTranslation(); - Vec3f Q = R.transposeTimes(new_s2.n); + Vec3f Q = R.transpose() * new_s2.n; Vec3f A(Q[0] * s1.side[0], Q[1] * s1.side[1], Q[2] * s1.side[2]); - Vec3f B = abs(A); + Vec3f B = A.cwiseAbs(); FCL_REAL signed_dist = new_s2.signedDistance(T); FCL_REAL depth = 0.5 * (B[0] + B[1] + B[2]) - std::abs(signed_dist); @@ -2081,7 +2081,7 @@ bool cylinderPlaneIntersect(const Cylinder& s1, const Transform3f& tf1, const Matrix3f& R = tf1.getRotation(); const Vec3f& T = tf1.getTranslation(); - Vec3f Q = R.transposeTimes(new_s2.n); + Vec3f Q = R.transpose() * new_s2.n; FCL_REAL term = std::abs(Q[2]) * s1.lz + s1.radius * std::sqrt(Q[0] * Q[0] + Q[1] * Q[1]); FCL_REAL dist = new_s2.distance(T); diff --git a/src/shape/geometric_shapes_utility.cpp b/src/shape/geometric_shapes_utility.cpp index 3cc09c6e3586f686f16cfdbf08ab7b86ba976d79..058a59088c3f1834e9ffa4acb4c444f6b407e642 100644 --- a/src/shape/geometric_shapes_utility.cpp +++ b/src/shape/geometric_shapes_utility.cpp @@ -945,33 +945,33 @@ void constructBox(const AABB& bv, Box& box, Transform3f& tf) void constructBox(const OBB& bv, Box& box, Transform3f& tf) { box = Box(bv.extent * 2); - tf = Transform3f(Matrix3f(bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], - bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], - bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]), bv.To); + tf = Transform3f((Matrix3f() << bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], + bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], + bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]).finished(), bv.To); } void constructBox(const OBBRSS& bv, Box& box, Transform3f& tf) { box = Box(bv.obb.extent * 2); - tf = Transform3f(Matrix3f(bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], + tf = Transform3f((Matrix3f() << bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], bv.obb.axis[0][1], bv.obb.axis[1][1], bv.obb.axis[2][1], - bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]), bv.obb.To); + bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]).finished(), bv.obb.To); } void constructBox(const kIOS& bv, Box& box, Transform3f& tf) { box = Box(bv.obb.extent * 2); - tf = Transform3f(Matrix3f(bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], + tf = Transform3f((Matrix3f() << bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], bv.obb.axis[0][1], bv.obb.axis[1][1], bv.obb.axis[2][1], - bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]), bv.obb.To); + bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]).finished(), bv.obb.To); } void constructBox(const RSS& bv, Box& box, Transform3f& tf) { box = Box(bv.width(), bv.height(), bv.depth()); - tf = Transform3f(Matrix3f(bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], + tf = Transform3f((Matrix3f() << bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], - bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]), bv.Tr); + bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]).finished(), bv.Tr); } void constructBox(const KDOP<16>& bv, Box& box, Transform3f& tf) @@ -1003,33 +1003,33 @@ void constructBox(const AABB& bv, const Transform3f& tf_bv, Box& box, Transform3 void constructBox(const OBB& bv, const Transform3f& tf_bv, Box& box, Transform3f& tf) { box = Box(bv.extent * 2); - tf = tf_bv *Transform3f(Matrix3f(bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], + tf = tf_bv *Transform3f((Matrix3f() << bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], - bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]), bv.To); + bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]).finished(), bv.To); } void constructBox(const OBBRSS& bv, const Transform3f& tf_bv, Box& box, Transform3f& tf) { box = Box(bv.obb.extent * 2); - tf = tf_bv * Transform3f(Matrix3f(bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], + tf = tf_bv * Transform3f((Matrix3f() << bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], bv.obb.axis[0][1], bv.obb.axis[1][1], bv.obb.axis[2][1], - bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]), bv.obb.To); + bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]).finished(), bv.obb.To); } void constructBox(const kIOS& bv, const Transform3f& tf_bv, Box& box, Transform3f& tf) { box = Box(bv.obb.extent * 2); - tf = tf_bv * Transform3f(Matrix3f(bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], + tf = tf_bv * Transform3f((Matrix3f() << bv.obb.axis[0][0], bv.obb.axis[1][0], bv.obb.axis[2][0], bv.obb.axis[0][1], bv.obb.axis[1][1], bv.obb.axis[2][1], - bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]), bv.obb.To); + bv.obb.axis[0][2], bv.obb.axis[1][2], bv.obb.axis[2][2]).finished(), bv.obb.To); } void constructBox(const RSS& bv, const Transform3f& tf_bv, Box& box, Transform3f& tf) { box = Box(bv.width(), bv.height(), bv.depth()); - tf = tf_bv * Transform3f(Matrix3f(bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], + tf = tf_bv * Transform3f((Matrix3f() << bv.axis[0][0], bv.axis[1][0], bv.axis[2][0], bv.axis[0][1], bv.axis[1][1], bv.axis[2][1], - bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]), bv.Tr); + bv.axis[0][2], bv.axis[1][2], bv.axis[2][2]).finished(), bv.Tr); } void constructBox(const KDOP<16>& bv, const Transform3f& tf_bv, Box& box, Transform3f& tf) diff --git a/test/test_fcl_geometric_shapes.cpp b/test/test_fcl_geometric_shapes.cpp index 7199b307c626a24e7273fd035ac030553d5194af..fa8bae9e89088c7b303fbbc06fdfacc32c37098f 100644 --- a/test/test_fcl_geometric_shapes.cpp +++ b/test/test_fcl_geometric_shapes.cpp @@ -179,7 +179,7 @@ void compareContact(const S1& s1, const Transform3f& tf1, { if (expected_point) { - bool contact_equal = contact.equal(*expected_point, tol); + bool contact_equal = isEqual(contact, *expected_point, tol); BOOST_CHECK(contact_equal); if (!contact_equal) printComparisonError("contact", s1, tf1, s2, tf2, solver_type, contact, *expected_point, false, tol); @@ -195,10 +195,10 @@ void compareContact(const S1& s1, const Transform3f& tf1, if (expected_normal) { - bool normal_equal = normal.equal(*expected_normal, tol); + bool normal_equal = isEqual(normal, *expected_normal, tol); if (!normal_equal && check_opposite_normal) - normal_equal = normal.equal(-(*expected_normal), tol); + normal_equal = isEqual(normal, -(*expected_normal), tol); BOOST_CHECK(normal_equal); if (!normal_equal) @@ -407,7 +407,7 @@ void testBoxBoxContactPoints(const Matrix3f& R) std::sort(vertices.begin(), vertices.end(), compareContactPoints); // The lowest vertex along z-axis should be the contact point - BOOST_CHECK(vertices[0].equal(point)); + BOOST_CHECK(isEqual(vertices[0], point)); } BOOST_AUTO_TEST_CASE(shapeIntersection_boxbox) @@ -761,11 +761,11 @@ BOOST_AUTO_TEST_CASE(shapeIntersection_spheretriangle) res = solver1.shapeTriangleIntersect(s, Transform3f(), t[0], t[1], t[2], NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, Vec3f(1, 0, 0), 1e-9)); res = solver1.shapeTriangleIntersect(s, transform, t[0], t[1], t[2], transform, NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); } BOOST_AUTO_TEST_CASE(shapeIntersection_halfspacetriangle) @@ -802,11 +802,11 @@ BOOST_AUTO_TEST_CASE(shapeIntersection_halfspacetriangle) res = solver1.shapeTriangleIntersect(hs, Transform3f(), t[0], t[1], t[2], Transform3f(), NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, Vec3f(1, 0, 0), 1e-9)); res = solver1.shapeTriangleIntersect(hs, transform, t[0], t[1], t[2], transform, NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); } BOOST_AUTO_TEST_CASE(shapeIntersection_planetriangle) @@ -843,11 +843,11 @@ BOOST_AUTO_TEST_CASE(shapeIntersection_planetriangle) res = solver1.shapeTriangleIntersect(hs, Transform3f(), t[0], t[1], t[2], Transform3f(), NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, Vec3f(1, 0, 0), 1e-9)); res = solver1.shapeTriangleIntersect(hs, transform, t[0], t[1], t[2], transform, NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); } BOOST_AUTO_TEST_CASE(shapeIntersection_halfspacesphere) @@ -3126,11 +3126,11 @@ BOOST_AUTO_TEST_CASE(shapeIntersectionGJK_spheretriangle) res = solver2.shapeTriangleIntersect(s, Transform3f(), t[0], t[1], t[2], NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, Vec3f(1, 0, 0), 1e-9)); res = solver2.shapeTriangleIntersect(s, transform, t[0], t[1], t[2], transform, NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); } BOOST_AUTO_TEST_CASE(shapeIntersectionGJK_halfspacetriangle) @@ -3167,11 +3167,11 @@ BOOST_AUTO_TEST_CASE(shapeIntersectionGJK_halfspacetriangle) res = solver2.shapeTriangleIntersect(hs, Transform3f(), t[0], t[1], t[2], Transform3f(), NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, Vec3f(1, 0, 0), 1e-9)); res = solver2.shapeTriangleIntersect(hs, transform, t[0], t[1], t[2], transform, NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); } BOOST_AUTO_TEST_CASE(shapeIntersectionGJK_planetriangle) @@ -3208,11 +3208,11 @@ BOOST_AUTO_TEST_CASE(shapeIntersectionGJK_planetriangle) res = solver2.shapeTriangleIntersect(hs, Transform3f(), t[0], t[1], t[2], Transform3f(), NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, Vec3f(1, 0, 0), 1e-9)); res = solver2.shapeTriangleIntersect(hs, transform, t[0], t[1], t[2], transform, NULL, NULL, &normal); BOOST_CHECK(res); - BOOST_CHECK(normal.equal(transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); + BOOST_CHECK(isEqual(normal, transform.getRotation() * Vec3f(1, 0, 0), 1e-9)); } @@ -3451,8 +3451,8 @@ void testReversibleShapeIntersection(const S1& s1, const S2& s2, FCL_REAL distan BOOST_CHECK(resA); BOOST_CHECK(resB); - BOOST_CHECK(contactA.equal(contactB, tol)); // contact point should be same - BOOST_CHECK(normalA.equal(-normalB, tol)); // normal should be opposite + BOOST_CHECK(isEqual(contactA, contactB, tol)); // contact point should be same + BOOST_CHECK(isEqual(normalA, -normalB, tol)); // normal should be opposite BOOST_CHECK_CLOSE(depthA, depthB, tol); // penetration depth should be same resA = solver2.shapeIntersect(s1, tf1, s2, tf2, &contactA, &depthA, &normalA); @@ -3460,8 +3460,8 @@ void testReversibleShapeIntersection(const S1& s1, const S2& s2, FCL_REAL distan BOOST_CHECK(resA); BOOST_CHECK(resB); - BOOST_CHECK(contactA.equal(contactB, tol)); - BOOST_CHECK(normalA.equal(-normalB, tol)); + BOOST_CHECK(isEqual(contactA, contactB, tol)); + BOOST_CHECK(isEqual(normalA, -normalB, tol)); BOOST_CHECK_CLOSE(depthA, depthB, tol); } @@ -3539,8 +3539,8 @@ void testReversibleShapeDistance(const S1& s1, const S2& s2, FCL_REAL distance) BOOST_CHECK(resA); BOOST_CHECK(resB); BOOST_CHECK_CLOSE(distA, distB, tol); // distances should be same - BOOST_CHECK(p1A.equal(p2B, tol)); // closest points should in reverse order - BOOST_CHECK(p2A.equal(p1B, tol)); + BOOST_CHECK(isEqual(p1A, p2B, tol)); // closest points should in reverse order + BOOST_CHECK(isEqual(p2A, p1B, tol)); resA = solver2.shapeDistance(s1, tf1, s2, tf2, &distA, &p1A, &p2A); resB = solver2.shapeDistance(s2, tf2, s1, tf1, &distB, &p1B, &p2B); @@ -3548,8 +3548,8 @@ void testReversibleShapeDistance(const S1& s1, const S2& s2, FCL_REAL distance) BOOST_CHECK(resA); BOOST_CHECK(resB); BOOST_CHECK_CLOSE(distA, distB, tol); - BOOST_CHECK(p1A.equal(p2B, tol)); - BOOST_CHECK(p2A.equal(p1B, tol)); + BOOST_CHECK(isEqual(p1A, p2B, tol)); + BOOST_CHECK(isEqual(p2A, p1B, tol)); } BOOST_AUTO_TEST_CASE(reversibleShapeDistance_allshapes) diff --git a/test/test_fcl_sphere_capsule.cpp b/test/test_fcl_sphere_capsule.cpp index 4e0763aafff2b247e3635b5263c617b16f995b51..c74e41a91aab3bf930fe143a0197afe1c8c3ddad 100644 --- a/test/test_fcl_sphere_capsule.cpp +++ b/test/test_fcl_sphere_capsule.cpp @@ -100,7 +100,8 @@ BOOST_AUTO_TEST_CASE(Sphere_Capsule_Intersect_test_separated_capsule_rotated) Capsule capsule (50, 200.); Matrix3f rotation; - rotation.setEulerZYX (M_PI * 0.5, 0., 0.); + // rotation.setEulerZYX (M_PI * 0.5, 0., 0.); + setEulerZYX(rotation, M_PI * 0.5, 0., 0.); Transform3f capsule_transform (rotation, Vec3f (150., 0., 0.)); BOOST_CHECK (!solver.shapeIntersect(sphere1, sphere1_transform, capsule, capsule_transform, NULL, NULL, NULL)); @@ -125,8 +126,8 @@ BOOST_AUTO_TEST_CASE(Sphere_Capsule_Intersect_test_penetration_z) BOOST_CHECK (is_intersecting); BOOST_CHECK (penetration == 25.); - BOOST_CHECK (Vec3f (0., 0., 1.).equal(normal)); - BOOST_CHECK (Vec3f (0., 0., 0.).equal(contact_point)); + BOOST_CHECK (isEqual(Vec3f (0., 0., 1.), normal)); + BOOST_CHECK (isEqual(Vec3f (0., 0., 0.), contact_point)); } BOOST_AUTO_TEST_CASE(Sphere_Capsule_Intersect_test_penetration_z_rotated) @@ -139,7 +140,8 @@ BOOST_AUTO_TEST_CASE(Sphere_Capsule_Intersect_test_penetration_z_rotated) Capsule capsule (50, 200.); Matrix3f rotation; - rotation.setEulerZYX (M_PI * 0.5, 0., 0.); + // rotation.setEulerZYX (M_PI * 0.5, 0., 0.); + setEulerZYX(rotation, M_PI * 0.5, 0., 0.); Transform3f capsule_transform (rotation, Vec3f (0., 50., 75)); FCL_REAL penetration = 0.; @@ -150,8 +152,8 @@ BOOST_AUTO_TEST_CASE(Sphere_Capsule_Intersect_test_penetration_z_rotated) BOOST_CHECK (is_intersecting); BOOST_CHECK_CLOSE (25, penetration, solver.collision_tolerance); - BOOST_CHECK (Vec3f (0., 0., 1.).equal(normal)); - BOOST_CHECK (Vec3f (0., 0., 50.).equal(contact_point, solver.collision_tolerance)); + BOOST_CHECK (isEqual(Vec3f (0., 0., 1.),normal)); + BOOST_CHECK (isEqual(Vec3f (0., 0., 50.), contact_point, solver.collision_tolerance)); } BOOST_AUTO_TEST_CASE(Sphere_Capsule_Distance_test_collision)