/*  carray.c  28 April 1999 */


/* 
  Tools for the processing of Cell Arrays 
*/

/* Copyright Phil Andrews, Pittsburgh Supercomputing Center */
/* all rights reserved */

#include "defs.h"
#include <stdio.h>

/* a few globals */

#define byte_size 8

unsigned char *allocate_mem();

/* need some macros to pick up encoded values */

/* apologise for the "magic numbers" in this macro, trying to optimise */

static int temp;

/* grab an unsigned  cgm integer at arbitrary (legal) precision */

#define mcr_aint(ptr, precision, out, done) switch (precision){\
case 32: out = (((unsigned int)ptr[0]) << 24) + (((unsigned int)ptr[1]) << 16)\
 + (((unsigned int)ptr[2]) << 8) + ptr[3]; ptr += 4; break;\
case 24:out = (((unsigned int)ptr[0]) << 16) + (((unsigned int)ptr[1]) << 8) +\
 ptr[2]; ptr += 3; break;\
case 16:out = (((unsigned int)ptr[0]) << 8) + ptr[1]; ptr += 2; break;\
case 8: out = *ptr++; break;\
case 4:if (!(temp=(((done + 1) * precision) % byte_size)))out = *ptr++ & 15;\
else out = (*ptr >> temp) & 15; break;\
case 2:if (!(temp=(((done + 1) * precision) % byte_size)))out = *ptr++ & 3;\
else out = (*ptr >> (byte_size - temp)) & 3; break;\
case 1:if (!(temp=(((done + 1) * precision) % byte_size)))out = *ptr++ & 1;\
else out = (*ptr >> (byte_size - temp)) & 1; break;}

/******************************************************************************/

/*  PROCEDURE DECLARATIONS  */

/******************************************************************************/

int c_expand ( unsigned char *dat_ptr, unsigned int nx, unsigned int ny, 
  unsigned row_size, unsigned int c_size, unsigned char *new_ptr );

static do_line ( unsigned char *inptr, int nxin, int in_rbytes, int no_rows, 
  int cpin, int csin, float *inctab, int in_floats, unsigned char *outptr, 
  int nxout, int out_rbytes, int cpout, int csout, float *outctab, 
  int out_floats, int *my_inds, struct c_v_struct *cvext, 
  unsigned int *int_array );

/******************************************************************************/

int c_expand ( unsigned char *dat_ptr, unsigned int nx, unsigned int ny, 
  unsigned row_size, unsigned int c_size, unsigned char *new_ptr ) {

/******************************************************************************/

/* 

  Routine to expand out the compressed raster and dump in supplied memory.

  Returns the number of bytes the original pointer shouldbe stepped forward.

  Run length encoding consists of each row being represented by a 
   series of list items, each list item consists of 2-byte integer followed 
   by a colour value (r,g,b or index). Each row begins on a word boundary.
*/

  int i, j, no_repeats, no_bits, bits_out;
  unsigned char *row_ptr, *old_ptr, *start_ptr;
  void c_dup();
  unsigned int count, bit_offset, bit_out;
  static unsigned char bit1_mask[byte_size] = {255, 127, 63, 31, 15, 7, 3, 1};
  
  /* initialise pointers */
  row_ptr = new_ptr;
  start_ptr = dat_ptr;

/* loop thru each requested row */
  for (i=0; i<ny; ++i) {   /* each row */
      j = 0;    /* no of bytes written in this row */
      new_ptr = row_ptr;  /* start of a new output row */
      old_ptr = dat_ptr;  /* start of a new input row */
      no_repeats = 0;
      bit_offset = 0;
      bit_out = 0;
      while (j < nx) {
    if (!bit_offset) {
        count = (*dat_ptr << byte_size) + *(dat_ptr + 1);
    } else {
        count = ((*dat_ptr & bit1_mask[bit_offset])
      << (bit_offset + byte_size)) +
      (*(dat_ptr + 1) << bit_offset) +
      (*(dat_ptr + 2) >> (byte_size - bit_offset));
    }
    if ((j + count) > nx) {  /* trouble */
        fprintf(stderr, "too many pixels/row ! %d > %d\n",
      count + j, nx);
        count = nx - j;
        return(0);
    }
    dat_ptr += 2;
    /* now need to replicate the next c_size bits count times */
    c_dup(new_ptr, bit_out, dat_ptr,  bit_offset, c_size, count);
    j += count;
    ++no_repeats;
    no_bits = no_repeats * (16 + c_size);
    dat_ptr = old_ptr + no_bits / byte_size;
    bit_offset = no_bits % byte_size;
    bits_out = j * c_size;
    new_ptr = row_ptr + bits_out / byte_size;
    bit_out = bits_out % byte_size;
      }
      /* now see if we have to skip a padded row end */
      if (bit_offset) {  /* incomplete byte */
    if (((no_bits + byte_size - 1) / byte_size) % 2) dat_ptr += 2;
    else dat_ptr +=1;
      } else {
    if (((no_bits + byte_size - 1) / byte_size) % 2) dat_ptr += 1;
      }
      row_ptr += row_size;
  }
  return(dat_ptr - start_ptr);
}
/* take the next size bits from in_ptr plus bit_offset and replicate them 
   count times beginning at out_ptr plus bit_out */
/* the posible values for size are (1, 2, 4, 8, 16, 24, 32) or 3 * those
   values and cannot change within a row, so the largest non-byte-oriented
   value occupies 12 bits, but could easily expand program to handle more
   simply by allowing multiple middle bytes */
void c_dup(out_ptr, bit_out, in_ptr, bit_offset, size, count)
unsigned int bit_offset, size, count, bit_out;
unsigned char *out_ptr, *in_ptr;
{
int i, j, byte_count, no_bytes, bits_left;
unsigned int inval, outval;
static unsigned char bit1_mask[byte_size] = {255, 127, 63, 31, 15, 7, 3, 1};
static unsigned char bit2_mask[byte_size] = 
  {128, 192, 224, 240, 248, 252, 254, 255};

  /* do the easy one first, everything byte oriented */
  if (!((bit_offset) || (size % byte_size) || (bit_out))) {
      byte_count = size / byte_size;
      for (j=0; j<count; ++j)
    for (i=0; i<byte_count; ++i)
        *(out_ptr++) = in_ptr[i];
      return;
  }
  /* now have to deal with tougher stuff, cannot be more than 12 bits */

  if (size > 12) (void) fprintf(stderr, "big size ! = %d\n", size);

  /* get the first byte's worth */
  bits_left = size - (byte_size - bit_offset);  /* what will be left */
  if (bits_left >= 0) 
      inval = (*in_ptr & bit1_mask[bit_offset]) << bits_left;
  else inval = (*in_ptr & bit1_mask[bit_offset]) >> - bits_left;

  /* is there a middle byte ? */
  if (bits_left >= byte_size) {
      ++in_ptr;
      inval += *in_ptr << (bits_left - byte_size);
      bits_left -= byte_size;
  }
  
  /* now the final byte, if necessary */
  if (bits_left > 0) {
      ++in_ptr;
      inval += 
    *in_ptr & bit2_mask[bits_left - 1];
  }
/* now put it out */

  for (j=0; j<count; ++j) {
      bits_left = size - (byte_size - bit_out);  /* what will be left */
      /* put the first byte's worth */
      if (bits_left >= 0) {
        *out_ptr |= (255 & (inval >> bits_left));
    ++out_ptr;
      }
      else *out_ptr |= (255 & (inval << - bits_left));

      /* is there a middle byte ? */
      if (bits_left >= byte_size) {
    *out_ptr = 255 & (inval >> (bits_left - byte_size));
        ++out_ptr;
    bits_left -= byte_size;
      }
  
      /* now the final byte, if necessary */
      if (bits_left > 0) {
    *out_ptr |= (255 & (inval << (byte_size - bits_left)));
      }
      bit_out = (bit_out + size) % byte_size;
  }
  return;
}
/* routine to convert a cell array to a grey-scale array */
/* first devices like postscript that can represent grey levels/colours */
/* if components is 1 then grey only */
/* macro to convert from r, g, b to grey_scale, NTSC standard */
#define col_grey(r, g, b)  (1.0 - (0.3 * r + 0.59 * g + 0.11 * b))

grey_array(dat_ptr, nx, ny, new_ptr, c_s_mode, col_prec,
  ctab, bits_sample, row_size, new_size, components)
unsigned char *dat_ptr;      /* input data */
float *ctab;       /* colour table */
unsigned char *new_ptr;      /* output grey values */
int nx, ny, col_prec, bits_sample, row_size, new_size, components;
enum cs_enum c_s_mode;  
{
unsigned char *start_ptr;
unsigned int index;
int i, j, ir, ig, ib, max_val;
float  r, g, b, rmax_col;

  if (bits_sample != 8) 
      fprintf(stderr, "only 8 bits/sample implemented\n");
  max_val = (1 << bits_sample) - 1;
  start_ptr = dat_ptr;
  switch (c_s_mode) {
case d_c_mode:  /* direct colour */
      rmax_col = (1 << col_prec) - 1;
      for (i=0; i < ny; ++i){
    dat_ptr = start_ptr;
    for (j=0; j < nx; ++j){
        mcr_aint(dat_ptr, col_prec, ir, j)
        r = ir / rmax_col;
        mcr_aint(dat_ptr, col_prec, ig, j)
        g = ig / rmax_col;
        mcr_aint(dat_ptr, col_prec, ib, j)
        b = ib / rmax_col;
        if (components == 1) {
                  new_ptr[j] = col_grey(r, g, b) * max_val;
                    }
                    else if (components == 3) {
                        new_ptr[j * components] = (1 - r) * max_val;
            new_ptr[j * components + 1] = (1 - g) * max_val;
            new_ptr[j * components + 2] = (1 - b) * max_val;
        }
    }   
    start_ptr += row_size;
    new_ptr += new_size;
      }
      break;
case i_c_mode:  /* indexed colour */
      for (i=0; i < ny; ++i){
    dat_ptr = start_ptr;
        for (j=0; j < nx; ++j){
        mcr_aint(dat_ptr, col_prec, index, j);
        r = *(ctab + index * 3);
        g = *(ctab + index * 3 + 1);
        b = *(ctab + index * 3 + 2);
        if (components == 1) {
                  new_ptr[j] = col_grey(r, g, b) * max_val;
                    }
                    else if (components == 3) {
                        new_ptr[j * components] = (1 - r) * max_val;
            new_ptr[j * components + 1] = (1 - g) * max_val;
            new_ptr[j * components + 2] = (1 - b) * max_val;
        }
    }
    start_ptr += row_size;
    new_ptr += new_size;
      }
      break;
  }
  return(1);
}

/* routine to convert a arbitrary cell array to an indexed, precision 1
   cellarray, i.e., a boolean array for monochome devices who cannot
   handle grey scales (e.g., simple tektronix) The output format will be
   identical to that in cell array */
bit_array(dat_ptr, nx, ny, row_size, new_ptr, c_s_mode, col_prec,
  ctab, new_size)
unsigned char *dat_ptr;      /* input data */
float *ctab;         /* colour table */
unsigned char *new_ptr;      /* output bit values */
int nx, ny, row_size, col_prec, new_size;
enum cs_enum c_s_mode;  
{
int i, j, ir, ig, ib, flag;
unsigned char *start_ptr;
unsigned int index;
float  r, g, b, rmax_col;
  start_ptr = dat_ptr;
  switch (c_s_mode) {
case d_c_mode:  /* direct colour */
      rmax_col = (1 << col_prec) - 1;
      for (i=0; i < ny; ++i){
    dat_ptr = start_ptr;
    for (j=0; j < nx; ++j){
        mcr_aint(dat_ptr, col_prec, ir, j)
        r = ir / rmax_col;
        mcr_aint(dat_ptr, col_prec, ig, j)
        g = ig / rmax_col;
        mcr_aint(dat_ptr, col_prec, ib, j)
        b = ib / rmax_col;
        if (col_grey(r, g, b) < 0.5)  /* is it on ? */
          *(new_ptr + (j / byte_size)) |= 
      (1 << (byte_size - 1 - (j % byte_size)));
    }
    new_ptr += new_size;
    start_ptr += row_size;
      }
      break;
case i_c_mode:  /* indexed colour */
      for (i=0; i < ny; ++i){
    dat_ptr = start_ptr;
        for (j=0; j < nx; ++j){
        mcr_aint(dat_ptr, col_prec, index, j);
        r = *(ctab + index * 3);
        g = *(ctab + index * 3 + 1);
        b = *(ctab + index * 3 + 2);
        if (col_grey(r, g, b) < 0.5)  /* is it on ? */
          *(new_ptr + (j / byte_size)) |= 
      (1 << (byte_size - 1 - (j % byte_size)));
    }
    new_ptr += new_size;
    start_ptr += row_size;
      }
      break;
  }
  return(1);
}
/* routine to expand a precision 1 cell array into a coloured one */
bit_col_carray(inptr, old_size, x_bits, no_rows, outptr, new_size, loc_prec,
val, old_val)
char *inptr;  /* input ptr */
char *outptr;  /* output ptr */
int x_bits;  /* no. of columns */
int no_rows;  /* no. of rows */
int old_size;  /* bytes per old row */
int new_size;  /* bytes per new row */
int loc_prec;  /* bits of precision */
int val;  /* value to expand to */
int old_val;  /* old index value */
{
int i, j, k;
char *inrow, *outrow;

/* note we are guarranteed that loc_prec >=8 will align on byte boundaries */

  /* check for the quickest case first */
  if ((loc_prec == sizeof(char)) && old_val) {
      for (i=0; i<(no_rows * new_size); ++i) {
    if (*outptr) *outptr = val;
    ++outptr;
      }
      return(1);
  }

  for (i=0; i<no_rows; ++i) {
      inrow = inptr + i * old_size;
      outrow = outptr + i * new_size;
      for (j = 0; j<x_bits; ++j) {
    if (inrow[j / byte_size] & 
        (1 << (byte_size - 1 - (j % byte_size)))) {
        if (loc_prec >= byte_size) 
      for (k=0; k < (loc_prec/byte_size); ++k)
          outrow[k + (j * loc_prec) / byte_size] = val;
        else outrow[(j * loc_prec) / byte_size] |=
        val << (byte_size - loc_prec - ((j*loc_prec) % byte_size));
    }
      }
  }

  return(1);
}
/* general routine to change the precision of a pixel array */
/* we might as well make this work for arbitrary precisions/padding */
/* with the restriction that if the precision is > byte size, then values */
/* start on a byte boundary */
change_prec(inptr, old_size, old_prec, x_pxls, no_rows, outptr, 
new_size, new_prec, new_col, c_s_mode)
char *inptr;  /* input ptr */
char *outptr;  /* output ptr */
int x_pxls;  /* no. of columns */
int no_rows;  /* no. of rows */
int old_size;  /* bytes per old row */
int new_size;  /* bytes per new row */
int old_prec;  /* bits of precision, old */
int new_prec;  /* bits of precision, new */
enum cs_enum c_s_mode;  /* the colour selection mode, indexed or direct */
struct rgbi_struct *new_col;  /* optional new colour */
{
int i, j, k, old_bits, old_bytes, in_byte, bits_left, new_bytes, out_byte,
bits_togo;
char *inrow, *outrow;
struct rgbi_struct old_col;  /* old colour */
static char bit_mask[byte_size] = {0, 128, 192, 224, 240, 248, 252, 254};
#define old_max 256
char old_array[old_max], *old_ptr = old_array;

  if ((new_size * byte_size < x_pxls * new_prec) ||
      (old_size * byte_size < x_pxls * old_prec)) {
      fprintf(stderr, "illegal row size in change_prec\n");
      return(0);
  }  /* note not quite sufficient to catch all illegal cases */

  old_bytes = old_prec / byte_size;
  new_bytes = new_prec / byte_size;
  if (((old_prec + byte_size - 1) / byte_size) > old_max)
      old_ptr = (char *) allocate_mem(old_bytes, 0);
/* now make the conversion */
switch (c_s_mode) {
case i_c_mode:
  for (i=0; i<no_rows; ++i) {
      inrow = inptr + i * old_size;
      outrow = outptr + i * new_size;
      for (j = 0; j < x_pxls; ++j) {
    in_byte = (j * old_prec) / byte_size;
    out_byte = (j * new_prec) / byte_size;
    bits_left = (j * old_prec) % byte_size;

    for (k=0; k < old_bytes; ++k)  /* slightly tricky */
        old_ptr[k] = inrow[k + in_byte];
    old_ptr[k] = inrow[k + in_byte] & bit_mask[bits_left];

    for (k=0; (k < old_bytes) && (k < new_bytes); ++k)
        outrow[out_byte + new_bytes - 1 - k] = 
      old_ptr[old_bytes - 1 - k];

/*
  FIX THIS!
*/

/*
     outrow[k] |= old_ptr[k] <<
*/

        if (new_prec >= byte_size) 
      for (k=0; k < (new_prec/byte_size); ++k)
          outrow[k + (j * new_prec) / byte_size] 
        = new_col->ind;
        else outrow[(j * new_prec) / byte_size] |= new_col->ind
        << (byte_size - new_prec - ((j*new_prec) % byte_size));
      }
  }
  break;
case d_c_mode:
  break;
  }  /* end switch */
  if (old_ptr != old_array) (void) free(old_ptr);
  return(1);
}
/* this is the general cell array modification routine */
/* it takes in a square cell array and writes out a parallelogram */

new_carray(host_name, nxin, nyin, cpin, rmin, rbin, inptr, in_bytes, csin,
    inctab,in_floats, px, py, qx, qy, rx, ry, nxout, nyout, 
    cpout, rmout, rbout, 
    outptr, out_bytes, csout, outctab, out_floats, in_flag, out_flag,
    first_index, my_inds, cvext)

char *host_name;  /* may be distributing this */
/* first the input data */
int nxin;  /* no. of input row pxls */
int nyin;  /* no. of input column pxls */
int cpin;  /* input colour precision */
int rmin;  /* input representation mode */
int rbin;  /* input no of bits per row */
unsigned char *inptr;  /* input data */
int in_bytes;  /* size of input data */
int csin;  /* the colour selection mode, indexed or direct */
float *inctab;  /* input colour table */
int in_floats;  /* size of input colour table (no. of floats) */
/* now the output data */
int px, py, qx, qy, rx, ry;  /* the output parallelogram corners */
int nxout, nyout;   /* output row and column size */
int cpout;  /* output colour precision */
int rmout;  /* output representation mode */
int rbout;  /* output no of bits per row */
unsigned char *outptr;  /* output data */
int out_bytes;  /* size of output data */
int csout;  /* the colour selection mode, indexed or direct */
float *outctab;  /* output colour table */
int out_floats;  /* size of output colour table (no. of floats) */
int in_flag, *out_flag;  /* input and output info flags */
int first_index;  /* first index in the output colour table */
int *my_inds;  /* pivot indices for direct colour table */
struct c_v_struct *cvext;  /* colour value extent */
{
/* local variables */
unsigned char *my_ptr, *row_ptr = NULL;  /* for convenience */
int row_out, row_in, row_bytes, pxl_along, x, y, xs, ys, col_in;
int skip_bytes, c_size;
int index_in, in_entries;
double rin, gin, bin;
int yins, out_rbytes, i;
unsigned int *int_array = NULL;
/* get down to business */

  if (rmout != 1) {  /* only packed list output for now */
      (void) fprintf(stderr, 
    "unimplemented output rep_mode (%d) in new_carray\n", rmout);
      return(0);
  }

  xs = (rx >= px) ? 1 : -1;
  ys = (qy >= py) ? 1 : -1;

  /* how many bytes for an input row ? */
  c_size = (csin == (int) d_c_mode) ? 3 * cpin : cpin;
  row_bytes = (nxin * c_size + byte_size - 1) / byte_size;
  if (row_bytes % 2) ++row_bytes;

  /* how many bytes for an output row ? */
  out_rbytes = (rbout + byte_size - 1) / byte_size;
  if ((out_rbytes * nyout) > out_bytes) {
      (void) fprintf(stderr, 
    "not enough memory in new_carray, %d > %d\n",
    out_rbytes * nyout, out_bytes);
      return(2);
  }
  /* how many input rows per output row ? */
  yins = (nyin > nyout) ? (nyin + nyout - 1) / nyout : 1;
  if (yins < 1) yins = 1;    /* for safety */
  if (yins > nyin) yins = nyin;  /* for safety */

  /* and make sure we have enough memory */  
  if (rmin == 0) {  /* run length encoded */
      /* allocate enough memory for a row */
      if (!(row_ptr = (unsigned char *) 
    allocate_mem(row_bytes * yins, 1))) {
        (void) fprintf(stderr, "no memory for new_carray !\n");
        return(0);
      }
  }
  /* get some memory for our sampling block */
  in_entries = (csin == (int) d_c_mode) ? 3 * yins * nxin :
      yins * nxin;
  /* now get it */
  if (!(int_array = (unsigned int *) 
      allocate_mem(in_entries * sizeof(int), 0))) {
      (void) fprintf(stderr, 
    "couldn't get enough memory (%d) in new_carray\n",
    in_entries * sizeof(int));
      return(2);
  }
/* initialise pointers */
  if (rmin == 0) {  /* run length encoded */
      my_ptr = row_ptr;
      skip_bytes = c_expand(inptr, nxin, yins, row_bytes,
    c_size, my_ptr);
      if (skip_bytes) inptr += skip_bytes;
      else return(0);
  } else my_ptr = inptr;  /* packed list */
  row_in = yins;  /* got some */

  /* loop thru all of the rows */

  if (ys < 0) outptr += (nyout + ys) * out_rbytes; /* going down */

  for (row_out = 1; row_out <= nyout; ++row_out) {
      /* do we need an input block ? */
      if ((row_out * nyin) > (row_in * nyout)) {
    /* may be at last step, need to reduce step */
    if ((nyin - row_in) < yins) yins = nyin - row_in;
    if (yins <= 0) continue;
    if (rmin == 0) {  /* run length encoded */
        my_ptr = row_ptr;
        for (i=0; i<row_bytes * yins; ++i) my_ptr[i] = 0;
        skip_bytes = c_expand(inptr, nxin, yins, row_bytes,
      c_size, my_ptr);
        if (skip_bytes) inptr += skip_bytes;
        else return(0);
    } else my_ptr += row_bytes * yins;
    row_in += yins;
      }

      /* the block we want is in my_ptr now */
      do_line(my_ptr, nxin, row_bytes, yins, cpin, csin, inctab, 
    in_floats,
    outptr, nxout, out_rbytes, cpout, csout, outctab, out_floats,
    my_inds, cvext, int_array);
      outptr += ys * out_rbytes;
  }

  /* free up the memory */
  if (int_array) free(int_array);
  if (row_ptr ) free(row_ptr);

  return(1);
}
/* do one output line */

static do_line ( unsigned char *inptr, int nxin, int in_rbytes, int no_rows, 
  int cpin, int csin, float *inctab, int in_floats, unsigned char *outptr, 
  int nxout, int out_rbytes, int cpout, int csout, float *outctab, 
  int out_floats, int *my_inds, struct c_v_struct *cvext, 
  unsigned int *int_array )
{
/* we have an input block of no_rows rows, with nxin pixels per row */
/* from this we want to produce one output row */
#define R_SIZE 10
static unsigned int yrarray[R_SIZE] = {2348, 947, 894, 4948, 8889, 63790, 449,
    844, 4783, 4940};  /* scientifically produced random integers */
static unsigned int xrarray[R_SIZE] = {2384, 378, 4903, 340, 9389, 9893, 783,
    8895, 3049, 909};
int i, j, k, ir, ig, ib;
unsigned char *my_ptr;
int ysample, xsample, bytes_cell;
unsigned int index;
unsigned int new_index;
unsigned int *int_ptr;
unsigned int itemp;
unsigned int get_index();
float r;
float g;
float b;

  /* no of bytes per output cell */
  bytes_cell = cpout / byte_size;
/* first we'll read in the input block */
  
  /* handle an input indexed colour first */

  if (csin == (int) i_c_mode) {

      /* fill out the entries */

      for (i=0; i<no_rows; ++i) {
        my_ptr = inptr + i * in_rbytes;
        int_ptr = int_array + i * nxin;
        for (j=0; j<nxin; ++j) {
        mcr_aint(my_ptr, cpin, itemp, j);
        int_ptr[j] = (itemp < in_floats) ? itemp :
      itemp % in_floats;
          }
      }
  } else if (csin == (int) d_c_mode) {  /* direct colour mode */
      /* fill out the entries */
      for (i=0; i<no_rows; ++i) {
        my_ptr = inptr + i * in_rbytes;
        int_ptr = int_array + i * nxin * 3;
        for (j=0; j<nxin; ++j) {
        for (k=0; k<3; ++k) {
            mcr_aint(my_ptr, cpin, itemp, (j * 3 + k));
      int_ptr[j * 3 + k] = itemp;
        }
          }
      }
  } else {
      (void) fprintf(stderr, "illegal csin in new_carray, %d\n", csin);
      return(2);
  }


/* now write the output block, four possibilities */
/* the possible values for cpout are (1, 2, 4, 8, 16, 24, 32) */
  
  ysample = 0;        /* input row we will use */
  if (csout == (int) i_c_mode) {
      /* now fill out the entries */
      if (csin == (int) i_c_mode) {
    /* indexed input, indexed output */

    for (i=0; i<nxout; ++i) {
        /* get a single (possibly random) input pixel */
        if (no_rows > 1) ysample =  yrarray[i % R_SIZE] % no_rows;

        if (nxin > nxout) xsample = (i * nxin) / nxout +
      (xrarray[i % R_SIZE] % (nxin / nxout));
        else if (nxout > nxin) xsample = (i * nxin) / nxout;
        else xsample = i;
        index = int_array[ysample * nxin + xsample];      

        /* may need to reduce precision */
        new_index = index;
/*        new_index = get_index(index, inctab, in_floats, 
      outctab, out_floats); */  /* try later */
        new_index &= ~(~0 << cpout);  /* for safety */

        /* now write out this index */
        if (bytes_cell) for(j=0; j<(bytes_cell); ++j){
      outptr[bytes_cell * (i + 1) - 1 - j] |= 
          (new_index & 255);
      new_index >>= byte_size;
        } else outptr[(i * cpout) / byte_size] |= new_index <<
      (byte_size - cpout - ((i * cpout) % byte_size));
        }
      } else 

/* direct input, indexed output */

      if (csin == (int) d_c_mode) {  /* direct input */
    for (i=0; i<nxout; ++i) {
        /* get a single (possibly random) input pixel */
        if (no_rows > 1) ysample =  yrarray[i % R_SIZE] % no_rows;

        if (nxin > nxout) xsample = (i * nxin) / nxout +
      (xrarray[i % R_SIZE] % (nxin / nxout));
        else if (nxout > nxin) xsample = (i * nxin) / nxout;
        else xsample = i;
        /* now get the r, g, b value */
        ir = int_array[ysample * nxin * 3 + xsample * 3];
        ig = int_array[ysample * nxin * 3 + xsample * 3 + 1];
        ib = int_array[ysample * nxin * 3 + xsample * 3 + 2];
        /* and convert it to an index */
        /* first the floating point values */
        r = dcind(0, ir, (*cvext));
        g = dcind(1, ig, (*cvext));
        b = dcind(2, ib, (*cvext));
        /* now the index */
        new_index = int_rgb(r, g, b, my_inds);
        new_index &= ~(~0 << cpout);  /* for safety */

        /* now write out this index */
        if (bytes_cell >= 1) for(j=0; j<bytes_cell; ++j){
      outptr[(i + 1) * bytes_cell - j - 1] |= 
          new_index & 255;
      new_index >>= byte_size;
        } else outptr[(i * cpout) / byte_size] |= new_index <<
      (byte_size - cpout - ((i * cpout) % byte_size));
        }
      }
/* direct output, indexed input */
  } else if (csout == (int) d_c_mode) {
      /* now fill out the entries */
      if (csin == (int) i_c_mode) {
    for (i=0; i<nxout; ++i) {
        }
      } else 
/* direct input, direct output */

      if (csin == (int) d_c_mode) {  /* direct input */
    for (i=0; i<nxout; ++i) {
        }
      }
/* none of these ? */
  } else {
      (void) fprintf(stderr, "illegal csout in new_carray, %d\n", csout);
      return(2);
  }

  return(1);
}
#undef R_SIZE
/* function to return the best index in the 2nd colour table */
/* do a better job later */
unsigned int get_index(index, inctab, in_floats, outctab, out_floats)
unsigned int index;
float *inctab, *outctab;
int in_floats, out_floats;
{
unsigned int ret, i, best_index;
float best_dif, temp_dif;
#define EPS 0.001
#define CDIF(j, k) ((inctab[3*j] - outctab[3*k])* (inctab[3*j] - outctab[3*k])\
 + (inctab[3*j + 1] - outctab[3*k + 1]) * (inctab[3*j + 1] - outctab[3*k + 1])\
 + (inctab[3*j + 2] - outctab[3*k + 2]) * (inctab[3*j + 2] - outctab[3*k + 2]))

  /* first check to see it we have a good match in colour tables */
  if ((index < out_floats) && (CDIF(index, index) < EPS))
      return(index);

  /* now get the closest match */
  best_index = 0;
  best_dif = CDIF(index, best_index);

  for (i=1; i<out_floats; ++i) {
      if ((temp_dif = CDIF(index, i)) < best_dif) {
    best_index = i;
    best_dif = temp_dif;
      }
  }

  return(best_index);
}
#undef EPS