 /*
  * Khoros: $Id: draw_3D.c,v 1.4 1992/03/20 22:46:28 dkhoros Exp $
  */

#if !defined(lint) && !defined(SABER)
static char rcsid[] = "Khoros: $Id: draw_3D.c,v 1.4 1992/03/20 22:46:28 dkhoros Exp $";
#endif

 /*
  * $Log: draw_3D.c,v $
 * Revision 1.4  1992/03/20  22:46:28  dkhoros
 * VirtualPatch5
 *
  */ 


/*
 *----------------------------------------------------------------------
 *
 * Copyright 1990, University of New Mexico.  All rights reserved.
 *
 * Permission to copy and modify this software and its documen-
 * tation only for internal use in your organization is hereby
 * granted, provided that this notice is retained thereon and
 * on all copies.  UNM makes no representations as too the sui-
 * tability and operability of this software for any purpose.
 * It is provided "as is" without express or implied warranty.
 * 
 * UNM DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
 * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FIT-
 * NESS.  IN NO EVENT SHALL UNM BE LIABLE FOR ANY SPECIAL,
 * INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY OTHER DAMAGES WHAT-
 * SOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
 * IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS
 * ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PER-
 * FORMANCE OF THIS SOFTWARE.
 * 
 * No other rights, including for example, the right to redis-
 * tribute this software and its documentation or the right to
 * prepare derivative works, are granted unless specifically
 * provided in a separate license agreement.
 *---------------------------------------------------------------------
 */

#include "unmcopyright.h"	 /* Copyright 1990 by UNM */
#include "X3D.h"	

#define FindColor(indx, colors) \
	{ while(colors[indx].flags == NULL && indx > 0) indx--; }

static void fill_colorplane(), fill_colorplane_elev();
static int find_level(), draw_colorline();

/* >>>>>>>>>>>>>>>>>>>>>>>>>>>>> <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
   >>>>                                                       <<<<
   >>>>	    file name: draw_3D.c			      <<<<
   >>>>               			      		      <<<<
   >>>>   description: 3D Drawing Utilities		      <<<<
   >>>>               			      		      <<<<
   >>>>      routines:					      <<<<
   >>>>               			      		      <<<<
   >>>>                 3D Drawing Utilities   	      	      <<<<
   >>>>                                                       <<<<
   >>>>			X3D_draw_polyline()		      <<<<
   >>>>			X3D_draw_lineplot()		      <<<<
   >>>>			X3D_draw_wireframe()		      <<<<
   >>>>			X3D_draw_bezier()		      <<<<
   >>>>                                                       <<<<
   >>>>			_X3D_blend_patch()		      <<<<
   >>>>			_X3D_draw_patch()		      <<<<
   >>>>			_X3D_fill_patch()		      <<<<
   >>>>                                                       <<<<
   >>>>			X3D_bezier_surface()		      <<<<
   >>>>			X3D_draw_impulse()		      <<<<
   >>>>			X3D_draw_polymarker()		      <<<<
   >>>>                                                       <<<<
   >>>>			draw_markers()			      <<<<
   >>>>                                                       <<<<
   >>>>			X3D_draw_polygon()		      <<<<
   >>>>			X3D_draw_segments()		      <<<<
   >>>>			X3D_draw_line()			      <<<<
   >>>>			X3D_draw_rectangle()		      <<<<
   >>>>                                                       <<<<
   >>>>			find_corner()			      <<<<
   >>>>                                                       <<<<
   >>>>			X3D_draw_mesh()	   		      <<<<
   >>>>			X3D_draw_horizon()   		      <<<<
   >>>>			X3D_draw_surface()   		      <<<<
   >>>>                                                       <<<<
   >>>>			fill_colorplane()		      <<<<
   >>>>                                                       <<<<
   >>>>			X3D_draw_colormesh()   		      <<<<
   >>>>			X3D_draw_colorline()   		      <<<<
   >>>>                                                       <<<<
   >>>>			fill_colorline()		      <<<<
   >>>>			find_level()		  	      <<<<
   >>>>                                                       <<<<
   >>>>			X3D_draw_contour()   		      <<<<
   >>>>			X3D_get_contour_levels()	      <<<<
   >>>>                                                       <<<<
   >>>>			compute_contour()	  	      <<<<
   >>>>			display_contour()	  	      <<<<
   >>>>                                                       <<<<
   >>>> modifications:
   >>>>                                                       <<<<
   >>>>>>>>>>>>>>>>>>>>>>>>>>>>> <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< */




/************************************************************
*
*  MODULE NAME:  X3D_draw_polyline
*
*      PURPOSE:  Uses X11 graphics routines to connect a series of points
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a line (or a set of lines) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/


X3D_draw_polyline (id, coords, size, fg)

int	id, size;
Coord	*coords;
XColor  *fg;
{
	X3DGraphics	*graphics;
	XPoint		*points;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_polyline:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_polyline:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

        /*
	 *  Plot the desired polyline.
	 */
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	X3D_draw[graphics->device].lines(graphics, points,size,CoordModeOrigin);
        free(points);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_lineplot
*
*      PURPOSE:  Uses X11 graphics routines to connect a series
*		 of polylines to make a 3D plot.
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		row_size  - number of points in each polyline
*		size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a line (or a set of lines) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/


X3D_draw_lineplot(id, coords, row_size, size, fg)

int	id, row_size, size;
Coord	*coords;
XColor  *fg;
{
	int		i, j, col_size;
	X3DGraphics	*graphics;
	XPoint		*points;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_lineplot:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_lineplot:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

        /*
	 *  see how many segments are to be plotted.  If there is only one
	 *  segment, then plot the whole thing.
	 */
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	if (row_size == 1)
	{
	   X3D_draw[graphics->device].lines(graphics, points, size,
			CoordModeOrigin);
	}
	else
	{
	   col_size = size/row_size;

           /* for each segment, plot a separate curve */
	   for (i = 0; i < col_size; i++)
           {
	      j = i * row_size;
	      X3D_draw[graphics->device].lines(graphics, &points[j], row_size,
	      		CoordModeOrigin);
	   }
	}
        free(points);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_wireframe
*
*      PURPOSE:  Uses X11 graphics routines to connect a series
*		 of polylines to make a 3D plot in a mesh pattern.
*		 This is done by drawing both lines along both
*		 the columns and rows.
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		row_size  - number of points in each polyline
*		size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a line (or a set of lines) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/


X3D_draw_wireframe(id, coords, row_size, size, fg)

int	id, row_size, size;
Coord	*coords;
XColor  *fg;
{
	X3DGraphics	*graphics;
	int		i, j, l, col_size;
	XPoint		*points, *columns;

	col_size = size/row_size;
	if (col_size == 1 || row_size == 1)
	{
	   X3D_draw_lineplot(id, coords, row_size, size, fg);
	   return;
	}

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_wireframe:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_wireframe:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if (!(columns = (XPoint *) malloc ((unsigned)col_size*sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_wireframe:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

        /*
	 *  Draw the rows first
	 */
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	for (i = 0; i < col_size; i++)
        {
	    j = i * row_size;
	    X3D_draw[graphics->device].lines(graphics, &points[j], row_size,
		      CoordModeOrigin);
	}

        /*
	 *  Draw the cols last
	 */
	for (i = 0; i < row_size; i++)
        {
	    for (j = 0; j < col_size; j++)
	    {
	        l = j * row_size + i;
		columns[j] = points[l];
	    }
	    X3D_draw[graphics->device].lines(graphics, columns, col_size,
		         CoordModeOrigin);
	}
        free(points);
        free(columns);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_bezier
*
*      PURPOSE:  Uses X11 graphics routines to connect a series
*		 of polylines to make a 3D plot in a mesh pattern.
*		 This is done by drawing both lines along both
*		 the columns and rows.
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		row_size  - number of points in each polyline
*		size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a line (or a set of lines) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/

int bezier_step;

X3D_draw_bezier(id, pts, row_size, size, fg)

int	id, row_size, size;
Coord	*pts;
XColor  *fg;
{
	X3DGraphics	*graphics;
	Coord		control[4];
	int		i, col_size;
	float		unit_step, step = 0.0;
	void		_X3D_blend_patch(), _X3D_draw_patch();

	col_size = size/row_size;
	if (col_size != 4 && row_size != 1)
	{
	   fprintf (stderr,"X3D_draw_bezier:\n");
	   fprintf (stderr,"\t Only able to plot 4x4 plots (a patch)\n");
	   return;
	}

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_bezier:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	unit_step = 1.0/bezier_step;
	for (i = 0; i <= bezier_step; i++)
	{
	   _X3D_blend_patch(pts[0],pts[1],pts[2],pts[3],step, &control[0]);
	   _X3D_blend_patch(pts[4],pts[5],pts[6],pts[7],step, &control[1]);
	   _X3D_blend_patch(pts[8],pts[9],pts[10],pts[11],step, &control[2]);
	   _X3D_blend_patch(pts[12],pts[13],pts[14],pts[15],step, &control[3]);
	   _X3D_draw_patch(graphics, control[0], control[1], control[2],
			     control[3], bezier_step);

	   _X3D_blend_patch(pts[0],pts[4],pts[8],pts[12],step, &control[0]);
	   _X3D_blend_patch(pts[1],pts[5],pts[9],pts[13],step, &control[1]);
	   _X3D_blend_patch(pts[2],pts[6],pts[10],pts[14],step, &control[2]);
	   _X3D_blend_patch(pts[3],pts[7],pts[11],pts[15],step, &control[3]);
	   _X3D_draw_patch(graphics, control[0], control[1], control[2],
			     control[3], bezier_step);

	   step += unit_step;
	}
}



static void _X3D_blend_patch(pt0, pt1, pt2, pt3, step, control)

register float	step;
Coord	 pt0, pt1, pt2, pt3, *control;
{
	register float	t, t0, t1, t2, t3;

	t   = (1 - step);
	t0  = t*t*t;
	t1  = 3*step*t*t;
	t2  = 3*step*step*t;
	t3  = step*step*step;

	control->x = pt0.x*t0 + pt1.x*t1 + pt2.x*t2 + pt3.x*t3;
	control->y = pt0.y*t0 + pt1.y*t1 + pt2.y*t2 + pt3.y*t3;
	control->z = pt0.z*t0 + pt1.z*t1 + pt2.z*t2 + pt3.z*t3;
}



static void _X3D_draw_patch(graphics, pt0, pt1, pt2, pt3, bezier_step)

int		bezier_step;
X3DGraphics	*graphics;
Coord		pt0, pt1, pt2, pt3;
{
	int	i, size;
	Coord	*coords;
	XPoint	*points;
	float	unit_step, step = 0.0;


	size = bezier_step +1;
	if (!(coords = (Coord *) malloc(size * sizeof (Coord))))
	{
	   fprintf (stderr,"X3D_display_patch:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if (!(points = (XPoint *) malloc(size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_display_patch:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	unit_step = 1.0/bezier_step;
	for (i = 0; i <= bezier_step; i++)
	{
	   _X3D_blend_patch(pt0, pt1, pt2, pt3, step, &coords[i]);
	   step += unit_step;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(coords);
	   free(points);
	   return;
	}
	X3D_draw[graphics->device].lines(graphics, points, size,
			CoordModeOrigin);
	
	free(coords);
	free(points);
}



/************************************************************
*
*  MODULE NAME:  X3D_bezier_surface
*
*      PURPOSE:  Uses X11 graphics routines to connect a series
*		 of polylines to make a 3D plot in a mesh pattern.
*		 This is done by drawing both lines along both
*		 the columns and rows.
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		row_size  - number of points in each polyline
*		size      - number of points
*
*       OUTPUT: displays a line (or a set of lines) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/


X3D_bezier_surface(id, pts, row_size, size, levels, colors, num_levels,
		   wcmin, wcmax)

int	 id; 
Coord 	 *pts, wcmin, wcmax;
int      row_size, size, num_levels;
Real	 *levels;
XColor	 *colors;
{
	X3DGraphics	*graphics;
	Coord		control[4];
	int		i, col_size;
	float		unit_step, step = 0.0;
	void		_X3D_blend_patch(), _X3D_fill_patch();

	col_size = size/row_size;
	if (col_size != 4 && row_size != 1)
	{
	   fprintf (stderr,"X3D_bezier_surface:\n");
	   fprintf (stderr,"\t Only able to plot 4x4 plots (a patch)\n");
	   return;
	}

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_bezier_surface:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	unit_step = 1.0/bezier_step;
	_X3D_blend_patch(pts[0],pts[1],pts[2],pts[3],step, &control[0]);
	_X3D_blend_patch(pts[4],pts[5],pts[6],pts[7],step, &control[1]);
	_X3D_blend_patch(pts[8],pts[9],pts[10],pts[11],step, &control[2]);
	_X3D_blend_patch(pts[12],pts[13],pts[14],pts[15],step, &control[3]);

	_X3D_blend_patch(pts[0],pts[4],pts[8],pts[12],step, &control[4]);
	_X3D_blend_patch(pts[1],pts[5],pts[9],pts[13],step, &control[5]);
	_X3D_blend_patch(pts[2],pts[6],pts[10],pts[14],step, &control[6]);
	_X3D_blend_patch(pts[3],pts[7],pts[11],pts[15],step, &control[7]);
	step += unit_step;
	for (i = 1; i <= bezier_step; i++)
	{
	   _X3D_blend_patch(pts[0],pts[1],pts[2],pts[3],step, &control[8]);
	   _X3D_blend_patch(pts[4],pts[5],pts[6],pts[7],step, &control[9]);
	   _X3D_blend_patch(pts[8],pts[9],pts[10],pts[11],step, &control[10]);
	   _X3D_blend_patch(pts[12],pts[13],pts[14],pts[15],step, &control[11]);

	   _X3D_blend_patch(pts[0],pts[4],pts[8],pts[12],step, &control[12]);
	   _X3D_blend_patch(pts[1],pts[5],pts[9],pts[13],step, &control[13]);
	   _X3D_blend_patch(pts[2],pts[6],pts[10],pts[14],step, &control[14]);
	   _X3D_blend_patch(pts[3],pts[7],pts[11],pts[15],step, &control[15]);
	   _X3D_fill_patch(graphics, control, bezier_step, levels, colors,
			  num_levels, wcmin, wcmax);

	   for (i = 0; i < 8; i++)
	      control[i] = control[i + 8];

	   step += unit_step;
	}
}



static void _X3D_fill_patch(graphics, pts, bezier_step, levels, colors,
			   num_levels, wcmin, wcmax)


X3DGraphics	*graphics;
Coord 	 *pts, wcmin, wcmax;
int      bezier_step, num_levels;
Real	 *levels;
XColor   *colors;
{
	int	i, j, m, indx, size;
	Coord	*coords;
	float	unit_step, step;


	size = bezier_step +1;
	if (!(coords = (Coord *) malloc(size * size * sizeof (Coord))))
	{
	   fprintf (stderr,"X3D_display_patch:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	m = 0;
	unit_step = 1.0/bezier_step;
	for (i = 0; i <= bezier_step; i++)
	{
	   step = 0.0;
	   for (j = 0; j <= bezier_step; j++)
	   {
	       indx = i*size + j;
	       _X3D_blend_patch(pts[indx],pts[indx+1],pts[indx+2],pts[indx+3],
			       step, &coords[m++]);
	       step += unit_step;
	   }
	}
	X3D_draw_surface(graphics->id, coords, 4, 16, levels, colors,
			 num_levels, wcmin, wcmax, FALSE, TRUE);
	free(coords);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_impulse
*
*      PURPOSE:  Uses X11 graphics routines to connect a series of points
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a discrete line plot to the graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/

X3D_draw_impulse (id, coords, size, fg)

int	id, size;
Coord	*coords;
XColor  *fg;
{
	X3DGraphics	*graphics;
        XSegment        *segments;
	Coord		*impulse;

	int		i, j, num;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_discrete:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(segments = (XSegment *) malloc ((unsigned) size * 
			 sizeof(XSegment))))
	{
	   fprintf (stderr,"X3D_draw_discrete:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	num = size * 2;
	if (!(impulse = (Coord *) malloc ((unsigned) num * sizeof(Coord))))
	{
	   fprintf (stderr,"X3D_draw_discrete:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   free(segments); return;
	}
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);

	i = j = 0;
	while (i < size)
	{
	   impulse[j] = impulse[j+1] = coords[i++];
	   if(graphics->wc_min.z < 0.0 && graphics->wc_max.z > 0.0)
		impulse[j+1].z = 0.0;
	   else
	   	impulse[j+1].z = graphics->wc_min.z;
	   j += 2;
	}

	if (!(_X3D_convert_wc_to_dc_seg (graphics, impulse, segments, num)))
	{
	   free(impulse);
	   free(segments);
	   return;
	}

	X3D_draw[graphics->device].segments(graphics, segments, size);
	free(impulse);
	free(segments);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_polymarker
*
*      PURPOSE:  Uses X11 graphics routines to draw a series of points
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a point (or a set of points) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/

X3D_draw_polymarker (id, coords, size, marker_type, fg)

int	id, size, marker_type;
Coord	*coords;
XColor  *fg;
{
	X3DGraphics	*graphics;
	XPoint		*points;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_polymarker:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_polymarker:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}
	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	draw_markers(id, points, NULL, size, marker_type, NULL, 0);

        free(points);
}



draw_markers(id, points, coords, size, marker_type, colors, num_colors)

int	id, size, num_colors, marker_type;
Coord   *coords;
XPoint  *points;
XColor  *colors;
{
	register int	i, psize, scale;
	X3DGraphics	*graphics;
	XPoint		*draw;

	XPoint  triangle[4];
	XPoint  diamond[6];

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"draw_marker:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (graphics->device == X11)
	   scale = 1;
	else if (graphics->device == POSTSCR)
	   scale = -6;
	else if (graphics->device == IMPRESS)
	   scale = 4;
	else if (graphics->device == HPGL)
	   scale = 50;

	switch (marker_type)
	{
	   case MarkerPoint:	/* marker star or astrix */
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
				 points[i].x, points[i].y +3*scale,
				 points[i].x, points[i].y -3*scale);
		    X3D_draw[graphics->device].line(graphics,
				 points[i].x +3*scale, points[i].y +2*scale,
				 points[i].x -3*scale, points[i].y -2*scale);
		    X3D_draw[graphics->device].line(graphics,
				 points[i].x +3*scale, points[i].y -2*scale,
				 points[i].x -3*scale, points[i].y +2*scale);
		}
		break;

	   case MarkerSquare:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    if (graphics->device == POSTSCR || graphics->device == HPGL)
		       X3D_draw[graphics->device].rectangle(graphics,
			      points[i].x-2*scale, points[i].y + 2*scale,
			      4*scale, 4*scale);
		    else
		       X3D_draw[graphics->device].rectangle(graphics,
			      points[i].x-2*scale, points[i].y -2*scale,
			      4*scale, 4*scale);
		}
		break;

	   case MarkerTriangle:
		if (graphics->device == POSTSCR || graphics->device == HPGL)
		   scale = -scale;
	        triangle[1].x = triangle[3].x = -2 * scale;
		triangle[1].y = -4 * scale;
	        triangle[2].y = 0;
	        triangle[3].y = triangle[2].x =  4 * scale;

		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    triangle[0].x = points[i].x +2*scale;
		    triangle[0].y = points[i].y +2*scale;
		    X3D_draw[graphics->device].lines(graphics,
				  triangle, 4, CoordModePrevious);
		}
		break;

	   case MarkerDiamond:
		diamond[0].x = diamond[0].y = diamond[1].x = 0;
		diamond[2].y = diamond[3].x = diamond[3].y =
		diamond[4].x =  2*scale;
		diamond[1].y = diamond[2].x = diamond[4].y = diamond[5].x =
		diamond[5].y = -2*scale;
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    diamond[0] = points[i];
		    X3D_draw[graphics->device].lines(graphics, diamond, 6,
			          CoordModePrevious);
		}
		break;

	   case MarkerDagger:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
				 points[i].x -2*scale, points[i].y,
				 points[i].x +2*scale,points[i].y);
		    X3D_draw[graphics->device].line(graphics,
				 points[i].x, points[i].y +2*scale,
				 points[i].x, points[i].y -4*scale);
		}
		break;

	   case MarkerCross:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
				 points[i].x -3*scale, points[i].y,
				 points[i].x+3*scale, points[i].y);
		    X3D_draw[graphics->device].line(graphics,
				 points[i].x, points[i].y+3*scale,
				 points[i].x, points[i].y-3*scale);
		}
		break;

	   case MarkerX:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
				 points[i].x -2*scale, points[i].y -2*scale,
				 points[i].x +2*scale, points[i].y +2*scale);
		    X3D_draw[graphics->device].line(graphics,
				 points[i].x -2*scale, points[i].y +2*scale,
				 points[i].x +2*scale, points[i].y -2*scale);
		}
		break;

	   case MarkerV:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
			      points[i].x, points[i].y +2*scale,
			      points[i].x -2*scale, points[i].y -2*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x, points[i].y +2*scale,
			      points[i].x +2*scale, points[i].y -2*scale);
		}
		break;

	   case MarkerCaret:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
			      points[i].x, points[i].y -2*scale,
			      points[i].x -2*scale, points[i].y +2*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x, points[i].y -2*scale,
			      points[i].x +2*scale, points[i].y +2*scale);
		}
		break;

	   case MarkerBowTie:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -3*scale, points[i].y -3*scale,
			      points[i].x +3*scale, points[i].y +3*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -3*scale, points[i].y +3*scale,
			      points[i].x +3*scale, points[i].y -3*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -3*scale, points[i].y +3*scale,
			      points[i].x -3*scale, points[i].y -3*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x +3*scale, points[i].y +3*scale,
			      points[i].x +3*scale, points[i].y -3*scale);
		}
		break;

	   case MarkerBox:
		if (graphics->device == POSTSCR || graphics->device == HPGL)
		   scale = -scale;
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    if (graphics->device == POSTSCR || graphics->device == HPGL)
		       X3D_draw[graphics->device].rectangle(graphics,
			      points[i].x-2*scale, points[i].y + 2*scale,
			      4*scale, 4*scale);
		    else
		       X3D_draw[graphics->device].rectangle(graphics,
			      points[i].x-2*scale, points[i].y -2*scale,
			      4*scale, 4*scale);

		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -3*scale, points[i].y -3*scale,
			      points[i].x +3*scale, points[i].y +3*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -3*scale, points[i].y +3*scale,
			      points[i].x +3*scale, points[i].y -3*scale);
		}
		break;

	   case MarkerCircle:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    if (graphics->device == POSTSCR || graphics->device == HPGL)
		       X3D_draw[graphics->device].arc(graphics,
			     points[i].x-2*scale, points[i].y,
			     4*scale, 4*scale, 0, 360);
		    else
		       X3D_draw[graphics->device].arc(graphics,
			     points[i].x-2*scale, points[i].y -2*scale,
			     4*scale, 4*scale, 0, 360);
		}
		break;

	   case MarkerArc:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    if (graphics->device == POSTSCR || graphics->device == HPGL)
		       X3D_draw[graphics->device].arc(graphics,
			     points[i].x -3*scale, points[i].y,
			     6*scale, 6*scale, 0, 180);
		    else
		       X3D_draw[graphics->device].arc(graphics,
			     points[i].x -3*scale, points[i].y -3*scale,
			     6*scale, 6*scale, 0, 180);
		}
		break;

	   case MarkerHexagon:
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].line(graphics,
			      points[i].x, points[i].y +3*scale,
			      points[i].x +2*scale, points[i].y +1*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x +2*scale, points[i].y +1*scale,
			      points[i].x +2*scale, points[i].y -1*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x +2*scale, points[i].y -1*scale,
			      points[i].x, points[i].y -3*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x, points[i].y -3*scale,
			      points[i].x -2*scale, points[i].y -1*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -2*scale, points[i].y -1*scale,
			      points[i].x -2*scale, points[i].y +1*scale);
		    X3D_draw[graphics->device].line(graphics,
			      points[i].x -2*scale, points[i].y +1*scale,
			      points[i].x, points[i].y +3*scale);
		}
		break;

           case MarkerDot:
		if (graphics->device == POSTSCR || graphics->device == HPGL)
		   scale = 2;
		else
		   scale = 1;
		for (i = 0; i < size; i++)
		{
		    if (colors != NULL && coords != NULL)
		       X3D_draw[graphics->device].draw_color(graphics,
					&colors[coords[i].d], FALSE);

		    X3D_draw[graphics->device].arc(graphics,
			     points[i].x - scale, points[i].y,
			     2*scale, 2*scale, 0, 360);
		}
		break;
	}
}




/************************************************************
*
*  MODULE NAME:  X3D_draw_polygon
*
*      PURPOSE:  Uses X11 graphics routines to draw a polygon
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		size      - number of points
*		fg	  - foreground color
*		bg	  - background color
*
*       OUTPUT: displays a polygon to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/

X3D_draw_polygon (id, coords, size, fg, bg)

int	id, size;
Coord	*coords;
XColor  *fg, *bg;
{
	X3DGraphics	*graphics;
	XPoint		*points;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_polyline:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned)(size +1)*sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_polygon:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

        /*
	 *  Fill the desired polygon.
	 */
	X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
	X3D_draw[graphics->device].fill_polygon(graphics,points, size, Complex, 
		        CoordModeOrigin);

	if (points[0].x != points[size-1].x || points[0].y != points[size-1].y)
	   points[size++] = points[0];

        /*
	 *  Draw the desired polyline outline.
	 */
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	X3D_draw[graphics->device].lines(graphics, points, size,
			CoordModeOrigin);
        free(points);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_segments
*
*      PURPOSE:  Uses graphics routines to draw a series of segments
*
*        INPUT: id        - X3D graphics structure ID
*               coords    - array of world coordinate points
*               size      - number of points
*		fg	  - foreground color
*
*       OUTPUT: displays a set of lines to the graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/

X3D_draw_segments (id, coords, size, fg)

int     id, size;
Coord   *coords;
XColor  *fg;
{
        int             num;
        X3DGraphics     *graphics;
        XSegment        *segments;

        if (!(graphics = _X3D_get_graphics(id)))
        {
           fprintf (stderr,"X3D_draw_segments:\n");
           fprintf (stderr,"\t unknown graphics id\n");
           return;
        }

        num = size/2;
        if (!(segments = (XSegment *) malloc ((unsigned) num *
                         sizeof (XSegment))))
        {
           fprintf (stderr,"X3D_draw_segments:\n");
           fprintf (stderr,"\t not enough memory\n");
           return;
        }

        if(!(_X3D_convert_wc_to_dc_seg(graphics, coords, segments, size)))
	{
	   free(segments);
	   return;
	}
	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
        X3D_draw[graphics->device].segments(graphics, segments, num);
        free(segments);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_line
*
*      PURPOSE:  Uses X11 graphics routines to connect a two points
*
*        INPUT: id        - X3D graphics structure ID
*               coord1    - the beginning world coordinate point
*               coord2    - the end world coordinate point
*		fg	  - foreground color
*
*       OUTPUT: displays a line (or a set of lines) to graphics workstation
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/

X3D_draw_line (id, coord1, coord2, fg)

int     id;
Coord   coord1, coord2;
XColor  *fg;
{
        X3DGraphics     *graphics;
        XPoint          pt1, pt2;

        if (!(graphics = _X3D_get_graphics(id)))
        {
           fprintf (stderr,"X3D_draw_line:\n");
           fprintf (stderr,"\t unknown graphics id\n");
           return;
        }

        if(!(_X3D_convert_wc_to_dc(graphics, &coord1, &pt1, 1)))
	   return;
	
        if(!(_X3D_convert_wc_to_dc(graphics, &coord2, &pt2, 1)))
	   return;

	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
        X3D_draw[graphics->device].line(graphics, pt1.x, pt1.y, pt2.x, pt2.y);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_rectangle
*
*      PURPOSE:  Uses X11 graphics routines to create a
*		 rectangle in 3D space.  Therefore the user
*		 specifies a starting point, width, height,
*		 and depth of that plane.
*
*        INPUT: id        - X3D graphics structure ID
*		coord1    - the beginning world coordinate point
*		width,
*		height,
*		depth     - this is used to describe the vector
*			    of the desired plane to be displayed.
*		fg	  - foreground color
*
*       OUTPUT: displays a plane (as a set of lines) to the
*		graphics workspace
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young
*
*************************************************************/


X3D_draw_rectangle (id, coord1, width, height, depth, fg)

int	id;
Coord	coord1;
Real	width, height, depth;
XColor  *fg;
{
	X3DGraphics	*graphics;
	Coord		coord[5];
	XPoint		point[5];


	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X2D_draw_rectangle:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	coord[0] = coord[4] = coord1;

	coord[1] = coord[0];
	coord[1].x += width;
	coord[1].y += depth;

	coord[2] = coord[1];
	coord[2].z += height;

	coord[3] = coord[1];
	coord[3].z += height;

	if(!(_X2D_convert_wc_to_dc(graphics, coord, point, 5)))
	   return;

	X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	X3D_draw[graphics->device].lines(graphics, point, 5, CoordModeOrigin);
}



/************************************************************
*
*  MODULE NAME:  find_corner
*
*      PURPOSE:  This routine finds the rear most corner to start
*		 drawing a plot.  Horizon & Mesh plots use painter's
*		 alorithm in order to draw the plot, thus it is
*		 important to find the opposite corner in which to
*		 start drawing from.
*
*	INPUT :  id - id of the gwin
*
*       OUTPUT:  short --  the corner found.
*
*    CALLED BY:  X3D_draw_mesh
*		 X3D_draw_horizon
*		 X3D_draw_surface
*		 X3D_draw_colormesh
*
*   WRITTEN BY:  Mark Young
*
*
*	SUMMARY:  BL - Back Left  BR - Back Right  
*		  FL - Front Left  FR - Front Right
*
*			  t  h  e  t  a
*
* 		       | 0 - 90 |  90 - 180  180 - 270  270 - 360
*	   -------------------------------------------------------
*	a   0 - 90     |  BL    |    FL     |    FR   |    BR
*	   -------------------------------------------------------
*	l   90 - 180   |  BL    |    FL     |    FR   |    BR
*	   -------------------------------------------------------
*	p   180 - 270  |  FL    |    FR     |    BR   |    BL
*	   -------------------------------------------------------
*	h   270 - 360  |  FL    |    FR     |    BR   |    BL
*	   -------------------------------------------------------
*	a
*
******************************************************************/


#define BACKLEFT	1
#define BACKRIGHT	2
#define FRONTLEFT	3
#define FRONTRIGHT	4


static short find_corner (id)

int	id;
{
	Real   alpha, theta, eye_dist;

	X3D_inquire_viewpoint (id, &alpha, &theta, &eye_dist);

	if ((theta >= 0) && (theta <= 90)) 
	{
	     if ((alpha >= 0) && (alpha <= 180)) 
		 return (BACKLEFT);
	     else return (FRONTLEFT);
	}
	else if ((theta >= 90) && (theta <= 180)) 
	{
	     if ((alpha >= 0) && (alpha <= 180)) 
		 return (FRONTLEFT);
	     else return (FRONTRIGHT);
	}
	else if ((theta >= 180) && (theta <= 270)) 
	{
	     if ((alpha >= 0) && (alpha <= 180)) 
		 return (FRONTRIGHT);
	     else return (BACKRIGHT);
	}
	else if ((theta >= 270) && (theta <= 360)) 
	{
	     if ((alpha >= 0) && (alpha <= 180)) 
		 return (BACKRIGHT);
	     else return (BACKLEFT);
	}
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_mesh
*
*      PURPOSE:  
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		row_size  - number of points in each polyline
*		size      - number of points
*		fg	  - foreground color
*		bg	  - background color
*
*       OUTPUT:  mesh plot displayed in current graphics window
*
*    CALLED BY:  application program
*
*   WRITTEN BY:  Mark Young
*
*
*************************************************************/


X3D_draw_mesh(id, coords, row_size, size, fg, bg)

int	id; 
Coord 	*coords;
int     row_size, size;
XColor  *fg, *bg;
{
        X3DGraphics     *graphics;
	XPoint		*points, polygon[5];

	register int indx;
	short	 corner, find_corner();
	int	 row, col, num_row, num_col, i, j;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_mesh:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}
	
	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_mesh:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

	num_col = row_size;
	num_row = size/num_col;
	row = num_row -1;
	col = num_col -1;

	/*  Compute which corner to start drawing from */
	corner = find_corner(id);

	if (corner == BACKLEFT)
	{
	    for (i = 0; i < col; i++)
	    {
	        for (j = row; j > 0; j--)
		{
		    indx = j * num_col + i;
		    polygon[0] =
		    polygon[4] = points[indx];

		    polygon[1] = points[++indx];

		    indx -= num_col;
		    polygon[2] = points[indx];

		    polygon[3] = points[--indx];

		    /*
		     *  Fill the desired polygon outline
		     */
		    X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		    X3D_draw[graphics->device].fill_polygon(graphics, polygon,
				4, Nonconvex, CoordModeOrigin);

		    /*
		     *  Draw the desired polyline outline.
		     */
		   X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
		   X3D_draw[graphics->device].lines(graphics, polygon, 5,
			CoordModeOrigin);
		}
	    }
	}
	else if (corner == BACKRIGHT)
	{
	    for (j = row; j > 0; j--)
	    {
	        for (i = col; i > 0; i--)
		{
	            indx = j * num_col + i;
		    polygon[0] =
		    polygon[4] = points[indx];

		    polygon[1] = points[--indx];

		    indx -= num_col;
		    polygon[2] = points[indx];

		    polygon[3] = points[++indx];

		    /*
		     *  Fill the desired polygon outline
		     */
		    X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		    X3D_draw[graphics->device].fill_polygon(graphics, polygon,
				4, Nonconvex, CoordModeOrigin);

		    /*
		     *  Draw the desired polyline outline.
		     */
		   X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
		   X3D_draw[graphics->device].lines(graphics, polygon, 5,
			CoordModeOrigin);
		}
	    }
	}
	else if (corner == FRONTRIGHT)
	{
	    for (i = col; i > 0; i--)
	    {
	        for (j = 0; j < row; j++)
		{
	            indx = j * num_col + i;
		    polygon[0] =
		    polygon[4] = points[indx];

		    polygon[1] = points[--indx];

		    indx += num_col;
		    polygon[2] = points[indx];

		    polygon[3] = points[++indx];

		    /*
		     *  Fill the desired polygon outline
		     */
		    X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		    X3D_draw[graphics->device].fill_polygon(graphics, polygon,
				4, Nonconvex, CoordModeOrigin);

		    /*
		     *  Draw the desired polyline outline.
		     */
		   X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
		   X3D_draw[graphics->device].lines(graphics, polygon, 5,
			CoordModeOrigin);
		}
	    }
	}
	else if (corner == FRONTLEFT)
	{
	    for (j = 0; j < row; j++)
	    {
	        for (i = 0; i < col; i++)
		{
	            indx = j * num_col + i;
		    polygon[0] =
		    polygon[4] = points[indx];

		    polygon[1] = points[++indx];

		    indx += num_col;
		    polygon[2] = points[indx];

		    polygon[3] = points[--indx];

		    /*
		     *  Fill the desired polygon outline
		     */
		    X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		    X3D_draw[graphics->device].fill_polygon(graphics, polygon,
				4, Nonconvex, CoordModeOrigin);

		    /*
		     *  Draw the desired polyline outline.
		     */
		   X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
		   X3D_draw[graphics->device].lines(graphics, polygon, 5,
			CoordModeOrigin);
		}
	    }
	}
	free(points);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_horizon
*
*      PURPOSE:  Draws a horizon plot.  This is done by finding
*		 opposite side of the plot and drawing each row
*		 from back to front.  Painter's algorithm is used
*		 to fill the row, thus obscuring the previous row.
*		 This is repeated until all rows have been drawn.
*
*		 Note:  We say row, but we also mean column.  We draw
*			rows or columns depending on which side of the
*			plot is furthest from the viewer.  Thus if the
*			back left corner is furthest then we would draw
*			each row, but if the front left corner is the
*			furthest then we draw the each column.
*			
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		row_size  - number of points in each polyline
*		size      - number of points
*		fg	  - foreground color
*		bg	  - background color
*		wcmin	  - the minimum extent of the plot
*		wcmax	  - the maximum extent of the plot
*
*       OUTPUT:  horizon plot displayed in current graphics window
*
*    CALLED BY:  application program
*
*   WRITTEN BY:  Mark Young
*
*
*************************************************************/


X3D_draw_horizon(id, coords, rows, size, fg, bg, wcmin, wcmax)

int	id; 
Coord 	*coords, wcmin, wcmax;
int     rows, size;
XColor  *fg, *bg;
{
        X3DGraphics     *graphics;
	register XPoint	*points;
	register Coord	*horizon;

	register int indx, i, j, m;
	short	 corner, find_corner();
	int	 cols, rsize = rows +2;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_horizon:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}
	
	if (!(points = (XPoint *) malloc ((unsigned) rsize * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_horizon:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if (!(horizon = (Coord *) malloc ((unsigned) rsize * sizeof(Coord))))
	{
	   fprintf (stderr,"X3D_draw_horizon:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   free(points); return;
	}


	/*  Compute which corner is closest to you */
	corner = find_corner(id);
	cols = size/rows;

	if (corner == BACKLEFT)
	{
	   horizon[0] = horizon[rows +1] = wcmin;
	   horizon[rows +1].x = wcmax.x;

	   for (i = 0; i < cols; i++)
	   {
	      m = (cols -i -1) * rows;
	      horizon[0].y = coords[m].y;

	      for (j = 1; j <= rows; j++)
		 horizon[j] = coords[m++];
	      horizon[rows +1].y = horizon[rows].y;

	      if (!(_X3D_convert_wc_to_dc(graphics, horizon, points, rsize)))
	      {
	         free(horizon); free(points); return;
	      }

	      /*
	       *  Fill the desired polygon outline
	       */
	      X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
	      X3D_draw[graphics->device].fill_polygon(graphics, points, rsize,
				Nonconvex, CoordModeOrigin);

	      /*
	       *  Draw the desired polyline outline.
	       */
	      X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	      X3D_draw[graphics->device].lines(graphics, points, rsize,
			CoordModeOrigin);
	   }
	}
	else if (corner == FRONTLEFT)
	{
	   horizon[0] = horizon[rows +1] = wcmin;
	   horizon[rows +1].x = wcmax.x;

	   for (i = 0; i < cols; i++)
	   {
	      m = i * rows;
	      horizon[0].y = coords[m].y;

	      for (j = 1; j <= rows; j++)
		 horizon[j] = coords[m++];
	      horizon[rows +1].y = horizon[rows].y;

	      if (!(_X3D_convert_wc_to_dc(graphics, horizon, points, rsize)))
	      {
	         free(horizon); free(points); return;
	      }

	      /*
	       *  Fill the desired polygon outline
	       */
	      X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
	      X3D_draw[graphics->device].fill_polygon(graphics, points, rsize,
				Nonconvex, CoordModeOrigin);

	      /*
	       *  Draw the desired polyline outline.
	       */
	      X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	      X3D_draw[graphics->device].lines(graphics, points, rsize,
			CoordModeOrigin);
	   }
	}
	else if (corner == BACKRIGHT)
	{
	   horizon[0] = horizon[rows +1] = wcmin;
	   horizon[0].x = wcmax.x;

	   for (i = 0; i < cols; i++)
	   {
	      m = (cols -i) * rows -1;
	      horizon[0].y = coords[m].y;

	      for (j = 1; j <= rows; j++)
		 horizon[j] = coords[m--];
	      horizon[rows +1].y = horizon[rows].y;

	      if (!(_X3D_convert_wc_to_dc(graphics, horizon, points, rsize)))
	      {
	         free(horizon); free(points); return;
	      }

	      /*
	       *  Fill the desired polygon outline
	       */
	      X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
	      X3D_draw[graphics->device].fill_polygon(graphics, points, rsize,
				Nonconvex, CoordModeOrigin);

	      /*
	       *  Draw the desired polyline outline.
	       */
	      X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	      X3D_draw[graphics->device].lines(graphics, points, rsize,
			CoordModeOrigin);
	   }
	}
	else if (corner == FRONTRIGHT)
	{
	   horizon[0] = horizon[rows +1] = wcmin;
	   horizon[0].x = wcmax.x;

	   for (i = 0; i < cols; i++)
	   {
	      m = (i +1) * rows -1;
	      horizon[0].y = coords[m].y;

	      for (j = 1; j <= rows; j++)
		 horizon[j] = coords[m--];
	      horizon[rows +1].y = horizon[rows].y;

	      if (!(_X3D_convert_wc_to_dc(graphics, horizon, points, rsize)))
	      {
	         free(horizon); free(points); return;
	      }

	      /*
	       *  Fill the desired polygon outline
	       */
	      X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
	      X3D_draw[graphics->device].fill_polygon(graphics, points, rsize,
				Nonconvex, CoordModeOrigin);

	      /*
	       *  Draw the desired polyline outline.
	       */
	      X3D_draw[graphics->device].draw_color(graphics, fg, FALSE);
	      X3D_draw[graphics->device].lines(graphics, points, rsize,
			CoordModeOrigin);
	   }
	}
	free(points); free(horizon);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_surface
*
*      PURPOSE:  Draws a mesh plot using contour levels for
*		 the different color levels.  Surfaces are built
*		 exactly the same way that X3D_draw_mesh() does,
*		 accept the surfaces are left in 3D coordinate
*		 system.  The surface is then converted to device
*		 coordinates and the area is filled using the gc_fill
*		 color.  The surface is finally passed to the
*		 draw_colorline() routine, which draws the color
*		 border around the filled area.
*
*        INPUT:
*
*       OUTPUT:  surface plot displayed in current graphics window
*
*    CALLED BY:  application program
*
*   WRITTEN BY:  Mark Young
*
*
*************************************************************/

static int nump = 0;

X3D_draw_surface(id, coords, row_size, size, levels, colors, num_levels, wcmin,
		 wcmax, use_image, interp_image)

int	 id; 
Coord 	 *coords, wcmin, wcmax;
int      row_size, size, num_levels;
Real	 *levels;
XColor   *colors;
Boolean  use_image, interp_image;
{
        X3DGraphics     *graphics;
	Coord		polygon[4];
	XPoint		*points, pts[4];
	Real		*X3D_get_contour_levels();

	register int indx;
	int	 free_levels = FALSE;
	short	 corner, find_corner();
	int	 row, col, num_row, num_col, i, j;
/*
printf("entering X3D_draw_surface\n");
nump = 0;
 */
	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_surface:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_surface:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if (!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

	num_col = row_size;
	num_row = size/num_col;
	row = num_row -1;
	col = num_col -1;

	/*  Compute which corner to start drawing from */
	corner = find_corner(id);

	if (levels == NULL)
	{
	   free_levels = TRUE;
	   levels = X3D_get_contour_levels(wcmin.z, wcmax.z, 0.0, num_levels);
	   if (levels == NULL)
	   {
	      free(points);
	      return;
	   }
	}


	if (corner == BACKLEFT)
	{
	   for (i = 0; i < col; i++)
	   {
	      for (j = row; j > 0; j--)
	      {
		 indx = j * num_col + i;
		 polygon[0] = coords[indx];
		 pts[0]	    = points[indx];
		 polygon[1] = coords[++indx];
		 pts[1]	    = points[indx];

		 indx -= num_col;
		 polygon[2] = coords[indx];
		 pts[2]	    = points[indx];
		 polygon[3] = coords[--indx];
		 pts[3]	    = points[indx];

		 if (use_image)
		 {
	           fill_colorplane(graphics, polygon[0], polygon[1], polygon[3],
				 pts[0], pts[1], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
	           fill_colorplane(graphics, polygon[1], polygon[2], polygon[3],
				 pts[1], pts[2], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
		 }
		 else
		 {
	           fill_colorplane_elev(graphics, polygon[0], polygon[1],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
	           fill_colorplane_elev(graphics, polygon[1], polygon[2],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
		 }
	      }
	   }
	}
	else if (corner == BACKRIGHT)
	{
	   for (j = row; j > 0; j--)
	   {
	      for (i = col; i > 0; i--)
	      {
	         indx = j * num_col + i;
		 polygon[0] = coords[indx];
		 pts[0]	    = points[indx];
		 polygon[1] = coords[--indx];
		 pts[1]	    = points[indx];

		 indx -= num_col;
		 polygon[2] = coords[indx];
		 pts[2]	    = points[indx];
		 polygon[3] = coords[++indx];
		 pts[3]	    = points[indx];

		 if (use_image)
		 {
	           fill_colorplane(graphics, polygon[0], polygon[1], polygon[3],
				 pts[0], pts[1], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
	           fill_colorplane(graphics, polygon[1], polygon[2], polygon[3],
				 pts[1], pts[2], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
		 }
		 else
		 {
	           fill_colorplane_elev(graphics, polygon[0], polygon[1],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
	           fill_colorplane_elev(graphics, polygon[1], polygon[2],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
		 }
	      }
	   }
	}
	else if (corner == FRONTRIGHT)
	{
	   for (i = col; i > 0; i--)
	   {
	      for (j = 0; j < row; j++)
	      {
	         indx = j * num_col + i;
		 polygon[0] = coords[indx];
		 pts[0]	    = points[indx];
		 polygon[1] = coords[--indx];
		 pts[1]	    = points[indx];

		 indx += num_col;
		 polygon[2] = coords[indx];
		 pts[2]	    = points[indx];
		 polygon[3] = coords[++indx];
		 pts[3]	    = points[indx];

		 if (use_image)
		 {
	           fill_colorplane(graphics, polygon[0], polygon[1], polygon[3],
				 pts[0], pts[1], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
	           fill_colorplane(graphics, polygon[1], polygon[2], polygon[3],
				 pts[1], pts[2], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
		 }
		 else
		 {
	           fill_colorplane_elev(graphics, polygon[0], polygon[1],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
	           fill_colorplane_elev(graphics, polygon[1], polygon[2],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
		 }
	      }
	   }
	}
	else if (corner == FRONTLEFT)
	{
	   for (j = 0; j < row; j++)
	   {
	      for (i = 0; i < col; i++)
	      {
	         indx = j * num_col + i;
		 polygon[0] = coords[indx];
		 pts[0]	    = points[indx];
		 polygon[1] = coords[++indx];
		 pts[1]	    = points[indx];

		 indx += num_col;
		 polygon[2] = coords[indx];
		 pts[2]	    = points[indx];
		 polygon[3] = coords[--indx];
		 pts[3]	    = points[indx];


		 if (use_image)
		 {
	           fill_colorplane(graphics, polygon[0], polygon[1], polygon[3],
				 pts[0], pts[1], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
	           fill_colorplane(graphics, polygon[1], polygon[2], polygon[3],
				 pts[1], pts[2], pts[3], levels, colors,
				 num_levels, use_image, interp_image);
		 }
		 else
		 {
	           fill_colorplane_elev(graphics, polygon[0], polygon[1],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
	           fill_colorplane_elev(graphics, polygon[1], polygon[2],
				 polygon[3], levels, colors, num_levels,
				 interp_image);
		 }
	      }
	   }
	}

	/*
	 *  Cleanup time.
	 */
	free(points);
	if (free_levels) free(levels);
/*
printf("leaving X3D_draw_surface (num = %d)\n", nump);
 */
}



/************************************************************
*
*  MODULE NAME:  fill_colorplane()
*
*      PURPOSE:  Draws a multi-colored plane according to the
*		 different color levels.
*
*		 fill_colorplane is the real routine that computes
*		 the different colored line segments and draws them.
*		 This is done in a seperate routine in order to speed
*		 up the graphics processing by avoiding the costly
*		 setup time required for X3D_color_polygon().
*
*		 fill_colorplane fills an arbitary three dimension
*		 plane which is described by three points.  The
*		 three points describe a triangle, which may or may
*		 not intersect many different levels.  This routine
*		 is used in association with X3D_color_polygon() and
*		 X3D_draw_surface(),  routines which wish to draw one
*		 or more triangular color plane.
*
*        INPUT:
*
*       OUTPUT:  
*
*    CALLED BY:  internal graphics routines
*
*   WRITTEN BY:  Mark Young and Mike Lang
*
*
*************************************************************/



static void fill_colorplane(graphics, poly1, poly2, poly3, pt1, pt2, pt3,
		       levels, colors, num, use_image, interp_image)

X3DGraphics	*graphics;
Coord		poly1, poly2, poly3;
XPoint		pt1, pt2, pt3;
XColor          *colors;
Real		*levels;
int		num;
Boolean		use_image, interp_image;
{
	Coord	 poly[4];
	XPoint   pts[4];
	int	 i, incr, indx;
	register int index1, index2, indx3;


	/*
	 *  Make a quick check to see if the triangle is one color or if we
	 *  don't want to interpolate the plane, if so then we can just fill
	 *  it with one color and return.
	 */
	if ((poly1.d == poly2.d && poly2.d == poly3.d && use_image == TRUE) ||
	    (interp_image == FALSE))
	{
	   if (use_image)
	      index1 = poly1.d;
	   else index1 = find_level(poly1.z, levels, num);

	   pts[0] = pt1; pts[1] = pt2; pts[2] = pt3;
	   X3D_draw[graphics->device].fill_color(graphics, &colors[index1],
				TRUE);
	   X3D_draw[graphics->device].fill_polygon(graphics, pts, 3,
				Convex, CoordModeOrigin);
	   return;
	}

	/*
	 *  "interp_image" is set to true, so we need to be interpolating
	 *  the image.  Therefore, we need to find the indx for each point
	 *  depending if we are suppose to use the imagery or elevation "z"
	 *  values.
	 */
	if (use_image)
	{
	   index1 = poly1.d;
	   index2 = poly2.d;
	   indx3 = poly3.d;
	}
	else
	{
	   index1 = find_level(poly1.z, levels, num);
	   index2 = find_level(poly2.z, levels, num);
	   indx3 = find_level(poly3.z, levels, num);
	}

	/*
	 *  Now that we have the three values we are going to interpolate
	 *  the data from each of the corners.
	 */

	/*
	 *  Interpolate from index1 & index2 to indx3
	 */
	if ((index2 < index1 && indx3 < index1) ||
	    (index2 > index1 && indx3 > index1))
	{
	   poly[0] = poly[1] = poly1;
	   incr = ((index1 < indx3) ? 1 : -1);
	   for (i = index1; i != indx3 && i != index2; i+= incr)
	   {
	      poly[2].x = CINTRP(poly1.x, poly2.x, index1, index2, i+incr);
	      poly[2].y = CINTRP(poly1.y, poly2.y, index1, index2, i+incr);
	      poly[2].z = CINTRP(poly1.z, poly2.z, index1, index2, i+incr);
	      poly[3].x = CINTRP(poly1.x, poly3.x, index1, indx3, i+incr);
	      poly[3].y = CINTRP(poly1.y, poly3.y, index1, indx3, i+incr);
	      poly[3].z = CINTRP(poly1.z, poly3.z, index1, indx3, i+incr);

	      if (_X3D_convert_wc_to_dc(graphics, poly, pts, 4))
	      {
	         indx = i;
	         FindColor(indx, colors);
	         X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
	         X3D_draw[graphics->device].fill_polygon(graphics, pts, 4,
				Convex, CoordModeOrigin);
	      }
	      poly[0] = poly[3];
	      poly[1] = poly[2];
	   }
	}

	/*
	 *  Interpolate from index1 & indx3 to index2
	 */
	if ((index1 < indx3 && index2 < indx3) ||
	    (index1 > indx3 && index2 > indx3))
	{
	   poly[0] = poly[1] = poly3;
	   incr = ((indx3 < index2) ? 1 : -1);
	   for (i = indx3; i != index2 && i != index1; i += incr)
	   {
	      poly[2].x = CINTRP(poly1.x, poly3.x, index1, indx3, i+incr);
	      poly[2].y = CINTRP(poly1.y, poly3.y, index1, indx3, i+incr);
	      poly[2].z = CINTRP(poly1.z, poly3.z, index1, indx3, i+incr);
	      poly[3].x = CINTRP(poly2.x, poly3.x, index2, indx3, i+incr);
	      poly[3].y = CINTRP(poly2.y, poly3.y, index2, indx3, i+incr);
	      poly[3].z = CINTRP(poly2.z, poly3.z, index2, indx3, i+incr);

	      if (_X3D_convert_wc_to_dc(graphics, poly, pts, 4))
	      {
	         indx = i;
	         FindColor(indx, colors);
	         X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
	         X3D_draw[graphics->device].fill_polygon(graphics, pts, 4,
				Convex, CoordModeOrigin);
	      }
	      poly[0] = poly[3];
	      poly[1] = poly[2];
	   }
	}

	/*
	 *  Interpolate from index2 & indx3 to index1
	 */
	if ((index1 < index2 && indx3 < index2) ||
	    (index1 > index2 && indx3 > index2))
	{
	   poly[0] = poly[1] = poly2;
	   incr = ((index2 < index1) ? 1 : -1);
	   for (i = index2; i != index1 && i != indx3; i += incr)
	   {
	      poly[2].x = CINTRP(poly2.x, poly3.x, index2, indx3, i+incr);
	      poly[2].y = CINTRP(poly2.y, poly3.y, index2, indx3, i+incr);
	      poly[2].z = CINTRP(poly2.z, poly3.z, index2, indx3, i+incr);
	      poly[3].x = CINTRP(poly1.x, poly2.x, index1, index2, i+incr);
	      poly[3].y = CINTRP(poly1.y, poly2.y, index1, index2, i+incr);
	      poly[3].z = CINTRP(poly1.z, poly2.z, index1, index2, i+incr);

	      indx = i;
	      FindColor(indx, colors);
	      if (_X3D_convert_wc_to_dc(graphics, poly, pts, 4))
	      {
	         X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
	         X3D_draw[graphics->device].fill_polygon(graphics, pts, 4,
				Convex, CoordModeOrigin);
	      }
	      poly[0] = poly[3];
	      poly[1] = poly[2];
	   }
	}
}



/************************************************************
*
*  MODULE NAME:  fill_colorplane_elev()
*
*      PURPOSE:  Draws a multi-colored plane according to the
*		 different color levels.
*
*		 fill_colorplane is the real routine that computes
*		 the different colored line segments and draws them.
*		 This is done in a seperate routine in order to speed
*		 up the graphics processing by avoiding the costly
*		 setup time required for X3D_color_polygon().
*
*		 fill_colorplane fills an arbitary three dimension
*		 plane which is described by three points.  The
*		 three points describe a triangle, which may or may
*		 not intersect many different levels.  This routine
*		 is used in association with X3D_color_polygon() and
*		 X3D_draw_surface(),  routines which wish to draw one
*		 or more triangular color plane.
*
*        INPUT:
*
*       OUTPUT:  
*
*    CALLED BY:  internal graphics routines
*
*   WRITTEN BY:  Mark Young and Mike Lang
*
*
*************************************************************/


static void fill_colorplane_elev(graphics, pt1, pt2, pt3, levels, colors, num,
		interp_image)

X3DGraphics	*graphics;
Coord		pt1, pt2, pt3;
Real		*levels;
XColor		*colors;
int		num;
Boolean		interp_image;
{
	XPoint	points[4];
	Coord	min, mid, max, pts[4];
	int	indx, index1, index2, indx3; 

	if (pt1.z < pt2.z)
	{
	   if (pt3.z < pt1.z)
	   {
	      min = pt3;
	      mid = pt1;
	      max = pt2;
	   }
	   else if (pt3.z > pt2.z)
	   {
	      min = pt1;
	      mid = pt2;
	      max = pt3;
	   }
	   else
	   {
	      min = pt1;
	      mid = pt3;
	      max = pt2;
	   }
	}
	else
	{
	   if (pt3.z < pt2.z)
	   {
	      min = pt3;
	      mid = pt2;
	      max = pt1;
	   }
	   else if (pt3.z > pt1.z)
	   {
	      min = pt2;
	      mid = pt1;
	      max = pt3;
	   }
	   else
	   {
	      min = pt2;
	      mid = pt3;
	      max = pt1;
	   }
	}

	index1 = find_level(min.z, levels, num);
	index2 = find_level(mid.z, levels, num);
	indx3 = find_level(max.z, levels, num);

	if ((index1 == index2 && index2 == indx3) || interp_image == FALSE)
	{
	   pts[0] = min; pts[1] = mid; pts[2] = max;
	   if (_X3D_convert_wc_to_dc(graphics, pts, points, 3))
	   {
	      FindColor(index1, colors);
	      X3D_draw[graphics->device].fill_color(graphics, &colors[index1],
				TRUE);
	      X3D_draw[graphics->device].fill_polygon(graphics, points, 3,
				Convex, CoordModeOrigin);
	   }
	}
	else
	{
	   pts[0] = pts[1] = min;
	   while (index1 != index2)
	   {
	      /*
	       *  Compute the end point for 
	       */
	      pts[3].x = CINTRP(pts[0].x,mid.x,pts[0].z,mid.z,levels[index1+1]);
	      pts[3].y = CINTRP(pts[0].y,mid.y,pts[0].z,mid.z,levels[index1+1]);
	      pts[3].z = levels[index1 +1];

	      pts[2].x = CINTRP(pts[1].x,max.x,pts[1].z,max.z,levels[index1+1]);
	      pts[2].y = CINTRP(pts[1].y,max.y,pts[1].z,max.z,levels[index1+1]);
	      pts[2].z = levels[index1 +1];

	      if (_X3D_convert_wc_to_dc(graphics, pts, points, 4))
	      {
		 indx = index1;
	         FindColor(indx, colors);
		 X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
		 X3D_draw[graphics->device].fill_polygon(graphics, points, 4,
				Convex, CoordModeOrigin);
	      }
	      pts[0] = pts[3];
	      pts[1] = pts[2];
	      index1++;
	   }

	   if (pts[0].z < mid.z)
	   {
	      pts[2] = mid;
	      if (_X3D_convert_wc_to_dc(graphics, pts, points, 3))
	      {
		 indx = index2;
	         FindColor(indx, colors);
		 X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
	         X3D_draw[graphics->device].fill_polygon(graphics, points, 3,
			        Convex, CoordModeOrigin);
	      }
	   }
	   pts[0] = mid;

	   while (index2 != indx3)
	   {
	      /*
	       *  Compute the end point for 
	       */
	      pts[3].x = CINTRP(pts[0].x,max.x,pts[0].z,max.z,levels[index2+1]);
	      pts[3].y = CINTRP(pts[0].y,max.y,pts[0].z,max.z,levels[index2+1]);
	      pts[3].z = levels[index2 +1];

	      pts[2].x = CINTRP(pts[1].x,max.x,pts[1].z,max.z,levels[index2+1]);
	      pts[2].y = CINTRP(pts[1].y,max.y,pts[1].z,max.z,levels[index2+1]);
	      pts[2].z = levels[index2 +1];

	      if (_X3D_convert_wc_to_dc(graphics, pts, points, 4))
	      {
		 indx = index2;
	         FindColor(indx, colors);
		 X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
		 X3D_draw[graphics->device].fill_polygon(graphics, points, 4,
				Convex, CoordModeOrigin);
	      }
	      pts[0] = pts[3];
	      pts[1] = pts[2];
	      index2++;
	   }

	   /*
	    *  Draw the last part of the triangle in case the max
	    *  point is greater than indx3.
	    */
	   if (levels[indx3] < max.z && levels[indx3] < mid.z)
	   {
	      pts[2] = max;
	      pts[3] = mid;
	      if (_X3D_convert_wc_to_dc(graphics, pts, points, 4))
	      {
		 indx = indx3;
	         FindColor(indx, colors);
		 X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
	         X3D_draw[graphics->device].fill_polygon(graphics, points, 4,
			        Convex, CoordModeOrigin);
	      }
	   }
	   else if (levels[indx3] < max.z)
	   {
	      pts[2] = max;
	      if (_X3D_convert_wc_to_dc(graphics, pts, points, 3))
	      {
		 indx = indx3;
	         FindColor(indx, colors);
		 X3D_draw[graphics->device].fill_color(graphics, &colors[indx],
				TRUE);
	         X3D_draw[graphics->device].fill_polygon(graphics, points, 3,
			        Convex, CoordModeOrigin);
	      }
	   }
	}
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_colormesh
*
*      PURPOSE:  Draws a mesh plot using contour levels for
*		 the different color levels.  Surfaces are built
*		 exactly the same way that X3D_draw_mesh() does,
*		 accept the surfaces are left in 3D coordinate
*		 system.  The surface is then converted to device
*		 coordinates and the area is filled using the gc_fill
*		 color.  The surface is finally passed to the
*		 draw_colorline() routine, which draws the color
*		 border around the filled area.
*
*        INPUT:
*
*       OUTPUT:  surface plot displayed in current graphics window
*
*    CALLED BY:  application program
*
*   WRITTEN BY:  Mark Young
*
*
*************************************************************/


X3D_draw_colormesh(id, coords, row_size, size, levels, colors, num_levels, bg,
		   wcmin, wcmax)

int	 id; 
Coord 	 *coords, wcmin, wcmax;
int      row_size, size, num_levels;
Real	 *levels;
XColor   *colors, *bg;
{
        X3DGraphics     *graphics;
	Coord		polygon[5];
	XPoint		points[5];
	Real		*X3D_get_contour_levels();

	register int indx;
	int	 free_levels = FALSE;
	short	 corner, find_corner();
	int	 row, col, num_row, num_col, i, j;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_colormesh:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}
	X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);

	num_col = row_size;
	num_row = size/num_col;
	row = num_row -1;
	col = num_col -1;

	/*  Compute which corner to start drawing from */
	corner = find_corner(id);

	if (levels == NULL)
	{
	   free_levels = TRUE;
	   levels = X3D_get_contour_levels(wcmin.z, wcmax.z, 0.0, num_levels);
	   if (levels == NULL)
	      return;
	}


	if (corner == BACKLEFT)
	{
	   for (i = 0; i < col; i++)
	   {
	      for (j = row; j > 0; j--)
	      {
		 indx = j * num_col + i;
		 polygon[0] =
		 polygon[4] = coords[indx];
		 polygon[1] = coords[++indx];

		 indx -= num_col;
		 polygon[2] = coords[indx];
		 polygon[3] = coords[--indx];

	         if (!(_X3D_convert_wc_to_dc(graphics, polygon, points, 5)))
		    return;

		 X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		 X3D_draw[graphics->device].fill_polygon(graphics, points, 5,
				Nonconvex, CoordModeOrigin);
		 draw_colorline(graphics, polygon, points, 5, levels, colors,
				num_levels);
	      }
	   }
	}
	else if (corner == BACKRIGHT)
	{
	   for (j = row; j > 0; j--)
	   {
	      for (i = col; i > 0; i--)
	      {
	         indx = j * num_col + i;
		 polygon[0] =
		 polygon[4] = coords[indx];
		 polygon[1] = coords[--indx];

		 indx -= num_col;
		 polygon[2] = coords[indx];
		 polygon[3] = coords[++indx];

	         if (!(_X3D_convert_wc_to_dc(graphics, polygon, points, 5)))
		    return;

		 X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		 X3D_draw[graphics->device].fill_polygon(graphics, points, 5,
				Nonconvex, CoordModeOrigin);
		 draw_colorline(graphics, polygon, points, 5, levels, colors,
				num_levels);
	      }
	   }
	}
	else if (corner == FRONTRIGHT)
	{
	   for (i = col; i > 0; i--)
	   {
	      for (j = 0; j < row; j++)
	      {
	         indx = j * num_col + i;
		 polygon[0] =
		 polygon[4] = coords[indx];
		 polygon[1] = coords[--indx];

		 indx += num_col;
		 polygon[2] = coords[indx];
		 polygon[3] = coords[++indx];

	         if (!(_X3D_convert_wc_to_dc(graphics, polygon, points, 5)))
		    return;

		 X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		 X3D_draw[graphics->device].fill_polygon(graphics, points, 5,
				Nonconvex, CoordModeOrigin);
		 draw_colorline(graphics, polygon, points, 5, levels, colors,
				num_levels);
	      }
	   }
	}
	else if (corner == FRONTLEFT)
	{
	   for (j = 0; j < row; j++)
	   {
	      for (i = 0; i < col; i++)
	      {
	         indx = j * num_col + i;
		 polygon[0] =
		 polygon[4] = coords[indx];
		 polygon[1] = coords[++indx];

		 indx += num_col;
		 polygon[2] = coords[indx];
		 polygon[3] = coords[--indx];

	         if (!(_X3D_convert_wc_to_dc(graphics, polygon, points, 5)))
		    return;

		 X3D_draw[graphics->device].fill_color(graphics, bg, FALSE);
		 X3D_draw[graphics->device].fill_polygon(graphics, points, 5,
				Nonconvex, CoordModeOrigin);
		 draw_colorline(graphics, polygon, points, 5, levels, colors,
				num_levels);
	      }
	   }
	}

	/*
	 *  Cleanup time.
	 */
	if (free_levels) free(levels);
}



/************************************************************
*
*  MODULE NAME:  X3D_draw_colorline
*		 static draw_colorline()
*
*      PURPOSE:  Draws a multi-colored line according to the
*		 different color levels.  X3D_draw_colorline()
*		 is a frontend for the real routine static 
*		 draw_colorline().
*
*		 X3D_draw_colorline gets the graphics structure,
*		 perpares the contour levels by making sure they
*		 are sorted in ascending order, converts the world
*		 coordinates to device coordinates, then calls
*		 draw_colorline to actually draw the multi-colored
*		 line.
*
*		 draw_colorline is the real routine that computes
*		 the different colored line segments and draws them.
*
*        INPUT:
*
*       OUTPUT:  a multi-colored line
*
*    CALLED BY:  application program
*
*   WRITTEN BY:  Mark Young
*
*
*************************************************************/


X3D_draw_colorline(id, coords, size, levels, colors, num_levels)

int	 id;
Coord	 *coords;
int      size, num_levels;
Real     *levels;
XColor   *colors;
{
	X3DGraphics	*graphics;
	XPoint		*points;

	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_colorline:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	if (!(points = (XPoint *) malloc ((unsigned) size * sizeof (XPoint))))
	{
	   fprintf (stderr,"X3D_draw_colorline:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if(!(_X3D_convert_wc_to_dc(graphics, coords, points, size)))
	{
	   free(points);
	   return;
	}

	/*
	 *  Sort contours in ascending order.
	 *
	 sort_contour(levels, &contours, colors, &colors, num_levels);
	 */

	/*
	 *  Call draw_colorline() to compute the different colored line
	 *  segments and draw them.
	 */
	draw_colorline(graphics, coords, points, size, levels, colors,
		       num_levels);
	free(points);
}



/************************************************************
*
*  MODULE NAME:  draw_colorline()
*
*      PURPOSE:  Draws a multi-colored line according to the
*		 different color levels.
*
*		 draw_colorline is the real routine that computes
*		 the different colored line segments and draws them.
*		 This is done in a seperate routine in order to speed
*		 up the graphics processing by avoiding the costly
*		 setup time required for X3D_draw_colorline().
*		 draw_colorline runs thru the coordinates finding
*		 the contour level indxes for a given line segment.
*		 If the levels are the same then the line is drawn
*		 with that color, otherwise the line is divided into
*		 segements for each contour level for the given line.
*
*        INPUT:
*
*       OUTPUT:  
*
*    CALLED BY:  internal graphics routines
*
*   WRITTEN BY:  Mark Young
*
*
*************************************************************/


static int draw_colorline(graphics, coords, points, size, levels, colors, num)

X3DGraphics	*graphics;
Coord		*coords;
XPoint		*points;
int		size, num;
Real		*levels;
XColor		*colors;
{
	int	i, index1, index2; 
	XPoint	start, end;
	Coord   point;

	index1 = find_level(coords[0].z, levels, num);
	X3D_draw[graphics->device].draw_color(graphics, &colors[index1], TRUE);

	for (i = 1; i < size; i++)
	{
	    index2 = find_level(coords[i].z, levels, num);
	    if (index1 == index2)
	    {
	       X3D_draw[graphics->device].line(graphics, points[i -1].x,
			points[i -1].y, points[i].x, points[i].y);
	    }
	    else
	    {
	       start = points[i -1];

	       if (index2 < index1)
	       {
	          while (index1 != index2)
	          {
		     /*
		      *  Compute the end point for 
		      */
		     point.x = CINTRP(coords[i-1].x, coords[i].x, coords[i-1].z,
				      coords[i].z, levels[index1]);
		     point.y = CINTRP(coords[i-1].y, coords[i].y, coords[i-1].z,
				      coords[i].z, levels[index1]);
		     point.z = levels[index1];

		     if (_X3D_convert_wc_to_dc(graphics, &point, &end, 1))
		     {
			X3D_draw[graphics->device].draw_color(graphics,
				&colors[index1], TRUE);
	                X3D_draw[graphics->device].line(graphics, start.x,
			        start.y, end.x, end.y);
		     }

		     /*
		      *  Save the end of the segment as the start of the new
		      *  segment.
		      */
		     start = end;
		     index1--;
	          }
		  X3D_draw[graphics->device].draw_color(graphics,
				&colors[index2], TRUE);
	          X3D_draw[graphics->device].line(graphics, start.x,start.y,
			   points[i].x, points[i].y);
	       }
	       else
	       {
	          while (index1 != index2)
	          {
		     /*
		      *  Compute the end point for 
		      */
		     point.x = CINTRP(coords[i-1].x, coords[i].x, coords[i-1].z,
				      coords[i].z, levels[index1 +1]);
		     point.y = CINTRP(coords[i-1].y, coords[i].y, coords[i-1].z,
				      coords[i].z, levels[index1 +1]);
		     point.z = levels[index1 +1];

		     if (_X3D_convert_wc_to_dc(graphics, &point, &end, 1))
		     {
			X3D_draw[graphics->device].draw_color(graphics,
				&colors[index1], TRUE);
	                X3D_draw[graphics->device].line(graphics, start.x,
			         start.y, end.x, end.y);
		     }

		     /*
		      *  Save the end of the segment as the start of the new
		      *  segment.
		      */
		     start = end;
		     index1++;
	          }
		  X3D_draw[graphics->device].draw_color(graphics,
				&colors[index2], TRUE);
	          X3D_draw[graphics->device].line(graphics, start.x,start.y,
			   points[i].x, points[i].y);
	       }
	    }
	}
}



static int find_level(value, levels, num)

Real	value, *levels;
int	num;
{
	int	indx = num -1;

	while (value < levels[indx] && indx > 0)
	   indx--;

	return(indx);
}



/*
 *
 *  Routines for drawing contour plots.
 *
 *	X3D_draw_contour()
 *	X3D_get_contour_levels()
 *	static compute_contour()
 *	static display_contour()
 *	static sort_contour()
 */

/************************************************************
*
*  MODULE NAME:  X3D_draw_contour
*
*      PURPOSE:  Produces a contour map of an arbitrary function
*	 given in array points.  This routine uses a grid square
*	 scan algorithm where each grid square is  analyzed
*	 for a zero crossing on each edge of the square at the
*	 desired contour level.  For each crossing, an interpolation
*	 is performed to find the location of intersection with
*	 the grid square.  Two intersections implies a line.
*	 Four intersections implies a saddle point or ridge.
*	 The intersection coordinates are stored in coordinate array
*	 and are plotted according to the above criteria.
*
*        INPUT: id        - X3D graphics structure ID
*		coords    - array of world coordinate points
*		rows	  - number of points per row
*		size      - number of points
*		levels	  - array of Real number containing the desired
*			    contour levels.
*		colors	  - array of xcolors containing the color to be
*			    associated with each level
*	       num_levels - number of desired contour levels
*		base	  - a boolean flag indicating if the contours
*			    are to be projected to the bottom of the window
*			    or at the desired contour level.
*
*    CALLED BY: application program
*
*   WRITTEN BY: Mark Young & Scott Wilson
*
*************************************************************/

X3D_draw_contour(id, coords, rows, size, levels, colors, num_levels, base)

int	 id, rows, size, num_levels, base;
Coord 	 *coords;
Real	 *levels;
XColor	 *colors;
{
        X3DGraphics     *graphics;
	register Coord	*pt0, *pt1, *pt2, *pt3;
	register Real	*zmin, *zmax, *minz, *maxz;

	int	 i, j, cols, free_levels = FALSE;
	Real	 wc_zmin, wc_zmax, base_level, *X3D_get_contour_levels();


	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_draw_contour:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	/*
	 *   Compute the minimum and maximum values for each set of
	 *   contour points.  This will be used by X3D_compute_contour
	 *   for each contour level.  We use it as a quick check to see
	 *   if a contour exist for that given level. 
	 */
	zmin   = (Real *) malloc(sizeof(Real) * size);
	zmax   = (Real *) malloc(sizeof(Real) * size);
	minz   = zmin;
	maxz   = zmax;

	if ((!zmin) || (!zmax))
	{
	   fprintf (stderr,"X3D_draw_contour:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	cols = size/rows;
	pt0 = coords;
	pt1 = pt0 + 1;
	pt2 = coords + rows;
	pt3 = pt2 + 1;

	/*
	 *  Need to find the overall min & max value for Z.  This information
	 *  is used if levels is NULL.  We then create our own levels evenly
	 *  distributed between the wc_zmin and wc_zmax.
	 */
	wc_zmin = wc_zmax = MIN4(pt0->z, pt1->z, pt2->z, pt3->z);

	for (i = 0; i < cols -1; i++)
	{
           for (j = 0; j < rows -1; j++)
           {
	       *minz = MIN4(pt0->z, pt1->z, pt2->z, pt3->z);
	       *maxz = MAX4(pt0->z, pt1->z, pt2->z, pt3->z);

	       if (*minz < wc_zmin)  wc_zmin = *minz;
	       if (*maxz > wc_zmax)  wc_zmax = *maxz;

	       pt0++; pt1++; pt2++; pt3++;
	       minz++; maxz++;
	   }
	   pt0++; pt1++; pt2++; pt3++;
	   minz++; maxz++;
	}

	if (levels == NULL)
	{
	   free_levels = TRUE;
	   levels = X3D_get_contour_levels(wc_zmin, wc_zmax, 0.0, num_levels);
	   if (levels == NULL)
	   {
	      free(zmin); free(zmax);
	      return;
	   }
	}

	for (i = 0; i < num_levels; i++)
	{
	   if (base)
	      base_level = graphics->wc_min.z;
	   else
	      base_level = levels[i];

	   compute_contour(id, coords, rows, cols, levels[i], &colors[i],
			   zmin, zmax, base_level);
	}

	/*
	 *  Cleanup time.
	 */
	if (free_levels) free(levels);
	free(zmin); free(zmax);
}



/************************************************************
*
*  MODULE NAME: X3D_get_contour_levels
*
*      PURPOSE: Creates contour levels
*               
*        INPUT: min_level  -  The minimum contour level.
*		max_level  -  The maximum contour level.
*		offset	   -  The offset value to add the given contour
*			      levels.
*		num	   -  Number of contour levels to be generated.
*
*       OUTPUT: levels     -  A real array containing the defined 
*			      contour levels.
*
*    CALLED BY: application routine and xvgraphics routines
*
*   WRITTEN BY: Mark Young
*
*
*************************************************************/

Real	*X3D_get_contour_levels(min_level, max_level, offset, num)

Real	min_level, max_level, offset;
int	num;
{
	int	i;
	char	mesg[512];
	Real	factor, *levels;

	if (!(levels = (Real *) malloc(sizeof(Real) * num)))
	{
	   sprintf(mesg, "Could not alloc contour %d levels\n", num);
	   xvf_error_wait(mesg, "X3D_get_contour_levels", NULL);
	   return(NULL);
	}

	factor = (max_level - min_level)/(num -1);
	for (i = 0; i < num; i++)
	{
	    levels[i] = min_level + i * factor + offset;
	}
	return(levels);
}



/************************************************************
*
*  MODULE NAME: compute_contour
*
*      PURPOSE: 
*               
*        INPUT: none
*
*       OUTPUT: none
*
*    CALLED BY: X3D_draw_contour
*
*   WRITTEN BY: Mark Young and Scott Wilson
*
*
*************************************************************/


static compute_contour(id, coords, rows, cols, level, color, minz, maxz, base)

int	  id, rows, cols;
Coord	  *coords;
Real	  base, level;
register  Real *minz, *maxz;
XColor    *color;
{
	Coord	 box[4];
	int	 size, i, j, c, n, num;
	register Coord *pt0, *pt1, *pt2, *pt3, *contours;


	size = rows * cols;
	if (!(contours = (Coord *) malloc(sizeof(Coord) * size)))
	{
	   fprintf (stderr,"X3D_compute_contour:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	/*
	 *  Initialize the box Z value to the specified base.  We do this
	 *  now rather than for each and every contour point.
	 */
	for (i = 0; i < 4; i++)
	   box[i].z = base;

	pt0 = coords;
	pt1 = pt0 + 1;
	pt2 = coords + rows;
	pt3 = pt2 + 1;
	num = 0;

	for (i = 0; i < cols -1; i++)
	{
	   for (j = 0; j < rows -1; j++)
	   {
	      if (*minz == *maxz)
		 ;
	      else if (level >= *minz && level <=  *maxz)
	      {
                 c = 0;

                 /*
                  *  Check for intersections and find coords.
                  */
                 if ((pt0->z <= level && pt2->z >= level) ||
		     (pt0->z >= level && pt2->z <= level))
                 {
		    box[c].x   = CINTRP(pt0->x, pt2->x, pt0->z, pt2->z, level);
		    box[c++].y = CINTRP(pt0->y, pt2->y, pt0->z, pt2->z, level);
                 }

                 if ((pt1->z <= level && pt0->z >= level) ||
		     (pt1->z >= level && pt0->z <= level))
                 {
		    box[c].x   = CINTRP(pt0->x, pt1->x, pt0->z, pt1->z, level);
		    box[c++].y = CINTRP(pt0->y, pt1->y, pt0->z, pt1->z, level);
                 }

                 if ((pt3->z <= level && pt2->z >= level) || 
		     (pt3->z >= level && pt2->z <= level))
                 {
		    box[c].x   = CINTRP(pt2->x, pt3->x, pt2->z, pt3->z, level);
		    box[c++].y = CINTRP(pt2->y, pt3->y, pt2->z, pt3->z, level);
                 }

                 if ((pt3->z <= level && pt1->z >= level) ||
		     (pt3->z >= level && pt1->z <= level))
                 {
		    box[c].x   = CINTRP(pt1->x, pt3->x, pt1->z, pt3->z, level);
		    box[c++].y = CINTRP(pt1->y, pt3->y, pt1->z, pt3->z, level);
		 }

                 /*
                  *  Perform plotting chores for this box
                  */
                 if (c > 1)
                 {
                    for (n = 0; n < (c -1); n++) 
                    {
		       contours[num++] = box[n];
		       contours[num++] = box[n +1];
		    }  
 
		    if ((num + 6) >= size)
		    {
		       display_contour(id, contours, num, color);
		       num = 0;
		    }
		 }
	      }
	      pt0++; pt1++; pt2++; pt3++;
	      minz++; maxz++;
	   }
	   pt0++; pt1++; pt2++; pt3++;
	   minz++; maxz++;
	}
	display_contour(id, contours, num, color);
	free(contours);
}



/************************************************************
*
*  MODULE NAME: display_contour
*
*      PURPOSE: 
*               
*        INPUT: none
*
*       OUTPUT: none
*
*    CALLED BY: compute_contour
*
*   WRITTEN BY: Mark Young and Scott Wilson
*
*
*************************************************************/


static display_contour(id, contours, num, color)

int	 id, num;
Coord	 *contours;
XColor   *color;
{
        X3DGraphics     *graphics;
        XSegment	*segments;
	int		size;


	if (!(graphics = _X3D_get_graphics(id)))
	{
	   fprintf (stderr,"X3D_display_contour:\n");
	   fprintf (stderr,"\t unknown graphics id\n");
	   return;
	}

	size   = num/2;
	if (!(segments = (XSegment *) malloc ((unsigned) size * 
			 sizeof(XSegment))))
	{
	   fprintf (stderr,"X3D_display_contour:\n");
	   fprintf (stderr,"\t not enough memory\n");
	   return;
	}

	if (!(_X3D_convert_wc_to_dc_seg(graphics, contours, segments, num)))
	{
	   free(segments);
	   return;
	}
	X3D_draw[graphics->device].draw_color(graphics, color, FALSE);
	X3D_draw[graphics->device].segments(graphics, segments, size);
	free(segments);
}
