3

I’m trying (again) to find an automatic way to extract peaks and cols from a DTM raster (derived from LiDAR surveys).

By “peak” I mean simply the cell with the highest value of a neighbour of cells;

By “col” I mean the lowest value of a ridge that divides two water basins which “connects” two areas of higher elevations (that possibly culminate in their respective peak). This concept is related to the prominence of a peak as depicted in this figure:

enter image description here

I’ve searched some similar topics on this and most of the threads converge to r.param.scale algorithm and r.geomorphon but even spending some time trying several different combinations of their parameters none of them gave me a satisfying result.

As for the definitions I’ve given to peaks and cols, the r.param.scale algorithm on QGIS actually analyses the DTM to find peaks and cols (saddles) but It’s very far from being precise and consistent. If I leave all parameters as default I get as a result that many of the features are not correctly classified (omitted or misplaced). Here’s an example:

enter image description here

In this case, the peak on the lower-right corner of the screen is correctly classified as a peak, conversely, the other peak is omitted (and adjacent to some cells classified as “col”). Nearby there are some cells that seem correctly representing the col between these two peaks.

If I increase the threshold of the first parameter (slope tolerance) the omissions decrease, but many other cells are classified as peaks or cols incorrectly:

enter image description here

With r.geomorphon (QGIS) the matter is quite similar. It seems slightly more consistent in classifying peaks, but has no classification for cols:

enter image description here

Here, what I don’t understand is why, for example, the lower-right peak is classified with many adjacent cells, while conceptually the peak should be one and one only cell (except in cases with cells with the exact same DN, but I’ve confirmed it’s not the case).

The ideal solution for my case would be classifying peaks and cols as single, isolated cells, then converting DN values to shapefile points. Even more useful would be the possibility to set thresholds to limit the classification, for example, by prominence. This would avoid the over-classification of features due to surface roughness.

Since the main threads on this topic are already some years old, is there any new advanced tool/plugin/script that may extract this information? I mainly use QGIS.

6
  • Please Edit your Question to focus on one software stack. Specifying two is likely to get this post closed as lacking focus.
    – Vince
    Jan 15 at 12:01
  • You used r.param.scale but please take a look at the much more powerful grass.osgeo.org/grass-stable/manuals/r.geomorphon.html
    – markusN
    Jan 15 at 15:49
  • Impossible to answer this without your data
    – BERA
    Jan 15 at 16:05
  • @markusN I've used also r.geomorphon, it's written in the main post. It seems a little more accurate than r.param.scale, but it classifies peaks as clusters of cells even if the real peak should be only one cell
    – NorthSon
    Jan 15 at 16:06
  • @BERA it's a normal DTM raster. Here's the link: drive.google.com/file/d/1QXA62JEQh0xmLliZ74A6OYaWq-A3KciC/… The data is projected in EPSG 3003 and the peak is at coordinates 1614510, 4737730
    – NorthSon
    Jan 15 at 16:19
3

It is all about definitions. If you assume that your peak is any cell, that is higher than 8 of its neighbors, you'll pick any tiny bump in your terrain:

enter image description here

However if definition is any cell that is highest for neighborhood square 101*101, you'll get this:

enter image description here

To locate mountain passes you'll need to isolate ridges first. Picture below shows divides between watersheds (catchments) greater than 10,000 cells:

enter image description here

Use cells under divides as a mask for your analysis, i.e. limit focal statistics to the cells underneath mask. This time compute minimum for neighborhood square:

enter image description here

Unfortunately this discards passes on the ridges inside watersheds. To overcome that I'd compute divides between same size subcatchments to use as a mask raster.

Perhaps correct order is defining ridge cells first and proceed with peaks, cols after.

I'll keep peaks and cols zip for next couple of days.

Difference of peaks and mean value of neighbors will give you a 'prominence'.

7
  • It might be useful to incorporate Topographic Position Index and Slope Position classifications in the workflow to further define relief features jennessent.com/downloads/tpi-poster-tnc_18x22.pdf
    – Matt
    Jan 15 at 23:41
  • That looks good! I've compared your results with the data I've had already produced some time ago for that territory with a manual workflow and the peak extraction is really accurate. There are some that are missing but I guess is due to the neighbourhood square set to 101x101. I'd ask you, which tool did you use to scan the raster with arbitrary neighbourhood square dimensions? One more thing: Cols have been as well extracted with the same neighbourhood square, right?
    – NorthSon
    Jan 16 at 7:40
  • I used ArcGis tool 'Focal Statistics', similar should exist in Qgis, because it's basic numpy function. You can also specify circle instead of square neighborhood. Yes, for cols, see layer name.
    – FelixIP
    Jan 16 at 18:17
  • Most tricky part are ridges. Most certainly they are divides of smaller subcatchments if you're going to pick them inside basins shown, but you'll need to eliminate stream cells, because they will also pop up as cols. Making hundreds of 'same' size subcatchments is also not easy. I use python for that.
    – FelixIP
    Jan 16 at 18:24
  • 1
    This is raster calculator expression Con(dem >= FSTAT, dem). This is essentially IF statement meaning condition. Follow by raster to points. Use less or equal for passes.
    – FelixIP
    Jan 17 at 18:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.