Итеративно загружать блок изображения блоком, где блоки частично перекрываются

Question

Попытка обработать большое спутниковое изображение (~ 10 ГБ).Для эффективной обработки памяти кусок изображения (block / tile) загружается в память на каждой итерации.

Пример кода дляэто как показано ниже:

def process_image(src_img, dst_img, band_id=1):
    with rasterio.open(src_img) as src:
        kwargs = src.meta
        tiles = src.block_windows(band_id)
        with rasterio.open(dst_img, 'w', **kwargs) as dst:
            for idx, window in tiles:
                print("Processing Block: ", idx[0]+1, ", ", idx[1]+1)
                src_data = src.read(band_id, window=window)
                dst_data = src_data ** 2 # Do the Processing Here
                dst.write_band( band_id, dst_data, window=window)
    return 0

Однако, для любого вида обработки, который требует мудрой операции ядра (такой как любой фильтр convolve, такой как smoothing), это приводит к проблеме при обработке вблизи краевблоки.Чтобы решить эту проблему, каждый блок должен быть слегка перекрыт, как показано ниже:

Моя цель состоит в том, чтобы загрузить блоки следующим образом, где количество перекрытийможно отрегулировать в соответствии с необходимостью.

До сих пор я не нашел прямого пути для достижения этой цели.Буду признателен за любую помощь в этом отношении.

Answer 1

Вы можете просто расширить свои окна:

import rasterio as rio
from rasterio import windows

def overlapping_blocks(src, overlap=0, band=1):
    nols, nrows = src.meta['width'], src.meta['height']
    big_window = windows.Window(col_off=0, row_off=0, width=nols, height=nrows)
    for ji, window in src.block_windows(band):

        if overlap == 0:
            yield ji, window

        else:
            col_off = window.col_off - overlap
            row_off = window.row_off - overlap
            width = window.width + overlap * 2
            height = window.height + overlap * 2
            yield ji, windows.Window(col_off, row_off, width, height).intersection(big_window)

И вы бы использовали его так:

def process_image(src_img, dst_img, band_id=1):
    with rasterio.open(src_img) as src:
        kwargs = src.meta
        with rasterio.open(dst_img, 'w', **kwargs) as dst:
            for idx, window in overlapping_block_windows(src, overlap=1, band=band_id):
                print("Processing Block: ", idx[0]+1, ", ", idx[1]+1)
                src_data = src.read(band_id, window=window)
                dst_data = src_data ** 2 # Do the Processing Here
                dst.write_band( band_id, dst_data, window=window)
    return 0

И вот способ перекрывать как блочные окна, так и неблокированныеокна, с дополнительным аргументом boundless, указывающим, следует ли расширять окна за пределы растрового экстента, полезно для безграничное чтение :

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


def overlapping_windows(src, overlap, width, height, boundless=False):
    """"width & height not including overlap i.e requesting a 256x256 window with 
        1px overlap will return a 258x258 window (for non edge windows)"""
    offsets = product(range(0, src.meta['width'], width), range(0, src.meta['height'], height))
    big_window = windows.Window(col_off=0, row_off=0, width=src.meta['width'], height=src.meta['height'])
    for col_off, row_off in offsets:

        window = windows.Window(
            col_off=col_off - overlap,
            row_off=row_off - overlap,
            width=width + overlap * 2,
            height=height + overlap * 2)

        if boundless:
            yield window
        else:
            yield window.intersection(big_window)


def overlapping_blocks(src, overlap=0, band=1, boundless=False):

    big_window = windows.Window(col_off=0, row_off=0, width=src.meta['width'], height=src.meta['height'])
    for ji, window in src.block_windows(band):

        if overlap == 0:
            yield window

        else:
            window = windows.Window(
                col_off=window.col_off - overlap,
                row_off=window.row_off - overlap,
                width=window.width + overlap * 2,
                height=window.height + overlap * 2)
            if boundless:
                yield window
            else:
                yield window.intersection(big_window)

Окна вывода: