Resit for the assignments

About point clouds, DTMs, and visualisation

Deadline is Friday 19 April 2024 at 8:59 (as the exam resit starts)

Late submission? 10% will be removed for each day that you are late.

This single assignment (with several steps) is the resit for the 4 assignments combined. It is worth 40% of the total mark of the course.

This is an individual assignment.



Overview

In this assignment you need to:

  1. automatically classify a 300mX300m region of the AHN4 point cloud into ground and non-ground using the CSF algorithm
  2. make a raster DTM using both kriging and Laplace interpolation
  3. compare your DTM rasters with the official AHN4 DTM raster
  4. to help with visualisation you’ll have to implement a contouring algorithm

You have to use the following Python libraries (others are not accepted, or ask me first): numpy, startinpy, pyinterpolate, rasterio, scipy, laspy, fiona, shapely.

The dataset to use

You have to create the 0.50cm-resolution DTM of this 300mX300m area in the Netherlands:

min = (190250.0, 313225.0) max = (190550.0, 313525.0)

This is located in tile 69EZ1_21 (GeoTiles).

Step 0: Register for the resit

If you plan on doing this resit: send me an email (h.ledoux@tudelft.nl) please.

Step 1: Ground filtering + DTM creation with Laplace

For this step, you need to write a Python program that performs all the steps automatically, and it should be useable for other areas in the Netherlands so the parameters should be exposed (see below).

You need to implement the Cloth Simulation Filter algorithm (CSF) as explained in the terrainbook.

You need to classify the points as ground and use those to generate the DTM with Laplace interpolation. You are not allowed to use startinpy’s interpolate() function, you need to implement it yourself (but you can use startinpy triangulation and other functions).

Your code should be a file called step1.py and it should have the following input parameters:

name step1.py
inputfile LAZ
minx minx of the bbox for which we want to obtain the DTM
miny miny of the bbox for which we want to obtain the DTM
maxx maxx of the bbox for which we want to obtain the DTM
maxy maxy of the bbox for which we want to obtain the DTM
res DTM resolution in meters
csf_res resolution in meters for the CSF grid to use
epsilon threshold in meters to classify the ground

It should output 2 files in the same folder as `step1.py:

  1. a file called dtm.tiff representing the 50cm-resolution DTM of the area created with Laplace
  2. the ground points in a file called ground.laz

For instance,

python step1.py mypc.laz 190000 303600 190100 303700 1.0 5.0 0.2

Use this in Python:

import argparse
parser = argparse.ArgumentParser()
parser.add_argument("ifile", help="input file in LAZ")
parser.add_argument("param1", help="param1 used for this and that")
parser.add_argument("param2", help="param2 used for this and that")
args = parser.parse_args()
print("ifile={}; param1={}; param2={})".format(ifile, param1, param2))

I will test your code with other areas/files.

Step 2: Create the DTM of your region with Ordinary kriging (OK)

Perform OK on the ground.laz created in the step above, and create a 0.5mX0.5m DTM.

For this step, you should use Pyinterpolate and document in the report the steps you took and show the results of the intermediate steps (eg the experimental variogram and the nugget/range/sill/function you used).

For this step, the parameters you find by modelling the dataset can be hardcoded, but these should be documented in the report and you still have to submit the code you code (which is only for this dataset). The code could be a Jupyter notebook, this would be fine.

Step 3: Compare your two created DTMs

  1. with each other (Laplace versus OK),
  2. and with the official AHN4 file (0.5m DTM).

How to compare them is left to you (you can use existing libraries and software), but highlight their differences and try to explain why they are different in the report.

For this step I don’t need the code, simply your analysis in the report.

Step 4: Visualisation with iso-contours

Implement the ideas described in “Conversion to isolines” in the terrainbook to produce a GeoPackage file of the contours lines at every meter (“full meters”, eg 1.0/2.0/3.0/etc).

Show me the results in the report and submit the code.

The Python code step4.py should also be submitted and have the following behaviour

name step4.py
inputfile GeoTIFF

and it should output a GeoPackage file called isocontours.gpkg where the lines have the elevation as an attribute.

You can create GeoPackage easily with “sister package” to rasterio called “fiona” (which uses the excellent “shapely”).

Marking

The following will be evaluated, each part refers to both the code (we will try to run it) and to the report (quality, completeness, analysis).

Criterion Points
ground filtering 3.0
DTM Laplace 1.0
kriging 3.0
iso-contours 1.0
comparison DTMs 1.0
quality of report 1.0

What to submit and how to submit it

You have to submit:

  1. All the code you wrote with the names as defined above
  2. a report (about ~20-25 pages) in PDF (no Word file please) where you elaborate on:
    • for each step, how did you perform it (in theory, I can read the code)
    • present the main issues you had and how you fixed them
    • give a reason why you picked certain parameters for the CSF algorithm
    • the differences between your results and the official DTM of AHN4

Those files need to be zipped and the name of the file is your studentID, for example 5015199.zip.

email me (h.ledoux@tudelft.nl) the ZIP before the deadline.

[last updated: 2024-03-18 08:42]