如何在 GDAL Python 中获得波段阵列的正确比例？

Question

我是在 Python 中使用 GDAL 的新手，我正在尝试使用它从 Sentinel 2 SAFE 产品中检索波段数据。我设法提取波段数组，但无法正确缩放它。

这将提取波段 4 的未缩放数组：[=20=]

import gdal

product_path = "S2B_MSIL2A_20200124T101219_N0213_R022_T33UUU_20200124T121752.SAFE"
dataset = gdal.Open(product_path + "MTD_MSIL2A.xml")
bands10m_path = dataset.GetSubDatasets()[0][0]
bands10m_dataset = gdal.Open(bands10m_path)
b4_band = bands10m_dataset.GetRasterBand(1)
b4_array = b4_band.ReadArray()

到目前为止一切顺利，但数组的数据类型为 uint16，值范围从 0 到 16896。

b4_band.GetMinimum() 和 b4_band.GetMaximum() 都 return None.

和b4_band.GetStatistics(True,True) returns [0.0, 2829.0, 347.05880000000104, 334.8397839901348]（作为最小值、最大值、平均值、标准偏差）。

这是否以某种方式帮助我提取正确的比例？我一无所知...

Answer 1

请注意，即使在文件中指定了 scale/offset，GDAL 也不会自动应用它们。

在 Sentinel 2 的情况下，它们没有在文件中指定，而是在元数据中指定 (xml)。您可以使用文本编辑器查看示例中使用的 xml 文件。并按照@Val 的建议搜索 "QUANTIFICATION_VALUE"。

它也可以从 GDAL 解析的元数据中检索。这可以用 dataset.GetMetadata() 来完成，这将 return 一个字典。您也可以调用 gdal.Info 实用程序，两种方法如下所示。

import gdal

archive = 'S2A_MSIL2A_20200126T082221_N0213_R121_T34HDH_20200126T110103.zip'

# Use a dataset
ds = gdal.Open(archive)
meta = ds.GetMetadata()
ds = None

# Alternatively use gdal.Info
r = gdal.Info(archive, format='json')
meta = r['metadata']['']

您可以过滤掉相关值，并使用类似以下内容将它们从字符串转换为浮点数：

{k: float(v) for k,v in meta.items() if k.endswith('_QUANTIFICATION_VALUE')}

# Result in:
{'AOT_QUANTIFICATION_VALUE': 1000.0,
 'BOA_QUANTIFICATION_VALUE': 10000.0,
 'WVP_QUANTIFICATION_VALUE': 1000.0}

如何在 GDAL Python 中获得波段阵列的正确比例？

How do I get the correct scale of a band array in GDAL Python?

python

gdal

sentinel2