解释Numpy数组中的特定位

Question

我正在尝试解释HDF文件，或更准确地解释文件中的Rasterband。有一个用于质量评估的频段，它将位信息表示为相应的整数（例如01000000为64）。

取决于具体位，我想得到一个计数，例如位置5上有多少个像素为1。如果以前的计数考虑了像素，则不应再次计数。现在，我将根据自己的优先级列表更改每个单元格中的条目。这需要花很多时间，而且我真的可以确定有更快的方法，因为我以前从未使用过钻头。

这是我的代码：

from osgeo import gdal
import numpy as np
QA_Band = gdal.Open(hdf.GetSubDatasets()[B][0],gdal.GA_ReadOnly)
QA = QA_Band.ReadAsArray()

# Calculate Bit-Representation of QA-Band
bin_vec = np.vectorize(np.binary_repr)
QAB = bin_vec(QA, width = 8)

# Filter by Bit-Values
QAB = np.where(QAB == '11111111', "OutOfSwath", QAB)

for i in range(QAB.shape[0]):
    for j in range(QAB.shape[1]):
        if QAB[i,j][6] == '1':
            QAB[i,j] = "Cloud"
            Cloud = (QAB == "Cloud").sum()
        elif QAB[i,j][4] == '1':
            QAB[i,j] = "Cloud Shadow"      
            Shadow = (QAB == "Cloud Shadow").sum()
        elif QAB[i,j][5] == '1':
            QAB[i,j] = "Adjacent Cloud"
            AC = (QAB == "Adjacent Cloud").sum()

        elif QAB[i,j][7] == '1':
            QAB[i,j] = "Cirrus"
            Cirrus = (QAB == "Cirrus").sum()

        elif QAB[i,j][3] == '1':
            QAB[i,j] = "Snow/Ice"
            SnowC = (QAB == "Snow/Ice").sum()
        elif QAB[i,j][2] == '1':
            QAB[i,j] = "Water"
            WaterC = (QAB == "Water").sum()

        elif QAB[i,j][0:1] == '11':
            QAB[i,j] = "High Aerosol"        
        elif QAB[i,j][0:1] == '10':
            QAB[i,j] = "Avrg. Aerosol"        
        elif QAB[i,j][0:1] == '01':
            QAB[i,j] = "Low Aerosol"
        elif QAB[i,j][0:1] == '00':
            QAB[i,j] = "Aerosol Climatology"

这将导致代表不同事物的字符串数组，但需要花费一定的时间，如前所述。

任何有关如何访问位表示的帮助都将是有用的：）

Answer 1

您可以使用numpy函数unpackbits处理位。使用numpy时，我们更喜欢使用numpy方法和函数而不是python进行循环-通常更快。因此，您可以将每个数字拆箱到第三个轴，然后将QAB[i, j][5] == '1'之类的条件转换为result[bits[5]]。我将您的elif子句的顺序颠倒过来，就好像您QAB[i, j][6] == '1'一样，然后将其设置为"Cloud"，并且从不覆盖所有子句，因此，如果我们运行每个条件，这种情况下应该是最后一次重写覆盖。同样，您最后的情况，例如QAB[i,j][0:1] == '11'从未触发，因为它的左侧始终为1。所以我用您的意思是QAB[i,j][0:2] == ...进行重写。

"""
>>> print(bits)
[[[False False False]
  [False False  True]]
<BLANKLINE>
 [[False False False]
  [False False  True]]
<BLANKLINE>
 [[False False False]
  [False False  True]]
<BLANKLINE>
 [[False False False]
  [False False  True]]
<BLANKLINE>
 [[False False False]
  [False  True  True]]
<BLANKLINE>
 [[False False False]
  [False False  True]]
<BLANKLINE>
 [[False  True  True]
  [ True False  True]]
<BLANKLINE>
 [[ True False  True]
  [ True  True  True]]]
>>> print(result)
[['Cirrus' 'Cloud' 'Cloud']
 ['Cloud' 'Cloud Shadow' 'OutOfSwath']]
>>> Cloud
3
"""
import numpy as np
QA = [[1, 2, 3], [3, 9, 255]]

bits = np.unpackbits(np.array(QA, dtype=np.uint8, ndmin=3), axis=0).astype(bool)
result = np.array([[""] * bits.shape[2] for _ in range(bits.shape[1])], dtype=object)
result[~bits[0] & ~bits[1]] = "Aerosol Climatology"
result[~bits[0] & bits[1]] = "Low Aerosol"
result[bits[0] & ~bits[1]] = "Avrg. Aerosol"
result[bits[0] & bits[1]] = "High Aerosol"
result[bits[2]] = "Water"
result[bits[3]] = "Snow/Ice"
result[bits[7]] = "Cirrus"
result[bits[5]] = "Adjacent Cloud"
result[bits[4]] = "Cloud Shadow"
result[bits[6]] = "Cloud"
result[bits.all(axis=0)] = "OutOfSwath"
Cloud = (result == "Cloud").sum()
Shadow = (result == "Cloud Shadow").sum()
AC = (result == "Adjacent Cloud").sum()
Cirrus = (result == "Cirrus").sum()
SnowC = (result == "Snow/Ice").sum()
WaterC = (result == "Water").sum()

Answer 2

您程序的逻辑实际上精确地映射到numpy.select上，该C0根据条件（布尔）数组的列表从数组列表中逐个元素地选择，第一个匹配获胜。因此，您可以紧凑地写类似

的内容

conditions = QAB&(128>>np.arange(8))[:,None,None]
values = ["OutOfSwath","Cloud","Cloud Shadow","Adjacent Cloud","Cirrus","Snow/Ice",
          "Water","High Aerosol","Avrg. Aerosol","Low Aerosol","Aerosol Climatology"]

np.select([QAB==255,*conditions[[6,4,5,7,3,2]],*(QAB==np.arange(0,256,64)[::-1,None,None])],values)