Bilinear interpolation.

- Compute src delta for each rotated steps.
master
Fabien Freling 2014-06-28 18:23:43 +02:00
parent 258ddce814
commit 0a0fbabeac
1 changed files with 89 additions and 57 deletions

View File

@ -424,6 +424,18 @@ void draw_outline(Image const& input, unsigned int degrees, string const& name)
//
//
// Math approximation
//
void round_if_very_small(double& d)
{
if (abs(d) < 1.0e-10)
d = 0.0;
}
//
//
// Image rotation
@ -440,26 +452,41 @@ DPoint get_mapped_point(Image const& src, Point const& p, double const rotation)
void rotate_pixel(Image const& src, Image& rotated,
DPoint const& src_rotated_point, Point const& rot_point,
unsigned int const src_limit, unsigned int const rot_limit,
Point const& previous)
unsigned int const src_limit, unsigned int const rot_limit)
{
// TODO: Interpolation would be done here
// For now we take the nearest point
// Out-of-bounds check
unsigned int src_index = ((int) src_rotated_point.y * src.width + (int) src_rotated_point.x) * 3;
unsigned int rot_index = (rot_point.y * rotated.width + rot_point.x) * 3;
if (src_index >= src_limit
|| rot_index >= rot_limit)
return;
memcpy(&rotated.buffer[rot_index], &src.buffer[src_index], 3 * sizeof (uint8_t));
static double interpolation_bottom[3];
static double interpolation_top[3];
static uint8_t interpolated[3];
// Fill missing points, created by interpolation
if (previous.x != rot_point.x && previous.y != rot_point.y)
// Bilinear interpolation
unsigned int src_index_1 = (floor(src_rotated_point.y) * src.width + floor(src_rotated_point.x)) * 3;
unsigned int src_index_2 = (floor(src_rotated_point.y) * src.width + ceil(src_rotated_point.x)) * 3;
unsigned int src_index_3 = (ceil(src_rotated_point.y) * src.width + floor(src_rotated_point.x)) * 3;
unsigned int src_index_4 = (ceil(src_rotated_point.y) * src.width + ceil(src_rotated_point.x)) * 3;
if (src_index_1 >= src_limit || src_index_1 >= src_limit)
return;
double x_delta = src_rotated_point.x - floor(src_rotated_point.x);
round_if_very_small(x_delta);
double y_delta = src_rotated_point.y - floor(src_rotated_point.y);
round_if_very_small(y_delta);
// Interpolation
// SIMD?
for (unsigned int i = 0; i < 3; ++i)
{
src_index -= 3;
unsigned int previous_index = (previous.y * rotated.width + rot_point.x) * 3;
memcpy(&rotated.buffer[previous_index], &src.buffer[src_index], 3 * sizeof (uint8_t));
interpolation_bottom[i] = src.buffer[src_index_1 + i] * (1 - x_delta) + src.buffer[src_index_2 + i] * x_delta;
interpolation_top[i] = src.buffer[src_index_3 + i] * (1 - x_delta) + src.buffer[src_index_4 + i] * x_delta;
interpolated[i] = interpolation_bottom[i] * (1 - y_delta) + interpolation_top[i] * y_delta;
}
memcpy(&rotated.buffer[rot_index], interpolated, 3 * sizeof (uint8_t));
}
Image rotate(Image const& src, double angle)
@ -471,33 +498,42 @@ Image rotate(Image const& src, double angle)
Image rotated(w, h);
// corner points in rotated image
DPoint tl_grid = get_mapped_point(src, Point(0, 0), rotation);
Point tl = convert_img_coord(rotated, tl_grid);
DPoint tr_grid = get_mapped_point(src, Point(src.width - 1, 0), rotation);
Point tr = convert_img_coord(rotated, tr_grid);
DPoint bl_grid = get_mapped_point(src, Point(0, src.height - 1), rotation);
Point bl = convert_img_coord(rotated, bl_grid);
// TODO: add one ligne for smooth border
DPoint const tl_grid = get_mapped_point(src, Point(0, 0), rotation);
Point const tl = convert_img_coord(rotated, tl_grid);
DPoint const tr_grid = get_mapped_point(src, Point(src.width - 1, 0), rotation);
Point const tr = convert_img_coord(rotated, tr_grid);
DPoint const bl_grid = get_mapped_point(src, Point(0, src.height - 1), rotation);
Point const bl = convert_img_coord(rotated, bl_grid);
// corner points in source image
DPoint src_tl = get_mapped_point(rotated, tl, -rotation);
src_tl = convert_img_coord_precision(src, src_tl);
DPoint src_tr = get_mapped_point(rotated, tr, -rotation);
src_tr = convert_img_coord_precision(src, src_tr);
DPoint src_bl = get_mapped_point(rotated, bl, -rotation);
src_bl = convert_img_coord_precision(src, src_bl);
// steps for first column in source image
DPoint src_origin = get_mapped_point(rotated, Point(0, 0), -rotation);
DPoint src_delta_x = get_mapped_point(rotated, Point(1, 0), -rotation);
DPoint src_delta_y = get_mapped_point(rotated, Point(0, 1), -rotation);
src_delta_x.x = src_delta_x.x - src_origin.x;
src_delta_x.y = src_delta_x.y - src_origin.y;
round_if_very_small(src_delta_x.x);
round_if_very_small(src_delta_x.y);
src_delta_y.x = src_delta_y.x - src_origin.x;
src_delta_y.y = src_delta_y.y - src_origin.y;
round_if_very_small(src_delta_y.x);
round_if_very_small(src_delta_y.y);
// // steps for first column in source image (y)
int origin_nb_steps = max(abs(bl.x - tl.x), abs(bl.y - tl.y));
DPoint origin_step((src_bl.x - src_tl.x) / origin_nb_steps, (src_bl.y - src_tl.y) / origin_nb_steps);
// steps for line in source image
// // steps for line in source image (x)
int line_nb_steps = max(abs(tr.x - tl.x), abs(tr.y - tl.y));
DPoint line_step((src_tr.x - src_tl.x) / line_nb_steps, (src_tr.y - src_tl.y) / line_nb_steps);
// steps for first column in rotated image
// steps for first column in rotated image (y)
DPoint rotated_step((bl.x - tl.x) / (float) origin_nb_steps, (bl.y - tl.y) / (float) origin_nb_steps);
// steps for line in rotated image
// steps for line in rotated image (x)
DPoint bresenham((tr.x - tl.x) / (float) line_nb_steps, (tr.y - tl.y) / (float) line_nb_steps);
unsigned int const src_limit = src.width * src.height * 3;
@ -506,19 +542,34 @@ Image rotate(Image const& src, double angle)
for (int y_i = 0; y_i <= (int) origin_nb_steps; ++y_i)
{
// first column origin
DPoint const src_origin(src_tl.x + y_i * origin_step.x, src_tl.y + y_i * origin_step.y);
Point const rot_origin(tl.x + y_i * rotated_step.x, tl.y + y_i * rotated_step.y);
Point previous = rot_origin;
for (int x_i = 0; x_i <= (int) line_nb_steps; ++x_i)
{
DPoint const src_rotated_point(src_origin.x + x_i * line_step.x, src_origin.y + x_i * line_step.y);
Point const rot_point(rot_origin.x + x_i * bresenham.x, rot_origin.y + x_i * bresenham.y);
// cout << "src rot point: " << src_rotated_point << endl;
// cout << "rot point: " << rot_point << endl;
Point rot_point(rot_origin.x + x_i * bresenham.x, rot_origin.y + x_i * bresenham.y);
Point const delta(rot_point.x - tl.x, rot_point.y - tl.y);
DPoint src_rotated_point(src_tl.x + delta.x * src_delta_x.x + delta.y * src_delta_y.x,
src_tl.y + delta.x * src_delta_x.y + delta.y * src_delta_y.y);
rotate_pixel(src, rotated, src_rotated_point, rot_point, src_limit, rot_limit);
if (previous.x != rot_point.x && previous.y != rot_point.y)
{
int y_slope = rot_point.y > previous.y ? 1 : -1;
int tmp_y = rot_point.y;
rot_point.y = previous.y;
src_rotated_point.x -= y_slope * src_delta_y.x;
src_rotated_point.y -= y_slope * src_delta_y.y;
rotate_pixel(src, rotated, src_rotated_point, rot_point, src_limit, rot_limit);
rot_point.y = tmp_y;
}
rotate_pixel(src, rotated, src_rotated_point, rot_point, src_limit, rot_limit, previous);
previous = rot_point;
}
}
@ -632,25 +683,6 @@ bool check_trigo()
return true;
}
void check_lines()
{
Image const square(500, 500);
draw_outline(square, 5, "square");
draw_outline(square, 12, "square");
draw_outline(square, 22, "square");
draw_outline(square, 33, "square");
draw_outline(square, 45, "square");
draw_outline(square, 60, "square");
draw_outline(square, 75, "square");
draw_outline(square, 90, "square");
Image const rect1(640, 480);
draw_outline(rect1, 22, "rect1");
draw_outline(rect1, 33, "rect1");
draw_outline(rect1, 45, "rect1");
draw_outline(rect1, 90, "rect1");
}
bool check_90(string const& path)
{
Image const src(path);
@ -702,19 +734,18 @@ int main(int argc, char* argv[])
if (!check_trigo())
return 1;
//check_lines();
if (!check_90(argv[1]))
{
cerr << __LINE__ << " | 90 degrees check failed" << endl;
return 1;
// return 1;
}
}
Image img(argv[1]);
for (double rotation = 2; rotation < 360; rotation += 400)
for (double rotation = 0; rotation < 360; rotation += 15)
{
// double rotation = 45;
auto const before = chrono::high_resolution_clock::now();
Image rotated = rotate(img, rotation);
@ -724,7 +755,8 @@ int main(int argc, char* argv[])
cout << "rotate(" << rotation << "): " << duration_ms.count() << " ms" << endl;
stringstream filename;
filename << "/tmp/rotated_";
// filename << "/tmp/";
filename << "rotated_";
if (rotation < 100)
filename << "0";
if (rotation < 10)