Short answer

Two rasters are aligned when their cells share the same size, the same orientation, and an origin offset that is an exact whole number of cells. Equal resolution alone is not enough. If grid A starts at x = 500000.0 and grid B starts at x = 500015.0 while both have 30 m cells, the cell boundaries never coincide, and any map-algebra or zonal operation must resample one layer, smearing values across a half-cell shift.

A snap raster is the cure: you nominate one grid as the reference and force every other raster to its cell size and origin. In GDAL you reproduce this with gdalwarp -tr (target resolution) plus -tap (target aligned pixels) and, where needed, -te (target extent). In ArcGIS it is the snapRaster environment setting.

The affine geotransform: what "alignment" really means

Every GeoTIFF stores a six-element affine geotransform that converts a pixel column/row to a map coordinate:

X = GT0 + col * GT1 + row * GT2
Y = GT3 + col * GT4 + row * GT5
  • GT0, GT3 - map coordinates of the top-left corner of the top-left pixel (the origin).
  • GT1 - pixel width (cell size in x); GT5 - pixel height, normally negative because rows increase downward.
  • GT2, GT4 - rotation terms, almost always 0 for a north-up raster.

Run gdalinfo dem.tif and you see these as "Origin", "Pixel Size", and the corner coordinates. Two rasters are aligned when GT1 and GT5 are identical, both are unrotated, and the difference between their origins is an integer multiple of the cell size:

(GT0_a - GT0_b) / GT1   is a whole number
(GT3_a - GT3_b) / GT5   is a whole number

If those quotients have a fractional part, the grids are offset and will not stack cleanly.

Why misalignment happens

  • Independent reprojection. Reprojecting two layers to the same CRS without a shared snap target lets gdalwarp choose each output extent independently, so origins end up on arbitrary coordinates.
  • Different source tiling. A DEM tile starting at a 1-degree boundary and a land-cover raster starting at a UTM-rounded coordinate will not share an origin even after both are in UTM.
  • Resampling at non-integer ratios. Going from 30 m to 25 m without snapping creates a fresh origin that matches nothing.
  • Clipping by a vector extent. Clipping to a polygon bounding box snaps the origin to the polygon, not to the raster grid.

The symptom is subtle: overlays look fine on screen, but a raster calculator expression like (slope > 30) AND (landcover == 5) either errors on mismatched grids or silently resamples and shifts the result by part of a cell, putting "steep + forest" pixels in slightly the wrong place.

Worked example: snapping to a reference DEM

Take a reference DEM ref.tif in EPSG:32631, 30 m cells, origin (500000, 5400000). You want to bring a Sentinel-derived raster s2.tif onto exactly this grid.

  1. Read the reference grid: 
    gdalinfo ref.tif
# Pixel Size = (30.0, -30.0); Origin = (500000.0, 5400000.0)
  2. Warp the second raster to match resolution and aligned pixels: 
    gdalwarp -t_srs EPSG:32631 -tr 30 30 -tap \
  -r bilinear s2.tif s2_aligned.tif
    -tap snaps the output extent to whole multiples of 30, guaranteeing the origin lands on a clean grid coordinate, so the offset from ref.tif is an integer number of cells.
  3. If you also need identical extent, supply the reference bounds explicitly: 
    gdalwarp -t_srs EPSG:32631 -tr 30 30 -tap \
  -te 500000 5100000 600000 5400000 \
  -r bilinear s2.tif s2_clip_aligned.tif
    With matching -tr, -tap, and -te, the two rasters now share cell size, origin, and extent.

In QGIS, the GDAL > Raster projections > Warp (reproject) dialog exposes the same controls, and the Raster calculator and Align Rasters tool (gdal:alignrasters) can snap a batch of layers to a reference in one pass. In ArcGIS Pro, set arcpy.env.snapRaster = "ref.tif" and arcpy.env.cellSize = "ref.tif" before running any Spatial Analyst tool.

Choosing the resampler when you resample

Snapping usually involves a small shift that triggers resampling, so the method matters:

  • bilinear / cubic for continuous data (elevation, temperature, reflectance).
  • nearest for categorical data (land cover, geology class) so codes are preserved exactly.
  • mode for downsampling categorical data when you want the majority class per output cell.

Using nearest on a DEM produces stair-step artefacts; using bilinear on a class raster invents non-existent classes between codes. Match the resampler to the data type, not to the default.

Common pitfalls and why they happen

  • Off-by-one-pixel zonal statistics. The zone raster and value raster are not snapped, so each zone borrows a fringe of neighbouring cells. Align both to a common snap raster first.
  • NoData creeping in at edges after warping. The target extent is slightly larger than the source; fill or clip to the valid extent.
  • Map algebra runs but results drift. Two grids with the same resolution but different origins were resampled on the fly. Pre-align rather than trusting implicit resampling.
  • -tap without a fixed -tr. -tap only makes sense with an explicit -tr; without it the resolution is inferred and alignment is not guaranteed.

Validation

  • Run gdalinfo on every layer and confirm identical Pixel Size and corner coordinates, or origins that differ by an exact multiple of the cell size.
  • Load all layers in QGIS and use the Raster calculator on a simple expression; if it runs without a "rasters do not align" warning, the grids match.
  • For categorical layers, check the value table before and after snapping to confirm no classes were invented.
  • Overlay at maximum zoom and verify cell corners coincide pixel-for-pixel.

Bathyl perspective

We treat the snap raster as a project-level decision made once, up front: one reference grid, one CRS, one cell size, recorded in the project metadata. Every derivative inherits that grid, so slope, aspect, land cover, and hazard layers stack cell-for-cell and map algebra means what it says. Alignment is invisible when it is right and corrosive when it is wrong.

Related reading

Sources