Реализация би-кубического изменения размера
Я пытался закодировать би-кубический алгоритм изменения размера для растровых изображений в памяти. Я знаком с тем, как работает бикубическая интерполяция, и я использовал как статью Википедии, так и существующие реализации в качестве руководства по написанию своей собственной версии.
Следующее - моя простая реализация. Вот, bmap это vector содержащий растровые данные, и get_subpixel это просто функция, которая обрабатывает растровое изображение как трехмерный массив, состоящий из X x Y x Channel пикселей и возвращает один субпиксель в указанных координатах.
std::vector<unsigned char> bicubic_resize(
std::vector<unsigned char>& bmap, std::size_t bmap_width, std::size_t bmap_height,
std::size_t channels, std::size_t dest_width, std::size_t dest_height)
{
std::vector<unsigned char> out(dest_width * dest_height * 3);
const double tx = double(bmap_width) / dest_width;
const double ty = double(bmap_height) / dest_height;
const std::size_t row_stride = dest_width * channels;
unsigned char C[5] = { 0 };
for (unsigned i = 0; i < dest_height; ++i)
{
for (unsigned j = 0; j < dest_width; ++j)
{
const int x = int(tx * j);
const int y = int(ty * i);
const double dx = tx * j - x;
const double dy = ty * i - y;
for (int k = 0; k < 3; ++k)
{
for (int jj = 0; jj < 4; ++jj)
{
const int idx = y - 1 + jj;
unsigned char a0 = get_subpixel(bmap, idx, x, k);
unsigned char d0 = get_subpixel(bmap, idx, x - 1, k) - a0;
unsigned char d2 = get_subpixel(bmap, idx, x + 1, k) - a0;
unsigned char d3 = get_subpixel(bmap, idx, x + 2, k) - a0;
unsigned char a1 = -1.0 / 3 * d0 + d2 - 1.0 / 6 * d3;
unsigned char a2 = 1.0 / 2 * d0 + 1.0 / 2 * d2;
unsigned char a3 = -1.0 / 6 * d0 - 1.0 / 2 * d2 + 1.0 / 6 * d3;
C[jj] = a0 + a1 * dx + a2 * dx * dx + a3 * dx * dx * dx;
d0 = C[0] - C[1];
d2 = C[2] - C[1];
d3 = C[3] - C[1];
a0 = C[1];
a1 = -1.0 / 3 * d0 + d2 -1.0 / 6 * d3;
a2 = 1.0 / 2 * d0 + 1.0 / 2 * d2;
a3 = -1.0 / 6 * d0 - 1.0 / 2 * d2 + 1.0 / 6 * d3;
out[i * row_stride + j * channels + k] = a0 + a1 * dy + a2 * dy * dy + a3 * dy * dy * dy;
}
}
}
}
return out;
}
Этот код прекрасно работает для определенных размеров назначения. Например, если исходное растровое изображение 500 X 366 и размер назначения 250 x 183 Алгоритм работает отлично:
Оригинал: 
Изменен: 
Тем не менее, для некоторых других размеров назначения, таких как 100 x 73 изображение получателя искажается: 
Я просматривал код интерполяции и не вижу, что делаю неправильно.
Буду признателен за любые советы, предложения или ответы.
1 ответ
Помимо смешения с плавающей запятой и целочисленной арифметики, я подозреваю, что вы получаете переполнение / переполнение чисел с некоторыми из ваших промежуточных значений.
Одно простое решение - просто быть последовательным и использовать с плавающей точкой всюду. Прямо сейчас у вас есть:
unsigned char C[5] = { 0 };
for (unsigned i = 0; i < dest_height; ++i)
{
for (unsigned j = 0; j < dest_width; ++j)
{
const int x = int(tx * j);
const int y = int(ty * i);
const double dx = tx * j - x;
const double dy = ty * i - y;
for (int k = 0; k < 3; ++k)
{
for (int jj = 0; jj < 4; ++jj)
{
const int idx = y - 1 + jj;
unsigned char a0 = get_subpixel(bmap, idx, x, k);
unsigned char d0 = get_subpixel(bmap, idx, x - 1, k) - a0;
unsigned char d2 = get_subpixel(bmap, idx, x + 1, k) - a0;
unsigned char d3 = get_subpixel(bmap, idx, x + 2, k) - a0;
unsigned char a1 = -1.0 / 3 * d0 + d2 - 1.0 / 6 * d3;
unsigned char a2 = 1.0 / 2 * d0 + 1.0 / 2 * d2;
unsigned char a3 = -1.0 / 6 * d0 - 1.0 / 2 * d2 + 1.0 / 6 * d3;
C[jj] = a0 + a1 * dx + a2 * dx * dx + a3 * dx * dx * dx;
d0 = C[0] - C[1];
d2 = C[2] - C[1];
d3 = C[3] - C[1];
a0 = C[1];
a1 = -1.0 / 3 * d0 + d2 -1.0 / 6 * d3;
a2 = 1.0 / 2 * d0 + 1.0 / 2 * d2;
a3 = -1.0 / 6 * d0 - 1.0 / 2 * d2 + 1.0 / 6 * d3;
out[i * row_stride + j * channels + k] = a0 + a1 * dy + a2 * dy * dy + a3 * dy * dy * dy;
}
}
}
}
У вас есть смесь unsigned char, int а также double, Каждый из тех 1.0 / 3 преобразует ваши 8-битные данные в число с плавающей запятой двойной точности, а затем присваивает их обратно.
Вместо этого, почему бы просто не использовать float на протяжении?
float C[5] = { 0 };
for (unsigned i = 0; i < dest_height; ++i)
{
for (unsigned j = 0; j < dest_width; ++j)
{
const float x = float(tx * j);
const float y = float(ty * i);
const float dx = tx * j - x;
const float dy = ty * i - y;
for (int k = 0; k < 3; ++k)
{
for (int jj = 0; jj < 4; ++jj)
{
const int idx = y - 1 + jj;
float a0 = get_subpixel(bmap, idx, x, k);
float d0 = get_subpixel(bmap, idx, x - 1, k) - a0;
float d2 = get_subpixel(bmap, idx, x + 1, k) - a0;
float d3 = get_subpixel(bmap, idx, x + 2, k) - a0;
float a1 = -(1.0f / 3.0f) * d0 + d2 - (1.0f / 6.0f) * d3;
float a2 = 0.5f * d0 + 0.5f * d2;
float a3 = -(1.0f / 6.0f) * d0 - 0.5f * d2 + (1.0f / 6.0f) * d3;
C[jj] = a0 + a1 * dx + a2 * dx * dx + a3 * dx * dx * dx;
d0 = C[0] - C[1];
d2 = C[2] - C[1];
d3 = C[3] - C[1];
a0 = C[1];
a1 = -(1.0f / 3.0f) * d0 + d2 -(1.0f / 6.0f) * d3;
a2 = 0.5f * d0 + 0.5f * d2;
a3 = -(1.0f / 6.0f) * d0 - 0.5f * d2 + (1.0f / 6.0f) * d3;
out[i * row_stride + j * channels + k] = saturate( a0 + a1 * dy + a2 * dy * dy + a3 * dy * dy * dy );
}
}
}
}
Затем определите функцию saturate это делает это:
inline unsigned char saturate( float x )
{
return x > 255.0f ? 255
: x < 0.0f ? 0
: unsigned char(x);
}
Это позволит устранить проблемы с переполнением, повысить точность и, вероятно, повысить производительность.
Если вам нужно еще больше улучшить производительность, вы должны изучить арифметику с фиксированной запятой. Но сейчас я думаю, что вышеописанная реализация лучше.
Кроме того, еще одна мысль: вы должны быть в состоянии получить дополнительную эффективность, предварительно вычислив dx * dx, dx * dx * dx, и так далее:
float C[5] = { 0 };
for (unsigned i = 0; i < dest_height; ++i)
{
for (unsigned j = 0; j < dest_width; ++j)
{
const float x = float(tx * j);
const float y = float(ty * i);
const float dx = tx * j - x, dx2 = dx * dx, dx3 = dx2 * dx;
const float dy = ty * i - y, dy2 = dy * dy, dy3 = dy2 * dy;
for (int k = 0; k < 3; ++k)
{
for (int jj = 0; jj < 4; ++jj)
{
const int idx = y - 1 + jj;
float a0 = get_subpixel(bmap, idx, x, k);
float d0 = get_subpixel(bmap, idx, x - 1, k) - a0;
float d2 = get_subpixel(bmap, idx, x + 1, k) - a0;
float d3 = get_subpixel(bmap, idx, x + 2, k) - a0;
float a1 = -(1.0f / 3.0f) * d0 + d2 - (1.0f / 6.0f) * d3;
float a2 = 0.5f * d0 + 0.5f * d2;
float a3 = -(1.0f / 6.0f) * d0 - 0.5f * d2 + (1.0f / 6.0f) * d3;
C[jj] = a0 + a1 * dx + a2 * dx2 + a3 * dx3;
d0 = C[0] - C[1];
d2 = C[2] - C[1];
d3 = C[3] - C[1];
a0 = C[1];
a1 = -(1.0f / 3.0f) * d0 + d2 -(1.0f / 6.0f) * d3;
a2 = 0.5f * d0 + 0.5f * d2;
a3 = -(1.0f / 6.0f) * d0 - 0.5f * d2 + (1.0f / 6.0f) * d3;
out[i * row_stride + j * channels + k] = saturate( a0 + a1 * dy + a2 * dy2 + a3 * dy3 );
}
}
}
}