¿Convertir archivo LAS a matriz numpy?

15

He comenzado a aprender a manipular los datos LAS en python y quería ver cómo otros manejan los archivos LAS. Me gustaría leer los puntos (estoy usando una matriz numpy), y filtrar las clases 1 y 2 (sin clasificar y suelo) a una matriz separada. Tengo el siguiente código, pero parece que no puedo filtrar los puntos.

# Import modules
from liblas import file
import numpy as np

if __name__=="__main__":
    '''Read LAS file and create an array to hold X, Y, Z values'''
    # Get file
    las_file = r"E:\Testing\ground_filtered.las"
    # Read file
    f = file.File(las_file, mode='r')
    # Get number of points from header
    num_points = int(f.__len__())
    # Create empty numpy array
    PointsXYZIC = np.empty(shape=(num_points, 5))
    # Load all LAS points into numpy array
    counter = 0
    for p in f:
        newrow = [p.x, p.y, p.z, p.intensity, p.classification]
        PointsXYZIC[counter] = newrow
        counter += 1

He visto arcpy.da.featureClassToNumpyArray, pero no quería importar arcpy ni tener que convertir a shapefile.

¿De qué otra forma puedo filtrar / leer datos LAS en una matriz numpy?

Barbarroja
fuente
¿Cuál es el mensaje de error (si lo hay)?
til_b
No hay error. Simplemente no sabía cómo filtrar, y no estaba seguro de si había una mejor manera de poner LAS en la matriz.
Barbarroja

Respuestas:

14

Tu PointsXYZICahora es una matriz numpy. Lo que significa que puede usar indexación numpy para filtrar los datos que le interesan. Por ejemplo, puede usar un índice de booleanos para determinar qué puntos tomar.

#the values we're classifying against
unclassified = 1
ground = 2

#create an array of booleans
filter_array = np.any(
    [
        PointsXYZIC[:, 4] == unclassified, #The final column to index against
        PointsXYZIC[:, 4] == ground,
    ],
    axis=0
)

#use the booleans to index the original array
filtered_rows = PointsXYZIC[filter_array]

Ahora debe tener una matriz numpy con todos los valores en los que los datos no están clasificados o son tierra. Para obtener los valores que se han clasificado, puede usar:

filter_array = np.all(
    [
        PointsXYZIC[:, 4] != unclassified, #The final column to index against
        PointsXYZIC[:, 4] != ground,
    ],
    axis=0
)
om_henners
fuente
El filtro parece funcionar pero solo escribe 5 registros. Traté de filtrar solo las clases 1 y 2, y luego intenté filtrar todas menos 1 y 2, ambas me dieron solo 5 resultados. ¿Algunas ideas?
Barbarroja
Estos 5 registros están en una matriz 1-d.
Barbarroja
Lo sentimos, se actualizó el código anterior, ya que requiere la especificación del eje para realizar cualquier cálculo (sin que realice todas las dimensiones de la matriz).
om_henners
5

Use laspy para leer archivos LAS y devolver fácilmente los datos como matrices numpy con las que puede interactuar. Laspy es Python puro, es casi tan rápido como libLAS, tiene más funciones que los enlaces libLAS Python y es mucho más fácil de implementar.

Howard Butler
fuente
0

Pido disculpas si ya sabe de esto, pero LASTools es una fantástica herramienta de código abierto que ahora se integra con ArcGIS y QGIS 2.0: tiene opciones para filtrar los datos de la manera que está viendo.

Nicholas Duggan
fuente
Gracias @Nicholas, estoy usando la biblioteca liblas python, que está estrechamente vinculada a LASTools
Barbarossa