/*
 * Copyright (C) 1997-2003, R3vis Corporation.
 *
 * 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-1307, USA,
 * or visit http://www.gnu.org/copyleft/gpl.html.
 *
 * Contributor(s):
 *   Wes Bethel, R3vis Corporation, Marin County, California
 *
 * The OpenRM project is located at http://openrm.sourceforge.net/.
 */
/* 
 * $Id: libdio.c,v 1.6 2003/10/15 05:47:57 wes Exp $
 * $Revision: 1.6 $
 * $Name: OpenRM-1-5-2-RC1 $
 * $Log: libdio.c,v $
 * Revision 1.6  2003/10/15 05:47:57  wes
 * Increased the size of the char buffer shared between reader routines.
 * This increase fixes a problem whereby long lines of data weren't
 * being read in properly. While not a final fix, the size increase should
 * be adequate to accommodate most data files. The line length limit for
 * input ASCII files is now on the order of 64K.
 *
 * Revision 1.5  2003/04/13 18:13:23  wes
 * Updated copyright dates.
 *
 * Revision 1.4  2003/01/16 22:22:45  wes
 * Updated all source files to reflect new organization of header files: all
 * headers that were formerly located in include/rmaux, include/rmv and
 * include/rmi are now located in include/rm.
 *
 * Revision 1.3  2002/06/17 00:39:15  wes
 * updated copyright line.
 *
 * Revision 1.2  2001/03/31 16:55:18  wes
 * Added procmode.h, which defines an RMpipe processing mode used in
 * most demonstration programs. The default processing mode is
 * RM_PIPE_MULTISTAGE_VIEW_PARALLEL.
 *
 * Revision 1.1.1.1  2000/02/28 21:55:30  wes
 * OpenRM 1.2 Release
 *
 * Revision 1.8  2000/02/28 17:21:55  wes
 * RM 1.2, pre-OpenRM
 */

#include <stdio.h>

#ifndef RM_WIN
#include <unistd.h>
#endif

#include <rm/rm.h>
#include "libdio.h"


/* image handling stuff */

/*
 * use a static pointer here to verify that the copy/don't copy stuff
 * for images works ok.
 */
static  unsigned char *static_image_data=NULL;
static  float *static_image_data_float=NULL;

#define IMAGE_TYPE RM_UNSIGNED_BYTE 
/*#define IMAGE_TYPE RM_FLOAT */

#define BSIZE (8192 << 2)

void
dioImageDataFreeFunc(void *d)
{
    free(d);
}

RMimage
*dioReadAVSImage(char *fname)
{
    RMimage *t;
    int w,h;
    int i,j,indx, indx2;
    FILE *f;

    f = fopen(fname,"r");
    if (f == NULL)
    {
	fprintf(stderr," can't open the AVS image file named <%s> \n",fname);
	return(NULL);
    }

#if 0
    /* works in unix, gives bogus results in windows */
    fread((void *)&w,4,1,f);
    fread((void *)&h,4,1,f);
#endif
    {
        unsigned char buf[4];
	read(fileno(f),(void *)buf,4);
	w = buf[3] | (buf[2] << 8);
	read(fileno(f),(void *)buf,4);
	h = buf[3] | (buf[2] << 8);
    }

    /* sanity check w,h */

    if ((w < 0) || (w > 4096) ||
	(h < 0) || (h > 4096))
    {
	fprintf(stderr," dubious image dimensions: %d, %d, aborting image read. \n",w,h);
	fclose(f);
	return(NULL);
    }

    static_image_data = (unsigned char *)malloc(sizeof(unsigned char)*w*h*4);

    /* this reads everything at once */
    read(fileno(f),(void *)static_image_data,w*h*4);

    fclose(f);

#define ONE_OVER_255 0.0039215686
    
    /* hack to change from ARGB to RGBA */
    indx = 0;
    indx2 = 3;
    for (i=0;i<w*h;i++)
    {
        unsigned char t;

	t = static_image_data[indx];
	for (j=0;j<3;j++)
	    static_image_data[indx+j] = static_image_data[indx+j+1];
	static_image_data[indx+3] = 0xFF;

	indx += 4;
    }

    if (IMAGE_TYPE == RM_FLOAT)
    {
	static_image_data_float = (float *)malloc(sizeof(float)*w*h*4);
	for (i=0;i<w*h*4;i++)
	    static_image_data_float[i] = (float)(static_image_data[i])*
		ONE_OVER_255;
    }


#if 0
    /* try RGBA, but it's wrong because our bytes are in ARGB  */
    /* scrunch the data from argb to rgb */
    hack_data_ptr = (unsigned char *)malloc(sizeof(unsigned char)*w*h*3);
    /*    hack_data_ptr = data;*/
    indx = 0;
    indx2 = 1;
    for (i=0;i<w*h;i++)
    {
	for (j=0;j<3;j++)
	{
	    hack_data_ptr[indx] = data[indx2];
	    indx++;
	    indx2++;
	}
	indx2++;
    }
#endif

    /* 
     * problem in win-32 with RM_COPY_DATA: not sure what's up, top
     * part of image gets chewed off.
     *
     * problem w/RGB mode in glDrawPixels . it looks like it wants all
     * the image data to be 4-byte word aligned.
     */


    t = rmImageNew(2,w,h,1,RM_IMAGE_RGBA,
		   IMAGE_TYPE,
		   RM_COPY_DATA);  
/*    t = rmImageNew(2,w,h,1,RM_IMAGE_RGBA,RM_UNSIGNED_BYTE,RM_COPY_DATA);  */
/*    t = rmImageNew(2,w,h,1,RM_IMAGE_RGBA,RM_UNSIGNED_BYTE,
      RM_DONT_COPY_DATA); */

    if (IMAGE_TYPE == RM_FLOAT)
    {
	rmImageSetPixelData(t, (void *)static_image_data_float,
			    RM_COPY_DATA,
			    NULL);
    }
    else
    {
	rmImageSetPixelData(t, (void *)static_image_data,
			    RM_COPY_DATA,
			    NULL);
    }

    rmImageMirror(t,RM_IMAGE_MIRROR_HEIGHT);

    /*    free((void *)data);*/
    
    return(t);
}

int
dioWriteAVSImage(const RMimage *img,
		 char *fname)
{
    int w,h;
    unsigned char *data;
    FILE *f;
    int indx,i;

    f = fopen(fname,"w");
    if (f == NULL)
    {
	fprintf(stderr," can't open the AVS image file named <%s> \n",fname);
	return(0);
    }

    rmImageGetImageSize(img,NULL,&w,&h,NULL,NULL,NULL);

    {
        unsigned char buf[4];

	buf[0] = buf[1] = 0;
	buf[2] = (w & 0x0ff00) >> 8;
	buf[3] = w & 0xff;
	
	write(fileno(f),(void *)buf,1);
	write(fileno(f),(void *)(buf+1),1);
	write(fileno(f),(void *)(buf+2),1);
	write(fileno(f),(void *)(buf+3),1);
	
	buf[0] = buf[1] = 0;
	buf[3] = h & 0xff;
	buf[2] = (h & 0x0ff00) >> 8;
	write(fileno(f),(void *)buf,1);
	write(fileno(f),(void *)(buf+1),1);
	write(fileno(f),(void *)(buf+2),1);
	write(fileno(f),(void *)(buf+3),1);

    }

    data = rmImageGetPixelData(img);

    /* assume image data is rgb */
    
    /* expand from rgb to argb */
    
    indx = 0;
    for (i=0;i<w*h;i++)
    {
	unsigned char d[4];
	
	d[0] = 0;
	
	memcpy(d+1,data+indx,sizeof(char)*3);
	fwrite(d,sizeof(unsigned char),4,f);

	indx += 3;
    }

    fclose(f);
    return(1);
}

dioDataObject
*dioDataObjectNew()
{
    dioDataObject *t;
    t = (dioDataObject *)malloc(sizeof(dioDataObject));
    memset(t,0,sizeof(dioDataObject));

    t->datamin = RM_MAXFLOAT;
    t->datamax = RM_MINFLOAT;

    return(t);
}


/* regular data handling stuff */

static char dioWSDelims[]={" \t"};
static char dioWSNDelims[]={" \t\n"};


int
private_dioReadBinaryUbyteSource(float *d, /* dest array */
				 int npts, /* total number of grid points */
				 char *fname, /* filename to read */
				 int normalize, /* 1 or 0 */
				 int skip_bytes) /* at start */
{
    /*
     * reads raw byte data from a file with one fread.
     * each input byte is converted to a float by:
     * if normalize == 1:
     *    bytes in range 0..255 are normalized to 0..1
     * if normalize == 0:
     *    bytes in range 0..255 map to floats in range 0..255
     */
    FILE *f;
    unsigned char *buf;
    int i;
    float divisor;
    unsigned char junk;
    
    f = fopen(fname,"r");
    
    if (f == NULL)
    {
	fprintf(stderr,"dio error: unable to open binary file named <%s> \n",fname);
	exit(1);
    }

    buf = (unsigned char *)malloc(sizeof(unsigned char)*npts);
    memset(buf,0,sizeof(unsigned char)*npts);

    /*
     * read raw byte data
     */
    for (i=0;i<skip_bytes;i++)
	fread((void *)&junk,sizeof(unsigned char),1,f);
    
    fread((void *)buf,sizeof(unsigned char),npts,f);

    fclose(f);

    if (normalize == 0)
	divisor = 1.0F;
    else
	divisor = (float)(1.0/255.0);

    for (i=0;i<npts;i++)
	d[i] = (float)(buf[i])*divisor;

    free((void *)buf);

    return(RM_CHILL);
}

int
private_dioReadFloats(float *d,
		      int npts,
		      FILE *f,
		      int *line_no,
		      char *buf)
{
    char *c;
    int k;
    
    c = strtok(NULL,dioWSDelims);
    for (k=0;k<npts;)
    {
	while (c != NULL)
	{
	    sscanf(c,"%f",d+k);
	    c = strtok(NULL,dioWSNDelims);
	    k++;
	}
	fgets(buf,BSIZE,f);
	*line_no += 1;
	c = strtok(buf,dioWSDelims);
    }
    return(k);
}

		    
dioDataObject
*dioReadDataObject(char *fname)
{
    char buf[BSIZE];
    dioDataObject *t=NULL;
    int line_no=0;
    int skip_bytes = 0;
    
    FILE *f;
    
    f = fopen(fname,"r");
    if (f == NULL)
    {
	fprintf(stderr,"dioReadDataObject(): error opening input file <%s> \n", fname);
	return(NULL);
    }

    t = dioDataObjectNew();

    fgets(buf,BSIZE,f);
    line_no++;
    while (!feof(f))
    {
	if (buf[0] == '#')
	{
	    /* a comment line, skip it*/
	    fgets(buf,BSIZE,f);

	    line_no++;
	}
	else			/* look for a keyword */
	{
	    char *c;
	    c = strtok(buf,dioWSDelims);
	    
	    if (strcmp(c,"width") == 0)
	    {
		c = strtok(NULL,dioWSDelims);
		sscanf(c,"%d",&(t->width));
	    }
	    else if (strcmp(c,"height") == 0)
	    {
		c = strtok(NULL,dioWSDelims);
		sscanf(c,"%d",&(t->height));
	    }
	    else if (strcmp(c,"depth") == 0)
	    {
		c = strtok(NULL,dioWSDelims);
		sscanf(c,"%d",&(t->depth));
	    }
	    else if (strcmp(c,"veclen") == 0)
	    {
		c = strtok(NULL,dioWSDelims);
		sscanf(c,"%d",&(t->veclen));
	    }
	    else if (strcmp(c,"binary-data-skip") == 0)
	    {
		c = strtok(NULL,dioWSDelims);
		sscanf(c,"%d",&skip_bytes);
	    }
	    else if (strncmp(c,"binary-data-file-ubyte-source",strlen("binary-data-file-ubyte-source")) == 0)
	    {
		int npts;
		c = strtok(NULL, dioWSDelims);
		
		npts = t->width * t->height * t->depth * t->veclen;
		t->rawdata = (float *)malloc(sizeof(float)*npts);
		memset(t->rawdata,0,sizeof(float));

		/* read raw ubytes, and convert to floats */
		private_dioReadBinaryUbyteSource(t->rawdata,npts,c,1,skip_bytes);

	    }
	    else if (strncmp(c,"data",strlen("data")) == 0)
	    {
				/* read the data. */
		int npts,nread;

		npts = t->width * t->height * t->depth * t->veclen;
		t->rawdata = (float *)malloc(sizeof(float)*npts);
		memset(t->rawdata,0,sizeof(float));

		nread = private_dioReadFloats(t->rawdata,npts,f,&line_no,buf);

		continue;
	    }
	    else if (strncmp(c,"grid-corners",strlen("grid-corners")) == 0)
	    {

		t->corners = rmVertex3DNew(2);
		private_dioReadFloats(&(t->corners->x),6,f,&line_no,buf);

		continue;
	    }
	    else if (strncmp(c,"grid-xy",strlen("grid-xy")) == 0)
	    {
		float *junk;
		int npts;
		int i,indx;
		
		npts = t->width * t->height * t->depth;
		junk = (float *)malloc(sizeof(float)*npts*2);
		private_dioReadFloats(junk,npts*2,f,&line_no,buf);

		t->xcoords = (float *)malloc(sizeof(float)*npts);
		t->ycoords = (float *)malloc(sizeof(float)*npts);

		indx = 0;
		for (i=0;i<npts;i++)
		{
		    t->xcoords[i] = junk[indx++];
		    t->ycoords[i] = junk[indx++];
		}

		free((void *)junk);
		
	    }
	    else
	    {
		fprintf(stderr," unrecognized keyword on line %d: %s \n",line_no,c);
	    }
	    
	    fgets(buf,BSIZE,f);
	    line_no++;
	} /* not a comment */
	
	    
    } /* while !feof */

    fclose(f);
    return(t);
}


void
dioCompute1DGrid(float *xcoords,
		 float x,
		 float dx,
		 float *ycoords,
		 float y,
		 float dy,
		 float *zcoords,
		 float z,
		 float dz,
		 int npts)
{
    int i;

    for (i=0;i<npts;i++)
    {
	xcoords[i] = x;  x += dx;
	ycoords[i] = y;  y += dy;
	zcoords[i] = z;  z += dz;
    }
}

void
dioCompute2DGrid(float *xcoords,
		 float x,
		 float dx,
		 float *ycoords,
		 float y,
		 float dy,
		 float *zcoords,
		 float z,
		 float dz,
		 int *dims)
{
    int i,j;
    int usize,vsize;

    float *fastptr, fast, dfast, savefast;
    float *slowptr, slow, dslow;
    float *nochangeptr, nochange;

    if (dims[0] != 1)		/* x varies fastest */
    {
	fastptr = xcoords;
	fast = savefast = x;
	dfast = dx;
	usize = dims[0];

	if (dims[1] != 1)	/* y varies next fastest */
	{
	    slowptr = ycoords;
	    slow = y;
	    dslow = dy;
	    vsize = dims[1];

	    nochangeptr = zcoords;
	    nochange = z;
	}
	else			/* assume zsize != 1, y is constant */
	{
	    slowptr = zcoords;
	    slow = z;
	    dslow = dz;

	    vsize = dims[2];

	    nochangeptr = ycoords;
	    nochange = y;
	}
    }
    else			/* xsize == 1, we'll assume the other
				 two dimensions are != 1*/
    {
	usize = dims[1];
	fastptr = ycoords;
	fast = savefast = y;
	dfast = dy;

	vsize = dims[2];
	slowptr = zcoords;
	slow = z;
	dslow = dz;

	nochangeptr = xcoords;
	nochange = x;
    }

    for (j=0;j<vsize;j++)
    {
	fast = savefast;
	for (i=0;i<usize;i++)
	{
	    *fastptr++ = fast;   fast += dfast;
	    *slowptr++ = slow;
	    *nochangeptr++ = nochange;
	}
	slow += dslow;
    }
}

void
dioCompute3DGrid(float *xcoords,
		 float x,
		 float dx,
		 float *ycoords,
		 float y,
		 float dy,
		 float *zcoords,
		 float z,
		 float dz,
		 int *dims)
{
    /*
     * assume:
     * 1. x varies fastest, y varies second fastest, z varies slowest
     * 2. dims[0] > 0, dims[1] > 0, dims[2] > 0
     */
    int i,j,k;
    float tx,ty,tz;
    int indx=0;

    tz = z;
    for (k=0;k<dims[2];k++, tz+=dz)
    {
	ty = y;
	for (j=0;j<dims[1];j++, ty+=dy)
	{
	    tx = x;
	    for (i=0;i<dims[0];i++, tx+=dx)
	    {
		xcoords[indx] = tx;
		ycoords[indx] = ty;
		zcoords[indx] = tz;
		indx++;
	    }
	}
    }
}

void
dioObjectConditioner(dioDataObject *d)
{
    int npts;
    /*
     * this routine "conditions" a dioDataObject for use by rmv vis
     * routines.
     *
     * the RMV routines expect separate arrays for the x,y,z coords.
     */

    npts = d->width * d->height * d->depth;
    if ((d->xcoords == NULL) && (d->ycoords == NULL) && (d->zcoords == NULL))
    {
	/* no input grid specified. create a default grid. */
	int i;
	float x,dx,y,dy,z,dz;
	int ndim=0;
	int dims[3]={1,1,1};

	dims[0] = d->width;
	dims[1] = d->height;
	dims[2] = d->depth;

	for (i=0;i<3;i++)
	    d->dims[i] = 0;
	
	if (d->width > 1)
	{
	    d->dims[ndim] = d->width;
	    ndim++;
	}

	if (d->height > 1)
	{
	    d->dims[ndim] = d->height;
	    ndim++;
	}
	
	if (d->depth > 1)
        {
	    d->dims[ndim] = d->depth;
	    ndim++;
	}

	if (d->corners == NULL)
        {
	    x = 0.0;
	    if (d->width != 1)
		dx = 1.0;
	    else
		dx = 0.0;
	}
	else
        {
	    x = d->corners[0].x;
	    if (d->width == 1)
		dx = 0.0;
	    else
		dx = (d->corners[1].x - d->corners[0].x)/(d->width-1);
	}

	if (d->corners == NULL)
	{
	    y = 0.0;
	    if (d->height != 1)
		dy = 1.0;
	    else
		dy = 0.0;
	}
	else
	{
	    y = d->corners[0].y;
	    if (d->height == 1)
		dy = 0.0;
	    else
		dy = (d->corners[1].y - d->corners[0].y)/(d->height-1);
	}

	if (d->corners == NULL)
	{
	    z = 0.0;
	    if (d->depth != 1)
		dz = 1.0;
	    else
		dz = 0.;
	}
	else
	{
	    z = d->corners[0].z;
	    if (d->depth == 1)
		dz = 0.0;
	    else
		dz = (d->corners[1].z - d->corners[0].z)/(d->depth-1);
	}

	d->xcoords = (float *)malloc(sizeof(float)*npts);
	d->ycoords = (float *)malloc(sizeof(float)*npts);
	d->zcoords = (float *)malloc(sizeof(float)*npts);

	if (ndim == 1)
	    dioCompute1DGrid(d->xcoords,x,dx,d->ycoords,y,dy,d->zcoords,z,dz,npts);
	else if (ndim == 2)
	    dioCompute2DGrid(d->xcoords,x,dx,d->ycoords,y,dy,d->zcoords,z,dz,dims);
	else			/* assume ndim ==3  */
	    dioCompute3DGrid(d->xcoords,x,dx,
			     d->ycoords,y,dy,
			     d->zcoords,z,dz,
			     dims);
	
    } /* if no input coords defined */

    else			/* there are some input coords, compute
				 corners from whatever's there */
    {
	int i;
	d->corners = rmVertex3DNew(2);

	if (d->xcoords)
	{
	    d->corners[0].x = RM_MAXFLOAT;
	    d->corners[1].x = RM_MINFLOAT;
	}
	else
	    d->corners[0].x = d->corners[1].x = 0.0F;
	
	if (d->ycoords)
        {
	    d->corners[0].y = RM_MAXFLOAT;
	    d->corners[1].y = RM_MINFLOAT;
	}
	else
	    d->corners[0].y = d->corners[1].y = 0.0F;

	if (d->zcoords)
	{
	    d->corners[0].z = RM_MAXFLOAT;
	    d->corners[1].z = RM_MINFLOAT;
	}
	else
	    d->corners[0].z = d->corners[1].z = 0.0F;
	    
	for (i=0;i<npts;i++)
	{
	    if (d->xcoords)
	    {
		if (d->xcoords[i] < d->corners[0].x)
		    d->corners[0].x = d->xcoords[i];
		if (d->xcoords[i] > d->corners[1].x)
		    d->corners[1].x = d->xcoords[i];
	    }
	    
	    if (d->ycoords)
	    {
		if (d->ycoords[i] < d->corners[0].y)
		    d->corners[0].y = d->ycoords[i];
		if (d->ycoords[i] > d->corners[1].y)
		    d->corners[1].y = d->ycoords[i];
	    }
	    
	    if (d->zcoords)
	    {
		if (d->zcoords[i] < d->corners[0].z)
		    d->corners[0].z = d->zcoords[i];
		if (d->zcoords[i] > d->corners[1].z)
		    d->corners[1].z = d->zcoords[i];
	    }
	}
    }

    if (d->datamin == RM_MAXFLOAT) /* find min in data */
    {
	int i;
	for (i=0;i<npts;i++)
	    if (d->rawdata[i] < d->datamin)
		d->datamin = d->rawdata[i];
    }

    if (d->datamax == RM_MINFLOAT) /* find max in data */
    {
	int i;
	for (i=0;i<npts;i++)
	    if (d->rawdata[i] > d->datamax)
		d->datamax = d->rawdata[i];
    }

}

void
dioObjectScaleData(dioDataObject *d,
		   float s)
{
    int i,npts;
    float *data;
    
    npts = d->width * d->height * d->depth * d->veclen;
    data = d->rawdata;
    
    for (i=0;i<npts;i++)
	data[i] = data[i] * s;
}

void
dioDeleteDataObject(dioDataObject *d)
{
    if (d == NULL)
	return;

    free((void *)(d->rawdata));
    free((void *)(d->xcoords));
    free((void *)(d->ycoords));
    free((void *)(d->zcoords));
    free((void *)(d->corners));

    free((void *)d);
}