diff --git a/.gitignore b/.gitignore index 8e94d97..7a1eb5c 100644 --- a/.gitignore +++ b/.gitignore @@ -9,3 +9,6 @@ __pycache__ *.egg-info *.idea + +# test output +tmp diff --git a/README.md b/README.md index 8b1e758..5bf9a35 100755 --- a/README.md +++ b/README.md @@ -31,9 +31,8 @@ The code is structured as: ├── macro # Python module with ready-to-use filters combinations │   ├── __init__.py │   ├── macro.py -│   └── version.py -├── scripts -│   ├── *.py # Example scripts to use the plugin filters + the filters combinations contained in `macro` +│   ├── version.py +│   └── *.py # Example scripts to use the plugin filters + the filters combinations contained in `macro` ├── test ├── CMakeLists.txt ├── environment*.yml diff --git a/environment.yml b/environment.yml index 64da390..a15c196 100755 --- a/environment.yml +++ b/environment.yml @@ -18,4 +18,4 @@ dependencies: # --------- pip & pip libraries --------- # - pip - pip: - - ign-pdal-tools + - ign-pdal-tools==1.7.1 diff --git a/environment_docker.yml b/environment_docker.yml index 6d00e85..4dab5e9 100755 --- a/environment_docker.yml +++ b/environment_docker.yml @@ -11,5 +11,5 @@ dependencies: # --------- pip & pip libraries --------- # - pip - pip: - - ign-pdal-tools - + - ign-pdal-tools==1.7.1 + diff --git a/pdal_ign_macro/mark_points_to_use_for_digital_models_with_new_dimension.py b/pdal_ign_macro/mark_points_to_use_for_digital_models_with_new_dimension.py index e455489..09f8e3f 100755 --- a/pdal_ign_macro/mark_points_to_use_for_digital_models_with_new_dimension.py +++ b/pdal_ign_macro/mark_points_to_use_for_digital_models_with_new_dimension.py @@ -1,6 +1,10 @@ import argparse +import shutil +import tempfile import pdal +from pdaltools.las_add_buffer import run_on_buffered_las +from pdaltools.las_remove_dimensions import remove_dimensions_from_las from pdal_ign_macro import macro @@ -39,23 +43,55 @@ def parse_args(): parser.add_argument( "--output_dtm", "-t", type=str, required=False, default="", help="Output dtm tiff file" ) + parser.add_argument( + "--keep_temporary_dims", + "-k", + action="store_true", + help="If set, do not delete temporary dimensions", + ) + parser.add_argument( + "--skip_buffer", + "-s", + action="store_true", + help="If set, skip adding a buffer from the neighbor tiles based on their name", + ) + parser.add_argument( + "--buffer_width", + type=float, + default=25, + help="width of the border to add to the tile (in meters)", + ) + parser.add_argument( + "--spatial_ref", + type=str, + default="EPSG:2154", + help="spatial reference for the writer (required when running with a buffer)", + ) + parser.add_argument( + "--tile_width", + type=int, + default=1000, + action="store_true", + help="width of tiles in meters (required when running with a buffer)", + ) + parser.add_argument( + "--tile_coord_scale", + type=int, + default=1000, + action="store_true", + help="scale used in the filename to describe coordinates in meters (required when running with a buffer)", + ) + return parser.parse_args() -def mark_points_to_use_for_digital_models_with_new_dimension( - input_las, output_las, dsm_dimension, dtm_dimension, output_dsm, output_dtm -): +def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension): pipeline = pdal.Pipeline() | pdal.Reader.las(input_las) # 0 - ajout de dimensions temporaires et de sortie - added_dimensions = [ - dtm_dimension, - dsm_dimension, - "PT_VEG_DSM", - "PT_ON_BRIDGE", - "PT_ON_BUILDING", - "PT_ON_VEGET", - ] + temporary_dimensions = ["PT_VEG_DSM", "PT_ON_BRIDGE", "PT_ON_BUILDING", "PT_ON_VEGET"] + added_dimensions = [dtm_dimension, dsm_dimension] + temporary_dimensions + pipeline |= pdal.Filter.ferry(dimensions="=>" + ", =>".join(added_dimensions)) # 1 - recherche des points max de végétation (4,5) sur une grille régulière, avec prise en @@ -220,31 +256,99 @@ def mark_points_to_use_for_digital_models_with_new_dimension( ) # ERREUR EN 4!###############################################################################################! # 5 - export du nuage et des DSM - # TODO: n'ajouter que les dimensions de sortie utiles ! pipeline |= pdal.Writer.las(extra_dims="all", forward="all", filename=output_las) - if output_dtm: - pipeline |= pdal.Writer.gdal( - gdaldriver="GTiff", - output_type="max", - resolution=0.5, - filename=output_dtm, - where=f"{dtm_dimension}==1", + return pipeline, temporary_dimensions + + +def mark_points_to_use_for_digital_models_with_new_dimension( + input_las, + output_las, + dsm_dimension, + dtm_dimension, + output_dsm, + output_dtm, + keep_temporary_dimensions=False, +): + with tempfile.NamedTemporaryFile(suffix="_with_temporary_dims.las", dir=".") as tmp_las: + pipeline, temporary_dimensions = define_marking_pipeline( + input_las, + tmp_las.name, + dsm_dimension, + dtm_dimension, ) - if output_dsm: - pipeline |= pdal.Writer.gdal( - gdaldriver="GTiff", - output_type="max", - resolution=0.5, - filename=output_dsm, - where=f"{dsm_dimension}==1", + if output_dtm: + pipeline |= pdal.Writer.gdal( + gdaldriver="GTiff", + output_type="max", + resolution=0.5, + filename=output_dtm, + where=f"{dtm_dimension}==1", + ) + + if output_dsm: + pipeline |= pdal.Writer.gdal( + gdaldriver="GTiff", + output_type="max", + resolution=0.5, + filename=output_dsm, + where=f"{dsm_dimension}==1", + ) + + pipeline.execute() + + if keep_temporary_dimensions: + shutil.copy(tmp_las.name, output_las) + else: + remove_dimensions_from_las( + tmp_las.name, + temporary_dimensions + ["SRC_DOMAIN", "REF_DOMAIN", "radius_search"], + output_las, + ) + + +def main( + input_las, + output_las, + dsm_dimension, + dtm_dimension, + output_dsm, + output_dtm, + keep_temporary_dimensions=False, + skip_buffer=False, + buffer_width=25, + spatial_ref="EPSG:2154", + tile_width=1000, + tile_coord_scale=1000, +): + if skip_buffer: + mark_points_to_use_for_digital_models_with_new_dimension( + input_las, + output_las, + dsm_dimension, + dtm_dimension, + output_dsm, + output_dtm, + keep_temporary_dimensions, ) + else: + mark_with_buffer = run_on_buffered_las( + buffer_width, spatial_ref, tile_width, tile_coord_scale + )(mark_points_to_use_for_digital_models_with_new_dimension) - pipeline.execute() + mark_with_buffer( + input_las, + output_las, + dsm_dimension, + dtm_dimension, + output_dsm, + output_dtm, + keep_temporary_dimensions, + ) if __name__ == "__main__": args = parse_args() - mark_points_to_use_for_digital_models_with_new_dimension(**vars(args)) + main(**vars(args)) diff --git a/test/data/buffer/test_data_77055_627755_LA93_IGN69.laz b/test/data/buffer/test_data_77055_627755_LA93_IGN69.laz new file mode 100644 index 0000000..6d7cc9c Binary files /dev/null and b/test/data/buffer/test_data_77055_627755_LA93_IGN69.laz differ diff --git a/test/data/buffer/test_data_77055_627760_LA93_IGN69.laz b/test/data/buffer/test_data_77055_627760_LA93_IGN69.laz new file mode 100644 index 0000000..13129bd Binary files /dev/null and b/test/data/buffer/test_data_77055_627760_LA93_IGN69.laz differ diff --git a/test/data/buffer/test_data_77060_627755_LA93_IGN69.laz b/test/data/buffer/test_data_77060_627755_LA93_IGN69.laz new file mode 100644 index 0000000..70e61b9 Binary files /dev/null and b/test/data/buffer/test_data_77060_627755_LA93_IGN69.laz differ diff --git a/test/pdal_ign_macro/test_mark_points_to_use_for_digital_models_with_new_dimension.py b/test/pdal_ign_macro/test_mark_points_to_use_for_digital_models_with_new_dimension.py index d66a811..f551d01 100644 --- a/test/pdal_ign_macro/test_mark_points_to_use_for_digital_models_with_new_dimension.py +++ b/test/pdal_ign_macro/test_mark_points_to_use_for_digital_models_with_new_dimension.py @@ -4,11 +4,12 @@ import pdal from pdal_ign_macro.mark_points_to_use_for_digital_models_with_new_dimension import ( + main, mark_points_to_use_for_digital_models_with_new_dimension, ) -def test_main(): +def test_mark_points_to_use_for_digital_models_with_new_dimension(): ini_las = "test/data/4_6.las" dsm_dimension = "dsm_marker" dtm_dimension = "dtm_marker" @@ -17,9 +18,102 @@ def test_main(): ini_las, las_output.name, dsm_dimension, dtm_dimension, "", "" ) pipeline = pdal.Pipeline() + pipeline |= pdal.Reader.las(ini_las) + input_dimensions = set(pipeline.quickinfo["readers.las"]["dimensions"].split(", ")) + pipeline = pdal.Pipeline() + pipeline |= pdal.Reader.las(las_output.name) + output_dimensions = set(pipeline.quickinfo["readers.las"]["dimensions"].split(", ")) + assert output_dimensions == input_dimensions.union([dsm_dimension, dtm_dimension]) + + pipeline.execute() + arr = pipeline.arrays[0] + assert np.any(arr[dsm_dimension] == 1) + assert np.any(arr[dtm_dimension] == 1) + + +def test_mark_points_to_use_for_digital_models_with_new_dimension_keep_dimensions(): + ini_las = "test/data/4_6.las" + dsm_dimension = "dsm_marker" + dtm_dimension = "dtm_marker" + with tempfile.NamedTemporaryFile(suffix="_mark_points_output.las") as las_output: + mark_points_to_use_for_digital_models_with_new_dimension( + ini_las, + las_output.name, + dsm_dimension, + dtm_dimension, + "", + "", + keep_temporary_dimensions=True, + ) + pipeline = pdal.Pipeline() + pipeline |= pdal.Reader.las(las_output.name) + output_dimensions = set(pipeline.quickinfo["readers.las"]["dimensions"].split(", ")) + assert dsm_dimension in output_dimensions + assert dtm_dimension in output_dimensions + + assert all( + [ + dim in output_dimensions + for dim in ["PT_VEG_DSM", "PT_ON_BRIDGE", "PT_ON_BUILDING", "PT_ON_VEGET"] + ] + ) + + pipeline.execute() + arr = pipeline.arrays[0] + assert np.any(arr[dsm_dimension] == 1) + assert np.any(arr[dtm_dimension] == 1) + + +def test_main_no_buffer(): + ini_las = "test/data/4_6.las" + dsm_dimension = "dsm_marker" + dtm_dimension = "dtm_marker" + with tempfile.NamedTemporaryFile(suffix="_mark_points_output.las") as las_output: + main( + ini_las, + las_output.name, + dsm_dimension, + dtm_dimension, + "", + "", + keep_temporary_dimensions=False, + skip_buffer=True, + ) + pipeline = pdal.Pipeline() + pipeline |= pdal.Reader.las(las_output.name) + output_dimensions = pipeline.quickinfo["readers.las"]["dimensions"].split(", ") + assert dsm_dimension in output_dimensions + assert dtm_dimension in output_dimensions + + pipeline.execute() + arr = pipeline.arrays[0] + assert np.any(arr[dsm_dimension] == 1) + assert np.any(arr[dtm_dimension] == 1) + + +def test_main_with_buffer(): + ini_las = "test/data/buffer/test_data_77055_627755_LA93_IGN69.laz" + dsm_dimension = "dsm_marker" + dtm_dimension = "dtm_marker" + with tempfile.NamedTemporaryFile(suffix="_mark_points_output.las") as las_output: + main( + ini_las, + las_output.name, + dsm_dimension, + dtm_dimension, + "", + "", + keep_temporary_dimensions=False, + skip_buffer=False, + buffer_width=10, + tile_width=50, + tile_coord_scale=10, + ) + pipeline = pdal.Pipeline() pipeline |= pdal.Reader.las(las_output.name) - assert dsm_dimension in pipeline.quickinfo["readers.las"]["dimensions"].split(", ") - assert dtm_dimension in pipeline.quickinfo["readers.las"]["dimensions"].split(", ") + output_dimensions = pipeline.quickinfo["readers.las"]["dimensions"].split(", ") + assert dsm_dimension in output_dimensions + assert dtm_dimension in output_dimensions pipeline.execute() arr = pipeline.arrays[0]