Update fast crossing

build.py

@@ -6,7 +6,7 @@
 from setuptools import setup, find_packages
 
 name = "cubao_headers"
-version = "0.0.5"
+version = "0.0.6"
 
 pattern = ["*"]
 

include/cubao/polyline_in_polygon.hpp

@@ -26,44 +26,50 @@ polyline_in_polygon(const RowVectors &polyline, //
                     const FastCrossing &fc)
 {
     auto intersections = fc.intersections(polyline);
-    // pt, (t, s), cur_label=(poly1, seg1), tree_label=(poly2, seg2)
     auto ruler = PolylineRuler(polyline, fc.is_wgs84());
-    // 0.0, [r1, r2, ..., length]
-    const int N = intersections.size() + 1;
-    Eigen::VectorXd ranges(N);
-    int idx = -1;
-    for (auto &inter : intersections) {
-        int seg_idx = std::get<2>(inter)[1];
-        double t = std::get<1>(inter)[0];
-        double r = ruler.range(seg_idx, t);
-        ranges[++idx] = r;
+    if (intersections.empty()) {
+        int inside = point_in_polygon(polyline.block(0, 0, 1, 2), polygon)[0];
+        if (!inside) {
+            return {};
+        }
+        return PolylineChunks{
+            {{0, 0.0, 0.0, ruler.N() - 2, 1.0, ruler.length()}, polyline}};
     }
-    ranges[++idx] = ruler.length();
-    RowVectorsNx2 midpoints(N, 3);
+    // pt, (t, s), cur_label=(poly1, seg1), tree_label=(poly2, seg2)
+    const int N = intersections.size() + 2;
+    // init ranges
+    Eigen::VectorXd ranges(N);
     {
-        idx = 0;
-        double r = 0.0;
-        while (idx < N) {
-            double rr = ranges[idx];
-            midpoints.row(idx) = ruler.at((r + rr) / 2.0).head(2);
-            r = rr;
-            ++idx;
+        int idx = -1;
+        ranges[++idx] = 0.0;
+        for (auto &inter : intersections) {
+            int seg_idx = std::get<2>(inter)[1];
+            double t = std::get<1>(inter)[0];
+            double r = ruler.range(seg_idx, t);
+            ranges[++idx] = r;
         }
+        ranges[++idx] = ruler.length();
+    }
+    // ranges o------o--------o-----------------o
+    // midpts     ^       ^             ^
+    RowVectorsNx2 midpoints(N - 1, 2);
+    for (int i = 0; i < N - 1; ++i) {
+        double rr = (ranges[i] + ranges[i + 1]) / 2.0;
+        midpoints.row(i) = ruler.along(rr).head(2);
     }
     auto mask = point_in_polygon(midpoints, polygon);
     PolylineChunks ret;
     {
-        idx = 0;
-        double r = 0.0;
-        while (idx < N) {
-            if (mask[idx] && ranges[idx] > r) {
-                auto [seg1, t1] = ruler.segment_index_t(r);
-                auto [seg2, t2] = ruler.segment_index_t(ranges[idx]);
-                ret.emplace(std::make_tuple(seg1, t1, r, seg2, t2, ranges[idx]),
-                            ruler.lineSliceAlong(r, ranges[idx]));
+        for (int i = 0; i < N - 1; ++i) {
+            double r1 = ranges[i];
+            double r2 = ranges[i + 1];
+            if (r2 <= r1 || mask[i] == 0) {
+                continue;
             }
-            r = ranges[idx];
-            ++idx;
+            auto [seg1, t1] = ruler.segment_index_t(r1);
+            auto [seg2, t2] = ruler.segment_index_t(r2);
+            ret.emplace(std::make_tuple(seg1, t1, r1, seg2, t2, r2),
+                        ruler.lineSliceAlong(r1, r2));
         }
     }
     return ret;