Векторные инструкции. Часть 2: Векторизация
В качестве входных данных во всех примерах выступят блоки пикселей некоторого изображения. Для простоты возьмем монохромное изображение, а значения пикселей пусть лежат в диапазоне 0...2^bitdepth — 1, где bitdepth — разрядность пикселей. Изображение представлено одномерным массивом; в примерах значение разрядности либо равно 8 бит, либо лежит в диапазоне 9–16 бит включительно. Следовательно, изображение представлено либо массивом байтов, либо массивом беззнаковых 16-битных целочисленных переменных.
Программный код намеренно упрощён: во всех примерах блоки изображения квадратные, и их размер равен 4, 8 или 16 пикселям. Каждая функция работает с блоками одного выбранного размера и выбранной разрядности пикселей.
Во всех примерах приняты следующие обозначения: src и dst — указатели на левый верхний угол блока-источника и блока-приёмника соответственно; stride — действительная длина одной строки изображения в пикселях, которая может быть равной видимой ширине изображения, либо может превышать её (например, для выравнивания, чтобы адрес пикселя в начале любой строки был кратен 16); bitdepth — разрядность пикселей изображения в том случае, если она отличается от 8 бит.
Копирование блоков
Одна из самых простых задач обработки изображения — копирование фрагмента изображения. Ниже приведена реализация функции копирования блока 8-битного изображения, оформленная в виде шаблона C++.
Пример 1. Копирование данных из src в dst
template void copy_mb(const uint8_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { for(int i = 0; i < SIZE; i++) { for(int j = 0; j < SIZE; j++) { dst[ j ] = src[ j ]; } src += src_stride; dst += dst_stride; } } Копирование происходит во внутреннем цикле, который можно заменить, например, вызовом стандартной функции memset. Но вместо вызова этой функции воспользуемся векторными инструкциями, которые считывают данные из оперативной памяти в регистр, а затем записывают обратно в оперативную память. Для копирования блоков данных в наборе инструкций SSE2 имеются три пары инструкций: _mm_cvtsi32_si128 и _mm_cvtsi128_si32 для блоков длиной 4 байта, _mm_loadl_epi64 и _mm_storel_epi64 для блоков длиной 8 байт, _mm_loadu_si128 и _mm_storeu_si128 для блоков длиной 16 байт. В последнем случае также можно воспользоваться парой _mm_load_si128 и _mm_store_si128 cпри условии, что адрес чтения и записи выровнен по границе (кратен) 16 байтам. Получаем три реализации функции копирования для блоков разных размеров. Пример 2. Копирование данных из src в dst (с векторизацией) #include void copy_mb_4(const uint8_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { __m128i x0; for(int i = 0; i < 4; i++) { x0 = _mm_cvtsi32_si128(*(int32_t*) src); *(int32_t*) dst = _mm_cvtsi128_si32(x0); src += src_stride; dst += dst_stride; } }// copy_mb_4 void copy_mb_8(const uint8_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { __m128i x0; for(int i = 0; i < 8; i++) { x0 = _mm_loadl_epi64((__m128i*)src); _mm_storel_epi64((_m128i*)dst, x0); src += src_stride; dst += dst_stride; } }// copy_mb_8 void copy_mb_16(const uint8_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { __m128i x0; for(int i = 0; i < 16; i++) { x0 = _mm_loadu_si128((__m128i*)src); _mm_storeu_si128((_m128i*)dst, x0); src += src_stride; dst += dst_stride; } }// copy_mb_16 Можно избавиться и от внешнего цикла, как показано ниже в примере 3. Такое «развёртывание цикла» — известный приём, полезный в тех случаях, когда количество итераций невелико и заранее известно, а тело цикла сравнительно небольшое. Пример 3. Копирование данных из src в dst (с векторизацией, без циклов) #include void copy_mb_4(const uint8_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { __m128i x0, x1, x2, x3; x0 = _mm_cvtsi32_si128(*(int32_t*)(src + 0 * src_stride)); x1 = _mm_cvtsi32_si128(*(int32_t*)(src + 1 * src_stride)); x2 = _mm_cvtsi32_si128(*(int32_t*)(src + 2 * src_stride)); x3 = _mm_cvtsi32_si128(*(int32_t*)(src + 3 * src_stride)); *(int32_t*)(dst + 0 * dst_stride) = _mm_cvtsi128_si32(x0); *(int32_t*)(dst + 1 * dst_stride) = _mm_cvtsi128_si32(x1); *(int32_t*)(dst + 2 * dst_stride) = _mm_cvtsi128_si32(x2); *(int32_t*)(dst + 3 * dst_stride) = _mm_cvtsi128_si32(x3); } Функции копирования блоков разрядностью 9–16 бит реализуются аналогично тому, как это сделано для восьмибитных блоков. Нужно просто воспользоваться инструкцией, копирующей большее количество байтов (например, _mm_loadl_epi64 вместо _mm_cvtsi32_si128), либо одной и той же инструкцией несколько раз (лучше всего для блоков 16×16 пикселей). Рассмотрим сложение пикселей двух блоков. Это неотъемлемая часть декодирования видео, в ходе которого один блок пикселей (в нашем случае dst) вычисляется путём интерполяции из пикселей текущего или предыдущих кадров, а другой блок (в нашем случае src) вычисляется с помощью обратного дискретно-косинусного преобразования из декодированных коэффициентов. Вычисление блоков src и dst в примерах не показано, а приведена только заключительная компенсация, то есть сложение значений пикселей указанных блоков и запись результата в блок dst. В простейшем случае (одинаковый тип переменных блоков пикселей и отсутствие целочисленного переполнения) компенсация реализуется следующим образом. Пример 4. Компенсация пикселей (16 бит) template void compensate_mb(const uint16_t * src, uint16_t * dst, size_t src_stride, size_t dst_stride) { for(int i = 0; i < SIZE; i++) { for(int j = 0; j < SIZE; j++) { dst[ j ] = dst[ j ] + src[ j ]; } src += src_stride; dst += dst_stride; } } Инструкции копирования блоков данных мы уже рассматривали ранее. Дополним их инструкцией сложения (в нашем случае это_mm_add_epi16). Тогда векторный вариант кода выглядит как в следующем примере. Для примера взят блок 8×8 пикселей. Пример 5. Компенсация пикселей (16 бит), блок 8×8, с векторизацией void compensate_8(const uint16_t * src uint16_t * dst, size_t src_stride, size_t dst_stride) { __m128i x0, x1; for(int i = 0; i < 8; i++) { x0 = _mm_loadu_si128((__m128i*)src); // 8 pixels x1 = _mm_loadu_si128((__m128i*)dst); x0 = _mm_add_epi16( x0, x1); _mm_storeu_si128((_m128i*)dst, x0); src += src_stride; dst += dst_stride; } } На практике, при декодировании видео разрядность блока, полученного путём обратного дискретно-косинусного преобразования (src), превышает разрядность блока, полученного путём интерполяции (dst). Кроме того, элементы блока src могут принимать и отрицательные значения. В этом случае более реалистичный вариант компенсации выглядит следующим образом. Пример 6. «Реалистичная» компенсация пикселей (16 и 8 бит) template void compensate_mb(const int16_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { for(int i = 0; i < SIZE; i++) { for(int j = 0; j < SIZE; j++) { int tmp = dst[ j ] + src[ j ]; if(tmp > MAX(uint8_t)) dst[ j ] = MAX(uint8_t); else if (tmp < 0) dst[ j ] = 0; else dst[ j ] = tmp; } src += src_stride; dst += dst_stride; } } Здесь значение переменной tmp не должно выходить за пределы допустимых значений пикселей dst (в данном случае 0...255), оно ограничивается верхним либо нижним пределом. Векторная реализация по сравнению с Примером 5 усложняется: необходимо преобразовать значения dst из 8-битных в 16-битные со знаком или без знака. Необходимо также реализовать обратное преобразование в 8-битное значение без знака с ограничениями сверху и снизу. Преобразование беззнаковых целых чисел в целые с вдвое большей разрядностью легче всего осуществляется при помощи инструкций _mm_unpacklo_epiX и_mm_unpackhi_epiX (X = 8, 16, 32, or 64). Например, вектор 8-битных беззнаковых целых чисел преобразуется в два вектора 16-битных целых следующим образом: zero = _mm_setzero_si128(); x1 = x0; x0 = _mm_unpacklo_epi8(x0, zero); x1 = _mm_unpackhi_epi8(x1, zero); Аналогично преобразуются величины большей разрядности. Чтобы преобразовать 16-битные величины в беззнаковые 8-битные, используется инструкция _mm_packus_epi16 , которая упаковывает содержимое двух векторных регистров в один. Она также производит усечение по границам диапазона 0...255, если значение какого-нибудь 16-битного элемента выходит за указанные границы. Тогда векторная реализация для блока 8×8 выглядит, как показано в Примере 7. Пример 7. «Реалистичная» компенсация пикселей (8 бит), блок 8×8, с векторизацией void compensate_8(const int16_t * src, uint8_t * dst, size_t src_stride, size_t dst_stride) { __m128i x0, x1, zero; zero = _mm_setzero_si128(); for(int i = 0; i < 8; i++) { x0 = _mm_loadu_si128((__m128i*)src); // 8 pixels x1 = _mm_loadl_epi64((__m128i*)dst); // 8 bit ! x1 = _mm_unpacklo_epi8(x1, zero); // from 8 to 16 bit x0 = _mm_add_epi16( x0, x1); x0 = _mm_packus_epi16(x0, x0); // back to 8 bit _mm_storel_epi64((_m128i*)dst, x0); src += src_stride; dst += dst_stride; } } Пусть теперь разрядность пикселей изображения (dst) равна некоторому значению из диапазона 9–16 бит, а блок src, вычисленный путём обратного дискретно-косинусного преобразования представлен 32-битными значениями. Преобразовать 32-битную величину в величину из диапазона 0...2bitdepth−1 не так просто, так как нет инструкции, аналогичной _mm_packus_epi16, но для произвольной разрядности данных. Преобразование выполняется в два шага: сначала с помощью инструкции _mm_packs_epi32 два вектора 32-битных целых со знаком упаковываются в один вектор, содержащий 16-битные целые со знаком (диапазон — 32768..32767). Затем, при помощи инструкций _mm_min_epi16 и _mm_max_epi16 (сравнивают пару соответствующих элементов двух регистров и выбирают минимальную и максимальную соответственно) производится усечение значений, лежащих за пределами диапазона. (Обратите внимание, если воспользоваться инструкцией_mm_packus_epi32 iиз набора SSE4.1, которая преобразовывает 32-битные величины в беззнаковые 16-битные, то достаточно всего одного сравнения). Полностью код для блока 4×4 выглядит приведён далее (предполагается, что переполнения при 32-битном сложении не происходит). Пример 8. «Реалистичная» компенсация пикселей (9-16 бит), блок 4×4, с векторизацией void compensate_4(const int32_t * src, uint16_t *dst, size_t src_stride, size_t dst_stride, int bitdepth) { __m128i x0, x1, zero, max_val; zero = _mm_setzero_si128(); max_val = _mm_set1_epi16((1 << bitdepth) — 1); for(int i = 0; i < 4; i++) { x0 = _mm_loadu_si128((__m128i*)src); // 4×32 x1 = _mm_loadl_epi64((__m128i*)dst); // 4×16 x1 = _mm_unpacklo_epi16(x1, zero); // from 16 bit to 32 bit x0 = _mm_add_epi32(x0, x1); x0 = _mm_packs_epi32( x0, x0 ); // from 32 bit to 16 bit /* if x0[k] < max_val, then x0[k]. else max_val */ x0 = _mm_min_epi16(x0, max_val); x0 = _mm_max_epi16(x0, zero); _mm_storel_epi64((_m128i*)dst, x0); src += src_stride; dst += dst_stride; } } В функции используется _mm_set1_epi16, которая устанавливает значения 16-битных элементов векторного регистра равными одной и той же величине. В действительности _mm_set1_epi16 является псевдоинструкцией (как и другие, аналогичные ей), реализация которой зависит от компилятора. Когда требуется высокая производительность, следует избегать подобных псевдоинструкций. Чтобы определить, насколько отличаются изображения, применяют метрики — выражения, использующие значения пикселей обоих изображений. С помощью метрик можно определить оптимальный способ кодирования. Рассмотрим вычисление двух метрик: суммы разностей абсолютных значений (SAD) и средне-квадратичной ошибки (MSE). Для двух изображений одинакового размера они определяются следующими формулами:Компенсация
Вычисление метрик