Point in a polygon

An important, yet basic problem in computational geometry is the following task: given a polygon (a sequence of points) and a point, determine if the point is inside, or outside the polygon.

To a human eye, this problem is trivial. After all, it seems obvious that point A is outside the polygon, whereas point B is inside it.

Computers, on the other hand, lack this simple intuition. It is not entirely obvious how to go about creating an algorithm that correctly determines this.

Triangles containing the origin

Let us first focus on an easier version of the problem. Instead of any arbitrary polygon, let our polygon be limited to three points, that is, a triangle. Instead of any point in space, let us determine if the triangle contains the origin.

We observe that we can draw a ray from our point outwards. Any direction will do. If the point is inside the triangle, then the ray will eventually hit a side of it.

On the other hand, if the point is outside the triangle, the ray may hit two sides of the triangle, or it may hit none. However, it will never hit exactly one side:

Aha; this is the approach. Now we can transpose our diagram to a Cartesian plane. For simplicity, let our ray be a line directly upwards from the origin, that is, the entire y-axis above the origin.

Our algorithm now takes the three line segments of the triangle, and counts how many of them intersect the y-axis above the origin.

There are three cases. In the first case, the triangle doesn’t intersect the y-intercept at all above the origin, thus the triangle doesn’t contain the origin:

In the second scenario, the triangle intersects the y-axis twice above the origin, in which case the triangle still does not contain the origin:

In the last case, the triangle intersects the y-axis once above the origin, in which case the origin is in the triangle:

So the algorithm is, for every line segment (pair of coordinates) in the triangle:

  1. Calculate the y-intercept, say, b
  2. If b<0 this line segment doesn’t cross the y-axis above the origin so throw it away. Otherwise, continue.
  3. Make sure that the line segment actually crosses the y-axis at all. Or, check that the two points are on opposite sides of the y-axis.

Here’s my implementation of the algorithm in C++:


// Slope of a line
double slope(double x1, double y1, double x2, double y2){
  return (y1-y2) / (x1-x2);
}

// Y intercept given slope and point
double yint(double px, double py, double m){
  return py - m*px;
}

// Determine if a line crosses the y axis
bool cyax(double ax, double bx){
  if(ax<0 && bx>0) return true;
  if(bx<0 && ax>0) return true;
  return false;
}

// Contains origin, or not?
// Input is (A_x, A_y, B_x, B_y, C_x, C_y).
bool cto(double ax, double ay, double bx, double by, double cx, double cy){

  // Find slopes
  double mAB = slope(ax,ay,bx,by);
  double mBC = slope(bx,by,cx,cy);
  double mAC = slope(ax,ay,cx,cy);

  // Find y intercepts
  double yAB = yint(ax,ay,mAB);
  double yBC = yint(bx,by,mBC);
  double yAC = yint(cx,cy,mAC);

  // How many times it crosses the y intercept above 0
  int c = 0;

  if(yAB>0 && cyax(ax,bx)) c++;
  if(yBC>0 && cyax(bx,cx)) c++;
  if(yAC>0 && cyax(ax,cx)) c++;

  return c==1;
}

Generalizing the algorithm

We can proceed to generalize the algorithm to not only triangles containing the origin. First, for the problem of determining if an arbitrary point is in a triangle, we can create a bijection:

Translate all the points so that the test point lies on the origin. The position of this point relative to the triangle is exactly the same, so we can proceed with the y-intercept algorithm.

For a polygon with more than 3 points, the same algorithm applies. Sort of. With more complicated polygons, we may encounter this edge case:

Here the ray crosses the polygon three times.

There’s an easy fix to this problem, however. It turns out that if the ray crosses the polygon an odd number of times, then the point is inside the polygon.

2 Responses to Point in a polygon

  1. ternarym01 says:

    Why don’t you give the implementation for the generalized algorithm?

    There is another way, that is to break up the polygon into triangles, and testing for each triangle. I don’t know if this is better or a bad idea, but there is a question left: how do we break up the polygon into triangles?

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

Follow

Get every new post delivered to your Inbox.

Join 61 other followers

%d bloggers like this: