8

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

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

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

  2. 与步骤 1 类似,但尝试使用 rasterio.open() in 打开 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.

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

4

1 回答 1

-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 个波段的数据集,因此它也应该适用于您的数据。如果没有,请告诉我。

于 2020-04-02T09:31:05.200 回答