lcsim-contrib/src/main/java/org/lcsim/contrib/sATLAS/UTOPIA7
diff -u -r1.1 -r1.2
--- TrackAnalysisDriver.java 20 Aug 2009 05:06:59 -0000 1.1
+++ TrackAnalysisDriver.java 20 Aug 2009 17:30:41 -0000 1.2
@@ -7,12 +7,11 @@
import java.io.IOException;
import java.util.logging.Level;
import java.util.logging.Logger;
-import org.lcsim.contrib.mgraham.sATLASDigi.*;
+//import org.lcsim.contrib.mgraham.sATLASDigi.*;
import hep.physics.matrix.SymmetricMatrix;
import hep.physics.vec.BasicHep3Vector;
import hep.physics.vec.Hep3Vector;
-import hep.physics.vec.VecOp;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
@@ -20,9 +19,9 @@
import hep.aida.*;
import java.util.Set;
-//import org.lcsim.contrib.Partridge.TrackingTest.FindableTrack.Ignore;
-//import org.lcsim.contrib.Partridge.TrackingTest.TrackAnalysis;
-import org.lcsim.detector.IDetectorElement;
+
+import org.lcsim.contrib.sATLAS.FindableTrack;
+import org.lcsim.contrib.sATLAS.TrackAnalysis;
import org.lcsim.event.EventHeader;
import org.lcsim.event.LCRelation;
import org.lcsim.event.MCParticle;
@@ -31,19 +30,12 @@
import org.lcsim.event.RawTrackerHit;
import org.lcsim.event.SimTrackerHit;
import org.lcsim.event.base.BaseRelationalTable;
-import org.lcsim.fit.helicaltrack.HelicalTrackCross;
import org.lcsim.fit.helicaltrack.HelicalTrackFit;
import org.lcsim.fit.helicaltrack.HelicalTrackHit;
-import org.lcsim.fit.helicaltrack.HelicalTrackStrip;
import org.lcsim.fit.helicaltrack.HelixParamCalculator;
-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;
import org.lcsim.recon.tracking.seedtracker.StrategyXMLUtils;
import org.lcsim.util.Driver;
import org.lcsim.util.aida.AIDA;
@@ -55,80 +47,70 @@
public class TrackAnalysisDriver extends Driver {
private AIDA aida = AIDA.defaultInstance();
- private IHistogram1D pTeff1;
- private IHistogram1D pTeff2;
- private IHistogram1D thetaeff;
- private IHistogram1D ctheff;
- private IHistogram1D etaeff;
- private IHistogram1D d0eff1;
- 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 IProfile1D pTeff1;
+ private IProfile1D pTeff2;
+ private IProfile1D thetaeff;
+ private IProfile1D ctheff;
+ private IProfile1D etaeff;
+ private IProfile1D d0eff1;
+ private IProfile1D z0eff1;
+ private IProfile1D z0eff2;
+ private IProfile1D pTeff1Findable;
+ private IProfile1D pTeff2Findable;
+ private IProfile1D thetaeffFindable;
+ private IProfile1D ctheffFindable;
+ private IProfile1D etaeffFindable;
+ private IProfile1D d0eff1Findable;
+ private IProfile1D d0eff2Findable;
+ private IProfile1D z0eff1Findable;
+ private IProfile1D z0eff2Findable;
private IHistogram1D fakes;
private IHistogram1D nfakes;
- private IHistogram1D etafake;
+ private IProfile1D etafake;
public String outputPlots = "myplots.aida";
Map<String, IProfile1D> clsizeMap = new HashMap<String, IProfile1D>();
// String[] detNames = {"VtxPixelBarrel", "VtxPixelEndcap", "SCTShortBarrel", "SCTLongBarrel", "SCTShortEndcap", "SCTLongEndcap"};
- String[] detNames = {"VtxPixelBarrel", "VtxPixelEndcap", "SCTShortBarrel", "SCTLongBarrel", "SCTEndcap", };
+ String[] detNames = {"VtxPixelBarrel", "VtxPixelEndcap", "SCTShortBarrel", "SCTLongBarrel", "SCTEndcap",};
// Integer[] nlayers = {4, 6, 3, 2, 5, 5};
- Integer[] nlayers = {4, 6, 3, 2, 5};
+ Integer[] nlayers = {4, 6, 3, 2, 5};
int trk_count = 0;
int nevt = 0;
int _nmcTrk = 0;
double _nrecTrk = 0;
- IAnalysisFactory af = IAnalysisFactory.create();
- ITree tree = af.createTreeFactory().create();
- ITupleFactory tf = af.createTupleFactory(tree);
- String columnString = "int ievt=0 , nMCPart=0; ITuple MCInfo={double pXMC=0, pYMC=0, pXMC=0};";
- ITuple tuple = tf.create("tuple", "MyNtpule", columnString);
+ double ptTrkCut = 1.0; //GeV
+ double d0TrkCut = 2.0; //mm
+ double z0TrkCut = 200.0; //mm
+ double etaTrkCut = 2.5;
+ double etaBarrel = 0.97;
+ double etaBarrPixECStrips=1.5;
public TrackAnalysisDriver() {
// Define the efficiency histograms
IHistogramFactory hf = aida.histogramFactory();
- pTeff1 = hf.createHistogram1D("Efficiency vs pT", "", 100, 0., 5., "type=efficiency");
- 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.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");
+ pTeff1 = hf.createProfile1D("Efficiency vs pT", "", 100, 0., 5.);
+ pTeff2 = hf.createProfile1D("Efficiency vs pT full", "", 100, 0., 50.);
+ thetaeff = hf.createProfile1D("Efficiency vs theta", "", 72, 0., 180.);
+ ctheff = hf.createProfile1D("Efficiency vs cos(theta)", "", 50, -1., 1.);
+ etaeff = hf.createProfile1D("Efficiency vs eta", "", 50, -2.5, 2.5);
+ d0eff1 = hf.createProfile1D("Efficiency vs d0", "", 50, -2., 2.);
+
+ z0eff1 = hf.createProfile1D("Efficiency vs z0", "", 50, -50., 50.);
+ z0eff2 = hf.createProfile1D("Efficiency vs z0 full", "", 50, -200., 200.);
+
+ pTeff1Findable = hf.createProfile1D("Findable Efficiency vs pT", "", 100, 0., 5.);
+ pTeff2Findable = hf.createProfile1D("Findable Efficiency vs pT full", "", 100, 0., 50.);
+ thetaeffFindable = hf.createProfile1D("Findable Efficiency vs theta", "", 72, 0., 180.);
+ ctheffFindable = hf.createProfile1D("Findable Efficiency vs cos(theta)", "", 50, -1., 1.);
+ etaeffFindable = hf.createProfile1D("Findable Efficiency vs eta", "", 50, -2.5, 2.5);
+ d0eff1Findable = hf.createProfile1D("Findable Efficiency vs d0", "", 50, -0.5, 0.5);
+ d0eff2Findable = hf.createProfile1D("Findable Efficiency vs d0 full", "", 50, -5., 5.);
+ z0eff1Findable = hf.createProfile1D("Findable Efficiency vs z0", "", 50, -50., 50.);
+ z0eff2Findable = hf.createProfile1D("Findable Efficiency vs z0 full", "", 50, -200., 200.);
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");
- int i, j;
-// for (i = 0; i < 6; i++) {
- for (i = 0; i < 5; i++) {
- for (j = 0; j < nlayers[i]; j++) {
- int laynum = j + 1;
- String profname = detNames[i] + "_layer" + laynum + " cluster size vs eta";
- String key = detNames[i] + "_layer" + laynum;
- clsizeMap.put(key, hf.createProfile1D(profname, 25, -2.5, 2.5));
- }
- }
-
+ etafake = hf.createProfile1D("Fake rate vs eta", "", 50, -2.5, 2.5);
}
@Override
@@ -148,7 +130,7 @@
// dump SThit information
// String[] input_hit_collections = {"VtxBarrHits", "VtxEndcapHits", "SCTShortBarrHits", "SCTLongBarrHits", "SCTShortEndcapHits", "SCTLongEndcapHits"};
- String[] input_hit_collections = {"VtxBarrHits", "VtxEndcapHits", "SCTShortBarrHits", "SCTLongBarrHits", "SCTEndcapHits"};
+ String[] input_hit_collections = {"VtxBarrHits", "VtxEndcapHits", "SCTShortBarrHits", "SCTLongBarrHits", "SCTEndcapHits"};
for (String input : input_hit_collections) {
List<SimTrackerHit> sthits = event.getSimTrackerHits(input);
int[] nhits = {0, 0, 0, 0, 0, 0, 0};
@@ -183,101 +165,6 @@
List<RawTrackerHit> rawHits = event.get(RawTrackerHit.class, "RawTrackerHitMaker_RawTrackerHits");
List<HelicalTrackHit> hthits = event.get(HelicalTrackHit.class, "HelicalTrackHits");
- 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(occDir + keys + " # of hits").fill(occupancyMap.get(keys));
- }
- int detNum = 0;
- int layerNum = 0;
-
- for (SiTrackerHitPixel stripCluster : pixelHits) {
-
- //System.out.println("Number of pixels: "+stripCluster.getReadoutElectrodes().getNCells());
- stripCluster.getSensor().getReadout().getHits(RawTrackerHit.class); //gets all hits on sensor
- Hep3Vector strCluPos = stripCluster.getPositionAsVector();
- double rHit = Math.sqrt(strCluPos.x() * strCluPos.x() + strCluPos.y() * strCluPos.y());
- double zHit = strCluPos.z();
- double etaHit = -Math.log(Math.tan(Math.atan2(rHit, zHit) / 2));
- List<RawTrackerHit> rthList = stripCluster.getRawHits();
- int nhits = rthList.size();
- //System.out.println("Number of rawhits: "+nhits);
- String detlayer = "Foobar";
- for (RawTrackerHit rth : rthList) {
- IDetectorElement rhDetE = rth.getDetectorElement();
- String rhDetName = rhDetE.getName();
-// System.out.println(rhDetName);
- int rhLayer = rth.getLayerNumber();
- String[] detLayerName = rhDetName.split("_module");
-// System.out.println(detLayerName[0]);
- for (String myname : detNames) {
- if (rhDetName.contains(myname)) {
- detlayer = myname + "_layer" + rhLayer;
- }
- }
- }
-// System.out.println(detlayer);
- clsizeMap.get(detlayer).fill(etaHit, nhits);
- aida.cloud1D(occDir + detlayer + " cluster size").fill(nhits);
-// ClusterSize[detNum][layerNum].fill(etaHit, nhits);
-
- }
-
- for (SiTrackerHitStrip1D stripCluster : stripHits) {
- Hep3Vector strCluPos = stripCluster.getPositionAsVector();
- double rHit = Math.sqrt(strCluPos.x() * strCluPos.x() + strCluPos.y() * strCluPos.y());
- double zHit = strCluPos.z();
- double etaHit = -Math.log(Math.tan(Math.atan2(rHit, zHit) / 2));
- List<RawTrackerHit> rthList = stripCluster.getRawHits();
- int nhits = rthList.size();
- String detlayer = "Foobar";
-
- for (RawTrackerHit rth : rthList) {
-
- IDetectorElement rhDetE = rth.getDetectorElement();
- String rhDetName = rhDetE.getName();
- int rhLayer = rth.getLayerNumber();
- for (String myname : detNames) {
- if (rhDetName.contains(myname)) {
- detlayer = myname + "_layer" + rhLayer;
- }
- }
- }
-// System.out.println(detlayer);
- clsizeMap.get(detlayer).fill(etaHit, nhits);
- aida.cloud1D(occDir + detlayer + " cluster size").fill(nhits);
- }
-
- for (HelicalTrackHit HelTrHit : hthits) {
- }
-
- // Get the list of strategies being used
-// String sfile = "autogen_ttbar_sid02_vs.xml";
-// List<SeedStrategy> slist = StrategyXMLUtils.getStrategyListFromResource(
-// StrategyXMLUtils.getDefaultStrategiesPrefix() + sfile);
-
-// String sfile = StrategyXMLUtils.getDefaultStrategiesPrefix() + "sATLASFull-JeffMarch26.xml";
-// String sfile = StrategyXMLUtils.getDefaultStrategiesPrefix() + "sATLASFull-Simple.xml";
String sfile = StrategyXMLUtils.getDefaultStrategiesPrefix() + "sATLASFull-UTOPIA7.xml";
List<SeedStrategy> slist = StrategyXMLUtils.getStrategyListFromResource(sfile);
@@ -317,135 +204,6 @@
double d0 = track.getTrackParameter(HelicalTrackFit.dcaIndex);
double z0 = track.getTrackParameter(HelicalTrackFit.z0Index);
- SeedTrack st = (SeedTrack) track;
-
- SeedCandidate seed = st.getSeedCandidate();
-
- HelicalTrackFit helixTrack = seed.getHelix();
-
- double[] chisq = helixTrack.chisq();
- double nhchisq = helixTrack.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(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(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();
- double etaHit = -Math.log(Math.tan(Math.atan2(rHit, zHit) / 2));
- // System.out.println("Looping over hit in " + hit.Detector());
- double hitchisq = hit.chisq();
- double s = helixTrack.PathMap().get(hit);
- Hep3Vector posonhelix = HelixUtils.PointOnHelix(helixTrack, s);
-
- if (hit instanceof HelicalTrackCross) {
- HelicalTrackCross cross = (HelicalTrackCross) hit;
- TrackDirection trkdir = HelixUtils.CalculateTrackDirection(helixTrack, s);
- cross.setTrackDirection(trkdir, helixTrack.covariance());
- List<HelicalTrackStrip> clusterlist = cross.getStrips();
- double du_stereo = 0;
- double du_axial = 0;
- for (HelicalTrackStrip cluster : clusterlist) {
- int nstrips = cluster.rawhits().size();
- aida.cloud1D(debugDir + hit.Detector() + " nStrips-per-layer").fill(nstrips);
- Hep3Vector corigin = cluster.origin();
- Hep3Vector u = cluster.u();
- List<RawTrackerHit> rawhits = cluster.rawhits();
- double umc = -999999;
- double stenergy = -999999;
- String stripdir = "axial";
- double umeas = cluster.umeas();
- double umeasErr = cluster.du();
-
- double charge = cluster.dEdx();
- double layer = cluster.layer();
- for (RawTrackerHit rhit : rawhits) {
- String deName = rhit.getDetectorElement().getName();
- if (deName.contains("sensor1")) {
- stripdir = "stereo";
- }
- // System.out.println("Layer number " + rhit.getLayerNumber() + " " + deName);
- List<SimTrackerHit> sthits = rhit.getSimTrackerHits();
- int nsthits = sthits.size();
- 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);
- stenergy = sthits.get(0).getdEdx();
- Hep3Vector vdiff = VecOp.sub(sthit, corigin);
- umc = VecOp.dot(vdiff, u);
- }
- }
-
-
- // System.out.println("filling...");
- if (umc != -999999) {
- 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(debugDir + hit.Detector() + "layer=" + stripdir + " delta(u)").fill(umeas - umc);
- aida.cloud1D(debugDir + hit.Detector() + " delta(u)").fill(umeas - umc);
- aida.cloud1D(debugDir + hit.Detector() + " delta(u) Error").fill(umeasErr);
- if (nstrips == 1) {
- 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);
- aida.cloud1D(debugDir + hit.Detector() + " delta(u) Error--1 strip").fill(umeasErr);
- }
- if (nstrips == 2) {
- 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);
- aida.cloud1D(debugDir + hit.Detector() + " delta(u) Error--2 strip").fill(umeasErr);
- }
- if (nstrips == 3) {
- 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.cloud1D(debugDir + hit.Detector() + " delta(u) Error--3 strip").fill(umeasErr);
- }
- }
-
- }
- 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();
- double msz = ms.dz();
- Hep3Vector hitpos = hit.getCorrectedPosition();
- SymmetricMatrix cov = hit.getCorrectedCovMatrix();
- // double rHit = getr(hitpos.x(), hitpos.y());
- // double rHel = getr(posonhelix.x(), posonhelix.y());
- // double phiHet = getphi(hitpos.x(), hitpos.y());
- // double phiHel = getr(posonhelix.x(), posonhelix.y());
- double dxdy = Math.sqrt(Math.pow(hitpos.x() - posonhelix.x(), 2) + Math.pow(hitpos.y() - posonhelix.y(), 2));
- double dxdyErr = getdxdyErr(hitpos, posonhelix, cov);
- double dz = posonhelix.z() - hitpos.z();
-
- double dzErr = Math.sqrt(cov.e(2, 2));
- 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);
- aida.cloud2D(resDir + hit.Detector() + " dz vs hit z").fill(hitpos.z(), dz);
- aida.cloud2D(resDir + hit.Detector() + " dz Pull vs hit z").fill(hitpos.z(), dz / dzErr);
- if (Math.abs(dz) > 4) {
- aida.cloud1D(debugDir + hit.Detector() + "Bad dz--nHits").fill(nhits);
- }
- aida.cloud1D(debugDir + "NH Chi2 for Hits on Track").fill(hitchisq);
- }
-
// Analyze the hits on the track
TrackAnalysis tkanal = new TrackAnalysis(track, hittomc);
@@ -653,10 +411,7 @@
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;
lcsim-contrib/src/main/java/org/lcsim/contrib/sATLAS/UTOPIA7
diff -u -r1.1 -r1.2
--- OccupancyDriver.java 20 Aug 2009 05:06:59 -0000 1.1
+++ OccupancyDriver.java 20 Aug 2009 17:30:41 -0000 1.2
@@ -56,9 +56,9 @@
_process_paths.add("SCTShortBarrel");
_process_paths.add("SCTLongBarrel");
_process_paths.add("VtxPixelEndcap");
- _process_paths.add("SCTShortEndcap");
- _process_paths.add("SCTLongEndcap");
-
+// _process_paths.add("SCTShortEndcap");
+// _process_paths.add("SCTLongEndcap");
+ _process_paths.add("SCTEndcap");
// Define the efficiency histograms
_hf = aida.histogramFactory();
}
@@ -70,7 +70,7 @@
// Increment the event counter
nevt++;
- System.out.println("B: "+event.getDetector().getFieldMap().getField(new BasicHep3Vector(0., 0.,0.)));
+// System.out.println("B: " + event.getDetector().getFieldMap().getField(new BasicHep3Vector(0., 0., 0.)));
ReadoutChip chip = new GenericReadoutChip();
@@ -78,10 +78,10 @@
for (SiSensor sensor : _process_sensors) {
SiSensorElectrodes electrodes = sensor.getReadoutElectrodes(ChargeCarrier.HOLE);
- if(sensor.hasPixels()){
+ if (sensor.hasPixels()) {
// System.out.println("Found some pixels");
electrodes = sensor.getReadoutElectrodes(ChargeCarrier.ELECTRON);
- } else{
+ } else {
// System.out.println("...got strips");
}
int nchan = electrodes.getNCells();
@@ -92,20 +92,20 @@
int nhits = raw_hits.size();
hittot += nhits;
Set<SimTrackerHit> simhits = new HashSet<SimTrackerHit>();
-/*
+ /*
for (RawTrackerHit hit : raw_hits) {
- double charge = chip.decodeCharge(hit);
- aida.histogram1D("Hit charge", 100, 0., 100000.).fill(charge);
- aida.cloud1D("Number of SimTrackerHits per raw hit").fill(hit.getSimTrackerHits().size());
- double energy = 0.;
- for (SimTrackerHit shit : hit.getSimTrackerHits()) {
- simhits.add(shit);
- energy += 1000 * shit.getdEdx();
- }
- aida.histogram1D("SimTrackerHit energy", 100, 0., 10.).fill(energy);
-// aida.cloud2D("charge vs SimTrackerHit energy").fill(charge, energy);
+ double charge = chip.decodeCharge(hit);
+ aida.histogram1D("Hit charge", 100, 0., 100000.).fill(charge);
+ aida.cloud1D("Number of SimTrackerHits per raw hit").fill(hit.getSimTrackerHits().size());
+ double energy = 0.;
+ for (SimTrackerHit shit : hit.getSimTrackerHits()) {
+ simhits.add(shit);
+ energy += 1000 * shit.getdEdx();
+ }
+ aida.histogram1D("SimTrackerHit energy", 100, 0., 10.).fill(energy);
+ // aida.cloud2D("charge vs SimTrackerHit energy").fill(charge, energy);
}
-*/
+ */
double occ = 4. * ((double) nhits) / ((double) nchan);
// double occ = ((double) nhits) / ((double) nchan);
@@ -121,75 +121,77 @@
int zi = 2 * ((int) Math.round(Math.abs(z) / 20));
boolean barrel = !sensor.getName().toUpperCase().contains("ENDCAP");
- if (!barrel) continue;
+ if (!barrel) {
+ continue;
+ }
String identifier;
- String identifier2;
+ String identifier2;
if (barrel) {
identifier = "Barrel at r = " + ri;
- identifier2= "Barrel vs r";
+ identifier2 = "Barrel vs r";
} else {
identifier = "Endcap at z = " + zi;
identifier2 = "Endcap vs z";
}
if (!occMap.containsKey(identifier)) {
- occMap.put(identifier, _hf.createProfile1D(identifier, 50, 0., 120.));
+ occMap.put(identifier, _hf.createProfile1D(identifier, 50, 0., 120.));
+ }
+ if (!occMap.containsKey(identifier2)) {
+ occMap.put(identifier2, _hf.createProfile1D(identifier2, 50, 0., 120.));
}
- if (!occMap.containsKey(identifier2)) {
- occMap.put(identifier2, _hf.createProfile1D(identifier2, 50, 0., 120.));
- }
if (barrel) {
occMap.get(identifier).fill(zi, occ);
occMap.get(identifier2).fill(ri, occ);
} else {
occMap.get(identifier).fill(ri, occ);
- occMap.get(identifier2).fill(zi, occ);
+ occMap.get(identifier2).fill(zi, occ);
}
}
System.out.println("Total number of hit channels: " + hittot);
-/*
+ /*
List<List<SimTrackerHit>> simall =
- (List<List<SimTrackerHit>>) event.get(SimTrackerHit.class);
+ (List<List<SimTrackerHit>>) event.get(SimTrackerHit.class);
for (List<SimTrackerHit> simcol : simall) {
- for (SimTrackerHit hit : simcol) {
- IDetectorElement de = hit.getDetectorElement();
- Hep3Vector pos = hit.getPositionVec();
- double x = pos.x();
- double y = pos.y();
- double z = pos.z();
- double r = Math.sqrt(x * x + y * y);
- Hep3Vector p = hit.getMCParticle().getMomentum();
- double px = p.x();
- double py = p.y();
- double pz = p.z();
- double pr = Math.sqrt(px * px + py * py);
- double eta = -Math.log(Math.tan(0.5 * Math.atan2(pr, pz)));
- String detname = hit.getSubdetector().getName();
- int layer = hit.getLayer();
- aida.histogram1D("MC eta for Layer " + detname + layer, 50, -2.5, 2.5).fill(eta);
- aida.histogram1D("SimTracker Hi z for Layer " + detname + layer, 50, -120., 120.).fill(z / 10);
-// de.
- }
+ for (SimTrackerHit hit : simcol) {
+ IDetectorElement de = hit.getDetectorElement();
+ Hep3Vector pos = hit.getPositionVec();
+ double x = pos.x();
+ double y = pos.y();
+ double z = pos.z();
+ double r = Math.sqrt(x * x + y * y);
+ Hep3Vector p = hit.getMCParticle().getMomentum();
+ double px = p.x();
+ double py = p.y();
+ double pz = p.z();
+ double pr = Math.sqrt(px * px + py * py);
+ double eta = -Math.log(Math.tan(0.5 * Math.atan2(pr, pz)));
+ String detname = hit.getSubdetector().getName();
+ int layer = hit.getLayer();
+ aida.histogram1D("MC eta for Layer " + detname + layer, 50, -2.5, 2.5).fill(eta);
+ aida.histogram1D("SimTracker Hi z for Layer " + detname + layer, 50, -120., 120.).fill(z / 10);
+ // de.
+ }
}
for (MCParticle mcp : event.getMCParticles()) {
- if (mcp.getCharge() == 0) continue;
- Hep3Vector p = mcp.getMomentum();
- double px = p.x();
- double py = p.y();
- double pz = p.z();
- double pr = Math.sqrt(px * px + py * py);
- double eta = -Math.log(Math.tan(0.5 * Math.atan2(pr, pz)));
- if (mcp.getGeneratorStatus() == mcp.FINAL_STATE) {
- aida.histogram1D("MCP eta - final state",50,-2.5, 2.5).fill(eta);
- }
- if (mcp.getSimulatorStatus().isCreatedInSimulation()) {
- aida.histogram1D("MCP eta - created in simulation",50, -2.5, 2.5).fill(eta);
- }
+ if (mcp.getCharge() == 0) continue;
+ Hep3Vector p = mcp.getMomentum();
+ double px = p.x();
+ double py = p.y();
+ double pz = p.z();
+ double pr = Math.sqrt(px * px + py * py);
+ double eta = -Math.log(Math.tan(0.5 * Math.atan2(pr, pz)));
+ if (mcp.getGeneratorStatus() == mcp.FINAL_STATE) {
+ aida.histogram1D("MCP eta - final state",50,-2.5, 2.5).fill(eta);
+ }
+ if (mcp.getSimulatorStatus().isCreatedInSimulation()) {
+ aida.histogram1D("MCP eta - created in simulation",50, -2.5, 2.5).fill(eta);
+ }
}
-*/
+ */
return;
}
@@ -208,12 +210,12 @@
}
public void detectorChanged(Detector detector) {
- System.out.println("In Occupancy Driver : "+detector.getName());
+ System.out.println("In Occupancy Driver : " + detector.getName());
super.detectorChanged(detector);
// Process detectors specified by path, otherwise process entire detector
IDetectorElement detector_de = detector.getDetectorElement();
- System.out.println("In Occupancy Driver : detector_de = "+detector_de.getName());
+ System.out.println("In Occupancy Driver : detector_de = " + detector_de.getName());
for (String de_path : _process_paths) {
_process_de.add(detector_de.findDetectorElement(de_path));
}
