improve sharpen for screen work

it was a bit granular (since it was designed for offset work), make sharpen
easier to adjust for smaller effects

- new 'sigma' parameter
- greater parameter range
- more self-tests

ChangeLog

@@ -7,6 +7,9 @@
 - added @page param to magickload 
 - matload is more specific (thanks bithive)
 - lower mem use for progressive jpg decode
+- sharpen has a new @sigma param, @radius is deprecated
+- sharpen allows a much greater range of parameters
+- better handling of deprecated args in python
 
 27/1/16 started 8.2.3
 - fix a crash with SPARC byte-order labq vips images

TODO

@@ -1,3 +1,17 @@
+try
+
+	$ vips cast babe.jpg x.v ushort
+	$ vips colourspace x.v x2.v labs
+	$ vips avg x2.v
+	27.403674
+
+	$ vips cast babe.jpg x.v short
+	$ vips colourspace x.v x2.v labs
+	$ vips avg x2.v
+	10356.775608
+
+  how odd! cause of break in test suite
+
 - could load pdf thumbnails?
 
 - new vips_reduce:

libvips/convolution/sharpen.c

@@ -1,11 +1,4 @@
 /* Cored sharpen of LABQ image.
- * 
- * Usage:
- *
- *   	int im_sharpen( IMAGE *in, IMAGE *out, 
- *		int mask_size, 
- *		int x1, int x2,
- *		double m1, double m2 )
  * 
  * Returns 0 on success and -1 on error
  *
@@ -42,6 +35,9 @@
  * 13/11/13
  * 	- redo as a class
  * 	- does any type, any number of bands
+ * 24/2/16
+ * 	- swap "radius" for "sigma", allows finer control
+ * 	- allow a much greater range of parameters
  */
 
 /*
@@ -92,29 +88,21 @@ typedef struct _VipsSharpen {
 	VipsImage *in;
 	VipsImage *out;
 
-	int radius; 
+	double sigma; 
 	double x1;
 	double y2;
 	double y3;
 	double m1;
 	double m2;
 
-	double x2;
-	double x3;
-
-	/* Parameters scaled up to int.
-	 *
-	 * We need indexes in the range [-x3,x2], so add x3 to 
-	 * indexes before starting to index table.
-	 */
-	int ix1;
-	int ix2;
-	int ix3;
-
 	/* The lut we build.
 	 */
 	int *lut;		
 
+	/* We used to have a radius control.
+	 */
+	int radius;
+
 } VipsSharpen;
 
 typedef VipsOperationClass VipsSharpenClass;
@@ -128,8 +116,6 @@ vips_sharpen_generate( VipsRegion *or,
 	VipsRegion **in = (VipsRegion **) vseq;
 	VipsSharpen *sharpen = (VipsSharpen *) b;
 	VipsRect *r = &or->valid;
-	int ix3 = sharpen->ix3;
-	int range = sharpen->ix2 + sharpen->ix3;
 	int *lut = sharpen->lut;
 
 	int x, y; 
@@ -141,43 +127,33 @@ vips_sharpen_generate( VipsRegion *or,
 	VIPS_GATE_START( "vips_sharpen_generate: work" ); 
 
 	for( y = 0; y < r->height; y++ ) {
-		short *p1 = (short *) 
+		short *p1 = (short * restrict ) 
 			VIPS_REGION_ADDR( in[0], r->left, r->top + y );
-		short *p2 = (short *) 
+		short *p2 = (short * restrict ) 
 			VIPS_REGION_ADDR( in[1], r->left, r->top + y );
-		short *q = (short *) 
+		short *q = (short * restrict ) 
 			VIPS_REGION_ADDR( or, r->left, r->top + y );
 
 		for( x = 0; x < r->width; x++ ) {
 			int v1 = p1[x];
 			int v2 = p2[x];
 
-			/* v2 is the area average. If this is zero, then we 
+			/* Our LUT is -32768 - 32767. For the v1, v2
-			 * pass the original image through unaltered.
+			 * difference to be in this range, both must be 0 -
+			 * 32767.
 			 */
-			if( v2 == 0 ) 
+			int d = (v1 & 0x8fff) - (v2 & 0x8fff);
-				q[x] = v1;
-			else {
-				/* Find difference. Offset by x3 to get the 
-				 * expected range of values.
-				 */
-				int s1 = ix3 + (v1 - v2);
-				int s2;
 
-				if( s1 < 0 )
+			int out;
-					s1 = 0;
-				else if( s1 > range )
-					s1 = range;
 
-				s2 = v1 + lut[s1];
+			out = v1 + lut[d + 32768];
 
-				if( s2 < 0 ) 
+			if( out < 0 )
-					s2 = 0;
+				out = 0;
-				else if( s2 > 32767 ) 
+			if( out > 32767 )
-					s2 = 32767;
+				out = 32767;
 
-				q[x] = s2;
+			q[x] = out;
-			}
 		}
 	}
 
@@ -195,7 +171,6 @@ vips_sharpen_build( VipsObject *object )
 	VipsImage **args = (VipsImage **) vips_object_local_array( object, 2 );
 
 	VipsImage *in;
-	int ix1, ix2, ix3;
 	int i;
 
 	VIPS_GATE_START( "vips_sharpen_build: build" ); 
@@ -203,12 +178,12 @@ vips_sharpen_build( VipsObject *object )
 	if( VIPS_OBJECT_CLASS( vips_sharpen_parent_class )->build( object ) )
 		return( -1 );
 
-	/* Turn y parameters into xs.
+	/* We used to have a radius control. If that's set but sigma isn't,
+	 * use it to set a reasonable value for sigma.
 	 */
-	sharpen->x2 = (sharpen->y2 - 
+	if( !vips_object_argument_isset( object, "sigma" )  &&
-		sharpen->x1 * (sharpen->m1 - sharpen->m2)) / sharpen->m2;
+		vips_object_argument_isset( object, "radius" ) )
-	sharpen->x3 = (sharpen->y3 - 
+		sharpen->sigma = 1 + sharpen->radius / 2;
-		sharpen->x1 * (sharpen->m1 - sharpen->m2)) / sharpen->m2;
 
 	in = sharpen->in; 
 
@@ -221,20 +196,10 @@ vips_sharpen_build( VipsObject *object )
 		vips_check_format( class->nickname, in, VIPS_FORMAT_SHORT ) )
   		return( -1 );
 
-	if( sharpen->x1 < 0 || sharpen->x1 > 99 || 
+	/* Stop at 20% of max ... a bit mean. We always sharpen a short, 
-		sharpen->x2 < 0 || sharpen->x2 > 99 || 
+	 * so there's no point using a float mask. 
-		sharpen->x1 > sharpen->x2 ||
-		sharpen->x3 < 0 || sharpen->x3 > 99 || 
-		sharpen->x1 > sharpen->x3 ) {
-		vips_error( class->nickname, 
-			"%s", _( "parameters out of range" ) );
-		return( -1 );
-	}
-
-	/* Stop at 20% of max ... bit mean, but means mask radius is roughly
-	 * right. We always sharpen a short, so no point using a float mask. 
 	 */
-	if( vips_gaussmat( &t[1], 1 + sharpen->radius / 2, 0.2, 
+	if( vips_gaussmat( &t[1], sharpen->sigma, 0.2, 
 		"separable", TRUE,
 		"precision", VIPS_PRECISION_INTEGER,
 		NULL ) )
@@ -245,25 +210,51 @@ vips_sharpen_build( VipsObject *object )
 	vips_matrixprint( t[1], NULL ); 
 #endif /*DEBUG*/
 
-	/* Build the int lut.
+	/* Index with the signed difference between two 0 - 32767 images.
 	 */
-	sharpen->ix1 = ix1 = sharpen->x1 * 327.67;
+	if( !(sharpen->lut = VIPS_ARRAY( object, 65536, int )) )
-	sharpen->ix2 = ix2 = sharpen->x2 * 327.67;
-	sharpen->ix3 = ix3 = sharpen->x3 * 327.67;
-
-	if( !(sharpen->lut = VIPS_ARRAY( object, ix2 + ix3 + 1, int )) )
 		return( -1 );
 
-	for( i = 0; i < ix1; i++ ) {
+	for( i = 0; i < 65536; i++ ) { 
-		sharpen->lut[ix3 + i] = i * sharpen->m1;
+		/* Rescale to +/- 100.
-		sharpen->lut[ix3 - i] = -i * sharpen->m1;
+		 */
+		double v = (i - 32767) / 327.67;
+		double y;
+
+		if( v < -sharpen->x1 ) 
+			/* Left of -x1.
+			 */
+			y = (v - sharpen->x1) * sharpen->m2 +
+				-sharpen->x1 * sharpen->m1;
+		else if( v < sharpen->x1 ) 
+			/* Centre section.
+			 */
+			y = v * sharpen->m1;
+		else 
+			/* Right of x1.
+			 */
+			y = (v - sharpen->x1) * sharpen->m2 +
+				sharpen->x1 * sharpen->m1;
+
+		if( y < -sharpen->y3 )
+			y = -sharpen->y3;
+		if( y > sharpen->y2 )
+			y = sharpen->y2;
+
+		sharpen->lut[i] = VIPS_RINT( y * 327.67 );
 	}
-	for( i = ix1; i <= ix2; i++ ) 
+
-		sharpen->lut[ix3 + i] = 
+#ifdef DEBUG
-			ix1 * sharpen->m1 + (i - ix1) * sharpen->m2; 
+{
-	for( i = ix1; i <= ix3; i++ )
+	VipsImage *mat;
-		sharpen->lut[ix3 - i] = 
+
-			-(ix1 * sharpen->m1 + (i - ix1) * sharpen->m2);
+	mat = vips_image_new_matrix( 65536, 1 );
+	for( i = 0; i < 65536; i++ )
+		*VIPS_MATRIX( mat, i, 0 ) = sharpen->lut[i];
+	vips_image_write_to_file( mat, "x.v", NULL ); 
+	g_object_unref( mat );
+}
+#endif /*DEBUG*/
 
 	/* Extract L and the rest, convolve L.
 	 */
@@ -326,33 +317,33 @@ vips_sharpen_class_init( VipsSharpenClass *class )
 		VIPS_ARGUMENT_REQUIRED_OUTPUT, 
 		G_STRUCT_OFFSET( VipsSharpen, out ) );
 
-	VIPS_ARG_INT( class, "radius", 3, 
+	VIPS_ARG_DOUBLE( class, "sigma", 3, 
-		_( "Radius" ), 
+		_( "Sigma" ), 
-		_( "Mask radius" ),
+		_( "Sigma of Gaussian" ),
 		VIPS_ARGUMENT_OPTIONAL_INPUT,
-		G_STRUCT_OFFSET( VipsSharpen, radius ),
+		G_STRUCT_OFFSET( VipsSharpen, sigma ),
-		1, 1000000, 3 );
+		0.000001, 10000.0, 1.0 );
 
 	VIPS_ARG_DOUBLE( class, "x1", 5, 
 		_( "x1" ), 
 		_( "Flat/jaggy threshold" ),
 		VIPS_ARGUMENT_OPTIONAL_INPUT,
 		G_STRUCT_OFFSET( VipsSharpen, x1 ),
-		1, 1000000, 1.5 );
+		0, 1000000, 1.5 );
 
 	VIPS_ARG_DOUBLE( class, "y2", 6, 
 		_( "y2" ), 
 		_( "Maximum brightening" ),
 		VIPS_ARGUMENT_OPTIONAL_INPUT,
 		G_STRUCT_OFFSET( VipsSharpen, y2 ),
-		1, 1000000, 20 );
+		0, 1000000, 20 );
 
 	VIPS_ARG_DOUBLE( class, "y3", 7, 
 		_( "y3" ), 
 		_( "Maximum darkening" ),
 		VIPS_ARGUMENT_OPTIONAL_INPUT,
 		G_STRUCT_OFFSET( VipsSharpen, y3 ),
-		1, 1000000, 50 );
+		0, 1000000, 50 );
 
 	VIPS_ARG_DOUBLE( class, "m1", 8, 
 		_( "m1" ), 
@@ -368,12 +359,21 @@ vips_sharpen_class_init( VipsSharpenClass *class )
 		G_STRUCT_OFFSET( VipsSharpen, m2 ),
 		0, 1000000, 2 );
 
+	/* We used to have a radius control.
+	 */
+	VIPS_ARG_INT( class, "radius", 3, 
+		_( "Radius" ), 
+		_( "radius of Gaussian" ),
+		VIPS_ARGUMENT_OPTIONAL_INPUT | VIPS_ARGUMENT_DEPRECATED,
+		G_STRUCT_OFFSET( VipsSharpen, radius ),
+		1, 100, 1 );
+
 }
 
 static void
 vips_sharpen_init( VipsSharpen *sharpen )
 {
-	sharpen->radius = 3; 
+	sharpen->sigma = 1.0; 
 	sharpen->x1 = 1.5; 
 	sharpen->y2 = 20; 
 	sharpen->y3 = 50; 
@@ -389,7 +389,7 @@ vips_sharpen_init( VipsSharpen *sharpen )
  *
  * Optional arguments:
  *
- * @radius: how large a mask to use
+ * @sigma: sigma of gaussian
  * @x1: flat/jaggy threshold
  * @y2: maximum amount of brightening
  * @y3: maximum amount of darkening
@@ -429,7 +429,7 @@ vips_sharpen_init( VipsSharpen *sharpen )
  * For printing, we recommend the following settings (the defaults):
  *
  * |[
- *   radius == 3
+ *   sigma == 1.0
  *   x1 == 1.5
  *   y2 == 20         (don't brighten by more than 20 L*)
  *   y3 == 50         (can darken by up to 50 L*)
@@ -441,10 +441,10 @@ vips_sharpen_init( VipsSharpen *sharpen )
  * If you want more or less sharpening, we suggest you just change the m1 
  * and m2 parameters. 
  *
- * The @radius parameter changes the width of the fringe and can be 
+ * The @sigma parameter changes the width of the fringe and can be 
  * adjusted according to the output printing resolution. As an approximate 
- * guideline, use 1 for 4 pixels/mm (CRT display resolution), 2 for 8 
+ * guideline, use 0.5 for 4 pixels/mm (CRT display resolution), 
- * pixels/mm, 3 for 12 pixels/mm and 4 for 16 pixels/mm (300 dpi == 12 
+ * 1.0 for 12 pixels/mm and 1.5 for 16 pixels/mm (300 dpi == 12 
  * pixels/mm). These figures refer to the image raster, not the half-tone 
  * resolution.
  *

libvips/deprecated/vips7compat.c

@@ -2518,7 +2518,7 @@ im_sharpen( IMAGE *in, IMAGE *out,
 	 * upon this behaviour. 
 	 */
 	if( vips_call( "sharpen", in, &t[0], 
-		"radius", mask_size / 2,
+		"sigma", mask_size / 4.0,
 		"x1", x1,
 		"y2", y2,
 		"y3", y3,

test/test_convolution.py

@@ -218,15 +218,28 @@ class TestConvolution(unittest.TestCase):
     def test_sharpen(self):
         for im in self.all_images:
             for fmt in noncomplex_formats:
-                for radius in range(1, 7):
+                # old vipses used "radius", check that that still works
+                sharp = im.sharpen(radius = 5)
+
+                for sigma in [0.5, 1, 1.5, 2]:
                     im = im.cast(fmt)
                     if im.bands == 3:
                         im = im.copy(interpretation = Vips.Interpretation.SRGB)
-                    sharp = im.sharpen(radius = radius)
+                    elif im.bands == 1:
+                        im = im.copy(interpretation = Vips.Interpretation.B_W)
+                    sharp = im.sharpen(sigma = sigma)
 
                     # hard to test much more than this
                     self.assertEqual(im.width, sharp.width)
                     self.assertEqual(im.height, sharp.height)
 
+                    # if m1 and m2 are zero, sharpen should do nothing
+                    sharp = im.sharpen(sigma = sigma, m1 = 0, m2 = 0)
+                    sharp = sharp.colourspace(im.interpretation)
+                    #print("testing sig = %g" % sigma)
+                    #print("testing fmt = %s" % fmt)
+                    #print("max diff = %g" % (im - sharp).abs().max())
+                    self.assertEqual((im - sharp).abs().max(), 0)
+
 if __name__ == '__main__':
     unittest.main()