ReebGraphConstruction.h

Go to the documentation of this file.

#ifndef PMPL_REEB_GRAPH_CONSTRUCTION_H_
#define PMPL_REEB_GRAPH_CONSTRUCTION_H_
 
#include <vector>
#include <set>
#include <unordered_set>
using namespace std;
 
#include <boost/functional/hash.hpp>
 
#include <containers/sequential/graph/directed_preds_graph.h>
 
#include <Vector.h>
using namespace mathtool;
 
#include "Utilities/IOUtils.h"
 
class Environment;
class WorkspaceDecomposition;
class WorkspaceSkeleton;
class XMLNode;
 
 
class ReebGraphConstruction {
 
  public:
 
 
    struct Parameters {
      std::string filename;  
      bool write{false};     
      bool read{false};     
      bool debug{false};     
    };
 
    static Parameters m_defaultParams; 
 
    typedef tuple<size_t, size_t, size_t> Triangle; 
 
    struct ReebNode {
 
      ReebNode(size_t _vIndx = -1, const Vector3d& _v = Vector3d(),
          double _w = numeric_limits<double>::max()) :
          m_vertexIndex(_vIndx), m_vertex(_v), m_w(_w) {}
 
      friend istream& operator>>(istream& _is, ReebNode& _rn);
      friend ostream& operator<<(ostream& _os, const ReebNode& _rn);
 
      size_t m_vertexIndex; 
      Vector3d m_vertex;    
      double m_w;           
      size_t m_order{0};    
      size_t m_tetra{0};    
    };
 
    struct ReebNodeComp {
 
      bool operator()(const ReebNode& _a, const ReebNode& _b) {
 
        auto vertComp = [](const Vector3d _a, const Vector3d _b) {
          return _a[0] - _b[0] > ReebTolerance() ||
            (fabs(_a[0] - _b[0]) < ReebTolerance() && (_a[1] - _b[1] > ReebTolerance() ||
            (fabs(_a[1] - _b[1]) < ReebTolerance() && _a[2] - _b[2] > ReebTolerance())));
        };
 
        const Vector3d& a = _a.m_vertex;
        const Vector3d& b = _b.m_vertex;
        bool ret = _b.m_vertexIndex != _a.m_vertexIndex &&
          (_a.m_w - _b.m_w > ReebTolerance() ||
           (fabs(_a.m_w - _b.m_w) < ReebTolerance() &&
            (vertComp(a, b) || (a == b &&
            _a.m_vertexIndex < _b.m_vertexIndex))));
        return ret;
      }
 
      static constexpr double ReebTolerance() {return 0.000001;}
    };
 
    struct MeshEdge;
 
    struct ReebArc {
 
      ReebArc(size_t _s = -1, size_t _t = -1, MeshEdge* _m = nullptr) :
        m_source(_s), m_target(_t) {
          if(_m)
            m_edges.insert(_m);
        }
 
      friend istream& operator>>(istream& _is, ReebArc& _ra);
      friend ostream& operator<<(ostream& _os, const ReebArc& _ra);
 
      size_t m_source;                  
      size_t m_target;                  
      unordered_set<MeshEdge*> m_edges; 
      unordered_set<size_t> m_tetra;    
      vector<Vector3d> m_path;          
    };
 
    typedef stapl::sequential::directed_preds_graph<
      stapl::MULTIEDGES, ReebNode, ReebArc>
      ReebGraph; 
    typedef ReebGraph::edge_descriptor RGEID; 
 
    struct ReebArcComp {
 
      ReebArcComp(ReebGraph* _rg) : m_rg(_rg) {}
 
      bool operator()(const RGEID& _a, const RGEID& _b) {
        bool ret = false;
        ReebNode& s0 = m_rg->find_vertex(_a.source())->property();
        ReebNode& s1 = m_rg->find_vertex(_b.source())->property();
        if(s0.m_order < s1.m_order)
          ret = true;
        else if(s0.m_order == s1.m_order) {
          ReebNode& t0 = m_rg->find_vertex(_a.target())->property();
          ReebNode& t1 = m_rg->find_vertex(_b.target())->property();
          if(t0.m_order < t1.m_order)
            ret = true;
          else if(t0.m_order == t1.m_order)
            if(_a.id() < _b.id())
              ret = true;
        }
        return ret;
      }
 
      ReebGraph* m_rg; 
    };
 
    typedef set<RGEID, ReebArcComp> ArcSet; 
 
    struct MeshEdge {
 
      MeshEdge(size_t _s, size_t _t, ReebGraph* _rg) :
        m_source(_s), m_target(_t), m_numTris(0), m_arcs(ReebArcComp(_rg)) {
        }
 
      size_t m_source;  
      size_t m_target;  
      size_t m_numTris; 
      ArcSet m_arcs;    
    };
 
    struct MeshEdgeHash {
 
      size_t operator()(const MeshEdge* const _m) const {
        return h(make_pair(_m->m_source, _m->m_target));
      }
 
      boost::hash<pair<size_t, size_t>> h; 
    };
 
    struct MeshEdgeEq {
 
      bool operator()(const MeshEdge* const _a, const MeshEdge* const _b) const {
        return _a->m_source == _b->m_source && _a->m_target == _b->m_target;
      }
    };
 
    typedef stapl::sequential::directed_preds_graph<
      stapl::MULTIEDGES, Vector3d, vector<Vector3d>
        > FlowGraph; 
 
 
    ReebGraphConstruction();
 
    static void SetDefaultParameters(XMLNode& _node);
 
 
    void Construct(const WorkspaceDecomposition* _decomposition);
 
    WorkspaceSkeleton GetSkeleton();
 
 
    ReebGraph& GetReebGraph() {return m_reebGraph;}
 
 
    void Read(const string& _filename);
 
    void Write(const string& _filename);
 
 
  private:
 
 
    void Initialize(const WorkspaceDecomposition* _decomposition);
 
    void Construct();
 
    void CreateNode(size_t _i, const Vector3d& _v, double _w);
 
    MeshEdge* CreateArc(size_t _s, size_t _t,
        const unordered_set<size_t>& _tetra = unordered_set<size_t>());
 
    void MergePaths(MeshEdge* _e0, MeshEdge* _e1, MeshEdge* _e2);
 
    void GlueByMergeSorting(ArcSet& _a0, MeshEdge* _e0,
        ArcSet& _a1, MeshEdge* _e1);
 
    void MergeArcs(RGEID _a0, RGEID _a1);
 
    ReebArc& GetReebArc(RGEID _a);
 
    void DeleteMeshEdge(MeshEdge* _e);
 
    void Remove2Nodes();
 
    void Embed(const WorkspaceDecomposition* _tetrahedralization);
 
    void Uniqueify();
 
 
    Parameters m_params; 
 
    vector<Vector3d> m_vertices; 
 
    unordered_set<MeshEdge*, MeshEdgeHash, MeshEdgeEq> m_edges;
 
    vector<pair<Triangle, unordered_set<size_t>>> m_triangles;
 
    ReebGraph m_reebGraph;       
 
};
 
#endif