/*
 * xvsmooth.c - smoothing/color dither routines for XV
 * 
 *  Author:    John Bradley, University of Pennsylvania
 *                (bradley@cis.upenn.edu)
 *
 *  Contains:
 *            void  Smooth()
 *            byte *Smooth24()
 *     static int   SmoothX(pic24)
 *     static int   SmoothY(pic24)
 *     static int   SmoothXY(pic24)
 *            void  ColorDither(pic24, w, h)
 */


/*
 * Copyright 1989, 1990, 1991, 1992 by John Bradley and
 *                       The University of Pennsylvania
 *
 * Permission to use, copy, and distribute for non-commercial purposes,
 * is hereby granted without fee, providing that the above copyright
 * notice appear in all copies and that both the copyright notice and this
 * permission notice appear in supporting documentation.
 *
 * The software may be modified for your own purposes, but modified versions
 * may not be distributed.
 *
 * This software is provided "as is" without any expressed or implied warranty.
 *
 * The author may be contacted via:
 *    US Mail:   John Bradley
 *               GRASP Lab, Room 301C
 *               3401 Walnut St.
 *               Philadelphia, PA  19104
 *
 *    Phone:     (215) 898-8813
 *    EMail:     bradley@cis.upenn.edu
 */


#include "xv.h"

#ifdef __STDC__
static int SmoothX(byte *);
static int SmoothY(byte *);
static int SmoothXY(byte *);
#else
static int  SmoothX(), SmoothY(), SmoothXY();
#endif


/***************************************************/
void Smooth()
{
  /* does a SMOOTH resize from cpic into a eWIDE*eHIGH epic.  
     Generates a smooth 24bit eWIDE*eHIGH image, calls ColorDither()
     to map it onto the screen, then throws away 24-bit version */

  byte *pic24;

  pic24 = Smooth24();
  if (pic24) {
    epicmode = EM_SMOOTH;   SetEpicMode();
    ColorDither(pic24, eWIDE, eHIGH);
    free(pic24);
  }
}




/***************************************************/
byte *Smooth24()
{
  /* does a SMOOTH resize from cpic into a eWIDE*eHIGH epic.  
     Generates a smooth 24bit eWIDE*eHIGH image.  Returns NULL if
     it can't (failed malloc, no doubt) */

  byte *pic24, *pp;
  int *cxtab, *pxtab;
  int y1cWIDE, cycWIDE;
  int ex, ey, cx, cy, px, py, apx, apy, x1, y1;
  int cA, cB, cC, cD;
  int pA, pB, pC, pD;
  int retval;

  pp = pic24 = (byte *) malloc(eWIDE * eHIGH * 3);
  if (!pic24) {
    fprintf(stderr,"unable to malloc in 'Smooth24()'\n");
    return NULL;
  }

  retval = 0;

  /* decide which smoothing routine to use based on type of expansion */
  if      (eWIDE <  cWIDE && eHIGH <  cHIGH) retval = SmoothXY(pic24);
  else if (eWIDE <  cWIDE && eHIGH >= cHIGH) retval = SmoothX(pic24);
  else if (eWIDE >= cWIDE && eHIGH <  cHIGH) retval = SmoothY(pic24);
  else {
    /* eWIDE >= cWIDE && eHIGH >= cHIGH */
    /* cx,cy = original pixel in cpic.  px,py = relative position 
       of pixel ex,ey inside of cx,cy as percentages +-50%, +-50%.  
       0,0 = middle of pixel */

    /* we can save a lot of time by precomputing cxtab[] and pxtab[], both
       arrays of eWIDE ints that contain values for the equations:
         cx = (ex * cWIDE) / eWIDE;
         px = ((ex * cWIDE * 100) / eWIDE) - (cx * 100) - 50; */

    cxtab = (int *) malloc(eWIDE * sizeof(int));
    if (!cxtab) { free(pic24);  return NULL; }

    pxtab = (int *) malloc(eWIDE * sizeof(int));
    if (!pxtab) { free(pic24);  free(cxtab);  return NULL; }

    for (ex=0; ex<eWIDE; ex++) {
      cxtab[ex] = (ex * cWIDE) / eWIDE;
      pxtab[ex] = (((ex * cWIDE)* 100) / eWIDE) 
	           - (cxtab[ex] * 100) - 50;
    }
    
    for (ey=0; ey<eHIGH; ey++) {
      cy = (ey * cHIGH) / eHIGH;
      py = (((ey * cHIGH) * 100) / eHIGH) - (cy * 100) - 50;
      if (py<0) { y1 = cy-1;  if (y1<0) y1=0; }
           else { y1 = cy+1;  if (y1>cHIGH-1) y1=cHIGH-1; }

      y1cWIDE = y1 * cWIDE;
      cycWIDE = cy * cWIDE;

      if ((ey&15) == 0) WaitCursor();

      for (ex=0; ex<eWIDE; ex++) {
	cx = cxtab[ex];
	px = pxtab[ex];

	if (px<0) { x1 = cx-1;  if (x1<0) x1=0; }
	     else { x1 = cx+1;  if (x1>cWIDE-1) x1=cWIDE-1; }

	cA = cpic[y1cWIDE + x1];
	cB = cpic[y1cWIDE + cx];
	cC = cpic[cycWIDE + x1];
	cD = cpic[cycWIDE + cx];

	/* quick check */
	if (cA == cB && cB == cC && cC == cD) {
	  /* set this pixel to the same color as in cpic */
	  *pp++ = r[cD];  *pp++ = g[cD];  *pp++ = b[cD];
	}

	else {
	  /* compute weighting factors */
	  apx = abs(px);  apy = abs(py);
	  pA = (apx * apy) / 100;
	  pB = (apy * (100 - apx)) / 100;
	  pC = (apx * (100 - apy)) / 100;
	  pD = 100 - (pA + pB + pC);

	  *pp++ = (pA * r[cA])/100 + (pB * r[cB])/100 + 
	          (pC * r[cC])/100 + (pD * r[cD])/100;

	  *pp++ = (pA * g[cA])/100 + (pB * g[cB])/100 + 
	          (pC * g[cC])/100 + (pD * g[cD])/100;

	  *pp++ = (pA * b[cA])/100 + (pB * b[cB])/100 + 
	          (pC * b[cC])/100 + (pD * b[cD])/100;
	}
      }
    }

  free(cxtab);  
  free(pxtab);
  }

  if (retval) {  /* failed */
    free(pic24);
    pic24 = (byte *) NULL;
  }

  return pic24;
}




/***************************************************/
static int SmoothX(pic24)
byte *pic24;
{
  byte *cptr, *cptr1;
  int  i, j;
  int  *lbufR, *lbufG, *lbufB;
  int  pixR, pixG, pixB;
  int  pcnt0, pcnt1, lastpix, pixcnt, thisline, ypcnt;
  int  *pixarr, *paptr;

  /* returns '0' if okay, '1' if failed (malloc) */

  /* for case where cpic is shrunk horizontally and stretched vertically
     shrinks cpic into an eWIDExeHIGH 24-bit picture.  Only works correctly
     when cWIDE>=eWIDE and cHIGH<=eHIGH */

  lbufR = (int *) calloc(cWIDE, sizeof(int));
  if (!lbufR) return 1;
  lbufG = (int *) calloc(cWIDE, sizeof(int));
  if (!lbufG) { free(lbufR); return 1; }
  lbufB = (int *) calloc(cWIDE, sizeof(int));
  if (!lbufB) { free(lbufG);  free(lbufR); return 1; }

  pixarr = (int *) calloc(cWIDE+1, sizeof(int));
  if (!pixarr) { free(lbufB);  free(lbufG);  free(lbufR); return 1; }
  
  for (j=0; j<=cWIDE; j++) 
    pixarr[j] = (j*eWIDE + (15*cWIDE)/16) / cWIDE;

  cptr = cpic;  cptr1 = cptr + cWIDE;

  for (i=0; i<eHIGH; i++) {
    if ((i&15) == 0) WaitCursor();

    ypcnt = (((i*cHIGH)<<6) / eHIGH) - 32;
    if (ypcnt<0) ypcnt = 0;

    pcnt1 = ypcnt & 0x3f;                     /* 64ths of NEXT line to use */
    pcnt0 = 64 - pcnt1;                       /* 64ths of THIS line to use */

    thisline = ypcnt>>6;

    cptr  = cpic + thisline * cWIDE;
    if (thisline+1 < cHIGH) cptr1 = cptr + cWIDE;
    else cptr1 = cptr;

    for (j=0; j<cWIDE; j++, cptr++, cptr1++) {
      lbufR[j] = ((r[*cptr] * pcnt0) + (r[*cptr1] * pcnt1)) >> 6;
      lbufG[j] = ((g[*cptr] * pcnt0) + (g[*cptr1] * pcnt1)) >> 6;
      lbufB[j] = ((b[*cptr] * pcnt0) + (b[*cptr1] * pcnt1)) >> 6;
    }

    pixR = pixG = pixB = pixcnt = lastpix = 0;

    for (j=0, paptr=pixarr; j<=cWIDE; j++,paptr++) {
      if (*paptr != lastpix) {   /* write a pixel to pic24 */
	*pic24++ = pixR / pixcnt;
	*pic24++ = pixG / pixcnt;
	*pic24++ = pixB / pixcnt;
	lastpix = *paptr;
	pixR = pixG = pixB = pixcnt = 0;
      }

      if (j<cWIDE) {
	pixR += lbufR[j];
	pixG += lbufG[j];
	pixB += lbufB[j];
	pixcnt++;
      }
    }
  }

  free(lbufR);  free(lbufG);  free(lbufB);  free(pixarr);
  return 0;
}

	
      



/***************************************************/
static int SmoothY(pic24)
byte *pic24;
{
  byte *clptr, *cptr, *cptr1;
  int  i, j;
  int  *lbufR, *lbufG, *lbufB, *pct0, *pct1, *cxarr, *cxptr;
  int  lastline, thisline, linecnt;
  int  retval;


  /* returns '0' if okay, '1' if failed (malloc) */

  /* for case where cpic is shrunk vertically and stretched horizontally
     shrinks cpic into an eWIDExeHIGH 24-bit picture.  Only works correctly
     when cWIDE<=eWIDE and cHIGH>=eHIGH */

  retval = 0;   /* no probs, yet... */

  lbufR = lbufG = lbufB = pct0 = pct1 = cxarr = NULL;
  lbufR = (int *) calloc(eWIDE, sizeof(int));
  lbufG = (int *) calloc(eWIDE, sizeof(int));
  lbufB = (int *) calloc(eWIDE, sizeof(int));
  pct0  = (int *) calloc(eWIDE, sizeof(int));
  pct1  = (int *) calloc(eWIDE, sizeof(int));
  cxarr = (int *) calloc(eWIDE, sizeof(int));

  if (!lbufR || !lbufG || !lbufB || !pct0 || ! pct1 || !cxarr)
    { retval = 1;  goto smyexit; }



  for (i=0; i<eWIDE; i++) {                /* precompute some handy tables */
    int cx64;
    cx64 = (((i * cWIDE) << 6) / eWIDE) - 32;
    if (cx64<0) cx64 = 0;
    pct1[i] = cx64 & 0x3f;
    pct0[i] = 64 - pct1[i];
    cxarr[i] = cx64 >> 6;
  }


  lastline = linecnt = 0;

  for (i=0, clptr=cpic; i<=cHIGH; i++, clptr+=cWIDE) {
    if ((i&15) == 0) WaitCursor();

    thisline = (i * eHIGH + (15*cHIGH)/16) / cHIGH;
    if (thisline != lastline) {  /* copy a line to pic24 */

      for (j=0; j<eWIDE; j++) {
	*pic24++ = lbufR[j] / linecnt;
	*pic24++ = lbufG[j] / linecnt;
	*pic24++ = lbufB[j] / linecnt;
      }

      bzero( (char *) lbufR, eWIDE * sizeof(int));  /* clear out line bufs */
      bzero( (char *) lbufG, eWIDE * sizeof(int));
      bzero( (char *) lbufB, eWIDE * sizeof(int));
      linecnt = 0;  lastline = thisline;
    }

    for (j=0, cxptr=cxarr; j<eWIDE; j++, cxptr++) {
      cptr  = clptr + *cxptr;
      if (*cxptr < cWIDE-1) cptr1 = cptr + 1;
                       else cptr1 = cptr;

      lbufR[j] += (((r[*cptr] * pct0[j]) + (r[*cptr1] * pct1[j])) >> 6);
      lbufG[j] += (((g[*cptr] * pct0[j]) + (g[*cptr1] * pct1[j])) >> 6);
      lbufB[j] += (((b[*cptr] * pct0[j]) + (b[*cptr1] * pct1[j])) >> 6);
    }
   
    linecnt++;
  }


 smyexit:
  if (lbufR) free(lbufR);
  if (lbufG) free(lbufG);
  if (lbufB) free(lbufB);
  if (pct0)  free(pct0);
  if (pct1)  free(pct1);
  if (cxarr) free(cxarr);

  return retval;
}

	
      



/***************************************************/
static int SmoothXY(pic24)
byte *pic24;
{
  byte *cptr;
  int  i,j;
  int  *lbufR, *lbufG, *lbufB;
  int  pixR, pixG, pixB;
  int  lastline, thisline, lastpix, linecnt, pixcnt;
  int  *pixarr, *paptr;


  /* returns '0' if okay, '1' if failed (malloc) */

  /* shrinks cpic into an eWIDExeHIGH 24-bit picture.  Only works correctly
     when cWIDE>=eWIDE and cHIGH>=eHIGH (ie, the picture is shrunk on both
     axes) */

  lbufR = (int *) calloc(cWIDE, sizeof(int));
  if (!lbufR) return 1;
  lbufG = (int *) calloc(cWIDE, sizeof(int));
  if (!lbufG) { free(lbufR); return 1; }
  lbufB = (int *) calloc(cWIDE, sizeof(int));
  if (!lbufB) { free(lbufG);  free(lbufR); return 1; }

  pixarr = (int *) calloc(cWIDE+1, sizeof(int));
  if (!pixarr) { free(lbufB);  free(lbufG);  free(lbufR); return 1; }
  
  for (j=0; j<=cWIDE; j++) 
    pixarr[j] = (j*eWIDE + (15*cWIDE)/16) / cWIDE;

  lastline = linecnt = pixR = pixG = pixB = 0;
  cptr = cpic;

  for (i=0; i<=cHIGH; i++) {
    if ((i&15) == 0) WaitCursor();

    thisline = (i * eHIGH + (15*cHIGH)/16 ) / cHIGH;

    if ((thisline != lastline)) {      /* copy a line to pic24 */
      pixR = pixG = pixB = pixcnt = lastpix = 0;

      for (j=0, paptr=pixarr; j<=cWIDE; j++,paptr++) {
	if (*paptr != lastpix) {                 /* write a pixel to pic24 */
	  *pic24++ = (pixR/linecnt) / pixcnt;
	  *pic24++ = (pixG/linecnt) / pixcnt;
	  *pic24++ = (pixB/linecnt) / pixcnt;
	  lastpix = *paptr;
	  pixR = pixG = pixB = pixcnt = 0;
	}

	if (j<cWIDE) {
	  pixR += lbufR[j];
	  pixG += lbufG[j];
	  pixB += lbufB[j];
	  pixcnt++;
	}
      }
      lastline = thisline;
      bzero( (char *) lbufR, cWIDE * sizeof(int));  /* clear out line bufs */
      bzero( (char *) lbufG, cWIDE * sizeof(int));
      bzero( (char *) lbufB, cWIDE * sizeof(int));
      linecnt = 0;
    }

    if (i<cHIGH) {
      for (j=0; j<cWIDE; j++, cptr++) {
	lbufR[j] += r[*cptr];
	lbufG[j] += g[*cptr];
	lbufB[j] += b[*cptr];
      }
      linecnt++;
    }
  }

  free(lbufR);  free(lbufG);  free(lbufB);  free(pixarr);
  return 0;
}

	
      



/********************************************/
void ColorDither(pic24,w,h)
byte *pic24;
int   w,h;
{
  /* takes a 24 bit picture, of size w*h, dithers with the colors that 
     have already been allocated, and sticks it's output in 'epic'.
     Calls CreateXImage() and DrawWindow to display results

     if pic24 is NULL, uses the existing 'epic' (an 8-bit image) as
     the source */

  byte *np, *ep, *newpic; 
  short *cache;
  int r2, g2, b2;
  int *thisline, *nextline, *thisptr, *nextptr, *tmpptr;
  int  i, j, rerr, gerr, berr, pwide3;
  int  imax, jmax;
  XColor ctab[256], *cptr;
  int key;
  long cnt1, cnt2;

  if (epicmode != EM_SMOOTH) {  /* don't set EM_DITH if called from Smooth() */
    epicmode = EM_DITH;   SetEpicMode();
  }


  SetISTR(ISTR_EXPAND, "%.3g x %.3g  (%d x %d)",
	  ((float) w) / cWIDE, ((float) h) / cHIGH, w, h);

  /* if mono and we're given a 24-bit image, monoify it before attempting
     to dither */
  if ((mono || ncols==0) && pic24) {
    for (i=w*h,np=pic24; i; i--, np+=3) {
      j = MONO(np[0],np[1],np[2]);
      np[0] = np[1] = np[2] = j;
    }
  }



  cnt1 = cnt2 = 0;
  cache = (short *) calloc(2<<14, sizeof(short));
  if (!cache) return;

  
  /* load up ctab[] with the colors we'll be dithering with */

  if (theVisual->class == TrueColor || theVisual->class == DirectColor) {  
    /* if truecolor, we have exactly the colors in the pic's colormap */
    for (i=0; i<numcols; i++) {
      ctab[i].red   = r[i];
      ctab[i].green = g[i];
      ctab[i].blue  = b[i];
    }
  }
    
  else if (ncols>0) {
    /* read the colors currently being used in the colormap */
    for (i=0; i<numcols; i++) ctab[i].pixel = cols[i];
    XQueryColors(theDisp, (LocalCmap) ? LocalCmap : theCmap, ctab, numcols);

    /* downshift color defs to 8 bit R,G,B vals instead of 16 bit */
    for (i=0; i<numcols; i++) 
      { ctab[i].red >>= 8;  ctab[i].green >>= 8;  ctab[i].blue >>=8; }
  }

  else {
    /* ncols==0 case */
    for (i=0; i<numcols; i++) {
      ctab[i].red = ctab[i].green = ctab[i].blue = 
	MONO(r[i],g[i],b[i]);
    }
  }

  if (!pic24) { w = eWIDE;  h = eHIGH;  ep = epic; }
  else ep = pic24;

  newpic = (byte *) malloc(w * h);         /* create output 8-bit pic */
  if (!newpic) { free(cache); return; }

  np = newpic;  pwide3 = w*3;  imax = h-1;  jmax = w-1;

  thisline = (int *) malloc(pwide3 * sizeof(int));
  nextline = (int *) malloc(pwide3 * sizeof(int));
  if (!thisline || !nextline)
    FatalError("unable to allocate memory in ColorDither()");

  /* get first line of picture */
  if (!pic24)
    for (j=w, tmpptr=nextline; j; j--, ep++) {
      *tmpptr++ = (int) r[*ep];
      *tmpptr++ = (int) g[*ep];
      *tmpptr++ = (int) b[*ep];
    }
  else
    for (j=pwide3, tmpptr=nextline; j; j--, ep++) *tmpptr++ = (int) *ep;


  for (i=0; i<h; i++) {
    if ((i&15) == 0) WaitCursor();

    tmpptr = thisline;  thisline = nextline;  nextline = tmpptr;   /* swap */

    if (i!=imax) {  /* get next line */
      if (!pic24)
	for (j=w, tmpptr=nextline; j; j--, ep++) {
	  *tmpptr++ = (int) r[*ep];
	  *tmpptr++ = (int) g[*ep];
	  *tmpptr++ = (int) b[*ep];
	}
      else
	for (j=pwide3, tmpptr=nextline; j; j--, ep++) *tmpptr++ = (int) *ep;
    }

    /* dither a line */
    for (j=0, thisptr=thisline, nextptr=nextline; j<w; j++,np++) {
      int k, d, mind, closest;

      r2 = *thisptr++;  g2 = *thisptr++;  b2 = *thisptr++;

      /* map r2,g2,b2 components (could be outside 0..255 range) 
	 into 0..255 range */
      
      if (r2<0 || g2<0 || b2<0) {   /* are there any negatives in RGB? */
	if (r2<g2) { if (r2<b2) k = 0; else k = 2; }
	else { if (g2<b2) k = 1; else k = 2; }

	switch (k) {
	case 0:  g2 -= r2;  b2 -= r2;  d = (abs(r2) * 3) / 2;    /* RED */
	         r2 = 0;  
	         g2 = (g2>d) ? g2 - d : 0;
	         b2 = (b2>d) ? b2 - d : 0;
	         break;

	case 1:  r2 -= g2;  b2 -= g2;  d = (abs(g2) * 3) / 2;    /* GREEN */
	         r2 = (r2>d) ? r2 - d : 0;
	         g2 = 0;  
	         b2 = (b2>d) ? b2 - d : 0;
	         break;

	case 2:  r2 -= b2;  g2 -= b2;  d = (abs(b2) * 3) / 2;    /* BLUE */
	         r2 = (r2>d) ? r2 - d : 0;
	         g2 = (g2>d) ? g2 - d : 0;
	         b2 = 0;  
	         break;
	}
      }

      if (r2>255 || g2>255 || b2>255) {   /* any overflows in RGB */
	if (r2>g2) { if (r2>b2) k = 0; else k = 2; }
              else { if (g2>b2) k = 1; else k = 2; }

	switch (k) {
	case 0:   g2 = (g2*255)/r2;  b2 = (b2*255)/r2;  r2=255;  break;
	case 1:   r2 = (r2*255)/g2;  b2 = (b2*255)/g2;  g2=255;  break;
	case 2:   r2 = (r2*255)/b2;  g2 = (g2*255)/b2;  b2=255;  break;
	}
      }

      key = ((r2&0xf8)<<6) | ((g2&0xf8)<<1) | (b2>>4);

      if (cache[key]) { *np = (byte) (cache[key] - 1);  cnt1++;	}
      else {
	/* not in cache, have to search the colortable */
	cnt2++;

        mind = 10000;
	for (k=closest=0, cptr=ctab; k<numcols && mind>7; k++,cptr++) {
	  d = abs(r2 - cptr->red)
	    + abs(g2 - cptr->green) 
	    + abs(b2 - cptr->blue);
	  if (d<mind) { mind = d;  closest = k; }
	}
	cache[key] = closest + 1;
	*np = closest;
      }


      /* propogate the error */
      rerr = r2 - ctab[*np].red;
      gerr = g2 - ctab[*np].green;
      berr = b2 - ctab[*np].blue;

      if (j!=jmax) {  /* adjust RIGHT pixel */
	thisptr[0] += (rerr*7)/16;
	thisptr[1] += (gerr*7)/16;
	thisptr[2] += (berr*7)/16;
      }
      
      if (i!=imax) {	/* do BOTTOM pixel */
	nextptr[0] += (rerr*5)/16;
	nextptr[1] += (gerr*5)/16;
	nextptr[2] += (berr*5)/16;

	if (j>0) {  /* do BOTTOM LEFT pixel */
	  nextptr[-3] += (rerr*3)/16;
	  nextptr[-2] += (gerr*3)/16;
	  nextptr[-1] += (berr*3)/16;
	}

	if (j!=jmax) {  /* do BOTTOM RIGHT pixel */
	  nextptr[3] += rerr/16;
	  nextptr[4] += gerr/16;
	  nextptr[5] += berr/16;
	}
	nextptr += 3;
      }
    }
  }

  if (epic!=cpic && epic != NULL) free(epic);
  epic = newpic;  eWIDE = w;  eHIGH = h;

  free(cache);

  CreateXImage();
  SetCursors(-1);
}


