Skip to content

Commit e899240

Browse files
author
cfarmer
committedMar 10, 2010
Uses spatial index to select intersecting features (used selections before). Tools should both run faster and with fewer map canvas glitches.
git-svn-id: http://svn.osgeo.org/qgis/trunk/qgis@13035 c8812cc2-4d05-0410-92ff-de0c093fc19c

File tree

2 files changed

+231
-272
lines changed

2 files changed

+231
-272
lines changed
 

‎python/plugins/fTools/tools/doIntersectLines.py

Lines changed: 125 additions & 149 deletions
Original file line numberDiff line numberDiff line change
@@ -39,157 +39,133 @@
3939

4040
class Dialog(QDialog, Ui_Dialog):
4141

42-
def __init__(self, iface):
43-
QDialog.__init__(self)
44-
self.iface = iface
45-
# Set up the user interface from Designer.
46-
self.setupUi(self)
47-
QObject.connect(self.toolOut, SIGNAL("clicked()"), self.outFile)
48-
QObject.connect(self.inLine1, SIGNAL("currentIndexChanged(QString)"), self.update1)
49-
QObject.connect(self.inLine2, SIGNAL("currentIndexChanged(QString)"), self.update2)
50-
self.setWindowTitle( self.tr("Line intersections") )
51-
# populate layer list
52-
self.progressBar.setValue(0)
53-
mapCanvas = self.iface.mapCanvas()
54-
for i in range(mapCanvas.layerCount()):
55-
layer = mapCanvas.layer(i)
56-
if layer.type() == layer.VectorLayer:
57-
if layer.geometryType() == QGis.Line:
58-
self.inLine1.addItem(layer.name())
59-
self.inLine2.addItem(layer.name())
42+
def __init__(self, iface):
43+
QDialog.__init__(self)
44+
self.iface = iface
45+
# Set up the user interface from Designer.
46+
self.setupUi(self)
47+
QObject.connect(self.toolOut, SIGNAL("clicked()"), self.outFile)
48+
QObject.connect(self.inLine1, SIGNAL("currentIndexChanged(QString)"), self.update1)
49+
QObject.connect(self.inLine2, SIGNAL("currentIndexChanged(QString)"), self.update2)
50+
self.setWindowTitle( self.tr("Line intersections") )
51+
# populate layer list
52+
self.progressBar.setValue(0)
53+
mapCanvas = self.iface.mapCanvas()
54+
layers = ftools_utils.getLayerNames([QGis.Line])
55+
self.inLine1.addItems(layers)
56+
self.inLine2.addItems(layers)
6057

61-
def update1(self, inputLayer):
62-
self.inField1.clear()
63-
changedLayer = self.getVectorLayerByName(unicode(inputLayer))
64-
changedField = self.getFieldList(changedLayer)
65-
for i in changedField:
66-
self.inField1.addItem(unicode(changedField[i].name()))
58+
def update1(self, inputLayer):
59+
self.inField1.clear()
60+
changedLayer = ftools_utils.getVectorLayerByName(unicode(inputLayer))
61+
changedField = ftools_utils.getFieldList(changedLayer)
62+
for i in changedField:
63+
self.inField1.addItem(unicode(changedField[i].name()))
6764

68-
def update2(self, inputLayer):
69-
self.inField2.clear()
70-
changedLayer = self.getVectorLayerByName(unicode(inputLayer))
71-
changedField = self.getFieldList(changedLayer)
72-
for i in changedField:
73-
self.inField2.addItem(unicode(changedField[i].name()))
65+
def update2(self, inputLayer):
66+
self.inField2.clear()
67+
changedLayer = ftools_utils.getVectorLayerByName(unicode(inputLayer))
68+
changedField = ftools_utils.getFieldList(changedLayer)
69+
for i in changedField:
70+
self.inField2.addItem(unicode(changedField[i].name()))
7471

75-
def accept(self):
76-
if self.inLine1.currentText() == "":
77-
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify input line layer") )
78-
elif self.outShape.text() == "":
79-
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify output shapefile") )
80-
elif self.inLine2.currentText() == "":
81-
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify line intersect layer") )
82-
elif self.inField1.currentText() == "":
83-
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify input unique ID field") )
84-
elif self.inField2.currentText() == "":
85-
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify intersect unique ID field") )
86-
else:
87-
line1 = self.inLine1.currentText()
88-
line2 = self.inLine2.currentText()
89-
field1 = self.inField1.currentText()
90-
field2 = self.inField2.currentText()
91-
outPath = self.outShape.text()
92-
if outPath.contains("\\"):
93-
outName = outPath.right((outPath.length() - outPath.lastIndexOf("\\")) - 1)
94-
else:
95-
outName = outPath.right((outPath.length() - outPath.lastIndexOf("/")) - 1)
96-
if outName.endsWith(".shp"):
97-
outName = outName.left(outName.length() - 4)
98-
self.outShape.clear()
99-
self.compute(line1, line2, field1, field2, outPath, self.progressBar)
100-
self.progressBar.setValue(100)
101-
addToTOC = QMessageBox.question(self, self.tr("Generate Centroids"), self.tr("Created output point shapefile:\n%1\n\nWould you like to add the new layer to the TOC?").arg( outPath ), QMessageBox.Yes, QMessageBox.No, QMessageBox.NoButton)
102-
if addToTOC == QMessageBox.Yes:
103-
self.vlayer = QgsVectorLayer(outPath, unicode(outName), "ogr")
104-
QgsMapLayerRegistry.instance().addMapLayer(self.vlayer)
105-
self.progressBar.setValue(0)
72+
def accept(self):
73+
if self.inLine1.currentText() == "":
74+
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify input line layer") )
75+
elif self.outShape.text() == "":
76+
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify output shapefile") )
77+
elif self.inLine2.currentText() == "":
78+
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify line intersect layer") )
79+
elif self.inField1.currentText() == "":
80+
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify input unique ID field") )
81+
elif self.inField2.currentText() == "":
82+
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Please specify intersect unique ID field") )
83+
else:
84+
line1 = self.inLine1.currentText()
85+
line2 = self.inLine2.currentText()
86+
field1 = self.inField1.currentText()
87+
field2 = self.inField2.currentText()
88+
outPath = self.outShape.text()
89+
self.outShape.clear()
90+
self.compute(line1, line2, field1, field2, outPath, self.progressBar)
91+
self.progressBar.setValue(100)
92+
addToTOC = QMessageBox.question(self, self.tr("Generate Centroids"), self.tr("Created output point shapefile:\n%1\n\nWould you like to add the new layer to the TOC?").arg( outPath ), QMessageBox.Yes, QMessageBox.No, QMessageBox.NoButton)
93+
if addToTOC == QMessageBox.Yes:
94+
if not ftools_utils.addShapeToCanvas( unicode( outPath ) ):
95+
QMessageBox.warning( self, self.tr("Geoprocessing"), self.tr( "Error loading output shapefile:\n%1" ).arg( unicode( outPath ) ))
96+
self.progressBar.setValue(0)
10697

107-
def outFile(self):
108-
self.outShape.clear()
109-
( self.shapefileName, self.encoding ) = ftools_utils.saveDialog( self )
110-
if self.shapefileName is None or self.encoding is None:
111-
return
112-
self.outShape.setText( QString( self.shapefileName ) )
98+
def outFile(self):
99+
self.outShape.clear()
100+
( self.shapefileName, self.encoding ) = ftools_utils.saveDialog( self )
101+
if self.shapefileName is None or self.encoding is None:
102+
return
103+
self.outShape.setText( QString( self.shapefileName ) )
113104

114-
def compute(self, line1, line2, field1, field2, outPath, progressBar):
115-
layer1 = self.getVectorLayerByName(line1)
116-
layer2 = self.getVectorLayerByName(line2)
117-
provider1 = layer1.dataProvider()
118-
provider2 = layer2.dataProvider()
119-
allAttrs = provider1.attributeIndexes()
120-
provider1.select(allAttrs)
121-
allAttrs = provider2.attributeIndexes()
122-
provider2.select(allAttrs)
123-
fieldList = self.getFieldList(layer1)
124-
index1 = provider1.fieldNameIndex(field1)
125-
field1 = fieldList[index1]
126-
field1.setName(unicode(field1.name()) + "_1")
127-
fieldList = self.getFieldList(layer2)
128-
index2 = provider2.fieldNameIndex(field2)
129-
field2 = fieldList[index2]
130-
field2.setName(unicode(field2.name()) + "_2")
131-
fieldList = {0:field1, 1:field2}
132-
sRs = provider1.crs()
133-
check = QFile(self.shapefileName)
134-
if check.exists():
135-
if not QgsVectorFileWriter.deleteShapeFile(self.shapefileName):
136-
return
137-
writer = QgsVectorFileWriter(self.shapefileName, self.encoding, fieldList, QGis.WKBPoint, sRs)
138-
#writer = QgsVectorFileWriter(outPath, "UTF-8", fieldList, QGis.WKBPoint, sRs)
139-
inFeat = QgsFeature()
140-
outFeat = QgsFeature()
141-
inGeom = QgsGeometry()
142-
tempGeom = QgsGeometry()
143-
start = 15.00
144-
add = 85.00 / layer1.featureCount()
145-
while provider1.nextFeature(inFeat):
146-
inGeom = inFeat.geometry()
147-
lineList = []
148-
#(check, lineList) = layer2.featuresInRectangle(inGeom.boundingBox(), True, True)
149-
# Below is a work-around for featuresInRectangle
150-
# Which appears to have been removed for QGIS version 1.0
151-
layer2.select(inGeom.boundingBox(), False)
152-
lineList = layer2.selectedFeatures()
153-
if len(lineList) > 0: check = 0
154-
else: check = 1
155-
if check == 0:
156-
for i in lineList:
157-
tempGeom = i.geometry()
158-
tempList = []
159-
atMap1 = inFeat.attributeMap()
160-
atMap2 = i.attributeMap()
161-
if inGeom.intersects(tempGeom):
162-
tempGeom = inGeom.intersection(tempGeom)
163-
if tempGeom.type() == QGis.Point:
164-
if tempGeom.isMultipart():
165-
tempList = tempGeom.asMultiPoint()
166-
else:
167-
tempList.append(tempGeom.asPoint())
168-
for j in tempList:
169-
outFeat.setGeometry(tempGeom.fromPoint(j))
170-
outFeat.addAttribute(0, atMap1[index1])
171-
outFeat.addAttribute(1, atMap2[index2])
172-
writer.addFeature(outFeat)
173-
start = start + add
174-
progressBar.setValue(start)
175-
del writer
176-
177-
def getVectorLayerByName(self, myName):
178-
mc = self.iface.mapCanvas()
179-
nLayers = mc.layerCount()
180-
for l in range(nLayers):
181-
layer = mc.layer(l)
182-
if layer.name() == unicode(myName):
183-
vlayer = QgsVectorLayer(unicode(layer.source()), unicode(myName), unicode(layer.dataProvider().name()))
184-
if vlayer.isValid():
185-
return vlayer
186-
else:
187-
QMessageBox.information(self, self.tr("Locate Line Intersections"), self.tr("Vector layer is not valid") )
188-
189-
def getFieldList(self, vlayer):
190-
fProvider = vlayer.dataProvider()
191-
feat = QgsFeature()
192-
allAttrs = fProvider.attributeIndexes()
193-
fProvider.select(allAttrs)
194-
myFields = fProvider.fields()
195-
return myFields
105+
def compute(self, line1, line2, field1, field2, outPath, progressBar):
106+
layer1 = ftools_utils.getVectorLayerByName(line1)
107+
layer2 = ftools_utils.getVectorLayerByName(line2)
108+
provider1 = layer1.dataProvider()
109+
provider2 = layer2.dataProvider()
110+
allAttrs = provider1.attributeIndexes()
111+
provider1.select(allAttrs)
112+
allAttrs = provider2.attributeIndexes()
113+
provider2.select(allAttrs)
114+
fieldList = ftools_utils.getFieldList(layer1)
115+
index1 = provider1.fieldNameIndex(field1)
116+
field1 = fieldList[index1]
117+
field1.setName(unicode(field1.name()) + "_1")
118+
fieldList = ftools_utils.getFieldList(layer2)
119+
index2 = provider2.fieldNameIndex(field2)
120+
field2 = fieldList[index2]
121+
field2.setName(unicode(field2.name()) + "_2")
122+
fieldList = {0:field1, 1:field2}
123+
sRs = provider1.crs()
124+
check = QFile(self.shapefileName)
125+
if check.exists():
126+
if not QgsVectorFileWriter.deleteShapeFile(self.shapefileName):
127+
return
128+
writer = QgsVectorFileWriter(self.shapefileName, self.encoding, fieldList, QGis.WKBPoint, sRs)
129+
#writer = QgsVectorFileWriter(outPath, "UTF-8", fieldList, QGis.WKBPoint, sRs)
130+
inFeat = QgsFeature()
131+
inFeatB = QgsFeature()
132+
outFeat = QgsFeature()
133+
inGeom = QgsGeometry()
134+
tempGeom = QgsGeometry()
135+
start = 15.00
136+
add = 85.00 / layer1.featureCount()
137+
index = ftools_utils.createIndex( provider2 )
138+
while provider1.nextFeature(inFeat):
139+
inGeom = inFeat.geometry()
140+
lineList = []
141+
#(check, lineList) = layer2.featuresInRectangle(inGeom.boundingBox(), True, True)
142+
# Below is a work-around for featuresInRectangle
143+
# Which appears to have been removed for QGIS version 1.0
144+
#layer2.select(inGeom.boundingBox(), False)
145+
#lineList = layer2.selectedFeatures()
146+
lineList = index.intersects( inGeom.boundingBox() )
147+
if len(lineList) > 0: check = 0
148+
else: check = 1
149+
if check == 0:
150+
for i in lineList:
151+
provider2.featureAtId( int( i ), inFeatB , True, allAttrs )
152+
tmpGeom = QgsGeometry( inFeatB.geometry() )
153+
#tempGeom = i.geometry()
154+
tempList = []
155+
atMap1 = inFeat.attributeMap()
156+
atMap2 = inFeatB.attributeMap()
157+
if inGeom.intersects(tmpGeom):
158+
tempGeom = inGeom.intersection(tmpGeom)
159+
if tempGeom.type() == QGis.Point:
160+
if tempGeom.isMultipart():
161+
tempList = tempGeom.asMultiPoint()
162+
else:
163+
tempList.append(tempGeom.asPoint())
164+
for j in tempList:
165+
outFeat.setGeometry(tempGeom.fromPoint(j))
166+
outFeat.addAttribute(0, atMap1[index1])
167+
outFeat.addAttribute(1, atMap2[index2])
168+
writer.addFeature(outFeat)
169+
start = start + add
170+
progressBar.setValue(start)
171+
del writer
Lines changed: 106 additions & 123 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,5 @@
1-
#-----------------------------------------------------------
1+
# -*- coding: utf-8 -*-
2+
#-----------------------------------------------------------
23
#
34
# Points in Polygon
45
#
@@ -38,127 +39,109 @@
3839

3940
class Dialog(QDialog, Ui_Dialog):
4041

41-
def __init__(self, iface):
42-
QDialog.__init__(self)
43-
self.iface = iface
44-
# Set up the user interface from Designer.
45-
self.setupUi(self)
46-
QObject.connect(self.toolOut, SIGNAL("clicked()"), self.outFile)
47-
self.setWindowTitle(self.tr("Count Points in Polygon"))
48-
# populate layer list
49-
self.progressBar.setValue(0)
50-
mapCanvas = self.iface.mapCanvas()
51-
for i in range(mapCanvas.layerCount()):
52-
layer = mapCanvas.layer(i)
53-
if layer.type() == layer.VectorLayer:
54-
if layer.geometryType() == QGis.Polygon:
55-
self.inPolygon.addItem(layer.name())
56-
elif layer.geometryType() == QGis.Point:
57-
self.inPoint.addItem(layer.name())
58-
59-
def accept(self):
60-
if self.inPolygon.currentText() == "":
61-
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify input polygon vector layer"))
62-
elif self.outShape.text() == "":
63-
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify output shapefile"))
64-
elif self.inPoint.currentText() == "":
65-
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify input point vector layer"))
66-
elif self.lnField.text() == "":
67-
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify output count field"))
68-
else:
69-
inPoly = self.inPolygon.currentText()
70-
inPts = self.inPoint.currentText()
71-
inField = self.lnField.text()
72-
outPath = self.outShape.text()
73-
if outPath.contains("\\"):
74-
outName = outPath.right((outPath.length() - outPath.lastIndexOf("\\")) - 1)
75-
else:
76-
outName = outPath.right((outPath.length() - outPath.lastIndexOf("/")) - 1)
77-
if outName.endsWith(".shp"):
78-
outName = outName.left(outName.length() - 4)
79-
self.compute(inPoly, inPts, inField, outPath, self.progressBar)
80-
self.outShape.clear()
81-
addToTOC = QMessageBox.question(self, self.tr("Count Points in Polygon"), self.tr("Created output shapefile:\n%1\n\nWould you like to add the new layer to the TOC?").arg(unicode(outPath)), QMessageBox.Yes, QMessageBox.No, QMessageBox.NoButton)
82-
if addToTOC == QMessageBox.Yes:
83-
self.vlayer = QgsVectorLayer(outPath, unicode(outName), "ogr")
84-
QgsMapLayerRegistry.instance().addMapLayer(self.vlayer)
85-
self.progressBar.setValue(0)
86-
87-
def outFile(self):
88-
self.outShape.clear()
89-
( self.shapefileName, self.encoding ) = ftools_utils.saveDialog( self )
90-
if self.shapefileName is None or self.encoding is None:
91-
return
92-
self.outShape.setText( QString( self.shapefileName ) )
42+
def __init__(self, iface):
43+
QDialog.__init__(self)
44+
self.iface = iface
45+
# Set up the user interface from Designer.
46+
self.setupUi(self)
47+
QObject.connect(self.toolOut, SIGNAL("clicked()"), self.outFile)
48+
self.setWindowTitle(self.tr("Count Points in Polygon"))
49+
# populate layer list
50+
self.progressBar.setValue(0)
51+
mapCanvas = self.iface.mapCanvas()
52+
layers = ftools_utils.getLayerNames([QGis.Polygon])
53+
self.inPolygon.addItems(layers)
54+
layers = ftools_utils.getLayerNames([QGis.Point])
55+
self.inPoint.addItems(layers)
56+
57+
def accept(self):
58+
if self.inPolygon.currentText() == "":
59+
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify input polygon vector layer"))
60+
elif self.outShape.text() == "":
61+
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify output shapefile"))
62+
elif self.inPoint.currentText() == "":
63+
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify input point vector layer"))
64+
elif self.lnField.text() == "":
65+
QMessageBox.information(self, self.tr("Count Points In Polygon"), self.tr("Please specify output count field"))
66+
else:
67+
inPoly = self.inPolygon.currentText()
68+
inPts = self.inPoint.currentText()
69+
inField = self.lnField.text()
70+
outPath = self.outShape.text()
71+
if outPath.contains("\\"):
72+
outName = outPath.right((outPath.length() - outPath.lastIndexOf("\\")) - 1)
73+
else:
74+
outName = outPath.right((outPath.length() - outPath.lastIndexOf("/")) - 1)
75+
if outName.endsWith(".shp"):
76+
outName = outName.left(outName.length() - 4)
77+
self.compute(inPoly, inPts, inField, outPath, self.progressBar)
78+
self.outShape.clear()
79+
addToTOC = QMessageBox.question(self, self.tr("Count Points in Polygon"), self.tr("Created output shapefile:\n%1\n\nWould you like to add the new layer to the TOC?").arg(unicode(outPath)), QMessageBox.Yes, QMessageBox.No, QMessageBox.NoButton)
80+
if addToTOC == QMessageBox.Yes:
81+
self.vlayer = QgsVectorLayer(outPath, unicode(outName), "ogr")
82+
QgsMapLayerRegistry.instance().addMapLayer(self.vlayer)
83+
self.progressBar.setValue(0)
9384

94-
def compute(self, inPoly, inPts, inField, outPath, progressBar):
95-
polyLayer = self.getVectorLayerByName(inPoly)
96-
pointLayer = self.getVectorLayerByName(inPts)
97-
polyProvider = polyLayer.dataProvider()
98-
pointProvider = pointLayer.dataProvider()
99-
if polyProvider.crs() <> pointProvider.crs():
100-
QMessageBox.warning(self, self.tr("CRS warning!"), self.tr("Warning: Input layers have non-matching CRS.\nThis may cause unexpected results."))
101-
allAttrs = polyProvider.attributeIndexes()
102-
polyProvider.select(allAttrs)
103-
allAttrs = pointProvider.attributeIndexes()
104-
pointProvider.select(allAttrs)
105-
fieldList = self.getFieldList(polyLayer)
106-
index = polyProvider.fieldNameIndex(unicode(inField))
107-
if index == -1:
108-
index = polyProvider.fieldCount()
109-
field = QgsField(unicode(inField), QVariant.Int, "real", 24, 15, self.tr("point count field"))
110-
fieldList[index] = field
111-
sRs = polyProvider.crs()
112-
check = QFile(self.shapefileName)
113-
if check.exists():
114-
if not QgsVectorFileWriter.deleteShapeFile(self.shapefileName):
115-
return
116-
writer = QgsVectorFileWriter(self.shapefileName, self.encoding, fieldList, polyProvider.geometryType(), sRs)
117-
#writer = QgsVectorFileWriter(outPath, "UTF-8", fieldList, polyProvider.geometryType(), sRs)
118-
inFeat = QgsFeature()
119-
outFeat = QgsFeature()
120-
inGeom = QgsGeometry()
121-
start = 15.00
122-
add = 85.00 / polyProvider.featureCount()
123-
while polyProvider.nextFeature(inFeat):
124-
inGeom = inFeat.geometry()
125-
atMap = inFeat.attributeMap()
126-
outFeat.setAttributeMap(atMap)
127-
outFeat.setGeometry(inGeom)
128-
pointList = []
129-
count = 0
130-
#(check, pointList) = pointLayer.featuresInRectangle(inGeom.boundingBox(), True, True)
131-
pointLayer.select(inGeom.boundingBox(), False)
132-
pointList = pointLayer.selectedFeatures()
133-
if len(pointList) > 0: check = 0
134-
else: check = 1
135-
if check == 0:
136-
for i in pointList:
137-
if inGeom.contains(i.geometry().asPoint()):
138-
count = count + 1
139-
outFeat.setAttributeMap(atMap)
140-
outFeat.addAttribute(index, QVariant(count))
141-
writer.addFeature(outFeat)
142-
start = start + 1
143-
progressBar.setValue(start)
144-
del writer
145-
146-
def getVectorLayerByName(self, myName):
147-
mc = self.iface.mapCanvas()
148-
nLayers = mc.layerCount()
149-
for l in range(nLayers):
150-
layer = mc.layer(l)
151-
if layer.name() == unicode(myName):
152-
vlayer = QgsVectorLayer(unicode(layer.source()), unicode(myName), unicode(layer.dataProvider().name()))
153-
if vlayer.isValid():
154-
return vlayer
155-
else:
156-
QMessageBox.information(self, self.tr("Counts Points In Polygon"), self.tr("Vector layer is not valid"))
85+
def outFile(self):
86+
self.outShape.clear()
87+
( self.shapefileName, self.encoding ) = ftools_utils.saveDialog( self )
88+
if self.shapefileName is None or self.encoding is None:
89+
return
90+
self.outShape.setText( QString( self.shapefileName ) )
15791

158-
def getFieldList(self, vlayer):
159-
fProvider = vlayer.dataProvider()
160-
feat = QgsFeature()
161-
allAttrs = fProvider.attributeIndexes()
162-
fProvider.select(allAttrs)
163-
myFields = fProvider.fields()
164-
return myFields
92+
def compute(self, inPoly, inPts, inField, outPath, progressBar):
93+
polyLayer = ftools_utils.getVectorLayerByName(inPoly)
94+
pointLayer = ftools_utils.getVectorLayerByName(inPts)
95+
polyProvider = polyLayer.dataProvider()
96+
pointProvider = pointLayer.dataProvider()
97+
if polyProvider.crs() <> pointProvider.crs():
98+
QMessageBox.warning(self, self.tr("CRS warning!"), self.tr("Warning: Input layers have non-matching CRS.\nThis may cause unexpected results."))
99+
allAttrs = polyProvider.attributeIndexes()
100+
polyProvider.select(allAttrs)
101+
allAttrs = pointProvider.attributeIndexes()
102+
pointProvider.select(allAttrs)
103+
fieldList = ftools_utils.getFieldList(polyLayer)
104+
index = polyProvider.fieldNameIndex(unicode(inField))
105+
if index == -1:
106+
index = polyProvider.fieldCount()
107+
field = QgsField(unicode(inField), QVariant.Int, "real", 24, 15, self.tr("point count field"))
108+
fieldList[index] = field
109+
sRs = polyProvider.crs()
110+
check = QFile(self.shapefileName)
111+
if check.exists():
112+
if not QgsVectorFileWriter.deleteShapeFile(self.shapefileName):
113+
return
114+
writer = QgsVectorFileWriter(self.shapefileName, self.encoding, fieldList, polyProvider.geometryType(), sRs)
115+
#writer = QgsVectorFileWriter(outPath, "UTF-8", fieldList, polyProvider.geometryType(), sRs)
116+
inFeat = QgsFeature()
117+
inFeatB = QgsFeature()
118+
outFeat = QgsFeature()
119+
inGeom = QgsGeometry()
120+
start = 15.00
121+
add = 85.00 / polyProvider.featureCount()
122+
spatialIndex = ftools_utils.createIndex( pointProvider )
123+
while polyProvider.nextFeature(inFeat):
124+
inGeom = inFeat.geometry()
125+
atMap = inFeat.attributeMap()
126+
outFeat.setAttributeMap(atMap)
127+
outFeat.setGeometry(inGeom)
128+
pointList = []
129+
count = 0
130+
#(check, pointList) = pointLayer.featuresInRectangle(inGeom.boundingBox(), True, True)
131+
#pointLayer.select(inGeom.boundingBox(), False)
132+
#pointList = pointLayer.selectedFeatures()
133+
pointList = spatialIndex.intersects(inGeom.boundingBox())
134+
if len(pointList) > 0: check = 0
135+
else: check = 1
136+
if check == 0:
137+
for i in pointList:
138+
pointProvider.featureAtId( int( i ), inFeatB , True, allAttrs )
139+
tmpGeom = QgsGeometry( inFeatB.geometry() )
140+
if inGeom.contains(tmpGeom.asPoint()):
141+
count = count + 1
142+
outFeat.setAttributeMap(atMap)
143+
outFeat.addAttribute(index, QVariant(count))
144+
writer.addFeature(outFeat)
145+
start = start + 1
146+
progressBar.setValue(start)
147+
del writer

0 commit comments

Comments
 (0)
Please sign in to comment.