Dense photogrammetric point cloud classification

611 views
Skip to first unread message

Muriel Lavy

unread,
Jun 20, 2017, 3:30:40 AM6/20/17
to LAStools - efficient tools for LiDAR processing
Dear all,

I got an aerial images dataset, collected from an airplane flight (2500 ground ft) over a flooded area. The dataset was acquired with a Panasonic DMC-GX7 camera (695 images)
I obtained a LAS file of the point cloud with Agisoft Photoscan (which looks quite good but not great)

I would like to generate a good lastools processing pipeline for the classification of the photogrammetric point cloud which could be used in future for different datasets with a minimum change of the control parameters.
I aim to obtain a proper classification of the las file, and then get a good DTM and DSM.

The point cloud has a lot of noise points and outlines (above and below terrain) generated from the Agisoft process, I would like to automatically remove those point with no manual interaction.

I took some ideas from several "photogrammetric" posts in Rapidlasso blog and the lastools google group and I prepared the script attached. The workflow was good on a little sample part of the area (following a step by step process) but I do not get a good result on the full dataset with the scripting.

Do you have any suggestion to properly modify the script in order to obtain a good classification of the point cloud?
Here the script:


:: specify input file
set INPUT_POINTS=C:\LAS\20161127_Pancalieri_UTM.las

:: specify file name of output classified point cloud
set OUTPUT_POINTS=C:\LAS\20161127_Pancalieri_classified.las

:: specify temporary directory for tiles
set TEMP_DIR=C:\LAS\temp

:: specify size of tile and edge-artifact-avoiding buffer
set TILE_SIZE=500
set BUFFER=25

:: specify number of cores
set CORES=11

:: create or empty temporary directory
rmdir %TEMP_DIR% /s /q
mkdir %TEMP_DIR%


:: create a compressed ZIP las
laszip  -i %INPUT_POINTS% ^
-rescale 0.01 0.01 0.01 ^
-olaz

:: create buffered tiling
lastile -i *.laz ^
        -tile_size %TILE_SIZE% -buffer %BUFFER% ^
        -o %TEMP_DIR%\tile.laz -olaz

:: run lassort on all tiles
lassort -i %TEMP_DIR%\tile*.laz ^
        -odix _s -olaz ^
        -cores %CORES%

:: run lasnoise on all tiles
lasnoise    -i %TEMP_DIR%\tile*_s.laz ^
-step_xy 4 ^
-step_z 1^
-isolated 15 ^
-odix n -olaz^
-cores %CORES%

:: run lasground on all tiles
lasground -i %TEMP_DIR%\tile*_sn.laz ^
              -all_returns ^
              -step 20 -bulge 0.5 -spike 0.5 -extra_fine ^
      -compute_height ^
              -odix g -olaz ^
              -cores %CORES%
 
:: run lascanopy + lasheight to remove point below terrain and above sky
lasheight -i %TEMP_DIR%\tile*_sng.laz ^
              -replace_z ^
              -odix h -olaz ^
      -cores %CORES%

lascanopy -i %TEMP_DIR%\tile*_sngh.laz ^
          -height_cutoff -1000 -step 5 ^
          -p 10 ^
          -obil
   
lasheight -i %TEMP_DIR%\tile*_sng.laz ^
          -ground_points *p10.bil ^
          -classify_below -1 7 ^
          -odix d -olaz ^
 -cores %CORES%
 
:: run lasclassify on all tiles
lasclassify -i %TEMP_DIR%\tile*_sngd.laz ^
-odix c -olaz ^
-cores %CORES%
:: run lastile to remove buffers on all tiles
lastile -i %TEMP_DIR%\tile*_sngdc.laz ^
-ignore_class 7 ^
-remove_buffer ^
-odix buf -olaz ^
-cores %CORES%
:: merge classified tiles into a single output
lasmerge -i %TEMP_DIR%\tile*_sngdcbuf.laz ^
-o %OUTPUT_POINTS%

:: remove the temporary files and directory
rmdir %TEMP_DIR% /s /q


here the link of the dataset:


Thanks for your Help :)

Susana Gonzalez

unread,
Jun 21, 2017, 2:00:16 AM6/21/17
to LAStools - efficient tools for LiDAR processing

Hi Muriel,

Could I have access to the raw data, the original images? And I would try to run it in Pix4D to see if the amount of noise is less.

Cheers

Susana


From: last...@googlegroups.com <last...@googlegroups.com> on behalf of Muriel Lavy <muri...@gmail.com>
Sent: Tuesday, 20 June 2017 7:29:04 p.m.
To: LAStools - efficient tools for LiDAR processing
Subject: [LAStools] Dense photogrammetric point cloud classification
 

Muriel Lavy

unread,
Jun 21, 2017, 10:23:07 AM6/21/17
to LAStools - efficient tools for LiDAR processing, susana....@interpine.nz
Thanks Susana,

here the dataset with the original images:

Tobias K Kohoutek

unread,
Jun 21, 2017, 6:58:14 PM6/21/17
to LAStools - efficient tools for LiDAR processing
Dear Muriel,

I'll paste my pipeline here. Of course also not the best and there's always space to improve. I do some things differently than you do, so check it out. (comments are in Spanish but I guess you'll figure)


@echo off
cls
echo #################################################
echo "batch para procesar NDP generado en base de imagenes"
echo "                        by Dr. Tobias K. Kohoutek                     "
echo "                                   20.06.2017                                "
echo #################################################


REM incluir LAStools en Path
set PATH=%PATH%;..;C:\lastools\bin;


REM ===================================================
REM = DESDE AQUI UNO TIENE QUE CAMBIAR LOS PARAMETROS =
REM ===================================================

REM la ruta de la NDP
set NDP_IMG=D:\...\NDP

REM el formato de los archivos originales
set NDP_IMG_FORMAT=laz

REM especificar la projeccion UTM
set PROJECTION=utm ...south

REM la ruta de los archivos temporarios
set TEMP_FILES=D:\...\temp

REM especificar la ruta de los resultados
set OUTPUT_FILES=D:\...\NDP_clasificada

REM especificar altura (m) por clasificacion ruido alto (clasa 18)
set ABOVE=100

REM especificar altura (m) por clasificacion ruido bajo (clasa 7)
set BELOW=-3

REM especificar los nucleos de computador para usar
set NUM_CORES
set /A NUM_CORES=%NUMBER_OF_PROCESSORS% - 1

REM especificar el tamaño de los tiles
set TILE_SIZE=1000

REM especificar el buffer de los tiles
set TILE_BUFFER=50

REM especificar la version LAS
set LAS_V=1.2

REM especificar hasta cual altura relativa de suele objectos se clasifican como suelo
set GROUND=0.3

REM especificar hasta cual altura relativa de suele objectos se clasifican como vegetacion baja
set VEG_L=1.3

REM especificar hasta cual altura relativa de suele objectos se clasifican como vegetacion media
set VEG_M=4

REM ===============================================
REM = DESDE AQUI TODO ESTA LISTO PARA EL PROCESSO =
REM ===============================================

REM eliminar carpeta temporales
rmdir %TEMP_FILES% /s /q
REM crear carpeta temporaria y resultada
mkdir %TEMP_FILES% 
mkdir %OUTPUT_FILES%

echo START Automatic Pointcloud Classification %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt


REM copiar archivos de lineas de vuelo con las2las desde odin a ruta local para procesar
REM Todos los puntos se clasifican en clasa 1 '-set_classification 1'. La version LAS se 
REM cambia con '-set_version %LAS_V%'. La  proyeccion '-%PROJECTION%'. Los coordinadas 
REM se escaladan a una precision de centimetro '-rescale 0.01 0.01 0.01' (see: las2las_README.txt)
mkdir %TEMP_FILES%\reprojection
las2las -i %NDP_IMG%\*.%NDP_IMG_FORMAT% ^
        -set_classification 1 -set_version %LAS_V% ^
-%PROJECTION% -rescale 0.01 0.01 0.01 ^
        -odir %TEMP_FILES%\reprojection -olaz 
echo las2las END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo las2las END %TIME% %DATE%

REM crear tiles con un buffer con '-tile_size %TILE_SIZE%' para especificar el tamaño de los tiles con 
REM '-buffer %TILE_BUFFER%' en metros. El Buffer es necesario para reducir artefactos en los bordes de
REM los tiles. Los resultados se guardan en el formato LAZ con '-olaz' (see: lastile_README.txt)
mkdir %TEMP_FILES%\tiles_raw
lasindex -i %TEMP_FILES%\reprojection\*.laz
echo lasindex END %TIME% %DATE% \n >>%OUTPUT_FILES%\NDP_log.txt
echo lasindex END %TIME% %DATE%
lastile -i %TEMP_FILES%\reprojection\*.laz ^
        -tile_size %TILE_SIZE% -buffer %TILE_BUFFER% ^
        -odir %TEMP_FILES%\tiles_raw -olaz
echo lastile END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lastile END %TIME% %DATE%

REM busca puntos de traslape (overage) entré lineas de vuelo y eliminarlos
REM y usar todos los nucleos definidos con '-cores %NUM_CORES%'. (see: lasoverage_README.txt)
mkdir %TEMP_FILES%\sort 
lassort -i %TEMP_FILES%\tiles_raw\*.laz ^
           -odir %TEMP_FILES%\sort -olaz ^
  -cores %NUM_CORES%
echo lassort END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lassort END %TIME% %DATE%

REM usar lasground_new para detectar puntos del suelo en todos los tiles con opcion '-city' (incremento 25m)
REM y '-ultra_fine' (see: lasground_README.txt). '-odir tiles_ground -olaz' parametros definen a guardar los
REM resultados compresados en la ruta 'tiles_ground'. el proceso corre en paralelo en caso de varios tiles
mkdir %TEMP_FILES%\tiles_ground
lasground_new -i %TEMP_FILES%\sort\*.laz ^
-step 3 -sub 8 -ultra_fine -compute_height ^
-odir %TEMP_FILES%\tiles_ground -olaz ^
-cores %NUM_CORES%
echo lasground END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lasground END %TIME% %DATE%

REM detectar construcciones en base de altura arriba suelo
mkdir %TEMP_FILES%\tiles_classify
lasclassify -i %TEMP_FILES%\tiles_ground\*.laz ^
-odir %TEMP_FILES%\tiles_classify -olaz ^
-cores %NUM_CORES%
echo lasclassify END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lasclassify END %TIME% %DATE%

REM clasificar suelo de nuevo con opcion de ignorar pre-detectados construcciones
mkdir %TEMP_FILES%\tiles_ground2
lasground_new -i %TEMP_FILES%\tiles_classify\*.laz ^
-wilderness -ignore_class 6 -compute_height00 ^
-odir %TEMP_FILES%\tiles_ground2 -olaz ^
-cores %NUM_CORES%
echo lasground2 END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lasground2 END %TIME% %DATE%

REM clasificar todos no-suelo/no-construccion ecos como vegetacion
mkdir %TEMP_FILES%\tiles_vegetation
lasheight -i %TEMP_FILES%\tiles_ground2\*.laz ^
-ignore_class 6 ^
-drop_below -3 -drop_above 100 ^
-classify_between 0.3 1.3 3 -classify_between 1.3 3.5 4 ^
-classify_between 3.5 100 5 ^
-odir %TEMP_FILES%\tiles_vegetation -olaz ^
-cores %NUM_CORES%
echo lasheight END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lasheight END %TIME% %DATE%

REM detectar construcciones de pequeño tamaño
mkdir %OUTPUT_FILES%\NDP
lasclassify -i %TEMP_FILES%\tiles_vegetation\*.laz ^
-odir %OUTPUT_FILES%\NDP -olaz ^
-cores %NUM_CORES%
echo lasclassify2 END %TIME% %DATE% >>%OUTPUT_FILES%\NDP_log.txt
echo lasclassify2 END %TIME% %DATE%

REM generar DTM y curvas de nivel desde los ecos de clasificacion suelo
mkdir %OUTPUT_FILES%
blast2dem -i %OUTPUT_FILES%\NDP\*.laz ^
-keep_class 2 ^
-step 1 -use_tile_bb ^
-odir %OUTPUT_FILES% -oasc ^
-cores %NUM_CORES%
lasgrid -i %OUTPUT_FILES%\NDP\*.laz - merged ^
-step 1 -average ^
-o %OUTPUT_FILES%\DTM.tif
blast2iso %OUTPUT_FILES%\NDP\*.laz - merged ^
 -iso_every 5.0 ^
 -smooth 2 -simplify_length 0.5 -simplify_area 0.5 -clean ^
 -o %OUTPUT_FILES%\CN.shp

REM eliminar carpetas temporales
rmdir %TEMP_FILES% /s /q

echo End Automatic Pointcloud Classification %TIME% %DATE% \n >>%OUTPUT_FILES%\NDP_log.txt

Cheers,
Tobias

Susana Gonzalez

unread,
Jun 23, 2017, 2:00:46 AM6/23/17
to Muriel Lavy, LAStools - efficient tools for LiDAR processing

Hi Muriel,

I have tried to run your model with different settings and I am not getting a better point cloud than your one.

The aerial images don’t seem to have enough overlap between them to create a good point cloud (areas where the flight lines are further apart do have a few gaps, such as the gap to the left-hand side of the orthomosaic) and the water bodies in the area are making a lot of noise (which is to be expected).

Cheers

Susana

Martin Isenburg

unread,
Jul 4, 2017, 12:48:06 PM7/4/17
to LAStools - efficient command line tools for LIDAR processing
Hello,

thank you Muriel for sharing your data. That allowed me to have a crack at it and design a new LAStools pipeline for attacking dense-matching points with the kind of excessive low noise as your data set is exhibiting it along the river area. Folks may download the data (and the latest version of LAStools) to follow along step by step:


Regards,

Martin @rapidlasso

--
muriel_pancalieri_processing_raw.jpg
muriel_pancalieri_processing_denoised_marked.jpg
Reply all
Reply to author
Forward
0 new messages