/*----------------  INCLUDE FILES  ----------------*/

#include "common.h"


/*----------------  External Variables	---------------*/

extern struct GfxBase *GfxBase;
extern struct Window *win;
extern struct RastPort *rastport;
extern struct ViewPort *viewport;
extern AllPlotData *apd;

extern DOUBLE Evaluate(DOUBLE, DOUBLE);


/*----------------  Global User Parameters  ------------*/

LONG   plottype;	   /* How to plot graph, see menuselect.h */

static
LONG   xcenter = 160,
       ycenter = 100,	   /* Center of screen */

       xp, yp,
       vp, zplane,

       xp1, xp2, xp3, xp4, yp1, yp2, yp3, yp4; /* for hidden and filled */

static
DOUBLE aspectratio = .927,	 /* to correct y-axis */

       fartherest, nextfartherest,
       xbegin, xend, ybegin, yend,
       xspc, xinc,
       yspc, yinc,
       xminymin, xminymax, xmaxymin, xmaxymax,
       xyminmax[4],

       dtoplane,
       dmin, dmax,
       colorspread,

       RotMatrix[3][3],
       TempMatrix[3][3];

static
BOOL   plottinggraph,
       onscreen,
       backflag;


/*------------	Global Color Palette Parameters  -------------*/

LONG numcolors = 24,
     pennumbase = 4,
     pennummax = 27,

     bkgpen  = 0,
     borderpen = 1,
     axespen = 2,
     textpen = 3,
     highlightpen1 = 30,
     highlightpen2 = 31,

     bkgred = 0,
     bkggreen = 0,
     bkgblue = 0,	  /* Black background */
     axesred = 8,
     axesgreen = 0,
     axesblue = 0,	  /* Red axes */
     borderred = 10,
     bordergreen = 10,
     borderblue = 10,	  /* Medium gray border */
     textred = 8,
     textgreen = 0,
     textblue = 8,	  /* Purple text */
     graphfarred = 0,
     graphfargreen = 0,
     graphfarblue = 7,	 /* Blue graph far color */
     graphnearred = 15,
     graphneargreen = 15,
     graphnearblue = 15;  /* White graph near color */


/*--------------  Temp Raster Global Parameters  ------------*/

static struct TmpRas tmpras;
static struct AreaInfo areainfo;
static PLANEPTR myplane;
static UWORD areabuffer[15];


/*------------	IsStopPlot Function  ---------------*/

BOOL IsStopPlot(message)
  struct IntuiMessage *message;
{
  ULONG mclass;
  UWORD menunum;

  mclass = message->Class;
  menunum = message->Code;
  ReplyMsg((struct Message *)message);

  if(mclass == MENUPICK) {
     if( (MENUNUM(menunum) == MENUFUNCTION) &&
	 (ITEMNUM(menunum) == FUNCTIONSTOPPLOT) )
       return TRUE;
  } /* if */

  return FALSE;

} /* IsStopPlot */


/*------------	SetScreenColors Function  --------------*/

VOID SetScreenColors()
{
SetRGB4(viewport,bkgpen,bkgred,bkggreen,bkgblue);
SetRGB4(viewport,borderpen,borderred,bordergreen,borderblue);
SetRGB4(viewport,axespen,axesred,axesgreen,axesblue);
SetRGB4(viewport,textpen,textred,textgreen,textblue);
SetRGB4(viewport,highlightpen1,0,4,5);
SetRGB4(viewport,highlightpen2,12,12,0);

}  /* SetScreenColors */


/*------------	SetGraphColors Function  --------------*/

VOID SetGraphColors()
{
  LONG num,
       redvalue,
       greenvalue,
       bluevalue;

  for (num = pennumbase ; num <= pennummax ; ++num) {
    redvalue = graphfarred + (num-pennumbase) * (graphnearred-graphfarred) / numcolors;
    greenvalue = graphfargreen + (num-pennumbase) * (graphneargreen-graphfargreen) / numcolors;
    bluevalue = graphfarblue + (num-pennumbase) * (graphnearblue-graphfarblue) / numcolors;
    SetRGB4(viewport,num,redvalue,greenvalue,bluevalue);
  }
}   /* SetGraphColors */


/*-------------  SetLoresColors Function  --------------*/

void SetLoresColors()
{
  SetScreenColors();

  /* Set fartherest colors as group */
  numcolors = 8;
  pennumbase = 4;
  pennummax = 11;
  graphfarred = 0;
  graphfargreen = 0;
  graphfarblue = 7;
  graphnearred = 0;
  graphneargreen = 0;
  graphnearblue = 15;

  SetGraphColors();

  /* Set nearest colors as group */
  numcolors = 16;
  pennumbase = 12;
  pennummax = 27;
  graphfarred = 0;
  graphfargreen = 0;
  graphfarblue = 15;
  graphnearred = 15;
  graphneargreen = 15;
  graphnearblue = 15;

  SetGraphColors();

  /* Set values to actual range for calculations */
  numcolors = 24;
  pennumbase = 4;

} /* SetLoresColors */


/*-------------  SetHiresColors Function  --------------*/

void SetHiresColors()
{
  SetScreenColors();

  numcolors = 12;
  pennumbase = 4;
  pennummax = 15;
  graphfarred = 0;
  graphfargreen = 0;
  graphfarblue = 15;
  graphnearred = 15;
  graphneargreen = 15;
  graphnearblue = 15;

  SetGraphColors();

} /* SetHiresColors */


/*-------------  ChangePlotParameters  Function  -------------*/

void ChangePlotParameters(scr)
  struct Screen *scr;
{
  WORD scrW, scrH;

  Forbid();
     scrW = scr->Width;
     scrH = scr->Height;
  Permit();

  if(scrW == 320) {
     apd->dd.OriginX = xcenter = 160;

     if(scrH == 200) {
	aspectratio = .927;
	apd->dd.OriginY = ycenter = 100;
     }
     else {
	aspectratio = 1.854;
	apd->dd.OriginY = ycenter = 200;
     }

     SetLoresColors();

  } /* if */

  else {  /* scrH == 640 */
     apd->dd.OriginX = xcenter = 320;

     if(scrH == 200) {
	aspectratio = .4635;
	apd->dd.OriginY = ycenter = 100;
     }
     else {
	aspectratio = .927;
	apd->dd.OriginY = ycenter = 200;
     }

     SetHiresColors();

  } /* else */

} /* ChangePlotParameters */


/*-----------------  calc Function  ------------------*/

static DOUBLE calc(x,y)
  DOUBLE x, y;
{
  return (Evaluate(x, y));

}    /* calc */


/*----------------  xcalc, ycalc, zcalc Functions  ----------------*/

static DOUBLE xcalc(x,y,z)
  DOUBLE x, y, z;
{
  return (RotMatrix[0][0] * x + RotMatrix[0][1] * y + RotMatrix[0][2] * z);

}  /* xcalc */


static DOUBLE ycalc(x,y,z)
  DOUBLE x, y, z;
{
  return (RotMatrix[1][0] * x + RotMatrix[1][1] * y + RotMatrix[1][2] * z);

}  /* ycalc */


static DOUBLE zcalc(x,y,z)
  DOUBLE x, y, z;
{
  return (RotMatrix[2][0] * x + RotMatrix[2][1] * y + RotMatrix[2][2] * z);

} /* zcalc */


/*--------------------	Calculate Function  -------------------*/

static VOID Calculate(x, y, z)
  DOUBLE x, y, z;
{
  DOUBLE xpixels, ypixels, zpixels,
	 zdiff;
  DOUBLE tforline;

  xpixels = (xcalc(x,y,z) * apd->dd.Scale) + (DOUBLE)(apd->dd.OriginX - xcenter);
  ypixels = (-ycalc(x,y,z) * apd->dd.Scale) + (DOUBLE)(apd->dd.OriginY - ycenter);
  zpixels = (zcalc(x,y,z) * apd->dd.Scale);

  zdiff = (DOUBLE)(vp + zplane) - zpixels; /* negative if behind viewer */

  if (zdiff > 1) {
    tforline = ((DOUBLE) vp) / zdiff;
    xp = (LONG)(xpixels * tforline + .5) + xcenter;
    yp = (LONG)(ypixels * tforline * aspectratio + .5) + ycenter;
    onscreen = TRUE;
  }
  else {
    xp = -1;
    yp = -1;
    onscreen = FALSE;
  }
}   /* Calculate */


/*------------------  CalcAndPlot Function  ----------------*/

static VOID CalcAndPlot(x, y)
  DOUBLE x, y;
{
  LONG	 pennum;
  DOUBLE z;

  z = calc(x,y);
  Calculate(x, y, z);

  if (onscreen == TRUE) {
    dtoplane = zcalc(x,y,z);
    pennum = (LONG)((dtoplane - dmin)/colorspread) + pennumbase;
    if (pennum > pennummax)  pennum = pennummax;
    if (pennum < pennumbase)  pennum = pennumbase;

      SetAPen(rastport,pennum);

      if (backflag == TRUE) {
	Move(rastport,xp,yp);
	WritePixel(rastport,xp,yp);
	backflag = FALSE;
      }
      else {
	Draw(rastport,xp,yp);
      }

  }  /* if */
}   /* CalcAndPlot */


/*------------------  CalcAndPlot2 Function  ----------------*/

static VOID CalcAndPlot2(x1, y1, flag)
  DOUBLE x1, y1;
  BOOL	 flag;
{
  static LONG	pennum;
  static DOUBLE x2, y2, z1;

  /* Assign last 2 points to first 2 points */
  xp1 = xp4;
  xp2 = xp3;
  yp1 = yp4;
  yp2 = yp3;

  x2 = x1 + xspc;  /* used if flag = FALSE, plot along Y-axis */
  y2 = y1 + yspc;  /* used if flag = TRUE , plot along X-axis */

  if(flag)
     Calculate(x1, y2, calc(x1, y2));
  else
     Calculate(x2, y1, calc(x2, y1));

  xp3 = xp;
  yp3 = yp;

  z1 = calc(x1, y1);

  Calculate(x1, y1, z1);
  xp4 = xp;
  yp4 = yp;

  if (onscreen == TRUE) {
    dtoplane = zcalc(x1, y1, z1);
    pennum = (LONG)((dtoplane - dmin)/colorspread) + pennumbase;
    if (pennum > pennummax)  pennum = pennummax;
    if (pennum < pennumbase)  pennum = pennumbase;

      SetOPen(rastport,pennum);

      if(plottype == PLOTHIDDEN)
	SetAPen(rastport, bkgpen);
      else /* PLOTFILLED */
	SetAPen(rastport, pennum);

      AreaMove(rastport, xp1, yp1);
      AreaDraw(rastport, xp2, yp2);
      AreaDraw(rastport, xp3, yp3);
      AreaDraw(rastport, xp4, yp4);
      AreaEnd(rastport);

  }  /* if */
}   /* CalcAndPlot2 */


/*------------------  AxesCalcAndPlot Function	----------------*/

static VOID AxesCalcAndPlot(x, y, z)
  DOUBLE x, y, z;
{
  LONG pixelcolor, potencolor;

  Calculate(x, y, z);

  if (onscreen == TRUE) {
    pixelcolor = ReadPixel(rastport,xp,yp);
    dtoplane = zcalc(x,y,z);
    potencolor = (LONG)((dtoplane - dmin)/colorspread) + pennumbase;

    if (potencolor > pixelcolor) {
      SetAPen(rastport,axespen);
      WritePixel(rastport,xp,yp);
    }
  }
}  /* AxesCalcAndPlot */


/*-----------------  RotateAboutX Function  ----------------*/

void RotateAboutX(angle)
  DOUBLE angle;
{
  DOUBLE sinx = sin(angle);
  DOUBLE cosx = cos(angle);

  RotMatrix[0][0] = TempMatrix[0][0];
  RotMatrix[0][1] = TempMatrix[0][1]*cosx+TempMatrix[0][2]*sinx;
  RotMatrix[0][2] = TempMatrix[0][1]*(-sinx)+TempMatrix[0][2]*cosx;
  RotMatrix[1][0] = TempMatrix[1][0];
  RotMatrix[1][1] = TempMatrix[1][1]*cosx+TempMatrix[1][2]*sinx;
  RotMatrix[1][2] = TempMatrix[1][1]*(-sinx)+TempMatrix[1][2]*cosx;
  RotMatrix[2][0] = TempMatrix[2][0];
  RotMatrix[2][1] = TempMatrix[2][1]*cosx+TempMatrix[2][2]*sinx;
  RotMatrix[2][2] = TempMatrix[2][1]*(-sinx)+TempMatrix[2][2]*cosx;
} /* RotateAboutX */


/*-----------------  RotateAboutY Function  ----------------*/

void RotateAboutY(angle)
  DOUBLE angle;
{
  DOUBLE siny = sin(angle);
  DOUBLE cosy = cos(angle);

  RotMatrix[0][0] = TempMatrix[0][0]*cosy+TempMatrix[0][2]*(-siny);
  RotMatrix[0][1] = TempMatrix[0][1];
  RotMatrix[0][2] = TempMatrix[0][0]*siny+TempMatrix[0][2]*cosy;
  RotMatrix[1][0] = TempMatrix[1][0]*cosy+TempMatrix[1][2]*(-siny);
  RotMatrix[1][1] = TempMatrix[1][1];
  RotMatrix[1][2] = TempMatrix[1][0]*siny+TempMatrix[1][2]*cosy;
  RotMatrix[2][0] = TempMatrix[2][0]*cosy+TempMatrix[2][2]*(-siny);
  RotMatrix[2][1] = TempMatrix[2][1];
  RotMatrix[2][2] = TempMatrix[2][0]*siny+TempMatrix[2][2]*cosy;
} /* RotateAboutY */


/*-----------------  RotateAboutZ Function  ----------------*/

void RotateAboutZ(angle)
  DOUBLE angle;
{
  DOUBLE sinz = sin(angle);
  DOUBLE cosz = cos(angle);

  RotMatrix[0][0] = TempMatrix[0][0]*cosz+TempMatrix[0][1]*sinz;
  RotMatrix[0][1] = TempMatrix[0][0]*(-sinz)+TempMatrix[0][1]*cosz;
  RotMatrix[0][2] = TempMatrix[0][2];
  RotMatrix[1][0] = TempMatrix[1][0]*cosz+TempMatrix[1][1]*sinz;
  RotMatrix[1][1] = TempMatrix[1][0]*(-sinz)+TempMatrix[1][1]*cosz;
  RotMatrix[1][2] = TempMatrix[1][2];
  RotMatrix[2][0] = TempMatrix[2][0]*cosz+TempMatrix[2][1]*sinz;
  RotMatrix[2][1] = TempMatrix[2][0]*(-sinz)+TempMatrix[2][1]*cosz;
  RotMatrix[2][2] = TempMatrix[2][2];
} /* RotateAboutZ */


/*-----------------  CopyRotToTemp Function  ----------------*/

void CopyRotToTemp()
{
  LONG i, j;

  for(i=0 ; i<=2 ; i++) {
     for(j=0 ; j<=2 ; j++) {
	TempMatrix[i][j] = RotMatrix[i][j];
     }
  }
}  /* CopytRotToTemp */


/*-----------------  PutIdentityInTemp Function  ----------------*/

void PutIdentityInTemp()
{
  TempMatrix[0][0] = 1;
  TempMatrix[0][1] = 0;
  TempMatrix[0][2] = 0;
  TempMatrix[1][0] = 0;
  TempMatrix[1][1] = 1;
  TempMatrix[1][2] = 0;
  TempMatrix[2][0] = 0;
  TempMatrix[2][1] = 0;
  TempMatrix[2][2] = 1;
}  /* PutIdentityInTemp */


/*--------------------	DrawGraph Function  -----------------*/

VOID DrawGraph()
{
  struct IntuiMessage *message;

  DOUBLE calcspacex, calcspacey, calcincx, calcincy,
	 tx, ty, tz,
	 yminend, ymaxend, ybeginspc,
	 xminend, xmaxend, xbeginspc,
	 zradians, yradians, xradians;
  BOOL	 alongxflag,
	 alongyflag;


  plottinggraph = TRUE;

  SetRast(rastport,0);

  /* Calculate Equation Parameter */

  vp = (LONG)(apd->dd.ViewDist * apd->dd.Scale); /* Viewer distance in pixels from projection plane */

  zplane = (LONG)(apd->dd.ProjPlane * apd->dd.Scale); /* Distance in pixels of projection plane from origin */

  calcspacex = apd->dd.LineSpacingX;
  calcspacey = apd->dd.LineSpacingY;

  if(plottype == PLOTNORMAL) {
     calcincx = 5*apd->dd.PlotPrecisionX;
     calcincy = 5*apd->dd.PlotPrecisionY;
  }
  else {
     calcincx = apd->dd.LineSpacingX;
     calcincy = apd->dd.LineSpacingY;
  }	 /* Values to speed up min and max calcs to get colorspread */

  /* Calculate Equation Coefficients */

  zradians = ((DOUBLE) apd->dd.RotationZ) * PI / 180.;
  yradians = ((DOUBLE) apd->dd.RotationY) * PI / 180.;
  xradians = ((DOUBLE) apd->dd.RotationX) * PI / 180.;

  PutIdentityInTemp();
  if((apd->pi.RotationOrder == ROTATEXYZ) || (apd->pi.RotationOrder == ROTATEXZY)) {
     RotateAboutX(xradians);
     CopyRotToTemp();
     if(apd->pi.RotationOrder == ROTATEXYZ) {
	RotateAboutY(yradians);
	CopyRotToTemp();
	RotateAboutZ(zradians);
     }
     else {
	RotateAboutZ(zradians);
	CopyRotToTemp();
	RotateAboutY(yradians);
     }
  }
  else if((apd->pi.RotationOrder == ROTATEYXZ) || (apd->pi.RotationOrder == ROTATEYZX)) {
     RotateAboutY(yradians);
     CopyRotToTemp();
     if(apd->pi.RotationOrder == ROTATEYXZ) {
	RotateAboutX(xradians);
	CopyRotToTemp();
	RotateAboutZ(zradians);
     }
     else {
	RotateAboutZ(zradians);
	CopyRotToTemp();
	RotateAboutX(xradians);
     }
  }
  else {
     RotateAboutZ(zradians);
     CopyRotToTemp();
     if(apd->pi.RotationOrder == ROTATEZXY) {
	RotateAboutX(xradians);
	CopyRotToTemp();
	RotateAboutY(yradians);
     }
     else {
	RotateAboutY(yradians);
	CopyRotToTemp();
	RotateAboutX(xradians);
     }
  }


  /* CALCULATE MINIMUM & MAXIMUM DISTANCES FROM USER */

  /* Initialize minimum & maximum */

  tz = calc(apd->dd.PlotXmin,apd->dd.PlotYmin);
  dtoplane = zcalc(apd->dd.PlotXmin,apd->dd.PlotYmin,tz);
  dmax = dtoplane;
  dmin = dtoplane;

  /* Loop on X at constant Y's */

  if ( (apd->pi.Surface == XONLY) || (apd->pi.Surface == XANDY) ) {

    for (ty = apd->dd.PlotYmin ; ty <= apd->dd.PlotYmax ; ty += calcspacey) {

	      /* Check for Stop Plot from user */
      if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
      } /* if */

      for (tx = apd->dd.PlotXmin ; tx <= apd->dd.PlotXmax ; tx += calcincx) {
	tz = calc(tx,ty);
	dtoplane = zcalc(tx,ty,tz);
	dmax = MAX(dmax,dtoplane);
	dmin = MIN(dmin,dtoplane);
      }
    }

  }

  /* Loop on Y at constant X's */

  if ( (apd->pi.Surface == YONLY) || (apd->pi.Surface == XANDY) ) {

    for (tx = apd->dd.PlotXmin ; tx <= apd->dd.PlotXmax ; tx += calcspacex) {

	     /* Check for Stop Plot from user */
      if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
      } /* if */

      for (ty = apd->dd.PlotYmin ; ty <= apd->dd.PlotYmax ; ty += calcincy) {
	tz = calc(tx,ty);
	dtoplane = zcalc(tx,ty,tz);
	dmax = MAX(dmax,dtoplane);
	dmin = MIN(dmin,dtoplane);
      }
    }

  }

  /* Calculate depth of each color range */

  colorspread = (dmax - dmin)/numcolors;
  if (colorspread == 0) colorspread = 1;  /* in case all points equidistant */

  /* Calculate two farthest corners from viewer */

  xyminmax[0] = xminymin = zcalc(apd->dd.PlotXmin, apd->dd.PlotYmin, 0);
  xyminmax[1] = xminymax = zcalc(apd->dd.PlotXmin, apd->dd.PlotYmax, 0);
  xyminmax[2] = xmaxymin = zcalc(apd->dd.PlotXmax, apd->dd.PlotYmin, 0);
  xyminmax[3] = xmaxymax = zcalc(apd->dd.PlotXmax, apd->dd.PlotYmax, 0);

  dqsort(xyminmax, 4);       /* Sort the distances */
  fartherest = xyminmax[0];
  nextfartherest = xyminmax[1];

  /* Set plotting order along X-axis */

  if( (fartherest == xminymin) || (fartherest == xminymax) ) {
     xbegin = apd->dd.PlotXmin;
     xend = apd->dd.PlotXmax;
     xspc = apd->dd.LineSpacingX;
     xinc = apd->dd.PlotPrecisionX;

     if( (nextfartherest == xmaxymin) || (nextfartherest == xmaxymax) )
	alongxflag = TRUE;
     else
	alongxflag = FALSE;

  }  /* if */

  else {
     xbegin = apd->dd.PlotXmax;
     xend = apd->dd.PlotXmin;
     xspc = -apd->dd.LineSpacingX;
     xinc = -apd->dd.PlotPrecisionX;

     if( (nextfartherest == xminymin) || (nextfartherest == xminymax) )
	alongxflag = TRUE;
     else
	alongxflag = FALSE;

  }  /* else */

  /* Set plotting order along Y-axis */

  if( (fartherest == xminymin) || (fartherest == xmaxymin) ) {
     ybegin = apd->dd.PlotYmin;
     yend = apd->dd.PlotYmax;
     yspc = apd->dd.LineSpacingY;
     yinc = apd->dd.PlotPrecisionY;
  }
  else {
     ybegin = apd->dd.PlotYmax;
     yend = apd->dd.PlotYmin;
     yspc = -apd->dd.LineSpacingY;
     yinc = -apd->dd.PlotPrecisionY;
  }

  /* PLOT THE GRAPH */

  switch(plottype) {

     case PLOTNORMAL:

	 /* Plot X's at constant Y's for specified range of Y */

  if (apd->pi.Surface == XONLY) {   /* Only if user specified */
    for (ty = ybegin ; (ty <= apd->dd.PlotYmax) && (ty >= apd->dd.PlotYmin) ; ty += yspc) {

		/* Check for Stop Plot from user */
      if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
      } /* if */

      backflag = TRUE;	/* Prevent last point from connecting with the next first point */

      for (tx = xbegin ; (tx <= apd->dd.PlotXmax) && (tx >= apd->dd.PlotXmin) ; tx += xinc)
	 CalcAndPlot(tx, ty);
    } /* for */
  }  /* if */

       /* Plot Y's at constant X for specified range of X */

  else if (apd->pi.Surface == YONLY) {      /* Do only if user specified */
    for (tx = xbegin ; (tx <= apd->dd.PlotXmax) && (tx >= apd->dd.PlotXmin) ; tx += xspc) {

	       /* Check for Stop Plot from user */
      if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
      } /* if */

      backflag = TRUE;	/* Prevent last point from connecting with the next first point */

      for (ty = ybegin ; (ty <= apd->dd.PlotYmax) && (ty >= apd->dd.PlotYmin) ; ty += yinc)
	 CalcAndPlot(tx, ty);
    } /* for */
  }  /* else if */

  else {   /* Plot crosshatched */
     DOUBLE ty2 = ybegin,
	    tx2 = xbegin;

     alongxflag = TRUE;
     alongyflag = TRUE;

     do {
	if(alongxflag) {
	   backflag = TRUE;
	   for(tx = xbegin ; (tx <= apd->dd.PlotXmax) && (tx >= apd->dd.PlotXmin) ; tx += xinc)
	      CalcAndPlot(tx, ty2);

	   ty2 += yspc;

	   if( (ty2 < apd->dd.PlotYmin) || (ty2 > apd->dd.PlotYmax) ) alongxflag = FALSE;

		   /* Check for Stop Plot from user */
	   if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	      if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
	   } /* if */

	} /* if */

	if(alongyflag) {
	   backflag = TRUE;
	   for(ty = ybegin ; (ty <= apd->dd.PlotYmax) && (ty >= apd->dd.PlotYmin) ; ty += yinc)
	      CalcAndPlot(tx2, ty);

	   tx2 += xspc;

	   if( (tx2 < apd->dd.PlotXmin) || (tx2 > apd->dd.PlotXmax) ) alongyflag = FALSE;

		     /* Check for Stop Plot from user */
	   if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	      if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
	   } /* if */

	} /* if */

     } while (alongxflag || alongyflag);

  } /* else */

  break;

     case PLOTHIDDEN:
     case PLOTFILLED:

#if DEBUG
  printf("Pre-InitArea\n");
#endif

	InitArea(&areainfo, &areabuffer[0], 6);

#if DEBUG
  printf("InitArea\n");
#endif

	rastport->AreaInfo = &areainfo;

#if DEBUG
  printf("rastport->AreaInfo\n");
#endif

	myplane = (PLANEPTR)AllocRaster(320,200);

#if DEBUG
  printf("AllocRaster\n");
#endif

	if(myplane == NULL) break;

#if DEBUG
  printf("myplane == NULL\n");
#endif

	rastport->TmpRas = (struct TmpRas *)InitTmpRas(&tmpras, myplane, RASSIZE(320,200));

#if DEBUG
  printf("Allocated tmpras\n");
#endif

     if(alongxflag) {  /* Plot along X-axis */
	yminend = apd->dd.PlotYmin;
	ymaxend = apd->dd.PlotYmax;

		/* Prevent going one yspc to far at end of plot */
	if(ybegin < yend)
	   ymaxend = apd->dd.PlotYmax - yspc;
	else
	   yminend = apd->dd.PlotYmin - yspc;

	for(ty = ybegin ; (ty <= ymaxend) && (ty >= yminend) ; ty += yspc) {

		    /* Check for Stop Plot from user */
	   if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	      if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
	   } /* if */

	   /* Initialize starting plot points */

	   ybeginspc = ty + yspc;

	   Calculate(xbegin, ty, calc(xbegin, ty));
	      xp4 = xp;
	      yp4 = yp;

	   Calculate(xbegin, ybeginspc, calc(xbegin, ybeginspc));
	      xp3 = xp;
	      yp3 = yp;

	   /* Plot along X-axis */

	   for(tx = xbegin + xspc ; (tx <= apd->dd.PlotXmax) && (tx >= apd->dd.PlotXmin) ; tx += xspc)
	      CalcAndPlot2(tx, ty, alongxflag);

	} /* for */
     } /* if */

     else {  /* Plot along Y-axis */
	xminend = apd->dd.PlotXmin;
	xmaxend = apd->dd.PlotXmax;

		/* Prevent going one xspc to far at end of plot */
	if(xbegin < xend)
	   xmaxend = apd->dd.PlotXmax - xspc;
	else
	   xminend = apd->dd.PlotXmin - xspc;

	for(tx = xbegin ; (tx <= xmaxend) && (tx >= xminend) ; tx += xspc) {

		    /* Check for Stop Plot from user */
	   if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	      if(IsStopPlot(message) == TRUE) goto EndDrawGraph;
	   } /* if */

	   /* Initialize starting plot points */

	   xbeginspc = tx + xspc;

	   Calculate(tx, ybegin, calc(tx, ybegin));
	      xp4 = xp;
	      yp4 = yp;

	   Calculate(xbeginspc, ybegin, calc(xbeginspc, ybegin));
	      xp3 = xp;
	      yp3 = yp;

	   /* Plot along X-axis */

	   for(ty = ybegin + yspc ; (ty <= apd->dd.PlotYmax) && (ty >= apd->dd.PlotYmin) ; ty += yspc)
	      CalcAndPlot2(tx, ty, alongxflag);

	} /* for */
     } /* else */

	FreeRaster(myplane, 320, 200);

	break;

     default:
	break;
  }  /* switch */


  /* DRAW AXES */

  switch(apd->pi.AxesType) {

     case AXESTYPENONE:
	break;

     case AXESTYPESTAR:

  /* X-axis */

  if (apd->ad.AxesXmax > apd->ad.AxesXmin) {   /* Do not plot if equal */
    ty = 0;
    tz = 0;
    for (tx = apd->ad.AxesXmin ; tx <= apd->ad.AxesXmax ; tx += apd->ad.AxesPrecision)  AxesCalcAndPlot(tx, ty, tz);
  }

  /* Y-axis */

  if (apd->ad.AxesYmax > apd->ad.AxesYmin) {
    tx = 0;
    tz = 0;
    for (ty = apd->ad.AxesYmin ; ty <= apd->ad.AxesYmax ; ty += apd->ad.AxesPrecision)  AxesCalcAndPlot(tx, ty, tz);
  }

  /* Z-axis */

  if (apd->ad.AxesZmax > apd->ad.AxesZmin) {
    tx = 0;
    ty = 0;
    for (tz = apd->ad.AxesZmin ; tz <= apd->ad.AxesZmax ; tz += apd->ad.AxesPrecision)  AxesCalcAndPlot(tx, ty, tz);
  }

  break;

     case AXESTYPEBOX:

  /* X-axis */

  if (apd->ad.AxesXmax > apd->ad.AxesXmin) {   /* Do not plot if equal */
    for (tx = apd->ad.AxesXmin ; tx <= apd->ad.AxesXmax ; tx += apd->ad.AxesPrecision) {
       AxesCalcAndPlot(tx, apd->ad.AxesYmin, apd->ad.AxesZmin);
       AxesCalcAndPlot(tx, apd->ad.AxesYmin, apd->ad.AxesZmax);
       AxesCalcAndPlot(tx, apd->ad.AxesYmax, apd->ad.AxesZmin);
       AxesCalcAndPlot(tx, apd->ad.AxesYmax, apd->ad.AxesZmax);
    } /* for */
  }   /* if */

  /* Y-axis */

  if (apd->ad.AxesYmax > apd->ad.AxesYmin) {
    for (ty = apd->ad.AxesYmin ; ty <= apd->ad.AxesYmax ; ty += apd->ad.AxesPrecision) {
       AxesCalcAndPlot(apd->ad.AxesXmin, ty, apd->ad.AxesZmin);
       AxesCalcAndPlot(apd->ad.AxesXmin, ty, apd->ad.AxesZmax);
       AxesCalcAndPlot(apd->ad.AxesXmax, ty, apd->ad.AxesZmin);
       AxesCalcAndPlot(apd->ad.AxesXmax, ty, apd->ad.AxesZmax);
    } /* for */
  }  /* if */

  /* Z-axis */

  if (apd->ad.AxesZmax > apd->ad.AxesZmin) {
    for (tz = apd->ad.AxesZmin ; tz <= apd->ad.AxesZmax ; tz += apd->ad.AxesPrecision) {
       AxesCalcAndPlot(apd->ad.AxesXmin, apd->ad.AxesYmin, tz);
       AxesCalcAndPlot(apd->ad.AxesXmin, apd->ad.AxesYmax, tz);
       AxesCalcAndPlot(apd->ad.AxesXmax, apd->ad.AxesYmin, tz);
       AxesCalcAndPlot(apd->ad.AxesXmax, apd->ad.AxesYmax, tz);
    } /* for */
  }   /* if */

  break;

     default:
	break;

  } /* switch */

EndDrawGraph:
   plottinggraph = FALSE;

}  /* DrawGraph */


VOID DrawContour(VOID)
{
   struct IntuiMessage *message;
   WORD lowX, lowY;
   WORD highX, highY;
   LONG   pennum;
   DOUBLE tx, ty, tz;
   DOUBLE stepvalX, stepvalY;
   DOUBLE calcspacex, calcspacey;

   plottinggraph = TRUE;

   /* Set up ranges */
   lowX = 0;
   lowY = 0;
   highX = win->Width - 1;
   highY = win->Height - 1;
   stepvalX = (apd->dd.PlotXmax - apd->dd.PlotXmin) / highX;
   stepvalY = (apd->dd.PlotYmax - apd->dd.PlotYmin) / highY;
   calcspacex = apd->dd.LineSpacingX;
   calcspacey = apd->dd.LineSpacingY;

   /* CALCULATE MINIMUM & MAXIMUM DISTANCES FROM USER */

   /* Initialize minimum & maximum */
   tz = calc(apd->dd.PlotXmin,apd->dd.PlotYmin);
   dmax = tz;
   dmin = tz;

   /* Loop over entire screen */
   for (ty = apd->dd.PlotYmin ; ty <= apd->dd.PlotYmax ; ty += calcspacey) {
	      /* Check for Stop Plot from user */
      if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	if(IsStopPlot(message) == TRUE) goto EndDrawContour;
      } /* if */
      for (tx = apd->dd.PlotXmin ; tx <= apd->dd.PlotXmax ; tx += calcspacex) {
	tz = calc(tx,ty);
	dmax = MAX(dmax,tz);
	dmin = MIN(dmin,tz);
      } /* for */
   } /* for */

   /* Calculate depth of each color range */
   colorspread = (dmax - dmin)/numcolors;
   if (colorspread == 0) colorspread = 1;  /* in case all points equidistant */

   /* Plot over entire screen */
   for (ty = apd->dd.PlotYmax, yp = lowY ;
	yp <= highY ; ty -= stepvalY, yp += 1) {
	      /* Check for Stop Plot from user */
      if(message = (struct IntuiMessage *) GetMsg(win->UserPort)) {
	if(IsStopPlot(message) == TRUE) goto EndDrawContour;
      } /* if */

      for (tx = apd->dd.PlotXmin, xp = lowX ;
	   xp <= highX ; tx += stepvalX, xp += 1)    {
	 tz = calc(tx,ty);
	 pennum = (LONG)((tz - dmin)/colorspread) + pennumbase;
	 if (pennum > pennummax)  pennum = pennummax;
	 if (pennum < pennumbase)  pennum = pennumbase;
	 SetAPen(rastport,pennum);
	 WritePixel(rastport,xp,yp);
      } /* for */
   } /* for */

EndDrawContour:
   plottinggraph = FALSE;

} /* DrawContour */