2

I'm trying to get this script to work by @gene but I'm not sure how to load the results into QGIS. I tried to make a memory layer but not sure how to add the geometries of the point, line and the intersection point:

point_layer = QgsVectorLayer("Point", "point_layer", "memory")
pr = point_layer.dataProvider()
fet = QgsFeature()
fet.setGeometry(QgsGeometry.fromPoint(QgsPoint(30,120)))

line_layer = QgsVectorLayer("LineString", "line_layer", "memory")
lr = line_layer.dataProvider()
seg = QgsFeature()
line_start = QgsPoint(50,50)
line_end = QgsPoint(100,150)
line = QgsGeometry.fromPolyline([line_start,line_end])

def intersect_point_to_line(point, line_start, line_end):
    magnitude2 = line_start.sqrDist(line_end)
    u = ((point.x() - line_start.x()) * (line_end.x() - line_start.x()) + (point.y() - line_start.y()) * (line_end.y() - line_start.y()))/(magnitude2)
    ix = line_start.x() + u * (line_end.x() - line_start.x())
    iy = line_start.y() + u * (line_end.y() - line_start.y())
    return QgsPoint(ix,iy) 

line = QgsGeometry.fromPolyline([point,intersect_point_to_line(point, line_start, line_end)])

I am using QGIS 2.10.

Community
  • 1
user35594
  • 555
  • 1
  • 5
  • 22

1 Answers1

2

You're getting there! Some things to note:

  • Point layer

    You defined the coordinate of the point feature but you still need to add it to the point memory layer. Also, the script uses point to calculate the intersection point which you would need to define. Therefore, replace fet.setGeometry(QgsGeometry.fromPoint(QgsPoint(30,120))) with this:

    point = QgsPoint(30,120)
fet.setGeometry(QgsGeometry.fromPoint(point))
pr.addFeatures([fet])

  • Line layer

    Again, you need to add the line feature to the line memory layer so replace line = QgsGeometry.fromPolyline([line_start,line_end]) with this:

    seg.setGeometry(QgsGeometry.fromPolyline([line_start, line_end]))
lr.addFeatures([seg])

  • Intersect layer

    Personally, I like my results to be stored in another shapefile so you could create another line memory layer to show the intersection line. You already created one so you can use the same method to create another and instead of adding coordinates, you would need to add the intersection line from the def intersect_point_to_line function:

    intersect_layer = QgsVectorLayer("LineString?crs=epsg:4326", "intersect_layer", "memory")
ir = intersect_layer.dataProvider()
seg = QgsFeature()
seg.setGeometry(QgsGeometry.fromPolyline([point,intersect_point_to_line(point, line_start, line_end)]))
ir.addFeatures([seg])

Here is the complete script which you should be able to directly copy/paste into your Python Console:

point_layer = QgsVectorLayer("Point?crs=epsg:4326", "point_layer", "memory")
pr = point_layer.dataProvider()
fet = QgsFeature()
point = QgsPoint(30,120)
fet.setGeometry(QgsGeometry.fromPoint(point))
pr.addFeatures([fet])

line_layer = QgsVectorLayer("LineString?crs=epsg:4326", "line_layer", "memory")
lr = line_layer.dataProvider()
seg = QgsFeature()
line_start = QgsPoint(50,50)
line_end = QgsPoint(100,150)
seg.setGeometry(QgsGeometry.fromPolyline([line_start, line_end]))
lr.addFeatures([seg])

def intersect_point_to_line(point, line_start, line_end):
    magnitude2 = line_start.sqrDist(line_end) 
    u = ((point.x() - line_start.x()) * (line_end.x() - line_start.x()) + (point.y() - line_start.y()) * (line_end.y() - line_start.y()))/(magnitude2)
    ix = line_start.x() + u * (line_end.x() - line_start.x())
    iy = line_start.y() + u * (line_end.y() - line_start.y())
    return QgsPoint(ix,iy)
    print QgsPoint(ix,iy)

intersect_layer = QgsVectorLayer("LineString?crs=epsg:4326", "intersect_layer", "memory")
ir = intersect_layer.dataProvider()
seg = QgsFeature()
seg.setGeometry(QgsGeometry.fromPolyline([point,intersect_point_to_line(point, line_start, line_end)]))
ir.addFeatures([seg])

QgsMapLayerRegistry.instance().addMapLayers([point_layer, line_layer, intersect_layer])

Result

Joseph
  • 75,746
  • 7
  • 171
  • 282