Какой эффективный способ выполнить 2D матричную фильтрацию в Python?

Question

Я пытаюсь реализовать матричную фильтрацию в Python, и пока что реализация кажется очень медленной и неэффективной.Интересно, существует ли эффективный способ выполнения такой фильтрации.

При наличии большой матрицы A и матрицы фильтрации M функция должна возвращать «повторно смешанную» матрицу R, которая получается умножением каждого элемента (i,j) из A с помощью M, то результат накладывается / вставляется в R в позиции (i, j).Ниже приведен код, который, как ожидается, сделает это.

Пример ниже занимает около 68 секунд (!) На моем компьютере, что выглядит очень неэффективно.

Я был бы очень признателен, если вымог бы порекомендовать способ ускорить эту функцию.Заранее большое спасибо!

import numpy as np
import time

nx = ny = 1500
n_mix = 50

# matrix to be filtered
A = np.random.random_sample( (nx, ny) )

# filter to be applied to each point:
M = np.random.random_sample( (2*n_mix+1, 2*n_mix+1) )

# the result is stored in "remix":
remix = np.zeros_like(A)

start = time.time()

for i in range(n_mix, nx-n_mix):
    for j in range(n_mix, ny-n_mix):
        remix[i - n_mix:i + n_mix + 1, j - n_mix:j + n_mix + 1 ] += M * A[i,j]

print remix

duration = time.time() - start
print(round(duration))

ОБНОВЛЕНИЕ

На самом деле пакет ndimage в scipy имеет общую функцию свертки, которая выполняет эту работу.Я публикую ниже 3 варианта выполнения фильтрации с соблюдением сроков.Самый быстрый - с ndimage.convolution (24 секунды против 56 и 68 другими методами).Тем не менее, это все еще кажется довольно медленным ...

import numpy as np
from scipy import ndimage
import time
import sys


def remix_function(A, M):
    n = (np.shape(M)[0]-1)/2
    R = np.zeros_like(A)

    for k in range(-n, n+1):
        for l in range(-n, n+1):
            # Ak  = np.roll(A, -k, axis = 0)
            # Akl = np.roll(Ak, -l, axis = 1)
            R += np.roll(A, (-k,-l), axis = (0,1) ) * M[n-k, n-l]
    return R

if __name__ == '__main__':
    np.set_printoptions(precision=2)

    nx = ny = 1500
    n_mix = 50
    nb = 2*n_mix+1

    # matrix to be filtered
    A = np.random.random_sample( (nx, ny) )
    # filter to be applied to each point:
    M = np.random.random_sample( (nb, nb) )


    # the result is stored in "remix":
    remix1 = np.zeros_like(A)
    remix2 = np.zeros_like(A)
    remix3 = np.zeros_like(A)


#------------------------------------------------------------------------------
# var 1
#------------------------------------------------------------------------------
    start = time.time()

    remix1 = remix_function(A, M)

    duration = time.time() - start
    print('time for var1 =', round(duration))

#------------------------------------------------------------------------------
# var 2
#------------------------------------------------------------------------------
    start = time.time()

    for i in range(n_mix, nx-n_mix):
        for j in range(n_mix, ny-n_mix):
            remix2[i - n_mix:i + n_mix + 1, j - n_mix:j + n_mix + 1 ] += M * A[i,j]


    duration = time.time() - start
    print('time for var2  =', round(duration))

#------------------------------------------------------------------------------
# var 3
#------------------------------------------------------------------------------
    start = time.time()

    remix3 = ndimage.convolve(A, M)

    duration = time.time() - start
    print('time for var3 (convolution) =', round(duration))

Answer 1

Я пока не могу комментировать посты, но проблема в двойном цикле for. Вы пытались определить функцию и затем использовать np.vectorize?