add search radius plugin

.gitignore

@@ -1,4 +1,5 @@
 xcode
 install
 __pycache__
-test/__pycache__
+test/__pycache_
+test/.idea

.pre-commit-config.yaml

@@ -0,0 +1,16 @@
+repos:
+  - repo: https://github.com/ambv/black
+    rev: 23.12.0
+    hooks:
+      - id: black
+        language_version: python3.11
+  - repo: https://github.com/pycqa/flake8
+    rev: 6.1.0
+    hooks:
+    -   id: flake8
+  - repo: https://github.com/pycqa/isort
+    rev: 5.13.2
+    hooks:
+      - id: isort
+        name: isort (python)
+        args: ["--profile", "black"]

CMakeLists.txt

@@ -3,6 +3,7 @@ cmake_minimum_required( VERSION 3.5 )
 project(MY_READER LANGUAGES CXX)
 
 set(CMAKE_PREFIX_PATH ${CONDA_PREFIX})
+set(CMAKE_XCODE_ATTRIBUTE_OTHER_CODE_SIGN_FLAGS "-o linker-signed")
 
 find_package(PDAL REQUIRED)
 
@@ -11,4 +12,5 @@ set(CMAKE_DEBUG_POSTFIX d)
 
 ## add plugin
 add_subdirectory(src/filter_grid_decimation)
+add_subdirectory(src/filter_radius_assign)
 

Dockerfile

@@ -1,6 +1,6 @@
 FROM mambaorg/micromamba:bullseye-slim as build
 
-COPY environment.yml /environment_docker.yml
+COPY environment_docker.yml /environment_docker.yml
 
 USER root
 RUN micromamba env create -f /environment_docker.yml
@@ -9,15 +9,19 @@ RUN apt-get update && apt-get install --no-install-recommends -y cmake make buil
 
 COPY src src
 COPY CMakeLists.txt CMakeLists.txt
+COPY macro macro
 
-RUN cmake -G"Unix Makefiles" -DCONDA_PREFIX=$CONDA_PREFIX -DCMAKE_BUILD_TYPE=Release  
+RUN cmake -G"Unix Makefiles" -DCONDA_PREFIX=$CONDA_PREFIX -DCMAKE_BUILD_TYPE=Release
 RUN make -j4 install
 
 FROM debian:bullseye-slim
 
-COPY --from=build /opt/conda/envs/pdal_ign_plugin /opt/conda/envs/pdal_ign_plugin
-COPY --from=build /tmp/install/lib /tmp/install/lib
-
+COPY --from=build /opt/conda/envs/pdal_ign_plugin /opt/conda/envs/pdal_ign_plugin
+RUN mkdir -p /pdal_ign_plugin
+COPY --from=build /tmp/install/lib /pdal_ign_plugin/install/lib
+COPY --from=build /tmp/macro /macro 
+
 ENV PATH=$PATH:/opt/conda/envs/pdal_ign_plugin/bin/
 ENV PROJ_LIB=/opt/conda/envs/pdal_ign_plugin/share/proj/
-ENV PDAL_DRIVER_PATH=/tmp/install/lib
+ENV PDAL_DRIVER_PATH=/pdal_ign_plugin/install/lib
+

README.md

@@ -45,6 +45,8 @@ python -m pytest -s
 
 [grid decimation](./doc/grid_decimation.md)
 
+[radius assign](./doc/radius_assign.md)
+
 ## Adding a filter
 
 In order to add a filter, you have to add a new folder in the src directory : 

ci/build.sh

@@ -24,6 +24,7 @@ fi
 conda activate pdal_ign_plugin
 
 export CONDA_PREFIX=$CONDA_PREFIX
+echo conda is $CONDA_PREFIX 
 
 mkdir build
 cd build

doc/grid_decimation.md

@@ -1,9 +1,11 @@
 # filter grid decimation
 
+**Deprecated** : *better use the gridDecimation filter of pdal > 2.7*
+
 Purpose
 ---------------------------------------------------------------------------------------------------------
 
-The **grid decimation filter** transform only one point in each cells of a grid calculated from the points cloud and a resolution therm. The transformation is done by the value information. The selected point could be the highest or the lowest point on the cell. It can be used, for exemple, to quickly filter vegetation points in order to keep only the canopy points. A new attribut is created with the value '1' for the grid, and '0' for the other points.
+The **grid decimation filter** transform only one point in each cells of a grid calculated from the points cloud and a resolution therm. The transformation is done by the value information. The selected point could be the highest or the lowest point on the cell. It can be used, for exemple, to quickly filter vegetation points in order to keep only the canopy points. A new dimension is created with the value '1' for the grid, and '0' for the other points.
 
 
 Example
@@ -18,7 +20,7 @@ This example transform highest points of classification 5 in classification 9, o
     {
         "type": "filters.gridDecimation",
 	"output_type":"max",
-        "output_name_attribut": "grid",
+        "output_dimension": "grid",
 	"output_wkt":"file-output.wkt"
     },
     {
@@ -37,6 +39,6 @@ Options
 **resolution** :
   The resolution of the cells in meter. [Default: 1.]
 
-**output_name_attribut**: The name of the new attribut. [Default: grid]
+**output_dimension**: The name of the new dimension. [Default: grid]
 
 **output_wkt**: the name of the export grid file as wkt polygon. If none, no export [Default:""]

doc/radius_assign.md

@@ -0,0 +1,51 @@
+# filter radius assign
+
+Purpose
+---------------------------------------------------------------------------------------------------------
+
+The **radius assign filter**  overwrites the output_dimension_ dimension with boolean values:
+* 1 if the point has any neighbor with a distance lower than radius_ that belongs to the domain reference_domain_
+* 0 otherwise.
+
+
+Example
+---------------------------------------------------------------------------------------------------------
+
+This pipeline updates the Keypoint dimension of all points with classification 1 to 2 (unclassified and ground) that are closer than 1 meter from a point with classification 6 (building)
+
+
+```
+  [
+     "file-input.las",
+      {
+          "type" : "filters.radius_assign",
+          "src_domain" : "Classification[1:2]",
+          "reference_domain" : "Classification[6:6]",
+          "radius" : 1,
+	  	  "output_dimension": "radius",
+	  	  "is3d": True
+      },
+      "output.las"
+  ]
+```
+
+Options
+---------------------------------------------------------------------------------------------------------------------------------------------------------------------
+
+**src_domain** :
+  A :ref:`range <ranges>` which selects points to be processed by the filter. Can be specified multiple times.  Points satisfying any range will be processed
+
+**reference_domain** :
+  A :ref:`range <ranges>` which selects points that can are considered as potential neighbors. Can be specified multiple times.
+
+**radius** :
+  An positive float which specifies the radius for the neighbors search.
+
+**output_dimension**: The name of the new dimension'. [Default: radius]
+
+**is3d**: Search in 3d (as a ball). [Default: false]
+
+**max2d_above**: If search in 2d : upward maximum distance in Z for potential neighbors (corresponds to a search in a cylinder with a height = max2d_above above the source point). Default (0) = infinite height [Default: 0.]
+
+**max2d_below**: If search in 2d : upward maximum distance in Z for potential neighbors (corresponds to a search in a cylinder with a height = max2d_below below the source point). Default (0) = infinite height [Default: 0.]
+

environment.yml

@@ -1,14 +1,19 @@
 name: pdal_ign_plugin
 channels:
   - conda-forge
+  - anaconda
 dependencies:
   - pdal
   - python-pdal
+  - gdal
+# --------- dev dep --------- #
+  - cmake
+  - pre-commit # hooks for applying linters on commit
+  - black # code formatting
+  - isort # import sorting
+  - flake8 # code analysis
   - pytest
-  - black
-  - isort
-  - shapely
+# --------- pip & pip librairies --------- #
+  - pip
   - pip:
     - ign-pdal-tools
-
-

environment_docker.yml

@@ -1,6 +1,13 @@
 name: pdal_ign_plugin
 channels:
   - conda-forge
+  - anaconda
 dependencies:
   - pdal
+  - python-pdal
+  - gdal
+# --------- pip & pip librairies --------- #
+  - pip
+  - pip:
+    - ign-pdal-tools
 

macro/ex_filtering_points.py

@@ -0,0 +1,94 @@
+import argparse
+import pdal
+import macro
+
+"""
+This tool shows how to use functions of macro in a pdal pipeline
+"""
+
+def parse_args():
+    parser = argparse.ArgumentParser("Tool to apply pdal pipelines for DSM and DTM calculation")
+    parser.add_argument("--input", "-i", type=str, required=True, help="Input las file")
+    parser.add_argument("--output_las", "-o", type=str, required=True, help="Output cloud las file")
+    parser.add_argument("--output_dsm", "-s", type=str, required=True, help="Output dsm tiff file")
+    parser.add_argument("--output_dtm", "-t", type=str, required=True, help="Output dtm tiff file")
+    return parser.parse_args()
+
+
+if __name__ == "__main__":
+    args = parse_args()
+
+    pipeline = pdal.Reader.las(args.input)
+
+    ## 1 - recherche des points max de végétation (4,5) sur une grille régulière, avec prise en compte des points sol (2) et basse
+    ##     vegetation (3) proche de la végétation : on les affecte en 100
+
+    # bouche trou : assigne les points sol en 102 à l'intérieur de la veget (4,5)
+    pipeline = macro.add_radius_assign(pipeline, 1, False, condition_src="Classification==2", condition_ref=macro.build_condition("Classification", [4,5]), condition_out="Classification=102")
+    pipeline = macro.add_radius_assign(pipeline, 1, False, condition_src="Classification==102", condition_ref="Classification==2", condition_out="Classification=2")
+
+    # selection des points de veget basse proche de la veget haute : assigne 103
+    pipeline = macro.add_radius_assign(pipeline, 1, False, condition_src="Classification==3", condition_ref="Classification==5", condition_out="Classification=103")
+
+    # max des points de veget (et surement veget - 102,103) sur une grille régulière : assigne 100
+    pipeline |= pdal.Filter.gridDecimation(resolution=0.75, value="Classification=100", output_type="max", where=macro.build_condition("Classification", [4,5,102,103]))
+
+    # remise à zero des codes 102 et 103
+    pipeline |= pdal.Filter.assign(value="Classification=2", where="Classification==102")
+    pipeline |= pdal.Filter.assign(value="Classification=3", where="Classification==103")
+
+    ## 2 - sélection des points pour DTM et DSM
+
+    # selection de points sol (max) sur une grille régulière
+    pipeline |= pdal.Filter.gridDecimation(resolution=0.5, value="Classification=102", output_type="max", where="Classification==2")
+
+    # selection de points DSM (max) sur une grille régulière
+    pipeline |= pdal.Filter.gridDecimation(resolution=0.5, value="Classification=200", output_type="max", where=macro.build_condition("Classification", [2,3,4,5,6,9,17,64,100]))
+
+    # assigne des points sol sélectionnés (102) en 100 : les points proches de la végaétation, des ponts, de l'eau et 64
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==102",
+                                       condition_ref=macro.build_condition("Classification", [4,5,6,9,17,64,100]), condition_out="Classification=100")
+
+    # remise à zero du code 102
+    pipeline |= pdal.Filter.assign(value="Classification=2", where="Classification==102")
+
+    ## 3 - gestion des ponts
+
+
+
+    # bouche trou : on élimine les points sol (2) au milieu du pont en les mettant à 102
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==2", condition_ref="Classification==17", condition_out="Classification=102")
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==102",
+                                       condition_ref=macro.build_condition("Classification", [2,3,4,5]), condition_out="Classification=2")
+
+    # bouche trou : on élimine les points basse végétation (3) au milieu du pont en les mettant à 103
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==3", condition_ref="Classification==17", condition_out="Classification=103")
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==103",
+                                       condition_ref=macro.build_condition("Classification", [2,3,4,5]), condition_out="Classification=3")
+
+    # bouche trou : on élimine les points moyenne végétation (4) au milieu du pont en les mettant à 104
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==4", condition_ref="Classification==17", condition_out="Classification=104")
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==104",
+                                       condition_ref=macro.build_condition("Classification", [2,3,4,5]), condition_out="Classification=4")
+
+    # bouche trou : on élimine les points haute végétation (5) au milieu du pont en les mettant à 105
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==5", condition_ref="Classification==17", condition_out="Classification=105")
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==105",
+                                       condition_ref=macro.build_condition("Classification", [2,3,4,5]), condition_out="Classification=5")
+
+    # bouche trou : on élimine les points eau (9) au milieu du pont en les mettant à 109
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==9", condition_ref="Classification==17", condition_out="Classification=109")
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="Classification==109",
+                                       condition_ref="Classification==9", condition_out="Classification=9")
+
+    # step 15 et supression des points ??
+
+    # 4 - export du nuage
+    pipeline |= pdal.Writer.las(extra_dims="all",forward="all",filename=args.output_las)
+
+    # export des DSM/DTM
+    pipeline |= pdal.Writer.gdal(gdaldriver="GTiff", output_type="max", resolution=2.0, filename=args.output_dtm, where=macro.build_condition("Classification", [2,66]))
+    pipeline |= pdal.Writer.gdal(gdaldriver="GTiff", output_type="max", resolution=2.0, filename=args.output_dsm, where=macro.build_condition("Classification", [2,3,4,5,17,64]))
+
+    pipeline.execute()
+

macro/ex_filtering_points_with_add_dimensions.py

@@ -0,0 +1,78 @@
+import argparse
+import pdal
+import macro
+
+"""
+This tool shows how to use functions of macro in a pdal pipeline
+"""
+
+def parse_args():
+    parser = argparse.ArgumentParser("Tool to apply pdal pipelines for DSM and DTM calculation (with add dimensions for the concerned points)")
+    parser.add_argument("--input", "-i", type=str, required=True, help="Input las file")
+    parser.add_argument("--output_las", "-o", type=str, required=True, help="Output cloud las file")
+    parser.add_argument("--output_dsm", "-s", type=str, required=True, help="Output dsm tiff file")
+    parser.add_argument("--output_dtm", "-t", type=str, required=True, help="Output dtm tiff file")
+    return parser.parse_args()
+
+
+if __name__ == "__main__":
+    args = parse_args()
+
+    pipeline = pdal.Reader.las(args.input)
+
+    # 0 - ajout de dimensions temporaires
+    pipeline |= pdal.Filter.ferry(dimensions=f"=>PT_GRID_DSM, =>PT_VEG_DSM, =>PT_GRID_DTM, =>PT_ON_BRIDGE")
+
+
+    ## 1 - recherche des points max de végétation (4,5) sur une grille régulière, avec prise en compte des points sol (2) et basse
+    ##     vegetation (3) proche de la végétation
+    ##     pour le calcul du DSM
+
+    pipeline |= pdal.Filter.assign(value=["PT_VEG_DSM = 1 WHERE " + macro.build_condition("Classification", [4,5])])
+
+    # bouche trou : assigne les points sol à l'intérieur de la veget (4,5)
+    pipeline = macro.add_radius_assign(pipeline, 1, False, condition_src="Classification==2", condition_ref=macro.build_condition("Classification", [4,5]), condition_out="PT_VEG_DSM=1")
+    pipeline = macro.add_radius_assign(pipeline, 1, False, condition_src="PT_VEG_DSM==1 && Classification==2", condition_ref="Classification==2", condition_out="PT_VEG_DSM=0")
+
+    # selection des points de veget basse proche de la veget haute
+    pipeline = macro.add_radius_assign(pipeline, 1, False, condition_src="Classification==3", condition_ref="Classification==5", condition_out="PT_VEG_DSM=1")
+
+    # max des points de veget (PT_VEG_DSM==1) sur une grille régulière :
+    pipeline |= pdal.Filter.gridDecimation(resolution=0.75, value="PT_GRID_DSM=1", output_type="max", where="PT_VEG_DSM==1")
+
+
+    ## 2 - sélection des points pour DTM et DSM
+
+    # selection de points DTM (max) sur une grille régulière
+    pipeline |= pdal.Filter.gridDecimation(resolution=0.5, value="PT_GRID_DTM=1", output_type="max", where="Classification==2")
+
+    # selection de points DSM (max) sur une grille régulière
+    pipeline |= pdal.Filter.gridDecimation(resolution=0.5, value="PT_GRID_DSM=1", output_type="max",
+                                           where="(" + macro.build_condition("Classification", [6,9,17,64]) + ") || PT_GRID_DSM==1")
+
+    # assigne des points sol sélectionnés : les points proches de la végétation, des ponts, de l'eau, 64
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="PT_GRID_DTM==1",
+                                       condition_ref=macro.build_condition("Classification", [4,5,6,9,17,64]),
+                                       condition_out="PT_GRID_DSM=1")
+
+
+    ## 3 - gestion des ponts
+    #      bouche trou : on filtre les points (2,3,4,5,9) au milieu du pont en les mettant à PT_ON_BRIDGE=1
+
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src=macro.build_condition("Classification", [2,3,4,5,9]), condition_ref="Classification==17", condition_out="PT_ON_BRIDGE=1")
+    pipeline = macro.add_radius_assign(pipeline, 1.5, False, condition_src="PT_ON_BRIDGE==1",
+                                       condition_ref=macro.build_condition("Classification", [2,3,4,5]), condition_out="PT_ON_BRIDGE=0")
+    pipeline |= pdal.Filter.assign(value=["PT_GRID_DSM=0 WHERE PT_ON_BRIDGE==1"])
+
+
+    ## 4 - point pour DTM servent au DSM également
+    pipeline |= pdal.Filter.assign(value=["PT_GRID_DSM=1 WHERE PT_GRID_DTM==1"])
+
+    ## 5 - export du nuage et des DSM
+
+    pipeline |= pdal.Writer.las(extra_dims="all", forward="all", filename=args.output_las)
+    pipeline |= pdal.Writer.gdal(gdaldriver="GTiff", output_type="max", resolution=2.0, filename=args.output_dtm, where="PT_GRID_DTM==1")
+    pipeline |= pdal.Writer.gdal(gdaldriver="GTiff", output_type="max", resolution=2.0, filename=args.output_dsm, where="PT_GRID_DSM==1")
+
+    pipeline.execute()
+