开窗读写

从Rasterio 0.3开始，您可以读取和写入栅格文件的“窗口”。此功能允许您处理比计算机RAM大的栅格或并行处理大块的栅格。

Windows

窗口是栅格数据集矩形子集上的视图，在栅格中按列和行的偏移量、宽度和高度（以像素为单位）进行描述。这些可以是int或float。

from rasterio.windows import Window

Window(col_off, row_off, width, height)

窗口也可以由numpy数组索引元组或切片对象构造。在这些情况下，只允许使用int值。

Window.from_slices((row_start, row_stop), (col_start, col_stop))
Window.from_slices(slice(row_start, row_stop), slice(col_start, col_stop))

如果将height和width关键字参数传递给 from_slices 可以使用相对切片和开放切片。

Window.from_slices(slice(None), slice(None), height=100, width=100)
# Window(col_off=0.0, row_off=0.0, width=100.0, height=100.0)

Window.from_slices(slice(10, -10), slice(10, -10), height=100, width=100)
# Window(col_off=10, row_off=10, width=80, height=80)

阅读

下面是读取Rasterio测试文件的256行x 512列子集的示例。

>>> import rasterio
>>> with rasterio.open('tests/data/RGB.byte.tif') as src:
...     w = src.read(1, window=Window(0, 0, 512, 256))
...
>>> print(w.shape)
(256, 512)

写作

写作也有类似的效果。下面创建一个空白的500列x 300行geotiff，并将值为127的37500像素从geotiff左上角向下拖动30像素，向右拖动50像素到窗口中。

image = numpy.ones((150, 250), dtype=rasterio.ubyte) * 127

with rasterio.open(
        '/tmp/example.tif', 'w',
        driver='GTiff', width=500, height=300, count=1,
        dtype=image.dtype) as dst:
    dst.write(image, window=Window(50, 30, 250, 150), indexes=1)

结果：

抽取

如果写入窗口小于数据，数据将被抽取。在下面，窗口被缩放到源图像的三分之一。

with rasterio.open('tests/data/RGB.byte.tif') as src:
    b, g, r = (src.read(k) for k in (1, 2, 3))
# src.height = 718, src.width = 791

write_window = Window.from_slices((30, 269), (50, 313))
# write_window.height = 239, write_window.width = 263

with rasterio.open(
        '/tmp/example.tif', 'w',
        driver='GTiff', width=500, height=300, count=3,
        dtype=r.dtype) as dst:
    for k, arr in [(1, b), (2, g), (3, r)]:
        dst.write(arr, indexes=k, window=write_window)

结果是：

数据窗口

有时需要裁剪数据集周围节点数据值的外部边界：

from rasterio.windows import get_data_window

with rasterio.open('tests/data/RGB.byte.tif') as src:
    window = get_data_window(src.read(1, masked=True))
    # window = Window(col_off=13, row_off=3, width=757, height=711)

    kwargs = src.meta.copy()
    kwargs.update({
        'height': window.height,
        'width': window.width,
        'transform': rasterio.windows.transform(window, src.transform)})

    with rasterio.open('/tmp/cropped.tif', 'w', **kwargs) as dst:
        dst.write(src.read(window=window))

窗口转换

可以使用数据集的 window_transform 方法：

>>> import rasterio
>>> from rasterio.windows import Window
>>> win = Window(256, 256, 128, 128)
>>> with rasterio.open('tests/data/RGB.byte.tif') as src:
...     src_transform = src.transform
...     win_transform = src.window_transform(win)
...
>>> print(src_transform)
| 300.04, 0.00, 101985.00|
| 0.00,-300.04, 2826915.00|
| 0.00, 0.00, 1.00|
>>> print(win_transform)
| 300.04, 0.00, 178794.71|
| 0.00,-300.04, 2750104.30|
| 0.00, 0.00, 1.00|

窗口实用程序

基本的联合和交叉操作可用于Windows，以简化在动态创建的窗口中对一系列带区或数据集的操作，使其具有相同的完整范围。

>>> from rasterio import windows
>>> # Full window is ((0, 1000), (0, 500))
>>> window1 = Window(10, 100, 490, 400)
>>> window2 = Window(50, 10, 200, 140)
>>> windows.union(window1, window2)
Window(col_off=10, row_off=10, width=490, height=490)
>>> windows.intersection(window1, window2)
Window(col_off=50, row_off=100, width=200, height=50)

阻碍

栅格数据集通常由多个数据块组成，当窗口与数据集自身的块结构匹配时，窗口读写效率最高。打开文件进行读取时，任何带区的块形状都可以从Block_Shapes属性中获得。

>>> with rasterio.open('tests/data/RGB.byte.tif') as src:
...     for i, shape in enumerate(src.block_shapes, 1):
...         print((i, shape))
...
(1, (3, 791))
(2, (3, 791))
(3, (3, 791))

块窗口本身可以从块窗口功能中获得。

>>> with rasterio.open('tests/data/RGB.byte.tif') as src:
...     for ji, window in src.block_windows(1):
...         print((ji, window))
...
((0, 0), ((0, 3), (0, 791)))
((1, 0), ((3, 6), (0, 791)))
...

此函数返回一个迭代器，该迭代器生成一对值。第二个是一个窗口元组，可用于调用 read 或 write . 第一个是数据集所有块中该块的行索引和列索引对。

您可以像这样从一个文件块读取数据窗口。

>>> with rasterio.open('tests/data/RGB.byte.tif') as src:
...     for ji, window in src.block_windows(1):
...         r = src.read(1, window=window)
...         print(r.shape)
...         break
...
(3, 791)

培育良好的文件具有相同的阻塞带，但gdal允许使用其他方法，最好在代码中测试这种假设。

>>> with rasterio.open('tests/data/RGB.byte.tif') as src:
...     assert len(set(src.block_shapes)) == 1
...     for ji, window in src.block_windows(1):
...         b, g, r = (src.read(k, window=window) for k in (1, 2, 3))
...         print((ji, r.shape, g.shape, b.shape))
...         break
...
((0, 0), (3, 791), (3, 791), (3, 791))

Block_Shapes属性是块形状和 set(src.block_shapes) 为您提供一组独特的形状。断言集合中只有一个项实际上与断言所有带区都具有相同的块结构相同。如果有，您可以对每个窗口使用相同的窗口。