hps-java/src/main/java/org/lcsim/HPSVertexing
diff -u -r1.4 -r1.5
--- VertexFitter.java 12 May 2010 21:51:59 -0000 1.4
+++ VertexFitter.java 14 Nov 2010 00:49:11 -0000 1.5
@@ -73,7 +73,7 @@
// Extract the error matrix and it's inverse for the measured variables
Matrix GInv = FillCovariance(sltlist);
- Matrix GG = MatrixInvert(GInv);
+ Matrix GG = MatrixOp.inverse(GInv);
// The inverse error matrix for the unmeasured variables is all 0's
Matrix GStar = new BasicMatrix(3, 3);
@@ -85,8 +85,8 @@
for (int iter = 0; iter < _maxIterations; iter++) {
// Update the fitted and unmeasured variables
- Matrix xx = MatrixAdd(mm, cc);
- Matrix xstar = MatrixAdd(mstar, cstar);
+ Matrix xx = MatrixOp.add(mm, cc);
+ Matrix xstar = MatrixOp.add(mstar, cstar);
// Calculate the errors in the constraint equations
Matrix FF = CalculateFF(ntrks, xx, xstar);
@@ -101,24 +101,24 @@
// Calculate the constant terms in the system of linear equations
MutableMatrix YY = new BasicMatrix(AA.getNRows(), 1);
- Matrix RR = MatrixSub(MatrixAdd(MatrixOp.mult(MatrixTranspose(BStar), cstar), MatrixOp.mult(MatrixTranspose(BB), cc)), FF);
- setSubMatrix(YY, RR, 0, 0);
+ Matrix RR = MatrixOp.sub(MatrixOp.add(MatrixOp.mult(MatrixOp.transposed(BStar), cstar), MatrixOp.mult(MatrixOp.transposed(BB), cc)), FF);
+ MatrixOp.setSubMatrix(YY, RR, 0, 0);
// Solve the system of linear equations
- Matrix AAInv = MatrixInvert(AA);
+ Matrix AAInv = MatrixOp.inverse(AA);
Matrix XX = MatrixOp.mult(AAInv, YY);
// Get the Lagrange multipliers
- Matrix alpha = getSubMatrix(XX, 0, 0, FF.getNRows(), 1);
+ Matrix alpha = MatrixOp.getSubMatrix(XX, 0, 0, FF.getNRows(), 1);
if (iter > 0 && Converged(FF)) {
// Fit converged - extract the fit quantities
// Calculate the fit chi square
- double chisq = MatrixOp.mult(MatrixTranspose(cc), MatrixOp.mult(GG, cc)).e(0, 0);
+ double chisq = MatrixOp.mult(MatrixOp.transposed(cc), MatrixOp.mult(GG, cc)).e(0, 0);
// Extract the covariance matrix for the vertex position from the matrix AAInv
- Matrix vtxcov = getSubMatrix(AAInv, FF.getNRows(), FF.getNRows(), GStar.getNRows(), GStar.getNColumns());
+ Matrix vtxcov = MatrixOp.getSubMatrix(AAInv, FF.getNRows(), FF.getNRows(), GStar.getNRows(), GStar.getNColumns());
// Find the fitted track directions and save them in a map
Map<StraightLineTrack, Hep3Vector> dirmap = new HashMap<StraightLineTrack, Hep3Vector>();
@@ -128,13 +128,13 @@
}
// Create a new VertexFit and return with a success flag
- _vfit = new VertexFit(getVector(xstar), new SymmetricMatrix(vtxcov), chisq, dirmap);
+ _vfit = new VertexFit(MatrixOp.as3Vector(xstar), new SymmetricMatrix(vtxcov), chisq, dirmap);
return true;
}
// Update the change in the fitted and unmeasured variables
- cc = MatrixMultiply(-1., MatrixOp.mult(GInv, MatrixOp.mult(BB, alpha)));
- cstar = getSubMatrix(XX, FF.getNRows(), 0, GStar.getNRows(), 1);
+ cc = MatrixOp.mult(-1., MatrixOp.mult(GInv, MatrixOp.mult(BB, alpha)));
+ cstar = MatrixOp.getSubMatrix(XX, FF.getNRows(), 0, GStar.getNRows(), 1);
}
return false;
@@ -315,13 +315,13 @@
*/
private Matrix CalculateAA(int ntrks, Matrix BB, Matrix BStar, Matrix GInv, Matrix GStar) {
- Matrix HH = MatrixOp.mult(MatrixTranspose(BB), MatrixOp.mult(GInv, BB));
- Matrix BStarT = MatrixTranspose(BStar);
+ Matrix HH = MatrixOp.mult(MatrixOp.transposed(BB), MatrixOp.mult(GInv, BB));
+ Matrix BStarT = MatrixOp.transposed(BStar);
MutableMatrix AA = new BasicMatrix(HH.getNRows() + BStar.getNRows(), HH.getNColumns() + BStarT.getNColumns());
- setSubMatrix(AA, MatrixMultiply(-1., HH), 0, 0);
- setSubMatrix(AA, BStarT, 0, HH.getNColumns());
- setSubMatrix(AA, BStar, HH.getNRows(), 0);
- setSubMatrix(AA, GStar, HH.getNRows(), HH.getNColumns());
+ MatrixOp.setSubMatrix(AA, MatrixOp.mult(-1., HH), 0, 0);
+ MatrixOp.setSubMatrix(AA, BStarT, 0, HH.getNColumns());
+ MatrixOp.setSubMatrix(AA, BStar, HH.getNRows(), 0);
+ MatrixOp.setSubMatrix(AA, GStar, HH.getNRows(), HH.getNColumns());
return AA;
}
@@ -371,165 +371,4 @@
}
return true;
}
-
- /**
- * Returns the transpose of the matrix (inexplicably not handled by
- * the matrix package for non-square matrices).
- *
- * @param m matrix to be transposed
- * @return transposed matrix
- */
- private Matrix MatrixTranspose(Matrix m) {
-
- MutableMatrix mt = new BasicMatrix(m.getNColumns(), m.getNRows());
- for (int i = 0; i < m.getNRows(); i++) {
- for (int j = 0; j < m.getNColumns(); j++) {
- mt.setElement(j, i, m.e(i, j));
- }
- }
- return mt;
- }
-
- /**
- * Add two matrices together. Matrices must have the same dimensions.
- *
- * @param a first matrix
- * @param b second matrix
- * @return sum of the two matrices
- */
- private Matrix MatrixAdd(Matrix a, Matrix b) {
-
- if (a.getNRows() != b.getNRows() ||
- a.getNColumns() != b.getNColumns()) {
- throw new RuntimeException("Trying to add two matrices with different dimensions");
- }
- MutableMatrix c = new BasicMatrix(a.getNRows(), a.getNColumns());
- for (int i = 0; i < a.getNRows(); i++) {
- for (int j = 0; j < a.getNColumns(); j++) {
- c.setElement(i, j, a.e(i, j) + b.e(i, j));
- }
- }
- return c;
- }
-
- /**
- * Subtract matrix b from matrix a. Matrices must have the same dimensions.
- *
- * @param a starting matrix
- * @param b matrix to be subtracted
- * @return difference (a - b)
- */
- private Matrix MatrixSub(Matrix a, Matrix b) {
-
- if (a.getNRows() != b.getNRows() ||
- a.getNColumns() != b.getNColumns()) {
- throw new RuntimeException("Trying to add two matrices with different dimensions");
- }
- MutableMatrix c = new BasicMatrix(a.getNRows(), a.getNColumns());
- for (int i = 0; i < a.getNRows(); i++) {
- for (int j = 0; j < a.getNColumns(); j++) {
- c.setElement(i, j, a.e(i, j) - b.e(i, j));
- }
- }
- return c;
- }
-
- /**
- * Multiply a matrix by a scaler constant.
- *
- * @param c constant that will mutliply the matrix
- * @param m matrix to be scaled
- * @return scaled matrix (c * m)
- */
- private Matrix MatrixMultiply(double c, Matrix m) {
-
- MutableMatrix b = new BasicMatrix(m.getNRows(), m.getNColumns());
- for (int i = 0; i < m.getNRows(); i++) {
- for (int j = 0; j < m.getNColumns(); j++) {
- b.setElement(i, j, c * m.e(i, j));
- }
- }
- return b;
- }
-
- /**
- * Convenience method to invert a matrix in one step.
- *
- * @param m matrix to be inverted
- * @return inverted matrix
- */
- private Matrix MatrixInvert(Matrix m) {
-
- MutableMatrix minv = new BasicMatrix(m.getNRows(), m.getNColumns());
- MatrixOp.inverse(m, minv);
- return minv;
- }
-
- /**
- * Fill in part of a matrix with the contents of another matrix.
- *
- * @param mat matrix to be modified
- * @param sub submatrix to be inserted
- * @param row row where insertion is to take place
- * @param col column where insertion is to take place
- */
- private void setSubMatrix(MutableMatrix mat, Matrix sub, int row, int col) {
-
- // First check that the submatrix fits in the matrix
- if (row < 0 ||
- col < 0 ||
- sub.getNRows() + row > mat.getNRows() ||
- sub.getNColumns() + col > mat.getNColumns()) {
- throw new RuntimeException("Invalid attempt to insert a submatrix into a matrix");
- }
-
- // Loop over rows & columns to insert matrix
- for (int i = 0; i < sub.getNRows(); i++) {
- for (int j = 0; j < sub.getNColumns(); j++) {
- mat.setElement(i + row, j + col, sub.e(i, j));
- }
- }
- }
-
- /**
- * Extract part of a matrix.
- *
- * @param m matrix containing the desired submatrix
- * @param row row where the submatrix is located
- * @param col column where the submatrix is located
- * @param nrow number of rows in the submatrix
- * @param ncol number of columns in the submatrix
- * @return submatrix
- */
- private Matrix getSubMatrix(Matrix m, int row, int col, int nrow, int ncol) {
-
- // First check that the submatrix is a valid size
- if (row < 0 ||
- row + nrow > m.getNRows() ||
- col < 0 ||
- col + ncol > m.getNColumns()) {
- throw new RuntimeException("Invalid attempt to get a submatrix from a matrix");
- }
-
- // Loop over rows & columns and get submatrix
- MutableMatrix sm = new BasicMatrix(nrow, ncol);
- for (int i = 0; i < nrow; i++) {
- for (int j = 0; j < ncol; j++) {
- sm.setElement(i, j, m.e(i + row, j + col));
- }
- }
- return sm;
- }
-
- /**
- * Convenience method to turn a 3 row column matrix into a vector.
- * @param m column matrix
- * @return vector containing contents of column matrix
- */
- private Hep3Vector getVector(Matrix m) {
-
- // First check that we have a 3x1 matrix
- if (m.getNRows() != 3 || m.getNColumns() != 1) throw new RuntimeException("Invalid attempt to form a vector from a matrix");
- return new BasicHep3Vector(m.e(0, 0), m.e(1, 0), m.e(2, 0));
- }
}
\ No newline at end of file
