Merge pull request #9802 from Nickolays:Fix#9797

modules/imgproc/src/shapedescr.cpp

@@ -42,11 +42,7 @@
 namespace cv
 {
 
-// inner product
-static float innerProduct(Point2f &v1, Point2f &v2)
-{
-    return v1.x * v2.y - v1.y * v2.x;
-}
+const float EPS = 1.0e-4f;
 
 static void findCircle3pts(Point2f *pts, Point2f &center, float &radius)
 {
@@ -54,72 +50,21 @@ static void findCircle3pts(Point2f *pts, Point2f &center, float &radius)
     Point2f v1 = pts[1] - pts[0];
     Point2f v2 = pts[2] - pts[0];
 
-    if (innerProduct(v1, v2) == 0.0f)
-    {
-        // v1, v2 colineation, can not determine a unique circle
-        // find the longtest distance as diameter line
-        float d1 = (float)norm(pts[0] - pts[1]);
-        float d2 = (float)norm(pts[0] - pts[2]);
-        float d3 = (float)norm(pts[1] - pts[2]);
-        if (d1 >= d2 && d1 >= d3)
-        {
-            center = (pts[0] + pts[1]) / 2.0f;
-            radius = (d1 / 2.0f);
-        }
-        else if (d2 >= d1 && d2 >= d3)
-        {
-            center = (pts[0] + pts[2]) / 2.0f;
-            radius = (d2 / 2.0f);
-        }
-        else if (d3 >= d1 && d3 >= d2)
-        {
-            center = (pts[1] + pts[2]) / 2.0f;
-            radius = (d3 / 2.0f);
-        }
-    }
-    else
-    {
-        // center is intersection of midperpendicular lines of the two edges v1, v2
-        // a1*x + b1*y = c1 where a1 = v1.x, b1 = v1.y
-        // a2*x + b2*y = c2 where a2 = v2.x, b2 = v2.y
-        Point2f midPoint1 = (pts[0] + pts[1]) / 2.0f;
-        float c1 = midPoint1.x * v1.x + midPoint1.y * v1.y;
-        Point2f midPoint2 = (pts[0] + pts[2]) / 2.0f;
-        float c2 = midPoint2.x * v2.x + midPoint2.y * v2.y;
-        float det = v1.x * v2.y - v1.y * v2.x;
-        float cx = (c1 * v2.y - c2 * v1.y) / det;
-        float cy = (v1.x * c2 - v2.x * c1) / det;
-        center.x = (float)cx;
-        center.y = (float)cy;
-        cx -= pts[0].x;
-        cy -= pts[0].y;
-        radius = (float)(std::sqrt(cx *cx + cy * cy));
-    }
-}
-
-const float EPS = 1.0e-4f;
-
-static void findEnclosingCircle3pts_orLess_32f(Point2f *pts, int count, Point2f &center, float &radius)
-{
-    switch (count)
-    {
-    case 1:
-        center = pts[0];
-        radius = 0.0f;
-        break;
-    case 2:
-        center.x = (pts[0].x + pts[1].x) / 2.0f;
-        center.y = (pts[0].y + pts[1].y) / 2.0f;
-        radius = (float)(norm(pts[0] - pts[1]) / 2.0);
-        break;
-    case 3:
-        findCircle3pts(pts, center, radius);
-        break;
-    default:
-        break;
-    }
-
-    radius += EPS;
+    // center is intersection of midperpendicular lines of the two edges v1, v2
+    // a1*x + b1*y = c1 where a1 = v1.x, b1 = v1.y
+    // a2*x + b2*y = c2 where a2 = v2.x, b2 = v2.y
+    Point2f midPoint1 = (pts[0] + pts[1]) / 2.0f;
+    float c1 = midPoint1.x * v1.x + midPoint1.y * v1.y;
+    Point2f midPoint2 = (pts[0] + pts[2]) / 2.0f;
+    float c2 = midPoint2.x * v2.x + midPoint2.y * v2.y;
+    float det = v1.x * v2.y - v1.y * v2.x;
+    float cx = (c1 * v2.y - c2 * v1.y) / det;
+    float cy = (v1.x * c2 - v2.x * c1) / det;
+    center.x = (float)cx;
+    center.y = (float)cy;
+    cx -= pts[0].x;
+    cy -= pts[0].y;
+    radius = (float)(std::sqrt(cx *cx + cy * cy)) + EPS;
 }
 
 template<typename PT>
@@ -145,7 +90,7 @@ static void findThirdPoint(const PT *pts, int i, int j, Point2f &center, float &
             ptsf[0] = (Point2f)pts[i];
             ptsf[1] = (Point2f)pts[j];
             ptsf[2] = (Point2f)pts[k];
-            findEnclosingCircle3pts_orLess_32f(ptsf, 3, center, radius);
+            findCircle3pts(ptsf, center, radius);
         }
     }
 }
@@ -210,8 +155,6 @@ void cv::minEnclosingCircle( InputArray _points, Point2f& _center, float& _radiu
     Mat points = _points.getMat();
     int count = points.checkVector(2);
     int depth = points.depth();
-    Point2f center;
-    float radius = 0.f;
     CV_Assert(count >= 0 && (depth == CV_32F || depth == CV_32S));
 
     _center.x = _center.y = 0.f;
@@ -224,52 +167,62 @@ void cv::minEnclosingCircle( InputArray _points, Point2f& _center, float& _radiu
     const Point* ptsi = points.ptr<Point>();
     const Point2f* ptsf = points.ptr<Point2f>();
 
-    // point count <= 3
-    if (count <= 3)
+    switch (count)
     {
-        Point2f ptsf3[3];
-        for (int i = 0; i < count; ++i)
+        case 1:
         {
-            ptsf3[i] = (is_float) ? ptsf[i] : Point2f((float)ptsi[i].x, (float)ptsi[i].y);
+            _center = (is_float) ? ptsf[0] : Point2f((float)ptsi[0].x, (float)ptsi[0].y);
+            _radius = EPS;
+            break;
         }
-        findEnclosingCircle3pts_orLess_32f(ptsf3, count, center, radius);
-        _center = center;
-        _radius = radius;
-        return;
-    }
-
-    if (is_float)
-    {
-        findMinEnclosingCircle<Point2f>(ptsf, count, center, radius);
-#if 0
-        for (size_t m = 0; m < count; ++m)
+        case 2:
         {
-            float d = (float)norm(ptsf[m] - center);
-            if (d > radius)
-            {
-                printf("error!\n");
-            }
+            Point2f p1 = (is_float) ? ptsf[0] : Point2f((float)ptsi[0].x, (float)ptsi[0].y);
+            Point2f p2 = (is_float) ? ptsf[1] : Point2f((float)ptsi[1].x, (float)ptsi[1].y);
+            _center.x = (p1.x + p2.x) / 2.0f;
+            _center.y = (p1.y + p2.y) / 2.0f;
+            _radius = (float)(norm(p1 - p2) / 2.0) + EPS;
+            break;
         }
-#endif
-    }
-    else
-    {
-        findMinEnclosingCircle<Point>(ptsi, count, center, radius);
-#if 0
-        for (size_t m = 0; m < count; ++m)
+        default:
         {
-            double dx = ptsi[m].x - center.x;
-            double dy = ptsi[m].y - center.y;
-            double d = std::sqrt(dx * dx + dy * dy);
-            if (d > radius)
+            Point2f center;
+            float radius = 0.f;
+            if (is_float)
+            {
+                findMinEnclosingCircle<Point2f>(ptsf, count, center, radius);
+                #if 0
+                    for (size_t m = 0; m < count; ++m)
+                    {
+                        float d = (float)norm(ptsf[m] - center);
+                        if (d > radius)
+                        {
+                            printf("error!\n");
+                        }
+                    }
+                #endif
+            }
+            else
             {
-                printf("error!\n");
+                findMinEnclosingCircle<Point>(ptsi, count, center, radius);
+                #if 0
+                    for (size_t m = 0; m < count; ++m)
+                    {
+                        double dx = ptsi[m].x - center.x;
+                        double dy = ptsi[m].y - center.y;
+                        double d = std::sqrt(dx * dx + dy * dy);
+                        if (d > radius)
+                        {
+                            printf("error!\n");
+                        }
+                    }
+                #endif
             }
+            _center = center;
+            _radius = radius;
+            break;
         }
-#endif
     }
-    _center = center;
-    _radius = radius;
 }
 
 

modules/imgproc/test/test_convhull.cpp

@@ -1017,6 +1017,62 @@ _exit_:
     return code;
 }
 
+/****************************************************************************************\
+*                                 MinEnclosingCircle Test 2                              *
+\****************************************************************************************/
+
+class CV_MinCircleTest2 : public CV_BaseShapeDescrTest
+{
+public:
+    CV_MinCircleTest2();
+protected:
+    RNG rng;
+    void run_func(void);
+    int validate_test_results( int test_case_idx );
+    float delta;
+};
+
+
+CV_MinCircleTest2::CV_MinCircleTest2()
+{
+    rng = ts->get_rng();
+}
+
+
+void CV_MinCircleTest2::run_func()
+{
+    Point2f center = Point2f(rng.uniform(0.0f, 1000.0f), rng.uniform(0.0f, 1000.0f));;
+    float radius = rng.uniform(0.0f, 500.0f);
+    float angle = (float)rng.uniform(0.0f, (float)(CV_2PI));
+    vector<Point2f> pts;
+    pts.push_back(center + Point2f(radius * cos(angle), radius * sin(angle)));
+    angle += (float)CV_PI;
+    pts.push_back(center + Point2f(radius * cos(angle), radius * sin(angle)));
+    float radius2 = radius * radius;
+    float x = rng.uniform(center.x - radius, center.x + radius);
+    float deltaX = x - center.x;
+    float upperBoundY = sqrt(radius2 - deltaX * deltaX);
+    float y = rng.uniform(center.y - upperBoundY, center.y + upperBoundY);
+    pts.push_back(Point2f(x, y));
+    // Find the minimum area enclosing circle
+    Point2f calcCenter;
+    float calcRadius;
+    minEnclosingCircle(pts, calcCenter, calcRadius);
+    delta = (float)norm(calcCenter - center) + abs(calcRadius - radius);
+}
+
+int CV_MinCircleTest2::validate_test_results( int test_case_idx )
+{
+    float eps = 1.0F;
+    int code = CV_BaseShapeDescrTest::validate_test_results( test_case_idx );
+    if (delta > eps)
+    {
+        ts->printf( cvtest::TS::LOG, "Delta center and calcCenter > %f\n", eps );
+        code = cvtest::TS::FAIL_BAD_ACCURACY;
+        ts->set_failed_test_info( code );
+    }
+    return code;
+}
 
 /****************************************************************************************\
 *                                   Perimeter Test                                     *
@@ -1905,6 +1961,7 @@ TEST(Imgproc_ConvexHull, accuracy) { CV_ConvHullTest test; test.safe_run(); }
 TEST(Imgproc_MinAreaRect, accuracy) { CV_MinAreaRectTest test; test.safe_run(); }
 TEST(Imgproc_MinTriangle, accuracy) { CV_MinTriangleTest test; test.safe_run(); }
 TEST(Imgproc_MinCircle, accuracy) { CV_MinCircleTest test; test.safe_run(); }
+TEST(Imgproc_MinCircle2, accuracy) { CV_MinCircleTest2 test; test.safe_run(); }
 TEST(Imgproc_ContourPerimeter, accuracy) { CV_PerimeterTest test; test.safe_run(); }
 TEST(Imgproc_FitEllipse, accuracy) { CV_FitEllipseTest test; test.safe_run(); }
 TEST(Imgproc_FitEllipse, parallel) { CV_FitEllipseParallelTest test; test.safe_run(); }