Commit 1962d2da authored by Téo Mazars's avatar Téo Mazars

operations: attenpt to fix lens-distortion.c

- Never allocate a too large buffer
- Give it GeglOperationAreaFilter's flavors
- Improve the use of the allocated buffer
- Format fixes and some renaming for clarity

Still not perfect, get_required_for_output may not
be 100% bullet-proof with corner values
parent becdbd5d
......@@ -13,15 +13,15 @@
* You should have received a copy of the GNU Lesser General Public
* License along with GEGL; if not, see <http://www.gnu.org/licenses/>.
*
* Lens plug-in - adjust for lens distortion
*
* Copyright (C) 1995 Spencer Kimball and Peter Mattis
* Copyright (C) 2001-2005 David Hodson <hodsond@acm.org>
* Copyright (C) 2011 Robert Sasu <sasu.robert@gmail.com>
* Copyright (C) 2013 Téo Mazars <teo.mazars@ensimag.fr>
*
* Lens plug-in - adjust for lens distortion
* Copyright (C) 2001-2005 David Hodson hodsond@acm.org
* Many thanks for Lars Clausen for the original inspiration,
* useful discussion, optimisation and improvements.
* Framework borrowed from many similar plugins...
*
* Copyright (C) 2011 Robert Sasu <sasu.robert@gmail.com>
*/
#include "config.h"
......@@ -36,6 +36,7 @@ gegl_chant_double (main, _("Main"),
gegl_chant_double (zoom, _("Zoom"),
-100.0, 100.0, 0.0,
_("Main value of distortion"))
gegl_chant_double (edge, _("Edge"),
-100.0, 100.0, 0.0,
_("Edge value of distortion"))
......@@ -59,6 +60,17 @@ gegl_chant_double (y_shift, _("Y shift"),
#include "gegl-chant.h"
#include <math.h>
#include <stdio.h>
#define SQR(x) ((x)*(x))
#define MIN3(x,y,z) (MIN (MIN ((x),(y)),(z)))
#define MAX3(x,y,z) (MAX (MAX ((x),(y)),(z)))
#define MAX_WH 1024
#define CHUNK_SIZE 512
typedef struct
{
......@@ -69,76 +81,186 @@ typedef struct
gdouble rescale;
gdouble brighten;
gdouble norm;
} LensDistortion;
} LensValues;
static void
lens_setup_calc (GeglChantO *o,
GeglRectangle boundary,
LensDistortion *old)
reorder (gdouble *low,
gdouble *high)
{
old->norm = 4.0 / (boundary.width * boundary.width +
boundary.height * boundary.height);
old->centre_x = boundary.width * (100.0 + o->x_shift) / 200.0;
old->centre_y = boundary.height * (100.0 + o->y_shift) / 200.0;
old->mult_sq = o->main / 200.0;
old->mult_qd = o->edge / 200.0;
old->rescale = pow (2.0, -o->zoom / 100.0);
old->brighten = -o->brighten / 10.0;
gdouble temp;
if (*low < *high) return;
temp = *low;
*low = *high;
*high = temp;
}
static GeglRectangle
get_bounding_box (GeglOperation *operation)
static LensValues
lens_setup_calc (GeglChantO *o,
GeglRectangle boundary)
{
GeglRectangle result = {0,0,0,0};
GeglRectangle *in_rect;
LensValues lens;
in_rect = gegl_operation_source_get_bounding_box (operation, "input");
if (!in_rect)
return result;
lens.norm = 4.0 / (SQR (boundary.width) + SQR (boundary.height));
lens.centre_x = boundary.width * (100.0 + o->x_shift) / 200.0;
lens.centre_y = boundary.height * (100.0 + o->y_shift) / 200.0;
lens.mult_sq = o->main / 200.0;
lens.mult_qd = o->edge / 200.0;
lens.rescale = pow (2.0, - o->zoom / 100.0);
lens.brighten = - o->brighten / 10.0;
return *in_rect;
return lens;
}
static GeglRectangle
get_required_for_output (GeglOperation *operation,
const gchar *input_pad,
const GeglRectangle *roi)
static void
lens_get_source_coord (gdouble i,
gdouble j,
gdouble *x,
gdouble *y,
gdouble *mag,
LensValues *lens)
{
return get_bounding_box (operation);
gdouble radius_sq, off_x, off_y, radius_mult;
off_x = i - lens->centre_x;
off_y = j - lens->centre_y;
radius_sq = SQR (off_x) + SQR (off_y);
radius_sq *= lens->norm;
radius_mult = radius_sq * lens->mult_sq + SQR (radius_sq) * lens->mult_qd;
*mag = radius_mult;
radius_mult = lens->rescale * (1.0 + radius_mult);
*x = lens->centre_x + radius_mult * off_x;
*y = lens->centre_y + radius_mult * off_y;
}
static void
prepare (GeglOperation *operation)
/* FIXME: not 100% bullet proof */
static GeglRectangle
get_required (GeglRectangle *boundary,
const GeglRectangle *roi,
GeglOperation *operation)
{
gegl_operation_set_format (operation, "input", babl_format ("RGBA float"));
gegl_operation_set_format (operation, "output", babl_format ("RGBA float"));
GeglChantO *o;
GeglRectangle area;
LensValues lens;
gdouble x1, y1, x2, y2, x3, y3, x4, y4, mag;
gint x, y, width, height;
o = GEGL_CHANT_PROPERTIES (operation);
lens = lens_setup_calc (o, *boundary);
x = roi->x;
y = roi->y;
width = roi->width;
height = roi->height;
lens_get_source_coord (x, y, &x1, &y1, &mag, &lens);
lens_get_source_coord (x + width, y, &x2, &y2, &mag, &lens);
lens_get_source_coord (x, y + height, &x3, &y3, &mag, &lens);
lens_get_source_coord (x + width, y + height, &x4, &y4, &mag, &lens);
/* This is ugly, and happens
* with a crazy set of parameters */
reorder (&x1, &x2);
reorder (&x3, &x4);
reorder (&y1, &y3);
reorder (&y2, &y4);
if (lens.centre_y > y && lens.centre_y < y + height)
{
gdouble x5, y5, x6, y6;
lens_get_source_coord (x, lens.centre_y, &x5, &y5, &mag, &lens);
lens_get_source_coord (x + width, lens.centre_y, &x6, &y6, &mag, &lens);
reorder (&x5, &x6);
area.x = floor (MIN3 (x1, x3, x5)) - 1;
area.width = ceil (MAX3 (x2, x4, x6)) + 3 - area.x;
}
else
{
area.x = floor (MIN (x1, x3)) - 1;
area.width = ceil (MAX (x2, x4)) + 3 - area.x;
}
if (lens.centre_x > x && lens.centre_x < x + width)
{
gdouble x5, y5, x6, y6;
lens_get_source_coord (lens.centre_x, y, &x5, &y5, &mag, &lens);
lens_get_source_coord (lens.centre_x, y + height, &x6, &y6, &mag, &lens);
reorder (&y5, &y6);
area.y = floor (MIN3 (y1, y2, y5)) - 1;
area.height = ceil (MAX3 (y3, y4, y6)) + 3 - area.y;
}
else
{
area.y = floor (MIN (y1, y2)) - 1;
area.height = ceil (MAX (y3, y4)) + 3 - area.y;
}
return area;
}
static void
lens_get_source_coord (gdouble i,
gdouble j,
gdouble *x,
gdouble *y,
gdouble *mag,
LensDistortion *o)
clamp_area (GeglRectangle *area,
gdouble center_x,
gdouble center_y)
{
gdouble radius_sq, off_x, off_y, radius_mult;
if (center_x <= area->x || area->width < 1)
{
area->width = CLAMP (area->width, 1, MAX_WH);
}
else
{
area->x += area->width - CLAMP (area->width, 1, MAX_WH);
area->width = CLAMP (area->width, 1, MAX_WH);
}
off_x = i - o->centre_x;
off_y = j - o->centre_y;
radius_sq = (off_x * off_x) + (off_y * off_y);
if (center_y <= area->y || area->height < 1)
{
area->height = CLAMP (area->height, 1, MAX_WH);
}
else
{
area->y += area->height - CLAMP (area->height, 1, MAX_WH);
area->height = CLAMP (area->height, 1, MAX_WH);
}
}
radius_sq *= o->norm;
static GeglRectangle
get_required_for_output (GeglOperation *operation,
const gchar *input_pad,
const GeglRectangle *roi)
{
GeglRectangle *boundary;
radius_mult = radius_sq * o->mult_sq + radius_sq * radius_sq * o->mult_qd;
boundary = gegl_operation_source_get_bounding_box (operation, "input");
*mag = radius_mult;
radius_mult = o->rescale * (1.0 + radius_mult);
if (strcmp (input_pad, "input") || !boundary)
return *GEGL_RECTANGLE (0, 0, 0, 0);
*x = o->centre_x + radius_mult * off_x;
*y = o->centre_y + radius_mult * off_y;
return get_required (boundary, roi, operation);
}
static void
prepare (GeglOperation *operation)
{
gegl_operation_set_format (operation, "input", babl_format ("RGBA float"));
gegl_operation_set_format (operation, "output", babl_format ("RGBA float"));
}
/*
......@@ -176,26 +298,24 @@ lens_cubic_interpolate (gfloat *src,
vp1 = ((-1.5 * dy + 2.0) * dy + 0.5) * dy;
vp2 = (0.5 * dy - 0.5) * dy * dy;
/* Note: if dst_depth < src_depth, we calculate unneeded pixels here */
/* later - select or create index array */
for (c = 0; c < 4 * 4; ++c)
{
verts[c] = vm1 * src[c] + v * src[c+row_stride] +
vp1 * src[c+row_stride*2] + vp2 * src[c+row_stride*3];
verts[c] = vm1 * src[c] + v * src[c + row_stride] +
vp1 * src[c + row_stride * 2] + vp2 * src[ c + row_stride * 3];
}
for (c = 0; c < 4; ++c)
{
gfloat result;
result = um1 * verts[c] + u * verts[c+4] +
up1 * verts[c+4*2] + up2 * verts[c+4*3];
result = um1 * verts[c] + u * verts[c + 4] +
up1 * verts[c + 4 * 2] + up2 * verts[c + 4 * 3];
result *= brighten;
if (c != 3)
result *= brighten;
dst[c] = CLAMP (result, 0.0, 1.0);
}
}
static void
......@@ -204,7 +324,7 @@ lens_distort_func (gfloat *src_buf,
const GeglRectangle *extended,
const GeglRectangle *result,
const GeglRectangle *boundary,
LensDistortion old,
LensValues *lens,
gint xx,
gint yy,
GeglBuffer *input)
......@@ -217,9 +337,10 @@ lens_distort_func (gfloat *src_buf,
temp[0] = temp[1] = temp[2] = temp[3] = 0.0;
lens_get_source_coord ((gdouble)xx, (gdouble)yy, &sx, &sy, &mag, &old);
lens_get_source_coord ((gdouble) xx, (gdouble) yy, &sx, &sy, &mag, lens);
brighten = 1.0 + mag * old.brighten;
/* pseudo gamma transformation, since the input is scRGB */
brighten = pow (MAX (1.0 + mag * lens->brighten, 0.0), 2.4);
x_int = floor (sx);
dx = sx - x_int;
......@@ -233,31 +354,31 @@ lens_distort_func (gfloat *src_buf,
{
gint b;
if (x >= extended->x && x<(extended->x + extended->width)
if (x >= extended->x && x< (extended->x + extended->width)
&& y >= extended->y && y < (extended->y + extended->height))
{
gint src_off;
src_off = (y - extended->y) * extended->width * 4 +
(x - extended->x) * 4;
for (b=0; b<4; b++)
temp[b] = src_buf[src_off++];
for (b = 0; b < 4; b++)
temp[b] = src_buf[src_off++];
}
else if (x >= boundary->x && x < boundary->x + boundary->width &&
y >= boundary->y && y < boundary->y + boundary->height)
else
{
gegl_buffer_sample (input, x, y, NULL, temp,
babl_format ("RGBA float"),
GEGL_SAMPLER_CUBIC,
GEGL_ABYSS_NONE);
GEGL_SAMPLER_LINEAR,
GEGL_ABYSS_CLAMP);
}
else
if (x < boundary->x || x >= (boundary->x + boundary->width) ||
y < boundary->y || y >= (boundary->y + boundary->height))
{
for (b=0; b<4; b++)
temp[b] = 0.0;
temp[3] = 0.0;
}
for (b=0; b<4; b++)
for (b = 0; b < 4; b++)
pixel_buffer[offset++] = temp[b];
}
}
......@@ -265,7 +386,8 @@ lens_distort_func (gfloat *src_buf,
lens_cubic_interpolate (pixel_buffer, temp, dx, dy, brighten);
offset = (yy - result->y) * result->width * 4 + (xx - result->x) * 4;
for (x=0; x<4; x++)
for (x = 0; x < 4; x++)
dst_buf[offset++] = temp[x];
}
......@@ -276,32 +398,51 @@ process (GeglOperation *operation,
const GeglRectangle *result,
gint level)
{
GeglChantO *o = GEGL_CHANT_PROPERTIES (operation);
LensDistortion old_lens;
GeglRectangle boundary = *gegl_operation_source_get_bounding_box
(operation, "input");
GeglChantO *o = GEGL_CHANT_PROPERTIES (operation);
LensValues lens;
GeglRectangle boundary;
gint i, j;
gfloat *src_buf, *dst_buf;
boundary = *gegl_operation_source_get_bounding_box (operation, "input");
lens = lens_setup_calc (o, boundary);
gint x, y;
gfloat *src_buf, *dst_buf;
src_buf = g_new0 (gfloat, SQR (MAX_WH) * 4);
dst_buf = g_new0 (gfloat, SQR (CHUNK_SIZE) * 4);
src_buf = g_new0 (gfloat, result->width * result->height * 4);
dst_buf = g_new0 (gfloat, result->width * result->height * 4);
for (j = 0; (j-1) * CHUNK_SIZE < result->height; j++)
for (i = 0; (i-1) * CHUNK_SIZE < result->width; i++)
{
GeglRectangle chunked_result;
GeglRectangle area;
gint x, y;
lens_setup_calc (o, boundary, &old_lens);
chunked_result = *GEGL_RECTANGLE (result->x + i * CHUNK_SIZE,
result->y + j * CHUNK_SIZE,
CHUNK_SIZE, CHUNK_SIZE);
gegl_buffer_get (input, result, 1.0, babl_format ("RGBA float"), src_buf,
GEGL_AUTO_ROWSTRIDE, GEGL_ABYSS_NONE);
gegl_rectangle_intersect (&chunked_result, &chunked_result, result);
for (y = result->y; y < result->y + result->height; y++)
for (x = result->x; x < result->x + result->width; x++)
{
lens_distort_func (src_buf, dst_buf, result, result, &boundary,
old_lens, x, y, input);
}
if (chunked_result.width < 1 || chunked_result.height < 1)
continue;
area = get_required (&boundary, &chunked_result, operation);
gegl_buffer_set (output, result, 0, babl_format ("RGBA float"),
dst_buf, GEGL_AUTO_ROWSTRIDE);
clamp_area (&area, lens.centre_x, lens.centre_y);
gegl_buffer_get (input, &area, 1.0, babl_format ("RGBA float"), src_buf,
GEGL_AUTO_ROWSTRIDE, GEGL_ABYSS_CLAMP);
for (y = chunked_result.y; y < chunked_result.y + chunked_result.height; y++)
for (x = chunked_result.x; x < chunked_result.x + chunked_result.width; x++)
{
lens_distort_func (src_buf, dst_buf, &area, &chunked_result, &boundary,
&lens, x, y, input);
}
gegl_buffer_set (output, &chunked_result, 0, babl_format ("RGBA float"),
dst_buf, GEGL_AUTO_ROWSTRIDE);
}
g_free (dst_buf);
g_free (src_buf);
......@@ -337,7 +478,6 @@ gegl_chant_class_init (GeglChantClass *klass)
filter_class = GEGL_OPERATION_FILTER_CLASS (klass);
operation_class->prepare = prepare;
operation_class->get_bounding_box = get_bounding_box;
operation_class->get_required_for_output = get_required_for_output;
filter_class->process = process;
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment