From ac358ff4b8b4b0783f9c9f35c6f58f3771becfce Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:01:34 +0200 Subject: [PATCH 1/7] Fix the pixel shift within reduce (#703) --- libvips/resample/reduce.c | 36 +++++-------- libvips/resample/reduceh.cpp | 80 ++++++++++++++++++----------- libvips/resample/reducev.cpp | 87 +++++++++++++++++++++----------- libvips/resample/resize.c | 2 - test/test-suite/test_resample.py | 5 +- 5 files changed, 124 insertions(+), 86 deletions(-) diff --git a/libvips/resample/reduce.c b/libvips/resample/reduce.c index b71ea36d..c5e1bed4 100644 --- a/libvips/resample/reduce.c +++ b/libvips/resample/reduce.c @@ -6,6 +6,8 @@ * - rename xshrink -> hshrink for greater consistency * 9/9/16 * - add @centre option + * 6/6/20 kleisauke + * - deprecate @centre option, it's now always on */ /* @@ -78,14 +80,10 @@ * Optional arguments: * * * @kernel: #VipsKernel to use to interpolate (default: lanczos3) - * * @centre: %gboolean use centre rather than corner sampling convention * * Reduce @in vertically by a float factor. The pixels in @out are * interpolated with a 1D mask generated by @kernel. * - * Set @centre to use centre rather than corner sampling convention. Centre - * convention can be useful to match the behaviour of other systems. - * * This is a very low-level operation: see vips_resize() for a more * convenient way to resize images. * @@ -107,14 +105,10 @@ * Optional arguments: * * * @kernel: #VipsKernel to use to interpolate (default: lanczos3) - * * @centre: %gboolean use centre rather than corner sampling convention * * Reduce @in horizontally by a float factor. The pixels in @out are * interpolated with a 1D mask generated by @kernel. * - * Set @centre to use centre rather than corner sampling convention. Centre - * convention can be useful to match the behaviour of other systems. - * * This is a very low-level operation: see vips_resize() for a more * convenient way to resize images. * @@ -136,7 +130,7 @@ typedef struct _VipsReduce { */ VipsKernel kernel; - /* Use centre rather than corner sampling convention. + /* Deprecated. */ gboolean centre; @@ -152,18 +146,16 @@ vips_reduce_build( VipsObject *object ) VipsResample *resample = VIPS_RESAMPLE( object ); VipsReduce *reduce = (VipsReduce *) object; VipsImage **t = (VipsImage **) - vips_object_local_array( object, 3 ); + vips_object_local_array( object, 2 ); if( VIPS_OBJECT_CLASS( vips_reduce_parent_class )->build( object ) ) return( -1 ); if( vips_reducev( resample->in, &t[0], reduce->vshrink, "kernel", reduce->kernel, - "centre", reduce->centre, NULL ) || vips_reduceh( t[0], &t[1], reduce->hshrink, "kernel", reduce->kernel, - "centre", reduce->centre, NULL ) || vips_image_write( t[1], resample->out ) ) return( -1 ); @@ -210,13 +202,6 @@ vips_reduce_class_init( VipsReduceClass *class ) G_STRUCT_OFFSET( VipsReduce, kernel ), VIPS_TYPE_KERNEL, VIPS_KERNEL_LANCZOS3 ); - VIPS_ARG_BOOL( class, "centre", 7, - _( "Centre" ), - _( "Use centre sampling convention" ), - VIPS_ARGUMENT_OPTIONAL_INPUT, - G_STRUCT_OFFSET( VipsReduce, centre ), - FALSE ); - /* The old names .. now use h and v everywhere. */ VIPS_ARG_DOUBLE( class, "xshrink", 8, @@ -233,6 +218,15 @@ vips_reduce_class_init( VipsReduceClass *class ) G_STRUCT_OFFSET( VipsReduce, vshrink ), 1.0, 1000000.0, 1.0 ); + /* We used to let people pick centre or corner, but it's automatic now. + */ + VIPS_ARG_BOOL( class, "centre", 7, + _( "Centre" ), + _( "Use centre sampling convention" ), + VIPS_ARGUMENT_OPTIONAL_INPUT | VIPS_ARGUMENT_DEPRECATED, + G_STRUCT_OFFSET( VipsReduce, centre ), + FALSE ); + } static void @@ -252,14 +246,10 @@ vips_reduce_init( VipsReduce *reduce ) * Optional arguments: * * * @kernel: #VipsKernel to use to interpolate (default: lanczos3) - * * @centre: %gboolean use centre rather than corner sampling convention * * Reduce @in by a pair of factors with a pair of 1D kernels. This * will not work well for shrink factors greater than three. * - * Set @centre to use centre rather than corner sampling convention. Centre - * convention can be useful to match the behaviour of other systems. - * * This is a very low-level operation: see vips_resize() for a more * convenient way to resize images. * diff --git a/libvips/resample/reduceh.cpp b/libvips/resample/reduceh.cpp index acb7bb31..1a6486cc 100644 --- a/libvips/resample/reduceh.cpp +++ b/libvips/resample/reduceh.cpp @@ -8,6 +8,9 @@ * - rename xshrink as hshrink for consistency * 9/9/16 * - add @centre option + * 6/6/20 kleisauke + * - deprecate @centre option, it's now always on + * - fix pixel shift */ /* @@ -66,14 +69,14 @@ typedef struct _VipsReduceh { */ VipsKernel kernel; - /* Use centre rather than corner sampling convention. - */ - gboolean centre; - /* Number of points in kernel. */ int n_point; + /* Horizontal displacement. + */ + double hoffset; + /* Precalculated interpolation matrices. int (used for pel * sizes up to short), and double (for all others). We go to * scale + 1 so we can round-to-nearest safely. @@ -81,6 +84,10 @@ typedef struct _VipsReduceh { int *matrixi[VIPS_TRANSFORM_SCALE + 1]; double *matrixf[VIPS_TRANSFORM_SCALE + 1]; + /* Deprecated. + */ + gboolean centre; + } VipsReduceh; typedef VipsResampleClass VipsReducehClass; @@ -315,12 +322,10 @@ vips_reduceh_gen( VipsRegion *out_region, void *seq, r->width, r->height, r->left, r->top ); #endif /*DEBUG*/ - s.left = r->left * reduceh->hshrink; + s.left = r->left * reduceh->hshrink - reduceh->hoffset; s.top = r->top; s.width = r->width * reduceh->hshrink + reduceh->n_point; s.height = r->height; - if( reduceh->centre ) - s.width += 1; if( vips_region_prepare( ir, &s ) ) return( -1 ); @@ -334,9 +339,8 @@ vips_reduceh_gen( VipsRegion *out_region, void *seq, q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); - X = r->left * reduceh->hshrink; - if( reduceh->centre ) - X += 0.5; + X = (r->left + 0.5) * reduceh->hshrink - 0.5 - + reduceh->hoffset; /* We want p0 to be the start (ie. x == 0) of the input * scanline we are reading from. We can then calculate the p we @@ -351,7 +355,7 @@ vips_reduceh_gen( VipsRegion *out_region, void *seq, ir->valid.left * ps; for( int x = 0; x < r->width; x++ ) { - int ix = (int) X; + const int ix = (int) X; VipsPel *p = p0 + ix * ps; const int sx = X * VIPS_TRANSFORM_SCALE * 2; const int six = sx & (VIPS_TRANSFORM_SCALE * 2 - 1); @@ -441,7 +445,7 @@ vips_reduceh_build( VipsObject *object ) vips_object_local_array( object, 2 ); VipsImage *in; - int width; + double width, extra_pixels; if( VIPS_OBJECT_CLASS( vips_reduceh_parent_class )->build( object ) ) return( -1 ); @@ -457,8 +461,6 @@ vips_reduceh_build( VipsObject *object ) if( reduceh->hshrink == 1 ) return( vips_image_write( in, resample->out ) ); - /* Build the tables of pre-computed coefficients. - */ reduceh->n_point = vips_reduce_get_points( reduceh->kernel, reduceh->hshrink ); g_info( "reduceh: %d point mask", reduceh->n_point ); @@ -467,6 +469,28 @@ vips_reduceh_build( VipsObject *object ) "%s", _( "reduce factor too large" ) ); return( -1 ); } + + /* Output size. We need to always round to nearest, so round(), not + * rint(). + */ + width = VIPS_ROUND_UINT( + (double) resample->in->Xsize / reduceh->hshrink ); + + /* How many pixels we are inventing in the input, -ve for + * discarding. + */ + extra_pixels = + width * reduceh->hshrink - resample->in->Xsize; + + /* If we are rounding down, we are not using some input + * pixels. We need to move the origin *inside* the input image + * by half that distance so that we discard pixels equally + * from left and right. + */ + reduceh->hoffset = (1 + extra_pixels) / 2.0; + + /* Build the tables of pre-computed coefficients. + */ for( int x = 0; x < VIPS_TRANSFORM_SCALE + 1; x++ ) { reduceh->matrixf[x] = VIPS_ARRAY( object, reduceh->n_point, double ); @@ -499,15 +523,10 @@ vips_reduceh_build( VipsObject *object ) in = t[0]; /* Add new pixels around the input so we can interpolate at the edges. - * In centre mode, we read 0.5 pixels more to the right, so we must - * enlarge a little further. */ - width = in->Xsize + reduceh->n_point - 1; - if( reduceh->centre ) - width += 1; if( vips_embed( in, &t[1], reduceh->n_point / 2 - 1, 0, - width, in->Ysize, + in->Xsize + reduceh->n_point, in->Ysize, "extend", VIPS_EXTEND_COPY, (void *) NULL ) ) return( -1 ); @@ -524,8 +543,7 @@ vips_reduceh_build( VipsObject *object ) * example, vipsthumbnail knows the true reduce factor (including the * fractional part), we just see the integer part here. */ - resample->out->Xsize = VIPS_ROUND_UINT( - resample->in->Xsize / reduceh->hshrink ); + resample->out->Xsize = width; if( resample->out->Xsize <= 0 ) { vips_error( object_class->nickname, "%s", _( "image has shrunk to nothing" ) ); @@ -574,20 +592,13 @@ vips_reduceh_class_init( VipsReducehClass *reduceh_class ) G_STRUCT_OFFSET( VipsReduceh, hshrink ), 1, 1000000, 1 ); - VIPS_ARG_ENUM( reduceh_class, "kernel", 3, + VIPS_ARG_ENUM( reduceh_class, "kernel", 4, _( "Kernel" ), _( "Resampling kernel" ), VIPS_ARGUMENT_OPTIONAL_INPUT, G_STRUCT_OFFSET( VipsReduceh, kernel ), VIPS_TYPE_KERNEL, VIPS_KERNEL_LANCZOS3 ); - VIPS_ARG_BOOL( reduceh_class, "centre", 7, - _( "Centre" ), - _( "Use centre sampling convention" ), - VIPS_ARGUMENT_OPTIONAL_INPUT, - G_STRUCT_OFFSET( VipsReduceh, centre ), - FALSE ); - /* Old name. */ VIPS_ARG_DOUBLE( reduceh_class, "xshrink", 3, @@ -597,6 +608,15 @@ vips_reduceh_class_init( VipsReducehClass *reduceh_class ) G_STRUCT_OFFSET( VipsReduceh, hshrink ), 1, 1000000, 1 ); + /* We used to let people pick centre or corner, but it's automatic now. + */ + VIPS_ARG_BOOL( reduceh_class, "centre", 7, + _( "Centre" ), + _( "Use centre sampling convention" ), + VIPS_ARGUMENT_OPTIONAL_INPUT | VIPS_ARGUMENT_DEPRECATED, + G_STRUCT_OFFSET( VipsReduceh, centre ), + FALSE ); + } static void diff --git a/libvips/resample/reducev.cpp b/libvips/resample/reducev.cpp index a64a052e..965d2a06 100644 --- a/libvips/resample/reducev.cpp +++ b/libvips/resample/reducev.cpp @@ -17,6 +17,9 @@ * - add @centre option * 7/3/17 * - add a seq line cache + * 6/6/20 kleisauke + * - deprecate @centre option, it's now always on + * - fix pixel shift */ /* @@ -104,14 +107,14 @@ typedef struct _VipsReducev { */ VipsKernel kernel; - /* Use centre rather than corner sampling convention. - */ - gboolean centre; - /* Number of points in kernel. */ int n_point; + /* Vertical displacement. + */ + double voffset; + /* Precalculated interpolation matrices. int (used for pel * sizes up to short), and double (for all others). We go to * scale + 1 so we can round-to-nearest safely. @@ -128,6 +131,10 @@ typedef struct _VipsReducev { int n_pass; Pass pass[MAX_PASS]; + /* Deprecated. + */ + gboolean centre; + } VipsReducev; typedef VipsResampleClass VipsReducevClass; @@ -531,22 +538,22 @@ vips_reducev_gen( VipsRegion *out_region, void *vseq, #endif /*DEBUG*/ s.left = r->left; - s.top = r->top * reducev->vshrink; + s.top = r->top * reducev->vshrink - reducev->voffset; s.width = r->width; s.height = r->height * reducev->vshrink + reducev->n_point; - if( reducev->centre ) - s.height += 1; if( vips_region_prepare( ir, &s ) ) return( -1 ); VIPS_GATE_START( "vips_reducev_gen: work" ); + double Y = (r->top + 0.5) * reducev->vshrink - 0.5 - + reducev->voffset; + for( int y = 0; y < r->height; y ++ ) { VipsPel *q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); - const double Y = (r->top + y) * reducev->vshrink + - (reducev->centre ? 0.5 : 0.0); - VipsPel *p = VIPS_REGION_ADDR( ir, r->left, (int) Y ); + const int py = (int) Y; + VipsPel *p = VIPS_REGION_ADDR( ir, r->left, py ); const int sy = Y * VIPS_TRANSFORM_SCALE * 2; const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); const int ty = (siy + 1) >> 1; @@ -606,13 +613,15 @@ vips_reducev_gen( VipsRegion *out_region, void *vseq, case VIPS_FORMAT_DPCOMPLEX: case VIPS_FORMAT_DOUBLE: reducev_notab( reducev, - q, p, ne, lskip, Y - (int) Y ); + q, p, ne, lskip, Y - py ); break; default: g_assert_not_reached(); break; } + + Y += reducev->vshrink; } VIPS_GATE_STOP( "vips_reducev_gen: work" ); @@ -646,9 +655,8 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, s.left = r->left; s.top = r->top * reducev->vshrink; s.width = r->width; - s.height = r->height * reducev->vshrink + reducev->n_point; - if( reducev->centre ) - s.height += 1; + s.height = r->height * reducev->vshrink + reducev->n_point - + reducev->voffset; if( vips_region_prepare( ir, &s ) ) return( -1 ); @@ -663,11 +671,12 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, VIPS_GATE_START( "vips_reducev_vector_gen: work" ); + double Y = (r->top + 0.5) * reducev->vshrink - 0.5 - + reducev->voffset; + for( int y = 0; y < r->height; y ++ ) { VipsPel *q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); - const double Y = (r->top + y) * reducev->vshrink + - (reducev->centre ? 0.5 : 0.0); const int py = (int) Y; const int sy = Y * VIPS_TRANSFORM_SCALE * 2; const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); @@ -682,7 +691,7 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, printf( "first column of pixel values:\n" ); for( int i = 0; i < reducev->n_point; i++ ) printf( "\t%d - %d\n", i, - *VIPS_REGION_ADDR( ir, r->left, r->top + y + i ) ); + *VIPS_REGION_ADDR( ir, r->left, py ) ); #endif /*DEBUG_PIXELS*/ /* We run our n passes to generate this scanline. @@ -709,6 +718,8 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, printf( "pixel result:\n" ); printf( "\t%d\n", *q ); #endif /*DEBUG_PIXELS*/ + + Y += reducev->vshrink; } VIPS_GATE_STOP( "vips_reducev_vector_gen: work" ); @@ -830,7 +841,7 @@ vips_reducev_build( VipsObject *object ) VipsImage **t = (VipsImage **) vips_object_local_array( object, 4 ); VipsImage *in; - int height; + double height, extra_pixels; if( VIPS_OBJECT_CLASS( vips_reducev_parent_class )->build( object ) ) return( -1 ); @@ -855,6 +866,25 @@ vips_reducev_build( VipsObject *object ) return( -1 ); } + /* Output size. We need to always round to nearest, so round(), not + * rint(). + */ + height = VIPS_ROUND_UINT( + (double) resample->in->Ysize / reducev->vshrink ); + + /* How many pixels we are inventing in the input, -ve for + * discarding. + */ + extra_pixels = + height * reducev->vshrink - resample->in->Ysize; + + /* If we are rounding down, we are not using some input + * pixels. We need to move the origin *inside* the input image + * by half that distance so that we discard pixels equally + * from left and right. + */ + reducev->voffset = (1 + extra_pixels) / 2.0; + /* Unpack for processing. */ if( vips_image_decode( in, &t[0] ) ) @@ -863,12 +893,9 @@ vips_reducev_build( VipsObject *object ) /* Add new pixels around the input so we can interpolate at the edges. */ - height = in->Ysize + reducev->n_point - 1; - if( reducev->centre ) - height += 1; if( vips_embed( in, &t[1], 0, reducev->n_point / 2 - 1, - in->Xsize, height, + in->Xsize, in->Ysize + reducev->n_point, "extend", VIPS_EXTEND_COPY, (void *) NULL ) ) return( -1 ); @@ -939,13 +966,6 @@ vips_reducev_class_init( VipsReducevClass *reducev_class ) G_STRUCT_OFFSET( VipsReducev, kernel ), VIPS_TYPE_KERNEL, VIPS_KERNEL_LANCZOS3 ); - VIPS_ARG_BOOL( reducev_class, "centre", 7, - _( "Centre" ), - _( "Use centre sampling convention" ), - VIPS_ARGUMENT_OPTIONAL_INPUT, - G_STRUCT_OFFSET( VipsReducev, centre ), - FALSE ); - /* Old name. */ VIPS_ARG_DOUBLE( reducev_class, "yshrink", 3, @@ -955,6 +975,15 @@ vips_reducev_class_init( VipsReducevClass *reducev_class ) G_STRUCT_OFFSET( VipsReducev, vshrink ), 1, 1000000, 1 ); + /* We used to let people pick centre or corner, but it's automatic now. + */ + VIPS_ARG_BOOL( reducev_class, "centre", 7, + _( "Centre" ), + _( "Use centre sampling convention" ), + VIPS_ARGUMENT_OPTIONAL_INPUT | VIPS_ARGUMENT_DEPRECATED, + G_STRUCT_OFFSET( VipsReducev, centre ), + FALSE ); + } static void diff --git a/libvips/resample/resize.c b/libvips/resample/resize.c index 28235351..3048c682 100644 --- a/libvips/resample/resize.c +++ b/libvips/resample/resize.c @@ -232,7 +232,6 @@ vips_resize_build( VipsObject *object ) g_info( "residual reducev by %g", vscale ); if( vips_reducev( in, &t[2], 1.0 / vscale, "kernel", resize->kernel, - "centre", TRUE, NULL ) ) return( -1 ); in = t[2]; @@ -243,7 +242,6 @@ vips_resize_build( VipsObject *object ) hscale ); if( vips_reduceh( in, &t[3], 1.0 / hscale, "kernel", resize->kernel, - "centre", TRUE, NULL ) ) return( -1 ); in = t[3]; diff --git a/test/test-suite/test_resample.py b/test/test-suite/test_resample.py index 04c81580..6453d62e 100644 --- a/test/test-suite/test_resample.py +++ b/test/test-suite/test_resample.py @@ -82,8 +82,9 @@ class TestResample: for fac in [1, 1.1, 1.5, 1.999]: for fmt in all_formats: - for kernel in ["nearest", "linear", - "cubic", "lanczos2", "lanczos3"]: + # TODO: Add nearest kernel when https://github.com/libvips/libvips/issues/1518 is done. + # (running the test suite with VIPS_NOVECTOR=1 should also work) + for kernel in ["linear", "cubic", "lanczos2", "lanczos3"]: x = im.cast(fmt) r = x.reduce(fac, fac, kernel=kernel) d = abs(r.avg() - im.avg()) From b6e4e9e74b210b16edac72ce3cd9e2b55a2e4ea0 Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:12:37 +0200 Subject: [PATCH 2/7] Speed up the mask construction for uchar/ushort images By not calling vips_vector_to_fixed_point repeatedly. --- libvips/resample/reducev.cpp | 64 ++++++++++++++++-------------------- 1 file changed, 29 insertions(+), 35 deletions(-) diff --git a/libvips/resample/reducev.cpp b/libvips/resample/reducev.cpp index 965d2a06..a5137c76 100644 --- a/libvips/resample/reducev.cpp +++ b/libvips/resample/reducev.cpp @@ -20,6 +20,7 @@ * 6/6/20 kleisauke * - deprecate @centre option, it's now always on * - fix pixel shift + * - speed up the mask construction for uchar/ushort images */ /* @@ -737,41 +738,7 @@ vips_reducev_raw( VipsReducev *reducev, VipsImage *in, VipsImage **out ) VipsGenerateFn generate; - /* Build masks. - */ - for( int y = 0; y < VIPS_TRANSFORM_SCALE + 1; y++ ) { - reducev->matrixf[y] = - VIPS_ARRAY( NULL, reducev->n_point, double ); - if( !reducev->matrixf[y] ) - return( -1 ); - - vips_reduce_make_mask( reducev->matrixf[y], - reducev->kernel, reducev->vshrink, - (float) y / VIPS_TRANSFORM_SCALE ); - -#ifdef DEBUG - printf( "%6.2g", (double) y / VIPS_TRANSFORM_SCALE ); - for( int i = 0; i < reducev->n_point; i++ ) - printf( ", %6.2g", reducev->matrixf[y][i] ); - printf( "\n" ); -#endif /*DEBUG*/ - } - - /* uchar and ushort need an int version of the masks. - */ - if( VIPS_IMAGE_SIZEOF_ELEMENT( in ) <= 2 ) - for( int y = 0; y < VIPS_TRANSFORM_SCALE + 1; y++ ) { - reducev->matrixi[y] = - VIPS_ARRAY( NULL, reducev->n_point, int ); - if( !reducev->matrixi[y] ) - return( -1 ); - - vips_vector_to_fixed_point( - reducev->matrixf[y], reducev->matrixi[y], - reducev->n_point, VIPS_INTERPOLATE_SCALE ); - } - - /* And we need an 2.6 version if we will use the vector path. + /* We need an 2.6 version if we will use the vector path. */ if( in->BandFmt == VIPS_FORMAT_UCHAR && vips_vector_isenabled() ) @@ -885,6 +852,33 @@ vips_reducev_build( VipsObject *object ) */ reducev->voffset = (1 + extra_pixels) / 2.0; + /* Build the tables of pre-computed coefficients. + */ + for( int y = 0; y < VIPS_TRANSFORM_SCALE + 1; y++ ) { + reducev->matrixf[y] = + VIPS_ARRAY( NULL, reducev->n_point, double ); + reducev->matrixi[y] = + VIPS_ARRAY( NULL, reducev->n_point, int ); + if( !reducev->matrixf[y] || + !reducev->matrixi[y] ) + return( -1 ); + + vips_reduce_make_mask( reducev->matrixf[y], + reducev->kernel, reducev->vshrink, + (float) y / VIPS_TRANSFORM_SCALE ); + + for( int i = 0; i < reducev->n_point; i++ ) + reducev->matrixi[y][i] = reducev->matrixf[y][i] * + VIPS_INTERPOLATE_SCALE; + +#ifdef DEBUG + printf( "vips_reducev_build: mask %d\n ", y ); + for( int i = 0; i < reducev->n_point; i++ ) + printf( "%d ", reducev->matrixi[y][i] ); + printf( "\n" ); +#endif /*DEBUG*/ + } + /* Unpack for processing. */ if( vips_image_decode( in, &t[0] ) ) From dfdf899c92c05a4caf7f4f0ef225129a92e41d8d Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:16:24 +0200 Subject: [PATCH 3/7] Ensure reducev is THINSTRIP In line with reduceh. --- libvips/resample/reducev.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/libvips/resample/reducev.cpp b/libvips/resample/reducev.cpp index a5137c76..a704f739 100644 --- a/libvips/resample/reducev.cpp +++ b/libvips/resample/reducev.cpp @@ -765,7 +765,7 @@ vips_reducev_raw( VipsReducev *reducev, VipsImage *in, VipsImage **out ) *out = vips_image_new(); if( vips_image_pipelinev( *out, - VIPS_DEMAND_STYLE_FATSTRIP, in, (void *) NULL ) ) + VIPS_DEMAND_STYLE_THINSTRIP, in, (void *) NULL ) ) return( -1 ); /* Size output. We need to always round to nearest, so round(), not From 369b09809670f9744262d6e708e417796891f73a Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:21:15 +0200 Subject: [PATCH 4/7] Prefer an immediate calculation where possible In line with reducev. --- libvips/resample/reduceh.cpp | 27 ++++++++------------------- 1 file changed, 8 insertions(+), 19 deletions(-) diff --git a/libvips/resample/reduceh.cpp b/libvips/resample/reduceh.cpp index 1a6486cc..73ddab38 100644 --- a/libvips/resample/reduceh.cpp +++ b/libvips/resample/reduceh.cpp @@ -201,13 +201,11 @@ reduceh_signed_int_tab( VipsReduceh *reduceh, for( int z = 0; z < bands; z++ ) { int sum; - sum = reduce_sum( in, bands, cx, n ); + sum = reduce_sum( in + z, bands, cx, n ); sum = signed_fixed_round( sum ); sum = VIPS_CLIP( min_value, sum, max_value ); out[z] = sum; - - in += 1; } } @@ -223,10 +221,8 @@ reduceh_float_tab( VipsReduceh *reduceh, const T* restrict in = (T *) pin; const int n = reduceh->n_point; - for( int z = 0; z < bands; z++ ) { - out[z] = reduce_sum( in, bands, cx, n ); - in += 1; - } + for( int z = 0; z < bands; z++ ) + out[z] = reduce_sum( in + z, bands, cx, n ); } /* 32-bit int output needs a double intermediate. @@ -245,10 +241,8 @@ reduceh_unsigned_int32_tab( VipsReduceh *reduceh, for( int z = 0; z < bands; z++ ) { double sum; - sum = reduce_sum( in, bands, cx, n ); - out[z] = VIPS_CLIP( 0, sum, max_value ); - - in += 1; + sum = reduce_sum( in + z, bands, cx, n ); + out[z] = VIPS_CLIP( 0, sum, max_value ); } } @@ -265,11 +259,9 @@ reduceh_signed_int32_tab( VipsReduceh *reduceh, for( int z = 0; z < bands; z++ ) { double sum; - sum = reduce_sum( in, bands, cx, n ); + sum = reduce_sum( in + z, bands, cx, n ); sum = VIPS_CLIP( min_value, sum, max_value ); out[z] = sum; - - in += 1; } } @@ -289,11 +281,8 @@ reduceh_notab( VipsReduceh *reduceh, vips_reduce_make_mask( cx, reduceh->kernel, reduceh->hshrink, x ); - for( int z = 0; z < bands; z++ ) { - out[z] = reduce_sum( in, bands, cx, n ); - - in += 1; - } + for( int z = 0; z < bands; z++ ) + out[z] = reduce_sum( in + z, bands, cx, n ); } /* Tried a vector path (see reducev) but it was slower. The vectors for From c0ed106079fe1d2989530eae780c3469d6599e3b Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:25:03 +0200 Subject: [PATCH 5/7] Formatting and whitespace changes --- libvips/resample/reduceh.cpp | 2 +- libvips/resample/reducev.cpp | 8 ++++---- libvips/resample/templates.h | 8 ++++---- 3 files changed, 9 insertions(+), 9 deletions(-) diff --git a/libvips/resample/reduceh.cpp b/libvips/resample/reduceh.cpp index 73ddab38..07ddff32 100644 --- a/libvips/resample/reduceh.cpp +++ b/libvips/resample/reduceh.cpp @@ -498,7 +498,7 @@ vips_reduceh_build( VipsObject *object ) VIPS_INTERPOLATE_SCALE; #ifdef DEBUG - printf( "vips_reduceh_build: mask %d\n ", x ); + printf( "vips_reduceh_build: mask %d\n ", x ); for( int i = 0; i < reduceh->n_point; i++ ) printf( "%d ", reduceh->matrixi[x][i] ); printf( "\n" ); diff --git a/libvips/resample/reducev.cpp b/libvips/resample/reducev.cpp index a704f739..7f876433 100644 --- a/libvips/resample/reducev.cpp +++ b/libvips/resample/reducev.cpp @@ -553,8 +553,8 @@ vips_reducev_gen( VipsRegion *out_region, void *vseq, for( int y = 0; y < r->height; y ++ ) { VipsPel *q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); - const int py = (int) Y; - VipsPel *p = VIPS_REGION_ADDR( ir, r->left, py ); + const int py = (int) Y; + VipsPel *p = VIPS_REGION_ADDR( ir, r->left, py ); const int sy = Y * VIPS_TRANSFORM_SCALE * 2; const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); const int ty = (siy + 1) >> 1; @@ -678,7 +678,7 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, for( int y = 0; y < r->height; y ++ ) { VipsPel *q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); - const int py = (int) Y; + const int py = (int) Y; const int sy = Y * VIPS_TRANSFORM_SCALE * 2; const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); const int ty = (siy + 1) >> 1; @@ -872,7 +872,7 @@ vips_reducev_build( VipsObject *object ) VIPS_INTERPOLATE_SCALE; #ifdef DEBUG - printf( "vips_reducev_build: mask %d\n ", y ); + printf( "vips_reducev_build: mask %d\n ", y ); for( int i = 0; i < reducev->n_point; i++ ) printf( "%d ", reducev->matrixi[y][i] ); printf( "\n" ); diff --git a/libvips/resample/templates.h b/libvips/resample/templates.h index b64737df..d60144ba 100644 --- a/libvips/resample/templates.h +++ b/libvips/resample/templates.h @@ -322,7 +322,7 @@ calculate_coefficients_triangle( double *c, /* Needs to be in sync with vips_reduce_get_points(). */ const int n_points = 2 * rint( shrink ) + 1; - const double half = x + n_points / 2.0 - 1; + const double half = x + n_points / 2.0 - 1; int i; double sum; @@ -360,7 +360,7 @@ calculate_coefficients_cubic( double *c, /* Needs to be in sync with vips_reduce_get_points(). */ const int n_points = 2 * rint( 2 * shrink ) + 1; - const double half = x + n_points / 2.0 - 1; + const double half = x + n_points / 2.0 - 1; int i; double sum; @@ -409,14 +409,14 @@ calculate_coefficients_lanczos( double *c, /* Needs to be in sync with vips_reduce_get_points(). */ const int n_points = 2 * rint( a * shrink ) + 1; - const double half = x + n_points / 2.0 - 1; + const double half = x + n_points / 2.0 - 1; int i; double sum; sum = 0; for( i = 0; i < n_points; i++ ) { - double xp = (i - half) / shrink; + const double xp = (i - half) / shrink; double l; From ac30bad6953da8b4850bb17ef51bbedf7cfd69c3 Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:29:57 +0200 Subject: [PATCH 6/7] Remove round-to-nearest behaviour It seems that it generates the same image, with or without this change. Tested with https://github.com/kleisauke/vips-issue-703. --- libvips/resample/bicubic.cpp | 25 +++++++++---------------- libvips/resample/reduceh.cpp | 17 ++++++++--------- libvips/resample/reducev.cpp | 30 ++++++++++++++---------------- 3 files changed, 31 insertions(+), 41 deletions(-) diff --git a/libvips/resample/bicubic.cpp b/libvips/resample/bicubic.cpp index 97e4a403..7388ab06 100644 --- a/libvips/resample/bicubic.cpp +++ b/libvips/resample/bicubic.cpp @@ -79,16 +79,15 @@ typedef VipsInterpolate VipsInterpolateBicubic; typedef VipsInterpolateClass VipsInterpolateBicubicClass; /* Precalculated interpolation matrices. int (used for pel - * sizes up to short), and double (for all others). We go to - * scale + 1 so we can round-to-nearest safely. + * sizes up to short), and double (for all others). */ /* We could keep a large set of 2d 4x4 matricies, but this actually * works out slower since for many resizes the thing will no longer * fit in L1. */ -static int vips_bicubic_matrixi[VIPS_TRANSFORM_SCALE + 1][4]; -static double vips_bicubic_matrixf[VIPS_TRANSFORM_SCALE + 1][4]; +static int vips_bicubic_matrixi[VIPS_TRANSFORM_SCALE][4]; +static double vips_bicubic_matrixf[VIPS_TRANSFORM_SCALE][4]; /* We need C linkage for this. */ @@ -498,19 +497,13 @@ static void vips_interpolate_bicubic_interpolate( VipsInterpolate *interpolate, void *out, VipsRegion *in, double x, double y ) { - /* Find the mask index. We round-to-nearest, so we need to generate - * indexes in 0 to VIPS_TRANSFORM_SCALE, 2^n + 1 values. We multiply - * by 2 more than we need to, add one, mask, then shift down again to - * get the extra range. + /* Find the mask index. */ - const int sx = x * VIPS_TRANSFORM_SCALE * 2; - const int sy = y * VIPS_TRANSFORM_SCALE * 2; + const int sx = x * VIPS_TRANSFORM_SCALE; + const int sy = y * VIPS_TRANSFORM_SCALE; - const int six = sx & (VIPS_TRANSFORM_SCALE * 2 - 1); - const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); - - const int tx = (six + 1) >> 1; - const int ty = (siy + 1) >> 1; + const int tx = sx & (VIPS_TRANSFORM_SCALE - 1); + const int ty = sy & (VIPS_TRANSFORM_SCALE - 1); /* We know x/y are always positive, so we can just (int) them. */ @@ -643,7 +636,7 @@ vips_interpolate_bicubic_class_init( VipsInterpolateBicubicClass *iclass ) /* Build the tables of pre-computed coefficients. */ - for( int x = 0; x < VIPS_TRANSFORM_SCALE + 1; x++ ) { + for( int x = 0; x < VIPS_TRANSFORM_SCALE; x++ ) { calculate_coefficients_catmull( vips_bicubic_matrixf[x], (float) x / VIPS_TRANSFORM_SCALE ); diff --git a/libvips/resample/reduceh.cpp b/libvips/resample/reduceh.cpp index 07ddff32..cc829c01 100644 --- a/libvips/resample/reduceh.cpp +++ b/libvips/resample/reduceh.cpp @@ -11,6 +11,7 @@ * 6/6/20 kleisauke * - deprecate @centre option, it's now always on * - fix pixel shift + * - remove unnecessary round-to-nearest behaviour */ /* @@ -78,11 +79,10 @@ typedef struct _VipsReduceh { double hoffset; /* Precalculated interpolation matrices. int (used for pel - * sizes up to short), and double (for all others). We go to - * scale + 1 so we can round-to-nearest safely. + * sizes up to short), and double (for all others). */ - int *matrixi[VIPS_TRANSFORM_SCALE + 1]; - double *matrixf[VIPS_TRANSFORM_SCALE + 1]; + int *matrixi[VIPS_TRANSFORM_SCALE]; + double *matrixf[VIPS_TRANSFORM_SCALE]; /* Deprecated. */ @@ -320,7 +320,7 @@ vips_reduceh_gen( VipsRegion *out_region, void *seq, VIPS_GATE_START( "vips_reduceh_gen: work" ); - for( int y = 0; y < r->height; y ++ ) { + for( int y = 0; y < r->height; y++ ) { VipsPel *p0; VipsPel *q; @@ -346,9 +346,8 @@ vips_reduceh_gen( VipsRegion *out_region, void *seq, for( int x = 0; x < r->width; x++ ) { const int ix = (int) X; VipsPel *p = p0 + ix * ps; - const int sx = X * VIPS_TRANSFORM_SCALE * 2; - const int six = sx & (VIPS_TRANSFORM_SCALE * 2 - 1); - const int tx = (six + 1) >> 1; + const int sx = X * VIPS_TRANSFORM_SCALE; + const int tx = sx & (VIPS_TRANSFORM_SCALE - 1); const int *cxi = reduceh->matrixi[tx]; const double *cxf = reduceh->matrixf[tx]; @@ -480,7 +479,7 @@ vips_reduceh_build( VipsObject *object ) /* Build the tables of pre-computed coefficients. */ - for( int x = 0; x < VIPS_TRANSFORM_SCALE + 1; x++ ) { + for( int x = 0; x < VIPS_TRANSFORM_SCALE; x++ ) { reduceh->matrixf[x] = VIPS_ARRAY( object, reduceh->n_point, double ); reduceh->matrixi[x] = diff --git a/libvips/resample/reducev.cpp b/libvips/resample/reducev.cpp index 7f876433..fa965134 100644 --- a/libvips/resample/reducev.cpp +++ b/libvips/resample/reducev.cpp @@ -21,6 +21,7 @@ * - deprecate @centre option, it's now always on * - fix pixel shift * - speed up the mask construction for uchar/ushort images + * - remove unnecessary round-to-nearest behaviour */ /* @@ -117,15 +118,14 @@ typedef struct _VipsReducev { double voffset; /* Precalculated interpolation matrices. int (used for pel - * sizes up to short), and double (for all others). We go to - * scale + 1 so we can round-to-nearest safely. + * sizes up to short), and double (for all others). */ - int *matrixi[VIPS_TRANSFORM_SCALE + 1]; - double *matrixf[VIPS_TRANSFORM_SCALE + 1]; + int *matrixi[VIPS_TRANSFORM_SCALE]; + double *matrixf[VIPS_TRANSFORM_SCALE]; /* And another set for orc: we want 2.6 precision. */ - int *matrixo[VIPS_TRANSFORM_SCALE + 1]; + int *matrixo[VIPS_TRANSFORM_SCALE]; /* The passes we generate for this mask. */ @@ -154,7 +154,7 @@ vips_reducev_finalize( GObject *gobject ) for( int i = 0; i < reducev->n_pass; i++ ) VIPS_FREEF( vips_vector_free, reducev->pass[i].vector ); reducev->n_pass = 0; - for( int i = 0; i < VIPS_TRANSFORM_SCALE + 1; i++ ) { + for( int i = 0; i < VIPS_TRANSFORM_SCALE; i++ ) { VIPS_FREE( reducev->matrixf[i] ); VIPS_FREE( reducev->matrixi[i] ); VIPS_FREE( reducev->matrixo[i] ); @@ -550,14 +550,13 @@ vips_reducev_gen( VipsRegion *out_region, void *vseq, double Y = (r->top + 0.5) * reducev->vshrink - 0.5 - reducev->voffset; - for( int y = 0; y < r->height; y ++ ) { + for( int y = 0; y < r->height; y++ ) { VipsPel *q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); const int py = (int) Y; VipsPel *p = VIPS_REGION_ADDR( ir, r->left, py ); - const int sy = Y * VIPS_TRANSFORM_SCALE * 2; - const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); - const int ty = (siy + 1) >> 1; + const int sy = Y * VIPS_TRANSFORM_SCALE; + const int ty = sy & (VIPS_TRANSFORM_SCALE - 1); const int *cyi = reducev->matrixi[ty]; const double *cyf = reducev->matrixf[ty]; const int lskip = VIPS_REGION_LSKIP( ir ); @@ -675,13 +674,12 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, double Y = (r->top + 0.5) * reducev->vshrink - 0.5 - reducev->voffset; - for( int y = 0; y < r->height; y ++ ) { + for( int y = 0; y < r->height; y++ ) { VipsPel *q = VIPS_REGION_ADDR( out_region, r->left, r->top + y ); const int py = (int) Y; - const int sy = Y * VIPS_TRANSFORM_SCALE * 2; - const int siy = sy & (VIPS_TRANSFORM_SCALE * 2 - 1); - const int ty = (siy + 1) >> 1; + const int sy = Y * VIPS_TRANSFORM_SCALE; + const int ty = sy & (VIPS_TRANSFORM_SCALE - 1); const int *cyo = reducev->matrixo[ty]; #ifdef DEBUG_PIXELS @@ -742,7 +740,7 @@ vips_reducev_raw( VipsReducev *reducev, VipsImage *in, VipsImage **out ) */ if( in->BandFmt == VIPS_FORMAT_UCHAR && vips_vector_isenabled() ) - for( int y = 0; y < VIPS_TRANSFORM_SCALE + 1; y++ ) { + for( int y = 0; y < VIPS_TRANSFORM_SCALE; y++ ) { reducev->matrixo[y] = VIPS_ARRAY( NULL, reducev->n_point, int ); if( !reducev->matrixo[y] ) @@ -854,7 +852,7 @@ vips_reducev_build( VipsObject *object ) /* Build the tables of pre-computed coefficients. */ - for( int y = 0; y < VIPS_TRANSFORM_SCALE + 1; y++ ) { + for( int y = 0; y < VIPS_TRANSFORM_SCALE; y++ ) { reducev->matrixf[y] = VIPS_ARRAY( NULL, reducev->n_point, double ); reducev->matrixi[y] = From d7a735400aa5fa6dc30da9cb891570f036850255 Mon Sep 17 00:00:00 2001 From: Kleis Auke Wolthuizen Date: Sat, 6 Jun 2020 14:50:32 +0200 Subject: [PATCH 7/7] reducev: Fix undefined-behaviour within the vector path Found by UBSan. --- libvips/resample/reducev.cpp | 5 ++--- test/test-suite/test_resample.py | 5 ++--- 2 files changed, 4 insertions(+), 6 deletions(-) diff --git a/libvips/resample/reducev.cpp b/libvips/resample/reducev.cpp index fa965134..396069f2 100644 --- a/libvips/resample/reducev.cpp +++ b/libvips/resample/reducev.cpp @@ -653,10 +653,9 @@ vips_reducev_vector_gen( VipsRegion *out_region, void *vseq, #endif /*DEBUG_PIXELS*/ s.left = r->left; - s.top = r->top * reducev->vshrink; + s.top = r->top * reducev->vshrink - reducev->voffset; s.width = r->width; - s.height = r->height * reducev->vshrink + reducev->n_point - - reducev->voffset; + s.height = r->height * reducev->vshrink + reducev->n_point; if( vips_region_prepare( ir, &s ) ) return( -1 ); diff --git a/test/test-suite/test_resample.py b/test/test-suite/test_resample.py index 6453d62e..04c81580 100644 --- a/test/test-suite/test_resample.py +++ b/test/test-suite/test_resample.py @@ -82,9 +82,8 @@ class TestResample: for fac in [1, 1.1, 1.5, 1.999]: for fmt in all_formats: - # TODO: Add nearest kernel when https://github.com/libvips/libvips/issues/1518 is done. - # (running the test suite with VIPS_NOVECTOR=1 should also work) - for kernel in ["linear", "cubic", "lanczos2", "lanczos3"]: + for kernel in ["nearest", "linear", + "cubic", "lanczos2", "lanczos3"]: x = im.cast(fmt) r = x.reduce(fac, fac, kernel=kernel) d = abs(r.avg() - im.avg())