Calling get_GeodesicArea from ogr2ogr

Question

I am trying to find the area of some geometries in France using ogr2ogr. The spread of these geometries span multiple UTM zones (see Réunion, Guyane).

I have seen get_GeodesicArea, which seems to do exactly what I need:

Get the area of the (closed) curve, considered as a surface on the underlying ellipsoid of the SRS attached to the geometry.

It also says:

Note that geometries with circular arcs will be linearized in their original coordinate space first, so the resulting geodesic area will be an approximation.

This is ok for my usecase.

I have also found a reference to OGR2SQLITE_ST_GeodesicArea.

(I am aware that to calculate the area of geometries in different UTM zones, I could just calculate the UTM zone for each as described here, reproject using the appropriate EPSG code, and then use ST_Area. I am hoping to find a less convoluted solution.)

Reproducible setup

Here is franceTest.geojson:

{
"type": "FeatureCollection",
"name": "franceTest",
"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:OGC:1.3:CRS84" } },
"features": [
{ "type": "Feature", "properties": { }, "geometry": { "type": "Polygon", "coordinates": [ [ [ -53.27917836, 4.34657742 ], [ -53.65825646, 3.86953154 ], [ -52.93198534, 3.29672455 ], [ -53.27917836, 4.34657742 ] ] ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Polygon", "coordinates": [ [ [ 55.52739964, -21.00728513 ], [ 55.44252555, -21.09369598 ], [ 55.65702552, -21.13256448 ], [ 55.52739964, -21.00728513 ] ] ] } }
]
}

Attempts

I tried the following:

ogr2ogr franceOut.geojson franceTest.geojson -sql "SELECT ST_GeodesicArea(geometry) as area FROM franceTest"

Output: ERROR 1: Undefined function 'ST_GeodesicArea' used.

Then I tried to use SQLite dialect:

ogr2ogr franceOut.geojson franceTest.geojson -dialect SQLite -sql "SELECT ST_GeodesicArea(geometry) as area FROM franceTest"

Output: ERROR 1: In ExecuteSQL(): sqlite3_prepare_v2(SELECT ST_GeodesicArea(geometry) as area FROM franceTest): no such function: ST_GeodesicArea

Then I tried to use OGR SQL dialect:

ogr2ogr franceOut.geojson franceTest.geojson -dialect OGRSQL -sql "SELECT ST_GeodesicArea(geometry) as area FROM franceTest"

Output: ERROR 1: Undefined function 'ST_GeodesicArea' used.

Then I tried, without being aware of any reasonable justification:

ogr2ogr franceOut.geojson franceTest.geojson -sql "SELECT GeodesicArea(geometry) as area FROM franceTest"

Output: ERROR 1: Undefined function 'GeodesicArea' used.

How can I call get_GeodesicArea, or equivalent function using ogr2ogr?

Answer 1

The new feature was added with this commit https://github.com/OSGeo/gdal/pull/9552.

By the comment three options for using the use_ellipsoid option were implemented:

• Add OGR_G_GeodesicArea() / OGRSurface::get_GeodesicArea()
• SQLite SQL dialect: implement ST_Area(geom, use_ellipsoid)
• GPKG: implement ST_Area(geom, use_ellipsoid)

The SQLite dialect option works for all data sources. With GeoPackage source it is not necessary to select the SQL dialect but the function can be used directly. Usage example with your data:

ogrinfo -dialect SQLite -sql "select ST_Area(geometry,1) from franceTest" franceTest.geojson
INFO: Open of `franceTest.geojson'
      using driver `GeoJSON' successful.

Layer name: SELECT
Geometry: None
Feature Count: 2
Layer SRS WKT:
(unknown)
ST_Area(geometry,1): Real (0.0)
OGRFeature(SELECT):0
  ST_Area(geometry,1) (Real) = 3461416131.06739

OGRFeature(SELECT):1
  ST_Area(geometry,1) (Real) = 125614718.574345

Usage example with a GeoPackage and data in EPSG:4326:

ogrinfo -sql "select ST_Area(geom) from countries limit 1" natural_earth.gpkg
Layer name: SELECT
Geometry: None
Feature Count: 1
Layer SRS WKT:
(unknown)
ST_Area(geom): Real (0.0)
OGRFeature(SELECT):0
  ST_Area(geom) (Real) = 33.3893521404266


ogrinfo -sql "select ST_Area(geom,1) from countries limit 1" natural_earth.gpkg

Layer name: SELECT
Geometry: None
Feature Count: 1
Layer SRS WKT:
(unknown)
ST_Area(geom,1): Real (0.0)
OGRFeature(SELECT):0
  ST_Area(geom,1) (Real) = 388994003560.589