[GIS] Iterating over several layers and running process using PyQGIS

intersectionpyqgisqgis-processing

I am trying to find road intersections within each of my memory Polyline layers (which are named road_layer_i for i in (0-99)) so this is what I tried:

import processing

self=qgis.utils
layers = self.iface.mapCanvas().layers()

for layer in layers:
   for i in range(0,100):
      if layer.name()=='road_layer_%s' % (i) :
        processing.runalg('qgis:lineintersections', layer, layer, 'id', 'id', 'C:/Phoenix/intersection_%s.shp' % (i))
        IntersectionLayer=QgsVectorLayer('C:/Phoenix/intersection_%s.shp' % (i) , 'intersection_%s' % (i), "ogr")
        QgsMapLayerRegistry.instance().addMapLayers([IntersectionLayer])

I don't get any error message but each layer created contains points that are the road intersections for the layer 'road_layer_0'.

I really don't understand why it does that so I would be grateful if someone could tell me what I did wrong.

Also, whenever I simply try 

import processing

self=qgis.utils
layers = self.iface.mapCanvas().layers()

for layer in layers:
  if layer.name()=='road_layer_33' :
    processing.runalg('qgis:lineintersections', layer, layer, 'id', 'id', 'C:/Phoenix/intersection_33.shp')
    IntersectionLayer=QgsVectorLayer('C:/Phoenix/intersection_33.shp' , 'intersection', "ogr")    
    QgsMapLayerRegistry.instance().addMapLayers([IntersectionLayer])

I still get the intersections for the road_layer_0…

After answer by @gcarrillo it is still not working. I indeed obtain 100 layers labbeled interesection_i (i in 0-99) but they all contain only the intersections for road_layer_0.
Here is the code I used to create my road layers :

from fTools.fTools import fToolsPlugin
import shapely
from qgis.core import QgsRasterLayer, QgsCoordinateReferenceSystem, QgsGeometry


RoadLayer=QgsVectorLayer("C:/Phoenix/Project/Phoenix/Phoenix_roads.shp" ,    "Phoenix_roads", "ogr")
RoadLayer.setCrs( QgsCoordinateReferenceSystem(3857, QgsCoordinateReferenceSystem.EpsgCrsId) )
#QgsMapLayerRegistry.instance().addMapLayers([RoadLayer])
#Get the data from this layer
provider = RoadLayer.dataProvider()
fields = RoadLayer.pendingFields()
field_names = [field.name() for field in fields]


# create list of points used to create rectangles 
partit=partition(10,c,d)
for i in range(0,100):
#    create new layer in which we will put the roads that are inside the rectangle
   layer =  QgsVectorLayer('LineString', 'road_layer_%s' %(i) , "memory")
   pr = layer.dataProvider()
   attrib_layer=input_layer_attrib_names = RoadLayer.dataProvider().fields()
   oldattributeList = RoadLayer.dataProvider().fields().toList()
#Generate attributes in destination layer
   newattributeList=[]
     for attrib in oldattributeList:
        if layer.fieldNameIndex(attrib.name())==-1:
           newattributeList.append(QgsField(attrib.name(),attrib.type()))
     pr.addAttributes(newattributeList)
  #    add empty attribute fields to new layer
     layer.updateFields()
  #create points left-bottom and right-top of rectangle we want to select
     x=partit[2*i]
     pointA = QgsPoint(x[0],x[1])
     y=partit[2*i+1]
     pointB = QgsPoint(y[0],y[1])
   #   create rectangle
     Rect=QgsRectangle(x[0],x[1],y[0],y[1])  
   #   create request
     request = QgsFeatureRequest().setFilterRect(Rect)
     for f in RoadLayer.getFeatures(request):
    #        add geometry from road to new road
         geom=f.geometry()
         coord=geom.asPolyline()
         pt = QgsFeature()
         pt.setGeometry(QgsGeometry.fromPolyline(coord))
    #       add attributes
         atr = dict(zip(field_names, f.attributes()))
         pt.setAttributes(f.attributes())
     #        add new feature to layer
         pr.addFeatures([pt])
         layer.updateExtents()
    #        add layer to map canva
      QgsMapLayerRegistry.instance().addMapLayers([layer])

Is there anything wrong with this method that makes the use of processing.runalg('qgis:lineintersections', … ) impossible to iterate over the layers thus created ?

Best Answer

I assume your road_layer_i (for i in [0-99]) layers are not the same layer, that is, they do not have exactly the same lines.

Having said that, I've modified your code snippet a bit, making these changes:

  • Use QgsMapLayerRegistry.instance().mapLayers().values() instead of self.iface.mapCanvas().layers(). This ensures you'll iterate on all layers and not only on those that are visible.

  • Instead of checking matching names, extract the i from the layer name.

  • Load all resulting intersection layers at once instead of one by one.

This is the resulting code:

import processing
layers = QgsMapLayerRegistry.instance().mapLayers().values()
outputLayers=[]

for layer in layers:
    i = layer.name().split("_")[-1]
    processing.runalg('qgis:lineintersections', layer, layer, 'id', 'id', 'C:/Phoenix/intersection_%s.shp' % (i))
    outputLayers.append(QgsVectorLayer('C:/Phoenix/intersection_%s.shp' % (i) , 'intersection_%s' % (i), "ogr"))

QgsMapLayerRegistry.instance().addMapLayers( outputLayers )

I've checked it on GNU/Linux, QGIS 2.6.1.

Related Question