Usando indice de Qgis con Python

Si es que alguna vez tienen que asignar un lat long a una id de un polygono, una forma de hacerlo procesando por linea es usando la API de Qgis

from qgis.core import *

QgsApplication.setPrefixPath("/usr", True)
QgsApplication.initQgis()


vlayer = QgsVectorLayer( "polygons.shp", "temp_layer", "ogr" )
provider = vlayer.dataProvider()
#main dict to get precalculated values from layer
id2feat = {}
#spatial index to optimize lookup
index = QgsSpatialIndex()
#retrive just the first field from the map layer
provider.select([ 0 ])

#temporal feature
feat = QgsFeature()
while provider.nextFeature( feat ):
    index.insertFeature( feat )
    id2feat[ feat.id() ] = str( feat.attributeMap()[ 0 ].toString() )

for line in sys.stdin:
    fields = line.strip().split(',')
    #read longitude and latitude values
    if len( fields[ 0 ] ) and len( fields[ 1 ] ):
        longitude = float( fields[ 0 ] )
        latitude = float( fields[ 1 ] )
    else:
        print ','.join( fields + [ '' ] )
        continue
    #query geometries
    rect = QgsRectangle(longitude,latitude,longitude,latitude)
    point = QgsGeometry.fromPoint( QgsPoint( longitude, latitude ) )
    #query the index for bbox intersections
    intersects = index.intersects( rect )
    found = False

    for i in intersects:
        provider.featureAtId( i, feat )
        if point.intersects( feat.geometry() ):
            print ','.join( fields + [ id2feat[ i ] ] )
            ##IN CASE WE CARE ABOUT JUST ONE INTERSECTION
            found = True
            break
    if not found:
        print ','.join( fields + [ '' ] )