File indexing completed on 2021-02-14 13:30:40
0001 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0002 #include "FWCore/Framework/interface/Run.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/EventSetup.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/ServiceRegistry/interface/Service.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Utilities/interface/EDGetToken.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "FWCore/Utilities/interface/Exception.h"
0012 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0013 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0014 #include "SimTracker/TrackTriggerAssociation/interface/TTClusterAssociationMap.h"
0015 #include "SimTracker/Common/interface/TrackingParticleSelector.h"
0016 #include "DataFormats/DetId/interface/DetId.h"
0017 #include "DataFormats/Common/interface/Ptr.h"
0018 #include "DataFormats/Common/interface/Handle.h"
0019 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0020 #include "DataFormats/L1TrackTrigger/interface/TTDTC.h"
0021 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0022 #include "DataFormats/GeometrySurface/interface/Plane.h"
0023 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0024
0025 #include "L1Trigger/TrackerDTC/interface/Setup.h"
0026
0027 #include <TProfile.h>
0028 #include <TProfile2D.h>
0029 #include <TH1F.h>
0030 #include <TH2F.h>
0031 #include <TEfficiency.h>
0032
0033 #include <vector>
0034 #include <map>
0035 #include <utility>
0036 #include <set>
0037 #include <algorithm>
0038 #include <cmath>
0039 #include <numeric>
0040 #include <array>
0041 #include <initializer_list>
0042 #include <sstream>
0043
0044 using namespace std;
0045 using namespace edm;
0046
0047 namespace trackerDTC {
0048
0049
0050 typedef TTClusterAssociationMap<Ref_Phase2TrackerDigi_> TTClusterAssMap;
0051 typedef edm::Ptr<TrackingParticle> TPPtr;
0052
0053 enum Resolution { R, Phi, Z, NumResolution };
0054 constexpr initializer_list<Resolution> AllResolution = {R, Phi, Z};
0055 constexpr auto NameResolution = {"R", "Phi", "Z"};
0056 inline string name(Resolution r) { return string(*(NameResolution.begin() + r)); }
0057
0058 enum Efficiency { Phi0, Pt, InvPt, D0, Z0, Eta, NumEfficiency };
0059 constexpr initializer_list<Efficiency> AllEfficiency = {Phi0, Pt, InvPt, D0, Z0, Eta};
0060 constexpr auto NameEfficiency = {"Phi0", "Pt", "InvPt", "D0", "Z0", "Eta"};
0061 inline string name(Efficiency e) { return string(*(NameEfficiency.begin() + e)); }
0062
0063
0064
0065
0066
0067
0068 class Analyzer : public one::EDAnalyzer<one::WatchRuns, one::SharedResources> {
0069 public:
0070 Analyzer(const ParameterSet& iConfig);
0071 void beginJob() override {}
0072 void beginRun(const Run& iEvent, const EventSetup& iSetup) override;
0073 void analyze(const Event& iEvent, const EventSetup& iSetup) override;
0074 void endRun(const Run& iEvent, const EventSetup& iSetup) override {}
0075 void endJob() override;
0076
0077 private:
0078
0079 void configTPSelector();
0080
0081 void bookHistograms();
0082
0083 void assoc(const Handle<TTStubDetSetVec>&, const Handle<TTClusterAssMap>&, map<TPPtr, set<TTStubRef>>&);
0084
0085 void convert(const map<TPPtr, set<TTStubRef>>&, map<TTStubRef, set<TPPtr>>&);
0086
0087 bool reconstructable(const set<TTStubRef>& ttStubRefs) const;
0088
0089 bool select(const TrackingParticle& tp) const;
0090
0091 void fill(const TPPtr& tpPtr, const vector<TH1F*> th1fs) const;
0092
0093 void analyzeStubs(const TTDTC*, const TTDTC*, const map<TTStubRef, set<TPPtr>>&, map<TPPtr, set<TTStubRef>>&);
0094
0095 void analyzeStream(const TTDTC::Stream& stream, int region, int channel, int& sum, TH2F* th2f);
0096
0097 int layerId(const TTStubRef& ttStubRef) const;
0098
0099 void analyzeTPs(const map<TPPtr, set<TTStubRef>>& mapTPsStubs);
0100
0101 void endJobMC();
0102
0103 void endJobDTC();
0104
0105
0106 EDGetTokenT<TTDTC> getTokenTTDTCAccepted_;
0107
0108 EDGetTokenT<TTDTC> getTokenTTDTCLost_;
0109
0110 EDGetTokenT<TTStubDetSetVec> getTokenTTStubDetSetVec_;
0111
0112 EDGetTokenT<TTClusterAssMap> getTokenTTClusterAssMap_;
0113
0114 ESGetToken<Setup, SetupRcd> esGetToken_;
0115
0116 Setup setup_;
0117
0118 TrackingParticleSelector tpSelector_;
0119
0120 bool useMCTruth_;
0121
0122 bool hybrid_;
0123
0124
0125
0126 TProfile* profMC_;
0127 TProfile* profDTC_;
0128 TProfile* profChannel_;
0129 TH1F* hisChannel_;
0130 TH2F* hisRZStubs_;
0131 TH2F* hisRZStubsLost_;
0132 TH2F* hisRZStubsEff_;
0133 vector<TH1F*> hisResolution_;
0134 vector<TProfile2D*> profResolution_;
0135 vector<TH1F*> hisEff_;
0136 vector<TH1F*> hisEffMC_;
0137 vector<TEfficiency*> eff_;
0138
0139
0140 stringstream log_;
0141 };
0142
0143 Analyzer::Analyzer(const ParameterSet& iConfig)
0144 : useMCTruth_(iConfig.getParameter<bool>("UseMCTruth")), hybrid_(iConfig.getParameter<bool>("UseHybrid")) {
0145 usesResource("TFileService");
0146
0147 const auto& inputTagAccepted = iConfig.getParameter<InputTag>("InputTagAccepted");
0148 const auto& inputTagLost = iConfig.getParameter<InputTag>("InputTagLost");
0149 getTokenTTDTCAccepted_ = consumes<TTDTC>(inputTagAccepted);
0150 getTokenTTDTCLost_ = consumes<TTDTC>(inputTagLost);
0151 if (useMCTruth_) {
0152 const auto& inputTagTTStubDetSetVec = iConfig.getParameter<InputTag>("InputTagTTStubDetSetVec");
0153 const auto& inputTagTTClusterAssMap = iConfig.getParameter<InputTag>("InputTagTTClusterAssMap");
0154 getTokenTTStubDetSetVec_ = consumes<TTStubDetSetVec>(inputTagTTStubDetSetVec);
0155 getTokenTTClusterAssMap_ = consumes<TTClusterAssMap>(inputTagTTClusterAssMap);
0156 }
0157
0158 esGetToken_ = esConsumes<Setup, SetupRcd, Transition::BeginRun>();
0159
0160 log_.setf(ios::fixed, ios::floatfield);
0161 log_.precision(4);
0162 }
0163
0164 void Analyzer::beginRun(const Run& iEvent, const EventSetup& iSetup) {
0165
0166 setup_ = iSetup.getData(esGetToken_);
0167
0168 configTPSelector();
0169
0170 bookHistograms();
0171 }
0172
0173 void Analyzer::analyze(const Event& iEvent, const EventSetup& iSetup) {
0174
0175 map<TTStubRef, set<TPPtr>> mapAllStubsTPs;
0176 if (useMCTruth_) {
0177 Handle<TTStubDetSetVec> handleTTStubDetSetVec;
0178 iEvent.getByToken<TTStubDetSetVec>(getTokenTTStubDetSetVec_, handleTTStubDetSetVec);
0179 Handle<TTClusterAssMap> handleTTClusterAssMap;
0180 iEvent.getByToken<TTClusterAssMap>(getTokenTTClusterAssMap_, handleTTClusterAssMap);
0181
0182 map<TPPtr, set<TTStubRef>> mapAllTPsAllStubs;
0183 assoc(handleTTStubDetSetVec, handleTTClusterAssMap, mapAllTPsAllStubs);
0184
0185 convert(mapAllTPsAllStubs, mapAllStubsTPs);
0186 }
0187
0188 Handle<TTDTC> handleTTDTCAccepted;
0189 iEvent.getByToken<TTDTC>(getTokenTTDTCAccepted_, handleTTDTCAccepted);
0190 Handle<TTDTC> handleTTDTCLost;
0191 iEvent.getByToken<TTDTC>(getTokenTTDTCLost_, handleTTDTCLost);
0192 map<TPPtr, set<TTStubRef>> mapTPsTTStubs;
0193
0194 analyzeStubs(handleTTDTCAccepted.product(), handleTTDTCLost.product(), mapAllStubsTPs, mapTPsTTStubs);
0195
0196 analyzeTPs(mapTPsTTStubs);
0197 }
0198
0199 void Analyzer::endJob() {
0200
0201 TH2F th2f("", ";;", 400, -300, 300., 400, 0., 120.);
0202 th2f.Add(hisRZStubsLost_);
0203 th2f.Add(hisRZStubs_);
0204 hisRZStubsEff_->Add(hisRZStubsLost_);
0205 hisRZStubsEff_->Divide(&th2f);
0206
0207 if (useMCTruth_) {
0208 for (Efficiency e : AllEfficiency) {
0209 eff_[e]->SetPassedHistogram(*hisEff_[e], "f");
0210 eff_[e]->SetTotalHistogram(*hisEffMC_[e], "f");
0211 }
0212 }
0213
0214 endJobMC();
0215
0216 endJobDTC();
0217 log_ << "=============================================================" << endl;
0218 LogPrint("L1Trigger/TrackerDTC") << log_.str();
0219 }
0220
0221
0222 void Analyzer::assoc(const Handle<TTStubDetSetVec>& handleTTStubDetSetVec,
0223 const Handle<TTClusterAssMap>& handleTTClusterAssMap,
0224 map<TPPtr, set<TTStubRef>>& mapTPsStubs) {
0225 int nStubs(0);
0226 int nStubsMatched(0);
0227 for (TTStubDetSetVec::const_iterator ttModule = handleTTStubDetSetVec->begin();
0228 ttModule != handleTTStubDetSetVec->end();
0229 ttModule++) {
0230 nStubs += ttModule->size();
0231 for (TTStubDetSet::const_iterator ttStub = ttModule->begin(); ttStub != ttModule->end(); ttStub++) {
0232 set<TPPtr> tpPtrs;
0233 for (unsigned int iClus = 0; iClus < 2; iClus++) {
0234 const vector<TPPtr>& assocPtrs = handleTTClusterAssMap->findTrackingParticlePtrs(ttStub->clusterRef(iClus));
0235 copy_if(assocPtrs.begin(), assocPtrs.end(), inserter(tpPtrs, tpPtrs.begin()), [](const TPPtr& tpPtr) {
0236 return tpPtr.isNonnull();
0237 });
0238 }
0239 for (const TPPtr& tpPtr : tpPtrs)
0240 mapTPsStubs[tpPtr].emplace(makeRefTo(handleTTStubDetSetVec, ttStub));
0241 if (!tpPtrs.empty())
0242 nStubsMatched++;
0243 }
0244 }
0245 profMC_->Fill(1, nStubs / (double)setup_.numRegions());
0246 profMC_->Fill(2, nStubsMatched / (double)setup_.numRegions());
0247 }
0248
0249
0250 void Analyzer::convert(const map<TPPtr, set<TTStubRef>>& mapTPsStubs, map<TTStubRef, set<TPPtr>>& mapStubsTPs) {
0251 int nTPsReco(0);
0252 int nTPsEff(0);
0253 for (const auto& mapTPStubs : mapTPsStubs) {
0254 if (!reconstructable(mapTPStubs.second))
0255 continue;
0256 nTPsReco++;
0257 const bool useForAlgEff = select(*mapTPStubs.first.get());
0258 if (useForAlgEff) {
0259 nTPsEff++;
0260 fill(mapTPStubs.first, hisEffMC_);
0261 for (const TTStubRef& ttStubRef : mapTPStubs.second)
0262 mapStubsTPs[ttStubRef].insert(mapTPStubs.first);
0263 }
0264 }
0265 profMC_->Fill(3, nTPsReco);
0266 profMC_->Fill(4, nTPsEff);
0267 }
0268
0269
0270 bool Analyzer::reconstructable(const set<TTStubRef>& ttStubRefs) const {
0271 const TrackerGeometry* trackerGeometry = setup_.trackerGeometry();
0272 const TrackerTopology* trackerTopology = setup_.trackerTopology();
0273 set<int> hitPattern;
0274 set<int> hitPatternPS;
0275 for (const TTStubRef& ttStubRef : ttStubRefs) {
0276 const DetId detId = ttStubRef->getDetId();
0277 const bool barrel = detId.subdetId() == StripSubdetector::TOB;
0278 const bool psModule = trackerGeometry->getDetectorType(detId) == TrackerGeometry::ModuleType::Ph2PSP;
0279 const int layerId = barrel ? trackerTopology->layer(detId) : trackerTopology->tidWheel(detId) + 10;
0280 hitPattern.insert(layerId);
0281 if (psModule)
0282 hitPatternPS.insert(layerId);
0283 }
0284 return (int)hitPattern.size() >= setup_.tpMinLayers() && (int)hitPatternPS.size() >= setup_.tpMinLayersPS();
0285 }
0286
0287
0288 bool Analyzer::select(const TrackingParticle& tp) const {
0289 const bool selected = tpSelector_(tp);
0290 const double cot = sinh(tp.eta());
0291 const double s = sin(tp.phi());
0292 const double c = cos(tp.phi());
0293 const TrackingParticle::Point& v = tp.vertex();
0294 const double z0 = v.z() - (v.x() * c + v.y() * s) * cot;
0295 const double d0 = v.x() * s - v.y() * c;
0296 return selected && (fabs(d0) < setup_.tpMaxD0()) && (fabs(z0) < setup_.tpMaxVertZ());
0297 }
0298
0299
0300 void Analyzer::fill(const TPPtr& tpPtr, const vector<TH1F*> th1fs) const {
0301 const double s = sin(tpPtr->phi());
0302 const double c = cos(tpPtr->phi());
0303 const TrackingParticle::Point& v = tpPtr->vertex();
0304 const vector<double> x = {tpPtr->phi(),
0305 tpPtr->pt(),
0306 tpPtr->charge() / tpPtr->pt(),
0307 v.x() * s - v.y() * c,
0308 v.z() - (v.x() * c + v.y() * s) * sinh(tpPtr->eta()),
0309 tpPtr->eta()};
0310 for (Efficiency e : AllEfficiency)
0311 th1fs[e]->Fill(x[e]);
0312 }
0313
0314
0315 void Analyzer::analyzeStubs(const TTDTC* accepted,
0316 const TTDTC* lost,
0317 const map<TTStubRef, set<TPPtr>>& mapStubsTPs,
0318 map<TPPtr, set<TTStubRef>>& mapTPsStubs) {
0319 for (int region = 0; region < setup_.numRegions(); region++) {
0320 int nStubs(0);
0321 int nLost(0);
0322 for (int channel = 0; channel < setup_.numDTCsPerTFP(); channel++) {
0323 const TTDTC::Stream& stream = accepted->stream(region, channel);
0324 hisChannel_->Fill(stream.size());
0325 profChannel_->Fill(region * setup_.numDTCsPerTFP() + channel, stream.size());
0326 for (const TTDTC::Frame& frame : stream) {
0327 if (frame.first.isNull())
0328 continue;
0329 const auto it = mapStubsTPs.find(frame.first);
0330 if (it == mapStubsTPs.end())
0331 continue;
0332 for (const TPPtr& tp : it->second)
0333 mapTPsStubs[tp].insert(frame.first);
0334 }
0335 analyzeStream(stream, region, channel, nStubs, hisRZStubs_);
0336 analyzeStream(lost->stream(region, channel), region, channel, nLost, hisRZStubsLost_);
0337 }
0338 profDTC_->Fill(1, nStubs);
0339 profDTC_->Fill(2, nLost);
0340 }
0341 }
0342
0343
0344 void Analyzer::analyzeStream(const TTDTC::Stream& stream, int region, int channel, int& sum, TH2F* th2f) {
0345 for (const TTDTC::Frame& frame : stream) {
0346 if (frame.first.isNull())
0347 continue;
0348 sum++;
0349 const GlobalPoint& pos = setup_.stubPos(hybrid_, frame, region, channel);
0350 const GlobalPoint& ttPos = setup_.stubPos(frame.first);
0351 const vector<double> resolutions = {
0352 ttPos.perp() - pos.perp(), deltaPhi(ttPos.phi() - pos.phi()), ttPos.z() - pos.z()};
0353 for (Resolution r : AllResolution) {
0354 hisResolution_[r]->Fill(resolutions[r]);
0355 profResolution_[r]->Fill(ttPos.z(), ttPos.perp(), abs(resolutions[r]));
0356 }
0357 th2f->Fill(ttPos.z(), ttPos.perp());
0358
0359 if (!hybrid_)
0360 continue;
0361 const vector<int>& encodingLayerId = setup_.encodingLayerId(channel);
0362 const auto it = find(encodingLayerId.begin(), encodingLayerId.end(), layerId(frame.first));
0363 if (it == encodingLayerId.end())
0364 throw cms::Exception("LogicError") << "Stub send from a DTC which is not connected to stub's layer.";
0365 }
0366 }
0367
0368
0369 int Analyzer::layerId(const TTStubRef& ttStubRef) const {
0370 const TrackerTopology* trackerTopology = setup_.trackerTopology();
0371 const DetId detId = ttStubRef->getDetId() + setup_.offsetDetIdDSV();
0372 const bool barrel = detId.subdetId() == StripSubdetector::TOB;
0373 return barrel ? trackerTopology->layer(detId) : trackerTopology->tidWheel(detId) + setup_.offsetLayerDisks();
0374 }
0375
0376
0377 void Analyzer::analyzeTPs(const map<TPPtr, set<TTStubRef>>& mapTPsStubs) {
0378 int nTPs(0);
0379 for (const auto& mapTPStubs : mapTPsStubs) {
0380 if (!reconstructable(mapTPStubs.second))
0381 continue;
0382 nTPs++;
0383 fill(mapTPStubs.first, hisEff_);
0384 }
0385 profDTC_->Fill(3, nTPs);
0386 }
0387
0388
0389 void Analyzer::endJobMC() {
0390 const double numStubs = profMC_->GetBinContent(1);
0391 const double numStubsMatched = profMC_->GetBinContent(2);
0392 const double numTPsReco = profMC_->GetBinContent(3);
0393 const double numTPsEff = profMC_->GetBinContent(4);
0394 const double errStubs = profMC_->GetBinError(1);
0395 const double errStubsMatched = profMC_->GetBinError(2);
0396 const double errTPsReco = profMC_->GetBinError(3);
0397 const double errTPsEff = profMC_->GetBinError(4);
0398 const vector<double> nums = {numStubs, numStubsMatched, numTPsReco, numTPsEff};
0399 const vector<double> errs = {errStubs, errStubsMatched, errTPsReco, errTPsEff};
0400 const int wNums = ceil(log10(*max_element(nums.begin(), nums.end()))) + 5;
0401 const int wErrs = ceil(log10(*max_element(errs.begin(), errs.end()))) + 5;
0402 log_ << "=============================================================" << endl;
0403 log_ << " MC SUMMARY " << endl;
0404 log_ << "number of stubs per TFP = " << setw(wNums) << numStubs << " +- " << setw(wErrs) << errStubs
0405 << endl;
0406 log_ << "number of matched stubs per TFP = " << setw(wNums) << numStubsMatched << " +- " << setw(wErrs)
0407 << errStubsMatched << endl;
0408 log_ << "number of TPs per TFP = " << setw(wNums) << numTPsReco << " +- " << setw(wErrs) << errTPsReco
0409 << endl;
0410 log_ << "number of TPs for eff per TFP = " << setw(wNums) << numTPsEff << " +- " << setw(wErrs) << errTPsEff
0411 << endl;
0412 }
0413
0414
0415 void Analyzer::endJobDTC() {
0416 const double numStubs = profDTC_->GetBinContent(1);
0417 const double numStubsLost = profDTC_->GetBinContent(2);
0418 const double numTPs = profDTC_->GetBinContent(3);
0419 const double errStubs = profDTC_->GetBinError(1);
0420 const double errStubsLost = profDTC_->GetBinError(2);
0421 const double totalTPs = profMC_->GetBinContent(4);
0422 const double eff = numTPs / totalTPs;
0423 const double errEff = sqrt(eff * (1. - eff) / totalTPs);
0424 const vector<double> nums = {numStubs, numStubsLost};
0425 const vector<double> errs = {errStubs, errStubsLost};
0426 const int wNums = ceil(log10(*max_element(nums.begin(), nums.end()))) + 5;
0427 const int wErrs = ceil(log10(*max_element(errs.begin(), errs.end()))) + 5;
0428 log_ << "=============================================================" << endl;
0429 log_ << " DTC SUMMARY " << endl;
0430 log_ << "number of stubs per TFP = " << setw(wNums) << numStubs << " +- " << setw(wErrs) << errStubs << endl;
0431 log_ << "number of lost stubs per TFP = " << setw(wNums) << numStubsLost << " +- " << setw(wErrs) << errStubsLost
0432 << endl;
0433 log_ << " max tracking efficiency = " << setw(wNums) << eff << " +- " << setw(wErrs) << errEff << endl;
0434 }
0435
0436
0437 void Analyzer::configTPSelector() {
0438 const double ptMin = hybrid_ ? setup_.hybridMinPt() : setup_.minPt();
0439 constexpr double ptMax = 9999999999.;
0440 const double etaMax = setup_.tpMaxEta();
0441 const double tip = setup_.tpMaxVertR();
0442 const double lip = setup_.tpMaxVertZ();
0443 constexpr int minHit = 0;
0444 constexpr bool signalOnly = true;
0445 constexpr bool intimeOnly = true;
0446 constexpr bool chargedOnly = true;
0447 constexpr bool stableOnly = false;
0448 tpSelector_ = TrackingParticleSelector(
0449 ptMin, ptMax, -etaMax, etaMax, tip, lip, minHit, signalOnly, intimeOnly, chargedOnly, stableOnly);
0450 }
0451
0452
0453 void Analyzer::bookHistograms() {
0454 Service<TFileService> fs;
0455 TFileDirectory dir;
0456
0457 dir = fs->mkdir("MC");
0458 profMC_ = dir.make<TProfile>("Counts", ";", 4, 0.5, 4.5);
0459 profMC_->GetXaxis()->SetBinLabel(1, "Stubs");
0460 profMC_->GetXaxis()->SetBinLabel(2, "Matched Stubs");
0461 profMC_->GetXaxis()->SetBinLabel(3, "reco TPs");
0462 profMC_->GetXaxis()->SetBinLabel(4, "eff TPs");
0463 constexpr array<int, NumEfficiency> binsEff{{9 * 8, 10, 16, 10, 30, 24}};
0464 constexpr array<pair<double, double>, NumEfficiency> rangesEff{
0465 {{-M_PI, M_PI}, {0., 100.}, {-1. / 3., 1. / 3.}, {-5., 5.}, {-15., 15.}, {-2.4, 2.4}}};
0466 if (useMCTruth_) {
0467 hisEffMC_.reserve(NumEfficiency);
0468 for (Efficiency e : AllEfficiency)
0469 hisEffMC_.emplace_back(
0470 dir.make<TH1F>(("HisTP" + name(e)).c_str(), ";", binsEff[e], rangesEff[e].first, rangesEff[e].second));
0471 }
0472
0473 dir = fs->mkdir("DTC");
0474 profDTC_ = dir.make<TProfile>("Counts", ";", 3, 0.5, 3.5);
0475 profDTC_->GetXaxis()->SetBinLabel(1, "Stubs");
0476 profDTC_->GetXaxis()->SetBinLabel(2, "Lost Stubs");
0477 profDTC_->GetXaxis()->SetBinLabel(3, "TPs");
0478
0479 constexpr int maxOcc = 180;
0480 const int numChannels = setup_.numDTCs() * setup_.numOverlappingRegions();
0481 hisChannel_ = dir.make<TH1F>("His Channel Occupancy", ";", maxOcc, -.5, maxOcc - .5);
0482 profChannel_ = dir.make<TProfile>("Prof Channel Occupancy", ";", numChannels, -.5, numChannels - .5);
0483
0484 if (useMCTruth_) {
0485 dir = fs->mkdir("DTC/Effi");
0486 hisEff_.reserve(NumEfficiency);
0487 for (Efficiency e : AllEfficiency)
0488 hisEff_.emplace_back(
0489 dir.make<TH1F>(("HisTP" + name(e)).c_str(), ";", binsEff[e], rangesEff[e].first, rangesEff[e].second));
0490 eff_.reserve(NumEfficiency);
0491 for (Efficiency e : AllEfficiency)
0492 eff_.emplace_back(
0493 dir.make<TEfficiency>(("Eff" + name(e)).c_str(), ";", binsEff[e], rangesEff[e].first, rangesEff[e].second));
0494 }
0495
0496 dir = fs->mkdir("DTC/Loss");
0497 constexpr int bins = 400;
0498 constexpr double maxZ = 300.;
0499 constexpr double maxR = 120.;
0500 hisRZStubs_ = dir.make<TH2F>("RZ Stubs", ";;", bins, -maxZ, maxZ, bins, 0., maxR);
0501 hisRZStubsLost_ = dir.make<TH2F>("RZ Stubs Lost", ";;", bins, -maxZ, maxZ, bins, 0., maxR);
0502 hisRZStubsEff_ = dir.make<TH2F>("RZ Stubs Eff", ";;", bins, -maxZ, maxZ, bins, 0., maxR);
0503
0504 dir = fs->mkdir("DTC/Res");
0505 constexpr array<double, NumResolution> ranges{{.2, .0001, .5}};
0506 constexpr int binsHis = 100;
0507 hisResolution_.reserve(NumResolution);
0508 profResolution_.reserve(NumResolution);
0509 for (Resolution r : AllResolution) {
0510 hisResolution_.emplace_back(dir.make<TH1F>(("HisRes" + name(r)).c_str(), ";", binsHis, -ranges[r], ranges[r]));
0511 profResolution_.emplace_back(
0512 dir.make<TProfile2D>(("ProfRes" + name(r)).c_str(), ";;", bins, -maxZ, maxZ, bins, 0., maxR));
0513 }
0514 }
0515
0516 }
0517
0518 DEFINE_FWK_MODULE(trackerDTC::Analyzer);