OpenCV中两个旋转矩形（RotatedRect）的交集计算优化技巧

namespace cv
{

int rotatedRectangleIntersection( const RotatedRect& rect1, const RotatedRect& rect2, OutputArray intersectingRegion )
{
    CV_INSTRUMENT_REGION();

    // 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; intersection.reserve(24);

    rect1.points(pts1);
    rect2.points(pts2);

    int ret = INTERSECT_FULL;

    // Specical case of rect1 == rect2
    {
        bool same = true;

        for( int i = 0; i < 4; i++ )
        {
            if( fabs(pts1[i].x - pts2[i].x) > samePointEps || (fabs(pts1[i].y - pts2[i].y) > samePointEps) )
            {
                same = false;
                break;
            }
        }

        if(same)
        {
            intersection.resize(4);

            for( int i = 0; i < 4; i++ )
            {
                intersection[i] = pts1[i];
            }

            Mat(intersection).copyTo(intersectingRegion);

            return INTERSECT_FULL;
        }
    }

    // Line vector
    // A line from p1 to p2 is: p1 + (p2-p1)*t, t=[0,1]
    for( int i = 0; i < 4; i++ )
    {
        vec1[i].x = pts1[(i+1)%4].x - pts1[i].x;
        vec1[i].y = pts1[(i+1)%4].y - pts1[i].y;

        vec2[i].x = pts2[(i+1)%4].x - pts2[i].x;
        vec2[i].y = pts2[(i+1)%4].y - pts2[i].y;
    }

    // Line test - test all line combos for intersection
    for( int i = 0; i < 4; i++ )
    {
        for( int j = 0; j < 4; j++ )
        {
            // Solve for 2x2 Ax=b
            float x21 = pts2[j].x - pts1[i].x;
            float y21 = pts2[j].y - pts1[i].y;

            float vx1 = vec1[i].x;
            float vy1 = vec1[i].y;

            float vx2 = vec2[j].x;
            float vy2 = vec2[j].y;

            float det = vx2*vy1 - vx1*vy2;

            float t1 = (vx2*y21 - vy2*x21) / det;
            float t2 = (vx1*y21 - vy1*x21) / det;

            // This takes care of parallel lines
            if( cvIsInf(t1) || cvIsInf(t2) || cvIsNaN(t1) || cvIsNaN(t2) )
            {
                continue;
            }

            if( t1 >= 0.0f && t1 <= 1.0f && t2 >= 0.0f && t2 <= 1.0f )
            {
                float xi = pts1[i].x + vec1[i].x*t1;
                float yi = pts1[i].y + vec1[i].y*t1;

                intersection.push_back(Point2f(xi,yi));
            }
        }
    }

    if( !intersection.empty() )
    {
        ret = INTERSECT_PARTIAL;
    }

    // Check for vertices from rect1 inside recct2
    for( int i = 0; i < 4; i++ )
    {
        // We do a sign test to see which side the point lies.
        // If the point all lie on the same sign for all 4 sides of the rect,
        // then there's an intersection
        int posSign = 0;
        int negSign = 0;

        float x = pts1[i].x;
        float y = pts1[i].y;

        for( int j = 0; j < 4; j++ )
        {
            // line equation: Ax + By + C = 0
            // see which side of the line this point is at
            float A = -vec2[j].y;
            float B = vec2[j].x;
            float C = -(A*pts2[j].x + B*pts2[j].y);

            float s = A*x+ B*y+ C;

            if( s >= 0 )
            {
                posSign++;
            }
            else
            {
                negSign++;
            }
        }

        if( posSign == 4 || negSign == 4 )
        {
            intersection.push_back(pts1[i]);
        }
    }

    // Reverse the check - check for vertices from rect2 inside recct1
    for( int i = 0; i < 4; i++ )
    {
        // We do a sign test to see which side the point lies.
        // If the point all lie on the same sign for all 4 sides of the rect,
        // then there's an intersection
        int posSign = 0;
        int negSign = 0;

        float x = pts2[i].x;
        float y = pts2[i].y;

        for( int j = 0; j < 4; j++ )
        {
            // line equation: Ax + By + C = 0
            // see which side of the line this point is at
            float A = -vec1[j].y;
            float B = vec1[j].x;
            float C = -(A*pts1[j].x + B*pts1[j].y);

            float s = A*x + B*y + C;

            if( s >= 0 )
            {
                posSign++;
            }
            else
            {
                negSign++;
            }
        }

        if( posSign == 4 || negSign == 4 )
        {
            intersection.push_back(pts2[i]);
        }
    }

    int N = (int)intersection.size();
    if (N == 0)
    {
        return INTERSECT_NONE;
    }

    // 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; )
        {
            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;
            }
            distPt[i*Nstride + j] = d2;
            ++j;
        }
    }
    while (N > 8) // we still have duplicate points after samePointEps threshold (eliminate closest points)
    {
        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;
                }
            }
        }
        CV_Assert(fabs(normL2Sqr<float>(intersection[minI] - intersection[minJ]) - minD) < 1e-6);  // ptDistRemap is not corrupted
        // drop minJ point
        if (minJ < N - 1)
        {
            intersection[minJ] =  intersection[N - 1];
            ptDistRemap[minJ] = ptDistRemap[N - 1];
        }
        N--;
    }

    // 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;
}

} // end namespace

这段代码实现了OpenCV中两个旋转矩形（RotatedRect）的交集计算。以下是对代码的分析和一些优化技巧：

代码分析

​​基本思路​​：
计算两个旋转矩形的所有边之间的交点
检查一个矩形的顶点是否在另一个矩形内部
合并所有交点并去除重复点
对交点进行排序以形成凸多边形
​​

主要步骤​​：

处理完全相同的矩形特殊情况
计算所有边之间的交点（16种组合）
检查矩形顶点是否在另一个矩形内
去除重复点
对点进行排序

优化技巧

​​数学计算优化​​：
使用normL2Sqr代替norm避免平方根计算
提前计算并重用向量和线方程系数
​​内存优化​​：
使用AutoBuffer而不是动态分配内存
预分配空间intersection.reserve(24)
​​数值稳定性​​：
使用samePointEps处理浮点精度问题
检查cvIsInf和cvIsNaN避免数值不稳定
​​算法优化​​：
提前终止检查（如发现完全相同的矩形）
使用符号测试判断点是否在矩形内

可能的改进

​​更高效的凸包算法​​：
当前的点排序算法复杂度较高，可以考虑Andrew’s monotone chain算法
​​并行化​​：
边交点计算可以并行化（4x4组合）
​​更精确的数值处理​​：
使用更高精度的中间计算减少误差累积
​​提前裁剪测试​​：
在详细计算前先进行简单的包围盒测试
​​更智能的重复点处理​​：
当前方法在极端情况下可能需要多次遍历

使用建议

对于性能敏感的应用，可以考虑：
如果矩形大多是轴对齐的，使用更简单的算法
添加快速拒绝测试（如分离轴测试）
对于高精度需求：
可以增加samePointEps的值
使用双精度计算关键部分
这段代码整体设计良好，考虑了数值稳定性和性能，但在极端情况下（如非常接近的矩形）可能会有性能波动。