Векторизация функции Python, которая размывает изображения

Question

Я сделал функцию, используя чистый Python для размытия любого изображения, которое я даю.Теперь я хочу использовать NumPy и удалить все циклы for.

def blur_image(src, dst):
    (h, w, c) = src.shape

    for x in range(h-1):
        for y in range(w-1):
            for z in range(c):
                if x != 0 or y != 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x-1, y, z]
                    + src[x+1, y, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x-1, y-1, z]
                    + src[x-1, y+1, z]
                    + src[x+1, y-1, z]
                    + src[x+1, y+1, z]) / 9
                if x == 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x, y, z]
                    + src[x+1, y, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x+1, y-1, z]
                    + src[x+1, y+1, z]) / 9
                if y == 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x-1, y, z]
                    + src[x+1, y, z]
                    + src[x, y, z]
                    + src[x, y+1, z]
                    + src[x-1, y, z]
                    + src[x-1, y+1, z]
                    + src[x+1, y, z]
                    + src[x+1, y+1, z]) / 9
    return dst

Как я могу улучшить этот код с помощью NumPy?Лучший вариант - удалить все петли for?Любые советы?

РЕДАКТИРОВАТЬ: Я пытаюсь решить эту проблему без каких-либо библиотек, чтобы помочь мне размыть мое изображение.Я размываю его, используя чистый python, и теперь ищу подсказки, чтобы все вычислительно тяжелые биты использовали массивы numpy.Не имею большого опыта использования numpy: //

Answer 1

Если я правильно читаю ваш код, вы описываете фильтрацию с помощью фильтра 3 × 3.Если вы заботитесь о производительности, вероятно, для этого стоит использовать существующую библиотеку, поскольку она будет гораздо более оптимизирована, чем все, что вы можете легко реализовать.

Например, то, что вы делаете, может быть достигнутосо сверткой, которая доступна в scipy.signal.Таким образом, вы можете сделать это для каждого канала, а затем сложить получившиеся массивы:

import numpy as np
import scipy.signal as ss

arr = np.random.random((100, 100))  # Some fake data.
kernel = np.ones((3, 3)) / 9  # The 'boxcar'.

ss.convolve(arr, kernel, mode='same')

Размытие также доступно в PIL:

import numpy as np
from PIL import Image, ImageFilter

arr = np.random.randint(0, 256, size=((100, 100, 3)), dtype=np.uint8)
img = Image.fromarray(arr)

img.filter(ImageFilter.BoxBlur(radius=3))

Естьмножество других способов сделать это, в том числе в пространстве Фурье ( один пример ).

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

Answer 2

Просто скомпилируйте

Я предполагаю, что вы хотите ускорить вашу функцию.Самый простой и, вероятно, самый быстрый способ сделать это - использовать компилятор, такой как Cython или Numba.

Проблема также смущающе параллельна.

Пример

import numpy as np
import numba as nb

@nb.njit(error_model="numpy",parallel=True)
def blur_image_nb(src, dst):
    (h, w, c) = src.shape
    for x in nb.prange(h-1):
        for y in range(w-1):
            for z in range(c):
                if x != 0 or y != 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x-1, y, z]
                    + src[x+1, y, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x-1, y-1, z]
                    + src[x-1, y+1, z]
                    + src[x+1, y-1, z]
                    + src[x+1, y+1, z]) / 9
                if x == 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x, y, z]
                    + src[x+1, y, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x+1, y-1, z]
                    + src[x+1, y+1, z]) / 9
                if y == 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x-1, y, z]
                    + src[x+1, y, z]
                    + src[x, y, z]
                    + src[x, y+1, z]
                    + src[x-1, y, z]
                    + src[x-1, y+1, z]
                    + src[x+1, y, z]
                    + src[x+1, y+1, z]) / 9
    return dst

Редактировать

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

#You can achieve the same with fastmath=True
#Divisions are quite costly
@nb.njit(error_model="numpy",parallel=True)
def blur_image_nb_2(src, dst):
    (h, w, c) = src.shape
    fact=1./9.
    for x in nb.prange(h-1):
        for y in range(w-1):
            for z in range(c):
                if x != 0 or y != 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x-1, y, z]
                    + src[x+1, y, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x-1, y-1, z]
                    + src[x-1, y+1, z]
                    + src[x+1, y-1, z]
                    + src[x+1, y+1, z]) *fact
                if x == 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x, y, z]
                    + src[x+1, y, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x, y-1, z]
                    + src[x, y+1, z]
                    + src[x+1, y-1, z]
                    + src[x+1, y+1, z]) *fact
                if y == 0:
                    dst[x, y, z] = (src[x, y, z]
                    + src[x-1, y, z]
                    + src[x+1, y, z]
                    + src[x, y, z]
                    + src[x, y+1, z]
                    + src[x-1, y, z]
                    + src[x-1, y+1, z]
                    + src[x+1, y, z]
                    + src[x+1, y+1, z]) *fact
    return dst

Сроки

src=np.random.rand(1024,1024,3)
dst=np.empty((1024,1024,3))

%timeit blur_image(src, dst)
8.08 s ± 52.1 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
%timeit blur_image_nb(src, dst)
5.19 ms ± 954 µs per loop (mean ± std. dev. of 7 runs, 1 loop each)
%timeit blur_image_nb_2(src, dst)
3.13 ms ± 167 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

#Just for comparison
%timeit res=np.copy(src)
11.1 ms ± 28.6 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
%timeit np.copyto(dst, src)
2.44 ms ± 53.2 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

import numpy as np
from PIL import Image, ImageFilter
import scipy.signal as ss

src=(np.random.rand(1024,1024,3)*255).astype(np.uint8)
img = Image.fromarray(src.astype(np.uint8))

%timeit img.filter(ImageFilter.BoxBlur(radius=3))
16.1 ms ± 65 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

kernel = np.ones((3, 3)) / 9  # The 'boxcar'.
#3x in a example where src.shape[2]==3
%timeit ss.convolve(src[:,:,0], kernel, mode='same')
76.6 ms ± 16.6 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

Вывод

Обычно вы просто используете функции, которые уже доступны, как @kwinkunks показал в своем ответе, но простое выполнение вашей функции показывает лучшую производительность, чем PIL (которая также имеет некоторые ограничения в отношении типа данных), и более общую свертку.