// Cutmesh - Copyright (C) 2010-2018 T. Mouton, E. Bechet 
//
// See the LICENSE file for license information and contributions.
// bugs and problems to <thibaud.mouton@gmail.com>.

#ifndef _MESH_CUTMESH_H_
#define _MESH_CUTMESH_H_

#include "common.h"
#include "LTetra.h"
#include "LVertex.h"
#include "hashstruct.h"
#include "discreteSurface.h"

#if defined(HAVE_ANN)
#include "ANN/ANN.h"
#endif

class IBrep;

class Bucket;

class Shell
{
  private :
    LVertex * _v[2];
    std::vector<LTetra *> _t;
    std::vector<LVertex *> _o;
    std::vector<int> _e;
    int _lock;

  public :
    Shell()
    {
      _v[0] = _v[1] = NULL;
      _lock = 0;
    }
    void set_edge ( LVertex *v0, LVertex *v1 );
    void set_edge ( LTetra *t, int edge );
    void build();
    inline int size()
    {
      return _t.size();
    }
    void get_tetra ( int i, LTetra **t, int *edge );
    LVertex *get_vertex ( int i );
    inline int vertex_size()
    {
      return _o.size();
    }
    inline void lock()
    {
      _lock = 1;
    }
    inline void unlock()
    {
      _lock = 0;
    }
    inline int locked()
    {
      return _lock;
    }
    void getNext ( LTetra *t, int opp, LVertex *common, LTetra **tn, int *oppn );
    void dump();
    void dump_vtk();
};

class Ball
{
  private :
    LVertex * _v;
    std::vector<LTetra *> _t;
    std::vector<int> _o;
    int _lock;

  public :
    Ball()
    {
      _v = NULL;
      _lock = 0;
    }
    void set_vertex ( LVertex *v );
    void build();
    inline int size()
    {
      return _t.size();
    }
    void get_tetra ( int i, LTetra **t, int *vertex );
    inline LVertex *get_vertex ( int i );

    inline void lock()
    {
      _lock = 1;
    }
    inline void unlock()
    {
      _lock = 0;
    }
    inline int locked()
    {
      return _lock;
    }
    void getNext ( LTetra *t, int opp, LVertex *common, LTetra **tn, int *oppn );
    void dump();
    void dump_vtk();
};

class Merge
{
  private :
    LTetra *_t;
    int _o;
    int _s;
    int _imposed;
    LVertex *_vlist[6];
    LVertex *_n;
    LVertex *_c;
    Merge *_ref;
    int flag;

  public :
    Merge ( LTetra *t/*, int opp*/, int surf );

    inline void set_ref ( Merge *m )
    {
      _ref = m;
    }

    int compare ( Merge *m, double dotmin, double distmax );
    void merge_with ( Merge *m );
    inline LTetra *get_tetra()
    {
      return _t;
    }
    inline int get_opp()
    {
      return _o;
    }
    inline int get_surf()
    {
      return _s;
    }
    inline LVertex *get_normal()
    {
      return _n;
    }
    inline LVertex *get_centroid()
    {
      return _c;
    }
    inline Merge *get_ref()
    {
      return _ref;
    }
    inline void set_flag ( int i )
    {
      flag = i;
    }
    inline int get_flag()
    {
      return flag;
    }
    void set ( LVertex *v0, LVertex *v1, LVertex *v2 );
    void set ( LVertex *v0, LVertex *v1, LVertex *v2, LVertex *v3 );
    void build();
    int apply_merge();
};

class EntityManager
{
  private:
    /////////////////////// volumes ///////////////////////
    std::vector<int> _id_vol;
    std::multimap<int, int> _associated_id_vol;
    int _max_vol_id;
    std::map<int, int> _lookup_id_vol;

    std::set<int> _realv; // volume reel
    std::set<int> _restv; // volume restore

    /////////////////////// faces ///////////////////////
    std::vector<int> _id_face;
    std::multimap<int, int> _associated_id_face;
    int _max_face_id;
    std::map<int, int> _lookup_id_face;

    std::set<int> _realf; // surface reelle
    std::set<int> _virtf; // surface virtuelle
    std::set<int> _restf; // surface restoree

    std::map<int, std::set<int> > _map_tang;
    std::multimap<std::pair<int, int>, int> _inv_map_tang;

    /////////////////////// edges ///////////////////////
    std::vector<int> _id_edge;
    std::multimap<int, int> _associated_id_edge;
    int _max_edge_id;
    std::map<int, int> _lookup_id_edge;

    std::set<int> _reale; // surface reelle
    std::set<int> _reste; // surface restoree

  public:
    EntityManager()
    {
      _max_edge_id = _max_face_id = _max_vol_id = 0;
    }
    inline int add_face_id ( int f )
    {
//         printf("%p : add surf %d\n", this, f);
      std::map <int, int >::iterator it = _lookup_id_face.find ( f );
      if ( it != _lookup_id_face.end() )
      {
        dbgprintf ( 5, "/!\\ --> surf %d is already used !\n", f );
        return it->second;
      }
      int index = _id_face.size();
      _id_face.push_back ( f );
      if ( _max_face_id < f )
        _max_face_id = f;
      _lookup_id_face.insert ( std::pair<int, int> ( f, index ) );
      return index;
    }
    inline int add_face_id()
    {
//         printf("%p : add new surf %d\n", this, _max_face_id+1);
      int index = _id_face.size();
      _max_face_id++;
      _id_face.push_back ( _max_face_id );
      _lookup_id_face.insert ( std::pair<int, int> ( _max_face_id, index ) );
      return index;
    }
    inline int add_cut_face_id ( int f, int s1, int s2 )
    {
//      _id_face.push_back(f);
//      _lookup_id_face[f] = _id_face.size()-1;
      int lid = add_face_id ( f );
      int ls1 = get_lookup_face_id ( s1 );
      int ls2 = get_lookup_face_id ( s2 );
      _map_tang[lid].insert ( ls1 );
      _map_tang[lid].insert ( ls2 );
      assert ( ls1 <= ls2 );
      _inv_map_tang.insert ( std::pair<std::pair<int, int>, int> ( std::pair<int, int> ( ls1, ls2 ), lid ) );
      return lid;
    }
    inline int get_face_id ( int i )
    {
      return _id_face[i];
    }
    inline int get_lookup_face_id ( int f )
    {
//         printf("%p : lookup id %d\n", this, f);
//  for (std::map<int, int>::iterator it= _lookup_id_face.begin(); it != _lookup_id_face.end(); it++)
//    printf("%d <--> %d\n", it->first, it->second);
      std::map<int, int>::iterator it = _lookup_id_face.find ( f );
      if ( it == _lookup_id_face.end() )
        return -1;
      assert ( it != _lookup_id_face.end() );
      return it->second;
    }
    inline void get_cutted_face_id ( int i, int *s1, int *s2 )
    {
      std::map<int, std::set<int> >::iterator it = _map_tang.find ( i );
      if ( it != _map_tang.end() )
      {
        *s1 = *it->second.begin();
        *s2 = *++it->second.begin();
      }
      else
        *s1 = *s2 = -1;
    }
    inline void get_virtual_face_id ( int s1, int s2, std::vector<int> &ts )
    {
      assert ( s1 < s2 );
      std::multimap<std::pair<int, int>, int>::iterator it;
      std::pair<std::multimap<std::pair<int, int>, int>::iterator, std::multimap<std::pair<int, int>, int>::iterator> ret;
      ret = _inv_map_tang.equal_range ( std::pair<int, int> ( s1, s2 ) );
      for ( it = ret.first; it != ret.second; it++ )
      {
        ts.push_back ( it->second );
      }
    }
    inline int get_max_face_id()
    {
      return _max_face_id;
    }
    inline int face_size()
    {
      return _id_face.size();
    }
    inline void set_real_face ( int surf )
    {
      _realf.insert ( surf );
    }
    inline void set_virtual_face ( int surf )
    {
      _virtf.insert ( surf );
    }
    inline void set_restore_face ( int surf )
    {
//         printf("set %d as restore face\n", surf);
      _restf.insert ( surf );
    }

    inline int is_real_face ( int surf )
    {
      return _realf.find ( surf ) != _realf.end();
    }
    inline int is_virtual_face ( int surf )
    {
      return _virtf.find ( surf ) != _virtf.end();
    }
    inline int is_restore_face ( int surf )
    {
//         printf("surf %d is restore face ? : %d\n", surf, _restf.find(surf) != _restf.end());
      return _restf.find ( surf ) != _restf.end();
    }

    inline void associate_face ( int surf, int asurf )
    {
      dbgprintf ( 2, "associate face %d to face %d\n", asurf, surf );
      _associated_id_face.insert ( std::pair<int, int> ( surf, asurf ) );
    }

    inline std::vector<int> get_associated_faces ( int surf )
    {
      std::vector<int> list;
      std::pair<std::multimap<int, int>::iterator, std::multimap<int, int>::iterator >  ret;
      ret = _associated_id_face.equal_range ( surf );
      std::multimap<int, int>::iterator it;
      for ( it = ret.first; it != ret.second; it++ )
        list.push_back ( it->second );
      return list;
    }

    inline int associated_faces_size ( int surf )
    {
      return _associated_id_face.count ( surf );
    }


    inline int add_vol_id ( int v )
    {
      std::map <int, int >::iterator it = _lookup_id_vol.find ( v );
      if ( it != _lookup_id_vol.end() )
      {
//        tprintf("/!\\ --> vol %d is already used !\n", v);
        return it->second;
      }
      int index = _id_vol.size();
      _id_vol.push_back ( v );
      if ( _max_vol_id < v )
        _max_vol_id = v;
      _lookup_id_vol.insert ( std::pair<int, int> ( v, _id_vol.size()-1 ) );
      return _id_vol.size()-1;
    }
    inline int add_vol_id()
    {
      int index = _id_vol.size();
      _max_vol_id++;
      _id_vol.push_back ( _max_vol_id );
      _lookup_id_vol.insert ( std::pair<int, int> ( _max_vol_id, index ) );
      return index;
    }
    inline int get_vol_id ( int i )
    {
      assert ( i < _id_vol.size() );
      return _id_vol[i];
    }
    inline int get_max_vol_id()
    {
      return _max_vol_id;
    }
    inline int vol_size()
    {
      return _id_vol.size();
    }
    inline void set_real_vol ( int vol )
    {
      _realv.insert ( vol );
    }
    inline void set_restore_vol ( int vol )
    {
//         printf("set %d as restore face\n", surf);
      _restv.insert ( vol );
    }
    inline int is_real_vol ( int vol )
    {
      return _realv.find ( vol ) != _realv.end();
    }
    inline int is_restore_vol ( int vol )
    {
//         printf("surf %d is restore face ? : %d\n", surf, _restf.find(surf) != _restf.end());
      return _restv.find ( vol ) != _restv.end();
    }
    inline void associate_vol ( int vol, int avol )
    {
      dbgprintf ( 2, "associate vol %d to vol %d\n", avol, vol );
      _associated_id_vol.insert ( std::pair<int, int> ( vol, avol ) );
    }

    inline std::vector<int> get_associated_vols ( int vol )
    {
      std::vector<int> list;
      std::pair<std::multimap<int, int>::iterator, std::multimap<int, int>::iterator >  ret;
      ret = _associated_id_vol.equal_range ( vol );
      std::multimap<int, int>::iterator it;
      for ( it = ret.first; it != ret.second; it++ )
        list.push_back ( it->second );
      return list;
    }

    inline int associated_vols_size ( int vol )
    {
      return _associated_id_vol.count ( vol );
    }

    inline int add_edge_id ( int e )
    {
      std::map <int, int >::iterator it = _lookup_id_edge.find ( e );
      if ( it != _lookup_id_edge.end() )
      {
        dbgprintf ( 5, "/!\\ --> edge %d is already used !\n", e );
        return it->second;
      }
      int index = _id_edge.size();
      _id_edge.push_back ( e );
      if ( _max_edge_id < e )
        _max_edge_id = e;
      _lookup_id_edge.insert ( std::pair<int, int> ( e, index ) );
      return index;
    }
    inline int add_edge_id()
    {
      int index = _id_edge.size();
      _max_edge_id++;
      _id_edge.push_back ( _max_edge_id );
      _lookup_id_edge.insert ( std::pair<int, int> ( _max_edge_id, index ) );
      return index;
    }
    inline int get_edge_id ( int i )
    {
      return _id_edge[i];
    }
    inline int get_lookup_edge_id ( int f )
    {
      std::map<int, int>::iterator it = _lookup_id_edge.find ( f );
      if ( it == _lookup_id_edge.end() )
        return -1;
      assert ( it != _lookup_id_edge.end() );
      return it->second;
    }
    inline int get_max_edge_id()
    {
      return _max_edge_id;
    }
    inline int edge_size()
    {
      return _id_edge.size();
    }
    inline void set_real_edge ( int edge )
    {
      _reale.insert ( edge );
    }
    inline void set_restore_edge ( int edge )
    {
//         printf("set %d as restore face\n", surf);
      _reste.insert ( edge );
    }

    inline int is_real_edge ( int edge )
    {
      return _reale.find ( edge ) != _reale.end();
    }
    inline int is_restore_edge ( int edge )
    {
//         printf("surf %d is restore face ? : %d\n", surf, _restf.find(surf) != _restf.end());
      return _reste.find ( edge ) != _reste.end();
    }
    inline void associate_edge ( int edge, int aedge )
    {
      dbgprintf ( 2, "associate edge %d to edge %d\n", aedge, edge );
      _associated_id_edge.insert ( std::pair<int, int> ( edge, aedge ) );
    }

    inline std::vector<int> get_associated_edges ( int edge )
    {
      std::vector<int> list;
      std::pair<std::multimap<int, int>::iterator, std::multimap<int, int>::iterator >  ret;
      ret = _associated_id_edge.equal_range ( edge );
      std::multimap<int, int>::iterator it;
      for ( it = ret.first; it != ret.second; it++ )
        list.push_back ( it->second );
      return list;
    }

    inline int associated_edges_size ( int edge )
    {
      return _associated_id_edge.count ( edge );
    }
};

class Mesh
{
  private :
    int global_tetra_index;
    int global_vertex_index;

    hashFaceStruct *_hfs;
    hashEdgeStruct *_hes;
    EntityManager *_em;
    Shell *_sot;
    Ball *_bot;
    Bucket *_bs;

    // variables globales tetra, face et edge
    std::queue<LTetra *> _gtq;
    std::queue<int> _gfq;
    std::queue<int> _geq;

    std::vector<LTetra *> _t;
    std::vector<LVertex *> _v;
    std::set<int> free_t; // index des emplacements libres de tt
    std::set<int> free_v; // index des emplacements libres de tt
    std::vector<int> _cs; // indice des surface commune a plusieurs sommets

    std::queue<LTetra *> _sc; // pile de coloration

    options *_opts;

    double _le;

    std::map<int, int> _lookup_vertex_id;
    std::map<int, int> _lookup_tetra_id;

    std::map<int, int>::iterator _lvit;
    std::map<int, int>::iterator _ltit;

#if defined(HAVE_ANN)
    ANNkd_tree *_kdtree;
#endif

    // refinement
    std::queue<LTetra *> _new_tetra_stack;
    std::queue<LTetra *> _refining_stack;
    std::queue<LTetra *> _conforming_stack;
    std::queue<LTetra *> _refined_stack;
    std::queue<LTetra *> _father_stack;
    std::queue<LTetra *> _residual_stack;
    int _max_generation;

    int _start_index;

  public :
    Mesh ( options *opt );
    ~Mesh();
//     inline size_t memory()
//     {
//       return sizeof(*this) + hfs->memory() + hes->memory() + _t.size()*sizeof(LTetra) + _t.size()*sizeof(LTetra*) + _v.size()*sizeof(LVertex) + _v.size()*sizeof(LVertex*);
//     }

    void allocate_vertices ( int size )
    {
      _v.resize ( size );
    }

    void allocate_tetra ( int size )
    {
      _t.resize ( size );
    }

    void set_start_index ( int pos, int index )
    {
      _start_index = index - pos;
    }

    LVertex * new_vertex ( double x, double y, double z )
    {
      return new LVertex ( global_vertex_index++, x, y, z );
    }

    LTetra * new_tetra ( LVertex *v0, LVertex *v1, LVertex *v2, LVertex *v3, LTetra *parent = NULL )
    {
      LTetra *t = new LTetra ( v0, v1, v2, v3, parent );
      t->setId ( global_tetra_index++ );
      return t;
    }

    inline LVertex * get_vertex ( int i )
    {
      assert ( i<_v.size() );
      return _v[i];
    }
    inline LTetra * get_tetra ( int i )
    {
      assert ( i<_t.size() );
      return _t[i];
    }
    inline LVertex * get_vertex_from_lookup ( int id )
    {
      _lvit = _lookup_vertex_id.find ( id );
      if ( _lvit == _lookup_vertex_id.end() )
        return NULL;
      return get_vertex ( _lvit->second );
    }
    inline LTetra * get_tetra_from_lookup ( int id )
    {
      _ltit = _lookup_tetra_id.find ( id );
      if ( _ltit == _lookup_tetra_id.end() )
        return NULL;
      return get_tetra ( _ltit->second );
    }
    inline int vertex_size()
    {
      return _v.size();
    }
    inline int tetra_size()
    {
      return _t.size();
    }
    inline void add_tetra_to_queue ( LTetra *t )
    {
      _gtq.push ( t );
    }
    inline void add_face_to_queue ( int face )
    {
      _gfq.push ( face );
    }
    inline void add_edge_to_queue ( int edge )
    {
      _geq.push ( edge );
    }
    inline LTetra *tetra_from_queue()
    {
      return _gtq.front();
    }
    inline int face_from_queue()
    {
      return _gfq.front();
    }
    inline int edge_from_queue()
    {
      return _geq.front();
    }
    inline void tetra_dequeue()
    {
      _gtq.pop();
    }
    inline void face_dequeue()
    {
      _gfq.pop();
    }
    inline void edge_dequeue()
    {
      _gtq.pop();
      _gfq.pop();
      _geq.pop();
    }
    inline int common_surf_size()
    {
      return _cs.size();
    }
    inline int get_common_surf ( int i )
    {
      assert ( i<_cs.size() );
      return _cs[i];
    }


    void compute_common_surf ( LVertex *v0, LVertex *v1 );
    void compute_common_zero_surf ( LVertex *v0, LVertex *v1 );
    void compute_common_surf ( std::vector<LVertex*> &vx );
    void compute_common_zero_surf ( std::vector<LVertex*> &vx );
    void compute_common_surf ( LTetra *t );
    void compute_common_zero_surf ( LTetra *t );
    int same_surf_zero ( LVertex *v0, LVertex *v1 );
    int same_root_parents ( LVertex *v0, LVertex *v1 );
    int same_surf_zero ( std::vector<LVertex*> &vx );
    int same_root_parents ( std::vector<LVertex*> &vx );
    int get_generation ( LTetra *t );
    LVertex * intersection ( LTetra *t, int edge, int s );
    int cross ( LVertex *v0, LVertex *v1, int s );
    int cross ( LTetra *t, int edge, int s );
    int crossing_surf ( LVertex *v0, LVertex *v1 );
    int crossing_surf ( LTetra *t, int edge );
    void process_edge();
    void split_edge ( LVertex *newv, int update );
    void split_face ( LVertex *newv, int update );
    void split_tetra ( LVertex *newv, int update );
    LTetra *locate_vertex ( LVertex *v, double b[4] );
    LVertex *insert_topo_vertex ( LVertex *v, LTetra *t, double b[4], std::vector<int> surf );
    LTetra *insert_topo_edge ( LVertex *v0, LVertex *v1, LTetra *t, std::vector<int> surf );
    void set_tetra ( int i, LTetra *t );
    void add_tetra ( LTetra *t );
    void add_tetra ( LTetra *t, int id );
    void set_vertex ( int i, LVertex *v );
    void add_vertex ( LVertex *v );
    void add_vertex ( LVertex *v, int id );
    void remove_tetra ( LTetra *t );
    void remove_vertex ( LVertex *v );
    void shrink_tetra();
    void shrink_vertex();
    void bucketize();
    void compute_adjacencies();
    void check_adjacence();
    void check_indexing();
    void check_surfaces ( int s );
    void readMeshFile ( const char *name );
    void ImportIbrep ( IBrep *ib ); // import ibrep = readmesh + readlevelset
    void readLevelsetType0 ( FILE *fp );
    void readLevelsetType1 ( FILE *fp );
    void addLevelset ( int ls, LVertex *v, double dist );
    void addLevelset ( int ls, LTetra *t, double dist[4] );
    void collect_edges();
    void readEdges ( FILE *fp );
    void readLevelsetFile ( const char *name );
    void load_existing_entities ( const char *filename );
    void build_iso_zero_tri ( LTetra *t, std::vector<Merge*> &mergelist );
//     int merged_surfaces_enforced(LTetra *t);
    int merged_surfaces ( LTetra *t, double dotmin, double distmax );
    void surface_merging ( double dotmin, double distmax );
//     void surface_merging_enforced();
    void writeMeshVTK ( const char *name );
    void writeSupportVTK ( int local_id_surf );
    void writeMesh ( /*const char *name*/ );
    void debugMesh ( LTetra *t );
    void debugMesh ( std::vector<LTetra*> tl );
    void writeMeshByVolumeId ( const char *name );
    void writeShell ( std::vector<LTetra *> shell );
    void cut_mesh();
    void buildRootParents();
    int delete_degenerated_tetra ( LTetra *t, double threshold, std::vector<LTetra*> &modified );
    int delete_extra_tetra ( LTetra *t, double threshold, std::vector<LTetra*> &modified );
    void clean_invalid_tetras ( double threshold );
    void remove_extra_tetras ( double threshold );
    void optimize();
    int merge_vertices ( std::vector<LVertex*> list, LVertex *mv, int move_allowed, double threshold, std::vector<LTetra*> &modified ); // merge list to mv

    inline void compute_longuest_edge()
    {
      LVertex d;
      for ( int i = 0; i < tetra_size(); i++ )
      {
        LTetra *t = get_tetra ( i );
        for ( int j = 0; j < 6; j++ )
        {
          LVertex *a = t->getEdgeVertex ( j, 0 );
          LVertex *b = t->getEdgeVertex ( j, 1 );
          d.setVector ( a, b );
          double l = d.norm();
          if ( l > _le )
            _le = l;
        }
      }
    }

    inline double longuest_edge()
    {
      return _le;
    }

    inline void set_tetra_size ( int i )
    {
      _t.resize ( i, NULL );
    }

    inline void seed ( LTetra *t )
    {
      _sc.push ( t );
    }
    int color ( int vol, int unfilled_vol );
    inline void set_all_vertex_flag ( int vf )
    {
      for ( int i = 0; i < vertex_size(); i++ )
        if ( get_vertex ( i ) != NULL )
          get_vertex ( i )->setFlag ( vf );
    }
    inline void set_all_tetra_flag ( int tf )
    {
      for ( int i = 0; i < tetra_size(); i++ )
        if ( get_tetra ( i ) != NULL )
          get_tetra ( i )->setFlag ( tf );
    }

    inline void clear_all_tetra_child()
    {
      for ( int i = 0; i < tetra_size(); i++ )
      {
        LTetra *at = get_tetra ( i )->getAncestor();
        if ( at )
          at->clear_child();
        assert ( get_tetra ( i ) != NULL );
        get_tetra ( i )->clear_child();
      }
    }

    inline EntityManager* get_entity_manager()
    {
      return _em;
    }

    // refinement
    int normalize_curv ( double rcurv, double seglenght, double factor );
    double evaluate_refinement_params();
    LTetra* find_tetra_containing ( LVertex *pt );
    void add_included ( DiscreteSurface *surf );
    void refinement();
    void refine_triangulation();
    void compute_refinement_adjacencies();
    int adapt ( double max_curv_allowed );
    void renumber();

    int valid_scheme ( LTetra *t );
    int split ( LTetra *t );
    void collect_neighbors ( LTetra *t );
    void subdivide ( LTetra *t );
    void init_stack();
    void refining();
    void check_adjacencies();
    void equilibrate_refinement();
};


class Bucket
{
  private:
    LVertex _min;
    LVertex _max;
    int _nx;
    int _ny;
    int _nz;
    int _maxn;;
    double _ref_size;
    std::vector<std::vector<LVertex *> > _bucket;
    std::vector<int> _flag;
    int _request_id;
  public:
    Bucket ( LVertex min, LVertex max, double enlarge, int max_size );
//     inline size_t memory()
//     {
//         size_t total = 0;
//         for (int i = 0; i < _bucket.size(); i++)
//             total += _bucket[i].size()*sizeof(LVertex*);
//         return sizeof(*this) + _flag.size()*sizeof(int) + total;
//     }
    void get_coord ( LVertex *pt, int &x, int &y, int &z );
    int get_index ( int x, int y, int z );
    void get_indices ( int index, int &x, int &y, int &z );
    LVertex *get_closest_in_bucket ( int index, LVertex *pt, double *sq_dist );
    void add_point ( LVertex *pt );
    LVertex *get_closest ( LVertex *pt );
    int contain ( LVertex *pt );
    void dump();
};

#endif
