Point in a polygon algorithm

Question

I am implementing a Point in polygon algorithm.

Inputs:

• M, N: size of the matrix
• poly: a list of tuples that represent the points of polygon.

Output:

A mask matrix which is ones everywhere, except the point in the polygon should be 0.

import numpy as np
import cv2
def getABC(x1, y1, x2, y2):
 A = y2 - y1
 B = x1 - x2
 C = A*x1 + B*y1
 return (A, B, C)
def polygon(M, N, poly):
 out = np.ones((M, N))*255
 n = len(poly)
 for i in range(M):
 intersection_x = i
 intersection_y = []
 # check through all edges
 for edge in range(n + 1):
 v1_x, v1_y = poly[edge % n]
 v2_x, v2_y = poly[(edge + 1) % n]
 A1, B1, C1 = getABC(v1_x, v1_y, v2_x, v2_y)
 A2 = 1
 B2 = 0
 C2 = i
 # find intersection
 det = A1*B2 - A2*B1
 if (det != 0):
 tmp = (A1 * C2 - A2 * C1)/det
 if tmp >= min(v1_y, v2_y) and tmp <= max(v1_y, v2_y):
 intersection_y.append(tmp)
 intersection_y = list(set(intersection_y))
 print intersection_y
 if len(intersection_y) == 1:
 intersection_y.append(intersection_y[0])
 for k in range(1, len(intersection_y), 2):
 out[intersection_x, intersection_y[k - 1]:intersection_y[k]] = 0
 return out
poly = [(10,20), (10,40), (30,20), (30,40)]
out = polygon(100,100, poly)
cv2.imwrite("out.png", out)

Is this the optimal algorithm? I use a scanning line, find the intersection with edges and make all the point between the intersections to be zero.

Answer 1

When you say "points of polygon", I am assuming you are referring to vertices.

I don't think what you have is the fastest. A simple improvement to this could be to divide your matrix in to a grid of p x p cells, where p is a parameter, and classify each grid-cell as completely inside or completely outside of the polygon. This can be done by checking if all the vertices of the grid-cell are inside the polygon. This way, large regions inside grid cells can be directly labeled as inside or outside, without needing to test every point inside the cell against all the edges of the polygon. For the grid-cells that are not completely inside the polygon, you just go by what you are doing. The optimal choice of p depends on the size of polygons you may encounter.

More details about the grid method here (with comparisons to other strategies): http://erich.realtimerendering.com/ptinpoly/

There could be better ones. I just suggested a simple improvement, which works well in practice.

• Soonts: "This can be done by checking if all the vertices of the grid-cell are inside the polygon" – no sir, it can’t be done your way: const.me/tmp/cell-vertices-inside-polygon.png
• user1669710: I was assuming convex polygons, somehow, although the OP did not mention it. You are right, this doesn't work in the non-convex case.