计算几何（一）——基础算法

待续 《计算几何资料》为提纲

1.

(1) 有向线段 ： 如果一条线段的端点有次序之分，我们把这种线段称为有向线段(directed segment)。

(2) 矢量 ：如果有向线段p1p2的起点p1在坐标原点，则称之为矢量(vector) p2

(3) (二维)矢量叉乘 ： P=(x1,y1),  Q=(x2,y2) ==> P×Q

    随便的定义： 大小|x1  y1|， 方向右手法则

                                  |x2  y2|

    严格定义：(    将二维扩展到三维 P=(x1,y1,0), Q=(x2,y2,0)    )

                   P×Q=|i     j     k|=k*|x1  y1|

                              |x1  y1  0|      |x2  y2|

                              |x2  y2  0|

    P×Q>0，则P在Q的顺时针方向（指：Q沿顺时针方向相比沿逆时针方向更早到达P）

    P×Q<0，则P在Q的逆时针方向

    P×Q=0，则P和Q共线（同向或反向）

 

2. 折线段拐向 ：利用二维矢量叉乘性质

 

3. “点”与“线段”关系

    设点Q，线段P1P2.

    点Q在线段上：共线(Q-P1)×(P2-P1) 且在P1,P2为对角顶点矩形内

 

2、3 代码：

#include<iostream>
using namespace std;

struct Point{
       double x;
       double y;
       Point(double a=0, double b=0){x=a;y=b;}
};

struct LineSeg{
       Point s;
       Point e;
       LineSeg(){}
       LineSeg(Point a, Point b){s=a;e=b;}
};

/*有向线段 dl<op,sp>×dl<op,ep>*/
double multiply(Point sp,Point ep,Point op) { 
  return((sp.x-op.x)*(ep.y-op.y)-(ep.x-op.x)*(sp.y-op.y)); 
}

/*点在线段上*/
bool ptOnLine(LineSeg l, Point p){
  return ( multiply(l.e, p, l.s)==0 &&
         (p.x-l.s.x)*(p.x-l.e.x)<=0 &&
         (p.y-l.s.y)*(p.y-l.e.y)<=0 );
}

int main(void){
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(2,2)), *new Point(1.5,1.5) )<<endl; //1
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(2,2)), *new Point(3,3) )<<endl;     //0
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(2,2)), *new Point(2,3) )<<endl;     //0
    
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(2,1)), *new Point(1.5,1) )<<endl;   //1
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(2,1)), *new Point(3,1) )<<endl;     //0
    
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(1,2)), *new Point(1,1.5) )<<endl;   //1
    cout<<ptOnLine( LineSeg(*new Point(1,1), *new Point(1,2)), *new Point(1,3) )<<endl;     //0
    
    system("pause");
    return 0;
}

 

12. “点”与“多边形”关系

(1)射线法——适用于任何多边形（凹凸），C:\Users\Administrator\Desktop\Test\计算几何\PointInPoly.cpp

    过此点向任意角度发一条射线，若与多边形的各条边交点个数之和为偶数，则此点在 多边形之外，否则在多边形之内。
    若有交点为多边形顶点则要另选一条射线重算

射线法代码：

#include<iostream>
#include<cmath>
using namespace std;

class CPoint{
public:
       CPoint(){this->x=x;this->y=y;}
       CPoint(double x, double y){
         this->x=x;
         this->y=y;
       }
       double x;
       double y;
};

//向右射线与线段关系：相交返回1，不相交返回0，在边上返回-1(-1时认为点在多边形内)
int IsIntersectant( CPoint ptStart, CPoint ptEnd, CPoint pd ){
  double tempx = 0;
  double tempy = 0;

  double startx = ptStart.x;   
  double starty = ptStart.y;
  double endx = ptEnd.x;   
  double endy = ptEnd.y;
 
  //横标范围[minX, maxX] 
  double maxX = startx;
  double minX = endx;
  if (startx<endx){
    minX = startx;
    maxX = endx;
  }
 
  //保证开始点在下(y坐标大)，作为下端点
  if(starty<endy){
    double temp=starty;
    starty=endy;
    endy=temp;
   
    temp=startx;
    startx=endx;
    endx=temp;
  }
 
  //1. 没有交点 （特例一） 
  if((pd.y>starty && pd.y>endy) || (pd.y<starty && pd.y<endy)){
    //水平射线在该直线段的上下两端之外
    return 0;
  }
  if(pd.x>startx && pd.x>endx){
    //水平射线的起点在该直线段的右边
    return 0;
  }
 
  //2. 对于水平线段 : 在线段上返回-1，否则返回0 （特例二） 
  if (starty ==endy)
  {
    if (pd.x<=maxX&&pd.x>=minX)
      return -1;
    return 0;
  }
 
  //3. 非水平线段 
  //   3.1 (向右射线端点)pd在线段上: 通过三角形三边关系判断 
  tempx = pd.x - startx;
  tempy = pd.y - starty;
  double dStartToX = sqrt(tempx*tempx+tempy*tempy);
 
  tempx = pd.x - endx;
  tempy = pd.y - endy;
  double dXToEnd = sqrt(tempx*tempx+tempy*tempy);
 
  tempx = startx - endx;
  tempy = starty - endy;
  double dStarToEnd = sqrt(tempx*tempx+tempy*tempy);
 
  if (dStarToEnd == (dStartToX + dXToEnd)){
    return -1;
  }
  
  //   3.2 向右射线和线段不相交 
  //h_x记录点pd(x,y)引水平线与直线段所在直线的交点x坐标
  //注：直线两点式公式 (y-y1)/(y2-y1)=(x-x1)/(x2-x1),其中 x1≠x2, y1≠y2
  //          也可写作 y-y1=(y2-y1)/(x2-x1)*(x-x1),其中x1≠x2 
  double h_x=(pd.y-starty)*(endx-startx)/(endy-starty)+startx; //只需保证y1≠y2  
  if(pd.x>h_x){//pd点在交点的右边，则射线与线段不相交
    return 0;
  }
  
  //   3.3 向右射线和线段相交
  return 1;
}
 
/**点是否在多边形内 (射线法)*/
bool PtInPolygon( CPoint *ptList, int ptCount, CPoint pd ){
  if (ptCount<3){
    cout<<"注：多边形顶点至少有3个"<<endl; 
    return false;
  }
  
  int cross_num = 0;
  int iFlag = 0;
  //从点pd，向右引水平射线
 
  //除了末端点与首端点连接线，其他各条边的线段 -->与水平射线的交点和 cross_num
  for(int i=0;i<ptCount-1;i++)
  {
    iFlag = IsIntersectant(ptList[i], ptList[i+1], pd);
    if (iFlag < 0){
      return true;
    }else{
     cross_num += iFlag;
    }
  }
  //末端点与首端点连接线 --> 与水平射线的交点和 cross_num
  iFlag = IsIntersectant(ptList[ptCount-1], ptList[0], pd);
  if(iFlag < 0){
    return true;
  }else{
    cross_num += iFlag;
  }
 
  if(cross_num%2==1){
    //1. 交点个数为奇，则点在多边形内 
    return true;
  }else if(cross_num%2==0){
    //2. 交点个数为偶，则点不在多边形内 
    return false;
  }
  return false;
}

int main(void){   
    //如果先声明 CPoint ptList[5]; 必须提供无参构造函数 
    CPoint ptList[5]={CPoint(0,0),CPoint(2,0),CPoint(2,2),CPoint(1,1),CPoint(0,2)};
    int ptNum=5;
    
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(1,1))<<endl;     //1
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(0.5,1))<<endl;   //1
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(1.5,1.1))<<endl; //1
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(1,0))<<endl;     //1
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(1,1.5))<<endl;   //0
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(2.5,3))<<endl;   //0
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(-1,-1))<<endl;   //0
    cout<<PtInPolygon(ptList,ptNum,*new CPoint(2,2))<<endl;     //1
    
    system("pause");
    return 0;
}

 

(2)叉乘判别法——适用于凸多边形，见 http://www.cppblog.com/w2001/archive/2011/06/17/31694.html

     想象一个凸多边形，其每一个边都将整个2D屏幕划分成为左右两边，连接每一边的第一个端点和要测试的点得到一 个矢量v，将两个2维矢量扩展成3维的，然后将该边与v叉乘，判断结果3维矢量中Z分量的符号是否发生变化，进而推导出点是否处于凸多边形内外。这里要注 意的是，多边形顶点究竟是左手序还是右手序，这对具体判断方式有影响。

 

(3)累积角度法——适用于任何多边形（凹凸）：
    过此点连接多边形的每一顶点，各相邻边角度（±）之和为360度 ，则此点在多边形内。否则为0度 ，在多边形外部。