java/trunk/recon/src/main/java/org/hps/recon/particle
--- java/trunk/recon/src/main/java/org/hps/recon/particle/ReconParticleDriver.java 2014-05-20 20:20:51 UTC (rev 610)
+++ java/trunk/recon/src/main/java/org/hps/recon/particle/ReconParticleDriver.java 2014-05-20 22:57:16 UTC (rev 611)
@@ -25,64 +25,62 @@
*
* @author Mathew Graham <[log in to unmask]>
* @author Omar Moreno <[log in to unmask]>
- * @version $Id$
+ * @version $Id$
*/
public abstract class ReconParticleDriver extends Driver {
-
- // Flags
- boolean debug = false;
-
- // Reconstructed particle collections
- List<ReconstructedParticle> finalStateParticles;
+ // Flags
+ boolean debug = false;
+
+ // Reconstructed particle collections
+ List<ReconstructedParticle> finalStateParticles;
List<ReconstructedParticle> unconstrainedV0Candidates;
List<ReconstructedParticle> beamConV0Candidates;
List<ReconstructedParticle> targetConV0Candidates;
List<ReconstructedParticle> electrons;
List<ReconstructedParticle> positrons;
-
+
// Collections
- String ecalClustersCollectionName = "EcalClusters";
- String tracksCollectionName = "MatchedTracks";
- String finalStateParticlesColName = "FinalStateParticles";
- String unconstrainedV0CandidatesColName = null;
- String beamConV0CandidatesColName = null;
- String targetV0ConCandidatesColName = null;
- String vertexCandidatesColName = null;
- String vertexBeamConsCandidatesName = null;
-
+ String ecalClustersCollectionName = "EcalClusters";
+ String tracksCollectionName = "MatchedTracks";
+ String finalStateParticlesColName = "FinalStateParticles";
+ String unconstrainedV0CandidatesColName = null;
+ String beamConV0CandidatesColName = null;
+ String targetV0ConCandidatesColName = null;
+ String vertexCandidatesColName = null;
+ String vertexBeamConsCandidatesName = null;
+
// The beamsize array is in the tracking frame
/* TODO mg-May 14, 2014: the the beam size from the conditions db...also beam position! */
double[] beamsize = {0.001, 0.2, 0.02};
double maxTrackClusterDistance = 10000; // [mm]
- double bField;
-
-
+ double bField;
+
// 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
+ // so we set the B-field in +iveZ and flip signs of fitted tracks
//
// Note: This should be -1 for test run configurations and +1 for
// prop-2014 configurations
- int flipSign = 1;
+ int flipSign = 1;
public ReconParticleDriver() {
}
- public void setDebug(boolean debug){
- this.debug = debug;
+ public void setDebug(boolean debug) {
+ this.debug = debug;
}
-
- public void setMaxTrackClusterDistance(double maxTrackClusterDistance){
- this.maxTrackClusterDistance = maxTrackClusterDistance;
+
+ public void setMaxTrackClusterDistance(double maxTrackClusterDistance) {
+ this.maxTrackClusterDistance = maxTrackClusterDistance;
}
-
+
public void setBeamSigmaX(double sigma_x) {
- beamsize[1] = sigma_x;
+ beamsize[1] = sigma_x;
}
public void setBeamSigmaY(double sigma_y) {
- beamsize[2] = sigma_y;
+ beamsize[2] = sigma_y;
}
public void setEcalClusterCollectionName(String ecalClustersCollectionName) {
@@ -91,39 +89,41 @@
public void setTrackCollectoinName(String tracksCollectionName) {
this.tracksCollectionName = tracksCollectionName;
- }
-
- @Override
- protected void detectorChanged(Detector detector){
-
- Hep3Vector ip = new BasicHep3Vector(0., 0., 1.);
- bField = detector.getFieldMap().getField(ip).y();
- if(bField < 0) flipSign = -1;
-
- }
+ }
+ @Override
+ protected void detectorChanged(Detector detector) {
+
+ Hep3Vector ip = new BasicHep3Vector(0., 0., 1.);
+ bField = detector.getFieldMap().getField(ip).y();
+ if (bField < 0)
+ flipSign = -1;
+
+ }
+
public void process(EventHeader event) {
-
+
// All events should have a collection of Ecal clusters. If the event
- // doesn't have one, skip the event.
- if (!event.hasCollection(HPSEcalCluster.class, ecalClustersCollectionName)) return;
-
+ // doesn't have one, skip the event.
+ if (!event.hasCollection(HPSEcalCluster.class, ecalClustersCollectionName))
+ return;
+
// Get the collection of Ecal clusters from the event. A triggered
// event should have Ecal clusters. If it doesn't, skip the event.
List<HPSEcalCluster> clusters = event.get(HPSEcalCluster.class, ecalClustersCollectionName);
//if(clusters.isEmpty()) return;
this.printDebug("Number of Ecal clusters: " + clusters.size());
-
+
// Get the collection of tracks from the event
List<Track> tracks = event.get(Track.class, tracksCollectionName);
- this.printDebug("Number of Tracks: " + tracks.size());
+ this.printDebug("Number of Tracks: " + tracks.size());
- finalStateParticles = new ArrayList<ReconstructedParticle>();
- electrons = new ArrayList<ReconstructedParticle>();
- positrons = new ArrayList<ReconstructedParticle>();
- unconstrainedV0Candidates = new ArrayList<ReconstructedParticle>();
- beamConV0Candidates = new ArrayList<ReconstructedParticle>();
- targetConV0Candidates = new ArrayList<ReconstructedParticle>();
+ finalStateParticles = new ArrayList<ReconstructedParticle>();
+ electrons = new ArrayList<ReconstructedParticle>();
+ positrons = new ArrayList<ReconstructedParticle>();
+ unconstrainedV0Candidates = new ArrayList<ReconstructedParticle>();
+ beamConV0Candidates = new ArrayList<ReconstructedParticle>();
+ targetConV0Candidates = new ArrayList<ReconstructedParticle>();
//
finalStateParticles = this.makeReconstructedParticles(clusters, tracks);
@@ -132,150 +132,138 @@
// Put all the reconstructed particles in the event
event.put(finalStateParticlesColName, finalStateParticles, ReconstructedParticle.class, 0);
-
// Loop through the list of final state particles and separate the
// charged particles to either electrons or positrons. These lists
// will be used for vertexing purposes.
- for(ReconstructedParticle finalStateParticle : finalStateParticles){
- if(finalStateParticle.getCharge() > 0) positrons.add(finalStateParticle);
- else if(finalStateParticle.getCharge() < 0) electrons.add(finalStateParticle);
+ for (ReconstructedParticle finalStateParticle : finalStateParticles) {
+ this.printDebug("\t\tThis Charge = " + finalStateParticle.getCharge());
+ if (finalStateParticle.getCharge() > 0)
+ positrons.add(finalStateParticle);
+ else if (finalStateParticle.getCharge() < 0)
+ electrons.add(finalStateParticle);
}
-
+ this.printDebug("Number of Electrons: " + electrons.size());
+ this.printDebug("Number of Positrons: " + positrons.size());
// Vertex electron and positron candidates
vertexParticles(electrons, positrons);
-
+
// If the list exist, put the vertexed candidates into the event
- if(unconstrainedV0CandidatesColName != null)
- this.printDebug("Total number of unconstrained V0 candidates: " + unconstrainedV0Candidates.size());
- event.put(unconstrainedV0CandidatesColName, unconstrainedV0Candidates, ReconstructedParticle.class, 0);
- if(beamConV0CandidatesColName != null)
- this.printDebug("Total number of beam constrained V0 candidates: " + unconstrainedV0Candidates.size());
- event.put(beamConV0CandidatesColName, beamConV0Candidates, ReconstructedParticle.class, 0);
- if(targetV0ConCandidatesColName != null)
- this.printDebug("Total number of target constrained V0 candidates: " + unconstrainedV0Candidates.size());
- event.put(targetV0ConCandidatesColName, targetConV0Candidates, ReconstructedParticle.class, 0);
+ if (unconstrainedV0CandidatesColName != null)
+ this.printDebug("Total number of unconstrained V0 candidates: " + unconstrainedV0Candidates.size());
+ event.put(unconstrainedV0CandidatesColName, unconstrainedV0Candidates, ReconstructedParticle.class, 0);
+ if (beamConV0CandidatesColName != null)
+ this.printDebug("Total number of beam constrained V0 candidates: " + unconstrainedV0Candidates.size());
+ event.put(beamConV0CandidatesColName, beamConV0Candidates, ReconstructedParticle.class, 0);
+ if (targetV0ConCandidatesColName != null)
+ this.printDebug("Total number of target constrained V0 candidates: " + unconstrainedV0Candidates.size());
+ event.put(targetV0ConCandidatesColName, targetConV0Candidates, ReconstructedParticle.class, 0);
}
/**
- *
+ *
*/
abstract void vertexParticles(List<ReconstructedParticle> electrons, List<ReconstructedParticle> positrons);
-
+
+
/**
- *
+ * make the final state particles from clusters & tracks
+ * loop over the tracks first and try to match with clusters
*/
- protected List<ReconstructedParticle> makeReconstructedParticles(List<HPSEcalCluster> clusters, List<Track> tracks){
-
- // Instantiate the list of reconstructed particles
- List<ReconstructedParticle> particles = new ArrayList<ReconstructedParticle>();
-
- // Instantiate the list of matched tracks and clusters. This list
- // will be used to check if a track has already been previously matched
- List<Track> unmatchedTracks = new ArrayList<Track>(tracks);
-
- Track matchedTrack = null;
- // Loop over all of the Ecal clusters and pair them with tracks
- for(HPSEcalCluster cluster : clusters){
-
- // Get the position of the Ecal cluster
- Hep3Vector clusterPosition = new BasicHep3Vector(cluster.getPosition());
- this.printDebug("Ecal cluster position: " + clusterPosition.toString());
-
- double rMax = Double.MAX_VALUE;
- // Loop over all of the tracks in the event
- for(Track track : tracks){
-
- // Check if the Ecal cluster and track are within the same
- // detector volume i.e. both top or bottom
- // TODO: mg-May 14, 2014 does getTrackStates().get(0).getZ0()
- // really get the z-pos (y in detector frame) @ the ECAL or anywhere
- // other than the POCA??? Where is this calculated
- if(clusterPosition.y()*track.getTrackStates().get(0).getZ0() < 0){
- this.printDebug("Track and Ecal cluster are in opposite volumes. Track Z0 = " + track.getTrackStates().get(0).getZ0());
- continue;
- }
-
- // Extrapolate the track to the Ecal cluster position
- Hep3Vector trackPosAtEcal = TrackUtils.extrapolateTrack(track, clusterPosition.z());
-
- // Check if any of the extrapolated values are invalid.
- // TODO: There are some track whose extrapolated coordinates
- // are NaN. The problem seems to be that the y-coordinate
- // of the extrapolated helix is found to be non-real. This
- // needs to be fixed.
- if(Double.isNaN(trackPosAtEcal.x()) || Double.isNaN(trackPosAtEcal.y())) continue;
- this.printDebug("Track position at shower max: " + trackPosAtEcal.toString());
-
- double r = VecOp.sub(trackPosAtEcal, clusterPosition).magnitude();
- this.printDebug("Distance between Ecal cluster and track position: " + r + " mm");
-
+ protected List<ReconstructedParticle> makeReconstructedParticles(List<HPSEcalCluster> clusters, List<Track> tracks) {
+
+ // Instantiate the list of reconstructed particles
+ List<ReconstructedParticle> particles = new ArrayList<ReconstructedParticle>();
+
+ // Instantiate the list of unmatched clusters. Remove if we find track match
+ List<HPSEcalCluster> unmatchedClusters = new ArrayList<HPSEcalCluster>(clusters);
+ HPSEcalCluster matchedCluster = null;
+
+ for (Track track : tracks) {
+ //make the containers for the reconstructed particle
+ ReconstructedParticle particle = new BaseReconstructedParticle();
+ HepLorentzVector fourVector = new BasicHepLorentzVector(0, 0, 0, 0);
+ //add the track information
+ particle.addTrack(track);
+ Hep3Vector momentum = new BasicHep3Vector(track.getTrackStates().get(0).getMomentum());
+ this.printDebug("Momentum in tracking frame: " + momentum.toString());
+ momentum = CoordinateTransformations.transformVectorToDetector(momentum);
+ this.printDebug("Momentum in detector frame: " + momentum.toString());
+ ((BasicHepLorentzVector) fourVector).setV3(fourVector.t(), momentum);
+ ((BaseReconstructedParticle) particle).setCharge(track.getCharge() * flipSign);
+ if (particle.getCharge() > 0)
+ ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(-11, 0, 0, 0));
+ else if (particle.getCharge() < 0)
+ ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(11, 0, 0, 0));
+
+ for (HPSEcalCluster cluster : unmatchedClusters) {
+
+ // Get the position of the Ecal cluster
+ Hep3Vector clusterPosition = new BasicHep3Vector(cluster.getPosition());
+ // Extrapolate the track to the Ecal cluster position
+ Hep3Vector trackPosAtEcal = TrackUtils.extrapolateTrack(track, clusterPosition.z());
+ this.printDebug("Ecal cluster position: " + clusterPosition.toString());
+
+ double rMax = Double.MAX_VALUE;
+
+ // Check if any of the extrapolated values are invalid.
+ // TODO: There are some track whose extrapolated coordinates
+ // are NaN. The problem seems to be that the y-coordinate
+ // of the extrapolated helix is found to be non-real. This
+ // needs to be fixed.
+ if (Double.isNaN(trackPosAtEcal.x()) || Double.isNaN(trackPosAtEcal.y()))
+ continue;
+ this.printDebug("Track position at shower max: " + trackPosAtEcal.toString());
+
+// double r = VecOp.sub(trackPosAtEcal, clusterPosition).magnitude();
+ //don't trust extrapolation...just do y-difference for now
+ double r = Math.abs(clusterPosition.y() - trackPosAtEcal.z());
+ this.printDebug("Distance between Ecal cluster and track position: " + r + " mm");
+
+ // Check if the Ecal cluster and track are within the same
+ // detector volume i.e. both top or bottom
+ if (clusterPosition.y() * trackPosAtEcal.z() < 0) {
+ this.printDebug("Track and Ecal cluster are in opposite volumes. Track Y @ ECAL = " + trackPosAtEcal.z());
+ continue;
+ }
+
if (r < rMax && r <= maxTrackClusterDistance) {
- rMax = r;
- matchedTrack = track;
+ rMax = r;
+ matchedCluster = cluster;
}
- }
-
- // Create a reconstructed particle and add it to the collection of particles
- ReconstructedParticle particle = new BaseReconstructedParticle();
- HepLorentzVector fourVector = new BasicHepLorentzVector(0, 0, 0, 0);
- particle.addCluster(cluster);
- ((BasicHepLorentzVector) fourVector).setT(cluster.getEnergy());
- if(matchedTrack != null){
- particle.addTrack(matchedTrack);
- Hep3Vector momentum = new BasicHep3Vector(matchedTrack.getTrackStates().get(0).getMomentum());
- this.printDebug("Momentum in tracking frame: " + momentum.toString());
- momentum = CoordinateTransformations.transformVectorToDetector(momentum);
- this.printDebug("Momentum in detector frame: " + momentum.toString());
- ((BasicHepLorentzVector) fourVector).setV3(fourVector.t(), momentum);
- ((BaseReconstructedParticle) particle).setCharge(matchedTrack.getCharge()*flipSign);
- if(particle.getCharge() > 0){
- ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(-11, 0, 0, 0));
- } else if(particle.getCharge() < 0){
- ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(11, 0, 0, 0));
- }
- if(unmatchedTracks.contains(matchedTrack)) unmatchedTracks.remove(matchedTrack);
- } else {
- ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(22, 0, 0, 0));
- }
- ((BaseReconstructedParticle) particle).set4Vector(fourVector);
-
- particles.add(particle);
- }
-
- if(!unmatchedTracks.isEmpty()){
- for(Track unmatchedTrack : unmatchedTracks){
-
- // Create a reconstructed particle and add it to the
- // collection of particles
- ReconstructedParticle particle = new BaseReconstructedParticle();
- HepLorentzVector fourVector = new BasicHepLorentzVector(0, 0, 0, 0);
-
- particle.addTrack(unmatchedTrack);
- Hep3Vector momentum = new BasicHep3Vector(unmatchedTrack.getTrackStates().get(0).getMomentum());
- this.printDebug("Momentum in tracking frame: " + momentum.toString());
- momentum = CoordinateTransformations.transformVectorToDetector(momentum);
- this.printDebug("Momentum in detector frame: " + momentum.toString());
- ((BasicHepLorentzVector) fourVector).setV3(momentum.magnitude(), momentum);
- ((BaseReconstructedParticle) particle).setCharge(unmatchedTrack.getCharge()*flipSign);
- ((BaseReconstructedParticle) particle).set4Vector(fourVector);
- if(particle.getCharge() > 0){
- ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(-11, 0, 0, 0));
- } else if(particle.getCharge() < 0){
- ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(11, 0, 0, 0));
- }
- particles.add(particle);
- }
- }
-
- return particles;
+ }
+ if (matchedCluster != null) {
+ particle.addCluster(matchedCluster);
+ ((BasicHepLorentzVector) fourVector).setT(matchedCluster.getEnergy());
+ unmatchedClusters.remove(matchedCluster);
+ }
+ ((BaseReconstructedParticle) particle).set4Vector(fourVector);
+ particles.add(particle);
+ }
+
+ if (!unmatchedClusters.isEmpty())
+ for (HPSEcalCluster unmatchedCluster : unmatchedClusters) {
+ // Create a reconstructed particle and add it to the
+ // collection of particles
+ ReconstructedParticle particle = new BaseReconstructedParticle();
+ HepLorentzVector fourVector = new BasicHepLorentzVector(0, 0, 0, 0);
+
+ particle.addCluster(unmatchedCluster);
+ ((BasicHepLorentzVector) fourVector).setT(unmatchedCluster.getEnergy());
+ ((BaseReconstructedParticle) particle).setCharge(0);
+ ((BaseReconstructedParticle) particle).set4Vector(fourVector);
+ ((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(22, 0, 0, 0));
+ particles.add(particle);
+ }
+
+ return particles;
}
-
+
/**
- *
+ *
* @param debugMessage
*/
- protected void printDebug(String debugMessage){
- if(debug)
- System.out.println(this.getClass().getSimpleName() + ": " + debugMessage);
+ protected void printDebug(String debugMessage) {
+ if (debug)
+ System.out.println(this.getClass().getSimpleName() + ": " + debugMessage);
}
}