使用np.select()或if else子句加速

Question

Python GIS开发人员的问题。我有一个高密度的三维点云，每个像素有73个观察值。用SFM获得云量。我正在将云中的DEM作为新的波段融合在HSI图像中，以供将来分析。

我有一个numpy数组，表示8个波段(n点，min z，maz z，sum z...)。我需要基于点云更改此数组的值。因此，我一次只能使用一个像素，这使得numpy向量化不是一个好主意。

变量为: map_z =平均海平面以上单个云点高度的浮点数组是8波段numpy数组

array = np.full((nbands, 2000, 3000), np.nan, dtype=OUTPUT_RASTER_DTYPE)

我已经尝试了两个版本的代码。一个矢量化版本和一个Else If版本。

slice_arr = array[0:4, pixel_y, pixel_x]
slice_arr[0] = np.where(slice_arr[0] == 0, 1, slice_arr[0] + 1)
slice_arr[1] = np.select([np.isnan(slice_arr[1]) == True, slice_arr[1] > map_z], [map_z, map_z], slice_arr[1])
slice_arr[2] = np.select([np.isnan(slice_arr[2]) == True, slice_arr[2] < map_z], [map_z, map_z], slice_arr[2])
slice_arr[3] = np.where(np.isnan(slice_arr[3]) == True, map_z, slice_arr[3] + map_z)
array[0:4, pixel_y, pixel_x] = slice_arr

或

slice_arr = array[0:4, pixel_y, pixel_x] # Take a slice of the array
if slice_arr[0] == 0:
    slice_arr[0] = 1
else:
    slice_arr[0] = slice_arr[0] + 1
if np.isnan(slice_arr[1]) == True:
    slice_arr[1] = map_z
elif slice_arr[1] > map_z:
    slice_arr[1] = map_z
else:
    pass
if np.isnan(slice_arr[2]) == True:  
    slice_arr[2] = map_z
elif slice_arr[2] < map_z:
    slice_arr[2] = map_z
else:
    pass
if np.isnan(slice_arr[3]) == True:  
    slice_arr[3] = map_z
else:
    slice_arr[3] = slice_arr[3] + map_z
array[0:4, pixel_y, pixel_x] = slice_arr

显然，If else版本的速度要快4倍，因为我一次只操作一个像素。这是因为map_z是一个点云，我一次只能处理一个点，而这些点只落在我的数组中的一个像素上。向量化应该在更大的数组上操作。

有没有办法使用numpy来加速这部分代码？如果数组(8,2000,3000)比(8,2000,3000)大得多，对数组(slice_arr)进行分片是否具有速度优势？

Answer 1

最后，我提高了If Else子句的速度，但把它放在一个函数中并使用numba。这将时间从2.33分钟缩短。为.42 min.！

from numba import jit
@jit(nopython=True)
def query_speed_up(slice_arr, map_z):

    if slice_arr[0] == 0:
        slice_arr[0] = 1
    else:
        slice_arr[0] = slice_arr[0] + 1

    # min z
    if np.isnan(slice_arr[1]) == True:  # not is
        slice_arr[1] = map_z
    elif slice_arr[1] > map_z:
        slice_arr[1] = map_z
    else:
        pass

    # max z
    if np.isnan(slice_arr[2]) == True:  # not is
        slice_arr[2] = map_z
    elif slice_arr[2] < map_z:
        slice_arr[2] = map_z
    else:
        pass

    # sum z
    if np.isnan(slice_arr[3]) == True:  # not is
        slice_arr[3] = map_z
    else:
        slice_arr[3] = slice_arr[3] + map_z

    return slice_arr

slice_arr = array[0:4, pixel_y, pixel_x]  # Take a slice of the array
jittest = query_speed_up(slice_arr, map_z)
array[0:4, pixel_y, pixel_x] = jittest