TangencyComputer.h

 
#pragma once
 
#if defined(TangencyComputer_RECURSES)
#error Recursive header files inclusion detected in TangencyComputer.h
#else // defined(TangencyComputer_RECURSES)
#define TangencyComputer_RECURSES
 
#if !defined TangencyComputer_h
#define TangencyComputer_h
 
// Inclusions
#include <iostream>
#include <vector>
#include <string>
#include <limits>
#include <unordered_set>
#include "DGtal/base/Common.h"
#include "DGtal/base/Clone.h"
#include "DGtal/kernel/domains/HyperRectDomain.h"
#include "DGtal/topology/CCellularGridSpaceND.h"
#include "DGtal/kernel/LatticeSetByIntervals.h"
#include "DGtal/geometry/volumes/CellGeometry.h"
#include "DGtal/geometry/volumes/DigitalConvexity.h"
 
namespace DGtal
{
 
  // template class TangencyComputer
  template < typename TKSpace >
  class TangencyComputer
  {
    BOOST_CONCEPT_ASSERT(( concepts::CCellularGridSpaceND< TKSpace > ));
 
  public:
    typedef TangencyComputer< TKSpace > Self;
    typedef TKSpace                     KSpace;
    typedef typename KSpace::Space      Space;
    typedef typename KSpace::Point      Point;
    typedef typename KSpace::Vector     Vector;
    typedef HyperRectDomain< Space >    Domain;
    typedef std::size_t                 Index;
    typedef std::size_t                 Size;
    typedef std::vector< Index >        Path;
    typedef CellGeometry< KSpace >      CellCover;
    typedef LatticeSetByIntervals< Space > LatticeCellCover;
    
    // ------------------------- Shortest path services --------------------------------
  public:
 
    struct ShortestPaths {
      typedef std::tuple< Index, Index, double > Node;
 
      struct Comparator {
        bool operator()  ( const Node& p1,
                           const Node& p2 ) const
        {
          return std::get<2>( p1 ) > std::get<2>( p2 );
        }
      };
 
      ShortestPaths()
        : myTgcyComputer( nullptr ), mySecure( 0.0 )
      {}
 
      ShortestPaths( const ShortestPaths& other ) = default;
      ShortestPaths( ShortestPaths&& other ) = default;
      ShortestPaths& operator=( const ShortestPaths& other ) = default;
      ShortestPaths& operator=( ShortestPaths&& other ) = default;
      
      ShortestPaths( ConstAlias< TangencyComputer > tgcy_computer,
                     double secure = sqrt( KSpace::dimension ) )
        : myTgcyComputer( &tgcy_computer ),
          mySecure( std::max( secure, 0.0 ) )
      {
        clear();
      }
 
      const TangencyComputer* tangencyComputerPtr() const
      {
        return myTgcyComputer;
      }
      
      Index size() const
      {
        return myTgcyComputer->size();
      }
      
      void clear()
      {
        const auto nb = size();
        myAncestor = std::vector< Index > ( nb, nb );
        myDistance = std::vector< double >( nb, std::numeric_limits<double>::infinity() );
        myVisited  = std::vector< bool >  ( nb, false );
        myQ        = std::priority_queue< Node, std::vector< Node >, Comparator >();
  }
 
      void init( Index i )
      {
        ASSERT( i < size() );
        myQ.push( std::make_tuple( i, i, 0.0 ) );
        myAncestor[ i ] = i;
        myDistance[ i ] = 0.0;
        myVisited [ i ] = true;
      }
                     
      template < typename IndexFwdIterator >
      void init( IndexFwdIterator it, IndexFwdIterator itE )
      {
        for ( ; it != itE; ++it )
          {
            const auto i = *it;
            ASSERT( i < size() );
            myQ.push( std::make_tuple( i, i, 0.0 ) );
          }
        const auto elem = myQ.top();
        const auto i    = std::get<0>( elem );
        myAncestor[ i ] = i;
        myDistance[ i ] = 0.0;
        myVisited [ i ] = true;
      }
 
      bool finished() const
      {
        return myQ.empty();
      }
 
      const Node& current() const
      {
        ASSERT( ! finished() );
        return myQ.top();
      }
 
      void expand();
 
      
      bool isValid() const
      {
        return myTgcyComputer != nullptr;
      }
 
      const Point& point( Index i ) const
      {
        ASSERT( i < size() );
        return myTgcyComputer->point( i );
      }
 
      Index ancestor( Index i ) const
      {
        ASSERT( i < size() );
        return myAncestor[ i ];
      }
 
      double distance( Index i ) const
      {
        ASSERT( i < size() );
        return myDistance[ i ];
      }
 
      bool isVisited( Index i ) const
      {
        return ancestor( i ) < size();
      }
 
      static double infinity()
      {
        return std::numeric_limits<double>::infinity();
      }
 
      Path pathToSource( Index i ) const
      {
        Path P;
        if ( ! isVisited( i ) ) return P;
        P.push_back( i );
        while ( ancestor( i ) != i )
          {
            i = ancestor( i );
            P.push_back( i );
          }
        return P;
      }
 
      const std::vector< Index >& ancestors() const
      { return myAncestor; }
      
      const std::vector< double >& distances() const
      { return myDistance; }
      
      const std::vector< bool >& visitedPoints() const
      { return myVisited; }
      
    protected:
      const TangencyComputer* myTgcyComputer;
      double                  mySecure;
      std::vector< Index >    myAncestor;
      std::vector< double >   myDistance;
      std::vector< bool >     myVisited;
      std::priority_queue< Node, std::vector< Node >, Comparator > myQ;
 
    protected:
 
      void propagate( Index current );
      
      std::vector< Index >
      getCotangentPoints( Index i ) const;
 
    };
 
    friend struct ShortestPaths;
    
    // ------------------------- Standard services --------------------------------
  public:
 
    TangencyComputer() = default;
 
    TangencyComputer( const Self& other ) = default;
 
    TangencyComputer( Self&& other ) = default;
 
    Self& operator=( const Self& other ) = default;
 
    Self& operator=( Self&& other ) = default;
    
    TangencyComputer( Clone< KSpace > aK );
 
    template < typename PointIterator >
    void init( PointIterator itB, PointIterator itE,
               bool use_lattice_cell_cover = false );
 
 
    // ------------------------- Accessors services --------------------------------
  public:
    
    const KSpace& space() const
    { return myK; }
 
    Size size() const
    { return myX.size(); }
    
    const std::vector< Point >& points() const
      { return myX; }
      
    const Point& point( Index i ) const
    { return myX[ i ]; }
      
    Size index( const Point& a ) const
    {
      const auto p = myPt2Index.find( a );
      return p == myPt2Index.cend() ? size() : p->second;
    }
      
    const CellCover& cellCover() const
    { return myCellCover; }
 
    const LatticeCellCover& latticeCellCover() const
    { return myLatticeCellCover; }
 
    double length( const Path& path ) const
    {
      auto eucl_d = [] ( const Point& p, const Point& q )
      { return ( p - q ).norm(); };
      double l = 0.0;
      for ( auto i = 1; i < path.size(); i++ )
        l += eucl_d( point( path[ i-1 ] ), point( path[ i ] ) );
      return l;
    }
    
 
    // ------------------------- Tangency services --------------------------------
  public:
    
    bool arePointsCotangent( const Point& a, const Point& b ) const;
 
    bool arePointsCotangent( const Point& a, const Point& b, const Point& c ) const;
 
    std::vector< Index >
    getCotangentPoints( const Point& a ) const;
    
    std::vector< Index >
    getCotangentPoints( const Point& a,
                        const std::vector< bool > & to_avoid ) const;
    
    
    // ------------------------- Shortest paths services --------------------------------
  public:
    
    ShortestPaths
    makeShortestPaths( double secure = sqrt( KSpace::dimension ) ) const;
        
    std::vector< Path >
    shortestPaths( const std::vector< Index >& sources,
                   const std::vector< Index >& targets,
                   double secure = sqrt( KSpace::dimension ),
                   bool verbose = false ) const;
    
    Path
    shortestPath( Index source, Index target,
                  double secure = sqrt( KSpace::dimension ),
                  bool verbose = false ) const;
    
    
    // ------------------------- Protected Datas ------------------------------
  protected:
    
    KSpace myK;
    DigitalConvexity< KSpace > myDConv;
    std::vector< Vector > myN;
    std::vector< double > myDN;
    std::vector< Point >  myX;
    bool myUseLatticeCellCover;
    CellCover myCellCover;
    LatticeCellCover myLatticeCellCover;
    
    std::unordered_map< Point, Index > myPt2Index;
    
    // ------------------------- Private Datas --------------------------------
  private:
 
 
    // ------------------------- Internals ------------------------------------
  private:
 
    void setUp();
    
  }; // end of class TangencyComputer
 
 
  template <typename TKSpace>
  std::ostream&
  operator<< ( std::ostream & out,
               const TangencyComputer<TKSpace> & object );
 
 
} // namespace DGtal
 
 
// Includes inline functions.
#include "TangencyComputer.ih"
 
//                                                                           //
 
#endif // !defined TangencyComputer_h
 
#undef TangencyComputer_RECURSES
#endif // else defined(TangencyComputer_RECURSES)