[Qgis-developer] Why do I get an empty attribute map when using featureAtId()

Anita Graser anitagraser at gmx.at
Wed Jan 2 06:59:46 PST 2013


Hi,

I'm working on an extended version of Sextante Ftools "Sum lines"
where I want to calculate line length sum for different classes.
When trying to get the classification value from the line features, I
run into the problem that the attribute map is empty. Is this because
of my usage of featureAtId() method? Or am I missing something else
here?

Thanks
Anita

Script here:


import os.path

from PyQt4 import QtGui
from PyQt4.QtCore import *

from qgis.core import *

from sextante.core.GeoAlgorithm import GeoAlgorithm
from sextante.core.QGisLayers import QGisLayers
from sextante.core.SextanteLog import SextanteLog

from sextante.parameters.ParameterVector import ParameterVector
from sextante.parameters.ParameterString import ParameterString
from sextante.parameters.ParameterTableField import ParameterTableField
from sextante.outputs.OutputVector import OutputVector

from sextante.ftools import FToolsUtils as utils

class SumLinesClassified(GeoAlgorithm):

    LINES = "LINES"
    POLYGONS = "POLYGONS"
    CLASS_FIELD = "CLASS_FIELD"
    LEN_FIELD = "LEN_FIELD"
    COUNT_FIELD = "COUNT_FIELD"
    OUTPUT = "OUTPUT"

    def getIcon(self):
        return QtGui.QIcon(os.path.dirname(__file__) + "/icons/sum_lines.png")

    def defineCharacteristics(self):
        self.name = "Sum line lengths by class"
        self.group = "Analysis tools"

        self.addParameter(ParameterVector(self.LINES, "Lines",
ParameterVector.VECTOR_TYPE_LINE))
        self.addParameter(ParameterVector(self.POLYGONS, "Polygons",
ParameterVector.VECTOR_TYPE_POLYGON))
        self.addParameter(ParameterTableField(self.CLASS_FIELD,
"Unique ID field",self.LINES ))
        #self.addParameter(ParameterString(self.LEN_FIELD, "Lines
length field name", "LENGTH"))
        #self.addParameter(ParameterString(self.COUNT_FIELD, "Lines
count field name", "COUNT"))

        self.addOutput(OutputVector(self.OUTPUT, "Result"))

    def processAlgorithm(self, progress):
        lineLayer =
QGisLayers.getObjectFromUri(self.getParameterValue(self.LINES))
        polyLayer =
QGisLayers.getObjectFromUri(self.getParameterValue(self.POLYGONS))
        classFieldName = self.getParameterValue(self.CLASS_FIELD)

        #lengthFieldName = self.getParameterValue(self.LEN_FIELD)
        #countFieldName = self.getParameterValue(self.COUNT_FIELD)

        output = self.getOutputValue(self.OUTPUT)

        polyProvider = polyLayer.dataProvider()
        lineProvider = lineLayer.dataProvider()
        allLineAttrs = lineProvider.attributeIndexes()

        classFieldIndex = lineProvider.fieldNameIndex(classFieldName)
        if polyProvider.crs() != lineProvider.crs():
            SextanteLog.addToLog(SextanteLog.LOG_WARNING,
                                 "CRS warning: Input layers have
non-matching CRS. This may cause unexpected results.")

        # get unique classification values and create a field for each
one of them
        classValues = lineLayer.uniqueValues(classFieldIndex)
        fieldIndizes = []
        lengthValues = []
        fieldList = polyLayer.pendingFields()
        for value in classValues:
            index, fieldList = utils.findOrCreateField(polyLayer,
fieldList, 'LEN_'+value.toString())
            fieldIndizes.append(index)
            lengthValues.append(0)

        writer =
self.getOutputFromName(self.OUTPUT).getVectorWriter(fieldList,
polyProvider.geometryType(), polyProvider.crs())

        spatialIndex = utils.createSpatialIndex(lineProvider)



        allAttrs = polyLayer.pendingAllAttributesList()
        polyLayer.select(allAttrs)

        ftLine = QgsFeature()
        ftPoly = QgsFeature()
        outFeat = QgsFeature()
        inGeom = QgsGeometry()
        outGeom = QgsGeometry()
        distArea = QgsDistanceArea()

        current = 0
        total = 100.0 / float(polyProvider.featureCount())
        hasIntersections = False

        while polyLayer.nextFeature(ftPoly):
            inGeom = QgsGeometry(ftPoly.geometry())
            atMap = ftPoly.attributeMap()
            #count = 0
            #length = 0
            hasIntersections = False
            lines = spatialIndex.intersects(inGeom.boundingBox())

            lineProvider.rewind()
            lineProvider.select( allLineAttrs )

            if len(lines) > 0:
                hasIntersections = True

            if hasIntersections:
                # sum up line lengths per classification value
                for i in lines:
                    lineProvider.featureAtId(int(i), ftLine)
                    SextanteLog.addToLog(SextanteLog.LOG_INFO, 'len of
attribute map: '+str(len(ftLine.attributeMap())))

                    tmpGeom = QgsGeometry(ftLine.geometry())
                    if inGeom.intersects(tmpGeom):
                        outGeom = inGeom.intersection(tmpGeom)
                        classValue =
ftLine.attributeMap()[classFieldIndex].toString()
                        lengthValues[classValues.index(classValue)] +=
distArea.measure(outGeom)
                        #count += 1

            outFeat.setGeometry(inGeom)
            outFeat.setAttributeMap(atMap)
            for newFieldIndex in fieldIndizes:
                outFeat.addAttribute(newFieldIndex,
QVariant(lengthValues[fieldIndizes.index(newFieldIndex)]))
            #outFeat.addAttribute(idxLength, QVariant(length))
            #outFeat.addAttribute(idxCount, QVariant(count))
            writer.addFeature(outFeat)

            current += 1
            progress.setPercentage(int(current * total))

        del writer


More information about the Qgis-developer mailing list