#define GRAPHDEF extern
#include "graph.h"
#include "color.h"

int fitbez(double **pxv, double **pyv, int **pmv,int *pnp, int multivalued);
int var_clear_local(void);
int big_open(char *s);
int big_vec(int nomiss);
int windownorm(void);
int windowdn(int d);
int big_mark(int mnum, double msize);
double fnx(double vv);
double fnx(double vv)
{
	if (xx[1].log==true) {
		return (log10(vv)-log10(wxmin))/(log10(wxmax)-log10(wxmin))
			* xlength + xbl;
	} else {
		return (((vv-wxmin)/(wxmax-wxmin)) * xlength + xbl);
	}
}
double fny(double vv)
{
	if (xx[2].log==true) {
		return (log10(vv)-log10(wymin))/(log10(wymax)-log10(wymin))
			* ylength + ybl;
	} else {
		return (((vv-wymin)/(wymax-wymin)) * ylength + ybl);
	}
}

int bigalli;
static int bigcol1,bigcol2,bigally;
static long bigrecord;
FILE *fptr;
big_open(char *infile2)
{
	char infile[90];
	char *s1,*s2;
	strcpy(infile,infile2);
	if (infile[strlen(infile)-1] == '$') {
		{ int idx,typ;
		var_find(infile,&idx,&typ) ;
		if (idx>=0) var_getstr(idx,infile);
	}}
	s1 = strchr(infile,',');
	bigrecord = 0;
	bigally = bigalli = 0;
	bigcol1 = 1; bigcol2 = 2;
	if (s1!=NULL) {
		s2 = strchr(s1+1,',');
		bigcol1 = atoi(s1+1);
		if (s2!=NULL) {
			if (*(s2+1) == '*') {
				bigally = true;
			} else {
				bigcol2 = atoi(s2+1);
				if (bigcol2==0) {
					gprint("Expecting \"file.name,xcoloumn,ycolumn\" found [%s] \n",infile);
				}
			}
			*s1 = 0;
		}
	}
	if (infile[strlen(infile)-1] == '$') {
		{ int idx,typ;
		var_find(infile,&idx,&typ) ;
		if (idx>=0) var_getstr(idx,infile);
	}}
	strlwr(infile);
	fptr = fopen(infile,"r");
	if (s1!=NULL) *s1 = ',';
	if (fptr==NULL) {
		gprint("Unable to open data file {%s} \n",infile);
		return;
	}
}
int big_line(double *x1, double *y1, int *m);
big_line(double *x1, double *y1, int *m)
{
	static char *bigee;
	int bg1;
	static int ntk;
	if (fptr==NULL) return;
next_line:;
	if (bigally) {
		if (bigalli==0 || bigalli>=ntk) {
			if (bigee!=NULL) { myfree(bigee); bigee = 0;}
			if (feof(fptr)) return;
			fgets(ebuff,390,fptr);
			bigee = sdup(ebuff);
			bigalli = 0;
		}
		strcpy(ebuff,bigee);
		token_data(ebuff,tk,&ntk,tkbuff);
		*x1 = ++bigrecord;
		*y1 = atof(tk[++bigalli]);
		*m = false;
		if (feof(fptr) || abort_key()) {*m = true; return;}
		return;
	}
	if (feof(fptr)) return;
	fgets(ebuff,390,fptr);
	token_data(ebuff,tk,&ntk,tkbuff);
	if (ntk<bigcol1 || ntk<bigcol2) goto next_line;
	bg1 = bigcol1;
	if (bg1==0) bg1 = bigcol2;
	if (*tk[bg1] != '*' && *tk[bigcol2] != '*') {
		if (bigcol1==0) *x1 = ++bigrecord;
		else *x1 = atof(tk[bigcol1]);
		*y1 = atof(tk[bigcol2]);
		*m = false;
	} else { *m = true;}
	if (feof(fptr) || abort_key()) {*m = true; return;}
}
setrange(double x, double y, int m)
{
	if (m) return;
	if (x>range_x2) range_x2 = x;
	if (x<range_x1) range_x1 = x;
	if (y>range_y2) range_y2 = y;
	if (y<range_y1) range_y1 = y;
}
preview_big()
{
	int dn;
	double x1,y1,x2,y2;
	int m1,m2;

	for (dn=1;dn<=ndata;dn++) {
	  if (dp[dn]!=NULL) {
	   if (dp[dn]->bigfile!=NULL && dp[dn]->autoscale) {
		big_open(dp[dn]->bigfile);
		if (fptr==NULL) goto endbigm;
		if (feof(fptr))  goto xxx;
		for (;;) {
aaa2:			if (feof(fptr) || abort_key())  goto xxx;
			big_line(&x1,&y1,&m1);
			setrange(x1,y1,m1);
		}
xxx:		big_close();
endbigm:	;
	   }
	  }
	}	  
}

request()
{
	static char s[20];
	gprint("Press return to continue \n");
	gets(s);
}
draw_bars()
{
	int b,nb,d,i,di,df,dt;
	long sidecolor,topcolor;
	double x3d,y3d;
	int notop;
	double *px,*py,xzz,*xt,*yt,*xf,*yf,bwid,whole_wid,bdis,yfval;
	int *m;
	double x,y,w;

    for (b=1;b<=g_nbar;b++) {
	done_line = true;
	if (br[b]==0) { gprint("Error, bars struct zero \n"); return;}
	nb = br[b]->ngrp;	/* number of bars in each group */
	d = br[b]->to[0]; 	/* dataset for first bar */
	if (d==0 || dp[d]==NULL)  {gprint("Error, bars zero dataset \n");return;}
	/* find the smallest x distance in dataset 1 */
	px = dp[d]->xv;
	py = dp[d]->yv;
	if (px==NULL || nb==0) { gprint("error in bar data dn=%d  ngrp=%d\n",d,nb); return;}
	m = dp[d]->miss;
	xzz = 1e30;
	for (i=1;i<dp[d]->np;i++,px++) {
		w = *(px+1) - *(px);
		if (w>0) if (w<xzz) xzz = w;
	}
	if (br[b]->width==0) br[b]->width = xzz/(nb*2);
	if (br[b]->dist==0) br[b]->dist = br[b]->width*1.4 ;
	bwid = br[b]->width;
	bdis = br[b]->dist;

	windowdn(d);
	g_gsave();
	for (di=0;di<nb;di++) {
		df = br[b]->from[di];
		dt = br[b]->to[di];
		  	whole_wid = (nb-1)*bdis+bwid;
		g_set_line_width(br[b]->lwidth[di]);
		g_set_line_style(&br[b]->lstyle[di][0]);
		if (br[b]->color[di] == 0) br[b]->color[di] = COLOR_BLACK;
		g_set_color(br[b]->color[di]);
		g_set_fill(br[b]->fill[di]);

		x3d = br[b]->x3d;
		y3d = br[b]->y3d;
		topcolor = br[b]->top[di];
		sidecolor = br[b]->side[di];
		notop = br[b]->notop;

		df = br[b]->from[di];
		dt = br[b]->to[di];
		if (dp[df]==NULL || dp[dt]==NULL) {
			gprint("No data in bargraph datasets\n");
			goto bar_exit;
		}
		yf = dp[df]->yv;
		xt = dp[dt]->xv;
		yt = dp[dt]->yv;
		m = dp[dt]->miss;
		if (yt==NULL) {
			gprint("No data in bargraph dataset. d(%d) \n",dt);
			goto bar_exit;
		}

		windowdn(d);
		for (i=0;i<dp[d]->np;i++,xt++,yt++,m++) {
		  if (yf!=NULL) yfval = *yf++;
			else 	yfval = 0;
		  if (!*m) {
			draw_bar((*xt)-whole_wid/2+di*bdis,yfval,*yt,bwid
				,x3d,y3d,sidecolor,topcolor,notop);
		  }
		}
		windownorm();
	}
bar_exit:;
	g_grestore();
    }
}
/*--------------------------------------------------------------------------*/
int box3d(double x1, double y1, double x2, double y2,double x3d,double y3d,long sidecolor, long topcolor, int notop);
draw_bar(double x, double yf, double yt, double wd,double x3d,double y3d,long sidecolor, long topcolor, int notop)
{
	double x1,y1,x2,y2;

	x = x + wd/2;
	/* ! draw a bar, wd wide, CENTREed at x , from yf, to yt, */

	x1 = x - wd/2;
	y1 = yf;
	x2 = x + wd/2;
	y2 = yt;
	box_clip(&x1,&y1,wxmin,wymin,wxmax,wymax);
	box_clip(&x2,&y2,wxmin,wymin,wxmax,wymax);

	if (x3d!=0) {
		box3d(fnx(x1),fny(y1),fnx(x2),fny(y2),x3d,y3d,sidecolor,topcolor,notop);
	}
	g_box_fill(fnx(x1),fny(y1),fnx(x2),fny(y2));
	g_box_stroke(fnx(x1),fny(y1),fnx(x2),fny(y2));
}
box3d(double x1, double y1, double x2, double y2,double x3d,double y3d,long sidecolor, long topcolor, int notop)
{
	/* assuming x3d is positive for the moment */
	double xx;

	if (x1>x2) { xx = x1; x1 = x2; x2 = xx;}
	if (y1>y2) { xx = y1; y1 = y2; y2 = xx;}

	x3d = x3d*(x2-x1);
	y3d = y3d*(x2-x1);

	if (x3d<0) {
		xx = x1; x1 = x2; x2 = xx;
	}
	g_gsave();
	g_set_path(true);
	g_set_line_join(1);	/* use rounded lines to avoid ucky peeks */
	g_newpath();
	g_move(x2,y1);
	g_line(x2+x3d,y1+y3d);
	g_line(x2+x3d,y2+y3d);
	g_line(x2,y2);
	g_line(x2,y1);
	if (topcolor!=0) {
	   g_set_fill(sidecolor);
	   g_fill();
	}
	g_stroke();
	g_newpath();

	if (!notop) {	/* now draw top of bar  */
	 g_move(x2,y2);
	 g_line(x2+x3d,y2+y3d);
	 g_line(x1+x3d,y2+y3d);
	 g_line(x1,y2);
	 g_line(x2,y2);
	 if (topcolor!=0) {
		g_set_fill(topcolor);
		g_fill();
	 }
	 g_stroke();
	}
	g_newpath();
	g_set_path(false);
	g_newpath();
	g_grestore();
}
/*------------------------------------------------------------------*/
/*
/*do_smooth()
/*{
/*	int i,j;
/*	double *x,*y;
/*	for (i=1;i<=ndata;i++) {
/*		if (dp[i]!=NULL) if (dp[i]->smooth && dp[i]->np>3) {
/*			gr_nomiss(i);	/* throw away missing points */
/*			fitbez(&dp[i]->xv,&dp[i]->yv,&dp[i]->miss,&dp[i]->np);
/*		}
/*	}
/*}
*/

/*------------------------------------------------------------------*/
gr_thrownomiss()
{
	int i,j;
	for (i=1;i<=ndata;i++) {
		if (dp[i]!=NULL) if (dp[i]->nomiss && dp[i]->np>0) {
			gr_nomiss(i);	/* throw away missing points */
		}
	}
}
gr_nomiss(int i)
{
	double *nx,*ny,*xt,*yt,oldlwidth;
	int *m,*nm;
	int dt,dn,k,npnts;

	for (dn=1;dn<=ndata;dn++) {
	  if (dp[dn]!=NULL) if (dp[dn]->xv!=NULL) {
		k = 0;
		yt = dp[dn]->yv;
		xt = dp[dn]->xv;
		m = dp[dn]->miss;
		nx = xt; ny = yt; nm = m;
		npnts = dp[dn]->np;
		for (i=0;i<(npnts);i++,m++,xt++,yt++) {
			if (!*m) {
				*nx++ = *xt;
				*ny++ = *yt;
				*nm++ = *m;
				k++;
			}
		}
		dp[dn]->np = k;
	  }
	}
}
roundrange(double *gmin,double *gmax)
{
	double delta,st,expnt,n,dticks;
	int ni;

	if (*gmax < *gmin) {
		*gmax = 100; *gmin = 10;
		return;
	}
	delta = *gmax - *gmin;
	if (delta==0) return;
	st = delta/10;
	expnt = floor(log10(st));
	n = st/pow(10,expnt);
	if (n>5)
		ni = 10;
	else if (n>2)
		ni = 5;
	else if (n>1)
		ni = 2;
	else
		ni = 1;

	dticks = ni * pow(10,expnt);
	if (*gmin ==  floor( *gmin/ dticks) * dticks)
		*gmin = *gmin - dticks;
	else
		*gmin = (floor(*gmin/ dticks) * dticks );

	if (( floor( *gmax/ dticks) * dticks) != (*gmax / dticks) * dticks )
		*gmax = floor(*gmax / dticks ) * dticks + dticks;
	else
		*gmax = *gmax  + dticks;
}

window_set()
{
	wxmin = range_x1;
	wxmax = range_x2;
	wymin = range_y1;
	wymax = range_y2;
	if (!xx[1].log)  roundrange(&wxmin,&wxmax);
	if (!xx[2].log)  roundrange(&wymin,&wymax);

	if (xx[1].minset) wxmin = xx[1].min;
	xx[1].min = wxmin;
	if (!xx[3].minset) xx[3].min = wxmin;

	if (xx[1].maxset) wxmax = xx[1].max;
	xx[1].max = wxmax;
	if (!xx[3].maxset) xx[3].max = wxmax;

	if (xx[2].minset) wymin = xx[2].min;
	xx[2].min = wymin;
	if (!xx[4].minset) xx[4].min = wymin;

	if (xx[2].maxset) wymax = xx[2].max;
	xx[2].max = wymax;
	if (!xx[4].maxset) xx[4].max = wymax;

	/* should go thru datasets to find max,min */
	if (wxmin>=wxmax || wymin>=wymax) {
		gprint("Warning no data (or no spread) and no min,max so making up range\n");
		wxmin = 0;
		wxmax = 10;
		wymin = 0;
		wymax = 10;
	}
}
draw_fills()
{
	double x1,y1,x2,y2,ymx,lastx,lasty;
	double *ssxt,*ssyt,*xt,*yt;
	int dt,dn,i,n,*m,npnts,*ssm,free_smoothed=0;
	struct fill_data *ff;

	/* set clipping region */

	for (n=1;n<=nfd;n++) {
	  if (fd[n]->type==0) return;
	  done_line = true;
	  ff = fd[n];
	  if (wxmin > ff->xmin) ff->xmin = wxmin;
	  if (wxmax < ff->xmax) ff->xmax = wxmax;
	  if (wymin > ff->ymin) ff->ymin = wymin;
	  if (wymax < ff->ymax) ff->ymax = wymax;

	g_beginclip();		/* saves current clipping */
	g_set_path(true);
	g_newpath();
	g_move(fnx(ff->xmin),fny(ff->ymin));
	g_box_stroke(fnx(ff->xmin),fny(ff->ymin),fnx(ff->xmax),fny(ff->ymax));
	g_clip();		/* sets current path to be clipping */

	  gntmp = 0;
	  dn = ff->da;
	  if (dp[dn]==NULL) { gprint("No data in fill dataset at all \n");return;}
	  free_smoothed = false;
	  yt = dp[dn]->yv;
	  xt = dp[dn]->xv;
	  m = dp[dn]->miss;
	  npnts = dp[dn]->np;

		if (dp[dn]->smooth && npnts>3 && npnts<190) {
			gr_nomiss(dn);	/* throw away missing points */
			fitbez(&xt,&yt,&m,&npnts,dp[dn]->smoothm);
			ssxt = xt; ssyt = yt; ssm = m;
			free_smoothed = true;
		}

	  ymx = ff->ymax;
	  if (xt==NULL) { gprint("No data in fill dataset \n"); return;}
	  x1 = *xt;
	  y1 = *yt;
	  switch(ff->type) {
	    case 1: /* x1,d1 */
		ymx = ff->ymin;
	    case 2: /* d1,x2 */
		fill_vec(*xt,ymx,*xt,*yt);
		for (i=0;i<npnts-1;i++,xt++,yt++) {
			fill_vec(*xt,*yt,*(xt+1),*(yt+1));
		}
		fill_vec(*xt,*yt,*xt,ymx);
		fill_vec(*xt,ymx,*dp[dn]->xv,ymx);
		break;
	    case 3: /* d1,d2 */
		for (i=0;i<npnts-1;i++,m++,xt++,yt++) {
			fill_vec(*xt,*yt,*(xt+1),*(yt+1));
			x2 = *(xt+1); y2 = *(yt+1);
		}
		dn = ff->db;
		yt = dp[dn]->yv;
		xt = dp[dn]->xv;
		m = dp[dn]->miss;
		npnts = dp[dn]->np;

		if (free_smoothed) {myfrees(ssxt,"Fill1"); myfrees(ssyt,"x"); myfrees(ssm,"y");}
		free_smoothed = false;
		if (dp[dn]->smooth && npnts>3 && npnts<190) {
			gr_nomiss(dn);	/* throw away missing points */
			npnts = dp[dn]->np;
			fitbez(&xt,&yt,&m,&npnts,dp[dn]->smoothm);
			ssxt = xt; ssyt = yt; ssm = m;
			free_smoothed = true;
		}

		xt += npnts - 1;
		yt += npnts - 1;
		fill_vec(x2,y2,*xt,*yt);
		for (i=0;i<npnts-1;i++,m--,xt--,yt--)
			fill_vec(*xt,*yt,*(xt-1),*(yt-1));
		fill_vec(*xt,*yt,x1,y1);
		break;
	    case 4: /* d1 */
		for (i=0;i<npnts-1;i++,m++,xt++,yt++) {
			if (!*m && !*(m+1))
				fill_vec(*xt,*yt,*(xt+1),*(yt+1));
		}
		fill_vec(*xt,*yt,x1,y1);
		break;
	  }

	if (free_smoothed) {myfrees(ssxt,"Fill2"); myfrees(ssyt,"f4"); myfrees(ssm,"f5");}
	/* Paint the region defined */
	g_set_fill(ff->color);
	g_newpath();
	g_move(fnx(tmpf[1]),fny(tmpf[2]));
	lastx = tmpf[1]; lasty=tmpf[2];
	for (i=1;i<=gntmp;i+=4) {
		if (lastx!=tmpf[i] || lasty!=tmpf[i+1])  {
			g_closepath();
			g_move(fnx(tmpf[i]),fny(tmpf[i+1]));
		}
		g_line(fnx(tmpf[i+2]),fny(tmpf[i+3]));
		lastx = tmpf[i+2]; lasty = tmpf[i+3];
	}
	g_closepath();
	g_fill();
	g_set_path(false);
	g_endclip();
	}
}
/*------------------------------------------------------------------*/
fill_vec(double x1, double y1, double x2, double y2)
{
	alloc_temp(gntmp);
	if (tmpf == NULL) {gprint("Ran out of tmp storage for fill \n"); gle_abort("tmp data \n"); }
	tmpf[++gntmp] = x1;
	tmpf[++gntmp] = y1;
	tmpf[++gntmp] = x2;
	tmpf[++gntmp] = y2;
}
int setupdown(char *s, int *doup, int *upd, int *upp, double *upval);
setupdown(char *s, int *doup, int *upd, int *upp, double *upval)
{
	static char vbuf[9];
	int j=0;

	*upd = 0;
	*doup = true;
	*upp = false;
	if (*s==0) { *doup=false; return;}
	if (*s=='D') { *upd = atoi(s+1); return; }
	if (strstr(s,"%")!=NULL) { j=1; *upp = true; }
	strncpy(vbuf,s,strlen(s)-j+1);
	*upval = atof(vbuf);
}
/*------------------------------------------------------------------*/
/* 	d3 err .1			*/
/* 	d3 err 10%			*/
/* 	d3 errup 10% errdown d2		*/
/* 	d3 err d1 errwidth .2		*/
dataset_null(int i)
{
	if (dp[i]==NULL) {
		gprint("Dataset %d doesn't exist at all\n",i);
		return true;
	}
	if (dp[i]->yv == NULL) {
		gprint("Dataset %d doesn't exist\n",i);
		return true;
	}
	return false;
}
#ifdef HERROR		/* because it won't fit on the PC */
draw_herr()
{
	int errup,doup,upd,upp;
	double upval,eup,edown,ewid;
	int errdown,dodown,downd,downp;
	double downval;
	double *xt,*yt,x;
	int *m,mup,mdown;
	int dt,dn,i;

	g_gsave();
	for (dn=1;dn<=ndata;dn++) {
	  dpp = dp[dn];
	  if (dpp!=NULL) if (dpp->herrup[0]!=0 || dpp->herrdown[0]!=0) {
		done_line = true;
		g_get_hei(&x);
		if (dpp->herrwidth==0) dpp->herrwidth = x/3;
		ewid = dpp->herrwidth;
		setupdown(dpp->herrup,&doup,&upd,&upp,&upval);
		setupdown(dpp->herrdown,&dodown,&downd,&downp,&downval);
		g_set_color(dpp->color);
		g_set_line_width(dpp->lwidth);
		windowdn(dn);
		yt = dpp->yv;
		xt = dpp->xv;
		m = dpp->miss;
		if (upd) if (dataset_null(upd)) return;
		if (downd) if (dataset_null(downd)) return;
		for (i=0;i<dpp->np;i++,m++,xt++,yt++) {
			mup = false; mdown = false;
			if (upd) {
				eup = *(dp[upd]->yv+i);
				mup = *(dp[upd]->miss+i);
			} else {
				eup = upval;
				if (upp) eup = (eup * *xt)/100;
			}
			if (downd) {
				edown = *(dp[downd]->yv+i);
				mdown = *(dp[downd]->miss+i);
			} else {
				edown = downval;
				if (downp) edown = (edown * *xt)/100;
			}
			/* if (eup!=-999) { */
			 if (doup) if (!*m) if (!mup) draw_herrbar(*xt,*yt,eup,ewid);
			 if (dodown) if (!*m) if (!mdown) draw_herrbar(*xt,*yt,-edown,ewid);
			/* } */
		}
		windownorm();
	  }
	}
exit_err:;
	g_grestore();
}
draw_herrbar(double x, double y, double eup, double ewid)
{
	if (x>=wxmin && x<=wxmax && y>=wymin && y<=wymax ) {
		g_move(fnx(x),fny(y));
		g_line(fnx(x-eup),fny(y));
		g_move(fnx(x-eup),-ewid/2+fny(y));
		g_line(fnx(x-eup),ewid/2+fny(y));
	}
}
#endif
draw_err()
{
	int errup,doup,upd,upp;
	double upval,eup,edown,ewid;
	int errdown,dodown,downd,downp,big;
	double downval;
	double *xt,*yt,x;
	int *m,mup,mdown;
	int dt,dn,i;

	g_gsave();
	for (dn=1;dn<=ndata;dn++) {
	  big = false;
	  dpp = dp[dn];
	  if (dpp!=NULL) if (dpp->errup[0]!=0 || dpp->errdown[0]!=0) {
		if (dpp->bigfile!=NULL) {
			big = true;
			big_open(dpp->bigfile);
	   	}
		done_line = true;
		g_get_hei(&x);
		if (dpp->errwidth==0) dpp->errwidth = x/3;
		ewid = dpp->errwidth;
		setupdown(dpp->errup,&doup,&upd,&upp,&upval);
		setupdown(dpp->errdown,&dodown,&downd,&downp,&downval);
		g_set_color(dpp->color);
		g_set_line_width(dpp->lwidth);
		windowdn(dn);
		yt = dpp->yv;
		xt = dpp->xv;
		m = dpp->miss;
		if (upd) if (dataset_null(upd)) return;
		if (downd) if (dataset_null(downd)) return;
		if (big) {{
		 	double x1,y1,x2,y2;
		 	int m1,m2;

			if (fptr==NULL) goto endbigm;
			if (feof(fptr))  goto xxx;
			for (;;) {
aaa2:
				if (feof(fptr) || abort_key())  goto xxx;
				big_line(&x1,&y1,&m2);
				if (m2) goto aaa2;
				mup = false; mdown = false;
				if (upd) {
					eup = *(dp[upd]->yv+i);
					mup = *(dp[upd]->miss+i);
				} else {
					eup = upval;
					if (upp) eup = (eup * y1)/100;
				}
				if (downd) {
					edown = *(dp[downd]->yv+i);
					mdown = *(dp[downd]->miss+i);
				} else {
					edown = downval;
					if (downp) edown = (edown * y1)/100;
				}
				if (doup) if (!mup) draw_errbar(x1,y1,eup,ewid);
				if (dodown) if (!mdown) draw_errbar(x1,y1,-edown,ewid);
			}
xxx:			fclose(fptr); fptr=NULL;
endbigm:		;
		}} else {
		  for (i=0;i<dpp->np;i++,m++,xt++,yt++) {
			mup = false; mdown = false;
			if (upd) {
				eup = *(dp[upd]->yv+i);
				mup = *(dp[upd]->miss+i);
			} else {
				eup = upval;
				if (upp) eup = (eup * *yt)/100;
			}
			if (downd) {
				edown = *(dp[downd]->yv+i);
				mdown = *(dp[downd]->miss+i);
			} else {
				edown = downval;
				if (downp) edown = (edown * *yt)/100;
			}
			/* if (eup!=-999) { */
			 if (doup) if (!*m) if (!mup) draw_errbar(*xt,*yt,eup,ewid);
			 if (dodown) if (!*m) if (!mdown) draw_errbar(*xt,*yt,-edown,ewid);
			/* } */
		  }
		}
		windownorm();
	  }
	}
exit_err:;
	g_grestore();
#ifdef HERROR
	draw_herr();
#endif
}
draw_errbar(double x, double y, double eup, double ewid)
{
	if (x>=wxmin && x<=wxmax && y>=wymin && y<=wymax ) {
		g_move(fnx(x),fny(y));
		g_line(fnx(x),fny(y+eup));
		g_move(fnx(x)-ewid/2,fny(y+eup));
		g_line(fnx(x)+ewid/2,fny(y+eup));
	}
}
double dwx1,dwy1,dwx2,dwy2;
windowdn(int d)
{
	dwx1 = wxmin; dwy1 = wymin; dwx2 = wxmax; dwy2 = wymax;
	if (dp[d]==NULL) return;
	if (dp[d]->xmax > dp[d]->xmin) {
		wxmin = dp[d]->xmin; wxmax = dp[d]->xmax;
	}
	if (dp[d]->ymax > dp[d]->ymin) {
		wymin = dp[d]->ymin; wymax = dp[d]->ymax;
	}
}
windownorm()
{
	wxmin = dwx1 ; wymin = dwy1 ; wxmax = dwx2 ; wymax = dwy2 ;
}
draw_lines()
{
	double *ssxt,*ssyt,*xt,*yt,oldlwidth;
	int *m,*ssm;
	char oldlstyle[10];
	int dt,dn,i,npnts,free_smoothed,big;

	g_gsave();
	g_get_line_style(oldlstyle);
	g_get_line_width(&oldlwidth);
	for (dn=1;dn<=ndata;dn++) {
	  big = false;
	  last_vecx = 1e10;
	  last_vecy = 1e10;
	  if (dp[dn]!=NULL) {
	   if (dp[dn]->bigfile!=NULL) {
		big=true;
		big_open(dp[dn]->bigfile);
	   }
	   if (dp[dn]->xv!=NULL || big)
	    if (dp[dn]->line==true || dp[dn]->lstyle[0]!=0) {
		free_smoothed = false;
		g_set_line_style(oldlstyle); /* use defualts for each */
		g_set_line_width(oldlwidth); /* use defaults for each */
		g_set_line_style(dp[dn]->lstyle);
		g_set_color(dp[dn]->color);
		g_set_line_width(dp[dn]->lwidth);
		windowdn(dn);
		done_line = true;
		yt = dp[dn]->yv;
		xt = dp[dn]->xv;
		m = dp[dn]->miss;
		npnts = dp[dn]->np;

		if (dp[dn]->smooth && npnts>3 && npnts<190) {
			/* (BUG, gr_nomiss has not been written HAS NOW !!! */
			gr_nomiss(dn);	/* throw away missing points */
			npnts = dp[dn]->np;
			fitbez(&xt,&yt,&m,&npnts,dp[dn]->smoothm);
			ssxt = xt; ssyt = yt; ssm = m;
			free_smoothed = true;
		}

		if (big) big_vec(dp[dn]->nomiss);
		else {
		  for (i=0;i<(npnts-1);i++,m++,xt++,yt++) {
			if (!*m && !*(m+1))
				draw_vec(*xt,*yt,*(xt+1),*(yt+1));
		  }
		}
		if (free_smoothed) {myfrees(ssxt,"Line1"); myfrees(ssyt,"l2"); myfrees(ssm,"l3");}
		windownorm();
	  }
	  big_close();
	 }
	}
	g_grestore();
}
big_vec(int nomiss)
{
	double x1,y1,x2,y2;
	int m1,m2;


	if (fptr==NULL) return;
aaa:
	if (feof(fptr))  goto xxx;
	big_line(&x1,&y1,&m1);
	if (m1) goto aaa;
	for (;;) {
aaa2:
		if (feof(fptr) || abort_key())  break;
		big_line(&x2,&y2,&m2);
		if (nomiss && m2) goto aaa2;
		if (!m1 && !m2) draw_vec(x1,y1,x2,y2);
		x1 = x2; y1 = y2; m1 = m2;
	}
xxx:	fclose(fptr); fptr = NULL;
}
big_close()
{
	if (fptr!=NULL) fclose(fptr);
	fptr = NULL;
}
big_mark(int mnum, double msize)
{
	double x1,y1,x2,y2;
	int m1,m2;


	if (fptr==NULL) return;
	if (feof(fptr))  goto xxx;
	for (;;) {
aaa2:
		if (feof(fptr) || abort_key())  break;
		big_line(&x1,&y1,&m2);
		if (m2) goto aaa2;
		draw_mark(x1,y1,mnum,msize,1.0);
	}
xxx:	fclose(fptr); fptr=NULL;
}
draw_vec(double x1,double y1,double x2,double y2)
{
	int invis;
	if (x1>=wxmin && x1<=wxmax && y1>=wymin && y1<=wymax
 	  && x2 >= wxmin && x2<=wxmax && y2>=wymin && y2<=wymax);
	else { 		/* ok one or both are outside our box */
		gclip(&x1,&y1,&x2,&y2,wxmin,wymin,wxmax,wymax,&invis);
		if (invis) return;
	}
	if (x1!=last_vecx || y1!=last_vecy) {
		g_move(fnx(x1),fny(y1));
	}
	g_line(fnx(x2),fny(y2));
	last_vecx = x2;
	last_vecy = y2;
}
/*---------------------------------------------------------------------------*/
draw_markers()
{
	double msize;
	double *yd,dval,*xt,*yt,oldlwidth;
	char oldlstyle[10];
	int *m;
	int dt,dn,i,big;

	g_gsave();
	g_get_line_style(oldlstyle);
	g_get_line_width(&oldlwidth);
	for (dn=1;dn<=ndata;dn++) {
	  big = false;
	  last_vecx = 1e10;
	  last_vecy = 1e10;
	  if (dp[dn]!=NULL) if (dp[dn]->marker!=0) {
		if (dp[dn]->bigfile!=NULL) {
			big=true;
			big_open(dp[dn]->bigfile);
		} else big = false;
		g_set_line_width(oldlwidth); /* use defaults for each */
		g_set_color(dp[dn]->color);
		g_set_line_width(dp[dn]->lwidth);
		windowdn(dn);
		done_line = true;
		yt = dp[dn]->yv;
		xt = dp[dn]->xv;
		m = dp[dn]->miss;
		msize = dp[dn]->msize;
		if (msize==0) msize = g_fontsz*.6;
		if (dp[dn]->mscale != 0) msize = msize* dp[dn]->mscale;
		if (big) big_mark(dp[dn]->marker,msize);
		else {
		 if (dp[dn]->mdata!=0) yd = dp[dp[dn]->mdata]->yv;
		 for (i=0;i<dp[dn]->np;i++,m++,xt++,yt++,yd++) {
			if (!*m) {
				if (dp[dn]->mdata==0) dval = 1; else dval = *yd;
				draw_mark(*xt,*yt,dp[dn]->marker,msize,dval);
			}
		 }
		}
		windownorm();
	  }
	}
	g_grestore();
}
/*---------------------------------------------------------------------------*/
draw_mark(double x, double y,int mrk, double msize, double dval)
{
	/* globals   *** xmin, xmax, ymin, ymax,  ox,oy,  xbl, ybl */
	if (x >= wxmin && x<=wxmax && y>=wymin && y<=wymax) {
		g_move(fnx(x),fny(y));
		g_marker2(mrk,msize,dval);
	}
}
/*---------------------------------------------------------------------------*/
box_clip(double *x,double *y,double wxmin,double wymin,double wxmax,double wymax)
{
	if (*x > wxmax) *x = wxmax;
	if (*y > wymax) *y = wymax;
	if (*x < wxmin) *x = wxmin;
	if (*y < wymin) *y = wymin;
}
gclip(double *x1,double *y1,double *x2,double *y2
		,double xmin,double ymin,double xmax,double ymax,int *invis)
{
	double dx,dy;

	*invis = true;
	if (*x2 > xmax) {
		if (*x1 > xmax) return;
		dx = *x2 - *x1;
		dy = *y2 - *y1;
		if (dx==0) return;
		*y2 = *y1 + dy*(xmax-*x1)/dx;
		*x2 = xmax;
	}
	if (*x1 > xmax) {
		dx = *x1 - *x2;
		dy = *y1 - *y2;
		if (dx==0) return;
		*y1 = *y2 + dy*(xmax-*x2)/dx;
		*x1 = xmax;
	}
	if (*y2 > ymax) {
		if (*y1 > ymax) return;
		dy = *y2 - *y1;
		dx = *x2 - *x1;
		if (dy==0) return;
		*x2 = *x1 + dx*(ymax-*y1)/dy;
		*y2 = ymax;
	}
	if (*y1 > ymax) {
		dx = *x1 - *x2;
		dy = *y1 - *y2;
		if (dy==0) return;
		*x1 = *x2 + dx*(ymax-*y2)/dy;
		*y1 = ymax ;
	}
	if (*x2 < xmin) {
		if (*x1 < xmin) return;
		dx = *x2 - *x1;
		dy = *y2 - *y1;
		if (dx==0) return;
		*y2 = *y1 + dy*(xmin-*x1)/dx;
		*x2 = xmin;
	}
	if (*x1 < xmin) {
		dx = *x1 - *x2 ;
		dy = *y1 - *y2 ;
		if (dx==0) return;
		*y1 = *y2 + dy*(xmin-*x2)/dx;
		*x1 = xmin;
	}
	if (*y2 < ymin) {
		if (*y1 < ymin) return;
		dy = *y2 - *y1;
		dx = *x2 - *x1;
		if (dy==0) return;
		*x2 = *x1 + dx*(ymin-*y1)/dy;
		*y2 = ymin;
	}
	if (*y1 < ymin) {
		dy = *y1 - *y2;
		dx = *x1 - *x2;
		if (dy==0) return;
		*x1 = *x2 + dx*(ymin-*y2)/dy;
		*y1 = ymin;
	}
	*invis = false;
}
axis_type(char *s)
{
	if (strncmp(s,"X2",2)==0) return 3;
	if (strncmp(s,"Y2",2)==0) return 4;
	if (strncmp(s,"X",1)==0) return 1;
	if (strncmp(s,"Y",1)==0) return 2;
	gprint("axis type did not match {%s} assuming xaxis \n",s);
	return 1;
}
vinit_axis(int i)
{
	long c;
	memset(&xx[i],0,sizeof(xx[i]));
	g_get_color(&c);
	xx[i].color = c;
	xx[i].side_color = c;
	xx[i].ticks_color = c;
	xx[i].label_color = c;
	xx[i].subticks_color = c;
 	xx[i].side_lwidth = -1;
 	xx[i].ticks_lwidth = -1;
 	xx[i].subticks_lwidth = -1;
 	xx[i].label_lwidth = -1;
	xx[i].type = i;
	if (i>2) xx[i].label_off = true;
}
doskip(char *s,int *ct)
{
	if (*s==' ') (*ct)++;
}
int polish_eval(char *exp, double *x);
double get_next_exp(char (*(*tk)[]),int *ntok,int *curtok)
{
	static int elen,etype,cp,i;
	static double x;

	(*curtok)++;
	cp = 0;
	elen = 0;
	etype = 1;
	dbg for (i=1;i<=*ntok;i++)  gprint("{%s} ",(*tk)[i]);
	dbg gprint("\n");
	dbg gprint("**get exp token ct %d  {%s} \n",*curtok,(*tk)[*curtok]);
	polish_eval((*tk)[*curtok],&x);
	return x;
}
polish_eval(char *exp, double *x)
{
	static long ebuff[300];
	static char ostr[20];
	static int elen,etype,cp,otyp;
	cp = 0;
	elen = 0;
	etype = 1;
	polish(exp,(char *) &ebuff,&elen,&etype);
	eval(ebuff,&cp,x,ostr,&otyp);
}
checktok(char *t, char *want)
{
	if (strcmp(t,want)==0) return true;
	else {
		gprint("Found token {%s} Wanted {%s} \n",t,want);
		return false;
	}
}
/*  LET d2 = exp(x) [ FROM exp TO exp [ STEP exp ] ]	*/
/*  LET d3 = exp(dn,d...)				*/
do_let(char *letcmd)
{
    	static char *tk[200];
	static int ndn,dn_idx[11],dn_var[11];
    	static char tkbuff[200];
	static int ntk,ct,i,elen,etype,dn,cp,savect,tt,dd,np,j;
	double letfrom,letto,letstep,ix;
	int varx,npnts,*newmiss;
	double *xt,*yt,*xf,logstep,*newx,*newy;
	int *mt,*mm;
	static char ostr[20];
	int gotfrom,vartype,otyp;
    	double x,y,ox,oy;
    	static char space_str[] = " ";

	if (tk[199]==NULL) for (i=0;i<200;i++) tk[i] = space_str;
	token(letcmd,(TOKENS) tk,&ntk,tkbuff);
	ct = 2;

	dd = atoi(tk[ct]+1);
	ct++;
	if (strcmp(tk[ct],"=")!=0) {
		gprint("%s\n",letcmd);
		gprint("Found {%s} When expecting {=} at token %d \n",tk[ct]
			,ct);
		return;
	}
	ct++;
	savect = ct;
	gotfrom = false;
	for (ct=savect;ct<=ntk;ct++) {
		if (strcmp(tk[ct],"FROM")==0) {
			gotfrom = true;
			letfrom = get_next_exp((TOKENS) tk,&ntk,&ct);
			ct++;
			if (!checktok(tk[ct],"TO")) return;
			letto = get_next_exp((TOKENS) tk,&ntk,&ct);
			ct++;
			if (strcmp(tk[ct],"STEP")==0) {
				letstep = get_next_exp((TOKENS) tk,&ntk,&ct);
				ct++;
			} else
				letstep = (letto-letfrom)/100;
			ct = ntk;
		}
	}

	var_set_local();
	var_alloc_local();
	ct = savect;
	elen = 0;
	etype = 1;
	polish(tk[ct],(char *) &ebuff,&elen,&etype);

	ndn = 0;
	if (!gotfrom) {
		dn_var[0]=12; dn_idx[0]=111;
		var_find_dn(dn_idx, dn_var, &ndn);
		if (ndn==0) {
			gprint("Expecting    FROM exp TO exp [ STEP exp] \n");
			gprint("or expression containing Dataset d1,d2 etc\n");
			gprint("But only got {%s}\n",letcmd);
		}
	}


	dn = dn_var[0];
	if (ndn>0) 	np = dp[dn]->np;
	else {
		np = (letto-letfrom)/letstep + 1;
		if (xx[1].log==true) {
		 if (letstep<9 || letstep>300) {
		   letstep=60;
		   gprint("Due to LOG xaxis stepsize is taken as number of steps (set to 60)\n");
		 }
		 np = letstep*2;
		 letstep = letstep/3;
		 logstep = 1+ (log10(letto)-log10(letfrom)) /((double) letstep);
		}
	}
/*	dbg for (i=ndata;i<=dd;i++) printf("dat apointer [%d] %p \n",i,dp[i]);
*/

	if (ndata<dd) ndata = dd;
	if (dp[dd]==NULL) {
		dp[dd] = myallocz(sizeof(*dp[dd]));
		copy_default(dd);
	}
	/* copy default dataset settings */
	if (dp[dd]==NULL) gprint("Memory allocation error, graph dataset \n");
	np++;
	if (dp[dd]->xv != NULL) gprint("Over writing dataset %d \n",dd);
	newx = myallocz(sizeof(x)*np);
	newy = myallocz(sizeof(x)*np);
	newmiss = myallocz(sizeof(i)*np);
	np--;
	dp[dd]->np = np;
	xt = newx ; yt = newy; mt = newmiss;

	npnts = 0;
	var_find("X",&varx,&vartype);
	dn = dn_var[0];
	if (ndn>0) {
		xf = dp[dn]->xv;
		mm = dp[dn]->miss;
		for (i=0;i<dp[dn]->np;i++,xf++,mm++) {
			if (varx>=0) var_set(varx,*xf);
			for (j=0;j<ndn;j++) {
				if (dp[dn_var[j]] == NULL) {
					gprint("Dataset not defined {%d} \n",dn_var[j]);
				} else {
					var_set(dn_idx[j],dp[dn_var[j]]->yv[i]);
				}
			}
			cp = 0;
			eval((long *) ebuff,&cp,&x,ostr,&otyp);
			if (! *mm) {
				npnts++;
				*(xt++) = *xf;
				*(yt++) = x;
				*(mt++) = 0;
				setrange(*xf,x,0);
			}
	 	}
	} else {
		for (ix=letfrom;ix<=letto  && npnts<np ;) {
			if (varx>=0) var_set(varx,ix);
			cp = 0;
			eval((long *) ebuff,&cp,&x,ostr,&otyp);
			npnts++;
			*(xt++) = ix;
			*(yt++) = x;
			*(mt++) = 0;
			setrange(ix,x,0);
			if (xx[1].log==true) {
				ix = ix*logstep;
				if (ix>letto) ix = letto;
			} else ix+=letstep;
		}
	}
	dp[dd]->np = npnts;
	dp[dd]->miss = newmiss;
	dp[dd]->xv = newx;
	dp[dd]->yv = newy;
	var_free_local();
	var_clear_local();
	var_set_global();
}
struct key_struct {
	char lstyle[9];
	long color,fill;
	int marker;
	double msize,lwidth;
	char *descrip;
};
extern struct key_struct *kd[100];
int draw_key(int nkd, double koffsetx, double koffsety, char *kpos
	,double khei, int knobox);
gdraw_key(double koffsetx, double koffsety, char *kpos, double khei, int knobox)
{
	int dt,dn,i,nkd;

	nkd = 0;
	for (dn=1;dn<=ndata;dn++) {
	  if (dp[dn]!=NULL) if (dp[dn]->key_name!=NULL) {
		kd[++nkd] = myallocz(sizeof(*kd[1]));
		kd[nkd]->fill = dp[dn]->key_fill;
		kd[nkd]->color = dp[dn]->color;
		kd[nkd]->lwidth = dp[dn]->lwidth;
		kd[nkd]->marker = dp[dn]->marker;
		kd[nkd]->msize = dp[dn]->msize;
		strcpy(kd[nkd]->lstyle,dp[dn]->lstyle);
		if (kd[nkd]->lstyle[0]==0) if (dp[dn]->line==true)
			kd[nkd]->lstyle[0]='1';
		mystrcpy(&kd[nkd]->descrip,dp[dn]->key_name);
	  }
	}
	draw_key(nkd,koffsetx,koffsety,kpos,khei,knobox);
	for (i=1;i<=nkd;i++) {
		if (kd[i]->descrip != NULL) myfrees(kd[i]->descrip,"GKEY");
		myfrees(kd[i],"GKEY1");
	}
}
/*----------------------------------------------------------------*/
#define next_dn  (ct+=1,skipspace,atoi(tk[ct]+1))
do_dataset(int d)
{
	int ct=2;
	while (ct<=ntk)  {
	     kw("LINE") 	dp[d]->line = true;
	else kw("LSTYLE") 	next_str(dp[d]->lstyle);
	else kw("LWIDTH") 	dp[d]->lwidth = next_exp;
	else kw("MARKER") 	dp[d]->marker = next_marker;
	else kw("MDATA") 	dp[d]->mdata = next_dn;
	else kw("COLOR") 	dp[d]->color = next_color;
	else kw("KEYFILL") 	dp[d]->key_fill = next_color;
	else kw("MSIZE") 	dp[d]->msize = next_exp;
	else kw("MSCALE") 	dp[d]->mscale = next_exp;
	else kw("KEY") 		next_vquote(dp[d]->key_name);
	else kw("AUTOSCALE") 	dp[d]->autoscale = true;
	else kw("AUTO") 	dp[d]->autoscale = true;
	else kw("NOMISS") 	dp[d]->nomiss = true;
	else kw("NOMISSING") 	dp[d]->nomiss = true;
	else kw("BIGFILE") 	next_vquote(dp[d]->bigfile);
	else kw("SMOOTH") 	{
		dp[d]->smoothm = false;
		dp[d]->smooth = true;
		dp[d]->line = true;
	}
	else kw("SMOOTHM") 	{
		dp[d]->smoothm = true;
		dp[d]->smooth = true;
		dp[d]->line = true;
	}
	else kw("XMIN") 	dp[d]->xmin = next_exp;
	else kw("XMAX") 	dp[d]->xmax = next_exp;
	else kw("YMIN") 	dp[d]->ymin = next_exp;
	else kw("YMAX") 	dp[d]->ymax = next_exp;
#ifdef HERROR 		/* because it won't fit on the PC */
	else kw("HERR") {
		next_str(dp[d]->herrup);
		strcpy(dp[d]->herrdown, dp[d]->herrup);
	}
	else kw("HERRLEFT")	next_str(dp[d]->herrup);
	else kw("HERRRIGHT")	next_str(dp[d]->herrdown);
	else kw("HERRWIDTH")	dp[d]->herrwidth = next_exp;
#endif
	else kw("ERR") {
		next_str(dp[d]->errup);
		strcpy(dp[d]->errdown, dp[d]->errup);
	}
	else kw("ERRUP")	next_str(dp[d]->errup);
	else kw("ERRDOWN")	next_str(dp[d]->errdown);
	else kw("ERRWIDTH")	dp[d]->errwidth = next_exp;
	else gprint("Unrecognised GRAPH DN sub command {%s} \n ",tk[ct]);
	ct++;
	}
}
double graph_xgraph(double v)
{
	double vvv;
	if (graph_xmax == graph_xmin) return 0.0;
	if (xx[1].log==true) {
		vvv = graph_x1 +
			 (log10(v)-log10(graph_xmin))/(log10(graph_xmax)
			-log10(graph_xmin))
			* (graph_x2-graph_x1) ;
	} else {
		vvv = graph_x1 +  ((v-graph_xmin)/(graph_xmax-graph_xmin))
				* (graph_x2-graph_x1);
	}
	return vvv;
}
double graph_ygraph(double v)
{
	double vvv;
	if (graph_ymax == graph_ymin) return 0.0;
	if (xx[2].log==true) {
		vvv = graph_y1 +
			 (log10(v)-log10(graph_ymin))/(log10(graph_ymax)
			-log10(graph_ymin))
			* (graph_y2-graph_y1) ;
		return vvv;
	} else {
 		return graph_y1 +  ((v-graph_ymin)/(graph_ymax-graph_ymin)) * (graph_y2-graph_y1);
	}
}