First time posting here. I'm working on my first remote sensing project to estimate forest structure attributes (e.g. forest height) using satellite images with deep learning. I've struggled for a while constructing the datasets for learning. Specifially, making sure the reference forest maps overlap the satellite images. QGIS actually shows the correct overlap after preprocessing (see picture 1). However, when I load the image in Python using rasterio each index (pixel) does not correspond the same pixel in the other images (see picture 2).

Why does QGIS show nearly perfect alignment while Python show a clear misalignment? Is there something wrong in the way I've preprocessed the image?


Reference map (white to red gradient) overlaps background optical image nearly perfectly.

enter image description here


The same images (left: optical, mid: sar, right: reference map) but accessing by the same indices (pixels) show different areas (by up to 50 pixels)

enter image description here

Preprocessing steps

  • Reprojecting - EPSG:4326 -> EPSG 32632
    • gdalwarp -t_srs EPSG:32632 -r lanczos -wo SOURCE_EXTRA=1000 -co COMPRESS=LZW opt.tif opt_reprojected.tif
  • Clipping - Clip all images according to a polygon that is an intersection between all images. See clip_img_to_aoi()
  • Transforming - ​Manually fix the last misalignments between imgs
    • gdalwarp -ct "+proj=affine +xoff=16 +yoff=0" opt.tif transformed_opt.tif
  • Resizing - Convert image to correct resolution and size (i.e. numpy shape). See resize_img()

Post-processing Rasterio Stats

Ground Truth (Reference map)
(1, 8904, 11801)
| 0.50, 0.00, 290779.00|
| 0.00,-0.50, 6809696.50|
| 0.00, 0.00, 1.00|

(4, 8904, 11801)
| 0.50, 0.00, 290778.15|
| 0.00,-0.50, 6809677.65|
| 0.00, 0.00, 1.00|

Clip image to area of interest
def clip_img_to_aoi(img_path, aoi_path='data/aoi.geojson'):
    with rasterio.open(img_path) as img_ds:
        with fiona.open(aoi_path) as geo:
            aoi = [feature["geometry"] for feature in geo]
        img_clipped, img_transform = rasterio.mask.mask(img_ds, aoi, crop=True)
        meta = img_ds.meta.copy()
        meta.update({"transform": img_transform,
        with rasterio.open(img_ds.name[:-4] + '_clipped.tif', 'w', **meta) as f:
Resize image to appropriate resolution and shape
def resize_img(SAT_img, from_res, to_res):
    factor = to_res / from_res
    shape = (np.round(SAT_img.shape[1] / factor).astype(int) , np.round(SAT_img.shape[2] / factor).astype(int))
    SAT_image_resized = np.zeros((SAT_img.shape[0], shape[0], shape[1]))
    for channel in range(0, SAT_img.shape[0]):
        print("Resizing channel ", channel)
        SAT_image_resized[channel, :, :] = resize(SAT_img[channel, :, :], (shape[0], shape[1]), anti_aliasing=True)
    return SAT_image_resized

def change_meta_data(meta_data, factor):
    meta_dataNew = meta_data
    meta_dataT = meta_dataNew['transform']
    meta_dataNew['transform'] = Affine(meta_dataT[0] * factor, meta_dataT[1], meta_dataT[2], \
                                       meta_dataT[3], meta_dataT[4] * factor, meta_dataT[5])
    meta_dataNew['width'] = np.round(meta_data['width'] / factor).astype(int)
    meta_dataNew['height'] = np.round(meta_data['height'] / factor).astype(int)
    return meta_dataNew

1 Answer 1


I managed to solve this myself. The culprit seemed to be the polygon (area of interest). There was a small edge that was hard to notice which extended beyond the optical image a few pixels. Recreating the shape to a true intersection of all the images fixed the issue.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.