lcsim-contrib/src/main/java/org/lcsim/contrib/sATLAS
diff -u -r1.1 -r1.2
--- TrackAnalysisDriver.java 28 Apr 2009 21:18:08 -0000 1.1
+++ TrackAnalysisDriver.java 6 May 2009 01:00:17 -0000 1.2
@@ -4,6 +4,9 @@
*/
package org.lcsim.contrib.sATLAS;
+import java.io.IOException;
+import java.util.logging.Level;
+import java.util.logging.Logger;
import org.lcsim.contrib.mgraham.sATLASDigi.*;
import hep.aida.IHistogram1D;
import hep.aida.IHistogramFactory;
@@ -13,12 +16,17 @@
import hep.physics.vec.VecOp;
import java.util.ArrayList;
+import java.util.HashMap;
import java.util.List;
+import java.util.Map;
+
+
import java.util.Set;
//import org.lcsim.contrib.Partridge.TrackingTest.FindableTrack.Ignore;
//import org.lcsim.contrib.Partridge.TrackingTest.TrackAnalysis;
import org.lcsim.contrib.mgraham.sATLASDigi.FindableTrack.Ignore;
+import org.lcsim.detector.IDetectorElement;
import org.lcsim.event.EventHeader;
import org.lcsim.event.LCRelation;
import org.lcsim.event.MCParticle;
@@ -35,6 +43,8 @@
import org.lcsim.fit.helicaltrack.HelixUtils;
import org.lcsim.fit.helicaltrack.MultipleScatter;
import org.lcsim.fit.helicaltrack.TrackDirection;
+import org.lcsim.recon.tracking.digitization.sisim.SiTrackerHitPixel;
+import org.lcsim.recon.tracking.digitization.sisim.SiTrackerHitStrip1D;
import org.lcsim.recon.tracking.seedtracker.SeedCandidate;
import org.lcsim.recon.tracking.seedtracker.SeedStrategy;
import org.lcsim.recon.tracking.seedtracker.SeedTrack;
@@ -54,15 +64,25 @@
private IHistogram1D thetaeff;
private IHistogram1D ctheff;
private IHistogram1D etaeff;
- private IHistogram1D etafake;
private IHistogram1D d0eff1;
- private IHistogram1D d0eff2;
private IHistogram1D z0eff1;
private IHistogram1D z0eff2;
+ private IHistogram1D pTeff1Findable;
+ private IHistogram1D pTeff2Findable;
+ private IHistogram1D thetaeffFindable;
+ private IHistogram1D ctheffFindable;
+ private IHistogram1D etaeffFindable;
+ private IHistogram1D d0eff1Findable;
+ private IHistogram1D d0eff2Findable;
+ private IHistogram1D z0eff1Findable;
+ private IHistogram1D z0eff2Findable;
private IHistogram1D fakes;
private IHistogram1D nfakes;
+ private IHistogram1D etafake;
int trk_count = 0;
int nevt = 0;
+ int _nmcTrk = 0;
+ double _nrecTrk = 0;
public TrackAnalysisDriver() {
@@ -72,15 +92,25 @@
pTeff2 = hf.createHistogram1D("Efficiency vs pT full", "", 100, 0., 50., "type=efficiency");
thetaeff = hf.createHistogram1D("Efficiency vs theta", "", 72, 0., 180., "type=efficiency");
ctheff = hf.createHistogram1D("Efficiency vs cos(theta)", "", 50, -1., 1., "type=efficiency");
- etaeff = hf.createHistogram1D("Efficiency vs eta", "", 50, -2., 2., "type=efficiency");
- etafake = hf.createHistogram1D("Fake rate vs eta", "", 50, -2., 2., "type=efficiency");
- d0eff1 = hf.createHistogram1D("Efficiency vs d0", "", 50, -5., 5., "type=efficiency");
- d0eff2 = hf.createHistogram1D("Efficiency vs d0 full", "", 24, -12., 12., "type=efficiency");
- z0eff1 = hf.createHistogram1D("Efficiency vs z0", "", 50, -5., 5., "type=efficiency");
- z0eff2 = hf.createHistogram1D("Efficiency vs z0 full", "", 24, -12., 12., "type=efficiency");
+ etaeff = hf.createHistogram1D("Efficiency vs eta", "", 50, -2.5, 2.5, "type=efficiency");
+ d0eff1 = hf.createHistogram1D("Efficiency vs d0", "", 50, -2., 2., "type=efficiency");
+
+ z0eff1 = hf.createHistogram1D("Efficiency vs z0", "", 50, -50., 50., "type=efficiency");
+ z0eff2 = hf.createHistogram1D("Efficiency vs z0 full", "", 50, -200., 200., "type=efficiency");
+
+ pTeff1Findable = hf.createHistogram1D("Findable Efficiency vs pT", "", 100, 0., 5., "type=efficiency");
+ pTeff2Findable = hf.createHistogram1D("Findable Efficiency vs pT full", "", 100, 0., 50., "type=efficiency");
+ thetaeffFindable = hf.createHistogram1D("Findable Efficiency vs theta", "", 72, 0., 180., "type=efficiency");
+ ctheffFindable = hf.createHistogram1D("Findable Efficiency vs cos(theta)", "", 50, -1., 1., "type=efficiency");
+ etaeffFindable = hf.createHistogram1D("Findable Efficiency vs eta", "", 50, -2.5, 2.5, "type=efficiency");
+ d0eff1Findable = hf.createHistogram1D("Findable Efficiency vs d0", "", 50, -0.5, 0.5, "type=efficiency");
+ d0eff2Findable = hf.createHistogram1D("Findable Efficiency vs d0 full", "", 50, -5., 5., "type=efficiency");
+ z0eff1Findable = hf.createHistogram1D("Findable Efficiency vs z0", "", 50, -50., 50., "type=efficiency");
+ z0eff2Findable = hf.createHistogram1D("Findable Efficiency vs z0 full", "", 50, -200., 200., "type=efficiency");
+
fakes = hf.createHistogram1D("Number of mis-matched hits (unnormalized)", "", 10, 0., 10.);
nfakes = hf.createHistogram1D("Number of mis-matched hits (normalized)", "", 10, 0., 10.);
-
+ etafake = hf.createHistogram1D("Fake rate vs eta", "", 50, -2.5, 2.5, "type=efficiency");
}
@Override
@@ -88,16 +118,19 @@
// Increment the event counter
nevt++;
-
+ String resDir = "residualsPlots/";
+ String simDir = "STHitPlots/";
+ String debugDir = "debugPlots/";
// Get the magnetic field
Hep3Vector IP = new BasicHep3Vector(0., 0., 0.);
double bfield = event.getDetector().getFieldMap().getField(IP).z();
-// dump SThit information
+ String[] detNames = {"VtxPixelBarrel", "VtxPixelEndcap", "SCTShortBarrel", "SCTLongBarrel", "SCTShortEndcap", "SCTLongEndcap"};
+ // dump SThit information
String[] input_hit_collections = {"VtxBarrHits", "VtxEndcapHits", "SCTShortBarrHits", "SCTLongBarrHits", "SCTShortEndcapHits", "SCTLongEndcapHits"};
for (String input : input_hit_collections) {
List<SimTrackerHit> sthits = event.getSimTrackerHits(input);
- int[] nhits = {0, 0, 0, 0, 0, 0};
+ int[] nhits = {0, 0, 0, 0, 0, 0, 0};
for (SimTrackerHit st : sthits) {
String detector = st.getDetectorElement().getName();
int layer = st.getLayerNumber();
@@ -110,21 +143,55 @@
double phi = Math.atan2(hp[1], hp[0]);
// System.out.println("r= " + r + " theta = "+theta+" eta = " + eta+ " phi=" + phi);
nhits[layer]++;
- aida.cloud1D(input + " layer " + layer + " STHit eta").fill(eta);
- aida.cloud1D(input + " layer " + layer + " STHit phi").fill(phi);
- aida.cloud2D(input + " layer " + layer + " STHit phi vs eta").fill(eta, phi);
- aida.histogram2D(input + " layer " + layer + " STHit phi vs eta occupancy", 100, -2.5, 2.5, 100, -3.2, 3.2).fill(eta, phi);
+ aida.cloud1D(simDir + input + " layer " + layer + " STHit eta").fill(eta);
+ aida.cloud1D(simDir + input + " layer " + layer + " STHit phi").fill(phi);
+ aida.cloud2D(simDir + input + " layer " + layer + " STHit phi vs eta").fill(eta, phi);
+ aida.histogram2D(simDir + input + " layer " + layer + " STHit phi vs eta occupancy", 100, -2.5, 2.5, 100, -3.2, 3.2).fill(eta, phi);
}
int i = 0;
- while (i < 6) {
+ while (i < 7) {
if (nhits[i] > 0) {
- aida.cloud1D(input + "layer " + i + " number of ST hits").fill(nhits[i]);
+ aida.cloud1D(simDir + input + "layer " + i + " number of ST hits").fill(nhits[i]);
}
i++;
}
}
-
- List<HelicalTrackHit> hthits = event.get(HelicalTrackHit.class, "HelicalTrackHits");
+
+ List<SiTrackerHitStrip1D> stripHits = event.get(SiTrackerHitStrip1D.class, "StripClusterer_SiTrackerHitStrip1D");
+ List<SiTrackerHitPixel> pixelHits = event.get(SiTrackerHitPixel.class, "PixelClusterer_SiTrackerHitPixel");
+ List<RawTrackerHit> rawHits = event.get(RawTrackerHit.class, "RawTrackerHitMaker_RawTrackerHits");
+ List<HelicalTrackHit> hthits = event.get(HelicalTrackHit.class, "HelicalTrackHits");
+
+// int<String> occupancy;
+// Map occupancyMap;
+
+ Map<String, Integer> occupancyMap = new HashMap<String, Integer>();
+ for (RawTrackerHit rh : rawHits) {
+ IDetectorElement rhDetE = rh.getDetectorElement();
+
+ String rhDetName = rhDetE.getName();
+// System.out.println(rhDetName);
+ int rhLayer = rh.getLayerNumber();
+// String[] shortrhDetName=rhDetName.split("^[A-Z]+_layer[0-9]");
+
+ for (String myname : detNames) {
+ if (rhDetName.contains(myname)) {
+ String detlayer = myname + "_" + rhLayer;
+ Integer myint = occupancyMap.get(detlayer);
+ if (myint == null) {
+ myint = 1;
+ }
+ myint++;
+ occupancyMap.put(detlayer, myint);
+ }
+ }
+ }
+
+ Set<String> mykeyset = (Set<String>) occupancyMap.keySet();
+ for (String keys : mykeyset) {
+ aida.cloud1D("occupancyPlots/" + keys + " # of hits").fill(occupancyMap.get(keys));
+ }
+
for (HelicalTrackHit HelTrHit : hthits) {
}
@@ -133,7 +200,7 @@
// List<SeedStrategy> slist = StrategyXMLUtils.getStrategyListFromResource(
// StrategyXMLUtils.getDefaultStrategiesPrefix() + sfile);
- String sfile = StrategyXMLUtils.getDefaultStrategiesPrefix() + "sATLASBarrel-SM08.xml";
+ String sfile = StrategyXMLUtils.getDefaultStrategiesPrefix() + "sATLASFull-JeffMarch26.xml";
List<SeedStrategy> slist = StrategyXMLUtils.getStrategyListFromResource(sfile);
// Find the minimum pT among the strategies
@@ -179,18 +246,18 @@
HelicalTrackFit helixTrack = seed.getHelix();
double[] chisq = helixTrack.chisq();
double nhchisq = helixTrack.nhchisq();
- aida.cloud1D("Track Chi2-Circle Fit").fill(chisq[0]);
- aida.cloud1D("Track Chi2-RZ Fit").fill(chisq[1]);
- aida.cloud1D("NH Track Chi2").fill(nhchisq);
+ aida.cloud1D(debugDir + "Track Chi2-Circle Fit").fill(chisq[0]);
+ aida.cloud1D(debugDir + "Track Chi2-RZ Fit").fill(chisq[1]);
+ aida.cloud1D(debugDir + "NH Track Chi2").fill(nhchisq);
if (nhchisq != 0) {
- aida.cloud1D("NH!=0 Track Chi2-Circle Fit").fill(chisq[0]);
- aida.cloud1D("NH!=0 Track Chi2-RZ Fit").fill(chisq[1]);
+ aida.cloud1D(debugDir + "NH!=0 Track Chi2-Circle Fit").fill(chisq[0]);
+ aida.cloud1D(debugDir + "NH!=0 Track Chi2-RZ Fit").fill(chisq[1]);
}
List<HelicalTrackHit> hitlist = seed.getHits();
for (HelicalTrackHit hit : hitlist) {
int nhits = hit.getRawHits().size();
- aida.cloud1D(hit.Detector() + " nHits").fill(nhits);
+ aida.cloud1D(debugDir + hit.Detector() + " nHits").fill(nhits);
Hep3Vector HTHPos = hit.getCorrectedPosition();
double rHit = Math.sqrt(HTHPos.x() * HTHPos.x() + HTHPos.y() * HTHPos.y());
double zHit = HTHPos.z();
@@ -209,7 +276,7 @@
double du_axial = 0;
for (HelicalTrackStrip cluster : clusterlist) {
int nstrips = cluster.rawhits().size();
- aida.cloud1D(hit.Detector() + " nStrips-per-layer").fill(nstrips);
+ aida.cloud1D(debugDir + hit.Detector() + " nStrips-per-layer").fill(nstrips);
Hep3Vector corigin = cluster.origin();
Hep3Vector u = cluster.u();
List<RawTrackerHit> rawhits = cluster.rawhits();
@@ -227,8 +294,8 @@
// System.out.println("Layer number " + rhit.getLayerNumber() + " " + deName);
List<SimTrackerHit> sthits = rhit.getSimTrackerHits();
int nsthits = sthits.size();
- aida.cloud1D(hit.Detector() + " associated ST hits").fill(nsthits);
- aida.cloud1D(hit.Detector() + " layer" + stripdir + " associated ST hits").fill(nsthits);
+ aida.cloud1D(debugDir + hit.Detector() + " associated ST hits").fill(nsthits);
+ aida.cloud1D(debugDir + hit.Detector() + " layer" + stripdir + " associated ST hits").fill(nsthits);
if (nsthits == 1) {
double[] sthitD = sthits.get(0).getPoint();
BasicHep3Vector sthit = new BasicHep3Vector(sthitD);
@@ -239,35 +306,35 @@
}
- aida.cloud2D(hit.Detector() + "clusterSize vs eta").fill(etaHit, nstrips);
+ aida.cloud2D(debugDir + hit.Detector() + "clusterSize vs eta").fill(etaHit, nstrips);
// System.out.println("filling...");
if (umc != -999999) {
- aida.cloud2D(hit.Detector() + "cluster vs STHit dedx").fill(stenergy, charge);
- aida.cloud2D(hit.Detector() + "cluster dedx vs delte(u)").fill(umeas - umc, charge);
+ aida.cloud2D(debugDir + hit.Detector() + "cluster vs STHit dedx").fill(stenergy, charge);
+ aida.cloud2D(debugDir + hit.Detector() + "cluster dedx vs delte(u)").fill(umeas - umc, charge);
if (stripdir.contains("stereo")) {
du_stereo = umeas - umc;
}
if (stripdir.contains("axial")) {
du_axial = umeas - umc;
}
- aida.cloud1D(hit.Detector() + "layer=" + stripdir + " delta(u)").fill(umeas - umc);
- aida.cloud1D(hit.Detector() + " delta(u)").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + "layer=" + stripdir + " delta(u)").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + " delta(u)").fill(umeas - umc);
if (nstrips == 1) {
- aida.cloud1D(hit.Detector() + "layer=" + stripdir + " delta(u)--1 strip").fill(umeas - umc);
- aida.cloud1D(hit.Detector() + " delta(u)--1 strip").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + "layer=" + stripdir + " delta(u)--1 strip").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + " delta(u)--1 strip").fill(umeas - umc);
}
if (nstrips == 2) {
- aida.cloud1D(hit.Detector() + "layer=" + stripdir + " delta(u)--2 strip").fill(umeas - umc);
- aida.cloud1D(hit.Detector() + " delta(u)--2 strip").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + "layer=" + stripdir + " delta(u)--2 strip").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + " delta(u)--2 strip").fill(umeas - umc);
}
if (nstrips == 3) {
- aida.cloud1D(hit.Detector() + "layer=" + stripdir + " delta(u)--3 strip").fill(umeas - umc);
- aida.cloud1D(hit.Detector() + " delta(u)--3 strip").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + "layer=" + stripdir + " delta(u)--3 strip").fill(umeas - umc);
+ aida.cloud1D(debugDir + hit.Detector() + " delta(u)--3 strip").fill(umeas - umc);
}
}
}
- aida.cloud2D(hit.Detector() + " delta(u) stereo v axial").fill(du_stereo, du_axial);
+ aida.cloud2D(debugDir + hit.Detector() + " delta(u) stereo v axial").fill(du_stereo, du_axial);
}
MultipleScatter ms = seed.getMSMap().get(hit);
double msphi = ms.drphi();
@@ -283,14 +350,14 @@
double dz = posonhelix.z() - hitpos.z();
double dzErr = Math.sqrt(cov.e(2, 2));
- aida.cloud1D(hit.Detector() + " dxdy").fill(dxdy);
- aida.cloud1D(hit.Detector() + " dz").fill(dz);
- aida.cloud1D(hit.Detector() + " dxdy Pull").fill(dxdy / dxdyErr);
- aida.cloud1D(hit.Detector() + " dz Pull").fill(dz / dzErr);
+ aida.cloud1D(resDir + hit.Detector() + " dxdy").fill(dxdy);
+ aida.cloud1D(resDir + hit.Detector() + " dz").fill(dz);
+ aida.cloud1D(resDir + hit.Detector() + " dxdy Pull").fill(dxdy / dxdyErr);
+ aida.cloud1D(resDir + hit.Detector() + " dz Pull").fill(dz / dzErr);
if (Math.abs(dz) > 4) {
- aida.cloud1D(hit.Detector() + "Bad dz--nHits").fill(nhits);
+ aida.cloud1D(debugDir + hit.Detector() + "Bad dz--nHits").fill(nhits);
}
- aida.cloud1D("NH Chi2 for Hits on Track").fill(hitchisq);
+ aida.cloud1D(debugDir + "NH Chi2 for Hits on Track").fill(hitchisq);
}
// Analyze the hits on the track
@@ -317,7 +384,7 @@
aida.histogram1D("d0 for fake tracks", 50, -10., 10.).fill(d0);
aida.histogram1D("z0 for fake tracks", 50, -10., 10.).fill(z0);
aida.histogram1D("eta for fake tracks", 100, -2., 2.).fill(eta);
- etafake.fill(eta, 1.0);
+ etafake.fill(eta, 1.0);
} else {
aida.histogram1D("Hits for non-fake tracks", 20, 0., 20.).fill(nhits);
aida.histogram1D("pT for non-fake tracks", 100, 0., 10.).fill(pt);
@@ -325,7 +392,7 @@
aida.histogram1D("d0 for non-fake tracks", 50, -10., 10.).fill(d0);
aida.histogram1D("z0 for non-fake tracks", 50, -10., 10.).fill(z0);
aida.histogram1D("eta for non-fake tracks", 100, -2., 2.).fill(eta);
- etafake.fill(eta, 0.0);
+ etafake.fill(eta, 0.0);
}
aida.histogram1D("Hits for all tracks", 20, 0., 20.).fill(nhits);
aida.histogram1D("pT for all tracks", 100, 0., 10.).fill(pt);
@@ -364,33 +431,36 @@
double d0resid = (d0tk - d0mc);
// Plot the pt and d0 pulls for various purity intervals
if (nbad == 0) {
- aida.histogram2D("pT MC vs pT Reco for 0 Bad Hits",
+ aida.histogram2D(resDir + "pT MC vs pT Reco for 0 Bad Hits",
100, 0., 5., 100, 0., 5.).fill(ptmc, pttk);
- aida.histogram2D("d0 MC vs d0 Reco for 0 Bad Hits",
+ aida.histogram2D(resDir + "d0 MC vs d0 Reco for 0 Bad Hits",
100, -0.2, 0.2, 100, -0.2, 0.2).fill(d0mc, d0tk);
- aida.histogram1D("pT Pull for 0 Bad Hits", 100, -10., 10.).fill(ptpull);
- aida.histogram1D("d0 pull for 0 Bad Hits", 100, -10., 10.).fill(d0pull);
- aida.cloud1D("pT Residual for 0 Bad Hits").fill(ptresid);
- aida.cloud1D("d0 Residual for 0 Bad Hits").fill(d0resid);
+ aida.histogram1D(resDir + "pT Pull for 0 Bad Hits", 100, -10., 10.).fill(ptpull);
+ aida.histogram1D(resDir + "d0 pull for 0 Bad Hits", 100, -10., 10.).fill(d0pull);
+ aida.cloud1D(resDir + "pT Residual for 0 Bad Hits").fill(ptresid);
+ aida.cloud1D(resDir + "d0 Residual for 0 Bad Hits").fill(d0resid);
+ aida.cloud1D(resDir + "1/pT for 0 Bad Hits").fill(1 / pttk);
} else if (purity > 0.5) {
- aida.histogram2D("pT MC vs pT Reco for 0.5 < purity < 1",
+ aida.histogram2D(resDir + "pT MC vs pT Reco for 0.5 < purity < 1",
100, 0., 5., 100, 0., 5.).fill(ptmc, pttk);
- aida.histogram2D("d0 MC vs d0 Reco for 0.5 < purity < 1",
+ aida.histogram2D(resDir + "d0 MC vs d0 Reco for 0.5 < purity < 1",
100, -0.2, 0.2, 100, -0.2, 0.2).fill(d0mc, d0tk);
- aida.histogram1D("pT Pull for 0.5 < purity < 1", 100, -10., 10.).fill(ptpull);
- aida.histogram1D("d0 pull for 0.5 < purity < 1", 100, -10., 10.).fill(d0pull);
- aida.cloud1D("pT Residual for 0.5 < purity < 1").fill(ptresid);
- aida.cloud1D("d0 Residual for 0.5 < purity < 1").fill(d0resid);
+ aida.histogram1D(resDir + "pT Pull for 0.5 < purity < 1", 100, -10., 10.).fill(ptpull);
+ aida.histogram1D(resDir + "d0 pull for 0.5 < purity < 1", 100, -10., 10.).fill(d0pull);
+ aida.cloud1D(resDir + "pT Residual for 0.5 < purity < 1").fill(ptresid);
+ aida.cloud1D(resDir + "d0 Residual for 0.5 < purity < 1").fill(d0resid);
+ aida.cloud1D(resDir + "1/pT for 0.5 < purity < 1").fill(1 / pttk);
} else if (purity < 0.5) {
- aida.histogram2D("pT MC vs pT Reco for purity <= 0.5",
+ aida.histogram2D(resDir + "pT MC vs pT Reco for purity <= 0.5",
100, 0., 5., 100, 0., 5.).fill(ptmc, pttk);
- aida.histogram2D("d0 MC vs d0 Reco for purity <= 0.5",
+ aida.histogram2D(resDir + "d0 MC vs d0 Reco for purity <= 0.5",
100, -0.2, 0.2, 100, -0.2, 0.2).fill(d0mc, d0tk);
- aida.histogram1D("pT Pull for purity <= 0.5", 100, -10., 10.).fill(ptpull);
- aida.histogram1D("d0 pull for purity <= 0.5", 100, -10., 10.).fill(d0pull);
- aida.cloud1D("pT Residual for purity <= 0.5").fill(ptresid);
- aida.cloud1D("d0 Residial for purity <= 0.5").fill(d0resid);
+ aida.histogram1D(resDir + "pT Pull for purity <= 0.5", 100, -10., 10.).fill(ptpull);
+ aida.histogram1D(resDir + "d0 pull for purity <= 0.5", 100, -10., 10.).fill(d0pull);
+ aida.cloud1D(resDir + "pT Residual for purity <= 0.5").fill(ptresid);
+ aida.cloud1D(resDir + "d0 Residial for purity <= 0.5").fill(d0resid);
+ aida.cloud1D(resDir + "1/pT for purity <= 0.5").fill(1 / pttk);
}
}
}
@@ -453,8 +523,8 @@
if (ntrk > 0) {
wgt = 1.;
}
- pTeff1.fill(pt, wgt);
- pTeff2.fill(pt, wgt);
+ pTeff1Findable.fill(pt, wgt);
+ pTeff2Findable.fill(pt, wgt);
}
// Make angular efficiency plot
@@ -462,10 +532,13 @@
double wgt = 0.;
if (ntrk > 0) {
wgt = 1.;
+ } else {
+ System.out.println("Findable Track Not Found! eta=" + eta);
}
- thetaeff.fill(theta, wgt);
- ctheff.fill(cth, wgt);
- etaeff.fill(eta, wgt);
+ thetaeffFindable.fill(theta, wgt);
+ ctheffFindable.fill(cth, wgt);
+ etaeffFindable.fill(eta, wgt);
+
}
// Make d0 efficiency plot
@@ -474,8 +547,8 @@
if (ntrk > 0) {
wgt = 1.;
}
- d0eff1.fill(d0, wgt);
- d0eff2.fill(d0, wgt);
+ d0eff1Findable.fill(d0, wgt);
+ d0eff2Findable.fill(d0, wgt);
}
// Make z0 efficiency plot
@@ -484,45 +557,73 @@
if (ntrk > 0) {
wgt = 1.;
}
- z0eff1.fill(z0, wgt);
- z0eff2.fill(z0, wgt);
+ z0eff1Findable.fill(z0, wgt);
+ z0eff2Findable.fill(z0, wgt);
}
// Select charged MC particles
if (mcp.getCharge() == 0) {
continue;
}
+// make the true efficiency plots
+ double ptTrkCut = 1.0; //GeV
+ double d0TrkCut = 2.0; //mm
+ double z0TrkCut = 200.0; //mm
+ double etaTrkCut = 2.5;
+// System.out.println("Final Stat Part? "+mcp.FINAL_STATE+"; pt = "+pt+"; d0 = "+d0);
+ if (pt > ptTrkCut && mcp.getGeneratorStatus() == mcp.FINAL_STATE && Math.abs(d0) < d0TrkCut && Math.abs(eta) < etaTrkCut && Math.abs(z0) < z0TrkCut) {
+ double wgt = 0.0;
+ if (ntrk > 0) {
+ wgt = 1.0;
+// System.out.println("found track!");
+ System.out.println("Found this track! eta = " + eta + "; pT = " + pt + "; z0 = " + z0 + "; d0 = " + d0);
+ } else {
+ System.out.println("Missed this track! eta = " + eta + "; pT = " + pt + "; z0 = " + z0 + "; d0 = " + d0);
+ }
+ pTeff1.fill(pt, wgt);
+ pTeff2.fill(pt, wgt);
+ ctheff.fill(cth, wgt);
+ thetaeff.fill(theta, wgt);
+ etaeff.fill(eta, wgt);
+ d0eff1.fill(d0, wgt);
+ z0eff1.fill(z0, wgt);
+ z0eff2.fill(z0, wgt);
+ if (eta < etaTrkCut) {
+ _nmcTrk++;
+ _nrecTrk += wgt;
+ }
+ }
// Select mcp that fail the final state requirement
if (mcp.getGeneratorStatus() != mcp.FINAL_STATE) {
- aida.histogram1D("Hits for non-final state particles", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for non-final state particles", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for non-final state particles", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for non-final state particles", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for non-final state particles", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for non-final state particles", 100, -100., 100.).fill(z0);
- aida.cloud2D("Hits vs eta for non-final state particles").fill(eta, nhits);
+ aida.cloud1D("findable/Hits for non-final state particles").fill(nhits);
+ aida.cloud1D("findable/pT for non-final state particles").fill(pt);
+ aida.cloud1D("findable/cos(theta) for non-final state particles").fill(cth);
+ aida.cloud1D("findable/eta for non-final state particles").fill(eta);
+ aida.cloud1D("findable/d0 for non-final state particles").fill(d0);
+ aida.cloud1D("findable/z0 for non-final state particles").fill(z0);
+ aida.cloud2D("findable/Hits vs eta for non-final state particles").fill(eta, nhits);
continue;
}
// Make plots for the base sample
- aida.histogram1D("Hits for base MC selection", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for base MC selection", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for base MC selection", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for base MC selection", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for base MC selection", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for base MC selection", 100, -100., 100.).fill(z0);
- aida.cloud2D("Hits vs eta for base MC selection").fill(eta, nhits);
+ aida.cloud1D("findable/Hits for base MC selection").fill(nhits);
+ aida.cloud1D("findable/pT for base MC selection").fill(pt);
+ aida.cloud1D("findable/cos(theta) for base MC selection").fill(cth);
+ aida.cloud1D("findable/eta for base MC selection").fill(eta);
+ aida.cloud1D("findable/d0 for base MC selection").fill(d0);
+ aida.cloud1D("findable/z0 for base MC selection").fill(z0);
+ aida.cloud2D("findable/Hits vs eta for base MC selection").fill(eta, nhits);
// Make plots for findable tracks
if (findable.isFindable(mcp, slist)) {
- aida.histogram1D("Hits for findable tracks", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for findable tracks", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for findable tracks", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for findable tracks", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for findable tracks", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for findable tracks", 100, -100., 100.).fill(z0);
- aida.cloud2D("Hits vs eta for findable tracks").fill(eta, nhits);
+ aida.cloud1D("findable/Hits for findable tracks").fill(nhits);
+ aida.cloud1D("findable/pT for findable tracks").fill(pt);
+ aida.cloud1D("findable/cos(theta) for findable tracks").fill(cth);
+ aida.cloud1D("findable/eta for findable tracks").fill(eta);
+ aida.cloud1D("findable/d0 for findable tracks").fill(d0);
+ aida.cloud1D("findable/z0 for findable tracks").fill(z0);
+ aida.cloud2D("findable/Hits vs eta for findable tracks").fill(eta, nhits);
continue;
}
@@ -533,80 +634,90 @@
ignores.add(Ignore.NoZ0Cut);
if (findable.isFindable(mcp, slist, ignores)) {
- aida.histogram1D("Hits for z0 check failures", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for z0 check failures", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for z0 check failures", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for z0 check failures", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for z0 check failures", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for z0 check failures", 100, -100., 100.).fill(z0);
+ aida.cloud1D("findable/Hits for z0 check failures").fill(nhits);
+ aida.cloud1D("findable/pT for z0 check failures").fill(pt);
+ aida.cloud1D("findable/cos(theta) for z0 check failures").fill(cth);
+ aida.cloud1D("findable/eta for z0 check failures").fill(eta);
+ aida.cloud1D("findable/d0 for z0 check failures").fill(d0);
+ aida.cloud1D("findable/z0 for z0 check failures").fill(z0);
continue;
}
// Select mc particles that fail on the d0 cut
ignores.add(Ignore.NoDCACut);
if (findable.isFindable(mcp, slist, ignores)) {
- aida.histogram1D("Hits for d0 check failures", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for d0 check failures", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for d0 check failures", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for d0 check failures", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for d0 check failures", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for d0 check failures", 100, -100., 100.).fill(z0);
+ aida.cloud1D("findable/Hits for d0 check failures").fill(nhits);
+ aida.cloud1D("findable/pT for d0 check failures").fill(pt);
+ aida.cloud1D("findable/cos(theta) for d0 check failures").fill(cth);
+ aida.cloud1D("findable/eta for d0 check failures").fill(eta);
+ aida.cloud1D("findable/d0 for d0 check failures").fill(d0);
+ aida.cloud1D("findable/z0 for d0 check failures").fill(z0);
continue;
}
// select mc particles that fail the confirm check
ignores.add(Ignore.NoConfirmCheck);
if (findable.isFindable(mcp, slist, ignores)) {
- aida.histogram1D("Hits for confirm check failures", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for confir check failures", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for confirm check failures", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for confirm check failures", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for seed confirm failures", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for seed confirm failures", 100, -100., 100.).fill(z0);
+ aida.cloud1D("findable/Hits for confirm check failures").fill(nhits);
+ aida.cloud1D("findable/pT for confir check failures").fill(pt);
+ aida.cloud1D("findable/cos(theta) for confirm check failures").fill(cth);
+ aida.cloud1D("findable/eta for confirm check failures").fill(eta);
+ aida.cloud1D("findable/d0 for seed confirm failures").fill(d0);
+ aida.cloud1D("findable/z0 for seed confirm failures").fill(z0);
continue;
}
// select mc particles that fail on the seed check
ignores.add(Ignore.NoSeedCheck);
if (findable.isFindable(mcp, slist, ignores)) {
- aida.histogram1D("Hits for seed check failures", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for seed check failures", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for seed check failures", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for seed check failures", 100, -1., 1.).fill(eta);
- aida.histogram1D("d0 for seed check failures", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for seed check failures", 100, -100., 100.).fill(z0);
+ aida.cloud1D("findable/Hits for seed check failures").fill(nhits);
+ aida.cloud1D("findable/pT for seed check failures").fill(pt);
+ aida.cloud1D("findable/cos(theta) for seed check failures").fill(cth);
+ aida.cloud1D("findable/eta for seed check failures").fill(eta);
+ aida.cloud1D("findable/d0 for seed check failures").fill(d0);
+ aida.cloud1D("findable/z0 for seed check failures").fill(z0);
continue;
}
// Select mc particles that fail the number of hit cut
ignores.add(Ignore.NoMinHitCut);
if (findable.isFindable(mcp, slist, ignores)) {
- aida.histogram1D("Hits for nhit check failures", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for nhit check failures", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for nhit check failures", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for nhit check failures", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for nhit check failures", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for nhit check failures", 100, -100., 100.).fill(z0);
+ aida.cloud1D("findable/Hits for nhit check failures").fill(nhits);
+ aida.cloud1D("findable/pT for nhit check failures").fill(pt);
+ aida.cloud1D("findable/cos(theta) for nhit check failures").fill(cth);
+ aida.cloud1D("findable/eta for nhit check failures").fill(eta);
+ aida.cloud1D("findable/d0 for nhit check failures").fill(d0);
+ aida.cloud1D("findable/z0 for nhit check failures").fill(z0);
continue;
}
// Select mc particles that fail on the pT cut
ignores.add(Ignore.NoPTCut);
if (findable.isFindable(mcp, slist, ignores)) {
- aida.histogram1D("Hits for pT check failures", 20, 0., 20.).fill(nhits);
- aida.histogram1D("pT for pT check failures", 100, 0., 10.).fill(pt);
- aida.histogram1D("cos(theta) for pT check failures", 100, -1., 1.).fill(cth);
- aida.histogram1D("eta for pT check failures", 100, -2., 2.).fill(eta);
- aida.histogram1D("d0 for pT check failures", 100, -100., 100.).fill(d0);
- aida.histogram1D("z0 for pT check failures", 100, -100., 100.).fill(z0);
+ aida.cloud1D("findable/Hits for pT check failures").fill(nhits);
+ aida.cloud1D("findable/pT for pT check failures").fill(pt);
+ aida.cloud1D("findable/cos(theta) for pT check failures").fill(cth);
+ aida.cloud1D("findable/eta for pT check failures").fill(eta);
+ aida.cloud1D("findable/d0 for pT check failures").fill(d0);
+ aida.cloud1D("findable/z0 for pT check failures").fill(z0);
} else {
System.out.println("MC Particle is not findable with all ignores set!!");
}
-
+
}
return;
}
+ @Override
+ public void endOfData() {
+ try {
+ aida.saveAs("myplots.aida");
+ } catch (IOException ex) {
+ Logger.getLogger(TrackAnalysisDriver.class.getName()).log(Level.SEVERE, null, ex);
+ }
+ System.out.println("# of reco tracks = " + _nrecTrk + "; # of MC tracks = " + _nmcTrk + "; Efficiency = " + _nrecTrk / _nmcTrk);
+ }
+
private double getr(double x, double y) {
return Math.sqrt(x * x + y * y);
}