3

Someone asked a nearly identical question about selecting features in a layer by location here Selecting features within polygon from another layer using R or Python?, but the answers are not for general gdal/ogr python users. I'm wondering if anyone has a pythonic solution to this problem that does not require qgis modules. In particular, I'd like to be able to specify the options qgis offers (include input features that touch and/or overlap/cross and/or are completely within and/or are selected).

I've been trying all day to remove polygons from a shapefile that share a border with (but are contained completely within) a polygon in another shapefile. This task takes <5min in QGIS, but after 5+ hours I still haven't found a way to "select by location" and specify "options" using just ogr (I've tried using layer.SetSpatialFilter(polygon_geometry), but the docs specify that this only works with GetNextFeature(), which I am not using to loop through my features, so I was not terribly surprised when it did not produce the desired results. Also, I suspect that features sharing only a border with the selection layer are included rather than excluded). I've also tried ogr2ogr

ogr2ogr -f 'ESRI Shapefile' -dialect sqlite -sql "SELECT * FROM bergs WHERE ST_Intersects(bergs.geometry, AOI_buffered100.geometry)" -overwrite bergs_testborderrem.shp bergs.shp

but cannot get a working SQL statement (I keep getting no such column errors, but all the ogr-sql examples I've found use .geometry, so I have no idea how to query on the geometry otherwise).

Does anyone have any suggestions?

2

Maybe something like

ogr2ogr -f 'ESRI Shapefile' -dialect sqlite -sql "SELECT a.* FROM bergs a, 
AOI_buffered100 b WHERE ST_Intersects(a.geometry, b.geometry)" -overwrite 
bergs_testborderrem.shp sub_directory

where the shapefiles are in sub_directory. Also you may want ST_Touches or something instead of ST_Intersects.

The key seems to be the specification of the directory as the source, though relative directories appear to be possible, at least for the second shapefile, as long as there are shapefiles in sub_directory. Finally, make sure that the name of the geometry columns are correct.

Also, see here.

  • This solution looks promising (thanks!), but I'm struggling to implement it because my shapefiles are organized such that they're not all in one directory. AOI_buffered100.shp is a generic shapefile in one location and the bergs.shp is in another sub directory. This is part of a larger processing and I don't want to have to copy AOI_buffered100.shp into each sub directory to run ogr2ogr. Is there a way to direct ogr to the multiple input files that I'm just not seeing? I'm also getting an "unrecognized token" error on "bergs" despite having copied the layer name directly from ogrinfo output – Jessica Dec 18 '17 at 15:26
  • Figured out the "unrecognized token" error. I forgot about the fact that shapefiles can't start with a number (the file has a different name on my system). – Jessica Dec 18 '17 at 15:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.