// static GLfloat radius;          // camera radius
static GLfloat near;            
static GLfloat far;

static twTriple VPN;            // View Plane Normal
static twTriple VUP;            // the ``up'' vector
static twTriple VRP;        // the View Reference Point (the eye location)

static GLfloat FieldOfView = DEFAULT_FOVY; // for square frustum
static GLfloat FieldOfViewY = DEFAULT_FOVY; // arg for gluPerspective

static const int DEPTH_BITS_TO_LOSE = 5; // see doc.tex for more info

// spinning stuff.  We can spin around any of the three major axes
typedef enum {X, Y, Z} axis_t;
axis_t axis = X;                // default is to spin around X
static GLfloat twSpinAngle = 0.5; // default angle per frame

void twNearFarSet(GLfloat& myNear, GLfloat& myFar) {
    myNear = near;
    myFar = far;
}

void twFovySet(GLfloat& myFovy) {
    myFovy = FieldOfViewY;
}

void twSpin () {
  GLfloat turn = twSpinAngle;
  twTriple xAxis={1,0,0};
  twTriple yAxis={0,1,0};
  twTriple zAxis={0,0,1};

  switch (axis) {
  case X: twRotateViewpoint(turn,xAxis); break;
  case Y: twRotateViewpoint(turn,yAxis); break;
  case Z: twRotateViewpoint(turn,zAxis); break;
  }
  glutPostRedisplay();
}

void startSpinning (axis_t a) {
  axis = a;
  twIdleFunc(twSpin);
}

void twSpinCommand (unsigned char key, int x, int y) {
  switch (key) {
  case '-': twSpinAngle *= 0.5; printf("twSpinAngle=%f\n",twSpinAngle); break; // halves the speed
  case '+': twSpinAngle *= 2.0; printf("twSpinAngle=%f\n",twSpinAngle); break; // doubles the speed
  case 'x': startSpinning(X); break; /* moves image about the x axis */
  case 'y': startSpinning(Y); break; /* moves image about the y axis */
  case 'z': startSpinning(Z); break; /* moves image about the z axis */
  }
}

// gets the modelview matrix and stores it into the given matrix.
void twGetModelView(GLdouble MV[16]) {
    glGetDoublev(GL_MODELVIEW_MATRIX,MV);
}

// gets the projection matrix and stores it into the given matrix.
void twGetProjection(GLdouble P[16]) {
    glGetDoublev(GL_PROJECTION_MATRIX,P);
}

// gets the viewport array and stores it into the given array
void twGetViewport(GLint V[4]) {
    glGetIntegerv(GL_VIEWPORT,V);
}

// computes the window projection of the vertex v and stores it in w.
// That is, it computes where on the screen the point v will project to.
void twProject(twTriple w, twTriple v) {
  GLdouble v0, v1, v2;
    GLdouble modelview[16];
    GLdouble proj[16];
    GLint view[4];
    GLdouble winx, winy, winz;

    winx=winy=winz=-2;
    twGetModelView(modelview);
    //twPrintMatrix4x4("MV",modelview);
    twGetProjection(proj);
    //twPrintMatrix4x4("P",proj);
    twGetViewport(view);
    //printf("viewport = %d %d %d %d\n",view[0],view[1],view[2],view[3]);
    v0=v[0];
    v1=v[1];
    v2=v[2];
    if(gluProject(v0,v1,v2,modelview,proj,view,&winx,&winy,&winz)!=GL_TRUE) {
        printf("Error in gluProject for world location (%f,%f,%f)\n",
               v[0],v[1],v[2]);
        exit(1);
    }
    twTripleInit(w,winx,winy,winz);
    // the following doesn't work.  
#ifdef TW_DEBUG
    // The following is a hand check of gluProject, to see how it works.
    printf("%f %f %f is the proj of V=%f %f %f\n",w[0],w[1],w[2],v[0],v[1],v[2]);
    mult3d(v,modelview,v);
    mult3d(v,proj,v);
    v[0]=view[0]+view[2]*(v[0]+1)/2;
    v[1]=view[1]+view[3]*(v[1]+1)/2;
    v[2]=(v[2]+1)/2;
    printf("%f %f %f is the hand double check\n\n",v[0],v[1],v[2]);
#endif    
}

// computes v, a world vertex coordinate corresponding to window location w
void twUnProject(twTriple v, twTriple w) {
    GLdouble modelview[16];
    GLdouble proj[16];
    GLint view[4];
    GLdouble winx=w[0],winy=w[1],winz=w[2];
    GLdouble worldx,worldy,worldz;

    twGetModelView(modelview);
    //    twPrintMatrix4x4("MV",modelview);
    twGetProjection(proj);
    //twPrintMatrix4x4("P",proj);
    twGetViewport(view);
    if(gluUnProject(winx,winy,winz,modelview,proj,view,&worldx,&worldy,&worldz)!=GL_TRUE) {
        printf("Error in gluUnProject for window location (%f,%f,%f)\n",winx,winy,winz);
        twPrintMatrix4x4("modelview",modelview);
        twPrintMatrix4x4("proj",proj);
        exit(1);
    }
    twTripleInit(v,worldx,worldy,worldz);
}

// ================================================================

// Stores into M a rotation matrix of angle degrees around vector (x,y,z).
// Normalizes the vector.  If the vector is zero length, returns
// immediately without modifying M.  From the GL reference manual for
// glRotate.
void rotationMatrix(GLfloat* M, GLfloat angle, GLfloat x, GLfloat y, GLfloat z) {
    GLfloat r=angle*M_PI/180.0;
    GLfloat c=cos(r);
    GLfloat s=sin(r);
    GLfloat l=x*x+y*y+z*z;
    
    if(l==0.0) return;

    x=x/l;
    y=y/l;
    z=z/l;

    M[0]=x*x*(1-c)+c;
    M[1]=y*x*(1-c)+z*s;
    M[2]=x*z*(1-c)-y*s;
    M[3]=0;
    M[4]=x*y*(1-c)-z*s;
    M[5]=y*y*(1-c)+c;
    M[6]=y*z*(1-c)+x*s;
    M[7]=0;
    M[8]=x*z*(1-c)+y*s;
    M[9]=y*z*(1-c)-x*s;
    M[10]=z*z*(1-c)+c;
    M[11]=0;
    M[12]=0;
    M[13]=0;
    M[14]=0;
    M[15]=1;
}
	
// ================================================================

// Rotates the viewpoint (changing the VRP, VPN and VUP so that it still
// looks at the center, but from a different place) by the given angle
// around the given rotation vector.  Returns if angle is zero.
void twRotateViewpoint(GLfloat angle, twTriple r) {
    GLfloat M[16];              // temporary for rotation matrix

    if(0==angle) return;
    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    twVectorNormalize(r);
    // printf("rotate by %f around %f %f %f\n",angle, r[0],r[1],r[2]);
    if( !Toggles[IMMERSE] ) 
        twVector(VPN,BBCenter,VRP);
    rotationMatrix(M,angle,r[0],r[1],r[2]);
    mult3(VPN,M,VPN);
    mult3(VUP,M,VUP);
    if( !Toggles[IMMERSE] ) {
        // recalculate VRP
        twTriple VPN_reverse;
        twVectorScale(VPN_reverse,VPN,-1);
        twPoint(VRP,BBCenter,VPN_reverse);
    }
    glutPostRedisplay();    
}    

// Rotates the view plane normal (VPN), keeping the VRP fixed, so that we
// are tilting or panning the camera. Rotates by the given angle around
// the given rotation vector.  Returns if angle is zero.  Leaves VUP
// unchanged.  If the VPN would end up outside the original frustum, it is
// limited to that angle.  This is done by limiting to 45 degrees from the
// vector from the VPN to the Center.
void twRotateVPN(GLfloat angle, twTriple r) {
    GLfloat M[16];              // temporary for rotation matrix
    twTriple newVPN;

    if(0==angle) return;
    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    twVectorNormalize(r);
    // printf("tilt/pan by %f around %f %f %f\n",angle, r[0],r[1],r[2]);
    rotationMatrix(M,angle,r[0],r[1],r[2]);
    mult3(newVPN,M,VPN);

    // Code to check limits on rotation
    twTriple originalVPN;
    twVector(originalVPN,BBCenter,VRP);
    GLfloat cosTotalAngle = twCosAngle(originalVPN,newVPN);
    if(cosTotalAngle < cos(M_PI/4.0)) {
        // printf("limiting rotation to 45 degrees\n");
        rotationMatrix(M,45,r[0],r[1],r[2]);
        mult3(newVPN,M,originalVPN);
    }
    twTripleCopy(VPN,newVPN);

    glutPostRedisplay();    
}    

// Rotates the viewpoint by an angle defined by two window locations:
// (Ax,Ay) and (Bx,By).  The y coordinate of the window locations should
// already have been subtracted from the window height.  This function
// retrieves the modelview matrix and sets originalView to false
void twTrackballOrientation_old(int Ax, int Ay, int Bx, int By) {
    twTriple A, B, C;
    twTriple v,w,n;
    GLfloat angle;
    GLfloat M[16];              // temporary for rotation matrix

    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    twTriple winA = {Ax, Ay, 0};
    twTriple winB = {Bx, By, 0};
    twUnProject(A,winA);
    twUnProject(B,winB);
    twVector(v,A,BBCenter);
    twVector(w,B,BBCenter);
#ifdef TW_DEBUG
    twTriplePrint("A",A);
    twTriplePrint("B",B);
    twTriplePrint("C",BBCenter);
    printf("AB = %f, AC = %f BC = %f\n",
           twPointDistance(A,B),
           twPointDistance(A,BBCenter),
           twPointDistance(B,BBCenter));
    twTriplePrint("v",v);
    twTriplePrint("w",w);
#endif
    angle = acos(twCosAngle(v,w))/M_PI*180.0;
    twCrossProduct(n,v,w);
    glTranslatef(BBCenter[0],BBCenter[1],BBCenter[2]);
    glRotated(angle,n[0],n[1],n[2]);
    glTranslatef(-BBCenter[0],-BBCenter[1],-BBCenter[2]);
}    

// Rotates the viewpoint by an angle defined by two window locations:
// (Ax,Ay) and (Bx,By).  The y coordinate of the window locations should
// already have been subtracted from the window height.  
void twTrackballOrientation(int Ax, int Ay, int Bx, int By) {
    twTriple A, B, C;
    twTriple v,w,n;
    GLfloat angle;

    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    twTriple winA = {Ax, Ay, 0};
    twTriple winB = {Bx, By, 0};
    twUnProject(A,winA);
    twUnProject(B,winB);
    twVector(v,A,BBCenter);
    twVector(w,B,BBCenter);
    angle = acos(twCosAngle(v,w))/M_PI*180.0;
    twCrossProduct(n,v,w);
    twRotateViewpoint(-angle,n);
}    

// ================================================================

// Modify the VPN so that the (x,y) pixel coordinates are in the center of
// the screen.  Ultimately calls RotateVPN.
void twOrientVPN(GLint x, GLint y) {
    twTriple C;
    twTriple win = {x, y, 0};
    twTriple v,n;
    GLfloat angle;
    GLfloat M[16];              // temporary for rotation matrix
    GLfloat w[4];

    twUnProject(C,win);
    twVector(v,C,VRP);          // v is from VRP to C
    angle = acos(twCosAngle(v,VPN))/M_PI*180.0;
    twCrossProduct(n,VPN,v);
    twVectorNormalize(n);
    //printf("change VPN by %f degrees around (%f,%f,%f)\n", -angle,n[0],n[1],n[2]);
    twRotateVPN(angle,n);
    return;
}
    
void twPerspective(GLfloat fovy, GLfloat ar,char* kind) {
    GLfloat fw, fh;

    twMessage(TW_CAMERA,"gluPerspective: fovy = %f ar = %f near = %f far = %f\n",
              fovy,ar,near,far);
    gluPerspective(fovy,ar,near,far);
    FieldOfViewY = fovy;        // store into global
    // yes, I know this recomputes things that are already computed, but I
    // think the modularity is best this way.
    fh = tan(fovy/2)*near*2;
    fw = fh*ar;
    twMessage(TW_CAMERA,"frustum WxH = %f x %f window %4d x %4d\n",
              fw, fh, windowWidth, windowHeight);
}


// This converts FOVX into FOVY, both in degrees
GLfloat twFOVX2FOVY(GLfloat fovx) {
    GLfloat frustumWidth, frustumHeight;
    
    // For our code, this is a common special case.  tan(90)=1
    if(90.0==fovx) 
        frustumWidth=near*2;
    else {
        // convert to radians
        frustumWidth=near*2*tan(fovx/2*M_PI/180.0);
    }
    frustumHeight=frustumWidth/aspectRatio;
    GLfloat fovy = 2*atan(frustumHeight/2/near)*180/M_PI;
    if(fovy<0) {
        printf("ack! fovy<0.  fovx=%f fovy=%f\n", fovx,fovy);
        printf("frustumWidth=%f frustumHeight=%f AR=%f\n",
               frustumWidth,frustumHeight,aspectRatio);
    }
    return fovy;
}

// function to set up the camera shape, based on global variables
// twAspectRatio, twfrustumMode, near.  If the argument is non-zero, the
// frustum shape is printed.

void twCameraShape() {
    glMatrixMode(GL_PROJECTION);
    glLoadIdentity();
    
    switch (frustumMode) {
    case DISTORT:
        twPerspective(FieldOfView,1,"DISTORT");
        break;

    case LETTERBOX:
        // To letterbox, the smaller angle of the frustum is determined by
        // the radius of the bounding sphere, thereby guaranteeing that
        // everything is inside it.  The larger angle will capture empty
        // space.
        if(aspectRatio<1) {
            // portrait, so the vertical angle is larger, and we have to
            // calculate it from the smaller angle, which is 90.
            twPerspective(twFOVX2FOVY(FieldOfView),aspectRatio,"LETTERBOX");
        } else {
            // landscape
            twPerspective(FieldOfView,aspectRatio,"LETTERBOX");
        }
        break;

    case CLIP:
        // To clip, the larger angle of the frustum is determined by the
        // radius of the bounding sphere, so that the other angle is
        // smaller than 90, cutting off part of the sphere.
        if(aspectRatio<1) {
            // portrait, so the vertical angle is 90 degrees
            twPerspective(FieldOfView,aspectRatio,"CLIP");
        } else {
            // landscape, so the vertical angle is determined from the
            // smaller, which is 90 degrees.
            twPerspective(twFOVX2FOVY(FieldOfView),aspectRatio,"LETTERBOX");
        }
        break;
    }
}

void twAxes() {
    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    glPushAttrib(GL_ALL_ATTRIB_BITS);
    glDisable(GL_LIGHTING);
    glPushMatrix();
    twTranslate(BBCenter);
    glEnable(GL_LINE_STIPPLE);
    glLineStipple(3, 0xcccc);
    glBegin(GL_LINES);
    glColor3f(1,0,0); // red for x axis
    glVertex3f(OuterRadius,0,0);
    glVertex3f(0,0,0);
    glColor3f(0,1,0); // green for y axis
    glVertex3f(0,OuterRadius,0);
    glVertex3f(0,0,0);
    glColor3f(0,0,1); // blue for z axis
    glVertex3f(0, 0,OuterRadius);
    glVertex3f(0, 0,0);
    glEnd();
    glDisable(GL_LINE_STIPPLE);
    glPopMatrix();
    glPopAttrib();
}

void twOriginalView() {
    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    FieldOfView = FieldOfViewY = DEFAULT_FOVY;
    GLfloat eyeRadius = OuterRadius*M_SQRT2;
    far  = eyeRadius+OuterRadius;
    twTripleCopy(VRP,BBCenter);
    if( Toggles[IMMERSE] ) {
        // InnerRadius/M_SQRT2 is the largest possible value of near.
        // Smaller means less stuff is cut off, at the price of, maybe
        // some problems with the depth buffer.
        // near = InnerRadius/M_SQRT2/2;
        near = far/pow(2,DEPTH_BITS_TO_LOSE);
        twMessage(TW_CAMERA,
                  "twCameraPosition: INSIDE InnerRadius = %f, near = %f far = %f\n",
                  InnerRadius, near, far);
    } else {
        VRP[2] += eyeRadius;
        near = eyeRadius-OuterRadius;
        twMessage(TW_CAMERA,
                  "twCameraPosition: OUTSIDE OuterRadius = %f, near = %f far = %f\n",
                  OuterRadius, near, far);
    }
    twTripleInit(VPN,0,0,-1);
    twTripleInit(VUP,0,1,0);
}

void twInitView(int axis) {
    if(!BoundingBoxInitialized) {
        printf("Bounding Box not initialized\n");
        return;
    }
    FieldOfView = FieldOfViewY = DEFAULT_FOVY;
    GLfloat eyeRadius = OuterRadius*M_SQRT2;
    far  = eyeRadius+OuterRadius;
    twTripleCopy(VRP,BBCenter);
    if( Toggles[IMMERSE] ) {
        // InnerRadius/M_SQRT2 is the largest possible value of near.
        // Smaller means less stuff is cut off, at the price of, maybe
        // some problems with the depth buffer.
        // near = InnerRadius/M_SQRT2/2;
        near = far/pow(2,DEPTH_BITS_TO_LOSE);
        twMessage(TW_CAMERA,
                  "twCameraPosition: INSIDE InnerRadius = %f, near = %f far = %f\n",
                  InnerRadius, near, far);
    } else {
        VRP[axis] += eyeRadius;
        near = eyeRadius-OuterRadius;
        twMessage(TW_CAMERA,
                  "twCameraPosition: OUTSIDE OuterRadius = %f, near = %f far = %f\n",
                  OuterRadius, near, far);
    }
}

void twZview() {
    twInitView(2);
    twTripleInit(VPN,0,0,-1);
    twTripleInit(VUP,0,1,0);
}

void twYview() {
    twInitView(1);
    twTripleInit(VPN,0,-1,0);
    twTripleInit(VUP,0,0,-1);
}

void twXview() {
    twInitView(0);
    twTripleInit(VPN,-1,0,0);
    twTripleInit(VUP,0,1,0);
}

void twViewCommand (unsigned char key, int x, int y) {
  switch (key) {
  case 'X': twXview(); break;
  case 'Y': twYview(); break;
  case 'Z': twZview(); break;
  }
  glutPostRedisplay();
}


void twCameraPosition() {
    glMatrixMode(GL_MODELVIEW);
    glLoadIdentity();
    if(twMessageKinds & TW_CAMERA) {
        twTriplePrint("VRP",VRP);
        twTriplePrint("VPN",VPN);
        twTriplePrint("VUP",VUP);
    }
    gluLookAt(VRP[0],VRP[1],VRP[2],
              VRP[0]+VPN[0],VRP[1]+VPN[1],VRP[2]+VPN[2],
              VUP[0],VUP[1],VUP[2]);
}  

// Sets the frustum mode to LETTERBOX, DISTORT, or CLIP
void twFrustumMode(frustumMode_t mode) {
    frustumMode = mode;
    if (mode == LETTERBOX) printf("Letterbox mode\n");
    else if (mode == CLIP) printf("clipping mode\n");
    else printf("distorting mode\n");
}

void twCamera() {  
    twError();
    glViewport(0,0,windowWidth,windowHeight);
    if( Toggles[IMMERSE] ) {
        twMessage(TW_CAMERA,"immerse:  InnerRadius = %f, near=%f far=%f\n",
                  InnerRadius, near, far);
    } else {
        twMessage(TW_CAMERA,"back off:  OuterRadius = %f, near=%f far=%f\n",
                  OuterRadius, near, far);
    }
    twError();
    twCameraShape();
    twError();
    twCameraPosition();
    twError();
    if( Toggles[AXES] ) twAxes();
    if( Toggles[BB] ) twDrawBoundingBox();
    twError();
}