Now it should work with overlapping polylines

Polyline_simplification_2/examples/Polyline_simplification_2/simplify_overlapping_polyline.cpp

@@ -0,0 +1,51 @@
+#include <iostream>
+#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+#include <CGAL/Polygon_2.h>
+#include <CGAL/Constrained_Delaunay_triangulation_2.h>
+#include <CGAL/Constrained_triangulation_plus_2.h>
+#include <CGAL/Polyline_simplification_2/simplify.h>
+#include <CGAL/Polyline_simplification_2/Squared_distance_cost.h>
+
+namespace PS = CGAL::Polyline_simplification_2;
+
+typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
+
+typedef PS::Vertex_base_2<K> Vb;
+typedef CGAL::Constrained_triangulation_face_base_2<K> Fb;
+typedef CGAL::Triangulation_data_structure_2<Vb, Fb> TDS;
+typedef CGAL::Constrained_Delaunay_triangulation_2<K, TDS, CGAL::Exact_predicates_tag> CDT;
+typedef CGAL::Constrained_triangulation_plus_2<CDT>     CT;
+
+typedef CT::Point                           Point;
+typedef CT::Constraint_iterator             Constraint_iterator;
+typedef CT::Vertices_in_constraint_iterator Vertices_in_constraint_iterator;
+typedef CT::Points_in_constraint_iterator   Points_in_constraint_iterator;
+typedef PS::Stop_below_count_ratio_threshold Stop;
+typedef PS::Squared_distance_cost Cost;
+
+int main()
+{
+  CT ct;
+
+  std::vector<Point> P = { Point(0,1), Point(1,1), Point(2,2), Point(3,1), Point(4,1), Point(5,1) };
+  std::vector<Point> Q = { Point(5,0), Point(4,1), Point(3,1), Point(2,2), Point(1,1), Point(0,0) };
+  std::vector<Point> R = { Point(3,1), Point(4,1) };
+
+  ct.insert_constraint(P);
+  ct.insert_constraint(Q);
+  ct.insert_constraint(R);
+
+  PS::simplify(ct, Cost(), Stop(0.5));
+
+  for(Constraint_iterator cit = ct.constraints_begin();
+      cit != ct.constraints_end();
+      ++cit) {
+    std::cout << "simplified polyline" << std::endl;
+    for(Points_in_constraint_iterator vit =
+          ct.points_in_constraint_begin(*cit);
+        vit != ct.points_in_constraint_end(*cit);
+        ++vit)
+      std::cout << *vit << std::endl;
+  }
+  return 0;
+}

Polyline_simplification_2/include/CGAL/Polyline_simplification_2/simplify.h

@@ -101,7 +101,6 @@ public:
     typedef boost::readable_property_map_tag category;
     typedef std::size_t                      value_type;
     typedef value_type                       reference;
-    // typedef Vertices_in_constraint_iterator  key_type;
     typedef Vertex_handle  key_type;
 
     Vertex_index_map* vertex_index_map;
@@ -112,12 +111,10 @@ public:
 
     reference operator[] ( key_type const& x ) const
     {
-      //  return x.base()->id ;
       return (*vertex_index_map)[x];
     }
   } ;
 
-  //typedef CGAL::Modifiable_priority_queue<Vertices_in_constraint_iterator,Compare_cost,Id_map> MPQ ;
   typedef CGAL::Modifiable_priority_queue<Vertex_handle,Compare_cost,Id_map> MPQ ;
 
   MPQ* mpq;
@@ -189,10 +186,20 @@ public:
           it != pct.vertices_in_constraint_end(cid);
           ++it){
         if((*it)->is_removable()){
+          typename std::map<Vertex_handle, std::list<Vertices_in_constraint_iterator> >::iterator lit;
+          lit = vertex_to_iterator.find(*it);
+
+          if(lit != vertex_to_iterator.end()){
+            std::list<Vertices_in_constraint_iterator>& ilist = lit->second;
+            if(std::find(ilist.begin(),ilist.end(),it) == ilist.end()){
+              ilist.push_back(it);
+            }
+          }else{
             vertex_to_iterator[*it].push_back(it);
           }
         }
       }
+    }
 
     /*
     //  debug output
@@ -202,26 +209,6 @@ public:
       }
     }
     */
-
-    /*
-    // the previous code does mark vertices with more than one constraint on it as unremovable
-    std::set<Vertex_handle> vertices;
-    Constraint_iterator cit = pct.constraints_begin(), e = pct.constraints_end();
-    for(; cit!=e; ++cit){
-      Constraint_id cid = *cit;
-      Vertices_in_constraint_iterator it = pct.vertices_in_constraint_begin(cid);
-      (*it)->set_removable(false);
-      for(; it != pct.vertices_in_constraint_end(cid); ++it){
-        if(vertices.find(*it) != vertices.end()){
-          (*it)->set_removable(false);
-        } else {
-          vertices.insert(*it);
-        }
-      }
-      it = boost::prior(it);
-      (*it)->set_removable(false);
-    }
-    */
   }
 
   // For all polyline constraints we compute the cost of all unremovable and not removed vertices
@@ -236,7 +223,9 @@ public:
         boost::optional<FT> dist = cost(pct, it);
         if(dist){
           (*it)->set_cost(*dist);
+          if(! (*mpq).contains(*it)){
               (*mpq).push(*it);
+            }
           ++n;
         } else {
           // no need to set the costs as this vertex is not in the priority queue
@@ -342,6 +331,7 @@ operator()()
   if(stop(pct, v, v->cost(), pct_initial_number_of_vertices, pct.number_of_vertices())){
     return false;
   }
+
   Vertices_in_constraint_iterator vit = vertex_to_iterator[v].front();
   if(is_removable(vit)){
     Vertices_in_constraint_iterator u = boost::prior(vit), w = boost::next(vit);

Triangulation_2/include/CGAL/Triangulation_2/internal/Polyline_constraint_hierarchy_2.h

@@ -46,14 +46,13 @@ private:
   class Node {
   public:
     explicit Node(Vertex_handle vh, bool input = false)
-      : vertex_(vh), /* id(-1),*/ input(input)
+      : vertex_(vh), input(input)
     {}
     const Point& point() const { return vertex_->point(); }
     Vertex_handle vertex() const { return vertex_; }
   private:
     Vertex_handle vertex_;
   public:
-    // int id;
     bool input;
   };
 
@@ -585,45 +584,58 @@ void Polyline_constraint_hierarchy_2<T,Compare,Point>::simplify(Vertex_it uc,
                                                        Vertex_it wc)
 
 {
+  // TODO: How do we (want to) deal with u == w ???
   Vertex_handle u = *uc, v = *vc, w = *wc;
   typename Sc_to_c_map::iterator uv_sc_iter = sc_to_c_map.find(make_edge(u, v));
-  CGAL_assertion_msg( uv_sc_iter != sc_to_c_map.end(), "not a subconstraint" );
+  typename Sc_to_c_map::iterator vw_sc_iter = sc_to_c_map.find(make_edge(v, w));
   Context_list*  uv_hcl = uv_sc_iter->second;
-  CGAL_assertion_msg((u == w) || (uv_hcl->size() == 1), "more than one constraint passing through the subconstraint" );
-
-  if(*(uv_hcl->front().current()) != u) {
-    std::swap(u,w);
-    uv_sc_iter = sc_to_c_map.find(make_edge(u, v));
-    CGAL_assertion_msg( uv_sc_iter != sc_to_c_map.end(), "not a subconstraint" );
-    uv_hcl = (*uv_sc_iter).second;
-    CGAL_assertion_msg((u == w) || (uv_hcl->size() == 1), "more than one constraint passing through the subconstraint" );
-  }
-  // now u,v, and w are ordered along the polyline constraint
+  Context_list*  vw_hcl = vw_sc_iter->second;
+  // AF:  what is input() about???
   if(vc.input()){
     uc.input() = true;
     wc.input() = true;
   }
-  typename Sc_to_c_map::iterator vw_sc_iter = sc_to_c_map.find(make_edge(v, w));
-  CGAL_assertion_msg( vw_sc_iter != sc_to_c_map.end(), "not a subconstraint" );
-  Context_list*  vw_hcl = vw_sc_iter->second;
-    CGAL_assertion_msg((u == w) || (vw_hcl->size() == 1), "more than one constraint passing through the subconstraint" );
 
-  Vertex_list* vertex_list = uv_hcl->front().id().vl_ptr();
-  CGAL_assertion_msg(vertex_list  == vw_hcl->front().id().vl_ptr(), "subconstraints from different polyline constraints" );
+  // Take contexts from the two context lists depending on the orientation of the constraints
+  // These are the contexts where current is either u or w
+  // remove from uv_hcl the contexts where current is not u
+  // remove from vw_hcl the contexts where current is not w
+  // splice into uv_hcl
+  typename Context_list::iterator it = uv_hcl->begin();
+  while(it != uv_hcl->end()){
+    if((*it->current()) != u){
+      it = uv_hcl->erase(it);
+    }else{
       // Remove the list item which points to v
-  vertex_list->skip(vc.base());
+      Vertex_list* vertex_list = it->id().vl_ptr();
+      Vertex_it vc_in_context = it->current();
+      vc_in_context = boost::next(vc_in_context);
+      vertex_list->skip(vc_in_context.base());
+      ++it;
+    }
+  }
+  it = vw_hcl->begin();
+  while(it != vw_hcl->end()){
+    if((*it->current()) != w){
+      it = vw_hcl->erase(it);
+    }else{
+      // Remove the list item which points to v
+      Vertex_list* vertex_list = it->id().vl_ptr();
+      Vertex_it vc_in_context = it->current();
+      vc_in_context = boost::next(vc_in_context);
+      vertex_list->skip(vc_in_context.base());
+      ++it;
+    }
+  }
+
+  uv_hcl->splice(uv_hcl->end(),*vw_hcl);
+  delete vw_hcl;
 
-  if(u != w){
-    // Remove the entries for [u,v] and [v,w]
   sc_to_c_map.erase(uv_sc_iter);
   sc_to_c_map.erase(vw_sc_iter);
-    delete vw_hcl;
+
   // reuse other context list
   sc_to_c_map[make_edge(u,w)] = uv_hcl;
-  }else{
-    sc_to_c_map.erase(uv_sc_iter);
-    delete vw_hcl;
-  }
 }