[GRASS-dev] area calculations in several GIS
Markus Metz
markus.metz.giswork at gmail.com
Tue Oct 2 10:55:28 PDT 2018
On Tue, Oct 2, 2018 at 7:30 PM Charles Karney <charles at karney.com> wrote:
>
> Sorry, I'm coming late to this conversation. The area of the polygon
> in the Baltic Sea posted by Kristian is
>
> 14737935340.098511 m^2
>
> assuming the WGS84 ellipsoid. This is accurate to 1 square-mm and was
> obtained with the MPFR-enabled version of GeographicLib's Planimeter.
>
> This is consistent with the result reported by the PROJ.4
> geod_polygonarea() and it agrees with the result I get with the regular
> Planimeter utility. I don't why someone else is getting a different
> result from Planimeter.
I can confirm this, the results of Planimeter and PROJ 5.2.0 are identical.
Some more confusion:
I created simple boxes for the test polygon, one a bit larger, one a bit
smaller than the test polygon
GRASS: 13,222.778
Planimeter: 13,221.965
geod_polygonarea(): 13,221.965
14.62569 55.36254
14.6256944444444 54.0786
17.1177 54.0786111111111
17.1177777777778 55.3625
GRASS: 22,950.510
Planimeter: 22,946.901
geod_polygonarea(): 22,946.901
GRASS native geodesic area is not too far off from GeographicLib/PROJ
With the test polygon
Planimeter/PROJ: 14,737.935 km^2
GRASS GIS: 14.718.098 km^2
GRASS is much farther off, and here smaller instead of larger than
Planimeter/PROJ.
I will have a look.
Markus M
> It's probably operator error -- but who knows?
>
> --Charles
>
> On 10/2/18 10:08 AM, Kristian Evers wrote:
> > Sorry, I don’t know. It’s possible there’s a bug somewhere.
> >
> > Charles, do you have any insights here?
> >
> > /Kristian
> >
> > *Fra:*Markus Metz <markus.metz.giswork at gmail.com>
> > *Sendt:* 2. oktober 2018 15:01
> > *Til:* Kristian Evers <kreve at sdfe.dk>
> > *Cc:* GRASS developers list <grass-dev at lists.osgeo.org>; Helmut
> > Kudrnovsky <hellik at web.de>
> > *Emne:* Re: [GRASS-dev] area calculations in several GIS
> >
> >
> >
> > On Tue, Oct 2, 2018 at 12:08 PM Kristian Evers <kreve at sdfe.dk
> > <mailto:kreve at sdfe.dk>> wrote:
> > >
> >
> > > Use GeographicLib/PROJ as the reference.
> >
> > Which one? The results of Planimeter from GeographicLib and
> > geod_polygonarea() from PROJ are different. If in doubt, use Planimeter
> > as reference?
> >
> > Markus M
> >
> > > I would say that this polygon is large, yes. If you want to get a
> > better sense for what constitutes a large polygon in this sense, try
> > creating a series of polygons of varying size (e.g. polygons of your
> > house, your city, your county, your country and your continent) and
> > calculate the area with both GRASS and GeographicLib. At some point the
> > calculated areas should start to diverge significantly. I suspect it
> > happens somewhere between county and country.
> >
> > >
> > >
> > >
> > > /Kristian
> > >
> > >
> > >
> > > Fra: Markus Metz <markus.metz.giswork at gmail.com
> > <mailto:markus.metz.giswork at gmail.com>>
> > > Sendt: 2. oktober 2018 11:51
> > > Til: Kristian Evers <kreve at sdfe.dk <mailto:kreve at sdfe.dk>>
> > > Cc: GRASS developers list <grass-dev at lists.osgeo.org
> > <mailto:grass-dev at lists.osgeo.org>>; Helmut Kudrnovsky <hellik at web.de
> > <mailto:hellik at web.de>>
> > > Emne: Re: [GRASS-dev] area calculations in several GIS
> > >
> > >
> > >
> > >
> > >
> > > On Tue, Oct 2, 2018 at 8:47 AM Kristian Evers <kreve at sdfe.dk
> > <mailto:kreve at sdfe.dk>> wrote:
> > > >
> > > > Markus,
> > > >
> > > >
> > > > Thanks for extending the list. I do wonder why Planimeter gives
> > different results than geod_polygonarea(). I’ll run that by Charles
> > Karney when I get a chance.
> > > >
> > >
> > > for completeness, I used geod_polygonarea() from proj-5.2.0
> > >
> > >
> > >
> > > >> Regarding areas in projected space: UTM is not area true, laea is
> > but is instead scewing angles. You can check that stuff by running proj
> > in very verbose mode.
> > > >
> > >
> > > [...]
> > >
> > > >
> > >
> > > > So yes, the LAEA is the better choice of the two but it is never
> > going to represent the true area, especially for large polygons, and I
> > would not advice using it as reference for ellipsoidal area
calculations.
> > >
> > >
> > >
> > > The question is, what can be used as reference? And is the test
> > polygon a "large" polygon, causing "large" deviations from the true area
> > when measured in LAEA?
> > >
> > >
> > >
> > > Markus M
> > >
> > > >
> > >
> > > > /Kristian
> > > >
> > > >
> > > >
> > > >
> > > >
> > > > Fra: Markus Metz <markus.metz.giswork at gmail.com
> > <mailto:markus.metz.giswork at gmail.com>>
> > > > Sendt: 1. oktober 2018 23:22
> > > > Til: GRASS developers list <grass-dev at lists.osgeo.org
> > <mailto:grass-dev at lists.osgeo.org>>; Kristian Evers <kreve at sdfe.dk
> > <mailto:kreve at sdfe.dk>>
> > > > Cc: Helmut Kudrnovsky <hellik at web.de <mailto:hellik at web.de>>
> > > > Emne: Re: [GRASS-dev] area calculations in several GIS
> > > >
> > > >
> > > >
> > > > Updated list with area calculations for
> > > >
> > > > https://gist.github.com/kbevers/207b5bcb9be20e7554abe5f56742ec2c
> > > >
> > > >
> > > >
> > > > PROJ [1]: 14,737.935 km^2
> > > > Caris LOTS: 14,737 km^2
> > > > ArcMap: 14,727.446 km^2
> > > > MapInfo: 14,727.352 km^2
> > > > GeoMedia: 14,726.443 km^2
> > > > Planimeter: 14,722.522 km^2
> > > > GRASS GIS: 14.718.098 km^2
> > > >
> > > > EU LAEA [2]: 14,718.098 km^2
> > > > UTM 33 N: 14,707.742 km^2
> > > >
> > > > QGIS 3.2: 14,652.181 km^2
> > > > QGIS 2.8: 14,652.181 km^2
> > > >
> > > > [1] geodesic.h:geod_polygonarea()
> > > > [2] EPSG:3035
> > > >
> > > >
> > > >
> > > > In this case, GRASS GIS provides the best match of geodesic area to
> > metric area.
> > > >
> > > >
> > > >
> > > > @Kristian: are the metric area measurements in "EU LAEA" and "UTM
> > 33 N" suitable as reference?
> > > >
> > > >
> > > >
> > > > Discussion started on
> > > >
> > > >
> >
http://osgeo-org.1560.x6.nabble.com/Re-Qgis-user-New-Features-in-Shape-Tools-3-2-0-td5378898.html
> > > >
> > > >
> > > >
> > > > Markus M
> > > >
> > > >
> > > >
> > > >
> > > >
> > > > On Mon, Oct 1, 2018 at 2:56 PM Markus Metz
> > <markus.metz.giswork at gmail.com <mailto:markus.metz.giswork at gmail.com>>
> > wrote:
> > > >
> > > >
> > > >
> > > > On Tue, Sep 25, 2018 at 7:38 PM Helmut Kudrnovsky <hellik at web.de
> > <mailto:hellik at web.de>> wrote:
> > > > >
> > > > > fyi see
> > > > >
> >
https://lists.osgeo.org/pipermail/qgis-developer/2018-September/054644.html
> > > > >
> > > > > with GRASS mentioned
> > > > > ------------------
> > > > > Kristian Evers:
> > > > >
> > > > > Right, here are the calculated areas as returned by a number of
> > different
> > > > > GIS applications and the planimeter app of GeographicLib for
> > reference:
> > > > >
> > > > > Caris LOTS: 14.737 km^2
> > > > > ArcMap: 14.727,446 km^2
> > > > > MapInfo: 14.727,352 km^2
> > > > > GeoMedia: 14.726,443 km^2
> > > > > Planimeter: 14.722,522 km^2
> > > > > QGIS 3.2: 14.652,181 km^2
> > > > > QGIS 2.8: 14.652,181 km^2
> > > >
> > > >
> > > >
> > > > adding to the confusion:
> > > >
> > > >
> > > >
> > > > I used the geographiclib API as included in PROJ 5.2.0 following
> > the example for geod_polygonarea() in geodesic.h and get
> > > >
> > > > geographiclib: 14,737.935 km^2
> > > >
> > > > quite different from
> > > >
> > > > Planimeter: 14,722.522 km^2
> > > >
> > > >
> > > >
> > > > GRASS native gives 14,718.097679
> > > >
> > > > as also reported by Helmut and Stefan
> > > >
> > > >
> > > >
> > > > Markus M
> > > >
> > > > >
> > > >
> > > > > The polygon that I have used to get the numbers above can be
> > found here:
> > > > > https://gist.github.com/kbevers/207b5bcb9be20e7554abe5f56742ec2c
> > > > >
> > > > > I am quite confident that GeographicLib delivers the most
> > accurate result
> > > > > (if you have doubts, this reference [0] should convince you). As
> > can be seen
> > > > > from the table above all but QGIS come fairly close. I expect
> > some variation
> > > > > in the results as these are numerical approximations, although I
> > think QGIS
> > > > > is too far of the mark. My suspicion is that the geodesic
> > algorithm used by
> > > > > QGIS (and apparently GRASS) is to blame here.
> > > > >
> > > > > /Kristian
> > > > >
> > > > > [0] https://arxiv.org/pdf/1102.1215.pdf
> > > > > -----------------
> > > > >
> > > > >
> > > > >
> > > > >
> > > > >
> > > > > -----
> > > > > best regards
> > > > > Helmut
> > > > > --
> > > > > Sent from:
> > http://osgeo-org.1560.x6.nabble.com/Grass-Dev-f3991897.html
> > > > > _______________________________________________
> > > > > grass-dev mailing list
> > > > > grass-dev at lists.osgeo.org <mailto:grass-dev at lists.osgeo.org>
> > > > > https://lists.osgeo.org/mailman/listinfo/grass-dev
> >
>
> --
> Charles Karney <charles at karney.com>
> Princeton, NJ
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20181002/48354f0e/attachment.html>
More information about the grass-dev
mailing list