0

I have downloaded some .las files from a LiDAR project listed on the USGS website. I am hoping to process them in Python. I was able to get the set of X,Y, and Z points using the following snippit in the documentation.

import laspy as lp
import numpy as np

las_file = r'lidar.las'
inFile = lp.file.File(las_file)

dataset = np.vstack([inFile.X, inFile.Y, inFile.Z]).transpose() #X, Y, Z data

The problem I was running into is that I didn't understand how to get a subset of the data using Latitude and Longitude coordinates. I have explored the data in the header and I see that the data is scaled by scaling factors but I still was confused how to scale this to a coordinate system that I would expect to see on a map.

>>>inFile.header.min
[504110.62, 4828500.0, 2728.0]
>>>inFile.header.max
[505499.99, 4829999.99, 3142.9900000000002]
>>>inFile.header.scale
[0.01, 0.01, 0.001]

I looked at the information sheet that was also downloaded with the .las file and it included the following information. I think this is trying to tell me how to convert the X,Y,Z data to Latitude,Longitude, and elevation but I didn't understand how to make this transformation in python.

Data that came with the .las file

Is there a package or library or custom function that I can implement that can assist me to make the conversion to latitude and longitude coordinates so that I will be able to select points from this LiDAR point cloud that are inside of a polygon defined as a series of latitude/longitude points?

3
  • 1
    Your file CRS seems to be UTM zone 12. The only way to extract point within a polygon is to project the polygon in the same CRS than the point cloud. – JRR Mar 1 at 12:41
  • Thanks for the comment. I have been reading up on the different coordinate systems and your remark seems to have pointed me in the right direction. It seems I need to do something like this – sTr8_Struggin Mar 1 at 12:53
  • I'm not in python but I'm sure it exists geospatial libraries to perform CRS transformation in python. But you should not do the math yourself for sure. Can also transform your polygon in a GIS software. – JRR Mar 1 at 12:58
0

I am sharing what I ended up doing here in case it helps someone else.

I realized that I needed a coordinate transformation tool and I found one that seemed quite lightweight here

To use it on my dataset I did the following:

lat_lon = utm.to_latlon(grid[:,0]*0.01, grid[:,1]*0.01, 12, nothern=True)

This worked well enough for my needs.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.