Finding which side of closest polyline point lies on in QGIS?

Question

I would like to find on which side of the line (left or right) my points are located with respect to the direction of the line (like knowing the rivers bank along a meandering river). But, if there are several line, I am only interested in the closest line.

I know that I could use one sided buffers like it is suggested in this post, but this doesn't work well when several lines are close to each other near the point.

For instance, there are sometimes two (or more) lines close to my point of interest, which then create two (or more) overlapping polygons covering this point. In this case, the point could be on the left of one line and on the right of the other line.

I only want to know the side of the closest road to the point (avoiding multiple buffers and selections because the database is very large. A single query maybe)

Answer 1

Not a single query but the following code can be executed in the Python Console. The code essentially does the following:

• Iterates each line feature and calculates the distance from this to a selected point feature;
• Stores these in a dictionary;
• A proximity is used to limit iterating through each line feature (i.e. a value of 100 means only lines within 100m of the point will be included);
• Sorts the dictionary and selects the line feature with the minimum distance from the point;
• Gets the coordinates of the start and end vertices;
• Determines if the point lies on the left- or right-hand side of the line (assuming the line is a single straight line).

---

pointLayer = QgsProject.instance().mapLayersByName('point')[0]
lineLayer = QgsProject.instance().mapLayersByName('line')[0]

# Set proximity limit
proximity = 1
distance_dict = {}

for point in pointLayer.selectedFeatures():
    start_point = point.geometry().asPoint().x()
    end_point = point.geometry().asPoint().y()
    for line in lineLayer.getFeatures():
        distance = point.geometry().distance(line.geometry())
        if distance < proximity:
            distance_dict[line.id()] = distance
    closestFeature = sorted(distance_dict.keys(), key=(lambda key:  distance_dict[key]))
    lineLayer.selectByIds([closestFeature[0]])
    for line in lineLayer.selectedFeatures():
        start_geom = line.geometry().asPolyline()[0]
        end_geom = line.geometry().asPolyline()[1]
    d = ((start_point - start_geom.x()) * (end_geom.y() - start_geom.y())) - (((end_point - start_geom.y()) * (end_geom.x() - start_geom.x())))
    if d < 0: 
        print('Left')
    elif d > 0:
        print('Right')
    else:
        print('On line')

---

Example of selected points (highlighted in yellow):

• Right of line -

  

• Left of line -

  

Answer 2

@Joseph's answer is great for segment and has implemented an algorithm. However, I've found that clostsSegmentWithContext() could be used to know if the point is on the left (result[3]=-1) or on the right (result[3]=1). Note that I first screen out the possible closest lines with spatial index, which returns the closts lines based on boxes. Therefore I a a loop on the 4 closest lines, more lines could be necessary in more complex situations (but of cours it would take longer). If there are many point but only a few lines, the loop can be done systematically on each line (without "preselecting" with spatial index)

pointlayer = QgsProject.instance().mapLayersByName('p')[0] #iface.mapCanvas().layer(0)
lineLayer = QgsProject.instance().mapLayersByName('lines')[0] # iface.mapCanvas().layer(1)

testIndex = QgsSpatialIndex(lineLayer)

pointlayer.startEditing() #edit is necessary to be able to update the points

feats=[]

for p in pointlayer.getFeatures():
    mindist = 10000. #initialize a large mindist value
    near_ids = testIndex.nearestNeighbor(p.geometry().asPoint(),4)
    features = lineLayer.getFeatures(QgsFeatureRequest().setFilterFids(near_ids))
    for tline in features: #loop on 4 closest line candidtates
        closeSegResult = tline.geometry().closestSegmentWithContext(p.geometry().asPoint())
        if mindist > closeSegResult[0]:
            mindist = closeSegResult[0]
            p['leftright'] = closeSegResult[3] #update the side of line in 'leftRight' field for the closest segment 
    pointlayer.updateFeature(p)

pointlayer.commitChanges() #close edit session to flush changes