/*
 * Fibracion de Hopf, proyeccion de circulos de Hopf del hiperespacio CxC a CxR
 * esto es para mostrarselo a Ximena :) 
 * compila 
 * gcc -lm -lGL -lglut -lGLU hopf.c -o hopf 
 *
 * beck@math.co.ro
 */

#include <GL/glut.h>
#include <complex.h>
#include <math.h>
#include <stdio.h>
float rot = 0.0f;

/* Vector en R4 */
typedef struct r4v
{
  double x1, x2, x3, x4;
} r4v_t;

/* Vector en R3 */
typedef struct r3v
{
  double x1, x2, x3;
} r3v_t;

/* Vector en CxR */
typedef struct crv
{
  complex z;
  double x;
} crv_t;

/* Vector en CxC */
typedef struct c2v
{
  complex z0, z1;
} c2v_t;

/* Identificamos RxRxRxR = R4 con el hiperespacio complejo CxC=C2 */

c2v_t
r4v_a_c2v (r4v_t u)
{
  c2v_t v;
  v.z0 = u.x1 + I * u.x2;
  v.z1 = u.x3 + I * u.x4;
  return v;
}

/* Identificamos RxRxR=R3 con el espacio de dimension 3 CxR */
crv_t
r3v_a_crv (r3v_t u)
{
  crv_t v;
  v.z = u.x1 + I * u.x2;
  v.x = u.x3;
  return v;
}

/* Fibracion de Hopf rho:CxC |-> CxR */
crv_t
rho (c2v_t u)
{
  crv_t v;
  v.z = 2 * u.z0 * conj (u.z1);
  v.x = (cabs (u.z0) * cabs (u.z0)) - (cabs (u.z1) * cabs (u.z1));
  return v;
}

/* parametrizacion de la hiperesfera S3 psi y theta en [0,pi], phi en [0,2pi] */
r4v_t
punto_s3 (double psi, double theta, double phi)
{
  r4v_t z;
  z.x1 = cos (psi);
  z.x2 = sin (psi) * cos (theta);
  z.x3 = sin (psi) * sin (theta) * cos (phi);
  z.x4 = sin (psi) * sin (theta) * sin (phi);
  return z;
}

/* Imprime vector de R4 */
void
print_r4v (r4v_t v)
{
  printf ("(%.2f,%.2f,%.2f,%.2f)\n", v.x1, v.x2, v.x3, v.x4);
  return;
}

/* Imprime vector de C2 */
void
print_c2v (c2v_t v)
{
  printf ("(%.2f+%.2fi,%.2f+%.2fi)\n", creal (v.z0), cimag (v.z0),
	  creal (v.z1), cimag (v.z1));
  return;
}

/* Imprime vector de CxR */
void
print_crv (crv_t v)
{
  printf ("(%.2f+%.2fi,%.2f)\n", creal (v.z), cimag (v.z), v.x);
  return;
}



void
display (void)
{
  double psi,theta,phi;
  r4v_t p; 
  glLoadIdentity ();
  glTranslatef (1.5f, 0.0f, -7.0f);
  glRotatef (rot, 1.0f, 1.0f, 1.0f);
  glClear (GL_COLOR_BUFFER_BIT);
  glBegin (GL_LINE_STRIP);

  for(psi=0;psi<=M_PI;psi+=0.1)
   for(theta=0;theta<=M_PI;theta+=0.1)
    for(phi=0;phi<=2*M_PI;phi+=0.1) {
 	p=punto_s3(psi,theta,phi);
		print_r4v(p);
	} 
			


  glEnd ();
  rot += 0.05;			// velocidad
  glutSwapBuffers ();
}

void
reshape (int w, int h)
{
  glViewport (0, 0, (GLsizei) w, (GLsizei) h);
  glMatrixMode (GL_PROJECTION);
  glLoadIdentity ();
  gluPerspective (50, (GLfloat) w / (GLfloat) h, 1.0, 20.0);
  glMatrixMode (GL_MODELVIEW);
  glLoadIdentity ();
  gluLookAt (0, 5, 5, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0);
}

int
main (int argc, char **argv)
{
  glutInit (&argc, argv);
  glutInitDisplayMode (GLUT_DOUBLE | GLUT_RGB);
  glutInitWindowSize (600, 400);
  glutInitWindowPosition (10, 10);
  glutCreateWindow ("test");
  glClearColor (0, 0, 0, 0);
  glShadeModel (GL_FLAT);
  glMatrixMode (GL_PROJECTION);
  glutReshapeFunc (reshape);
  glutDisplayFunc (display);
  glutIdleFunc (display);
  glutMainLoop ();
  return 0;
}