3

I have a vector tile layer in QGIS, loaded from: https://cdn.kadastr.live/tiles/maps/kadastr/land_polygons/{z}/{x}/{y}.pbf.

How to select a single feature of this layer by coordinates, in PyQGIS?

I am trying to use selectByGeometry method, but it is not working.

latitude = 48.23117391
longitude = 22.74110528

layer = iface.activeLayer()


iface.mapCanvas().refresh()
geom = QgsGeometry.fromPointXY(QgsPointXY(longitude, latitude))


canvas = iface.mapCanvas()
context = QgsSelectionContext()
context.setScale(10000)


print(layer.selectedFeatureCount())
layer.selectByGeometry(geom,context)#,context,behavior= Qgis.SelectBehavior.AddToSelection,relationship= Qgis.SelectGeometryRelationship.Intersect,flags=Qgis.SelectionFlags(2))
print(layer.selectedFeatureCount())
Taras
  • 32,823
  • 4
  • 66
  • 137
Стьопа Брич
  • 33
  • 3

2 Answers2

4

Your vector tile layer is in EPSG:3857 but your lat and long are in geographic degrees (EPSG:4326). You need to transform the selection geometry to the layer crs.

Additionally, I found it necessary to use a polygon selection geometry. I created a QgsGeometry from the bounding box of a very small buffer applied to the transformed point. To be honest I'm not sure why this is necessary - looking at the source code it seems that it should work with a single point selection geometry, though the tests for QgsVectorTileLayer use a polygon selection geometry.

This worked for me:

from qgis.utils import iface
from qgis.core import (QgsGeometry, QgsPointXY, QgsCoordinateReferenceSystem,
                       QgsCoordinateTransform, QgsProject, QgsSelectionContext)

latitude = 48.23117391
longitude = 22.74110528


layer = iface.activeLayer()


geom = QgsGeometry.fromPointXY(QgsPointXY(longitude, latitude))


def transformed_geom(g):
    pt_geom = QgsGeometry().fromWkt(g.asWkt())
    xform = QgsCoordinateTransform(QgsCoordinateReferenceSystem('epsg:4326'), layer.crs(), QgsProject.instance())
    pt_geom.transform(xform)
    bb_geom = QgsGeometry().fromRect(pt_geom.buffer(0.001, 1).boundingBox())
    return bb_geom


canvas = iface.mapCanvas()
context = QgsSelectionContext()
context.setScale(canvas.scale())


print(layer.selectedFeatureCount())


layer.selectByGeometry(transformed_geom(geom), context)


print(layer.selectedFeatureCount())

enter image description here

And at higher zoom level:

enter image description here

Taras
  • 32,823
  • 4
  • 66
  • 137
Ben W
  • 21,426
  • 3
  • 15
  • 39
  • Thank you! This works fine for me. Could you please explain why do you use: pt_geom = QgsGeometry().fromWkt(g.asWkt()) ? – Стьопа Брич Sep 25 '23 at 17:53
  • @СтьопаБрич, it's just copying the input geometry so that the original geom is not modified. – Ben W Sep 26 '23 at 10:35
  • 2
    I still do not know why the creation of a tiny buffer around a point is necessary ... it precisely says bool isPointOrRectangle – Taras Nov 22 '23 at 09:52
  • QgsGeometry(g) also creates a copy of a geometry and will be more efficient because it saves on a roundtrip through WKT. – bugmenot123 Nov 22 '23 at 11:29
-1

This is all you need to select all the features of a layer at a certain point:

latitude = 48.23117391
longitude = 22.74110528
layer = iface.activeLayer()

point = QgsPointXY(longitude, latitude)
rect = QgsRectangle(point, point)


layer.selectByRect(rect)
Vsevolod Tsukanov
  • 1,047
  • 8
  • 13
  • 1
    This solution is not applicable for the QgsVectorTileLayer, otherwise you will get this: "AttributeError: 'QgsVectorTileLayer' object has no attribute 'selectByRect'" – Taras Nov 22 '23 at 08:33