Skip to content

Instantly share code, notes, and snippets.

@dvgodoy
Last active November 24, 2022 16:54
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save dvgodoy/605a6bb1f1dadba1ebef45346fd382a5 to your computer and use it in GitHub Desktop.
Save dvgodoy/605a6bb1f1dadba1ebef45346fd382a5 to your computer and use it in GitHub Desktop.
#!pip install area
from area import area
from functools import partial
def make_polygon(coords, size):
left_long, bottom_lat = coords
return np.array([[[left_long, bottom_lat],
[left_long, bottom_lat + size],
[left_long + size, bottom_lat + size],
[left_long + size, bottom_lat],
[left_long, bottom_lat]]], dtype=np.float32)
def surface_polygons(grid_size=5):
long, lat = np.meshgrid(np.arange(-180, 180, grid_size),
np.arange(90-grid_size, -90-grid_size, -grid_size))
coords = np.dstack([long, lat])
polygons = np.apply_along_axis(func1d=partial(make_polygon, size=grid_size),
axis=2, arr=coords)
return polygons
def surface_area(polygons, perc=True):
meridian_surface_area = np.array(list(map(area,
map(lambda c: {'type': 'Polygon', 'coordinates': c.tolist()},
polygons[:, 0]))))/1e6
surface_area = np.tile(meridian_surface_area.reshape(-1, 1), polygons.shape[1])
surface_perc = (surface_area / surface_area.sum())
return surface_perc if perc else surface_area
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment