8 modified files
lcsim/src/org/lcsim/contrib/CosminDeaconu/EfficiencyPlotter
diff -u -r1.1 -r1.2
--- EfficiencyPlotManager.java 15 Nov 2008 02:54:25 -0000 1.1
+++ EfficiencyPlotManager.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -6,6 +6,7 @@
package org.lcsim.contrib.CosminDeaconu.EfficiencyPlotter;
import java.util.ArrayList;
+import java.util.Arrays;
import java.util.Collection;
import java.util.Collections;
import org.lcsim.contrib.CosminDeaconu.*;
@@ -33,7 +34,8 @@
public EfficiencyPlotManager(AIDA instance){
this.instance = instance;
}
-
+
+ //WARNING: remembering entries uses a lot of memory...
public EfficiencyPlotManager(AIDA instance, boolean rememberEntries) {
this.instance = instance;
remember = rememberEntries;
@@ -59,6 +61,18 @@
defaultCuts = Collections.EMPTY_LIST;
}
+
+ public void addPlots(EfficiencyPlot[] plots) {
+ addPlots(Arrays.asList(plots));
+ }
+
+ public void addPlots(Collection<EfficiencyPlot> plots) {
+ for (EfficiencyPlot plot : plots) {
+ addPlot(plot);
+ }
+
+ }
+
public void addPlot(EfficiencyPlot p) {
p.setupHistogram(instance);
plots.add(p);
lcsim/src/org/lcsim/contrib/CosminDeaconu/EfficiencyPlotter
diff -u -r1.1 -r1.2
--- EfficiencyPlotPresets.java 15 Nov 2008 02:54:25 -0000 1.1
+++ EfficiencyPlotPresets.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -20,10 +20,27 @@
public class EfficiencyPlotPresets {
+ public static int[] PARTICLE_TYPES = new int[] {
+ 11,//"e-",
+ -11,//"e+",
+ 2212,//"p",
+ -2212,//pbar
+ 211,//"pi+",
+ -211,//"pi-",
+ 321,//"K+",
+ -321,//"K-",
+ -13,//"u+",
+ 13//"u-"
+ };
+
+ public static double PT_CUT = 0.2;
+ public static double COS_THETA_CUT=0.985;
+ public static double MAX_RADIUS = 10;
+
private static List<EfficiencyCut> commonCuts = Arrays.asList(new EfficiencyCut[] {
new SimpleEfficiencyCut(EfficiencyCutOperator.NOTEQUALS, EfficiencyParameter.CHARGE, 0.0),
new SimpleEfficiencyCut(EfficiencyCutOperator.EQUALS, EfficiencyParameter.GENERATOR_STATUS, Particle.FINAL_STATE),
- new SimpleEfficiencyCut(EfficiencyCutOperator.MAXIMUM, EfficiencyParameter.RADIUS, 10.0)
+ new SimpleEfficiencyCut(EfficiencyCutOperator.MAXIMUM, EfficiencyParameter.RADIUS, MAX_RADIUS)
});
private static List<EfficiencyCut> totalCuts = new ArrayList<EfficiencyCut>();
@@ -34,9 +51,9 @@
private static void setup() {
- thetaCuts.add(new SimpleEfficiencyCut(EfficiencyCutOperator.MINIMUM, EfficiencyParameter.PT, 1.1));
+ thetaCuts.add(new SimpleEfficiencyCut(EfficiencyCutOperator.MINIMUM, EfficiencyParameter.PT,PT_CUT));
thetaCuts.addAll(commonCuts);
- ptCuts.add(new SimpleEfficiencyCut(EfficiencyCutOperator.MAXIMUM, EfficiencyParameter.COS_THETA, 0.985));
+ ptCuts.add(new SimpleEfficiencyCut(EfficiencyCutOperator.MAXIMUM, EfficiencyParameter.COS_THETA, COS_THETA_CUT));
ptCuts.addAll(commonCuts);
totalCuts.addAll(thetaCuts);
totalCuts.addAll(ptCuts);
@@ -49,12 +66,46 @@
new FindableTrackCut(strategy));
}
+ public static List<EfficiencyPlot> getEfficiencyVsParticlePlots() {
+ if (!ready) setup();
+ List<EfficiencyPlot> plotlist = new ArrayList<EfficiencyPlot>();
+ for (int i : PARTICLE_TYPES) {
+ List<EfficiencyCut> theseCuts = new ArrayList<EfficiencyCut>();
+ theseCuts.addAll(totalCuts);
+ theseCuts.add(new SimpleEfficiencyCut(EfficiencyCutOperator.EQUALS, EfficiencyParameter.PARTICLE_ID, i));
+ plotlist.add(new EfficiencyPlot(null, "Eff vs. Particle "+i, 20, 0, 1.1, theseCuts));
+
+ }
+
+ return plotlist;
+
+ }
+
public static EfficiencyPlot getEfficiencyVsPTPlotUsingFindable(List<SeedStrategy> strategy) {
return new EfficiencyPlot(EfficiencyParameter.PT, "Findable: Eff vs. Pt", 100, 0, 50,
new FindableTrackCut(strategy, Ignore.NoPTCut));
}
-
+ public static EfficiencyPlot getEfficiencyPlotUsingFindable(List<SeedStrategy> strategy) {
+ return new EfficiencyPlot(null, "Findable: Efficiency",10,0,1.2, new FindableTrackCut(strategy));
+ }
+
+ public static EfficiencyPlot getEfficiencyVsNumHitsPlot(){
+ if(!ready) setup();
+ return new EfficiencyPlot(EfficiencyParameter.NUM_HITS, "EPM: Eff vs. NumHits", 15, 0, 15, totalCuts);
+ }
+
+ public static EfficiencyPlot getEfficiencyVsCharge(){
+ if(!ready) setup();
+ return new EfficiencyPlot(EfficiencyParameter.CHARGE, "EPM: Eff vs. Charge", 10, -1.2,1.2, totalCuts);
+ }
+
+
+ public static ParameterPlot getMCNumHitsPlot(){
+ if(!ready) setup();
+ return new ParameterPlot(EfficiencyParameter.NUM_HITS, "EPM: MC NumHits", 15, 0, 15, totalCuts);
+ }
+
public static EfficiencyPlot getEfficiencyVsThetaPlot(){
if(!ready) setup();
return new EfficiencyPlot(EfficiencyParameter.THETA_DEGREES, "EPM: Eff vs. Theta", 180, 0, 180, thetaCuts);
lcsim/src/org/lcsim/contrib/CosminDeaconu/EfficiencyPlotter
diff -u -r1.1 -r1.2
--- EfficiencyParameter.java 15 Nov 2008 02:54:25 -0000 1.1
+++ EfficiencyParameter.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -24,5 +24,6 @@
MASS,
PRODUCTION_TIME,
GENERATOR_STATUS,
- RADIUS;
+ RADIUS,
+ NUM_HITS;
}
lcsim/src/org/lcsim/contrib/CosminDeaconu/EfficiencyPlotter
diff -u -r1.1 -r1.2
--- EfficiencyPlot.java 15 Nov 2008 02:54:25 -0000 1.1
+++ EfficiencyPlot.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -63,25 +63,26 @@
if (var!=null) {
hist = instance.histogramFactory().createHistogram1D(name,"",bins,min,max,"type=efficiency");
- hist.setTitle("EFFICIENCY vs. "+var.toString());
+// hist.setTitle("EFFICIENCY vs. "+var.toString());
} else {
hist = instance.histogramFactory().createHistogram1D(name,bins,min,max);
- hist.setTitle("EFFICIENCY");
+// hist.setTitle("EFFICIENCY");
}
-
ready = true;
-
}
public IHistogram1D getHistogram(){
if(!ready) {
throw new RuntimeException("Efficiency Plot uninitialized. Must be initialized by calling setupHistogram() before you can call getHistogram().");
}
+
return hist;
}
+
+
public void addEntry(EfficiencyPlotEntry e) {
if(!ready) {
lcsim/src/org/lcsim/contrib/CosminDeaconu/EfficiencyPlotter
diff -u -r1.1 -r1.2
--- EfficiencyPlotEntry.java 15 Nov 2008 02:54:25 -0000 1.1
+++ EfficiencyPlotEntry.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -1,13 +1,21 @@
-/*
+/*
* To change this template, choose Tools | Templates
* and open the template in the editor.
*/
package org.lcsim.contrib.CosminDeaconu.EfficiencyPlotter;
+import java.util.ArrayList;
+import java.util.HashMap;
+import java.util.HashSet;
+import java.util.List;
+import java.util.Map;
+import java.util.Set;
import org.lcsim.event.EventHeader;
import org.lcsim.event.MCParticle;
+import org.lcsim.event.SimTrackerHit;
import org.lcsim.fit.helicaltrack.HelixParamCalculator;
+import org.lcsim.fit.helicaltrack.HitIdentifier;
/**
*
@@ -23,6 +31,11 @@
//Cached calculations for stuff not in calc
private double cth;
+ //allow number of hits calculations
+ private static Map<MCParticle, Integer> hitmap = null;
+ private static EventHeader hitmapEvent = null;
+ private static HitIdentifier id = null;
+
public EfficiencyPlotEntry(MCParticle mcp, boolean found, EventHeader event) {
this.mcp = mcp;
this.found = found;
@@ -61,7 +74,50 @@
case GENERATOR_STATUS: return mcp.getGeneratorStatus();
case RADIUS: return mcp.getOrigin().magnitude();
case CHARGE: return mcp.getCharge();
+ case NUM_HITS: return getNumHits(mcp);
default: throw new UnsupportedOperationException();
}
}
+
+ public static Map<MCParticle,Integer> generateHitNumberMap(EventHeader event) {
+
+ if (id==null) id = new HitIdentifier(HitIdentifier.getSpecialLayers(event.getDetectorName()));
+
+ Map<MCParticle,Set<String>> identmap = new HashMap<MCParticle,Set<String>>();
+
+ //get SimHits
+ List<List<SimTrackerHit>> simhitlistlist = event.get(SimTrackerHit.class);
+ List<SimTrackerHit> simhits = new ArrayList<SimTrackerHit>();
+ for (List<SimTrackerHit> l : simhitlistlist) {
+ simhits.addAll(l);
+ }
+
+ for (SimTrackerHit h : simhits) {
+ MCParticle p = h.getMCParticle();
+ if (!identmap.containsKey(p)) {
+ identmap.put(p, new HashSet<String>());
+ }
+ identmap.get(p).add(id.Identifier(h.getDetectorElement()));
+ }
+
+ Map<MCParticle,Integer> returnme = new HashMap<MCParticle,Integer>();
+ for (MCParticle p : identmap.keySet()) {
+ returnme.put(p, identmap.get(p).size());
+ }
+
+ return returnme;
+ }
+
+ private int getNumHits(MCParticle mcp) {
+ if (getEvent()!=hitmapEvent) {
+ hitmapEvent = getEvent();
+ hitmap = generateHitNumberMap(getEvent());
+ }
+ if (!hitmap.containsKey(mcp)) {
+ // System.out.println("Can't Find MCParticle w/ p="+mcp.getMomentum().magnitude());
+ return 0;
+ }
+ return hitmap.get(mcp);
+ }
+
}
\ No newline at end of file
lcsim/src/org/lcsim/fit/helicaltrack
diff -u -r1.2 -r1.3
--- HitIdentifier.java 10 Nov 2008 19:46:01 -0000 1.2
+++ HitIdentifier.java 5 Dec 2008 20:25:30 -0000 1.3
@@ -30,7 +30,7 @@
}
public HitIdentifier(List<String> special) {
- super();
+ this();
_special.addAll(special);
}
@@ -100,4 +100,18 @@
_special.add(name);
return;
}
+
+ public static List<String> getSpecialLayers(String detectorName) {
+
+ List<String> returnMe = new ArrayList<String>();
+ if (detectorName.equals("sid02")) {
+ returnMe.add("TrackerEndcap");
+ } else if (detectorName.equals("sid01")) {
+ returnMe.add("TrackerForward");
+ returnMe.add("TrackerEndcap");
+ }
+
+ return returnMe;
+ }
+
}
\ No newline at end of file
lcsim/src/org/lcsim/contrib/CosminDeaconu/TrackingAnalysis
diff -u -r1.1 -r1.2
--- HelicalTrackHitReconstructionDriver.java 15 Nov 2008 02:54:25 -0000 1.1
+++ HelicalTrackHitReconstructionDriver.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -12,7 +12,9 @@
import java.util.HashMap;
import java.util.List;
import java.util.Map;
+import org.lcsim.digisim.MyLCRelation;
import org.lcsim.event.EventHeader;
+import org.lcsim.event.LCRelation;
import org.lcsim.event.Track;
import org.lcsim.event.TrackerHit;
import org.lcsim.recon.tracking.seedtracker.ReconTracking.SiD02ReconHybridHitMaker;
@@ -26,10 +28,10 @@
public class HelicalTrackHitReconstructionDriver extends Driver{
- private String remadeHitsName = "ReconstructedHelicalTrackHits";
- private String outputSeedCandidatesName = "ReconstructedSeedCandidates";
+ private static String remadeHitsName = "ReconstructedHelicalTrackHits";
+ private static String outputSeedCandidatesName = "ReconstructedSeedCandidates";
private String hitCollectionName = "HelicalTrackHits";
-
+ private static String outputTrackSeedCandidateRelationsName="SeedCandidateTrackMapping";
public HelicalTrackHitReconstructionDriver(){
@@ -72,8 +74,9 @@
List<SeedCandidate> seedCandidates = new ArrayList<SeedCandidate>();
Collections.sort(inputHits, hitcmp);
Collections.sort(outputHits, hitcmp);
-
- Map<TrackerHit, HelicalTrackHit> hitmap = new HashMap<TrackerHit,HelicalTrackHit>();
+
+ Map<TrackerHit, HelicalTrackHit> hitmap = new HashMap<TrackerHit,HelicalTrackHit>();
+ List<LCRelation> trackSeedCandidateRelations = new ArrayList<LCRelation>();
for (int i =0; i < inputHits.size(); i++) {
hitmap.put(inputHits.get(i), outputHits.get(i));
@@ -85,14 +88,24 @@
for (TrackerHit h : t.getTrackerHits()){
hits.add(hitmap.get(h));
}
- seedCandidates.add(new SeedCandidate(hits, null));
-
+ SeedCandidate sc = new SeedCandidate(hits, null);
+ seedCandidates.add(sc);
+ trackSeedCandidateRelations.add(new MyLCRelation(sc, t));
}
-
event.put(outputSeedCandidatesName, seedCandidates, SeedCandidate.class, 0);
+ event.put(outputTrackSeedCandidateRelationsName,trackSeedCandidateRelations, LCRelation.class, 0);
}
+ public static String getSeedCandidateName() {
+ return outputSeedCandidatesName;
+ }
+ public static String getTrackSeedCandidateRelationsName(){
+ return outputTrackSeedCandidateRelationsName;
+ }
+ public static String getReconstructedHelicalTrackHitsName(){
+ return remadeHitsName;
+ }
}
lcsim/src/org/lcsim/contrib/CosminDeaconu/TrackingAnalysis
diff -u -r1.1 -r1.2
--- ReconstructedAnalysisDriver.java 15 Nov 2008 02:54:25 -0000 1.1
+++ ReconstructedAnalysisDriver.java 5 Dec 2008 20:25:30 -0000 1.2
@@ -5,6 +5,8 @@
package org.lcsim.contrib.CosminDeaconu.TrackingAnalysis;
+import hep.physics.particle.properties.ParticleType;
+import java.util.Collection;
import java.util.HashMap;
import java.util.HashSet;
import java.util.List;
@@ -16,7 +18,7 @@
import org.lcsim.event.EventHeader;
import org.lcsim.event.MCParticle;
import org.lcsim.fit.helicaltrack.HelicalTrackHit;
-import org.lcsim.contrib.CosminDeaconu.TrackingAnalysis.HelicalTrackHitReconstructionDriver;
+import org.lcsim.fit.helicaltrack.HelixParamCalculator;
import org.lcsim.recon.tracking.seedtracker.SeedCandidate;
import org.lcsim.recon.tracking.seedtracker.SeedStrategy;
import org.lcsim.recon.tracking.seedtracker.StrategyXMLUtils;
@@ -31,12 +33,13 @@
private EfficiencyPlotManager epm;
-
+ AIDA aida = AIDA.defaultInstance();
+
public ReconstructedAnalysisDriver(){
add(new HelicalTrackHitReconstructionDriver());
- epm = new EfficiencyPlotManager(AIDA.defaultInstance());
+ epm = new EfficiencyPlotManager(aida);
// Get the list of strategies being used
@@ -46,22 +49,29 @@
epm.addPlot(EfficiencyPlotPresets.getEfficiencyVsPTPlotUsingFindable(slist));
epm.addPlot(EfficiencyPlotPresets.getEfficiencyVsThetaPlotUsingFindable(slist));
-
+ epm.addPlot(EfficiencyPlotPresets.getEfficiencyPlotUsingFindable(slist));
+ epm.addPlot(EfficiencyPlotPresets.getEfficiencyVsNumHitsPlot());
epm.addPlot(EfficiencyPlotPresets.getEfficiencyPlot());
epm.addPlot(EfficiencyPlotPresets.getEfficiencyVsPTPlot());
epm.addPlot(EfficiencyPlotPresets.getEfficiencyVsThetaPlot());
epm.addPlot(EfficiencyPlotPresets.getMCPTPlot());
epm.addPlot(EfficiencyPlotPresets.getMCThetaplot());
+ epm.addPlot(EfficiencyPlotPresets.getMCNumHitsPlot());
+ epm.addPlots(EfficiencyPlotPresets.getEfficiencyVsParticlePlots());
+ epm.addPlot(EfficiencyPlotPresets.getEfficiencyVsCharge());
}
+ @Override
protected void process(EventHeader event) {
super.process(event);
- List<SeedCandidate> candidates = event.get(SeedCandidate.class, "ReconstructedSeedCandidates");
+ List<SeedCandidate> candidates = event.get(SeedCandidate.class, HelicalTrackHitReconstructionDriver.getSeedCandidateName());
Set<MCParticle> found = new HashSet<MCParticle>();
+ Map<SeedCandidate,MCParticle> scmcmap = new HashMap<SeedCandidate,MCParticle>();
for (SeedCandidate c : candidates) {
+
Map<MCParticle,Integer> mcmap = new HashMap<MCParticle,Integer>();
- for (HelicalTrackHit h : c.getHits()){
+ for (HelicalTrackHit h : c.getHits()) {
for (MCParticle p : h.getMCParticles()){
Integer i = mcmap.get(p);
if (i==null) mcmap.put(p,1);
@@ -78,10 +88,74 @@
}
}
found.add(max);
+ scmcmap.put(c,max);
}
-
+ double foundEnergy = 0;
+ double totalEnergy = 0;
for (MCParticle mcp : event.getMCParticles()){
- epm.addEntry(new EfficiencyPlotEntry(mcp, found.contains(mcp), event));
+ boolean f = found.contains(mcp);
+ epm.addEntry(new EfficiencyPlotEntry(mcp, f, event));
+ if(f) foundEnergy+=mcp.getEnergy();
+ if (mcp.getGeneratorStatus()==mcp.FINAL_STATE) totalEnergy+=mcp.getEnergy();
}
+
+ //Calculate purity
+ Map<SeedCandidate,Double> purityMap = makePurityMap(candidates);
+
+ //Way to get Track from SeedCandidate... this is useful if we want to know about fitting parameters
+// List<LCRelation> scTrackRelations = event.get(LCRelation.class, HelicalTrackHitReconstructionDriver.getTrackSeedCandidateRelationsName());
+// RelationalTable scToTrack = new BaseRelationalTable(RelationalTable.Mode.ONE_TO_ONE, RelationalTable.Weighting.UNWEIGHTED);
+// for (LCRelation lcr : scTrackRelations) {
+// scToTrack.add(lcr.getFrom(), lcr.getTo());
+// }
+
+ for (SeedCandidate s: purityMap.keySet()){
+ double purity = purityMap.get(s);
+ aida.cloud1D("purity").fill(purity);
+ aida.cloud2D("purity vs. num hits").fill(purity, s.getHits().size());
+ MCParticle mcp = scmcmap.get(s);
+ HelixParamCalculator calc = new HelixParamCalculator(mcp,event);
+ aida.cloud2D("purity vs. theta").fill(purity, calc.getTheta());
+ }
+
+
+ //System.out.println(foundEnergy+"GeV found");
+ aida.cloud1D("foundEnergy").fill(foundEnergy/totalEnergy);
+
+ }
+
+ public static Map<SeedCandidate,Double> makePurityMap(Collection<SeedCandidate> coll) {
+
+ Map<SeedCandidate,Double> purityMap = new HashMap<SeedCandidate,Double>();
+ for (SeedCandidate sc : coll) {
+ Map<MCParticle,Integer> mcCounter = new HashMap<MCParticle,Integer>();
+ for (HelicalTrackHit h : sc.getHits()) {
+ for (MCParticle p : h.getMCParticles()) {
+ if (mcCounter.containsKey(p)) {
+ mcCounter.put(p,mcCounter.get(p)+1);
+ }
+ else mcCounter.put(p,1);
+ }
+ }
+ MCParticle maxParticle = null;
+ int max = 0;
+ for (MCParticle p : mcCounter.keySet()) {
+ if (mcCounter.get(p) > max) {
+ max = mcCounter.get(p);
+ maxParticle = p;
+ }
+ }
+ double total = sc.getHits().size();
+ int numFromMax = 0;
+ for (HelicalTrackHit h : sc.getHits()) {
+ if (h.getMCParticles().contains(maxParticle))
+ numFromMax++;
+ }
+ double purity = numFromMax/total;
+ purityMap.put(sc, purity);
+ }
+
+ return purityMap;
+
}
}
CVSspam 0.2.8