/* C Sociedad Protectora de Diplodocus soft. nov. 1991, Jun 93*/
/* sector de isoaxis unitario */

#include <gl/gl.h>

#include <math.h>
#include <stdio.h>
#include "points.h"
#include "fmath.h"

#define _X [0]
#define _Y [1]
#define _Z [2]
#define SQR3 1.732050807568877293527446341505872366942805
#define FEQ(a,b)	(fabs((double)(a)-(double)(b))<TOL)

/*#define	LOG*/

double		_beta,_gamma;

short int	fError;
extern int	cheat_flag;

short int	sector_isoaxis(double delta,double alfai, Point3D *aux);

short int sector_isoaxis(double delta,double alfa, Point3D *aux)
{
	static int npolig=10,nvert[10];
	int i;
	double ang1,ang2;
	Point3D tmp;
	short int TheResult=0;
	Float ag1, ag2, ag3, ag4;
	
#ifdef LOG
	FILE *arclog;
	static char	fa=0;
	
	if( fa==0 )
	{
		fa=1;
		arclog=fopen("isoaxis.out","w");
		fprintf(arclog,"delta\talfa\tag1\tag2\ty0\ty1\ty1\ty3\ty4\ty5\ty6\n");
	}
	else
	{
		arclog=fopen("isoaxis.out","a");
	}
	fprintf(arclog,"%lf\t%lf\t", DEGS(delta), alfa );
	
#endif /*LOG*/
		
	
	/* Asignacion  valores a las coordenadas dependientes
		de alfa y de delta	*/

	aux[0]_X = (Float)sin(alfa);
	aux[0]_Y = (Float)SQR3 * aux[0]_X;
	aux[0]_Z = (Float)(- cos(alfa)*sin(delta));
	
	aux[4]_X = 0.0;
	aux[4]_Y = (Float)(aux[0]_Y + cos(alfa)*cos(delta) - sin(delta));
	aux[4]_Z = (Float)cos(delta);
	
	aux[5]_X = 0.0;
	aux[5]_Y = aux[0]_Y + (Float)(cos(alfa)*cos(delta) + sin(delta));
	aux[5]_Z = -aux[4]_Z;
	
	aux[3]_X = aux[6]_X = 0.0;
	/*___________________________________________________________________ */
	tmp[0]=tmp[2]=0.0;
	tmp[1]=(aux[4]_Y+aux[5]_Y)/2.0;

	if((fError=fesf(aux[1],aux[5],aux[0],tmp))!=NO_ERR)
	{
		fprintf(stderr,"aborto fesf(5) =%04X, alfa incorrecto\n",fError);
		TheResult=-1;
		goto TheExit;
	}
	if((fError=ffsen(	(double)(aux[1]_Y-aux[5]_Y),
				(double)(aux[5]_Z-aux[1]_Z),
				&ang1,&ang2))!=NO_ERR)
	{
		fprintf(stderr,"aborto _gamma ffsen=%04X, error eleccion x fcuad\n",fError);
		TheResult=-1;
		goto TheExit;
	}
  	{
	Point3D vjm, vji;
	
		aux[6]_Y = aux[5]_Y + (Float)sin(ang1);
		aux[6]_Z = aux[5]_Z + (Float)cos(ang1);
		aux[6]_X = 0.0;
		
		ResV( aux[6], aux[5], vjm );
		ResV( aux[4], aux[5], vji );
		ag1=AngrV(vjm,vji);
		
		aux[6]_Y = aux[5]_Y + (Float)sin(ang2);
		aux[6]_Z = aux[5]_Z + (Float)cos(ang2);
		aux[6]_X = 0.0;
		
		ResV( aux[6], aux[5], vjm );
		ResV( aux[4], aux[5], vji );
		ag2=AngrV(vjm,vji);
		
		if((ag1<-TOL)&&(ag2<-TOL))
		{
			fprintf(stderr,"aborto _gamma error raro\n");
			TheResult=-2;
			goto TheExit;
		}
		if (ag1< 0 ) _gamma = ang2;
		else if ( ag2 < 0 ) _gamma = ang1;
		else
			_gamma = ( ag1 < ag2 ? ang1:ang2 );
	}
	/*___________________________________________________________________ */
	
	if((fError=fesf(aux[2],aux[4],aux[0],tmp))!=NO_ERR)
	{
		fprintf(stderr,"aborto fesf(4) =%04X, alfa incorrecto\n",fError);
		TheResult=-1;
		goto TheExit;
	}
	if((fError=ffsen(	(double)(aux[2]_Y-aux[4]_Y),
				(double)(aux[2]_Z-aux[4]_Z),
				&ang1,&ang2))!=NO_ERR)
	{
		fprintf(stderr,"aborto _beta ffsen=%04X, error eleccion x fcuad\n",fError);
		TheResult=-1;
		goto TheExit;
	}
  	{
		Point3D vjm, vji;
		
			aux[3]_Y = aux[4]_Y + (Float)sin(ang1);
			aux[3]_Z = aux[4]_Z - (Float)cos(ang1);
			aux[3]_X = 0.0;
			
			ResV( aux[3], aux[4], vjm );
			ResV( aux[5], aux[4], vji );
			ag3=AngrV(vjm,vji);
			
			aux[3]_Y = aux[4]_Y + (Float)sin(ang2);
			aux[3]_Z = aux[4]_Z - (Float)cos(ang2);
			aux[3]_X = 0.0;
			
			ResV( aux[3], aux[4], vjm );
			ResV( aux[5], aux[4], vji );
			ag4=AngrV(vjm,vji);
		
		if((ag3<-TOL)&&(ag4<-TOL))
		{
			fprintf(stderr,"aborto _beta error raro\n");
			TheResult=-2;
			goto TheExit;
		}
		if (ag3<0) _beta = ang2;
		else if ( ag4 < 0 ) _beta = ang1;
		else
			_beta = ( ag3 < ag4 ? ang1:ang2 );
	}
	
	/*___________________________________________________________________ */
	/*  verificacion y trampa	*/
for(i=0;i<7;i++)
{
	if( cheat_flag )
	{
		if ( aux[i]_Y < 0.0  ) aux[i]_Y=0.0;
		if ( aux[i]_X < 0.0  ) aux[i]_X=0.0;
	}
#ifdef DEBUG
	else
		if ( aux[i]_Y < -TOL  )
			fprintf(stderr, "Error Y%d=%f<0 \n",i,aux[i]_Y);
#endif
}
	/*___________________________________________________________________ */
	/* asignacion los valores de las coodenadas que dependen de _beta y de
	_gamma, asi como los simetricos de 0,1,2
	*/

	aux[3]_Y = aux[4]_Y + (Float)sin(_beta);
	aux[3]_Z = aux[4]_Z - (Float)cos(_beta);
	
	aux[6]_Y = aux[5]_Y + (Float)sin(_gamma);
	aux[6]_Z = aux[5]_Z + (Float)cos(_gamma);
	
	aux[9]_X = -aux[0]_X;	aux[9]_Y = aux[0]_Y;	aux[9]_Z = aux[0]_Z;
	aux[8]_X = -aux[2]_X;	aux[8]_Y = aux[2]_Y;	aux[8]_Z = aux[2]_Z;
	aux[7]_X = -aux[1]_X;	aux[7]_Y = aux[1]_Y;	aux[7]_Z = aux[1]_Z;
	/*___________________________________________________________________ */
	
	
	for(i=0;i<npolig;i++) nvert[i]=3;
	
		
#ifdef DEBUG

if ( !FEQ(DistV( aux[0], aux[1] ),2.0 ) )
	fprintf(stderr,"d 0-1: %lf ( %lf )\n", DistV( aux[0], aux[1] ), 2.0);
if ( !FEQ(DistV( aux[0], aux[2] ),2.0 ) )
	fprintf(stderr,"d 0-2: %lf ( %lf )\n", DistV( aux[0], aux[2] ), 2.0);
if ( !FEQ(DistV( aux[4], aux[5] ),2.0 ) )
	fprintf(stderr,"d 4-5: %lf ( %lf )\n", DistV( aux[4], aux[5] ), 2.0);

if ( !FEQ(DistV( aux[0], aux[4] ),sqrt(2.0) ) )
	fprintf(stderr,"d 0-4: %lf ( %lf )\n", DistV( aux[0], aux[4] ), sqrt(2.0));
if ( !FEQ(DistV( aux[0], aux[5] ),sqrt(2.0) ) )
	fprintf(stderr,"d 0-5: %lf ( %lf )\n", DistV( aux[0], aux[5] ), sqrt(2.0));
if ( !FEQ(DistV( aux[1], aux[5] ),sqrt(2.0) ) )
	fprintf(stderr,"d 1-5: %lf ( %lf )\n", DistV( aux[1], aux[5] ), sqrt(2.0));
if ( !FEQ(DistV( aux[2], aux[4] ),sqrt(2.0) ) )
	fprintf(stderr,"d 2-4: %lf ( %lf )\n", DistV( aux[2], aux[4] ), sqrt(2.0));
	
if ( !FEQ(DistV( aux[1], aux[6] ),1.0 ) )
	fprintf(stderr,"d 1-6: %lf ( %lf )\n", DistV( aux[1], aux[6] ), 1.0);
if ( !FEQ(DistV( aux[2], aux[3] ),1.0 ) )
	fprintf(stderr,"d 2-3: %lf ( %lf )\n", DistV( aux[2], aux[3] ), 1.0);
if ( !FEQ(DistV( aux[3], aux[4] ),1.0 ) )
	fprintf(stderr,"d 3-4: %lf ( %lf )\n", DistV( aux[3], aux[4] ), 1.0);

#endif
	
TheExit:
	if(TheResult!=0) aux[0][1]=-1.0;

#ifdef LOG
/*	fprintf(arclog,"%lf\t%lf\t%lf\t%lf", DEGS(delta), DEGS(alfa), DEGS(_beta), DEGS(_gamma));*/
	fprintf(arclog,"%f\t%f", ag1, ag2 );
	for(i=0;i<=6;i++)
		fprintf(arclog,"\t%lf", aux[i]_Y );
	fprintf(arclog,"\n");
	fclose(arclog);
#endif /*LOG*/
	return(TheResult);
} /* De la funcion */