Here is Python solution (for QGIS Python console) which is mostly created with the help of the following links:
- Create Intersecting lines and remove dangles
- How do I find vector line bearing in QGIS or GRASS?
- How to create points in a specified distance along the line in QGIS?
- Redirecting stdout to “nothing” in python
- QGIS 2.18 Network analysis library
Following steps are involved:
- Find point on nearest road to given point and create a linestring.
- Extend the linestring by a small distance to make sure the linestring is intersecting the road. (I recognized that there are less intersections without using the extended linestring.)
- Finding the shortest path from the intersection point to the linestring start point using
dijkstra()
method from the QGIS network analysis library.
Test data:
Road layer with 7220 features:

Random points layer with 7220 features (using the Toolbox: Bounding boxes
from road layer, Random points in layer bounds
(with minimum distance = 10
))

Python code:
import math
import os
import sys
from contextlib import contextmanager
from PyQt4.QtCore import *
from PyQt4.QtGui import *
from qgis.core import *
from qgis.gui import *
from qgis.networkanalysis import *
# define input layer points and roads
p_lyr = QgsMapLayerRegistry.instance().mapLayersByName('points')[0]
l_lyr = QgsMapLayerRegistry.instance().mapLayersByName('roads')[0]
# function to create the azimuth of line segments
def azimuth(point1, point2):
return point1.azimuth(point2)
# project a point in azimuth direction using direction cosines
def cosdir_azim(azim):
azim = math.radians(azim)
cosa = math.sin(azim)
cosb = math.cos(azim)
return cosa,cosb
# used to get line segments
def pair(list):
for i in range(1, len(list)):
yield list[i-1], list[i]
# create function to hide QGIS output ("True") while creating QgsSpatialIndex for every feature (takes a while for x*1000 features)
@contextmanager
def silence_stdout():
new_target = open(os.devnull, "w")
old_target, sys.stdout = sys.stdout, new_target
try:
yield new_target
finally:
sys.stdout = old_target
# creating index silently
with silence_stdout():
lines = [feature for feature in l_lyr.getFeatures()]
lines_spIndex = QgsSpatialIndex()
for elem in lines:
lines_spIndex.insertFeature(elem)
# set up memory layer for the normal vector
d_lyr = QgsVectorLayer('LineString', 'normalVector', 'memory')
QgsMapLayerRegistry.instance().addMapLayer(d_lyr)
prov = d_lyr.dataProvider()
# adding three attributes (holding point_id, road_id and the distance)
prov.addAttributes( [ QgsField("point_id", QVariant.Int), QgsField("road_id", QVariant.Int), QgsField("distance",QVariant.Int)])
# function to create normal vector (geometry) and distance
def normalVector(distance):
vect = []
feat = []
for points in p_lyr.getFeatures():
# find closest point to line
minDistPoint = min([l.geometry().closestSegmentWithContext(QgsPoint(points.geometry().asPoint())) for l in lines])[1]
feat = QgsFeature()
# create line from point to minDistPoint
feat.setGeometry(QgsGeometry.fromPolyline([QgsPoint(points.geometry().asPoint()), QgsPoint(minDistPoint[0], minDistPoint[1])]))
line = feat.geometry().asPolyline()
for seg_start, seg_end in pair(line):
cosa, cosb = cosdir_azim(azimuth(seg_start, seg_end))
lenght = distance
feat = QgsFeature()
# extend line for intersection with road feature
feat.setGeometry(QgsGeometry.fromPolyline([seg_start,QgsPoint(seg_end.x()+(lenght*cosa), seg_end.y()+(lenght*cosb))]))
for id in lines_spIndex.intersects(feat.geometry().boundingBox()):
sp_geom = lines[id].geometry().asPolyline()
start_point = QgsPoint(sp_geom[0])
sp_x, sp_y = start_point[0], start_point[1]
if feat.geometry().intersects(lines[id].geometry()):
intersection_point = feat.geometry().intersection(lines[id].geometry()).asPoint()
end_point = QgsPoint(intersection_point)
ep_x, ep_y = end_point[0], end_point[1]
# building graph and calculating shortest path
director = QgsLineVectorLayerDirector(l_lyr, -1, '', '', '', 3)
properter = QgsDistanceArcProperter()
director.addProperter(properter)
crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
builder = QgsGraphBuilder(crs)
pStart = QgsPoint(ep_x, ep_y )
pStop = QgsPoint(sp_x, sp_y)
tiedPoints = director.makeGraph(builder, [pStart, pStop])
graph = builder.graph()
tStart = tiedPoints[0]
tStop = tiedPoints[1]
idStart = graph.findVertex(tStart)
idStop = graph.findVertex(tStop)
(tree, cost) = QgsGraphAnalyzer.dijkstra(graph, idStart, 0)
if tree[idStop] == -1:
print "Path not found"
else:
p = []
curPos = idStop
while curPos != idStart:
p.append(graph.vertex(graph.arc(tree[curPos]).inVertex()).point())
curPos = graph.arc(tree[curPos]).outVertex()
p.append(tStart)
vect = QgsFeature()
vect.setGeometry(QgsGeometry.fromPolyline(p))
vect.setAttributes([int(points["id"]), int(lines[id]["id"]), vect.geometry().length()])
prov.addFeatures([vect])
d_lyr.updateExtents()
d_lyr.triggerRepaint()
d_lyr.updateFields()
# run function
normalVector(0.1)
Results:
For the amount of data as described above it takes about 16 minutes of processing time. For 11 times I got Path not found
.

QGIS Processing Toolbox | GRASS 7 | Vector
. (1) You will obtain newNetwork
line layer. (2) Then you can select the first segment which represents Start - Intersection interval. (3) Save them as new layer. (4) Calculate$length
.