ErmIg3 дек 2013 в 13:26

Оптимизация обработки изображений на C++ с использованием SIMD. Медианный фильтр

8 мин

31K

C++ * Алгоритмы * Обработка изображений *

Туториал

+20

Комментарии 28

1dash 3 дек 2013 в 13:45

Не сортировать, а находить сразу медиану не пробовали? e-maxx.ru/algo/kth_order_statistics

ErmIg 3 дек 2013 в 13:53

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

leventov 3 дек 2013 в 14:08

За счет чего получилось ускорение 40-60 раз, если за операцию обрабатывается 16 или 32 точки. Разве это не будет теоретическим пределом ускорения?

1dash 3 дек 2013 в 15:04

Не будет, в векторной версии есть замечательные быстрые инструкции mm256_min_*, mm256_max_*

leventov 3 дек 2013 в 15:11

Как-то странно, сделать инструкцию (причем очень полезную) в векторном варианте, но не сделать в скалярном.

encyclopedist 3 дек 2013 в 19:55

В скалярном варианте есть conditional move, которые тоже позволяют избавиться от условных переходов, хоть и не так красиво конечно.

mejedi 3 дек 2013 в 15:22

В бенчмарке не хватает версии на OpenCL.

encyclopedist 3 дек 2013 в 19:42

На всякий случай, вдруг кто не знает — в стандартной библиотеке C++ есть функция std::nth_element для частичной сортировки массива. (Конечно, оне не для этого конкретного случая)

tp7 3 дек 2013 в 20:14

Спасибо за статью!
Интересно было бы реализацию для фильтра с произвольным радиусом. Понятно, что оптимальные сети сортировки для маленьких радиусов известны и могут быть легко реализованы по-отдельности, а что делать в SIMD с большими областями?

tp7 3 дек 2013 в 20:39

Ну и насчет переключения режима работы процессора — есть такая инструкция, как vzeroupper, которая позволяет этого гигантского штрафа от переключения избежать. Конечно, писать вмеремешку SSE2 и AVX код всё равно не стоит, но и ужасных >100 циклов вы не потеряете.

Gumanoid 3 дек 2013 в 21:33

Алгоритмы были скомпилированы с максимальной оптимизацией

Значит результаты 1 и 2 — с автовекторизацией? Или для них был выключен векторизатор?

encyclopedist 3 дек 2013 в 21:48

Здесь использовался MSVC, и он видимо не очень справился с оптимизацией. На GCC результаты будут гораздо лучше.

encyclopedist 3 дек 2013 в 23:35

Я провёл стравнение на разных компиляторах и разных типах данных (автор почему-то выбрал int для промежуточных данных, что мешает автоматической векторизации). Методика аналогична авторской. Также добывил nth_element и вариант и испозованием std::min/max (T t = a; a = std::min(t, b); b = std::max(t, b);)

1   std::sort
2   std::nth_element
3   if
4   std::min std::max
5   bits
6   sse2

Компилятор  T           1       2       3       4       5       6
GCC native  int         260     305      11.3    11.3   18      2.82
            uint16_t    -       -         5.3     5.3    7.7    2.82
            uint8_t     -       -         2.79    2.82  -       2.83
GCC         int         -       -        21      21.2   17.7    2.79
            uint16_t    -       -         8.36    8.56   7.95   2.8
            uint8_t     -       -         2.85    2.85  -       2.8
clang       int         234     309     155      40.9   82.2    2.92
            uint16_t    -       -       155      40.9   79.5    2.92
            uint8_t     -       -       155     184     -       2.92

CXXFLAGS: -O3 -std=c++11

GCC:     gcc version 4.8.1 (Ubuntu/Linaro 4.8.1-10ubuntu9)
native:  -march=native
clang:   Ubuntu clang version 3.4-1ubuntu1 (trunk) (based on LLVM 3.4)

CPU:     Core2 Duo P8600 2.4 GHz
OS:      Ubuntu 13.10 64-bit

Выводы:
— Результаты очень сильно зависят от компилятора
— GCC способен векторизовывать такие вычисления
— Способность к векторизации сильно зависит от типа данных. Не следует использовать слишком длинные типы данных.
— Наилучшие результаты получаются при T=uint8_t с ключом -march=native. В этом случае автоматическая векторизация даёт ту же производительность, что и ручная.
— GCC для кода с if даёт код примерно той-же производительности, что и для битовых манипуляций. Условные переходы он заменяет на условные операции (cmov)
— clang c такими опциями с векторизацией справиться пока не способен
— В данной задаче nth_element медленнее чем sort, и оба они на порядки медленнее по сранению с другими методами на GCC (но не на clang).

ErmIg 4 дек 2013 в 07:22

Спасибо за такой развернутый комментарий. Честно говоря, я не знал что уже GCC дорос до автовекторизации подобных выражений. Буду экспереметировать и курить доки.

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

ErmIg 4 дек 2013 в 08:32

P.S. Проверил. Пока автовекторизация в большинстве перечисленных случаев не работает. Возможно в будущем это изменится.

Gumanoid 4 дек 2013 в 13:09

Осталось проверить эти примеры на ICC :)

DmitryBabokin 4 дек 2013 в 15:13

Для автовекторизации лучше использовать Intel C/C++ Compiler, его автовекторизатор на порядок продвинутей MSVS/gcc/clang.

encyclopedist 5 дек 2013 в 00:29

Я тоже так думал. Но сегодня я специально прверил и ICC:

ICC         int         -       -       132     240     66.7    2.88
            uint16_t    -       -       156     260     74.4    2.79
            uint8_t     -       -       158     268       -     2.97
ICC native  int         -       -       132     111     66.7    2.8
            uint16_t    -       -       156     109     74.5    2.8
            uint8_t     -       -       158     115       -     2.8

icc: icc version 14.0.1 (gcc version 4.8.0 compatibility)

То есть ICC с векторизацией не справился.

Вот что он выдал с ключом -vec-report2 (для программы с if):

scalar-if.cpp(46): (col. 9) remark: loop was not vectorized: existence of vector dependence
scalar-if.cpp(44): (col. 5) remark: loop was not vectorized: not inner loop

а вот что GCC c ключом -ftree-vectorizer-verbose=1:

Analyzing loop at scalar-if.cpp:44

Analyzing loop at scalar-if.cpp:46

Vectorizing loop at scalar-if.cpp:46

scalar-if.cpp:46: note: create runtime check for data references MEM[(const uint8_t *)_28 + -1B] and *_17
scalar-if.cpp:46: note: create runtime check for data references *_28 and *_17
scalar-if.cpp:46: note: create runtime check for data references MEM[(const uint8_t *)_28 + 1B] and *_17
scalar-if.cpp:46: note: create runtime check for data references MEM[(const uint8_t *)_13 + -1B] and *_17
scalar-if.cpp:46: note: create runtime check for data references *_13 and *_17
scalar-if.cpp:46: note: create runtime check for data references MEM[(const uint8_t *)_13 + 1B] and *_17
scalar-if.cpp:46: note: create runtime check for data references MEM[(const uint8_t *)_41 + -1B] and *_17
scalar-if.cpp:46: note: create runtime check for data references *_41 and *_17
scalar-if.cpp:46: note: create runtime check for data references MEM[(const uint8_t *)_41 + 1B] and *_17
scalar-if.cpp:46: note: created 9 versioning for alias checks.

scalar-if.cpp:46: note: === vect_do_peeling_for_loop_bound ===Setting upper bound of nb iterations for epilogue loop to 14

scalar-if.cpp:46: note: LOOP VECTORIZED.
scalar-if.cpp:41: note: vectorized 1 loops in function.

То есть ICC увиел какие-то зависимости и успокоился, а GCC вставил рантайм проверки (и видимо создал несколько вариантов цикла?)

DmitryBabokin 5 дек 2013 в 08:47

#pragma ivdep
перед внутренним циклом спасёт отца русской демократии :)

encyclopedist 5 дек 2013 в 23:08

Она ускоряет в полтора раза примерно, но всё равно сильно медленнее чем GCC.

AndreyDmitriev 3 дек 2013 в 22:22

В конце прошлого века помнится, реализовывал медианный фильтр на MMX. Положил кучу времени и сил, а результат по скорости был ненамного лучше того, что получалось на процессоре. Сейчас, конечно, прогресс ушёл далеко вперёд — на современных инструкциях работать куда как приятнее и легче. Ещё несколько лет назад понадобился мне медианный фильтр с очень большим ядром — 15х15 и больше, до 31х31. Наиболее оптимальным оказался алгоритм Хуанга. Ну или вот ещё — Median Filtering in Constant Time.

ErmIg 4 дек 2013 в 07:33

Реализовывал медианный фильтр по алгоритму из последней вашей ссылки. Да он позволяет обрабатывать изображения с одинаковой скоростью не зависимо от размера окна. Но эта постоянная скорость все же очень низкая (около секунды для HD разрешения, на сколько я помню). В области где я работаю (видеоаналитика), такая скорость неприемлема. Да и не требуются для нас медианные фильтры большого радиуса. максимум 5x5.

AndreyDmitriev 4 дек 2013 в 07:59

Да, он имеет выигрыш только на больших ядрах (где-то начиная от 9x9 и выше), кроме того довольно чувствителен к сильным перепадам яркости, поскольку ему по гистограмме слишком много бегать придётся.