I need to tell whether geometries from a shapefile are inside other geometries or not.
For example, I want to test if there are trees (represented as points in `trees.shp`) inside the urban areas (reprensented as polygons in `area.shp`)

What's the easiest way to proceed?  

I loop over the points, and within this loop, I have a nested loop over the polygons.
Then, I test if my point geometry is inside the polygon with the `Within(Geometry geo)` method from OGR.
It's very long, not optimized at all.

What should I do to fasten the computation?

    def PointsInPolygons(self, pointsLayer, polysLayer):
        nbInside = 0
        total = pointsLayer.GetFeatureCount()
        for i in range(total):
            pointFeat = pointsLayer.GetFeature(i)
            pointGeo = pointFeat.GetGeometryRef()
            
            for j in range(polysLayer.GetFeatureCount()):
                polyFeat = polysLayer.GetFeature(j)
                polyGeo = polyFeat.GetGeometryRef()
                
                if pointGeo.Within(polyGeo):
                    nbInside+=1
                    break