xarray 如何使用 open_rasterio 加载和索引大型 GeoTIFF 文件？答案

【问题标题】：How does xarray load and index large GeoTIFF files with open_rasterio?xarray 如何使用 open_rasterio 加载和索引大型 GeoTIFF 文件？
【发布时间】：2020-06-09 16:59:43
【问题描述】：

我正在使用 xarray 包来加载和访问大型 GeoTIFF 文件 (>50GB) 的数据，它运行良好。

import xarray as xr
img = xr.open_rasterio("path/to/large_geo_tiff.tif")
pixel_value = img[0,1225, 4321]
print("The pixel value is: ", pixel_value.values.item())

但是，我想知道 xarray 是如何实际加载大型 GeoTIFF 文件的。显然它不会将整个文件加载到内存中，因为它不适合但会执行某种延迟加载。我只习惯 dask 和 Dask Arrays 将数据拆分成块，并通过将相应的块加载到内存中来轻松进行像素访问。但是函数 open_rasterio 的签名如下所示

xarray.open_rasterio(filename, parse_coordinates=None, chunks=None, cache=None, lock=None)

由于我没有定义 chunks，所以 img 不应该被分块。所以我的问题是调用pixel_value = img[0,1225, 4321] 时会发生什么，xarray 怎么能这么快地访问给定位置的像素值？

我期待任何反馈。

【问题讨论】：

标签： python-xarray

【解决方案1】：

与 xarray

相比，dask_image 似乎在读取 tiff 文件 方面并不那么聪明

提供的 GeoTIFF 不是平铺，但我只是检查了 tiff 文件格式，如下所示

TIFF File: '/Users/jan/Desktop/ms-project/processed/metashape/export/untiled.tif'
                       Mode: 'r'
    Current Image Directory: 1
           Number Of Strips: 1336
                SubFileType: Tiff.SubFileType.Default
                Photometric: Tiff.Photometric.RGB
                ImageLength: 42746
                 ImageWidth: 35366
               RowsPerStrip: 32
              BitsPerSample: 8
                Compression: Tiff.Compression.LZW
               SampleFormat: Tiff.SampleFormat.UInt
            SamplesPerPixel: 3
        PlanarConfiguration: Tiff.PlanarConfiguration.Chunky
                Orientation: Tiff.Orientation.TopLeft

重要的部分是 tif 是 stripped，其中每个条带都是连续的位图图像数据行的单独集合，这使得随机访问更加容易，并且xarray 似乎通过加载所需的数据来利用 stripped 图像 - 我只是通过加载 GeoTIFF 调用使用的后端模块来检查这一点由 xarray.open_rasterio 称为 rasterio

>>> import rasterio
>>> ortho = rasterio.open("./orthomosaic.tif")
>>> ortho.is_tiled
False
>>> ortho.block_shapes
[(32, 248489), (32, 248489), (32, 248489)]

另一方面，

dask_image 并没有利用这些条纹，而是将整个图像加载到一个 chunk 中，因此内存不足——即使你写了一个 tiled tif 文件 dask_image 不会考虑这个

>>> import dask_image.imread
>>> ortho = dask_image.imread.imread("./orthomosaic.tif")
>>> ortho
dask.array<from-value, shape=(1, 104700, 248489, 3), dtype=uint8, chunksize=(1, 104700, 248489, 3), chunktype=numpy.ndarray>

【讨论】：