Using OGR_GEOM_AREA to add area attribute to shapefile?

Question

I've been trying to batch convert raster files to polygons then delete areas <100m2 before clipping the raster to the filtered polygons.

I have referred to both the below posts:

• Is there any way to add an area column to new shapefile created from raster with gdal's polygonize in the same step?
• https://stackoverflow.com/questions/42355370/sqlite-query-on-path-as-directory-ogr2ogr

I have tried to use the desired command in its simplest form:

ogr2ogr ouput.shp input.shp -sql "SELECT *, OGR_GEOM_AREA AS area FROM input"

Despite this I still get the error "ERROR 1: SELECT from table input failed, no such table/featureclass."

I know that the input.shp exists and is in the same directory as my batch file, so I can only assume that the OGR_GEOM_AREA can not be found.

I notice that the gdal_polygonize step I use prior writes the output in GML format: "Creating output input.shp of format GML."

Could this be the source of the issue?

If so I have not been able to find how to change the polygonize output format to something more ogr friendly.

I did some tests and can confirm that the issue is definitely the GML format being output by gdal_polygonize

Answer 1

Your issue is indeed with the input data. You have created a GML file that has the extension .shp in the name but it does not mean that the examples using real shapefiles as input would work out of the box.

Shapefiles have always just one layer and GDAL shapefile driver is using the basename as a layer name. Thus the one layer in "input.shp" is "input". GML file can have any number of layers and they can have whatever names. You can't guess what the name is but you must check that with ogrinfo. I did a test with gdal_polygonize and here is what ogrinfo reports:

ogrinfo output.shp
INFO: Open of `output.shp'
      using driver `GML' successful.
1: out (Polygon)

Now this ogr2ogr command works with the layer name "out"

ogr2ogr output_with_area.shp input.shp -sql "SELECT *, OGR_GEOM_AREA AS area FROM out"

The usage of gdal_polygonize is documented in https://gdal.org/programs/gdal_polygonize.html. Use "-f" for changing the output format, for example

gdal_polygonize input.tif -f "ESRI Shapefile" output.shp

My GDAL 3 version seems to be more clever than your version and it is guessing the outputformat by the filename extension.

gdal_polygonize input.tif output.shp
Creating output output.shp of format ESRI Shapefile.
0...10...20...30...40...50...60...70...80...90...100 - done.