-
Notifications
You must be signed in to change notification settings - Fork 560
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Polygon area calculation doesn't match other tools like google earth engine or geojson.io #726
Comments
As far as I understand, shapely assumes all points are cartesian, and calculates area in whatever coordinates it is given; in the example above wgs coordinates are converted to web mercator If we substitute conversion to an equal area coordinate system such as world mollweide (https://spatialreference.org/ref/esri/54009/) , the area calculation is performed as expected: coords = [(-97.59238135821987, 43.47456565304017),
(-97.59244690469288, 43.47962399877412),
(-97.59191951546768, 43.47962728271748),
(-97.59185396090983, 43.47456565304017),
(-97.59238135821987, 43.47456565304017)]
projection = partial(pyproj.transform, pyproj.Proj(init='epsg:4326'), pyproj.Proj(init='esri:54009'))
shapely.ops.transform(projection, shapely.geometry.Polygon(coords)).area Which gives We use this over in https://github.com/mapbox/fio-area. Hope this helps! |
Shapely only does Cartesian calculations, so this issue is mostly off topic. I've just fixed the broken link in There are several different ways to calculate areas from geodetic coordinates, and they will each give a slightly different result. To add another method, use Planimeter on WGS84 for the example to get 23988.9 m^2 with an error about 0.2 m^2, so I'd call this the "best answer". This method is also used in PostGIS on geography types. |
I appreciate the input from everyone. I'm pretty new to the GIS world so it's all helpful, thank you! |
@dnomadb when trying the method suggested I got the following error:
Have you come across this? |
@brentonmallen1 it looks like you have |
That was a typo on my part trying different things to see if I could get it to work. I still get the error with projection = partial(pyproj.transform, pyproj.Proj(init='epsg:4326'), pyproj.Proj(init='esri:54009'))
|
In case anyone else comes across the above issue w.r.t. re-projecting to ESRI:54009, here is a solution pyproj4/pyproj#341 |
Polygon area calculation doesn't seem to match up with other similar tools such as Geojson.io or Google Earth Engine. For example:
yields an area of
45573.993884405005
m^2 while Google Earth Engine and Geojson.io report23944.14737277293
and23997.77
, respectively.The area in shapely is calculated here quoting a method whose link is broken, while the Geojson.io is calculated here via a method described on starting on page 3 of this document.
I've converted this briefly to python with the argument expecting the same coordinate scheme as above (list of coordinates):
This method returns an area of
23997.769250450423
m^2.My main question is: Is the current method preferred for some reason? If so, can the link to the cited be fixed (i'm not sure where it's supposed to go). If not, is it possible to conform (or add the option) to the method described int he JPL doc?
I'd be happy to contribute to this if the option for this other method would like to be added but I didn't see any contributing guidelines.
The text was updated successfully, but these errors were encountered: