CGAL users discussion list

[MrVH <>]

I'm trying to understand your mapping code and to modify it. But I've got
some compilation errors and I don't catch some lines. I really have to use
this code cause to avoid to recode AlphaShape Weighted from scratch in
Python :(
I don't see where the Delaunay function is launched. 

In main I'm running:
Alpha_shape_3  as;
build_triangulation_with_indices(lwp.begin(),lwp.end(),as);

The build_triangulation_with_indices function take now an Alpha_Shape_3
object:
//function inserting points into a triangulation
//and setting the info field to the order in the input list.
template <class Alpha_Shape_3,class Point_iterator>
void
build_triangulation_with_indices(
        Point_iterator begin,
        Point_iterator end,
        Alpha_Shape_3& T)
{
// The vector is now a Weighted_point vector cause it is regular
triangulation
  std::vector<std::pair<const typename Gt::Weighted_point*,int> > points;
  int index=-1;
  for (Point_iterator it=begin;it!=end;++it){
    points.push_back(std::make_pair(&(*it),++index));
  }
  std::random_shuffle (points.begin(), points.end()); // shuffle point but
why ???
  spatial_sort (points.begin(),points.end(),
Traits_for_spatial_sort<Alpha_Shape_3>()); // build a structure struct
Traits_for_spatial_sort:public Alpha_Shape_3::Geom_traits

// from here I do not understand. Where Trigulation is launched ???
  typename Alpha_Shape_3::Cell_handle hint;
  for (typename std::vector< std::pair<const typename
Gt::Weighted_point*,int> >::const_iterator
        p = points.begin();p != points.end(); ++p)
  {
    typename Alpha_Shape_3::Locate_type lt;
    typename Alpha_Shape_3::Cell_handle c;
    int li, lj;
    c = T.locate (*(p->first), lt, li, lj, hint);

    typename Alpha_Shape_3::Vertex_handle v = T.insert (*(p->first), lt, c,
li, lj);
    if ( v==typename Alpha_Shape_3::Vertex_handle() )
      hint=c;
    else{
      v->info() = p->second ;
      hint=v->cell();
    }
  }
}  

the structure
//typedef CGAL::Exact_predicates_inexact_constructions_kernel K;

typedef CGAL::Triangulation_vertex_base_with_info_3<int, K> Wb; // avoid
conflict so Wb instead of Vb cause casting is different than the previous
<int, K>
typedef CGAL::Triangulation_data_structure_3<Wb>                    Tds2;
typedef CGAL::Delaunay_triangulation_3<K, Tds2>                     
Delaunay; // avoir conflict

typedef Delaunay::Point   Point;

//spatial sort traits to use with a pair of point pointers and integer.
template<class Alpha_Shape_3>
struct Traits_for_spatial_sort:public Alpha_Shape_3::Geom_traits{
  typedef typename Alpha_Shape_3::Geom_traits Gt;
  typedef std::pair<const typename Gt::Weighted_point*,int> Point_3;
  
  struct Less_x_3{
    bool operator()(const Point_3& p,const Point_3& q) const {
      return typename Gt::Less_x_3()(*(p.first),*(q.first));
    }
  };

  struct Less_y_3{
    bool operator()(const Point_3& p,const Point_3& q) const {
      return typename Gt::Less_y_3()(*(p.first),*(q.first));
    }
  };

  struct Less_z_3{
    bool operator()(const Point_3& p,const Point_3& q) const {
      return typename Gt::Less_z_3()(*(p.first),*(q.first));
    }
  };  
  
  Less_x_3      less_x_3_object () const {return Less_x_3();}
  Less_y_3      less_y_3_object () const {return Less_y_3();}
  Less_z_3      less_z_3_object () const {return Less_z_3();}
};


Thank you for your help and your explanations,




############################################################
The full code:
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Regular_triangulation_euclidean_traits_3.h>
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Alpha_shape_3.h>
#include <fstream>
#include <iostream>
#include <cstdio>
#include <cstring>
#include <string>
#include <map>
#include <set>
#include <vector>
#include <cmath>
#include <stdio.h>
#include <stdlib.h>


//#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_3.h>
#include <CGAL/Triangulation_vertex_base_with_info_3.h>
//#include <iostream>


typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Regular_triangulation_euclidean_traits_3<K> Gt;
typedef CGAL::Alpha_shape_vertex_base_3<Gt>         Vb;
typedef CGAL::Alpha_shape_cell_base_3<Gt>           Fb;
typedef CGAL::Triangulation_data_structure_3<Vb,Fb> Tds;
typedef CGAL::Regular_triangulation_3<Gt,Tds>       Triangulation_3;
typedef CGAL::Alpha_shape_3<Triangulation_3>        Alpha_shape_3;

typedef Alpha_shape_3::Cell_handle          Cell_handle;
typedef Alpha_shape_3::Vertex_handle        Vertex_handle;
typedef Alpha_shape_3::Facet                Facet;
typedef Alpha_shape_3::Edge                 Edge;

typedef Gt::Weighted_point                  Weighted_point;
typedef Gt::Bare_point                      Bare_point;












//typedef CGAL::Exact_predicates_inexact_constructions_kernel K;

typedef CGAL::Triangulation_vertex_base_with_info_3<int, K> Wb;
typedef CGAL::Triangulation_data_structure_3<Wb>                    Tds2;
typedef CGAL::Delaunay_triangulation_3<K, Tds2>                     
Delaunay;

typedef Delaunay::Point   Point;

//spatial sort traits to use with a pair of point pointers and integer.
template<class Alpha_Shape_3>
struct Traits_for_spatial_sort:public Alpha_Shape_3::Geom_traits{
  typedef typename Alpha_Shape_3::Geom_traits Gt;
  typedef std::pair<const typename Gt::Weighted_point*,int> Point_3;
  
  struct Less_x_3{
    bool operator()(const Point_3& p,const Point_3& q) const {
      return typename Gt::Less_x_3()(*(p.first),*(q.first));
    }
  };

  struct Less_y_3{
    bool operator()(const Point_3& p,const Point_3& q) const {
      return typename Gt::Less_y_3()(*(p.first),*(q.first));
    }
  };

  struct Less_z_3{
    bool operator()(const Point_3& p,const Point_3& q) const {
      return typename Gt::Less_z_3()(*(p.first),*(q.first));
    }
  };  
  
  Less_x_3      less_x_3_object () const {return Less_x_3();}
  Less_y_3      less_y_3_object () const {return Less_y_3();}
  Less_z_3      less_z_3_object () const {return Less_z_3();}
};


//function inserting points into a triangulation
//and setting the info field to the order in the input list.
template <class Alpha_Shape_3,class Point_iterator>
void
build_triangulation_with_indices(
        Point_iterator begin,
        Point_iterator end,
        Alpha_Shape_3& T)
{
  std::vector<std::pair<const typename Gt::Weighted_point*,int> > points;
  int index=-1;
  for (Point_iterator it=begin;it!=end;++it){
    points.push_back(std::make_pair(&(*it),++index));
  }
  std::random_shuffle (points.begin(), points.end());
  spatial_sort (points.begin(),points.end(),
Traits_for_spatial_sort<Alpha_Shape_3>());

  typename Alpha_Shape_3::Cell_handle hint;
  for (typename std::vector< std::pair<const typename
Gt::Weighted_point*,int> >::const_iterator
        p = points.begin();p != points.end(); ++p)
  {
    typename Alpha_Shape_3::Locate_type lt;
    typename Alpha_Shape_3::Cell_handle c;
    int li, lj;
    c = T.locate (*(p->first), lt, li, lj, hint);

    typename Alpha_Shape_3::Vertex_handle v = T.insert (*(p->first), lt, c,
li, lj);
    if ( v==typename Alpha_Shape_3::Vertex_handle() )
      hint=c;
    else{
      v->info() = p->second ;
      hint=v->cell();
    }
  }
}  













int main(int argc, char *argv[]) 
{
 
  std::list<Weighted_point> lwp;
     
  /*Get arguments*/
  char *filename = argv[1];    // Format of file %8.3f%8.3f%8.3f%8.3f, x, y,
z, weight 
  float sizeas = atof(argv[2]);  // size of AlphaSphere
  /*----------------*/
  
    
  std::cout << "# File" << filename << std::endl;
  
  /*Read the file*/
  std::ifstream fcoords;
  fcoords.open(filename);
  std::string line;
  while (std::getline(fcoords, line)){
    
        /* Get coordinates*/
        std::string type(line.substr(0,1));
        if (type == "#" || type == "\n"){continue;} // comment line
        
        std::string xstr(line.substr(0,8));
        std::string ystr(line.substr(8,8));
        std::string zstr(line.substr(16,8));
        std::string weightstr(line.substr(24,8));
        //std::cout << xstr << ystr << zstr << weightstr << std::endl;

        float x = atof(xstr.c_str());
        float y = atof(ystr.c_str());
        float z = atof(zstr.c_str());
        float weight = atof(weightstr.c_str());
        
        /* fill each point with its associate weight */
        lwp.push_back(Weighted_point(Bare_point( x, y, z), weight));
       
  }
  fcoords.close();
  
  
  std::cout << "# read " << lwp.size() << " atoms" << std::endl;
  
  /* Launch AlphaShape*/
  //Alpha_shape_3  as(lwp.begin(), lwp.end(), sizeas,
Alpha_shape_3::GENERAL);
  Alpha_shape_3  as;
  build_triangulation_with_indices(lwp.begin(),lwp.end(),as);
  
  
  /* Create empty storage*/
  std::list<Cell_handle>       cells;
  std::list<Facet>       facets;
  std::list<Edge>        edges;
  
  /*Get information 
    Cell -> Tetraedron
    Facet -> Triangle
    Edge -> segment
  */
 
as.get_alpha_shape_cells(std::back_inserter(cells),Alpha_shape_3::INTERIOR);
 
as.get_alpha_shape_facets(std::back_inserter(facets),Alpha_shape_3::REGULAR);
 
as.get_alpha_shape_facets(std::back_inserter(facets),Alpha_shape_3::SINGULAR);
 
as.get_alpha_shape_edges(std::back_inserter(edges),Alpha_shape_3::SINGULAR);
  
  
  /* Enumeration */
  std::map<Alpha_shape_3::Vertex*,unsigned int> index;
  
  
  unsigned int i=0;
  for (Alpha_shape_3::Finite_vertices_iterator
vit=as.finite_vertices_begin();vit!=as.finite_vertices_end();++vit)
    index.insert(std::make_pair(&(*vit),i++));
  
  /* Tetrahedron*/
  for (std::list<Cell_handle>::iterator
cell=cells.begin();cell!=cells.end();++cell){
    std::cout << "cell ";
    //for (int k=0;k<4;++k) std::cout << index[ &( *((*cell)->vertex(k)) ) ]
<< " ";
    //std::cout <<std::endl;
    for (int k=0;k<4;++k) std::cout <<  (*cell)->vertex(k)->point() << " ";
    std::cout <<std::endl; 
  }
  
  /* Triangle */
  for (std::list<Facet>::iterator
facet=facets.begin();facet!=facets.end();++facet){
    std::cout << "facet ";
    for (int k=0;k<3;++k)
      std::cout << (facet->first->vertex( (facet->second+k)%4 ))->point() <<
" ";
      //std::cout << index[ &(*(facet->first->vertex( (facet->second+k)%4
)))] << " ";
    std::cout <<std::endl;
  }

  /* Edge */
  for (std::list<Edge>::iterator
edge=edges.begin();edge!=edges.end();++edge)
    {
    std::cout << "edge ";
    std::cout << (edge->first->vertex( edge->second ))->point() <<" ";
    std::cout << (edge->first->vertex( edge->third ))->point() << std::endl;
    //std::cout << index[ &(*(edge->first->vertex( edge->second )))] <<" ";
    //std::cout << index[ &(*(edge->first->vertex( edge->third )))]
<<std::endl;
    }
  
  std::cout << "# The 0-shape has : " << std::endl;
  std::cout << "# " << cells.size() << " interior tetrahedra" << std::endl;
  std::cout << "# " << facets.size() << " boundary facets" << std::endl;
  std::cout << "# " << edges.size()  << " singular edges" << std::endl;
  return 1;
 

}


-- 
View this message in context: 
http://cgal-discuss.949826.n4.nabble.com/Screen-easely-Cells-Facets-Edges-AlphaShape3D-Weighted-tp3297041p3300898.html
Sent from the cgal-discuss mailing list archive at Nabble.com.