laspy filter

laspy_filter1.py

import laspy
import numpy as np
import sys


las = laspy.read('/Users/hugo/teaching/hw04-marking/data/area_full.laz')

#-- fetch all points within 100m of the median xyz of the file
coords = np.vstack((las.x, las.y, las.z)).transpose()
median = np.median(coords, axis=0)
distances = np.sqrt(np.sum((coords - median) ** 2, axis=1))
mask1 = distances < 100
# mask = (las.x >= min_x) & (las.x <= max_x) & (las.y >= min_y) & (las.y <= max_y)
# return np.column_stack((las.x[mask], las.y[mask], las.z[mask])) 

#-- create a new file
new_file = laspy.create(point_format=las.header.point_format, file_version=las.header.version)
new_file.points = las.points[mask1]

#-- filter to keep only the ground
new_file.points = new_file.points[new_file.classification == 2]
new_file.write("good_points.laz")