Short answer

Slope is rise over run. In a geographic CRS such as WGS84 (EPSG:4326) the run — the horizontal spacing between pixels — is expressed in degrees of latitude and longitude, while the rise — elevation — is in metres. The slope algorithm computes elevation difference divided by horizontal distance, so it divides metres by degrees. That is not just imprecise, it is dimensionally incoherent, and because one degree spans roughly 111 km while elevation differences are tens of metres, the run is enormously under-stated and the slope comes out wildly exaggerated, often pinned near 90 degrees across the whole raster.

The fix is to make the horizontal and vertical units the same. Either reproject the DEM into a projected, metre-based CRS before running slope, or supply a z-factor that converts the units — and the z-factor route is fragile because the conversion depends on latitude.

What slope algorithms actually compute

Standard slope (the Horn 1981 method used by GDAL gdaldem and the QGIS/ArcGIS slope tools) fits a plane to a 3×3 window of cells and computes the gradient from the elevation differences across the window, divided by the cell spacing. In pseudo-form:

dz/dx = (elevation change east–west) / (2 × cell_width)
dz/dy = (elevation change north–south) / (2 × cell_height)
slope = atan( sqrt( (dz/dx)^2 + (dz/dy)^2 ) )

The algorithm trusts that cell_width, cell_height, and the elevation values are all in the same linear unit. In a 1-arc-second SRTM tile, the cell size reported in the raster is 0.0002777778 degrees, not 30 metres. Feed that into the formula against metre elevations and dz/dx is inflated by a factor of about 111,000 — the number of metres in a degree. The slope saturates.

The latitude trap

Even if you try to compensate, geographic coordinates hide a second problem: a degree of longitude is not a fixed ground distance. It shrinks with the cosine of latitude.

  • At the equator, 1 degree of longitude ≈ 111.3 km.
  • At 45 degrees latitude, ≈ 78.7 km.
  • At 60 degrees latitude, ≈ 55.7 km.

A degree of latitude stays close to 111 km everywhere. So in a geographic CRS, your pixels are not even square in ground terms — they are progressively wider (E–W) than tall (N–S) as you move away from the equator. Any single horizontal-unit correction is therefore only valid at one latitude. Over a DEM that spans several degrees of latitude, no constant z-factor is correct everywhere.

The two fixes

Fix 1 (recommended): reproject to a metric CRS

Reproject the DEM into a projected coordinate system whose units are metres and that fits your area — usually the local UTM zone, or a national grid. With GDAL:

gdalwarp -t_srs EPSG:32633 -r bilinear -tr 30 30 \
  dem_wgs84.tif dem_utm33n.tif

Here EPSG:32633 is WGS84 / UTM zone 33N, -r bilinear is appropriate resampling for continuous elevation (never nearest-neighbour for a DEM you will derive slope from — it introduces stair-step artefacts), and -tr 30 30 sets a clean 30 m grid. Then run slope on the reprojected raster:

gdaldem slope dem_utm33n.tif slope_deg.tif

Now rise and run are both metres, the z-factor defaults to 1, and the output is correct. In QGIS this is Raster → Projections → Warp (Reproject) followed by Raster terrain analysis → Slope, or the GDAL slope algorithm in the Processing toolbox.

Fix 2 (last resort): set a z-factor

If you genuinely cannot reproject — say the data must stay in 4326 for a tiling pipeline — you can pass a z-factor that scales vertical units into the horizontal degree units. For a DEM with elevations in metres, a workable approximation at a representative latitude φ is:

z-factor = 1 / (111320 × cos(φ))

ArcGIS Pro exposes this directly: the Slope tool, when the input is in a geographic CRS, can apply a latitude-aware z-factor, and Esri publishes a lookup table of z-factors by latitude. QGIS and gdaldem expose a -s/scale parameter that does the same thing (gdaldem slope -s 111120 ... is the common shorthand, valid only near the equator). The catch is everything above: this is correct only for a thin latitude band. Use it for small, low-relief, single-latitude tiles and document it. Reproject for anything else.

Worked example

You download a 1-arc-second Copernicus DEM tile centred at 60°N. Elevations are in metres, CRS is EPSG:4326.

  • Run slope as-is: nearly every cell reports slope above 89 degrees. Obviously wrong.
  • Apply a single equatorial scale of 111120: slope drops but is still too steep, because at 60°N a degree of longitude is only ~55.7 km, roughly half the equatorial value — your run is overstated, slope understated in the E–W direction and inconsistent with N–S.
  • Reproject to EPSG:32635 (UTM 35N) with gdalwarp, then run gdaldem slope: values now read like real terrain — valley floors near 0–3 degrees, hillslopes 15–35 degrees, cliffs above 45. Cross-check against contours or a known road grade to confirm.

QA and validation

  • Inspect the CRS first. gdalinfo dem.tif shows the coordinate system and, near the bottom, Pixel Size = (0.000277..., -0.000277...). Degrees in the pixel size is your warning sign.
  • Check the slope histogram. A correct slope raster has a smooth distribution concentrated below 45 degrees for most landscapes. A spike at the maximum means a units mismatch.
  • Compare to a control. Measure a slope you know — a quarry face, a graded road, a ski run — and confirm the raster agrees within a few degrees.
  • Re-check vertical units. A DEM in feet against a metric horizontal grid needs a z-factor of 0.3048. This is a separate axis of the same units problem and equally capable of producing wrong slope.

Common pitfalls and why they happen

  • Trusting on-the-fly reprojection. QGIS will display a 4326 DEM and a UTM basemap together, which fools people into thinking the DEM is "in metres." Display reprojection never changes the stored pixel spacing the slope tool reads.
  • Using Web Mercator (EPSG:3857) as the metric CRS. It is metres, but it distorts distance badly with latitude — the scale factor at 60°N is about 2. Slope from a Mercator DEM is wrong in the same family of ways. Use a conformal local projection (UTM, national grid), not Mercator, for terrain.
  • Nearest-neighbour resampling on reproject. It preserves exact values but creates blocky elevation steps that become spurious slope spikes. Use bilinear or cubic for DEMs.
  • One z-factor over a tall study area. Correct at the centre latitude, increasingly wrong toward the top and bottom edges.

Bathyl perspective

We treat the CRS of an elevation raster as part of the measurement, not a display setting, because slope, aspect, and curvature all silently inherit a units error from the grid they are computed on. The reliable habit is simple: reproject DEMs to a local metric CRS before any terrain derivative, and reserve geographic coordinates for storage and exchange, never for analysis.

Related Bathyl reading

Sources