[GIS] Creating new layer in PyQGIS with same attributes of another preexisting and copy some elements in it

pyqgis

I'm new in QGIS and I'm trying to create my first plugin. I have a layer containing lots of elements of type Point, Line and Polygon. Suppose I want to create two new layers with the polygons above and below 15 vertices. I want the new two layers have the same attributes as the original one (CRS, colours, etc.) The general structure can be seen as

#original layer
layer = self.iface.activeLayer()
#HERE I WANT TO CREATE TWO EMPTY LAYERS WITH THE SAME ATTRIBUTES AS layer
layerAbove15 = NEW_LAYER(layerAbove15Name,same attributes as layer)
layerBelow15 = NEW_LAYER(layerBelow15Name,same attributes as layer)
#get layer features
feat = layer.getFeatures()
#loop over all elements
for elem in feat:
    #extract geometry
    geom = elem.geometry()
    #check if the element is a polygon
    if geom.type()==QGis.Polygon:
        #element as polygon
        poly = geom.asPolygon()
        #check the number of vertices of the outer ring
        if len(poly[0])>=15:
            #HERE I WANT TO COPY elem TO THE layerAbove15 LAYER
            layerAbove15.ADD(elem)
        else:
            #HERE I WANT TO COPY elem TO THE layerBelow15 LAYER
            layerBelow15.ADD(elem)

Is my pseudocode correct?

Best Answer

I have a layer containing lots of elements of type Point, Line and Polygon

It is impossible in QGIS, look at how to open layers with multiple geometry types in qgis? ( only points layers, polylines layers or polygons layers and multi-parts geometries layers, MultiPoints, etc.).

So, no need of if geom.type()==QGis.Polygon, unless the geometry type is QGis.MultiPolygon.

The correct syntax is (Geometry Handling):

a) Before iterating over the layer elements

layer = qgis.utils.iface.activeLayer()
if layer.wkbType()==QGis.WKBPolygon:
    print 'Layer is a polygon layer'
if layer.wkbType()==QGis.WKBMultiPolygon:
    print 'Layer is a multi-polygon layer'

b) during the iteration

print geom.type() # and the result is a value from QGis.GeometryType enumeration
print geom.wkbType() # same

Example in the Python console (different in a Plugin) with a simple polygon layer :

enter image description here

1) create the two empty shapefiles, layerAbove15 and layerBelow15

layer = qgis.utils.iface.activeLayer()
# source of the layer
provider = layer.dataProvider()
# creation of the shapefiles:
writer1 = QgsVectorFileWriter("layerAbove15Name.shp", "CP1250", provider.fields(), provider.geometryType(), provider.crs(), "ESRI Shapefile")
writer2 = QgsVectorFileWriter("layerBelow15Name.shp", "CP1250", provider.fields(), provider.geometryType(), provider.crs(), "ESRI Shapefile")

2) populate the shapefiles depending on the condition:

outelem = QgsFeature()
# iterating over the input layer
for elem in layer.getFeatures():
    if len(elem.geometry().asPolygon()[0])-1 >=15:
             outelem.setGeometry(elem.geometry() )
             outelem.setAttributes(elem.attributes())
             writer1.addFeature(outelem)
    else:
             ....
             writer2.addFeature(outelem)

3) with a simple function to avoid repetition:

def write(out,inp):
    out.setGeometry(inp.geometry() )
    out.setAttributes(inp.attributes())
    return out

for elem in layer.getFeatures():
    if len(elem.geometry().asPolygon()[0])-1 >=15:
        writer1.addFeature(write(outelem,elem))
    else:
        writer2.addFeature(write(outelem,elem))

enter image description here

If the layer has multi-parts geometries (MultiPolygons), it is necessary to adapt the script consequently.

Related Question