计算几何总结

1、点积

a·b的几何意义为a在b上的投影长度乘以b的模, a·b=|a||b|cosθ,其中θ为a,b之间的夹角

vector operator + (vector a,vector b)  {return vector (a.x+b.x,a.y+b.y); }
vector operator - (vector a,vector b)  {return vector (a.x-b.x,a.y-b.y);}
vector operator * (vector a,double p)  {return vector (a.x*p,a.y*p);}
vector operator / (vector a.double p)  {return vector (a.x/p,a.y/p);}

（1）判断两个向量是否垂直 a⊥b <=> a·b=0
（2）求两个向量的夹角，点积<0为钝角，点积>0为锐角

double dot(vector a,vector b)
{
    return a.x*b.x+a.y*b.y;
}

double Angle(vector a,vector b)
{
    return acos(dot(a,b)/len(a)/len(b));
}

求模长

double len(vector a)
{ 
    return sqrt(dot(a,a));
}

vector normal(vector a)
{ 
    double l=len(a);
    return vector (-a.y/l,a.x/l);
}

2、二维叉积

设矢量P = （x1,y1） ，Q = (x2,y2)
则矢量叉积定义为： P × Q = x1y2 -x2y1 得到的是一个标量
显然有性质 P × Q = - ( Q × P) P × ( - Q ) = - ( P × Q )
如不加说明，下面所有的点都看作矢量，点的乘法看作矢量叉积；
叉乘的重要性质：
> 若 P × Q > 0 , 则P 在Q的顺时针方向
> 若 P × Q < 0 , 则P 在Q的逆时针方向
> 若 P × Q = 0 , 则P 与Q共线，但可能同向也可能反向

二维叉积是一个标量，a×b的几何意义为他们所形成的平行四边形的有向面积
直观理解，假如b在a的左边，则有向面积为正，假如在右边则为负。假如b,a共线，则叉积为0, 所以叉积可以用来判断平行。

image.png

double cross(vector a,vector b)
{
    return a.x*b.y-a.y*b.x;
}

三维叉积是一个向量，更为熟知的叫法是法向量，该向量垂直于a和b向量构成的平面。

3、点到直线的距离

利用叉积求面积，然后除以平行四边形的底边长，得到平行四边形的高即点到直线的距离

double distl(point p,point a,point b)
{
    vector v=p-a; vector u=b-a;
    return fabs(cross(v,u))/len(u);
}

4、点到线段的距离

比点到直线的距离稍微复杂。因为是线段，所以如果平行四边形的高在区域之外的话就不合理，这时候需要计算点到距离较近的端点的距离

double dists(point p,point a,point b)
{
    if (a==b) return len(p-a);
    vector v1=b-a,v2=p-a,v3=p-b;
    if (dcmp(dot(v1,v2))<0) return len(v2);
    else if (dcmp(dot(v1,v3))>0) return len(v3);
    return fabs(cross(v1,v2))/len(v1);
}

5、判断两个线段是否相交

跨立实验：判断一条线段的两端是否在另一条线段的两侧（两个端点与另一线段的叉积乘积为负）。需要正反判断两侧。

bool segment(point a,point b,point c,point d)
{
    double c1=cross(b-a,c-a),c2=cross(b-a,d-a);
    double d1=cross(d-c,a-c),d2=cross(d-c,b-c);
    return dcmp(c1)*dcmp(c2)<0&&dcmp(d1)*dcmp(d2)<0;
}

6、求两条直线的交点

有两种方法，比较常用的一种是用叉积的比值计算。但是这种方法的精度不是很高。

point glt(point a,point a1,point b,point b1)  
{  
    vector v=a1-a; vector w=b1-b;  
    vector u=a-b;  
    double t=cross(w,u)/cross(v,w);  
    return a+v*t;  
}

还有一种比较麻烦，不常用但是精度相对较好

point line_intersection(point a,point a0,point b,point b0)  
{  
    double a1,b1,c1,a2,b2,c2;  
    a1 = a.y - a0.y;  
    b1 = a0.x - a.x;  
    c1 = cross(a,a0);  
    a2 = b.y - b0.y;  
    b2 = b0.x - b.x;  
    c2 = cross(b,b0);  
    double d = a1 * b2 - a2 * b1;  
    return point((b1 * c2 - b2 * c1) / d,(c1 * a2 - c2 * a1) / d);  
} 

7、求垂足

def getFootPoint(point, line_p1, line_p2):
    """
    @point, line_p1, line_p2 : [x, y, z]
    """
    x0 = point[0]
    y0 = point[1]
    z0 = point[2]

    x1 = line_p1[0]
    y1 = line_p1[1]
    z1 = line_p1[2]

    x2 = line_p2[0]
    y2 = line_p2[1]
    z2 = line_p2[2]

    k = -((x1 - x0) * (x2 - x1) + (y1 - y0) * (y2 - y1) + (z1 - z0) * (z2 - z1)) / \
        ((x2 - x1) ** 2 + (y2 - y1) ** 2 + (z2 - z1) ** 2)*1.0

    xn = k * (x2 - x1) + x1
    yn = k * (y2 - y1) + y1
    zn = k * (z2 - z1) + z1

    return (xn, yn, zn)