/*
 *
 * Imview, the portable image analysis application
 * http://www.cmis.csiro.au/Hugues.Talbot/imview
 * ----------------------------------------------------------
 *
 *  Imview is an attempt to provide an image display application
 *  suitable for professional image analysis. It was started in
 *  1997 and is mostly the result of the efforts of Hugues Talbot,
 *  Image Analysis Project, CSIRO Mathematical and Information
 *  Sciences, with help from others (see the CREDITS files for
 *  more information)
 *
 *  Imview is Copyrighted (C) 1997-2001 by Hugues Talbot and was
 *  supported in parts by the Australian Commonwealth Science and 
 *  Industry Research Organisation. Please see the COPYRIGHT file 
 *  for full details. Imview also includes the contributions of 
 *  many others. Please see the CREDITS file for full details.
 *
 *  This program is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation; either version 2 of the License, or
 *  (at your option) any later version.
 *  
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *  
 *  You should have received a copy of the GNU General Public License
 *  along with this program; if not, write to the Free Software
 *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111, USA.
 * */

/*------------------------------------------------------------------------
 *
 * Read data from Matlab level 5 MAT-files
 * 
 * Teemu Ikonen 	 16 Aug 2004
 *
 *      
 *-----------------------------------------------------------------------*/

#include <fstream>
#include "imnmspc.hxx" // name space definition if needed 

#include "imunistd.h"

#include "imview.hxx"
#include "imageIO.hxx"
#include "machine.hxx" // portability stuff
#include "readmat5.hxx"
#include "byte-swap.h"

#define PAD(l) (((l)<=4)?4:(((l)+7)/8)*8)
#define TAGLENGTH(l) ((l)<=4?4:8)
#define X_CAST(T, E) (T) (E)

extern imageIO *IOBlackBox;

typedef uchar   uint8;
typedef char    int8;
typedef ushort  uint16;
typedef short   int16;
typedef uint    uint32;
typedef int     int32;
//typedef unsigned long uint64;
//typedef long    int64;
typedef float   single;
// typedef double  double;

enum mat5_data_type
{
    miINT8 = 1,			// 8 bit signed
    miUINT8,			// 8 bit unsigned
    miINT16,			// 16 bit signed
    miUINT16,			// 16 bit unsigned
    miINT32,			// 32 bit signed
    miUINT32,			// 32 bit unsigned
    miSINGLE,			// IEEE 754 single precision float
    miRESERVE1,
    miDOUBLE,			// IEEE 754 double precision float
    miRESERVE2,
    miRESERVE3,
    miINT64,			// 64 bit signed
    miUINT64,			// 64 bit unsigned
    miMATRIX			// MATLAB array
};

enum arrayclasstype
  {
    mxCELL_CLASS=1,		// cell array
    mxSTRUCT_CLASS,		// structure
    mxOBJECT_CLASS,		// object
    mxCHAR_CLASS,		// character array
    mxSPARSE_CLASS,		// sparse array
    mxDOUBLE_CLASS,		// double precision array
    mxSINGLE_CLASS,		// single precision floating point
    mxINT8_CLASS,		// 8 bit signed integer
    mxUINT8_CLASS,		// 8 bit unsigned integer
    mxINT16_CLASS,		// 16 bit signed integer
    mxUINT16_CLASS,		// 16 bit unsigned integer
    mxINT32_CLASS,		// 32 bit signed integer
    mxUINT32_CLASS		// 32 bit unsigned integer
  };

int read_mat5_binary_file_header (std::istream& is, bool& swap);
static int read_mat5_tag (std::istream& is, bool swap, int& type, int& bytes);
int read_mat5_binary_data (std::istream& is, void **&buf, int len,
			   int matcount, bool swap,
			   mat5_data_type type, pixtype &imtype);
void read_int32 (std::istream& is, bool swap_bytes, int32& val);


int mat5nbsubfiles(const char *ident)
{
    return 1; // we may start to support more frames in the future
}

int readMat5(const char *ident, int index)
{
    void **buf;
    bool swap;
    int res;
    int type, bytes, len;
    
    dbgprintf("Read Mat5 called\n");

    std::ifstream is(ident);
    if(!is) {
	errprintf("readmat5: Cannot open file %s\n", ident);
	return 10;
    }
	
    res = read_mat5_binary_file_header (is, swap);
    if(res)
	return res;

    dbgprintf("readmat5: byte-swap is %s.\n", swap ? "on" : "off");
   
    read_mat5_tag (is, swap, type, bytes);
    if(type != miMATRIX) {
	errprintf("Only matrix type MAT5-datafiles supported.\n");
	return 100;
    }

    // array flags subelement
    if (read_mat5_tag (is, swap, type, len) || type != miUINT32 || len != 8)
    {
      errprintf ("readmat5: invalid array flags subelement\n");
      return 100;
    }
    int32 flags;
    read_int32(is, swap, flags);
    bool imag = (flags & 0x0800) != 0;	// has an imaginary part?
    if(imag) {
      errprintf ("readmat5: complex arrays not supported\n");
      return 100;
    }
//    bool global = (flags & 0x0400) != 0; // global variable?
//    bool logicalvar = (flags & 0x0200) != 0; 
    enum arrayclasstype arrayclass = (arrayclasstype)(flags & 0xff);
    int32 junk;
    read_int32 (is, swap, junk);	// an "undefined" entry

    // dimensions array subelement
    int32 dim_len;
    
    if (read_mat5_tag (is, swap, type, dim_len) || type != miINT32) {
	errprintf ("readmat5: invalid dimensions array subelement");
	return 100;
    }
    
    int ndims = dim_len / 4;
    if(ndims != 2) {
	errprintf ("readmat5: only two-dimensional arrays supported");
	return 100;
    }
    int32 xsize;
    read_int32 (is, swap, xsize);
    int32 ysize;
    read_int32 (is, swap, ysize);    

    dbgprintf("readmat5: matrix has size %d x %d\n", xsize, ysize);
    
    std::streampos tmp_pos = is.tellg ();
    is.seekg (tmp_pos + static_cast<std::streamoff>
	      (PAD (dim_len) - dim_len));
    
    if (read_mat5_tag (is, swap, type, len) || type != miINT8) {
	errprintf ("readmat5: invalid array name subelement\n");
	return 100;
    }
        
    tmp_pos = is.tellg ();

    char name[len+1];    
    if (len) {
	if (! is.read (name, len ))
	    return 100;
	is.seekg (tmp_pos + static_cast<std::streamoff> (PAD (len)));
    }

    pixtype imtype = IM_INVALID;
    
    switch (arrayclass) {
	case mxDOUBLE_CLASS:
	case mxSINGLE_CLASS:
	case mxINT8_CLASS:
	case mxUINT8_CLASS:
	case mxINT16_CLASS:
	case mxUINT16_CLASS:
	case mxINT32_CLASS:
	case mxUINT32_CLASS:
	    // read real data subelement
	    if (read_mat5_tag (is, swap, type, len)) {
		errprintf ("readmat5: error reading matrix data.\n");
		return 100;
	    }
	    dbgprintf("readmat5: data has length %d\n", len);
	    	    
	    tmp_pos = is.tellg ();	    
	    res = read_mat5_binary_data(is, buf, len, xsize*ysize, swap,
					(mat5_data_type) type, imtype);
	    if(res)
		return res;

	    is.seekg (tmp_pos + static_cast<std::streamoff> (PAD (len)));
	    
	    break;
	    
	default:
	    errprintf ("readmat5: unsupported array class.\n");
	    return 100;
    }

    if (res == 0) {
	IOBlackBox->setCurrBuffp(buf);
	IOBlackBox->setCurrImgWidth(xsize);
	IOBlackBox->setCurrImgHeight(ysize);
	IOBlackBox->setCurrImgThickness(1);
	IOBlackBox->setXOffset(0);
	IOBlackBox->setYOffset(0);
	IOBlackBox->setZOffset(0);
	IOBlackBox->setCurrZPos(0); // first plane
	IOBlackBox->setCurrImgType(IM_SINGLE);
	IOBlackBox->setCurrPixType(imtype);
	IOBlackBox->setCurrImgNbComps(1);
	IOBlackBox->setCurrImgNbSamples(1);
	IOBlackBox->setImgDesc("MAT5");
    }
    
    return res;

}

// Code below is adapted from GNU Octave (www.octave.org)

int
read_mat5_binary_file_header (std::istream& is, bool& swap)
{
  uint16 version=0, magic=0;

  is.seekg (124, std::ios::beg);
  is.read ((char *) &version, 2);
  is.read ((char *) &magic, 2);
  
  if (magic == 0x4d49)
      swap = 0;
  else if (magic == 0x494d)
      swap = 1;
  else {
      errprintf ("readmat5: can't read binary file");
      return 10;
  }

  if (! swap)			// version number is inverse swapped!
    version = ((version >> 8) & 0xff) + ((version & 0xff) << 8);

  if (version != 1)
      dbgprintf ("readmat5: found version %d binary MAT file, but only prepared for version 1");

  return 0;
}

// Read one element tag from stream IS, 
// place the type code in TYPE and the byte count in BYTES
// return nonzero on error
static int read_mat5_tag (std::istream& is, bool swap, int& type, int& bytes)
{
  unsigned int upper;
  int32 temp;

  if (! is.read ((char *) &temp, 4 ))
    goto data_read_error;

  if (swap)
    swap_4_bytes ((char *)&temp);

  upper = (temp >> 16) & 0xffff;
  type = temp & 0xffff;

  if (upper)
    {
      // "compressed" format
      bytes = upper;
    }
  else
    {
      if (! is.read (X_CAST (char *, &temp), 4 ))
	goto data_read_error;
      if (swap)
	swap_4_bytes ((char *)&temp);
      bytes = temp;
    }

  return 0;

 data_read_error:
  return 1;
}

#define swap_1_bytes(x, y)

#define LS_DO_READ(swap, data, size, len, stream) \
  do \
    { \
      if (len > 0) \
	{ \
	  stream.read (X_CAST (char *, data), size * len); \
	  if (swap) \
	    swap_ ## size ## _bytes (data, len); \
	} \
    } \
  while (0)


// Return the data read from a mat5 array in buf[0] and the type of
// a pixel in imtype. If the expected array size given in matcount
// does not match the binary array length, return a non-zero errorcode.
int read_mat5_binary_data (std::istream& is, void **&buf, int len,
			   int matcount, bool swap,
			   mat5_data_type type, pixtype &imtype)
{
    buf = (void **)calloc(1, sizeof(void *));
    if(buf == NULL) {
	errprintf("readmat5: Out of memory\n");
	return 20;
    }

    void *data = malloc(len);
    if(data == NULL) {
	errprintf("readmat5: Out of memory\n");
	free(buf);
	return 20;
    }
    
    int count;
    
    switch (type) {
	case miINT8:
	    count = len / 1;
	    LS_DO_READ (swap, ((int8 *) data), 1, count, is);
	    imtype = IM_INT1;
	    break;
	    
	case miUINT8:
	    count = len / 1;	    
	    LS_DO_READ (swap, ((uint8 *)data), 1, count, is);
	    imtype = IM_UINT1;
	    break;
	    
	case miINT16:
	    count = len / 2;	    
	    LS_DO_READ (swap, ((int16 *)data), 2, count, is);
	    imtype = IM_INT2;	  
	    break;
	    
	case miUINT16:
	    count = len / 2;	    
	    LS_DO_READ (swap, ((uint16 *)data), 2, count, is);
	    imtype = IM_UINT2;	  
	    break;
	    
	case miINT32:
	    count = len / 4;	    
	    LS_DO_READ (swap, ((int32 *)data), 4, count, is);
	    imtype = IM_INT4;	  
	    break;
	    
	case miUINT32:
	    count = len / 4;	    
	    LS_DO_READ (swap, ((uint32 *)data), 4, count, is);
	    imtype = IM_UINT4;
	    break;
	    
	case miSINGLE:
	    count = len / 4;	    
	    is.read (X_CAST (char *, data), 4 * count);
	    if(swap)
	        swap_4_bytes(data, count);
	    imtype = IM_FLOAT;	      
	    break;
		
	case miDOUBLE:
	    count = len / 8;	    
	    is.read (X_CAST (char *, data), 8 * count);
	    if(swap)
	        swap_8_bytes(data, count);	  
	    imtype = IM_DOUBLE;	      	  
	    break;	  
	    
	case miRESERVE1:
	case miRESERVE2:
	case miRESERVE3:
	case miINT64:
	case miUINT64:
	case miMATRIX:
	default:
	    free(buf);	    
	    free(data);
	    return 100;
	    break;
    }

    if(count != matcount) {
	errprintf("readmat5: array sizes do not match.\n");
	free(buf);
	free(data);
	return 100;
    }

    *buf = data;
    
    return 0;
}


void read_int32 (std::istream& is, bool swap_bytes, int32& val)
{

    is.read (X_CAST (char *, &val), sizeof (int32));

    if (swap_bytes) {
	swap_4_bytes ((char *) &val);
    }
}