imgproc: fix rotatedRectangleIntersection()

modules/imgproc/src/intersection.cpp

@@ -51,12 +51,13 @@ int rotatedRectangleIntersection( const RotatedRect& rect1, const RotatedRect& r
 {
     CV_INSTRUMENT_REGION()
 
-    const float samePointEps = 0.00001f; // used to test if two points are the same
+    // L2 metric
+    const float samePointEps = std::max(1e-16f, 1e-6f * (float)std::max(rect1.size.area(), rect2.size.area()));
 
     Point2f vec1[4], vec2[4];
     Point2f pts1[4], pts2[4];
 
-    std::vector <Point2f> intersection;
+    std::vector <Point2f> intersection; intersection.reserve(24);
 
     rect1.points(pts1);
     rect2.points(pts2);
@@ -219,41 +220,80 @@ int rotatedRectangleIntersection( const RotatedRect& rect1, const RotatedRect& r
         }
     }
 
-    // Get rid of dupes and order points.
-    for( int i = 0; i < (int)intersection.size()-1; i++ )
+    int N = (int)intersection.size();
+    if (N == 0)
     {
-        float dx1 = intersection[i + 1].x - intersection[i].x;
-        float dy1 = intersection[i + 1].y - intersection[i].y;
-        for( size_t j = i+1; j < intersection.size(); j++ )
-        {
-            float dx = intersection[j].x - intersection[i].x;
-            float dy = intersection[j].y - intersection[i].y;
-            double d2 = dx*dx + dy*dy; // can be a really small number, need double here
+        return INTERSECT_NONE;
+    }
 
-            if( d2 < samePointEps*samePointEps )
+    // Get rid of duplicated points
+    int Nstride = N;
+    cv::AutoBuffer<float, 100> distPt(N * N);
+    cv::AutoBuffer<int> ptDistRemap(N);
+    for (int i = 0; i < N; ++i)
+    {
+        const Point2f pt0 = intersection[i];
+        ptDistRemap[i] = i;
+        for (int j = i + 1; j < N; )
         {
-                // Found a dupe, remove it
-                std::swap(intersection[j], intersection.back());
-                intersection.pop_back();
-                j--; // restart check
+            const Point2f pt1 = intersection[j];
+            float d2 = normL2Sqr<float>(pt1 - pt0);
+            if(d2 <= samePointEps)
+            {
+                if (j < N - 1)
+                    intersection[j] =  intersection[N - 1];
+                N--;
+                continue;
             }
-            else if (dx1 * dy - dy1 * dx < 0)
+            distPt[i*Nstride + j] = d2;
+            ++j;
+        }
+    }
+    while (N > 8) // we still have duplicate points after samePointEps threshold (eliminate closest points)
     {
-                std::swap(intersection[i + 1], intersection[j]);
-                dx1 = dx;
-                dy1 = dy;
+        int minI = 0;
+        int minJ = 1;
+        float minD = distPt[1];
+        for (int i = 0; i < N - 1; ++i)
+        {
+            float* pDist = distPt.data() + Nstride * ptDistRemap[i];
+            for (int j = i + 1; j < N; ++j)
+            {
+                float d = pDist[ptDistRemap[j]];
+                if (d < minD)
+                {
+                    minD = d;
+                    minI = i;
+                    minJ = j;
                 }
             }
         }
-
-    if( intersection.empty() )
+        CV_Assert(fabs(normL2Sqr<float>(intersection[minI] - intersection[minJ]) - minD) < 1e-6);  // ptDistRemap is not corrupted
+        // drop minJ point
+        if (minJ < N - 1)
         {
-        return INTERSECT_NONE ;
+            intersection[minJ] =  intersection[N - 1];
+            ptDistRemap[minJ] = ptDistRemap[N - 1];
+        }
+        N--;
     }
 
-    // If this check fails then it means we're getting dupes, increase samePointEps
-    CV_Assert( intersection.size() <= 8 );
+    // order points
+    for (int i = 0; i < N - 1; ++i)
+    {
+        Point2f diffI = intersection[i + 1] - intersection[i];
+        for (int j = i + 2; j < N; ++j)
+        {
+            Point2f diffJ = intersection[j] - intersection[i];
+            if (diffI.cross(diffJ) < 0)
+            {
+                std::swap(intersection[i + 1], intersection[j]);
+                diffI = diffJ;
+            }
+        }
+    }
 
+    intersection.resize(N);
     Mat(intersection).copyTo(intersectingRegion);
 
     return ret;

modules/imgproc/test/test_intersection.cpp

@@ -118,8 +118,6 @@ void CV_RotatedRectangleIntersectionTest::run(int)
     test12();
     test13();
     test14();
-    test15();
-    test16();
 }
 
 void CV_RotatedRectangleIntersectionTest::test1()
@@ -380,7 +378,7 @@ void CV_RotatedRectangleIntersectionTest::test14()
     }
 }
 
-void CV_RotatedRectangleIntersectionTest::test15()
+TEST(Imgproc_RotatedRectangleIntersection, regression_12221_1)
 {
     RotatedRect r1(
         Point2f(259.65081787109375, 51.58895492553711),
@@ -393,11 +391,11 @@ void CV_RotatedRectangleIntersectionTest::test15()
 
     std::vector<Point2f> intersections;
     int interType = cv::rotatedRectangleIntersection(r1, r2, intersections);
-    EXPECT_TRUE(interType == INTERSECT_PARTIAL);
-    EXPECT_TRUE(intersections.size() == 8);
+    EXPECT_EQ(INTERSECT_PARTIAL, interType);
+    EXPECT_LE(intersections.size(), (size_t)8);
 }
 
-void CV_RotatedRectangleIntersectionTest::test16()
+TEST(Imgproc_RotatedRectangleIntersection, regression_12221_2)
 {
     RotatedRect r1(
         Point2f(239.78500366210938, 515.72021484375),
@@ -410,8 +408,8 @@ void CV_RotatedRectangleIntersectionTest::test16()
 
     std::vector<Point2f> intersections;
     int interType = cv::rotatedRectangleIntersection(r1, r2, intersections);
-    EXPECT_TRUE(interType == INTERSECT_PARTIAL);
-    EXPECT_TRUE(intersections.size() == 8);
+    EXPECT_EQ(INTERSECT_PARTIAL, interType);
+    EXPECT_LE(intersections.size(), (size_t)8);
 }
 
 TEST (Imgproc_RotatedRectangleIntersection, accuracy) { CV_RotatedRectangleIntersectionTest test; test.safe_run(); }