Python – Georeferencing an Unreferenced Image Using Rasterio

Question

I have 6-band unregistered/unreferenced TIFF images that I have been trying to export as georeferenced GeoTIFFs. The GPS coordinates of the corners of each image is known and each image is 1024 height and 1280 width. I attempted the GDAL solution in https://stackoverflow.com/questions/55681995/how-to-georeference-an-unreferenced-aerial-image-using-ground-control-points-in but was receiving item has no get_attr_ errors at the setGCP command.

Attempt 1:

import rasterio as rio
from rasterio import Affine
from rasterio.warp import reproject, Resampling, calculate_default_transform
from rasterio.transform import from_origin
from rasterio.control import GroundControlPoint

im = rio.open('path to a tif')
arr = im.read()
dst_fn1 = outdir + 'out.tif'

tl = GroundControlPoint(0, 0, -83.70113410256013, 42.307951446432604)
bl = GroundControlPoint(1024, 0, -83.69940501521428, 42.307603183805234)
br = GroundControlPoint(1024, 1280, -83.698829074736, 42.3091785425499)
tr = GroundControlPoint(0, 1280, -83.70055820297041, 42.309526812647555)
gcps = [tl, bl, br, tr]

src_crs = im._get_crs()  ## This is a blank object. No crs 
dst_crs = {'init': 'epsg:4326'}
src_tran = im.transform
dst_transform, dst_width, dst_height = calculate_default_transform(a, crs, 1280, 1024, gcps=gcps)
dest_ds = np.empty(shape=(arr.shape[0], dst_width, dst_height))
reproject(sand, dest_ds, src_crs=dst_crs, dst_crs=dst_crs, 
      src_transform=src_tran, dst_transform=dst_transform, 
      resampling=Resampling.nearest)

This method does not throw any errors but the resulting array is only zeroes. Also exported array is not rotated to match where the corners should be. Perhaps this is the wrong expectation but I thought that the data would be inside the yellow box and the outside would be nodata values to square it off.

Attempt 2:

( This is essentially the method found here https://www.earthdatascience.org/courses/earth-analytics-python/lidar-raster-data/reproject-raster/)

with rio.open('path to a tif') as src:
    transform, width, height = calculate_default_transform(src.crs, dst_crs, src.width, src.height, gcps=gcps)
    kwargs = src.meta.copy()
    kwargs.update({'crs': dst_crs, 'transform': transform, 'width': width, 'height': height})
        with rio.open('path/output.tif', 'w', **kwargs) as dst:
            for i in range(1, src.count + 1):
                reproject(
                    source=rio.band(src, i),
                    destination=rio.band(dst, i),
                    src_transform=src.transform,
                    src_crs=src.crs,
                    dst_transform=transform, dst_crs=dst_crs,
                    resampling=Resampling.nearest)

Results in this error:

rasterio._err.CPLE_AppDefinedError: The transformation is already "north up" or a transformation between pixel/line and georeferenced coordinates cannot be computed for TIFFPATH. There is no affine transformation and no GCPs. Specify transformation option SRC_METHOD=NO_GEOTRANSFORM to bypass this check.

Answer 1

Use rasterio.transform.from_gcps, then open the raster in append r+ mode and write the CRS and transform.

import rasterio as rio
from rasterio.transform import from_gcps
from rasterio.control import GroundControlPoint


tl = GroundControlPoint(0, 0, -83.70113410256013, 42.307951446432604)
bl = GroundControlPoint(1024, 0, -83.69940501521428, 42.307603183805234)
br = GroundControlPoint(1024, 1280, -83.698829074736, 42.3091785425499)
tr = GroundControlPoint(0, 1280, -83.70055820297041, 42.309526812647555)
gcps = [tl, bl, br, tr]

transform = from_gcps(gcps)
crs = 'epsg:4326'

with rio.open(filepath, 'r+') as ds:
    ds.crs = crs
    ds.transform = transform