/*
* ximage.c  Neil Gershenfeld  9/11/93
* demonstrates drawing X images by animating sin(k*x)/k*x
*/

#include <X11/Xlib.h>
#include <X11/Xutil.h>
#include <math.h>
#include <string.h>
#include <stdio.h>
#include <stdlib.h>

#define WIDTH 200
#define HEIGHT 200
#define NPTS 200
#define NCOLORS 256
#define DEPTH 8

Display  *D;
Visual   *V;
Colormap C;
int      S,Loop,Point;
Window   W;
GC       Gc,GcRev;
XImage   *I;
XColor   Colors[NCOLORS];
XSetWindowAttributes Attr;
XVisualInfo Info;
double     Data[NPTS][NPTS];
unsigned char     realData[NPTS][NPTS];
unsigned char     Boundary[NPTS][NPTS];
float    x,y,r,z;
int      i,j;
double alpha;
char readme[100];
char * sp;
char *sp2;
int val;

main(int argc, char **argv) {
   D = XOpenDisplay("");
   S = DefaultScreen(D);
   if (XMatchVisualInfo(D,S,DEPTH,PseudoColor,&Info) == 0) {
      printf("Display can not handle %d bit pseudo color\n",DEPTH);
      return;
   }
   V = Info.visual;
	fprintf(stderr, "alpha: %s\n", argv[1]);
   alpha = strtod(argv[1], NULL);
	fprintf(stderr, "alpha2: %f\n", alpha);
   for (i = 0; i < NCOLORS; ++i) {
      Colors[i].pixel = i;
      Colors[i].red = 0;
      Colors[i].green = 0;
      Colors[i].blue = 0;
	if(i < NCOLORS/3){
		Colors[i].green = i*255;
	}
	else if(i < (2*NCOLORS/3)){
		Colors[i].red = i*255;
	}
	else {
		Colors[i].blue = i*255;
	}

/*     Colors[i].red = (int) 4000*sqrt(i);
      Colors[i].green = (int) 256*i;
      Colors[i].blue = (int) 4000*sqrt(i);
*/
	fprintf(stderr, "Color %d is (%d, %d, %d)\n", i, Colors[i].red,
				Colors[i].blue, Colors[i].green);
      Colors[i].flags = DoRed | DoGreen | DoBlue; 
   }
/*	XAllocColor(D, C, Colors); */

   for (i = 0; i < NPTS; ++i)
     for (j = 0; j < NPTS; ++j) {
       Data[i][j] = 0;
       Boundary[i][j] = 0;
     }

   while(fgets(readme, 100, stdin)){
	sp = index(readme, (int) ' ');
	*sp = (char) '\0';
	sp2=index((char *) sp+1, (int) ' ');
	*sp2 = '\0';
	i = atoi(readme);
	j = atoi(sp+1);
	val = atoi(sp2+1);
	Boundary[i][j] = (unsigned char) val;
   }

   for (j = 0; j < NPTS; ++j) {
     Boundary[0][j] = 1;
     Boundary[NPTS-1][j] = 1;
   }
   for (i = 0; i < NPTS; ++i) {
     Boundary[i][0] = 1;
     Boundary[i][NPTS-1] = 1;
   }

/*
   for(i = 50; i< 100; i++)
     for(j = 50; j< 100; j++)
       if( (((i-50)*(i-50))+((j-50)*(j-50))) > 2500 )
	 Boundary[i][j] = (i*j)/100;
*/

   C = XCreateColormap(D,RootWindow(D,S),V,AllocAll);
   XStoreColors(D,C,Colors,NCOLORS); 
   Attr.colormap = C;
   Attr.background_pixel = WhitePixel(D,S);
   Attr.border_pixel = BlackPixel(D,S);
   W = XCreateWindow(D,DefaultRootWindow (D),0,0,WIDTH,HEIGHT,100,
      DEPTH,InputOutput,V,CWColormap|CWBackPixel|CWBorderPixel,&Attr);
   XStoreName(D, W, "ximage output");
   XMapRaised(D, W);
   Gc = XCreateGC (D, W, 0L, (XGCValues *) 0);
   I = XCreateImage(D,V,DEPTH,ZPixmap,0,realData,WIDTH,HEIGHT,DEPTH,0);
   while(1){
      for (i = 1; i < NPTS-1; ++i)
         for (j = 1; j < NPTS-1; ++j) {
	   if(!Boundary[i][j]) {
	     Data[i][j] = (((1.0-alpha)*((double) Data[i][j])) +
					   alpha*
					   (
					    (
					     ((double) Data[i-1][j])+
					     ((double) Data[i][j-1])+
					     ((double) Data[i+1][j])+
					     ((double) Data[i][j+1])
					     )/
					    4.0)
					   );
/*		if(Data[i][j] != 0){
			fprintf(stderr, "Data at %d %d is %g\n", i, j, Data[i][j]);
		} */
		}
	   else
	     Data[i][j] = (double) Boundary[i][j];
	   realData[i][j] = (unsigned char) Data[i][j];
/*	if(realData[i][j] != 0){
		fprintf(stderr, "realdata %d, %d = %d\n", i, j, realData[i][j]);
	}
*/
         }
      XPutImage(D,W,Gc,I,0,0,0,0,NPTS,NPTS);
   }
}