2

I have a large polygon (A) with a specific area. I also have another polygon (B) with area less than polygon A.

I want to buffer the polygon B so that its area becomes equal to the area of A.

from shapely.geometry import Polygon
ulx,uly= -105.645292, 28.094511
lrx, lry = -101.985052, 24.885167
poly_A = Polygon([(ulx,uly),(ulx,lry),(lrx,lry),(lrx,uly)])
print (poly_A.area)


ulx,uly= -110.398408,  58.136267
lrx, lry = -108.689382,  57.137692
poly_B = Polygon([(ulx,uly),(ulx,lry),(lrx,lry),(lrx,uly)])
print (poly_B.area)


buffer_size = (poly_A.area - poly_B.area) / poly_B.length
print (buffer_size)

poly_C = poly_B.buffer(buffer_size)
print (poly_C.area)

After knowing buffer_size, I can get bufferred polygon C with area equals to polygon A.

Note. All polygons dealt are rectangles. Both polygons are in latlon coordinates

6
3
  1. The .area in shapely is plane area, so we need to transform longitude, latitude (4326) to project spatial reference (e.g. 3857).
  2. For rectangle, if the original size is a * b, area is A, buffer distance is X, the buffer area B then:
B = A + 2 * (a + b) * X + pi * X * X   ----- (1)

2 * (a + b) = c = perimeter = poly_B.length

(1) => X = sqrt(c * c - 4 * (B - A) * PI) / (2 * PI)
import pyproj
import math
from functools import partial
from shapely.geometry import Polygon
from shapely.geometry import shape
from shapely.ops import transform

proj = partial(pyproj.transform, pyproj.Proj(init='epsg:4326'), \
               pyproj.Proj(init='epsg:8858'))

ulx, uly= -105.645292, 28.094511
lrx, lry = -101.985052, 24.885167
poly_A = Polygon([(ulx, uly),(ulx, lry),(lrx, lry),(lrx, uly)])
print (poly_A.area)  # 11.746969282560011


poly_A_proj = transform(proj, poly_A)
poly_A_proj_area = poly_A_proj.area
print(poly_A_proj_area) # 129724395664.97488


ulx, uly= -110.398408, 58.136267
lrx, lry = -108.689382, 57.137692
poly_B = Polygon([(ulx,uly),(ulx,lry),(lrx,lry),(lrx,uly)])
print (poly_B.area)  # 39511382337.41971


poly_B_proj = transform(proj, poly_B)
poly_B_proj_area = poly_B_proj.area


buffer_size = (math.sqrt(poly_B.length * poly_B.length - 4*(poly_B.area - poly_A.area)*math.pi)-poly_B.length) / (2 * math.pi)
print (buffer_size)  # 1.122771717045914


poly_C = poly_B.buffer(buffer_size)
print (poly_C.area)  # 11.740610528262806



12
  • I tried according to you and edited the question. Unfortunately, I did not get expected result – thunder Aug 9 '20 at 3:22
  • Sorry. I think you can project the polygon first and then calculate the area. – Leo Aug 9 '20 at 3:24
  • how can i project? – thunder Aug 9 '20 at 3:27
  • Both polygons are in latlon coordinates – thunder Aug 9 '20 at 3:33
  • And both polygons are rectangles, so a different trick, rather than .buffer, should work, as output is enlarged polygon – thunder Aug 9 '20 at 4:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.