Short answer

Quality control on a digital elevation model means verifying the header (CRS, vertical datum, units, pixel size, NoData), then hunting for the four artifact families that wreck terrain derivatives: voids, spikes and pits, striping, and tile seams. Run these checks before you compute slope, aspect, or hillshade, because every derivative operates on a 3x3 cell neighborhood and turns a single bad cell into a visible ring of garbage. The fastest first pass is gdalinfo -stats plus a hillshade preview; the durable pass is the nine-step checklist below.

Why DEM QC has to come first

A DEM is just a grid of elevation values, and almost every terrain product you build reads each cell together with its eight neighbors. Horn's slope algorithm (the default in gdaldem, ArcGIS, and QGIS) fits a plane to that 3x3 window. So one cell that is 200 m too high does not produce one bad slope value; it produces eight, arranged as a ring of near-vertical slope around the spike. After you render that as hillshade, it looks like a real knoll, and a reviewer three weeks later cannot tell whether it is terrain or an artifact. That asymmetry, cheap to catch now and expensive to catch later, is the whole reason QC belongs at the top of the workflow.

Step 1: Read the header with gdalinfo

gdalinfo -stats dem.tif

This single command answers most of the structural questions. Confirm:

  • Coordinate system. A projected CRS (e.g. UTM, EPSG:32633) with metres on both axes is what you want for slope. If it reports a geographic CRS like EPSG:4326, the pixel size is in degrees and slope will be wrong unless you reproject or set -s_srs/scale correctly.
  • Pixel Size. Look at the Pixel Size = (x, y) line. For a UTM grid you should see something like (30.0, -30.0). A negative y is normal (north-up). Wildly non-square pixels are a flag.
  • NoData Value. Note it. Common values are -9999, -32768, or 3.402823e+38 (float max).
  • Min/Max. STATISTICS_MINIMUM and STATISTICS_MAXIMUM. If max is 9000 m for a coastal site, you have spikes or NoData leaking in.

Step 2: Confirm units and vertical datum

Horizontal and vertical units must agree for slope to be meaningful. If horizontal is metres and elevation is feet, every slope is overstated by the 3.28 ratio. Read the band description or accompanying metadata; if the source is a US product, suspect NAVD88 in feet. Cross-check the elevation range against a known summit or sea level. Datum (ellipsoidal vs orthometric/geoid) matters less for slope and hillshade but is critical the moment you compare two DEMs or quote absolute heights, so record it now. See the dedicated checks in the related reading below.

Step 3: Verify the NoData value actually matches the fill

The single most common silent failure: the GeoTIFF declares NoData as -9999 but the void cells were written as 0 or 3.4e38. Then 0 is treated as sea-level terrain (a cliff at every void edge) or the float-max value blows up the statistics. Confirm with:

gdalinfo -stats dem.tif | grep -i nodata
gdal_calc.py -A dem.tif --calc="A==0" --outfile=zeros.tif

If the declared NoData is wrong, fix it before anything else:

gdalwarp -srcnodata 3.402823466e+38 -dstnodata -9999 dem.tif dem_fixed.tif

Step 4: Find voids visually

Voids appear in SRTM (radar shadow in steep terrain), in photogrammetric DEMs (water, deep shadow), and in LiDAR (specular returns over water). Load the DEM in QGIS, set the NoData rendering to a loud colour, and zoom the whole extent. For a quick raster of where data is missing:

gdal_calc.py -A dem.tif --calc="A==-9999" --outfile=voids.tif --NoDataValue=0

Decide per project whether to fill (gdal_fillnodata.py with a sensible -md max distance, e.g. -md 10) or to mask the affected area out of the analysis. Never fill silently across large gaps and then quote slope inside them.

Step 5: Hunt spikes and pits

Single-cell spikes and pits come from sensor noise and interpolation. They dominate slope because a 50 m jump over a 10 m cell is a 79-degree slope. Detect them by computing slope first as a diagnostic and looking for isolated maxima:

gdaldem slope dem.tif slope_check.tif -compute_edges

Histogram the slope. Anything above roughly 60-70 degrees over open terrain is suspect. To remove isolated outliers, a 3x3 median filter (QGIS Processing > SAGA > Simple Filter, or a focal statistics median in ArcGIS) flattens single-cell spikes while preserving real edges. Apply it sparingly; over-smoothing erases legitimate microtopography.

Step 6: Detect striping

Acquisition striping shows as regular north-south, east-west, or diagonal banding, common in older SRTM tiles and some satellite-derived products. It is invisible in a plain elevation colour ramp but jumps out under hillshade with a low sun angle:

gdaldem hillshade dem.tif hs.tif -az 315 -alt 25

Try two azimuths (e.g. 315 and 45) so you do not mistake real slopes aligned with one light direction for artifacts. Striping that persists across azimuths and follows a perfectly regular grid spacing is an artifact, not terrain.

Step 7: Check tile seams and edges

Mosaicked DEMs often have a vertical offset at tile boundaries because adjacent tiles used different vertical references or processing batches. A 1-2 m step at a straight, axis-aligned line is a seam. The hillshade from Step 6 reveals these as sharp linear features that ignore the actual topography. Use -compute_edges on all gdaldem calls so the outermost row and column are not left as a NoData fringe.

Step 8: Compare against a control

Where you have any trusted reference (a survey benchmark, a GNSS point, a higher-grade DEM), sample the DEM at those locations:

gdallocationinfo -valonly -geoloc dem.tif 512345 6712345

A consistent bias (every point off by +30 m) points to a datum mismatch; scattered random differences point to noise or resolution limits. Even three or four control points catch gross errors.

Step 9: Document the limits, then approve

Record CRS and EPSG code, vertical datum and units, native resolution, source and acquisition date, NoData handling, any filling or filtering applied, and the residual against control. This is what lets the next analyst reuse the DEM instead of re-deriving everything.

Common pitfalls and why they happen

  • Filling voids before checking the declared NoData. If NoData is wrong, the fill routine has nothing to fill and instead smooths real values. Fix the NoData tag first (Step 3).
  • Judging artifacts from a single hillshade azimuth. Real slopes that face away from the sun read as flat; artifacts aligned with the light hide. Always use two azimuths.
  • Median-filtering the whole DEM to kill a few spikes. This destroys legitimate sharp features like scarps and channel banks. Target the spikes, do not blanket-smooth.
  • Trusting reported Min/Max from a cached statistics file. If .aux.xml stats are stale, regenerate with gdalinfo -stats or -approx_stats so the numbers reflect the current data.

QA validation summary

Before the DEM is allowed into a derivative pipeline, confirm: header read and CRS projected; units and datum recorded; NoData verified against actual fill; voids mapped and handled explicitly; spike/pit histogram clean; two-azimuth hillshade shows no striping or seams; at least one control comparison done; limitations written down.

Bathyl perspective

We treat the DEM header and the first hillshade as a gate, not a formality. Most "mysterious" terrain features in a deliverable trace back to a NoData or seam problem that a two-minute gdalinfo and a low-sun hillshade would have caught. Documenting the QC outcome alongside the raster is what makes the elevation layer auditable and reusable instead of a black box.

Related reading

Sources