Updating dimensions of raster with rasterio without scaling data

Question

I have a raster with a width of 2979 and a height of 1867 pixels. I need to split the raster into tiles of equal size (e.g. 5 kilometres x 5 kilometres). To do so, I need to update the height and width of the raster to fit an exact (i.e. integer) number of tiles. If the chosen dimension is 5km x 5km and the spatial resolution of the raster is 10 metres, the width of the raster would need to be increased to 3000 pixels and the height to 2000 pixels. This would create 6 (3000 / 500) x 4 (3000 / 500) = 24 tiles. I am using rasterio to work with the raster. This is my current code:

# Determine padding
padding_x = int(new_raster_width - raster.width)
padding_y = int(new_raster_height - raster.height)

# Update values based on padding
out_meta.update({
    'width': new_raster_width,
    'height': new_raster_height
})

# Loop through raster bands
for band_no in range(1, raster.count + 1):
    # Read band
    raster_band = raster.read(band_no)

    # Add padding
    padded_raster_band = np.pad(
        raster_band, 
        pad_width = (padding_x, padding_y), 
        mode = 'constant',
        constant_values = 0)

    # Append padded raster band to list
    padded_bands.append(padded_raster_band)

with rio.open(raster_out, 'w', **out_meta) as dest:
    for band_nr, src in enumerate(padded_bands, start=1):
        dest.write(src, band_nr)

I am trying to do the following: Calculate by how much pixels I need to increase the raster, update the width and height, loop through each band of the raster and update the values and, finally, write these bands into a new dataset.

This works. I do have, however, the following problem: The raster is also scaled, which results in an incorrect 'projection'. I have also included two images of the raster before and after increasing the size (although it might be hard to tell from the images).

How do I add a padding/increase the raster size without changing the scale of the raster?

Answer 1

You don't need to pad at all. Just read each tile, using boundless=True and let rasterio pad for you.

from itertools import product
import rasterio as rio
from rasterio.windows import Window

raster = '/path/to/input.tif'
width, height = 500, 500  # size of tiles in pixels
fill_value = 0

with rio.open(raster) as src:
    offsets = product(range(0, src.meta['width'], width), range(0, src.meta['height'], height))
    for col_off, row_off in offsets:
        window = Window(col_off=col_off, row_off=row_off, width=width, height=height)
        data = src.read(boundless=True, window=window, fill_value=fill_value)
        print(data.shape)

Answer 2

If you change the origin of the raster (typically the upper left) by adding the padded values, you need to also change the affine transformation (src.transform). But if you only add padding to the bottom and right of your original extent, everything should theoretically stay where it should be.

However as far as I can tell you don't exactly pad your array the right way. The way you have it right now I think it pads both axis the same amount. Try it this way and see if it works:

# Add padding
padded_raster_band = np.pad(
    raster_band, 
    pad_width = ((0, padding_y), (0, padding_x)), 
    mode = 'constant',
    constant_values = 0)

The first tubple of pad_width is the padding of the y axis, with the first value being the padding before the array and second after the array. The same is the case for the second tuple only for the x axis.