wip18-09-11

11 files changed, 681 insertions, 144 deletions

diff --git a/lib/Slic3r/Print.pm b/lib/Slic3r/Print.pm
index a9763a824..d5e6b2a9d 100644
--- a/lib/Slic3r/Print.pm
+++ b/lib/Slic3r/Print.pm
@@ -34,6 +34,7 @@ sub size {
 # Slicing process, running at a background thread.
 sub process {
     my ($self) = @_;
+    print("perl process start\n");
     
     Slic3r::trace(3, "Staring the slicing process.");
     $_->make_perimeters for @{$self->objects};
@@ -61,7 +62,8 @@ sub process {
         eval "use Slic3r::Test::SectionCut";
         Slic3r::Test::SectionCut->new(print => $self)->export_svg("section_cut.svg");
     }
-    Slic3r::trace(3, "Slicing process finished.")
+    Slic3r::trace(3, "Slicing process finished.");
+    print("perl process end\n");
 }
 
 # G-code export process, running at a background thread.
@@ -213,6 +215,7 @@ EOF
 sub make_skirt {
     my $self = shift;
     
+    print("perl make_skirt start\n");
     # prerequisites
     $_->make_perimeters for @{$self->objects};
     $_->infill for @{$self->objects};
@@ -227,6 +230,7 @@ sub make_skirt {
         $self->_make_skirt();
     }
     $self->set_step_done(STEP_SKIRT);
+    print("perl make_skirt end\n");
 }
 
 sub make_brim {
diff --git a/lib/Slic3r/Print/Object.pm b/lib/Slic3r/Print/Object.pm
index e275fdde5..6194adc54 100644
--- a/lib/Slic3r/Print/Object.pm
+++ b/lib/Slic3r/Print/Object.pm
@@ -57,7 +57,7 @@ sub slice {
 # 3) Generates perimeters, gap fills and fill regions (fill regions of type stInternal).
 sub make_perimeters {
     my ($self) = @_;
-    
+    print("perl make_perimeters\n");
     # prerequisites
     $self->slice;
 
@@ -65,29 +65,40 @@ sub make_perimeters {
         $self->print->status_cb->(20, "Generating perimeters");
         $self->_make_perimeters;
     }
+    print("perl end make_perimeters\n");
 }
 
 sub prepare_infill {
     my ($self) = @_;
+    print("perl start prepare_infill\n");
     
     # prerequisites
     $self->make_perimeters;
     
-    return if $self->step_done(STEP_PREPARE_INFILL);
+    print("perl middle prepare_infill\n");
+	my $resul = $self->step_done(STEP_PREPARE_INFILL);
+    print("perl middle step_done?".$resul."\n");
+    return if $resul;
+    print("perl middleok prepare_infill\n");
     $self->set_step_started(STEP_PREPARE_INFILL);
     $self->print->status_cb->(30, "Preparing infill");
     
+    print("perl real prepare_infill\n");
     $self->_prepare_infill;
 
     $self->set_step_done(STEP_PREPARE_INFILL);
+    print("perl end prepare_infill\n");
 }
 
 sub infill {
     my ($self) = @_;
     
+    print("perl infill start\n");
     # prerequisites
     $self->prepare_infill;
+    print("perl infill middle\n");
     $self->_infill;
+    print("perl infill end\n");
 }
 
 sub generate_support_material {
diff --git a/xs/src/libslic3r/ExPolygon.cpp b/xs/src/libslic3r/ExPolygon.cpp
index 7c912b7c1..dc7cfe1c6 100644
--- a/xs/src/libslic3r/ExPolygon.cpp
+++ b/xs/src/libslic3r/ExPolygon.cpp
@@ -210,25 +210,44 @@ void ExPolygon::simplify(double tolerance, ExPolygons* expolygons) const
 
 /// remove point that are at SCALED_EPSILON * 2 distance.
 void
-remove_point_too_near(ThickPolyline* to_reduce)
-{
-        const coord_t smallest = SCALED_EPSILON * 2;
-        size_t id = 1;
-        while (id < to_reduce->points.size() - 1) {
-            size_t newdist = min(to_reduce->points[id].distance_to(to_reduce->points[id - 1])
-                , to_reduce->points[id].distance_to(to_reduce->points[id + 1]));
-            if (newdist < smallest) {
-                to_reduce->points.erase(to_reduce->points.begin() + id);
-                to_reduce->width.erase(to_reduce->width.begin() + id);
-                newdist = to_reduce->points[id].distance_to(to_reduce->points[id - 1]);
-            }
-            //go to next one
-            //if you removed a point, it check if the next one isn't too near from the previous one.
-            // if not, it byepass it.
-            if (newdist > smallest) {
-                ++id;
-            }
+remove_point_too_near(ThickPolyline& to_reduce, const coord_t smallest = SCALED_EPSILON * 2) {
+    size_t id = 1;
+    while (id < to_reduce.points.size() - 1) {
+        size_t newdist = min(to_reduce.points[id].distance_to(to_reduce.points[id - 1])
+            , to_reduce.points[id].distance_to(to_reduce.points[id + 1]));
+        if (newdist < smallest) {
+            to_reduce.points.erase(to_reduce.points.begin() + id);
+            to_reduce.width.erase(to_reduce.width.begin() + id);
+            to_reduce.curved.erase(to_reduce.curved.begin() + id);
+            newdist = to_reduce.points[id].distance_to(to_reduce.points[id - 1]);
         }
+        //go to next one
+        //if you removed a point, it check if the next one isn't too near from the previous one.
+        // if not, it byepass it.
+        if (newdist > smallest) {
+            ++id;
+        }
+    }
+}
+
+/// remove point that are at SCALED_EPSILON * 2 distance.
+void
+remove_point_too_near(ExPolygon& to_reduce, const coord_t smallest = SCALED_EPSILON * 2) {
+    size_t id = 1;
+    while (id < to_reduce.contour.points.size() - 1) {
+        size_t newdist = min(to_reduce.contour.points[id].distance_to(to_reduce.contour.points[id - 1])
+            , to_reduce.contour.points[id].distance_to(to_reduce.contour.points[id + 1]));
+        if (newdist < smallest) {
+            to_reduce.contour.points.erase(to_reduce.contour.points.begin() + id);
+            newdist = to_reduce.contour.points[id].distance_to(to_reduce.contour.points[id - 1]);
+        }
+        //go to next one
+        //if you removed a point, it check if the next one isn't too near from the previous one.
+        // if not, it byepass it.
+        if (newdist > smallest) {
+            ++id;
+        }
+    }
 }
 
 /// add points  from pattern to to_modify at the same % of the length
@@ -362,40 +381,12 @@ interpolate(double percent, Polyline& poly)
         }
     }
 }
-int id=0;
+
+int id = 0;
 void
-ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_width, ThickPolylines* polylines, double height) const
+fusion(const ExPolygon &bounds, double max_width, double min_width, ThickPolylines& pp, ExPolygon& polygon)
 {
-    id++;
-    std::cout << "\n============== me " << id << "\n";
-    // init helper object
-    Slic3r::Geometry::MedialAxis ma(max_width, min_width, this);
-    ma.lines = this->lines();
-    
-    // compute the Voronoi diagram and extract medial axis polylines
-    ThickPolylines pp;
-    ma.build(&pp);
-    
-
-    {
-        stringstream stri;
-        stri << "medial_axis" << id << ".svg";
-        SVG svg(stri.str());
-        svg.draw(bounds);
-        svg.draw(*this);
-        svg.draw(pp);
-        svg.Close();
-    }
-    
-    /* Find the maximum width returned; we're going to use this for validating and 
-       filtering the output segments. */
-    double max_w = 0;
-    for (ThickPolylines::const_iterator it = pp.begin(); it != pp.end(); ++it)
-        max_w = fmaxf(max_w, *std::max_element(it->width.begin(), it->width.end()));
 
-    concatThickPolylines(pp);
-    //reoder pp by length (ascending) It's really important to do that to avoid building the line from the width insteand of the length
-    std::sort(pp.begin(), pp.end(), [](const ThickPolyline & a, const ThickPolyline & b) { return a.length() < b.length(); });
 
     // Aligned fusion: Fusion the bits at the end of lines by "increasing thickness"
     // For that, we have to find other lines,
@@ -409,26 +400,27 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
         stri << "medial_axis_fus_" << id << "_" << id_f << "_" << i << ".svg";
         SVG svg(stri.str());
         svg.draw(bounds);
-        svg.draw(*this);
+        svg.draw(polygon);
         svg.draw(pp[i]);
         svg.Close();
     }
+    map<Point, double> coeff_angle_cache;
     while (changes) {
         changes = false;
         for (size_t i = 0; i < pp.size(); ++i) {
             ThickPolyline& polyline = pp[i];
-            std::cout << "next i : " << i << ", l=" << unscale(polyline.length()) << "<"<<EPSILON<<"\n";
-            
+            std::cout << "next i : " << i << ", l=" << unscale(polyline.length()) << "<" << EPSILON << "\n";
+
             //simple check to see if i can be fusionned
             if (!polyline.endpoints.first && !polyline.endpoints.second) continue;
-            
+
 
             ThickPolyline* best_candidate = nullptr;
             float best_dot = -1;
             size_t best_idx = 0;
             double dot_poly_branch = 0;
             double dot_candidate_branch = 0;
-            
+
             // find another polyline starting here
             for (size_t j = i + 1; j < pp.size(); ++j) {
                 std::cout << "next j : " << j << ", " << id << "\n";
@@ -444,15 +436,15 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 } else {
                     continue;
                 }
-                std::cout << " try : " << i << ":" << j << " : " << 
+                std::cout << " try : " << i << ":" << j << " : " <<
                     (polyline.points.size() < 2 && other.points.size() < 2) <<
                     (!polyline.endpoints.second || !other.endpoints.second) <<
                     ((polyline.points.back().distance_to(other.points.back())
                     + (polyline.width.back() + other.width.back()) / 4)
                     > max_width*1.05) <<
-                    (abs(polyline.length() - other.length()) > max_width / 2) << "\n";
+                    (abs(polyline.length() - other.length()) > max_width) << "\n";
                 std::cout << " polyline.length() : " << unscale(polyline.length()) << ", other.length()=" << unscale(other.length())
-                    << ",  diff=" << unscale(abs(polyline.length() - other.length())) << ", max= " << unscale(max_width / 2) << "\n";
+                    << ",  diff=" << unscale(abs(polyline.length() - other.length())) << ", max= " << unscale(max_width) << "\n";
                 std::cout << " polyline.dist() : " << unscale(polyline.points.back().distance_to(other.points.back())) << ", max=" << unscale(max_width*1.05) << "\n";
 
                 //// mergeable tests
@@ -462,15 +454,18 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 //note that this isn't the real calcul. It's just to avoid merging lines too far apart.
                 if (
                     ((polyline.points.back().distance_to(other.points.back())
-                        + (polyline.width.back() + other.width.back()) / 4)
-                        > max_width*1.05))
-                        continue;
+                    + (polyline.width.back() + other.width.back()) / 4)
+                > max_width*1.05))
+                continue;
                 // test if the lines are not too different in length.
                 if (abs(polyline.length() - other.length()) > max_width) continue;
 
 
                 //convex test
                 Point garbage_point;
+                std::cout << " other.length() > min_width: " << (other.length() > min_width)
+                    << ", garbage_point?" << bounds.contour.intersection(Line(interpolate(0.9, polyline), interpolate(0.9, other)), &garbage_point)
+                    << "\n";
                 if (other.length() > min_width &&
                     bounds.contour.intersection(Line(interpolate(0.9, polyline), interpolate(0.9, other)), &garbage_point)
                     ) continue;
@@ -523,8 +518,8 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     dot_poly_branch = 0.707;
                     dot_candidate_branch = 0.707;
                     std::cout << "no main branch... impossible!!\n";
-                } else if (!find_main_branch && 
-                    (pp[biggest_main_branch_id].length() < polyline.length() || pp[biggest_main_branch_id].length() < other.length()) ){
+                } else if (!find_main_branch &&
+                    (pp[biggest_main_branch_id].length() < polyline.length() || pp[biggest_main_branch_id].length() < other.length())) {
                     //the main branch should have no endpoint or be bigger!
                     //here, it have an endpoint, and is not the biggest -> bad!
                     std::cout << "no, the main is smaller than us\n";
@@ -536,7 +531,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                         stri << "medial_axis_fus_" << id << "_" << id_f << "_" << biggest_main_branch_id << "_main.svg";
                         SVG svg(stri.str());
                         svg.draw(bounds);
-                        svg.draw(*this);
+                        svg.draw(polygon);
                         svg.draw(pp[biggest_main_branch_id]);
                         svg.Close();
                     }
@@ -572,8 +567,8 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 std::cout << "nb points: < " << polyline.points.size() << ", " << best_candidate->points.size() << "\n";
 
                 // delete very near points
-                remove_point_too_near(&polyline);
-                remove_point_too_near(best_candidate);
+                remove_point_too_near(polyline);
+                remove_point_too_near(*best_candidate);
 
                 // add point at the same pos than the other line to have a nicer fusion
                 add_point_same_percent(&polyline, best_candidate);
@@ -583,7 +578,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     stri << "medial_axis_fus_" << id << "_" << id_f << "_" << i << "_if.svg";
                     SVG svg(stri.str());
                     svg.draw(bounds);
-                    svg.draw(*this);
+                    svg.draw(polygon);
                     svg.draw(polyline);
                     svg.Close();
                 }
@@ -592,7 +587,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     stri << "medial_axis_fus_" << id << "_" << id_f << "_" << best_idx << "_jf.svg";
                     SVG svg(stri.str());
                     svg.draw(bounds);
-                    svg.draw(*this);
+                    svg.draw(polygon);
                     svg.draw(*best_candidate);
                     svg.Close();
                 }
@@ -602,28 +597,46 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 //get the angle of the nearest points of the contour to see : _| (good) \_ (average) __(bad)
                 //sqrt because the result are nicer this way: don't over-penalize /_ angles
                 //TODO: try if we can achieve a better result if we use a different algo if the angle is <90°
-                const double coeff_angle_poly = (get_coeff_from_angle_countour(polyline.points.back(), *this, min(min_width, polyline.length() / 2)));
-                const double coeff_angle_candi = (get_coeff_from_angle_countour(best_candidate->points.back(), *this, min(min_width, best_candidate->length() / 2)));
+                double coeff_angle_poly = (coeff_angle_cache.find(polyline.points.back()) != coeff_angle_cache.end())
+                    ? coeff_angle_cache[polyline.points.back()]
+                    : (get_coeff_from_angle_countour(polyline.points.back(), polygon, min(min_width, polyline.length() / 2)));
+                double coeff_angle_candi = (coeff_angle_cache.find(best_candidate->points.back()) != coeff_angle_cache.end())
+                    ? coeff_angle_cache[best_candidate->points.back()]
+                    : (get_coeff_from_angle_countour(best_candidate->points.back(), polygon, min(min_width, best_candidate->length() / 2)));
 
                 //this will encourage to follow the curve, a little, because it's shorter near the center
                 //without that, it tends to go to the outter rim.
-                double weight_poly = 1 - polyline.length() / max(polyline.length(), best_candidate->length());
-                double weight_candi = 1 - best_candidate->length() / max(polyline.length(), best_candidate->length());
-                weight_poly *= 0;
-                weight_candi *= 0;
-                weight_poly += 1;
-                weight_candi += 1;
-                weight_poly = 0.5+coeff_angle_poly;
-                weight_candi = 0.5+coeff_angle_candi;
-                const double coeff_poly = (dot_poly_branch * weight_poly) / (dot_poly_branch * weight_poly + dot_candidate_branch * weight_candi);
+                double weight_poly = 0;
+                double weight_candi = 0;
+
+                std::cout
+                    << ", polylength=" << unscale(polyline.length())
+                    << ", weightDist_poly=" << (1 - polyline.length() / max(polyline.length(), best_candidate->length()))
+                    << ", dot_poly_branch=" << dot_poly_branch
+                    << ", coeff_angle_poly=" << coeff_angle_poly
+                    << "\n";
+                std::cout
+                    << ", candilengt=" << unscale(best_candidate->length())
+                    << ", weightDistcandi=" << (1 - best_candidate->length() / max(polyline.length(), best_candidate->length()))
+                    << ", dotcandi_branch=" << dot_candidate_branch
+                    << ", coeff_anglecandi=" << coeff_angle_candi
+                    << "\n";
+                weight_poly = .75;
+                weight_candi = .75;
+                weight_poly += (1 - polyline.length() / max(polyline.length(), best_candidate->length()));
+                weight_poly *= dot_poly_branch;
+                weight_poly *= coeff_angle_poly;
+                weight_candi += (1 - best_candidate->length() / max(polyline.length(), best_candidate->length()));
+                weight_candi *= dot_candidate_branch;
+                weight_candi *= coeff_angle_candi;
+                const double coeff_poly = weight_poly / (weight_poly + weight_candi); //(dot_poly_branch * weight_poly) / (dot_poly_branch * weight_poly + dot_candidate_branch * weight_candi);
                 const double coeff_candi = 1.0 - coeff_poly;
                 //iterate the points
-                // as voronoi should create symetric thing, we can iterate synchonously
+                // as voronoi should create symetric thing, we can iterate synchronously
                 size_t idx_point = 1;
                 while (idx_point < min(polyline.points.size(), best_candidate->points.size())) {
                     std::cout << " create point " << idx_point << " < " << polyline.points.size() << ", " << best_candidate->points.size() << "\n";
-                    std::cout << "dot_poly_branch=" << dot_poly_branch << ", weight_poly=" << weight_poly
-                        << ", dot_candidate_branch" << dot_candidate_branch << ", weight_candi=" << weight_candi << "\n";
+                    std::cout << "coeff_poly=" << coeff_poly << ", coeff_candi=" << coeff_candi << "\n";
                     //fusion
                     std::cout << "x= " << unscale(polyline.points[idx_point].x * coeff_poly + best_candidate->points[idx_point].x * coeff_candi)
                         << ", oldx1=" << unscale(polyline.points[idx_point].x)
@@ -635,7 +648,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                         << "\n";
                     polyline.points[idx_point].x = polyline.points[idx_point].x * coeff_poly + best_candidate->points[idx_point].x * coeff_candi;
                     polyline.points[idx_point].y = polyline.points[idx_point].y * coeff_poly + best_candidate->points[idx_point].y * coeff_candi;
-                   
+
                     // The width decrease with distance from the centerline.
                     // This formula is what works the best, even if it's not perfect (created empirically).  0->3% error on a gap fill on some tests.
                     //If someone find  an other formula based on the properties of the voronoi algorithm used here, and it works better, please use it.
@@ -680,12 +693,14 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 //remove points that are the same or too close each other, ie simplify
                 for (size_t idx_point = 1; idx_point < polyline.points.size(); ++idx_point) {
                     if (polyline.points[idx_point - 1].distance_to(polyline.points[idx_point]) < SCALED_EPSILON) {
-                        if (idx_point < polyline.points.size() -1) {
+                        if (idx_point < polyline.points.size() - 1) {
                             polyline.points.erase(polyline.points.begin() + idx_point);
                             polyline.width.erase(polyline.width.begin() + idx_point);
+                            polyline.curved.erase(polyline.curved.begin() + idx_point);
                         } else {
                             polyline.points.erase(polyline.points.begin() + idx_point - 1);
                             polyline.width.erase(polyline.width.begin() + idx_point - 1);
+                            polyline.curved.erase(polyline.curved.begin() + idx_point - 1);
                         }
                         --idx_point;
                     }
@@ -695,9 +710,14 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     if (!bounds.contains_b(polyline.points[idx_point])) {
                         polyline.points.erase(polyline.points.begin() + idx_point);
                         polyline.width.erase(polyline.width.begin() + idx_point);
+                        polyline.curved.erase(polyline.curved.begin() + idx_point);
                         --idx_point;
                     }
                 }
+
+                //update cache
+                coeff_angle_cache[polyline.points.back()] = coeff_angle_poly * coeff_poly + coeff_angle_candi * coeff_candi;
+
                 if (polyline.points.size() < 2) {
                     //remove self
                     pp.erase(pp.begin() + i);
@@ -717,7 +737,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 stri << "medial_axis_fus_" << id << "_" << id_f << ".svg";
                 SVG svg(stri.str());
                 svg.draw(bounds);
-                svg.draw(*this);
+                svg.draw(polygon);
                 svg.draw(pp);
                 svg.Close();
             }
@@ -734,25 +754,146 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 stri << "medial_axis_fus_" << id << "_" << id_f << "_" << i << ".svg";
                 SVG svg(stri.str());
                 svg.draw(bounds);
-                svg.draw(*this);
+                svg.draw(polygon);
                 svg.draw(pp[i]);
                 svg.Close();
             }
         }
     }
-    std::cout << "0";
+}
 
-    {
-        stringstream stri;
-        stri << "medial_axis_fusioned" << id << ".svg";
-        SVG svg(stri.str());
-        svg.draw(bounds);
-        svg.draw(*this);
-        svg.draw(pp);
-        svg.Close();
+
+void
+fusionv2(const ExPolygon &bounds, double max_width, double min_width, ThickPolylines& pp, ExPolygon& polygon) {
+
+    //fusion Y with only 1 '0' value => the "0" branch "pull" the cross-point
+
+    bool changes = false;
+    while (changes) {
+        changes = false;
+        for (size_t i = 0; i < pp.size(); ++i) {
+            ThickPolyline& polyline = pp[i];
+            // only consider polyline with no endpoints
+            if (polyline.endpoints.first || polyline.endpoints.second) continue;
+
+            //check both end
+            for (int j = 0; j < 2; j++) {
+                if (j == 1) polyline.reverse();
+
+                //find the others 
+
+
+            }
+
+            vector<size_t> crosspoint;
+            for (size_t j = 0; j < pp.size(); ++j) {
+                if (j == i) continue;
+                ThickPolyline& other = pp[j];
+                if (polyline.first_point().coincides_with(other.last_point())) {
+                    other.reverse();
+                    crosspoint.push_back(j);
+                } else if (polyline.first_point().coincides_with(other.first_point())) {
+                    crosspoint.push_back(j);
+                }
+            }
+            //check if it's a line that we can pull
+            if (crosspoint.size() != 2) continue;
+
+        }
+        if (changes) {
+            concatThickPolylines(pp);
+            ///reorder, in case of change
+            std::sort(pp.begin(), pp.end(), [](const ThickPolyline & a, const ThickPolyline & b) { return a.length() < b.length(); });
+        }
     }
+}
+
+
+void
+fusion_curve(ExPolygon& polygon, ThickPolylines &pp, double max_width, double min_width) {
+
     //fusion Y with only 1 '0' value => the "0" branch "pull" the cross-point
-    changes = false;
+    bool changes = false;
+    for (size_t i = 0; i < pp.size(); ++i) {
+        ThickPolyline& polyline = pp[i];
+        std::cout << "check curv : " << polyline.curved.front() << "-> " << polyline.curved.back()
+            << ", size: " << polyline.points.size()
+            << ", length: " << unscale(polyline.length())
+            << ", endf: " << polyline.endpoints.first
+            << ", endb: " << polyline.endpoints.second
+            << ", widthf: " << polyline.width.front()
+            << ", widthb: " << polyline.width.back()
+            << "\n";
+        // only consider 2-point polyline with endpoint
+        if (polyline.points.size() != 2) continue;
+        if (polyline.endpoints.first) polyline.reverse();
+        else if (!polyline.endpoints.second) continue;
+        if (polyline.width.back() > 0) continue;
+
+        //check my length is small
+        coord_t length = (coord_t)polyline.length();
+        if (length > max_width) continue;
+
+        std::cout << "polyline.curved.front() : " << polyline.curved.front() << ", get_coeff_from_angle_countour= " << get_coeff_from_angle_countour(polyline.points.back(), polygon, SCALED_RESOLUTION) << "\n";
+        //only consider "curved" segment or very shallow angle for contour
+        if (!polyline.curved.front() || 
+            get_coeff_from_angle_countour(polyline.points.back(), polygon, SCALED_RESOLUTION) < 0.2) continue;
+        std::cout << "CURV : " << polyline.curved.front() << "-> " << polyline.curved.back() << "\n";
+
+
+        // look if other end is a cross point with multiple other branch
+        vector<size_t> crosspoint;
+        for (size_t j = 0; j < pp.size(); ++j) {
+            if (j == i) continue;
+            ThickPolyline& other = pp[j];
+            if (polyline.first_point().coincides_with(other.last_point())) {
+                other.reverse();
+                crosspoint.push_back(j);
+            } else if (polyline.first_point().coincides_with(other.first_point())) {
+                crosspoint.push_back(j);
+            }
+        }
+        //check if it's a line that we can pull
+        if (crosspoint.size() != 2) continue;
+
+        // if true, pull it a bit, depends on my size, the dot?, and the coeff at my 0-end (~12% for a square, almost 0 for a gentle curve)
+        if (!polyline.curved.front()) {
+            coord_t length_pull = polyline.length();
+            length_pull *= 0.15 * get_coeff_from_angle_countour(polyline.points.back(), polygon, min(min_width, polyline.length() / 2));
+
+            //compute dir
+            Vectorf pull_direction(polyline.points[1].x - polyline.points[0].x, polyline.points[1].y - polyline.points[0].y);
+            pull_direction = normalize(pull_direction);
+            pull_direction.x *= length_pull;
+            pull_direction.y *= length_pull;
+
+            //pull the points
+            Point &p1 = pp[crosspoint[0]].points[0];
+            p1.x = p1.x + (coord_t)pull_direction.x;
+            p1.y = p1.y + (coord_t)pull_direction.y;
+
+            Point &p2 = pp[crosspoint[1]].points[0];
+            p2.x = p2.x + (coord_t)pull_direction.x;
+            p2.y = p2.y + (coord_t)pull_direction.y;
+        }
+
+        //delete the now unused polyline
+        pp.erase(pp.begin() + i);
+        --i;
+        changes = true;
+    }
+    if (changes) {
+        concatThickPolylines(pp);
+        ///reorder, in case of change
+        std::sort(pp.begin(), pp.end(), [](const ThickPolyline & a, const ThickPolyline & b) { return a.length() < b.length(); });
+    }
+}
+
+void
+fusion_corners(ExPolygon& polygon, ThickPolylines &pp, double max_width, double min_width) {
+
+    //fusion Y with only 1 '0' value => the "0" branch "pull" the cross-point
+    bool changes = false;
     for (size_t i = 0; i < pp.size(); ++i) {
         ThickPolyline& polyline = pp[i];
         // only consider polyline with 0-end
@@ -761,6 +902,10 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
         else if (!polyline.endpoints.second) continue;
         if (polyline.width.back() > 0) continue;
 
+        //check my length is small
+        coord_t length = polyline.length();
+        if (length > max_width) continue;
+
         // look if other end is a cross point with multiple other branch
         vector<size_t> crosspoint;
         for (size_t j = 0; j < pp.size(); ++j) {
@@ -773,24 +918,30 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                 crosspoint.push_back(j);
             }
         }
+        //check if it's a line that we can pull
         if (crosspoint.size() != 2) continue;
-        // if true, see if the dot from me and the other is <0
-        Line poly_direction(polyline.points[0], polyline.points[1]);
-        double dot1 = dot(poly_direction, Line(pp[crosspoint[0]].points[0], pp[crosspoint[0]].points[1]));
-        if (dot1 == 0) continue;
-        double dot2 = dot(poly_direction, Line(pp[crosspoint[1]].points[0], pp[crosspoint[1]].points[1]));
-        if (dot2 == 0) continue;
+
+        // check if i am at the external side of a curve
+        double angle1 = polyline.points[0].ccw_angle(polyline.points[1], pp[crosspoint[0]].points[1]);
+        if (angle1 >= PI) angle1 = 2 * PI - angle1;  // smaller angle
+        double angle2 = polyline.points[0].ccw_angle(polyline.points[1], pp[crosspoint[1]].points[1]);
+        if (angle2 >= PI) angle2 = 2 * PI - angle2;  // smaller angle
+        if (angle1 + angle2 < PI) continue;
+
+        //check if is smaller or the other ones are not endpoits
+        if (pp[crosspoint[0]].endpoints.second && length > pp[crosspoint[0]].length()) continue;
+        if (pp[crosspoint[1]].endpoints.second && length > pp[crosspoint[1]].length()) continue;
 
         // if true, pull it a bit, depends on my size, the dot?, and the coeff at my 0-end (~12% for a square, almost 0 for a gentle curve)
         coord_t length_pull = polyline.length();
-        length_pull *= 0.12 * get_coeff_from_angle_countour(polyline.points.back(), *this, min(min_width, polyline.length() / 2));
-        
+        length_pull *= 0.15 * get_coeff_from_angle_countour(polyline.points.back(), polygon, min(min_width, polyline.length() / 2));
+
         //compute dir
-        Vectorf pull_direction(poly_direction.b.x - poly_direction.a.x, poly_direction.b.y - poly_direction.a.y);
+        Vectorf pull_direction(polyline.points[1].x - polyline.points[0].x, polyline.points[1].y - polyline.points[0].y);
         pull_direction = normalize(pull_direction);
         pull_direction.x *= length_pull;
         pull_direction.y *= length_pull;
-        
+
         //pull the points
         Point &p1 = pp[crosspoint[0]].points[0];
         p1.x = p1.x + pull_direction.x;
@@ -799,7 +950,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
         Point &p2 = pp[crosspoint[1]].points[0];
         p2.x = p2.x + pull_direction.x;
         p2.y = p2.y + pull_direction.y;
-        
+
         //delete the now unused polyline
         pp.erase(pp.begin() + i);
         --i;
@@ -810,9 +961,248 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
         ///reorder, in case of change
         std::sort(pp.begin(), pp.end(), [](const ThickPolyline & a, const ThickPolyline & b) { return a.length() < b.length(); });
     }
+}
+
+void
+extends_line(ThickPolyline& polyline, const ExPolygon& inner_bound, const ExPolygon& outter_bounds, const ExPolygons& anchors, const coord_t max_width, const coord_t join_width) {
+
+    // 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
+    std::cout << " ok, polyline.endpoints.second:" << polyline.endpoints.second << ", bounds?" << outter_bounds.has_boundary_point(polyline.points.back()) << "\n";
+    if (polyline.endpoints.second && !outter_bounds.has_boundary_point(polyline.points.back())) {
+        Line line(*(polyline.points.end() - 2), polyline.points.back());
+        std::cout << " ok, line:" << unscale(line.a.x)<<"->"<<unscale(line.b.x) << " (x)\n";
+
+        // 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);
+        Point new_back;
+        if (inner_bound.contour.has_boundary_point(polyline.points.back())) {
+            std::cout << "inner boundary ok\n";
+            new_back = polyline.points.back();
+        } else {
+            Line width;
+            (void)inner_bound.contour.first_intersection(line, &new_back, &width);
+            polyline.points.push_back(new_back);
+            polyline.width.push_back(width.length());
+            polyline.curved.push_back(false);
+            std::cout << "extends inner boundary to" << unscale(new_back.x)<<"\n";
+        }
+        Point new_bound;
+        (void)outter_bounds.contour.first_intersection(line, &new_bound);
+        if (new_bound.coincides_with_epsilon(new_back)) {
+            std::cout << "outter boundary ok\n";
+            return;
+        }
+        //find anchor
+        Point best_anchor;
+        double shortest_dist = -1;
+        for (const ExPolygon& a : anchors) {
+            std::cout << "try anchor\n";
+            Point p_maybe_inside;
+            for (Point p : a.contour.points) {
+                std::cout << "add point " << unscale(p.x) << "\n";
+                p_maybe_inside.x += p.x;
+                p_maybe_inside.y += p.y;
+            }
+            p_maybe_inside.x /= a.contour.points.size();
+            p_maybe_inside.y /= a.contour.points.size();
+            std::cout << "total point " << unscale(p_maybe_inside.x) << "\n";
+            double test_dist = new_bound.distance_to(p_maybe_inside) + new_back.distance_to(p_maybe_inside);
+            std::cout << "test_dist= " << unscale(test_dist) << "\n";
+            if (test_dist < max_width / 2 && (test_dist < shortest_dist || shortest_dist < 0)) {
+                std::cout << " best one! " << "\n";
+                shortest_dist = test_dist;
+                best_anchor = p_maybe_inside;
+            }
+        }
+        if (new_bound.distance_to(new_back) < best_anchor.distance_to(new_back)) {
+            std::cout << "extends outter boundary to straith" << unscale(new_bound.x) << "\n";
+            polyline.points.push_back(new_bound);
+            polyline.width.push_back(join_width);
+            polyline.curved.push_back(false);
+        }
+        if (shortest_dist > 0) {
+            std::cout << "extends outter boundary to middle point" << unscale(best_anchor.x) << "\n";
+            //line.a = new_back;
+            //line.b = best_anchor;
+            //line.extend_end(max_width);
+            //(void)outter_bounds.contour.first_intersection(line, &new_bound);
+            polyline.points.push_back(best_anchor);
+            polyline.width.push_back(join_width);
+            polyline.curved.push_back(false);
+        }
+    }
+}
+
+void
+ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_width, ThickPolylines* polylines, double height) const
+{
+    id++;
+    std::cout << "\n============== me " << id << "\n";
+
+    //simplify the boundary between us and the bounds.
+    //int firstIdx = 0;
+    //while (firstIdx < contour.points.size() && bounds.contour.contains(contour.points[firstIdx])) firstIdx++;
+    ExPolygon simplified_poly = *this;
+    if (this != &bounds) {
+        bool need_intersect = false;
+        for (size_t i = 0; i < simplified_poly.contour.points.size(); i++) {
+            Point &p_check = simplified_poly.contour.points[i];
+            std::cout << "check point " << p_check.x << " : " << p_check.y << "\n";
+            //bool find = false;
+            //for (size_t j = 0; j<bounds.contour.points.size(); j++) {
+            //    const Point &p_iter = bounds.contour.points[j];
+            //    if (abs(p_iter.x - p_check.x)< SCALED_EPSILON
+            //        && abs(p_iter.y - p_check.y)< SCALED_EPSILON) {
+            //        std::cout << "find point " << "\n";
+            //        find = true;
+            //        break;
+            //    }
+            //    size_t next = (j + 1) % simplified_poly.contour.points.size();
+            //    double angle = p_check.ccw_angle(p_iter, bounds.contour.points[next]);
+            //    if (angle > PI) angle = PI * 2 - angle;
+            //    if (angle > PI*0.99) {
+            //        std::cout << "find point (almost inside a line) " << angle / PI << "\n";
+            //        find = true;
+            //        break;
+            //    }
+            //}
+            //if (!find) {
+            if (!bounds.has_boundary_point(p_check)) {
+                //check if we put it at a bound point instead of delete it
+                size_t prev_i = i == 0 ? simplified_poly.contour.points.size() - 1 : (i - 1);
+                size_t next_i = i == simplified_poly.contour.points.size() - 1 ? 0 : (i + 1);
+                const Point* closest = bounds.contour.closest_point(p_check);
+                if (closest != nullptr && closest->distance_to(p_check) + SCALED_EPSILON
+                    < min(p_check.distance_to(simplified_poly.contour.points[prev_i]), p_check.distance_to(simplified_poly.contour.points[next_i])) / 2) {
+                    p_check.x = closest->x;
+                    p_check.y = closest->y;
+                    std::cout << " move to " << closest->x << " : " << closest->x << "\n";
+                    need_intersect = true;
+                } else {
+                    std::cout << " erase " << "\n";
+                    simplified_poly.contour.points.erase(simplified_poly.contour.points.begin() + i);
+                    i--;
+                }
+            }
+        }
+        if (need_intersect) {
+            ExPolygons simplified_polys = intersection_ex(simplified_poly, bounds);
+            if (simplified_polys.size() == 1) {
+                simplified_poly = simplified_polys[0];
+            } else {
+                simplified_poly = *this;
+            }
+        }
+    }
+    std::cout << "SCALED_RESOLUTION=" << unscale(SCALED_RESOLUTION) << ", /10=" << unscale(SCALED_RESOLUTION / 10)
+        << ", SCALED_EPSILON=" << unscale(SCALED_EPSILON) << ", x5=" << unscale(SCALED_EPSILON * 10)
+        << "\n";
+    remove_point_too_near(simplified_poly, min(SCALED_RESOLUTION / 20, SCALED_EPSILON * 3));
+
+
+    // init helper object
+    Slic3r::Geometry::MedialAxis ma(max_width, min_width, &simplified_poly);
+    ma.lines = simplified_poly.lines();
+    
+    // compute the Voronoi diagram and extract medial axis polylines
+    ThickPolylines pp;
+    ma.build(&pp);
+
+    for (ThickPolyline& pc : pp) {
+        remove_point_too_near(pc, min(SCALED_RESOLUTION/10, SCALED_EPSILON*5));
+    }
+
+    {
+        stringstream stri;
+        stri << "medial_axis_0begin_" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(*this);
+        svg.draw(pp);
+        svg.Close();
+    }
+
+    {
+        stringstream stri;
+        stri << "medial_axis_0simpl" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(pp);
+        svg.Close();
+    }
+    {
+        ThickPolylines ppCurved = pp;
+        for (ThickPolyline& pc : ppCurved) {
+            for (int i = 0; i < pc.points.size(); i++) {
+                if (pc.curved[i]) {
+                    pc.points.erase(pc.points.begin() + i);
+                    pc.width.erase(pc.width.begin() + i);
+                    pc.curved.erase(pc.curved.begin() + i);
+                    i--;
+                }
+            }
+        }
+        stringstream stri;
+        stri << "medial_axis_0straith" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(ppCurved);
+        svg.Close();
+    }
+    
+    /* Find the maximum width returned; we're going to use this for validating and 
+       filtering the output segments. */
+    double max_width_voronoi = 0;
+    for (ThickPolylines::const_iterator it = pp.begin(); it != pp.end(); ++it)
+        max_width_voronoi = fmaxf(max_width_voronoi, *std::max_element(it->width.begin(), it->width.end()));
+
+    concatThickPolylines(pp);
+    //reoder pp by length (ascending) It's really important to do that to avoid building the line from the width insteand of the length
+    std::sort(pp.begin(), pp.end(), [](const ThickPolyline & a, const ThickPolyline & b) { return a.length() < b.length(); });
+
+
+    fusion_curve(simplified_poly, pp, max_width, min_width);
+    {
+        stringstream stri;
+        stri << "medial_axis_2_uncurved" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(pp);
+        svg.Close();
+    }
+    fusion(bounds, max_width, min_width, pp, simplified_poly);
+    std::cout << "0\n";
+
+    {
+        stringstream stri;
+        stri << "medial_axis_3fusioned" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(pp);
+        svg.Close();
+    }
+    fusion_corners(simplified_poly, pp, max_width, min_width);
+    {
+        stringstream stri;
+        stri << "medial_axis_4pulled" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(pp);
+        svg.Close();
+    }
 
     // remove too small extrusion at start & end of polylines
-    changes = false;
+    bool changes = false;
     for (size_t i = 0; i < pp.size(); ++i) {
         ThickPolyline& polyline = pp[i];
         // remove bits with too small extrusion
@@ -834,6 +1224,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
             }
             polyline.points.erase(polyline.points.begin());
             polyline.width.erase(polyline.width.begin());
+            polyline.curved.erase(polyline.curved.begin());
             changes = true;
         }
         while (polyline.points.size() > 1 && polyline.width.back() < min_width && polyline.endpoints.second) {
@@ -852,8 +1243,9 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     break;
                 }
             }
-            polyline.points.erase(polyline.points.end()-1);
+            polyline.points.erase(polyline.points.end() - 1);
             polyline.width.erase(polyline.width.end() - 1);
+            polyline.curved.erase(polyline.curved.end() - 1);
             changes = true;
         }
         if (polyline.points.size() < 2) {
@@ -864,42 +1256,143 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
     }
     if (changes) concatThickPolylines(pp);
 
+    {
+        stringstream stri;
+        stri << "medial_axis_5shortened" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(pp);
+        svg.Close();
+    }
+
+    ExPolygons anchors = diff_ex(bounds, *this);
+    {
+        stringstream stri;
+        stri << "medial_axis_5.5ANCHOR" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(anchors);
+        svg.draw(pp);
+        svg.Close();
+    }
     // Loop through all returned polylines in order to extend their endpoints to the 
     //   expolygon boundaries
     for (size_t i = 0; i < pp.size(); ++i) {
         ThickPolyline& polyline = pp[i];
+        std::cout << " extends poly " << i << " bounds\n";
+        extends_line(polyline, simplified_poly, bounds, anchors, max_width, min_width);
+        polyline.reverse();
+        std::cout << " extends poly " << i << " bounds (reverse)\n";
+        extends_line(polyline, simplified_poly, bounds, anchors, max_width, min_width);
 
         // 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.first && !bounds.has_boundary_point(new_front)) {
-            Line line(polyline.points[1], polyline.points.front());
-
-            // 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)bounds.contour.first_intersection(line, &new_front);
-        }
-        if (polyline.endpoints.second && !bounds.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)bounds.contour.first_intersection(line, &new_back);
-        }
-        polyline.points.front() = new_front;
-        polyline.points.back() = new_back;
+        //if (polyline.endpoints.first && !bounds.has_boundary_point(polyline.points.front())) {
+        //    Line line(polyline.points[1], polyline.points.front());
+
+        //    // 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);
+        //    Point new_front;
+        //    if (this->contour.has_boundary_point(polyline.points.front())) {
+        //        new_front = polyline.points.front();
+        //    } else {
+        //        Line width;
+        //        (void)this->contour.first_intersection(line, &new_front, &width);
+        //        polyline.points.insert(polyline.points.begin(), new_front);
+        //        polyline.width.insert(polyline.width.begin(), width.length());
+        //    }
+        //    Point new_bound;
+        //    (void)bounds.contour.first_intersection(line, &new_bound);
+        //    //find anchor
+        //    Point best_anchor;
+        //    double shortest_dist = -1;
+        //    for (ExPolygon& a : anchors) {
+        //        Point p_maybe_inside;
+        //        for (Point p : a.contour.points) {
+        //            p_maybe_inside.x += p.x;
+        //            p_maybe_inside.y += p.y;
+        //        }
+        //        p_maybe_inside.x /= a.contour.points.size();
+        //        p_maybe_inside.y /= a.contour.points.size();
+        //        double test_dist = new_bound.distance_to(p_maybe_inside) + new_front.distance_to(p_maybe_inside);
+        //        if(test_dist < max_width / 2 && (test_dist < shortest_dist || shortest_dist < 0)) {
+        //            shortest_dist = test_dist;
+        //            best_anchor = p_maybe_inside;
+        //        }
+        //    }
+        //    if (shortest_dist > 0) {
+        //        line.a = new_front;
+        //        line.b = best_anchor;
+        //        line.extend_end(max_width);
+        //        (void)bounds.contour.first_intersection(line, &new_bound);
+        //        polyline.points.insert(polyline.points.begin(), best_anchor);
+        //        polyline.width.insert(polyline.width.begin(), min_width);
+        //    }
+        //}
+        //if (polyline.endpoints.second && !bounds.has_boundary_point(polyline.points.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);
+
+        //    Point new_back;
+        //    if (this->contour.has_boundary_point(polyline.points.back())) {
+        //        new_back = polyline.points.back();
+        //    } else {
+        //        Line width;
+        //        (void)this->contour.first_intersection(line, &new_back, &width);
+        //        polyline.points.push_back(new_back);
+        //        polyline.width.push_back(width.length());
+        //    }
+        //    Point new_bound;
+        //    (void)bounds.contour.first_intersection(line, &new_bound);
+        //    //find anchor
+        //    Point best_anchor;
+        //    double shortest_dist = -1;
+        //    for (ExPolygon& a : anchors) {
+        //        Point p_maybe_inside;
+        //        for (Point p : a.contour.points) {
+        //            p_maybe_inside.x += p.x;
+        //            p_maybe_inside.y += p.y;
+        //        }
+        //        p_maybe_inside.x /= a.contour.points.size();
+        //        p_maybe_inside.y /= a.contour.points.size();
+        //        double test_dist = new_bound.distance_to(p_maybe_inside) + new_back.distance_to(p_maybe_inside);
+        //        if (test_dist < max_width / 2 && (test_dist < shortest_dist || shortest_dist < 0)) {
+        //            shortest_dist = test_dist;
+        //            best_anchor = p_maybe_inside;
+        //        }
+        //    }
+        //    if (shortest_dist > 0) {
+        //        line.a = new_back;
+        //        line.b = best_anchor;
+        //        line.extend_end(max_width);
+        //        (void)bounds.contour.first_intersection(line, &new_bound);
+        //        polyline.points.push_back(best_anchor);
+        //        polyline.width.push_back(min_width);
+        //    }
+        //}
 
     }
+    {
+        stringstream stri;
+        stri << "medial_axis_6extended" << id << ".svg";
+        SVG svg(stri.str());
+        svg.draw(bounds);
+        svg.draw(simplified_poly);
+        svg.draw(pp);
+        svg.Close();
+    }
 
 
     // concatenate, but even where multiple thickpolyline join, to create nice long strait polylines
@@ -971,24 +1464,30 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     //too thin at start
                     polyline.points.erase(polyline.points.begin());
                     polyline.width.erase(polyline.width.begin());
+                    polyline.curved.erase(polyline.curved.begin());
                     idx_point = 0;
                 } else if (idx_point == 1) {
                     //too thin at start
                     polyline.points.erase(polyline.points.begin());
                     polyline.width.erase(polyline.width.begin());
+                    polyline.curved.erase(polyline.curved.begin());
                     polyline.points.erase(polyline.points.begin());
                     polyline.width.erase(polyline.width.begin());
+                    polyline.curved.erase(polyline.curved.begin());
                     idx_point = 0;
                 } else if (idx_point == polyline.points.size() - 2) {
                     //too thin at (near) end
                     polyline.points.erase(polyline.points.end() - 1);
                     polyline.width.erase(polyline.width.end() - 1);
+                    polyline.curved.erase(polyline.curved.end() - 1);
                     polyline.points.erase(polyline.points.end() - 1);
                     polyline.width.erase(polyline.width.end() - 1);
+                    polyline.curved.erase(polyline.curved.end() - 1);
                 } else if (idx_point == polyline.points.size() - 1) {
                     //too thin at end
                     polyline.points.erase(polyline.points.end() - 1);
                     polyline.width.erase(polyline.width.end() - 1);
+                    polyline.curved.erase(polyline.curved.end() - 1);
                 } else {
                     //too thin in middle : split
                     pp.emplace_back();
@@ -997,6 +1496,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
                     newone.width.insert(newone.width.begin(), polyline.width.begin() + idx_point + 1, polyline.width.end());
                     polyline.points.erase(polyline.points.begin() + idx_point, polyline.points.end());
                     polyline.width.erase(polyline.width.begin() + idx_point, polyline.width.end());
+                    polyline.curved.erase(polyline.curved.begin() + idx_point, polyline.curved.end());
                 }
             } else idx_point++;
 
@@ -1014,7 +1514,7 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
     while (changes) {
         changes = false;
         
-        double shortest_size = max_w * 2;
+        double shortest_size = max_width_voronoi * 2;
         size_t shortest_idx = -1;
         for (size_t i = 0; i < pp.size(); ++i) {
             ThickPolyline& polyline = pp[i];
@@ -1076,10 +1576,11 @@ ExPolygon::medial_axis(const ExPolygon &bounds, double max_width, double min_wid
     }
     {
         stringstream stri;
-        stri << "medial_axis_end" << id << ".svg";
+        stri << "medial_axis_9end" << id << ".svg";
         SVG svg(stri.str());
         svg.draw(bounds);
         svg.draw(*this);
+        svg.draw(simplified_poly);
         svg.draw(pp);
         svg.Close();
     }
diff --git a/xs/src/libslic3r/Geometry.cpp b/xs/src/libslic3r/Geometry.cpp
index ed18f66db..c4b892599 100644
--- a/xs/src/libslic3r/Geometry.cpp
+++ b/xs/src/libslic3r/Geometry.cpp
@@ -905,6 +905,8 @@ MedialAxis::build(ThickPolylines* polylines)
         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);
+        polyline.curved.push_back(edge->is_curved());
+        polyline.curved.push_back(edge->is_curved());
         
         // remove this edge and its twin from the available edges
         (void)this->edges.erase(edge);
@@ -912,24 +914,25 @@ MedialAxis::build(ThickPolylines* polylines)
         
         // get next points
         this->process_edge_neighbors(edge, &polyline);
-        
+
         // get previous points
         {
             ThickPolyline rpolyline;
             this->process_edge_neighbors(edge->twin(), &rpolyline);
             polyline.points.insert(polyline.points.begin(), rpolyline.points.rbegin(), rpolyline.points.rend());
             polyline.width.insert(polyline.width.begin(), rpolyline.width.rbegin(), rpolyline.width.rend());
+            polyline.curved.insert(polyline.curved.begin(), rpolyline.curved.rbegin(), rpolyline.curved.rend());
             polyline.endpoints.first = rpolyline.endpoints.second;
         }
         
         assert(polyline.width.size() == polyline.points.size());
-        
+
         // prevent loop endpoints from being extended
         if (polyline.first_point().coincides_with(polyline.last_point())) {
             polyline.endpoints.first = false;
             polyline.endpoints.second = false;
         }
-        
+
         // append polyline to result
         polylines->push_back(polyline);
     }
@@ -984,6 +987,7 @@ MedialAxis::process_edge_neighbors(const VD::edge_type* edge, ThickPolyline* pol
             Point new_point(neighbor->vertex1()->x(), neighbor->vertex1()->y());
             polyline->points.push_back(new_point);
             polyline->width.push_back(this->thickness[neighbor].second);
+            polyline->curved.push_back(neighbor->is_curved());
             
             (void)this->edges.erase(neighbor);
             (void)this->edges.erase(neighbor->twin());
diff --git a/xs/src/libslic3r/Layer.cpp b/xs/src/libslic3r/Layer.cpp
index 652bcdaa0..a2582a9c6 100644
--- a/xs/src/libslic3r/Layer.cpp
+++ b/xs/src/libslic3r/Layer.cpp
@@ -75,6 +75,7 @@ void Layer::merge_slices()
 // The resulting fill surface is split back among the originating regions.
 void Layer::make_perimeters()
 {
+    std::cout << "make_perimeters\n";
     BOOST_LOG_TRIVIAL(trace) << "Generating perimeters for layer " << this->id();
     
     // keep track of regions whose perimeters we have already generated
@@ -84,6 +85,7 @@ void Layer::make_perimeters()
         size_t region_id = layerm - this->regions.begin();
         if (done.find(region_id) != done.end()) continue;
         BOOST_LOG_TRIVIAL(trace) << "Generating perimeters for layer " << this->id() << ", region " << region_id;
+        std::cout << "Generating perimeters for layer " << this->id() << ", region " << region_id<<"\n";
         done.insert(region_id);
         const PrintRegionConfig &config = (*layerm)->region()->config;
         
@@ -111,7 +113,9 @@ void Layer::make_perimeters()
         if (layerms.size() == 1) {  // optimization
             (*layerm)->fill_surfaces.surfaces.clear();
             (*layerm)->make_perimeters((*layerm)->slices, &(*layerm)->fill_surfaces);
+            std::cout << "fill_surfaces call\n";
             (*layerm)->fill_expolygons = to_expolygons((*layerm)->fill_surfaces.surfaces);
+            std::cout << "fill_surfaces call end\n";
         } else {
             SurfaceCollection new_slices;
             {
@@ -143,6 +147,7 @@ void Layer::make_perimeters()
         }
     }
     BOOST_LOG_TRIVIAL(trace) << "Generating perimeters for layer " << this->id() << " - Done";
+    std::cout << "end make_perimeters\n";
 }
 
 void Layer::make_fills()
diff --git a/xs/src/libslic3r/LayerRegion.cpp b/xs/src/libslic3r/LayerRegion.cpp
index a409f3885..a621caa33 100644
--- a/xs/src/libslic3r/LayerRegion.cpp
+++ b/xs/src/libslic3r/LayerRegion.cpp
@@ -54,6 +54,7 @@ void LayerRegion::slices_to_fill_surfaces_clipped()
 void
 LayerRegion::make_perimeters(const SurfaceCollection &slices, SurfaceCollection* fill_surfaces)
 {
+    std::cout << "make peri\n";
     this->perimeters.clear();
     this->thin_fills.clear();
     
@@ -86,6 +87,7 @@ LayerRegion::make_perimeters(const SurfaceCollection &slices, SurfaceCollection*
     g.process();
 
     this->fill_no_overlap_expolygons = g.fill_no_overlap;
+    std::cout << "end make peri\n";
 }
 
 //#define EXTERNAL_SURFACES_OFFSET_PARAMETERS ClipperLib::jtMiter, 3.
diff --git a/xs/src/libslic3r/MultiPoint.cpp b/xs/src/libslic3r/MultiPoint.cpp
index 2e65492cd..216c3e791 100644
--- a/xs/src/libslic3r/MultiPoint.cpp
+++ b/xs/src/libslic3r/MultiPoint.cpp
@@ -137,7 +137,7 @@ MultiPoint::intersection(const Line& line, Point* intersection) const
     return false;
 }
 
-bool MultiPoint::first_intersection(const Line& line, Point* intersection) const
+bool MultiPoint::first_intersection(const Line& line, Point* intersection, Line* line_out) const
 {
     bool   found = false;
     double dmin  = 0.;
@@ -148,11 +148,13 @@ bool MultiPoint::first_intersection(const Line& line, Point* intersection) const
                 found = true;
                 dmin = ip.distance_to(line.a);
                 *intersection = ip;
+                if (line_out != nullptr) *line_out = l;
             } else {
                 double d = ip.distance_to(line.a);
                 if (d < dmin) {
                     dmin = d;
                     *intersection = ip;
+                    if (line_out != nullptr) *line_out = l;
                 }
             }
         }
diff --git a/xs/src/libslic3r/MultiPoint.hpp b/xs/src/libslic3r/MultiPoint.hpp
index 0970e9a67..46fe6504a 100644
--- a/xs/src/libslic3r/MultiPoint.hpp
+++ b/xs/src/libslic3r/MultiPoint.hpp
@@ -74,7 +74,7 @@ public:
     }
 
     bool intersection(const Line& line, Point* intersection) const;
-    bool first_intersection(const Line& line, Point* intersection) const;
+    bool first_intersection(const Line& line, Point* intersection, Line* line_out = nullptr) const;
     std::string dump_perl() const;
     
     static Points _douglas_peucker(const Points &points, const double tolerance);
diff --git a/xs/src/libslic3r/PerimeterGenerator.cpp b/xs/src/libslic3r/PerimeterGenerator.cpp
index 0c2b0655a..8ec897352 100644
--- a/xs/src/libslic3r/PerimeterGenerator.cpp
+++ b/xs/src/libslic3r/PerimeterGenerator.cpp
@@ -251,6 +251,7 @@ void PerimeterGenerator::process()
                             no_thin_zone = diff_ex(last, offset_ex(expp, (float)(min_width / 2)), true);
                         }
                         // compute a bit of overlap to anchor thin walls inside the print.
+                        std::cout << "thin walls begin\n";
                         for (ExPolygon &ex : expp) {
                             //growing back the polygon
                             //a vary little bit of overlap can be created here with other thin polygon, but it's more useful than worisome.
@@ -264,13 +265,14 @@ void PerimeterGenerator::process()
                                     //be sure it's not too small to extrude reliably
                                     if (ex_bigger[0].area() > min_width*(ext_perimeter_width + ext_perimeter_spacing2)) {
                                         // the maximum thickness of our thin wall area is equal to the minimum thickness of a single loop
-                                        ex_bigger[0].medial_axis(bound, ext_perimeter_width + ext_perimeter_spacing2, min_width, 
+                                        ex_bigger[0].medial_axis(bound, ext_perimeter_width + ext_perimeter_spacing2, min_width,
                                             &thin_walls, this->layer_height);
                                     }
                                     break;
                                 }
                             }
                         }
+                        std::cout << "thin walls end\n";
                     }
                 } else {
                     //FIXME Is this offset correct if the line width of the inner perimeters differs
@@ -413,6 +415,7 @@ void PerimeterGenerator::process()
                 this->loops->append(entities);
         } // for each loop of an island
 
+        std::cout << "fill gaps begin\n";
         // fill gaps
         if (!gaps.empty()) {
             // collapse 
@@ -445,6 +448,7 @@ void PerimeterGenerator::process()
                 last = diff_ex(to_polygons(last), gap_fill.polygons_covered_by_width(10.f));
             }
         }
+        std::cout << "fill gaps end\n";
 
         // create one more offset to be used as boundary for fill
         // we offset by half the perimeter spacing (to get to the actual infill boundary)
@@ -484,6 +488,7 @@ void PerimeterGenerator::process()
             this->fill_no_overlap.insert(this->fill_no_overlap.end(), polyWithoutOverlap.begin(), polyWithoutOverlap.end());
         }
     } // for each island
+    std::cout << "perimeters end\n";
 }
 
 ExtrusionEntityCollection PerimeterGenerator::_traverse_loops(
diff --git a/xs/src/libslic3r/Polyline.hpp b/xs/src/libslic3r/Polyline.hpp
index 6f70942c4..52846ff2b 100644
--- a/xs/src/libslic3r/Polyline.hpp
+++ b/xs/src/libslic3r/Polyline.hpp
@@ -138,6 +138,7 @@ class ThickPolyline : public Polyline {
     public:
     /// width size must be == point size
     std::vector<coordf_t> width;
+    std::vector<bool> curved;
     /// if true => it's an endpoint, if false it join an other ThickPolyline. first is at front(), second is at back()
     std::pair<bool,bool> endpoints;
     ThickPolyline() : endpoints(std::make_pair(false, false)) {};
diff --git a/xs/src/libslic3r/PrintObject.cpp b/xs/src/libslic3r/PrintObject.cpp
index 91baa3392..4b96b35c5 100644
--- a/xs/src/libslic3r/PrintObject.cpp
+++ b/xs/src/libslic3r/PrintObject.cpp
@@ -1978,6 +1978,7 @@ void PrintObject::_make_perimeters()
         BOOST_LOG_TRIVIAL(debug) << "Generating extra perimeters for region " << region_id << " in parallel - end";
     }
 
+    std::cout << "Generating perimeters in parallel - start";
     BOOST_LOG_TRIVIAL(debug) << "Generating perimeters in parallel - start";
     tbb::parallel_for(
         tbb::blocked_range<size_t>(0, this->layers.size()),
@@ -1987,6 +1988,7 @@ void PrintObject::_make_perimeters()
         }
     );
     BOOST_LOG_TRIVIAL(debug) << "Generating perimeters in parallel - end";
+    std::cout << "Generating perimeters in parallel - end";
 
     /*
         simplify slices (both layer and region slices),