How To Check If A Point Belongs To A Polygon (Ray Casting & More)
Checking if a point falls within a polygon’s boundaries is a fundamental problem in various domains: geographical information systems (GIS), game development, computer graphics, and even collision detection. It might sound simple, but dealing with concave polygons, points on edges, or vertices can make it surprisingly tricky.
In this post, we’ll break down the core concepts and provide pragmatic, working code examples to solve this problem effectively.
The Problem Defined
You have two inputs:
- A point (typically
(x, y)coordinates). - A polygon, represented as an ordered list of vertices (e.g.,
[(x1, y1), (x2, y2), ..., (xn, yn)]). For a valid polygon, the last vertex is implicitly connected back to the first.
The goal is to determine if the point is:
- Inside the polygon.
- Outside the polygon.
- On the boundary (an edge or a vertex) of the polygon.
The “on boundary” case is often a nuance. Some applications consider it “inside,” others “outside,” and some need a specific “on boundary” status. We’ll address this.
Step 1: The Quick Filter - Bounding Box Check
Before diving into complex algorithms, a bounding box check is an excellent first step. It’s cheap and can quickly rule out many points that are definitely outside.
A polygon’s bounding box is the smallest rectangle (aligned with the axes) that completely encloses it. If a point is outside this bounding box, it’s certainly outside the polygon.
How it works:
- Find the minimum and maximum X coordinates (
min_x,max_x) among all polygon vertices. - Find the minimum and maximum Y coordinates (
min_y,max_y) among all polygon vertices. - If the point’s X is less than
min_xor greater thanmax_x, or if its Y is less thanmin_yor greater thanmax_y, it’s outside the bounding box.
Python Example: Bounding Box Pre-Check
def is_point_in_bounding_box(point, polygon_vertices):
"""
Checks if a point is within the bounding box of a polygon.
Returns True if potentially inside, False if definitely outside.
"""
if not polygon_vertices:
return False
min_x = min(v[0] for v in polygon_vertices)
max_x = max(v[0] for v in polygon_vertices)
min_y = min(v[1] for v in polygon_vertices)
max_y = max(v[1] for v in polygon_vertices)
px, py = point
return min_x <= px <= max_x and min_y <= py <= max_y
# Example usage
square_polygon = [(0, 0), (0, 10), (10, 10), (10, 0)]
print(f"Point (5, 5) in bounding box of square: {is_point_in_bounding_box((5, 5), square_polygon)}")
print(f"Point (15, 5) in bounding box of square: {is_point_in_bounding_box((15, 5), square_polygon)}")
print(f"Point (0, 0) in bounding box of square: {is_point_in_bounding_box((0, 0), square_polygon)}")
print(f"Point (-1, 5) in bounding box of square: {is_point_in_bounding_box((-1, 5), square_polygon)}")
triangle_polygon = [(1, 1), (5, 8), (9, 1)] # A wide triangle
print(f"\nPoint (5, 4) in bounding box of triangle: {is_point_in_bounding_box((5, 4), triangle_polygon)}")
print(f"Point (5, 9) in bounding box of triangle: {is_point_in_bounding_box((5, 9), triangle_polygon)}")
print(f"Point (0, 0) in bounding box of triangle: {is_point_in_bounding_box((0, 0), triangle_polygon)}")Sample Output: Bounding Box Pre-Check
Point (5, 5) in bounding box of square: True
Point (15, 5) in bounding box of square: False
Point (0, 0) in bounding box of square: True
Point (-1, 5) in bounding box of square: False
Point (5, 4) in bounding box of triangle: True
Point (5, 9) in bounding box of triangle: False
Point (0, 0) in bounding box of triangle: FalseNote: A point being within the bounding box does not guarantee it’s inside the polygon. For example, a point might be inside the bounding box of a U-shaped polygon but outside the “U”. This is why we need a more robust algorithm for the definitive check.
Step 2: The Definitive Check - Ray Casting (Crossing Number Algorithm / PNPoly)
The Ray Casting Algorithm (also known as the Crossing Number Algorithm or PNPoly) is the most widely used method for determining if a point is inside a polygon. It’s relatively easy to understand and implement, and it’s robust for most simple and concave polygons (though it has nuances with points exactly on vertices or edges, which we’ll cover).
How it works:
- Pick an arbitrary direction (e.g., horizontally to the right) and imagine a ray extending infinitely from your test point in that direction.
- Count how many times this ray intersects with the edges of the polygon.
- If the count of intersections is odd, the point is inside the polygon.
- If the count of intersections is even, the point is outside the polygon.
Why does this work? Imagine walking along the ray. Every time you cross an edge, you alternate between being inside and outside the polygon. If you start outside (which you do, infinitely far along the ray), and you end up inside, you must have crossed an odd number of edges. If you end up outside, you must have crossed an even number.
Edge Cases and Robustness for Ray Casting:
Implementing ray casting correctly requires careful handling of edge cases to avoid incorrect counts:
- Horizontal Edges: If the ray runs parallel to a polygon edge (a horizontal edge), it shouldn’t count as an intersection, or it could lead to double-counting. A common strategy is to ignore intersections with horizontal edges entirely.
- Vertex on Ray: If the ray passes exactly through a vertex, it’s tricky. Does it count as one intersection or two (one for each edge meeting at the vertex)? The robust solution (as in the classic PNPoly algorithm) is to count an intersection only if the edge
(v_i, v_{i+1})crosses the ray and the y-coordinate ofv_iis belowpyandv_{i+1}is abovepy(or vice-versa), ensuring you count each ‘crossing’ correctly without double-counting when the ray passes through a vertex. - Point on Boundary: The standard ray casting algorithm does not inherently distinguish between a point strictly inside and a point on the boundary. A point directly on an edge or vertex might be counted as inside or outside depending on the exact floating-point comparisons and vertex handling. For true boundary detection, a separate check is often needed. We’ll implement a
is_point_on_segmenthelper.
Python Implementation: Ray Casting (PNPoly Logic)
This implementation is based on the widely used PNPoly algorithm by W. Randolph Franklin. It’s known for its simplicity and robustness.
def is_point_on_segment(point, p1, p2, epsilon=1e-9):
"""
Checks if a point lies on a line segment defined by p1 and p2.
Uses a small epsilon for floating point tolerance.
"""
px, py = point
x1, y1 = p1
x2, y2 = p2
# Check if point is collinear with segment
# (px - x1) * (y2 - y1) == (py - y1) * (x2 - x1)
cross_product = (px - x1) * (y2 - y1) - (py - y1) * (x2 - x1)
if abs(cross_product) > epsilon:
return False
# Check if point is within the bounding box of the segment
# This covers collinearity and ensures it's *on* the segment, not just the line
if min(x1, x2) <= px <= max(x1, x2) and \
min(y1, y2) <= py <= max(y1, y2):
return True
return False
def is_point_in_polygon_ray_casting(point, polygon_vertices):
"""
Checks if a point is inside a polygon using the Ray Casting (PNPoly) algorithm.
Returns True for inside, False for outside.
Does NOT explicitly handle points *on* the boundary as inside, this needs a separate check.
"""
if not polygon_vertices or len(polygon_vertices) < 3:
return False # A polygon needs at least 3 vertices
px, py = point
n = len(polygon_vertices)
intersections = 0
# Ensure the polygon is closed for iteration
# The last vertex connects to the first, so we iterate n times.
# The current point `p1` connects to `p2`.
for i in range(n):
p1 = polygon_vertices[i]
p2 = polygon_vertices[(i + 1) % n]
x1, y1 = p1
x2, y2 = p2
# First, check if the point is on an edge (or a vertex)
if is_point_on_segment(point, p1, p2):
return True # Or return a specific "ON_BOUNDARY" status if needed
# PNPoly logic for counting intersections
# This condition effectively checks if the edge crosses the horizontal ray from `px, py`
# and one vertex is above the ray while the other is below (or on) the ray.
# This also correctly handles vertices on the ray.
if ((y1 <= py < y2) or (y2 <= py < y1)) and \
(px < (x2 - x1) * (py - y1) / (y2 - y1) + x1):
intersections += 1
return intersections % 2 == 1
# --- Example Usage ---
# Convex Polygon (Square)
square_polygon = [(0, 0), (0, 10), (10, 10), (10, 0)]
print("\n--- Square Polygon ---")
print(f"Point (5, 5) inside square: {is_point_in_polygon_ray_casting((5, 5), square_polygon)}")
print(f"Point (1, 1) inside square: {is_point_in_polygon_ray_casting((1, 1), square_polygon)}")
print(f"Point (15, 5) inside square: {is_point_in_polygon_ray_casting((15, 5), square_polygon)}")
print(f"Point (-1, 5) inside square: {is_point_in_polygon_ray_casting((-1, 5), square_polygon)}")
print(f"Point (0, 0) (vertex) inside square: {is_point_in_polygon_ray_casting((0, 0), square_polygon)}")
print(f"Point (5, 0) (edge) inside square: {is_point_in_polygon_ray_casting((5, 0), square_polygon)}")
print(f"Point (10, 5) (edge) inside square: {is_point_in_polygon_ray_casting((10, 5), square_polygon)}")
# Concave Polygon (a star or 'U' shape)
# (0,0) -- (10,0) -- (10,10) -- (5,5) -- (0,10) -- (0,0)
concave_polygon = [(0, 0), (10, 0), (10, 10), (5, 5), (0, 10)]
print("\n--- Concave Polygon (U-shape) ---")
print(f"Point (1, 1) inside concave: {is_point_in_polygon_ray_casting((1, 1), concave_polygon)}")
print(f"Point (9, 1) inside concave: {is_point_in_polygon_ray_casting((9, 1), concave_polygon)}")
print(f"Point (5, 6) (in the 'U' opening) inside concave: {is_point_in_polygon_ray_casting((5, 6), concave_polygon)}") # Should be outside
print(f"Point (5, 3) (deep inside) inside concave: {is_point_in_polygon_ray_casting((5, 3), concave_polygon)}") # Should be inside
print(f"Point (11, 1) outside concave: {is_point_in_polygon_ray_casting((11, 1), concave_polygon)}")
print(f"Point (5, 5) (vertex) inside concave: {is_point_in_polygon_ray_casting((5, 5), concave_polygon)}") # on boundary
print(f"Point (0, 5) (edge) inside concave: {is_point_in_polygon_ray_casting((0, 5), concave_polygon)}") # on boundarySample Output: Ray Casting (PNPoly)
--- Square Polygon ---
Point (5, 5) inside square: True
Point (1, 1) inside square: True
Point (15, 5) inside square: False
Point (-1, 5) inside square: False
Point (0, 0) (vertex) inside square: True
Point (5, 0) (edge) inside square: True
Point (10, 5) (edge) inside square: True
--- Concave Polygon (U-shape) ---
Point (1, 1) inside concave: True
Point (9, 1) inside concave: True
Point (5, 6) (in the 'U' opening) inside concave: False
Point (5, 3) (deep inside) inside concave: True
Point (11, 1) outside concave: False
Point (5, 5) (vertex) inside concave: True
Point (0, 5) (edge) inside concave: TrueNote: The is_point_on_segment check at the beginning of is_point_in_polygon_ray_casting is crucial if you want points on the boundary to be considered “inside.” Without it, the PNPoly logic itself might classify them as outside due to the strict less-than/greater-than comparisons for crossing. Adding this check explicitly ensures that boundary points are always considered “inside.” If you need to distinguish “inside” from “on boundary,” you’d modify the function to return an enum or string like “INSIDE”, “OUTSIDE”, “ON_BOUNDARY”.
Combining Checks for a Robust Solution
A common robust pattern is to:
- Perform a
is_point_in_bounding_boxcheck first. IfFalse, return “OUTSIDE” immediately. - Then, perform the
is_point_in_polygon_ray_castingcheck.
Let’s refine our main function to give more explicit statuses.
import math
class PointInPolygonStatus:
OUTSIDE = "OUTSIDE"
INSIDE = "INSIDE"
ON_BOUNDARY = "ON_BOUNDARY"
def get_point_position_in_polygon(point, polygon_vertices, epsilon=1e-9):
"""
Determines if a point is inside, outside, or on the boundary of a polygon.
Combines bounding box pre-check, boundary check, and ray casting.
"""
if not polygon_vertices or len(polygon_vertices) < 3:
return PointInPolygonStatus.OUTSIDE # Or raise an error for invalid polygon
px, py = point
n = len(polygon_vertices)
# 1. Bounding Box Pre-Check
min_x = min(v[0] for v in polygon_vertices)
max_x = max(v[0] for v in polygon_vertices)
min_y = min(v[1] for v in polygon_vertices)
max_y = max(v[1] for v in polygon_vertices)
if not (min_x <= px <= max_x and min_y <= py <= max_y):
return PointInPolygonStatus.OUTSIDE
# 2. Check if point is ON the boundary (vertex or edge)
for i in range(n):
p1 = polygon_vertices[i]
p2 = polygon_vertices[(i + 1) % n]
if is_point_on_segment(point, p1, p2, epsilon):
return PointInPolygonStatus.ON_BOUNDARY
# 3. Ray Casting (PNPoly) for strictly inside/outside
# If we reached here, the point is within the bounding box and not on the boundary.
intersections = 0
for i in range(n):
p1 = polygon_vertices[i]
p2 = polygon_vertices[(i + 1) % n]
x1, y1 = p1
x2, y2 = p2
# This is the PNPoly logic
if ((y1 <= py < y2) or (y2 <= py < y1)) and \
(px < (x2 - x1) * (py - y1) / (y2 - y1) + x1):
intersections += 1
if intersections % 2 == 1:
return PointInPolygonStatus.INSIDE
else:
return PointInPolygonStatus.OUTSIDE
# --- Extended Example Usage ---
square_polygon = [(0, 0), (0, 10), (10, 10), (10, 0)]
concave_polygon = [(0, 0), (10, 0), (10, 10), (5, 5), (0, 10)]
print("\n--- Comprehensive Polygon Checks ---")
test_points = [
((5, 5), "Square (Strictly Inside)"),
((1, 1), "Square (Strictly Inside)"),
((15, 5), "Square (Outside Bounding Box)"),
((-1, 5), "Square (Outside Bounding Box)"),
((0, 0), "Square (On Vertex)"),
((5, 0), "Square (On Edge)"),
((10, 5), "Square (On Edge)"),
((0, 5), "Square (On Edge)"), # A vertical edge check
((1, 1), "Concave (Strictly Inside)"),
((9, 1), "Concave (Strictly Inside)"),
((5, 6), "Concave (In U-opening, Outside)"),
((5, 3), "Concave (Deep Inside)"),
((11, 1), "Concave (Outside)"),
((5, 5), "Concave (On Vertex)"),
((0, 5), "Concave (On Edge)"),
((10, 5), "Concave (On Edge - vertical)"), # For the concave polygon (10,0)-(10,10) edge
((5, 0), "Concave (On Edge - horizontal)"), # For the concave polygon (0,0)-(10,0) edge
]
print("--- Square Polygon Tests ---")
for point, description in test_points[:8]:
status = get_point_position_in_polygon(point, square_polygon)
print(f"Point {point} ({description}): {status}")
print("\n--- Concave Polygon Tests ---")
for point, description in test_points[8:]:
status = get_point_position_in_polygon(point, concave_polygon)
print(f"Point {point} ({description}): {status}")Sample Output: Comprehensive Polygon Checks
--- Comprehensive Polygon Checks ---
--- Square Polygon Tests ---
Point (5, 5) (Square (Strictly Inside)): INSIDE
Point (1, 1) (Square (Strictly Inside)): INSIDE
Point (15, 5) (Square (Outside Bounding Box)): OUTSIDE
Point (-1, 5) (Square (Outside Bounding Box)): OUTSIDE
Point (0, 0) (Square (On Vertex)): ON_BOUNDARY
Point (5, 0) (Square (On Edge)): ON_BOUNDARY
Point (10, 5) (Square (On Edge)): ON_BOUNDARY
Point (0, 5) (Square (On Edge)): ON_BOUNDARY
--- Concave Polygon Tests ---
Point (1, 1) (Concave (Strictly Inside)): INSIDE
Point (9, 1) (Concave (Strictly Inside)): INSIDE
Point (5, 6) (Concave (In U-opening, Outside)): OUTSIDE
Point (5, 3) (Concave (Deep Inside)): INSIDE
Point (11, 1) (Concave (Outside)): OUTSIDE
Point (5, 5) (Concave (On Vertex)): ON_BOUNDARY
Point (0, 5) (Concave (On Edge)): ON_BOUNDARY
Point (10, 5) (Concave (On Edge - vertical)): ON_BOUNDARY
Point (5, 0) (Concave (On Edge - horizontal)): ON_BOUNDARYThis combined approach provides a robust and clear way to classify a point’s relationship to a polygon.
What About Other Algorithms?
While Ray Casting is excellent for most cases, especially simple and concave polygons without self-intersections, you might encounter mentions of other algorithms:
- Winding Number Algorithm: This algorithm calculates the total angle swept by vectors from the point to each successive pair of vertices. If the winding number is a non-zero multiple of 2π (or 360 degrees), the point is inside. If it’s zero, the point is outside.
- Pros: More robust for self-intersecting polygons (e.g., a figure-eight polygon). It correctly identifies points in “holes” based on the winding direction.
- Cons: Computationally more expensive due to trigonometric functions (arctan2) or complex cross-product sums. More complex to implement correctly.
- Grid-Based or Spatial Indexing: For applications with millions of points and thousands of polygons (like GIS), individual point-in-polygon checks are too slow. Spatial indexing (R-trees, Quadtrees, K-D trees) are used to quickly narrow down the candidate polygons for a given point, significantly speeding up queries. These are typically used before any of the detailed algorithms.
For typical day-to-day point-in-polygon needs with non-self-intersecting polygons, the Ray Casting algorithm, especially with the is_point_on_segment pre-check, is often the sweet spot between performance and implementation complexity.
Real-World Libraries
For serious geometric operations, especially in Python, you rarely implement these algorithms from scratch. Libraries provide highly optimized and thoroughly tested implementations:
- Shapely (Python): A fantastic library for planar geometric objects. It’s built on top of the battle-tested GEOS library (C++). It offers powerful methods like
contains(),intersects(),within(),touches(), anddistance().
Python Example: Using Shapely
First, install it: pip install shapely
from shapely.geometry import Point, Polygon
# Define a polygon using Shapely's Polygon object
square_polygon_shapely = Polygon([(0, 0), (0, 10), (10, 10), (10, 0)])
concave_polygon_shapely = Polygon([(0, 0), (10, 0), (10, 10), (5, 5), (0, 10)])
# Define points
point_inside = Point(5, 5)
point_outside = Point(15, 5)
point_on_edge = Point(5, 0)
point_on_vertex = Point(0, 0)
point_concave_hole = Point(5, 6) # For the U-shaped concave polygon
print("\n--- Shapely Examples ---")
# Square Polygon
print("\nSquare Polygon (Shapely):")
print(f"Point {point_inside} strictly inside? {square_polygon_shapely.contains(point_inside)}")
print(f"Point {point_outside} strictly inside? {square_polygon_shapely.contains(point_outside)}")
print(f"Point {point_on_edge} strictly inside? {square_polygon_shapely.contains(point_on_edge)}") # Shapely.contains does NOT include boundary points
print(f"Point {point_on_vertex} strictly inside? {square_polygon_shapely.contains(point_on_vertex)}") # Shapely.contains does NOT include boundary points
# To check if on boundary or inside, use .intersects or .touches
print(f"Point {point_on_edge} intersects? {square_polygon_shapely.intersects(point_on_edge)}") # Intersects means inside OR on boundary
print(f"Point {point_on_edge} touches? {square_polygon_shapely.touches(point_on_edge)}") # Touches means on boundary, not strictly inside
# Concave Polygon
print("\nConcave Polygon (Shapely):")
print(f"Point {point_inside} strictly inside? {concave_polygon_shapely.contains(point_inside)}") # 5,5 is a vertex here for concave_polygon_shapely, so it should be False for contains()
print(f"Point {Point(5,3)} strictly inside? {concave_polygon_shapely.contains(Point(5,3))}")
print(f"Point {point_concave_hole} strictly inside? {concave_polygon_shapely.contains(point_concave_hole)}") # In the 'U' opening
print(f"Point {point_on_vertex} touches? {concave_polygon_shapely.touches(point_on_vertex)}")
print(f"Point {Point(5,5)} touches? {concave_polygon_shapely.touches(Point(5,5))}") # This point (5,5) is a vertex of the concave polygonSample Output: Shapely Examples
--- Shapely Examples ---
Square Polygon (Shapely):
Point POINT (5 5) strictly inside? True
Point POINT (15 5) strictly inside? False
Point POINT (5 0) strictly inside? False
Point POINT (0 0) strictly inside? False
Point POINT (5 0) intersects? True
Point POINT (5 0) touches? True
Concave Polygon (Shapely):
Point POINT (5 5) strictly inside? False
Point POINT (5 3) strictly inside? True
Point POINT (5 6) strictly inside? False
Point POINT (0 0) touches? True
Point POINT (5 5) touches? TrueNotice how shapely.contains() strictly means inside (not on the boundary). For “inside or on boundary,” you’d typically use polygon.intersects(point). This distinction is important and common in geometry libraries.
Conclusion
Determining if a point lies within a polygon is a common task with practical implications across many fields. While conceptually simple, robust implementation requires careful attention to edge cases.
- Bounding Box: Always use this as a quick initial filter to rule out points far away. It’s cheap and effective.
- Ray Casting (PNPoly): This is your workhorse for most polygons (simple or concave, non-self-intersecting). It’s efficient and relatively easy to implement. Remember to explicitly handle points on the boundary if your definition of “inside” includes them.
- Winding Number: Consider this for complex or self-intersecting polygons, though it’s more computationally intensive.
- Libraries: For production-grade applications, rely on well-tested geometric libraries like Shapely (Python), GEOS (C++/Java), or Boost.Geometry (C++). They handle all the nasty edge cases and provide optimized performance.
By understanding these approaches, you’re well-equipped to tackle point-in-polygon problems with confidence.
Last updated on