这里讨论的是二维平面中的算法

基础常识

这里一改之前草率的作风，使用类来进行构造。

在这里默认把一个点看做还是一个以原点为起点的向量。

最开始的类

#include <bits/stdc++.h>
using namespace std;
class Point{
   
public:
    double x, y;
    Point(){
   
        x = 0, y = 0;
    }//默认构造函数
    Point(double x1, double y1){
   //重载构造函数
        x = x1;
        y = y1;
    }
    void disp()
    {
   
        printf("(%g, %g)", x, y);
    }
};
int main()
{
   
    Point p1(1, 2);
    p1.disp();   
    return 0;
}

实现向量的相关简单运算

//加运算
Point operator +(const Point &p1, const Point &p2)
{
   
    return Point(p1.x + p2.x, p1.y + p2.y);
}

//减运算
Point operator -(const Point &p1, const Point &p2)
{
   
    return Point(p1.x - p2.x, p1.y - p2.y);
}

//相等
bool operator ==(const Point &p1, const Point &p2)
{
   
    if(p1.x == p2.x && p1.y == p2.y)
        return true;
    else 
        return false;
}

//点积运算
double Dot(const Point &p1,const Point &p2)
{
   
    return p1.x*p2.x + p1.y * p2.y;
}
/* 通过点积可以判断两个向量的夹角 */

//求向量的模
double Length(const Point &p)
{
   
    return sqrt(Dot(p, p));
}

//对于具有相同起点的两个向量的夹角的计算
int Angle(const Point &p0, const Point &p1, const Point &p2)
{
   
    double t = Dot((p1-p0), (p2-p0));
    if(t > 1e-8) return 1;//锐角
    else if(t < -(1e-8)) return -1;//钝角
    else return 0;//直角
}

//pay attention to：
//如果是一个二维向量，那么叉积为一个标量
double Det(const Point &p1, const Point p2)
{
   
    return p1.x * p2.y - p1.y * p2.x;
}
//两点之间的距离（欧式距离）
double Distance(const Point &p1, const Point &p2)
{
   
    return Length(p2-p1);
}

叉积详解

叉积在上面已经提到过。叉积在运算中具有很重要的地位。

叉积的算法：p1.x * p2.y - p1.y * p2.x;

叉积的意义：

绝对值大小：表示(0, 0), p1, p2, p1+p2这四个点所组成的平行四边形的面积。
符号：（前提：Det(p1, p2)）
正：p1逆时针到p2
负：p1逆时针到p2
0：共线，但是可能同向，也可能反向

如果对于p0->p1, p0->p2

p0, p1, p2在右手螺旋方向上Det( (p1-p0), (p2-p0) ) > 0

p0, p1, p2在左手螺旋方向上Det( (p1-p0), (p2-p0) ) < 0

当仅仅使用叉积的大小时，别忘记加绝对值

点到线段的距离

有以下三种情况

double DistPtoSegment(const Point &p0, const Point &p1, const Point &p2)
{
   
    if(p1 == p2)
    {
   
        return Distance(p1, p0);
    }
    if(Dot((p0 - p1), (p2-p1)) < 0)
    {
   
        return Distance(p0, p1);
    }
    else if(Dot(p1-p2, p0-p2) < 0)
    {
   
        return Distance(p0, p2);
    }
    return (fabs(Det(p2-p1, p0-p1))/Distance(p1, p2));
}

判断一个点是不是在一个矩形里面

假设p1以及p2是矩形的左上角以及右下角

bool InRectangle(const Point &p0, const Point &p1, const Point &p2)
{
   
   return Dot(p2-p0, p1-p0) <= 0;
}

判断一个点是不是在一条线段上

注意：限制一个点的范围可以使用矩形框起来！！！

bool OnSegment(const Point &p0, const Point &p1, const Point &p2)
{
     
    return Det(p1-p0, p2-p0) == 0 && Dot(p1-p0, p2-p0) <= 0;
}

实际情况下，等于号应该是取一个精度

判断两条线段是否平行

bool Parallel(const Point &p1, const Point &p2, const Point &p3, const Point &p4)
{
   
    return Det(p2-p1, p4-p3) == 0;
}

判断两条线段是否相交

这里这一个思路比较巧妙，见代码

p1,p2在线段p3p4的两端，并且p3,p4在线段p1p2的两端
点在线上

bool SegIntersect(const Point &p1, const Point &p2, const Point &p3, const Point &p4)
{
   
    if(Det(p3-p1, p2-p1)*Det(p4-p1, p2-p1) < 0 && Det(p1-p4, p3-p4)*Det(p2-p4, p3-p4)<0)
        return true;
    else if(Det(p3-p1, p2-p1) == 0 && OnSegment(p3, p1, p2))
        return true;
    else if(Det(p4-p1, p2-p1) == 0 && OnSegment(p4, p1, p2))
        return true;
    else if(Det(p1-p4, p3-p4) == 0 && OnSegment(p1, p2, p4))
        return true;
    else if(Det(p2-p4, p3-p4) == 0 && OnSegment(p2, p3, p4))
        return true;
    else
        return false;
}

判断一个点是否在多边形内

简单多边形：所有的边不相交。（这里讨论的情况就是简单多边形）

要求：判断p0是否在a[0...n]中（a[0] == a[n]）[包含边界]

我们需要在这一个点做一条向右的射线，如果与多边形的边有奇数个交点，那么这一个点就在多边形的内部，否则就在外部。

教课书是一般不会错的，每一个异样处都有原因！

double x = (p0.y - p1.y) * (p2.x - p1.x) / (p2.y - p1.y) + p1.x;

这是算出p0向右发出的射线与多边形的这一条边的交点

bool PointInPolygon(Point p0, vector<Point> a)
{
   
    int cnt = 0;
    for(int i = 0; i < a.size()-1; i++)//我改了一下，感觉书上有越界现象
    {
   
        Point p1 = a[i], p2 = a[i+1];
        if(OnSegment(p0, p1, p2)) return true;//如果在多边形的边上
        if(p1.y == p2.y) continue;//跳过水平线
        if(p0.y < p1.y && p0.y < p2.y) continue;
        if(p0.y >= p1.y && p0.y >= p2.y) continue;//这里有玄机，等于号不要删除
        double x = (p0.y - p1.y) * (p2.x - p1.x) / (p2.y - p1.y) + p1.x;
        if(x > p0.x) cnt++;
    }
    return (cnt & 1);
}

对于“玄机”的解释：

先上一段代码：

bool PointInPolygon(Point p0, vector<Point> a)
{
    
    int cnt = 0;
    for(int i = 0; i < a.size()-1; i++)
    {
    
        Point p1 = a[i], p2 = a[i+1];
        if(OnSegment(p0, p1, p2)) return true;
        if(p1.y == p2.y) continue;
        if(p0.y < p1.y && p0.y < p2.y) continue;
        if(p0.y > p1.y && p0.y > p2.y) continue;//等于号已经删除！！！！！
        double x = (p0.y - p1.y) * (p2.x - p1.x) / (p2.y - p1.y) + p1.x;
        if(x > p0.x) cnt++;
    }
    return (cnt & 1);
}


int main()
{
    
    vector<Point> a;
    a.push_back(Point(-2, 0));
    a.push_back(Point(0, 2));
    a.push_back(Point(2, 0));
    a.push_back(Point(1, 0));
    a.push_back(Point(0, -1));//更改后面的坐标为-1和1，观察结果
    a.push_back(Point(-1, 0));
    a.push_back(Point(-2, 0));
    bool ret = PointInPolygon(Point(0, 0), a);
    cout << ret;
    return 0;
}

正是由于如果一个

三个点构成的三角形的面积

利用向量的叉积处理

double TriangleArea(const Point &p0, const Point &p1, const Point &p2)
{
   
    return fabs(Det(p1-p0, p2-p0))/2;
}

求一个多边形的面积

这里把这一个多边形分成多块三角形，然后分别计算面积。

注意：采用这一种方法可以有效地处理凹形的问题。

double PolyArea(vector<Point> p)
{
   
    double ans = 0.0;
    for(int i = 1; i < p.size()-1; i++)
    {
   
        ans += Det(p[i] - p[0], p[i+1] - p[0]);
        //不可以使用 TriangleArea(p[0], p[i], p[i+1]);，因为三角形的面积是带有绝对值的
    }
    return fabs(ans) / 2;
}

最终代码

#include <bits/stdc++.h>
using namespace std;
#include <bits/stdc++.h>
using namespace std;
class Point{
   
public:
    double x, y;
    Point(){
   
        x = 0, y = 0;
    }//默认构造函数
    Point(double x1, double y1){
   //重载构造函数
        x = x1;
        y = y1;
    }
    void disp()
    {
   
        printf("(%g, %g)", x, y);
    }
};
//加运算
Point operator +(const Point &p1, const Point &p2)
{
   
    return Point(p1.x + p2.x, p1.y + p2.y);
}

//减运算
Point operator -(const Point &p1, const Point &p2)
{
   
    return Point(p1.x - p2.x, p1.y - p2.y);
}

//相等
bool operator ==(const Point &p1, const Point &p2)
{
   
    if(p1.x == p2.x && p1.y == p2.y)
        return true;
    else 
        return false;
}

//点积运算
double Dot(const Point &p1,const Point &p2)
{
   
    return p1.x*p2.x + p1.y * p2.y;
}
/* 通过点积可以判断两个向量的夹角 */

//求向量的模
double Length(const Point &p)
{
   
    return sqrt(Dot(p, p));
}

//对于具有相同起点的两个向量的夹角的计算
int Angle(const Point &p0, const Point &p1, const Point &p2)
{
   
    double t = Dot((p1-p0), (p2-p0));
    if(t > 1e-8) return 1;//锐角
    else if(t < -(1e-8)) return -1;//钝角
    else return 0;//直角
}

//pay attention to：
//如果是一个二维向量，那么叉积为一个标量
double Det(const Point &p1, const Point p2)
{
   
    return p1.x * p2.y - p1.y * p2.x;
}

//两点之间的距离（欧式距离）
double Distance(const Point &p1, const Point &p2)
{
   
    return Length(p2-p1);
}

//点到线段的距离
double DistPtoSegment(const Point &p0, const Point &p1, const Point &p2)
{
   
    if(p1 == p2)
    {
   
        return Distance(p1, p0);
    }
    if(Dot((p0 - p1), (p2-p1)) < 0)
    {
   
        return Distance(p0, p1);
    }
    else if(Dot(p1-p2, p0-p2) < 0)
    {
   
        return Distance(p0, p2);
    }
    return (fabs(Det(p2-p1, p0-p1))/Distance(p1, p2));
}

//判断一个点是不是在一个矩形里面,假设`p1`以及`p2`是矩形的左上角以及右下角
bool InRectangle(const Point &p0, const Point &p1, const Point &p2)
{
   
   return Dot(p2-p0, p1-p0) <= 0;
}

//判断一个点是不是在一条线段上
bool OnSegment(const Point &p0, const Point &p1, const Point &p2)
{
     
    return Det(p1-p0, p2-p0) == 0 && Dot(p1-p0, p2-p0) <= 0;
}

//判断两条线段是否平行
bool Parallel(const Point &p1, const Point &p2, const Point &p3, const Point &p4)
{
   
    return Det(p2-p1, p4-p3) == 0;
}

//判断两条线段是否相交
bool SegIntersect(const Point &p1, const Point &p2, const Point &p3, const Point &p4)
{
   
    if(Det(p3-p1, p2-p1)*Det(p4-p1, p2-p1) < 0 && Det(p1-p4, p3-p4)*Det(p2-p4, p3-p4)<0)
        return true;
    else if(Det(p3-p1, p2-p1) == 0 && OnSegment(p3, p1, p2))
        return true;
    else if(Det(p4-p1, p2-p1) == 0 && OnSegment(p4, p1, p2))
        return true;
    else if(Det(p1-p4, p3-p4) == 0 && OnSegment(p1, p2, p4))
        return true;
    else if(Det(p2-p4, p3-p4) == 0 && OnSegment(p2, p3, p4))
        return true;
    else
        return false;
}

//判断一个点是否在多边形内
bool PointInPolygon(Point p0, vector<Point> a)
{
   
    int cnt = 0;
    for(int i = 0; i < a.size()-1; i++)//我改了一下，感觉书上有越界现象
    {
   
        Point p1 = a[i], p2 = a[i+1];
        if(OnSegment(p0, p1, p2)) return true;//如果在多边形的边上
        if(p1.y == p2.y) continue;//跳过水平线
        if(p0.y < p1.y && p0.y < p2.y) continue;
        if(p0.y >= p1.y && p0.y >= p2.y) continue;//这里有玄机，等于号不要删除
        double x = (p0.y - p1.y) * (p2.x - p1.x) / (p2.y - p1.y) + p1.x;
        if(x > p0.x) cnt++;
    }
    return (cnt & 1);
}

//三个点构成的三角形的面积
double TriangleArea(const Point &p0, const Point &p1, const Point &p2)
{
   
    return fabs(Det(p1-p0, p2-p0))/2;
}


//求一个多边形的面积
double PolyArea(vector<Point> p)
{
   
    double ans = 0.0;
    for(int i = 1; i < p.size()-1; i++)
    {
   
        ans += Det(p[i] - p[0], p[i+1] - p[0]);
        //不可以使用 TriangleArea(p[0], p[i], p[i+1]);，因为三角形的面积是带有绝对值的
    }
    return fabs(ans) / 2;
}
int main()
{
   
    vector<Point> a;
    a.push_back(Point(-2, 0));
    a.push_back(Point(0, 2));
    a.push_back(Point(2, 0));
    a.push_back(Point(1, 0));
    a.push_back(Point(0, 1));//
    a.push_back(Point(-1, 0));
    a.push_back(Point(-2, 0));
    double ret = TriangleArea(a[0], a[1], a[2]);
    cout << ret;
    return 0;
}

【算法设计与分析 李春葆】计算几何（一）