1    /*
2    Copyright 2004 by Ralph Hartley
3    This software is licenced under the terms of the
4    Gnu Public Licence
5    */
6    
7    public class Matrix2 implements java.io.Serializable {
8      
9      private static final long serialVersionUID = Version.getSUID();
10   
11     double[][] A = new double[2][2];
12   
13     public Matrix2(double a00, double a01, double a10,double a11) {
14       A[0][0] = a00;
15       A[0][1] = a01;
16       A[1][0] = a10;
17       A[1][1] = a11;
18     }
19   
20     public static Matrix2 I = new Matrix2(1,0,0,1);
21     public static Matrix2 Iminus = new Matrix2(-1,0,0,-1);
22   
23     public Matrix2(Matrix2 X,Matrix2 Y) {
24       for (int i=0;i<2;i++)
25         for (int j=0;j<2;j++)
26           for (int k=0;k<2;k++)
27           A[i][j] += X.A[i][k] * Y.A[k][j];
28     }
29   
30     public void add(Matrix2 term) {
31       for (int i=0;i<2;i++)
32         for (int j=0;j<2;j++)
33           A[i][j] += term.A[i][j];
34     }
35   
36     public void mul(double fact) {
37       for (int i=0;i<2;i++)
38         for (int j=0;j<2;j++)
39           A[i][j] *= fact;
40     }
41   
42     public Matrix2 transpose() {
43       return new Matrix2(A[0][0],A[1][0],A[0][1],A[1][1]);
44     }
45   
46     private double orthoStep() {
47       Matrix2 X = new Matrix2(transpose(),this);
48       X.add(Iminus);
49   
50       Matrix2 C = I;
51   
52       X.mul(-0.5);
53       C.add(X);
54   
55       Matrix2 Y = new Matrix2(X,X);
56       Y.mul(3.0/2.0);
57       C.add(Y);
58   
59       Y = new Matrix2(X,Y);
60       Y.mul(5.0/3.0);
61       C.add(Y);
62   
63       Y = new Matrix2(this,C);
64       A = Y.A;
65   
66       return(X.A[0][0]*X.A[1][1] - X.A[0][1]*X.A[1][0]);
67     }
68   
69     static double MAXERR = 0.0001;
70     static int MAXINTS = 5;
71   
72     public boolean orthogonalize() {
73       double err = 1;
74       for (int i=0;i<MAXINTS && err>MAXERR;i++)
75         err = orthoStep();
76       return(err<MAXERR);
77     }
78   
79     public double getScale() {
80       return(Math.sqrt(A[0][0]*A[1][1]-A[0][1]*A[1][0]));
81     }
82   
83     public double getAngle() {
84       return(Math.atan2(A[0][1],A[0][0]));
85     }
86   }
87