Rotar un vectorial tipo línea, dados el ángulo y el eje de rotación, con PyQGIS

En el post anterior se presentó un código que permitía rotar vectoriales tipo polígono dados el ángulo y el eje de rotación. El código de este post complementa al anterior y hace posible la rotación de un vectorial tipo línea utilizando los mismos eje y ángulo de rotación empleados para el vectorial tipo polígono. El código es el siguiente:

from math import sin, cos, pi

def rotatePoints(axis, r_angle, points):
    """Rotate points around an axis by a fixed angle. It uses linear algebra
    concepts to set a new reference system at axis point and rotation 
    equations for x,y plane around axis point. After rotation, converts 
    rotated points to (0,0) system reference newly. 
        
    :param axis: Rotation axis taken as (0,0) reference
    :type axis: QgsPoint
        
    :param r_angle: Rotation angle around axis in degrees
    :type r_angle: float
        
    :param points: Points to be rotated around the axis by r_angle
    :type points: List of no rotated QgsPoint
        
    :returns: Rotated points around the elected axis by r_angle
    :rtype: List of rotated QgsPoint
    """
        
    points_rot = []

    for point in points:
        #creating coordinate reference system at axis point 
        tr = QgsPoint( point.x() - axis.x(), point.y() - axis.y() )
        #rotation equations around axis point at x,y plane
        xrot = tr.x()*cos(r_angle*(pi/180)) - tr.y()*sin(r_angle*(pi/180))  
        yrot = tr.x()*sin(r_angle*(pi/180)) + tr.y()*cos(r_angle*(pi/180))
        #restauring coordinate reference system at origin point (0,0)
        xnew = xrot + axis.x()
        ynew = yrot + axis.y()
        #appending rotated point a points_rot list
        points_rot.append(QgsPoint(xnew,ynew))
        
    return points_rot

def createMemoryLayer(geom_idx, epsg, points):
    """Create a memory layer with rotated points taking in account
        its geometry type and each QgsFeature
        
    :param geom_idx: QGis.Point (0), QGis.Line (1) or QGis.Polygon (2) enumerate type
    :type geom_idx: Integer
        
    :param epsg: European Petroleum Survey Group Id (Coordinate Reference System)
    :type epsg: Integer
        
    :param points: Rotated points with "rotatePoints" function
    :type points: Nested List of QgsPoints according its geometry type (Point, Line or
        Polygon)
    """
    CRS = QgsCoordinateReferenceSystem(epsg)
        
    if geom_idx == QGis.Point:
        geom = "Point"
        label = "Points"
        
    if geom_idx == QGis.Line:
        geom = "LineString"
        label = "Line"

    if geom_idx == QGis.Polygon:
        geom = "Polygon"
        label = "Polygon"

    uri = geom + "?crs=epsg:" + str(epsg) + "&field=id:integer""&index=yes"

    mem_layer = QgsVectorLayer(uri,
                               label,
                               "memory")

    #add Map Layer to Registry
    QgsMapLayerRegistry.instance().addMapLayer(mem_layer)
                
    #Prepare mem_layer for editing
    mem_layer.startEditing()

    #Calculate number of features: points could be a nested list depending of its geometry 
    n = len(points)
          
    #Set feature
    feature = []
            
    for i in range(n):
        feat =QgsFeature()
        feature.append(feat)
        
    #set attributes values 
    if geom_idx == QGis.Point:
        for i in range(n):
            feature[i].setGeometry(QgsGeometry.fromPoint(points[i]))  #Set geometry
            feature[i].setAttributes([i])
            mem_layer.addFeature(feature[i], True)

    elif geom_idx == QGis.Line:
        #Set geometry
        for i in range(n):
            feature[i].setGeometry(QgsGeometry.fromPolyline(points[i]))
            #set attributes values
            feature[i].setAttributes([i])
            mem_layer.addFeature(feature[i], True)

    elif geom_idx == QGis.Polygon:
        #Set geometry
        for i in range(n):
            feature[i].setGeometry(QgsGeometry.fromPolygon(points[i]))
            #set attributes values
            feature[i].setAttributes([i])
            mem_layer.addFeature(feature[i], True)

    #stop editing and save changes
    mem_layer.commitChanges()

#Code starts here

layer = iface.activeLayer()

crs = layer.crs()
epsg = crs.postgisSrid()

features = layer.getFeatures()

geom = []

for feature in features:
    geom.append(feature.geometry().asPolyline())
    
points = geom

n = layer.featureCount()

line_rot = []

axis = QgsPoint(416047.568418, 4443411.85665)

r_angle = 70
    
for i in range(n):
    points_rot = rotatePoints(axis, r_angle, points[i])
    line_rot.append(points_rot)

createMemoryLayer(QGis.Line, epsg, line_rot)

Se ha usado el estilo “docstring” para la documentación y generalizado la función que crea la “memory layer” para que pueda ser aplicada a los vectoriales independientemente de su tipo de geometría (Punto, Línea o Polígono).

Cuando se aplica el recipe anterior a un vectorial de línea (arbitrario) con tres features se obtiene el resultado siguiente:

rotate_line

Después de aplicarlo 9 veces consecutivas adicionales a la memory layer resultante de cada ejecución (parece arte🙂 ):

rotate_line2

Después de codificado en un plugin, el código rota el vectorial original 180º con relación a un punto exterior situado a la izquierda de éste:

rotate_linea3

En nuestro curso de PyQGIS se contemplan con detalle el desarrollo de éste y otros plugins mediante el empleo de ‘Plugin Builder’.

PRÓXIMA SESIÓN EL 20 DE OCTUBRE

Esta entrada fue publicada en PyQGIS, QGIS, SIG, Software Libre. Guarda el enlace permanente.

2 respuestas a Rotar un vectorial tipo línea, dados el ángulo y el eje de rotación, con PyQGIS

  1. Pingback: Rotación de archivos vectoriales con método ‘rotate’ de QgsGeometry en PyQGIS | El Blog de José Guerrero

  2. Pingback: Newsletter Octubre 2015 | El Blog de José Guerrero

Responder

Por favor, inicia sesión con uno de estos métodos para publicar tu comentario:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s