[gdal-dev] proper way to get area of a lon/lat polygon

Alan Snow alansnow21 at gmail.com
Fri Sep 18 20:43:27 PDT 2020


EPSG 9835 is an conversion method which needs some parameters defined:
https://pyproj4.github.io/pyproj/stable/api/crs/coordinate_operation.html#lambertcylindricalequalareaconversion

You can use that conversion to build a projected CRS:
https://pyproj4.github.io/pyproj/stable/build_crs.html#projected-crs

That being said, you might also be interested in geodesic area if you have
geographic coordinates:
https://pyproj4.github.io/pyproj/stable/examples.html#geodesic-area

And you can use this to get the Geod class for EPSG:4326:
https://pyproj4.github.io/pyproj/stable/api/crs/crs.html#pyproj.crs.CRS.get_geod

Hopefully this helps,
Alan


On Fri, Sep 18, 2020, 9:31 PM Mikhail Kruk <meshko at gmail.com> wrote:

> Yes, you are correct, that was one of my first attempts and was wrong
> because I don't really understand these projection systems.
> I now have a better one which sort of works and produces areas in meters
> by setting first system to
> sr.SetWellKnownGeogCS("EPSG:4326")
> Now I am struggling with selecting the projection coordinate system.  Some
> suggested that I use cylindrical projection EPSG:9835 which should
> preserve the area anywhere on the planet but when I try
> sr_proj.ImportFromEPSG(9835)
> I get error:
> PROJ: proj_create_from_database: crs not found
>
>
> On Fri, Sep 18, 2020 at 9:50 PM Michael Patrick <geodesy99 at gmail.com>
> wrote:
>
>> ... but basically I came to the conclusion that when people talk about an
>>> area of something on Earth they really mean the area of a projection (am i
>>> right about that?) and so i need to specify a projection before getting a
>>> meaningful area.  So i ended up with this code:
>>>
>>> wkt = "POLYGON ((lon1 lat1, lon2 lat2, lon3 lat3))"
>>> poly = ogr.CreateGeometryFromWkt(wkt)
>>> sr = osr.SpatialReference()
>>> sr.SetWellKnownGeogCS("NAD83")
>>> sr.SetProjCS("SRS_PT_TRANSVERSE_MERCATOR")
>>> poly.AssignSpatialReference(sr)
>>> poly.GetArea() <https://lists.osgeo.org/mailman/listinfo/gdal-dev>
>>>
>>> I might be dense, but it seems you are assigning sr as a specific
>> transverse mercator, i.e.  UTM zone 18N, or a custom one centered on the
>> centroid of your polygon. i.e. you don't have sufficient arguments.
>>
>>
>>
>> <https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail&utm_term=icon> Virus-free.
>> www.avast.com
>> <https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail&utm_term=link>
>> <#m_-5507655572585127822_m_8065453534991225343_m_-2581199688950320765_m_6579958654730313444_DAB4FAD8-2DD7-40BB-A1B8-4E2AA1F9FDF2>
>> _______________________________________________
>> gdal-dev mailing list
>> gdal-dev at lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/gdal-dev
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20200918/6e162ca8/attachment.html>


More information about the gdal-dev mailing list