Search code examples
pythonmatlabmathgeometryarea

Python calculate surface area bounded by the parallels lat1 and lat2 and the meridians lon1 and lon2

I'm looking for an equivalent to this Matlab function areaquad in Python.

area = areaquad(lat1,lon1,lat2,lon2) returns the surface area bounded by the parallels lat1 and lat2 and the meridians lon1 and lon2. The output area is a fraction of the unit sphere's area of 4π, so the result ranges from 0 to 1.

Do you happen to know how I can calculate this using Python? Any tips are much appreciated!

Solution

  • It's not that hard to calculate it, if you need it only for spheres. Here is the formula for calculating total area of sphere surface between two latitudes. So:

    h = sin(lat2)-sin(lat1)
Az = 2 * pi * h

    Now, we can simply find the proportion of the region restricted between two longitudes:

    Aq = Az * (lon2-lon1)/(2*pi)

    And finally, to make the result a fraction of surface of unit sphere, divide it by 4*pi. Putting all together, with simplification and take angle units into account:

     A = (sind(lat2)-sind(lat1)) * deg2rad(lon2-lon1) / (4*pi);

    Hope you can translate it to python.

    enter image description here

    Edit: This is what I get in R2020b for your test case:

    lat1 = -90; lat2 = -89; lon1 = -180; lon2 = 180;
A = (sind(lat2)-sind(lat1)) * deg2rad(lon2-lon1) / (4*pi)

    A =
    7.6152e-05

    enter image description here

    Also, about Aq not being present in the final formula:

    h = sin(lat2)-sin(lat1)
Az = 2 * pi * h
Aq = Az * (lon2-lon1)/(2*pi)
   = 2*pi*h*(lon2-lon1)/(2*pi) // reducing 2*pi
   = h * (lon2-lon1)
   = (sin(lat2)-sin(lat1))*(lon2-lon1)
A  = Aq/(4*pi)
   = (sin(lat2)-sin(lat1))*(lon2-lon1)/(4*pi)