Variable-width gap fill. Yay! #2960

xs/src/libslic3r/ExPolygon.cpp

@@ -6,7 +6,6 @@
 #include "ClipperUtils.hpp"
 #include "polypartition.h"
 #include "poly2tri/poly2tri.h"
-
 #include <algorithm>
 #include <list>
 
@@ -171,10 +170,11 @@ ExPolygon::simplify(double tolerance, ExPolygons* expolygons) const
 }
 
 void
-ExPolygon::medial_axis(double max_width, double min_width, Polylines* polylines) const
+ExPolygon::medial_axis(double max_width, double min_width, ThickPolylines* polylines) const
 {
     // init helper object
     Slic3r::Geometry::MedialAxis ma(max_width, min_width);
+    ma.expolygon = this;
 
     // populate list of segments for the Voronoi diagram
     ma.lines = this->contour.lines();
@@ -184,41 +184,71 @@ ExPolygon::medial_axis(double max_width, double min_width, Polylines* polylines)
     }
 
     // compute the Voronoi diagram
-    Polylines pp;
+    ThickPolylines pp;
     ma.build(&pp);
 
-    // clip segments to our expolygon area
-    // (do this before extending endpoints as external segments coule be extended into
-    // expolygon, this leaving wrong things inside)
-    pp = intersection(pp, *this);
-
-    // extend initial and final segments of each polyline (they will be clipped)
-    // unless they represent closed loops
-    for (Polylines::iterator polyline = pp.begin(); polyline != pp.end(); ++polyline) {
-        if (polyline->points.front().distance_to(polyline->points.back()) < min_width) continue;
-        // TODO: we should *not* extend endpoints where other polylines start/end
-        // (such as T joints, which are returned as three polylines by MedialAxis)
-        polyline->extend_start(max_width);
-        polyline->extend_end(max_width);
-    }
+    /*
+    SVG svg("medial_axis.svg");
+    svg.draw(*this);
+    svg.draw(pp);
+    svg.Close();
+    */
 
-    // clip again after extending endpoints to prevent them from exceeding the expolygon boundaries
-    pp = intersection(pp, *this);
-
-    // remove too short polylines
-    // (we can't do this check before endpoints extension and clipping because we don't
-    // know how long will the endpoints be extended since it depends on polygon thickness
-    // which is variable - extension will be <= max_width/2 on each side)
     for (size_t i = 0; i < pp.size(); ++i) {
-        if (pp[i].length() < max_width) {
+        ThickPolyline& polyline = pp[i];
+
+        // extend initial and final segments of each polyline if they're actual endpoints
+        /* We assign new endpoints to temporary variables because in case of a single-line
+           polyline, after we extend the start point it will be caught by the intersection()
+           call, so we keep the inner point until we perform the second intersection() as well */
+        Point new_front = polyline.points.front();
+        Point new_back  = polyline.points.back();
+        if (polyline.endpoints.front() && !this->has_boundary_point(new_front)) {
+            Line line(polyline.points.front(), polyline.points[1]);
+
+            // prevent the line from touching on the other side, otherwise intersection() might return that solution
+            if (polyline.points.size() == 2) line.b = line.midpoint();
+
+            line.extend_start(max_width);
+            (void)this->contour.intersection(line, &new_front);
+        }
+        if (polyline.endpoints.back() && !this->has_boundary_point(new_back)) {
+            Line line(
+                *(polyline.points.end() - 2),
+                polyline.points.back()
+            );
+
+            // prevent the line from touching on the other side, otherwise intersection() might return that solution
+            if (polyline.points.size() == 2) line.a = line.midpoint();
+            line.extend_end(max_width);
+
+            (void)this->contour.intersection(line, &new_back);
+        }
+        polyline.points.front() = new_front;
+        polyline.points.back()  = new_back;
+
+        /*  remove too short polylines
+            (we can't do this check before endpoints extension and clipping because we don't
+            know how long will the endpoints be extended since it depends on polygon thickness
+            which is variable - extension will be <= max_width/2 on each side)  */
+        if (polyline.length() < max_width) {
             pp.erase(pp.begin() + i);
             --i;
+            continue;
         }
     }
 
     polylines->insert(polylines->end(), pp.begin(), pp.end());
 }
 
+void
+ExPolygon::medial_axis(double max_width, double min_width, Polylines* polylines) const
+{
+    ThickPolylines tp;
+    this->medial_axis(max_width, min_width, &tp);
+    polylines->insert(polylines->end(), tp.begin(), tp.end());
+}
+
 void
 ExPolygon::get_trapezoids(Polygons* polygons) const
 {

xs/src/libslic3r/ExPolygon.hpp

@@ -32,6 +32,7 @@ class ExPolygon
     Polygons simplify_p(double tolerance) const;
     ExPolygons simplify(double tolerance) const;
     void simplify(double tolerance, ExPolygons* expolygons) const;
+    void medial_axis(double max_width, double min_width, ThickPolylines* polylines) const;
     void medial_axis(double max_width, double min_width, Polylines* polylines) const;
     void get_trapezoids(Polygons* polygons) const;
     void get_trapezoids(Polygons* polygons, double angle) const;

xs/src/libslic3r/ExtrusionEntity.cpp

@@ -114,7 +114,7 @@ Polygons
 ExtrusionPath::grow() const
 {
     Polygons pp;
-    offset(this->polyline, &pp, +this->width/2);
+    offset(this->polyline, &pp, +scale_(this->width/2));
     return pp;
 }
 

xs/src/libslic3r/ExtrusionEntity.hpp

@@ -104,6 +104,10 @@ class ExtrusionLoop : public ExtrusionEntity
     ExtrusionLoop(ExtrusionLoopRole role = elrDefault) : role(role) {};
     ExtrusionLoop(const ExtrusionPaths &paths, ExtrusionLoopRole role = elrDefault)
         : paths(paths), role(role) {};
+    ExtrusionLoop(const ExtrusionPath &path, ExtrusionLoopRole role = elrDefault)
+        : role(role) {
+        this->paths.push_back(path);
+    };
     bool is_loop() const {
         return true;
     };

xs/src/libslic3r/ExtrusionEntityCollection.cpp

@@ -101,6 +101,20 @@ ExtrusionEntityCollection::append(const ExtrusionPaths &paths)
         this->append(*path);
 }
 
+void
+ExtrusionEntityCollection::replace(size_t i, const ExtrusionEntity &entity)
+{
+    delete this->entities[i];
+    this->entities[i] = entity.clone();
+}
+
+void
+ExtrusionEntityCollection::remove(size_t i)
+{
+    delete this->entities[i];
+    this->entities.erase(this->entities.begin() + i);
+}
+
 ExtrusionEntityCollection
 ExtrusionEntityCollection::chained_path(bool no_reverse, std::vector<size_t>* orig_indices) const
 {

xs/src/libslic3r/ExtrusionEntityCollection.hpp

@@ -36,6 +36,8 @@ class ExtrusionEntityCollection : public ExtrusionEntity
     void append(const ExtrusionEntity &entity);
     void append(const ExtrusionEntitiesPtr &entities);
     void append(const ExtrusionPaths &paths);
+    void replace(size_t i, const ExtrusionEntity &entity);
+    void remove(size_t i);
     ExtrusionEntityCollection chained_path(bool no_reverse = false, std::vector<size_t>* orig_indices = NULL) const;
     void chained_path(ExtrusionEntityCollection* retval, bool no_reverse = false, std::vector<size_t>* orig_indices = NULL) const;
     void chained_path_from(Point start_near, ExtrusionEntityCollection* retval, bool no_reverse = false, std::vector<size_t>* orig_indices = NULL) const;

xs/src/libslic3r/Geometry.cpp

@@ -5,12 +5,13 @@
 #include "PolylineCollection.hpp"
 #include "clipper.hpp"
 #include <algorithm>
+#include <cassert>
 #include <cmath>
 #include <list>
 #include <map>
 #include <set>
+#include <utility>
 #include <vector>
-#include <assert.h>
 
 #ifdef SLIC3R_DEBUG
 #include "SVG.hpp"
@@ -289,19 +290,8 @@ arrange(size_t total_parts, Pointf part, coordf_t dist, const BoundingBoxf* bb)
     return positions;
 }
 
-Line
-MedialAxis::edge_to_line(const VD::edge_type &edge) const
-{
-    Line line;
-    line.a.x = edge.vertex0()->x();
-    line.a.y = edge.vertex0()->y();
-    line.b.x = edge.vertex1()->x();
-    line.b.y = edge.vertex1()->y();
-    return line;
-}
-
 void
-MedialAxis::build(Polylines* polylines)
+MedialAxis::build(ThickPolylines* polylines)
 {
     /*
     // build bounding box (we use it for clipping infinite segments)
@@ -317,7 +307,7 @@ MedialAxis::build(Polylines* polylines)
         for (VD::const_edge_iterator edge = this->vd.edges().begin(); edge != this->vd.edges().end(); ++edge) {
             if (edge->is_infinite()) continue;
             
-            Polyline polyline;
+            ThickPolyline polyline;
             polyline.points.push_back(Point( edge->vertex0()->x(), edge->vertex0()->y() ));
             polyline.points.push_back(Point( edge->vertex1()->x(), edge->vertex1()->y() ));
             polylines->push_back(polyline);
@@ -373,7 +363,7 @@ MedialAxis::build(Polylines* polylines)
         assert(vertex_edges[v].size() == 1);
         edge_t* edge = *vertex_edges[v].begin();
 
-        if (!this->is_valid_edge(*edge)) {
+        if (!this->validate_edge(*edge)) {
             // if edge is not valid, erase it and its twin from edge list
             (void)this->edges.erase(edge);
             (void)this->edges.erase(edge->twin());
@@ -400,22 +390,36 @@ MedialAxis::build(Polylines* polylines)
         edge_t &edge = **this->edges.begin();
 
         // start a polyline
-        Polyline polyline;
+        ThickPolyline polyline;
         polyline.points.push_back(Point( edge.vertex0()->x(), edge.vertex0()->y() ));
         polyline.points.push_back(Point( edge.vertex1()->x(), edge.vertex1()->y() ));
+        polyline.width.push_back(this->thickness[&edge].first);
+        polyline.width.push_back(this->thickness[&edge].second);
 
         // remove this edge and its twin from the available edges
         (void)this->edges.erase(&edge);
         (void)this->edges.erase(edge.twin());
 
         // get next points
-        this->process_edge_neighbors(edge, &polyline.points);
+        this->process_edge_neighbors(edge, &polyline.points, &polyline.width, &polyline.endpoints);
 
         // get previous points
         {
             Points pp;
-            this->process_edge_neighbors(*edge.twin(), &pp);
+            std::vector<coordf_t> width;
+            std::vector<bool> endpoints;
+            this->process_edge_neighbors(*edge.twin(), &pp, &width, &endpoints);
             polyline.points.insert(polyline.points.begin(), pp.rbegin(), pp.rend());
+            polyline.width.insert(polyline.width.begin(), width.rbegin(), width.rend());
+            polyline.endpoints.insert(polyline.endpoints.begin(), endpoints.rbegin(), endpoints.rend());
+        }
+
+        assert(polyline.width.size() == polyline.points.size()*2 - 2);
+        assert(polyline.endpoints.size() == polyline.points.size());
+
+        if (polyline.first_point().coincides_with(polyline.last_point())) {
+            polyline.endpoints.front() = false;
+            polyline.endpoints.back() = false;
         }
 
         // append polyline to result
@@ -424,7 +428,16 @@ MedialAxis::build(Polylines* polylines)
 }
 
 void
-MedialAxis::process_edge_neighbors(const VD::edge_type& edge, Points* points)
+MedialAxis::build(Polylines* polylines)
+{
+    ThickPolylines tp;
+    this->build(&tp);
+    polylines->insert(polylines->end(), tp.begin(), tp.end());
+}
+
+void
+MedialAxis::process_edge_neighbors(const VD::edge_type& edge, Points* points,
+    std::vector<coordf_t>* width, std::vector<bool>* endpoints)
 {
     // Since rot_next() works on the edge starting point but we want
     // to find neighbors on the ending point, we just swap edge with
@@ -439,18 +452,26 @@ MedialAxis::process_edge_neighbors(const VD::edge_type& edge, Points* points)
 
     // if we have a single neighbor then we can continue recursively
     if (neighbors.size() == 1) {
+        endpoints->push_back(false);
         const VD::edge_type& neighbor = *neighbors.front();
         points->push_back(Point( neighbor.vertex1()->x(), neighbor.vertex1()->y() ));
+        width->push_back(this->thickness[&neighbor].first);
+        width->push_back(this->thickness[&neighbor].second);
         (void)this->edges.erase(&neighbor);
         (void)this->edges.erase(neighbor.twin());
-        this->process_edge_neighbors(neighbor, points);
+        this->process_edge_neighbors(neighbor, points, width, endpoints);
+    } else if (neighbors.size() == 0) {
+        endpoints->push_back(true);
+    } else {
+        // T-shaped or star-shaped joint
+        endpoints->push_back(false);
     }
 }
 
 bool
-MedialAxis::is_valid_edge(const VD::edge_type& edge) const
+MedialAxis::validate_edge(const VD::edge_type& edge)
 {
-    /* If the cells sharing this edge have a common vertex, we're not interested
+    /* If the cells sharing this edge have a common vertex, we're not (probably) interested
        in this edge. Why? Because it means that the edge lies on the bisector of
        two contiguous input lines and it was included in the Voronoi graph because
        it's the locus of centers of circles tangent to both vertices. Due to the 
@@ -481,31 +502,32 @@ MedialAxis::is_valid_edge(const VD::edge_type& edge) const
     // and we might need to skip the edge since it's not really part of
     // our skeleton
 
-    // get perpendicular distance of each edge vertex to the segment(s)
-    double dist0 = segment1.a.distance_to(segment2.b);
-    double dist1 = segment1.b.distance_to(segment2.a);
+    /*  Calculate perpendicular distance. We consider segment2 instead of segment1 
+        because our Voronoi edge is part of a CCW sequence going around its Voronoi cell
+        (segment).
+        This means that such segment is on the left on our edge, and goes backwards.
+        So we use the cell of the twin edge, which is located on the right of our edge
+        and goes in the same direction as it. This way we can map dist0 and dist1 
+        correctly.  */
+    const Line line(
+        Point( edge.vertex0()->x(), edge.vertex0()->y() ),
+        Point( edge.vertex1()->x(), edge.vertex1()->y() )
+    );
+    coordf_t dist0 = segment2.a.perp_distance_to(line)*2;
+    coordf_t dist1 = segment2.b.perp_distance_to(line)*2;
 
-    /*
-    Line line = this->edge_to_line(edge);
-    double diff = fabs(dist1 - dist0);
-    double dist_between_segments1 = segment1.a.distance_to(segment2);
-    double dist_between_segments2 = segment1.b.distance_to(segment2);
-    printf("w = %f/%f, dist0 = %f, dist1 = %f, diff = %f, seg1len = %f, seg2len = %f, edgelen = %f, s2s = %f / %f\n",
-        unscale(this->max_width), unscale(this->min_width),
-        unscale(dist0), unscale(dist1), unscale(diff),
-        unscale(segment1.length()), unscale(segment2.length()),
-        unscale(line.length()),
-        unscale(dist_between_segments1), unscale(dist_between_segments2)
-        );
-    */
-
     // if this edge is the centerline for a very thin area, we might want to skip it
     // in case the area is too thin
     if (dist0 < this->min_width && dist1 < this->min_width) {
         //printf(" => too thin, skipping\n");
         return false;
     }
 
+    if (this->expolygon != NULL && !this->expolygon->contains(line))
+        return false;
+
+    this->thickness[&edge] = std::make_pair(dist0, dist1);
+
     return true;
 }