AlignedBox.h
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr>
00005 //
00006 // This Source Code Form is subject to the terms of the Mozilla
00007 // Public License v. 2.0. If a copy of the MPL was not distributed
00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
00009 
00010 #ifndef EIGEN_ALIGNEDBOX_H
00011 #define EIGEN_ALIGNEDBOX_H
00012 
00013 namespace Eigen { 
00014 
00027 template <typename _Scalar, int _AmbientDim>
00028 class AlignedBox
00029 {
00030 public:
00031 EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar,_AmbientDim)
00032   enum { AmbientDimAtCompileTime = _AmbientDim };
00033   typedef _Scalar                                   Scalar;
00034   typedef NumTraits<Scalar>                         ScalarTraits;
00035   typedef DenseIndex                                Index;
00036   typedef typename ScalarTraits::Real               RealScalar;
00037   typedef typename ScalarTraits::NonInteger      NonInteger;
00038   typedef Matrix<Scalar,AmbientDimAtCompileTime,1>  VectorType;
00039 
00041   enum CornerType
00042   {
00044     Min=0, Max=1,
00045 
00047     BottomLeft=0, BottomRight=1,
00048     TopLeft=2, TopRight=3,
00049 
00051     BottomLeftFloor=0, BottomRightFloor=1,
00052     TopLeftFloor=2, TopRightFloor=3,
00053     BottomLeftCeil=4, BottomRightCeil=5,
00054     TopLeftCeil=6, TopRightCeil=7
00055   };
00056 
00057 
00059   inline AlignedBox()
00060   { if (AmbientDimAtCompileTime!=Dynamic) setEmpty(); }
00061 
00063   inline explicit AlignedBox(Index _dim) : m_min(_dim), m_max(_dim)
00064   { setEmpty(); }
00065 
00067   template<typename OtherVectorType1, typename OtherVectorType2>
00068   inline AlignedBox(const OtherVectorType1& _min, const OtherVectorType2& _max) : m_min(_min), m_max(_max) {}
00069 
00071   template<typename Derived>
00072   inline explicit AlignedBox(const MatrixBase<Derived>& a_p)
00073   {
00074     typename internal::nested<Derived,2>::type p(a_p.derived());
00075     m_min = p;
00076     m_max = p;
00077   }
00078 
00079   ~AlignedBox() {}
00080 
00082   inline Index dim() const { return AmbientDimAtCompileTime==Dynamic ? m_min.size() : Index(AmbientDimAtCompileTime); }
00083 
00085   inline bool isNull() const { return isEmpty(); }
00086 
00088   inline void setNull() { setEmpty(); }
00089 
00091   inline bool isEmpty() const { return (m_min.array() > m_max.array()).any(); }
00092 
00094   inline void setEmpty()
00095   {
00096     m_min.setConstant( ScalarTraits::highest() );
00097     m_max.setConstant( ScalarTraits::lowest() );
00098   }
00099 
00101   inline const VectorType& (min)() const { return m_min; }
00103   inline VectorType& (min)() { return m_min; }
00105   inline const VectorType& (max)() const { return m_max; }
00107   inline VectorType& (max)() { return m_max; }
00108 
00110   inline const CwiseUnaryOp<internal::scalar_quotient1_op<Scalar>,
00111                             const CwiseBinaryOp<internal::scalar_sum_op<Scalar>, const VectorType, const VectorType> >
00112   center() const
00113   { return (m_min+m_max)/2; }
00114 
00119   inline const CwiseBinaryOp< internal::scalar_difference_op<Scalar>, const VectorType, const VectorType> sizes() const
00120   { return m_max - m_min; }
00121 
00123   inline Scalar volume() const
00124   { return sizes().prod(); }
00125 
00130   inline CwiseBinaryOp< internal::scalar_difference_op<Scalar>, const VectorType, const VectorType> diagonal() const
00131   { return sizes(); }
00132 
00142   inline VectorType corner(CornerType corner) const
00143   {
00144     EIGEN_STATIC_ASSERT(_AmbientDim <= 3, THIS_METHOD_IS_ONLY_FOR_VECTORS_OF_A_SPECIFIC_SIZE);
00145 
00146     VectorType res;
00147 
00148     Index mult = 1;
00149     for(Index d=0; d<dim(); ++d)
00150     {
00151       if( mult & corner ) res[d] = m_max[d];
00152       else                res[d] = m_min[d];
00153       mult *= 2;
00154     }
00155     return res;
00156   }
00157 
00160   inline VectorType sample() const
00161   {
00162     VectorType r;
00163     for(Index d=0; d<dim(); ++d)
00164     {
00165       if(!ScalarTraits::IsInteger)
00166       {
00167         r[d] = m_min[d] + (m_max[d]-m_min[d])
00168              * internal::random<Scalar>(Scalar(0), Scalar(1));
00169       }
00170       else
00171         r[d] = internal::random(m_min[d], m_max[d]);
00172     }
00173     return r;
00174   }
00175 
00177   template<typename Derived>
00178   inline bool contains(const MatrixBase<Derived>& a_p) const
00179   {
00180     typename internal::nested<Derived,2>::type p(a_p.derived());
00181     return (m_min.array()<=p.array()).all() && (p.array()<=m_max.array()).all();
00182   }
00183 
00185   inline bool contains(const AlignedBox& b) const
00186   { return (m_min.array()<=(b.min)().array()).all() && ((b.max)().array()<=m_max.array()).all(); }
00187 
00189   template<typename Derived>
00190   inline AlignedBox& extend(const MatrixBase<Derived>& a_p)
00191   {
00192     typename internal::nested<Derived,2>::type p(a_p.derived());
00193     m_min = m_min.cwiseMin(p);
00194     m_max = m_max.cwiseMax(p);
00195     return *this;
00196   }
00197 
00199   inline AlignedBox& extend(const AlignedBox& b)
00200   {
00201     m_min = m_min.cwiseMin(b.m_min);
00202     m_max = m_max.cwiseMax(b.m_max);
00203     return *this;
00204   }
00205 
00207   inline AlignedBox& clamp(const AlignedBox& b)
00208   {
00209     m_min = m_min.cwiseMax(b.m_min);
00210     m_max = m_max.cwiseMin(b.m_max);
00211     return *this;
00212   }
00213 
00215   inline AlignedBox intersection(const AlignedBox& b) const
00216   {return AlignedBox(m_min.cwiseMax(b.m_min), m_max.cwiseMin(b.m_max)); }
00217 
00219   inline AlignedBox merged(const AlignedBox& b) const
00220   { return AlignedBox(m_min.cwiseMin(b.m_min), m_max.cwiseMax(b.m_max)); }
00221 
00223   template<typename Derived>
00224   inline AlignedBox& translate(const MatrixBase<Derived>& a_t)
00225   {
00226     const typename internal::nested<Derived,2>::type t(a_t.derived());
00227     m_min += t;
00228     m_max += t;
00229     return *this;
00230   }
00231 
00236   template<typename Derived>
00237   inline Scalar squaredExteriorDistance(const MatrixBase<Derived>& a_p) const;
00238 
00243   inline Scalar squaredExteriorDistance(const AlignedBox& b) const;
00244 
00249   template<typename Derived>
00250   inline NonInteger exteriorDistance(const MatrixBase<Derived>& p) const
00251   { using std::sqrt; return sqrt(NonInteger(squaredExteriorDistance(p))); }
00252 
00257   inline NonInteger exteriorDistance(const AlignedBox& b) const
00258   { using std::sqrt; return sqrt(NonInteger(squaredExteriorDistance(b))); }
00259 
00265   template<typename NewScalarType>
00266   inline typename internal::cast_return_type<AlignedBox,
00267            AlignedBox<NewScalarType,AmbientDimAtCompileTime> >::type cast() const
00268   {
00269     return typename internal::cast_return_type<AlignedBox,
00270                     AlignedBox<NewScalarType,AmbientDimAtCompileTime> >::type(*this);
00271   }
00272 
00274   template<typename OtherScalarType>
00275   inline explicit AlignedBox(const AlignedBox<OtherScalarType,AmbientDimAtCompileTime>& other)
00276   {
00277     m_min = (other.min)().template cast<Scalar>();
00278     m_max = (other.max)().template cast<Scalar>();
00279   }
00280 
00285   bool isApprox(const AlignedBox& other, const RealScalar& prec = ScalarTraits::dummy_precision()) const
00286   { return m_min.isApprox(other.m_min, prec) && m_max.isApprox(other.m_max, prec); }
00287 
00288 protected:
00289 
00290   VectorType m_min, m_max;
00291 };
00292 
00293 
00294 
00295 template<typename Scalar,int AmbientDim>
00296 template<typename Derived>
00297 inline Scalar AlignedBox<Scalar,AmbientDim>::squaredExteriorDistance(const MatrixBase<Derived>& a_p) const
00298 {
00299   typename internal::nested<Derived,2*AmbientDim>::type p(a_p.derived());
00300   Scalar dist2(0);
00301   Scalar aux;
00302   for (Index k=0; k<dim(); ++k)
00303   {
00304     if( m_min[k] > p[k] )
00305     {
00306       aux = m_min[k] - p[k];
00307       dist2 += aux*aux;
00308     }
00309     else if( p[k] > m_max[k] )
00310     {
00311       aux = p[k] - m_max[k];
00312       dist2 += aux*aux;
00313     }
00314   }
00315   return dist2;
00316 }
00317 
00318 template<typename Scalar,int AmbientDim>
00319 inline Scalar AlignedBox<Scalar,AmbientDim>::squaredExteriorDistance(const AlignedBox& b) const
00320 {
00321   Scalar dist2(0);
00322   Scalar aux;
00323   for (Index k=0; k<dim(); ++k)
00324   {
00325     if( m_min[k] > b.m_max[k] )
00326     {
00327       aux = m_min[k] - b.m_max[k];
00328       dist2 += aux*aux;
00329     }
00330     else if( b.m_min[k] > m_max[k] )
00331     {
00332       aux = b.m_min[k] - m_max[k];
00333       dist2 += aux*aux;
00334     }
00335   }
00336   return dist2;
00337 }
00338 
00355 #define EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, Size, SizeSuffix)    \
00356                                  \
00357 typedef AlignedBox<Type, Size>   AlignedBox##SizeSuffix##TypeSuffix;
00358 
00359 #define EIGEN_MAKE_TYPEDEFS_ALL_SIZES(Type, TypeSuffix) \
00360 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 1, 1) \
00361 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 2, 2) \
00362 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 3, 3) \
00363 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 4, 4) \
00364 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, Dynamic, X)
00365 
00366 EIGEN_MAKE_TYPEDEFS_ALL_SIZES(int,                  i)
00367 EIGEN_MAKE_TYPEDEFS_ALL_SIZES(float,                f)
00368 EIGEN_MAKE_TYPEDEFS_ALL_SIZES(double,               d)
00369 
00370 #undef EIGEN_MAKE_TYPEDEFS_ALL_SIZES
00371 #undef EIGEN_MAKE_TYPEDEFS
00372 
00373 } // end namespace Eigen
00374 
00375 #endif // EIGEN_ALIGNEDBOX_H


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Thu Aug 27 2015 11:57:48