是否有更像 Numpy 的方式来插值我的栅格？

Question

我写了一个小函数来插入不规则采样光栅图像的时间序列，以便它们在时间上均匀分布（如下）。它工作正常，但我只是通过查看它知道我缺少一些快捷方式。我正在寻找一个 Numpy 忍者来给我一些关于如何加强我的语法的专业提示，并且可能还会获得一点性能提升。

干杯！

import numpy as np

def interp_rasters(rasters, chrons, sampleFactor = 1):
    nFrames = round(len(chrons) * sampleFactor)
    interpChrons = np.linspace(np.min(chrons), np.max(chrons), nFrames)
    frames, rows, cols, channels = rasters.shape
    interpRasters = np.zeros((nFrames, rows, cols, channels), dtype = 'uint8')
    outs = []
    for row in range(rows):
        for col in range(cols):
            for channel in range(channels):
                pixelSeries = rasters[:, row, col, channel]
                interpRasters[:, row, col, channel] = np.interp(
                    interpChrons,
                    chrons,
                    pixelSeries
                    )
    return interpRasters

Answer 1

由于要查找的 y 值必须是 1d，所以我看不到不循环 np.arrays 的方法。如果 rasters 和 interpRasters 数组被重塑为如何在函数中使用一个循环，而无需显式索引。这使我编写的测试数据的速度提高了大约 10%。

import numpy as np

frames = 10
rows = 5
cols = 10
channels = 3

np.random.seed(1234)

rasters = np.random.randint(0,256, size=(frames, rows, cols, channels))
chrons = np.random.randint(0, 256, size  = 10 )


# The original function.
def interp_rasters(rasters, chrons, sampleFactor = 1):
    nFrames = round(len(chrons) * sampleFactor)
    interpChrons = np.linspace(np.min(chrons), np.max(chrons), nFrames)
    frames, rows, cols, channels = rasters.shape
    interpRasters = np.zeros((nFrames, rows, cols, channels), dtype = 'uint8')
    outs = []
    for row in range(rows):
        for col in range(cols):
            for channel in range(channels):
                pixelSeries = rasters[:, row, col, channel]
                interpRasters[:, row, col, channel] = np.interp(
                    interpChrons,
                    chrons,
                    pixelSeries
                    )
    return interpRasters

def interp_rasters2(rasters, chrons, sampleFactor = 1):
    nFrames = round(len(chrons) * sampleFactor)
    interpChrons = np.linspace(np.min(chrons), np.max(chrons), nFrames)
    frames, rows, cols, channels = rasters.shape
    interpRasters = np.zeros((nFrames, rows, cols, channels), dtype = 'uint8')

    # Create reshaped arrays pointing to the same data 
    dat_in = rasters.reshape(frames, rows*cols*channels).T  
    # shape (r*c*c, frames)

    dat_out = interpRasters.reshape(nFrames, rows*cols*channels).T  
    # shape (r*c*c, frames)

    for pixelseries, row_out in zip(dat_in, dat_out):
        # Loop through all data in one loop.
        row_out[:] = np.interp( interpChrons, chrons, pixelseries )
    return interpRasters  
    # As dat_out and interpRasters share the same data return interpRasters

print(np.isclose(interp_rasters(rasters, chrons), interp_rasters2(rasters, chrons)).all())
# True  # The results are the same from the two functions.

%timeit interp_rasters(rasters, chrons)
# 568 µs ± 2.55 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%timeit interp_rasters2(rasters, chrons)
# 520 µs ± 239 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

编辑还有np.apply_along_axis。这删除了任何显式 for 循环，减少了代码量，但比以前的解决方案慢。

def interp_rasters3(rasters, chrons, sampleFactor = 1):
    nFrames = round(len(chrons) * sampleFactor)
    interpChrons = np.linspace(np.min(chrons), np.max(chrons), nFrames)

    def func( arr ):  # Define the function to apply along the axis
        return np.interp( interpChrons, chrons, arr )

    return np.apply_along_axis( func, 0, rasters ).astype( np.uint8 )

print(np.isclose(interp_rasters(rasters, chrons), interp_rasters3(rasters, chrons)).all())
# True

如果速度不是关键的话，我认为在 6 个月内我会比版本 1 或 2 更好地理解版本 3。

HTH

是否有更像 Numpy 的方式来插值我的栅格？

Is there a more Numpy-esque way of interpolating my rasters?

python

numpy

raster