Author: [log in to unmask]
Date: Tue Dec 8 16:14:30 2015
New Revision: 4029
Log:
updates, including greatly loosening matching requirement while storing matching quality factor in ReconParticle
Modified:
java/branches/nathan-dev/recon/src/main/java/org/hps/recon/particle/ReconParticleDriver.java
java/branches/nathan-dev/recon/src/main/java/org/hps/recon/utils/TrackClusterMatcher.java
Modified: java/branches/nathan-dev/recon/src/main/java/org/hps/recon/particle/ReconParticleDriver.java
=============================================================================
--- java/branches/nathan-dev/recon/src/main/java/org/hps/recon/particle/ReconParticleDriver.java (original)
+++ java/branches/nathan-dev/recon/src/main/java/org/hps/recon/particle/ReconParticleDriver.java Tue Dec 8 16:14:30 2015
@@ -6,6 +6,8 @@
import hep.physics.vec.HepLorentzVector;
import hep.physics.vec.VecOp;
+import java.io.FileWriter;
+import java.io.IOException;
import java.util.ArrayList;
import java.util.HashMap;
import java.util.HashSet;
@@ -16,7 +18,6 @@
import org.hps.recon.tracking.CoordinateTransformations;
import org.hps.recon.tracking.TrackUtils;
import org.hps.recon.utils.TrackClusterMatcher;
-import org.lcsim.event.CalorimeterHit;
import org.lcsim.event.Cluster;
import org.lcsim.event.EventHeader;
import org.lcsim.event.ReconstructedParticle;
@@ -25,8 +26,8 @@
import org.lcsim.event.base.BaseCluster;
import org.lcsim.event.base.BaseReconstructedParticle;
import org.lcsim.geometry.Detector;
+import org.lcsim.geometry.subdetector.HPSEcal3;
import org.lcsim.util.Driver;
-import org.lcsim.geometry.subdetector.HPSEcal3;
/**
@@ -50,6 +51,9 @@
public static final int MOLLER_TOP = 0;
public static final int MOLLER_BOT = 1;
+ // normalized cluster-track distance required for qualifying as a match:
+ private double MAXNSIGMAPOSITIONMATCH=30.0;
+
HPSEcal3 ecal;
/**
@@ -180,6 +184,15 @@
this.trackCollectionNames = trackCollectionNames;
}
+ /**
+ * Set the requirement on cluster-track position matching in terms of N-sigma.
+ *
+ * @param nsigma
+ */
+ public void setNSigmaPositionMatch(double nsigma) {
+ MAXNSIGMAPOSITIONMATCH=nsigma;
+ }
+
/**
@@ -235,15 +248,16 @@
// Create a mapping of matched clusters to corresponding tracks.
HashMap<Cluster, Track> clusterToTrack = new HashMap<Cluster,Track>();
-
+
// Loop through all of the track collections and try to match every
// track to a cluster. Allow a cluster to be matched to multiple
// tracks and use a probability (to be coded later) to determine what
// the best match is.
// TODO: At some point, pull this out to it's own method
for (List<Track> tracks : trackCollections) {
+
for (Track track : tracks) {
-
+
// Create a reconstructed particle to represent the track.
ReconstructedParticle particle = new BaseReconstructedParticle();
@@ -257,6 +271,9 @@
// Derive the charge of the particle from the track.
((BaseReconstructedParticle) particle).setCharge(track.getCharge() * flipSign);
+
+ // initialize PID quality to a junk value:
+ ((BaseReconstructedParticle)particle).setGoodnessOfPid(-9999);
// Extrapolate the particle ID from the track. Positively
// charged particles are assumed to be positrons and those
@@ -267,9 +284,6 @@
((BaseReconstructedParticle) particle).setParticleIdUsed(new SimpleParticleID(11, 0, 0, 0));
}
- // normalized cluster-track distance required for qualifying as a match:
- final double maximumNSigma=5.0;
-
// normalized distance of the closest match:
double smallestNSigma=Double.MAX_VALUE;
@@ -278,33 +292,40 @@
for (Cluster cluster : clusters) {
// normalized distance between this cluster and track:
- final double thisNSigma=matcher.getNSigmaPosition(cluster, track);
-
- // ignore if distance doesn't make the cut:
- if (thisNSigma > maximumNSigma) continue;
-
- // the cluster with the smallest normalized distance will be the match:
- if (thisNSigma < smallestNSigma) {
-
- smallestNSigma = thisNSigma;
- matchedCluster = cluster;
-
- // prefer using GBL tracks to actually correct the clusters (later):
- if (track.getType() >= 32 || !clusterToTrack.containsKey(matchedCluster))
- clusterToTrack.put(matchedCluster,track);
- }
+ final double thisNSigma=matcher.getNSigmaPosition(cluster, particle);
+
+ // ignore if matching quality doesn't make the cut:
+ if (thisNSigma > MAXNSIGMAPOSITIONMATCH) continue;
+
+ // ignore if we already found a cluster that's a better match:
+ if (thisNSigma > smallestNSigma) continue;
+
+ // we found a new best cluster candidate for this track:
+ smallestNSigma = thisNSigma;
+ matchedCluster = cluster;
+
+ // prefer using GBL tracks to correct (later) the clusters, for some consistency:
+ if (track.getType() >= 32 || !clusterToTrack.containsKey(matchedCluster)) {
+ clusterToTrack.put(matchedCluster,track);
+ }
}
// If a cluster was found that matches the track...
if (matchedCluster != null) {
+ // add cluster to the particle:
particle.addCluster(matchedCluster);
- int pid = particle.getParticleIDUsed().getPDG();
+ // use pid quality to store track-cluster matching quality:
+ ((BaseReconstructedParticle)particle).setGoodnessOfPid(smallestNSigma);
+
+ // propogate pid to the cluster:
+ final int pid = particle.getParticleIDUsed().getPDG();
if (Math.abs(pid) == 11) {
((BaseCluster) matchedCluster).setParticleId(pid);
}
+ // unmatched clusters will (later) be used to create photon particles:
unmatchedClusters.remove(matchedCluster);
}
Modified: java/branches/nathan-dev/recon/src/main/java/org/hps/recon/utils/TrackClusterMatcher.java
=============================================================================
--- java/branches/nathan-dev/recon/src/main/java/org/hps/recon/utils/TrackClusterMatcher.java (original)
+++ java/branches/nathan-dev/recon/src/main/java/org/hps/recon/utils/TrackClusterMatcher.java Tue Dec 8 16:14:30 2015
@@ -1,8 +1,4 @@
package org.hps.recon.utils;
-
-import java.io.IOException;
-import java.util.HashMap;
-import java.util.Map;
import hep.aida.IAnalysisFactory;
import hep.aida.IHistogram1D;
@@ -13,12 +9,17 @@
import hep.physics.vec.BasicHep3Vector;
import hep.physics.vec.Hep3Vector;
+import java.io.IOException;
+import java.util.HashMap;
+import java.util.Map;
+
+import org.hps.recon.tracking.CoordinateTransformations;
+import org.hps.recon.tracking.TrackUtils;
import org.lcsim.event.Cluster;
+import org.lcsim.event.ReconstructedParticle;
import org.lcsim.event.Track;
import org.lcsim.event.TrackState;
import org.lcsim.geometry.FieldMap;
-import org.hps.recon.tracking.CoordinateTransformations;
-import org.hps.recon.tracking.TrackUtils;
/**
* Utility used to determine if a track and cluster are matched.
@@ -66,8 +67,11 @@
/**
* Rafo's parameterization of cluster-seed x/y position residuals as function of energy.
+ *
* Derived using GBL/seed tracks, non-analytic extrapolation, uncorrected cluster positions,
- * and EngRun2015-Nominal-v3-4-fieldmap detector
+ * and EngRun2015-Nominal-v3-4-fieldmap detector.
+ *
+ * f = p0+e*(p1+e*(p2+e*(p3+e*(p4+e*p5))))
*/
private static final double dxMeanTopPosiGBL[] = { 6.67414,-9.57296, 5.70647, 27.4523,-28.1103,-9.11424 };
private static final double dxSigmTopPosiGBL[] = { 52.6437,-478.805, 1896.73,-3761.48, 3676.77,-1408.31 };
@@ -255,8 +259,11 @@
*
* @return #sigma between cluster and track positions
*/
- public double getNSigmaPosition(Cluster cluster,Track track) {
-
+ public double getNSigmaPosition(Cluster cluster,ReconstructedParticle particle) {
+
+ if (particle.getTracks().size()<1) return Double.MAX_VALUE;
+ Track track=particle.getTracks().get(0);
+
if (this.useAnalyticExtrapolator)
throw new RuntimeException("This is to be used with non-analytic extrapolator only.");
@@ -276,10 +283,10 @@
// whether it's a GBL track:
final boolean isGBL = track.getType() >= 32;
-
+
// choose which parameterization of mean and sigma to use:
double dxMean[],dyMean[],dxSigm[],dySigm[];
- if (track.getCharge()==1) {
+ if (particle.getCharge()>0) {
if (isTopTrack) {
dxMean = isGBL ? dxMeanTopPosiGBL : dxMeanTopPosiSeed;
dxSigm = isGBL ? dxSigmTopPosiGBL : dxSigmTopPosiSeed;
@@ -293,7 +300,7 @@
dySigm = isGBL ? dySigmBotPosiGBL : dySigmBotPosiSeed;
}
}
- else if (track.getCharge()==-1) {
+ else if (particle.getCharge()<0) {
if (isTopTrack) {
dxMean = isGBL ? dxMeanTopElecGBL : dxMeanTopElecSeed;
dxSigm = isGBL ? dxSigmTopElecGBL : dxSigmTopElecSeed;
@@ -310,11 +317,12 @@
else return Double.MAX_VALUE;
// get particle energy:
- final double pp[] = track.getMomentum();
- double ee = Math.sqrt( pp[0]*pp[0] + pp[1]*pp[1] + pp[2]*pp[2] );
-
- // force the parameterizations to be constant above a certain energy:
- if (ee > 0.6) ee=0.6;
+ Hep3Vector p3 = new BasicHep3Vector(track.getTrackStates().get(0).getMomentum());
+ p3 = CoordinateTransformations.transformVectorToDetector(p3);
+ double ee = p3.magnitude();
+
+ // Rafo's parameterization isn't measured above 650 MeV/c but expected to be constant:
+ if (ee > 0.65) ee=0.65;
// calculate measured mean and sigma of deltaX and deltaY for this energy:
double aDxMean=0,aDxSigm=0,aDyMean=0,aDySigm=0;
@@ -327,6 +335,7 @@
final double nSigmaX = (cPos.x() - tPos.x() - aDxMean) / aDxSigm;
final double nSigmaY = (cPos.y() - tPos.y() - aDyMean) / aDySigm;
return Math.sqrt(nSigmaX*nSigmaX + nSigmaY*nSigmaY);
+ //return Math.sqrt( 1 / ( 1/nSigmaX/nSigmaX + 1/nSigmaY/nSigmaY ) );
}
@@ -520,4 +529,16 @@
e.printStackTrace();
}
}
+
+ /**
+ * Class to store track-cluster matching qualities.
+ */
+ public class TrackClusterMatch {
+ private double nSigmaPositionMatch=Double.MAX_VALUE;
+ public TrackClusterMatch(ReconstructedParticle pp, Cluster cc) {
+ nSigmaPositionMatch = getNSigmaPosition(cc,pp);
+ }
+ public double getNSigmaPositionMatch() { return nSigmaPositionMatch; }
+ }
+
}
|
