/* @(#)ag_gint.c	17.1.1.1 (ES0-DMD) 01/25/02 17:33:30 */
/*===========================================================================
  Copyright (C) 1995 European Southern Observatory (ESO)
 
  This program is free software; you can redistribute it and/or 
  modify it under the terms of the GNU General Public License as 
  published by the Free Software Foundation; either version 2 of 
  the License, or (at your option) any later version.
 
  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
 
  You should have received a copy of the GNU General Public 
  License along with this program; if not, write to the Free 
  Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, 
  MA 02139, USA.
 
  Corresponding concerning ESO-MIDAS should be addressed as follows:
	Internet e-mail: midas@eso.org
	Postal address: European Southern Observatory
			Data Management Division 
			Karl-Schwarzschild-Strasse 2
			D 85748 Garching bei Muenchen 
			GERMANY
===========================================================================*/

/* @(#)ag_gint.c	17.1.1.1  (OAA-ASTRONET) 01/25/02 17:33:30   */
/*
 * HEADER : ag_gint.c      - Vers 3.6.000  - Oct 1991 -  L. Fini, OAA
 *
 *
 *                           C INTERFACE MODULE
 */


#include <aglenvr.h>
#include <aglstat.h>


static void AG_npoints(x0,x1,y0,y1,nsteps,xstep,ystep)

/* this routine tries to guess the number of interpolation points needed to
 * Trace with suitable resolution a segment betveen two points.

* Algorithm:

*  1 - Find out an approximation of the length of the segment, summing 
*      up ten pieces of it (suitable if non-linearity is not too savage)

*  2 - Compute how many pixels will be traversed by each pieces both 
*      for X and Y directions

*  3 - take the sum of the two (X,Y) numbers as an extimate of the number 
*      of steps.
*/


#define NPIECES  10

double x0,x1;		/* X and X+1 values                                 */
double y0,y1;		/* Y and Y+1 values                                 */

			/* Output  values:                                  */
double *nsteps;		/* Number of step to interpolate between (x0,y0)    */
			/* and (x1,y1)   with suitable resolution           */
double *xstep;		/* X-axis interpolation step                        */
double *ystep;		/* Y-axis interpolation step                        */
{
double px[NPIECES+1],py[NPIECES+1];
double lx,ly;
double xst,yst;
double xd,yd;
double steps;
int k;

px[0]      = x0;
px[NPIECES]= x1;

py[0]      = y0;
py[NPIECES]= y1;

xd=(px[NPIECES]-px[0]);
yd=(py[NPIECES]-py[0]);
xst=xd/NPIECES;
yst=yd/NPIECES;

for(k=1;k<NPIECES;k++) {
	px[k] = px[k-1]+xst;
	py[k] = py[k-1]+yst;
}
for(k=0;k<=NPIECES;k++)
	AGL_cnvt(&px[k],&py[k]);

lx=0; ly=0;

for(k=0;k<NPIECES;k++) {
	lx += fabs(px[k+1]-px[k])/RSLX(AGL_status.devid);
	ly += fabs(py[k+1]-py[k])/RSLY(AGL_status.devid);
}

steps = floor(lx+ly);

*nsteps=steps;

steps += 1.0;

*xstep = xd / steps;
*ystep = yd / steps;
}

/*****************************************************************************/
/*++                               AG_GINT (C callable)                      */

/* AG_GINT    Draw a segments interpolating between points                   */


/* This module will trace a polyline interpolating in a suitable number of   */
/* intermediate points between the points defined in the call.               */
/* This module is of use when a non-linear transformation is active and a    */
/* line which is "straight" in the world coordinate system would not result  */
/* in a straight line on the graphic device.                                 */

/* This routine is less efficient than AG_GPLL() so when the input polyline  */
/* is known to contain points which are close enough to get proper resolution*/
/* AG_GPLL() should be called instead.                                       */

/* Line aspect (style, width, color) may be selected via AG_SSET.            */

/* The coordinate values are referred to the current graphic mode (either    */
/* NORMALIZED, USER or SPECIAL) and mapped accordingly.  Graphic mode is     */
/* affected by: AG_WDEF, AG_TRNS, and AG_SSET items "NORMAL","USER" and      */
/* "SPECIAL"                                                                 */


void AG_GINT (xv,yv,np)

 float xv[],yv[];             /* Coordinate arrays                           */
 int np;                      /* Array length                                */
                                                                         /*--*/
{
static char *modnam = "GINT";
int ks;
int k;

float xp[AUXPOLYBUFLNG],yp[AUXPOLYBUFLNG];
struct polyline my_poly;

if(AGL_DEBUG) 
	AG_DMSG(modnam,(char *)0);


#ifndef NO_METAFILE_SUPPORT
if ( DY(filact) == SOFTMETAFILE ) {
	struct polyline inpoly;
	enum METACODE cod = MFGINT;
	AGL_plmk(xv,yv,np,np,DY(modeflag),&inpoly);
	(void)fwrite((void *)&cod,sizeof(int),1,DY(metafile));
	AGL_mput(DY(metafile),&inpoly,FALSE);
	if(AGL_status.errstat != AGLNOERR) AGL_siger(modnam);
}
#endif

AGL_plmk(xp,yp,AUXPOLYBUFLNG,0,DY(modeflag),&my_poly);

my_poly.join=TRUE;

if(AGL_status.dashon) {
	AGL_spfls(AGL_dpll,&my_poly);
	AGL_dshz();
} else
	AGL_spfls(AGL_gpll,&my_poly);

np--;

for(k=0;k<np;k++) {
	double nsteps,xstep,ystep;
	double x0= *xv++,y0= *yv++;

	AG_npoints(x0,(double)*xv,y0,(double)*yv,&nsteps,&xstep,&ystep);
	for(ks=0;ks<=nsteps;ks++) {
		double xx,yy;
		xx=x0+ks*xstep;
		yy=y0+ks*ystep;
		AGL_pcoo(xx,yy,&my_poly);
	}
}

AGL_pflsh(&my_poly);

if(AGL_status.errstat != AGLNOERR) 
	AGL_siger(modnam);
}