эффективная реализация фильтра 2D-средних, минимизирующая избыточные нагрузки на память?

Question

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

Позвольте мне объяснить на практическом примере ... Предположим, чтоВы хотите выполнить медианный фильтр для большой двумерной матрицы с размером ядра (окна) 3x3.Первая позиция ядра (красная на изображении ниже) будет центрирована в (1,1), вторая позиция (зеленая) будет в центре (1,2).Обратите внимание на то, как желтая область перекрывается, и теперь эти значения необходимо перезагружать из памяти.

meanfilter http://luka.s3.amazonaws.com/meanfilter.png

Моя конкретная проблема - средний 3D-фильтр, поэтому перекрытие еще больше(3 ^ 3-3 ^ 2 = 18 для 3D против 3 ^ 2-3 = 6 для 2D).

Я уверен, что это общая проблема ... Кто-нибудь знает, как такие алгоритмы, как этотреализованы эффективно для устранения избыточных поисков памяти или для использования пространственной и временной локализации кэша ЦП на современных архитектурах (например, 2-сторонний ассоциативный кэш)?

Моя конкретная проблема в 3D требует только среднего значения от6 ближайших соседей (не диагональных) и реализованы в C следующим образом:

for( i = 0; i <= maxi; i++ ) {
    for( j = 0; j <= maxj; j++ ) {
        for( k = 0; k <= maxk; k++ ) {
            filteredData[ i ][ j ][ k ] = 
            ONE_SIXTH *
            ( 
             data[ i + 1 ][ j     ][ k     ] +
             data[ i - 1 ][ j     ][ k     ] +
             data[ i     ][ j + 1 ][ k     ] +
             data[ i     ][ j - 1 ][ k     ] +
             data[ i     ][ j     ][ k + 1 ] +
             data[ i     ][ j     ][ k - 1 ]
            );
        }
    }
}

Answer 1

То, что вы делаете, называется Свертка . Вы объединяете многомерные данные с меньшим ядром с таким же количеством измерений. Это очень распространенная задача, и для нее есть множество библиотек.

Быстрое решение ( в зависимости от размера ядра ) заключается в вычислении свертки в частотной области. Вы вычисляете (многомерное) БПФ как данных, так и ядра, умножаете их и вычисляете обратное БПФ. Вы найдете библиотеки, оптимизированные для этого, например. для Python есть scipy.ndimage.filters.convolve и scipy.signal.fftconvolve .

Черепица - это распространенный метод обработки изображений для оптимизации низкоуровневого доступа к памяти. Вы выделяете квадратные плитки (или кубы), которые хорошо вписываются в кэш процессора. При доступе к соседним пикселям они большую часть времени будут находиться близко друг к другу в памяти. Однако цикл по всему массиву становится немного сложнее.

Для дальнейшего чтения я рекомендую статью Почему современные процессоры голодают и что с этим можно сделать , в которой упоминается этот метод блокировки памяти и указываются на числовые библиотеки, которые его реализуют.

И, наконец, есть Интегральное изображение , которое позволяет вычислять среднее произвольного прямоугольника / кубоида с очень небольшим числом обращений к памяти.

Answer 2

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

for (i = 1; i < M - 1; ++i)
{
    // init first two column sums
    col0 = a[i - 1][0] + a[i][0] + a[i + 1][0];
    col1 = a[i - 1][1] + a[i][1] + a[i + 1][1];
    for (j = 1; j < N - 1; ++j)
    {
        // calc new col sum
        col2 = a[i - 1][j + 1] + a[i][j + 1] + a[i + 1][j + 1];
        // calc new mean
        mean[i][j] = (col0 + col1 + col2) / 9;
        // shuffle col sums
        col0 = col1;
        col1 = col2;
    }
}

Это приводит только к 3 нагрузкам на точку, а не к 9, как в простом случае, но все еще не совсем оптимально.

Вы можете оптимизировать это далее, обрабатывая две строки за итерацию и поддерживая перекрывающиеся суммы столбцов для строк i и i + 1.