GCompute

GCompute - 计算三角形对的交点（3) - 共面三角形相交

在上两篇文章中，我们介绍了使用不同的方法，计算空间中不共面的两个三角形的交点，这篇文章将介绍如何对空间中共面的两个三角形求交。

完整代码如下：

#include <iostream>
#include <cmath>
#include <algorithm>
#include <assert.h>
#include <vector>
#include <iterator>

#include <Eigen/Dense>

typedef Eigen::Vector2d Point2d;
typedef Eigen::Vector2d Vector2d;

struct Triangle2D
{
    Point2d m_pt[3];
};

struct Segment2D
{
    Point2d m_pt[2];
};

bool CalSegmentIntersectionPoint(const Segment2D& seg1, const Segment2D& seg2, Point2d& pt)
{
    // pt0, pt1 on segment
    // dir = pt1 - pt0
    // line(x,y) = pt0 + t * dir
    // line1 line2 intersection cal:
    // => 
    // 0 = pt0on1 - pt0on2 + ton1 * diron1 - ton2 * diron2 =>
    // pt0on2 - pt0on1 = ton1 * diron1 - ton2 * diron2  =>
    // mat(diron1.x, -diron2.x; diron1.y, -diron2.y) * mat(ton1; ton2) = mat(pt0on2.x - pt0on1.x; pt0on2.y - pt0on1.y) = mat(a, b)
    Vector2d dirOn1 = seg1.m_pt[1] - seg1.m_pt[0];
    Vector2d dirOn2 = seg2.m_pt[1] - seg2.m_pt[0];

    double a = seg2.m_pt[0].x() - seg1.m_pt[0].x();
    double b = seg2.m_pt[0].y() - seg1.m_pt[0].y();

    double denom = dirOn1.x() * (-dirOn2.y()) - (-dirOn2.x()) * (dirOn1.y());
    double denom_t1 = a * (-dirOn2.y()) - (-dirOn2.x()) * b;
    double denom_t2 = dirOn1.x() * b - a * (dirOn1.y());

    double t1 = denom_t1 / denom;
    double t2 = denom_t2 / denom;

    if (t1 >= 0 && t1 <= 1 && t2 >= 0 && t2 <= 1)
    {
        pt = seg1.m_pt[0] + t1 * dirOn1;
        return true;
    }

    return false;
}

void GetSegmentFromTriangle(const Triangle2D& tri, Segment2D (&seg)[3])
{
    seg[0].m_pt[0] = tri.m_pt[0];
    seg[0].m_pt[1] = tri.m_pt[1];

    seg[1].m_pt[0] = tri.m_pt[1];
    seg[1].m_pt[1] = tri.m_pt[2];

    seg[2].m_pt[0] = tri.m_pt[2];
    seg[2].m_pt[1] = tri.m_pt[0];
}

void GetTriangleIntersectionPoints(const Triangle2D& tri1, const Triangle2D& tri2, std::vector<Point2d>& pts)
{
    Segment2D seg1[3], seg2[3];
    GetSegmentFromTriangle(tri1, seg1);
    GetSegmentFromTriangle(tri2, seg2);

    Point2d pt;
    for (int i = 0; i < 3; ++i)
    {
        for (int j = 0; j < 3; ++j)
        {
            if (CalSegmentIntersectionPoint(seg1[i], seg2[j], pt))
            {
                pts.emplace_back(pt);
            }
        }
    }
}

void TriIntersectTestCase()
{
    Triangle2D tri1, tri2;
    tri1.m_pt[0] = Point2d(4, 3); tri1.m_pt[1] = Point2d(6, -2); tri1.m_pt[2] = Point2d(-2, 1);
    tri2.m_pt[0] = Point2d(1, -3); tri2.m_pt[1] = Point2d(2, 1); tri2.m_pt[2] = Point2d(13, -2);

    std::vector<Point2d> pts;
    GetTriangleIntersectionPoints(tri1, tri2, pts);
    std::cout << "Intersection points: 
";
    for (int i = 0; i < pts.size(); ++i)
    {
        std::cout << "=====" << "
";
        std::cout << pts[i] << "
";
    }
}

int main()
{
    TriIntersectTestCase();
    return 0;
}

得到结果如下：

Intersection points:
=====
 5.14286
0.142857
=====
  1.65714
-0.371429

示意图如下：