Hi,
Just by way of finality:
The CRS of the polygon SHP was not what it was purported to be.
I found that the .prj file was invalid when opening it in QGIS - maybe
next time I will first look at the maps I am given, rather than just
dive into python scripting.
All working as expected, now.
Thanks to those who gave ideas to try.
Kind regards,
Zoltan
On 2024/08/25 19:44, Zoltan Szecsei via QGIS-User wrote:
So as updated plan-of-attack as below.....
lyr_poly = QgsVectorLayer(poly_shp, '', 'Ogr') lyr_poly_si =
QgsSpatialIndex(lyr_poly.getFeatures(),
flags=QgsSpatialIndex.FlagStoreFeatureGeometries) # EPSG:4326 point =
QgsPointXY(lon, lat)
point_rect = QgsGeometry.fromPointXY(point).buffer(0.001,5).boundingBox()
# feature_ids = lyr_poly_si.intersects(point_rect)# intersects wants
rectangle requests = QgsFeatureRequest().setFilterFids(feature_ids)
# feature_idx = lyr_poly_si.nearestNeighbor(point,10,maxDistance=0.001)
requestx = QgsFeatureRequest().setFilterFids(feature_idx)
For the above code the lists feature_ids and feature_idx are always empty.
The Longitude and Latitude coords are correct, also for the box
created from the Point.
I'd be grateful for any insight on why nothing is found.
Thanks in advance,
Zoltan
On 2024/08/25 16:19, Zoltan Szecsei via QGIS-User wrote:
Hi,
I'm doing my "once a year need to do something in python with QGIS" :-/
Using QGIS 3.38 and PyCharm.
I have a SHP file with 28000 polygons.
*A SPREADSHEET file with 22000 Points (for which I am making a SHP
file and wanting attribute vaues from underlying polygons)*
Both in EPSG:4326
Without using QGIS processing, what is the quickest way to, for each
point 1 by 1, find the underlying polygon and read the attribute
field values of that polygon?
The code snippet below is far too inefficient - but I cannot find
code that selects at XY and benefits from any spatial index on the
polygons.
lyr_poly = QgsVectorLayer(polygon_shp, '', 'Ogr') point =
QgsPointXY(lon, lat) request =
QgsFeatureRequest().setFilterRect(lyr_poly.extent()) for
feature in lyr_poly.getFeatures(request): geom =
feature.geometry() if
geom.contains(QgsGeometry.fromPointXY(point)):
lyr_poly.select(feature.id()) print(f"Feature with
ID {feature.id()} is selected.") break
Thanks in advance,
Zoltan
_______________________________________________
QGIS-User mailing list
[email protected]
List info:https://lists.osgeo.org/mailman/listinfo/qgis-user
Unsubscribe:https://lists.osgeo.org/mailman/listinfo/qgis-user
_______________________________________________
QGIS-User mailing list
[email protected]
List info: https://lists.osgeo.org/mailman/listinfo/qgis-user
Unsubscribe: https://lists.osgeo.org/mailman/listinfo/qgis-user
