Point in convex polygon

Question

I'm trying to make a function that will return True if the given (x,y) point is inside a convex polygon. I'm trying to make it without numpy or any similar imports, just pure python code.

I've already found a sample solution, which seems OK at first sight, but it's not working correctly, and I can't figure out why. The code is as follows:

def point_in_poly(x,y,poly):

  n = len(poly)
  inside = False

  p1x,p1y = poly[0]
  for i in range(n+1):
      p2x,p2y = poly[i % n]
      if y > min(p1y,p2y):
          if y <= max(p1y,p2y):
              if x <= max(p1x,p2x):
                  if p1y != p2y:
                      xints = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x
                  if p1x == p2x or x <= xints:
                      inside = not inside
      p1x,p1y = p2x,p2y

  return inside

If I test it for (9,9), for the following polygon, it gives me True:

polygon = [(0,10),(10,10),(10,0),(0,0)]
point_x = 9
point_y = 9
print point_in_poly(point_x,point_y,polygon)

But when I change the order of the points of the polygon, for the same point, it gives me False:

polygon = [(0,0), (0,10), (10,0), (10,10)]
point_x = 9
point_y = 9
print point_in_poly(point_x,point_y,polygon)

Anybody knows the reason? Thanks!

Possibly relevant: http://stackoverflow.com/questions/217578/point-in-polygon-aka-hit-test?rq=1 — Hannele, May 01 '13 at 20:35

GeospatialPython.com · Accepted Answer · 2013-05-06T02:57:39.307

In the particular case you are having problems with is special: polygon = [(0,0), (0,10), (10,0), (10,10)]

Changing the order of points in a polygon can have significant impact on algorithms.

If you draw your polygon on a graph you'll see you have a horizontal hourglass shape. The polygon border overlaps itself. In geospatial analysis this overlap is not allowed because visually and logically you now have two closed polygons with a common intersection point. By the way most geospatial software doesn't deal well with triangles either.

In this case the point at 9,9 will trick the ray casting algorithm used in your method above because it can easily cross the doubled-over polygon boundary twice.

Please run the following code to see what is going on. (9,9) is on the line and this algorithm doesn't account for it. (5,8) is way outside:

import turtle as t
polygon = [(0,0), (0,100), (100,0), (100,100)]
t.goto(0,0)
fp = None
for p in polygon:
  t.goto(p)
  if not fp: fp=p
t.goto(fp)
t.up()
t.goto(90,90)
t.write("90,90")
t.dot(10)
t.goto(50,80)
t.write("50,80")
t.dot(10)
t.done()

This code handles the (9,9) edge case: http://geospatialpython.com/2011/08/point-in-polygon-2-on-line.html

The code above draws this image.

Well, you have a point. But I tried the 'hourglass'-shape with (5,8), it should be in the hourglass, because it's above the crossing of the two lines, inside the upper triangle. Result: also says "False"... :( — John Berry, May 05 '13 at 21:31
Wow, you know your job! :) Now I understand how it works... Thank you! And big thanks for the picture! — John Berry, May 08 '13 at 19:29
The link to code that handles this edge case appears broken. Please fix. — martineau, May 27 '16 at 14:17

score 0 · Answer 2 · answered May 01 '13 at 20:25

0

the point 9,0 is not inside the polygon [(0,10),(10,10),(10,0),(0,0)] its on the edge. Points exactly on the edge can be considered in or out depending on the specifics of your algorithm.

answered May 01 '13 at 20:25

cmd

5,754
16
30

OK, but if I try with (9,9), and modify the order of the polygon-points: FROM [(0,10),(10,10),(10,0),(0,0)] TO [(0,0),(0,10),(10,0),(10,10)] It's still not working, only with the given order.... – John Berry May 01 '13 at 20:47

score 0 · Answer 3 · answered May 01 '13 at 21:21

For polygons with lots of points it often pays to check if the point falls outside of the bounding box first:

def point_in_poly(x, y, poly):
    # Check bounding box first
    xl = [p[0] for p in poly]
    yl = [p[1] for p in poly]
    if x < min(xl) or x > max(xl) or y < min(yl) or y > max(yl):
       return False
   # Now check all points.

You can find a lot of algoriths dealing with polygons and meshes at one of Paul Bourke's webpages.

For a lot of these algorithms it pays to use numpy, though, since a lot of the steps have to be done for each point in the poly array.

Point in convex polygon

3 Answers3