diff --git a/opendm/dem/commands.py b/opendm/dem/commands.py index 56b80b05..47903a68 100755 --- a/opendm/dem/commands.py +++ b/opendm/dem/commands.py @@ -50,36 +50,13 @@ def rectify(lasFile, debug=False, reclassify_threshold=5, min_area=750, min_poin start = datetime.now() try: - # Currently, no Python 2 lib that supports reading and writing LAZ, so we will do it manually until ODM is migrated to Python 3 - # When migration is done, we can move to pylas and avoid using PDAL for conversion - tempLasFile = os.path.join(os.path.dirname(lasFile), 'tmp.las') - - # Convert LAZ to LAS - cmd = [ - 'pdal', - 'translate', - '-i %s' % lasFile, - '-o %s' % tempLasFile - ] - system.run(' '.join(cmd)) - log.ODM_INFO("Rectifying {} using with [reclassify threshold: {}, min area: {}, min points: {}]".format(lasFile, reclassify_threshold, min_area, min_points)) run_rectification( - input=tempLasFile, output=tempLasFile, debug=debug, \ + input=lasFile, output=lasFile, debug=debug, \ reclassify_plan='median', reclassify_threshold=reclassify_threshold, \ extend_plan='surrounding', extend_grid_distance=5, \ min_area=min_area, min_points=min_points) - # Convert LAS to LAZ - cmd = [ - 'pdal', - 'translate', - '-i %s' % tempLasFile, - '-o %s' % lasFile - ] - system.run(' '.join(cmd)) - os.remove(tempLasFile) - except Exception as e: raise Exception("Error rectifying ground in file %s: %s" % (lasFile, str(e))) diff --git a/opendm/dem/ground_rectification/extra_dimensions/distance_dimension.py b/opendm/dem/ground_rectification/extra_dimensions/distance_dimension.py index 33b56835..d2f72fdb 100755 --- a/opendm/dem/ground_rectification/extra_dimensions/distance_dimension.py +++ b/opendm/dem/ground_rectification/extra_dimensions/distance_dimension.py @@ -35,7 +35,7 @@ class DistanceDimension(Dimension): return 'distance_to_ground' def get_las_type(self): - return 10 + return 'float64' def __calculate_angle(self, model): "Calculate the angle between the estimated plane and the XY plane" diff --git a/opendm/dem/ground_rectification/extra_dimensions/extended_dimension.py b/opendm/dem/ground_rectification/extra_dimensions/extended_dimension.py index 741d041b..c371e83b 100755 --- a/opendm/dem/ground_rectification/extra_dimensions/extended_dimension.py +++ b/opendm/dem/ground_rectification/extra_dimensions/extended_dimension.py @@ -20,4 +20,4 @@ class ExtendedDimension(Dimension): return 'extended' def get_las_type(self): - return 3 + return 'uint16' diff --git a/opendm/dem/ground_rectification/extra_dimensions/partition_dimension.py b/opendm/dem/ground_rectification/extra_dimensions/partition_dimension.py index 8d39866b..b7fed6b6 100755 --- a/opendm/dem/ground_rectification/extra_dimensions/partition_dimension.py +++ b/opendm/dem/ground_rectification/extra_dimensions/partition_dimension.py @@ -22,4 +22,4 @@ class PartitionDimension(Dimension): return self.name def get_las_type(self): - return 5 + return 'uint32' diff --git a/opendm/dem/ground_rectification/io/las_io.py b/opendm/dem/ground_rectification/io/las_io.py index 49ce0361..07f50729 100755 --- a/opendm/dem/ground_rectification/io/las_io.py +++ b/opendm/dem/ground_rectification/io/las_io.py @@ -1,52 +1,39 @@ # TODO: Move to pylas when project migrates to python3 -from laspy.file import File -from laspy.header import Header +import laspy import numpy as np from ..point_cloud import PointCloud def read_cloud(point_cloud_path): # Open point cloud and read its properties - las_file = File(point_cloud_path, mode='r') - header = (las_file.header.copy(), las_file.header.scale, las_file.header.offset,las_file.header.evlrs, las_file.header.vlrs) - [x_scale, y_scale, z_scale] = las_file.header.scale - [x_offset, y_offset, z_offset] = las_file.header.offset + las_file = laspy.read(point_cloud_path) + header = las_file.header - # Calculate the real coordinates - x = las_file.X * x_scale + x_offset - y = las_file.Y * y_scale + y_offset - z = las_file.Z * z_scale + z_offset + x = las_file.x.scaled_array() + y = las_file.y.scaled_array() + z = las_file.z.scaled_array() - cloud = PointCloud.with_dimensions(x, y, z, las_file.Classification, las_file.red, las_file.green, las_file.blue) - - # Close the file - las_file.close() + cloud = PointCloud.with_dimensions(x, y, z, las_file.classification.array, las_file.red, las_file.green, las_file.blue) # Return the result return header, cloud def write_cloud(header, point_cloud, output_point_cloud_path, write_extra_dimensions=False): - (h, scale, offset, evlrs, vlrs) = header - # Open output file - output_las_file = File(output_point_cloud_path, mode='w', header=h, evlrs=evlrs, vlrs=vlrs) + output_las_file = laspy.LasData(header) if write_extra_dimensions: - # Create new dimensions - for name, dimension in point_cloud.extra_dimensions_metadata.items(): - output_las_file.define_new_dimension(name=name, data_type=dimension.get_las_type(), description="Dimension added by Ground Extend") - + extra_dims = [laspy.ExtraBytesParams(name=name, type=dimension.get_las_type(), description="Dimension added by Ground Extend") for name, dimension in point_cloud.extra_dimensions_metadata.items()] + output_las_file.add_extra_dims(extra_dims) # Assign dimension values for dimension_name, values in point_cloud.extra_dimensions.items(): setattr(output_las_file, dimension_name, values) # Adapt points to scale and offset - [x_scale, y_scale, z_scale] = scale - [x_offset, y_offset, z_offset] = offset [x, y] = np.hsplit(point_cloud.xy, 2) - output_las_file.X = (x.ravel() - x_offset) / x_scale - output_las_file.Y = (y.ravel() - y_offset) / y_scale - output_las_file.Z = (point_cloud.z - z_offset) / z_scale + output_las_file.x = x.ravel() + output_las_file.y = y.ravel() + output_las_file.z = point_cloud.z # Set color [red, green, blue] = np.hsplit(point_cloud.rgb, 3) @@ -55,11 +42,6 @@ def write_cloud(header, point_cloud, output_point_cloud_path, write_extra_dimens output_las_file.blue = blue.ravel() # Set classification - output_las_file.Classification = point_cloud.classification.astype(np.uint8) + output_las_file.classification = point_cloud.classification.astype(np.uint8) - # Set header - output_las_file.header.scale = scale - output_las_file.header.offset = offset - - # Close files - output_las_file.close() + output_las_file.write(output_point_cloud_path) \ No newline at end of file diff --git a/requirements.txt b/requirements.txt index c5441e9a..28d71347 100644 --- a/requirements.txt +++ b/requirements.txt @@ -7,7 +7,7 @@ ODMExifRead==3.0.4 Fiona==1.8.17 ; sys_platform == 'linux' or sys_platform == 'darwin' https://github.com/OpenDroneMap/windows-deps/raw/main/Fiona-1.8.19-cp38-cp38-win_amd64.whl ; sys_platform == 'win32' joblib==1.1.0 -laspy==1.7.0 +laspy[laszip]==2.3.0 lxml==4.6.1 matplotlib==3.3.3 networkx==2.5