5 added files
java/sandbox/tracking/src/main/java/org/hps/svt/alignment
--- java/sandbox/tracking/src/main/java/org/hps/svt/alignment/AlignmentParameters.java (rev 0)
+++ java/sandbox/tracking/src/main/java/org/hps/svt/alignment/AlignmentParameters.java 2014-03-26 02:23:40 UTC (rev 362)
@@ -0,0 +1,561 @@
+package org.hps.svt.alignment;
+
+import hep.physics.matrix.BasicMatrix;
+import hep.physics.matrix.MatrixOp;
+import hep.physics.vec.BasicHep3Matrix;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Matrix;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+
+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.logging.Level;
+import java.util.logging.Logger;
+
+import org.hps.recon.tracking.HPSTransformations;
+import org.lcsim.detector.IDetectorElement;
+import org.lcsim.detector.ITransform3D;
+import org.lcsim.detector.tracker.silicon.ChargeCarrier;
+import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.detector.tracker.silicon.SiSensorElectrodes;
+import org.lcsim.event.RawTrackerHit;
+import org.lcsim.event.Track;
+import org.lcsim.event.TrackerHit;
+import org.lcsim.fit.helicaltrack.HelicalTrackCross;
+import org.lcsim.fit.helicaltrack.HelicalTrackFit;
+import org.lcsim.fit.helicaltrack.HelicalTrackHit;
+import org.lcsim.fit.helicaltrack.HelicalTrackStrip;
+import org.lcsim.fit.helicaltrack.HelixUtils;
+import org.lcsim.fit.helicaltrack.MultipleScatter;
+import org.lcsim.fit.helicaltrack.TrackDirection;
+import org.lcsim.recon.tracking.seedtracker.SeedCandidate;
+import org.lcsim.recon.tracking.seedtracker.SeedTrack;
+
+/**
+ * Class to calculate and print the residuals and derivatives
+ * of the alignment parameters...used as input for MillePede
+ * Notation follows the MillePede manual:
+ * http://www.desy.de/~blobel/Mptwo.pdf
+ *
+ * the track is measured in the HelicalTrackFit frame
+ * and residuals are in the sensor frame (u,v,w)
+ *
+ * ordering of track parameters is
+ * double d0 = _trk.dca();
+ * double z0 = _trk.z0();
+ * double slope = _trk.slope();
+ * double phi0 = _trk.phi0();
+ * double R = _trk.R();
+ *
+ * @author mgraham
+ */
+public class AlignmentParameters {
+
+ private int _nlc = 5; //the five track parameters
+ private int _ngl = 1; //delta(u) and delta(gamma) for each plane
+ private BasicMatrix _dfdq;
+ private BasicMatrix _dfdp;
+ private HelicalTrackFit _trk;
+ private double[] _resid = new double[3];
+ private double[] _error = new double[3];
+ private int[] _globalLabel = new int[1];
+ FileWriter fWriter;
+ PrintWriter pWriter;
+ Set<SiSensor> _process_sensors = new HashSet<SiSensor>();
+ boolean _DEBUG = false;
+ double smax = 1e3;
+
+ public AlignmentParameters(String outfile) {
+ try {
+//open things up
+ fWriter = new FileWriter(outfile);
+ pWriter = new PrintWriter(fWriter);
+ } catch (IOException ex) {
+ Logger.getLogger(RunAlignment.class.getName()).log(Level.SEVERE, null, ex);
+ }
+
+ }
+
+ public void PrintResidualsAndDerivatives(Track track) {
+ SeedTrack st = (SeedTrack) track;
+ SeedCandidate seed = st.getSeedCandidate();
+ Map<HelicalTrackHit, MultipleScatter> msmap = seed.getMSMap();
+ _trk = seed.getHelix();
+ List<TrackerHit> hitsOnTrack = track.getTrackerHits();
+ for (TrackerHit hit : hitsOnTrack) {
+ HelicalTrackHit htc = (HelicalTrackHit) hit;
+ double msdrphi = msmap.get(htc).drphi();
+ double msdz = msmap.get(htc).dz();
+ double sHit = _trk.PathMap().get(htc);
+ HelicalTrackCross cross = (HelicalTrackCross) htc;
+ List<HelicalTrackStrip> clusterlist = cross.getStrips();
+ TrackDirection trkdir = HelixUtils.CalculateTrackDirection(_trk, sHit);
+ cross.setTrackDirection(trkdir, _trk.covariance());
+ for (HelicalTrackStrip cl : clusterlist) {
+ CalculateLocalDerivatives(cl);
+ CalculateGlobalDerivatives(cl);
+ CalculateResidual(cl, msdrphi, msdz);
+// CalculateResidual(cl, 0,0);
+ PrintStripResiduals(cl);
+ }
+ }
+ AddTarget(0.1, 0.02);
+ }
+
+ private void CalculateLocalDerivatives(HelicalTrackStrip strip) {
+ //get track parameters.
+ double d0 = _trk.dca();
+ double z0 = _trk.z0();
+ double slope = _trk.slope();
+ double phi0 = _trk.phi0();
+ double R = _trk.R();
+//strip origin is defined in the tracking coordinate system (x=beamline)
+ double xint = strip.origin().x();
+ double s = HelixUtils.PathToXPlane(_trk, xint, smax, _nlc).get(0);
+ double phi = s / R - phi0;
+ double[][] dfdq = new double[3][5];
+ //dx/dq
+ //these are wrong for X, but for now it doesn't matter
+ dfdq[0][0] = Math.sin(phi0);
+ dfdq[0][1] = 0;
+ dfdq[0][2] = 0;
+ dfdq[0][3] = d0 * Math.cos(phi0) + R * Math.sin(phi0) - s * Math.cos(phi0);
+ dfdq[0][4] = (phi - phi0) * Math.cos(phi0);
+ double[] mydydq = dydq(R, d0, phi0, xint, s);
+ double[] mydzdq = dzdq(R, d0, phi0, xint, slope, s);
+ for (int i = 0; i < 5; i++) {
+ dfdq[1][i] = mydydq[i];
+ dfdq[2][i] = mydzdq[i];
+ }
+
+ BasicMatrix dfdqGlobal = FillMatrix(dfdq, 3, 5);
+ Hep3Matrix trkToStrip = getTrackToStripRotation(strip);
+ _dfdq = (BasicMatrix) MatrixOp.mult(trkToStrip, dfdqGlobal);
+
+ if (_DEBUG) {
+ double[] trackpars = {d0, z0, slope, phi0, R, s, xint};
+ System.out.println("Strip Origin: ");
+ System.out.println(strip.origin());
+ System.out.println("trkToStrip Rotation:");
+ System.out.println(trkToStrip.toString());
+ printDerivatives(trackpars, dfdq);
+ }
+ }
+
+ private void CalculateGlobalDerivatives(HelicalTrackStrip strip) {
+ //1st index = alignment parameter (only u so far)
+ //2nd index = residual coordinate (on du so far)
+
+ double[][] dfdpLab = new double[3][1];
+ dfdpLab[0][0] = 0; //df/dx
+ dfdpLab[1][0] = 0; //df/dy
+ dfdpLab[2][0] = 1; //df/dz
+ 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\n", _dfdp.e(0, 0), _dfdp.e(1, 0), _dfdp.e(2, 0));
+ }
+ _globalLabel[0] = GetIdentifier(strip);
+// _globalLabel[0] = GetIdentifierModule(strip);
+
+ }
+
+ private void CalculateResidual(HelicalTrackStrip strip, double msdrdphi, double msdz) {
+
+ Hep3Vector u = strip.u();
+ Hep3Vector v = strip.v();
+ Hep3Vector w = strip.w();
+ Hep3Vector corigin = strip.origin();
+ double phi0 = _trk.phi0();
+ double R = _trk.R();
+ double xint = strip.origin().x();
+ double s = HelixUtils.PathToXPlane(_trk, xint, smax, _nlc).get(0);
+ double phi = s / R - phi0;
+ Hep3Vector trkpos = HelixUtils.PointOnHelix(_trk, s);
+
+ //System.out.println("trkpos = "+trkpos.toString());
+ //System.out.println("origin = "+corigin.toString());
+
+ Hep3Vector mserr = new BasicHep3Vector(msdrdphi * Math.sin(phi), msdrdphi * Math.sin(phi), msdz);
+ Hep3Vector vdiffTrk = VecOp.sub(trkpos, corigin);
+ Hep3Matrix trkToStrip = getTrackToStripRotation(strip);
+ Hep3Vector vdiff = VecOp.mult(trkToStrip, vdiffTrk);
+ double umc = vdiff.x();
+ double vmc = vdiff.y();
+ double wmc = vdiff.z();
+ double umeas = strip.umeas();
+ double uError = strip.du();
+ double msuError = VecOp.dot(mserr, u);
+ double vmeas = 0;
+ double vError = (strip.vmax() - strip.vmin()) / Math.sqrt(12);
+ double wmeas = 0;
+ double wError = 0.001;
+ //System.out.println("strip error="+uError+"; ms error ="+msuError);
+ _resid[0] = umeas - umc;
+ _error[0] = Math.sqrt(uError * uError + msuError * msuError);
+ _resid[1] = vmeas - vmc;
+ _error[1] = vError;
+ _resid[2] = wmeas - wmc;
+ _error[2] = wError;
+ if (_DEBUG) {
+ System.out.println("Strip Origin: ");
+ System.out.println(corigin.toString());
+ System.out.println("Position on Track:");
+ System.out.println(trkpos.toString());
+ System.out.println("vdiff :");
+ System.out.println(vdiff.toString());
+ System.out.println("u :");
+ System.out.println(u.toString());
+ System.out.println("umeas = " + umeas + "; umc = " + umc);
+ System.out.println("udiff = " + _resid[0] + " +/- " + _error[0]);
+
+ }
+
+ }
+
+ public double[] getResidual(Track track, int layer) {
+ double[] res = new double[7];
+ SeedTrack st = (SeedTrack) track;
+ SeedCandidate seed = st.getSeedCandidate();
+ Map<HelicalTrackHit, MultipleScatter> msmap = seed.getMSMap();
+ _trk = seed.getHelix();
+ List<TrackerHit> hitsOnTrack = track.getTrackerHits();
+ for (TrackerHit hit : hitsOnTrack) {
+ HelicalTrackHit htc = (HelicalTrackHit) hit;
+ double sHit = _trk.PathMap().get(htc);
+ HelicalTrackCross cross = (HelicalTrackCross) htc;
+ List<HelicalTrackStrip> clusterlist = cross.getStrips();
+ TrackDirection trkdir = HelixUtils.CalculateTrackDirection(_trk, sHit);
+ double msdrphi = msmap.get(htc).drphi();
+ double msdz = msmap.get(htc).dz();
+ cross.setTrackDirection(trkdir, _trk.covariance());
+ for (HelicalTrackStrip cl : clusterlist) {
+ if (cl.layer() == layer) {
+ CalculateResidual(cl, msdrphi, msdz);
+ res[0] = _resid[0];
+ res[1] = _resid[1];
+ res[2] = _resid[2];
+ res[3] = _error[0];
+ res[4] = _error[1];
+ res[5] = _error[2];
+ res[6] = layer;
+ if(hit.getPosition()[2]<0)res[6]=layer+10;
+ }
+ }
+ }
+ return res;
+
+ }
+
+ public void AddTarget(double beamdy, double beamdz) {
+ double[][] dfdp = new double[3][1];
+ double d0 = _trk.dca();
+ double z0 = _trk.z0();
+ double slope = _trk.slope();
+ double phi0 = _trk.phi0();
+ double R = _trk.R();
+ double xint = 0; //target
+ double s = HelixUtils.PathToXPlane(_trk, xint, smax, _nlc).get(0);
+ Hep3Vector ptAtTarget = HelixUtils.PointOnHelix(_trk, s);
+ double[] mydydq = dydq(R, d0, phi0, xint, s);
+ double[] mydzdq = dzdq(R, d0, phi0, xint, slope, s);
+ _resid[0] = ptAtTarget.z();
+ _resid[1] = ptAtTarget.y();
+ _resid[2] = ptAtTarget.x();
+ _error[0] = beamdz;
+ _error[1] = beamdy;
+ _error[2] = 666;
+ dfdp[0][0] = 1;
+ dfdp[1][0] = 0;
+ dfdp[2][0] = 0;
+ _dfdp = 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]);
+ pWriter.printf("%5.5e %5.5e %5.5e\n", _error[0], _error[1], _error[2]);
+ for (int i = 0; i < _nlc; i++) {
+ 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]);
+ }
+
+ }
+
+ private void PrintStripResiduals(HelicalTrackStrip strip) {
+ if (_DEBUG) {
+ System.out.printf("Strip Layer = %4d\n", strip.layer());
+ System.out.printf("Residuals (u,v,w) : %5.5e %5.5e %5.5e\n", _resid[0], _resid[1], _resid[2]);
+ System.out.printf("Errors (u,v,w) : %5.5e %5.5e %5.5e\n", _error[0], _error[1], _error[2]);
+ String[] q = {"d0", "z0", "slope", "phi0", "R"};
+ System.out.println("track parameter derivatives");
+ 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"};
+ 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]);
+ }
+
+ }
+ pWriter.printf("%4d\n", strip.layer());
+ pWriter.printf("%5.5e %5.5e %5.5e\n", _resid[0], _resid[1], _resid[2]);
+ pWriter.printf("%5.5e %5.5e %5.5e\n", _error[0], _error[1], _error[2]);
+ 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]);
+ }
+ }
+
+ private Hep3Matrix getTrackToStripRotation(HelicalTrackStrip strip) {
+ ITransform3D detToStrip = GetGlobalToLocal(strip);
+ Hep3Matrix detToStripMatrix = (BasicHep3Matrix) detToStrip.getRotation().getRotationMatrix();
+ Hep3Matrix detToTrackMatrix = (BasicHep3Matrix) HPSTransformations.getMatrix();
+
+ if (_DEBUG) {
+ System.out.println("gblToLoc translation:");
+ System.out.println(detToStrip.getTranslation().toString());
+ System.out.println("gblToLoc Rotation:");
+ System.out.println(detToStrip.getRotation().toString());
+ System.out.println("detToTrack Rotation:");
+ System.out.println(detToTrackMatrix.toString());
+ }
+
+ return (Hep3Matrix) VecOp.mult(detToStripMatrix, VecOp.inverse(detToTrackMatrix));
+ }
+
+ private ITransform3D GetGlobalToLocal(HelicalTrackStrip strip) {
+ RawTrackerHit rth = (RawTrackerHit) strip.rawhits().get(0);
+ IDetectorElement ide = rth.getDetectorElement();
+ SiSensor sensor = ide.findDescendants(SiSensor.class).get(0);
+ SiSensorElectrodes electrodes = sensor.getReadoutElectrodes(ChargeCarrier.HOLE);
+ return electrodes.getGlobalToLocal();
+ }
+
+ private int GetIdentifier(HelicalTrackStrip strip) {
+ RawTrackerHit rth = (RawTrackerHit) strip.rawhits().get(0);
+ IDetectorElement ide = rth.getDetectorElement();
+ SiSensor sensor = ide.findDescendants(SiSensor.class).get(0);
+ // return rth.getIdentifierFieldValue(sensor.getName());
+ return sensor.getSensorID(); //individual sensor positions
+// int sid=sensor.getSensorID();
+// int global=1;
+// if(sid>10)global=2;
+// return global; //return top/bottom plates
+ }
+
+ private int GetIdentifierModule(HelicalTrackStrip strip) {
+ RawTrackerHit rth = (RawTrackerHit) strip.rawhits().get(0);
+ IDetectorElement ide = rth.getDetectorElement();
+ SiSensor sensor = ide.findDescendants(SiSensor.class).get(0);
+ // return rth.getIdentifierFieldValue(sensor.getName());
+// return sensor.getSensorID(); //individual sensor positions
+ int sid = sensor.getSensorID();
+ int gid = -1;
+ switch (sid) {
+ case 1:
+ gid = 1; break;
+ case 2:
+ gid = 1;break;
+ case 3:
+ gid = 2;break;
+ case 4:
+ gid = 2;break;
+ case 5:
+ gid = 3;break;
+ case 6:
+ gid = 3;break;
+ case 7:
+ gid = 4;break;
+ case 8:
+ gid = 4;break;
+ case 9:
+ gid = 5;break;
+ case 10:
+ gid = 5;break;
+ case 11:
+ gid = 11;break;
+ case 12:
+ gid = 11;break;
+ case 13:
+ gid = 12;break;
+ case 14:
+ gid = 12;break;
+ case 15:
+ gid = 13;break;
+ case 16:
+ gid = 13;break;
+ case 17:
+ gid = 14;break;
+ case 18:
+ gid = 14;break;
+ case 19:
+ gid = 15;break;
+ case 20:
+ gid = 15;break;
+ }
+
+ return gid; //return top/bottom plates
+ }
+
+ private BasicMatrix FillMatrix(double[][] array, int nrow, int ncol) {
+ BasicMatrix retMat = new BasicMatrix(nrow, ncol);
+ for (int i = 0; i < nrow; i++) {
+ for (int j = 0; j < ncol; j++) {
+ retMat.setElement(i, j, array[i][j]);
+ }
+ }
+ return retMat;
+ }
+
+ public void closeFile() throws IOException {
+ pWriter.close();
+ fWriter.close();
+ }
+
+ private double dsdR(double R, double d0, double phi0, double xint) {
+ double sqrtTerm = Sqrt(R * R - Math.pow(((d0 - R) * Sin(phi0) + xint), 2));
+
+ double rsign = Math.signum(R);
+ double dsdr = (1 / sqrtTerm) * ((-rsign * xint) + (-rsign) * d0 * Sin(phi0)
+ + ArcTan(R * Cos(phi0), (-R) * Sin(phi0))
+ * sqrtTerm
+ - ArcTan(rsign * sqrtTerm, xint + (d0 - R) * Sin(phi0))
+ * sqrtTerm);
+
+
+ if (_DEBUG)
+ System.out.println("xint = " + xint + "; dsdr = " + dsdr);
+ return dsdr;
+
+ }
+
+ private double dsdphi(double R, double d0, double phi0, double xint) {
+ double sqrtTerm = Sqrt(R * R - Math.pow(((d0 - R) * Sin(phi0) + xint), 2));
+ double rsign = Math.signum(R);
+ double dsdphi = R * (sqrtTerm + rsign * d0 * Cos(phi0) - rsign * R * Cos(phi0)) / sqrtTerm;
+ if (_DEBUG)
+ System.out.println("xint = " + xint + "; dsdphi = " + dsdphi);
+ return dsdphi;
+ }
+
+ private double dsdd0(double R, double d0, double phi0, double xint) {
+ double sqrtTerm = Sqrt(R * R - Math.pow(((d0 - R) * Sin(phi0) + xint), 2));
+ double rsign = Math.signum(R);
+ double dsdd0 = rsign * (R * Sin(phi0)) / sqrtTerm;
+ if (_DEBUG)
+ System.out.println("xint = " + xint + "; dsdd0 = " + dsdd0);
+ return dsdd0;
+ }
+
+ private double[] dydq(double R, double d0, double phi0, double xint, double s) {
+ double[] dy = new double[5];
+// dy[0] = Cos(phi0) + Cot(phi0 - s / R) * Csc(phi0 - s / R) * dsdd0(R, d0, phi0, xint);
+ dy[0] = Cos(phi0) - Sec(phi0 - s / R) * Tan(phi0 - s / R) * dsdd0(R, d0, phi0, xint);
+ dy[1] = 0;
+ dy[2] = 0;
+// dy[3] = (-(d0 - R)) * Sin(phi0) - R * Cot(phi0 - s / R) * Csc(phi0 - s / R) * (1 - dsdphi(R, d0, phi0, xint) / R);
+ dy[3] = (-(d0 - R)) * Sin(phi0) + Sec(phi0 - s / R) * Tan(phi0 - s / R) * (R - dsdphi(R, d0, phi0, xint));
+ // dy[4] = -Cos(phi0) + Csc(phi0 - s / R) - R * Cot(phi0 - s / R) * Csc(phi0 - s / R) * (s / (R * R) - dsdR(R, d0, phi0, xint) / R);
+ dy[4] = -Cos(phi0) + Sec(phi0 - s / R) + (1 / R) * Sec(phi0 - s / R) * Tan(phi0 - s / R) * (s - R * dsdR(R, d0, phi0, xint));
+ return dy;
+ }
+
+ private double[] dzdq(double R, double d0, double phi0, double xint, double slope, double s) {
+ double[] dz = new double[5];
+ dz[0] = slope * dsdd0(R, d0, phi0, xint);
+ dz[1] = 1;
+ dz[2] = s;
+ dz[3] = slope * dsdphi(R, d0, phi0, xint);
+ dz[4] = slope * dsdR(R, d0, phi0, xint);
+ return dz;
+ }
+
+ private double Csc(double val) {
+ return 1 / Math.sin(val);
+ }
+
+ private double Cot(double val) {
+ return 1 / Math.tan(val);
+ }
+
+ private double Sec(double val) {
+ return 1 / Math.cos(val);
+ }
+
+ private double Sin(double val) {
+ return Math.sin(val);
+ }
+
+ private double Cos(double val) {
+ return Math.cos(val);
+ }
+
+ private double Tan(double val) {
+ return Math.tan(val);
+ }
+
+ private double ArcTan(double val1, double val2) {
+ return Math.atan2(val1, val2);
+ }
+
+ private double Sign(double val) {
+ return Math.signum(val);
+ }
+
+ private double Sqrt(double val) {
+ return Math.sqrt(val);
+ }
+
+ private void printDerivatives(double[] trackpars, double[][] dfdq) {
+ System.out.println("======================================================");
+ System.out.println("s xint");
+ System.out.printf("%5.5f %5.5f\n", trackpars[5], trackpars[6]);
+ System.out.println(" d0 z0 slope phi0 R");
+ System.out.printf("Values %5.5f %5.5f %5.5f %5.5f %5.5f\n", trackpars[0], trackpars[1], trackpars[2], trackpars[3], trackpars[4]);
+ System.out.printf("dzdq ");
+ for (int i = 0; i < 5; i++) {
+ System.out.printf("%5.3e ", dfdq[2][i]);
+ }
+ System.out.println();
+ System.out.printf("dudq ");
+ for (int i = 0; i < _nlc; i++) {
+ System.out.printf("%5.3e ", _dfdq.e(0, i));
+ }
+ System.out.println();
+ System.out.println();
+ System.out.printf("dydq ");
+ for (int i = 0; i < 5; i++) {
+ System.out.printf("%5.3e ", dfdq[1][i]);
+ }
+ System.out.println();
+ System.out.printf("dvdq ");
+ for (int i = 0; i < _nlc; i++) {
+ System.out.printf("%5.3e ", _dfdq.e(1, i));
+ }
+ System.out.println();
+ System.out.println();
+ System.out.printf("dxdq ");
+ for (int i = 0; i < 5; i++) {
+ System.out.printf("%5.3e ", dfdq[0][i]);
+ }
+ System.out.println();
+ System.out.printf("dwdq ");
+ for (int i = 0; i < _nlc; i++) {
+ System.out.printf("%5.3e ", _dfdq.e(2, i));
+ }
+ System.out.println();
+ // System.out.println( _trk.xc()+ "; "+_trk.yc());
+// System.out.println( _trk.x0()+ "; "+_trk.y0());
+ }
+}
java/sandbox/tracking/src/main/java/org/hps/svt/alignment
--- java/sandbox/tracking/src/main/java/org/hps/svt/alignment/BuildCompact.java (rev 0)
+++ java/sandbox/tracking/src/main/java/org/hps/svt/alignment/BuildCompact.java 2014-03-26 02:23:40 UTC (rev 362)
@@ -0,0 +1,645 @@
+package org.hps.svt.alignment;
+
+/**
+ * Class building a new compact.xml detector based on MillepedeII input corrections
+ * @author phansson
+ * created on 1/15/2014
+ */
+
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+
+import java.io.BufferedReader;
+import java.io.File;
+import java.io.FileInputStream;
+import java.io.FileNotFoundException;
+import java.io.FileWriter;
+import java.io.IOException;
+import java.io.InputStreamReader;
+import java.util.ArrayList;
+import java.util.HashMap;
+import java.util.List;
+import java.util.Map;
+import java.util.Map.Entry;
+
+import org.apache.commons.cli.CommandLine;
+import org.apache.commons.cli.CommandLineParser;
+import org.apache.commons.cli.HelpFormatter;
+import org.apache.commons.cli.Option;
+import org.apache.commons.cli.Options;
+import org.apache.commons.cli.ParseException;
+import org.apache.commons.cli.PosixParser;
+import org.apache.commons.lang.StringUtils;
+import org.hps.conditions.deprecated.HPSSVTSensorSetup;
+import org.hps.conditions.deprecated.SvtUtils;
+import org.hps.recon.tracking.HPSTransformations;
+import org.jdom.Document;
+import org.jdom.Element;
+import org.jdom.JDOMException;
+import org.jdom.input.SAXBuilder;
+import org.jdom.output.XMLOutputter;
+import org.lcsim.conditions.ConditionsManager;
+import org.lcsim.conditions.ConditionsManager.ConditionsNotFoundException;
+import org.lcsim.detector.IDetectorElement;
+import org.lcsim.detector.ITransform3D;
+import org.lcsim.detector.tracker.silicon.ChargeCarrier;
+import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.geometry.Detector;
+import org.lcsim.geometry.GeometryReader;
+import org.lcsim.util.xml.ElementFactory.ElementCreationException;
+
+
+public class BuildCompact {
+
+ private static int runNumber = -1; //1351;
+ private static String detectorName = ""; //"HPS-TestRun-v7";
+ private static ConditionsManager conditionsManager = null;
+ private static double corrScaleFactor = -1.;
+
+ private static Options createCmdLineOpts() {
+ Options options = new Options();
+ options.addOption(new Option("c",true,"The path to the compact xml file."));
+ options.addOption(new Option("o",true,"The name of the new compact xml file."));
+ options.addOption(new Option("d",true,"Detector name."));
+ options.addOption(new Option("r",true,"Run number."));
+
+ return options;
+ }
+
+ private static void printHelpAndExit(Options options) {
+ HelpFormatter help = new HelpFormatter();
+ help.printHelp(" ", options);
+ System.exit(1);
+ }
+
+ //private static buildDetector()
+
+
+
+
+
+
+ private static class MilleParameterSet {
+ private IDetectorElement _det = null;
+ List<MilleParameter> params = new ArrayList<MilleParameter>();
+ public MilleParameterSet(IDetectorElement d) {
+ setDetector(d);
+ }
+ public void setDetector(IDetectorElement d) {
+ _det = d;
+ }
+ public IDetectorElement getDetector() {
+ return _det;
+ }
+ public void add(MilleParameter par) {
+ params.add(par);
+ }
+ public Hep3Vector getLocalTranslation() {
+ Map<String,Double> m = new HashMap<String,Double>();
+ for(MilleParameter p : params) {
+ if (p.getType() == 1) {
+ if (p.getDim() == 1) m.put("u", p.getValue());
+ else if(p.getDim() == 2) m.put("v", p.getValue());
+ else m.put("w", p.getValue());
+ }
+ }
+ if(m.size() != 3) {
+ System.out.println("bad trans!!");
+ System.exit(1);
+ }
+ return new BasicHep3Vector(m.get("u"),m.get("v"),m.get("w"));
+
+ }
+ public Hep3Vector getLocalRotation() {
+ Map<String,Double> m = new HashMap<String,Double>();
+ for(MilleParameter p : params) {
+ if (p.getType() == 2) {
+ if (p.getDim() == 1) m.put("alpha",p.getValue());
+ else if(p.getDim() == 2) m.put("beta", p.getValue());
+ else m.put("gamma", p.getValue());
+ }
+ }
+ if(m.size() != 3) {
+ System.out.println("bad rot!!");
+ System.exit(1);
+ }
+ return new BasicHep3Vector(m.get("alpha"),m.get("beta"),m.get("gamma"));
+ }
+ public Hep3Vector getGlobalTranslation() {
+ ITransform3D localToGlobal = getLocalToGlobal();
+ return localToGlobal.getRotation().rotated(getLocalTranslation());
+ }
+ public double getGlobalTranslation(int d) {
+ return getGlobalTranslation().v()[d-1];
+ }
+ public Hep3Vector getGlobalRotation() {
+ ITransform3D localToGlobal = getLocalToGlobal();
+ return localToGlobal.getRotation().rotated(getLocalRotation());
+ }
+ public double getGlobalRotation(int d) {
+ return getGlobalRotation().v()[d-1];
+ }
+ public ITransform3D getLocalToGlobal() {
+ ITransform3D localToGlobal = ( (SiSensor) _det).getReadoutElectrodes(ChargeCarrier.HOLE).getLocalToGlobal();
+ return localToGlobal;
+ }
+
+ }
+
+ private static class MilleParameter {
+ private int id;
+ private double value;
+ private double presigma;
+ private static final Map<Integer,String> dMap;
+ private static final Map<Integer,String> tMap;
+ private static final Map<Integer,String> hMap;
+ static {
+ dMap = new HashMap<Integer,String>();
+ dMap.put(1, "x");dMap.put(2, "y"); dMap.put(3, "z");
+ tMap = new HashMap<Integer,String>();
+ tMap.put(1, "");tMap.put(2, "r");
+ hMap = new HashMap<Integer,String>();
+ hMap.put(1, "t");hMap.put(2, "b");
+ }
+
+ public MilleParameter(String line) {
+ String[] vals = StringUtils.split(line);// line.split("\\s+");
+ if(vals.length <3) {
+ System.out.println("this line is ill-formatted (" + vals.length + ")");
+ System.out.println(line);
+ System.exit(1);
+ }
+ try {
+ //for(String v : vals) System.out.println("\"" + v + "\"");
+ setId(Integer.parseInt(vals[0]));
+ setValue( corrScaleFactor * Double.parseDouble(vals[1]) );
+ setPresigma(Double.parseDouble(vals[2]));
+
+ } catch (NumberFormatException e) {
+ System.out.println(vals[0] + " " + vals[1] + " " + vals[2]);
+ throw new RuntimeException("problem parsing string ", e);
+ }
+ }
+
+ public String getXMLName() {
+ String d = dMap.get(getDim());
+ String t = tMap.get(getType());
+ String h = hMap.get(getHalf());
+ int s = getSensor();
+ return String.format("%s%s%d%s_align", t,d,s,h);
+
+ }
+
+ private int getDim() {
+ int h = (int) (getHalf() * 1e4);
+ int t = (int) (getType() * 1e3);
+ return (int) Math.floor((id- h -t)/1e2);
+ }
+
+ private int getSensor() {
+ int h = (int) (getHalf() * 1e4);
+ int t = (int) (getType() * 1e3);
+ int d = (int) (getDim() * 1e2);
+ return (id - h - t -d);
+ }
+
+ public int getType() {
+ int h = (int) (getHalf() * 1e4);
+ return (int) Math.floor((id -h)/1e3);
+ }
+
+ private int getHalf() {
+ return (int)Math.floor(id/1e4);
+ }
+
+ public int getId() {
+ return id;
+ }
+
+ public void setId(int id) {
+ this.id = id;
+ }
+
+ public double getValue() {
+ return value;
+ }
+
+ public void setValue(double value) {
+ this.value = value;
+ }
+
+ public double getPresigma() {
+ return presigma;
+ }
+
+ public void setPresigma(double presigma) {
+ this.presigma = presigma;
+ }
+ }
+
+
+ public static void main(String[] args) {
+
+ // Setup command line input
+ Options options = createCmdLineOpts();
+ if (args.length == 0) {
+ printHelpAndExit(options);
+ }
+
+ CommandLineParser parser = new PosixParser();
+ CommandLine cl = null;
+ try {
+ cl = parser.parse(options, args);
+ } catch (ParseException e) {
+ throw new RuntimeException("Problem parsing command line options.",e);
+ }
+
+ String compactFilename = null;
+ if(cl.hasOption("c")) {
+ compactFilename = cl.getOptionValue("c");
+ } else {
+ printHelpAndExit(options);
+ }
+
+ if(cl.hasOption("d")) {
+ detectorName = cl.getOptionValue("d");
+ } else {
+ printHelpAndExit(options);
+ }
+
+ if(cl.hasOption("r")) {
+ runNumber = Integer.parseInt(cl.getOptionValue("r"));
+ } else {
+ printHelpAndExit(options);
+ }
+
+ String compactFilenameNew = "compact_new.xml";
+ if(cl.hasOption("o")) {
+ compactFilenameNew = cl.getOptionValue("o");
+ }
+
+
+
+ File compactFile = new File(compactFilename);
+
+ // read XML
+ SAXBuilder builder = new SAXBuilder();
+ Document compact_document = null;
+ try {
+ compact_document = (Document) builder.build(compactFile);
+ } catch (JDOMException | IOException e1) {
+ throw new RuntimeException("problem with JDOM ", e1);
+ }
+
+ // read detector geometry
+ FileInputStream inCompact;
+ try {
+ inCompact = new FileInputStream(compactFile);
+ } catch (FileNotFoundException e) {
+ throw new RuntimeException("cannot open compact file",e);
+ }
+
+ GeometryReader reader = new GeometryReader();
+ Detector det;
+ try {
+ det = reader.read(inCompact);
+ } catch (IOException | JDOMException | ElementCreationException e) {
+ throw new RuntimeException("problem reading compact file",e);
+ }
+
+
+ // set conditions in order to be able to determine which sensors are where in the geometry
+ setConditions(detectorName,runNumber);
+
+
+ // setup sensors
+ HPSSVTSensorSetup sensorSetup = new HPSSVTSensorSetup();
+ sensorSetup.detectorChanged(det);
+
+ // Loop over all millepede input files and match parameters with detectors
+
+ List<MilleParameterSet> list_det = new ArrayList<MilleParameterSet>();
+
+ FileInputStream inMille = null;
+ BufferedReader brMille = null;
+ try {
+ for(String milleFilename : cl.getArgs()) {
+ inMille = new FileInputStream(milleFilename);
+ brMille = new BufferedReader(new InputStreamReader(inMille));
+ String line;
+ while((line = brMille.readLine()) != null) {
+ //System.out.printf("%s\n",line);
+ if(!line.contains("Parameter") && !line.contains("!")) {
+
+ MilleParameter par = new MilleParameter(line);
+ //System.out.println(par.getXMLName() + " " + par.getValue());
+
+ SiSensor detEle = getTrackerDetElement( det, par);
+ if(detEle == null) {
+ System.out.println("Couldn't find detector for param " + par.getId());
+ System.exit(1);
+ }
+ System.out.println("Found detector " + detEle.getName());
+ if(detEle.getClass().isInstance(SiSensor.class)) {
+ System.out.println("yeah");
+ }
+
+ // do we have it already?
+ MilleParameterSet useSet = null;
+ for(MilleParameterSet set : list_det) {
+ if(set.getDetector() == detEle) {
+ useSet = set;
+ }
+ }
+
+ if (useSet == null) {
+ useSet = new MilleParameterSet(detEle);
+ list_det.add(useSet);
+ }
+
+ //add the parameter
+ useSet.add(par);
+
+
+ }
+ }
+ brMille.close();
+ }
+ }
+ catch (IOException e) {
+ throw new RuntimeException("problem reading mille file",e);
+ }
+
+ for(MilleParameterSet set : list_det) {
+ System.out.println("Detector " + set.getDetector().getName());
+ List<MilleParameter> pars = set.params;
+ for(MilleParameter p : pars) {
+ System.out.println(p.getXMLName() + " " + p.getValue());
+ Element node = findXMLNode(compact_document,p.getXMLName());
+ double value = 0.0;
+ if(p.getType() == 1){
+ value = set.getGlobalTranslation(p.getDim());
+ } else if(p.getType() == 2){
+ value = set.getGlobalRotation(p.getDim());
+ } else {
+ System.out.println("This type is illdefnied " + p.getType());
+ System.exit(1);
+ }
+ node.setAttribute("value", String.format("%.6f",value));
+ }
+ Hep3Vector u = getMeasuredCoordinate( (SiSensor )set.getDetector());
+ System.out.println("u " + u.toString());
+ System.out.println("t (local) " + set.getLocalTranslation().toString());
+ System.out.println("t (global) " + set.getGlobalTranslation().toString());
+ System.out.println("r (local) " + set.getLocalRotation().toString());
+ System.out.println("r (global) " + set.getGlobalRotation().toString());
+
+
+ }
+
+ // Save new XML file
+
+ XMLOutputter xmlOutput = new XMLOutputter();
+ // display nice
+ //xmlOutput.setFormat(Format.getPrettyFormat());
+ try {
+ xmlOutput.output(compact_document, new FileWriter(compactFilenameNew));
+ } catch (IOException e) {
+ throw new RuntimeException("problem with xml output",e);
+ }
+
+
+
+
+
+ }
+
+ private static Element findXMLNode(Document document, String name) {
+ Element rootNode = document.getRootElement();
+ List list = rootNode.getChildren("define");
+ for(int i = 0; i < list.size(); ++i ) {
+ Element node = (Element) list.get(i);
+ List llist = node.getChildren("constant");
+ //System.out.println("length of list " + llist.size());
+ for(int ii = 0; ii < llist.size(); ++ii ) {
+ Element nnode = (Element) llist.get(ii);
+ //System.out.println("node name " + nnode.getAttributeValue("name") + " " + nnode.getAttributeValue("value") );
+ //if(nnode.getAttributeValue("name").contains(name)) {
+ if(nnode.getAttributeValue("name").compareTo(name) ==0 ) {
+ return nnode;
+ }
+ }
+ }
+
+ return null;
+ }
+
+ private static void setConditions(String detectorName, int run) {
+
+ try {
+ if(conditionsManager == null) {
+ conditionsManager = ConditionsManager.defaultInstance();
+ }
+ conditionsManager.setDetector(detectorName, run);
+
+ } catch (ConditionsNotFoundException e1) {
+ throw new RuntimeException("problem setting conditions",e1);
+ }
+
+ }
+
+ private static SiSensor getTrackerDetElement(Detector det, MilleParameter par) {
+ List<SiSensor> sensors = det.getSubdetector("Tracker").getDetectorElement().findDescendants(SiSensor.class);
+ //List<SiTrackerModule> modules = det.getDetectorElement().findDescendants(SiTrackerModule.class);
+ //System.out.printf("%d sensors\n",sensors.size());
+ for (SiSensor module: sensors) {
+ // Create DAQ Maps
+ if (!SvtUtils.getInstance().isSetup()) {
+ SvtUtils.getInstance().setup(det);
+ }
+ boolean isTop = SvtUtils.getInstance().isTopLayer(module);
+ int h = par.getHalf();
+ if ((isTop && h == 1) || (!isTop && h == 2)) {
+ int layer = SvtUtils.getInstance().getLayerNumber(module);
+ if (layer == par.getSensor()) {
+ //found match
+ return module;
+ }
+ }
+ }
+ return null;
+
+ }
+
+
+
+ private static class DetectorList<K> extends ArrayList<DetAlignConstants> {
+ //List<DetAlignConstants> _detectors = new ArrayList<DetAlignConstants>();
+ public DetectorList() {
+ }
+
+ public boolean contains(IDetectorElement detEle) {
+ return this.get(detEle) == null ? false : true;
+ }
+
+ public DetAlignConstants get(IDetectorElement detEle) {
+ for(DetAlignConstants d : this) {
+ if (d == detEle) {
+ return d;
+ }
+ }
+ return null;
+ }
+ public void print() {
+ System.out.println("==== " + this.size() + " detectors has alignment corrections ====");
+ for(DetAlignConstants det : this) {
+ det.print();
+ }
+ }
+
+ }
+
+
+
+ private static Hep3Vector getTrackingMeasuredCoordinate(SiSensor sensor)
+ {
+ // p-side unit vector
+ ITransform3D electrodes_to_global = sensor.getReadoutElectrodes(ChargeCarrier.HOLE).getLocalToGlobal();
+ Hep3Vector measuredCoordinate = sensor.getReadoutElectrodes(ChargeCarrier.HOLE).getMeasuredCoordinate(0);
+ measuredCoordinate = VecOp.mult(VecOp.mult(HPSTransformations.getMatrix(),electrodes_to_global.getRotation().getRotationMatrix()), measuredCoordinate);
+ return measuredCoordinate;
+ }
+
+ private static Hep3Vector getMeasuredCoordinate(SiSensor sensor)
+ {
+ // p-side unit vector
+ ITransform3D electrodes_to_global = sensor.getReadoutElectrodes(ChargeCarrier.HOLE).getLocalToGlobal();
+ Hep3Vector measuredCoordinate = sensor.getReadoutElectrodes(ChargeCarrier.HOLE).getMeasuredCoordinate(0);
+ measuredCoordinate = VecOp.mult(electrodes_to_global.getRotation().getRotationMatrix(), measuredCoordinate);
+ return measuredCoordinate;
+ }
+
+
+
+ private static class SiSensorDetAlignConstants extends DetAlignConstants {
+ public SiSensorDetAlignConstants(IDetectorElement det) {
+ super(det);
+ }
+ public void transform() {
+ ITransform3D localToGlobal = null;
+ if(_det.getClass().isInstance(SiSensor.class)) {
+ localToGlobal = ( (SiSensor) _det).getReadoutElectrodes(ChargeCarrier.HOLE).getLocalToGlobal();
+ }
+ //Translation
+ Hep3Vector t_local = _constants.getTranslationVector();
+ Hep3Vector t_global = localToGlobal.getRotation().rotated(t_local);
+ _constants.addGlobalTranslation(t_global);
+ //Rotation
+ Hep3Vector r_local = _constants.getRotationVector();
+ Hep3Vector r_global = localToGlobal.getRotation().rotated(r_local);
+ _constants.addGlobalRotation(r_global);
+ }
+
+ }
+
+ private static abstract class DetAlignConstants {
+ protected IDetectorElement _det = null;
+ protected AlignConstants<String,Double> _constants = new AlignConstants<String,Double>();
+ public DetAlignConstants(IDetectorElement det) {
+ _det = det;
+ }
+ public abstract void transform();
+ public void add(MilleParameter par) {
+ this._constants.add(par);
+ }
+ public void print() {
+ System.out.println(_det.getName());
+ for(Entry<String, Double> c : this._constants.entrySet()) {
+ System.out.println(c.getKey() + " " + c.getValue());
+ }
+ System.out.println("Local translation " + _constants.getTranslationVector().toString());
+ System.out.println("Global translation " + _constants.getTranslationVectorGlobal().toString());
+ System.out.println("Local rotation " + _constants.getRotationVector().toString());
+ System.out.println("Global rotation " + _constants.getRotationVectorGlobal().toString());
+
+
+ }
+
+
+ }
+
+
+
+ private static class AlignConstants<K,V> extends HashMap<String,Double> {
+ List<MilleParameter> _pars = new ArrayList<MilleParameter>();
+ public AlignConstants() {
+ super();
+ }
+ public void add(MilleParameter p) {
+ _pars.add(p);
+ if (p.getType() == 1) {
+ if (p.getDim() == 1) this.put("u", p.getValue());
+ else if(p.getDim() == 2) this.put("v", p.getValue());
+ else this.put("w", p.getValue());
+ }
+ else {
+ if (p.getDim() == 1) this.put("alpha", p.getValue());
+ else if(p.getDim() == 2) this.put("beta", p.getValue());
+ else this.put("gamma", p.getValue());
+ }
+ }
+ public void print() {
+ for(Entry<String,Double> e : this.entrySet()) {
+ System.out.println(e.getKey() + " " + e.getValue());
+ }
+ }
+ public Hep3Vector getTranslationVector() {
+ if(!this.containsKey("u") || !this.containsKey("v") || !this.containsKey("w")) {
+ System.out.println("missing pars for translation");
+ print();
+ System.exit(1);
+ }
+ return new BasicHep3Vector(this.get("u"),this.get("v"),this.get("w"));
+ }
+ public Hep3Vector getTranslationVectorGlobal() {
+ if(!this.containsKey("x") || !this.containsKey("y") || !this.containsKey("z")) {
+ System.out.println("missing pars for global translation");
+ print();
+ System.exit(1);
+ }
+ return new BasicHep3Vector(this.get("x"),this.get("y"),this.get("z"));
+ }
+ public Hep3Vector getRotationVector() {
+ if(!this.containsKey("alpha") || !this.containsKey("beta") || !this.containsKey("gamma")) {
+ System.out.println("missing pars for rotation");
+ print();
+ System.exit(1);
+ }
+ return new BasicHep3Vector(this.get("alpha"),this.get("beta"),this.get("gamma"));
+ }
+ public Hep3Vector getRotationVectorGlobal() {
+ if(!this.containsKey("rx") || !this.containsKey("ry") || !this.containsKey("rz")) {
+ System.out.println("missing pars for global rotation");
+ print();
+ System.exit(1);
+ }
+ return new BasicHep3Vector(this.get("rx"),this.get("ry"),this.get("rz"));
+ }
+ private void addGlobalTranslation(Hep3Vector t) {
+ this.put("x", t.x());
+ this.put("y", t.y());
+ this.put("z", t.z());
+ }
+ private void addGlobalRotation(Hep3Vector t) {
+ this.put("rx", t.x());
+ this.put("ry", t.y());
+ this.put("rz", t.z());
+ }
+
+
+ }
+
+
+
+}
java/sandbox/tracking/src/main/java/org/hps/svt/alignment
--- java/sandbox/tracking/src/main/java/org/hps/svt/alignment/HPSStrips.java (rev 0)
+++ java/sandbox/tracking/src/main/java/org/hps/svt/alignment/HPSStrips.java 2014-03-26 02:23:40 UTC (rev 362)
@@ -0,0 +1,483 @@
+package org.hps.svt.alignment;
+/*
+ * SiStrips.java
+ *
+ * Created on July 22, 2005, 4:07 PM
+ *
+ * To change this template, choose Tools | Options and locate the template under
+ * the Source Creation and Management node. Right-click the template and choose
+ * Open. You can then make changes to the template in the Source Editor.
+ */
+
+//import static org.lcsim.units.clhep.SystemOfUnits.*;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+
+import java.util.ArrayList;
+import java.util.HashSet;
+import java.util.List;
+import java.util.Set;
+import java.util.SortedMap;
+import java.util.TreeMap;
+
+import org.lcsim.detector.IDetectorElement;
+import org.lcsim.detector.ITransform3D;
+import org.lcsim.detector.Transform3D;
+import org.lcsim.detector.solids.GeomOp2D;
+import org.lcsim.detector.solids.GeomOp3D;
+import org.lcsim.detector.solids.Line3D;
+import org.lcsim.detector.solids.LineSegment3D;
+import org.lcsim.detector.solids.Point3D;
+import org.lcsim.detector.solids.Polygon3D;
+import org.lcsim.detector.tracker.silicon.ChargeCarrier;
+import org.lcsim.detector.tracker.silicon.ChargeDistribution;
+import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.detector.tracker.silicon.SiStrips;
+
+/**
+ *
+ * @author tknelson
+ */
+public class HPSStrips extends SiStrips
+{
+
+ // Fields
+
+ // Object definition
+ private ChargeCarrier _carrier; // charge carrier collected
+ private int _nstrips; // number of strips
+ private double _pitch; // sense pitch
+ private IDetectorElement _detector; // associated detector element
+ private ITransform3D _parent_to_local; // parent to local transform
+ private ITransform3D _local_to_global; // transformation to global coordinates
+ private ITransform3D _global_to_local; // transformation from global coordinates
+ private Polygon3D _geometry; // region in which strips are defined
+ private double _capacitance_intercept = 10.; // fixed capacitance independent of strip length
+ private double _capacitance_slope = 0.1; // capacitance per unit length of strip
+
+ // Cached for convenience
+ private double _strip_offset;
+
+
+
+ public HPSStrips(ChargeCarrier carrier, double pitch, IDetectorElement detector, ITransform3D parent_to_local,ITransform3D misalignment)
+ {
+
+// System.out.println("Plane of polygon in sensor coordinates has... ");
+// System.out.println(" normal: "+((SiSensor)detector).getBiasSurface(carrier).getNormal());
+// System.out.println(" distance: "+((SiSensor)detector).getBiasSurface(carrier).getDistance());
+
+ setCarrier(carrier);
+ setPitch(pitch);
+ setGeometry(((SiSensor)detector).getBiasSurface(carrier).transformed(parent_to_local));
+ setStripNumbering();
+ setDetectorElement(detector);
+ setParentToLocal(parent_to_local);
+ setGlobalToLocal(Transform3D.multiply(Transform3D.multiply(parent_to_local,detector.getGeometry().getGlobalToLocal()),misalignment));
+ setLocalToGlobal(getGlobalToLocal().inverse());
+ }
+
+ public HPSStrips(ChargeCarrier carrier, double pitch, int nstrips, IDetectorElement detector, ITransform3D parent_to_local,ITransform3D misalignment)
+ {
+ setCarrier(carrier);
+ setPitch(pitch);
+ setGeometry(((SiSensor)detector).getBiasSurface(carrier).transformed(parent_to_local));
+ setNStrips(nstrips);
+ setDetectorElement(detector);
+ setParentToLocal(parent_to_local);
+ setGlobalToLocal(Transform3D.multiply(Transform3D.multiply(parent_to_local,detector.getGeometry().getGlobalToLocal()),misalignment));
+ setLocalToGlobal(getGlobalToLocal().inverse());
+ }
+
+
+ // SiSensorElectrodes interface
+ //=============================
+
+ // Mechanical properties
+ public int getNAxes()
+ {
+ return 1;
+ }
+
+ public IDetectorElement getDetectorElement()
+ {
+ return _detector;
+ }
+
+ public ITransform3D getParentToLocal()
+ {
+ return _parent_to_local;
+ }
+
+ public ITransform3D getLocalToGlobal()
+ {
+ return _local_to_global;
+ }
+
+ public ITransform3D getGlobalToLocal()
+ {
+ return _global_to_local;
+ }
+
+ public Polygon3D getGeometry()
+ {
+ return _geometry;
+ }
+
+ public Hep3Vector getMeasuredCoordinate(int axis)
+ {
+ if (axis == 0) return new BasicHep3Vector(1.0,0.0,0.0);
+ else return null;
+ }
+
+ public Hep3Vector getUnmeasuredCoordinate(int axis)
+ {
+ if (axis == 0) return new BasicHep3Vector(0.0,1.0,0.0);
+ else return null;
+ }
+
+ public int getNeighborCell(int cell, int ncells_0, int ncells_1)
+ {
+ int neighbor_cell = cell + ncells_0;
+ if (isValidCell(neighbor_cell)) return neighbor_cell;
+ else return -1;
+ }
+
+ public Set<Integer> getNearestNeighborCells(int cell)
+ {
+ Set<Integer> neighbors = new HashSet<Integer>();
+ for (int ineigh = -1 ; ineigh <= 1; ineigh=ineigh+2)
+ {
+ int neighbor_cell = getNeighborCell(cell,ineigh,0);
+ if (isValidCell(neighbor_cell)) neighbors.add(neighbor_cell);
+ }
+ return neighbors;
+ }
+
+ public boolean isValidCell(int cell)
+ {
+ return (cell >= 0 && cell < getNCells());
+ }
+
+ public int getNCells()
+ {
+ return _nstrips;
+ }
+
+ public int getNCells(int axis)
+ {
+ if (axis == 0)
+ {
+ return _nstrips;
+ }
+ else return 1;
+ }
+
+ public double getPitch(int axis)
+ {
+ if (axis == 0)
+ {
+ return _pitch;
+ }
+ else return 0;
+ }
+ public int getCellID(Hep3Vector position)
+ {
+ return (int)Math.round((position.x()+_strip_offset)/_pitch);
+ }
+
+
+ public int getRowNumber(Hep3Vector position)
+ {
+ return 0;
+ }
+
+ public int getColumnNumber(Hep3Vector position)
+ {
+ return getCellID(position);
+ }
+
+ public int getCellID(int row_number, int column_number)
+ {
+ return column_number;
+ }
+
+ public int getRowNumber(int cell_id)
+ {
+ return 0;
+ }
+
+ public int getColumnNumber(int cell_id)
+ {
+ return cell_id;
+ }
+
+ public Hep3Vector getPositionInCell(Hep3Vector position)
+ {
+ return VecOp.sub(position,getCellPosition(getCellID(position)));
+ }
+
+ public Hep3Vector getCellPosition(int strip_number)
+ {
+ return new BasicHep3Vector(strip_number*_pitch-_strip_offset,0.0,0.0);
+ }
+
+ // Electrical properties
+
+ /**
+ * Capacitance intercept parameter. Units are pF.
+ *
+ * Capacitance is calculated as:
+ * C = capacitance_intercept + strip_length * capacitance slope
+ *
+ * @param capacitance_intercept
+ */
+ public void setCapacitanceIntercept(double capacitance_intercept) {
+ _capacitance_intercept = capacitance_intercept;
+ }
+
+ /**
+ * Capacitance per unit strip length. Units are pF / mm.
+ *
+ * @param capacitance_slope
+ */
+ public void setCapacitanceSlope(double capacitance_slope) {
+ _capacitance_slope = capacitance_slope;
+ }
+
+ public ChargeCarrier getChargeCarrier()
+ {
+ return _carrier;
+ }
+
+ /**
+ * Capacitance for a particular cell. Units are pF.
+ *
+ * @param cell_id
+ * @return
+ */
+ public double getCapacitance(int cell_id) // capacitance in pF
+ {
+ return _capacitance_intercept + _capacitance_slope*getStripLength(cell_id);
+ }
+
+ /**
+ * Nominal capacitance used for throwing random noise in the sensor.
+ * Calculated using middle strip. Units are pF.
+ *
+ * @return
+ */
+ public double getCapacitance() {
+ return getCapacitance(getNCells(0) / 2);
+ }
+
+ public SortedMap<Integer,Integer> computeElectrodeData(ChargeDistribution distribution)
+ {
+ SortedMap<Integer,Integer> electrode_data = new TreeMap<Integer,Integer>();
+
+ int base_strip = getCellID(distribution.getMean());
+
+ // put charge on strips in window 3-sigma strips on each side of base strip
+ int axis = 0;
+ int window_size = (int)Math.ceil(3.0*distribution.sigma1D(getMeasuredCoordinate(axis))/getPitch(axis));
+
+ double integral_lower = distribution.getNormalization();
+ double integral_upper = distribution.getNormalization();
+
+ for (int istrip = base_strip-window_size; istrip <= base_strip+window_size; istrip++)
+ {
+ double cell_edge_upper = getCellPosition(istrip).x() + getPitch(axis)/2.0;
+
+// System.out.println("cell_edge_upper: "+cell_edge_upper);
+
+ double integration_limit = cell_edge_upper; //cell_edge_upper-distribution.mean().x();
+
+// System.out.println("integration_limit: "+integration_limit);
+
+ integral_upper = distribution.upperIntegral1D(getMeasuredCoordinate(axis),integration_limit);
+
+// System.out.println("integral_upper: "+integral_upper);
+
+ if (integral_lower<integral_upper)
+ {
+ throw new RuntimeException("Error in integrating Gaussian charge distribution!");
+ }
+
+ int strip_charge = (int)Math.round(integral_lower-integral_upper);
+
+// System.out.println("strip_charge: "+strip_charge);
+
+ if (strip_charge != 0)
+ {
+ electrode_data.put(istrip,strip_charge);
+ }
+
+ integral_lower = integral_upper;
+ }
+
+ return electrode_data;
+
+ }
+
+ // Strip specific methods
+
+ // length of strip
+ public double getStripLength(int cell_id)
+ {
+// System.out.println("strip_length: "+getStrip(cell_id).getLength());
+ return getStrip(cell_id).getLength();
+ }
+
+ // center of strip
+ public Hep3Vector getStripCenter(int cell_id)
+ {
+ LineSegment3D strip = getStrip(cell_id);
+ return strip.getEndPoint(strip.getLength()/2);
+ }
+
+ // line segment for strip
+ public LineSegment3D getStrip(int cell_id)
+ {
+ Line3D strip_line = new Line3D(new Point3D(getCellPosition(cell_id)),getUnmeasuredCoordinate(0));
+
+// System.out.println("Number of strips: "+this._nstrips);
+// System.out.println("Strip offset: "+this._strip_offset);
+// System.out.println("Pitch: "+this._pitch);
+// System.out.println("cell_id: "+cell_id);
+// System.out.println("strip_line start point: "+strip_line.getStartPoint());
+// System.out.println("strip_line direction: "+strip_line.getDirection());
+
+ List<Point3D> intersections = new ArrayList<Point3D>();
+
+ // Get intersections between strip line and edges of electrode polygon
+ for (LineSegment3D edge : _geometry.getEdges())
+ {
+// System.out.println("edge start point: "+edge.getStartPoint());
+// System.out.println("edge end point: "+edge.getEndPoint());
+
+ if (GeomOp2D.intersects(strip_line,edge))
+ {
+ intersections.add(GeomOp3D.lineBetween(strip_line,edge).getStartPoint());
+ }
+ }
+
+ // Check for rare occurrence of duplicates (can happen at corners of polygon)
+ List<Point3D> strip_ends = new ArrayList<Point3D>(intersections);
+ if (intersections.size() > 2)
+ {
+ for (int ipoint1 = 0; ipoint1 < intersections.size(); ipoint1++)
+ {
+ Point3D point1 = intersections.get(ipoint1);
+ for (int ipoint2 = ipoint1+1; ipoint2 < intersections.size(); ipoint2++)
+ {
+ Point3D point2 = intersections.get(ipoint2);
+ if (GeomOp3D.intersects(point1,point2))
+ {
+ strip_ends.remove(point2);
+ if (strip_ends.size() == 2) break;
+ }
+ }
+ }
+ }
+
+ return new LineSegment3D(strip_ends.get(0),strip_ends.get(1));
+ }
+
+ // Private setters
+ //==================
+ public void setCarrier(ChargeCarrier carrier)
+ {
+ _carrier = carrier;
+ }
+
+ public void setGeometry(Polygon3D geometry)
+ {
+// System.out.println("Plane of polygon has... ");
+// System.out.println(" normal: "+geometry.getNormal());
+// System.out.println(" distance: "+geometry.getDistance());
+//
+// System.out.println("Working plane has... ");
+// System.out.println(" normal: "+GeomOp2D.PLANE.getNormal());
+// System.out.println(" distance: "+GeomOp2D.PLANE.getDistance());
+
+ if (GeomOp3D.equals(geometry.getPlane(),GeomOp2D.PLANE))
+ {
+ _geometry = geometry;
+ }
+ else
+ {
+ throw new RuntimeException("Electrode geometry must be defined in x-y plane!!");
+ }
+ }
+
+ private void setStripNumbering()
+ {
+ double xmin = Double.MAX_VALUE;
+ double xmax = Double.MIN_VALUE;
+ for (Point3D vertex : _geometry.getVertices())
+ {
+ xmin = Math.min(xmin,vertex.x());
+ xmax = Math.max(xmax,vertex.x());
+ }
+
+// System.out.println("xmin: " + xmin);
+// System.out.println("xmax: " + xmax);
+//
+//
+// System.out.println("# strips: " + (int)Math.ceil((xmax-xmin)/getPitch(0)) ) ;
+
+ setNStrips( (int)Math.ceil((xmax-xmin)/getPitch(0)) ) ;
+ }
+
+ private void setNStrips(int nstrips)
+ {
+ _nstrips = nstrips;
+ setStripOffset();
+// _strip_offset = (_nstrips-1)*_pitch/2.;
+ }
+
+ private void setStripOffset()
+ {
+ double xmin = Double.MAX_VALUE;
+ double xmax = Double.MIN_VALUE;
+ for (Point3D vertex : _geometry.getVertices())
+ {
+ xmin = Math.min(xmin,vertex.x());
+ xmax = Math.max(xmax,vertex.x());
+ }
+
+ double strips_center = (xmin+xmax)/2;
+
+ _strip_offset = ((_nstrips-1)*_pitch)/2 - strips_center;
+
+ }
+
+ private void setPitch(double pitch)
+ {
+ _pitch = pitch;
+ }
+
+ private void setDetectorElement(IDetectorElement detector)
+ {
+ _detector = detector;
+ }
+
+ private void setParentToLocal(ITransform3D parent_to_local)
+ {
+ _parent_to_local = parent_to_local;
+ }
+
+ private void setLocalToGlobal(ITransform3D local_to_global)
+ {
+ _local_to_global = local_to_global;
+ }
+
+ private void setGlobalToLocal(ITransform3D global_to_local)
+ {
+ _global_to_local = global_to_local;
+ }
+
+}
+
+
java/sandbox/tracking/src/main/java/org/hps/svt/alignment
--- java/sandbox/tracking/src/main/java/org/hps/svt/alignment/RunAlignment.java (rev 0)
+++ java/sandbox/tracking/src/main/java/org/hps/svt/alignment/RunAlignment.java 2014-03-26 02:23:40 UTC (rev 362)
@@ -0,0 +1,96 @@
+package org.hps.svt.alignment;
+
+import java.io.IOException;
+import java.util.List;
+import java.util.logging.Level;
+import java.util.logging.Logger;
+
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.Track;
+import org.lcsim.util.Driver;
+import org.lcsim.util.aida.AIDA;
+
+/**
+ *
+ * @author mgraham
+ */
+public class RunAlignment extends Driver {
+
+ private AIDA aida = AIDA.defaultInstance();
+ String[] detNames = {"Tracker"};
+ Integer _minLayers = 8;
+ Integer[] nlayers = {8};
+ int nevt = 0;
+ double[] beamsize = {0.001, 0.02, 0.02};
+ String _config = "1pt8";
+ AlignmentParameters ap;
+ int totalTracks=0;
+// flipSign is a kludge...
+// HelicalTrackFitter doesn't deal with B-fields in -ive Z correctly
+// so we set the B-field in +iveZ and flip signs of fitted tracks
+// note: this should be -1 for Test configurations and +1 for Full (v3.X and lower) configurations
+// this is set by the _config variable (detType in HeavyPhotonDriver)
+ int flipSign = 1;
+
+ public RunAlignment(int trackerLayers, int mintrkLayers, String config) {
+ nlayers[0] = trackerLayers;
+ _minLayers = mintrkLayers;
+ _config = config;
+ if (_config.contains("Test"))
+ flipSign = -1;
+ ap = new AlignmentParameters("/Users/mgraham/HPS/align.txt");
+
+ }
+
+ public void process(
+ EventHeader event) {
+
+
+ // Create a map between tracks and the associated MCParticle
+ List<Track> tracklist = event.get(Track.class, "MatchedTracks");
+// System.out.println("Number of Tracks = " + tracklist.size());
+ double duRange=0.1;
+ for (Track trk : tracklist) {
+ totalTracks++;
+ ap.PrintResidualsAndDerivatives(trk);
+
+ if(1==1){
+ aida.histogram1D("Track d0",50,-0.5,0.5).fill(trk.getTrackParameter(0));
+ aida.histogram1D("Track sin(phi0)",50,-0.5,0.5).fill(Math.sin(trk.getTrackParameter(1)));
+ aida.histogram1D("Track z0",50,-0.1,0.1).fill(Math.sin(trk.getTrackParameter(3)));
+ aida.histogram1D("Track chi^2",50,0,25).fill(trk.getChi2());
+ for (int i = 1; i < 11; i++) {
+ double[] res = ap.getResidual(trk, i);
+ int mylayer=(int)res[6];
+ if(mylayer<11){
+ aida.histogram1D("Track chi^2 Positive Side",50,0,25).fill(trk.getChi2());
+ }else{
+ aida.histogram1D("Track chi^2 Negative Side",50,0,25).fill(trk.getChi2());
+ }
+
+ aida.histogram1D("deltaU -- Layer " + mylayer,50,-duRange,duRange).fill(res[0]);
+ aida.histogram1D("deltaU Pull-- Layer " + mylayer,50,-3,3).fill(res[0]/res[3]);
+ if(i==3&&Math.sin(trk.getTrackParameter(1))>0){
+ aida.histogram1D("Positive phi0 deltaU -- Layer " + mylayer,50,-duRange,duRange).fill(res[0]);
+ aida.histogram1D("Positive phi0 deltaU Pull-- Layer " + mylayer,50,-3,3).fill(res[0]/res[3]);
+ }
+ if(i==3&&Math.sin(trk.getTrackParameter(1))<0){
+ aida.histogram1D("Negative phi0 deltaU -- Layer " + mylayer,50,-duRange,duRange).fill(res[0]);
+ aida.histogram1D("Negative phi0 deltaU Pull-- Layer " + mylayer,50,-3,3).fill(res[0]/res[3]);
+ }
+
+ }
+ }
+ }
+
+ }
+
+ public void endOfData() {
+ try {
+ System.out.println("Total Number of Tracks Found = "+totalTracks);
+ ap.closeFile();
+ } catch (IOException ex) {
+ Logger.getLogger(RunAlignment.class.getName()).log(Level.SEVERE, null, ex);
+ }
+ }
+}
java/sandbox/tracking/src/main/java/org/hps/svt/alignment
--- java/sandbox/tracking/src/main/java/org/hps/svt/alignment/SiTrackerSpectrometerSensorSetup.java (rev 0)
+++ java/sandbox/tracking/src/main/java/org/hps/svt/alignment/SiTrackerSpectrometerSensorSetup.java 2014-03-26 02:23:40 UTC (rev 362)
@@ -0,0 +1,265 @@
+package org.hps.svt.alignment;
+
+import hep.physics.matrix.BasicMatrix;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.VecOp;
+
+import java.util.List;
+
+import org.lcsim.detector.IDetectorElement;
+import org.lcsim.detector.IRotation3D;
+import org.lcsim.detector.ITranslation3D;
+import org.lcsim.detector.RotationPassiveXYZ;
+import org.lcsim.detector.Transform3D;
+import org.lcsim.detector.Translation3D;
+import org.lcsim.detector.solids.Polygon3D;
+import org.lcsim.detector.solids.Trd;
+import org.lcsim.detector.tracker.silicon.ChargeCarrier;
+import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.geometry.Detector;
+import org.lcsim.geometry.compact.Subdetector;
+import org.lcsim.geometry.subdetector.SiTrackerSpectrometer;
+import org.lcsim.util.Driver;
+
+public class SiTrackerSpectrometerSensorSetup extends Driver {
+
+ String subdetectorName;
+
+ public SiTrackerSpectrometerSensorSetup() {
+ }
+
+ public SiTrackerSpectrometerSensorSetup(String subdetectorName) {
+ this.subdetectorName = subdetectorName;
+ }
+
+ public void setSubdetectorName(String subdetectorName) {
+ this.subdetectorName = subdetectorName;
+ }
+
+ public void detectorChanged(Detector detector) {
+ if (subdetectorName == null) {
+ throw new RuntimeException("The subdetectorName was not set.");
+ }
+
+ Subdetector subdetector = detector.getSubdetector(subdetectorName);
+ if (subdetector instanceof SiTrackerSpectrometer) {
+ setupSensorDetectorElements(subdetector);
+ } else {
+ throw new RuntimeException("The subdetector " + subdetectorName + " is not an instance of SiTrackerSpectrometer.");
+ }
+ }
+
+ private void setupSensorDetectorElements(Subdetector subdet) {
+ System.out.println(this.getClass().getCanonicalName() + " - Setting up sensors for " + subdet.getName() + " ...");
+ int sensorId = 0;
+
+ for (IDetectorElement endcap : subdet.getDetectorElement().getChildren()) {
+ for (IDetectorElement layer : endcap.getChildren()) {
+ //int nwedges = layer.getChildren().size();
+ for (IDetectorElement wedge : layer.getChildren()) {
+ for (IDetectorElement module : wedge.getChildren()) {
+ List<SiSensor> sensors = module.findDescendants(SiSensor.class);
+
+ if (sensors.size() == 0) {
+ throw new RuntimeException("No sensors found in module.");
+ }
+
+ for (SiSensor sensor : sensors) {
+ Trd sensor_solid = (Trd) sensor.getGeometry().getLogicalVolume().getSolid();
+
+ Polygon3D n_side = sensor_solid.getFacesNormalTo(new BasicHep3Vector(0, -1, 0)).get(0);
+ Polygon3D p_side = sensor_solid.getFacesNormalTo(new BasicHep3Vector(0, 1, 0)).get(0);
+
+ // Bias the sensor
+// sensor.setBiasSurface(ChargeCarrier.ELECTRON, p_side);
+// sensor.setBiasSurface(ChargeCarrier.HOLE, n_side);
+
+ sensor.setBiasSurface(ChargeCarrier.HOLE, p_side);
+ sensor.setBiasSurface(ChargeCarrier.ELECTRON, n_side);
+
+// double strip_angle = Math.atan2(sensor_solid.getXHalfLength2() - sensor_solid.getXHalfLength1(), sensor_solid.getZHalfLength() * 2);
+ double strip_angle = 0.00;
+ ITranslation3D electrodes_position = new Translation3D(VecOp.mult(-p_side.getDistance(), new BasicHep3Vector(0, 0, 1))); // translate to outside of polygon
+ //ITranslation3D electrodes_position = new Translation3D(VecOp.mult(n_side.getDistance(), new BasicHep3Vector(0, 0, 1))); // translate to outside of polygon
+ // System.out.println("SensorID = " + sensorId + " " + electrodes_position.toString());
+ IRotation3D electrodes_rotation = new RotationPassiveXYZ(-Math.PI / 2, 0, strip_angle);
+ Transform3D electrodes_transform = new Transform3D(electrodes_position, electrodes_rotation);
+
+ // Free calculation of readout electrodes, sense electrodes determined thereon
+// SiStrips readout_electrodes = new SiStrips(ChargeCarrier.HOLE, 0.060, sensor, electrodes_transform);
+// SiStrips sense_electrodes = new SiStrips(ChargeCarrier.HOLE,0.030,(readout_electrodes.getNCells()*2-1),sensor,electrodes_transform);
+ ITranslation3D misalign_position;
+// System.out.println(layer.getName());
+// if (layer.getName().contains("3")) {
+// if (layer.getName().contains("3")&&layer.getName().contains("positive")) {
+/*
+ if (layer.getName().contains("positive")) {
+ System.out.println("Putting in a misalignment for layer "+layer.getName());
+ misalign_position = new Translation3D(0, 0.05, 0.0); // translate to outside of polygon
+ } else {
+ // misalign_position = new Translation3D(0, 0.0, 0.0);
+ misalign_position = new Translation3D(0, -0.05, 0.0);
+ }
+ */
+
+// if ((layer.getName().contains("3")||layer.getName().contains("4"))&&layer.getName().contains("positive")) {
+ // if (layer.getName().contains("positive")) {
+// System.out.println("Putting in a misalignment for layer "+layer.getName());
+// misalign_position = new Translation3D(0, 0.010, 0.0);
+// } else {
+// misalign_position = new Translation3D(0, 0.0, 0.0);
+ // misalign_position = new Translation3D(0, -0.05, 0.0);
+// }
+/*
+ if (layer.getName().contains("positive")) {
+ int gid = GetIdentifierModule(layer.getName());
+ misalign_position = new Translation3D(0, 0.0, 0.0);
+ if (gid == 1)
+ misalign_position = new Translation3D(0, -0.0144, 0.0);
+ if (gid == 2)
+ misalign_position = new Translation3D(0, 0.05-0.0297, 0.0);
+ if (gid == 3)
+ misalign_position = new Translation3D(0, -0.0253, 0.0);
+ if (gid == 4)
+ misalign_position = new Translation3D(0, -0.0346, 0.0);
+ if (gid == 5)
+ misalign_position = new Translation3D(0, -0.0433, 0.0);
+
+ } else {
+ misalign_position = new Translation3D(0, 0.0, 0.0);
+ }
+ */
+
+ int gid = GetIdentifierLayer(layer.getName());
+ misalign_position = new Translation3D(0, 0.0, 0.0);
+ /*
+ if (layer.getName().contains("positive")) {
+ if (gid == 1)
+ misalign_position = new Translation3D(0, -0.00144, 0.0);
+ if (gid == 2)
+ misalign_position = new Translation3D(0, 0.005, 0.0);
+ if (gid == 3)
+ misalign_position = new Translation3D(0, -0.00253, 0.0);
+ if (gid == 4)
+ misalign_position = new Translation3D(0, -0.00346, 0.0);
+ if (gid == 5)
+ misalign_position = new Translation3D(0, 0.00433, 0.0);
+ if (gid == 6)
+ misalign_position = new Translation3D(0, 0.0002, 0.0);
+ if (gid == 7)
+ misalign_position = new Translation3D(0, 0.002, 0.0);
+ if (gid == 8)
+ misalign_position = new Translation3D(0, -0.004, 0.0);
+ if (gid == 9)
+ misalign_position = new Translation3D(0, 0.006, 0.0);
+ if (gid == 10)
+ misalign_position = new Translation3D(0, -0.001, 0.0);
+ } else {
+ if (gid == 1)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 2)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 3)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 4)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 5)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 6)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 7)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 8)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 9)
+ misalign_position = new Translation3D(0, 0.00, 0.0);
+ if (gid == 10)
+ misalign_position = new Translation3D(0, 0.01, 0.0);
+ }
+
+*/
+ IRotation3D misalign_rotation = new RotationPassiveXYZ(0, 0, 0);
+ Transform3D misalign_transform = new Transform3D(misalign_position, misalign_rotation);
+
+ HPSStrips readout_electrodes = new HPSStrips(ChargeCarrier.HOLE, 0.060, sensor, electrodes_transform, misalign_transform);
+ HPSStrips sense_electrodes = new HPSStrips(ChargeCarrier.HOLE, 0.030, (readout_electrodes.getNCells() * 2 - 1), sensor, electrodes_transform, misalign_transform);
+
+
+// SiStrips readout_electrodes = new SiStrips(ChargeCarrier.ELECTRON, 0.060, sensor, electrodes_transform);
+// SiStrips sense_electrodes = new SiStrips(ChargeCarrier.ELECTRON, 0.030, (readout_electrodes.getNCells() * 2 - 1), sensor, electrodes_transform);
+
+ // SiSensorElectrodes sense_electrodes = new SiStrips(ChargeCarrier.HOLE, 0.060, sensor, electrodes_transform);
+
+//pristine conditions
+/*
+ readout_electrodes.setCapacitanceIntercept(0);
+ readout_electrodes.setCapacitanceSlope(0.12);
+ sense_electrodes.setCapacitanceIntercept(0);
+ sense_electrodes.setCapacitanceSlope(0.12);
+ */
+
+ readout_electrodes.setCapacitanceIntercept(0);
+ readout_electrodes.setCapacitanceSlope(0.16);
+ sense_electrodes.setCapacitanceIntercept(0);
+ sense_electrodes.setCapacitanceSlope(0.16);
+
+ sensor.setSenseElectrodes(sense_electrodes);
+ sensor.setReadoutElectrodes(readout_electrodes);
+//
+
+// double[][] transfer_efficiencies = {{1.0}};
+ double[][] transfer_efficiencies = {{0.986, 0.419}};
+ sensor.setTransferEfficiencies(ChargeCarrier.HOLE, new BasicMatrix(transfer_efficiencies));
+// sensor.setTransferEfficiencies(ChargeCarrier.ELECTRON, new BasicMatrix(transfer_efficiencies));
+ // here
+
+ sensor.setSensorID(++sensorId);
+ }
+ }
+ }
+ }
+ }
+ }
+
+ private int GetIdentifierModule(String mylayer) {
+ int gid = -1;
+ if (mylayer.contains("1") || mylayer.contains("2"))
+ gid = 1;
+ if (mylayer.contains("3") || mylayer.contains("4"))
+ gid = 2;
+ if (mylayer.contains("5") || mylayer.contains("6"))
+ gid = 3;
+ if (mylayer.contains("7") || mylayer.contains("8"))
+ gid = 4;
+ if (mylayer.contains("9") || mylayer.contains("10"))
+ gid = 5;
+
+ return gid; //return top/bottom plates
+ }
+
+ private int GetIdentifierLayer(String mylayer) {
+ int gid = -1;
+ if (mylayer.contains("1"))
+ gid = 1;
+ if (mylayer.contains("2"))
+ gid = 2;
+ if (mylayer.contains("3"))
+ gid = 3;
+ if (mylayer.contains("4"))
+ gid = 4;
+ if (mylayer.contains("5"))
+ gid = 5;
+ if (mylayer.contains("6"))
+ gid = 6;
+ if (mylayer.contains("7"))
+ gid = 7;
+ if (mylayer.contains("8"))
+ gid = 8;
+ if (mylayer.contains("9"))
+ gid = 9;
+ if (mylayer.contains("10"))
+ gid = 10;
+
+ return gid; //return top/bottom plates
+ }
+}
SVNspam 0.1