#include "DILevelsetTextile.h"
#include <assert.h>

#define CalculDistance 0

#if (CalculDistance==0)
double gLevelsetFils::operator()(const double x, const double y, const double z) const
{
  TVertex * initial;
  initial = new TVertex(x,y,z);

  TVertex * searched;
  searched = _lsTisse->getTisse()->getBucket(_numBucket)->get_closest(initial);

  TVertex vector1;
  vector1.setVector(initial, searched);

  Section * section =searched->getSection();
  assert(section!=NULL);
  TVertex * centre = section->getCentre();

  TVertex vector2;
  vector2.setVector(initial, centre);
  double signe = vector2.dotproduct(&vector1);
  if(signe<0)
    return -vector1.norm();
  else
    return vector1.norm();
}
#else
double gLevelsetFils::operator()(const double &x, const double &y, const double &z) const
{
  TVertex initial(x,y,z);

  TVertex * searched;
  searched = _lsTisse->getTisse()->getBucket(_numBucket)->get_closest(&initial);

  Section * section =searched->getSection();
  assert(section!=NULL);
  TVertex * centre = section->getCentre();
  int index=centre->getIndex();

  Section * sectionBefor =section->getSectionBefor();
  Section * sectionAfter =section->getSectionAfter();

  assert(section->getNbNodes()>2);

  // Interpolation PID
  std::vector<TVertex*> v;

  // Genere la liste des points interpolés connus
  if (sectionBefor!=NULL)
  {
    v.push_back( sectionBefor->getNode(index-1) );
    v.push_back( sectionBefor->getNode(index) );
    v.push_back( sectionBefor->getNode(index+1) );
  }

  v.push_back( section->getNode(index-1) );
  v.push_back( section->getNode(index) );
  v.push_back( section->getNode(index+1) );

  if (sectionAfter!=NULL)
  {
    v.push_back( sectionAfter->getNode(index-1) );
    v.push_back( sectionAfter->getNode(index) );
    v.push_back( sectionAfter->getNode(index+1) );
  }

  // calcul de la valeur interpolée
  double p=1, denom=0, w=0;
  TVertex numer, temp;
  for (int i=0; i<v.size(); i++)
  {
    temp.setVector(v[i], &initial);
    w=1/pow(temp.norm(),p);
    numer+=w*(*v[i]);
    denom+=w;
  }
  denom=1/denom;
  TVertex interp(denom*numer);

  TVertex vector1;
  vector1.setVector(&interp, &initial);

  TVertex vector2;
  vector2.setVector(centre, searched);
  double signe = vector2.dotproduct(&vector1);
  if(signe<0)
    return -vector1.norm();
  else
    return vector1.norm();
}
#endif


gLevelsetTextile::gLevelsetTextile (const double &epaisseur, const double &dimChaine, const double &dimTrame,
                                  const std::vector<double> &axesChaine, const std::vector<double> &axesTrame,
                                  const int &layer, const fullMatrix<int> &structure, const std::vector<int> &multipliciteChaine,
                                  const int &nbdx, const int &nbdyz,
                                  int &tag) : gLevelsetPrimitive(tag), _tisse(epaisseur, dimChaine, dimTrame, axesChaine, axesTrame)
{
  _tisse.analyseStructure(layer, structure, multipliciteChaine);
  _tisse.genereTisse(nbdx, nbdyz);

  _lsFils.resize(_tisse.getNbBuckets());
  for (int i=0; i<_lsFils.size(); i++)
    _lsFils[i] = new gLevelsetFils(this, i, tag+i);

  //_tisse.printTisseGeo();
  //_tisse.printTisseGeoSurf();
  //_tisse.printLsNodes();
  _tisse.printLsNodesByFils();
  
}

/*double gLevelsetTextile::operator()(const double &x, const double &y, const double &z) const
{
  TVertex * initial;
  initial = new TVertex(x,y,z);
  
  TVertex * searched;
  searched = _tisse.getBucket(0)->get_closest(initial);
  
  TVertex vector1;
  vector1.setVector(initial, searched);
  
  Section * section =searched->getSection();
  assert(section!=NULL);
  TVertex * centre = section->getCentre();
  
  TVertex vector2;
  vector2.setVector(initial, centre);
  double signe = vector2.dotproduct(&vector1);
  if(signe<0)
    return -vector1.norm();
  else
    return vector1.norm();
}

 // Interpolation : Inverse distance weighting (IDW)
double gLevelsetTextile::operator()(const double &x, const double &y, const double &z) const
{
  TVertex initial(x,y,z);
  
  TVertex * searched;
  searched = _tisse.getBucket()->get_closest(&initial);
  
  Section * section =searched->getSection();
  assert(section!=NULL);
  TVertex * centre = section->getCentre();
  int index=centre->getIndex();
  
  Section * sectionBefor =section->getSectionBefor();
  Section * sectionAfter =section->getSectionAfter();
  
  assert(section->getNbNodes()>2);
  
  // Interpolation PID
  std::vector<TVertex*> v;
  
  // Genere la liste des points interpolés connus
  if (sectionBefor!=NULL)
  {
    v.push_back( sectionBefor->getNode(index-1) );
    v.push_back( sectionBefor->getNode(index) );
    v.push_back( sectionBefor->getNode(index+1) );
  }
  
  v.push_back( section->getNode(index-1) );
  v.push_back( section->getNode(index) );
  v.push_back( section->getNode(index+1) );
  
  if (sectionAfter!=NULL)
  {
    v.push_back( sectionAfter->getNode(index-1) );
    v.push_back( sectionAfter->getNode(index) );
    v.push_back( sectionAfter->getNode(index+1) );
  }
  
  // calcul de la valeur interpolée
  double p=1, denom=0, w=0;
  TVertex numer, temp;
  for (int i=0; i<v.size(); i++)
  {
    temp.setVector(v[i], &initial);
    w=1/pow(temp.norm(),p);
    numer+=w*(*v[i]);
    denom+=w;
  }
  denom=1/denom;
  TVertex interp(denom*numer);
  
  TVertex vector1;
  vector1.setVector(&interp, &initial);
  
  TVertex vector2;
  vector2.setVector(centre, searched);
  double signe = vector2.dotproduct(&vector1);
  if(signe<0)
    return -vector1.norm();
  else
    return vector1.norm();
}
*/