revise affine, yet again

how can this be so hard

the tougher self-tests in nip2/test/extras now all pass
This commit is contained in:
John Cupitt 2014-08-01 16:34:09 +01:00
parent 9ddca0e99e
commit 274e6c1b2b
4 changed files with 124 additions and 111 deletions

View File

@ -3,6 +3,7 @@
- limit n_thr on tiny images
- don't exit() on memleak detected, just warn
- add "autocrop" option to openslide load
- argh fix affine, again
4/7/14 started 7.40.4
- fix vips_rawsave_fd(), thanks aferrero2707

View File

@ -709,10 +709,6 @@ vips_image_rewind( VipsObject *object )
/* Delayed save.
*/
/* If we write to (eg.) TIFF, actually do the write
* to a "p" and on "written" do im_vips2tiff() or whatever.
*/
/* From "written" callback: save to image->filename using VipsForeign.
*/
static void
@ -1102,10 +1098,10 @@ vips_image_class_init( VipsImageClass *class )
* can't :-(
*
* For example, a "p" image might be made with vips_image_new() and
* constructed, then passed to im_copy() of whatever to be written to.
* constructed, then passed to vips_copy() of whatever to be written to.
* That operation will then need to set width/height etc.
*
* We can't set_once either, since im_copy_set() etc. need to update
* We can't set_once either, since vips_copy() etc. need to update
* xoffset and friends on the way through.
*/
@ -2321,10 +2317,10 @@ vips_image_write_to_memory( VipsImage *in, void **buf_out, size_t *len_out )
len = VIPS_IMAGE_SIZEOF_IMAGE( in );
if( !(buf = g_try_malloc( len )) ) {
vips_error( "vips_image_write_to_buffer",
vips_error( "vips_image_write_to_memory",
_( "out of memory --- size == %dMB" ),
(int) (len / (1024.0 * 1024.0)) );
vips_warn( "vips_image_write_to_buffer",
vips_warn( "vips_image_write_to_memory",
_( "out of memory --- size == %dMB" ),
(int) (len / (1024.0*1024.0)) );
return( -1 );
@ -3025,7 +3021,7 @@ vips_image_pio_output( VipsImage *image )
case VIPS_IMAGE_PARTIAL:
if( image->generate_fn ) {
vips_error( "im_poutcheck",
vips_error( "vips_image_pio_output",
"%s", _( "image already written" ) );
return( -1 );
}

View File

@ -81,6 +81,9 @@
* - added input space translation
* 22/1/14
* - auto RAD decode
* 1/8/14
* - revise transform ... again
* - see new stress test in nip2/test/extras
*/
/*
@ -111,6 +114,7 @@
*/
/*
#define DEBUG_VERBOSE
#define DEBUG
*/
@ -177,11 +181,12 @@ G_DEFINE_TYPE( VipsAffine, vips_affine, VIPS_TYPE_RESAMPLE );
/* We have five (!!) coordinate systems. Working forward through them, these
* are:
*
* 1. The original input image
* 1. The original input image. iarea is defined on this image.
*
* 2. This is embedded in a larger image to provide borders for the
* interpolator. iarea->left/top give the offset. These are the coordinates we
* pass to VIPS_REGION_ADDR()/vips_region_prepare() for the input image.
* interpolator. window_offset and window_size control the embedding.
* These are the coordinates we pass to VIPS_REGION_ADDR()/
* vips_region_prepare() and the interpolator.
*
* The borders are sized by the interpolator's window_size property and offset
* by the interpolator's window_offset property. For example,
@ -194,7 +199,7 @@ G_DEFINE_TYPE( VipsAffine, vips_affine, VIPS_TYPE_RESAMPLE );
* 3. We need point (0, 0) in (1) to be at (0, 0) for the transformation. So
* shift everything up and left to make the displaced input image. This is the
* space that the transformation maps from, and can have negative pixels
* (up and left of the image, for interpolation).
* (up and left of the image, for interpolation). iarea works here too.
*
* 4. Output transform space. This is the where the transform maps to. Pixels
* can be negative, since a rotated image can go up and left of the origin.
@ -233,22 +238,23 @@ vips_affine_gen( VipsRegion *or, void *seq, void *a, void *b, gboolean *stop )
VipsRect image, want, need, clipped;
#ifdef DEBUG
#ifdef DEBUG_VERBOSE
printf( "vips_affine_gen: "
"generating left=%d, top=%d, width=%d, height=%d\n",
r->left,
r->top,
r->width,
r->height );
#endif /*DEBUG*/
#endif /*DEBUG_VERBOSE*/
/* We are generating this chunk of the transformed image.
/* We are generating this chunk of the transformed image. This takes
* us to space 4.
*/
want = *r;
want.left += oarea->left;
want.top += oarea->top;
/* Find the area of the input image we need.
/* Find the area of the input image we need. This takes us to space 3.
*/
vips__transform_invert_rect( &affine->trn, &want, &need );
@ -261,17 +267,18 @@ vips_affine_gen( VipsRegion *or, void *seq, void *a, void *b, gboolean *stop )
*/
vips_rect_marginadjust( &need, 1 );
/* Now go to space (2) above.
*/
need.left += iarea->left;
need.top += iarea->top;
/* Add a border for interpolation.
/* We need to fetch a larger area for the interpolator.
*/
need.left -= window_offset;
need.top -= window_offset;
need.width += window_size - 1;
need.height += window_size - 1;
need.left -= window_offset;
need.top -= window_offset;
/* Now go to space 2, the expanded input image. This is the one we
* read pixels from.
*/
need.left += window_offset;
need.top += window_offset;
/* Clip against the size of (2).
*/
@ -281,38 +288,34 @@ vips_affine_gen( VipsRegion *or, void *seq, void *a, void *b, gboolean *stop )
image.height = in->Ysize;
vips_rect_intersectrect( &need, &image, &clipped );
/* Outside input image? All black.
*/
if( vips_rect_isempty( &clipped ) ) {
vips_region_black( or );
return( 0 );
}
/* We do need some pixels from the input image to make our output -
* ask for them.
*/
if( vips_region_prepare( ir, &clipped ) )
return( -1 );
#ifdef DEBUG
printf( "affine: preparing left=%d, top=%d, width=%d, height=%d\n",
#ifdef DEBUG_VERBOSE
printf( "vips_affine_gen: "
"preparing left=%d, top=%d, width=%d, height=%d\n",
clipped.left,
clipped.top,
clipped.width,
clipped.height );
#endif /*DEBUG*/
#endif /*DEBUG_VERBOSE*/
if( vips_rect_isempty( &clipped ) ) {
vips_region_black( or );
return( 0 );
}
if( vips_region_prepare( ir, &clipped ) )
return( -1 );
VIPS_GATE_START( "vips_affine_gen: work" );
/* Resample! x/y loop over pixels in the output image (5).
*/
for( y = to; y < bo; y++ ) {
/* Input clipping rectangle.
/* Input clipping rectangle. We offset this so we can clip in
* space 2.
*/
const int ile = iarea->left;
const int ito = iarea->top;
const int iri = iarea->left + iarea->width;
const int ibo = iarea->top + iarea->height;
const int ile = iarea->left + window_offset;
const int ito = iarea->top + window_offset;
const int iri = ile + iarea->width;
const int ibo = ito + iarea->height;
/* Derivative of matrix.
*/
@ -335,16 +338,16 @@ vips_affine_gen( VipsRegion *or, void *seq, void *a, void *b, gboolean *stop )
ix = affine->trn.ia * ox + affine->trn.ib * oy;
iy = affine->trn.ic * ox + affine->trn.id * oy;
/* Now move to (2).
*/
ix += iarea->left;
iy += iarea->top;
/* And the input offset.
/* And the input offset in (3).
*/
ix -= affine->trn.idx;
iy -= affine->trn.idy;
/* Finally to 2.
*/
ix += window_offset;
iy += window_offset;
q = VIPS_REGION_ADDR( or, le, y );
for( x = le; x < ri; x++ ) {
@ -353,19 +356,31 @@ vips_affine_gen( VipsRegion *or, void *seq, void *a, void *b, gboolean *stop )
fx = FAST_PSEUDO_FLOOR( ix );
fy = FAST_PSEUDO_FLOOR( iy );
/* Clipping!
/* Clip against iarea.
*/
if( fx < ile ||
fx > iri ||
fy < ito ||
fy > ibo ) {
for( z = 0; z < ps; z++ )
q[z] = 0;
}
else {
if( fx >= ile &&
fx < iri &&
fy >= ito &&
fy < ibo ) {
/* Verify that we can read the whole stencil.
* With DEBUG on this will range-check.
*/
g_assert( VIPS_REGION_ADDR( ir,
(int) ix - window_offset,
(int) iy - window_offset ) );
g_assert( VIPS_REGION_ADDR( ir,
(int) ix - window_offset +
window_size - 1,
(int) iy - window_offset +
window_size - 1 ) );
interpolate( affine->interpolate,
q, ir, ix, iy );
}
else {
for( z = 0; z < ps; z++ )
q[z] = 0;
}
ix += ddx;
iy += ddy;
@ -417,8 +432,8 @@ vips_affine_build( VipsObject *object )
window_offset =
vips_interpolate_get_window_offset( affine->interpolate );
affine->trn.iarea.left = window_offset;
affine->trn.iarea.top = window_offset;
affine->trn.iarea.left = 0;
affine->trn.iarea.top = 0;
affine->trn.iarea.width = in->Xsize;
affine->trn.iarea.height = in->Ysize;
affine->trn.a = ((double *) affine->matrix->data)[0];
@ -430,6 +445,11 @@ vips_affine_build( VipsObject *object )
affine->trn.odx = 0;
affine->trn.ody = 0;
if( vips__transform_calc_inverse( &affine->trn ) )
return( -1 );
/* Set the default value for oarea.
*/
vips__transform_set_area( &affine->trn );
if( vips_object_argument_isset( object, "oarea" ) ) {
@ -449,18 +469,16 @@ vips_affine_build( VipsObject *object )
if( vips_object_argument_isset( object, "idy" ) )
affine->trn.idy = affine->idy;
if( vips__transform_calc_inverse( &affine->trn ) )
return( -1 );
#ifdef DEBUG
printf( "vips_affine_build: copy on identity transform disabled\n" );
#else /*!DEBUG*/
if( vips__transform_isidentity( &affine->trn ) &&
affine->trn.oarea.left == 0 &&
affine->trn.oarea.top == 0 &&
affine->trn.oarea.width == in->Xsize &&
affine->trn.oarea.height == in->Ysize )
return( vips_image_write( in, resample->out ) );
/*
printf( "vips_affine_build: identity shortcut disabled\n" );
*/
#endif /*!DEBUG*/
/* Check for coordinate overflow ... we want to be able to hold the
* output space inside INT_MAX / TRANSFORM_SCALE.
@ -483,13 +501,14 @@ vips_affine_build( VipsObject *object )
*/
if( vips_embed( in, &t[2],
window_offset, window_offset,
in->Xsize + window_size, in->Ysize + window_size,
in->Xsize + window_size - 1, in->Ysize + window_size - 1,
"extend", VIPS_EXTEND_COPY,
NULL ) )
return( -1 );
in = t[2];
/* Normally SMALLTILE ... except if this is a size up/down affine.
/* Normally SMALLTILE ... except if this is strictly a size
* up/down affine.
*/
if( affine->trn.b == 0.0 &&
affine->trn.c == 0.0 )
@ -504,21 +523,14 @@ vips_affine_build( VipsObject *object )
resample->out->Ysize = affine->trn.oarea.height;
#ifdef DEBUG
printf( "vips_affine_build: transform: "
"a = %g, b = %g, c = %g, d = %g\n",
affine->trn.a,
affine->trn.b,
affine->trn.c,
affine->trn.d );
#endif /*DEBUG*/
#ifdef DEBUG
printf( "vips_affine_build: output area: "
"left = %d, top = %d, width = %d, height = %d\n",
affine->trn.oarea.left,
affine->trn.oarea.top,
affine->trn.oarea.width,
affine->trn.oarea.height );
printf( "vips_affine_build: transform: " );
vips__transform_print( &affine->trn );
printf( " window_offset = %d, window_size = %d\n",
window_offset, window_size );
printf( " input image width = %d, height = %d\n",
in->Xsize, in->Ysize );
printf( " output image width = %d, height = %d\n",
resample->out->Xsize, resample->out->Ysize );
#endif /*DEBUG*/
/* Generate!

View File

@ -156,8 +156,8 @@ vips__transform_print( const VipsTransformation *trn )
*/
void
vips__transform_forward_point( const VipsTransformation *trn,
double x, double y, /* In input space */
double *ox, double *oy ) /* In output space */
double x, double y, /* In input space */
double *ox, double *oy )/* In output space */
{
x += trn->idx;
y += trn->idy;
@ -170,11 +170,11 @@ vips__transform_forward_point( const VipsTransformation *trn,
*/
void
vips__transform_invert_point( const VipsTransformation *trn,
double x, double y, /* In output space */
double *ox, double *oy ) /* In input space */
double x, double y, /* In output space */
double *ox, double *oy )/* In input space */
{
x -= trn->odx;
y -= trn->ody;
x -= trn->odx;
y -= trn->ody;
*ox = trn->ia * x + trn->ib * y - trn->idx;
*oy = trn->ic * x + trn->id * y - trn->idy;
@ -187,8 +187,8 @@ typedef void (*transform_fn)( const VipsTransformation *,
*/
static void
transform_rect( const VipsTransformation *trn, transform_fn transform,
const Rect *in, /* In input space */
Rect *out ) /* In output space */
const VipsRect *in, /* In input space */
VipsRect *out ) /* In output space */
{
double x1, y1; /* Map corners */
double x2, y2;
@ -196,25 +196,29 @@ transform_rect( const VipsTransformation *trn, transform_fn transform,
double x4, y4;
double left, right, top, bottom;
/* Map input Rect.
/* Map input VipsRect.
*/
transform( trn, in->left, in->top, &x1, &y1 );
transform( trn, in->left, IM_RECT_BOTTOM( in ), &x3, &y3 );
transform( trn, IM_RECT_RIGHT( in ), in->top, &x2, &y2 );
transform( trn, IM_RECT_RIGHT( in ), IM_RECT_BOTTOM( in ), &x4, &y4 );
transform( trn, in->left, in->top,
&x1, &y1 );
transform( trn, in->left, VIPS_RECT_BOTTOM( in ),
&x3, &y3 );
transform( trn, VIPS_RECT_RIGHT( in ), in->top,
&x2, &y2 );
transform( trn, VIPS_RECT_RIGHT( in ), VIPS_RECT_BOTTOM( in ),
&x4, &y4 );
/* Find bounding box for these four corners. Round-to-nearest to try
* to stop rounding errors growing images.
*/
left = IM_MIN( x1, IM_MIN( x2, IM_MIN( x3, x4 ) ) );
right = IM_MAX( x1, IM_MAX( x2, IM_MAX( x3, x4 ) ) );
top = IM_MIN( y1, IM_MIN( y2, IM_MIN( y3, y4 ) ) );
bottom = IM_MAX( y1, IM_MAX( y2, IM_MAX( y3, y4 ) ) );
left = VIPS_MIN( x1, VIPS_MIN( x2, VIPS_MIN( x3, x4 ) ) );
right = VIPS_MAX( x1, VIPS_MAX( x2, VIPS_MAX( x3, x4 ) ) );
top = VIPS_MIN( y1, VIPS_MIN( y2, VIPS_MIN( y3, y4 ) ) );
bottom = VIPS_MAX( y1, VIPS_MAX( y2, VIPS_MAX( y3, y4 ) ) );
out->left = IM_RINT( left );
out->top = IM_RINT( top );
out->width = IM_RINT( right - left );
out->height = IM_RINT( bottom - top );
out->left = VIPS_RINT( left );
out->top = VIPS_RINT( top );
out->width = VIPS_RINT( right - left );
out->height = VIPS_RINT( bottom - top );
}
/* Given an area in the input image, calculate the bounding box for those
@ -222,8 +226,8 @@ transform_rect( const VipsTransformation *trn, transform_fn transform,
*/
void
vips__transform_forward_rect( const VipsTransformation *trn,
const Rect *in, /* In input space */
Rect *out ) /* In output space */
const VipsRect *in, /* In input space */
VipsRect *out ) /* In output space */
{
transform_rect( trn, vips__transform_forward_point, in, out );
}
@ -233,8 +237,8 @@ vips__transform_forward_rect( const VipsTransformation *trn,
*/
void
vips__transform_invert_rect( const VipsTransformation *trn,
const Rect *in, /* In output space */
Rect *out ) /* In input space */
const VipsRect *in, /* In output space */
VipsRect *out ) /* In input space */
{
transform_rect( trn, vips__transform_invert_point, in, out );
}