Author: [log in to unmask]
Date: Mon Nov 16 10:00:18 2015
New Revision: 3960
Log:
Fix problem with MC truth and beamspot option.
Modified:
java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackUtils.java
java/trunk/tracking/src/main/java/org/hps/recon/tracking/gbl/GBLOutput.java
Modified: java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackUtils.java
=============================================================================
--- java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackUtils.java (original)
+++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/TrackUtils.java Mon Nov 16 10:00:18 2015
@@ -793,16 +793,14 @@
* @return helix object based on the MC particle
*/
public static HelicalTrackFit getHTF(MCParticle mcp, double Bz) {
- boolean debug = true;
- if (debug) {
- System.out.printf("getHTF\n");
- System.out.printf("mcp org %s mc p %s\n", mcp.getOrigin().toString(), mcp.getMomentum().toString());
- }
+ boolean debug = false;
+
+ if (debug) System.out.printf("getHTF\nmcp org %s mc p %s\n", mcp.getOrigin().toString(), mcp.getMomentum().toString());
+
Hep3Vector org = CoordinateTransformations.transformVectorToTracking(mcp.getOrigin());
Hep3Vector p = CoordinateTransformations.transformVectorToTracking(mcp.getMomentum());
- if (debug)
- System.out.printf("mcp org %s mc p %s (trans)\n", org.toString(), p.toString());
+ if (debug) System.out.printf("mcp org %s mc p %s (trans)\n", org.toString(), p.toString());
// Move to x=0 if needed
double targetX = BeamlineConstants.DIPOLE_EDGELOW_TESTRUN;
@@ -820,8 +818,7 @@
// old.toString(),p.toString(),org.toString());
}
- if (debug)
- System.out.printf("mcp org %s mc p %s (trans2)\n", org.toString(), p.toString());
+ if (debug) System.out.printf("mcp org %s mc p %s (trans2)\n", org.toString(), p.toString());
HelixParamCalculator helixParamCalculator = new HelixParamCalculator(p, org, -1 * ((int) mcp.getCharge()), Bz);
double par[] = new double[5];
@@ -831,8 +828,7 @@
par[HelicalTrackFit.curvatureIndex] = 1.0 / helixParamCalculator.getRadius();
par[HelicalTrackFit.z0Index] = helixParamCalculator.getZ0();
HelicalTrackFit htf = getHTF(par);
- System.out.printf("d0 %f z0 %f R %f phi %f lambda %s\n",
- htf.dca(), htf.z0(), htf.R(), htf.phi0(), htf.slope());
+ if(debug) 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;
}
Modified: java/trunk/tracking/src/main/java/org/hps/recon/tracking/gbl/GBLOutput.java
=============================================================================
--- java/trunk/tracking/src/main/java/org/hps/recon/tracking/gbl/GBLOutput.java (original)
+++ java/trunk/tracking/src/main/java/org/hps/recon/tracking/gbl/GBLOutput.java Mon Nov 16 10:00:18 2015
@@ -63,10 +63,12 @@
private boolean addBeamspot = false;
private double beamspotTiltZOverY = 0; //Math.PI/180* 15;
private double beamspotScatAngle = 0.000001;
- // beam spot with in tracking frame
+ // beam spot, in tracking frame
private double beamspotWidthZ = 0.05;
private double beamspotWidthY = 0.150;
- double beamspotPosition[] = {0,0,0};
+ private double beamspotPosition[] = {0,0,0};
+ // human readable ID for beam spot
+ private final int iBeamspotHit = -1;
/**
@@ -274,7 +276,6 @@
int istrip = 0;
- final int iBeamspotHit = -1; // human readable ID for beam spot
int beamSpotMillepedeId = 98; // just a random int number that I came up with
for (int ihit = -1; ihit != hits.size(); ++ihit) {
@@ -282,7 +283,7 @@
HelicalTrackHit hit = null;
HelicalTrackCross htc = null;
List<HelicalTrackStrip> strips;
- List<MCParticle> hitMCParticles = null;
+ List<MCParticle> hitMCParticles = new ArrayList<MCParticle>();
Hep3Vector correctedHitPosition = null;
// Add beamspot first
@@ -390,7 +391,7 @@
SimTrackerHit simHit = simHitsLayerMap.get(strip.layer());
if (isMC) {
- if (simHit == null) {
+ if (simHit == null && ihit != iBeamspotHit) {
System.out.printf("%s: no sim hit for strip hit at layer %d\n", this.getClass().getSimpleName(), strip.layer());
System.out.printf("%s: it as %d mc particles associated with it:\n", this.getClass().getSimpleName(), hitMCParticles.size());
for (MCParticle particle : hitMCParticles) {
@@ -684,6 +685,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(TrackUtils.getHTF(mcp_pair.get(0), -1.0 * this.bFieldVector.z()), TrackUtils.getHTF(mcp_pair.get(1), -1.0 * this.bFieldVector.z()));
+
System.out.printf("%s: invM = %f\n", this.getClass().getSimpleName(), invMassTruth);
System.out.printf("%s: invMTracks = %f\n", this.getClass().getSimpleName(), invMassTruthTrks);
}
|