import java.awt.*;
import java.awt.image.*;
import Raster;
import Point3D;

/* Matrix3D.java
 *
 * This class implements transformation matrices.  Bascially, instantiations
 *  of this class are matrices that, when applied to a vector representing a
 *  point in 3-dimensional space, transform the point in some way.  The types
 *  of transformations implemented include transaltions, rotations, scaling,
 *  shearing, projecting, etc.  These transformations are composed, meaning
 *  that if a transformation matrix to allow a translation has already been
 *  created, and the shearing method is called, then the matrix is transformed
 *  so as to implement both a translation and a shearing.  They are composed
 *  in such a way to create a transformation stack, meaning that the
 *  transformations are performed in the reverse order as they were composed
 *  on the matrix.  All transformation matrices, when not given other
 *  information, start out as the identity transformation.  3-dimensional
 *  transforms work by giveing a transformation matrix an array of points,
 *  with the output being the transformed array of points.  The class assumes
 *  that the points in 3-dimensional space are given as instantiations of the
 *  class Point3D, although no specific Point3D methods are used inside
 *  Matrix3D.
 *
 * Matthew Lennon
 * Version 1.0
 * November 5, 1998
 */
public class Matrix3D {

  // Members //
  private float m[][];

  // Constructors //

  /* This method initializes a new matrix with the identity transformation */
  public Matrix3D() {
    int i, j;

    m = new float[4][4];

    for(i = 0; i < 4; i++)
      for(j = 0; j < 4; j++)
	if(i == j)
	  m[i][j] = 1.0f;
	else
	  m[i][j] = 0.0f;
  }
  
  /* This method initializes a new matrix so that it is a copy of the matrix
   *  that is used as a source
   */
  public Matrix3D(Matrix3D copy) {
    int i, j;

    m = new float[4][4];

    for(i = 0; i < 4; i++)
      for(j = 0; j < 4; j++)
	m[i][j] = copy.get((i + 1), (j + 1));
  }

  /* This constructor creates a matrix that contains the screen-space
   *  transform appropriate to the raster which is given to it.  Note that
   *  since the matrix is a transformation stack, this screenspace transform
   *  is the last operation that is performed.
   */
  public Matrix3D(Raster r) {
    int i, j;
    float zmax;
    
    /* First, we intialize to the identity matrix just to make sure that the
     *  entries in the matrix is clean
     */
    m = new float[4][4];
    for(i = 0; i < 4; i++)
      for(j = 0; j < 4; j++)
	if(i == j)
	  m[i][j] = 1.0f;
	else
	  m[i][j] = 0.0f;

    /* Now we decide how big to make the z-buffer with the zmax variable */
    zmax = 1.0f;
        
    /* Now we fill in neccessary values */
    set(1, 1, ((float)r.width / 2.0f));
    set(1, 4, ((float)r.width / 2.0f));
    set(2, 2, ((float)r.height / 2.0f));
    set(2, 4, ((float)r.height / 2.0f));
    set(3, 3, (zmax / 2.0f));
    set(3, 4, (zmax / 2.0f));
  }
  
  // Accessors //

  /* This method sets the entry in row i and column j in the matrix to the
   *  value given by value
   */
  public void set(int i, int j, float value) {
    m[i - 1][j - 1] = value;
  }
  
  /* This method returns the value of the entry in row i and column j of the
   *  matrix
   */
  public float get(int i, int j) {
    return m[i - 1][j - 1];
  }

  // Methods //
  
  /* This method takes in two arrays of points, and makes the second array
   *  look like the first array under the transform captured by the matrix.
   *  The subset of points transformed begins at the start index and has the
   *  specified length.
   */
  public void transform(Point3D in[], Point3D out[], int start, int length) {
    int i;
    float w;

    /* We loop through our array, using the specified constraints, and
     *  transform each point
     */
    for(i = 0; i < length; i++) {
      out[start + i] = multVector(in[start + i]);

      /* Next, we need to scale each point to make sure that it is on the
       *  1-plane in 4-space.  This is especially true if the point underwent
       *  a perspective transform.
       */
      w = ((in[start + i].x * m[3][0]) +
	   (in[start + i].y * m[3][1]) +
	   (in[start + i].z * m[3][2]) +
	   m[3][3]);
      if(w != 1.0f) {

	/* In this case, we need to scale the out point to be on the
	 *  appropriate plane in 4-space.  This is simply a scalar division
	 *  because of the way that we defined our spaces.
	 */
	out[start + i].x /= w;
	out[start + i].y /= w;
	out[start + i].z /= w;
      }
    }
  }

  /* This method takes a matrix and adds the matrix given in A to it */
  public final void add(Matrix3D A) {
    int i, j;

    for(i = 0; i < 4; i++)
      for(j = 0; j < 4; j++)
	m[i][j] += A.get((i + 1), (j + 1));
  }

  /* This method takes a matrix and multiplies it by a scalar constant */
  public final void multScalar(float c) {
    int i, j;

    for(i = 0; i < 4; i++)
      for(j = 0; j < 4; j++)
	m[i][j] = c * m[i][j];
  }

  /* This method takes in a point in the form of a Point3D, and returns a new
   *  point which is the original point under the transformation of the matrix
   */
  public Point3D multVector(Point3D p) {
    Point3D q;
    float x, y, z;

    /* First we calculate the new coordinates of the transformed point */
    x = ((p.x * get(1, 1)) +
	 (p.y * get(1, 2)) +
	 (p.z * get(1, 3)) +
	 get(1, 4));
    y = ((p.x * get(2, 1)) +
	 (p.y * get(2, 2)) +
	 (p.z * get(2, 3)) +
	 get(2, 4));
    z = ((p.x * get(3, 1)) +
	 (p.y * get(3, 2)) +
	 (p.z * get(3, 3)) +
	 get(3, 4));

    /* Now we construct a point with these new coordinates */
    q = new Point3D(x, y, z);

    return q;
  }

  /* This method takes the matrix and multiplies it by the matrix src.  It is
   *  important to note that this multipication is right composed, meaning
   *  that this = this * src.
   */
  public final void compose(Matrix3D src) {
    int i, j;
    Matrix3D A;

    /* Because we will be changing the values in the matrix, but we will still
     *  need the old values at later times, we store the current values as
     *  another matrix.
     */
    A = new Matrix3D(this);

    /* Now we perform the actual multipication */
    for(i = 1; i <= 4; i++)
      for(j = 1; j <= 4; j++)
	m[i - 1][j - 1] = ((A.get(i, 1) * src.get(1, j)) +
			   (A.get(i, 2) * src.get(2, j)) +
			   (A.get(i, 3) * src.get(3, j)) +
			   (A.get(i, 4) * src.get(4, j)));
  }

  /* This method makes the matrix the identity matrix, which would make points
   *  undergo the identity transformation
   */
  public void loadIdentity() {
    int i, j;
    
    for(i = 0; i < 4; i++)
      for(j = 0; j < 4; j++)
	if(i == j)
	  m[i][j] = 1.0f;
	else
	  m[i][j] = 0.0f;
  }

  /* This method takes the matrix and makes it reflect a translation along the
   *  vector (tx, ty, tz)
   */
  public void translate(float tx, float ty, float tz) {
    Matrix3D A;
    
    /* First we make the translation matrix */
    A = new Matrix3D();
    A.set(1, 4, tx);
    A.set(2, 4, ty);
    A.set(3, 4, tz);

    /* Now we compose the translation with the already existing transformation
     *  matrix
     */
    compose(A);
  }
  
  /* This method takes the matrix and makes it represent a scaling of sx in the
   *  x direction, sy in the y direction, and sz in the z direction
   */
  public void scale(float sx, float sy, float sz) {
    Matrix3D A;

    /* First we make the scaling matrix */
    A = new Matrix3D();
    A.set(1, 1, sx);
    A.set(2, 2, sy);
    A.set(3, 3, sz);

    /* Now we compose the translation with the already existing transformation
     *  matrix
     */
    compose(A);
  }
  
  /* This method takes the matrix and makes it represent a shear with the
   *  quantization kxy, kxz, and kyz
   */
  public void skew(float kxy, float kxz, float kyz) {
    Matrix3D A;

    /* First we make the shearing matrix */
    A = new Matrix3D();
    A.set(1, 2, kxy);
    A.set(1, 3, kxz);
    A.set(2, 3, kyz);

    /* Now we compose this shearing matrix with the already existing
     *  transformation matrix
     */
    compose(A);
  }

  /* This method takes the matrix and makes it represent a rotation of angle
   *  around a rotation axis described by the vector (ax, ay, az)
   */
  public void rotate(float ax, float ay, float az, float angle) {
    Matrix3D A, temp;
    Point3D axis;
    float x, y, z, t;

    /* First, we make the rotation matrix.  The first step in doing that is
     *  normalize the rotation axis and to make it into a vector.  We can use
     *  the Point3D class to also represent a vector.
     */
    t = (float)Math.sqrt((double)((ax * ax) + (ay * ay) + (az * az)));
    x = ax / t;
    y = ay / t;
    z = ay / t;
    axis = new Point3D(x, y, z);

    /* Next we begin to construct the actual matrix, which can be done in three
     *  parts.  For the first part we take the symmetric of the axis multiplied
     *  by (1 - cos angle).
     */
    A = new Matrix3D();
    A.set(1, 1, (axis.x * axis.x));
    A.set(1, 2, (axis.x * axis.y));
    A.set(1, 3, (axis.x * axis.z));
    A.set(2, 1, (axis.y * axis.x));
    A.set(2, 2, (axis.y * axis.y));
    A.set(2, 3, (axis.y * axis.z));
    A.set(3, 1, (axis.z * axis.x));
    A.set(3, 2, (axis.z * axis.y));
    A.set(3, 3, (axis.z * axis.z));
    t = (1.0f - (float)Math.cos((double)angle));
    A.multScalar(t);

    /* The next part of the rotation matrix is made up of the skew of the axis
     *  multiplied by the sine of the angle.  After making the second part we
     *  add it to the previous part.
     */
    temp = new Matrix3D();
    temp.set(1, 1, 0.0f);
    temp.set(1, 2, (-1.0f * axis.z));
    temp.set(1, 3, axis.y);
    temp.set(2, 1, axis.z);
    temp.set(2, 2, 0.0f);
    temp.set(2, 3, (-1.0f * axis.x));
    temp.set(3, 1, (-1.0f * axis.y));
    temp.set(3, 2, axis.x);
    temp.set(3, 3, 0.0f);
    t = (float)Math.sin((double)angle);
    temp.multScalar(t);
    A.add(temp);

    /* The third part of the rotation matrix is given by the identity matrix
     *  multiplied by cos angle
     */
    temp.loadIdentity();
    t = (float)Math.cos((double)angle);
    temp.multScalar(t);
    A.add(temp);

    /* Now we make sure that the final row is (0, 0, 0, 1) as it should be */
    A.set(4, 1, 0.0f);
    A.set(4, 2, 0.0f);
    A.set(4, 3, 0.0f);
    A.set(4, 4, 1.0f);
    
    /* Now that we have the rotation transformation, which is represented by A,
     *  we can compose it with the already existing transforation matrix
     */
    compose(A);
  }

  /* This method takes a matrix and makes it represent a viewing transform.  In
   *  other words, it takes realworld coordinates and moves them to eye-space
   *  with the eye located at point (eyex, eyey, eyez), looking at point (atx,
   *  aty, atz), and with the up direction being defined as (upx, upy, upz).
   */
  public void lookAt(float eyex, float eyey, float eyez,
		     float atx,  float aty,  float atz,
		     float upx,  float upy,  float upz) {
    Matrix3D A;
    Point3D l, r, u;
    float x, y, z, t;

    /* First we need to make a matrix that will represent a transform from
     *  world space to the specified camera space.  To do this, we need to
     *  compute three special vectors, and use those to compute the matrix.
     *  We can use instantiations of the Point3D class to represent these
     *  vectors.  First we compute a vector in the lookat direction, and
     *  normalize that vector.
     */
    x = atx - eyex;
    y = aty - eyey;
    z = atz - eyez;
    t = (float)Math.sqrt((double)((x * x) + (y * y) + (z * z)));
    x /= t;
    y /= t;
    z /= t;
    l = new Point3D(x, y, z);

    /* Next we compute the vector in the right-hand direction and normalize
     *  that
     */
    x = ((l.y * upz) - (l.z * upy));
    y = ((l.z * upx) - (l.x * upz));
    z = ((l.x * upy) - (l.y * upx));
    t = (float)Math.sqrt((double)((x * x) + (y * y) + (z * z)));
    x /= t;
    y /= t;
    z /= t;
    r = new Point3D(x, y, z);
    
    /* Next, we compute the vector in the up direction (to make sure we have
     *  orthoginal space) and normalize that vector
     */
    x = ((r.y * l.z) - (r.z * l.y));
    y = ((r.z * l.x) - (r.x * l.z));
    z = ((r.x * l.y) - (r.y * l.x));
    t = (float)Math.sqrt((double)((x * x) + (y * y) + (z * z)));
    x /= t;
    y /= t;
    z /= t;
    u = new Point3D(x, y, z);

    /* Now we calculate a few dot products as added values for the matrix
     *  entries
     */
    x = (((-1.0f * r.x * eyex) + (-1.0f * r.y * eyey) + (-1.0f * r.z * eyez)));
    y = (((-1.0f * u.x * eyex) + (-1.0f * u.y * eyey) + (-1.0f * u.z * eyez)));
    z = ((l.x * eyex) + (l.y * eyey) + (l.z * eyez));

    /* Now we can finally construct the matrix */
    A = new Matrix3D();
    A.set(1, 1, r.x);
    A.set(1, 2, r.y);
    A.set(1, 3, r.z);
    A.set(1, 4, x);
    A.set(2, 1, u.x);
    A.set(2, 2, u.y);
    A.set(2, 3, u.z);
    A.set(2, 4, y);
    A.set(3, 1, (-1.0f * l.x));
    A.set(3, 2, (-1.0f * l.y));
    A.set(3, 3, (-1.0f * l.z));
    A.set(3, 4, z);

    /* Last, we compose this transformation with any other transformations that
     *  exist already in the matrix
     */
    compose(A);
  }
  
  /* This method makes the matrix represent a perspective transform into
   *  canonical space.  The viewing frustum is represented by the coordinates
   *  which are fed into the method.
   */
  public void perspective(float left, float right,
			  float bottom, float top,
			  float near, float far) {
    Matrix3D A;

    /* First we construct the matrix which accomplishes the transform from the
     *  viewing frustum into canonical space
     */
    A = new Matrix3D();
    A.set(1, 1, ((2.0f * near) / (right - left)));
    A.set(1, 3, ((-1.0f * (right + left)) / (right - left)));
    A.set(2, 2, ((2.0f * near) / (bottom - top)));
    A.set(2, 3, ((-1.0f * (bottom + top)) / (bottom - top)));
    A.set(3, 3, ((far + near) / (far - near)));
    A.set(3, 4, ((-2.0f * far * near) / (far - near)));
    A.set(4, 3, 1.0f);
    A.set(4, 4, 0.0f);	  

    /* Now we compose the viewing transform with the other transforms already
     *  represented by the matrix
     */
    compose(A);
  }

  /* This method makes the matrix represent an orthographic transform into
   *  canonical space.  The viewing frustum is represented by the coordinates
   *  which are fed into the method.
   */
  public void orthographic(float left, float right,
			   float bottom, float top,
			   float near, float far) {
    Matrix3D A;

    /* First we construct the matrix which accomplishes the transform from the
     *  viewing frustum into canonical space
     */
    A = new Matrix3D();
    A.set(1, 1, (2.0f / (right - left)));
    A.set(1, 4, ((-1.0f * (right + left)) / (right - left)));
    A.set(2, 2, (2.0f / (bottom - top)));
    A.set(2, 4, ((-1.0f * (bottom + top)) / (bottom - top)));
    A.set(3, 3, (2.0f / (far - near)));
    A.set(3, 4, ((-1.0f * (far + near)) / (far - near)));

    /* Now we compose the viewing transform with the other transforms already
     *  represented by the matrix
     */
    compose(A);
  }
}