rebase add-jp2ktiff on master

ChangeLog

@@ -21,7 +21,8 @@
 - png and gif load note background colour as metadata [781545872]
 - add vips_image_[set|get]_array_double()
 - add GIF load with libnsgif
-- add JPEG2000 load and save
+- add jp2kload, jp2ksave
+- add jp2k compression to tiff load and save
 - add JPEG-XL load and save
 - add black_point_compensation flag for icc transforms
 - add "rgba" flag to vips_text() to enable full colour text rendering

libvips/foreign/jp2kload.c

@@ -274,10 +274,6 @@ vips_foreign_load_jp2k_error_callback( const char *msg, void *client )
 
 	vips_error( class->nickname, "%s", msg ); 
 	jp2k->n_errors += 1;
-
-#ifdef DEBUG
-	printf( "%s: error %s",  class->nickname, msg );
-#endif /*DEBUG*/
 }
 
 /* The openjpeg info and warning callbacks are incredibly chatty.
@@ -586,7 +582,7 @@ vips_foreign_load_jp2k_header( VipsForeignLoad *load )
 	\
 	for( x = 0; x < length; x++ ) { \
 		for( i = 0; i < b; i++ ) { \
-			int dx = jp2k->image->comps[i].dx; \
+			int dx = image->comps[i].dx; \
 			int pixel = planes[i][x / dx]; \
 			\
 			tq[i] = pixel; \
@@ -596,28 +592,31 @@ vips_foreign_load_jp2k_header( VipsForeignLoad *load )
 	} \
 }
 
-/* Pack the set of openjpeg components into a libvips region. left/top are the
+/* Pack a line of openjpeg pixels into libvips format. left/top are the
- * offsets into the tile in pixel coordinates where we should start reading.
+ * offsets into the opj image in pixel coordinates where we should start 
+ * reading.
+ *
+ * Set upsample if any opj component is subsampled.
  */
 static void
-vips_foreign_load_jp2k_pack( VipsForeignLoadJp2k *jp2k, 
+vips_foreign_load_jp2k_pack( gboolean upsample, 
-	VipsImage *image, VipsPel *q, 
+	opj_image_t *image, VipsImage *im, 
-	int left, int top, int length )
+	VipsPel *q, int left, int top, int length )
 {
 	int *planes[MAX_BANDS];
-	int b = jp2k->image->numcomps;
+	int b = image->numcomps;
 
 	int x, i;
 
 	for( i = 0; i < b; i++ ) {
-		opj_image_comp_t *comp = &jp2k->image->comps[i];
+		opj_image_comp_t *comp = &image->comps[i];
 
 		planes[i] = comp->data + (top / comp->dy) * comp->w + 
 			(left / comp->dx);
 	}
 
-	if( jp2k->upsample ) 
+	if( upsample ) 
-		switch( image->BandFmt ) {
+		switch( im->BandFmt ) {
 		case VIPS_FORMAT_CHAR:
 		case VIPS_FORMAT_UCHAR:
 			PACK_UPSAMPLE( unsigned char );
@@ -640,7 +639,7 @@ vips_foreign_load_jp2k_pack( VipsForeignLoadJp2k *jp2k,
 	else 
 		/* Fast no-upsample path.
 		 */
-		switch( image->BandFmt ) {
+		switch( im->BandFmt ) {
 		case VIPS_FORMAT_CHAR:
 		case VIPS_FORMAT_UCHAR:
 			PACK( unsigned char );
@@ -692,17 +691,16 @@ vips_foreign_load_jp2k_pack( VipsForeignLoadJp2k *jp2k,
 /* YCC->RGB for a line of pels.
  */
 static void
-vips_foreign_load_jp2k_ycc_to_rgb( VipsForeignLoadJp2k *jp2k, 
+vips_foreign_load_jp2k_ycc_to_rgb( opj_image_t *image, VipsImage *im, 
 	VipsPel *q, int length )
 {
-	VipsForeignLoad *load = (VipsForeignLoad *) jp2k;
+	int prec = image->comps[0].prec;
-	int prec = jp2k->image->comps[0].prec;
 	int offset = 1 << (prec - 1);
 	int upb = (1 << prec) - 1;
 
 	int x;
 
-	switch( load->out->BandFmt ) {
+	switch( im->BandFmt ) {
 	case VIPS_FORMAT_CHAR:
 	case VIPS_FORMAT_UCHAR:
 		YCC_TO_RGB( unsigned char );
@@ -806,15 +804,16 @@ vips_foreign_load_jp2k_generate( VipsRegion *out,
 				VipsPel *q = VIPS_REGION_ADDR( out, 
 					hit.left, hit.top + z );
 
-				vips_foreign_load_jp2k_pack( jp2k,
+				vips_foreign_load_jp2k_pack( jp2k->upsample, 
-					out->im, q,
+					jp2k->image, out->im, q,
 					hit.left - tile.left,
 					hit.top - tile.top + z,
 					hit.width ); 
 
 				if( jp2k->ycc_to_rgb )
-					vips_foreign_load_jp2k_ycc_to_rgb( jp2k,
+					vips_foreign_load_jp2k_ycc_to_rgb( 
-						q, hit.width );
+						jp2k->image, out->im, q, 
+						hit.width );
 			}
 
 			x += hit.width;
@@ -1132,6 +1131,135 @@ vips_foreign_load_jp2k_source_init(
 {
 }
 
+static void 
+warning_callback( const char *msg G_GNUC_UNUSED, void *data G_GNUC_UNUSED ) 
+{
+	/* There are a lot of warnings ...
+	 */
+}
+
+static void 
+error_callback( const char *msg, void *data ) 
+{
+	printf( "OpenJPEG: %s", msg ); 
+	vips_error( "OpenJPEG", "%s", msg ); 
+}
+
+typedef struct _TileDecompress {
+	VipsSource *source;
+        opj_stream_t *stream;
+        opj_codec_t *codec;
+	opj_image_t *image;
+} TileDecompress;
+
+void
+vips__foreign_load_jp2k_decompress_free( TileDecompress *decompress )
+{
+	VIPS_FREEF( opj_destroy_codec, decompress->codec );
+	VIPS_FREEF( opj_image_destroy, decompress->image );
+	VIPS_FREEF( opj_stream_destroy, decompress->stream );
+	VIPS_UNREF( decompress->source );
+}
+
+/* Called from tiff2vips to decode a jp2k-compressed tile. 
+ *
+ * width/height is the tile size. If this is an edge tile, and smaller than 
+ * this, we still write a full-size tile and our caller will clip.
+ */
+int
+vips__foreign_load_jp2k_decompress( VipsImage *out, 
+	int width, int height, gboolean ycc_to_rgb, 
+	void *from, size_t from_length, 
+	void *to, size_t to_length )
+{
+	size_t pel_size = VIPS_IMAGE_SIZEOF_PEL( out );
+	size_t line_size = pel_size * width;
+
+	TileDecompress decompress = { 0 };
+	opj_dparameters_t parameters;
+	int i;
+	gboolean upsample;
+	VipsPel *q;
+	int y;
+
+	/* Our ycc->rgb only works for exactly 3 bands.
+	 */
+	ycc_to_rgb = ycc_to_rgb && out->Bands == 3;
+
+	decompress.codec = opj_create_decompress( OPJ_CODEC_J2K );
+	opj_set_default_decoder_parameters( &parameters );
+	opj_setup_decoder( decompress.codec, &parameters );
+	opj_set_warning_handler( decompress.codec, warning_callback, NULL );
+	opj_set_error_handler( decompress.codec, error_callback, NULL );
+
+	decompress.source = vips_source_new_from_memory( from, from_length );
+	decompress.stream = vips_foreign_load_jp2k_stream( decompress.source );
+	if( !opj_read_header( decompress.stream, 
+		decompress.codec, &decompress.image ) ) {
+		vips_error( "jp2kload", "%s", ( "header error" ) );
+		vips__foreign_load_jp2k_decompress_free( &decompress ); 
+		return( -1 );
+	}
+
+	if( decompress.image->x1 > width || 
+		decompress.image->y1 > height ||
+		line_size * height > to_length ) {
+		vips_error( "jp2kload", "%s", ( "bad dimensions" ) );
+		vips__foreign_load_jp2k_decompress_free( &decompress ); 
+    		return( -1 );
+	}
+
+	if( !opj_decode( decompress.codec, 
+		decompress.stream, decompress.image ) ) {
+		vips_error( "jp2kload", "%s", ( "decode error" ) );
+		vips__foreign_load_jp2k_decompress_free( &decompress ); 
+		return( -1 );
+	}
+
+	/* Do any components need upsampling?
+	 */
+	upsample = FALSE;
+	for( i = 0; i < decompress.image->numcomps; i++ ) {
+		opj_image_comp_t *this = &decompress.image->comps[i];
+
+		if( this->dx > 1 ||
+			this->dy > 1 )
+			upsample = TRUE;
+	}
+
+	/* Unpack hit pixels to buffer in vips layout. 
+	 */
+	q = to;
+	for( y = 0; y < height; y++ ) {
+		vips_foreign_load_jp2k_pack( upsample, 
+			decompress.image, out, q,
+			0, y, width ); 
+
+		if( ycc_to_rgb )
+			vips_foreign_load_jp2k_ycc_to_rgb( 
+				decompress.image, out, q, 
+				width );
+
+		q += line_size;
+	}
+
+	vips__foreign_load_jp2k_decompress_free( &decompress ); 
+
+	return( 0 );
+}
+
+#else /*!HAVE_LIBOPENJP2*/
+
+int
+vips__foreign_load_jp2k_decompress_buffer( void *data, size_t length, 
+	int width, int height, 
+	void *dest, size_t dest_length )
+{
+	vips_error( "tiff2vips", 
+		"%s", _( "libvips built without JPEG2000 support" ) );
+	return( -1 );
+}
+
 #endif /*HAVE_LIBOPENJP2*/
 
 /**

libvips/foreign/jp2ksave.c

@@ -42,6 +42,8 @@
  *
  * - could support png-like bitdepth parameter
  *
+ * - could support cp_comment field? not very useful
+ *
  */
 
 #ifdef HAVE_CONFIG_H
@@ -93,7 +95,6 @@ typedef struct _VipsForeignSaveJp2k {
 	opj_stream_t *stream;
 	opj_codec_t *codec;
 	opj_cparameters_t parameters;
-	opj_image_cmptparm_t comps[MAX_BANDS];
 	opj_image_t *image;
 
 	/* The line of tiles we are building, and the buffer we
@@ -102,9 +103,9 @@ typedef struct _VipsForeignSaveJp2k {
 	VipsRegion *strip;
 	VipsPel *tile_buffer;
 
-	/* If we need to downsample during unpacking.
+	/* If we need to subsample during unpacking.
 	 */
-	gboolean downsample;
+	gboolean subsample;
 
 	/* If we converto RGB to YCC during save.
 	 */
@@ -172,9 +173,7 @@ vips_foreign_save_jp2k_target( VipsTarget *target )
 static void 
 vips_foreign_save_jp2k_error_callback( const char *msg, void *client )
 {
-	VipsObjectClass *class = VIPS_OBJECT_GET_CLASS( client );
+	vips_error( "jp2ksave", "%s", msg ); 
-
-	vips_error( class->nickname, "%s", msg ); 
 }
 
 /* The openjpeg info and warning callbacks are incredibly chatty.
@@ -183,9 +182,7 @@ static void
 vips_foreign_save_jp2k_warning_callback( const char *msg, void *client )
 {
 #ifdef DEBUG
-	VipsObjectClass *class = VIPS_OBJECT_GET_CLASS( client );
+	g_warning( "jp2ksave: %s",  class->nickname, msg );
-
-	g_warning( "%s: %s",  class->nickname, msg );
 #endif /*DEBUG*/
 }
 
@@ -195,22 +192,19 @@ static void
 vips_foreign_save_jp2k_info_callback( const char *msg, void *client )
 {
 #ifdef DEBUG
-	VipsObjectClass *class = VIPS_OBJECT_GET_CLASS( client );
+	g_info( "jp2ksave: %s",  class->nickname, msg );
-
-	g_info( "%s: %s",  class->nickname, msg );
 #endif /*DEBUG*/
 }
 
 static void
-vips_foreign_save_jp2k_attach_handlers( VipsForeignSaveJp2k *jp2k,
+vips_foreign_save_jp2k_attach_handlers( opj_codec_t *codec )
-	opj_codec_t *codec )
 {
-	opj_set_info_handler( codec, 
+	opj_set_info_handler( codec,
-		vips_foreign_save_jp2k_info_callback, jp2k );
+		vips_foreign_save_jp2k_info_callback, NULL );
 	opj_set_warning_handler( codec, 
-		vips_foreign_save_jp2k_warning_callback, jp2k );
+		vips_foreign_save_jp2k_warning_callback, NULL );
 	opj_set_error_handler( codec, 
-		vips_foreign_save_jp2k_error_callback, jp2k );
+		vips_foreign_save_jp2k_error_callback, NULL );
 }
 
 /* See also https://en.wikipedia.org/wiki/YCbCr#JPEG_conversion
@@ -238,23 +232,25 @@ vips_foreign_save_jp2k_attach_handlers( VipsForeignSaveJp2k *jp2k,
 	} \
 }
 
-/* RGB->YCC for a line of pels.
+/* In-place RGB->YCC for a line of pels.
  */
 static void
-vips_foreign_save_jp2k_rgb_to_ycc( VipsForeignSaveJp2k *jp2k, VipsRect *tile )
+vips_foreign_save_jp2k_rgb_to_ycc( VipsRegion *region, 
+	VipsRect *tile, int prec ) 
 {
-	VipsForeignSave *save = (VipsForeignSave *) jp2k;
+	VipsImage *im = region->im;
-	int prec = jp2k->image->comps[0].prec;
 	int offset = 1 << (prec - 1);
 	int upb = (1 << prec) - 1;
 
 	int x, y;
 
+	g_assert( im->Bands == 3 );
+
 	for( y = 0; y < tile->height; y++ ) {
-		VipsPel *q = VIPS_REGION_ADDR( jp2k->strip, 
+		VipsPel *q = VIPS_REGION_ADDR( region, 
 			tile->left, tile->top + y );
 
-		switch( save->ready->BandFmt ) {
+		switch( im->BandFmt ) {
 		case VIPS_FORMAT_CHAR:
 		case VIPS_FORMAT_UCHAR:
 			RGB_TO_YCC( unsigned char );
@@ -282,9 +278,9 @@ vips_foreign_save_jp2k_rgb_to_ycc( VipsForeignSaveJp2k *jp2k, VipsRect *tile )
  *   2. add subsequent lines in comp.dy.
  *   3. horizontal average to output line
  */
-#define SHRINK( ACC_TYPE, PIXEL_TYPE ) { \
+#define SHRINK( OUTPUT_TYPE, ACC_TYPE, PIXEL_TYPE ) { \
-	ACC_TYPE *acc = (ACC_TYPE *) jp2k->accumulate; \
+	ACC_TYPE *acc = (ACC_TYPE *) accumulate; \
-	PIXEL_TYPE *tq = (PIXEL_TYPE *) q; \
+	OUTPUT_TYPE *tq = (OUTPUT_TYPE *) q; \
 	const int n_pels = comp->dx * comp->dy; \
 	\
 	PIXEL_TYPE *tp; \
@@ -318,23 +314,23 @@ vips_foreign_save_jp2k_rgb_to_ycc( VipsForeignSaveJp2k *jp2k, VipsRect *tile )
 }
 
 static void
-vips_foreign_save_jp2k_unpack_downsample( VipsForeignSaveJp2k *jp2k, 
+vips_foreign_save_jp2k_unpack_subsample( VipsRegion *region, VipsRect *tile,
-	VipsRect *tile )
+	opj_image_t *image, VipsPel *tile_buffer, VipsPel *accumulate )
 {
-	VipsForeignSave *save = (VipsForeignSave *) jp2k;
+	VipsImage *im = region->im;
-	size_t sizeof_element = VIPS_IMAGE_SIZEOF_ELEMENT( save->ready );
+	size_t sizeof_element = VIPS_REGION_SIZEOF_ELEMENT( region );
-	size_t lskip = VIPS_REGION_LSKIP( jp2k->strip );
+	size_t lskip = VIPS_REGION_LSKIP( region );
-	int n_bands = save->ready->Bands;
+	int n_bands = im->Bands;
 
 	VipsPel *q;
-	int x, y, z, band_index;
+	int x, y, z, i;
 
-	q = jp2k->tile_buffer;
+	q = tile_buffer;
-	for( band_index = 0; band_index < n_bands; band_index++ ) {
+	for( i = 0; i < n_bands; i++ ) {
-		opj_image_comp_t *comp = &jp2k->image->comps[band_index];
+		opj_image_comp_t *comp = &image->comps[i];
 
-		/* The number of pixels we write for this component. Round to
+		/* The number of pixels we write for this component. No
-		 * nearest, and we may have to write half-pixels at the edges.
+		 * padding.
 		 */
 		int output_width = VIPS_ROUND_UINT( 
 			(double) tile->width / comp->dx );
@@ -342,35 +338,35 @@ vips_foreign_save_jp2k_unpack_downsample( VipsForeignSaveJp2k *jp2k,
 			(double) tile->height / comp->dy );;
 
 		for( y = 0; y < output_height; y++ ) {
-			VipsPel *p = band_index * sizeof_element + 
+			VipsPel *p = i * sizeof_element + 
-				VIPS_REGION_ADDR( jp2k->strip, 
+				VIPS_REGION_ADDR( region, 
 					tile->left, tile->top + y * comp->dy );
 
 			/* Shrink a line of pels to q.
 			 */
-			switch( save->ready->BandFmt ) {
+			switch( im->BandFmt ) {
 			case VIPS_FORMAT_CHAR:
-				SHRINK( int, signed char );
+				SHRINK( signed char, int, signed char );
 				break;
 
 			case VIPS_FORMAT_UCHAR:
-				SHRINK( int, unsigned char );
+				SHRINK( unsigned char, int, unsigned char );
 				break;
 
 			case VIPS_FORMAT_SHORT:
-				SHRINK( int, signed short );
+				SHRINK( signed short, int, signed short );
 				break;
 
 			case VIPS_FORMAT_USHORT:
-				SHRINK( int, unsigned short );
+				SHRINK( unsigned short, int, unsigned short );
 				break;
 
 			case VIPS_FORMAT_INT:
-				SHRINK( gint64, signed int );
+				SHRINK( signed int, gint64, signed int );
 				break;
 
 			case VIPS_FORMAT_UINT:
-				SHRINK( gint64, unsigned int );
+				SHRINK( unsigned int, gint64, unsigned int );
 				break;
 
 			default:
@@ -383,60 +379,56 @@ vips_foreign_save_jp2k_unpack_downsample( VipsForeignSaveJp2k *jp2k,
 	}
 }
 
-#define UNPACK( TYPE ) { \
+#define UNPACK( OUT, IN ) { \
-	TYPE **tplanes = (TYPE **) planes; \
+	OUT *tq = (OUT *) q; \
-	TYPE *tp = (TYPE *) p; \
+	IN *tp = (IN *) p + i; \
 	\
-	for( i = 0; i < b; i++ ) { \
+	for( x = 0; x < tile->width; x++ ) { \
-		TYPE *q = tplanes[i]; \
+		tq[x] = *tp; \
-		TYPE *tp1 = tp + i; \
+		tp += b; \
-		\
-		for( x = 0; x < tile->width; x++ ) { \
-			q[x] = *tp1; \
-			tp1 += b; \
-		} \
-		\
-		tplanes[i] += tile->width; \
 	} \
 }
 
 static void
-vips_foreign_save_jp2k_unpack( VipsForeignSaveJp2k *jp2k, VipsRect *tile )
+vips_foreign_save_jp2k_unpack( VipsRegion *region, VipsRect *tile,
+	opj_image_t *image, VipsPel *tile_buffer )
 {
-	VipsForeignSave *save = (VipsForeignSave *) jp2k;
+	VipsImage *im = region->im;
-	size_t sizeof_element = VIPS_IMAGE_SIZEOF_ELEMENT( save->ready );
+	size_t sizeof_element = VIPS_REGION_SIZEOF_ELEMENT( region );
-	int b = save->ready->Bands;
+	size_t sizeof_line = sizeof_element * tile->width;
+	size_t sizeof_tile = sizeof_line * tile->height;
+	int b = im->Bands;
 
-	VipsPel *planes[MAX_BANDS];
 	int x, y, i;
 
-	for( i = 0; i < b; i++ )
-		planes[i] = jp2k->tile_buffer +
-			i * sizeof_element * tile->width * tile->height;
-
 	for( y = 0; y < tile->height; y++ ) {
-		VipsPel *p = VIPS_REGION_ADDR( jp2k->strip, 
+		VipsPel *p = VIPS_REGION_ADDR( region, 
 			tile->left, tile->top + y );
 
-		switch( save->ready->BandFmt ) {
+		for( i = 0; i < b; i++ ) {
-		case VIPS_FORMAT_CHAR:
+			VipsPel *q = tile_buffer + 
-		case VIPS_FORMAT_UCHAR:
+				i * sizeof_tile + y * sizeof_line;
-			UNPACK( unsigned char );
-			break;
 
-		case VIPS_FORMAT_SHORT:
+			switch( im->BandFmt ) {
-		case VIPS_FORMAT_USHORT:
+			case VIPS_FORMAT_CHAR:
-			UNPACK( unsigned short );
+			case VIPS_FORMAT_UCHAR:
-			break;
+				UNPACK( unsigned char, unsigned char );
+				break;
 
-		case VIPS_FORMAT_INT:
+			case VIPS_FORMAT_SHORT:
-		case VIPS_FORMAT_UINT:
+			case VIPS_FORMAT_USHORT:
-			UNPACK( unsigned int );
+				UNPACK( unsigned short, unsigned short );
-			break;
+				break;
 
-		default:
+			case VIPS_FORMAT_INT:
-			g_assert_not_reached();
+			case VIPS_FORMAT_UINT:
-			break;
+				UNPACK( unsigned int, unsigned int );
+				break;
+
+			default:
+				g_assert_not_reached();
+				break;
+			}
 		}
 	}
 }
@@ -472,13 +464,13 @@ static int
 vips_foreign_save_jp2k_write_tiles( VipsForeignSaveJp2k *jp2k )
 {
 	VipsForeignSave *save = (VipsForeignSave *) jp2k;
-	int tiles_across = 
+	VipsImage *im = save->ready;
-		VIPS_ROUND_UP( save->ready->Xsize, jp2k->tile_width ) /
+	int tiles_across = VIPS_ROUND_UP( im->Xsize, jp2k->tile_width ) /
-			jp2k->tile_width;
+		jp2k->tile_width;
 
 	int x;
 
-	for( x = 0; x < save->ready->Xsize; x += jp2k->tile_width ) {
+	for( x = 0; x < im->Xsize; x += jp2k->tile_width ) {
 		VipsRect tile;
 		size_t sizeof_tile;
 		int tile_index;
@@ -490,14 +482,20 @@ vips_foreign_save_jp2k_write_tiles( VipsForeignSaveJp2k *jp2k )
 		vips_rect_intersectrect( &tile, &jp2k->strip->valid, &tile );
 
 		if( jp2k->save_as_ycc ) 
-			vips_foreign_save_jp2k_rgb_to_ycc( jp2k, &tile );
+			vips_foreign_save_jp2k_rgb_to_ycc( jp2k->strip, 
+				&tile, jp2k->image->comps[0].prec ); 
 
-		if( jp2k->downsample )
+		if( jp2k->subsample )
-			vips_foreign_save_jp2k_unpack_downsample( jp2k, &tile );
+			vips_foreign_save_jp2k_unpack_subsample( jp2k->strip, 
+				&tile, jp2k->image, 
+				jp2k->tile_buffer, jp2k->accumulate ); 
 		else
-			vips_foreign_save_jp2k_unpack( jp2k, &tile );
+			vips_foreign_save_jp2k_unpack( jp2k->strip, 
+				&tile, jp2k->image, 
+				jp2k->tile_buffer );
 
-		sizeof_tile = vips_foreign_save_jp2k_sizeof_tile( jp2k, &tile );
+		sizeof_tile = 
+			vips_foreign_save_jp2k_sizeof_tile( jp2k, &tile );
 		tile_index = tiles_across * tile.top / jp2k->tile_height +
 			x / jp2k->tile_width;
 		if( !opj_write_tile( jp2k->codec, tile_index, 
@@ -550,7 +548,7 @@ vips_foreign_save_jp2k_write_block( VipsRegion *region, VipsRect *area,
 			break;
 
 		/* We have reached the bottom of the strip. Write this line of
-		 * pixels and ove the strip down.
+		 * pixels and move the strip down.
 		 */
 		if( vips_foreign_save_jp2k_write_tiles( jp2k ) )
 			return( -1 );
@@ -566,67 +564,90 @@ vips_foreign_save_jp2k_write_block( VipsRegion *region, VipsRect *area,
 	return( 0 );
 }
 
-static int
+/* We can't call opj_calloc on win, sadly.
-vips_foreign_save_jp2k_build( VipsObject *object )
+ */
-{
+#define VIPS_OPJ_CALLOC( N, TYPE ) \
-	VipsObjectClass *class = VIPS_OBJECT_GET_CLASS( object );
+	((TYPE *) calloc( (N), sizeof( TYPE ) ))
-	VipsForeignSave *save = (VipsForeignSave *) object;
-	VipsForeignSaveJp2k *jp2k = (VipsForeignSaveJp2k *) object;
 
+/* Allocate an openjpeg image structure. Openjpeg has opj_image_create(), but
+ * that always allocates memory for each channel, and we don't want that when
+ * we are doing tiled write.
+ */
+static opj_image_t * 
+vips_opj_image_create( OPJ_UINT32 numcmpts,
+	opj_image_cmptparm_t *cmptparms, OPJ_COLOR_SPACE clrspc,
+	gboolean allocate )
+{
+	OPJ_UINT32 compno;
+	opj_image_t *image = NULL;
+
+	if( !(image = VIPS_OPJ_CALLOC( 1, opj_image_t )) )
+		return( NULL );
+
+        image->color_space = clrspc;
+        image->numcomps = numcmpts;
+        image->comps = VIPS_OPJ_CALLOC( image->numcomps, opj_image_comp_t );
+        if( !image->comps ) {
+            opj_image_destroy( image );
+            return( NULL );
+        }
+
+        for( compno = 0; compno < numcmpts; compno++ ) {
+		opj_image_comp_t *comp = &image->comps[compno];
+
+		comp->dx = cmptparms[compno].dx;
+		comp->dy = cmptparms[compno].dy;
+		comp->w = cmptparms[compno].w;
+		comp->h = cmptparms[compno].h;
+		comp->x0 = cmptparms[compno].x0;
+		comp->y0 = cmptparms[compno].y0;
+		comp->prec = cmptparms[compno].prec;
+		comp->bpp = cmptparms[compno].bpp;
+		comp->sgnd = cmptparms[compno].sgnd;
+
+		if( comp->h != 0 &&
+			(OPJ_SIZE_T) comp->w > SIZE_MAX / comp->h / 
+				sizeof( OPJ_INT32 ) ) {
+			opj_image_destroy( image );
+			return( NULL );
+		}
+
+		/* Allocation is optional.
+		 */
+		if( allocate ) {
+			size_t bytes = (size_t) comp->w * comp->h * 
+                                sizeof( OPJ_INT32 );
+
+			comp->data = (OPJ_INT32*) opj_image_data_alloc( bytes );
+			if( !comp->data ) {
+				opj_image_destroy( image );
+				return NULL;
+			}
+			memset( comp->data, 0, bytes );
+		}
+	}
+
+	return( image );
+}
+
+static opj_image_t *
+vips_foreign_save_jp2k_new_image( VipsImage *im, 
+	int width, int height, 
+	gboolean subsample, gboolean save_as_ycc, gboolean allocate )
+{
 	OPJ_COLOR_SPACE color_space;
 	int expected_bands;
 	int bits_per_pixel;
+	opj_image_cmptparm_t comps[MAX_BANDS];
+	opj_image_t *image;
 	int i;
-	size_t sizeof_tile;
-	size_t sizeof_line;
-	VipsRect strip_position;
 
-	if( VIPS_OBJECT_CLASS( vips_foreign_save_jp2k_parent_class )->
+	if( im->Bands > MAX_BANDS )
-		build( object ) )
+		return( NULL );
-		return( -1 );
-
-	opj_set_default_encoder_parameters( &jp2k->parameters );
-
-	/* Analyze our arguments.
-	 */
-
-	if( !vips_band_format_isint( save->ready->BandFmt ) ) {
-		vips_error( class->nickname,
-			"%s", _( "not an integer format" ) );
-		return( -1 );
-	}
-
-	switch( jp2k->subsample_mode ) {
-	case VIPS_FOREIGN_SUBSAMPLE_AUTO:
-		jp2k->downsample =
-			!jp2k->lossless &&
-			jp2k->Q < 90 &&
-			save->ready->Xsize % 2 == 0 &&
-			save->ready->Ysize % 2 == 0 &&
-			(save->ready->Type == VIPS_INTERPRETATION_sRGB ||
-			 save->ready->Type == VIPS_INTERPRETATION_RGB16) &&
-			save->ready->Bands == 3;
-		break;
-
-	case VIPS_FOREIGN_SUBSAMPLE_ON:
-		jp2k->downsample = TRUE;
-		break;
-
-	case VIPS_FOREIGN_SUBSAMPLE_OFF:
-		jp2k->downsample = FALSE;
-		break;
-
-	default:
-		g_assert_not_reached();
-		break;
-	}
-
-	if( jp2k->downsample ) 
-		jp2k->save_as_ycc = TRUE;
 
 	/* CIELAB etc. do not seem to be well documented.
 	 */
-	switch( save->ready->Type ) {
+	switch( im->Type ) {
 	case VIPS_INTERPRETATION_B_W:
 	case VIPS_INTERPRETATION_GREY16:
 		color_space = OPJ_CLRSPC_GRAY;
@@ -635,8 +656,7 @@ vips_foreign_save_jp2k_build( VipsObject *object )
 
 	case VIPS_INTERPRETATION_sRGB:
 	case VIPS_INTERPRETATION_RGB16:
-		color_space = jp2k->save_as_ycc ? 
+		color_space = save_as_ycc ? OPJ_CLRSPC_SYCC : OPJ_CLRSPC_SRGB;
-			OPJ_CLRSPC_SYCC : OPJ_CLRSPC_SRGB;
 		expected_bands = 3;
 		break;
 
@@ -647,11 +667,11 @@ vips_foreign_save_jp2k_build( VipsObject *object )
 
 	default:
 		color_space = OPJ_CLRSPC_UNSPECIFIED;
-		expected_bands = save->ready->Bands;
+		expected_bands = im->Bands;
 		break;
 	}
 
-	switch( save->ready->BandFmt ) {
+	switch( im->BandFmt ) {
 	case VIPS_FORMAT_CHAR:
 	case VIPS_FORMAT_UCHAR:
 		bits_per_pixel = 8;
@@ -674,8 +694,128 @@ vips_foreign_save_jp2k_build( VipsObject *object )
 		break;
 	}
 
+	for( i = 0; i < im->Bands; i++ ) {
+		comps[i].dx = (subsample && i > 0) ? 2 : 1;
+		comps[i].dy = (subsample && i > 0) ? 2 : 1;
+		comps[i].w = width;
+		comps[i].h = height;
+		comps[i].x0 = 0;
+		comps[i].y0 = 0;
+		comps[i].prec = bits_per_pixel;
+		comps[i].bpp = bits_per_pixel;
+		comps[i].sgnd = !vips_band_format_isuint( im->BandFmt );
+	}
+
+	image = vips_opj_image_create( im->Bands, comps, color_space, 
+		allocate );
+	image->x1 = width;
+	image->y1 = height;
+
+	/* Tag alpha channels.
+	 */
+	for( i = 0; i < im->Bands; i++ )
+		image->comps[i].alpha = i >= expected_bands;
+
+	return( image );
+}
+
+/* Compression profile derived from the BM's recommenadations, see:
+ *
+ * https://purl.pt/24107/1/iPres2013_PDF/An%20Analysis%20of%20Contemporary%20JPEG2000%20Codecs%20for%20Image%20Format%20Migration.pdf
+ */
+static void
+vips_foreign_save_jp2k_set_profile( opj_cparameters_t *parameters, 
+	gboolean lossless, int Q )
+{
+	if( lossless )
+		parameters->irreversible = FALSE;
+	else {
+		int i;
+
+		/* Equivalent command-line flags:
+		 *
+		 *   -I -p RPCL -n 7 \
+		 *   	-c[256,256],[256,256],[256,256],[256,256],[256,256],[256,256],[256,256] \
+		 *   	-b 64,64
+		 */
+
+		parameters->irreversible = TRUE;
+		parameters->prog_order = OPJ_RPCL;
+		parameters->numresolution = 7;
+		parameters->cblockw_init = 64;
+		parameters->cblockh_init = 64;
+		parameters->cp_disto_alloc = 1;
+		parameters->cp_fixed_quality = TRUE;
+		parameters->tcp_numlayers = 1;
+
+		/* No idea what this does, but opj_compress sets it.
+		 */
+		parameters->csty = 1;
+
+		parameters->res_spec = 7;
+		for( i = 0; i < parameters->res_spec; i++ ) { 
+			parameters->prch_init[i] = 256;
+			parameters->prcw_init[i] = 256;
+			parameters->tcp_distoratio[i] = Q + 10 * i;
+		}
+	}
+}
+
+static int
+vips_foreign_save_jp2k_build( VipsObject *object )
+{
+	VipsObjectClass *class = VIPS_OBJECT_GET_CLASS( object );
+	VipsForeignSave *save = (VipsForeignSave *) object;
+	VipsForeignSaveJp2k *jp2k = (VipsForeignSaveJp2k *) object;
+
+	size_t sizeof_tile;
+	size_t sizeof_line;
+	VipsRect strip_position;
+
+	if( VIPS_OBJECT_CLASS( vips_foreign_save_jp2k_parent_class )->
+		build( object ) )
+		return( -1 );
+
+	/* Analyze our arguments.
+	 */
+
+	if( !vips_band_format_isint( save->ready->BandFmt ) ) {
+		vips_error( class->nickname,
+			"%s", _( "not an integer format" ) );
+		return( -1 );
+	}
+
+	switch( jp2k->subsample_mode ) {
+	case VIPS_FOREIGN_SUBSAMPLE_AUTO:
+		jp2k->subsample =
+			!jp2k->lossless &&
+			jp2k->Q < 90 &&
+			save->ready->Xsize % 2 == 0 &&
+			save->ready->Ysize % 2 == 0 &&
+			(save->ready->Type == VIPS_INTERPRETATION_sRGB ||
+			 save->ready->Type == VIPS_INTERPRETATION_RGB16) &&
+			save->ready->Bands == 3;
+		break;
+
+	case VIPS_FOREIGN_SUBSAMPLE_ON:
+		jp2k->subsample = TRUE;
+		break;
+
+	case VIPS_FOREIGN_SUBSAMPLE_OFF:
+		jp2k->subsample = FALSE;
+		break;
+
+	default:
+		g_assert_not_reached();
+		break;
+	}
+
+	if( jp2k->subsample ) 
+		jp2k->save_as_ycc = TRUE;
+
 	/* Set parameters for compressor.
 	 */ 
+	opj_set_default_encoder_parameters( &jp2k->parameters );
 
 	/* Always tile.
 	 */
@@ -694,55 +834,27 @@ vips_foreign_save_jp2k_build( VipsObject *object )
 		jp2k->parameters.numresolution );
 #endif /*DEBUG*/
 
-	for( i = 0; i < save->ready->Bands; i++ ) {
-		jp2k->comps[i].dx = (jp2k->downsample && i > 0) ? 2 : 1;
-		jp2k->comps[i].dy = (jp2k->downsample && i > 0) ? 2 : 1;
-		jp2k->comps[i].w = save->ready->Xsize;
-		jp2k->comps[i].h = save->ready->Ysize;
-		jp2k->comps[i].x0 = 0;
-		jp2k->comps[i].y0 = 0;
-		jp2k->comps[i].prec = bits_per_pixel;
-		jp2k->comps[i].bpp = bits_per_pixel;
-		jp2k->comps[i].sgnd = 
-			!vips_band_format_isuint( save->ready->BandFmt );
-	}
-
 	/* Makes three band images smaller, somehow.
 	 */
-	jp2k->parameters.tcp_mct = 
+	jp2k->parameters.tcp_mct = save->ready->Bands >= 3 ? 1 : 0;
-		(save->ready->Bands == 3 && !jp2k->downsample) ? 1 : 0;
 
-	/* Lossy mode.
+	/* Set compression profile.
 	 */
-	if( !jp2k->lossless ) {
+	vips_foreign_save_jp2k_set_profile( &jp2k->parameters, 
-		jp2k->parameters.irreversible = TRUE;
+		jp2k->lossless, jp2k->Q ); 
-
-		/* Map Q to allowed distortion.
-		 */
-		jp2k->parameters.cp_disto_alloc = 1;
-		jp2k->parameters.cp_fixed_quality = TRUE;
-		jp2k->parameters.tcp_distoratio[0] = jp2k->Q;
-		jp2k->parameters.tcp_numlayers = 1;
-	}
-
-	/* Create output image.
-	 */
-
-	jp2k->image = opj_image_create( save->ready->Bands, 
-		jp2k->comps, color_space );
-	jp2k->image->x1 = save->ready->Xsize;
-	jp2k->image->y1 = save->ready->Ysize;
-
-	/* Tag alpha channels.
-	 */
-	for( i = 0; i < save->ready->Bands; i++ )
-		jp2k->image->comps[i].alpha = i >= expected_bands;
 
 	/* Set up compressor.
 	 */
 
 	jp2k->codec = opj_create_compress( OPJ_CODEC_J2K );
-	vips_foreign_save_jp2k_attach_handlers( jp2k, jp2k->codec );
+	vips_foreign_save_jp2k_attach_handlers( jp2k->codec );
+	/* FALSE means don't alloc memory for image planes (we write in 
+	 * tiles, not whole images).
+	 */
+	if( !(jp2k->image = vips_foreign_save_jp2k_new_image( save->ready,
+		save->ready->Xsize, save->ready->Ysize, 
+		jp2k->subsample, jp2k->save_as_ycc, FALSE )) )
+		return( -1 );
         if( !opj_setup_encoder( jp2k->codec, &jp2k->parameters, jp2k->image ) ) 
 		return( -1 );
 
@@ -756,7 +868,7 @@ vips_foreign_save_jp2k_build( VipsObject *object )
 	if( !(jp2k->stream = vips_foreign_save_jp2k_target( jp2k->target )) )
 		return( -1 );
 
-	if( !opj_start_compress( jp2k->codec, jp2k->image,  jp2k->stream ) )
+	if( !opj_start_compress( jp2k->codec, jp2k->image, jp2k->stream ) )
 		return( -1 );
 
 	/* The buffer we repack tiles to for write. Large enough for one
@@ -813,7 +925,7 @@ vips_foreign_save_jp2k_class_init( VipsForeignSaveJp2kClass *class )
 	gobject_class->get_property = vips_object_get_property;
 
 	object_class->nickname = "jp2ksave_base";
-	object_class->description = _( "save image in HEIF format" );
+	object_class->description = _( "save image in JPEG2000 format" );
 	object_class->build = vips_foreign_save_jp2k_build;
 
 	foreign_class->suffs = vips__jp2k_suffs;
@@ -854,7 +966,7 @@ vips_foreign_save_jp2k_class_init( VipsForeignSaveJp2kClass *class )
 		_( "Q factor" ),
 		VIPS_ARGUMENT_OPTIONAL_INPUT,
 		G_STRUCT_OFFSET( VipsForeignSaveJp2k, Q ),
-		1, 100, 45 );
+		1, 100, 48 );
 
 }
 
@@ -864,9 +976,9 @@ vips_foreign_save_jp2k_init( VipsForeignSaveJp2k *jp2k )
 	jp2k->tile_width = 512;
 	jp2k->tile_height = 512;
 
-	/* 45 gives about the same filesize as default regular jpg.
+	/* Chosen to give about the same filesize as regular jpg Q75.
 	 */
-	jp2k->Q = 45;
+	jp2k->Q = 48;
 
 	jp2k->subsample_mode = VIPS_FOREIGN_SUBSAMPLE_AUTO;
 }
@@ -1048,6 +1160,241 @@ vips_foreign_save_jp2k_target_init( VipsForeignSaveJp2kTarget *target )
 {
 }
 
+/* Stuff we track during tile compress.
+ */
+typedef struct _TileCompress {
+        opj_codec_t *codec;
+	opj_image_t *image;
+	opj_stream_t *stream;
+	VipsPel *accumulate;
+} TileCompress;
+
+/* Unpack from @tile within @region to the int data pointers on @image with
+ * subsampling.
+ */
+static void
+vips_foreign_save_jp2k_unpack_subsample_image( VipsRegion *region, 
+	VipsRect *tile, opj_image_t *image, VipsPel *accumulate )
+{
+	VipsImage *im = region->im;
+	size_t sizeof_element = VIPS_REGION_SIZEOF_ELEMENT( region );
+	size_t lskip = VIPS_REGION_LSKIP( region );
+	int n_bands = im->Bands;
+
+	int x, y, z, i;
+
+	for( i = 0; i < n_bands; i++ ) {
+		opj_image_comp_t *comp = &image->comps[i];
+		int *q = comp->data;
+
+		/* The number of pixels we write for this component. Lines
+		 * align to scanlines on comp.
+		 */
+		int output_width = VIPS_ROUND_UINT( 
+			(double) comp->w / comp->dx );
+		int output_height = VIPS_ROUND_UINT( 
+			(double) comp->h / comp->dy );
+
+		for( y = 0; y < output_height; y++ ) {
+			VipsPel *p = i * sizeof_element + 
+				VIPS_REGION_ADDR( region, 
+					tile->left, tile->top + y * comp->dy );
+
+			/* Shrink a line of pels to q.
+			 */
+			switch( im->BandFmt ) {
+			case VIPS_FORMAT_CHAR:
+				SHRINK( int, int, signed char );
+				break;
+
+			case VIPS_FORMAT_UCHAR:
+				SHRINK( int, int, unsigned char );
+				break;
+
+			case VIPS_FORMAT_SHORT:
+				SHRINK( int, int, signed short );
+				break;
+
+			case VIPS_FORMAT_USHORT:
+				SHRINK( int, int, unsigned short );
+				break;
+
+			case VIPS_FORMAT_INT:
+				SHRINK( int, gint64, signed int );
+				break;
+
+			case VIPS_FORMAT_UINT:
+				SHRINK( int, gint64, unsigned int );
+				break;
+
+			default:
+				g_assert_not_reached();
+				break;
+			}
+
+			q += output_width;
+		}
+	}
+}
+
+/* Unpack from @tile within @region to the int data pointers on @image. No
+ * subsampling.
+ */
+static void
+vips_foreign_save_jp2k_unpack_image( VipsRegion *region, VipsRect *tile,
+	opj_image_t *image )
+{
+	VipsImage *im = region->im;
+	int b = im->Bands;
+
+	int x, y, i;
+
+	for( y = 0; y < tile->height; y++ ) {
+		VipsPel *p = VIPS_REGION_ADDR( region, 
+			tile->left, tile->top + y );
+
+		for( i = 0; i < b; i++ ) {
+			opj_image_comp_t *comp = &image->comps[i];
+                        int *q = comp->data + y * comp->w; 
+
+			switch( im->BandFmt ) {
+			case VIPS_FORMAT_CHAR:
+			case VIPS_FORMAT_UCHAR:
+				UNPACK( int, unsigned char );
+				break;
+
+			case VIPS_FORMAT_SHORT:
+			case VIPS_FORMAT_USHORT:
+				UNPACK( int, unsigned short );
+				break;
+
+			case VIPS_FORMAT_INT:
+			case VIPS_FORMAT_UINT:
+				UNPACK( int, unsigned int );
+				break;
+
+			default:
+				g_assert_not_reached();
+				break;
+			}
+		}
+	}
+}
+
+void
+vips__foreign_load_jp2k_compress_free( TileCompress *compress )
+{
+	VIPS_FREEF( opj_destroy_codec, compress->codec );
+	VIPS_FREEF( opj_image_destroy, compress->image );
+	VIPS_FREEF( opj_stream_destroy, compress->stream );
+	VIPS_FREE( compress->accumulate );
+}
+
+/* Compress area @tile within @region and write to @target as a @tile_width by
+ * @tile_height jp2k compressed image. This is called from eg. vips2tiff to 
+ * write jp2k-compressed tiles.
+ *
+ * You'd think we could reuse things like the encoder between calls but ...
+ * nope, openjpeg does not allow that.
+ */
+int
+vips__foreign_load_jp2k_compress( VipsRegion *region, 
+	VipsRect *tile, VipsTarget *target,
+	int tile_width, int tile_height,
+        gboolean save_as_ycc, gboolean subsample, gboolean lossless, int Q )
+{
+	TileCompress compress = { 0 };
+	opj_cparameters_t parameters;
+	size_t sizeof_line;
+
+	/* Our rgb->ycc only works for exactly 3 bands.
+	 */
+	save_as_ycc = save_as_ycc && region->im->Bands == 3;
+	subsample = subsample && save_as_ycc;
+
+	/* Set compression params.
+	 */
+	opj_set_default_encoder_parameters( &parameters );
+	parameters.numresolution = 1;
+
+	/* Makes three band images smaller, somehow.
+	 */
+	parameters.tcp_mct = (region->im->Bands >= 3 && !subsample) ? 1 : 0;
+
+	/* Set compression profile.
+	 */
+	vips_foreign_save_jp2k_set_profile( &parameters, lossless, Q ); 
+
+	/* Create output image. TRUE means we alloc memory for the image
+	 * planes.
+	 */
+	if( !(compress.image = vips_foreign_save_jp2k_new_image( region->im,
+		tile_width, tile_height, subsample, save_as_ycc, TRUE )) ) {
+		vips__foreign_load_jp2k_compress_free( &compress );
+		return( -1 );
+	}
+
+	/* We need a line of sums for chroma subsample. At worst, gint64.
+	 */
+	sizeof_line = sizeof( gint64 ) * tile->width;
+	if( !(compress.accumulate = 
+		VIPS_ARRAY( NULL, sizeof_line, VipsPel )) ) {
+		vips__foreign_load_jp2k_compress_free( &compress );
+		return( -1 );
+	}
+
+	compress.codec = opj_create_compress( OPJ_CODEC_J2K );
+	vips_foreign_save_jp2k_attach_handlers( compress.codec );
+        if( !opj_setup_encoder( compress.codec, 
+		&parameters, compress.image ) ) {
+		vips__foreign_load_jp2k_compress_free( &compress );
+		return( -1 );
+	}
+
+#ifdef HAVE_LIBOPENJP2_THREADING
+	/* Use eg. VIPS_CONCURRENCY etc. to set n-cpus, if this openjpeg has
+	 * stable support. 
+	 */
+	opj_codec_set_threads( compress.codec, vips_concurrency_get() );
+#endif /*HAVE_LIBOPENJP2_THREADING*/
+
+	if( save_as_ycc ) 
+		vips_foreign_save_jp2k_rgb_to_ycc( region, 
+			tile, compress.image->comps[0].prec );
+
+	/* we need to unpack to the int arrays on comps[i].data
+	 */
+	if( subsample )
+		vips_foreign_save_jp2k_unpack_subsample_image( region,
+			tile, compress.image, 
+			compress.accumulate ); 
+	else
+		vips_foreign_save_jp2k_unpack_image( region,
+			tile, compress.image ); 
+
+	if( !(compress.stream = vips_foreign_save_jp2k_target( target )) ) {
+		vips__foreign_load_jp2k_compress_free( &compress );
+		return( -1 );
+	}
+
+	if( !opj_start_compress( compress.codec, 
+		compress.image, compress.stream ) ) {
+		vips__foreign_load_jp2k_compress_free( &compress );
+		return( -1 );
+	}
+
+	if( !opj_encode( compress.codec, compress.stream ) ) {
+		vips__foreign_load_jp2k_compress_free( &compress );
+		return( -1 );
+	}
+
+	opj_end_compress( compress.codec, compress.stream );
+
+	vips__foreign_load_jp2k_compress_free( &compress );
+
+	return( 0 );
+}
+
 #endif /*HAVE_LIBOPENJP2*/
 
 /**
@@ -1067,9 +1414,9 @@ vips_foreign_save_jp2k_target_init( VipsForeignSaveJp2kTarget *target )
  * Write a VIPS image to a file in JPEG2000 format. 
  * The saver supports 8, 16 and 32-bit int pixel
  * values, signed and unsigned. It supports greyscale, RGB, CMYK and
- * multispectral images. 
+ * multispectral images.
  *
- * Use @Q to set the compression quality factor. The default value of 45
+ * Use @Q to set the compression quality factor. The default value
  * produces file with approximately the same size as regular JPEG Q 75.
  *
  * Set @lossless to enable lossless compresion.
@@ -1113,7 +1460,7 @@ vips_jp2ksave( VipsImage *in, const char *filename, ... )
  * * @tile_height: %gint for tile size
  * * @subsample_mode: #VipsForeignSubsample, chroma subsampling mode
  *
- * As vips_jp2ksave(), but save to a memory buffer.
+ * As vips_jp2ksave(), but save to a target.
  *
  * See also: vips_jp2ksave(), vips_image_write_to_target().
  *

libvips/foreign/pforeign.h

@@ -237,6 +237,14 @@ struct heif_image;
 void vips__heif_image_print( struct heif_image *img );
 
 extern const char *vips__jp2k_suffs[];
+int vips__foreign_load_jp2k_decompress( VipsImage *out,
+	int width, int height, gboolean ycc_to_rgb,
+	void *from, size_t from_length,
+	void *to, size_t to_length );
+int vips__foreign_load_jp2k_compress( VipsRegion *region,
+	VipsRect *tile, VipsTarget *target,
+	int tile_width, int tile_height,
+	gboolean save_as_ycc, gboolean subsample, gboolean lossless, int Q );
 
 extern const char *vips__jxl_suffs[];
 

libvips/foreign/tiff2vips.c

@@ -201,6 +201,8 @@
  * 	- add subifd
  * 6/6/20 MathemanFlo
  * 	- support 2 and 4 bit greyscale load
+ * 27/3/21
+ * 	- add jp2k decompresion
  */
 
 /*
@@ -253,6 +255,23 @@
 #include "pforeign.h"
 #include "tiff.h"
 
+/* Aperio TIFFs (svs) use these compression types for jp2k-compressed tiles.
+ */
+#define JP2K_YCC 33003
+#define JP2K_RGB 33005
+
+/* Bioformats uses this tag for jp2k compressed tiles.
+ */
+#define JP2K_LOSSY 33004
+
+/* Compression types we handle ourselves.
+ */
+static int rtiff_we_decompress[] = {
+	JP2K_YCC,
+	JP2K_RGB,
+	JP2K_LOSSY
+};
+
 /* What we read from the tiff dir to set our read strategy. For multipage
  * read, we need to read and compare lots of these, so it needs to be broken
  * out as a separate thing.
@@ -316,6 +335,11 @@ typedef struct _RtiffHeader {
 	 */
 	char *image_description;
 
+	/* TRUE if the compression type is not supported by libtiff directly
+	 * and we must read the raw data and decompress ourselves.
+	 */
+	gboolean we_decompress;
+
 } RtiffHeader;
 
 /* Scanline-type process function.
@@ -371,6 +395,12 @@ typedef struct _Rtiff {
 	 */
 	tdata_t contig_buf;
 
+	/* If we are decompressing, we need a buffer to read the raw tile to
+	 * before running the decompressor.
+	 */
+	tdata_t compressed_buf;
+	tsize_t compressed_buf_length;
+
 	/* The Y we are reading at. Used to verify strip read is sequential.
 	 */
 	int y_pos;
@@ -1720,12 +1750,47 @@ static int
 rtiff_read_tile( Rtiff *rtiff, tdata_t *buf, int x, int y )
 {
 #ifdef DEBUG_VERBOSE
-	printf( "rtiff_read_tile: x = %d, y = %d\n", x, y ); 
+	printf( "rtiff_read_tile: x = %d, y = %d, we_decompress = %d\n", 
+		x, y, rtiff->header.we_decompress ); 
 #endif /*DEBUG_VERBOSE*/
 
-	if( TIFFReadTile( rtiff->tiff, buf, x, y, 0, 0 ) < 0 ) { 
+	if( rtiff->header.we_decompress ) {
-		vips_foreign_load_invalidate( rtiff->out );
+		ttile_t tile_no = TIFFComputeTile( rtiff->tiff, x, y, 0, 0 );
-		return( -1 ); 
+
+		tsize_t size;
+
+		size = TIFFReadRawTile( rtiff->tiff, tile_no, 
+			rtiff->compressed_buf, rtiff->compressed_buf_length );
+		if( size <= 0 ) {
+			vips_foreign_load_invalidate( rtiff->out );
+			return( -1 ); 
+		}
+
+		switch( rtiff->header.compression ) {
+		case JP2K_YCC:
+		case JP2K_RGB:
+		case JP2K_LOSSY:
+			if( vips__foreign_load_jp2k_decompress( 
+				rtiff->out, 
+				rtiff->header.tile_width, 
+				rtiff->header.tile_height,
+				TRUE,
+				rtiff->compressed_buf, size,
+				buf, rtiff->header.tile_size ) ) 
+				return( -1 );
+			break;
+
+		default:
+			g_assert_not_reached();
+			break;
+		}
+
+	}
+	else {
+		if( TIFFReadTile( rtiff->tiff, buf, x, y, 0, 0 ) < 0 ) { 
+			vips_foreign_load_invalidate( rtiff->out );
+			return( -1 ); 
+		}
 	}
 
 	return( 0 ); 
@@ -1982,6 +2047,20 @@ rtiff_read_tilewise( Rtiff *rtiff, VipsImage *out )
 		return( -1 );
 	}
 
+	/* If we will be decompressing, we need a buffer large enough to hold
+	 * the largest compressed tile in any page.
+	 *
+	 * Allocate a buffer 2x the uncompressed tile size ... much simpler
+	 * than searching every page for the largest tile with
+	 * TIFFTAG_TILEBYTECOUNTS.
+	 */
+	if( rtiff->header.we_decompress ) {
+		rtiff->compressed_buf_length = 2 * rtiff->header.tile_size;
+		if( !(rtiff->compressed_buf = vips_malloc( VIPS_OBJECT( out ), 
+			rtiff->compressed_buf_length )) )
+			return( -1 );
+	}
+
 	/* Read to this image, then cache to out, see below.
 	 */
 	t[0] = vips_image_new(); 
@@ -2355,6 +2434,7 @@ rtiff_read_stripwise( Rtiff *rtiff, VipsImage *out )
 static int
 rtiff_header_read( Rtiff *rtiff, RtiffHeader *header )
 {
+	int i;
 	uint16 extra_samples_count;
 	uint16 *extra_samples_types;
 	toff_t *subifd_offsets;
@@ -2377,6 +2457,19 @@ rtiff_header_read( Rtiff *rtiff, RtiffHeader *header )
 	TIFFGetFieldDefaulted( rtiff->tiff, 
 		TIFFTAG_COMPRESSION, &header->compression );
 
+	/* One of the types we decompress?
+	 */
+	for( i = 0; i < VIPS_NUMBER( rtiff_we_decompress ); i++ )
+		if( header->compression == rtiff_we_decompress[i] ) {
+#ifdef DEBUG
+			printf( "rtiff_header_read: "
+				"compression %d handled by us\n", 
+				header->compression );
+#endif /*DEBUG*/
+			header->we_decompress = TRUE;
+			break;
+		}
+
 	/* We must set this here since it'll change the value of scanline_size.
 	 */
 	rtiff_set_decode_format( rtiff );
@@ -2471,6 +2564,8 @@ rtiff_header_read( Rtiff *rtiff, RtiffHeader *header )
 	 */
 	header->tiled = TIFFIsTiled( rtiff->tiff );
 
+
+
 #ifdef DEBUG
 	printf( "rtiff_header_read: header.width = %d\n", 
 		header->width );
@@ -2652,6 +2747,7 @@ rtiff_header_equal( RtiffHeader *h1, RtiffHeader *h2 )
 		h1->photometric_interpretation != 
 			h2->photometric_interpretation ||
 		h1->sample_format != h2->sample_format ||
+		h1->compression != h2->compression ||
 		h1->separate != h2->separate ||
 		h1->tiled != h2->tiled ||
 		h1->orientation != h2->orientation )

libvips/foreign/vips2tiff.c

@@ -270,6 +270,16 @@
  */
 #define MAX_ALPHA (64)
 
+/* Bioformats uses this tag for lossy jp2k compressed tiles.
+ */
+#define JP2K_LOSSY 33004
+
+/* Compression types we handle ourselves.
+ */
+static int wtiff_we_compress[] = {
+	JP2K_LOSSY
+};
+
 typedef struct _Layer Layer;
 typedef struct _Wtiff Wtiff;
 
@@ -351,7 +361,7 @@ struct _Wtiff {
 	int strip;			/* Don't write metadata */
 	VipsRegionShrink region_shrink; /* How to shrink regions */
 	int level;			/* zstd compression level */
-	gboolean lossless;		/* webp lossless mode */
+	gboolean lossless;		/* lossless mode */
 	VipsForeignDzDepth depth;	/* Pyr depth */
 	gboolean subifd;		/* Write pyr layers into subifds */
 
@@ -368,6 +378,16 @@ struct _Wtiff {
 	 * roll mode.
 	 */
 	int image_height;
+
+	/* TRUE if the compression type is not supported by libtiff directly
+	 * and we must compress ourselves. 
+	 */
+	gboolean we_compress;
+
+	/* If we are copying, we need a buffer to read the compressed tile to.
+	 */
+	tdata_t compressed_buf;
+	tsize_t compressed_buf_length;
 };
 
 /* Write an ICC Profile from a file into the JPEG stream.
@@ -636,6 +656,7 @@ wtiff_write_header( Wtiff *wtiff, Layer *layer )
 {
 	TIFF *tif = layer->tif;
 
+	int i;
 	int orientation; 
 
 #ifdef DEBUG
@@ -672,6 +693,12 @@ wtiff_write_header( Wtiff *wtiff, Layer *layer )
 		wtiff->predictor != VIPS_FOREIGN_TIFF_PREDICTOR_NONE ) 
 		TIFFSetField( tif, TIFFTAG_PREDICTOR, wtiff->predictor );
 
+	for( i = 0; i < VIPS_NUMBER( wtiff_we_compress ); i++ )
+		if( wtiff->compression == wtiff_we_compress[i] ) {
+			wtiff->we_compress = TRUE;
+			break;
+		}
+
 	/* Don't write mad resolutions (eg. zero), it confuses some programs.
 	 */
 	TIFFSetField( tif, TIFFTAG_RESOLUTIONUNIT, wtiff->resunit );
@@ -926,6 +953,20 @@ wtiff_allocate_layers( Wtiff *wtiff )
 			return( -1 );
 	}
 
+	/* If we will be copying layers we need a buffer large enough to hold
+	 * the largest compressed tile in any page.
+	 *
+	 * Allocate a buffer 2x the uncompressed tile size ... much simpler
+	 * than searching every page for the largest tile with
+	 * TIFFTAG_TILEBYTECOUNTS.
+	 */
+	if( wtiff->pyramid ) {
+		wtiff->compressed_buf_length = 2 * wtiff->tls * wtiff->tileh;
+		if( !(wtiff->compressed_buf = vips_malloc( NULL,
+			wtiff->compressed_buf_length )) )
+			return( -1 );
+	}
+
 	return( 0 );
 }
 
@@ -981,11 +1022,11 @@ static void
 wtiff_free( Wtiff *wtiff )
 {
 	wtiff_delete_temps( wtiff );
-
 	VIPS_UNREF( wtiff->ready );
 	VIPS_FREE( wtiff->tbuf );
 	VIPS_FREEF( layer_free_all, wtiff->layer );
 	VIPS_FREE( wtiff->filename );
+	VIPS_FREE( wtiff->compressed_buf );
 	VIPS_FREE( wtiff );
 }
 
@@ -1011,6 +1052,8 @@ get_compression( VipsForeignTiffCompression compression )
 	case VIPS_FOREIGN_TIFF_COMPRESSION_ZSTD:
 		return( COMPRESSION_ZSTD );
 #endif /*HAVE_TIFF_COMPRESSION_WEBP*/
+	case VIPS_FOREIGN_TIFF_COMPRESSION_JP2K:
+		return( JP2K_LOSSY );
 	
 	default:
 		return( COMPRESSION_NONE );
@@ -1237,21 +1280,6 @@ wtiff_new( VipsImage *input, const char *filename,
 		wtiff->miniswhite = FALSE;
 	}
 
-	/* lossless is for webp only.
-	 */
-#ifdef HAVE_TIFF_COMPRESSION_WEBP
-	if( wtiff->lossless ) {
-		if( wtiff->compression == COMPRESSION_NONE )
-			wtiff->compression = COMPRESSION_WEBP;
-
-		if( wtiff->compression != COMPRESSION_WEBP ) {
-			g_warning( "%s", 
-				_( "lossless is for WEBP compression only" ) );
-			wtiff->lossless = FALSE;
-		}
-	}
-#endif /*HAVE_TIFF_COMPRESSION_WEBP*/
-
 	/* Sizeof a line of bytes in the TIFF tile.
 	 */
 	if( wtiff->ready->Coding == VIPS_CODING_LABQ )
@@ -1489,7 +1517,7 @@ wtiff_pack2tiff( Wtiff *wtiff, Layer *layer,
 /* Write a set of tiles across the strip.
  */
 static int
-wtiff_layer_write_tile( Wtiff *wtiff, Layer *layer, VipsRegion *strip )
+wtiff_layer_write_tiles( Wtiff *wtiff, Layer *layer, VipsRegion *strip )
 {
 	VipsImage *im = layer->image;
 	VipsRect *area = &strip->valid;
@@ -1511,21 +1539,72 @@ wtiff_layer_write_tile( Wtiff *wtiff, Layer *layer, VipsRegion *strip )
 		tile.height = wtiff->tileh;
 		vips_rect_intersectrect( &tile, &image, &tile );
 
-		/* Have to repack pixels.
-		 */
-		wtiff_pack2tiff( wtiff, layer, strip, &tile, wtiff->tbuf );
-
 #ifdef DEBUG_VERBOSE
 		printf( "Writing %dx%d tile at position %dx%d to image %s\n",
 			tile.width, tile.height, tile.left, tile.top,
 			TIFFFileName( layer->tif ) );
 #endif /*DEBUG_VERBOSE*/
 
-		if( TIFFWriteTile( layer->tif, wtiff->tbuf, 
+		if( wtiff->we_compress ) {
-			tile.left, tile.top, 0, 0 ) < 0 ) {
+			ttile_t tile_no = TIFFComputeTile( layer->tif,
-			vips_error( "vips2tiff", 
+				tile.left, tile.top, 0, 0 );
-				"%s", _( "TIFF write tile failed" ) );
+
-			return( -1 );
+			VipsTarget *target;
+			int result;
+			unsigned char *buffer;
+			size_t length;
+
+			target = vips_target_new_to_memory();
+
+			switch( wtiff->compression ) {
+			case JP2K_LOSSY:
+				result = vips__foreign_load_jp2k_compress( 
+					strip, &tile, target,
+					wtiff->tilew, wtiff->tileh,
+					!wtiff->rgbjpeg,
+					!wtiff->rgbjpeg && wtiff->Q < 90,
+					wtiff->lossless, 
+					wtiff->Q );
+				break;
+
+			default:
+				result = -1;
+				g_assert_not_reached();
+				break;
+			}
+
+			if( result ) {
+				g_object_unref( target );
+				return( -1 );
+			}
+
+			buffer = vips_target_steal( target, &length );
+
+			g_object_unref( target );
+
+			result = TIFFWriteRawTile( layer->tif, tile_no, 
+				buffer, length );
+
+			g_free( buffer );
+		
+			if( result < 0 ) {
+				vips_error( "vips2tiff", 
+					"%s", _( "TIFF write tile failed" ) );
+				return( -1 );
+			}
+		}
+		else {
+			/* Have to repack pixels for libtiff.
+			 */
+			wtiff_pack2tiff( wtiff, 
+				layer, strip, &tile, wtiff->tbuf );
+
+			if( TIFFWriteTile( layer->tif, wtiff->tbuf, 
+				tile.left, tile.top, 0, 0 ) < 0 ) {
+				vips_error( "vips2tiff", 
+					"%s", _( "TIFF write tile failed" ) );
+				return( -1 );
+			}
 		}
 	}
 
@@ -1672,7 +1751,7 @@ layer_strip_arrived( Layer *layer )
 	VipsRect image_area;
 
 	if( wtiff->tile ) 
-		result = wtiff_layer_write_tile( wtiff, layer, layer->strip );
+		result = wtiff_layer_write_tiles( wtiff, layer, layer->strip );
 	else
 		result = wtiff_layer_write_strip( wtiff, layer, layer->strip );
 	if( result )
@@ -1803,7 +1882,7 @@ wtiff_copy_tiff( Wtiff *wtiff, TIFF *out, TIFF *in )
 	uint16 ui16_2;
 	float f;
 	tdata_t buf;
-	ttile_t tile;
+	ttile_t tile_no;
 	ttile_t n;
 	uint16 *a;
 
@@ -1891,19 +1970,15 @@ wtiff_copy_tiff( Wtiff *wtiff, TIFF *out, TIFF *in )
 
 	buf = vips_malloc( NULL, TIFFTileSize( in ) );
 	n = TIFFNumberOfTiles( in );
-	for( tile = 0; tile < n; tile++ ) {
+	for( tile_no = 0; tile_no < n; tile_no++ ) {
 		tsize_t len;
 
-		/* It'd be good to use TIFFReadRawTile()/TIFFWtiffRawTile() 
+		len = TIFFReadRawTile( in, tile_no, 
-		 * here to save compression/decompression, but sadly it seems
+			wtiff->compressed_buf, wtiff->compressed_buf_length );
-		 * not to work :-( investigate at some point.
+		if( len <= 0 ||
-		 */
+			TIFFWriteRawTile( out, tile_no, 
-		len = TIFFReadEncodedTile( in, tile, buf, -1 );
+				wtiff->compressed_buf, len ) < 0 )
-		if( len < 0 ||
-			TIFFWriteEncodedTile( out, tile, buf, len ) < 0 ) {
-			g_free( buf );
 			return( -1 );
-		}
 	}
 	g_free( buf );
 

libvips/include/vips/foreign.h

@@ -465,6 +465,7 @@ int vips_webpsave_mime( VipsImage *in, ... )
  * @VIPS_FOREIGN_TIFF_COMPRESSION_LZW: LZW compression
  * @VIPS_FOREIGN_TIFF_COMPRESSION_WEBP: WEBP compression
  * @VIPS_FOREIGN_TIFF_COMPRESSION_ZSTD: ZSTD compression
+ * @VIPS_FOREIGN_TIFF_COMPRESSION_JP2K: JP2K compression
  *
  * The compression types supported by the tiff writer.
  *
@@ -485,6 +486,7 @@ typedef enum {
 	VIPS_FOREIGN_TIFF_COMPRESSION_LZW,
 	VIPS_FOREIGN_TIFF_COMPRESSION_WEBP,
 	VIPS_FOREIGN_TIFF_COMPRESSION_ZSTD,
+	VIPS_FOREIGN_TIFF_COMPRESSION_JP2K,
 	VIPS_FOREIGN_TIFF_COMPRESSION_LAST
 } VipsForeignTiffCompression;
 

libvips/include/vips/region.h

@@ -165,18 +165,16 @@ void vips_region_invalidate( VipsRegion *reg );
 #define VIPS_COUNT_PIXELS( R, N ) 
 #endif /*DEBUG_LEAK*/
 
-/* Macros on VipsRegion.
- *	VIPS_REGION_LSKIP()		add to move down line
- *	VIPS_REGION_N_ELEMENTS()	number of elements across region
- *	VIPS_REGION_SIZEOF_LINE()	sizeof width of region
- *	VIPS_REGION_ADDR()		address of pixel in region
- */
 #define VIPS_REGION_LSKIP( R ) \
 	((size_t)((R)->bpl))
 #define VIPS_REGION_N_ELEMENTS( R ) \
 	((size_t)((R)->valid.width * (R)->im->Bands))
+#define VIPS_REGION_SIZEOF_ELEMENT( R ) \
+	(VIPS_IMAGE_SIZEOF_ELEMENT( (R)->im ))
+#define VIPS_REGION_SIZEOF_PEL( R ) \
+	(VIPS_IMAGE_SIZEOF_PEL( (R)->im ))
 #define VIPS_REGION_SIZEOF_LINE( R ) \
-	((size_t)((R)->valid.width * VIPS_IMAGE_SIZEOF_PEL( (R)->im) ))
+	((size_t)((R)->valid.width * VIPS_REGION_SIZEOF_PEL( R )))
 
 /* If DEBUG is defined, add bounds checking.
  */
@@ -184,7 +182,7 @@ void vips_region_invalidate( VipsRegion *reg );
 #define VIPS_REGION_ADDR( R, X, Y ) \
 	( (vips_rect_includespoint( &(R)->valid, (X), (Y) ))? \
 	  ((R)->data + ((Y) - (R)->valid.top) * VIPS_REGION_LSKIP(R) + \
-	  ((X) - (R)->valid.left) * VIPS_IMAGE_SIZEOF_PEL((R)->im)): \
+	  ((X) - (R)->valid.left) * VIPS_REGION_SIZEOF_PEL( R )): \
 	  (fprintf( stderr, \
 		"VIPS_REGION_ADDR: point out of bounds, " \
 		"file \"%s\", line %d\n" \
@@ -202,7 +200,7 @@ void vips_region_invalidate( VipsRegion *reg );
 #define VIPS_REGION_ADDR( R, X, Y ) \
 	((R)->data + \
 	((Y)-(R)->valid.top) * VIPS_REGION_LSKIP( R ) + \
-	((X)-(R)->valid.left) * VIPS_IMAGE_SIZEOF_PEL( (R)->im ))
+	((X)-(R)->valid.left) * VIPS_REGION_SIZEOF_PEL( R ))
 #endif /*DEBUG*/
 
 #define VIPS_REGION_ADDR_TOPLEFT( R ) ((R)->data)

libvips/iofuncs/enumtypes.c

@@ -574,6 +574,7 @@ vips_foreign_tiff_compression_get_type( void )
 			{VIPS_FOREIGN_TIFF_COMPRESSION_LZW, "VIPS_FOREIGN_TIFF_COMPRESSION_LZW", "lzw"},
 			{VIPS_FOREIGN_TIFF_COMPRESSION_WEBP, "VIPS_FOREIGN_TIFF_COMPRESSION_WEBP", "webp"},
 			{VIPS_FOREIGN_TIFF_COMPRESSION_ZSTD, "VIPS_FOREIGN_TIFF_COMPRESSION_ZSTD", "zstd"},
+			{VIPS_FOREIGN_TIFF_COMPRESSION_JP2K, "VIPS_FOREIGN_TIFF_COMPRESSION_JP2K", "jp2k"},
 			{VIPS_FOREIGN_TIFF_COMPRESSION_LAST, "VIPS_FOREIGN_TIFF_COMPRESSION_LAST", "last"},
 			{0, NULL, NULL}
 		};

test/test-suite/test_foreign.py

@@ -525,6 +525,16 @@ class TestForeign:
             z = y.hist_find(band=0)
             assert z(0, 0)[0] + z(255, 0)[0] == y.width * y.height
 
+    @skip_if_no("jp2kload")
+    @skip_if_no("tiffload")
+    def test_tiffjp2k(self):
+        self.save_load_file(".tif", "[tile,compression=jp2k]", self.colour, 80)
+        self.save_load_file(".tif",
+                            "[tile,pyramid,compression=jp2k]", self.colour, 80)
+        self.save_load_file(".tif",
+                            "[tile,pyramid,subifd,compression=jp2k]",
+                            self.colour, 80)
+
     @skip_if_no("magickload")
     def test_magickload(self):
         def bmp_valid(im):