如何在GDAL Python中获得带数组的正确比例？

Question

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

这将提取频带4的未缩放阵列：

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()都返回None

和b4_band.GetStatistics(True,True)返回[0.0, 2829.0, 347.05880000000104, 334.8397839901348]（作为min、max、mean、stddev）

这是否有助于我提取正确的刻度？我不知所措

Answer 1

请注意，即使在文件中指定了比例/偏移，GDAL也不会自动应用它们

对于Sentinel 2，它们不是在文件中指定的，而是在元数据（xml）中指定的。您可以使用文本编辑器查看示例中使用的xml文件。并搜索@Val建议的“量化值”

它也可以从GDAL解析的元数据中检索。这可以通过dataset.GetMetadata()完成，它将返回一个命令。您还可以调用gdal.Infoutilty，这两种方法如下所示

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}