import java.awt.*;
import java.awt.image.*;


// 4x4 matrix


public class Matrix3D {


  // instance variables
  float[][] matrix;


  //
  // Constructors
  //

public Matrix3D() {

  matrix = new float[4][4];
  // initialize with identity transform
  loadIdentity();

}


public Matrix3D(Matrix3D copy) {

  matrix = new float[4][4];
  // initialize with copty of source
  for (int i=0; i<4; i++)
    for (int j=0; j<4; j++)
      matrix[j][i] = copy.get(i, j);

}


public Matrix3D(Raster r) {

  this();

  int width = r.getWidth()/2;
  int height = r.getHeight()/2;

  // initialize with a mapping from canonical space to screen space
  matrix[0][0] = width;
  matrix[0][3] = width;
  matrix[1][1] = height;
  matrix[1][3] = height;

}


  //
  // General interface methods
  //


public void set(int i, int j, float value) {

  // set element [j][i] to value
  matrix[j][i] = value;

}


public float get(int i, int j) {

  // return element [j][i]
  return matrix[j][i];

}


  //
  // Helper methods
  //


public float dotProduct(Point3D a, Point3D b) {

  return a.x*b.x+a.y*b.y+a.z*b.z;

}
  

public Point3D crossProduct(Point3D a, Point3D b) {

  Matrix3D ma = new Matrix3D();

  // initialize ma with point a
  ma.set(0, 0, 0);
  ma.set(1, 0, -a.z);
  ma.set(2, 0, a.y);
  ma.set(0, 1, a.z);
  ma.set(1, 1, 0);
  ma.set(2, 1, -a.x);
  ma.set(0, 2, -a.y);
  ma.set(1, 2, a.x);
  ma.set(2, 2, 0);

  // multiply ma by b
  return ma.multiplyPoint(b);

}


public Point3D multiplyPoint(Point3D a) {

  float x, y, z, w;

  // multiply a 3x3 matrix by a point a
  x = matrix[0][0]*a.x + matrix[0][1]*a.y + matrix[0][2]*a.z + matrix[0][3];
  y = matrix[1][0]*a.x + matrix[1][1]*a.y + matrix[1][2]*a.z + matrix[1][3];
  z = matrix[2][0]*a.x + matrix[2][1]*a.y + matrix[2][2]*a.z + matrix[2][3];
  w = matrix[3][0]*a.x + matrix[3][1]*a.y + matrix[3][2]*a.z + matrix[3][3];

  x=x/w;
  y=y/w;
  z=z/w;

  Point3D out = new Point3D(x, y, z);
  return out;

}


public void multiplyScalar(float scalar) {

  // this = scalar * this
  for (int i=0; i<4; i++)
    for (int j=0; j<4; j++)
      matrix[j][i] *= scalar;

}


public void addMatrix(Matrix3D src) {

  // this = this + src
  for (int i=0; i<4; i++)
    for (int j=0; j<4; j++)
      matrix[j][i] += src.get(i, j);

}
    

  //
  // Transformation methods
  //


public void transform(Point3D in[], Point3D out[], int start, int length) {

  // transform points from in to out using current matrix
  for (int i=0; i<length; i++) {
    out[start+i] = multiplyPoint(in[start+i]);
  }
}


public final void compose(Matrix3D src) {


  // initialize new matrix with the current one
  Matrix3D temp = new Matrix3D(this);
  
  // matrix = temp*src
  for (int r=0; r<4; r++)
    for (int c=0; c<4; c++)
      matrix[r][c]= temp.get(0, r)*src.get(c, 0) +
	temp.get(1, r)*src.get(c, 1) +
	temp.get(2, r)*src.get(c, 2) +
	temp.get(3, r)*src.get(c, 3);

}


public void loadIdentity() {

  // this = identity matrix
  for (int r=0; r<4; r++)
    for (int c=0; c<4; c++)
      matrix[r][c] = (r==c) ? 1 : 0;

}


public void translate(float tx, float ty, float tz) {

  // set up a translation matrix
  Matrix3D t = new Matrix3D();
  t.set(3, 0, tx);
  t.set(3, 1, ty);
  t.set(3, 2, tz);

  // this = this * t
  compose(t);

}


public void scale(float sx, float sy, float sz) {

  // set up a scale matrix
  Matrix3D s = new Matrix3D();
  s.set(0, 0, sx);
  s.set(1, 1, sy);
  s.set(2, 2, sz);

  // this = this * s
  compose(s);

}


public void skew(float kxy, float kxz, float kyz) {

  // set up a skew matrix
  Matrix3D s = new Matrix3D();
  s.set(1, 0, kxy);
  s.set(2, 0, kxz);
  s.set(2, 1, kyz);

  // this = this * s
  compose(s);

}


public void rotate(float ax, float ay, float az, float angle) {

  // normalize the points
  float length = 1/(float)(Math.sqrt(ax*ax + ay*ay + az*az));
  ax *= length;
  ay *= length;
  az *= length;
  
  // symmetric matrix
  Matrix3D syma = new Matrix3D();
  syma.set(0, 0, ax*ax);
  syma.set(1, 0, ax*ay);
  syma.set(2, 0, ax*az); 
  syma.set(0, 1, ax*ay);
  syma.set(1, 1, ay*ay);
  syma.set(2, 1, ay*az);
  syma.set(0, 2, ax*az);
  syma.set(1, 2, ay*az);
  syma.set(2, 2, az*az);

  Matrix3D skewa = new Matrix3D();
  skewa.set(0, 0, 0);
  skewa.set(1, 1, 0);
  skewa.set(2, 2, 0);
  skewa.set(1, 0, -az);
  skewa.set(2, 0, ay);
  skewa.set(0, 1, az);
  skewa.set(2, 1, -ax);
  skewa.set(0, 2, -ay);
  skewa.set(1, 2, ax);

  // identity matrix
  Matrix3D i = new Matrix3D();

  // rotate = syma*(1-cos(angle)) + skewa*sin(angle) + i*cos(angle)
  Matrix3D rotate = new Matrix3D();
  rotate.set(0, 0, 0);
  rotate.set(1, 1, 0);
  rotate.set(2, 2, 0);
  
  // multiply syma by 1-cos(angle)
  float scalar = 1 - (float)(Math.cos(angle));
  syma.multiplyScalar(scalar);
  // multiply skewa by sin(angle)
  scalar = (float)(Math.sin(angle));
  skewa.multiplyScalar(scalar);
  // multiply i by cos(angle)
  scalar = (float)(Math.cos(angle));
  i.multiplyScalar(scalar);
  // add three matrices to get rotate
  rotate.addMatrix(i);
  rotate.addMatrix(skewa);
  rotate.addMatrix(syma);
  rotate.set(3, 3, 1);

  // this = this * rotate
  compose(rotate);

}


public void lookAt(float eyex, float eyey, float eyez,
		   float atx, float aty, float atz,
		   float upx, float upy, float upz) {

  // setup l
  float lx = atx - eyex;
  float ly = aty - eyey;
  float lz = atz - eyez;
  float length = (float)(1/Math.sqrt(lx*lx + ly*ly + lz*lz));
  lx *= length;
  ly *= length;
  lz *= length;
  Point3D l = new Point3D(lx, ly, lz);

  // setup u
  length = (float)(1/Math.sqrt(upx*upx + upy*upy + upz*upz));
  upx *= length;
  upy *= length;
  upz *= length;
  Point3D upminus = new Point3D(-upx, -upy, -upz);
  
  // setup r
  Point3D r = new Point3D();
  r = crossProduct( l, new Point3D(upx, upy, upz));
  float rx = r.x;
  float ry = r.y;
  float rz = r.z;
  Point3D rminus = new Point3D(-rx, -ry, -rz);

  // setup lookat
  Point3D eye = new Point3D(eyex, eyey, eyez);
  Matrix3D lookat = new Matrix3D();
  lookat.set(0, 0, rx);
  lookat.set(1, 0, ry);
  lookat.set(2, 0, rz);
  lookat.set(3, 0, dotProduct(rminus, eye));
  lookat.set(0, 1, upx);
  lookat.set(1, 1, upy);
  lookat.set(2, 1, upz);
  lookat.set(3, 1, dotProduct(upminus, eye));
  lookat.set(0, 2, -lx);
  lookat.set(1, 2, -ly);
  lookat.set(2, 2, -lz);
  lookat.set(3, 2, dotProduct(l, eye));  
  
  // this = this * lookat
  compose(lookat);

}
  

  //
  // Projection transformations
  //

  public void perspective(float left, float right,
			  float bottom, float top,
			  float near, float far) {

    // set up perspective matrix
    Matrix3D p = new Matrix3D();
    p.set(0, 0, 2*near/(right-left));
    p.set(2, 0, -(right+left)/(right-left));
    p.set(1, 1, 2*near/(bottom-top));
    p.set(2, 1, -(bottom+top)/(bottom-top));
    p.set(2, 2, (far+near)/(far-near));
    p.set(3, 2, -2*far*near/(far-near));
    p.set(2, 3, 1);
    p.set(3, 3, 0);

    // this = this * p
    compose(p);

  }


  public void orthographic(float left, float right,
			   float bottom, float top,
			   float near, float far) {

    // set up orthographic matrix
    Matrix3D o = new Matrix3D();
    o.set(0, 0, 2/(right-left));
    o.set(3, 0, -(right+left)/(right-left));
    o.set(1, 1, 2/(bottom-top));
    o.set(3, 1, -(bottom+top)/(bottom-top));
    o.set(2, 2, 2/(far-near));
    o.set(3, 2, -(far+near)/(far-near));

    // this = this * o
    compose(o);

  }


  public void print() {

    System.out.println();
    for (int i=0; i<4; i++) {
      for (int j=0; j<4; j++)
	System.out.print(matrix[i][j] + "  ");
      System.out.println();
    }
    System.out.println();

  }



}