lcsim-contrib/src/main/java/org/lcsim/contrib/sATLAS
diff -N FindableTrack.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ FindableTrack.java 20 Aug 2009 17:23:58 -0000 1.1
@@ -0,0 +1,286 @@
+/*
+ * FindableTrack.java
+ *
+ * Created on October 24, 2008, 9:50 PM
+ *
+ */
+package org.lcsim.contrib.sATLAS;
+
+import org.lcsim.contrib.mgraham.sATLASDigi.*;
+import hep.physics.vec.BasicHep3Vector;
+import hep.physics.vec.Hep3Vector;
+
+import java.util.ArrayList;
+import java.util.HashSet;
+import java.util.List;
+import java.util.Set;
+
+import org.lcsim.detector.DetectorElementStore;
+import org.lcsim.detector.IDetectorElement;
+import org.lcsim.detector.IDetectorElementContainer;
+import org.lcsim.detector.identifier.IExpandedIdentifier;
+import org.lcsim.detector.identifier.IIdentifier;
+import org.lcsim.detector.identifier.IIdentifierHelper;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+import org.lcsim.event.RelationalTable;
+import org.lcsim.event.SimTrackerHit;
+import org.lcsim.event.base.BaseRelationalTable;
+import org.lcsim.fit.helicaltrack.HelixParamCalculator;
+import org.lcsim.fit.helicaltrack.HitIdentifier;
+import org.lcsim.geometry.subdetector.BarrelEndcapFlag;
+import org.lcsim.recon.tracking.seedtracker.SeedLayer;
+import org.lcsim.recon.tracking.seedtracker.SeedLayer.SeedType;
+import org.lcsim.recon.tracking.seedtracker.SeedStrategy;
+
+/**
+ *
+ * @author Richard Partridge
+ * @version 1.0
+ */
+public class FindableTrack {
+
+ public enum Ignore {
+
+ NoPTCut, NoDCACut, NoZ0Cut, NoSeedCheck, NoConfirmCheck, NoMinHitCut
+ };
+ private double _bfield;
+ private RelationalTable _hittomc;
+ private HitIdentifier _ID;
+
+ /** Creates a new instance of FindableTrack */
+ public FindableTrack(EventHeader event) {
+
+ // Get the magnetic field
+ Hep3Vector IP = new BasicHep3Vector(0., 0., 1.);
+ _bfield = event.getDetector().getFieldMap().getField(IP).z();
+
+ // Instantiate the hit identifier class
+ _ID = new HitIdentifier();
+
+ // Create a relational table that maps SimTrackerHits to MCParticles
+ _hittomc = new BaseRelationalTable(RelationalTable.Mode.MANY_TO_MANY, RelationalTable.Weighting.UNWEIGHTED);
+
+ // Get the collections of SimTrackerHits
+ List<List<SimTrackerHit>> simcols = event.get(SimTrackerHit.class);
+
+ // Loop over the SimTrackerHits and fill in the relational table
+ for (List<SimTrackerHit> simlist : simcols) {
+ for (SimTrackerHit simhit : simlist) {
+ _hittomc.add(simhit, simhit.getMCParticle());
+ }
+ }
+ }
+
+ public boolean isFindable(MCParticle mcp, List<SeedStrategy> slist, Ignore ignore) {
+ List<Ignore> ignores = new ArrayList<Ignore>();
+ ignores.add(ignore);
+ return isFindable(mcp, slist, ignores);
+ }
+
+ public boolean isFindable(MCParticle mcp, List<SeedStrategy> slist) {
+ return isFindable(mcp, slist, new ArrayList<Ignore>());
+ }
+
+ public boolean isFindable(MCParticle mcp, List<SeedStrategy> slist, List<Ignore> ignores) {
+
+ // We can't find neutral particles'
+ if (mcp.getCharge() == 0) {
+ return false;
+ }
+
+ // Find the helix parameters in the L3 convention used by org.lcsim
+ HelixParamCalculator helix = new HelixParamCalculator(mcp, _bfield);
+
+ // We haven't yet determined the track is findable
+ boolean findable = false;
+
+ // Loop over strategies and check if the track is findable
+ for (SeedStrategy strat : slist) {
+
+ // Check the MC Particle's pT
+ if (!CheckPT(helix, ignores, strat)) {
+ continue;
+ }
+
+ // Check the MC Particle's DCA
+ if (!CheckDCA(helix, ignores, strat)) {
+ continue;
+ }
+
+ // Check the MC Particle's Z0
+ if (!CheckZ0(helix, ignores, strat)) {
+ continue;
+ }
+
+ // Check that we have hits on the seed layers
+ if (!CheckSeed(mcp, ignores, strat)) {
+ continue;
+ }
+
+ // Check that we have the required confirmation hits
+ if (!CheckConfirm(mcp, ignores, strat)) {
+ continue;
+ }
+
+ // Check for the minimum number of hits
+ if (!CheckMinHits(mcp, ignores, strat)) {
+ continue;
+ }
+
+ // Passed all the checks - track is findable
+ findable = true;
+ break;
+ }
+
+ return findable;
+ }
+
+ public int LayersHit(MCParticle mcp) {
+
+ // Get the list of hits associated with the MCParticle
+ Set<SimTrackerHit> hitlist = _hittomc.allTo(mcp);
+
+ // Create a set of the identifiers for the hit layers
+ Set<String> idset = new HashSet<String>();
+
+ // Create the set of identifiers
+ for (SimTrackerHit simhit : hitlist) {
+
+// String identifier = simhit.getDetectorElement();
+ String identifier_old = _ID.Identifier(getDetectorElement(simhit));
+ String identifier = _ID.Identifier(simhit);
+// if(!identifier.equals(identifier_old)) System.out.println("old identifier = "+identifier_old+ " new identifier = "+identifier);
+ if (!idset.contains(identifier)) {
+ idset.add(identifier);
+ }
+ }
+
+ return idset.size();
+ }
+
+ private boolean CheckPT(HelixParamCalculator helix, List<Ignore> ignores, SeedStrategy strat) {
+
+ // First see if we are skipping this check
+ if (ignores.contains(Ignore.NoPTCut)) {
+ return true;
+ }
+
+ return helix.getMCTransverseMomentum() >= strat.getMinPT();
+ }
+
+ private boolean CheckDCA(HelixParamCalculator helix, List<Ignore> ignores, SeedStrategy strat) {
+
+ // First see if we are skipping this check
+ if (ignores.contains(Ignore.NoDCACut)) {
+ return true;
+ }
+
+ return Math.abs(helix.getDCA()) <= strat.getMaxDCA();
+ }
+
+ private boolean CheckZ0(HelixParamCalculator helix, List<Ignore> ignores, SeedStrategy strat) {
+
+ // First see if we are skipping this check
+ if (ignores.contains(Ignore.NoZ0Cut)) {
+ return true;
+ }
+
+ return Math.abs(helix.getZ0()) <= strat.getMaxZ0();
+ }
+
+ private boolean CheckSeed(MCParticle mcp, List<Ignore> ignores, SeedStrategy strat) {
+
+ // First see if we are skipping this check
+ if (ignores.contains(Ignore.NoSeedCheck)) {
+ return true;
+ }
+
+ return HitCount(mcp, strat.getLayers(SeedType.Seed)) == 3;
+ }
+
+ private boolean CheckConfirm(MCParticle mcp, List<Ignore> ignores, SeedStrategy strat) {
+
+ // First see if we are skipping this check
+ if (ignores.contains(Ignore.NoConfirmCheck)) {
+ return true;
+ }
+
+ return HitCount(mcp, strat.getLayers(SeedType.Confirm)) >= strat.getMinConfirm();
+ }
+
+ private boolean CheckMinHits(MCParticle mcp, List<Ignore> ignores, SeedStrategy strat) {
+
+ // First see if we are skipping this check
+ if (ignores.contains(Ignore.NoMinHitCut)) {
+ return true;
+ }
+
+ return HitCount(mcp, strat.getLayerList()) >= strat.getMinHits();
+ }
+
+ private int HitCount(MCParticle mcp, List<SeedLayer> lyrlist) {
+
+ // Get the list of hits associated with the MCParticle
+ Set<SimTrackerHit> hitlist = _hittomc.allTo(mcp);
+
+ // Count the number of layers with hits in them
+ int hitcount = 0;
+ for (SeedLayer lyr : lyrlist) {
+
+ // Loop over the hits for this MCParticle
+ for (SimTrackerHit simhit : hitlist) {
+
+ // Get the detector element for this hit
+// IDetectorElement de = getDetectorElement(simhit);
+
+ // Check names
+// String detname_old = _ID.getName(de);
+ String detname_new = simhit.getSubdetector().getName();
+ // if (!detname_old.equals(detname_new)) {
+ // System.out.println("Detector name mismatch - old: "+detname_old+ " new: "+detname_new);
+ // }
+ // int layer_old = _ID.getLayer(de);
+ int layer_new = simhit.getLayer();
+ // if (layer_old != layer_new) {
+ // System.out.println("Layer number mismatch - old: "+layer_old+" new: "+layer_new);
+ // }
+// BarrelEndcapFlag be_old = _ID.getBarrelEndcapFlag(de);
+ BarrelEndcapFlag be_new = simhit.getBarrelEndcapFlag();
+ // if (!be_old.equals(be_new)) {
+ // System.out.println("BarrelEndcapFlag mismatch - old: "+be_old+" new: "+be_new);
+ // }
+
+ // See if this hit is on the layer we are checking
+// if (!lyr.getDetName().equals(_ID.getName(de))) continue;
+// if (lyr.getLayer() != _ID.getLayer(de)) continue;
+// if (!lyr.getBarrelEndcapFlag().equals(_ID.getBarrelEndcapFlag(de)))
+// continue;
+ if (!lyr.getDetName().equals(detname_new)) {
+ continue;
+ }
+ if (lyr.getLayer() != layer_new) {
+ continue;
+ }
+ if (!lyr.getBarrelEndcapFlag().equals(be_new)) {
+ continue;
+ }
+ hitcount++;
+ break;
+ }
+ }
+
+ return hitcount;
+ }
+
+ private IDetectorElement getDetectorElement(SimTrackerHit hit) {
+ IDetectorElementContainer cont = DetectorElementStore.getInstance().find(hit.getIdentifier());
+ IDetectorElement de;
+ if (cont.isEmpty()) {
+ throw new RuntimeException("Detector Container is empty!");
+ } else {
+ de = cont.get(0);
+ }
+ return de;
+ }
+}
\ No newline at end of file
