37

I'm looking for an algorithm, a high-level solution, or even a library which can help me determine if two polygons intersect, in Python.

I have the vertices of the two polygons (These are single-part polygons without any holes) in two different arrays. The polygons are 2D (i.e. just X and Y coordinates)

I'll like to make a function which will return a boolean indicating whether these two polygons intersect.

Please note that I cannot use arcpy, or any arcgis components in this.

Can you suggest an algorithm or library for doing this?

Taras
35.8k5 gold badges77 silver badges152 bronze badges
asked Mar 18, 2014 at 12:44

6 Answers 6

67

You could try shapely.

They describe spatial relationships and it works on Windows

The spatial data model is accompanied by a group of natural language relationships between geometric objects – contains, intersects, overlaps, touches, etc. – and a theoretical framework for understanding them using the 3x3 matrix of the mutual intersections of their component point sets

The following code shows how you can test for intersection:

from shapely.geometry import Polygon
p1 = Polygon([(0,0), (1,1), (1,0)])
p2 = Polygon([(0,1), (1,0), (1,1)])
print(p1.intersects(p2))
Taras
35.8k5 gold badges77 silver badges152 bronze badges
answered Mar 18, 2014 at 13:37
26

You can use the GDAL/OGR Python bindings for that.

from osgeo import ogr
wkt1 = "POLYGON ((1208064.271243039 624154.6783778917, 1208064.271243039 601260.9785661874, 1231345.9998651114 601260.9785661874, 1231345.9998651114 624154.6783778917, 1208064.271243039 624154.6783778917))"
wkt2 = "POLYGON ((1199915.6662253144 633079.3410163528, 1199915.6662253144 614453.958118695, 1219317.1067437078 614453.958118695, 1219317.1067437078 633079.3410163528, 1199915.6662253144 633079.3410163528)))"
poly1 = ogr.CreateGeometryFromWkt(wkt1)
poly2 = ogr.CreateGeometryFromWkt(wkt2)
intersection = poly1.Intersection(poly2)
print intersection.ExportToWkt()

It returns None if they don't intersect. If they intersect it returns the geometry where both intersect.

Also, you can find further info in the GDAL/OGR Cookbook.

Taras
35.8k5 gold badges77 silver badges152 bronze badges
answered Mar 18, 2014 at 13:00
3
  • I would love to use this, but I'm on windows, and on both the systems i've tried, I can't get the python bindings to work. I run into the problem described in this post: gis.stackexchange.com/questions/44958/… Commented Mar 18, 2014 at 13:13
  • 1
    Just in case someone else stumbles upon this, it is possible to use GDAL/OGR with Python in Windows (and within ArcGIS, no less): gis.stackexchange.com/questions/74524/… Commented Mar 19, 2014 at 12:16
  • you also can write intersection = poly1.Intersect(poly2) --- the value of intersection will be TRUE or FALSE depending on if the polygons intersect Commented Nov 19, 2015 at 10:01
4

I know this is an old question, but I've written a python library for handling collisions between concave and convex polygons, as well as circles.

It's pretty simple to use, here you go!

Example:

from collision import *
from collision import Vector as v
p0 = Concave_Poly(v(0,0), [v(-80,0), v(-20,20), v(0,80), v(20,20), v(80,0), v(20,-20), v(0,-80), v(-20,-20)])
p1 = Concave_Poly(v(20,20), [v(-80,0), v(-20,20), v(0,80), v(20,20), v(80,0), v(20,-20), v(0,-80), v(-20,-20)])
print(collide(p0,p1))

You can also have it generate a reponse, which includes:

overlap (how much they overlap)
overlap vector (when subtracted from second shapes position, the shapes will no longer be colliding)
overlap vector normalized (vector direction of collision)
a in b (whether the first shape is fully inside the second)
b in a (whether the second shape is fully inside the first)

https://github.com/QwekoDev/collision

answered Nov 14, 2018 at 18:19
4

You can also use PyQGIS.

from qgis.core import QgsGeometry
p1 = QgsGeometry.fromWkt('POLYGON ((0 0, 1 1, 1 0, 0 0))')
p2 = QgsGeometry.fromWkt('POLYGON ((0 1, 1 0, 1 1, 0 1))')
intersection = p1.intersection(p2)
print(intersection)
# OUT
# <QgsGeometry: Polygon ((1 1, 1 0, 0.5 0.5, 1 1))>
Taras
35.8k5 gold badges77 silver badges152 bronze badges
answered Jan 19, 2023 at 12:22
3

If you want to know the level you can use this. As an argument, you can give a list of polygons. And as a return value, you get a list of levels. In the list of levels, there are the polygons.

from shapely.geometry import Point
from shapely.geometry.polygon import Polygon
def isPolygonInPolygon(poly1,poly2):
 poly2 = Polygon(poly2)
 for poi in poly1:
 poi = Point(poi)
 if(poly2.contains(poi)):
 return True
def polygonTransformHierarchy(polygon_list):
 polygon_list_hierarchy = []
 for polygon1 in polygon_list:
 level = 0
 for polygon2 in polygon_list:
 if(isPolygonInPolygon(polygon1, polygon2)):
 level += 1
 if(level > len(polygon_list_hierarchy)-1):
 dif = (level+1)- len(polygon_list_hierarchy)
 for _ in range(dif):
 polygon_list_hierarchy.append([]) 
 polygon_list_hierarchy[level].append(polygon1)
 return polygon_list_hierarchy
Taras
35.8k5 gold badges77 silver badges152 bronze badges
answered Jul 18, 2019 at 14:36
1

If you know or are interested in learning R it has some useful spatial packages. http://cran.r-project.org/web/views/Spatial.html There is Python module to interact with R (RPy*)

answered Mar 19, 2014 at 11:50

Your Answer

Draft saved
Draft discarded

Sign up or log in

Sign up using Google
Sign up using Email and Password

Post as a guest

Required, but never shown

Post as a guest

Required, but never shown

By clicking "Post Your Answer", you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.