/Users/mourrain/Devel/mmx/shape/include/shape/ssi_dsearch.hpp

Go to the documentation of this file.

#ifndef SYNAPS_SHAPE_GAIA_DSEARCH_H
#define SYANSP_SHAPE_GAIA_DSEARCH_H
#include <shape/ssi_def.hpp>
#include <shape/ssi_qsegment.hpp>
#include <list>
#include <queue>
#include <shape/ssi_ordered_pair.hpp>

#define ParametricSurface shape::surface_parametric<double>

namespace mmx {
namespace shape_ssi {
  
  typedef std::list< point2 >  list2_t; 
  typedef std::list< point3 >  list3_t;
  typedef  list2_t::iterator p2ref_t;
  
  struct dcurve
  {      
    unsigned size;
    dcurve  *   owner; 
    list2_t     left,right;
    list3_t     inter;
    dcurve( const point2& l, const point2& r, const point3& i ) ;
    dcurve( const point2& l0, const point2& l1, const point2& r0, const point2& r1, const point3& i0, const point3& i1 );
  };

  template<class DPoint>
    struct dim_cmp   
    { 
      unsigned dim; 
      dim_cmp(unsigned d = 0): dim(d) {}; 
      dim_cmp next() const { return dim_cmp((dim+1)%4); };
      bool operator()(  DPoint*const& a, DPoint *const& b ) const { return (*a)[dim] < (*b)[dim]; };
      bool operator()( const DPoint& a, const DPoint&  b ) const { return    a[dim] < b[dim];    };
    };

  template<class DP>
    inline double distance( DP * pa, DP * pb )
    {
      double acc = 0.0;
      for ( unsigned i = 0; i < 4; i++ )
        acc += ((*pa)[i]-(*pb)[i])*((*pa)[i]-(*pb)[i]);
      return acc;
    };

  struct sdpoint_t
  {
    typedef double value_type;
    dcurve  *  curve;
    p2ref_t idl;
    p2ref_t idr;
  
    inline const double& operator[]( unsigned i ) const
    {
      if ( i < 2 ) return (*idl)[i];
      return (*idr)[i-2];
    };
  
    inline double& operator[]( unsigned i ) 
    {
      if ( i < 2 ) return (*idl)[i];
      return (*idr)[i-2];
    };
  };

  struct pdpoint_t
  {
    typedef double value_type;
    dcurve * curve;
    p2ref_t ref;
    pdpoint_t * a;

    inline value_type& operator[]( unsigned i ) 
    { 
      if ( i < 2 ) return (*ref)[i];
      return (*(a->ref))[i-2];
    };
    inline const value_type& operator[]( unsigned i ) const 
    { 
      if ( i < 2 ) return (*ref)[i];
      return (*(a->ref))[i-2];
    };
  };

  void append( dcurve *  a, dcurve   * b ) ;
  void prepend( dcurve * a, dcurve * b ) ;
  void reverse( dcurve * c ) ;


  struct sdknode
  {
    sdpoint_t   * data; 
    sdknode * l, * r; 
    sdknode( sdpoint_t * d, sdknode * _l, sdknode * _r ) : data(d), l(_l), r(_r){}; 
    ~sdknode() { if ( l ) delete l; if ( r ) delete r; };
  };

  struct pdknode
  {
    pdpoint_t   * data; 
    pdknode * l, * r; 
    pdknode( pdpoint_t * d, pdknode * _l, pdknode * _r ) : data(d), l(_l), r(_r){}; 
    ~pdknode() { if ( l ) delete l; if ( r ) delete r; };
  };



  inline  bool   empty( dcurve * c ) { return c->owner != c; };

  inline  pdpoint_t * sibble( pdpoint_t * p )  { return p->a; };
  void extremums( sdpoint_t* dst, dcurve* src );
  void extremums( pdpoint_t * dst, dcurve * c );
  dcurve * curve( pdpoint_t * p );
  dcurve *  curve( sdpoint_t * p );


  inline  bool isfirst ( sdpoint_t* p ) { return curve(p)->left.begin() == p->idl;     };
  inline  bool islast  ( sdpoint_t* p ) { return (--(curve(p)->left.end())) == p->idl; };

  inline bool isfirst ( pdpoint_t* p ) { return curve(p)->left.begin() == p->ref || curve(p)->right.begin() == p->ref; };
  inline bool islast ( pdpoint_t* p ) { return (--(curve(p)->left.end())) == p->ref ||(--(curve(p)->right.end())) == p->ref;};

  inline bool isextrem( pdpoint_t* p  ) { return  isfirst(p) || islast(p); };
  inline bool isextrem( sdpoint_t*  p ) { return  isfirst(p) || islast(p); };

  void link( sdpoint_t* a, sdpoint_t* b);
  void link( pdpoint_t* a, pdpoint_t* b);
  inline void swap( dcurve * c ) { c->left.swap( c->right ); };

  inline bool isright( pdpoint_t* p )
    { return p->ref == curve(p)->right.begin() || p->ref == --(curve(p)->right.end());};


  template<class DPoint>
    inline void _link(  DPoint* a, DPoint* b ) 
    {
      if ( islast(a)   ) 
        {
          if( islast(b) )
            reverse(curve(b));
          // last a + first b
          append(curve(a),curve(b));
        }
      else
        {
          if ( isfirst(b) )
            reverse(curve(b));
          // last b + first a
          prepend(curve(a),curve(b));
        };
    };


  sdknode * make( sdpoint_t * first, sdpoint_t * last, const dim_cmp<sdpoint_t>& f ) ;
  pdknode * make( pdpoint_t ** first, pdpoint_t ** last, const dim_cmp<pdpoint_t>& f ) ;


  template<class DPoint>
    struct assoc_t
    {
      std::pair<DPoint*,DPoint*>         pp;
      double                              d;
      // comparing association on pointers
      struct pair_cmp
      { bool operator()( const assoc_t& a1, const assoc_t& a2 ) const { return a1.pp < a2.pp;  }; };
      // comparing association on distances
      struct dist_cmp
      { bool operator()( const assoc_t& a1, const assoc_t& a2 ) const { return a1.d > a2.d;    }; };
      assoc_t(){};
      assoc_t( DPoint * a, DPoint * b, double _d ):
        pp(a,b), d(_d){};
    };

  typedef std::priority_queue< assoc_t<pdpoint_t>, 
    std::vector< assoc_t<pdpoint_t> >, 
    assoc_t< pdpoint_t>::dist_cmp > pheap_t;    
  typedef std::priority_queue< assoc_t<sdpoint_t>, 
    std::vector< assoc_t<sdpoint_t> >, 
    assoc_t< sdpoint_t>::dist_cmp > sheap_t;

  typedef std::set< assoc_t< pdpoint_t > ,  assoc_t<pdpoint_t>::pair_cmp > pset_t;

  typedef assoc_t<sdpoint_t> sdassc_t;
  typedef assoc_t<pdpoint_t> pdassc_t;

  void solve_sheap( sheap_t& h );
  void solve_pheap( pheap_t& h );

  typedef std::pair<pdpoint_t*,pdpoint_t*> ppair_t;
  typedef std::list< ppair_t > links_t;
  typedef std::map< dcurve*, links_t > curves_links_t;
  typedef ordered_pair<dcurve*> cpair_t;

  template <class A, class B>
    inline void container_add( A& h, const B& a ) { h.push(a); };
  inline void container_add( pset_t& s, const pdassc_t& a ) { s.insert(a); };
  inline cpair_t curves( const pdassc_t& ass ) { return cpair_t(curve(ass.pp.first),curve(ass.pp.second)); };
  inline pdassc_t sibbles( const pdassc_t& ass )
    {
      pdpoint_t * a = sibble(ass.pp.first);
      pdpoint_t * b = sibble(ass.pp.second);
      if ( a < b )
        return pdassc_t(a,b,0);
      else 
        return pdassc_t(b,a,0);
    };

  void satisfy_links(  curves_links_t::iterator it );
  void solve_pset( curves_links_t& links, pset_t& s );
  std::list<dcurve*> * reduce2( std::vector< dcurve * > * conflicts );
  typedef  pset_t::iterator pset_iterator;

  struct dsearch : qsegment
  {
    bool m_cnf;
  
  
    //------------------------------------------------------------------------------
    std::vector< dcurve * >  * curr_config;           // represent the current region/region conflict resultings segments
    std::list  < dcurve * >  ** _branches;            // represent the reduced set of curves for each region/region conflict
  
    void cnfpush( qtree f, qtree s );
    void push( qtree f, qtree s );
    inline  std::list  < dcurve * >  *& branches( unsigned i, unsigned j ) 
    { return _branches[i*this->regions.size()+j]; };
    // add a pair of segment to the segment soup
    void push (  const point2& l0,  const point2& l1, 
                 const point2& r0,  const point2& r1,
                 const point3& i0,  const point3& i1 )
    { 
      curr_config->push_back( new dcurve( l0,l1,r0, r1, i0, i1 ) ); 
    };
  
    // use of the qsegment_t bounding box hierarchy to retrieve unit quad that intersects
#define push_s push
#define push_f push
    // search specialisation s is known to be a leaf
    void search_f(  qtree f, qtree const s );
    // same thing here f is a leaf
    void search_s(   qtree const f, qtree s );
    // search generic case: both f and s are root node
    void search(qtree f, qtree s );
    
    // launch the search for region i vs region j
    void search( rid_t i, rid_t j );
  
    // try to reduce the number of double curves
    std::list<dcurve*> * reduce( std::vector< dcurve * > * conflicts );
    inline void find_branches( rid_t id0, rid_t id1 )
    {
      branches( id0, id1 ) = reduce( curr_config );
      delete curr_config;
    };
  
    // stats structure
    struct 
    {
      double search_time;
      double recons_time;
    } stats;
  
    std::list< dcurve * >  *result;
  
    point2 ** lefts;
    point2 ** rights;
  
    unsigned * sizes;
    unsigned nbcurves;
    unsigned nbpoints;
    unsigned conflict_count;
  
    void GnuPlotOutput( char * name );
  
    std::vector< point3 > cnfdata3;
    std::vector< point2 > cnfdata2;
  
  
    double _st;
    dsearch( const ParametricSurface * s , unsigned w, unsigned h, bool cnf = false  );
    ~dsearch();
  };
}

} //namespace mmx
# undef ParametricSurface
# endif