7

I'm using QGIS 3.4.4. I have a polygon layer and a points layer. The points have the attributes district and streetname.

How can a check if all the points within a polygon have the same district and streetname?

4
  • Have you tried this, in the Select by expression field: 'streetnames' == or = 'district' ?
    – pnz1337
    Jun 3, 2020 at 13:12
  • Please explain, my question is have all the pointobjects in the same polygon the same streetname or district. So some combination of location and expressions.
    – Jan
    Jun 3, 2020 at 14:14
  • I wrote a python script that selects each polygon if points in it have different district or street name. You OK with that or do you prefer different solution? Jun 3, 2020 at 14:33
  • Thanks, that's OK for me, please profide.
    – Jan
    Jun 3, 2020 at 15:00

2 Answers 2

7

You can use Python script. Copy the script to QGIS Python Editor, change the layer names.

This script finds points in polygon, creates a set containing (district, streetname) pairs from these points and checks length of the set. If not 1, it selects that polygon. Because, in this case, at least one streetname or district is different than others. If all pairs are the same, then length of the set will be 1.

I recommend you first to create spatial index using "Vector general > Create spatial index" tool for both point and polygon layers.

from qgis.PyQt.QtWidgets import qApp 

point_layer = QgsProject.instance().mapLayersByName("AdresSectieNN")[0]
polygon_layer = QgsProject.instance().mapLayersByName("PCPlus")[0]

qgis.utils.iface.messageBar().clearWidgets()
messageBar = qgis.utils.iface.messageBar()
progress = QProgressBar()
progress.setMaximum(polygon_layer.featureCount()) 
messageBar.pushWidget(progress)

for i, pol in enumerate(polygon_layer.getFeatures()):
    # create a set
    x = {(p["Sectie"], p["straatnaam"]) for p in point_layer.getFeatures()
                if p.geometry().within(pol.geometry())}    
    # if set (x) has 1 element, it means all points have the same district name and street name
    if len(x) != 1:
        polygon_layer.select(pol.id())

    # set current progress
    progress.setValue(i+1)
    qApp.processEvents()

iface.mapCanvas().refresh()    

enter image description here

8
  • Thank you Kadir, I need some time to study this.
    – Jan
    Jun 3, 2020 at 18:11
  • I was not succesful: I paste the following code in the script editor point_layer = QgsProject.instance().mapLayersByName("AdresSectieNN")[0] polygon_layer = QgsProject.instance().mapLayersByName("PCPlus")[0] for pol in polygon_layer.getFeatures(): # create a set x = {(p["Sectie"], p["straatnaam"]) for p in point_layer.getFeatures() if p.geometry().within(pol.geometry())} # if set (x) has 1 element, it means all points have the same district name and street name if len(x) != 1: polygon_layer.select(pol.id()) Responding error
    – Jan
    Jun 4, 2020 at 10:16
  • Sorry for poor layout, don't know how to do it better in 5 minutes edit time I stopped the script after 10 minutes Not Responding
    – Jan
    Jun 4, 2020 at 10:24
  • How many points are there in point layer? If you have large point dataset, try to run this script on small dataset. Jun 4, 2020 at 11:47
  • What is qApp ?
    – Taras
    Jun 4, 2020 at 12:20
3

You can use a virtual layer to create a summary table.

Go to the menu Layer > Add Layer > Add/Edit Virtual Layer... and enter the following query. Replace the layer name for yours

SELECT poly.id, pt.district, pt.streetname, count(*)
FROM my_polygon_layer poly
    JOIN my_point_layer pt ON st_intersects(poly.geometry, pt.geometry)
GROUP BY poly.id, pt.district, pt.streetname

If you only want to identify polygons having several districts/points, you can refine the query:

SELECT id, count(*) as cnt_diff
FROM (
    SELECT poly.id, pt.district, pt.streetname, count(*) as cnt_same
    FROM my_polygon_layer poly
        JOIN my_point_layer pt ON st_intersects(poly.geometry, pt.geometry)
    GROUP BY poly.id, pt.district, pt.streetname
)
GROUP BY id
HAVING COUNT(*) > 1
4
  • Thank you JGH, I need some time to study this.
    – Jan
    Jun 3, 2020 at 18:11
  • I was not succesful with your SQL code. I corrected the GROUP BY statement by removing ",count(*)" because I got an error. After correcting and press Test, QGis did not answer (Not Responding) after 10 minutes.
    – Jan
    Jun 4, 2020 at 9:47
  • @Jan good catch, I have edited the code accordingly. It should work though. If you are running it on a large dateset, or on very large polygons (country size), it might be slow (and since you have the same issue with Kadir's solution it might very well be the case). In this case you will want to use a database with indexed field/geometries and eventually on subdivided polygons
    – JGH
    Jun 4, 2020 at 11:48
  • The point layer has 6771 features and the polygon layer has 7403 features. Not a country size at all.
    – Jan
    Jun 4, 2020 at 13:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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