I have created a square grid and I would like to figure out the real world area of each polygon in a metric unit e.q. square metres via GeoDataFrame.area.

import geopandas as gpd
from shapely.geometry import Polygon

# Latitude (y)
lat_min = -90
lat_max = 90

# Longtidue (x)
lon_min = -180
lon_max = 180

# Edge length
side_length = 0.5

# List to which we will append the polygons
temp_polygons = []

for y in range(int(abs(lat_min-lat_max)/side_length)):

    for x in range(int(abs(lon_min-lon_max)/side_length)):
        x_start = lon_min + x * side_length
        y_start = lat_min + y * side_length
        bottom_left_x = x_start 
        bottom_left_y = y_start
        bottom_right_x = x_start + side_length
        bottom_right_y = y_start
        top_right_x = x_start + side_length
        top_right_y = y_start + side_length
        top_left_x = x_start
        top_left_y = y_start + side_length
        # Append polygon to list
            # Polygon
                # Bottom-left
                # Bottom-right
                # Top-left
                # Top-right

# Turn list of polygons into GeoSeries
temp_polygons = gpd.GeoSeries(temp_polygons)

# Turn GeoSeries into GeoDataFrame
raster = gpd.GeoDataFrame({'geometry': temp_polygons})

# Add coordinate reference system
raster.crs = {"init":"epsg:4326"}

# Change projection
# https://en.wikipedia.org/wiki/Web_Mercator_projection#EPSG:3857
raster['geometry'] = (
    raster['geometry'].to_crs({'init': 'EPSG:3857'})

# Get area
raster['area'] = raster.area

First, the geometries are as follows:

enter image description here

However, after applying ...

# Change projection
# https://en.wikipedia.org/wiki/Web_Mercator_projection#EPSG:3857
raster['geometry'] = (
    raster['geometry'].to_crs({'init': 'EPSG:3857'})

... there are a lot of inf values and hence the cells' area is nan:

enter image description here

What am I doing wrong? Did I select a wrong CRS? How can I determine the area of each grid cell in a metric unit?


This answer seems to work for me. However, with ~25 iterations/second it's rather slow (~3 hours for ~260,000 rows). It's not prohibitively slow but there must be a quicker solution.

  • 1
    This may help you: gis.stackexchange.com/questions/127607/…
    – Taras
    Oct 5, 2021 at 11:26
  • Thank you for your comment. This answer seems to work: gis.stackexchange.com/a/166421/97137; However, it's rather slow. ~25 iterations per second. That's too slow for my >250,000 rows.
    – Stücke
    Oct 5, 2021 at 11:41
  • 1
    Tangential comment: Do you remember why you thought that EPSG:3857 would give you actual meters? I would love to fix the web on that misleading thought so if it was a website or something, I'd be interested. Oct 5, 2021 at 16:34
  • 1
    While EPSG:3857 won't give you useful areal measurements, are you sure your axes are correct? I get no NaNs with pyproj==3.1.0 and geopandas==0.9.0. Oct 5, 2021 at 20:03
  • 1
    @bugmenot123 Thank you for your comments! Sorry but I cannot recall where I stumbled upon that information. It may have even been Stackexchange but I cannot find it anymore. Regarding the 2nd comment I am not a hundred percent sure what you mean. I provided all the code in the minimal working example.
    – Stücke
    Oct 5, 2021 at 20:15

1 Answer 1


Here is a way to calculate lat/long cell area in m^2 without any transformation as long as they are regular grid cells. In other words the polygons must be squares where each side is a latitudinal or longitudinal line. From the article Santini et al. 2010

It's quite fast. It did the ~260,000 cells in your raster data in ~10 seconds on my small desktop computer.

timeit.timeit(lambda: raster.geometry.apply(lat_lon_cell_area), number=1)

(259200, 2)

I confirmed the values by writing the raster grid to a shapefile and then opening it in qgis where I measured some cell areas manually. The two were approximately the same.

from math import radians, sin

def lat_lon_cell_area(lat_lon_grid_cell):
    Calculate the area of a cell, in meters^2, on a lat/lon grid.
    This applies the following equation from Santini et al. 2010.
    S = (λ_2 - λ_1)(sinφ_2 - sinφ_1)R^2
    S = surface area of cell on sphere
    λ_1, λ_2, = bands of longitude in radians
    φ_1, φ_2 = bands of latitude in radians
    R = radius of the sphere
    Santini, M., Taramelli, A., & Sorichetta, A. (2010). ASPHAA: A GIS‐Based 
    Algorithm to Calculate Cell Area on a Latitude‐Longitude (Geographic) 
    Regular Grid. Transactions in GIS, 14(3), 351-377.

        A shapely box with coordinates on the lat/lon grid

        The cell area in meters^2

    # mean earth radius - https://en.wikipedia.org/wiki/Earth_radius#Mean_radius
    west, south, east, north = lat_lon_grid_cell.bounds

    west = radians(west)
    east = radians(east)
    south = radians(south)
    north = radians(north)
    return (east - west) * (sin(north) - sin(south)) * (AVG_EARTH_RADIUS_METERS**2)

# Example using a 2 degree lat/lon grid

import numpy as np
import geopandas as gpd
from shapely.geometry import box

# Latitude (y)
lat_min = -90
lat_max = 90

# Longitude (x)
lon_min = -180
lon_max = 180

# Edge length
side_length = 2

all_lats, all_lons = np.meshgrid(
    np.arange(lat_min, lat_max, side_length),
    np.arange(lon_min, lon_max, side_length)

polygons = []
for lon, lat in zip(all_lons.flatten(), all_lats.flatten()):
        box(lon, lat, lon+side_length, lat+side_length)

raster = gpd.GeoDataFrame({'geometry': gpd.GeoSeries(polygons)})

raster['cell_area'] = raster.geometry.apply(lat_lon_cell_area)

  • 1
    Thank you for your superb answer! It seems to work perfectly fine. I checked raster['cell_area'].sum() and it deviates by only 1% (510063456035226.2 / 510065880972871.75 = 0.9999952458344382) from the result I cross-linked in my answer (I guess the deviation is caused by e.g. post-digit differences in e.g. the earth radius). But your approach is much much quicker! ... perfect answer! Thanks!
    – Stücke
    Oct 6, 2021 at 6:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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