Flattening of Water

[arjannv]

I've found an issue with the flattening of water in the DSM file.

This happend with the new AHN4 09DZ1, subtile 1 (with default 4x5 subtiling).

There are some 6m high ridges in the resulting geotif as you can see here:

This is not present in the DTM file, so also tried to create the DSM with Laplace interpolation instead of the QuadIDW:

With Laplace the results are different but still a 6m difference.

Maybe this is an result of the overhanging tree which is classified as water? However, the trees are filtered out the LAS for the DTM, which does not have this problem:

Could this be handled with filtering the outliers with different parameters?

Edit:
It seems that there's a hole in the data where these ridges appear.

Would it be possible to only use ground and water for the interpolation? And ignore the other classes for this?

[arjannv]

Ok, I think I found why this happened.
The flattening of water did not happen in that spot because of the polygon used for flattening.
The part with the ridges was not included in the polygon. See the yellow bit below.

After playing a bit with the options in the QuadIDW interpolator and the flatten.patch function, I think I found a nice looking solution. As this was the only thing I'm interested in. Creating nice looking, hole-free GeoTIFFs for 3D printing.

As I found the Flattening polygon rather course, I disabled it for better looks. This introduced the ridge problem all over the map though..

I tried to use the QuadIDW interpolator to fill the holes by just using data from ground, water and buildings. The first iteration of the QuadIDW interpolator takes direct data from the x,y point. The next iterations the search radius increases with 3 (m?). For the next iterations I've created an extra LAS array with another PDAL filter. This is with the following classification filter: Classification[2:2],Classification[6:6],Classification[9:17] To use this I increased the iterations to 4.

This didn't solve my problem though. The ridges where formed by the flatten.patch function. This takes the median of some surrounding points to patch the holes. raster[yi, xi] = np.median(vals)
I've replaced this with the lower 10% percentile raster[yi, xi] = np.percentile(vals, 10) to prefer the lowest measurements to make the ridges as small as possible.
This does seem to create some (rare) issues with the interpolation of holes in data from glaced roof tiles (hole is 1,5m high instead of 6-10m):

However, this problem was already there with the standard median interpolation. Only less noticeable:

When there is more (and better) surrounding data, the problem is less noticeable (the hole is like an inverted dakkapel):

The modifications above, however still resulted in some plateaus from the QuadIDW interpolator.

To improve this, I ditched the extra LAS array and used only 1 iteration, so it only took data from the exact x,y point.

Altough, the results are now much better, it's still not there entirely:

There are still some plateaus noticeable when there are overhanging trees. The plateaus are not that high, less then 1m above water level.

I still have to check how this works out along the coast.

Changes done for this result:

• Replace src/interpolation/flatten.py line 159 with: raster[yi, xi] = np.percentile(vals, 10)
• Disable flatten.water in src/interpolation/interpolation.py on line 171 - 180

[mdjong1]

Nice exploration! We're not actively maintaining this project, but I'm happy to think along if you need extra eyes