Преобразование файла LAS в массив NumPy?

15

Я начал изучать, как манипулировать данными LAS в python, и хотел посмотреть, как другие обрабатывают файлы LAS. Я хотел бы прочитать точки (я использую массив NumPy) и отфильтровать классы 1 и 2 (неклассифицированные и заземленные) в отдельный массив. У меня есть следующий код, но я не могу отфильтровать очки.

# 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

Я видел arcpy.da.featureClassToNumpyArray, но я не хотел импортировать arcpy или преобразовывать его в шейп-файл.

Как еще я могу отфильтровать / прочитать данные LAS в массив Numpy?

Barbarossa
источник
Что такое сообщение об ошибке (если есть)?
til_b
Нет ошибок. Я просто не знал, как фильтровать, и не был уверен, есть ли лучший способ получить LAS в массив.
Барбаросса

Ответы:

14

Вы PointsXYZICтеперь NumPy массив. Это означает, что вы можете использовать индексирование NumPy для фильтрации данных, которые вас интересуют. Например, вы можете использовать индекс логических значений, чтобы определить, какие точки нужно захватить.

#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]

Теперь у вас должен быть массив Numpy со всеми значениями, в которых данные не классифицированы или заземлены. Чтобы получить значения, которые были классифицированы, вы можете использовать:

filter_array = np.all(
    [
        PointsXYZIC[:, 4] != unclassified, #The final column to index against
        PointsXYZIC[:, 4] != ground,
    ],
    axis=0
)
om_henners
источник
Фильтр работает, но записывает только 5 записей. Я попытался отфильтровать только классы 1 и 2, а затем попытался отфильтровать все, кроме 1 и 2, и оба дали мне только 5 результатов. Есть идеи?
Барбаросса
Эти 5 записей находятся в 1-м массиве.
Барбаросса
Извините, обновленный код выше, так как он требует указания оси для выполнения любых вычислений (без него он выполняет или во всех измерениях массива).
om_henners
5

Используйте laspy для чтения файлов LAS и легко возвращайте данные в виде пустых массивов, с которыми вы можете взаимодействовать. laspy - это чистый Python, он почти такой же быстрый, как libLAS, имеет больше возможностей, чем привязки libLAS Python, и его гораздо проще развернуть.

Говард Батлер
источник
0

Я прошу прощения, если вы уже знаете об этом, но LASTools - это фантастический инструмент с открытым исходным кодом, который теперь интегрируется как с ArcGIS, так и с QGIS 2.0. У него есть опции для фильтрации данных так, как вы смотрите.

Николас Дагган
источник
Спасибо @Nicholas, я использую библиотеку liblas python, которая тесно связана с LASTools
Barbarossa