[gdal-dev] SetSpatialFilter not working?

Frank Warmerdam warmerdam at pobox.com
Thu Mar 4 19:11:53 EST 2010

Spencer Gardner wrote:
> I'm writing a script with 2 primary layers, a road layer and an intersection
> layer.  I have designed the script to cycle through the intersections and
> remove any that have less than 3 road segments connected to it.  I have set
> a buffer of 1 foot to account for any misalignment.  The script runs just
> fine (using almost 24,000 intersection nodes) but doesn't identify any
> intersections for removal.
> After troubleshooting, I have found that counting the road features after
> setting the spatial filter returns an unusually high number of road segments
> (usually between 9 and 12) for each intersection so that no intersection
> returns less than 3 segments.  The same phenomenon occurs if I remove the
> 1-foot buffer.  If I test the layers in ArcView by selecting roads that
> intersect a specific node, the result is correct.  Is there something in my
> code that could be causing this problem?  The documentation says the
> SetSpatialFilter method "may be inaccurately implemented" - would this be a
> symptom of this problem?  Is there another way of doing this type of
> analysis?
> Thanks for the help.  Please excuse my lack of knowledge as I'm fairly new
> to this.


If GEOS is linked into GDAL/OGR the spatial filtering should be accurate.
Otherwise a bounding box comparison is done and you can get significant
numbers of false positives depending on the structure of the data.

So, my suggestion is to rebuild GDAL with GEOS support.

Best regards,
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | Geospatial Programmer for Rent

More information about the gdal-dev mailing list