From cloud to landscape
Tutorial for processing LiDAR datasets and visualization (for Slovenian lidar)
Three hillforts above village Grahovo - near Cerkniško lake, Central Slovenia, Notranjska region
Tutorial for processing LiDAR datasets and visualization (for Slovenian lidar)
Three hillforts above village Grahovo - near Cerkniško lake, Central Slovenia, Notranjska region
Contents
2. Something about tutorial | 3. Software | 4. Slovenia - LiDAR portal | 5. From LiDAR pointcloud to ground points | 6. Further to DTM raster and 2D visualization | 7. Final word | 8. About author | 9. Literature1. Introduction
This tutorial is based on a workflow how to process Slovenian lidar data, but it can also be applied to other countries lidar data available in LAZ/LAS format. If you have lidar data available only (or already) classified and rasterized DTM see chapter 6.
A number of papers in Slovenian archaeological literature has been written about LiDAR and the use of LiDAR data in archaeology. If I only list the few: Štular 2011; Mlekuž and Budja 2008; Štular, Kokalj, Oštir, Nuninger 2012; Milk 2009, 2013, 2015, Mlekuž and Rutar 2013.
You probably thought when reading these (and others!) papers and admiring the maps that were made from the data of the LiDAR surface recording – how can I do this myself. Usually there is an issue right after downloading a file from the Slovenian environment agency portal in an unknown LAZ format. What to do next? How to open it?
Before we can use the LiDAR data point clouds for mapping, this data needs to be processed and put into the appropriate format. How to do this was my question few years ago when I started to get interested in LiDAR. I was a complete novice in the GIS environment. When I started, I haven’t found any integrated or slightly complete workflow on the internet about processing LiDAR that included the open-source GIS software. There were just some pieces on different GIS boards and forums.
This tutorial is product of experimenting and trying how to get things running from a classified cloud of lidar points to 2D visualization (maps) that is useful for archaeological research.
Data obtained through LiDAR recording is an exceptional set of information about the physical space (landscape) that can be used to explore archaeological sites and, in particular, landscapes. Even though LiDAR data is accessible online without restrictions, the value of these data is zero if we don’t know how to use it. The purpose of this tutorial is, in particular, to show one of the paths based on open-source and free software, how to bring this data to a level of readability and record that is useful for archaeological research.
2. Instructions
These workflow instructions here are an upgrade and a more detailed version of published paper in ARHEO no. X, journal of Slovenian archaeological society. In the future, there certainly will be some changes due to new programs and software updates. These changes will be clearly marked.
3. Software that we need
For our work we need to install programs and tools listed bellow. Most of them are open-source and all of them are free to use.
QGIS - One of the leading and most developed open source GIS programs. In addition to many features for creating maps and spatial functions, it also offers a wide range of plug-ins and languages. This is ArcGIS of open source scene….or even better. In this tutorial we are using 2.18.18 version, as version 3.0.xis not supported (yet!) by all plugins that we need. But you can have both installed to work with.
LAStools - A collection of tools for processing data in point clouds. Optimized for full data processing of LiDAR data. Some tools in the collection are open source and we'll use them in this workflow.
FUSION - Open source program with a full spectrum of algorithms for processing point clouds. It is being developed mainly for use in forestry. Please read install instructions, especially LASZip.dll file placement.
SAGA - GIS open source program, used and developed mainly for raster and vector analysis. It has so many tool, that experimenting is fun!
RVT - The main tool in our workflow for making DTM visualizations. Easy to use and simple, but powerful.
4. Source of LiDAR pointcloud data for Slovenian territory
For Slovenian territory all LiDAR data is available on a Slovenian enviroment agency portal. Here you can download data in different stages of data processing, formats and coordination systems.
4.1. Coordinate system
Lidar data is available in two coordinate systems (CRS): D48GK (EPSG: 3912) and D96TM (EPSG:3794). D48GK coordinate system. Most of the available raster and vector data managed and provided by the Surveying and Mapping Authority of the Republic of Slovenia (GURS) is written in the D48GK coordinate system. D96TM is a new and updated coordinate system in which all Slovenian GIS data will be recorded (or transformed) in the following years. If you have data from GURS that you have acquired over the years, I recommend using the D48GK. I also use D48GK for aesthetic reasons, because the edges are nice when overlaying different layers and it is easier to manage if we use same CRS - less room for errors. QGIS supports on-the-fly transformation of CRS within the project; there is nothing wrong with using D96TM, especially if you use data that was surveyed recently with GPS.
4.2. Formats of LiDAR data pointclouds
On the portal we have two different pointcloud formats available.
zLAS is a proprietary format developed by company ESRI. It is supported and it is native in ArcGIS program. Because it is not an open-source format we will not use it in this tutorial.
LAZ is open-source compressed pointcloud format and it is supported (is readable) by most open-source and free software. Uncompressed format is LAS. As LAZ is not 100% supported by software that we will use (SAGA), we will in first steps transform it to LAS format.
More information about LAZ/LAS format and its properties: LAS format
Something about "battle" between zLAS in LAS/LAZ formats: zLAS vs. LAS/LAZ
4.3. Levels of processed data
We have three different levels of processed data available.
5. From LiDAR pointcloud to ground points
At this point, we need to have installed all software that is listed above. At the beginning, we need to set some things to make the work easier. In order to be able to run the LAS and FUSION tools from the Command Prompt (Cmd), we need to set the path to the programs. Open the Control Panel and look for Edit the System Environment Variables and then click System Variables. In the System variables window, we will find Variable (PATH), mark it and click Edit. At the end of already entered paths to the folders in the Variable value field paste "C:\lastools\bin;C:\FUSION;". This paths are the folders where I have both programs installed. If you have folders with the programs Lastools and Fusion installed elsewhere, edit it accordingly. To check if you set it up correctly open a Command Prompt (Cmd) and type laszip.exe. If LAStools graphical interface opens, we are ready to go. If it does not work, check if you have entered correct path to the installation folders.
Within the QGIS in toolbox, we need to configure and enable LAStools and FUSION tools to work.
In QGIS we open Processing – Options. Inside Options window we enable under Providers Tools for LiDAR data both tools (FUSION in LAStools) with simply entering their install folder paths in blank fields and check box Activate. If we get error when we click OK saying that msys (wrong value) is missing, we need to set up under tab Providers – GRASS commands - line Msys folder like you see on screenshot below.
Both previous steps in chapter 5 where we set up cmd paths and QGIS tollbox need to be done ONLY FOR THE FIRST TIME :)
Now we can start working with *.LAZ format GKOT file(s), that we downloaded from eVode portal. First we will extract all points that are representing surface points. As I have mentioned before not all software that we use supports compressed LAZ format, so we have to decompress it.
We will do this with a tool called laszip.exe that we can run in few different ways. We can open it directly from a folder that we installed LAStools (C:LAStools\bin - default install folder). We can also run it through Command prompt (Cmd) or as a tool in QGIS. To work in Cmd, we need to run it and then go to folder where we have our LAZ files. When path is correct we enter laszip.exe and run it. The GUI of LAStools will open. Then we select all LAZ files that we want to decompress, mark output file as LAS and run decompress. If we want to run everything in cmd the command is laszip *****.LAZ or laszip -lof file_list.txt if we have more files. LAS file(s) will be in a folder where we have LAZ file(s)
We can also run LAStools GUI through QGIS. We find in toolbox laszip and with doubleclick run it. When QGIS GUI opens we find in input LAS/LAZ file field our LAZ file, then mark open LAStools GUI and press Run to open LAStools GUI. From here on it is the same as steps above.
As we usually have more than one LAZ file to work with, we can also do merging and decompressing in one go with lasmerge.exe tool. Procedure on how to run it it the same as laszip.exe. Here we select all file that we want to merge, file will be named by default merged.las, but we can change it. The only thing that we need to be careful is selecting LAS format before clicking RUN.
As LiDAR data for Slovenia is already classified we can use LAZ files with only ground points (OTR files, see the end of chapter 4). I usually use it when I am lazy. But it is better to make your own pointcloud of ground points. From my experience you gain some more info on potentially interesting features such as drywalls, stone mounds etc...
To get a pointcloud that represents surface we have to run GroundFilter tool that is part of FUSION software. We can run it through QGIS toolbox on the right side. We can find it under Tools for LiDAR data -> Fusion or just use search. In GroundFilter we enter Input LAS layer a *.LAS file that we created in previous step. Under Output LAS we enter a new *.LAS file where a new surface pointcloud will be stored. Under Cellsize for intermediate surfaces we need to enter 1 to get the most detailed surface pointcloud. It is important that if we have classified LAS file we can put under Additional modifiers a modifier /class:0,2,6. With this setting we will remove all vegetation and other points that are almost 100% above surface. We also gain some time with this as logarithm doesn't need to calculate all points. When we run this tool we will get a surface pointcloud.
Now a little about classified pointclouds. Classified pointcloud in LAS/LAZ format has 13 possible different classes. Not in every file we will have all of them. Classification is done with special tools that automatically define what a point in space represents in reality. When processing surface pointclouds I usually use class 0 (points that algorithm didnt classified), 2 (points classified as surface) and 6 (points classified as buildings). I don't recommend usage of other classes as they usually bring more noise into the surface pointcloud. This choice of classes is for Slovenian LiDAR data. For data from other countries try experimentig a little. :)
More about classes of classified pointcloud in LAS format you can find in FUSION manual, that you can find ithere.
6. Further to DTM raster and 2D visualization
For next step we will need SAGA software. First we have to open it and the on the left side inside toolbox menu we need to open Import/Export – LAS – Import LAS files tool, so we can import LAS file that we created in previous steps. Next we need to open in left toolbox menu a tool Shapes – Point clouds – orodje Point Cloud to Grid. In this tool we need to write into cellsize filed a value of 0.5 or 1. When this tool finish its work we will get a grid (raster) made out of pointcloud.
Now we have a grid made out of points, but space between them is empty. We need to fill it and this process is called interpolation. Again in (SAGA) left toolbox menu we runGrid – Gridding – Inverse Distance Weighted tool. In this tool we need to be careful to properly set 3 things. In Attribute field we need to choose value Z, then in Target Grid System field we need to choose grid or grid system and in Target Grid choose »create«. When we finish this process we have made a DTM raster.
I personally favour the tool we used (IDW) and another one that we find in SAGA toolbox (left side menu) Grid - Tools - CloseGaps tool. It is much faster tool than IDWV and it is easier set up, but sometimes it smoothens low profile features a little to much. Try it and decide it for yourself.
For all of you that wants to know more about different interpolation methods and what method produces best results in different terrains here is the link.
Before we start with visualization of created raster we need to export it in a readable file format (in our case GeoTIFF format). In left (SAGA) toolbox menu we open Import/Export – GDAL/OGR – orodje Export GeoTIFF tool. In tool under Grid field we choose raster, that we created last. It has if we created it with IDW method of interpolation a Inverse Distance Weighted name in it. With this we have finished working in SAGA software and we have created DTM raster of surface.
Now we need to make visualisation of DMR/DMT. We will do this in RVT tool, where we have all types of visualisation that we need in archaeology on our disposal. You can read more about different methods of visiualization in RVT manual. There is also a good article about different LiDAR visualizations that are avaliable in RVT software.
Output of visualization is in *.tif format and in two versions, 32 and 8 bit. The 8bit version is optimised for viewing in non-GIS software and it can be viewed in Microsoft Windows photo viewer and any other photo viewer (GIMP, Photoshop, IrfanView).
Now that we have made visualizations that we wanted we can open them in QGIS or any other GIS software. In QGIS go to Add Raster function to import it to project. It will ask us to select right coordination system. In our case we have to select MGI 1901/Slovene National Grid; EPSG:3912 if we used and downloaded files in D48GK coordination system. Slovenia 1996/Slovenian national grid EPSG:3794 coordinate system we need to use if we downloaded and used files in D96TM coordinate system.
We are ready to start to work on interpretation, creating maps or 3D models in QGIS. I have in plan to write some other tutorials about this themes in near future.
If you want to read more about QGIS - link. And Google finds almost anything :).
7. And for the end of tutorial
This tutorial is not set in stone and probably there are some other ways how to do it (maybe shorter, more easy?). I regularly check and test new things in GIS open source world. I will add new things to tutorial and update it. I'm open to feedback, opinions and new ideas. Please don't hesitate :)
8. Author
My name is Jošt Hobič and I am a archaeologist. I have worked for over 10 years as field arcaheologist in Slovenia. I have experience working in almost all environments (except under the water). I was also working for one year in UK for Cotswold archaeology - it was great!. As field archaeologist I have specialised in field photography. I am working as a volunteer raising awareness of local heritage and local archeological sites. I currently work at Ministry of Culture of Republic of Slovenia on a GIS/National heritage project.
What I love is researching landscapes and through this I started to work with GIS and LiDAR. Although I had some GIS lessons during my study I didn't get much out of it. That meant that i needed to learn almost everything anew (I am still learning).
Im regular user of Twitter where I usually publish different LiDAR visualizations of archaeological sites, most of them under #FridayLiDAR or #HillfortsWednesday.
I'm also publishing 3D models of LiDAR visualizations on a Sketchfab. You can check them out here.
On my blog I write in very irregular time intervals articles about archaeological sites in a landscape where I live and about LiDAR. There I maintain a list of List of freely accessible LiDAR data and digital terrain models from all over the World. You can find it here.
9. Literature
