二維平面上判斷點是否在三角形內

 

 

1. 已知三角形的三個頂點座標，判斷某個點是否在三角形中（在三角形的邊上，咱們也視做在三角形中），本文給出了三種方法。

算法1：利用面積法

如上圖所示，若是點P在三角形ABC的內部，則三個小三角形PAB, PBC, PAC的面積之和 = ABC的面積，反之則不相等。

已知三角形的三個頂點座標求其面積，能夠根據向量的叉乘

　　該算法詳見後面代碼中的函數：IsPointInTriangle1

算法2： 

　　首先看一下這個問題，如何判斷某兩個點在某條直線的同一側（代碼中函數：IsPointAtSameSideOfLine）？

　　

　　根據向量的叉乘以及右手螺旋定則，AB^AM （^表示叉乘，這裏向量省略了字母上面的箭頭符號）的方向爲向外指出屏幕，AB^AN也是向外指出屏幕，但AB^AO的方向是向內指向屏幕，所以M,N在直線AB的同側，M ,O在直線AB的兩側。實際計算時，只須要考慮叉積的數值正負。

假設以上各點座標爲A(0,0), B(4,0), M(1,2), N(3,4), O(3,-4)， 則：

AB^AM = (4,0)^(1,2) = 4*2 - 0*1 = 8

AB^AN = (4，0)^(3，4) = 4*4 – 0*3 = 16

AB^AO = (4,0)^(3,-4) = 4*-4 – 0*3 = –16

 　　由上面的數值可知，能夠根據數值的正負判斷叉乘後向量的方向。即，若是叉積AB^AM 和 AB^AN的結果同號，那麼M,N兩點就在直線的同側，不然不在同一側。特殊地，若是點M在直線AB上，則AB^AM的值爲0。（若是是在三維座標系中，求出的叉積是一個向量，能夠根據兩個向量的點積結果正負來判斷兩個向量的是否指向同一側）    

　　以上的問題解決了，就很容易的用來判斷某個點是否在三角形內，若是P在三角形ABC內部，則知足如下三個條件：P,A在BC的同側、P,B在AC的同側、PC在AB的同側。某一個不知足則表示P不在三角形內部。

該算法詳見後面代碼中的函數：IsPointInTriangle2

算法3： 

　　該方法也用到了向量。對於三角形ABC和一點P，能夠有以下的向量表示：

　　

p點在三角形內部的充分必要條件是：1 >= u >= 0,   1 >= v >= 0,   u+v <= 1。

 　　已知A,B,C,P四個點的座標，能夠求出u，v，把上面的式子分別點乘向量AC和向量AB

　　

解方程獲得：

　　

解出u，v後只須要看他們是否知足「1 >= u >= 0, 1 >= v >= 0, u+v <= 1」，如知足，則，p 在三角形內。

（u = 0時，p在AB上， v = 0時，p在AC上，二者均爲0時，p和A重合）

該算法詳見後面代碼中的函數：IsPointInTriangle3

 

#include<string>
 #include<sstream>
 #include<fstream>
 #include<cstdio>
 #include<cstring>
 #include<iostream>
 #include<vector>
 #include<stack>
 #include<queue>
 #include<list>
 #include<algorithm>
 #include<ctime>
 #include<unordered_map>
 #include<map>
 #include<typeinfo>
 #include<cmath>
 #include<ctime>
 #include<climits>
 using namespace std;
 
//類定義：二維向量
class Vector2d
{
public:
    double x_;
    double y_;
 
public:
    Vector2d(double x, double y):x_(x), y_(y){}
    Vector2d():x_(0), y_(0){}
 
    //二維向量叉乘, 叉乘的結果實際上是向量，方向垂直於兩個向量組成的平面，這裏咱們只須要其大小和方向
    double CrossProduct(const Vector2d vec)
    {
        return x_*vec.y_ - y_*vec.x_;
    }
 
    //二維向量點積
    double DotProduct(const Vector2d vec)
    {
        return x_ * vec.x_ + y_ * vec.y_;
    }
 
    //二維向量減法
    Vector2d Minus(const Vector2d vec) const
    {
        return Vector2d(x_ - vec.x_, y_ - vec.y_);
    }
 
    //判斷點M,N是否在直線AB的同一側
    static bool IsPointAtSameSideOfLine(const Vector2d &pointM, const Vector2d &pointN, 
                                        const Vector2d &pointA, const Vector2d &pointB)
    {
        Vector2d AB = pointB.Minus(pointA);
        Vector2d AM = pointM.Minus(pointA);
        Vector2d AN = pointN.Minus(pointA);
 
        //等於0時表示某個點在直線上
        return AB.CrossProduct(AM) * AB.CrossProduct(AN) >= 0;
    }
};
 
//三角形類
class Triangle
{
private:
    Vector2d pointA_, pointB_, pointC_;
 
public:
    Triangle(Vector2d point1, Vector2d point2, Vector2d point3)
        :pointA_(point1), pointB_(point2), pointC_(point3)
    {
        //todo 判斷三點是否共線
    }
 
    //計算三角形面積
    double ComputeTriangleArea()
    {
        //依據兩個向量的叉乘來計算，可參考http://blog.csdn.net/zxj1988/article/details/6260576
        Vector2d AB = pointB_.Minus(pointA_);
        Vector2d BC = pointC_.Minus(pointB_);
        return fabs(AB.CrossProduct(BC) / 2.0);
    }
 
    bool IsPointInTriangle1(const Vector2d pointP)
    {
        double area_ABC = ComputeTriangleArea();
        double area_PAB = Triangle(pointP, pointA_, pointB_).ComputeTriangleArea();
        double area_PAC = Triangle(pointP, pointA_, pointC_).ComputeTriangleArea();
        double area_PBC = Triangle(pointP, pointB_, pointC_).ComputeTriangleArea();
 
        if(fabs(area_PAB + area_PBC + area_PAC - area_ABC) < 0.000001)
            return true;
        else return false;
    }
 
    bool IsPointInTriangle2(const Vector2d pointP)
    {
        return Vector2d::IsPointAtSameSideOfLine(pointP, pointA_, pointB_, pointC_) &&
            Vector2d::IsPointAtSameSideOfLine(pointP, pointB_, pointA_, pointC_) &&
            Vector2d::IsPointAtSameSideOfLine(pointP, pointC_, pointA_, pointB_);
    }
 
    bool IsPointInTriangle3(const Vector2d pointP)
    {
        Vector2d AB = pointB_.Minus(pointA_);
        Vector2d AC = pointC_.Minus(pointA_);
        Vector2d AP = pointP.Minus(pointA_);
        double dot_ac_ac = AC.DotProduct(AC);
        double dot_ac_ab = AC.DotProduct(AB);
        double dot_ac_ap = AC.DotProduct(AP);
        double dot_ab_ab = AB.DotProduct(AB);
        double dot_ab_ap = AB.DotProduct(AP);
 
        double tmp = 1.0 / (dot_ac_ac * dot_ab_ab - dot_ac_ab * dot_ac_ab);
         
        double u = (dot_ab_ab * dot_ac_ap - dot_ac_ab * dot_ab_ap) * tmp;
        if(u < 0 || u > 1)
            return false;
        double v = (dot_ac_ac * dot_ab_ap - dot_ac_ab * dot_ac_ap) * tmp;
        if(v < 0 || v > 1)
            return false;
 
        return u + v <= 1;
    }
 
    bool IsPointInTriangle4(const Vector2d pointP)
    {
        Vector2d PA = pointA_.Minus(pointP);
        Vector2d PB = pointB_.Minus(pointP);
        Vector2d PC = pointC_.Minus(pointP);
        double t1 = PA.CrossProduct(PB);
        double t2 = PB.CrossProduct(PC);
        double t3 = PC.CrossProduct(PA);
        return t1*t2 >= 0 && t1*t3 >= 0;
    }
};
 
 
int main()
{
    Triangle tri(Vector2d(0,0), Vector2d(6,6), Vector2d(12,0));
    srand(time(0));
    for(int i = 0; i < 20; ++i)
    {
         Vector2d point((rand()*1.0 / RAND_MAX) * 12, (rand()*1.0 / RAND_MAX) * 6);
         cout<<point.x_<<" "<<point.y_<<":     ";
         cout<<tri.IsPointInTriangle1(point)<<" ";
         cout<<tri.IsPointInTriangle2(point)<<" ";
         cout<<tri.IsPointInTriangle3(point)<<" ";
         cout<<tri.IsPointInTriangle4(point)<<endl;
 
    }
}