hps-java/src/main/java/org/lcsim/hps/users/phansson
diff -u -r1.1 -r1.2
--- MPAlignmentParameters.java 23 May 2012 01:04:12 -0000 1.1
+++ MPAlignmentParameters.java 24 May 2012 01:30:31 -0000 1.2
@@ -11,10 +11,7 @@
import java.io.FileWriter;
import java.io.IOException;
import java.io.PrintWriter;
-import java.util.HashSet;
-import java.util.List;
-import java.util.Map;
-import java.util.Set;
+import java.util.*;
import java.util.logging.Level;
import java.util.logging.Logger;
import org.lcsim.detector.IDetectorElement;
@@ -58,8 +55,9 @@
private int _nlc = 5; //the five track parameters
private int _ngl = 1; //delta(u) and delta(gamma) for each plane
+ private GlobalParameters _glp;
private BasicMatrix _dfdq;
- private BasicMatrix _dfdp;
+ //private BasicMatrix _dfdp;
private HelicalTrackFit _trk;
private double[] _resid = new double[3];
private double[] _error = new double[3];
@@ -70,6 +68,7 @@
private HPSTransformations _detToTrk;
boolean _DEBUG = true;
double smax = 1e3;
+
public MPAlignmentParameters(String outfile) {
_detToTrk = new HPSTransformations();
@@ -80,7 +79,10 @@
} catch (IOException ex) {
Logger.getLogger(RunAlignment.class.getName()).log(Level.SEVERE, null, ex);
}
-
+
+ _glp = new GlobalParameters();
+ if(_DEBUG) _glp.print();
+
}
public void PrintResidualsAndDerivatives(Track track) {
@@ -151,6 +153,10 @@
}
private void CalculateGlobalDerivatives(HelicalTrackStrip strip) {
+ if(_DEBUG) {
+ System.out.println("Calculate global derivaties for this sensor hit in layer " + strip.layer());
+ }
+
//1st index = alignment parameter (only u so far)
//2nd index = residual coordinate (on du so far)
//Naming scheme:
@@ -167,30 +173,91 @@
// [Layer]
// 1-10
+
+ //go through all the global parameters defined and calculate dfdp
+
+ //The global derivatives is done in the tracking frame
+ //THis means that after they are calculated I need to transform them
+ //into the sensor frame because that is where the residuals are calculated
+
+ double[][] dfdpTRACK = new double[3][1];
+ dfdpTRACK[0][0] = 0; //df/dx
+ dfdpTRACK[1][0] = 0; //df/dy
+ dfdpTRACK[2][0] = 0; //df/dz
+
+ //Find out which global derivative this is
//align only layer 3 on top side!
int side = 10000;
- //if( strip.origin().z()>0) side = 10000;
- //else side = 20000;
- int l = 3;//strip.layer();
+ if( strip.origin().z()>0) side = 10000;
+ else side = 20000;
+ int l = strip.layer();
int type = 1000;
int dir = 300;
int label = side + type + dir + l;
- double[][] dfdpLab = new double[3][1];
- dfdpLab[0][0] = 0; //df/dx
- dfdpLab[1][0] = 0; //df/dy
- dfdpLab[2][0] = 0; //df/dz
- if(strip.origin().z()>0 && strip.layer()==3) dfdpLab[2][0] = 1; //df/dz
-
+ //if(strip.origin().z()>0 && strip.layer()==3) {
- BasicMatrix _dfdpLab = FillMatrix(dfdpLab, 3, 1);
- Hep3Matrix trkToStrip = getTrackToStripRotation(strip);
- _dfdp = (BasicMatrix) MatrixOp.mult(trkToStrip, _dfdpLab);
- if (_DEBUG) {
- System.out.printf("dfdz = %5.5f %5.5f %5.5f GL%d\n", _dfdp.e(0, 0), _dfdp.e(1, 0), _dfdp.e(2, 0), label);
+ //Go through each gl parameter and see if it has non-zero contribution
+ for (GlobalParameter gp : _glp.getList()) {
+ //if(gp.active()==false) continue;
+
+ dfdpTRACK[0][0] = 0;
+ dfdpTRACK[1][0] = 0;
+ dfdpTRACK[2][0] = 0;
- }
- _globalLabel[0] = label; //GetIdentifier(strip);
-// _globalLabel[0] = GetIdentifierModule(strip);
+ if(gp.getSide()==side) {
+
+
+ // correct side
+ if(gp.getLayer()==strip.layer()) {
+ //correct layer
+
+ //Calcuate dfdp derivatives depending on type of global parameter
+ if(gp.getType()==1000) {
+ //translation
+ if(gp.getDirection()==100) {
+ //x - tracking frame --> beamline direction
+ dfdpTRACK[0][0] = 1; //df/dx
+ dfdpTRACK[1][0] = 0; //df/dy
+ dfdpTRACK[2][0] = 0; //df/dz
+ }
+ else if(gp.getDirection()==200) {
+ //y - tracking frame --> almost unmeasured coordinate direction
+ dfdpTRACK[0][0] = 0; //df/dx
+ dfdpTRACK[1][0] = 1; //df/dy
+ dfdpTRACK[2][0] = 0; //df/dz
+
+ }
+ else {
+ //z - tracking frame --> almost measured coordinate direction
+ dfdpTRACK[0][0] = 0; //df/dx
+ dfdpTRACK[1][0] = 0; //df/dy
+ dfdpTRACK[2][0] = 1; //df/dz
+
+ }
+
+ }//type
+ }//layer
+
+ }//side
+
+ //Put it into a matrix to be able to transform easily
+ BasicMatrix _dfdpTRACK = FillMatrix(dfdpTRACK, 3, 1);
+ //Get transformation matrix from tracking frame to sensor frame where the residuals are calculated
+ Hep3Matrix trkToStrip = getTrackToStripRotation(strip);
+ //Transform derivatives to sensor frame!
+ BasicMatrix dfdp = (BasicMatrix) MatrixOp.mult(trkToStrip, _dfdpTRACK);
+ //Add it to the global parameter object
+ gp.setDfDp(dfdp);
+ if (_DEBUG) {
+ System.out.printf("dfdz = %5.5f %5.5f %5.5f GL%d name: %s\n", gp.dfdp(0), gp.dfdp(1), gp.dfdp(2), gp.getLabel(),gp.getName());
+
+ }
+ } //gps
+
+
+
+
+
}
@@ -280,7 +347,7 @@
return res;
}
-
+/*
public void AddTarget(double beamdy, double beamdz) {
double[][] dfdp = new double[3][1];
double d0 = _trk.dca();
@@ -302,7 +369,7 @@
dfdp[0][0] = 1;
dfdp[1][0] = 0;
dfdp[2][0] = 0;
- _dfdp = FillMatrix(dfdp, 3, 1);
+ _dfdp.add(FillMatrix(dfdp, 3, 1));
_globalLabel[0] = 666;
pWriter.printf("%4d\n", 666);
pWriter.printf("%5.5e %5.5e %5.5e\n", _resid[0], _resid[1], _resid[2]);
@@ -311,11 +378,11 @@
pWriter.printf("%5.5e %5.5e -1.0\n", mydzdq[i], mydydq[i]);
}
for (int j = 0; j < _ngl; j++) {
- pWriter.printf("%5.5e %5.5e %5.5e %5d\n", _dfdp.e(0, j), _dfdp.e(1, j), _dfdp.e(2, j), _globalLabel[j]);
+ pWriter.printf("%5.5e %5.5e %5.5e %5d\n", _dfdp.get(_dfdp.size()-1).e(0, j), _dfdp.get(_dfdp.size()-1).e(1, j), _dfdp.get(_dfdp.size()-1).e(2, j), _globalLabel.get(_globalLabel.size()-1)[j]);
}
}
-
+*/
private void PrintStripResiduals(HelicalTrackStrip strip) {
if (_DEBUG) {
System.out.printf("Strip Layer = %4d\n", strip.layer());
@@ -326,10 +393,11 @@
for (int i = 0; i < _nlc; i++) {
System.out.printf("%s %5.5e %5.5e %5.5e\n", q[i], _dfdq.e(0, i), _dfdq.e(1, i), _dfdq.e(2, i));
}
- String[] p = {"u-displacement"};
+ //String[] p = {"u-displacement"};
System.out.println("global parameter derivatives");
- for (int j = 0; j < _ngl; j++) {
- System.out.printf("%s %5.5e %5.5e %5.5e %5d\n", p[j], _dfdp.e(0, j), _dfdp.e(1, j), _dfdp.e(2, j), _globalLabel[j]);
+ for (GlobalParameter gl : _glp.getList()) {
+ System.out.printf("%s %5.5e %5.5e %5.5e %5d %10s\n", "", gl.dfdp(0), gl.dfdp(1), gl.dfdp(2), gl.getLabel(),gl.getName());
+ //System.out.printf("%s %5.5e %5.5e %5.5e %5d\n", p[j], _dfdp.e(0, j), _dfdp.e(1, j), _dfdp.e(2, j), _globalLabel[j]);
}
}
@@ -339,14 +407,21 @@
for (int i = 0; i < _nlc; i++) {
pWriter.printf("%5.5e %5.5e %5.5e\n", _dfdq.e(0, i), _dfdq.e(1, i), _dfdq.e(2, i));
}
- for (int j = 0; j < _ngl; j++) {
- pWriter.printf("%5.5e %5.5e %5.5e %5d\n", _dfdp.e(0, j), _dfdp.e(1, j), _dfdp.e(2, j), _globalLabel[j]);
+ for (GlobalParameter gl: _glp.getList()) {
+ if(gl.active()){
+ pWriter.printf("%5.5e %5.5e %5.5e %5d\n", gl.dfdp(0), gl.dfdp(1), gl.dfdp(2), gl.getLabel());
+ }
}
}
private Hep3Matrix getTrackToStripRotation(HelicalTrackStrip strip) {
+ //This function transforms the hit to the sensor coordinates
+
+ //Transform from JLab frame to sensor frame (done through the RawTrackerHit)
ITransform3D detToStrip = GetGlobalToLocal(strip);
+ //Get rotation matrix
Hep3Matrix detToStripMatrix = (BasicHep3Matrix) detToStrip.getRotation().getRotationMatrix();
+ //Transformation between the JLAB and tracking coordinate systems
Hep3Matrix detToTrackMatrix = (BasicHep3Matrix) _detToTrk.getMatrix();
if (_DEBUG) {
@@ -362,6 +437,7 @@
}
private ITransform3D GetGlobalToLocal(HelicalTrackStrip strip) {
+ //Transform from JLab frame (RawTrackerHit) to sensor frame (i.e. u,v,w)
RawTrackerHit rth = (RawTrackerHit) strip.rawhits().get(0);
IDetectorElement ide = rth.getDetectorElement();
SiSensor sensor = ide.findDescendants(SiSensor.class).get(0);