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 clearVisited( const std::vector<std::size_t>& visited )
      {
        const auto nb = size();
        for ( auto v : visited )
          {
            myAncestor[ v ] = nb;
            myDistance[ v ] = std::numeric_limits<double>::infinity();
            myVisited [ v ] = 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 ( size_t 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;

    std::vector< Index >
    getCotangentPoints( const Point& a, double max_discrete_distance ) 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 Data ------------------------------
  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 Data --------------------------------
  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)