0

I used the Zonal Statistics tool to extract the max values of a raster/DEM layer. The attribute table shows the result, but the problem now is to identify the coordinates of the said maximum elevation. I tried looking for solutions online and saw Find the coordinates of the max n values of a raster layer but it does not apply to QGIS.

How do I identify these max values?

0

If you know a little bit of using python in Qgis, you can easily use this code and get the coordinates. If you don't know how to use the Python console, follow the steps below:

  1. Press CTRL+ALT+P on keyboard. This will open the console.
  2. Click 'Show Editor' in console panel. Look picture.
  3. Paste the following code to the editor. And replace 'your_raster_name' variable with your own raster name.(Look to the code)
  4. And run it.(Icon: Green Triangle)
  5. After running, Coordinates and elevation will be displayed as text in the console
  6. Note: Raster Crs and Project Crs must be same.

raster = QgsProject.instance().mapLayersByName("your_raster_name")[0]
xmin = raster.extent().xMinimum()
ymin = raster.extent().yMinimum()
xunit = raster.rasterUnitsPerPixelX()
yunit = raster.rasterUnitsPerPixelY()
vals = list()
for i in range(0,int(raster.width())):
    for j in range(0,int(raster.height())):
        val, res = raster.dataProvider().sample(QgsPointXY(xmin+i*xunit+xunit/2,ymin+j*yunit+yunit/2), 1)
        vals.append([xmin+i*xunit+xunit/2,ymin+j*yunit+yunit/2,val,res])
result = max(vals, key=lambda x: x[2])
print(f'Max Elevation is {result[2]}')
print(f'Coordinates are {result[0]} by {result[1]}')

enter image description here

2
  • Thank you for the code. I tried it, but it does not seem to work in my case. – JAT86 Oct 5 '20 at 6:38
  • your project crs and raster crs must be same. Maybe that's why. – Xi Jin Oct 5 '20 at 15:24
0

So, as output of the zonal stat, you have the max value of the raster in the attribute table for each polygon.

Convert the polygons to raster (raster > conversion > rasterize) with the same pixel size and extend than you raster dataset (in the parameter "output extent (xmim,ymin,xmax,ymax)" press the button on the right and "use layer/canvas extent") and based on the value of your "_max" field in the parameter "field for use for a burn in value [optional]" .

Use raster calculator to get all values where your input_raster = converted_zonal_stat_raster

Convert this result to points, either with "Processing > SAGA > Raster values to points" (which gives you the XY coordinate) or with "vector creation > raster pixels to points" (faster, but then you need to compute the coordinate values in the attribute table)

1
  • Thank you, can you please clarify what "same pixel size and extend than you raster dataset and based on the value of your "_max" field" in your answer above means? – JAT86 Oct 5 '20 at 6:39

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.