/* demo_mat.c  study perspective and model view matrices */

#include <stdio.h>
#include <math.h>
#include <GL/glut.h>

static GLfloat roll = 30.0;      /* rotate about Z axis 'r' */
static GLfloat pitch = 30.0;    /* rotate about X axis 'p' */
static GLfloat heading = 30.0;  /* rotate about Y axis 'h' */
static GLfloat pos[3] = {1.0, 2.0, -20.0}; /* x, y, z translate */
static GLdouble mmatrix[16];    /* 4 by 4 model transformation matrix */
static GLdouble pmatrix[16];    /* 4 by 4 perspective transformation matrix */
static GLdouble imatrix[16]=    /* 4 by 4 identity matrix */
                                {1.0, 0.0, 0.0, 0.0,
                                 0.0, 1.0, 0.0, 0.0,
                                 0.0, 0.0, 1.0, 0.0,
                                 0.0, 0.0, 0.0, 1.0};
static int win_w, win_h;        /* window width and height in pixels */
static GLdouble size = 2.0;     /* -2.0 to +2.0 in X at near */
static GLdouble wx0, wy0, wz0;  /* window x,y,z at near */
static GLdouble wx1, wy1, wz1;  /* window x,y,z at far */
static GLint viewport[4];       /* screen viewport */
static GLint y = 100;           /* a point in pixel coordinates */
static GLint x = 150;
static GLint realy;

void display(void)
{
  glClearColor(1.0, 1.0, 1.0, 1.0);
  glClear(GL_COLOR_BUFFER_BIT); 
  glColor3f(0.0, 0.0, 0.0);
  
  glGetDoublev(GL_MODELVIEW_MATRIX, mmatrix); /* save transformation */
  printf("\n model view matrix in display \n");
  printf("mmatrix[ 0.. 3] %f, %f, %f, %f \n",
          mmatrix[0], mmatrix[1], mmatrix[2], mmatrix[3]);
  printf("mmatrix[ 4.. 7] %f, %f, %f, %f \n",
          mmatrix[4], mmatrix[5], mmatrix[6], mmatrix[7]);
  printf("mmatrix[ 8..11] %f, %f, %f, %f \n",
          mmatrix[8], mmatrix[9], mmatrix[10], mmatrix[11]);
  printf("mmatrix[12..15] %f, %f, %f, %f \n",
          mmatrix[12], mmatrix[13], mmatrix[14], mmatrix[15]);

  glGetIntegerv (GL_VIEWPORT, viewport);
  glGetDoublev(GL_PROJECTION_MATRIX, pmatrix);
  printf("\n projection matrix in display \n");
  printf("pmatrix[ 0.. 3] %f, %f, %f, %f \n",
          pmatrix[0], pmatrix[1], pmatrix[2], pmatrix[3]);
  printf("pmatrix[ 4.. 7] %f, %f, %f, %f \n",
          pmatrix[4], pmatrix[5], pmatrix[6], pmatrix[7]);
  printf("pmatrix[ 8..11] %f, %f, %f, %f \n",
          pmatrix[8], pmatrix[9], pmatrix[10], pmatrix[11]);
  printf("pmatrix[12..15] %f, %f, %f, %f \n",
          pmatrix[12], pmatrix[13], pmatrix[14], pmatrix[15]);

  /*  note viewport[3] is height of window in pixels  */
  realy = viewport[3] - (GLint)y - 1;
  gluUnProject((GLdouble) x, (GLdouble) realy, 0.0, /* near */
               mmatrix, pmatrix, viewport, &wx0, &wy0, &wz0); 
  gluUnProject((GLdouble) x, (GLdouble) realy, 1.0, /* far */
               mmatrix, pmatrix, viewport, &wx1, &wy1, &wz1); 

  glPointSize(5.0);
  glBegin(GL_POINTS);
    glVertex3f(0.5,  1.5, 2.5); /* point in world coordinates */
  glEnd();
  glFlush();
} /* end display */


void myReshape(int w, int h)
{
  win_w = w;
  win_h = h;
  glViewport(0, 0, win_w, win_h);
  glMatrixMode(GL_PROJECTION);
  glLoadIdentity();
  glFrustum(-size,                                /* xmin */
             size,                                /* xmax */
            -size*(GLfloat)win_h/(GLfloat)win_w,  /* ymin */
             size*(GLfloat)win_h/(GLfloat)win_w,  /* ymax */
             1.0*size,                            /* near */
             100.0*size                           /* far  */
            );
  glGetDoublev(GL_PROJECTION_MATRIX, pmatrix);
  printf("\n projection matrix \n");
  printf("pmatrix[ 0.. 3] %f, %f, %f, %f \n",
          pmatrix[0], pmatrix[1], pmatrix[2], pmatrix[3]);
  printf("pmatrix[ 4.. 7] %f, %f, %f, %f \n",
          pmatrix[4], pmatrix[5], pmatrix[6], pmatrix[7]);
  printf("pmatrix[ 8..11] %f, %f, %f, %f \n",
          pmatrix[8], pmatrix[9], pmatrix[10], pmatrix[11]);
  printf("pmatrix[12..15] %f, %f, %f, %f \n",
          pmatrix[12], pmatrix[13], pmatrix[14], pmatrix[15]);


  glMatrixMode(GL_MODELVIEW);
  glLoadIdentity();
  glGetDoublev(GL_MODELVIEW_MATRIX, mmatrix);
  printf("\n model view matrix, identity \n");
  printf("mmatrix[ 0.. 3] %f, %f, %f, %f \n",
          mmatrix[0], mmatrix[1], mmatrix[2], mmatrix[3]);
  printf("mmatrix[ 4.. 7] %f, %f, %f, %f \n",
          mmatrix[4], mmatrix[5], mmatrix[6], mmatrix[7]);
  printf("mmatrix[ 8..11] %f, %f, %f, %f \n",
          mmatrix[8], mmatrix[9], mmatrix[10], mmatrix[11]);
  printf("mmatrix[12..15] %f, %f, %f, %f \n",
          mmatrix[12], mmatrix[13], mmatrix[14], mmatrix[15]);

  glRotatef(pitch, 1.0, 0.0, 0.0);      /* rotate world about X axis    */
  glGetDoublev(GL_MODELVIEW_MATRIX, mmatrix);
  printf("\n model view matrix, after rotation about X axis \n");
  printf("mmatrix[ 0.. 3] %f, %f, %f, %f \n",
          mmatrix[0], mmatrix[1], mmatrix[2], mmatrix[3]);
  printf("mmatrix[ 4.. 7] %f, %f, %f, %f \n",
          mmatrix[4], mmatrix[5], mmatrix[6], mmatrix[7]);
  printf("mmatrix[ 8..11] %f, %f, %f, %f \n",
          mmatrix[8], mmatrix[9], mmatrix[10], mmatrix[11]);
  printf("mmatrix[12..15] %f, %f, %f, %f \n",
          mmatrix[12], mmatrix[13], mmatrix[14], mmatrix[15]);

  glRotatef(heading, 0.0, 1.0, 0.0);    /* rotate world about Y axis    */
  glGetDoublev(GL_MODELVIEW_MATRIX, mmatrix);
  printf("\n model view matrix, after rotation about X and Y axis \n");
  printf("mmatrix[ 0.. 3] %f, %f, %f, %f \n",
          mmatrix[0], mmatrix[1], mmatrix[2], mmatrix[3]);
  printf("mmatrix[ 4.. 7] %f, %f, %f, %f \n",
          mmatrix[4], mmatrix[5], mmatrix[6], mmatrix[7]);
  printf("mmatrix[ 8..11] %f, %f, %f, %f \n",
          mmatrix[8], mmatrix[9], mmatrix[10], mmatrix[11]);
  printf("mmatrix[12..15] %f, %f, %f, %f \n",
          mmatrix[12], mmatrix[13], mmatrix[14], mmatrix[15]);

  glRotatef(roll, 0.0, 0.0, 1.0);       /* rotate world about Z axis    */
  glGetDoublev(GL_MODELVIEW_MATRIX, mmatrix);
  printf("\n model view matrix, after rotation about X,Y,Z axis \n");
  printf("mmatrix[ 0.. 3] %f, %f, %f, %f \n",
          mmatrix[0], mmatrix[1], mmatrix[2], mmatrix[3]);
  printf("mmatrix[ 4.. 7] %f, %f, %f, %f \n",
          mmatrix[4], mmatrix[5], mmatrix[6], mmatrix[7]);
  printf("mmatrix[ 8..11] %f, %f, %f, %f \n",
          mmatrix[8], mmatrix[9], mmatrix[10], mmatrix[11]);
  printf("mmatrix[12..15] %f, %f, %f, %f \n",
          mmatrix[12], mmatrix[13], mmatrix[14], mmatrix[15]);

  glTranslatef(pos[0], pos[1], pos[2]); /* translate world to new x,y,z */
  glGetDoublev(GL_MODELVIEW_MATRIX, mmatrix);
  printf("\n model view matrix, after rotation and translation \n");
  printf("mmatrix[ 0.. 3] %f, %f, %f, %f \n",
          mmatrix[0], mmatrix[1], mmatrix[2], mmatrix[3]);
  printf("mmatrix[ 4.. 7] %f, %f, %f, %f \n",
          mmatrix[4], mmatrix[5], mmatrix[6], mmatrix[7]);
  printf("mmatrix[ 8..11] %f, %f, %f, %f \n",
          mmatrix[8], mmatrix[9], mmatrix[10], mmatrix[11]);
  printf("mmatrix[12..15] %f, %f, %f, %f \n",
          mmatrix[12], mmatrix[13], mmatrix[14], mmatrix[15]);

  glGetIntegerv (GL_VIEWPORT, viewport);
  gluProject(0.5, 1.5, 2.5, mmatrix, pmatrix, viewport, &wx0, &wy0, &wz0);
  printf("\n 0.5, 1.5, 2.5 at win x=%f, y=%f, z=%f\n", wx0, wy0, wz0);

  { /* not part of OpenGL */
    GLdouble xmin, xmax, ymin, ymax, near, far;
    xmin = -size;
    xmax = size;
    ymin = -size*(GLfloat)win_h/(GLfloat)win_w;
    ymax = size*(GLfloat)win_h/(GLfloat)win_w;
    near = 1.0*size;
    far = 100.0*size;
    gluProject(xmin, ymin, -near, imatrix, pmatrix, viewport, &wx0, &wy0, &wz0);
    printf(" proj %f, %f, %f to win x=%f, y=%f, z=%f\n", xmin, ymin, -near, wx0, wy0, wz0);
    gluProject(xmax, ymax, -near, imatrix, pmatrix, viewport, &wx0, &wy0, &wz0);
    printf(" proj %f, %f, %f to win x=%f, y=%f, z=%f\n", xmax, ymax, -near, wx0, wy0, wz0);
    gluProject(100.0*xmin, 100.0*ymax, -far, imatrix, pmatrix, viewport, &wx0, &wy0, &wz0);
    printf(" proj %f, %f, %f to win x=%f, y=%f, z=%f\n", 100.0*xmin, 100.0*ymax, -far, wx0, wy0, wz0);
    gluProject(100.0*xmax, 100.0*ymin, -far, imatrix, pmatrix, viewport, &wx0, &wy0, &wz0);
    printf(" proj %f, %f, %f to win x=%f, y=%f, z=%f\n", 100.0*xmax, 100.0*ymin, -far, wx0, wy0, wz0);
  }
}


int main(int argc, char *argv[])
{
  printf("output from demo_mat.c OpenGL\n");
  glutInit(&argc, argv);
  glutInitDisplayMode(GLUT_SINGLE | GLUT_RGB | GLUT_DEPTH);
  glutInitWindowSize(300, 400);
  glutCreateWindow(argv[0]);
  glutReshapeFunc(myReshape);
  glutDisplayFunc(display);
  glutMainLoop();
  return 0;
}