I was writing a plugin for QGIS when I wanted to copy all my data to memory layers and edit them instead of actual shapefiles to speed up the process. In a QgsMapTool-derived class, I am using a method call like this

def canvasPressEvent(self, e):
    mapCoords = self.toMapCoordinates(e.pos())
    features = self.__layer.getFeatures(QgsFeatureRequest(QgsRectangle(mapCoords, mapCoords)))

This worked just fine for the old, file-based layers, however, for the new in-memory ones, I'm getting horribly inaccurate results.

My code for copying the layer is pretty straightfoward:

with fiona.open(path) as f:
    geometry = f.schema["geometry"] # yeah, this is a bit hacky, sry

originalLayer = QgsVectorLayer(path, name + "_temp_layer", "ogr")

# build the creation string
props = geometry + "?"
props += "crs=WKT:" + originalLayer.crs().toWkt()
props += "&"
props += "index=yes"

# create layer and start filling it
layer = QgsVectorLayer(props, name, "memory")
#layer.setCrs(originalLayer.crs()) #this doesn't seem to change anything

# copy fields
for field in originalLayer.fields():
    layer.addAttribute(QgsField(field.name(), field.type()))

# copy features
newFeatures = []
for feature in originalLayer.getFeatures():
    newFeature = QgsFeature()
result = layer.dataProvider().addFeatures(newFeatures)

# commit

Now, this code seems to work just fine, the result looks just as it looked before, but when using the tool, I can click on a spot a few centimetres away from a feature and the tool detects it as if I clicked the feature itself. This also creates massive overlap between the features. These things worked just fine with file-based layers.

It feels like something is wrong with the CRS... Do I somehow have to apply the CRS to be taken into account during spatial querys?

I already tried manually creating spatial indices, but that yielded the same results.

  • Does it change something using "crs=" instead of "crs=WKT:"? – mgri Jan 22 '17 at 10:51

As no one could answer me, I opened this bug report: https://issues.qgis.org/issues/16068

It looks like you need to set the ExactIntersect flag in the QgsFeatureRequest. I still hope that this is getting fixed (or at least documented).

The proposed fix does not work, but you can refine your results by checking the geometry intersection directly:

features = [feature for feature in features 
    if feature.geometry().intersects(QgsRectangle(mapCoords, mapCoords))]

Note that feature.geometry().contains(mapCoords) method does not work. I'll probably open another bug report for this.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.