CGAL users discussion list

["Sebastien Loriot (GeometryFactory)" <>]

Do you expect the result to fit in a polyhedron?
In particular all edges should be adjacent to only 2 facets.

Sebastien.

On 11/14/2013 03:50 AM, costantinos wrote:
> I have a requirement to merge two Polyhedron_3 objects, one is closed
> geometry and the other is a surface. There are common vertices and edges
> between them. I have written a routine based on
> Polyhedron_incremental_builder_3, I am able to successfully combine the
> Polyhedron_3 objects (as checked by outputting results into OFF file),
> however assertion errors due to shared halfedges unfortunately inhibit the
> method. I am not even sure whether the approach I have taken works, I would
> be very grateful for any help/advice or a solution to my problem?
>
> Ideally, I would get to a point where I am able to combine the geometries
> by:
>
> CGAL::poly_merge<Polyhedron, HalfedgeDS> modifier(P1,P2);
> P3.delegate(modifier);
>
> or,
>
> CGAL::poly_merge<Polyhedron, HalfedgeDS> modifier(P1,P2);
>   P2.delegate(modifier);
>
> (I have checked orientation of surface but no effect + please note the code
> below has worked well so far for combining geometry but due to the sorting
> method there are exceptions!)
>
>
>
> int compare_Point_3(const Point_3* p, const Point_3* q)
> {
>      if ((fabs((*p)[0] - (*q)[0]) < GEOM_TOL) && (fabs((*p)[1] - (*q)[1]) <
> GEOM_TOL) && (fabs((*p)[2] - (*q)[2]) < GEOM_TOL))
>      {
>          return 0;
>      } else {
>          double delta = ((*p)[0] * (*p)[0] + (*p)[1] * (*p)[1] + (*p)[2] *
> (*p)[2]) - ((*q)[0] * (*q)[0] + (*q)[1] * (*q)[1] + (*q)[2] * (*q)[2]);
>          if (delta < 0) {
>              return -1;
>          } else {
>              return 1;
>          }
>      }
> }
>
> bool compare_lessthan(const Point_3* p, const Point_3* q)
> {
>      return compare_Point_3(p, q) == -1;
> }
>
>
>
> namespace CGAL {
>
>      template <class Poly, class HDS>
>      class poly_merge : public Modifier_base<HDS> {
>      public:
>          typedef Poly Polyhedron;
>          typedef HDS  Halfedge_data_structure;
>
>          poly_merge(const Polyhedron& poly_1, const Polyhedron& poly_2) :
> input_1(poly_1), input_2(poly_2) {}
>
>          void operator()(HDS& output);
>
>      protected:
>
>          const Poly& input_1;
>          const Poly& input_2;
>      };
>
>
>      template <class Poly, class HDS>
>      void poly_merge<Poly,HDS>:: operator()(HDS& output) {
>
>          typedef typename Poly::Vertex_const_iterator Vertex_const_iterator;
>          typedef typename Poly::Facet_const_iterator Facet_const_iterator;
>          typedef typename HDS::Vertex::Point Point;
>          typedef typename Poly::Point_3 Point_3;
>          typedef typename std::vector<const Point_3*> pPoints_3;
>          typename std::vector<const Point_3*>::iterator p_it;
>          typedef typename std::map<const Point_3*, int> Point_map;
>
>          pPoints_3 vert_input_1;
>          for (Vertex_const_iterator vi = input_1.vertices_begin(); vi !=
> input_1.vertices_end(); ++vi) vert_input_1.push_back(&(vi->point()));
>
>          std::sort(vert_input_1.begin(), vert_input_1.end(),
> compare_lessthan);
>
>
>          pPoints_3 vert_input_2;
>          for (Vertex_const_iterator vi = input_2.vertices_begin(); vi !=
> input_2.vertices_end(); ++vi) vert_input_2.push_back(&(vi->point()));
>
>          std::sort(vert_input_2.begin(), vert_input_2.end(),
> compare_lessthan);
>
>
>          Point_map pm;
>          pPoints_3 vert_input_2_diff, input_intersection;
>
>          int i = 0, j = 0;
>          while (i < vert_input_1.size() && j < vert_input_2.size()) {
>              std::cout << "i: " << i << "; j: " << j << std::endl;
>              int c = compare_Point_3(vert_input_1[i], vert_input_2[j]);
>              std::cout << "c: " << c << std::endl;
>              if (c == -1) {
>                  ++i;
>              } else if (c == 1) {
>                  vert_input_2_diff.push_back(vert_input_2[j]);
>                  ++j;
>              } else if (c == 0) {
>                  input_intersection.push_back(vert_input_1[i]);
>                  pm[vert_input_2[j]] = i;
>                 // std::cout << "Intersection: " << i << std::endl;
>                  ++i; ++j;
>              }
>          }
>
>          if (vert_input_2_diff.size() != vert_input_2.size())
> vert_input_2_diff.push_back(vert_input_2[vert_input_2.size() - 1]);
>
>
>          /*
>          std::cout << "Input 1 verts: " << vert_input_1.size() << "\n";
>          for (p_it = vert_input_1.begin(); p_it != vert_input_1.end();
> ++p_it) std::cout << *(*p_it) << std::endl;
>          std::cout << std::endl << std::endl;
>          std::cout << "Input 2 verts: " << vert_input_2.size() << "\n";
>          for (p_it = vert_input_2.begin(); p_it != vert_input_2.end();
> ++p_it) std::cout << *(*p_it) << std::endl;
>          std::cout << std::endl << std::endl;
>          std::cout << "Intersect points: " << input_intersection.size() <<
> "\n";
>         for (p_it = input_intersection.begin(); p_it !=
> input_intersection.end(); ++p_it) std::cout << *(*p_it) << std::endl;
>          std::cout << std::endl << std::endl;
>          std::cout << "Input 2 diff: " << vert_input_2_diff.size() << "\n";
>          for (p_it = vert_input_2_diff.begin(); p_it !=
> vert_input_2_diff.end(); ++p_it) std::cout << *(*p_it) << std::endl;
>          std::cout << std::endl << std::endl;
>          std::cout << "intersection map\n"
>          "Printing point and position from original vector through pointer"
> << std::endl;
>          for (auto it = pm.begin(); it != pm.end(); ++it) {
>              std::cout << "Position: " << it->second << std::endl;
>              const Point_3* p = it->first;
>              std::cout << *p << std::endl;
>          }
>          */
>
>          Polyhedron_incremental_builder_3<HDS> B(output, true);
>          B.begin_surface(vert_input_1.size() + vert_input_2_diff.size(),
> input_1.size_of_facets() + input_2.size_of_facets());
>          //std::cout << "OFF" << std::endl;
>          //std::cout << (vert_input_1.size() + vert_input_2_diff.size()) << "
> " << (input_1.size_of_facets() + input_2.size_of_facets())  << " 0" <<
> std::endl;
>          std::cout << std::endl;
>
>          for (i = 0; i < vert_input_1.size(); ++i)
>          {
>              B.add_vertex(Point_3(*(vert_input_1[i])));
>              pm[vert_input_1[i]] = i;
>
>              //std::cout << *(vert_input_1[i]) << std::endl;
>              //std::cout << "Index 1: " << pm[vert_input_1[i]] << std::endl;
>          }
>          std::cout << std::endl;
>
>          for (j = 0; j < vert_input_2_diff.size(); ++j, ++i)
>          {
>              B.add_vertex(Point_3(*(vert_input_2_diff[j])));
>              pm[vert_input_2_diff[j]] = i;
>
>              //std::cout << *(vert_input_2_diff[j]) << std::endl;
>              //std::cout << "Index 1: " << pm[vert_input_2_diff[j]] <<
> std::endl;
>          }
>
>          for (Facet_const_iterator fi = input_1.facets_begin(); fi !=
> input_1.facets_end(); ++fi)
>          {
>              //std::cout << "3 ";
>              B.begin_facet();
>              typedef typename Poly::Halfedge_around_facet_const_circulator
>              Halfedge_around_facet_const_circulator;
>              Halfedge_around_facet_const_circulator hc = fi->facet_begin();
>              Halfedge_around_facet_const_circulator hc_end = hc;
>              CGAL_assertion( hc != NULL);
>              do {
>                  B.add_vertex_to_facet(pm[&(hc->vertex()->point())]);
>                  //std::cout << "Point around face 1: " <<
> hc->vertex()->point() << std::endl;
>                  //std::cout << pm[&(hc->vertex()->point())] << " ";
>                  ++hc;
>              } while( hc != hc_end);
>              B.end_facet();
>
>              //std::cout << std::endl;
>          }
>
>          for (Facet_const_iterator fi = input_2.facets_begin(); fi !=
> input_2.facets_end(); ++fi)
>          {
>              //std::cout << "3 ";
>              B.begin_facet();
>              typedef typename Poly::Halfedge_around_facet_const_circulator
>              Halfedge_around_facet_const_circulator;
>              Halfedge_around_facet_const_circulator hc = fi->facet_begin();
>              Halfedge_around_facet_const_circulator hc_end = hc;
>              CGAL_assertion( hc != NULL);
>              do {
>                  B.add_vertex_to_facet(pm[&(hc->vertex()->point())]);
>                  //std::cout << "Point around face 2: " <<
> hc->vertex()->point() << std::endl;
>                  //std::cout << pm[&(hc->vertex()->point())] << " ";
>                  ++hc;
>              } while( hc != hc_end);
>              B.end_facet();
>              //std::cout << std::endl;
>          }
>
>          B.end_surface();
>          output.normalize_border();
>      }
>
> }
>
>
>
> --
> View this message in context: 
> http://cgal-discuss.949826.n4.nabble.com/Merging-Combining-two-Polyhedron-3-with-shared-vertices-and-edges-tp4658412.html
> Sent from the cgal-discuss mailing list archive at Nabble.com.