Finding points in LAS file by coordinates

Question

I am trying to find points in a LAS file (Lidar format):

Right now I am doing it the really slow way:

from laspy.file import File
import numpy as np

inFile = File('inputfile.las', mode='r')

coord = np.vstack((inFile.x, inFile.y, inFile.z)).transpose()

def find_pt(coord, x, y, z):
    found = []
    for i in coord:
        if(i[0] >= x and i[0] < x+1):
            if(i[1] >= y and i[1] < y+1):
                if(i[2] >= z and i[2] < z+1):
                    found.append(i)

    return found

Then I call it with:

find_pt(coord, 358397, 5280527, 550)

Which of course takes some time, especially when there are a lot of points in the file.

Is there a better/faster way? coords is of type numpy.ndarray

Answer 1

I don't know las or laspy, but if inFile.x & friends are separate numpy arrays (judging from the vstack call), you could simply use and compare them separately, then combine the results:

xgood = (inFile.x >= x) & (inFile.x < x+1)
ygood = (inFile.y >= y) & (inFile.y < y+1)
zgood = (inFile.z >= z) & (inFile.z < z+1)
good = xgood & ygood & zgood
found = (inFile.x[good], inFile.y[good], inFile.z[good])

The resulting found is slightly different, as it is a tuple of x, y and z coordinate arrays for valid point(s).