DGtal  1.3.beta
BoundedRationalPolytope.h
1 
17 #pragma once
18 
31 #if defined(BoundedRationalPolytope_RECURSES)
32 #error Recursive header files inclusion detected in BoundedRationalPolytope.h
33 #else // defined(BoundedRationalPolytope_RECURSES)
34 
35 #define BoundedRationalPolytope_RECURSES
36 
37 #if !defined BoundedRationalPolytope_h
38 
39 #define BoundedRationalPolytope_h
40 
42 // Inclusions
43 #include <iostream>
44 #include <list>
45 #include <vector>
46 #include <string>
47 #include "DGtal/base/Common.h"
48 #include "DGtal/kernel/CSpace.h"
49 #include "DGtal/kernel/domains/HyperRectDomain.h"
50 #include "DGtal/arithmetic/IntegerComputer.h"
51 #include "DGtal/arithmetic/ClosedIntegerHalfPlane.h"
53 
54 namespace DGtal
55 {
56 
58  // template class BoundedRationalPolytope
73  template < typename TSpace >
75  {
77 
78  public:
80  typedef TSpace Space;
81  typedef typename Space::Integer Integer;
82  typedef typename Space::Point Point;
83  typedef typename Space::Vector Vector;
84  typedef std::vector<Vector> InequalityMatrix;
85  typedef std::vector<Integer> InequalityVector;
88 #ifdef WITH_BIGINTEGER
90 #else
91  typedef DGtal::int64_t BigInteger;
92 #endif
94 
99  struct UnitSegment {
101  UnitSegment( Dimension d ) : k( d ) {}
102  };
103 
109  struct UnitCell {
110  std::vector<Dimension> dims;
111  UnitCell( std::initializer_list<Dimension> l )
112  : dims( l.begin(), l.end() ) {}
113 
120  friend std::ostream&
121  operator<< ( std::ostream & out,
122  const UnitCell & object )
123  {
124  out << "{";
125  for ( Dimension i = 0; i < object.dims.size(); ++i ) out << object.dims[ i ];
126  out << "}";
127  return out;
128  }
129  };
130 
131 
134  struct Rational {
135  Integer p; // numerator
136  Integer q; // denominator
140  inline Rational( Integer a, Integer b ) : p( a ), q( b ) {}
141  };
142 
145 
149  ~BoundedRationalPolytope() = default;
150 
155 
160  BoundedRationalPolytope ( const Self & other ) = default;
161 
162 
173  BoundedRationalPolytope( std::initializer_list<Point> l );
174 
188  template <typename PointIterator>
189  BoundedRationalPolytope( Integer d, PointIterator itB, PointIterator itE );
190 
217  template <typename HalfSpaceIterator>
219  const Domain& domain,
220  HalfSpaceIterator itB, HalfSpaceIterator itE,
221  bool valid_edge_constraints = false,
222  bool check_duplicate_constraints = false );
223 
250  template <typename HalfSpaceIterator>
251  void init( Integer d, const Domain& domain,
252  HalfSpaceIterator itB, HalfSpaceIterator itE,
253  bool valid_edge_constraints = false,
254  bool check_duplicate_constraints = false );
255 
273  template <typename PointIterator>
274  bool init( Integer d, PointIterator itB, PointIterator itE );
275 
281  Self & operator= ( const Self & other ) = default;
282 
284  void clear();
285 
287 
288  // ----------------------- Accessor services ------------------------------
289  public:
292 
294  const Domain& getDomain() const;
295 
298  const Domain& getLatticeDomain() const;
299 
302  const Domain& getRationalDomain() const;
303 
305  unsigned int nbHalfSpaces() const;
306 
310  Integer denominator() const;
311 
317  const Vector& getA( unsigned int i ) const;
318 
324  Integer getB( unsigned int i ) const;
325 
331  bool isLarge( unsigned int i ) const;
332 
334  const InequalityMatrix& getA() const;
335 
337  const InequalityVector& getB() const;
341  const std::vector<bool>& getI() const;
342 
347  bool canBeSummed() const;
348 
350 
351  // ----------------------- Check point services ------------------------------
352  public:
353 
356 
361  bool isInside( const Point& p ) const;
362 
369  bool isDomainPointInside( const Point& p ) const;
370 
375  bool isInterior( const Point& p ) const;
376 
381  bool isBoundary( const Point& p ) const;
382 
384 
385  // ----------------------- Modification services ------------------------------
386  public:
387 
390 
394 
407  unsigned int cut( Dimension k, bool pos, Integer b, bool large = true );
408 
426  unsigned int cut( const Vector& a, Integer b, bool large = true,
427  bool valid_edge_constraint = false );
428 
445  unsigned int cut( const HalfSpace & hs, bool large = true,
446  bool valid_edge_constraint = false );
447 
452  void swap( BoundedRationalPolytope & other );
453 
454 
460  Self& operator*=( Integer t );
461 
467  Self& operator*=( Rational r );
468 
475  Self& operator+=( UnitSegment s );
476 
483  Self& operator+=( UnitCell c );
484 
486 
487  // ----------------------- Enumeration services ------------------------------
488  public:
489 
492 
500  Integer count() const;
501 
512  Integer countInterior() const;
513 
524  Integer countBoundary() const;
525 
536  Integer countWithin( Point low, Point hi ) const;
537 
554  Integer countUpTo( Integer max ) const;
555 
564  void getPoints( std::vector<Point>& pts ) const;
565 
574  void getInteriorPoints( std::vector<Point>& pts ) const;
575 
584  void getBoundaryPoints( std::vector<Point>& pts ) const;
585 
596  template <typename PointSet>
597  void insertPoints( PointSet& pts_set ) const;
598 
600 
601 
602  // ----------------------- Interface --------------------------------------
603  public:
606 
611  void selfDisplay ( std::ostream & out ) const;
612 
619  bool isValid() const;
620 
624  std::string className() const;
625 
627 
628  // ------------------------- Protected Datas ------------------------------
629  protected:
630  // Denominator for constraints, i.e. \f$ q A x \le B \f$.
632  // The matrix A in the polytope representation \f$ q A x \le B \f$.
634  // The vector B in the polytope representation \f$ q A x \le B \f$.
636  // Tight bounded box
638  // Lattice bounded box (i.e. floor( rationalD/q ))
640  // Are inequalities large ?
641  std::vector<bool> I;
642  // Indicates if Minkowski sums with segments will be valid
644 
645  // ------------------------- Private Datas --------------------------------
646  private:
647 
648 
649  // ------------------------- Internals ------------------------------------
650  private:
658 
665 
672 
679  Domain computeLatticeDomain( const Domain& d );
680 
686  Domain computeRationalDomain( const Domain& d );
687 
688  }; // end of class BoundedRationalPolytope
689 
690  namespace detail {
699  template <DGtal::Dimension N, typename TInteger>
701  typedef TInteger Integer;
703  typedef typename Space::Point Point;
704  typedef typename Space::Vector Vector;
707 
720  static void
721  addEdgeConstraint( Polytope& , unsigned int , unsigned int ,
722  const std::vector<Point>& )
723  {
724  trace.error() << "[BoundedRationalPolytopeHelper::addEdgeConstraint]"
725  << " this method is only implemented in 3D." << std::endl;
726  }
727 
730  static
731  Vector crossProduct( const Vector& , const Vector& )
732  {
733  trace.error() << "[BoundedRationalPolytopeHelper::crossProduct]"
734  << " this method is only implemented in 3D." << std::endl;
735  return Vector::zero;
736  }
737  };
738 
746  template <typename TInteger>
748  typedef TInteger Integer;
750  typedef typename Space::Point Point;
751  typedef typename Space::Vector Vector;
754 
764  static void
765  addEdgeConstraint( Polytope& P, unsigned int i, unsigned int j,
766  const std::vector<Point>& pts )
767  {
768  Vector ab = pts[ i ] - pts[ j ];
769  for ( int s = 0; s < 2; s++ )
770  for ( Dimension k = 0; k < dimension; ++k )
771  {
772  Vector n = ab.crossProduct( Point::base( k, (s == 0) ? 1 : -1 ) );
773  Integer b = n.dot( pts[ i ] );
774  std::size_t nb_in = 0;
775  for ( auto p : pts ) {
776  Integer v = n.dot( p );
777  if ( v < b ) nb_in++;
778  }
779  if ( nb_in == dimension - 1 ) {
780  P.cut( n, b, true, true );
781  }
782  }
783  }
788  static
789  Vector crossProduct( const Vector& v1, const Vector& v2 )
790  {
791  return v1.crossProduct( v2 );
792  }
793  };
794  }
795 
798 
805  template <typename TSpace>
806  std::ostream&
807  operator<< ( std::ostream & out,
808  const BoundedRationalPolytope<TSpace> & object );
809 
810 
816  template <typename TSpace>
817  BoundedRationalPolytope<TSpace>
819  const BoundedRationalPolytope<TSpace> & P );
820 
826  template <typename TSpace>
827  BoundedRationalPolytope<TSpace>
828  operator* ( typename BoundedRationalPolytope<TSpace>::Rational r,
829  const BoundedRationalPolytope<TSpace> & P );
830 
831 
839  template <typename TSpace>
840  BoundedRationalPolytope<TSpace>
841  operator+ ( const BoundedRationalPolytope<TSpace> & P,
842  typename BoundedRationalPolytope<TSpace>::UnitSegment s );
843 
851  template <typename TSpace>
852  BoundedRationalPolytope<TSpace>
853  operator+ ( const BoundedRationalPolytope<TSpace> & P,
854  typename BoundedRationalPolytope<TSpace>::UnitCell c );
855 
857 
858 } // namespace DGtal
859 
860 
862 // Includes inline functions.
863 #include "BoundedRationalPolytope.ih"
864 
865 // //
867 
868 #endif // !defined BoundedRationalPolytope_h
869 
870 #undef BoundedRationalPolytope_RECURSES
871 #endif // else defined(BoundedRationalPolytope_RECURSES)
DGtal::BigInteger
mpz_class BigInteger
Multi-precision integer with GMP implementation.
Definition: BasicTypes.h:79
DGtal::BoundedRationalPolytope::isInterior
bool isInterior(const Point &p) const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Space
SpaceND< 3, Integer > Space
Definition: BoundedRationalPolytope.h:749
DGtal::PointVector< dim, Integer >::zero
static Self zero
Static const for zero PointVector.
Definition: PointVector.h:1595
DGtal::detail::BoundedRationalPolytopeSpecializer
Aim: It is just a helper class for BoundedRationalPolytope to add dimension specific static methods.
Definition: BoundedRationalPolytope.h:700
DGtal::BoundedRationalPolytope::Rational::p
Integer p
Definition: BoundedRationalPolytope.h:135
DGtal::BoundedRationalPolytope::Rational::Rational
Rational(Integer a, Integer b)
Definition: BoundedRationalPolytope.h:140
DGtal::HyperRectDomain< Space >
DGtal::BoundedRationalPolytope::isLarge
bool isLarge(unsigned int i) const
DGtal::concepts::CSpace
Aim: Defines the concept describing a digital space, ie a cartesian product of integer lines.
Definition: CSpace.h:105
DGtal::BoundedRationalPolytope::init
void init(Integer d, const Domain &domain, HalfSpaceIterator itB, HalfSpaceIterator itE, bool valid_edge_constraints=false, bool check_duplicate_constraints=false)
max
int max(int a, int b)
Definition: testArithmeticalDSS.cpp:1108
DGtal::BoundedRationalPolytope::operator=
Self & operator=(const Self &other)=default
DGtal::BoundedRationalPolytope::internalInitFromSegment3D
bool internalInitFromSegment3D(Point a, Point b)
DGtal::BoundedRationalPolytope::getInteriorPoints
void getInteriorPoints(std::vector< Point > &pts) const
DGtal::Trace::error
std::ostream & error()
DGtal::BoundedRationalPolytope::UnitCell::UnitCell
UnitCell(std::initializer_list< Dimension > l)
Definition: BoundedRationalPolytope.h:111
DGtal::BoundedRationalPolytope::isDomainPointInside
bool isDomainPointInside(const Point &p) const
DGtal::BoundedRationalPolytope::Vector
Space::Vector Vector
Definition: BoundedRationalPolytope.h:83
DGtal::BoundedRationalPolytope::latticeD
Domain latticeD
Definition: BoundedRationalPolytope.h:639
DGtal::BoundedRationalPolytope::dimension
static const Dimension dimension
Definition: BoundedRationalPolytope.h:93
DGtal::BoundedRationalPolytope::BOOST_CONCEPT_ASSERT
BOOST_CONCEPT_ASSERT((concepts::CSpace< TSpace >))
DGtal::BoundedRationalPolytope::countUpTo
Integer countUpTo(Integer max) const
DGtal::BoundedRationalPolytope::I
std::vector< bool > I
Definition: BoundedRationalPolytope.h:641
DGtal::trace
Trace trace
Definition: Common.h:154
DGtal::BoundedRationalPolytope::Space
TSpace Space
Definition: BoundedRationalPolytope.h:80
DGtal::BoundedRationalPolytope::countInterior
Integer countInterior() const
DGtal::Dimension
DGtal::uint32_t Dimension
Definition: Common.h:137
DGtal::BoundedRationalPolytope::getA
const InequalityMatrix & getA() const
DGtal::BoundedRationalPolytope::isInside
bool isInside(const Point &p) const
DGtal::BoundedRationalPolytope::canBeSummed
bool canBeSummed() const
DGtal::detail::BoundedRationalPolytopeSpecializer::addEdgeConstraint
static void addEdgeConstraint(Polytope &, unsigned int, unsigned int, const std::vector< Point > &)
Definition: BoundedRationalPolytope.h:721
DGtal::BoundedRationalPolytope::Self
BoundedRationalPolytope< TSpace > Self
Definition: BoundedRationalPolytope.h:79
DGtal::detail::BoundedRationalPolytopeSpecializer::dimension
static const Dimension dimension
Definition: BoundedRationalPolytope.h:706
DGtal::BoundedRationalPolytope::nbHalfSpaces
unsigned int nbHalfSpaces() const
DGtal::detail::BoundedRationalPolytopeSpecializer::Integer
TInteger Integer
Definition: BoundedRationalPolytope.h:701
DGtal::operator<<
std::ostream & operator<<(std::ostream &out, const ATu0v1< TKSpace, TLinearAlgebra > &object)
DGtal::BoundedRationalPolytope::BoundedRationalPolytope
BoundedRationalPolytope()
DGtal::BoundedRationalPolytope::Domain
HyperRectDomain< Space > Domain
Definition: BoundedRationalPolytope.h:86
DGtal::PointVector< dim, Integer >::base
static Self base(Dimension k, Component val=1)
DGtal::BoundedRationalPolytope::InequalityVector
std::vector< Integer > InequalityVector
Definition: BoundedRationalPolytope.h:85
DGtal::SpaceND
Definition: SpaceND.h:95
DGtal::PointVector::crossProduct
auto crossProduct(const PointVector< dim, OtherComponent, OtherStorage > &v) const -> decltype(DGtal::crossProduct(*this, v))
Cross product with a PointVector.
DGtal::BoundedRationalPolytope::UnitCell::operator<<
friend std::ostream & operator<<(std::ostream &out, const UnitCell &object)
Definition: BoundedRationalPolytope.h:121
DGtal::BoundedRationalPolytope::InequalityMatrix
std::vector< Vector > InequalityMatrix
Definition: BoundedRationalPolytope.h:84
DGtal::BoundedRationalPolytope::getLatticeDomain
const Domain & getLatticeDomain() const
DGtal::BoundedRationalPolytope::selfDisplay
void selfDisplay(std::ostream &out) const
DGtal::detail::BoundedRationalPolytopeSpecializer::crossProduct
static Vector crossProduct(const Vector &, const Vector &)
Definition: BoundedRationalPolytope.h:731
DGtal::BoundedRationalPolytope::computeLatticeDomain
Domain computeLatticeDomain(const Domain &d)
DGtal::operator+
Circulator< TIterator > operator+(typename IteratorCirculatorTraits< TIterator >::Difference d, Circulator< TIterator > &object)
Definition: Circulator.h:453
DGtal::detail::BoundedRationalPolytopeSpecializer::Space
SpaceND< N, Integer > Space
Definition: BoundedRationalPolytope.h:702
DGtal
DGtal is the top-level namespace which contains all DGtal functions and types.
DGtal::BoundedRationalPolytope::swap
void swap(BoundedRationalPolytope &other)
DGtal::detail::BoundedRationalPolytopeSpecializer::Polytope
BoundedRationalPolytope< Space > Polytope
Definition: BoundedRationalPolytope.h:705
Domain
HyperRectDomain< Space > Domain
Definition: testSimpleRandomAccessRangeFromPoint.cpp:44
DGtal::BoundedRationalPolytope::Rational::q
Integer q
Definition: BoundedRationalPolytope.h:136
DGtal::BoundedRationalPolytope::getBoundaryPoints
void getBoundaryPoints(std::vector< Point > &pts) const
DGtal::BoundedRationalPolytope::getDomain
const Domain & getDomain() const
DGtal::BoundedRationalPolytope::UnitSegment::k
Dimension k
Definition: BoundedRationalPolytope.h:100
DGtal::BoundedRationalPolytope::rationalD
Domain rationalD
Definition: BoundedRationalPolytope.h:637
DGtal::SpaceND::dimension
static const Dimension dimension
static constants to store the dimension.
Definition: SpaceND.h:132
DGtal::BoundedRationalPolytope::UnitCell::dims
std::vector< Dimension > dims
Definition: BoundedRationalPolytope.h:110
DGtal::BoundedRationalPolytope::Point
Space::Point Point
Definition: BoundedRationalPolytope.h:82
DGtal::BoundedRationalPolytope::UnitSegment::UnitSegment
UnitSegment(Dimension d)
Definition: BoundedRationalPolytope.h:101
DGtal::BoundedRationalPolytope::~BoundedRationalPolytope
~BoundedRationalPolytope()=default
DGtal::detail::BoundedRationalPolytopeSpecializer::Vector
Space::Vector Vector
Definition: BoundedRationalPolytope.h:704
DGtal::BoundedRationalPolytope::Integer
Space::Integer Integer
Definition: BoundedRationalPolytope.h:81
DGtal::BoundedRationalPolytope::B
InequalityVector B
Definition: BoundedRationalPolytope.h:635
DGtal::BoundedRationalPolytope::denominator
Integer denominator() const
Vector
FreemanChain< int >::Vector Vector
Definition: testCombinDSS.cpp:60
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::crossProduct
static Vector crossProduct(const Vector &v1, const Vector &v2)
Definition: BoundedRationalPolytope.h:789
DGtal::BoundedRationalPolytope::count
Integer count() const
DGtal::BoundedRationalPolytope::Rational
Definition: BoundedRationalPolytope.h:134
DGtal::BoundedRationalPolytope::insertPoints
void insertPoints(PointSet &pts_set) const
Integer
Point::Coordinate Integer
Definition: examplePlaneProbingParallelepipedEstimator.cpp:44
DGtal::ClosedIntegerHalfPlane
Aim: A half-space specified by a vector N and a constant c. The half-space is the set .
Definition: ClosedIntegerHalfPlane.h:63
DGtal::BoundedRationalPolytope
Aim: Represents an nD rational polytope, i.e. a convex polyhedron bounded by vertices with rational c...
Definition: BoundedRationalPolytope.h:74
DGtal::BoundedRationalPolytope::BigInteger
DGtal::BigInteger BigInteger
Definition: BoundedRationalPolytope.h:89
DGtal::BoundedRationalPolytope::isBoundary
bool isBoundary(const Point &p) const
DGtal::BoundedRationalPolytope::className
std::string className() const
DGtal::BoundedRationalPolytope::countBoundary
Integer countBoundary() const
DGtal::BoundedRationalPolytope::HalfSpace
ClosedIntegerHalfPlane< Space > HalfSpace
Definition: BoundedRationalPolytope.h:87
DGtal::BoundedRationalPolytope::operator*=
Self & operator*=(Integer t)
DGtal::detail::BoundedRationalPolytopeSpecializer::Point
Space::Point Point
Definition: BoundedRationalPolytope.h:703
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Integer
TInteger Integer
Definition: BoundedRationalPolytope.h:748
domain
Domain domain
Definition: testProjection.cpp:88
DGtal::PointVector< dim, Integer >
DGtal::BoundedRationalPolytope::internalInitFromTriangle3D
bool internalInitFromTriangle3D(Point a, Point b, Point c)
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Vector
Space::Vector Vector
Definition: BoundedRationalPolytope.h:751
DGtal::BoundedRationalPolytope::countWithin
Integer countWithin(Point low, Point hi) const
DGtal::BoundedRationalPolytope::getPoints
void getPoints(std::vector< Point > &pts) const
DGtal::BoundedRationalPolytope::operator+=
Self & operator+=(UnitSegment s)
DGtal::BoundedRationalPolytope::UnitCell
Definition: BoundedRationalPolytope.h:109
DGtal::operator*
KForm< Calculus, order, duality > operator*(const typename Calculus::Scalar &scalar, const KForm< Calculus, order, duality > &form)
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Point
Space::Point Point
Definition: BoundedRationalPolytope.h:750
DGtal::int64_t
boost::int64_t int64_t
signed 94-bit integer.
Definition: BasicTypes.h:74
DGtal::BoundedRationalPolytope::isValid
bool isValid() const
DGtal::BoundedRationalPolytope::q
Integer q
Definition: BoundedRationalPolytope.h:631
DGtal::BoundedRationalPolytope::A
InequalityMatrix A
Definition: BoundedRationalPolytope.h:633
DGtal::BoundedRationalPolytope::getI
const std::vector< bool > & getI() const
DGtal::BoundedRationalPolytope::internalInitFromSegment2D
bool internalInitFromSegment2D(Point a, Point b)
DGtal::BoundedRationalPolytope::getRationalDomain
const Domain & getRationalDomain() const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::Polytope
BoundedRationalPolytope< Space > Polytope
Definition: BoundedRationalPolytope.h:752
DGtal::BoundedRationalPolytope::cut
unsigned int cut(Dimension k, bool pos, Integer b, bool large=true)
DGtal::BoundedRationalPolytope::computeRationalDomain
Domain computeRationalDomain(const Domain &d)
Point
MyPointD Point
Definition: testClone2.cpp:383
DGtal::BoundedRationalPolytope::myValidEdgeConstraints
bool myValidEdgeConstraints
Definition: BoundedRationalPolytope.h:643
DGtal::BoundedRationalPolytope::UnitSegment
Definition: BoundedRationalPolytope.h:99
DGtal::BoundedRationalPolytope::interiorPolytope
BoundedRationalPolytope interiorPolytope() const
DGtal::SpaceND::Integer
TInteger Integer
Arithmetic ring induced by (+,-,*) and Integer numbers.
Definition: SpaceND.h:102
DGtal::BoundedRationalPolytope::clear
void clear()
Clears the polytope.
DGtal::BoundedRationalPolytope::getB
const InequalityVector & getB() const
DGtal::detail::BoundedRationalPolytopeSpecializer< 3, TInteger >::addEdgeConstraint
static void addEdgeConstraint(Polytope &P, unsigned int i, unsigned int j, const std::vector< Point > &pts)
Definition: BoundedRationalPolytope.h:765