search

OSGeo / gdal

GDAL is an open source MIT licensed translator library for raster and vector geospatial data formats.
https://gdal.org
Other
4.78k stars 2.5k forks source link

gdaldem: Slopes in geographical coordinates should use geodesic distances #3243

Open CarlosGrohmann opened 3 years ago

CarlosGrohmann commented 3 years ago

Expected behavior and actual behavior.

After some testing with different software, we noticed that GDAL has a limitation for computing slopes with DEMs in geographical coordinates. Since the distance between pixels is not calculated using geodesic formulas, the results are wrong for high-latitude data. In the histogram, we can see the difference in slopes calculated by GRASS, MicroDEM, ArcGIS (using the "geodetic" option), and GDAL. All produce similar results, but GDAL requires the scale parameter, which will then be applied to the entire dataset. The example is for a DEM in Alaska, at 60°N, so we can use s=55560 (from cos(60)*111120, but using 111120 directly also gives wrong results).

slope_gdal_bug

Expected: GDAL is everywhere today, so it should output correct results. The solution to this is to use geodetic distances to calculate the slope (GRASS uses the Haversine formula if I remember correctly).

The files can be accessed here:

DEM_latLong DEM_UTM

Steps to reproduce

gdaldem slope alaska_elevation_LL.tif alaska_slope_LL_GDAL.tif -of GTiff -b 1 -s 55560.0 -p

Operating system

Tested on Linux (Ubuntu 20.04 and Windows 10)

GDAL version and provenance

3.1.3 version ( I think that is from ubuntugis-unstable PPA)

jratike80 commented 3 years ago

It would be a nice new feature. I believe that geographical coordinates do not work well for other gdaldem processes either. Fortunately users are at least warned about this in https://gdal.org/programs/gdaldem.html

... For locations not near the equator, it would be best to reproject your grid using gdalwarp before using gdaldem.