fix matlab load

ChangeLog

@@ -14,6 +14,7 @@
 - check for SetImageOption() so we work with GraphicsMagick too
 - "header" sets a non-zero exit code if anything failed
 - add and use im_check_uncoded() and friends
+- matlab load handles column-major and plane-separated images (thanks Mikhail)
 
 25/3/09 started 7.18.0
 - revised version numbers

libsrc/format/matlab.c

@@ -1,4 +1,7 @@
 /* Read matlab save files with libmatio
+ *
+ * 4/8/09
+ *	- transpose on load, assemble planes into bands (thanks Mikhail)
  */
 
 /*
@@ -31,10 +34,6 @@
 
 	Remaining issues:
 
-+ it transposes images, I think, since Matlab seems to use column-major order
-
-+ will colour images work? no idea, needs testing
-
 + it will not do complex images
 
 + it will not handle sparse matricies
@@ -44,9 +43,6 @@
 
 + load only, no save
 
-+ could use much less memory --- we use Mat_VarReadDataAll() to read the 
-  whole variable into mem, then copy to a vips buffer, yuk
-
  */
 
 /*
@@ -226,6 +222,13 @@ static int
 mat2vips_get_data( mat_t *mat, matvar_t *var, IMAGE *im )
 {
 	int y;
+	PEL *buffer;
+	const int es = IM_IMAGE_SIZEOF_ELEMENT( im );
+
+	/* Matlab images are plane-separate, so we have to assemble bands in
+	 * image-size chunks.
+	 */
+	const int is = es * im->Xsize * im->Ysize;
 
 	if( Mat_VarReadDataAll( mat, var ) ) {
 		im_error( "mat2vips", "%s", _( "Mat_VarReadDataAll failed" ) );
@@ -235,10 +238,37 @@ mat2vips_get_data( mat_t *mat, matvar_t *var, IMAGE *im )
 		im_setupout( im ) )
 		return( -1 );
 
-	for( y = 0; y < im->Ysize; y++ )
+	/* Matlab images are in columns, so we have to transpose into
-		if( im_writeline( y, im, 
+	 * scanlines with this buffer.
-			var->data + y * IM_IMAGE_SIZEOF_LINE( im ) ) )
+	 */
+	if( !(buffer = IM_ARRAY( im, IM_IMAGE_SIZEOF_LINE( im ), PEL )) )
+		return( -1 );
+
+	for( y = 0; y < im->Ysize; y++ ) {
+		const PEL *p = var->data + y * es;
+		int x;
+		PEL *q;
+
+		q = buffer;
+		for( x = 0; x < im->Xsize; x++ ) {
+			int b;
+
+			for( b = 0; b < im->Bands; b++ ) {
+				const PEL *p2 = p + b * is;
+				int z;
+
+				for( z = 0; z < es; z++ )
+					q[z] = p2[z];
+
+				q += es;
+			}
+
+			p += es * im->Ysize;
+		}
+
+		if( im_writeline( y, im, buffer ) )
 			return( -1 );
+	}
 
 	return( 0 );
 }