AlignedBox.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // Eigen is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 3 of the License, or (at your option) any later version.
10 //
11 // Alternatively, you can redistribute it and/or
12 // modify it under the terms of the GNU General Public License as
13 // published by the Free Software Foundation; either version 2 of
14 // the License, or (at your option) any later version.
15 //
16 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
17 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU Lesser General Public
22 // License and a copy of the GNU General Public License along with
23 // Eigen. If not, see <http://www.gnu.org/licenses/>.
24 
25 #ifndef EIGEN_ALIGNEDBOX_H
26 #define EIGEN_ALIGNEDBOX_H
27 
28 namespace Eigen {
29 
42 template <typename _Scalar, int _AmbientDim>
44 {
45 public:
47  enum { AmbientDimAtCompileTime = _AmbientDim };
48  typedef _Scalar Scalar;
50  typedef DenseIndex Index;
51  typedef typename ScalarTraits::Real RealScalar;
52  typedef typename ScalarTraits::NonInteger NonInteger;
54 
57  {
59  Min=0, Max=1,
60 
64 
70  };
71 
72 
74  inline explicit AlignedBox()
76 
78  inline explicit AlignedBox(Index _dim) : m_min(_dim), m_max(_dim)
79  { setEmpty(); }
80 
82  template<typename OtherVectorType1, typename OtherVectorType2>
83  inline AlignedBox(const OtherVectorType1& _min, const OtherVectorType2& _max) : m_min(_min), m_max(_max) {}
84 
86  template<typename Derived>
87  inline explicit AlignedBox(const MatrixBase<Derived>& a_p)
88  {
89  const typename internal::nested<Derived,2>::type p(a_p.derived());
90  m_min = p;
91  m_max = p;
92  }
93 
95 
97  inline Index dim() const { return AmbientDimAtCompileTime==Dynamic ? m_min.size()-1 : Index(AmbientDimAtCompileTime); }
98 
100  inline bool isNull() const { return isEmpty(); }
101 
103  inline void setNull() { setEmpty(); }
104 
106  inline bool isEmpty() const { return (m_min.array() > m_max.array()).any(); }
107 
109  inline void setEmpty()
110  {
111  m_min.setConstant( ScalarTraits::highest() );
112  m_max.setConstant( ScalarTraits::lowest() );
113  }
114 
116  inline const VectorType& (min)() const { return m_min; }
118  inline VectorType& (min)() { return m_min; }
120  inline const VectorType& (max)() const { return m_max; }
122  inline VectorType& (max)() { return m_max; }
123 
127  center() const
128  { return (m_min+m_max)/2; }
129 
135  { return m_max - m_min; }
136 
138  inline Scalar volume() const
139  { return sizes().prod(); }
140 
146  { return sizes(); }
147 
158  {
159  EIGEN_STATIC_ASSERT(_AmbientDim <= 3, THIS_METHOD_IS_ONLY_FOR_VECTORS_OF_A_SPECIFIC_SIZE);
160 
161  VectorType res;
162 
163  Index mult = 1;
164  for(Index d=0; d<dim(); ++d)
165  {
166  if( mult & corner ) res[d] = m_max[d];
167  else res[d] = m_min[d];
168  mult *= 2;
169  }
170  return res;
171  }
172 
175  inline VectorType sample() const
176  {
177  VectorType r;
178  for(Index d=0; d<dim(); ++d)
179  {
181  {
182  r[d] = m_min[d] + (m_max[d]-m_min[d])
183  * internal::random<Scalar>(Scalar(0), Scalar(1));
184  }
185  else
186  r[d] = internal::random(m_min[d], m_max[d]);
187  }
188  return r;
189  }
190 
192  template<typename Derived>
193  inline bool contains(const MatrixBase<Derived>& a_p) const
194  {
195  typename internal::nested<Derived,2>::type p(a_p.derived());
196  return (m_min.array()<=p.array()).all() && (p.array()<=m_max.array()).all();
197  }
198 
200  inline bool contains(const AlignedBox& b) const
201  { return (m_min.array()<=(b.min)().array()).all() && ((b.max)().array()<=m_max.array()).all(); }
202 
204  template<typename Derived>
206  {
207  typename internal::nested<Derived,2>::type p(a_p.derived());
208  m_min = m_min.cwiseMin(p);
209  m_max = m_max.cwiseMax(p);
210  return *this;
211  }
212 
214  inline AlignedBox& extend(const AlignedBox& b)
215  {
216  m_min = m_min.cwiseMin(b.m_min);
217  m_max = m_max.cwiseMax(b.m_max);
218  return *this;
219  }
220 
222  inline AlignedBox& clamp(const AlignedBox& b)
223  {
224  m_min = m_min.cwiseMax(b.m_min);
225  m_max = m_max.cwiseMin(b.m_max);
226  return *this;
227  }
228 
230  inline AlignedBox intersection(const AlignedBox& b) const
231  {return AlignedBox(m_min.cwiseMax(b.m_min), m_max.cwiseMin(b.m_max)); }
232 
234  inline AlignedBox merged(const AlignedBox& b) const
235  { return AlignedBox(m_min.cwiseMin(b.m_min), m_max.cwiseMax(b.m_max)); }
236 
238  template<typename Derived>
240  {
241  const typename internal::nested<Derived,2>::type t(a_t.derived());
242  m_min += t;
243  m_max += t;
244  return *this;
245  }
246 
251  template<typename Derived>
252  inline Scalar squaredExteriorDistance(const MatrixBase<Derived>& a_p) const;
253 
258  inline Scalar squaredExteriorDistance(const AlignedBox& b) const;
259 
264  template<typename Derived>
267 
272  inline NonInteger exteriorDistance(const AlignedBox& b) const
274 
280  template<typename NewScalarType>
281  inline typename internal::cast_return_type<AlignedBox,
283  {
284  return typename internal::cast_return_type<AlignedBox,
286  }
287 
289  template<typename OtherScalarType>
291  {
292  m_min = (other.min)().template cast<Scalar>();
293  m_max = (other.max)().template cast<Scalar>();
294  }
295 
300  bool isApprox(const AlignedBox& other, RealScalar prec = ScalarTraits::dummy_precision()) const
301  { return m_min.isApprox(other.m_min, prec) && m_max.isApprox(other.m_max, prec); }
302 
303 protected:
304 
306 };
307 
308 
309 
310 template<typename Scalar,int AmbientDim>
311 template<typename Derived>
313 {
314  const typename internal::nested<Derived,2*AmbientDim>::type p(a_p.derived());
315  Scalar dist2(0);
316  Scalar aux;
317  for (Index k=0; k<dim(); ++k)
318  {
319  if( m_min[k] > p[k] )
320  {
321  aux = m_min[k] - p[k];
322  dist2 += aux*aux;
323  }
324  else if( p[k] > m_max[k] )
325  {
326  aux = p[k] - m_max[k];
327  dist2 += aux*aux;
328  }
329  }
330  return dist2;
331 }
332 
333 template<typename Scalar,int AmbientDim>
335 {
336  Scalar dist2(0);
337  Scalar aux;
338  for (Index k=0; k<dim(); ++k)
339  {
340  if( m_min[k] > b.m_max[k] )
341  {
342  aux = m_min[k] - b.m_max[k];
343  dist2 += aux*aux;
344  }
345  else if( b.m_min[k] > m_max[k] )
346  {
347  aux = b.m_min[k] - m_max[k];
348  dist2 += aux*aux;
349  }
350  }
351  return dist2;
352 }
353 
370 #define EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, Size, SizeSuffix) \
371  \
372 typedef AlignedBox<Type, Size> AlignedBox##SizeSuffix##TypeSuffix;
373 
374 #define EIGEN_MAKE_TYPEDEFS_ALL_SIZES(Type, TypeSuffix) \
375 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 1, 1) \
376 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 2, 2) \
377 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 3, 3) \
378 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, 4, 4) \
379 EIGEN_MAKE_TYPEDEFS(Type, TypeSuffix, Dynamic, X)
380 
384 
385 #undef EIGEN_MAKE_TYPEDEFS_ALL_SIZES
386 #undef EIGEN_MAKE_TYPEDEFS
387 
388 } // end namespace Eigen
389 
390 #endif // EIGEN_ALIGNEDBOX_H