6 added files
lcsim-contrib/src/main/java/org/lcsim/contrib/onoprien/crux/cheat
diff -N CheatClusteringDriver.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ CheatClusteringDriver.java 22 Jan 2009 21:12:42 -0000 1.1
@@ -0,0 +1,230 @@
+package org.lcsim.contrib.onoprien.crux.cheat;
+
+import java.util.*;
+import org.lcsim.contrib.onoprien.crux.CruxManager;
+import org.lcsim.contrib.onoprien.crux.geom.CalGeometry;
+
+import org.lcsim.event.CalorimeterHit;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+import org.lcsim.event.SimCalorimeterHit;
+import org.lcsim.recon.cluster.util.BasicCluster;
+
+import org.lcsim.contrib.onoprien.util.job.Driver;
+
+import org.lcsim.contrib.onoprien.crux.infrastructure.*;
+
+/**
+ * Cheater that attaches hits and tracks to reconstructed particles.
+ *
+ * @author D. Onoprienko
+ * @version $Id: CheatClusteringDriver.java,v 1.1 2009/01/22 21:12:42 onoprien Exp $
+ */
+public class CheatClusteringDriver extends Driver {
+
+// -- Constructors : ----------------------------------------------------------
+
+ public CheatClusteringDriver() {
+ _geom = CruxManager.defaultInstance().getGeometry();
+ }
+
+// -- Setters : ---------------------------------------------------------------
+
+ /**
+ * Set any parameter.
+ * The following parameters can be set with this method:
+ * <p><dl>
+ * <dt>"RECO_PARTICLE_LIST"</dt> <dd>Name of the reconstructed particle list this
+ * driver works with.
+ * <br>Default: "CruxParticles".</dd>
+ * <dt>"INPUT_HIT_MAP"</dt> <dd>Name of input map of <tt>CalorimeterHits</tt> to
+ * be used by this driver.
+ * <br>No default, must be set explicitly if this class is to be used as a driver.</dd>
+ * <dt>"OUTPUT_HIT_MAP"</dt> <dd>Name of output map of <tt>CalorimeterHits</tt> to
+ * be created and saved into event by this driver.
+ * <br>No default, must be set explicitly if this class is to be used as a driver.</dd>
+ * <dt>"INPUT_TRACK_LIST"</dt> <dd>Name of input list of <tt>Tracks</tt>. Tracks from this list will be attached
+ * to reconstructed particles based on their MC truth information.
+ * <br>Default: <tt>null</tt> (no tracks will be attached).</dd>
+ * <dt>"OUTPUT_TRACK_LIST"</dt> <dd>Name of output list of <tt>Tracks</tt>. If this name is set,
+ * a list of tracks from the input list that have not been attached to any reconstructed
+ * particles will be saved into the event record.
+ * <br>Default: <tt>null</tt>.</dd></dl>
+ *
+ * @param name Name of parameter to be set. Case is ignored.
+ * @param values List of values to be used for setting the parameter.
+ * @throws NoSuchParameterException Thrown if the supplied parameter name is unknown.
+ * @throws IllegalArgumentException Thrown if incorrect nunber of values, or a value
+ * of incorrect type is supplied.
+ */
+ public void set(String name, Object... values) {
+ Object value = values.length == 0 ? null : values[0];
+ try {
+ if (name.equalsIgnoreCase("RECO_PARTICLE_LIST")) {
+ _recoListName = (String) value;
+ } else if (name.equalsIgnoreCase("INPUT_HIT_MAP")) {
+ _inMapName = (String) value;
+ } else if (name.equalsIgnoreCase("OUTPUT_HIT_MAP")) {
+ _outMapName = (String) value;
+ } else if (name.equalsIgnoreCase("INPUT_TRACK_LIST")) {
+ _inTrackListName = (String) value;
+ } else if (name.equalsIgnoreCase("OUTPUT_TRACK_LIST")) {
+ _outTrackListName = (String) value;
+ } else {
+ super.set(name, values);
+ }
+ } catch (ClassCastException x) {
+ throw new IllegalArgumentException(ERR_VIT, x);
+ }
+ }
+
+// -- Event processing : ------------------------------------------------------
+
+ /** Called by the framework to process an event. */
+ public void process(EventHeader event) {
+
+ // Fetch reconstructed particle list
+
+ List<CruxParticle> recoList = null;
+ try {
+ recoList = (List<CruxParticle>) event.get(_recoListName);
+ } catch (IllegalArgumentException x) {
+ throw new RuntimeException("Failed to fetch the list of reco particles: " + _recoListName);
+ }
+
+ // Fetch track list and calorimeter hit map
+
+ List<CruxTrack> inTrackList = null;
+ try {
+ inTrackList = (List<CruxTrack>) event.get(_inTrackListName);
+ } catch (IllegalArgumentException x) {}
+
+ CruxHitMap inHitMap = null;
+ try {
+ inHitMap = (CruxHitMap) event.get(_inMapName);
+ } catch (IllegalArgumentException x) {}
+
+ // Attach hits and tracks to reconstructed particles and put unused into the event
+
+ if (inHitMap != null) {
+ CruxHitMap outMap = attachHits(recoList, inHitMap);
+ if (_outMapName != null) event.put(_outMapName, inHitMap);
+ }
+
+ if (inTrackList != null) {
+ List<CruxTrack> outTrackList = attachTracks(recoList, inTrackList);
+ event.put(_outTrackListName, outTrackList, CruxTrack.class, 0);
+ }
+ }
+
+// -- Attaching tracks and hits : ---------------------------------------------
+
+ /**
+ * Attach hits from the given map to reconstructed particles from the given list.
+ * Returns a map containing unused hits.
+ */
+ public CruxHitMap attachHits(List<CruxParticle> recoList, CruxHitMap inMap) {
+
+ ArrayList<CalorimeterHit> remainingHits = new ArrayList<CalorimeterHit>(inMap.size());
+ HashMap<CruxParticle, ArrayList<CalorimeterHit>> cMap
+ = new HashMap<CruxParticle, ArrayList<CalorimeterHit>>(recoList.size()*3);
+
+ for (CalorimeterHit hit : inMap.values()) {
+ MCParticle mc = getMCParticle(hit);
+ boolean notFound = true;
+ while ((mc != null) && notFound) {
+ for (CruxParticle part : recoList) {
+ if (part.getMCParticles().contains(mc)) {
+ ArrayList<CalorimeterHit> hits = cMap.get(part);
+ if (hits == null) {
+ hits = new ArrayList<CalorimeterHit>(50);
+ cMap.put(part, hits);
+ }
+ hits.add(hit);
+ notFound = false;
+ break;
+ }
+ }
+ List<MCParticle> mcList = mc.getParents();
+ mc = (mcList.size() == 1) ? mcList.get(0) : null;
+ }
+ if (notFound) remainingHits.add(hit);
+ }
+
+ for (Map.Entry<CruxParticle, ArrayList<CalorimeterHit>> entry : cMap.entrySet()) {
+ createClusters(entry.getKey(), entry.getValue());
+ }
+
+ CruxHitMap outMap = new CruxHitMap((int)(remainingHits.size()*1.34)+1);
+ outMap.addAll(remainingHits);
+ return outMap;
+ }
+
+ /**
+ * Attach tracks from the given list to reconstructed particles from the given list.
+ * Returns a list of unused tracks.
+ */
+ public ArrayList<CruxTrack> attachTracks(List<CruxParticle> particles, List<CruxTrack> tracks) {
+ ArrayList<CruxTrack> remainingTracks = new ArrayList<CruxTrack>(tracks.size());
+ for (CruxTrack track : tracks) {
+ boolean notFound = true;
+ MCParticle mc = track.getMCParticle();
+ if (mc == null) break;
+ for (CruxParticle particle : particles) {
+ if (particle.getMCParticles().contains(mc)) {
+ particle.setTrack(track);
+ notFound = false;
+ break;
+ }
+ }
+ if (notFound) remainingTracks.add(track);
+ }
+ remainingTracks.trimToSize();
+ return remainingTracks;
+ }
+
+// -- Helper methods : --------------------------------------------------------
+
+ /**
+ * Returns <tt>MCParticle</tt> associated with the given hit.
+ * Default implementation chooses <tt>MCParticle</tt> that made the biggest contribution
+ * to the hit's energy. Subclusses may override this method to use different association
+ * algorithm.
+ */
+ protected MCParticle getMCParticle(CalorimeterHit hit) {
+ SimCalorimeterHit sHit = (SimCalorimeterHit) hit;
+ MCParticle mc = null;
+ double maxEnergy = -1.;
+ int count = sHit.getMCParticleCount();
+ for (int i=0; i<count; i++) {
+ double energy = sHit.getContributedEnergy(i);
+ if (energy > maxEnergy) {
+ maxEnergy = energy;
+ mc = sHit.getMCParticle(i);
+ }
+ }
+ return mc;
+ }
+
+ /**
+ * Create <tt>Cluster</tt> from the list of hits and associate it with the given
+ * reconstructed particle.
+ */
+ protected void createClusters(CruxParticle particle, ArrayList<CalorimeterHit> hitList) {
+ BasicCluster cluster = new BasicCluster();
+ for (CalorimeterHit hit : hitList) {
+ cluster.addHit(hit);
+ }
+ particle.addCluster(cluster);
+ }
+
+// -- Private parts : ---------------------------------------------------------
+
+ private CalGeometry _geom;
+
+ protected String _recoListName = "CruxParticles";
+ protected String _inMapName;
+ protected String _outMapName;
+ private String _inTrackListName;
+ private String _outTrackListName;
+}
lcsim-contrib/src/main/java/org/lcsim/contrib/onoprien/crux/cheat
diff -N CheatRecoDefinition.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ CheatRecoDefinition.java 22 Jan 2009 21:12:42 -0000 1.1
@@ -0,0 +1,171 @@
+package org.lcsim.contrib.onoprien.crux.cheat;
+
+import java.util.*;
+
+import hep.physics.particle.properties.ParticleType;
+import hep.physics.vec.Hep3Vector;
+import org.lcsim.event.MCParticle;
+
+import org.lcsim.contrib.onoprien.util.NoSuchParameterException;
+
+/**
+ * Class that defines reconstructable MCParticle. Used by track finder cheaters and
+ * tracking performance analysis classes. The definition can be customized either by
+ * setting parameters through the {@link #setCut setCut(String,Object)} method (by default, all cuts are
+ * disabled, so every particle is reconstructable) or by overriding {@link #isReconstructable(MCParticle)}.
+ *
+ * @author D. Onoprienko
+ * @version $Id: CheatRecoDefinition.java,v 1.1 2009/01/22 21:12:42 onoprien Exp $
+ */
+public class CheatRecoDefinition implements ICheatRecoDefinition {
+
+// -- Constructors : ----------------------------------------------------------
+
+ public CheatRecoDefinition() {
+ this("");
+ }
+
+ public CheatRecoDefinition(String name) {
+
+ _name = name;
+
+ _cut_Charged = 0;
+ _cut_Pt = 0.;
+ _cut_ThetaMin = 0.;
+ _cut_ThetaMax = Math.PI;
+ _cut_MinTrackerSuperLayers = 0;
+ _cut_Types = null;
+ }
+
+
+// -- Setters : ---------------------------------------------------------------
+
+ /**
+ * Set cut parameters. Calling this method this a name of a cat but no values diables the cut.
+ * The following cuts can be configured with this method:
+ * <p><dl>
+ * <dt>"CHARGED"</dt> <dd>Boolean value: if <tt>true</tt>, only charged particles are selected.
+ * If <tt>false</tt>, only neutral particles are selected.</dd>
+ * <dt>"PT"</dt> <dd>Double value: minimum transverse momentum.</dd>
+ * <dt>"THETA"</dt> <dd>Two double values: mimimum and maximum theta. If only one value is specified,
+ * it is used as minimum, and symmetric theta range is selected.</dd>
+ * <dt>"MIN_TRACK_HIT_LAYERS"</dt> <dd>Integer value: minimum number of tracker superlayers with hits.
+ * {@link #getNumTrackerSuperLayers getNumTrackerSuperLayers(MCParticle)} method must be implemented.</dd>
+ * <dt>"PARTICLE_TYPES"</dt> <dd>One ore more <tt>ParticleType</tt> values: only particles
+ * of these types will be selected.</dd></dl>
+ *
+ * By default, all cuts are disabled.
+ *
+ * @param name Name of the cut. Case is ignored.
+ * @param values List of parameters to be used for configuring the cut.
+ * @throws NoSuchParameterException Thrown if the supplied parameter name is unknown.
+ * @throws IllegalArgumentException Thrown if incorrect number of values, or a value
+ * of incorrect type is supplied.
+ */
+ public void setCut(String name, Object... values) {
+ Object value = values.length == 0 ? null : values[0];
+ try {
+ if (name.equalsIgnoreCase("CHARGED")) {
+ if (values.length == 0) {
+ _cut_Charged = 0;
+ } else {
+ if ((Boolean)values[0]) {
+ _cut_Charged = 1;
+ } else {
+ _cut_Charged = 2;
+ }
+ }
+ } else if (name.equalsIgnoreCase("PT")) {
+ if (values.length == 0) _cut_Pt = 0.;
+ _cut_Pt = (Double) value;
+ } else if (name.equalsIgnoreCase("THETA")) {
+ if (values.length == 0) {
+ _cut_ThetaMin = 0.;
+ _cut_ThetaMax = Math.PI;
+ } else if (values.length == 1) {
+ _cut_ThetaMin = (Double) values[0];
+ _cut_ThetaMax = Math.PI - _cut_ThetaMin;
+ } else {
+ _cut_ThetaMin = (Double) values[0];
+ _cut_ThetaMax = (Double) values[1];
+ }
+ } else if (name.equalsIgnoreCase("MIN_TRACK_HIT_LAYERS")) {
+ if (values.length == 0) _cut_MinTrackerSuperLayers = 0;
+ _cut_MinTrackerSuperLayers = (Integer) value;
+ } else if (name.equalsIgnoreCase("PARTICLE_TYPES")) {
+ if (values.length == 0) {
+ _cut_Types = null;
+ } else {
+ _cut_Types = new HashSet<ParticleType>();
+ for (Object o : values) {
+ _cut_Types.add((ParticleType)o);
+ }
+ }
+ } else {
+ throw new NoSuchParameterException();
+ }
+ } catch (ClassCastException x) {
+ throw new IllegalArgumentException(x);
+ }
+ }
+
+
+// -- Implementing ICheatRecoDefinition : ------------------------------------------
+
+ public String getName() {return _name;}
+
+
+// -- Methods defining "reconstructable" and "reconstructed" : ----------------
+
+ public boolean isTrackFindable(MCParticle mcParticle) {
+
+ // Charge cut
+
+ double charge = mcParticle.getCharge();
+ if (_cut_Charged == 1) {
+ if (Math.abs(charge) < .9) return false;
+ } else if (_cut_Charged == 2) {
+ if (Math.abs(charge) > .1) return false;
+ }
+
+ // Pt cut
+
+ Hep3Vector p = mcParticle.getMomentum();
+ double pT = Math.hypot(p.x(), p.y());
+ if (pT < _cut_Pt) return false;
+
+ // Theta cut
+
+ double theta = Math.acos(p.z()/p.magnitude());
+ if (theta < _cut_ThetaMin || theta > _cut_ThetaMax) return false;
+
+ // Minimum number of hit layers cut
+
+ if (_cut_MinTrackerSuperLayers > 0 && _cut_MinTrackerSuperLayers > getNumTrackerSuperLayers(mcParticle)) return false;
+
+ // Particle type cut
+
+ if (_cut_Types != null && (! _cut_Types.contains(mcParticle.getType()))) return false;
+
+ // Passed all cuts !
+
+ return true;
+ }
+
+
+// -- Helper methods to be implemented by subclasses: -------------------------
+
+ protected int getNumTrackerSuperLayers(MCParticle mc) {
+ throw new UnsupportedOperationException();
+ }
+
+// -- Private parts : ---------------------------------------------------------
+
+ private String _name;
+
+ protected int _cut_Charged;
+ protected double _cut_Pt;
+ protected double _cut_ThetaMin, _cut_ThetaMax;
+ protected int _cut_MinTrackerSuperLayers;
+ protected Set<ParticleType> _cut_Types;
+}
lcsim-contrib/src/main/java/org/lcsim/contrib/onoprien/crux/cheat
diff -N CheatRecoParticleCreator.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ CheatRecoParticleCreator.java 22 Jan 2009 21:12:42 -0000 1.1
@@ -0,0 +1,344 @@
+package org.lcsim.contrib.onoprien.crux.cheat;
+
+import java.util.*;
+import static java.util.logging.Level.*;
+
+import hep.physics.particle.properties.ParticleType;
+import hep.physics.particle.properties.UnknownParticleIDException;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+
+import org.lcsim.contrib.onoprien.util.job.Driver;
+
+import org.lcsim.contrib.onoprien.crux.infrastructure.*;
+import org.lcsim.contrib.onoprien.crux.mctruth.MCTruthCrux;
+
+/**
+ * Driver that creates a tree of reconstructed particles and vertices based on
+ * Monte Carlo truth information, packs it into {@link CruxParticleList}, and
+ * puts it into the event record.
+ *
+ * @author D. Onoprienko
+ * @version $Id: CheatRecoParticleCreator.java,v 1.1 2009/01/22 21:12:42 onoprien Exp $
+ */
+public class CheatRecoParticleCreator extends Driver {
+
+// -- Constructors : ----------------------------------------------------------
+
+ public CheatRecoParticleCreator() {
+ //set("LOG_LEVEL", FINEST);
+ }
+
+// -- Setters : ---------------------------------------------------------------
+
+ /**
+ * Set any parameter.
+ * The following parameters can be set with this method:
+ * <p><dl>
+ * <dt>"RECO_PARTICLE_LIST"</dt> <dd>Name of output list of reconstructed particles.
+ * <br>Default: "CruxParticles".</dd></dl>
+ *
+ * @param name Name of parameter to be set. Case is ignored.
+ * @param values List of values to be used for setting the parameter.
+ * @throws NoSuchParameterException Thrown if the supplied parameter name is unknown.
+ * @throws IllegalArgumentException Thrown if incorrect nunber of values, or a value
+ * of incorrect type is supplied.
+ */
+ public void set(String name, Object... values) {
+ Object value = values.length == 0 ? null : values[0];
+ try {
+ if (name.equalsIgnoreCase("RECO_PARTICLE_LIST")) {
+ _recoListName = (String) value;
+ } else {
+ super.set(name, values);
+ }
+ } catch (ClassCastException x) {
+ throw new IllegalArgumentException(ERR_VIT, x);
+ }
+ }
+
+// -- Event processing : ------------------------------------------------------
+
+ /** Called by the framework to process an event. */
+ public void process(EventHeader event) {
+
+ log("Event " + event.getEventNumber(), FINE);
+
+ super.process(event);
+
+ List<MCParticle> mcList = event.getMCParticles();
+
+ // Fetch or create a list of reconstructed particles
+
+ try {
+ _recoList = (CruxParticleList) event.get(_recoListName);
+ _recoList.clear();
+ } catch (IllegalArgumentException x) {
+ _recoList = new CruxParticleList();
+ event.put(_recoListName, _recoList);
+ }
+
+ // Create CruxParticle for each FINAL_STATE MCParticle
+
+ ArrayList<CruxParticle> fsList = new ArrayList<CruxParticle>(mcList.size());
+
+ for (MCParticle mc : mcList) {
+ if (mc.getGeneratorStatus() == MCParticle.FINAL_STATE) {
+ CruxParticle part = new CruxParticle();
+ part.setMCParticle(mc);
+ fsList.add(part);
+ _recoList.add(part);
+ }
+ }
+
+ // Grow tree from FS particles up, skipping (multiple parent particles) and
+ // merging vertices if they have the same parent
+
+ ArrayList<CruxParticle> pList = new ArrayList<CruxParticle>(mcList.size());
+ ArrayList<CruxVertex> pvList = new ArrayList<CruxVertex>(mcList.size());
+
+ for (CruxParticle fsPart : fsList) {
+ CruxParticle rp = fsPart;
+ while (rp != null) rp = addParent(rp, pList, pvList);
+ }
+
+ // Create CruxParticles for all decendents of FS MCParticles
+
+ ArrayList<CruxParticle> dList = new ArrayList<CruxParticle>(mcList.size());
+ ArrayList<CruxVertex> dvList = new ArrayList<CruxVertex>(mcList.size());
+ for (CruxParticle fsPart : fsList) {
+ addDaughters(fsPart);
+ }
+
+ log(" ");
+ log("Reco List before simplifying chains");
+ if (isLogging(FINEST)) {
+ _recoList.print();
+ }
+
+ // Replace linear chains by single particles
+
+ List<CruxVertex> verList = _recoList.getOrphanVertexes();
+ log(" ");
+ log("Particles " + _recoList.size() + " Orphan vertexes " + verList.size());
+ for (CruxVertex v : verList) {
+ List<CruxParticle> parList = new ArrayList<CruxParticle>(v.getDaughters());
+ for (CruxParticle p : parList) {
+ simplifyChains(p, null);
+ }
+ }
+
+ log(" ");
+ log("Reco List after simplifying chains");
+ if (isLogging(FINEST)) {
+ _recoList.print();
+ }
+
+ // Compact the particle list
+
+ _recoList.trimToSize();
+ _recoList = null;
+ log("Finished");
+ }
+
+// -- Helper methods : --------------------------------------------------------
+
+ /**
+ * Add all decendents of the given particle.
+ */
+ private void addDaughters(CruxParticle rp) {
+ MCParticle mc = rp.getMCParticle();
+ CruxVertex endVertex = null;
+ for (MCParticle dmc : mc.getDaughters()) {
+ CruxParticle drp = new CruxParticle();
+ _recoList.add(drp);
+ drp.setMCParticle(dmc);
+ if (dmc.getSimulatorStatus().vertexIsNotEndpointOfParent()) {
+ rp.getDirectDaughters().add(drp);
+ } else {
+ if (endVertex == null) {
+ endVertex = new CruxVertex();
+ endVertex.setParent(rp);
+ rp.setEndVertex(endVertex);
+ }
+ endVertex.getDaughters().add(drp);
+ drp.setStartVertex(endVertex);
+ }
+ addDaughters(drp);
+ }
+ }
+
+ /**
+ * Add parent of the given particle.
+ */
+ private CruxParticle addParent(CruxParticle rp, ArrayList<CruxParticle> pList, ArrayList<CruxVertex> vList) {
+
+ MCParticle mc = rp.getMCParticle();
+ List<MCParticle> parents = mc.getParents();
+
+ // No parents - create a new vertex and return
+ if (parents == null || parents.isEmpty()) {
+ CruxVertex vertex = new CruxVertex();
+ rp.setStartVertex(vertex);
+ vertex.getDaughters().add(rp);
+ vList.add(vertex);
+ return null;
+ }
+
+ // Multiple parents - go to previous generation until converge to a single parent or reach the top of particle tree
+ while (parents.size() > 1) {
+ Set<MCParticle> temp = new HashSet<MCParticle>();
+ for (MCParticle mcp : parents) {
+ List<MCParticle> mcpParents = mcp.getParents();
+ if (mcpParents != null) temp.addAll(mcpParents);
+ }
+
+ if (temp.isEmpty()) { // Reached the top of the tree
+ CruxVertex vertex = null;
+ for (CruxVertex v : vList) {
+ if (parents.containsAll(v.getMCParticles())) {
+ vertex = v;
+ break;
+ }
+ }
+ if (vertex == null) {
+ vertex = new CruxVertex();
+ ArrayList<MCParticle> mcList = new ArrayList<MCParticle>(parents);
+ vertex.setMCParticles(mcList);
+ vList.add(vertex);
+ }
+ rp.setStartVertex(vertex);
+ vertex.getDaughters().add(rp);
+ return null;
+ }
+
+ parents = new ArrayList<MCParticle>(temp);
+ }
+
+ // Found single parent
+ MCParticle pmc = parents.get(0);
+ CruxParticle prp = null;
+ for (CruxParticle p : pList) {
+ if (p.getMCParticle() == pmc) {
+ prp = p;
+ break;
+ }
+ }
+ CruxParticle out = null;
+ if (prp == null) {
+ prp = new CruxParticle();
+ prp.setMCParticle(pmc);
+ pList.add(prp);
+ _recoList.add(prp);
+ out = prp;
+ }
+ rp.setParent(prp);
+ if (mc.getSimulatorStatus().vertexIsNotEndpointOfParent()) {
+ prp.getDirectDaughters().add(rp);
+ } else {
+ CruxVertex vertex = prp.getEndVertex();
+ if (vertex == null) {
+ vertex = new CruxVertex();
+ vertex.setMCParticles(new ArrayList<MCParticle>(parents));
+ vertex.setParent(prp);
+ prp.setEndVertex(vertex);
+ }
+ vertex.getDaughters().add(rp);
+ rp.setStartVertex(vertex);
+ }
+
+ return out;
+ }
+
+ private void simplifyChains(CruxParticle particle, List<CruxParticle> chain) {
+ List<CruxParticle> daughters = particle.getDaughters();
+ int nAll = daughters.size();
+ int nDir = particle.getDirectDaughters().size();
+ if (nAll == 1 && nDir == 0) {
+ if (chain == null) chain = new ArrayList<CruxParticle>(2);
+ chain.add(particle);
+ simplifyChains(daughters.get(0), chain);
+ } else {
+ if (chain != null) {
+ if (nDir == 0) chain.add(particle);
+ if (chain.size() > 1) reduceChain(chain);
+ }
+ for (CruxParticle p : daughters) {
+ simplifyChains(p, null);
+ }
+ }
+ }
+
+ private void reduceChain(List<CruxParticle> chain) {
+ log("--- Chain Before (DOC: "+ MCParticle.DOCUMENTATION +" INT: "+ MCParticle.INTERMEDIATE +") --");
+ for (CruxParticle p : chain) {
+ log(" "+ p.getMCParticle().getType().getName() +" "+ p.getMCParticle().getPDGID() +
+ " "+ p.getMCParticle().getType() +" "+ p.getMCParticle().getGeneratorStatus());
+ }
+ while (chain.size() > 1) {
+ CruxParticle candidate = null;
+ int rating = 0;
+ for (CruxParticle p : chain) {
+ int status = p.getMCParticle().getGeneratorStatus();
+ if (status != MCParticle.FINAL_STATE && status != 0) {
+ double charge = Double.NaN;
+ try {
+ charge = p.getMCParticle().getType().getCharge();
+ } catch (UnknownParticleIDException x) {
+ candidate = p;
+ break;
+ }
+ if (charge == Double.NaN) {
+ candidate = p;
+ break;
+ } else {
+ charge = Math.abs(charge);
+ if (charge > .1 && charge < .9) {
+ if (rating < 4) {
+ rating = 4;
+ candidate = p;
+ }
+ } else if (status == MCParticle.DOCUMENTATION) {
+ if (rating < 3) {
+ rating = 3;
+ candidate = p;
+ }
+ } else if (status == MCParticle.INTERMEDIATE) {
+ if (rating < 2) {
+ rating = 2;
+ candidate = p;
+ }
+ }
+ }
+ }
+ }
+ if (candidate == null) {
+ break;
+ } else {
+ int i = chain.indexOf(candidate);
+ chain.remove(i);
+ List<MCParticle> mcList = candidate.getMCParticles();
+ CruxParticle part = (chain.size() > i) ? chain.get(i) : chain.get(i-1);
+ part.getMCParticles().addAll(mcList);
+ _recoList.remove(candidate);
+ CruxVertex start = candidate.getStartVertex();
+ CruxVertex end = candidate.getEndVertex();
+ start.getDaughters().remove(candidate);
+ if (end != null) {
+ start.getDaughters().addAll(end.getDaughters());
+ }
+ }
+ }
+ log("--- Chain After ");
+ for (CruxParticle p : chain) {
+ log(" "+ p.getMCParticle().getType().getName() +" "+ p.getMCParticle().getPDGID() +
+ " "+ p.getMCParticle().getType() +" "+ p.getMCParticle().getGeneratorStatus());
+ }
+ }
+
+// -- Private parts : ---------------------------------------------------------
+
+ protected String _recoListName = "CruxParticles";
+
+ CruxParticleList _recoList;
+}
lcsim-contrib/src/main/java/org/lcsim/contrib/onoprien/crux/cheat
diff -N CheatTrackFinderDriver.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ CheatTrackFinderDriver.java 22 Jan 2009 21:12:42 -0000 1.1
@@ -0,0 +1,277 @@
+package org.lcsim.contrib.onoprien.crux.cheat;
+
+import java.util.*;
+
+import hep.physics.vec.Hep3Vector;
+import org.lcsim.geometry.Subdetector;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+import org.lcsim.event.SimTrackerHit;
+
+import org.lcsim.contrib.onoprien.util.NoSuchParameterException;
+import org.lcsim.contrib.onoprien.util.ConstHep3Vector;
+import org.lcsim.contrib.onoprien.util.job.Driver;
+import org.lcsim.contrib.onoprien.util.job.JobEvent;
+import org.lcsim.contrib.onoprien.util.job.JobEventListener;
+import org.lcsim.contrib.onoprien.util.job.JobManager;
+import org.lcsim.contrib.onoprien.util.swim.BField;
+import org.lcsim.contrib.onoprien.util.swim.Helix;
+
+import org.lcsim.contrib.onoprien.vsegment.hit.ITrackerHit;
+import org.lcsim.contrib.onoprien.vsegment.mctruth.MCTruthVS;
+
+import org.lcsim.contrib.onoprien.crux.infrastructure.*;
+
+/**
+ * Driver that uses MC truth info to produce a list of tracks.
+ *
+ * @author D. Onoprienko
+ * @version $Id: CheatTrackFinderDriver.java,v 1.1 2009/01/22 21:12:42 onoprien Exp $
+ */
+public class CheatTrackFinderDriver extends Driver implements JobEventListener {
+
+// -- Constructors and initialization : ---------------------------------------
+
+ public CheatTrackFinderDriver() {
+ JobManager.defaultInstance().addListener(this);
+ }
+
+ /**
+ * Called by the framework to perform detector-dependent initialization.
+ * @throws IllegalStateException if some of the required parameters have not been set.
+ */
+ public void detectorChanged(JobEvent jEvent) {
+ _bField = JobManager.defaultInstance().get(BField.class);
+ if (_hitMapNames == null || _trackListName == null || _def == null) {
+ throw new IllegalStateException(ERR_NS);
+ }
+ }
+
+// -- Setters : ---------------------------------------------------------------
+
+ /**
+ * Set any parameter.
+ * The following parameters can be set with this method:
+ * <p><dl>
+ * <dt>"HIT_COLLECTIONS"</dt> <dd>One or more names of input tracker hit collections.
+ * <tt>ITrackerHits</tt> from these collections will be attached to created tracks.
+ * <br>No default, must be set explicitly.</dd>
+ * <dt>"TRACK_LIST"</dt> <dd>Name of output list of <tt>Tracks</tt>.
+ * <br>No default, must be set explicitly.</dd>
+ * <dt>"DEFINITION"</dt> <dd>The value is an instance of {@link ICheatRecoDefinition} that will
+ * be used to determine which MCParticles produce reconstructable tracks.
+ * <br>No default, must be either set explicitly or configured after getting.</dd></dl>
+ *
+ * @param name Name of parameter to be set. Case is ignored.
+ * @param values List of values to be used for setting the parameter.
+ * @throws NoSuchParameterException Thrown if the supplied parameter name is unknown.
+ * @throws IllegalArgumentException Thrown if incorrect number of values, or a value
+ * of incorrect type is supplied.
+ */
+ public void set(String name, Object... values) {
+ Object value = (values.length == 0) ? null : values[0];
+ try {
+ if (name.equalsIgnoreCase("HIT_COLLECTIONS")) {
+ _hitMapNames = new ArrayList<String>(values.length);
+ for (Object o : values) _hitMapNames.add((String)o);
+ } else if (name.equalsIgnoreCase("TRACK_LIST")) {
+ _trackListName = (String) value;
+ } else if (name.equalsIgnoreCase("DEFINITION")) {
+ _def = (ICheatRecoDefinition) value;
+ } else {
+ super.set(name, values);
+ }
+ } catch (ClassCastException x) {
+ throw new IllegalArgumentException(ERR_VIT, x);
+ }
+ }
+
+// -- Getters : ---------------------------------------------------------------
+
+ /**
+ * Returns reconstructable particle definition used by this driver.
+ * If the definition is not an instance of <tt>CheatRecoDefinition</tt>, returns <tt>null</tt>.
+ */
+ public CheatRecoDefinition getDefinition() {
+ if (_def == null) {
+ _def = new CheatRecoDefinition() {
+ public int getNumTrackerSuperLayers(MCParticle mc) {
+ List<ITrackerHit> hits = _mc2hit.get(mc);
+ int nLayers = 0;
+ Subdetector sd = null;
+ int layer = -1;
+ for (ITrackerHit hit : hits) {
+ Subdetector sd1 = hit.getSubdetector();
+ int layer1 = hit.getLayer();
+ if (sd != sd1 || layer != layer1) {
+ sd = sd1;
+ layer = layer1;
+ if (hit.getType().nUnmeasured() == 0) {
+ nLayers +=2;
+ } else {
+ nLayers +=1;
+ }
+ }
+ }
+ nLayers /= 2;
+ return nLayers;
+ }
+ };
+ }
+ return (_def instanceof CheatRecoDefinition) ? (CheatRecoDefinition)_def : null ;
+ }
+
+// -- Event processing : ------------------------------------------------------
+
+ /** Called by the framework to process an event. */
+ public void process(EventHeader event) {
+
+ super.process(event);
+
+ // Fetch MCTruth and list of MCParticles
+
+ _mcTruth = (MCTruthVS) event.get(MCTruthVS.KEY);
+
+ // Create list of tracks and save it
+
+ ArrayList<CruxTrack> trackList = new ArrayList<CruxTrack>(100);
+ event.put(_trackListName, trackList, CruxTrack.class, 0);
+
+ // Compile set of hits
+
+ HashSet<ITrackerHit> hitSet = new HashSet<ITrackerHit>();
+ for (String name : _hitMapNames) {
+ try {
+ Collection<ITrackerHit> hitMap = (Collection<ITrackerHit>) event.get(name);
+ hitSet.addAll(hitMap);
+ } catch (IllegalArgumentException x) {
+ continue;
+ }
+ }
+
+ // Associate hits to MCParticles :
+
+ _mc2hit = new HashMap<MCParticle,List<ITrackerHit>>();
+ for (ITrackerHit hit : hitSet) {
+ List<MCParticle> mcList = _mcTruth.getMCParticles(hit);
+ for (MCParticle mc : mcList) {
+ List<ITrackerHit> hits = _mc2hit.get(mc);
+ if (hits == null) {
+ hits = new ArrayList<ITrackerHit>();
+ _mc2hit.put(mc, hits);
+ }
+ hits.add(hit);
+ }
+ }
+ for (Map.Entry<MCParticle, List<ITrackerHit>> entry : _mc2hit.entrySet()) {
+ _mc = entry.getKey();
+ Collections.sort(entry.getValue(), _hitTimeComparator);
+ _mc = null;
+ }
+
+ // Loop over MCParticle, creating Tracks for those deemed reconstructable
+
+ List<MCParticle> mcList = event.getMCParticles();
+ for (MCParticle mc : mcList) {
+ if (_def.isTrackFindable(mc)) {
+ CruxTrack track = makeTrack(mc);
+ trackList.add(track);
+ }
+ }
+
+ // Cleanup
+
+ trackList.trimToSize();
+ _mcTruth = null;
+ _mc2hit = null;
+ }
+
+
+// -- Helper methods : --------------------------------------------------------
+
+
+ /**
+ * Create <tt>CruxTrack</tt> object for the given <tt>MCParticle</tt>, and set its parameters.
+ * Default implementation only sets <tt>MCParticle</tt> associated with the track -
+ * subclasses can override.
+ */
+ private CruxTrack makeTrack(MCParticle mcParticle) {
+
+ // Get hits produced by MCParticle
+
+ List<ITrackerHit> hits = _mc2hit.get(mcParticle);
+
+ // Wrap hits into anchors :
+
+ ArrayList<ITrackAnchor> anchors = new ArrayList<ITrackAnchor>(hits.size());
+ for (ITrackerHit hit : hits) {
+ anchors.add(new CruxTrackAnchorHit(hit));
+ }
+
+ // Create track node at the particle origin :
+
+ Hep3Vector pos = mcParticle.getOrigin();
+ Hep3Vector momentum = mcParticle.getMomentum();
+ int charge = (int) mcParticle.getCharge();
+ Helix helix = _bField.makeHelix(pos, momentum, charge);
+ CruxTrackNode startNode = new CruxTrackNode(helix);
+
+ ArrayList<ITrackNode> nodes = new ArrayList<ITrackNode>(2);
+ nodes.add(startNode);
+
+ // Create node at last hit position :
+
+ SimTrackerHit simHit = _mcTruth.getSimTrackerHits(hits.get(hits.size()-1)).get(0);
+ pos = new ConstHep3Vector(simHit.getPoint());
+ momentum = new ConstHep3Vector(simHit.getMomentum());
+ helix = _bField.makeHelix(pos, momentum, charge);
+ CruxTrackNode lastHitNode = new CruxTrackNode(helix);
+
+ nodes.add(lastHitNode);
+
+ // Create track
+
+ CruxTrack track = new CruxTrack(anchors, nodes, charge);
+ track.setMCParticle(mcParticle);
+ track.setFunctionalPoint(ITrack.Point.START, startNode, 0.);
+ track.setFunctionalPoint(ITrack.Point.LAST_HIT, lastHitNode, 0.);
+
+ return track;
+ }
+
+// -- Private parts : ---------------------------------------------------------
+
+ private ArrayList<String> _hitMapNames;
+ private String _trackListName;
+
+ private ICheatRecoDefinition _def;
+
+ private MCTruthVS _mcTruth;
+ private BField _bField;
+
+ private HashMap<MCParticle,List<ITrackerHit>> _mc2hit;
+
+ private Comparator<ITrackerHit> _hitTimeComparator = new Comparator<ITrackerHit>() {
+ public int compare(ITrackerHit hit1, ITrackerHit hit2) {
+ SimTrackerHit sHit1 = null;
+ SimTrackerHit sHit2 = null;
+ List<SimTrackerHit> simList = _mcTruth.getSimTrackerHits(hit1);
+ for (SimTrackerHit sHit : simList) {
+ if (sHit.getMCParticle() == _mc) {
+ sHit1 = sHit;
+ break;
+ }
+ }
+ simList = _mcTruth.getSimTrackerHits(hit2);
+ for (SimTrackerHit sHit : simList) {
+ if (sHit.getMCParticle() == _mc) {
+ sHit2 = sHit;
+ break;
+ }
+ }
+ return (int) Math.signum(sHit1.getTime() - sHit2.getTime());
+ }
+ };
+ private MCParticle _mc;
+
+}
lcsim-contrib/src/main/java/org/lcsim/contrib/onoprien/crux/cheat
diff -N CheatVertexFinderDriver.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ CheatVertexFinderDriver.java 22 Jan 2009 21:12:42 -0000 1.1
@@ -0,0 +1,192 @@
+package org.lcsim.contrib.onoprien.crux.cheat;
+
+import java.util.*;
+
+import hep.physics.vec.Hep3Vector;
+import hep.physics.vec.VecOp;
+import org.lcsim.event.EventHeader;
+import org.lcsim.event.MCParticle;
+import org.lcsim.units.clhep.SystemOfUnits;
+
+import org.lcsim.contrib.onoprien.util.NoSuchParameterException;
+import org.lcsim.contrib.onoprien.util.job.Driver;
+
+import org.lcsim.contrib.onoprien.crux.infrastructure.*;
+import org.lcsim.contrib.onoprien.crux.mctruth.MCTruthCrux;
+
+/**
+ * Driver that uses MC truth info to produce a list of reconstructable vertices in the
+ * event, and to associate them with previously reconstructed tracks.
+ * <p>
+ * A vertex is created if more than one track from the given list originates, terminates,
+ * or flies by it.
+ *
+ * @author D. Onoprienko
+ * @version $Id: CheatVertexFinderDriver.java,v 1.1 2009/01/22 21:12:42 onoprien Exp $
+ */
+public class CheatVertexFinderDriver extends Driver {
+
+// -- Constructors : ----------------------------------------------------------
+
+ public CheatVertexFinderDriver() {
+ }
+
+// -- Setters : ---------------------------------------------------------------
+
+ /**
+ * Set any parameter.
+ * The following parameters can be set with this method:
+ * <p><dl>
+ * <dt>"TRACK_LIST"</dt> <dd>Name of input list of <tt>Tracks</tt>
+ * <br>No default, must be set explicitly if this class is to be used as a driver.</dd>
+ * <dt>"VERTEX_LIST"</dt> <dd>Name of output list of found <tt>Vertexes</tt>
+ * <br>Default: <tt>null</tt> (vertices will be associated with tracks but the list will not
+ * be explicitly saved into the event record.</dd>
+ * <dt>"MERGE_DISTANCE"</dt> <dd>Double value defining the minimal distance between vertices -
+ * closer vertices will be merged.
+ * <br>Default: 1 micrometer.</dd></dl>
+ *
+ * @param name Name of parameter to be set. Case is ignored.
+ * @param values List of values to be used for setting the parameter.
+ * @throws NoSuchParameterException Thrown if the supplied parameter name is unknown.
+ * @throws IllegalArgumentException Thrown if incorrect nunber of values, or a value
+ * of incorrect type is supplied.
+ */
+ public void set(String name, Object... values) {
+ Object value = (values.length == 0) ? null : values[0];
+ try {
+ if (name.equalsIgnoreCase("TRACK_LIST")) {
+ _trackListName = (String) value;
+ } else if (name.equalsIgnoreCase("VERTEX_LIST")) {
+ _vertexListName = (String) value;
+ } else if (name.equalsIgnoreCase("MERGE_DISTANCE")) {
+ _mergeDist2 = (Double) value;
+ _mergeDist2 *= _mergeDist2;
+ } else {
+ super.set(name, values);
+ }
+ } catch (ClassCastException x) {
+ throw new IllegalArgumentException("Value of incompatible type", x);
+ }
+ }
+
+// -- Event processing : ------------------------------------------------------
+
+ /** Called by the framework to process an event. */
+ public void process(EventHeader event) {
+
+ super.process(event);
+
+ _mcTruth = (MCTruthCrux) event.get("CruxMCTruth");
+ List<MCParticle> mcList = event.getMCParticles();
+ List<CruxTrack> trackList = event.get(CruxTrack.class, _trackListName);
+
+ ArrayList<CruxTrackVertex> vertexList = new ArrayList<CruxTrackVertex>();
+
+ // Create vertexes for start and end points of each track, merging those less than _mergeDist apart
+
+ for (CruxTrack track : trackList) {
+ MCParticle mc = track.getMCParticle();
+ Hep3Vector endPoint = mc.getEndPoint();
+ Hep3Vector startPoint = mc.getOrigin();
+ boolean start = true;
+ boolean end = true;
+ for (CruxTrackVertex vertex : vertexList) {
+ Hep3Vector pos = vertex.getPosition();
+ if ( start && (VecOp.sub(pos, startPoint).magnitudeSquared() < _mergeDist2) ) {
+ vertex.getDaughterTracks().add(track);
+ start = false;
+ } else if ( end && (VecOp.sub(pos, endPoint).magnitudeSquared() < _mergeDist2) ) {
+ vertex.setParentTrack(track);
+ end = false;
+ }
+ }
+ if (start) {
+ CruxTrackVertex vertex = new CruxTrackVertex();
+ vertex.setPosition(startPoint);
+ vertex.getDaughterTracks().add(track);
+ vertexList.add(vertex);
+ }
+ if (end) {
+ CruxTrackVertex vertex = new CruxTrackVertex();
+ vertex.setPosition(endPoint);
+ vertex.setParentTrack(track);
+ vertexList.add(vertex);
+ }
+ }
+
+ // Save vertexes with more than one track
+
+ ArrayList<CruxTrackVertex> finalList = new ArrayList<CruxTrackVertex>();
+ for (CruxTrackVertex vertex : vertexList) {
+ boolean isIntermediate = false;
+ List<CruxTrack> daughters = vertex.getDaughterTracks();
+ int nTracks = daughters.size();
+ if (vertex.getParentTrack() == null) {
+ for (CruxTrack track : daughters) {
+ MCParticle mc = track.getMCParticle();
+ CruxTrack parentTrack = getFlyByParent(mc);
+ if (parentTrack != null) {
+ vertex.setParentTrack(parentTrack);
+ isIntermediate = true;
+ nTracks++;
+ break;
+ }
+ }
+ } else {
+ nTracks++;
+ }
+ if (nTracks > 1) {
+ finalList.add(vertex);
+ for (CruxTrack track : daughters) {
+ track.setStartVertex(vertex);
+ }
+ CruxTrack parent = vertex.getParentTrack();
+ if (parent != null) {
+ if (isIntermediate) {
+ parent.getIntermediateVertexes().add(vertex);
+ } else {
+ parent.setEndVertex(vertex);
+ }
+ }
+ }
+ }
+
+ // Put list of vertexes into event record
+
+ if (_vertexListName != null) {
+ finalList.trimToSize();
+ event.put(_vertexListName, finalList, CruxTrackVertex.class, 0);
+ }
+
+ // Clean up
+
+ _mcTruth = null;
+ }
+
+ /**
+ * Try to find a track whose intermediate vertex the given MCParticle was produced at.
+ */
+ private CruxTrack getFlyByParent(MCParticle mcParticle) {
+ boolean bornInFlight = mcParticle.getSimulatorStatus().vertexIsNotEndpointOfParent();
+ List<MCParticle> pList = mcParticle.getParents();
+ for (MCParticle mc : pList) {
+ CruxTrack track = _mcTruth.getTrack(mc);
+ if (bornInFlight && track != null) {
+ return track;
+ } else if (VecOp.sub(mc.getOrigin(), mc.getEndPoint()).magnitudeSquared() < _mergeDist2) {
+ track = getFlyByParent(mc);
+ if (track != null) return track;
+ }
+ }
+ return null;
+ }
+
+// -- Private parts : ---------------------------------------------------------
+
+ protected String _trackListName;
+ protected String _vertexListName;
+ protected double _mergeDist2 = (1.*SystemOfUnits.micrometer) * (1.*SystemOfUnits.micrometer);
+
+ protected MCTruthCrux _mcTruth;
+}
lcsim-contrib/src/main/java/org/lcsim/contrib/onoprien/crux/cheat
diff -N ICheatRecoDefinition.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ ICheatRecoDefinition.java 22 Jan 2009 21:12:42 -0000 1.1
@@ -0,0 +1,25 @@
+package org.lcsim.contrib.onoprien.crux.cheat;
+
+import org.lcsim.event.MCParticle;
+
+/**
+ * Interface to be implemented by classes that define reconstructable MCParticle.
+ * Used by cheaters and performance analysis classes.
+ *
+ * @author D. Onoprienko
+ * @version $Id: ICheatRecoDefinition.java,v 1.1 2009/01/22 21:12:42 onoprien Exp $
+ */
+public interface ICheatRecoDefinition {
+
+ /**
+ * Returns the name of this definition.
+ */
+ String getName();
+
+ /**
+ * Returns <tt>true</tt> if the specified <tt>MCParticle</tt> produced a findable
+ * track, according to this definition.
+ */
+ boolean isTrackFindable(MCParticle mcParticle);
+
+}
CVSspam 0.2.8