Блог

Я пытаюсь реализовать матричную фильтрацию в 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))

  

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