如何使用 Rasterio 更改栅格的 dtype

Question

在栅格数据集上应用多边形掩码时，我在处理 Python 栅格数据包中没有数据值时遇到了问题。这个特定的栅格是具有 7 个波段的 Landsat uint8，并且没有固有地指定无数据值，因为 255 是无数据的保留值。但是，有时 uint8 数据是从 uint16 数据压缩而来的，而 255 值是一个有效的数据值，我不想将其视为“无数据”（数据是完整的位范围）。如果未指定此参数，则 rasterio 的掩码函数的默认值是将 0 视为“无数据”值，这同样存在问题，因为 0 有时被视为有效数据值。有没有办法覆盖“无数据”的元数据值？

我尝试了几种不同的方法来解决这个问题（详见下文），但都没有成功。

1. 使用 rasterio.open() 将 uint8 数据转换为 uint16 数据并将“256”指定为无数据值，因为它超出了任何 uint8 数据的范围，但在 uint16 数据范围内被接受。这就是某些软件程序（例如 ArcMap）有时会处理不分配数据值的方式。

2. 与步骤 1 类似，但尝试使用 rasterio.open() 打开 uint8 数据并在函数中设置“nodata=np.nan”。收到错误：“给定 nodata 值 nan，超出了其数据类型的有效范围。”尽管在文档中 nan 被列为“nodata”参数的有效条目。

3. 在使用 rasterio.mask() 的掩码过程中，指定 nodata=nan。收到错误“无法将 fill_value nan 转换为 dtype。”

import rasterio
import fiona
import numpy as np

fp_src = ''
fp_dst = ''
shape = ''

# get shapes
with fiona.open(shape, 'r') as shapefile:
    geoms = [feature['geometry'] for feature in shapefile]



# Method Number 1
# ~~~~~~~~~~~~~~~~~~~~~~~~~

# open original raster, copy meta & alter dtype
with rasterio.open(fp_src) as src_dataset:
    kwds = src_dataset.profile
    kwds['dtype'] = 'uint16'
    src_meta = src_dataset.meta

# write a new raster with the copied and altered meta
with rasterio.open(fp_dst, 'w', **kwds) as dst_dataset:
    dst_meta = dst_dataset.meta

src_dataset.close()
dst_dataset.close()

img = rasterio.open(fp_dst)

# mask img and set nodata to 256 (out of the uint8 range)
out_image, out_transform = mask(img, geoms, nodata=256)
# out_image output: values outside of the geoms are 256 & values inside are 0.



# Method Number 2
# ~~~~~~~~~~~~~~~~~~~~~~~~~

# open original raster, copy meta & alter dtype
with rasterio.open(fp_src) as src_dataset:
    kwds = src_dataset.profile
    kwds['nodata'] = np.nan
    kwds['dtype'] = 'uint16'
    src_meta = src_dataset.meta

# write a new raster with the copied and altered meta
with rasterio.open(fp_dst, 'w', **kwds) as dst_dataset:
    dst_meta = dst_dataset.meta

src_dataset.close()
dst_dataset.close()

img = rasterio.open(fp_dst)

# mask img and let the mask function default to the image's newly created nodata (np.nan from inside with rastario.open...)
out_image, out_transform = mask(img, geoms)
# out_image output: nodata value, nan, is beyond the valid range of its data type



# Method Number 3
# ~~~~~~~~~~~~~~~~~~~~~~~~~

# mask img and set nodata to nan
out_image, out_transform = mask(fp_src, geoms, nodata=np.nan)
# out_image output: Cannot convert fill_value nan to dtype.

我希望看到给定多边形之外的所有像素都转换为“无数据”条目，该条目不一定是有效范围的一部分，这样脚本就不会意外地将有效值视为无数据.

Answer 1

问题在于np.nan 是一个无法转换为整数的浮点数。以下是我将如何解决这个问题：

with rasterio.open(fp_src) as src_dataset:
    meta = src_dataset.meta
    meta.update(
        {
            "nodata": np.iinfo(src_dataset.dtypes[0]).max
        }
    )
    data = src_dataset.read()

    with rasterio.open(fp_dst, 'w', **meta) as dst_dataset:
        dst_dataset.write(data)

函数np.iinfo(dtype).max 查找给定整数类型的最大值。我将此方法用于具有 1 个波段的数据集，因此它也应该适用于您的数据。如果没有，请告诉我。