c 点在多边形内python 点多边形算法




点在多边形内python (10)

我看到下面的算法用于检查点是否在此link的给定多边形中:

int pnpoly(int nvert, float *vertx, float *verty, float testx, float testy)
{
  int i, j, c = 0;
  for (i = 0, j = nvert-1; i < nvert; j = i++) {
    if ( ((verty[i]>testy) != (verty[j]>testy)) &&
     (testx < (vertx[j]-vertx[i]) * (testy-verty[i]) / (verty[j]-verty[i]) + vertx[i]) )
       c = !c;
  }
  return c;
}

我试过这个算法,它实际上工作得很完美。 但遗憾的是,在花了一些时间试图了解它之后,我无法理解它。

因此,如果有人能够理解这个算法,请向我解释一下。

谢谢。


只要多边形的边不交叉,此算法就可以在任何闭合多边形中工作。 三角形,五边形,正方形,甚至是一个非常弯曲的分段线性橡皮筋,它不会穿过自身。

1)将多边形定义为有向矢量组。 这意味着多边形的每一边都是由从顶点a到顶点+ 1的矢量描述的。 向量是如此定向的,以便一个人的头部接触下一个的尾部,直到最后一个向量接触第一个向量的尾部。

2)选择要测试多边形内部或外部的点。

3)对于沿着多边形的周长的每个向量Vn,找到向量Dn,其在测试点上开始并且在Vn的尾部结束。 计算定义为DnXVn / DN * VN的向量Cn(X表示叉积; *表示点积)。 用名称Mn调用Cn的大小。

4)添加所有Mn并称此数量为K.

5)如果K为零,则该点在多边形之外。

6)如果K不为零,则该点在多边形内。

理论上,位于多边形边缘的点将产生未定义的结果。

K的几何意义是坐在我们的测试点上的跳蚤“看到”走在多边形行走边缘的蚂蚁行走到左边减去向右走的角度的总角度。 在闭路中,蚂蚁在它开始的地方结束。 在多边形之外,无论位置如何,答案都是零。
在多边形内部,无论位置如何,答案都是“一次点”。


要扩展@ chowlette的answer ,其中第二行检查点是否在该行的左侧,没有给出推导但这是我的计算出来的:首先它有助于想象2个基本情况:

  • 这一点就在线的左边. / . /
  • 这一点是正确的/ .

如果我们的观点是水平射出一条射线,它会撞击线段。 我们指向它的左侧还是右侧? 在里面还是外面? 我们知道它的y坐标,因为它的定义与点相同。 x坐标是什么?

采用传统的线公式y = mx + b 。 m是跑步的上升。 在这里,我们试图找到该线段上与我们的点具有相同高度(y)的点x坐标

所以我们求解x: x = (y - b)/mm在运行中上升,因此变为上升或(yj - yi)/(xj - xi)变为(xj - xi)/(yj - yi) 。 b是与原点的偏移量。 如果我们假设yi作为我们坐标系的基础,则b变为yi。 我们的点测试是我们的输入,减去yi将整个公式变成yi的偏移量。

我们现在有(xj - xi)/(yj - yi)1/m (testy - yi) y或(testy - yi) :( (xj - xi)(testy - yi)/(yj - yi)但是testx不是基于yi所以我们将它添加回来以便比较两者(或者也是零testx)


我更改了here ,使其更具可读性(这也使用了Eigen)。 算法是相同的。

// This uses the ray-casting algorithm to decide whether the point is inside
// the given polygon. See https://en.wikipedia.org/wiki/Point_in_polygon#Ray_casting_algorithm
bool pnpoly(const Eigen::MatrixX2d &poly, float x, float y)
{
    // If we never cross any lines we're inside.
    bool inside = false;

    // Loop through all the edges.
    for (int i = 0; i < poly.rows(); ++i)
    {
        // i is the index of the first vertex, j is the next one.
        // The original code uses a too-clever trick for this.
        int j = (i + 1) % poly.rows();

        // The vertices of the edge we are checking.
        double xp0 = poly(i, 0);
        double yp0 = poly(i, 1);
        double xp1 = poly(j, 0);
        double yp1 = poly(j, 1);

        // Check whether the edge intersects a line from (-inf,y) to (x,y).

        // First check if the line crosses the horizontal line at y in either direction.
        if ((yp0 <= y) && (yp1 > y) || (yp1 <= y) && (yp0 > y))
        {
            // If so, get the point where it crosses that line. This is a simple solution
            // to a linear equation. Note that we can't get a division by zero here -
            // if yp1 == yp0 then the above if be false.
            double cross = (xp1 - xp0) * (y - yp0) / (yp1 - yp0) + xp0;

            // Finally check if it crosses to the left of our test point. You could equally
            // do right and it should give the same result.
            if (cross < x)
                inside = !inside;
        }
    }
    return inside;
}

这可能与解释实际代码中的光线跟踪算法一样详细。 它可能没有被优化,但必须始终在完全掌握系统之后。

    //method to check if a Coordinate is located in a polygon
public boolean checkIsInPolygon(ArrayList<Coordinate> poly){
    //this method uses the ray tracing algorithm to determine if the point is in the polygon
    int nPoints=poly.size();
    int j=-999;
    int i=-999;
    boolean locatedInPolygon=false;
    for(i=0;i<(nPoints);i++){
        //repeat loop for all sets of points
        if(i==(nPoints-1)){
            //if i is the last vertex, let j be the first vertex
            j= 0;
        }else{
            //for all-else, let j=(i+1)th vertex
            j=i+1;
        }

        float vertY_i= (float)poly.get(i).getY();
        float vertX_i= (float)poly.get(i).getX();
        float vertY_j= (float)poly.get(j).getY();
        float vertX_j= (float)poly.get(j).getX();
        float testX  = (float)this.getX();
        float testY  = (float)this.getY();

        // following statement checks if testPoint.Y is below Y-coord of i-th vertex
        boolean belowLowY=vertY_i>testY;
        // following statement checks if testPoint.Y is below Y-coord of i+1-th vertex
        boolean belowHighY=vertY_j>testY;

        /* following statement is true if testPoint.Y satisfies either (only one is possible) 
        -->(i).Y < testPoint.Y < (i+1).Y        OR  
        -->(i).Y > testPoint.Y > (i+1).Y

        (Note)
        Both of the conditions indicate that a point is located within the edges of the Y-th coordinate
        of the (i)-th and the (i+1)- th vertices of the polygon. If neither of the above
        conditions is satisfied, then it is assured that a semi-infinite horizontal line draw 
        to the right from the testpoint will NOT cross the line that connects vertices i and i+1 
        of the polygon
        */
        boolean withinYsEdges= belowLowY != belowHighY;

        if( withinYsEdges){
            // this is the slope of the line that connects vertices i and i+1 of the polygon
            float slopeOfLine   = ( vertX_j-vertX_i )/ (vertY_j-vertY_i) ;

            // this looks up the x-coord of a point lying on the above line, given its y-coord
            float pointOnLine   = ( slopeOfLine* (testY - vertY_i) )+vertX_i;

            //checks to see if x-coord of testPoint is smaller than the point on the line with the same y-coord
            boolean isLeftToLine= testX < pointOnLine;

            if(isLeftToLine){
                //this statement changes true to false (and vice-versa)
                locatedInPolygon= !locatedInPolygon;
            }//end if (isLeftToLine)
        }//end if (withinYsEdges
    }

    return locatedInPolygon;
}

关于优化的一句话:最短(和/或最简洁)的代码实现最快是不正确的。 从数组中读取和存储元素并在代码块的执行中(可能)多次使用它比在每次需要时访问数组要快得多。 如果阵列非常大,这尤其重要。 在我看来,通过将一个数组的每个术语存储在一个命名良好的变量中,它也更容易评估其目的,从而形成更易读的代码。 只是我的两分钱......


这是一个php的实现:

<?php
class Point2D {

    public $x;
    public $y;

    function __construct($x, $y) {
        $this->x = $x;
        $this->y = $y;
    }

    function x() {
        return $this->x;
    }

    function y() {
        return $this->y;
    }

}

class Point {

    protected $vertices;

    function __construct($vertices) {

        $this->vertices = $vertices;
    }

    //Determines if the specified point is within the polygon. 
    function pointInPolygon($point) {
        /* @var $point Point2D */
    $poly_vertices = $this->vertices;
    $num_of_vertices = count($poly_vertices);

    $edge_error = 1.192092896e-07;
    $r = false;

    for ($i = 0, $j = $num_of_vertices - 1; $i < $num_of_vertices; $j = $i++) {
        /* @var $current_vertex_i Point2D */
        /* @var $current_vertex_j Point2D */
        $current_vertex_i = $poly_vertices[$i];
        $current_vertex_j = $poly_vertices[$j];

        if (abs($current_vertex_i->y - $current_vertex_j->y) <= $edge_error && abs($current_vertex_j->y - $point->y) <= $edge_error && ($current_vertex_i->x >= $point->x) != ($current_vertex_j->x >= $point->x)) {
            return true;
        }

        if ($current_vertex_i->y > $point->y != $current_vertex_j->y > $point->y) {
            $c = ($current_vertex_j->x - $current_vertex_i->x) * ($point->y - $current_vertex_i->y) / ($current_vertex_j->y - $current_vertex_i->y) + $current_vertex_i->x;

            if (abs($point->x - $c) <= $edge_error) {
                return true;
            }

            if ($point->x < $c) {
                $r = !$r;
            }
        }
    }

    return $r;
}

测试运行:

        <?php
        $vertices = array();

        array_push($vertices, new Point2D(120, 40));
        array_push($vertices, new Point2D(260, 40));
        array_push($vertices, new Point2D(45, 170));
        array_push($vertices, new Point2D(335, 170));
        array_push($vertices, new Point2D(120, 300));
        array_push($vertices, new Point2D(260, 300));


        $Point = new Point($vertices);
        $point_to_find = new Point2D(190, 170);
        $isPointInPolygon = $Point->pointInPolygon($point_to_find);
        echo $isPointInPolygon;
        var_dump($isPointInPolygon);

我修改了代码以检查点是否在多边形中,包括点在边缘上。

bool point_in_polygon_check_edge(const vec<double, 2>& v, vec<double, 2> polygon[], int point_count, double edge_error = 1.192092896e-07f)
{
    const static int x = 0;
    const static int y = 1;
    int i, j;
    bool r = false;
    for (i = 0, j = point_count - 1; i < point_count; j = i++)
    {
        const vec<double, 2>& pi = polygon[i);
        const vec<double, 2>& pj = polygon[j];
        if (fabs(pi[y] - pj[y]) <= edge_error && fabs(pj[y] - v[y]) <= edge_error && (pi[x] >= v[x]) != (pj[x] >= v[x]))
        {
            return true;
        }

        if ((pi[y] > v[y]) != (pj[y] > v[y]))
        {
            double c = (pj[x] - pi[x]) * (v[y] - pi[y]) / (pj[y] - pi[y]) + pi[x];
            if (fabs(v[x] - c) <= edge_error)
            {
                return true;
            }
            if (v[x] < c)
            {
                r = !r;
            }
        }
    }
    return r;
}

我认为基本思想是从点开始计算矢量,每个边缘一个。 如果矢量穿过一条边,则该点在多边形内。 通过凹多边形,如果它穿过奇数个边缘,它也在内部(免责声明:虽然不确定它是否适用于所有凹多边形)。


算法向右射线投射。 循环的每次迭代,都会根据多边形的一个边缘检查测试点。 如果点的y-coord在边缘范围内,则if测试的第一行成功。 第二行检查测试点是否在线的左侧(我想 - 我没有任何废纸可供检查)。 如果确实如此,则从测试点向右绘制的线穿过该边缘。

通过反复反转c的值,算法计算向右线穿过多边形的次数。 如果它经过了奇数次,则该点在内部; 如果是偶数,则该点在外面。

我会担心a)浮点运算的准确性,以及b)具有水平边缘或具有相同y坐标的测试点作为顶点的效果。


Chowlett在各方面,形状和形式上都是正确的。 该算法假定如果您的点位于多边形的线上,那么这是在外面 - 在某些情况下,这是错误的。 将两个'>'运算符更改为'> ='并将'<'更改为'<='将解决此问题。

bool PointInPolygon(Point point, Polygon polygon) {
  vector<Point> points = polygon.getPoints();
  int i, j, nvert = points.size();
  bool c = false;

  for(i = 0, j = nvert - 1; i < nvert; j = i++) {
    if( ( (points[i].y >= point.y ) != (points[j].y >= point.y) ) &&
        (point.x <= (points[j].x - points[i].x) * (point.y - points[i].y) / (points[j].y - points[i].y) + points[i].x)
      )
      c = !c;
  }

  return c;
}

该算法被分解为最必要的元素。 在开发和测试之后,所有不必要的东西都被删除了。 因此,您无法轻松应对,但它可以完成工作并且性能非常好。

我冒昧地把它翻译成ActionScript-3

// not optimized yet (nvert could be left out)
public static function pnpoly(nvert: int, vertx: Array, verty: Array, x: Number, y: Number): Boolean
{
    var i: int, j: int;
    var c: Boolean = false;
    for (i = 0, j = nvert - 1; i < nvert; j = i++)
    {
        if (((verty[i] > y) != (verty[j] > y)) && (x < (vertx[j] - vertx[i]) * (y - verty[i]) / (verty[j] - verty[i]) + vertx[i]))
            c = !c;
    }
    return c;
}






algorithm