[gdal-dev] OSR Python: nearest polygon distance from each point feature?

[Even Rouault]

Le samedi 22 décembre 2012 17:17:38, David Verbyla a écrit :
> > I am a newbie to OGR/OSR using Python 2.7.2, please excuse my ignorance:
> > I am trying to find the distance to the closest polygon for each point.
> > Is there a more efficient way than simply computing the distance to all
> 
> > polygons for each point?  For example:
> nPts = pointsLayer.GetFeatureCount()
> nPolys = polysLayer.GetFeatureCount()
> 
> for pt in range(0,nPts):
>       minDist = 1000000
>       pointFeature = pointsLayer.GetFeature(pt)
>       for poly in range(0,nPolys):
>             polyFeature = polysLayer.GetFeature(poly)
>             Dist = pointFeature.GetGeometryRef().Distance(
> polyFeature.GetGeometryRef() )
>             if (Dist < minDist):
>                   minDist = Dist
>       print (pt,minDist)

Several points to keep in mind : GetFeature() is not always efficiently 
implemented depending on the capability of the format and/or driver. For 
shapefiles, it is OK, but for some other formats, GetFeature() might be slow. 
If you need to iterate over all the features of a layer, using 
GetNextFeature() is recommanded.

If the layer you're querying has a spatial index ( PostGIS database with 
spatial index, shapefiles with .qix files), you could try to set a spatial filter 
around the target point. But you might need to have a strategy where you start 
with a small radius and increase it until the spatial query isn't empty.

> 
> Thank you.
> 
> Dave Verbyla