Downscaling a tiff image by aggregating 5x5 pixels based on stats

Question

How to downscale a tiff image of 10m resolution and create a new image of 50m where each pixel is stats from the first image?

The initial tiff image is a binary classification map - meaning each pixel (10m) belongs either to class "water" (value =0) or class "ice" (value=1). I would like to create a new image, where each pixel is the percentage of water in a 5 x 5 block of the initial map, meaning each pixel of the new image will have a 50 m resolution and represents the ratio or percentage of "water" pixel on every 5x5 pixel of the former map. You can see the example here: Example

Here is an image sample (can be downloaded from google drive): https://drive.google.com/uc?export=download&id=19hWQODERRsvoESiUZuL0GQHg4Mz4RbXj

user_2836621 · Answer

Your image is saved in a rather odd format, using a 32-bit float to represent just two classes of data which could be represented in a single bit, so I converted it to PNG with ImageMagick using:

magick YOURIMAGE.TIF -alpha off image.png

Many Python libraries will stutter on your actual TIFF so maybe think about using a different way of writing it.

Once that is done, the code might look something like this:

#!/usr/bin/env python3

from PIL import Image
import numpy as np

# Set size of tiles for subsampling
tileX = tileY = 5

# Open image and convert to greyscale and thence to Numpy array
im = Image.open('image.png').convert('L')
na = np.array(im)

# Round height and width down to next lower multiple of tile sizes
h = (na.shape[0] // tileY) * tileY
w = (na.shape[1] // tileX) * tileX

# Create empty output array to fill
res = np.empty((h//tileY,w//tileX), np.uint8)

pxPerTile = tileX * tileY
for yoffset in range(0,h,tileY):
    for xoffset in range(0,w,tileX):
        # Count ice pixels in this 5x5 tile
        nonZero = np.count_nonzero(na[yoffset:yoffset+tileY, xoffset:xoffset+tileX])
        percent = int((100.0 * (pxPerTile - nonZero))/pxPerTile)
        res[yoffset//tileY, xoffset//tileX] = percent

# Make Numpy array back into PIL Image and save
Image.fromarray(res.astype(np.uint8)).save('result.png')

On reflection, you can probably do it faster and more simply with cv2.resize() and a decimation of 0.2 on both axes and interpolation cv2.INTER_AREA

user_894763 · Answer

I did a version in pyvips:

#!/usr/bin/python3

import sys
import pyvips

image = pyvips.Image.new_from_file(sys.argv[1])

# label (0 == water, 1 == ice) is in the first band
label = image[0]

# average 5x5 areas
label = label.shrink(5, 5)

# turn into a percentage of water
water_percent = 100 * (1 - label)

# ... and save
water_percent.write_to_file(sys.argv[2])

I can run it on your test image like this:

$ ./average.py ~/pics/meltPondClassifiedAndS12.tif x.png

To make this (rather dark) output:

Downscaling a tiff image by aggregating 5x5 pixels based on stats

2 Answers