Commit in hps-java/src/main/java/org/lcsim/hps/recon/tracking on MAIN | |||
gbl/GBLOutput.java | +44 | -17 | 1.11 -> 1.12 |
/GBLOutputDriver.java | +3 | -1 | 1.6 -> 1.7 |
/TruthResiduals.java | +10 | -30 | 1.2 -> 1.3 |
TrackUtils.java | +46 | -1 | 1.24 -> 1.25 |
+103 | -49 |
Moved common function to util class.
diff -u -r1.11 -r1.12 --- GBLOutput.java 15 Sep 2013 19:24:49 -0000 1.11 +++ GBLOutput.java 17 Sep 2013 00:52:52 -0000 1.12 @@ -121,7 +121,7 @@
// Get track parameters from MC particle
- HelicalTrackFit htfTruth = getHTF(mcp);
+ HelicalTrackFit htfTruth = TrackUtils.getHTF(mcp,-1.0*this._B.z());
// Use the truth helix as the initial track for GBL? //htf = htfTruth;
@@ -363,7 +363,7 @@
if(_debug>0) { double invMassTruth = Math.sqrt( Math.pow(mcp_pair.get(0).getEnergy()+mcp_pair.get(1).getEnergy(),2) - VecOp.add(mcp_pair.get(0).getMomentum(), mcp_pair.get(1).getMomentum()).magnitudeSquared() );
- double invMassTruthTrks = getInvMassTracks(getHTF(mcp_pair.get(0)),getHTF(mcp_pair.get(1)));
+ double invMassTruthTrks = getInvMassTracks(TrackUtils.getHTF(mcp_pair.get(0),-1.0*this._B.z()),TrackUtils.getHTF(mcp_pair.get(1),-1.0*this._B.z()));
System.out.printf("%s: invM = %f\n",this.getClass().getSimpleName(),invMassTruth); System.out.printf("%s: invMTracks = %f\n",this.getClass().getSimpleName(),invMassTruthTrks); }
@@ -560,21 +560,48 @@
} }
- private HelicalTrackFit getHTF(MCParticle mcp) { - Hep3Vector org = this._hpstrans.transformVectorToTracking(mcp.getOrigin()); - Hep3Vector p = this._hpstrans.transformVectorToTracking(mcp.getMomentum()); - HelixParamCalculator helixParamCalculator = new HelixParamCalculator(p, org, -1*((int)mcp.getCharge()), -1.0*this._B.z()); - double par[] = new double[5]; - par[HelicalTrackFit.dcaIndex] = helixParamCalculator.getDCA(); - par[HelicalTrackFit.slopeIndex] = helixParamCalculator.getSlopeSZPlane(); - par[HelicalTrackFit.phi0Index] = helixParamCalculator.getPhi0(); - par[HelicalTrackFit.curvatureIndex] = 1.0/helixParamCalculator.getRadius(); - par[HelicalTrackFit.z0Index] = helixParamCalculator.getZ0(); - SymmetricMatrix cov = new SymmetricMatrix(5); - for(int i=0;i<cov.getNRows();++i) cov.setElement(i, i, 1.); - HelicalTrackFit htf = new HelicalTrackFit(par, cov, new double[2], new int[2], null, null); - return htf; - }
+ + + /** + * Transform MCParticle into a Helix object. + * Note that it produces the helix parameters at nominal x=0 and assumes that there is no field at x<0 + * + * @param mcp MC particle to be transformed + * @return helix object based on the MC particle + */ +// private HelicalTrackFit getHTF(MCParticle mcp) { +// Hep3Vector org = this._hpstrans.transformVectorToTracking(mcp.getOrigin()); +// Hep3Vector p = this._hpstrans.transformVectorToTracking(mcp.getMomentum()); +// // Move to x=0 if needed +// if(org.x() < 0.) { +// double dydx = p.y()/p.x(); +// double dzdx = p.z()/p.x(); +// double delta_x = -1. * org.x(); +// double y = delta_x * dydx; +// double z = delta_x * dzdx; +// double x = org.x() + delta_x; +// if( Math.abs(x) > 1e-8) throw new RuntimeException("Error: origin is not zero!"); +// Hep3Vector old = org; +// org = new BasicHep3Vector(x,y,z); +// System.out.printf("org %s p %s -> org %s\n", old.toString(),p.toString(),org.toString()); +// } else { +// org = this._hpstrans.transformVectorToTracking(mcp.getOrigin()); +// } +// +// +// +// HelixParamCalculator helixParamCalculator = new HelixParamCalculator(p, org, -1*((int)mcp.getCharge()), -1.0*this._B.z()); +// double par[] = new double[5]; +// par[HelicalTrackFit.dcaIndex] = helixParamCalculator.getDCA(); +// par[HelicalTrackFit.slopeIndex] = helixParamCalculator.getSlopeSZPlane(); +// par[HelicalTrackFit.phi0Index] = helixParamCalculator.getPhi0(); +// par[HelicalTrackFit.curvatureIndex] = 1.0/helixParamCalculator.getRadius(); +// par[HelicalTrackFit.z0Index] = helixParamCalculator.getZ0(); +// SymmetricMatrix cov = new SymmetricMatrix(5); +// for(int i=0;i<cov.getNRows();++i) cov.setElement(i, i, 1.); +// HelicalTrackFit htf = new HelicalTrackFit(par, cov, new double[2], new int[2], null, null); +// return htf; +// }
private double truthTrackFitChi2(PerigeeParams perPar, PerigeeParams perParTruth, SymmetricMatrix covariance) { //re-shuffle the param vector to match the covariance order of parameters
diff -u -r1.6 -r1.7 --- GBLOutputDriver.java 15 Sep 2013 19:24:49 -0000 1.6 +++ GBLOutputDriver.java 17 Sep 2013 00:52:52 -0000 1.7 @@ -65,6 +65,8 @@
gbl = new GBLOutput(gblFile,bfield); gbl.setDebug(_debug); gbl.buildModel(detector);
+ gbl.setAPrimeEventFlag(false); + gbl.setXPlaneFlag(false);
truthRes = new TruthResiduals(bfield); truthRes.setDebug(_debug); truthRes.setHideFrame(hideFrame);
@@ -94,7 +96,7 @@
List<SimTrackerHit> simTrackerHits = event.getSimTrackerHits("TrackerHits");
- //truthRes.processSim(mcParticles, simTrackerHits);
+ truthRes.processSim(mcParticles, simTrackerHits);
diff -u -r1.2 -r1.3 --- TruthResiduals.java 30 Aug 2013 01:23:22 -0000 1.2 +++ TruthResiduals.java 17 Sep 2013 00:52:52 -0000 1.3 @@ -7,10 +7,12 @@
import hep.aida.*; import hep.physics.vec.Hep3Vector; import hep.physics.vec.VecOp;
+
import java.util.ArrayList; import java.util.HashMap; import java.util.List; import java.util.Map;
+
import org.lcsim.detector.tracker.silicon.ChargeCarrier; import org.lcsim.detector.tracker.silicon.SiSensor; import org.lcsim.event.MCParticle;
@@ -20,6 +22,7 @@
import org.lcsim.fit.helicaltrack.HelixUtils; import org.lcsim.geometry.Detector; import org.lcsim.hps.event.HPSTransformations;
+import org.lcsim.hps.recon.tracking.TrackUtils;
import org.lcsim.hps.recon.tracking.TrackerHitUtils; import org.lcsim.util.aida.AIDA;
@@ -90,6 +93,10 @@
mcPartSimHitsMap.get(part).add(sh); }
+ + for(MCParticle mcp : mcPartSimHitsMap.keySet()) { + this.h_mcp_org.fill(mcp.getOriginX(), mcp.getOriginY()); + }
// Find the particle responsible for the hit in each layer and compute the residual
@@ -132,7 +139,8 @@
MCParticle mcp = simHit.getMCParticle(); // Get track parameters from MC particle
- HelicalTrackFit htfTruth = getHTF(mcp);
+ //HelicalTrackFit htfTruth = TrackUtils.getHTF(mcp, -1*this._B.z()); + HelicalTrackFit htfTruth = TrackUtils.getHTF(mcp, -1*this._B.z());
// Find pathlength to sim hit double s = HelixUtils.PathToXPlane(htfTruth, simHitPosTracking.x(), 0., 0).get(0);
@@ -172,27 +180,6 @@
} }
- - - - - - - - - private HelicalTrackFit getHTF(MCParticle mcp) { - Hep3Vector org = this._hpstrans.transformVectorToTracking(mcp.getOrigin()); - Hep3Vector p = this._hpstrans.transformVectorToTracking(mcp.getMomentum()); - HelixParamCalculator helixParamCalculator = new HelixParamCalculator(p, org, -1*((int)mcp.getCharge()), -1.0*this._B.z()); - double par[] = new double[5]; - par[HelicalTrackFit.dcaIndex] = helixParamCalculator.getDCA(); - par[HelicalTrackFit.slopeIndex] = helixParamCalculator.getSlopeSZPlane(); - par[HelicalTrackFit.phi0Index] = helixParamCalculator.getPhi0(); - par[HelicalTrackFit.curvatureIndex] = 1.0/helixParamCalculator.getRadius(); - par[HelicalTrackFit.z0Index] = helixParamCalculator.getZ0(); - HelicalTrackFit htf = new HelicalTrackFit(par, null, new double[2], new int[2], null, null); - return htf; - }
@@ -314,14 +301,7 @@
if(this._hideFrame) pl_org.hide(); else pl_org.show();
- this.h_mcp_orgy_resy = hf.createHistogram2D("MC particle res y vs origin y ", 50, -0.2,0.2,50,-0.1,0.1); - IPlotter pl_org_2 = af.createPlotterFactory().create("MC particle res y vs origin y"); - pl_org_2.createRegions(1, 1); - pl_org_2.region(0).plot(this.h_mcp_orgy_resy); - pl_org_2.region(0).style().setParameter("hist2DStyle", "colorMap"); - pl_org_2.region(0).style().dataStyle().fillStyle().setParameter("colorMapScheme", "rainbow"); - if(this._hideFrame) pl_org_2.hide(); - else pl_org_2.show();
+
}
diff -u -r1.24 -r1.25 --- TrackUtils.java 15 Sep 2013 19:24:49 -0000 1.24 +++ TrackUtils.java 17 Sep 2013 00:52:52 -0000 1.25 @@ -17,12 +17,14 @@
import org.lcsim.detector.solids.Point3D; import org.lcsim.detector.solids.Polygon3D; import org.lcsim.detector.tracker.silicon.SiSensor;
+import org.lcsim.event.MCParticle;
import org.lcsim.event.RawTrackerHit; //--- org.lcsim ---// import org.lcsim.event.Track; import org.lcsim.event.TrackerHit; import org.lcsim.fit.helicaltrack.*;
+import org.lcsim.hps.event.HPSTransformations;
import org.lcsim.hps.users.phansson.WTrack; import org.lcsim.recon.tracking.seedtracker.SeedCandidate; import org.lcsim.recon.tracking.seedtracker.SeedTrack;
@@ -30,7 +32,7 @@
/** * * @author Omar Moreno <[log in to unmask]>
- * @version $Id: TrackUtils.java,v 1.24 2013/09/15 19:24:49 phansson Exp $
+ * @version $Id: TrackUtils.java,v 1.25 2013/09/17 00:52:52 phansson Exp $
* TODO: Switch to JLab coordinates */
@@ -753,4 +755,47 @@
return !isTopTrack(htf); }
+ /** + * Transform MCParticle into a Helix object. + * Note that it produces the helix parameters at nominal x=0 and assumes that there is no field at x<0 + * + * @param mcp MC particle to be transformed + * @return helix object based on the MC particle + */ + public static HelicalTrackFit getHTF(MCParticle mcp, double Bz) { + HPSTransformations hpstrans = new HPSTransformations(); + Hep3Vector org = hpstrans.transformVectorToTracking(mcp.getOrigin()); + Hep3Vector p = hpstrans.transformVectorToTracking(mcp.getMomentum()); + // Move to x=0 if needed + if(org.x() < 0.) { + double dydx = p.y()/p.x(); + double dzdx = p.z()/p.x(); + double delta_x = -1. * org.x(); + double y = delta_x * dydx + org.y(); + double z = delta_x * dzdx + org.z(); + double x = org.x() + delta_x; + if( Math.abs(x) > 1e-8) throw new RuntimeException("Error: origin is not zero!"); + Hep3Vector old = org; + org = new BasicHep3Vector(x,y,z); + System.out.printf("org %s p %s -> org %s\n", old.toString(),p.toString(),org.toString()); + } + + + + HelixParamCalculator helixParamCalculator = new HelixParamCalculator(p, org, -1*((int)mcp.getCharge()), Bz); + double par[] = new double[5]; + par[HelicalTrackFit.dcaIndex] = helixParamCalculator.getDCA(); + par[HelicalTrackFit.slopeIndex] = helixParamCalculator.getSlopeSZPlane(); + par[HelicalTrackFit.phi0Index] = helixParamCalculator.getPhi0(); + par[HelicalTrackFit.curvatureIndex] = 1.0/helixParamCalculator.getRadius(); + par[HelicalTrackFit.z0Index] = helixParamCalculator.getZ0(); + SymmetricMatrix cov = new SymmetricMatrix(5); + for(int i=0;i<cov.getNRows();++i) cov.setElement(i, i, 1.); + HelicalTrackFit htf = new HelicalTrackFit(par, cov, new double[2], new int[2], null, null); + System.out.printf("d0 %f z0 %f R %f phi %f lambda %s\n", htf.dca(),htf.z0(),htf.R(),htf.phi0(),htf.slope() ); + return htf; + } +
}
+ +
Use REPLY-ALL to reply to list
To unsubscribe from the LCD-CVS list, click the following link:
https://listserv.slac.stanford.edu/cgi-bin/wa?SUBED1=LCD-CVS&A=1