1

I have two layers: point and line (postgis or shp) The line layer represent a hiking trail. The point layer the mountain passes. The points are exactly on the same position like the node of the line. Now I would like to rotate the points in orientation of the line. Is there a way to do that automatically?

enter image description here

0

1 Answer 1

8

You need to determine the azimuth of the polyline segment next to the point. It is possible to calculate this value and save it to an attribute, or calculate it as an expression which give the azimuth as result.

How to use the following Python script:

  1. Point layer and line layer are added to the legend (in the script named bridges and trails) and change the names in the code
  2. Add an attribut of type REAL to the point layer (in the script named angle), when using a different name, change the code in line 8
  3. Open Python Console -> Editor and copy & paste the code
  4. Run the script
from math import atan2

# get layers
bridges = QgsMapLayerRegistry.instance().mapLayersByName('bridges')[0]
trails = QgsMapLayerRegistry.instance().mapLayersByName('trails')[0]

# get fieldindex
fni = bridges.fieldNameIndex('angle')

# initializations
tolerance = 0.01  # search tolerance
bridges.startEditing()

# loop over all bridges
for bridge in bridges.getFeatures():
    x = bridge.geometry().asPoint().x()
    y = bridge.geometry().asPoint().y()

    # get the rectangular search area 
    searchRect = QgsRectangle(x - tolerance, y - tolerance,  x + tolerance, y + tolerance)

    # find trails 
    for trail in trails.getFeatures(QgsFeatureRequest().setFilterRect(searchRect)):
        # get the nearest vertex on trail and the one before and after
        pnt, v, b, a, d = trail.geometry().closestVertex(bridge.geometry().asPoint())
        p1 = trail.geometry().vertexAt(v)
        # when vertex before exists look back, otherwise look forward
        if v>-1 and b>-1:
            p2 = trail.geometry().vertexAt(b)
        elif v>-1 and a>-1:
            p2 = trail.geometry().vertexAt(a)

        # calculate azimuth
        angle = atan2(p2.x() - p1.x(), p2.y() - p1.y()) / 0.017453
        bridge[fni] = angle
        bridges.updateFeature(bridge)

# save changes and stop editing
bridges.commitChanges()

The value for variable tolerance depends on the units of length measurements. Tt may be necessary to choose the appropriate value to really pick the nearest vertex.

Use this attribut angle as variable for data defined rotation of the symbols:

rotated symbols

You may wish to complicate things a little bit and calculate azimuths to the vertex before AND after and take the average.

The dynamic expression version:

from qgis.core import *
from qgis.gui import *
from math import atan2

@qgsfunction(args='auto', group='Custom')
def get_azimuth(layername, tolerance, feature, parent):
    # initializations
    trails  = QgsMapLayerRegistry.instance().mapLayersByName(layername)[0]
    x = feature.geometry().asPoint().x()
    y = feature.geometry().asPoint().y()

    # get the rectangular search area 
    searchRect = QgsRectangle(x - tolerance, y - tolerance,  x + tolerance, y + tolerance)

    # find trails 
    for trail in trails.getFeatures(QgsFeatureRequest().setFilterRect(searchRect)):
        # get the nearest vertex on trail and the one before and after
        pnt, v, b, a, d = trail.geometry().closestVertex(feature.geometry().asPoint())
        if v == -1:
            continue
        p1 = trail.geometry().vertexAt(v)
        # when vertex before exists look back, otherwise look forward
        if v>-1 and b>-1:
            p2 = trail.geometry().vertexAt(b)
        elif v>-1 and a>-1:
            p2 = trail.geometry().vertexAt(a)
        # calculate azimuth
        angle = atan2(p2.x() - p1.x(), p2.y() - p1.y()) / 0.017453
        return angle

expression

field calculator

9
  • Looks pretty nice! But I'm sorry, I don't know how to run that script in qgis. I would like to get the result directly as expression in the Rotation-Filed.
    – MartinMap
    Commented Jan 29, 2016 at 19:18
  • 1
    @MAP I added the answer, hope you get it. Sometimes I had to wait several seconds for the symbols to get correct rotated. It named the layer bridges but you seem to be interested in passes, don't you? But the principle remains the same.
    – Detlev
    Commented Jan 29, 2016 at 20:21
  • Great, works like it should! You've saved me a lot of work! Thanks!
    – MartinMap
    Commented Jan 30, 2016 at 11:47
  • I've used this technique to build an expression, but it doesn't quite work correctly. The cross lines (point symbols) should be perpendicular to the stream, and some of them appear correct while others are closer to parallel. Any idea where my problem is? i.sstatic.net/0r8Wv.png
    – Nathan
    Commented Oct 30, 2017 at 4:15
  • 1
    The variable tolerance is used to search lines near to the points. Choose this value as small as possible, otherwise the script has to be changed to evaluate for the closest vertex in the search rectangle. With the data you mentioned (in NAD83) I get reasonable results with tolerance equal to 0.0001. Keep in mind that this value depends on the projection respectively length units used.
    – Detlev
    Commented Oct 30, 2017 at 17:17

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