TransWikia.com

Determining how many of my coordinates overlap with layer in QGIS

Geographic Information Systems Asked by Fleur on October 1, 2021

I am working on a project about a certain animal species. I have gotten coordinates from where the species has been observed.

Now I want to know what percentage of the observations are from gardens/residential area and what percentage from forest/ edges of forests.

I have found a layer for living areas so the animals found there would count as ‘found in garden/residential area’.

So my question is how do I find out how many of my coordinates overlap with this layer and how many do not?

2 Answers

Let's assume that there is a csv-file (in UTF-8) 'animals.csv', see image below

input_csv

Step 1. In QGIS create a new QGIS Project. Then proceed with Layer > Add Layer > Add Delimited Text Layer ....

step_1

Step 2. Add your area via Layer > Add Layer > Add Vector Layer .... After adding areas your QGIS working window may look as following.

step_2

Step 3. Apply 'Vector > Data Management Tools > Join attributes by location' as I demonstrated on the image below.

step_3

Step 4. In the attribute table of the result from a previous step apply 'Field calculator' with a small expression for

CASE
   WHEN "Name" IS NULL THEN 0
   ELSE 1
END

Alternative is: if("Description", 1, 0)

step_4

So, the updated Attribute table will now look as demonstrated on the image below

new_at

Step 5. On this last stage I will offer not the best solution in terms of performance for large data but I reckon it is simple to understand. Attaining a basic statistics by means of a "Virtual Layer" through Layer > Add Layer > Add/Edit Virtual Layer... apply this query

WITH o AS (
    SELECT COUNT() AS "overlap"
    FROM "Joined layer"
    WHERE "calc" = 1
),

no AS (
    SELECT COUNT() AS "not_overlap"
    FROM "Joined layer"
    WHERE "calc" = 0
)

SELECT "o"."overlap", "no"."not_overlap"
FROM "o", "no"

And your output is

output

Note: that the output does not include geometry.

Additionally I am providing you with some references that are necessary to read.


References:

Answered by Taras on October 1, 2021

A solution using PyQGIS.

Let's assume there is a project folder with two files: 'test_points.csv' and 'test_polygons.shp', see image below.

input

In QGIS they will look like

qgis

Proceed with Plugins > Python Console > Show Editor and copy&edit&paste the script below

# defining inputs
project_path = 'C:/Users/DUBRTARA/Desktop/test/project/'
csv_file_name = 'test_points.csv'
shp_file_name = 'test_polygons.shp'

# Step 1: reading a csv file with points and converting it into a shapefile
uri_csv_file = 'file:///' + project_path + csv_file_name + "?encoding={0}&delimiter={1}&xField={2}&yField={3}&crs={4}".format("UTF-8",",","x","y","epsg:31469")
# In .format("UTF-8",",","x","y","epsg:31469")
# "UTF-8" stands for encoding : data encoding (optional)
# "," is a delimiter inm the input file
# xField  is a column name for longitude value
# yField is a column name for latitude value
# crs is a Coordinate system in EPSG number
points = QgsVectorLayer(uri_csv_file, csv_file_name.split(".")[0], "delimitedtext")
if not points.isValid():
    print ("{} layer was not loaded".format(csv_file_name))

# Step 2: reading a shapefile with polygons
path_to_shp_file_name = project_path + shp_file_name
polygons = QgsVectorLayer(path_to_shp_file_name, shp_file_name.split(".")[0], "ogr")
if not polygons.isValid():
    print ("{} layer was not loaded".format(shp_file_name))

# Step 3: selecting point features that intersect polygon features
points_in = processing.run("qgis:selectbylocation",{
                                                    'INPUT':points,
                                                    'PREDICATE':0,#corresponds to 'intersects'
                                                    'INTERSECT':polygons,
                                                    'METHOD':0
                                                    })['OUTPUT']

# creating a dictionary for output
result = {}
# nesting results into the dictionary
result['points in polygons'] = points_in.selectedFeatureCount()
result['points not in polygons'] = points.featureCount() - points_in.selectedFeatureCount()

# printing the result
print(result)

# for more details run the following command >>> processing.algorithmHelp("qgis:selectbylocation")

Press Run script run script and get the output that will look like {'points in polygons': 4, 'points not in polygons': 1}, see image below

result


References:

Answered by Taras on October 1, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP