diff --git a/CMakeLists.txt b/CMakeLists.txt index 0f8e54f69f..ce6de98b99 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -57,7 +57,7 @@ else() endif() set(PDAL_UTIL_LIB_NAME pdal_util) set(PDAL_BOOST_LIB_NAME pdal_boost) -#set(PDAL_ARBITER_LIB_NAME pdal_arbiter) +set(PDAL_KAZHDAN_LIB_NAME pdal_kazhdan) set(CMAKE_INCLUDE_DIRECTORIES_PROJECT_BEFORE ON) @@ -172,6 +172,7 @@ endif() add_subdirectory(dimbuilder) add_subdirectory(vendor/pdalboost) add_subdirectory(vendor/arbiter) +add_subdirectory(vendor/kazhdan) if (NOT PDAL_HAVE_JSONCPP) add_subdirectory(vendor/jsoncpp/dist) @@ -247,6 +248,7 @@ target_link_libraries(${PDAL_BASE_LIB_NAME} ${PDAL_REEXPORT} ${PDAL_UTIL_LIB_NAME} ${PDAL_ARBITER_LIB_NAME} + ${PDAL_KAZHDAN_LIB_NAME} ${JSON_CPP_LINK_TYPE} ${PDAL_JSONCPP_LIB_NAME} INTERFACE @@ -285,6 +287,11 @@ install(DIRECTORY ${PDAL_KERNELS_DIR} FILES_MATCHING PATTERN "*.hpp" PATTERN "private" EXCLUDE ) +install(DIRECTORY ${PDAL_VENDOR_DIR}/nanoflann + DESTINATION include/pdal + FILES_MATCHING PATTERN "*.hpp" + PATTERN "private" EXCLUDE +) install(DIRECTORY ${PDAL_IO_DIR} DESTINATION include/pdal FILES_MATCHING PATTERN "*.hpp" diff --git a/cmake/modules/FindLASzip.cmake b/cmake/modules/FindLASzip.cmake index 320984a11c..08683256aa 100644 --- a/cmake/modules/FindLASzip.cmake +++ b/cmake/modules/FindLASzip.cmake @@ -54,6 +54,46 @@ FIND_LIBRARY(LASZIP_LIBRARY /usr/local/lib ${OSGEO4W_ROOT_DIR}/lib) +# Comment out laszip.hpp version info +#[[ +SET(LASZIP_VERSION_H "${LASZIP_INCLUDE_DIR}/laszip.hpp") +IF(LASZIP_INCLUDE_DIR AND EXISTS ${LASZIP_VERSION_H}) + SET(LASZIP_VERSION 0) + + message("Looking for: ${LASZIP_INCLUDE_DIR}/laszip.hpp") + + FILE(READ ${LASZIP_VERSION_H} LASZIP_VERSION_H_CONTENTS) + + IF (DEFINED LASZIP_VERSION_H_CONTENTS) + string(REGEX REPLACE ".*#define[ \t]LASZIP_VERSION_MAJOR[ \t]+([0-9]+).*" "\\1" LASZIP_VERSION_MAJOR "${LASZIP_VERSION_H_CONTENTS}") + string(REGEX REPLACE ".*#define[ \t]LASZIP_VERSION_MINOR[ \t]+([0-9]+).*" "\\1" LASZIP_VERSION_MINOR "${LASZIP_VERSION_H_CONTENTS}") + string(REGEX REPLACE ".*#define[ \t]LASZIP_VERSION_REVISION[ \t]+([0-9]+).*" "\\1" LASZIP_VERSION_REVISION "${LASZIP_VERSION_H_CONTENTS}") + + if(NOT ${LASZIP_VERSION_MAJOR} MATCHES "[0-9]+") + message(FATAL_ERROR "LASzip version parsing failed for LASZIP_VERSION_MAJOR!") + endif() + if(NOT ${LASZIP_VERSION_MINOR} MATCHES "[0-9]+") + message(FATAL_ERROR "LASzip version parsing failed for LASZIP_VERSION_MINOR!") + endif() + if(NOT ${LASZIP_VERSION_REVISION} MATCHES "[0-9]+") + message(FATAL_ERROR "LASzip version parsing failed for LASZIP_VERSION_REVISION!") + endif() + + + SET(LASZIP_VERSION "${LASZIP_VERSION_MAJOR}.${LASZIP_VERSION_MINOR}.${LASZIP_VERSION_REVISION}" + CACHE INTERNAL "The version string for LASzip library") + + IF (LASZIP_VERSION VERSION_LESS LASzip_FIND_VERSION) + MESSAGE(FATAL_ERROR "LASzip version check failed. Version ${LASZIP_VERSION} was found, at least version ${LASzip_FIND_VERSION} is required") + ENDIF() + ELSE() + MESSAGE(FATAL_ERROR "Failed to open ${LASZIP_VERSION_H} file") + ENDIF() +ELSE() + return() +ENDIF() +]] + # Handle the QUIETLY and REQUIRED arguments and set LASZIP_FOUND to TRUE # if all listed variables are TRUE INCLUDE(FindPackageHandleStandardArgs) diff --git a/doc/stages/filters.crop.rst b/doc/stages/filters.crop.rst index fbe9a10f62..0dab0b4571 100644 --- a/doc/stages/filters.crop.rst +++ b/doc/stages/filters.crop.rst @@ -8,6 +8,14 @@ box (2D), polygon, or point+distance. If more than one bounding region is specified, the filter will pass all input points through each bounding region, creating an output point set for each input crop region. +The provided bounding regions are assumed to have the same spatial reference +as the points unless the option `a_srs` provides an explicit spatial reference +for bounding regions. +If the point input consists of multiple point views with differing +spatial references, one is chosen at random and assumed to be the +spatial reference of the input bounding region. In this case a warning will +be logged. + Example ------- @@ -46,3 +54,8 @@ point distance Distance in units of common X, Y, and Z :ref:`dimensions` to crop circle or sphere in combination with ``point``. + +a_srs + Indicates the spatial reference of the bounding regions. If not provided, + it is assumed that the spatial reference of the bounding region matches + that of the points (if possible). diff --git a/doc/stages/filters.pmf.rst b/doc/stages/filters.pmf.rst index 6046843303..d49f97e14c 100644 --- a/doc/stages/filters.pmf.rst +++ b/doc/stages/filters.pmf.rst @@ -52,6 +52,12 @@ Notes extra iteration. This parameter can have a strongly negative impact on computation performance. +* ``exponential`` is used to control the rate of growth of morphological window + sizes toward ``max_window_size``. Linear growth preserves gradually changing + topographic features well, but demands considerable compute time. The default + behavior is to grow the window sizes exponentially, thus reducing the number + of iterations. + * This filter will mark all returns deemed to be ground returns with a classification value of 2 (per the LAS specification). To extract only these returns, users can add a :ref:`range filter` to the pipeline. @@ -71,27 +77,27 @@ Notes Options ------------------------------------------------------------------------------- -max_window_size - Maximum window size. [Default: **33**] - -slope - Slope. [Default: **1.0**] - -max_distance - Maximum distance. [Default: **2.5**] - -initial_distance - Initial distance. [Default: **0.15**] - cell_size Cell Size. [Default: **1**] -approximate - Use approximate algorithm? [Default: **false**] - +exponential + Use exponential growth for window sizes? [Defualt: **true**] + ignore Optional range of values to ignore. + +initial_distance + Initial distance. [Default: **0.15**] last Consider only last returns (when return information is available)? [Default: **true**] + +max_distance + Maximum distance. [Default: **2.5**] + +max_window_size + Maximum window size. [Default: **33**] + +slope + Slope. [Default: **1.0**] diff --git a/doc/stages/filters.poisson.rst b/doc/stages/filters.poisson.rst index fa99fa2f74..1735bfa977 100644 --- a/doc/stages/filters.poisson.rst +++ b/doc/stages/filters.poisson.rst @@ -4,14 +4,26 @@ filters.poisson =============================================================================== -The Poisson filter passes data through the Point Cloud Library (`PCL`_) Poisson -surface reconstruction algorithm. - -Poisson is an implementation of the method described in [Kazhdan2006]_. +The poisson filter passes data Mischa Kazhdan's poisson surface reconstruction +algorithm. [Kazhdan2006]_ It creates a watertight surface from the original +point set by creating an entirely new point set representing the imputed +isosurface. The algorithm requires normal vectors to each point in order +to run. If the x, y and z normal dimensions are present in the input point +set, they will be used by the algorithm. If they don't exist, the poisson +filter will invoke the PDAL normal filter to create them before running. + +The poisson algorithm will usually create a larger output point set +than the input point set. Because the algorithm constructs new points, data +associated with the original points set will be lost, as the algorithm has +limited ability to impute associated data. However, if color dimensions +(red, green and blue) are present in the input, colors will be reconstruced +in the output point set. .. [Kazhdan2006] Kazhdan, Michael, Matthew Bolitho, and Hugues Hoppe. "Poisson surface reconstruction." Proceedings of the fourth Eurographics symposium on Geometry processing. Vol. 7. 2006. -.. _`PCL`: http://www.pointclouds.org +This integration of the algorithm with PDAL only supports a limited set of +the options available to the implementation. If you need support for further +options, please let us know. Example ------------------------------------------------------------------------------- @@ -23,12 +35,10 @@ Example "dense.las", { "type":"filters.poisson", - "depth":"8", - "point_weight":"4" }, { - "type":"writers.las", - "filename":"thinned.las", + "type":"writers.ply", + "filename":"isosurface.ply", } ] } @@ -37,8 +47,12 @@ Example Options ------------------------------------------------------------------------------- +density + Write an estimate of neighborhood density for each point in the output + set. + depth - Maximum depth of the tree used for reconstruction. [Default: **8**] + Maximum depth of the tree used for reconstruction. The output is sentsitve + to this parameter. Increase if the results appear unsatisfactory. + [Default: **8**] -point_weight - Importance of interpolation of point samples in the screened Poisson equation. [Default: **4.0**] diff --git a/doc/tutorial/python-filtering.rst b/doc/tutorial/python-filtering.rst index 28d4c71b96..30492fcc7f 100644 --- a/doc/tutorial/python-filtering.rst +++ b/doc/tutorial/python-filtering.rst @@ -33,9 +33,8 @@ be a very quick way to prototype a tool that identified specific points we would like to filter. PDAL has three different ways to manipulate data with Python. The first is -:ref:`filters.python`, which we will be using in this tutorial. The -second is :ref:`filters.`, which allows you to keep or remove points -given a Python filtering operation. The third is the Python extension at +:ref:`filters.python`, which we will be using in this tutorial. +The second is the Python extension at https://pypi.python.org/pypi/PDAL that allows you to utilize PDAL processing operations in your own Python programs. @@ -75,12 +74,11 @@ running. Python Filter ------------------------------------------------------------------------------- -Through the use of the :ref:`filters.python` and :ref:`filters.` -filters, PDAL allows the use of Python and |NumPy| to process point cloud -data. This can be very useful in prototyping situations, where PDAL can provide -convenient data access and the processing logic of software that is still -taking shape can be constructed with the rapid prototyping tools that Python -can provide. +Through the use of the :ref:`filters.python`, PDAL allows the use of Python and +|NumPy| to process point cloud data. This can be very useful in prototyping +situations, where PDAL can provide convenient data access and the processing +logic of software that is still taking shape can be constructed with the rapid +prototyping tools that Python can provide. .. code-block:: python @@ -112,7 +110,7 @@ can provide. # Print our dict to stdout print output - # filters. must return True to tell + # filters.python must return True to tell # PDAL it successfully completed return True @@ -221,7 +219,7 @@ three standard deviations: # Print our dict to stdout print output - # filters. must return True to tell + # filters.python must return True to tell # PDAL it successfully completed return True diff --git a/doc/tutorial/writing-reader.rst b/doc/tutorial/writing-reader.rst index af5a92181a..be19a9b237 100644 --- a/doc/tutorial/writing-reader.rst +++ b/doc/tutorial/writing-reader.rst @@ -33,15 +33,7 @@ These methods are required to fulfill the specs for defining a new plugin. .. literalinclude:: ../../examples/writing-reader/MyReader.hpp :language: cpp - :lines: 20 - :linenos: - -``getDefaultDimensions`` returns a list of :ref:`dimensions` that the -reader provides. - -.. literalinclude:: ../../examples/writing-reader/MyReader.hpp - :language: cpp - :lines: 23-25 + :lines: 21-23 :linenos: ``m_stream`` is used to process the input, while ``m_index`` is used to track @@ -111,16 +103,7 @@ dimension MyData. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 31-40 - :linenos: - -This method returns the list of :ref:`dimensions` that the reader can -provide. - - -.. literalinclude:: ../../examples/writing-reader/MyReader.cpp - :language: cpp - :lines: 42-46 + :lines: 31-35 :linenos: This method is called when the Reader is ready for use. It will only be @@ -129,7 +112,7 @@ processed. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 49-62 + :lines: 37-50 :linenos: This is a helper function, which will convert a string value into the type @@ -139,7 +122,7 @@ strings to doubles when reading from the input stream. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 65 + :lines: 53 :linenos: This method is the main processing method for the reader. It takes a @@ -152,7 +135,7 @@ PointView object. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 74-76 + :lines: 62-64 :linenos: In preparation for reading the file, we prepare to skip some header lines. In @@ -160,7 +143,7 @@ our case, the header is only a single line. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 77-82 + :lines: 65-70 :linenos: Here we begin our main loop. In our example file, the first line is a header, @@ -171,7 +154,7 @@ sure we are skipping the header lines here before moving on. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 85-94 + :lines: 73-82 :linenos: Here we take the line we read in the for block header, split it, and make sure @@ -179,10 +162,9 @@ that we have the proper number of fields. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 96-109 + :lines: 84-97 :linenos: - Here we take the values we read and put them into the PointView object. The X and Y fields are simply converted from the file and put into the respective fields. MyData is done likewise with the custom dimension we defined. The Z @@ -195,7 +177,7 @@ each iteration of the loop), and the dimension value. .. literalinclude:: ../../examples/writing-reader/MyReader.cpp :language: cpp - :lines: 111-113 + :lines: 99-101 :linenos: Finally, we increment the nextId and make a call into the progress callback diff --git a/doc/tutorial/writing-writer.rst b/doc/tutorial/writing-writer.rst index a4c99db8f8..373c6a495f 100644 --- a/doc/tutorial/writing-writer.rst +++ b/doc/tutorial/writing-writer.rst @@ -58,7 +58,7 @@ that corresponds to the data field for easier lookup. As mentioned, there cen be additional configurations done as needed. -The header +The source ------------------------------------------------------------------------------- We will start with a full listing of the writer source. diff --git a/doc/workshop/agenda.rst b/doc/workshop/agenda.rst index 11db4c8986..b51993d10f 100644 --- a/doc/workshop/agenda.rst +++ b/doc/workshop/agenda.rst @@ -23,7 +23,7 @@ Materials Slides ................................................................................ -* `Slides `__ +* `Slides `__ Workshop Materials ................................................................................ @@ -31,8 +31,8 @@ Workshop Materials These materials are available at http://pdal.io/workshop/ as both a PDF and an HTML website. -* `PDF download `__ -* `Website `__ +* `PDF download `__ +* `Website `__ USB Example Data Drive ................................................................................ diff --git a/doc/workshop/exercises/analysis/boundary/boundary-command-boundary.txt b/doc/workshop/exercises/analysis/boundary/boundary-command-boundary.txt index 60ec810a42..d8cd870dad 100644 --- a/doc/workshop/exercises/analysis/boundary/boundary-command-boundary.txt +++ b/doc/workshop/exercises/analysis/boundary/boundary-command-boundary.txt @@ -1,4 +1,3 @@ pdal info ^ -c:/Users/hobu/PDAL/exercises/analysis/density/uncompahgre.laz ^ ---boundary - + c:/Users/hobu/PDAL/exercises/analysis/density/uncompahgre.laz ^ + --boundary diff --git a/doc/workshop/exercises/analysis/boundary/boundary-command-tindex.txt b/doc/workshop/exercises/analysis/boundary/boundary-command-tindex.txt index 25de72798d..109ad0cd0b 100644 --- a/doc/workshop/exercises/analysis/boundary/boundary-command-tindex.txt +++ b/doc/workshop/exercises/analysis/boundary/boundary-command-tindex.txt @@ -1,4 +1,4 @@ pdal tindex ^ ---tindex c:/Users/hobu/PDAL/exercises/analysis/boundary/boundary.sqlite ^ ---filespec c:/Users/hobu/PDAL/exercises/analysis/density/uncompahgre.laz ^ --f SQLite + --tindex c:/Users/hobu/PDAL/exercises/analysis/boundary/boundary.sqlite ^ + --filespec c:/Users/hobu/PDAL/exercises/analysis/density/uncompahgre.laz ^ + -f SQLite diff --git a/doc/workshop/exercises/analysis/clipping/attributes.vrt b/doc/workshop/exercises/analysis/clipping/attributes.vrt index 1362d722df..ca52240347 100644 --- a/doc/workshop/exercises/analysis/clipping/attributes.vrt +++ b/doc/workshop/exercises/analysis/clipping/attributes.vrt @@ -1,7 +1,7 @@ - /data/exercises/analysis/clipping/attributes.json + c:/Users/hobu/PDAL/exercises/analysis/clipping/attributes.json EPSG:4326 +proj=lcc +lat_1=43 +lat_2=45.5 +lat_0=41.75 +lon_0=-120.5 +x_0=399999.9999999999 +y_0=0 +ellps=GRS80 +units=ft +no_defs diff --git a/doc/workshop/exercises/analysis/clipping/clipping.rst b/doc/workshop/exercises/analysis/clipping/clipping.rst index af6e060028..ab9ce40fed 100644 --- a/doc/workshop/exercises/analysis/clipping/clipping.rst +++ b/doc/workshop/exercises/analysis/clipping/clipping.rst @@ -124,8 +124,8 @@ Visualization ................................................................................ Use one of the point cloud visualization tools you installed to take a look at -your ``c:\Users\Howard\PDAL\exercises\analysis\clipping\stadium.las`` output. -In the example below, we simply opened the file using the http://plas.io +your ``c:/Users/Howard/PDAL/exercises/analysis/clipping/stadium.las`` output. +In the example below, we opened the file to view it using the http://plas.io website. @@ -142,7 +142,6 @@ Notes 2. Points that are *on* the boundary are included. -.. _`NASA Airborne Snow Observatory`: http://aso.jpl.nasa.gov/ .. _`CloudCompare`: http://www.danielgm.net/cc/ .. _`ASPRS LAS`: http://www.asprs.org/Committee-General/LASer-LAS-File-Format-Exchange-Activities.html diff --git a/doc/workshop/exercises/analysis/denoising/denoising.rst b/doc/workshop/exercises/analysis/denoising/denoising.rst index aae0b90cbd..2e2e753f22 100644 --- a/doc/workshop/exercises/analysis/denoising/denoising.rst +++ b/doc/workshop/exercises/analysis/denoising/denoising.rst @@ -1,14 +1,13 @@ -.. _denoising: +.. _workshop-denoising: Removing noise ================================================================================ .. include:: ../../../includes/substitutions.rst -.. index:: Denoising, Filtering - -This exercise uses PDAL to remove unwanted noise in an ALS collection. +.. index:: Denoising, outliers +This exercise uses PDAL to remove unwanted noise in an airborne LiDAR collection. Exercise @@ -83,6 +82,8 @@ Both :ref:`ranges ` are passed as a comma-separated list to the "limits": "Classification![7:7],Z[-100:3000]" }, +.. index:: range filter, classifications + 4. :ref:`writers.las` ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ diff --git a/doc/workshop/exercises/analysis/dtm/dtm.rst b/doc/workshop/exercises/analysis/dtm/dtm.rst index f7dd0915b5..eeb7386b2c 100644 --- a/doc/workshop/exercises/analysis/dtm/dtm.rst +++ b/doc/workshop/exercises/analysis/dtm/dtm.rst @@ -1,4 +1,4 @@ -.. _dtm: +.. _workshop-dtm: Generating a DTM ================================================================================ diff --git a/doc/workshop/exercises/analysis/ground/ground.rst b/doc/workshop/exercises/analysis/ground/ground.rst index 4cdfbb239c..48a4c7cffc 100644 --- a/doc/workshop/exercises/analysis/ground/ground.rst +++ b/doc/workshop/exercises/analysis/ground/ground.rst @@ -1,4 +1,4 @@ -.. _ground: +.. _workshop-ground: Identifying ground ================================================================================ diff --git a/doc/workshop/exercises/analysis/thinning/thinning.rst b/doc/workshop/exercises/analysis/thinning/thinning.rst index 2479ff73f6..135401c11a 100644 --- a/doc/workshop/exercises/analysis/thinning/thinning.rst +++ b/doc/workshop/exercises/analysis/thinning/thinning.rst @@ -1,4 +1,4 @@ -.. _thinning: +.. _workshop-thinning: Thinning ================================================================================ diff --git a/doc/workshop/exercises/index.rst b/doc/workshop/exercises/index.rst index 1ec253f403..b4ca067d16 100644 --- a/doc/workshop/exercises/index.rst +++ b/doc/workshop/exercises/index.rst @@ -26,6 +26,7 @@ Translation translation/compression translation/reprojection + translation/greyhound .. _analysisa: @@ -45,6 +46,14 @@ Analysis analysis/ground/ground analysis/dtm/dtm +Python +-------------------------------------------------------------------------------- + +.. toctree:: + :maxdepth: 3 + + python/histogram + Georeferencing -------------------------------------------------------------------------------- diff --git a/doc/workshop/exercises/translation/compression-command.txt b/doc/workshop/exercises/translation/compression-command.txt index 215660b09f..4de5ed7f6f 100644 --- a/doc/workshop/exercises/translation/compression-command.txt +++ b/doc/workshop/exercises/translation/compression-command.txt @@ -1,2 +1,3 @@ -pdal translate c:/Users/hobu/pdal/exercises/translation/interesting.las ^ -c:/Users/hobu/pdal/exercises/translation/interesting.laz +pdal translate ^ + c:/Users/hobu/pdal/exercises/translation/interesting.las ^ + c:/Users/hobu/pdal/exercises/translation/interesting.laz diff --git a/doc/workshop/exercises/translation/reprojection-command-1.txt b/doc/workshop/exercises/translation/reprojection-command-1.txt index f2422d47ae..adef0e3ee9 100644 --- a/doc/workshop/exercises/translation/reprojection-command-1.txt +++ b/doc/workshop/exercises/translation/reprojection-command-1.txt @@ -1,4 +1,4 @@ -pdal translate c:/Users/hobu/PDAL/exercises/analysis/ground/CSite1_orig-utm.laz ^ +pdal translate ^ c:/Users/hobu/PDAL/exercises/translation/csite-dd.laz ^ reprojection ^ --filters.reprojection.out_srs="EPSG:4326" diff --git a/doc/workshop/images/agenda-usb-drive.jpg b/doc/workshop/images/agenda-usb-drive.jpg index acdd72ceab..3b6fd32de9 100644 Binary files a/doc/workshop/images/agenda-usb-drive.jpg and b/doc/workshop/images/agenda-usb-drive.jpg differ diff --git a/doc/workshop/includes/substitutions.rst b/doc/workshop/includes/substitutions.rst index 20b0989713..d0e9cf593b 100644 --- a/doc/workshop/includes/substitutions.rst +++ b/doc/workshop/includes/substitutions.rst @@ -24,7 +24,9 @@ .. |Hobu| replace:: `Hobu `__ .. |Optech| replace:: `Optech `__ .. |Riegl| replace:: `Riegl `__ +.. |Greyhound| replace:: `Greyhound `__ .. |NCALM| replace:: `NCALM `__ .. |UTM| replace:: `UTM `__ .. |WGS84| replace:: `WGS84 `__ .. |WellKnownText| replace:: `Well Known Text `__ +.. |Terminal| replace:: `OSGeo4W Shell` diff --git a/doc/workshop/index.rst b/doc/workshop/index.rst index 3a0556b329..3e72767fe3 100644 --- a/doc/workshop/index.rst +++ b/doc/workshop/index.rst @@ -7,7 +7,7 @@ Point Cloud Processing and Analysis with PDAL :Author: Pete Gadomski :Author: Dr. Craig Glennie :Contact: howard@hobu.co -:Date: 03/30/2016 +:Date: 07/25/2017 .. include:: ./includes/substitutions.rst @@ -16,6 +16,9 @@ Point Cloud Processing and Analysis with PDAL :maxdepth: 3 agenda + about lidar-introduction software exercises/index + capstone + notes diff --git a/doc/workshop/slides/source/basic_info.rst b/doc/workshop/slides/source/basic_info.rst index 9d11ee4c5a..3765e9c028 100644 --- a/doc/workshop/slides/source/basic_info.rst +++ b/doc/workshop/slides/source/basic_info.rst @@ -3,6 +3,9 @@ Basic Information ================================================================================ +.. include:: ../../includes/substitutions.rst + + * Inspection of file contents * Location determination * Investigate suitability @@ -10,46 +13,38 @@ Basic Information Exercises ================================================================================ -1. :ref:`single-point` -2. :ref:`metadata` -3. :ref:`near` +1. :ref:`pdal:workshop-single-point` +2. :ref:`pdal:workshop-metadata` +3. :ref:`pdal:near` Printing a Single Point ================================================================================ Purpose: - * Learn how to run PDAL via :ref:`docker` * Verify PDAL is working correctly * Learn about point cloud data types + * Familiarize yourself with terminal -:ref:`Workshop Materials ` +:ref:`Workshop Materials ` -Command (first point) +Start Terminal (first point) ================================================================================ -In your `Docker Quickstart Terminal`, issue the following: - -.. literalinclude:: ../../exercises/info/single-point-command.txt - :linenos: - -1. ``docker``: All our commands start with ``docker`` - -2. ``run``: Tells docker we're going to run an image +1. Start the |Terminal| -3. ``-v /c/Users/Howard/PDAL:/data``: Maps our workshop directory to a directory called - ``/data`` inside the container. +.. figure:: ./img/osgeo4w-shell.png Command (first point) ================================================================================ +In your |Terminal|, issue the following: + .. literalinclude:: ../../exercises/info/single-point-command.txt :linenos: -4. ``pdal/pdal``: Run the image with name ``pdal/pdal`` - -5. ``pdal``: Inside the ``pdal/pdal`` image, run the `pdal` command +1. ``pdal``: Inside the ``pdal/pdal`` image, run the `pdal` command -6. ``info``: `info` is a "kernel" in PDAL-speak. It is a unit +2. ``info``: `info` is a "kernel" in PDAL-speak. It is a unit of functionality driven by the command line. Command (first point) @@ -58,12 +53,12 @@ Command (first point) .. literalinclude:: ../../exercises/info/single-point-command.txt :linenos: -7. ``/data/exercises/info/interesting.las``: Our directory - is now mounted at ``/data`` +3. ``c:/Users/hobu/PDAL/exercises/info/interesting.las``: + File to run ``info`` command on -8. ``-p``: ``-p`` corresponds to "print a point". +4. ``-p``: argument corresponds to "print a point". -9. ``0`` means to print the first +5. ``0`` means to print the first (starting from 0) Run (first point) ================================================================================ @@ -87,12 +82,12 @@ Purpose: * Compute bounding box * Retrieve supporting information -:ref:`Metadata Workshop Materials ` +:ref:`Metadata Workshop Materials ` Command (metadata) ================================================================================ -In your `Docker Quickstart Terminal`, issue the following: +In your |Terminal|, issue the following: .. literalinclude:: ../../exercises/info/metadata-command.txt :linenos: @@ -121,7 +116,7 @@ Purpose: Command (near - metadata) ================================================================================ -In your `Docker Quickstart Terminal`, issue the following: +In your |Terminal|, issue the following: .. literalinclude:: ../../exercises/info/near-command-1.txt :linenos: @@ -133,7 +128,7 @@ In your `Docker Quickstart Terminal`, issue the following: Command (near - query) ================================================================================ -In your `Docker Quickstart Terminal`, issue the following: +In your |Terminal|, issue the following: .. literalinclude:: ../../exercises/info/near-command-2.txt :linenos: diff --git a/doc/workshop/slides/source/boundary.rst b/doc/workshop/slides/source/boundary.rst index c0a597d9b6..468fe5ee51 100644 --- a/doc/workshop/slides/source/boundary.rst +++ b/doc/workshop/slides/source/boundary.rst @@ -3,13 +3,16 @@ Boundary ================================================================================ + +.. include:: ../../includes/substitutions.rst + Purpose: * Estimate spatial footprint * Do it quickly * Balance speed/quality * Look with :ref:`qgis` -:ref:`Boundary Workshop Materials ` +:ref:`Boundary Workshop Materials ` Boundary (command) ================================================================================ diff --git a/doc/workshop/slides/source/clipping.rst b/doc/workshop/slides/source/clipping.rst index 432e670fef..bea7708de6 100644 --- a/doc/workshop/slides/source/clipping.rst +++ b/doc/workshop/slides/source/clipping.rst @@ -3,12 +3,16 @@ Cliping data with polygons ================================================================================ +.. include:: ../../includes/substitutions.rst + + + Purpose: * Subset data with a polygon * Assign attribute data * Use PDAL pipeline -:ref:`Clipping Workshop Materials ` +:ref:`Clipping Workshop Materials ` Clipping (Autzen) ================================================================================ @@ -94,6 +98,17 @@ Clipping (verify) .. image:: ../../images/clipping-stadium-clipped.png + +Other ways to clip +================================================================================ + +* Clip using multiple :ref:`filters.range` +* :ref:`filters.divider` or :ref:`filters.chipper` + + .. image:: ../../../stages/filters.chipper.img2.png + + + Next ================================================================================ diff --git a/doc/workshop/slides/source/colorization.rst b/doc/workshop/slides/source/colorization.rst index 89d1c2897f..a2a10a8c5e 100644 --- a/doc/workshop/slides/source/colorization.rst +++ b/doc/workshop/slides/source/colorization.rst @@ -3,12 +3,14 @@ Colorizing points with imagery ================================================================================ +.. include:: ../../includes/substitutions.rst + Purpose: * Top-down projection of RGB imagery * Improve visualization * Assign other GDAL data to points -:ref:`Colorizing Workshop Materials ` +:ref:`Colorizing Workshop Materials ` Colorization (setup) diff --git a/doc/workshop/slides/source/conf.py b/doc/workshop/slides/source/conf.py index f9f6aa8c0a..77052c5676 100644 --- a/doc/workshop/slides/source/conf.py +++ b/doc/workshop/slides/source/conf.py @@ -35,6 +35,7 @@ 'sphinx.ext.coverage', 'sphinx.ext.pngmath', 'sphinx.ext.ifconfig', + 'hieroglyph', ] # Add any paths that contain templates here, relative to this directory. @@ -53,7 +54,7 @@ # General information about the project. project = u'Point Cloud Processing with PDAL' -copyright = u'2016, Howard Butler' +copyright = u'2017, Howard Butler' author = u'Howard Butler' # The version info for the project you're documenting, acts as replacement for @@ -61,9 +62,9 @@ # built documents. # # The short X.Y version. -version = u'2016.03.23' +version = u'2017.07.25' # The full version, including alpha/beta/rc tags. -release = u'2016.03.23' +release = u'2017.07.25' # The language for content autogenerated by Sphinx. Refer to documentation # for a list of supported languages. @@ -314,6 +315,6 @@ # Example configuration for intersphinx: refer to the Python standard library. -intersphinx_mapping = {'pdal': ('http://www.pdal.io/', None), +intersphinx_mapping = {'pdal': ('https://www.pdal.io/', None), 'unavco': ('../../_build/html/../../../../_build/html/', '../../_build/html/objects.inv')} diff --git a/doc/workshop/slides/source/denoising.rst b/doc/workshop/slides/source/denoising.rst index c6cdf259bd..063561523f 100644 --- a/doc/workshop/slides/source/denoising.rst +++ b/doc/workshop/slides/source/denoising.rst @@ -3,11 +3,13 @@ Removing noise ================================================================================ +.. include:: ../../includes/substitutions.rst + Purpose: * Remove noise with a statistical method * Maintain point integrity -:ref:`Denoising Workshop Materials ` +:ref:`Denoising Workshop Materials ` Pipeline @@ -28,7 +30,6 @@ Denoising (pipeline) ================================================================================ .. literalinclude:: ../../exercises/analysis/denoising/denoise.json - :linenos: Denoising (command) ================================================================================ diff --git a/doc/workshop/slides/source/density.rst b/doc/workshop/slides/source/density.rst index c4a5a84b2f..3f26e99205 100644 --- a/doc/workshop/slides/source/density.rst +++ b/doc/workshop/slides/source/density.rst @@ -3,10 +3,13 @@ Visualizing acquisition density ================================================================================ + +.. include:: ../../includes/substitutions.rst + Purpose: * Characterize collection density -:ref:`Density Workshop Materials ` +:ref:`Density Workshop Materials ` Density (execution) diff --git a/doc/workshop/slides/source/dtm.rst b/doc/workshop/slides/source/dtm.rst index f1cb8e8006..bf18103308 100644 --- a/doc/workshop/slides/source/dtm.rst +++ b/doc/workshop/slides/source/dtm.rst @@ -3,20 +3,26 @@ Generating a DTM ================================================================================ + +.. include:: ../../includes/substitutions.rst + + + Purpose: * Generate a DTM from ground-filtered data * Data was created in :ref:`ground` * Visualize DTM using :ref:`qgis` -:ref:`DTM Workshop Materials ` +:ref:`DTM Workshop Materials ` writers.gdal ================================================================================ * :ref:`writers.gdal` -* Generated using points2grid (OpenTopography) +* Ported points2grid (OpenTopography) +* IDW available for fill-in * Write TIFF/ASCII raster * Control pixel size @@ -24,7 +30,6 @@ DTM (pipeline) ================================================================================ .. literalinclude:: ../../exercises/analysis/dtm/gdal.json - :linenos: DTM (execution) ================================================================================ diff --git a/doc/workshop/slides/source/ground.rst b/doc/workshop/slides/source/ground.rst index 2fb6415c71..30f22a0642 100644 --- a/doc/workshop/slides/source/ground.rst +++ b/doc/workshop/slides/source/ground.rst @@ -3,11 +3,13 @@ Identifying ground ================================================================================ +.. include:: ../../includes/substitutions.rst + + Purpose: * Filter and classify ground points - -:ref:`Ground Workshop Materials ` +:ref:`Ground Workshop Materials ` filters.ground diff --git a/doc/workshop/slides/source/introduction.rst b/doc/workshop/slides/source/introduction.rst index e455d01c64..c52492672e 100644 --- a/doc/workshop/slides/source/introduction.rst +++ b/doc/workshop/slides/source/introduction.rst @@ -3,11 +3,10 @@ Order of business ================================================================================ -* Materials https://hobu.co/workshops/unavco/ -* Follow along https://hobu.co/workshops/unavco/slides/ -* PDF https://hobu.co/workshops/unavco/PDAL-UNAVCO-ShortCourse.pdf +* Materials https://s3.amazonaws.com/pdal/workshop/PDAL.zip +* Follow along https://pdal.io/workshop/slides/ +* PDF http://pdal.io/PDAL.pdf * Links everywhere -* First time guinea pigs Guinea Pigs ================================================================================ @@ -28,8 +27,9 @@ Who are we ================================================================================ * Howard Butler -* Pete Gadomski -* Dr. Craig Glennie +* Connor Manning +* Rick Brown + Howard Butler ================================================================================ @@ -41,20 +41,20 @@ Howard Butler .. _`Hobu, Inc.`: http://hobu.co -Pete Gadomski +Connor Manning ================================================================================ -* Physical Scientist at CRREL -* Graduate student in Geosensing Systems Engineering and Sciences at UH -* Writes code, collects data -* https://github.com/gadomski +* Works @ `Hobu, Inc.`_ +* Develops Entwine https://entwine.io +* Develops Greyhound http://greyhound.io +* Tweets @ http://twitter.com/csmannin -Dr. Craig Glennie +Rick Brown ================================================================================ -* Co-PI National Center for Airborne Laser Mapping (NCALM) -* University of Houston -* Formerly VP Engineering, Terrapoint +* Works here! +* LiDAR data processing and management +* Entwine with Connor Manning Next ================================================================================ diff --git a/doc/workshop/slides/source/lidar_intro.rst b/doc/workshop/slides/source/lidar_intro.rst index ac3e05957c..122770b77b 100644 --- a/doc/workshop/slides/source/lidar_intro.rst +++ b/doc/workshop/slides/source/lidar_intro.rst @@ -35,7 +35,7 @@ How is LiDAR used? Discrete vs. full-waveform ============================ -- Outgoing puses are *not* instantaneous — they have a finite width and height (approximated by a Gaussian) +- Outgoing pulses are *not* instantaneous — they have a finite width and height (approximated by a Gaussian) - Return energy is not instantaneous either (see picture) - The "full waveform" output (see picture) is usually simplified down to single points (discrete-return) @@ -55,7 +55,7 @@ How is LiDAR processed? * Vendor-specific software (e.g. RiSCAN Pro and RiProcess, from |Riegl|) * Other commercial softwares (e.g. TerraStation, QT Modeler) -* Mixed-source software (e.g. LasTools) +* Mixed-source software (e.g. LASTools) * Open-source software (e.g. CloudCompare, PDAL, laspy) diff --git a/doc/workshop/slides/source/pdal_intro.rst b/doc/workshop/slides/source/pdal_intro.rst index 336d305551..1f8af8be29 100644 --- a/doc/workshop/slides/source/pdal_intro.rst +++ b/doc/workshop/slides/source/pdal_intro.rst @@ -3,6 +3,8 @@ Introduction to PDAL ================================================================================ +.. include:: ../../includes/substitutions.rst + * Point Data Abstraction Library * "GDAL for point cloud data" * Focus (in priority) @@ -19,11 +21,12 @@ Introduction to PDAL Open Source ================================================================================ -* `BSD`_ license +* https://pdal.io * https://github.com/PDAL/PDAL * `Tested with every commit`_ * 30+ `contributors`_ * Driven by real world use +* `BSD`_ license .. _`Tested with every commit`: http://travis-ci.org/PDAL/PDAL/builds/ .. _`BSD`: https://opensource.org/licenses/BSD-2-Clause @@ -70,7 +73,7 @@ Processing Pipeline Pipeline Architecture ================================================================================ -.. image:: ../../images/pdal-reader-writer.png +.. image:: ../../../images/las-crop-bpf-pipeline.png Pipeline Architecture ================================================================================ @@ -128,11 +131,11 @@ Filters Filters (cont) ================================================================================ -* :ref:`Ground classification` -* :ref:`Normalized heights ` +* :ref:`Ground classification` +* :ref:`Normalized heights ` * :ref:`Sorting ` * :ref:`Spatial curve sorting ` -* :ref:`Noise filtering ` +* :ref:`Noise filtering ` * :ref:`Merge` Filters (cont) @@ -181,6 +184,15 @@ Python (cont) outs['Mask'] = keep return True +Matlab +================================================================================ + +* :ref:`writers.matlab` and :ref:`readers.matlab` for i/o + +* Embed + + * Use Matlab scripts inline with :ref:`pipeline ` operations + Next ================================================================================ diff --git a/doc/workshop/slides/source/software.rst b/doc/workshop/slides/source/software.rst index f30f1c7ac8..7160da5b8f 100644 --- a/doc/workshop/slides/source/software.rst +++ b/doc/workshop/slides/source/software.rst @@ -3,29 +3,33 @@ Software installation ================================================================================ -1. :ref:`Docker ` +.. include:: ../../includes/substitutions.rst + +(Done for you already) + +.. 1. :ref:`Docker ` * Run PDAL and GDAL command line applications -2. :ref:`QGIS` +.. 2. :ref:`QGIS` * Visualize vector and raster data -3. `QT Reader`_ +.. 3. `QT Reader`_ * Visualize point cloud data -4. `Fugro Viewer`_ +.. 4. `Fugro Viewer`_ * Visualize point cloud data -5. http://plas.io +.. 5. http://plas.io * Visualize point cloud data -6. `CloudCompare`_ +.. 6. `CloudCompare`_ * Visualize point cloud data .. _`QT Reader`: http://appliedimagery.com/download/ .. _`Fugro Viewer`: http://www.fugroviewer.com/ .. _`CloudCompare`: http://www.danielgm.net/cc/ -File copying +Data installation ================================================================================ 1. Copy entire contents: @@ -34,85 +38,7 @@ File copying C:\Users\Howard\PDAL -2. Navigate to ``software`` directory: - -:: - - c:\Users\Howard\PDAL\software - - - -Docker Installation -================================================================================ - -.. image:: ../../images/docker-file-navigate.png - -1. :ref:`Docker ` -2. :ref:`QGIS` - - -Docker -================================================================================ - -* Containerization -* Not virtual machines -* Operating system isolation - -Virtual Machine -================================================================================ - -.. image:: ./img/vm-diagram.png - :target: https://www.docker.com/what-docker - :scale: 75% - -Docker -================================================================================ - -.. image:: ./img/docker-diagram.png - :target: https://www.docker.com/what-docker - :scale: 75% - -Why Docker? -================================================================================ - -* Rapid development -* Compiling is pain -* Windows development is pain -* Old software is pain - -Compile yourself? -================================================================================ - -.. image:: ./img/life-so-hard.gif - :scale: 50% - -Benefits -================================================================================ - -* Freshies - * Short bug fix cycle - -* Known configuration - * Developer configures - * *All* the features - -* Replicability of configuration - * Runs the same on all machines - * Inspectable - * Repeatable bugs - -How it works -================================================================================ - -* Linux -* On Windows, you actually run on a Linux VM - * Same on Mac - * We will use `VirtualBox`_ - * `Docker Toolbox`_ does all of it for us - * Both free and open source -.. _`VirtualBox`: https://www.virtualbox.org/wiki/Downloads -.. _`Docker Toolbox`: https://www.docker.com/products/docker-toolbox Next ================================================================================ diff --git a/doc/workshop/slides/source/thinning.rst b/doc/workshop/slides/source/thinning.rst index efa51886bb..fd24ad9aea 100644 --- a/doc/workshop/slides/source/thinning.rst +++ b/doc/workshop/slides/source/thinning.rst @@ -3,16 +3,18 @@ Thinning ================================================================================ +.. include:: ../../includes/substitutions.rst + Purpose: * Reduce collection density -:ref:`Thinning Workshop Materials ` +:ref:`Thinning Workshop Materials ` PDAL Thinning Techniques ================================================================================ -* :ref:`filters.dartsample` +* :ref:`filters.sample` * :ref:`filters.voxelgrid` * :ref:`filters.randomize` + :ref:`filters.decimation` diff --git a/doc/workshop/slides/source/translation.rst b/doc/workshop/slides/source/translation.rst index 53396f3533..c2d440fa76 100644 --- a/doc/workshop/slides/source/translation.rst +++ b/doc/workshop/slides/source/translation.rst @@ -3,6 +3,9 @@ Translation ================================================================================ +.. include:: ../../includes/substitutions.rst + + * Convert one format to another * Do stuff along the way @@ -10,13 +13,12 @@ Compression ================================================================================ Purpose: - * Output a compressed `LAZ`_ + * Output a compressed |LASzip| * Learn about point cloud metadata * Learn about |ASPRSLAS| -:ref:`Compression Workshop Materials ` +:ref:`Compression Workshop Materials ` -.. include:: ../../includes/substitutions.rst .. _`LAZ`: http://laszip.org LASzip @@ -54,15 +56,14 @@ Purpose: * Utilize driver options * Scale output data -:ref:`Reprojection Workshop Materials ` +:ref:`Reprojection Workshop Materials ` Command (reproject) ================================================================================ .. literalinclude:: ../../exercises/translation/reprojection-command-1.txt - :linenos: -* We tell ``filters.reprojection`` to output to ``EPSG:4326`` +* We tell :ref:`filters.reprojection` to output to ``EPSG:4326`` * Add ``reprojection`` filter directly to ``translate`` command * Define reprojection filter option via command line @@ -75,7 +76,6 @@ Command (scale) ================================================================================ .. literalinclude:: ../../exercises/translation/reprojection-command-2.txt - :linenos: * Scale set to ``1e-7`` * Offset to ``auto`` (PDAL calculates minimum) diff --git a/doc/workshop/software.rst b/doc/workshop/software.rst index c16f8f3f8c..22d66a37c7 100644 --- a/doc/workshop/software.rst +++ b/doc/workshop/software.rst @@ -9,8 +9,3 @@ Software Installation :maxdepth: 2 osgeo4w - qgis - - - - diff --git a/examples/writing-reader/MyReader.cpp b/examples/writing-reader/MyReader.cpp index b662ef588b..12961ab143 100644 --- a/examples/writing-reader/MyReader.cpp +++ b/examples/writing-reader/MyReader.cpp @@ -28,24 +28,12 @@ namespace pdal layout->registerOrAssignDim("MyData", Dimension::Type::Unsigned64); } - Dimension::IdList MyReader::getDefaultDimensions() - { - Dimension::IdList ids; - - ids.push_back(Dimension::Id::X); - ids.push_back(Dimension::Id::Y); - ids.push_back(Dimension::Id::Z); - - return ids; - } - void MyReader::ready(PointTableRef) { SpatialReference ref("EPSG:4385"); setSpatialReference(ref); } - template T convert(const StringList& s, const std::string& name, size_t fieldno) { diff --git a/examples/writing-reader/MyReader.hpp b/examples/writing-reader/MyReader.hpp index e7b8ff01e2..76057fae8f 100644 --- a/examples/writing-reader/MyReader.hpp +++ b/examples/writing-reader/MyReader.hpp @@ -17,8 +17,6 @@ namespace pdal static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions(); - private: std::unique_ptr m_stream; point_count_t m_index; diff --git a/filters/ChipperFilter.cpp b/filters/ChipperFilter.cpp index 868dcd4ba9..65d8d63f1f 100644 --- a/filters/ChipperFilter.cpp +++ b/filters/ChipperFilter.cpp @@ -251,13 +251,8 @@ void ChipperFilter::split(ChipRefList& wide, ChipRefList& narrow, ChipRefList& s } } - // Save away the direction so we know which array is X and which is Y - // so that when we emit, we can properly label the max/min points. - Direction dir = narrow.m_dir; - spare.m_dir = dir; decideSplit(wide, spare, narrow, pleft, pcenter); decideSplit(wide, spare, narrow, pcenter, pright); - narrow.m_dir = dir; } } diff --git a/filters/ChipperFilter.hpp b/filters/ChipperFilter.hpp index 88bb548900..74985a4695 100644 --- a/filters/ChipperFilter.hpp +++ b/filters/ChipperFilter.hpp @@ -57,14 +57,6 @@ class Stage; class PDAL_DLL ChipperFilter; -enum Direction -{ - DIR_X, - DIR_Y, - DIR_NONE -}; - - class PDAL_DLL ChipPtRef { friend class ChipRefList; @@ -89,10 +81,7 @@ class PDAL_DLL ChipRefList private: std::vector m_vec; - Direction m_dir; - ChipRefList(Direction dir = DIR_NONE) : m_dir(dir) - {} std::vector::size_type size() const { return m_vec.size(); @@ -121,25 +110,13 @@ class PDAL_DLL ChipRefList { return m_vec[pos]; } - std::string Dir() - { - if (m_dir == DIR_X) - return "X"; - else if (m_dir == DIR_Y) - return "Y"; - else - return "NONE"; - } }; class PDAL_DLL ChipperFilter : public pdal::Filter { public: - ChipperFilter() : Filter(), - m_xvec(DIR_X), m_yvec(DIR_Y), m_spare(DIR_NONE) - {} - + ChipperFilter() {} static void * create(); static int32_t destroy(void *); std::string getName() const; diff --git a/filters/CropFilter.cpp b/filters/CropFilter.cpp index 825ffb1ca0..7759cf215e 100644 --- a/filters/CropFilter.cpp +++ b/filters/CropFilter.cpp @@ -94,6 +94,11 @@ void CropFilter::initialize() m_geoms.push_back(poly); } } + + m_boxes.clear(); + for (auto& bound : m_bounds) + m_boxes.push_back(bound.to2d()); + m_distance2 = m_distance * m_distance; } @@ -102,7 +107,13 @@ void CropFilter::ready(PointTableRef table) { // If the user didn't provide an SRS, take one from the table. if (m_assignedSrs.empty()) + { m_assignedSrs = table.anySpatialReference(); + if (!table.spatialReferenceUnique()) + log()->get(LogLevel::Warning) << "Can't determine spatial " + "reference for provided bounds. Consider using 'a_srs' " + "option.\n"; + } for (auto& geom : m_geoms) geom.setSpatialReference(m_assignedSrs); } @@ -114,8 +125,8 @@ bool CropFilter::processOne(PointRef& point) if (!crop(point, geom)) return false; - for (auto& box : m_bounds) - if (!crop(point, box.to2d())) + for (auto& box : m_boxes) + if (!crop(point, box)) return false; for (auto& center: m_centers) @@ -153,17 +164,19 @@ void CropFilter::transform(const SpatialReference& srs) throwError("Unable to transform crop geometry to point " "coordinate system."); - for (auto& box : m_bounds) + for (auto& box : m_boxes) { - BOX3D b3d = box.to3d(); - gdal::reprojectBounds(b3d, m_assignedSrs.getWKT(), srs.getWKT()); - box = b3d; + if (!gdal::reprojectBounds(box, m_assignedSrs.getWKT(), srs.getWKT())) + throwError("Unable to reproject bounds."); } for (auto& point : m_centers) { - gdal::reprojectPoint(point.x, point.y, point.z, - m_assignedSrs.getWKT(), srs.getWKT()); + if (!gdal::reprojectPoint(point.x, point.y, point.z, + m_assignedSrs.getWKT(), srs.getWKT())) + throwError("Unable to reproject point center."); } + // Set the assigned SRS for the points/bounds to the one we've + // transformed to. m_assignedSrs = srs; } @@ -180,10 +193,10 @@ PointViewSet CropFilter::run(PointViewPtr view) viewSet.insert(outView); } - for (auto& box : m_bounds) + for (auto& box : m_boxes) { PointViewPtr outView = view->makeNew(); - crop(box.to2d(), *view, *outView); + crop(box, *view, *outView); viewSet.insert(outView); } diff --git a/filters/CropFilter.hpp b/filters/CropFilter.hpp index 197eb9e5ab..e9339d1285 100644 --- a/filters/CropFilter.hpp +++ b/filters/CropFilter.hpp @@ -71,6 +71,7 @@ class PDAL_DLL CropFilter : public Filter double m_distance2; std::vector m_centers; std::vector m_geoms; + std::vector m_boxes; void addArgs(ProgramArgs& args); virtual void initialize(); diff --git a/filters/PMFFilter.cpp b/filters/PMFFilter.cpp index 373181d89b..a3fa1c7e98 100644 --- a/filters/PMFFilter.cpp +++ b/filters/PMFFilter.cpp @@ -32,15 +32,18 @@ * OF SUCH DAMAGE. ****************************************************************************/ +// PDAL implementation of K. Zhang, S.-C. Chen, D. Whitman, M.-L. Shyu, J. Yan, +// and C. Zhang, “A progressive morphological filter for removing nonground +// measurements from airborne LIDAR data,” Geosci. Remote Sensing, IEEE Trans., +// vol. 41, no. 4, pp. 872–882, 2003. + #include "PMFFilter.hpp" #include #include -#include #include #include #include -#include #include "private/DimRange.hpp" @@ -60,14 +63,15 @@ std::string PMFFilter::getName() const void PMFFilter::addArgs(ProgramArgs& args) { - args.add("max_window_size", "Maximum window size", m_maxWindowSize, 33.0); - args.add("slope", "Slope", m_slope, 1.0); - args.add("max_distance", "Maximum distance", m_maxDistance, 2.5); - args.add("initial_distance", "Initial distance", m_initialDistance, 0.15); args.add("cell_size", "Cell size", m_cellSize, 1.0); - args.add("approximate", "Use approximate algorithm?", m_approximate); + args.add("exponential", "Exponential growth of window size?", m_exponential, + true); args.add("ignore", "Ignore values", m_ignored); + args.add("initial_distance", "Initial distance", m_initialDistance, 0.15); args.add("last", "Consider last returns only?", m_lastOnly, true); + args.add("max_distance", "Maximum distance", m_maxDistance, 2.5); + args.add("max_window_size", "Maximum window size", m_maxWindowSize, 33.0); + args.add("slope", "Slope", m_slope, 1.0); } void PMFFilter::addDimensions(PointLayoutPtr layout) @@ -95,203 +99,121 @@ void PMFFilter::prepared(PointTableRef table) } } -std::vector PMFFilter::morphOpen(PointViewPtr view, float radius) -{ - point_count_t np(view->size()); - - QuadIndex idx(*view); - - std::vector minZ(np), maxZ(np); - - // erode - for (PointId i = 0; i < np; ++i) - { - double x = view->getFieldAs(Dimension::Id::X, i); - double y = view->getFieldAs(Dimension::Id::Y, i); - - std::vector ids = - idx.getPoints(x - radius, y - radius, x + radius, y + radius); - - double localMin(std::numeric_limits::max()); - for (auto const& j : ids) - { - double z = view->getFieldAs(Dimension::Id::Z, j); - if (z < localMin) - localMin = z; - } - minZ[i] = localMin; - } - - // dilate - for (PointId i = 0; i < np; ++i) - { - double x = view->getFieldAs(Dimension::Id::X, i); - double y = view->getFieldAs(Dimension::Id::Y, i); - - std::vector ids = - idx.getPoints(x - radius, y - radius, x + radius, y + radius); - - double localMax(std::numeric_limits::lowest()); - for (auto const& j : ids) - { - double z = minZ[j]; - if (z > localMax) - localMax = z; - } - maxZ[i] = localMax; - } - - return maxZ; -} - -std::vector PMFFilter::processGround(PointViewPtr view) +PointViewSet PMFFilter::run(PointViewPtr input) { - // Compute the series of window sizes and height thresholds - std::vector htvec; - std::vector wsvec; - int iter = 0; - float ws = 0.0f; - float ht = 0.0f; - - while (ws < m_maxWindowSize) - { - // Determine the initial window size. - if (1) // exponential - ws = m_cellSize * (2.0f * std::pow(2, iter) + 1.0f); - else - ws = m_cellSize * (2.0f * (iter + 1) * 2 + 1.0f); - - // Calculate the height threshold to be used in the next iteration. - if (iter == 0) - ht = m_initialDistance; - else - ht = m_slope * (ws - wsvec[iter - 1]) * m_cellSize + - m_initialDistance; - - // Enforce max distance on height threshold - if (ht > m_maxDistance) - ht = m_maxDistance; - - wsvec.push_back(ws); - htvec.push_back(ht); - - iter++; - } - - std::vector groundIdx; - for (PointId i = 0; i < view->size(); ++i) - groundIdx.push_back(i); - - // Progressively filter ground returns using morphological open - for (size_t j = 0; j < wsvec.size(); ++j) - { - // Limit filtering to those points currently considered ground returns - PointViewPtr ground = view->makeNew(); - for (auto const& i : groundIdx) - ground->appendPoint(*view, i); + PointViewSet viewSet; + if (!input->size()) + return viewSet; - log()->get(LogLevel::Debug) - << "Iteration " << j << " (height threshold = " << htvec[j] - << ", window size = " << wsvec[j] << ")...\n"; + // Segment input view into ignored/kept views. + PointViewPtr ignoredView = input->makeNew(); + PointViewPtr keptView = input->makeNew(); + if (m_ignored.m_id == Dimension::Id::Unknown) + keptView->append(*input); + else + Segmentation::ignoreDimRange(m_ignored, input, keptView, ignoredView); - // Create new cloud to hold the filtered results. Apply the - // morphological opening operation at the current window size. - auto maxZ = morphOpen(ground, wsvec[j] * 0.5); + // Classify remaining points with value of 1. processGround will mark ground + // returns as 2. + for (PointId i = 0; i < keptView->size(); ++i) + keptView->setField(Dimension::Id::Classification, i, 1); - // Find indices of the points whose difference between the source and - // filtered point clouds is less than the current height threshold. - std::vector groundNewIdx; - for (PointId i = 0; i < ground->size(); ++i) - { - double z0 = ground->getFieldAs(Dimension::Id::Z, i); - float diff = z0 - maxZ[i]; - if (diff < htvec[j]) - groundNewIdx.push_back(groundIdx[i]); - } + // Segment kept view into last/other-than-last return views. + PointViewPtr lastView = keptView->makeNew(); + PointViewPtr nonlastView = keptView->makeNew(); + if (m_lastOnly) + Segmentation::segmentLastReturns(keptView, lastView, nonlastView); + else + lastView->append(*keptView); - groundIdx.swap(groundNewIdx); + // Run the actual PMF algorithm. + processGround(lastView); - log()->get(LogLevel::Debug) - << "Ground now has " << groundIdx.size() << " points.\n"; - } + // Prepare the output PointView. + PointViewPtr outView = input->makeNew(); + outView->append(*ignoredView); + outView->append(*nonlastView); + outView->append(*lastView); + viewSet.insert(outView); - return groundIdx; + return viewSet; } -std::vector PMFFilter::fillNearest(PointViewPtr view, size_t rows, - size_t cols, double cell_size, - BOX2D bounds) +void PMFFilter::processGround(PointViewPtr view) { - using namespace Dimension; + // initialize bounds, rows, columns, and surface + BOX2D bounds; + view->calculateBounds(bounds); + size_t cols = ((bounds.maxx - bounds.minx) / m_cellSize) + 1; + size_t rows = ((bounds.maxy - bounds.miny) / m_cellSize) + 1; + // initialize surface to NaN std::vector ZImin(rows * cols, std::numeric_limits::quiet_NaN()); + // loop through all points, identifying minimum Z value for each populated + // cell for (PointId i = 0; i < view->size(); ++i) { - double x = view->getFieldAs(Id::X, i); - double y = view->getFieldAs(Id::Y, i); - double z = view->getFieldAs(Id::Z, i); - - int c = static_cast(floor(x - bounds.minx) / cell_size); - int r = static_cast(floor(y - bounds.miny) / cell_size); - - if (z < ZImin[c * rows + r] || std::isnan(ZImin[c * rows + r])) - ZImin[c * rows + r] = z; + double x = view->getFieldAs(Dimension::Id::X, i); + double y = view->getFieldAs(Dimension::Id::Y, i); + double z = view->getFieldAs(Dimension::Id::Z, i); + int c = static_cast(floor(x - bounds.minx) / m_cellSize); + int r = static_cast(floor(y - bounds.miny) / m_cellSize); + size_t idx = c * rows + r; + if (z < ZImin[idx] || std::isnan(ZImin[idx])) + ZImin[idx] = z; } - // convert cz into PointView + // convert vector to PointView for indexing PointViewPtr temp = view->makeNew(); PointId i(0); for (size_t c = 0; c < cols; ++c) { for (size_t r = 0; r < rows; ++r) { - if (std::isnan(ZImin[c * rows + r])) + size_t idx = c * rows + r; + if (std::isnan(ZImin[idx])) continue; - - temp->setField(Id::X, i, bounds.minx + (c + 0.5) * cell_size); - temp->setField(Id::Y, i, bounds.miny + (r + 0.5) * cell_size); - temp->setField(Id::Z, i, ZImin[c * rows + r]); + double x = bounds.minx + (c + 0.5) * m_cellSize; + double y = bounds.miny + (r + 0.5) * m_cellSize; + temp->setField(Dimension::Id::X, i, x); + temp->setField(Dimension::Id::Y, i, y); + temp->setField(Dimension::Id::Z, i, ZImin[idx]); i++; } } - // make a 2D KDIndex + // build the 2D KD-tree KD2Index kdi(*temp); kdi.build(); + // loop through all cells, and for each NaN, replace with elevation of + // nearest neighbor std::vector out = ZImin; for (size_t c = 0; c < cols; ++c) { for (size_t r = 0; r < rows; ++r) { - if (!std::isnan(out[c * rows + r])) + size_t idx = c * rows + r; + if (!std::isnan(out[idx])) continue; - - // find k nearest points - double x = bounds.minx + (c + 0.5) * cell_size; - double y = bounds.miny + (r + 0.5) * cell_size; + double x = bounds.minx + (c + 0.5) * m_cellSize; + double y = bounds.miny + (r + 0.5) * m_cellSize; int k = 1; std::vector neighbors(k); std::vector sqr_dists(k); kdi.knnSearch(x, y, k, &neighbors, &sqr_dists); - - out[c * rows + r] = - temp->getFieldAs(Dimension::Id::Z, neighbors[0]); + out[idx] = temp->getFieldAs(Dimension::Id::Z, neighbors[0]); } } - return out; -}; - -std::vector PMFFilter::processGroundApprox(PointViewPtr view) -{ - BOX2D bounds; - view->calculateBounds(bounds); + ZImin.swap(out); - size_t cols = ((bounds.maxx - bounds.minx) / m_cellSize) + 1; - size_t rows = ((bounds.maxy - bounds.miny) / m_cellSize) + 1; + // initialize ground indices + std::vector groundIdx; + for (PointId i = 0; i < view->size(); ++i) + groundIdx.push_back(i); // Compute the series of window sizes and height thresholds std::vector htvec; @@ -300,10 +222,11 @@ std::vector PMFFilter::processGroundApprox(PointViewPtr view) float ws = 0.0f; float ht = 0.0f; + // pre-compute window sizes and height thresholds while (ws < m_maxWindowSize) { // Determine the initial window size. - if (1) // exponential + if (m_exponential) ws = m_cellSize * (2.0f * std::pow(2, iter) + 1.0f); else ws = m_cellSize * (2.0f * (iter + 1) * 2 + 1.0f); @@ -325,13 +248,6 @@ std::vector PMFFilter::processGroundApprox(PointViewPtr view) iter++; } - std::vector groundIdx; - for (PointId i = 0; i < view->size(); ++i) - groundIdx.push_back(i); - - std::vector ZImin = - fillNearest(view, rows, cols, m_cellSize, bounds); - // Progressively filter ground returns using morphological open for (size_t j = 0; j < wsvec.size(); ++j) { @@ -339,10 +255,10 @@ std::vector PMFFilter::processGroundApprox(PointViewPtr view) << "Iteration " << j << " (height threshold = " << htvec[j] << ", window size = " << wsvec[j] << ")...\n"; - std::vector me = - eigen::erodeDiamond(ZImin, rows, cols, 0.5 * (wsvec[j] - 1)); - std::vector mo = - eigen::dilateDiamond(me, rows, cols, 0.5 * (wsvec[j] - 1)); + int iters = 0.5 * (wsvec[j] - 1); + using namespace eigen; + std::vector me = erodeDiamond(ZImin, rows, cols, iters); + std::vector mo = dilateDiamond(me, rows, cols, iters); std::vector groundNewIdx; for (auto p_idx : groundIdx) @@ -365,73 +281,13 @@ std::vector PMFFilter::processGroundApprox(PointViewPtr view) << "Ground now has " << groundIdx.size() << " points.\n"; } - return groundIdx; -} - -PointViewSet PMFFilter::run(PointViewPtr input) -{ - PointViewSet viewSet; - if (!input->size()) - return viewSet; - - // Segment input view into ignored/kept views. - PointViewPtr ignoredView = input->makeNew(); - PointViewPtr keptView = input->makeNew(); - if (m_ignored.m_id == Dimension::Id::Unknown) - keptView->append(*input); - else - Segmentation::ignoreDimRange(m_ignored, input, keptView, ignoredView); - - // Segment kept view into last/other-than-last return views. - PointViewPtr lastView = keptView->makeNew(); - PointViewPtr nonlastView = keptView->makeNew(); - if (m_lastOnly) - Segmentation::segmentLastReturns(keptView, lastView, nonlastView); - else - lastView->append(*keptView); - - for (PointId i = 0; i < nonlastView->size(); ++i) - nonlastView->setField(Dimension::Id::Classification, i, 1); - - for (PointId i = 0; i < lastView->size(); ++i) - lastView->setField(Dimension::Id::Classification, i, 1); - - std::vector idx; - if (m_approximate) - idx = processGroundApprox(lastView); - else - idx = processGround(lastView); + log()->get(LogLevel::Debug2) + << "Labeled " << groundIdx.size() << " ground returns!\n"; - PointViewPtr outView = input->makeNew(); - if (!idx.empty()) - { - - log()->get(LogLevel::Debug2) - << "Labeled " << idx.size() << " ground returns!\n"; - - // set the classification label of ground returns as 2 - // (corresponding to ASPRS LAS specification) - for (const auto& i : idx) - { - lastView->setField(Dimension::Id::Classification, i, 2); - } - - outView->append(*ignoredView); - outView->append(*nonlastView); - outView->append(*lastView); - } - else - { - if (idx.empty()) - log()->get(LogLevel::Debug2) << "Filtered cloud has no " - "ground returns!\n"; - - // return the input buffer unchanged - outView->append(*input); - } - viewSet.insert(outView); - - return viewSet; + // set the classification label of ground returns as 2 + // (corresponding to ASPRS LAS specification) + for (const auto& i : groundIdx) + view->setField(Dimension::Id::Classification, i, 2); } } // namespace pdal diff --git a/filters/PMFFilter.hpp b/filters/PMFFilter.hpp index d2655d4568..d41b1eb052 100644 --- a/filters/PMFFilter.hpp +++ b/filters/PMFFilter.hpp @@ -39,19 +39,12 @@ #include "private/DimRange.hpp" -#include - extern "C" int32_t PMFFilter_ExitFunc(); extern "C" PF_ExitFunc PMFFilter_InitPlugin(); namespace pdal { -class Options; -class PointLayout; -class PointTable; -class PointView; - class PDAL_DLL PMFFilter : public Filter { public: @@ -64,25 +57,22 @@ class PDAL_DLL PMFFilter : public Filter std::string getName() const; private: - double m_maxWindowSize; - double m_slope; - double m_maxDistance; - double m_initialDistance; double m_cellSize; - bool m_approximate; + bool m_exponential; DimRange m_ignored; + double m_initialDistance; bool m_lastOnly; + double m_maxDistance; + double m_maxWindowSize; + double m_slope; virtual void addDimensions(PointLayoutPtr layout); virtual void addArgs(ProgramArgs& args); - std::vector fillNearest(PointViewPtr view, size_t rows, size_t cols, - double cell_size, BOX2D bounds); - std::vector morphOpen(PointViewPtr view, float radius); virtual void prepared(PointTableRef table); - std::vector processGround(PointViewPtr view); - std::vector processGroundApprox(PointViewPtr view); virtual PointViewSet run(PointViewPtr view); + void processGround(PointViewPtr view); + PMFFilter& operator=(const PMFFilter&); // not implemented PMFFilter(const PMFFilter&); // not implemented }; diff --git a/io/GDALReader.hpp b/io/GDALReader.hpp index d811a183ec..e625e85ffb 100644 --- a/io/GDALReader.hpp +++ b/io/GDALReader.hpp @@ -64,8 +64,6 @@ class PDAL_DLL GDALReader : public Reader GDALReader(); - static Dimension::IdList getDefaultDimensions(); - private: virtual void initialize(); virtual void addDimensions(PointLayoutPtr layout); diff --git a/io/GeotiffSupport.cpp b/io/GeotiffSupport.cpp index bf117d3828..7650547a21 100644 --- a/io/GeotiffSupport.cpp +++ b/io/GeotiffSupport.cpp @@ -103,17 +103,25 @@ GeotiffSrs::GeotiffSrs(const std::vector& directoryRec, size_t declaredSize = (header->numKeys + 1) * 4; if (directoryRec.size() < declaredSize) return; + + uint8_t *dirData = const_cast(directoryRec.data()); ST_SetKey(ctx.tiff, GEOTIFF_DIRECTORY_RECORD_ID, - (1 + header->numKeys) * 4, STT_SHORT, (void *)directoryRec.data()); + (1 + header->numKeys) * 4, STT_SHORT, (void *)dirData); if (doublesRec.size()) + { + uint8_t *doubleData = const_cast(doublesRec.data()); ST_SetKey(ctx.tiff, GEOTIFF_DOUBLES_RECORD_ID, doublesRec.size() / sizeof(double), STT_DOUBLE, - (void *)doublesRec.data()); + (void *)doubleData); + } if (asciiRec.size()) + { + uint8_t *asciiData = const_cast(asciiRec.data()); ST_SetKey(ctx.tiff, GEOTIFF_ASCII_RECORD_ID, - asciiRec.size(), STT_ASCII, (void *)asciiRec.data()); + asciiRec.size(), STT_ASCII, (void *)asciiData); + } ctx.gtiff = GTIFNewSimpleTags(ctx.tiff); diff --git a/io/Ilvis2Reader.cpp b/io/Ilvis2Reader.cpp index 175ba8f6a9..19dd164c7a 100644 --- a/io/Ilvis2Reader.cpp +++ b/io/Ilvis2Reader.cpp @@ -115,19 +115,6 @@ void Ilvis2Reader::addDimensions(PointLayoutPtr layout) } -Dimension::IdList Ilvis2Reader::getDefaultDimensions() -{ - using namespace pdal::Dimension; - Dimension::IdList ids; - - ids.push_back(Id::GpsTime); - ids.push_back(Id::Y); - ids.push_back(Id::X); - ids.push_back(Id::Z); - return ids; -} - - void Ilvis2Reader::initialize(PointTableRef) { if (!m_metadataFile.empty() && !FileUtils::fileExists(m_metadataFile)) diff --git a/io/Ilvis2Reader.hpp b/io/Ilvis2Reader.hpp index cf65be665f..ee0f64e5f0 100644 --- a/io/Ilvis2Reader.hpp +++ b/io/Ilvis2Reader.hpp @@ -80,8 +80,6 @@ class PDAL_DLL Ilvis2Reader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions(); - private: std::ifstream m_stream; IlvisMapping m_mapping; diff --git a/io/OptechReader.cpp b/io/OptechReader.cpp index 8c8fd7aa4e..411ced7ed7 100644 --- a/io/OptechReader.cpp +++ b/io/OptechReader.cpp @@ -82,22 +82,6 @@ OptechReader::OptechReader() } -Dimension::IdList OptechReader::getDefaultDimensions() -{ - Dimension::IdList dims; - dims.push_back(Dimension::Id::X); - dims.push_back(Dimension::Id::Y); - dims.push_back(Dimension::Id::Z); - dims.push_back(Dimension::Id::GpsTime); - dims.push_back(Dimension::Id::ReturnNumber); - dims.push_back(Dimension::Id::NumberOfReturns); - dims.push_back(Dimension::Id::EchoRange); - dims.push_back(Dimension::Id::Intensity); - dims.push_back(Dimension::Id::ScanAngleRank); - return dims; -} - - const CsdHeader& OptechReader::getHeader() const { return m_header; } @@ -136,10 +120,11 @@ void OptechReader::initialize() void OptechReader::addDimensions(PointLayoutPtr layout) { - for (auto it : getDefaultDimensions()) - { - layout->registerDim(it); - } + using namespace Dimension; + + layout->registerDims( { Id::X, Id::Y, Id::Z, Id::GpsTime, Id::ReturnNumber, + Id::NumberOfReturns, Id::EchoRange, Id::Intensity, + Id::ScanAngleRank } ); } diff --git a/io/OptechReader.hpp b/io/OptechReader.hpp index 38399c46ed..6cdfde1dc0 100644 --- a/io/OptechReader.hpp +++ b/io/OptechReader.hpp @@ -66,7 +66,6 @@ class PDAL_DLL OptechReader : public Reader OptechReader(); - static Dimension::IdList getDefaultDimensions(); const CsdHeader& getHeader() const; private: diff --git a/io/PlyReader.cpp b/io/PlyReader.cpp index b8a0867bb5..b1095c3af1 100644 --- a/io/PlyReader.cpp +++ b/io/PlyReader.cpp @@ -267,6 +267,7 @@ bool PlyReader::extractElement() void PlyReader::extractHeader() { + m_elements.clear(); extractMagic(); extractFormat(); while (extractElement()) diff --git a/io/QfitReader.cpp b/io/QfitReader.cpp index 99a863d096..d93e3301a0 100644 --- a/io/QfitReader.cpp +++ b/io/QfitReader.cpp @@ -401,30 +401,6 @@ point_count_t QfitReader::read(PointViewPtr data, point_count_t count) } -Dimension::IdList QfitReader::getDefaultDimensions() -{ - Dimension::IdList ids; - - ids.push_back(Dimension::Id::OffsetTime); - ids.push_back(Dimension::Id::Y); - ids.push_back(Dimension::Id::X); - ids.push_back(Dimension::Id::Z); - ids.push_back(Dimension::Id::StartPulse); - ids.push_back(Dimension::Id::ReflectedPulse); - ids.push_back(Dimension::Id::Azimuth); - ids.push_back(Dimension::Id::Pitch); - ids.push_back(Dimension::Id::Roll); - ids.push_back(Dimension::Id::Pdop); - ids.push_back(Dimension::Id::PulseWidth); - ids.push_back(Dimension::Id::PassiveSignal); - ids.push_back(Dimension::Id::PassiveY); - ids.push_back(Dimension::Id::PassiveX); - ids.push_back(Dimension::Id::PassiveZ); - - return ids; -} - - void QfitReader::done(PointTableRef) { m_istream.reset(); diff --git a/io/QfitReader.hpp b/io/QfitReader.hpp index ea31bdd55e..7b158c1597 100644 --- a/io/QfitReader.hpp +++ b/io/QfitReader.hpp @@ -66,8 +66,6 @@ class PDAL_DLL QfitReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions(); - private: QFIT_Format_Type m_format; std::ios::off_type m_point_bytes; diff --git a/io/SbetReader.cpp b/io/SbetReader.cpp index 44ada0db34..6eed8233ec 100644 --- a/io/SbetReader.cpp +++ b/io/SbetReader.cpp @@ -53,14 +53,14 @@ std::string SbetReader::getName() const { return s_info.name; } void SbetReader::addDimensions(PointLayoutPtr layout) { - layout->registerDims(getDefaultDimensions()); + layout->registerDims(fileDimensions()); } void SbetReader::ready(PointTableRef) { size_t fileSize = FileUtils::fileSize(m_filename); - size_t pointSize = getDefaultDimensions().size() * sizeof(double); + size_t pointSize = fileDimensions().size() * sizeof(double); if (fileSize % pointSize != 0) throwError("Invalid file size."); m_numPts = fileSize / pointSize; @@ -114,7 +114,7 @@ bool SbetReader::eof() void SbetReader::seek(PointId idx) { - m_stream->seek(idx * sizeof(double) * getDefaultDimensions().size()); + m_stream->seek(idx * sizeof(double) * fileDimensions().size()); } } // namespace pdal diff --git a/io/SbetReader.hpp b/io/SbetReader.hpp index 85ed2893d4..63526f4be8 100644 --- a/io/SbetReader.hpp +++ b/io/SbetReader.hpp @@ -57,9 +57,6 @@ class PDAL_DLL SbetReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions() - { return fileDimensions(); } - private: std::unique_ptr m_stream; // Number of points in the file. diff --git a/io/SbetWriter.cpp b/io/SbetWriter.cpp index 81f84a4644..9ed7367145 100644 --- a/io/SbetWriter.cpp +++ b/io/SbetWriter.cpp @@ -64,7 +64,7 @@ void SbetWriter::ready(PointTableRef) void SbetWriter::write(const PointViewPtr view) { - Dimension::IdList dims = getDefaultDimensions(); + Dimension::IdList dims = fileDimensions(); for (PointId idx = 0; idx < view->size(); ++idx) { for (auto di = dims.begin(); di != dims.end(); ++di) diff --git a/io/SbetWriter.hpp b/io/SbetWriter.hpp index 262209a3d4..d6ca9bf83d 100644 --- a/io/SbetWriter.hpp +++ b/io/SbetWriter.hpp @@ -53,9 +53,6 @@ class PDAL_DLL SbetWriter : public Writer static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions() - { return fileDimensions(); } - private: std::unique_ptr m_stream; std::string m_filename; diff --git a/io/TIndexReader.cpp b/io/TIndexReader.cpp index 10da1d525c..affde7a0cd 100644 --- a/io/TIndexReader.cpp +++ b/io/TIndexReader.cpp @@ -133,12 +133,6 @@ void TIndexReader::addDimensions(PointLayoutPtr layout) } -pdal::Dimension::IdList TIndexReader::getDefaultDimensions() -{ - return Dimension::IdList(); -} - - void TIndexReader::initialize() { if (!m_bounds.empty()) diff --git a/io/TIndexReader.hpp b/io/TIndexReader.hpp index 09d88ec7ab..20772b1509 100644 --- a/io/TIndexReader.hpp +++ b/io/TIndexReader.hpp @@ -73,8 +73,6 @@ class PDAL_DLL TIndexReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions(); - private: virtual void addDimensions(PointLayoutPtr layout); virtual void addArgs(ProgramArgs& args); diff --git a/io/TerrasolidReader.cpp b/io/TerrasolidReader.cpp index 989b1aebc7..35f99d4979 100644 --- a/io/TerrasolidReader.cpp +++ b/io/TerrasolidReader.cpp @@ -129,31 +129,6 @@ void TerrasolidReader::addDimensions(PointLayoutPtr layout) } -Dimension::IdList TerrasolidReader::getDefaultDimensions() -{ - using namespace Dimension; - - IdList dims; - - dims.push_back(Id::Classification); - dims.push_back(Id::PointSourceId); - dims.push_back(Id::ReturnNumber); - dims.push_back(Id::NumberOfReturns); - dims.push_back(Id::Flag); - dims.push_back(Id::Mark); - dims.push_back(Id::Intensity); - dims.push_back(Id::X); - dims.push_back(Id::Y); - dims.push_back(Id::Z); - dims.push_back(Id::Red); - dims.push_back(Id::Green); - dims.push_back(Id::Blue); - dims.push_back(Id::Alpha); - dims.push_back(Id::OffsetTime); - return dims; -} - - void TerrasolidReader::ready(PointTableRef) { m_istream.reset(new IStream(m_filename)); diff --git a/io/TerrasolidReader.hpp b/io/TerrasolidReader.hpp index 796029a207..5c45095735 100644 --- a/io/TerrasolidReader.hpp +++ b/io/TerrasolidReader.hpp @@ -96,8 +96,6 @@ class PDAL_DLL TerrasolidReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions(); - point_count_t getNumPoints() const { return m_header->PntCnt; } diff --git a/pdal/Dimension.json b/pdal/Dimension.json index 80b51bc062..2a3f1ff717 100644 --- a/pdal/Dimension.json +++ b/pdal/Dimension.json @@ -340,16 +340,19 @@ }, { "name": "NormalX", + "alt_names": "nx", "type": "double", "description": "X component of a vector normal to surface at this point" }, { "name": "NormalY", + "alt_names": "ny", "type": "double", "description": "Y component of a vector normal to surface at this point" }, { "name": "NormalZ", + "alt_names": "nz", "type": "double", "description": "Z component of a vector normal to surface at this point" }, @@ -357,5 +360,11 @@ "name": "Curvature", "type": "double", "description": "Curvature of surface at this point" + }, + { + "name": "Density", + "type": "double", + "description": "Estimate of point density" } + ] } diff --git a/pdal/GDALUtils.cpp b/pdal/GDALUtils.cpp index d7ca4116d0..dbe88800b0 100644 --- a/pdal/GDALUtils.cpp +++ b/pdal/GDALUtils.cpp @@ -146,6 +146,24 @@ bool reprojectBounds(BOX3D& box, const std::string& srcSrs, } +/** + Reproject a bounds box from a source projection to a destination. + \param box 2D Bounds box to be reprojected in-place. + \param srcSrs String in WKT or other suitable format of box coordinates. + \param dstSrs String in WKT or other suitable format to which + coordinates should be projected. + \return Whether the reprojection was successful or not. +*/ +bool reprojectBounds(BOX2D& box, const std::string& srcSrs, + const std::string& dstSrs) +{ + BOX3D b(box); + bool res = reprojectBounds(b, srcSrs, dstSrs); + box = b.to2d(); + return res; +} + + /** Reproject a point from a source projection to a destination. \param x X coordinate of point to be reprojected. diff --git a/pdal/GDALUtils.hpp b/pdal/GDALUtils.hpp index b0209ad211..119c20052f 100644 --- a/pdal/GDALUtils.hpp +++ b/pdal/GDALUtils.hpp @@ -65,6 +65,8 @@ PDAL_DLL void registerDrivers(); PDAL_DLL void unregisterDrivers(); PDAL_DLL bool reprojectBounds(BOX3D& box, const std::string& srcSrs, const std::string& dstSrs); +PDAL_DLL bool reprojectBounds(BOX2D& box, const std::string& srcSrs, + const std::string& dstSrs); PDAL_DLL bool reprojectPoint(double& x, double& y, double& z, const std::string& srcSrs, const std::string& dstSrs); PDAL_DLL std::string lastError(); diff --git a/pdal/StageFactory.cpp b/pdal/StageFactory.cpp index cac9c9747e..6d0ee39cf7 100644 --- a/pdal/StageFactory.cpp +++ b/pdal/StageFactory.cpp @@ -66,6 +66,7 @@ #include #include #include +#include #include #include #include @@ -276,6 +277,7 @@ StageFactory::StageFactory(bool no_plugins) PluginManager::initializePlugin(OutlierFilter_InitPlugin); PluginManager::initializePlugin(OverlayFilter_InitPlugin); PluginManager::initializePlugin(PMFFilter_InitPlugin); + PluginManager::initializePlugin(PoissonFilter_InitPlugin); PluginManager::initializePlugin(RadialDensityFilter_InitPlugin); PluginManager::initializePlugin(RandomizeFilter_InitPlugin); PluginManager::initializePlugin(RangeFilter_InitPlugin); diff --git a/pdal/util/Bounds.hpp b/pdal/util/Bounds.hpp index 45d93ca48b..2438a062d7 100644 --- a/pdal/util/Bounds.hpp +++ b/pdal/util/Bounds.hpp @@ -566,8 +566,8 @@ class PDAL_DLL Bounds Bounds() {} - Bounds(const BOX3D& box); - Bounds(const BOX2D& box); + explicit Bounds(const BOX3D& box); + explicit Bounds(const BOX2D& box); BOX3D to3d() const; BOX2D to2d() const; diff --git a/plugins/geowave/io/GeoWaveReader.cpp b/plugins/geowave/io/GeoWaveReader.cpp index d269ca78d4..1f2d7b03b9 100644 --- a/plugins/geowave/io/GeoWaveReader.cpp +++ b/plugins/geowave/io/GeoWaveReader.cpp @@ -191,7 +191,7 @@ std::string pdal::GeoWaveReader::getName() const { return s_info.name; } void GeoWaveReader::addDimensions(PointLayoutPtr layout) { - layout->registerDims(getDefaultDimensions()); + layout->registerDims( { Dimension::Id::X, Dimension::Id::Y } ); BasicAccumuloOperations accumuloOperations; try @@ -229,14 +229,6 @@ std::string pdal::GeoWaveReader::getName() const { return s_info.name; } } } - Dimension::IdList GeoWaveReader::getDefaultDimensions() - { - Dimension::IdList ids; - ids.push_back(Dimension::Id::X); - ids.push_back(Dimension::Id::Y); - return ids; - } - void GeoWaveReader::ready(PointTableRef table) { if (m_bounds.empty()) diff --git a/plugins/geowave/io/GeoWaveReader.hpp b/plugins/geowave/io/GeoWaveReader.hpp index 9b1bdadbfd..22d185c504 100644 --- a/plugins/geowave/io/GeoWaveReader.hpp +++ b/plugins/geowave/io/GeoWaveReader.hpp @@ -65,8 +65,6 @@ namespace pdal virtual void ready(PointTableRef table); virtual point_count_t read(PointViewPtr view, point_count_t count); virtual void done(PointTableRef table); - - Dimension::IdList getDefaultDimensions(); int createJvm(); std::string m_zookeeperUrl; diff --git a/plugins/icebridge/io/IcebridgeReader.cpp b/plugins/icebridge/io/IcebridgeReader.cpp index 454d8b2c92..2bd9d59f59 100644 --- a/plugins/icebridge/io/IcebridgeReader.cpp +++ b/plugins/icebridge/io/IcebridgeReader.cpp @@ -73,31 +73,23 @@ CREATE_SHARED_PLUGIN(1, 0, IcebridgeReader, Reader, s_info) std::string IcebridgeReader::getName() const { return s_info.name; } -Dimension::IdList IcebridgeReader::getDefaultDimensions() +namespace { - Dimension::IdList ids; +Dimension::IdList dimensions() +{ using namespace Dimension; - ids.push_back(Id::OffsetTime); - ids.push_back(Id::Y); - ids.push_back(Id::X); - ids.push_back(Id::Z); - ids.push_back(Id::StartPulse); - ids.push_back(Id::ReflectedPulse); - ids.push_back(Id::Azimuth); - ids.push_back(Id::Pitch); - ids.push_back(Id::Roll); - ids.push_back(Id::Pdop); - ids.push_back(Id::PulseWidth); - ids.push_back(Id::GpsTime); - return ids; + return IdList { Id::OffsetTime, Id::Y, Id::X, Id::Z, Id::StartPulse, + Id::ReflectedPulse, Id::Azimuth, Id::Pitch, Id::Roll, Id::Pdop, + Id::PulseWidth, Id::GpsTime }; } +} // unnamed namespace void IcebridgeReader::addDimensions(PointLayoutPtr layout) { - return layout->registerDims(getDefaultDimensions()); + layout->registerDims(dimensions()); } @@ -134,7 +126,7 @@ point_count_t IcebridgeReader::read(PointViewPtr view, point_count_t count) rawData(new unsigned char[count * sizeof(float)]); //Not loving the position-linked data, but fine for now. - Dimension::IdList dims = getDefaultDimensions(); + Dimension::IdList dims = dimensions(); auto di = dims.begin(); for (auto ci = hdf5Columns.begin(); ci != hdf5Columns.end(); ++ci, ++di) { diff --git a/plugins/icebridge/io/IcebridgeReader.hpp b/plugins/icebridge/io/IcebridgeReader.hpp index d19d50910d..dddc9db312 100644 --- a/plugins/icebridge/io/IcebridgeReader.hpp +++ b/plugins/icebridge/io/IcebridgeReader.hpp @@ -69,8 +69,6 @@ class PDAL_DLL IcebridgeReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions(); - private: Hdf5Handler m_hdf5Handler; point_count_t m_index; diff --git a/plugins/matlab/io/MatlabReader.cpp b/plugins/matlab/io/MatlabReader.cpp index faf28346b1..641196ce57 100644 --- a/plugins/matlab/io/MatlabReader.cpp +++ b/plugins/matlab/io/MatlabReader.cpp @@ -163,7 +163,7 @@ bool MatlabReader::processOne(PointRef& point) if (numElements >= m_pointIndex ) { size_t size = mxGetElementSize(f); - char* p = (char*)mxGetData(f) + size; + char* p = (char*)mxGetData(f) + (m_pointIndex*size); point.setField(d, t, (void*)p); } diff --git a/plugins/pcl/io/PcdReader.cpp b/plugins/pcl/io/PcdReader.cpp index 0f20604d6f..aa1f4e486b 100644 --- a/plugins/pcl/io/PcdReader.cpp +++ b/plugins/pcl/io/PcdReader.cpp @@ -68,7 +68,7 @@ void PcdReader::ready(PointTableRef table) void PcdReader::addDimensions(PointLayoutPtr layout) { - layout->registerDims(getDefaultDimensions()); + layout->registerDims(fileDimensions()); } diff --git a/plugins/pcl/io/PcdReader.hpp b/plugins/pcl/io/PcdReader.hpp index c824bb932d..7d9ebe7cd9 100644 --- a/plugins/pcl/io/PcdReader.hpp +++ b/plugins/pcl/io/PcdReader.hpp @@ -52,11 +52,6 @@ class PDAL_DLL PcdReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions() - { - return fileDimensions(); - }; - private: point_count_t m_numPts; diff --git a/plugins/rxp/io/RxpReader.hpp b/plugins/rxp/io/RxpReader.hpp index ae79caebc8..86258e5b2c 100644 --- a/plugins/rxp/io/RxpReader.hpp +++ b/plugins/rxp/io/RxpReader.hpp @@ -80,11 +80,6 @@ class PDAL_DLL RxpReader : public pdal::Reader static int32_t destroy(void *); std::string getName() const; - static Dimension::IdList getDefaultDimensions() - { - return getRxpDimensions(DEFAULT_SYNC_TO_PPS, DEFAULT_MINIMAL); - } - private: virtual void addArgs(ProgramArgs& args); virtual void initialize(); diff --git a/scripts/ci/script.sh b/scripts/ci/script.sh index 63e1f11c37..454ee475f6 100755 --- a/scripts/ci/script.sh +++ b/scripts/ci/script.sh @@ -68,7 +68,7 @@ if [ "${OPTIONAL_COMPONENT_SWITCH}" == "ON" ]; then python setup.py test # JNI tests - cd /pdal/java; PDAL_DEPEND_ON_NATIVE=false ./sbt -Djava.library.path=/pdal/_build/lib core/test +# cd /pdal/java; PDAL_DEPEND_ON_NATIVE=false ./sbt -Djava.library.path=/pdal/_build/lib core/test # Build all examples for EXAMPLE in writing writing-filter writing-kernel writing-reader writing-writer diff --git a/scripts/docker/Dockerfile b/scripts/docker/Dockerfile index 6767b04ad8..9b4407827e 100644 --- a/scripts/docker/Dockerfile +++ b/scripts/docker/Dockerfile @@ -11,7 +11,7 @@ RUN apt-get update && apt-get install -y --fix-missing --no-install-recommends \ libsqlite3-mod-spatialite \ && rm -rf /var/lib/apt/lists/* -RUN git clone --depth=1 https://github.com/PDAL/PDAL \ +RUN git clone https://github.com/PDAL/PDAL \ && cd PDAL \ && git checkout master \ && mkdir build \ diff --git a/test/unit/OldPCLBlockTest.cpp b/test/unit/OldPCLBlockTest.cpp index 4a284f283b..47231ba45c 100644 --- a/test/unit/OldPCLBlockTest.cpp +++ b/test/unit/OldPCLBlockTest.cpp @@ -193,101 +193,7 @@ TEST(OldPCLBlockTests, RadiusOutliers2) EXPECT_EQ(3u, view->size()); } -TEST(OldPCLBlockTests, PMF1) -{ - StageFactory f; - - Options ro; - ro.add("filename", Support::datapath("autzen/autzen-point-format-3.las")); - - Stage* r(f.createStage("readers.las")); - EXPECT_TRUE(r); - r->setOptions(ro); - - Options ao; - ao.add("assignment", "Classification[:]=0"); - - Stage* assign(f.createStage("filters.assign")); - EXPECT_TRUE(assign); - assign->setOptions(ao); - assign->setInput(*r); - - Options fo; - fo.add("max_window_size", 200); - fo.add("last", false); - - Stage* outlier(f.createStage("filters.pmf")); - EXPECT_TRUE(outlier); - outlier->setOptions(fo); - outlier->setInput(*assign); - - Options rangeOpts; - rangeOpts.add("limits", "Classification[2:2]"); - - Stage* range(f.createStage("filters.range")); - EXPECT_TRUE(range); - range->setOptions(rangeOpts); - range->setInput(*outlier); - - PointTable table; - range->prepare(table); - PointViewSet viewSet = range->execute(table); - - EXPECT_EQ(1u, viewSet.size()); - PointViewPtr view = *viewSet.begin(); - EXPECT_EQ(93u, view->size()); -} - -TEST(OldPCLBlockTests, PMF2) -{ - StageFactory f; - - Options ro; - ro.add("filename", Support::datapath("autzen/autzen-point-format-3.las")); - - Stage* r(f.createStage("readers.las")); - EXPECT_TRUE(r); - r->setOptions(ro); - - Options ao; - ao.add("assignment", "Classification[:]=0"); - - Stage* assign(f.createStage("filters.assign")); - EXPECT_TRUE(assign); - assign->setOptions(ao); - assign->setInput(*r); - - Options fo; - fo.add("max_window_size", 200); - fo.add("cell_size", 1.0); - fo.add("slope", 1.0); - fo.add("initial_distance", 0.05); - fo.add("max_distance", 3.0); - fo.add("last", false); - - Stage* outlier(f.createStage("filters.pmf")); - EXPECT_TRUE(outlier); - outlier->setOptions(fo); - outlier->setInput(*assign); - - Options rangeOpts; - rangeOpts.add("limits", "Classification[2:2]"); - - Stage* range(f.createStage("filters.range")); - EXPECT_TRUE(range); - range->setOptions(rangeOpts); - range->setInput(*outlier); - - PointTable table; - range->prepare(table); - PointViewSet viewSet = range->execute(table); - - EXPECT_EQ(1u, viewSet.size()); - PointViewPtr view = *viewSet.begin(); - EXPECT_EQ(94u, view->size()); -} - -TEST(OldPCLBlockTests, PMF3) +TEST(OldPCLBlockTests, PMF) { StageFactory f; @@ -335,38 +241,3 @@ TEST(OldPCLBlockTests, PMF3) PointViewPtr view = *viewSet.begin(); EXPECT_EQ(106u, view->size()); } - -/* -// Test these with "autzen/autzen-thin.las" -TEST(PCLBlockFilterTest, PCLBlockFilterTest_filter_PMF) -{ -#if RUN_SLOW_TESTS - // explicitly with all defaults - test_filter("filters/pcl/filter_PMF_1.json", 9223, true); - - // with CellSize=3 - test_filter("filters/pcl/filter_PMF_2.json", 8298, true); - - // with WindowSize=50 - test_filter("filters/pcl/filter_PMF_3.json", 7970, true); - - // with Slope=0.25 - test_filter("filters/pcl/filter_PMF_4.json", 9206, true); - - // with MaxDistance=5 - test_filter("filters/pcl/filter_PMF_5.json", 9373, true); - - // with InitialDistance=0.25 - test_filter("filters/pcl/filter_PMF_6.json", 9229, true); - - // with Base=3 - test_filter("filters/pcl/filter_PMF_7.json", 8298, true); - - // with Exponential=false - test_filter("filters/pcl/filter_PMF_8.json", 9138, true); - - // with Negative=true - test_filter("filters/pcl/filter_PMF_9.json", 1430, true); -#endif -} -*/