# Determine if a polygon can be translated to contain another

I would like to determine whether a polygon could or could not fit inside another polygon. This latter polygon can be translated arbitrarily. A Python example:

``````from shapely.geometry import Polygon
import geopandas as gpd

p = Polygon([[0, 0], [1, 0], [1, 1], [0, 1], [0, 0]])
q0 = Polygon([[1.5, 1.6], [2, 1.6], [2, 2],[1.5, 1.6]])
q1 = Polygon([[1.4, -0.45], [1.6, -0.45], [1.95, 1.4], [1.4, 1.4], [1.4, -0.2]])
q2 = Polygon([[1.1, -0.05], [1.2, -0.05], [1.2,1.05],[1.1,1.05],[1.1,-0.05]])

gpd.GeoSeries([p,q0,q1,q2]).plot(color=['r','g','b','y'],alpha=0.5)
``````

Clearly, the red square (`p0`) can be translated to contain the green triangle (`q0`). The blue trapezoid (`q1`) and the yellow rectangle (`q2`) cannot be contained, no matter how we translate the square. (Rotation is not allowed in this example.) Given GeoSeries:

``````qs = gpd.GeoSeries([q0,q1,q2])
``````

and polygon `p`, I am looking for function `does_it_fit`, which, when applied to `qs`:

``````qs.geometry.apply(lambda row: does_it_fit(row,p))
``````

returns:

``````0    True
1    False
2    False
dtype: bool
``````

ie:

``````def does_it_fit(q,p):
### what comes here? ###
return result
``````

Related but different questions:

This problem was inspired by the interesting task of finding parcels on a map which can contain a building with certain size & orientation. Notebook to play with the above example here.

• Do you want to try each combination e.g. `(p,q0)`, `(p,q1)`, `(p,q2)`, `(q0,q1)`, ... or each `q` against single `p`? Oct 2, 2022 at 8:51
• Each `q` against single `p`. I provided more `q`s to explain requirements better (green: fits, blue: too big, yellow: only fits if rotation is allowed, so does not fit in this case). Oct 2, 2022 at 9:35

If rotations are strictly forbidden, what you basically want to do is to translate all your polygons towards the origin. More specifically, by attaching their most lowest-left point to the origin.

After that, it will be trivial to see for each one, if it does overlap another one.

Let's dive into the code:

``````
import numpy as np
from shapely.geometry import Point, Polygon
from shapely.affinity import affine_transform
from shapely.ops import nearest_points
import geopandas as gpd
import matplotlib.pyplot as plt

ref = Point([-10,-10])
p = Polygon([[0, 0], [1, 0], [1, 1], [0, 1], [0, 0]])
q0 = Polygon([[1.5, 1.6], [2, 1.6], [2, 2],[1.5, 1.6]])
q1 = Polygon([[1.4, -0.45], [1.6, -0.45], [1.95, 1.4], [1.4, 1.4], [1.4, -0.2]])
q2 = Polygon([[1.1, -0.05], [1.2, -0.05], [1.2,1.05],[1.1,1.05],[1.1,-0.05]])

polys = [p, q0, q1, q2]

gpd.GeoSeries([p,q0,q1,q2]).plot(color=['r','g','b','y'],alpha=0.5)
``````

Here, I simply add some imports and defined a reference point `ref` which will serve us to find the lowest left point for each polygon. Of course, if you have other polygon, you may want to displace this point elsewhere but setting [-10,-10] is perfectly fine with the data your provided.

Then we define convenience functions:

``````def find_nearest_point(poly,ref):
"""Find the closest point of a polygon to a reference point."""
p0, _ = nearest_points(poly, ref)
p0 = tuple(p0.coords)[0]

return p0

def transform_affine(poly):
"""Translate the polygon so that it's lowest left point
is at the origin (0,0)."""
x,y = find_nearest_point(poly,ref)
m = [1,0,0,1,-x,-y]
poly_new = affine_transform(poly, m)

return poly_new
``````

We can now check that everything went fine by plotting them:

``````new_polys = [transform_affine(poly) for poly in polys]
g = gpd.GeoSeries(new_polys)
fig, ax = plt.subplots()
g.plot(color=['r','g','b','y'], alpha=0.24, ax=ax)
ax.grid()
ax.set_aspect('equal', 'box')
plt.show()
``````

Things become more obvious now...

But we have to verify it mathematically:

``````n = len(g)
fig, axs = plt.subplots(n,n, figsize=(15, 15))
for i, gi in enumerate(g):
gs = gpd.GeoSeries(np.roll(g, shift=i+1))
for j, gsi in enumerate(gs):
if gsi.contains(gi) and not gsi.equals_exact(gi,0):
print(f'iter {i}:\n', f'{gsi.wkt} contains {gi.wkt}')
xa,ya = gi.exterior.xy
xb,yb = gsi.exterior.xy
axs[i,j].plot(xa,ya,'-r', lw=4)
axs[i,j].plot(xb,yb,'-b', lw=2)
axs[i,j].grid()
``````

This code is purposefully expanded in order to make the logic clear enough.

Results is as follow:

``````iter 1:
POLYGON ((0 0, 1 0, 1 1, 0 1, 0 0)) contains POLYGON ((0 0, 0.5 0, 0.5 0.4, 0 0))
iter 3:
POLYGON ((0 0, 0.200 0, 0.55 1.85, 0 1.85, 0 0.25, 0 0)) contains POLYGON ((0 0, 0.1 0, 0.1 1.1, 0 1.1, 0 0))
``````

This states which polygon contain which other. I removed self containing polygons with the test `and not gsi.equals_exact(gi,0)`.

And the figure shows every occurrence when it happens, otherwise it's a white plot. You can figure out how to discard those by plotting each individual figure separately when conditions are met. But I like this kind of matrix representation:

This is a pairwise operation. As your problem is apparently a subset of that, I let you figure out the was to discard what's unnecessary to you by, e.g. removing a loop if needed. He who can do more can do less.

The mathematics behind: https://en.wikipedia.org/wiki/Transformation_matrix#Affine_transformations

• Imagine you wish to check for containment of a small pentagon in a square large enough to contain it. There is no vertex of the polygon that you can translate to any vertex of the square and have the pentagon contained in the square. Oct 3, 2022 at 3:41
• what is the "lower-left" vertex of a polygon? Oct 3, 2022 at 3:42
• Cool approach! I agree with @Llaves: I fail to see how this approach covers all possible cases. Or, to give another example: I think above method would give wrong results if the red polygon was `Polygon([[0,0],[10,0],[10,2],[8,2],[7,0.1],[0,0]])`. Still, nice things to learn from this answer! Oct 3, 2022 at 7:29