[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