Производительность путем расчета среднего значения растрового стека с использованием NumPy

Question

Я вычисляю среднее значение стека изображений. Я создал виртуальный растр, который представляет собой стек из 68 изображений часового-1. Я хочу вычислить среднее значение для каждого пикселя.

Вот мой код:

import os
import rasterio 
import numpy as np
import time 


outvrt = "stacked_s1.vrt"

in_ds = rasterio.open(outvrt)
print(in_ds.count)

array_mean = np.zeros(120560400, dtype=float)

for i in range(in_ds.count):
    start = time.time()
    index = i + 1
    band = in_ds.read(index)
    print(f"reading band index {i + 1}")
    flatten_band = band.flatten()
    for r in np.nditer(flatten_band):

        array_mean[r] = (array_mean[r] + r)/2
    end_time = time.time()
    hours, rem = divmod(end_time-start, 3600)
    minutes, seconds = divmod(rem,60)
    print("Finished visit of array in:")
    print(f"{int(hours)}:{int(minutes)}:{int(seconds)}")

print(len(array_mean))

Как вы можете видеть, я создал нулевой массив, который будет содержать результат вычислениясделано в итерации. Затем я перебрал каждую полосу в стеке и прочитал полосу как пустой массив. После этого я сгладил массив и использовал np.nditer для итерации по каждому значению сглаженного массива. Затем я обновил значения array_mean, добавив значения массива и разделив его на 2. Это ошибка для первой итерации, поскольку начальное значение в array_mean равно 0. Поэтому я вставил несколько операторов if, чтобы проверить это условие:

for r in np.nditer(flatten_band):
        if array_mean[r] == 0:
            array_mean[r] = r
        else:
            array_mean[r] = (array_mean[r] + r)/2

В первом коде у меня было время завершения ~ 1,40 минуты на группу. Благодаря введению этой дополнительной проверки условий у меня была производительность ~ 2,30 минут на полосу.

Должен ли я избавиться от вычисления среднего значения во время итерации и просто добавить новые значения для создания массива накопительной суммы и затем применитьnumpy.divide в конечном массиве по количеству полос (68)?

Этот подход помогает мне вычислять среднее значение временного ряда данных в моем маленьком ноутбуке, и не впадает в ошибки памятидоступ к значениям в плотном временном ряду.

Отказ от ответственности: я не написал код для изменения формы и выписал окончательный результат. Этот вопрос был задан мной также в канале GeoGraphic Stackoverflow, но пока не получил никакого ответа. Я немного помог бы с благодарностью.