第33章 ：计算几何学

一：判定给定的两条线段是否相交：

当两条线段的四个点p1,p2,p3,p4给定时，如下的代码能够判断线段p1p2与线段p3p4是否相交，如果相交则返回true，否则返回false;

//pair<double,double>第一个元素存储坐标点的x坐标，第二个元素存储坐标点的y坐标
using point=pair<double,double>;
// 计算矢量(pk-pi)与矢量(pj-pi)的叉乘
inline int direction(const point& pi,const point& pj,const point& pk)
{
        return (pk.first-pi.first)*(pj.second-pi.second)-(pj.first-pi.first)*(pk.second-pi.second);

}

//在矢量(pk-pi)与矢量(pj-pi)的叉乘为0的基础上，判断点pk是否在线段pipj上，如果在的话返回true，否则返回false;
inline bool on_segment(const point& pi,const point& pj,const point& pk)
{
        return min(pi.first,pj.first)<=pk.first&&pk.first<=max(pi.first,pj.first)&&min(pi.second,pj.second)<=pk.second&&pk.second<=max(pi.second,pj.second);
}

//判断线段p1p2是否与线段p3p4相交,相交返回true,不相交返回false;
bool segment_intersect(const point& p1,const point& p2,const point& p3,const point& p4)
{
        int d1=direction(p3,p4,p1);
        int d2=direction(p3,p4,p2);
        int d3=direction(p1,p2,p3);
        int d4=direction(p1,p2,p4);

        if(((d1>0&&d2<0)||(d1<0&&d2>0))&&((d3>0&&d4<0)||(d3<0&&d4>0)))
                return true;
        else if(d1==0&&on_segment(p3,p4,p1))
                return true;
        else if(d2==0&&on_segment(p3,p4,p2))
                return true;
        else if(d3==0&&on_segment(p1,p2,p3))
                return true;
        else if(d4==0&&on_segment(p1,p2,p4))
                return true;
        else
                return false;
}

二：寻找最近点对：

假设有一个xy平面上的点列P，点与点之间的距离定义为欧几里得距离，现在要求出一对最近的点。

1：表示点的坐标的代码如下：

//
struct point{
        double xval;
        double yval;
};

2： 计算点与点之间欧几里得距离：

inline double distance(point p1,point p2)
{
        return sqrt((p1.xval-p2.xval)*(p1.xval-p2.xval)+(p1.yval-p2.yval)*(p1.yval-p2.yval));
}

3：用暴力搜索法求出点对的最近距离：

double shortestDistanceEnumeration(vector<point>& xvalSort)
{
        double MinDistance=DBL_MAX;

        for(int i=0;i!=xvalSort.size();++i)
                for(int j=i+1;j!=xvalSort.size();++j)
                        if(distance(xvalSort[j],xvalSort[i])<MinDistance)
                                MinDistance=distance(xvalSort[j],xvalSort[i]);

        return MinDistance;
}

4：用分治算法求出点对的最近距离，输入xvalSort表示的是一系列以x坐标值从小到大排序的点，yvalSort表示的是一些列以y坐标值从小到大排序的点，当然这两个集合中点集是一样的，只是将它们排序的方法不一样，一个是按x值排序，另外一个是按y值排序：

double shortestDistance(vector<point>& xvalSort,vector<point>& yvalSort)
{
        if(xvalSort.size()<=10)
                return shortestDistanceEnumeration(xvalSort);

        vector<point> xvalSortLeft;
        vector<point> xvalSortRight;
        vector<point> yvalSortLeft;
        vector<point> yvalSortRight;

        double mediumXval=(xvalSort[0].xval+xvalSort[xvalSort.size()-1].xval)/2;
        for(int i=0;i!=xvalSort.size();++i)
        {
                if(xvalSort[i].xval<=mediumXval)
                        xvalSortLeft.push_back(xvalSort[i]);
                else
                        xvalSortRight.push_back(xvalSort[i]);
        }

        for(int i=0;i!=yvalSort.size();++i)
        {
                if(yvalSort[i].xval<=mediumXval)
                        yvalSortLeft.push_back(yvalSort[i]);
                else
                        yvalSortRight.push_back(yvalSort[i]);
        }

        if(xvalSortRight.size()==0)
                return shortestDistanceEnumeration(xvalSortLeft);
        else if(xvalSortRight.size()==1){
                xvalSortLeft.push_back(xvalSortRight[0]);
                yvalSortLeft.push_back(xvalSortRight[0]);
                return shortestDistanceEnumeration(xvalSortLeft);
        }

        double leftMin=shortestDistance(xvalSortLeft,yvalSortLeft);
        double rightMin=shortestDistance(xvalSortRight,yvalSortRight);

        double Min=leftMin<=rightMin?leftMin:rightMin;
        vector<point> pointsInStrip;
        for(int i=0;i!=yvalSort.size();++i)
                if((mediumXval-Min)<=yvalSort[i].xval&&yvalSort[i].xval<=(mediumXval+Min))

        const int lowestPointNumInStrip=2;
        if(pointsInStrip.size()>=lowestPointNumInStrip){
                for(int i=0;i!=pointsInStrip.size();++i)
                        for(int j=i+1;j!=pointsInStrip.size();++j)
                                if(Min<pointsInStrip[j].yval-pointsInStrip[i].yval)
                                        break;
                                else{
                                        if(distance(pointsInStrip[j],pointsInStrip[i])<Min)
                                                Min=distance(pointsInStrip[j],pointsInStrip[i]);
                                }
        }

        return Min;
}