Example of calculating an all-ones kernel given the pixel area (assuming pixel width/height in degrees)

kernel-by-pixel-area.py

import math
import sys

import numpy
from osgeo import gdal


def make_ones_kernel(radius_in_pixels):
    side_length = radius_in_pixels * 2
    if side_length % 2 == 0:
        side_length += 1
    return numpy.ones((side_length, side_length), dtype=numpy.uint8)


def create_kernel_from_global_raster_pixel_size(raster_path, kernel_radius_km):
    raster = gdal.Open(raster_path)

    _, pixel_width, _, _, pixel_height, _ = raster.GetGeoTransform()

    # convert pixel width/height from degrees to km
    mean_pixel_dim_degrees = (abs(pixel_width) + abs(pixel_height)) / 2
    earth_radius_km = 6371
    pixel_width_km = (mean_pixel_dim_degrees * 2 * math.pi * earth_radius_km) / 360.

    # convert from km to pixels
    print(f'Requested kernel radius in km: {kernel_radius_km}')
    print(f'Pixel width in km: {pixel_width_km}')
    radius_in_pixels = math.ceil(kernel_radius_km / pixel_width_km)

    print(f'Kernel radius in pixels: {radius_in_pixels}')
    kernel = make_ones_kernel(radius_in_pixels)

    actual_kernel_radius_in_km = ((kernel.shape[0]-1)/2)*pixel_width_km
    print(f'Computed kernel radius in km: {actual_kernel_radius_in_km}')
    return kernel


if __name__ == '__main__':
    print(create_kernel_from_global_raster_pixel_size(sys.argv[1],
                                                      float(sys.argv[2])))