#ifndef SIMPLEX_H
#define SIMPLEX_H

#include <stdlib.h>
#include <stdio.h>
#include <vector>
#include <math.h>
#include <assert.h>

#define NOT_REFINED 0
#define TO_REFINE 1
#define SPLITTED 2
#define REFINED 3

#define FULL_SPLIT 63

// #define REFINED_1_T0_2 1
// #define REFINED_1_T0_4 2
// #define REFINED_1_T0_8 3
// #define TO_REFINE_1_T0_2 4
// #define TO_REFINE_1_T0_4 5
// #define TO_REFINE_1_T0_8 6
// #define EDGES_SPLITTED 7
#define TO_ERASE 8
#define ERASED 100
#define UNDEFINED 1000


class Tetra;
class Vertex;

class Tetra
{
private :
    Vertex *_v[4];
    Tetra *_a[4];
    int _o[4];
    int id;
    int index;
    int durty_bit;
    
    int flag;
    char bit_edge;
    int _state;
    int _refdepth;
    Tetra *_father;
    Tetra *_child[8];
    int _generation;  
    std::vector<Vertex *> _inc;
    double _rl;

public :
    Tetra(Vertex *v0, Vertex *v1, Vertex *v2, Vertex *v3);
    Tetra(int i, Vertex *v0, Vertex *v1, Vertex *v2, Vertex *v3);
    inline void setAdj(int i, Tetra *t, int opp)
    {
        _a[i] = t;
        _o[i] = opp;
    }
    void setMutualAdj(int opp1, Tetra *t2, int opp2);
    inline Tetra *getAdj(int i)
    {
        assert(i>=0 && i<4);
        return _a[i];
    }
    int getOpp(int i)
    {
        assert(i>=0 && i<4);
        return _o[i];
    }
    inline void getAdjOpp(int i, Tetra **t, int *opp)
    {
        assert(i>=0 && i<4);
        *t = _a[i];
        *opp = _o[i];
    }
    inline Vertex *getVertex(int i)
    {
        assert(i>=0 && i<4);
        return _v[i];
    }
    inline Vertex *setVertex(int i, Vertex *v)
    {
        assert(i>=0 && i<4);
        assert(v != NULL);
        _v[i] = v;
    }
    int getEdgeIndex(Vertex *v0, Vertex *v1);
    int getFaceEdgeIndex(int face, int edge);
    int getEdgeVertexIndex(int edge, int vertex);
    int getEdgeCode(int edge);    
    Vertex *getEdgeVertex(int edge, int vertex);
    int getOppEdgeVertexIndex(int edge, int vertex);
    Vertex *getOppEdgeVertex(int edge, int vertex);
    int getFaceVertexIndex(int face, int vertex);
    Vertex *getFaceVertex(int face, int vertex);
    int getFaceCode(int face); 
    int getEdgeIndex(int face, int edge);
    int getFaceEdgeVertexIndex(int face, int edge, int vertex);
    Vertex *getFaceEdgeVertex(int face, int edge, int vertex);
    int getOppEdgeOppVertexIndex(int edge, int vertex);
    Vertex *getOppEdgeOppVertex(int edge, int vertex);
    inline void setFlag ( int f )
    {
        flag = f;
    }
    inline void incFlag ()
    {
        flag++;
    }
    inline void decFlag ()
    {
        flag--;
    }
    inline int getFlag()
    {
        return flag;
    }
    inline void setIndex ( int i )
    {
        index = i;
    }
    inline int getIndex()
    {
        return index;
    }
    inline int getId()
    {
        return id;
    }
    
    void split_edge(int edge);
    int edge_state(int edge);
    int edge_split_count();
    
    
    inline int edge_code()
    {
      return (int)bit_edge;
    }
    
    inline void setState(int state)
    {
        _state = state;
    }
//     inline void setCurrentState(int state)
//     {
//         _current_state = state;
//     }
    inline int getState()
    {
        return _state;
    }
//     inline int getCurrentState()
//     {
//         return _current_state;
//     }
    inline void setRefDepth(int i)
    {
      _refdepth = i;
    }
    inline int getRefDepth()
    {
        return _refdepth;
    }
    inline void setFather(Tetra* father)
    {
        _father = father;
	setGen(father->getGen()+1);
	setRefDepth(father->getRefDepth()-1);
// 	if (getRefDepth()<0)
// 	{
// 	  setRefDepth(0);
// 	  Tetra *p = getFather();
// 	  do
// 	  {
// 	    p->setRefDepth(p->getChild(0)->getRefDepth()+1);
// 	    p = p->getFather();
// 	  }
// 	  while (p != NULL);
// 	}
// 	assert(getRefDepth()>=0);
    }
    inline Tetra* getFather()
    {
        return _father;
    }   
    inline void setChild(int i, Tetra *c)
    {
	assert(i>=0 && i<8);
        _child[i] = c;
    }
    inline Tetra* getChild(int i)
    {
	assert(i>=0 && i<8);
        return _child[i];
    }
    inline void setGen(int g)
    {
        _generation = g;
    }
    inline int getGen()
    {
        return _generation;
    }
    
    inline void addInc(Vertex *iv)
    {
        _inc.push_back(iv);
    }

    inline Vertex *getInc(int i)
    {
        return _inc[i];
    }

    inline int inc_size()
    {
        return _inc.size();
    }
    
    inline double getlenght()
    {
        return _rl;
    }
    
    void reset();    
    
    int barycentric(int face, Vertex *c, double b[4]);
    int barycentric(Vertex *c, double b[4]);
//    int visibility(int face, Vertex *v0, Vertex *v1, double b[4], Vertex **inter, int *inside, int *vertex, int *edge);
    int contain(Vertex *v);
};

class Vertex
{
private :
    double _p[3];
    Tetra *_t;
    int id;
    int index;
    int surf;
    Vertex *_n;
    double _c;
    double dist;
    double projdist;
//     std::vector<double> _value;
//     std::vector<int> _surf;
//     std::vector<int> _zero;
    int flag;

public :
    Vertex()
    {
        _p[0] = _p[1] = _p[2] = 0.0;
        _t = NULL;
        _n = NULL;
	_c = 0.0;
        id = 0;
        index = 0;
        flag = 0;
    }
    Vertex(double x, double y, double z)
    {
        _p[0] = x;
        _p[1] = y;
        _p[2] = z;
        id = 0;
        _t = NULL;
        _n = NULL;
	_c = 0.0;
        index = 0;
        flag = 0;
    }
    Vertex(int i, double x, double y, double z)
    {
        _p[0] = x;
        _p[1] = y;
        _p[2] = z;
        id = i;
        _t = NULL;
        _n = NULL;
	_c = 0.0;
        index = 0;
        flag = 0;
    }

    inline void setPosition(double x, double y, double z)
    {
      _p[0]=x; _p[1]=y; _p[2]=z;
    }
    inline void setVector(Vertex *v0, Vertex *v1)
    {
      _p[0]=v1->x()-v0->x(); _p[1]=v1->y()-v0->y(); _p[2]=v1->z()-v0->z();
//       printf("vector --> %lf %lf %lf\n", _p[0], _p[1], _p[2]);
    }
    inline void setLinComb(Vertex *v0, Vertex *v1, double t)
    {
      _p[0]=v0->x()+t*(v1->x()-v0->x()); _p[1]=v0->y()+t*(v1->y()-v0->y()); _p[2]=v0->z()+t*(v1->z()-v0->z());
//       printf("vector --> %lf %lf %lf\n", _p[0], _p[1], _p[2]);
    }
    inline double x()
    {
        return _p[0];
    }
    inline double y()
    {
        return _p[1];
    }
    inline double z()
    {
        return _p[2];
    }
    inline void setTetra(Tetra *t)
    {
        _t = t;
    }
    inline Tetra *getTetra()
    {
        return _t;
    }
    inline void setNormal(Vertex *n)
    {
        _n = n;
    }
    inline void setNormal(double x, double y, double z)
    {
	_n = new Vertex(x, y, z);
    }
    inline Vertex *getNormal()
    {
//       printf("n --> %lf %lf %lf\n", _n->x(), _n->y(), _n->z());
        return _n;
    }
    inline void setCurvature(double c)
    {
	_c = c;
    }
    inline double getCurvature()
    {
        return _c;
    }
    inline void setFlag(int f)
    {
        flag = f;
    }
    inline void incFlag ()
    {
        flag++;
    }
    inline void decFlag ()
    {
        flag--;
    }
    inline int getFlag()
    {
        return flag;
    }
    inline void setIndex(int i)
    {
        index = i;
    }
    inline int getIndex()
    {
        return index;
    }
    inline void setSurf(int i)
    {
        surf = i;
    }
    inline int getSurf()
    {
        return surf;
    }
    inline int getId()
    {
        return id;
    }
    inline void crossproduct(Vertex *v1, Vertex *r)
    {
        r->setPosition(_p[1]*v1->z()-_p[2]*v1->y(), _p[2]*v1->x()-_p[0]*v1->z(), _p[0]*v1->y()-_p[1]*v1->x());
// 	r[0] = _p[1]*v1->z()-_p[2]*v1->y();
//         r[1] = _p[2]*v1->x()-_p[0]*v1->z();
//         r[2] = _p[0]*v1->y()-_p[1]*v1->x();
    }
    inline double dotproduct(Vertex *v1)
    {
        return _p[0]*v1->x() + _p[1]*v1->y() + _p[2]*v1->z();
    }
    inline double norm()
    {
      return sqrt(_p[0]*_p[0] + _p[1]*_p[1] + _p[2]*_p[2]);
    }
    inline double sq_norm()
    {
      return _p[0]*_p[0] + _p[1]*_p[1] + _p[2]*_p[2];
    }
    inline void normalize()
    {
      double l = norm();
      _p[0] /= l;
      _p[1] /= l;
      _p[2] /= l;
    }
    inline double get_projdist()
    {
      return projdist;      
    }
    inline double set_projdist(double d)
    {
      projdist = d;      
    }
    inline double get_dist()
    {
      return dist;      
    }
    inline double set_dist(double d)
    {
      dist = d;      
    }
    
};

#endif