Skip to content

Commit

Permalink
review and cleanup Convex Hull tool, also move some common routines to
Browse files Browse the repository at this point in the history 
utils module and update other algs accordingly
  • Loading branch information
@alexbruy
alexbruy committed Mar 21, 2013
1 parent 11d567c commit 48e9f42
Show file tree
Hide file tree
Showing 3 changed files with 151 additions and 141 deletions.
213 changes: 114 additions & 99 deletions python/plugins/sextante/algs/ftools/ConvexHull.py
View file
Expand Up @@ -23,144 +23,159 @@
# This will get replaced with a git SHA1 when you do a git archive
__revision__ = '$Format:%H$'

from sextante.core.GeoAlgorithm import GeoAlgorithm
from PyQt4.QtCore import *
from PyQt4.QtGui import *

from qgis.core import *
from sextante.parameters.ParameterVector import ParameterVector

from sextante.core.GeoAlgorithm import GeoAlgorithm
from sextante.core.QGisLayers import QGisLayers
from sextante.outputs.OutputVector import OutputVector
from sextante.algs.ftools import FToolsUtils as utils
from sextante.core.SextanteLog import SextanteLog

from sextante.parameters.ParameterVector import ParameterVector
from sextante.parameters.ParameterTableField import ParameterTableField
from sextante.parameters.ParameterSelection import ParameterSelection

from sextante.outputs.OutputVector import OutputVector

from sextante.algs.ftools import FToolsUtils as utils

class ConvexHull(GeoAlgorithm):

INPUT = "INPUT"
OUTPUT = "OUTPUT"
FIELD = "FIELD"
METHOD = "METHOD"
METHODS = ["Create single minimum convex hull", "Create convex hulls based on field"]
METHODS = ["Create single minimum convex hull",
"Create convex hulls based on field"
]

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

def defineCharacteristics(self):
self.name = "Convex hull"
self.group = "Vector geometry tools"
self.addParameter(ParameterVector(ConvexHull.INPUT, "Input layer", ParameterVector.VECTOR_TYPE_ANY))
self.addParameter(ParameterTableField(ConvexHull.FIELD, "Field", ConvexHull.INPUT))
self.addParameter(ParameterSelection(ConvexHull.METHOD, "Method", ConvexHull.METHODS))
self.addOutput(OutputVector(ConvexHull.OUTPUT, "Convex hull"))

def processAlgorithm(self, progress):
useField = (self.getParameterValue(ConvexHull.METHOD) == 1)
field = self.getParameterValue(ConvexHull.FIELD)
vlayerA = QGisLayers.getObjectFromUri(self.getParameterValue(ConvexHull.INPUT))
fieldName = self.getParameterValue(ConvexHull.FIELD)
layer = QGisLayers.getObjectFromUri(self.getParameterValue(ConvexHull.INPUT))

GEOS_EXCEPT = True
FEATURE_EXCEPT = True
vproviderA = vlayerA.dataProvider()
#allAttrsA = vproviderA.attributeIndexes()
#vproviderA.select(allAttrsA)
fields = [QgsField("ID", QVariant.Int),
QgsField("Area", QVariant.Double),
QgsField("Perim", QVariant.Double)]
writer = self.getOutputFromName(ConvexHull.OUTPUT).getVectorWriter(fields, QGis.WKBPolygon, vproviderA.crs())
inFeat = QgsFeature()

index = layer.fieldNameIndex(fieldName)
fType = layer.pendingFields()[index].type()
f = QgsField("value")
f.setType(QVariant.String)
f.setLength(255)
if useField:
if fType == QVariant.Int:
f.setType(QVariant.Int)
f.setLength(20)
elif fType == QVariant.Double:
f.setType(QVariant.Double)
f.setLength(20)
f.setPrecision(6)
else:
f.setType(QVariant.String)
f.setLength(255)

fields = [QgsField("id", QVariant.Int, "", 20),
f,
QgsField("area", QVariant.Double, "", 20, 6),
QgsField("perim", QVariant.Double, "", 20, 6)
]

writer = self.getOutputFromName(ConvexHull.OUTPUT).getVectorWriter(fields, QGis.WKBPolygon, layer.dataProvider().crs())

outFeat = QgsFeature()
inGeom = QgsGeometry()
outGeom = QgsGeometry()
nElement = 0
index = vproviderA.fieldNameIndex(field)

features = QGisLayers.features(vlayerA)
nFeat = len(features)
current = 0

fid = 0
val = ""
if useField:
unique = utils.getUniqueValues(vproviderA, index)
nFeat = nFeat * len(unique)
for i in unique:
hull = []
first = True
outID = 0
#vproviderA.select(allAttrsA)
features = QGisLayers.features(vlayerA)
for inFeat in features:
atMap = inFeat.attributes()
idVar = atMap[ index ]
if idVar.toString().trimmed() == i.toString().trimmed():
if first:
outID = idVar
first = False
inGeom = QgsGeometry(inFeat.geometry())
points = utils.extractPoints(inGeom)
hull.extend(points)
nElement += 1
progress.setPercentage(int(nElement / nFeat * 100))
if len(hull) >= 3:
tmpGeom = QgsGeometry(outGeom.fromMultiPoint(hull))
try:
outGeom = tmpGeom.convexHull()
outFeat.setGeometry(outGeom)
(area, perim) = self.simpleMeasure(outGeom)
outFeat.addAttribute(0, QVariant(outID))
outFeat.addAttribute(1, QVariant(area))
outFeat.addAttribute(2, QVariant(perim))
writer.addFeature(outFeat)
except:
GEOS_EXCEPT = False
continue
unique = layer.uniqueValues(index)
total = 100.0 / float(layer.featureCount() * len (unique))

for i in unique:
hull = []
first = True
features = QGisLayers.features(layer)
for f in features:
idVar = f.attribute(fieldName)
if idVar.toString().trimmed() == i.toString().trimmed():
if first:
val = idVar
first = False
inGeom = QgsGeometry(f.geometry())
points = utils.extractPoints(inGeom)
hull.extend(points)
current += 1
progress.setPercentage(int(current * total))

if len(hull) >= 3:
tmpGeom = QgsGeometry(outGeom.fromMultiPoint(hull))
try:
outGeom = tmpGeom.convexHull()
(area, perim) = utils.simpleMeasure(outGeom)
outFeat.setGeometry(outGeom)
outFeat.setAttributes([QVariant(fid),
QVariant(val),
QVariant(area),
QVariant(perim)
])
#~ outFeat.setAttribute("id", QVariant(fid))
#~ outFeat.setAttribute("value", QVariant(val))
#~ outFeat.setAttribute("area", QVariant(area))
#~ outFeat.setAttribute("perim", QVariant(perim))
writer.addFeature(outFeat)
except:
GEOS_EXCEPT = False
continue
fid += 1
else:
hull = []
#vproviderA.select(allAttrsA)
features = QGisLayers.features(vlayerA)
for inFeat in features:
inGeom = QgsGeometry(inFeat.geometry())
points = utils.extractPoints(inGeom)
hull.extend(points)
nElement += 1
progress.setPercentage(int(nElement / nFeat * 100))
total = 100.0 / float(layer.featureCount())
features = QGisLayers.features(layer)
for f in features:
inGeom = QgsGeometry(f.geometry())
points = utils.extractPoints(inGeom)
hull.extend(points)
current += 1
progress.setPercentage(int(current * total))

tmpGeom = QgsGeometry(outGeom.fromMultiPoint(hull))
try:
outGeom = tmpGeom.convexHull()
outFeat.setGeometry(outGeom)
(area, perim) = self.simpleMeasure(outGeom)
#outFeat.addAttribute(0, QVariant("1"))
#outFeat.addAttribute(1, QVariant(area))
#outFeat.addAttribute(2, QVariant(perim))
writer.addFeature(outFeat)
outGeom = tmpGeom.convexHull()
(area, perim) = utils.simpleMeasure(outGeom)
outFeat.setGeometry(outGeom)
outFeat.setAttributes([QVariant(0),
QVariant("all"),
QVariant(area),
QVariant(perim)
])
#print outFeat.setAttribute("id", QVariant(0))
#print outFeat.setAttribute("value", QVariant("all"))
#print outFeat.setAttribute("area", QVariant(area))
#print outFeat.setAttribute("perim", QVariant(perim))
writer.addFeature(outFeat)
except:
GEOS_EXCEPT = False
GEOS_EXCEPT = False

del writer

if not GEOS_EXCEPT:
SextanteLog.addToLog(SextanteLog.LOG_WARNING, "Geometry exception while computing convex hull")
if not FEATURE_EXCEPT:
SextanteLog.addToLog(SextanteLog.LOG_WARNING, "Feature exception while computing convex hull")

def simpleMeasure(self, inGeom):
measure = QgsDistanceArea()
attr1 = measure.measure(inGeom)
if inGeom.type() == QGis.Polygon:
attr2 = self.perimMeasure(inGeom, measure)
else:
attr2 = attr1
return (attr1, attr2)

def perimMeasure(self, inGeom, measure):
value = 0.00
if inGeom.isMultipart():
poly = inGeom.asMultiPolygon()
for k in poly:
for j in k:
value = value + measure.measureLine(j)
else:
poly = inGeom.asPolygon()
for k in poly:
value = value + measure.measureLine(k)
return value

def defineCharacteristics(self):
self.name = "Convex hull"
self.group = "Vector geometry tools"
self.addParameter(ParameterVector(ConvexHull.INPUT, "Input layer", ParameterVector.VECTOR_TYPE_ANY))
self.addParameter(ParameterTableField(ConvexHull.FIELD, "Field", ConvexHull.INPUT))
self.addParameter(ParameterSelection(ConvexHull.METHOD, "Method", ConvexHull.METHODS))
self.addOutput(OutputVector(ConvexHull.OUTPUT, "Convex hull"))
#=========================================================
50 changes: 8 additions & 42 deletions python/plugins/sextante/algs/ftools/ExportGeometryInfo.py
View file
Expand Up @@ -24,13 +24,19 @@
__revision__ = '$Format:%H$'

from PyQt4.QtCore import *

from qgis.core import *

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

from sextante.parameters.ParameterVector import ParameterVector
from sextante.parameters.ParameterSelection import ParameterSelection

from sextante.outputs.OutputVector import OutputVector

from sextante.algs.ftools import FToolsUtils as utils

class ExportGeometryInfo(GeoAlgorithm):

INPUT = "INPUT"
Expand All @@ -56,7 +62,6 @@ def defineCharacteristics(self):

self.addOutput(OutputVector(self.OUTPUT, "Output layer"))


def processAlgorithm(self, progress):
layer = QGisLayers.getObjectFromUri(self.getParameterValue(self.INPUT))
method = self.getParameterValue(self.METHOD)
Expand Down Expand Up @@ -90,7 +95,7 @@ def processAlgorithm(self, progress):
# 2 - ellipsoidal
if method == 2:
settings = QSettings()
ellips = settings.value("/qgis/measure/ellipsoid", "WGS84").toString()
ellips = QgsProject.instance().readEntry("Measure", "/Ellipsoid", GEO_NONE)[0]
crs = layer.crs().srsid()
elif method == 1:
mapCRS = QGisLayers.iface.mapCanvas().mapRenderer().destinationCrs()
Expand All @@ -109,7 +114,7 @@ def processAlgorithm(self, progress):
if method == 1:
inGeom.transform(coordTransform)

(attr1, attr2) = self.simpleMeasure(inGeom, method, ellips, crs)
(attr1, attr2) = utils.simpleMeasure(inGeom, method, ellips, crs)

outFeat.setGeometry(inGeom)
atMap = inFeat.attributes()
Expand All @@ -123,42 +128,3 @@ def processAlgorithm(self, progress):
progress.setPercentage(int(current * total))

del writer

def simpleMeasure(self, geom, method, ellips, crs):
if geom.wkbType() in [QGis.WKBPoint, QGis.WKBPoint25D]:
pt = geom.asPoint()
attr1 = pt.x()
attr2 = pt.y()
elif geom.wkbType() in [QGis.WKBMultiPoint, QGis.WKBMultiPoint25D]:
pt = inGeom.asMultiPoint()
attr1 = pt[0].x()
attr2 = pt[0].y()
else:
measure = QgsDistanceArea()

if method == 2:
measure.setSourceCrs(crs)
measure.setEllipsoid(ellips)
measure.setEllipsoidalMode(True)

attr1 = measure.measure(geom)
if geom.type() == QGis.Polygon:
attr2 = self.perimMeasure(geom, measure)
else:
attr2 = None

return (attr1, attr2)

def perimMeasure(self, geom, measure):
value = 0.0
if geom.isMultipart():
polygons = geom.asMultiPolygon()
for p in polygons:
for line in p:
value += measure.measureLine(line)
else:
poly = geom.asPolygon()
for r in poly:
value += measure.measureLine(r)

return value
29 changes: 29 additions & 0 deletions python/plugins/sextante/algs/ftools/FToolsUtils.py
View file
Expand Up @@ -96,6 +96,35 @@ def extractPoints( geom ):

return points

def simpleMeasure(geom, method=0, ellips=None, crs=None):
# method defines calculation type:
# 0 - layer CRS
# 1 - project CRS
# 2 - ellipsoidal
if geom.wkbType() in [QGis.WKBPoint, QGis.WKBPoint25D]:
pt = geom.asPoint()
attr1 = pt.x()
attr2 = pt.y()
elif geom.wkbType() in [QGis.WKBMultiPoint, QGis.WKBMultiPoint25D]:
pt = inGeom.asMultiPoint()
attr1 = pt[0].x()
attr2 = pt[0].y()
else:
measure = QgsDistanceArea()

if method == 2:
measure.setSourceCrs(crs)
measure.setEllipsoid(ellips)
measure.setEllipsoidalMode(True)

attr1 = measure.measure(geom)
if geom.type() == QGis.Polygon:
attr2 = measure.measurePerimeter(geom)
else:
attr2 = attr1

return (attr1, attr2)

def getUniqueValues(layer, fieldIndex):
values = []
layer.select([fieldIndex], QgsRectangle(), False)
Expand Down

0 comments on commit 48e9f42

Please sign in to comment.