/
githubmirror
/
cmssw
Обзор
Документация
Войти
/
githubmirror
/
cmssw
Код
Запросы
0
Пакеты
0
Релизы
0
Аналитика
Безопасность
master
L1Trigger/TrackerTFP/test/AnalyzerHT.cc
322 строки
15 KB
tschuh
quinnanm 1
16 апр 2025, 23:34
16 апр 2025, 23:34
e51e86a
Код
Авторство
О чём код?
#include "FWCore/Framework/interface/one/EDAnalyzer.h" #include "FWCore/Framework/interface/Run.h" #include "FWCore/Framework/interface/Event.h" #include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/Framework/interface/MakerMacros.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/ServiceRegistry/interface/Service.h" #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "FWCore/Utilities/interface/EDGetToken.h" #include "FWCore/Utilities/interface/InputTag.h" #include "FWCore/Utilities/interface/Exception.h" #include "CommonTools/UtilAlgos/interface/TFileService.h" #include "DataFormats/Common/interface/Handle.h" #include "SimTracker/TrackTriggerAssociation/interface/StubAssociation.h" #include "L1Trigger/TrackTrigger/interface/Setup.h" #include "L1Trigger/TrackerTFP/interface/DataFormats.h" #include <TProfile.h> #include <TH1F.h> #include <TEfficiency.h> #include <vector> #include <deque> #include <set> #include <cmath> #include <numeric> #include <sstream> namespace trackerTFP { /*! \class trackerTFP::AnalyzerHT * \brief Class to analyze hardware like structured TTStub Collection generated by Geometric Processor * \author Thomas Schuh * \date 2020, Apr */ class AnalyzerHT : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> { public: AnalyzerHT(const edm::ParameterSet& iConfig); void beginJob() override {} void beginRun(const edm::Run& iEvent, const edm::EventSetup& iSetup) override; void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override; void endRun(const edm::Run& iEvent, const edm::EventSetup& iSetup) override {} void endJob() override; private: // void formTracks(const tt::StreamStub& stream, std::vector<std::vector<TTStubRef>>& tracks) const; // void associate(const std::vector<std::vector<TTStubRef>>& tracks, const tt::StubAssociation* ass, std::set<TPPtr>& tps, int& sum) const; // ED input token of stubs edm::EDGetTokenT<tt::StreamsStub> edGetToken_; // ED input token of TTStubRef to TPPtr association for tracking efficiency edm::EDGetTokenT<tt::StubAssociation> edGetTokenSelection_; // ED input token of TTStubRef to recontructable TPPtr association edm::EDGetTokenT<tt::StubAssociation> edGetTokenReconstructable_; // Setup token edm::ESGetToken<tt::Setup, tt::SetupRcd> esGetTokenSetup_; // DataFormats token edm::ESGetToken<DataFormats, DataFormatsRcd> esGetTokenDataFormats_; // stores, calculates and provides run-time constants const tt::Setup* setup_ = nullptr; // helper class to extract structured data from tt::Frames const DataFormats* dataFormats_ = nullptr; // enables analyze of TPs bool useMCTruth_; // int nEvents_ = 0; // Histograms TProfile* prof_; TProfile* profChannel_; TH1F* hisChannel_; TH1F* hisLayers_; TH1F* hisNumLayers_; TProfile* profNumLayers_; TH1F* hisEffEta_; TH1F* hisEffEtaTotal_; TEfficiency* effEta_; TH1F* hisEffZT_; TH1F* hisEffZTTotal_; TEfficiency* effZT_; TH1F* hisEffInv2R_; TH1F* hisEffInv2RTotal_; TEfficiency* effInv2R_; TH1F* hisEffPT_; TH1F* hisEffPTTotal_; TEfficiency* effPT_; // printout std::stringstream log_; }; AnalyzerHT::AnalyzerHT(const edm::ParameterSet& iConfig) : useMCTruth_(iConfig.getParameter<bool>("UseMCTruth")) { usesResource("TFileService"); // book in- and output ED products const std::string& label = iConfig.getParameter<std::string>("OutputLabelHT"); const std::string& branch = iConfig.getParameter<std::string>("BranchStubs"); edGetToken_ = consumes<tt::StreamsStub>(edm::InputTag(label, branch)); if (useMCTruth_) { const auto& inputTagSelecttion = iConfig.getParameter<edm::InputTag>("InputTagSelection"); const auto& inputTagReconstructable = iConfig.getParameter<edm::InputTag>("InputTagReconstructable"); edGetTokenSelection_ = consumes<tt::StubAssociation>(inputTagSelecttion); edGetTokenReconstructable_ = consumes<tt::StubAssociation>(inputTagReconstructable); } // book ES products esGetTokenSetup_ = esConsumes<edm::Transition::BeginRun>(); esGetTokenDataFormats_ = esConsumes<edm::Transition::BeginRun>(); // log config log_.setf(std::ios::fixed, std::ios::floatfield); log_.precision(4); } void AnalyzerHT::beginRun(const edm::Run& iEvent, const edm::EventSetup& iSetup) { // helper class to store configurations setup_ = &iSetup.getData(esGetTokenSetup_); // helper class to extract structured data from tt::Frames dataFormats_ = &iSetup.getData(esGetTokenDataFormats_); // book histograms edm::Service<TFileService> fs; TFileDirectory dir; dir = fs->mkdir("HT"); prof_ = dir.make<TProfile>("Counts", ";", 9, 0.5, 9.5); prof_->GetXaxis()->SetBinLabel(1, "Stubs"); prof_->GetXaxis()->SetBinLabel(2, "Tracks"); prof_->GetXaxis()->SetBinLabel(3, "Truncated Tracks"); prof_->GetXaxis()->SetBinLabel(4, "Matched Tracks"); prof_->GetXaxis()->SetBinLabel(5, "All Tracks"); prof_->GetXaxis()->SetBinLabel(6, "Found TPs"); prof_->GetXaxis()->SetBinLabel(7, "Found selected TPs"); prof_->GetXaxis()->SetBinLabel(8, "Truncated TPs"); prof_->GetXaxis()->SetBinLabel(9, "All TPs"); // Efficiencies hisEffEtaTotal_ = dir.make<TH1F>("HisTPEtaTotal", ";", 48, -2.4, 2.4); hisEffEta_ = dir.make<TH1F>("HisTPEta", ";", 48, -2.4, 2.4); effEta_ = dir.make<TEfficiency>("EffEta", ";", 48, -2.4, 2.4); const double rangeZT = dataFormats_->format(Variable::zT, Process::dr).range(); const int zTBins = setup_->gpNumBinsZT(); hisEffZTTotal_ = dir.make<TH1F>("HisTPZTTotal", ";", zTBins, -rangeZT / 2, rangeZT / 2); hisEffZT_ = dir.make<TH1F>("HisTPZT", ";", zTBins, -rangeZT / 2, rangeZT / 2); effZT_ = dir.make<TEfficiency>("EffZT", ";", zTBins, -rangeZT / 2, rangeZT / 2); const double rangeInv2R = dataFormats_->format(Variable::inv2R, Process::dr).range(); const int inv2RBins = (setup_->htNumBinsInv2R() + 2) * 2; hisEffInv2R_ = dir.make<TH1F>("HisTPInv2R", ";", inv2RBins, -rangeInv2R / 2., rangeInv2R / 2.); hisEffInv2RTotal_ = dir.make<TH1F>("HisTPInv2RTotal", ";", inv2RBins, -rangeInv2R / 2., rangeInv2R / 2.); effInv2R_ = dir.make<TEfficiency>("EffInv2R", ";", inv2RBins, -rangeInv2R / 2., rangeInv2R / 2.); hisEffPT_ = dir.make<TH1F>("HisTPPT", ";", 100, 0, 100); hisEffPTTotal_ = dir.make<TH1F>("HisTPPTTotal", ";", 100, 0, 100); effPT_ = dir.make<TEfficiency>("EffPT", ";", 100, 0, 100); // binInv2R occupancy constexpr int maxOcc = 180; const int numChannel = dataFormats_->numChannel(Process::ht); hisChannel_ = dir.make<TH1F>("His binInv2R Occupancy", ";", maxOcc, -.5, maxOcc - .5); profChannel_ = dir.make<TProfile>("Prof binInv2R Occupancy", ";", numChannel, -.5, numChannel - .5); // layers hisLayers_ = dir.make<TH1F>("HisLayers", ";", 8, 0, 8); hisNumLayers_ = dir.make<TH1F>("HisNumLayers", ";", 9, 0, 9); profNumLayers_ = dir.make<TProfile>("Prof NumLayers", ";", 32, 0, 2.4); } void AnalyzerHT::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) { auto fill = [this](const TPPtr& tpPtr, TH1F* hisEta, TH1F* hisZT, TH1F* hisInv2R, TH1F* hisPT) { const double tpPhi0 = tpPtr->phi(); const double tpCot = sinh(tpPtr->eta()); const math::XYZPointD& v = tpPtr->vertex(); const double tpZ0 = v.z() - tpCot * (v.x() * cos(tpPhi0) + v.y() * sin(tpPhi0)); const double tpZT = tpZ0 + tpCot * setup_->chosenRofZ(); hisEta->Fill(tpPtr->eta()); hisZT->Fill(tpZT); hisInv2R->Fill(tpPtr->charge() / tpPtr->pt() * setup_->invPtToDphi()); hisPT->Fill(tpPtr->pt()); }; // read in ht products edm::Handle<tt::StreamsStub> handleStubs; iEvent.getByToken<tt::StreamsStub>(edGetToken_, handleStubs); // read in MCTruth const tt::StubAssociation* selection = nullptr; const tt::StubAssociation* reconstructable = nullptr; if (useMCTruth_) { edm::Handle<tt::StubAssociation> handleSelection; iEvent.getByToken<tt::StubAssociation>(edGetTokenSelection_, handleSelection); selection = handleSelection.product(); prof_->Fill(9, selection->numTPs()); edm::Handle<tt::StubAssociation> handleReconstructable; iEvent.getByToken<tt::StubAssociation>(edGetTokenReconstructable_, handleReconstructable); reconstructable = handleReconstructable.product(); for (const auto& p : selection->getTrackingParticleToTTStubsMap()) fill(p.first, hisEffEtaTotal_, hisEffZTTotal_, hisEffInv2RTotal_, hisEffPTTotal_); } // analyze ht products and associate found tracks with reconstrucable TrackingParticles std::set<TPPtr> tpPtrs; std::set<TPPtr> tpPtrsSelection; int allMatched(0); int allTracks(0); for (int region = 0; region < setup_->numRegions(); region++) { int nStubs(0); int nTracks(0); for (int channel = 0; channel < dataFormats_->numChannel(Process::ht); channel++) { const int index = region * dataFormats_->numChannel(Process::ht) + channel; const tt::StreamStub& accepted = handleStubs->at(index); hisChannel_->Fill(accepted.size()); profChannel_->Fill(channel, accepted.size()); nStubs += accepted.size(); std::vector<std::vector<TTStubRef>> tracks; formTracks(accepted, tracks); nTracks += tracks.size(); allTracks += tracks.size(); if (!useMCTruth_) continue; int tmp(0); associate(tracks, selection, tpPtrsSelection, tmp); associate(tracks, reconstructable, tpPtrs, allMatched); } prof_->Fill(1, nStubs); prof_->Fill(2, nTracks); } for (const TPPtr& tpPtr : tpPtrsSelection) fill(tpPtr, hisEffEta_, hisEffZT_, hisEffInv2R_, hisEffPT_); prof_->Fill(4, allMatched); prof_->Fill(5, allTracks); prof_->Fill(6, tpPtrs.size()); prof_->Fill(7, tpPtrsSelection.size()); nEvents_++; } void AnalyzerHT::endJob() { if (nEvents_ == 0) return; // effi effEta_->SetPassedHistogram(*hisEffEta_, "f"); effEta_->SetTotalHistogram(*hisEffEtaTotal_, "f"); effZT_->SetPassedHistogram(*hisEffZT_, "f"); effZT_->SetTotalHistogram(*hisEffZTTotal_, "f"); effInv2R_->SetPassedHistogram(*hisEffInv2R_, "f"); effInv2R_->SetTotalHistogram(*hisEffInv2RTotal_, "f"); effPT_->SetPassedHistogram(*hisEffPT_, "f"); effPT_->SetTotalHistogram(*hisEffPTTotal_, "f"); // printout HT summary const double totalTPs = prof_->GetBinContent(9); const double numStubs = prof_->GetBinContent(1); const double numTracks = prof_->GetBinContent(2); const double totalTracks = prof_->GetBinContent(5); const double numTracksMatched = prof_->GetBinContent(4); const double numTPsAll = prof_->GetBinContent(6); const double numTPsEff = prof_->GetBinContent(7); const double errStubs = prof_->GetBinError(1); const double errTracks = prof_->GetBinError(2); const double fracFake = (totalTracks - numTracksMatched) / totalTracks; const double fracDup = (numTracksMatched - numTPsAll) / totalTracks; const double eff = numTPsEff / totalTPs; const double errEff = sqrt(eff * (1. - eff) / totalTPs / nEvents_); const std::vector<double> nums = {numStubs, numTracks}; const std::vector<double> errs = {errStubs, errTracks}; const int wNums = std::ceil(std::log10(*std::max_element(nums.begin(), nums.end()))) + 5; const int wErrs = std::ceil(std::log10(*std::max_element(errs.begin(), errs.end()))) + 5; log_ << " HT SUMMARY " << std::endl; log_ << "number of stubs per TFP = " << std::setw(wNums) << numStubs << " +- " << std::setw(wErrs) << errStubs << std::endl; log_ << "number of tracks per TFP = " << std::setw(wNums) << numTracks << " +- " << std::setw(wErrs) << errTracks << std::endl; log_ << " max tracking efficiency = " << std::setw(wNums) << eff << " +- " << std::setw(wErrs) << errEff << std::endl; log_ << " fake rate = " << std::setw(wNums) << fracFake << std::endl; log_ << " duplicate rate = " << std::setw(wNums) << fracDup << std::endl; log_ << "============================================================="; edm::LogPrint(moduleDescription().moduleName()) << log_.str(); } // void AnalyzerHT::formTracks(const tt::StreamStub& stream, std::vector<std::vector<TTStubRef>>& tracks) const { static const DataFormat& layer = dataFormats_->format(Variable::layer, Process::ctb); auto toTrkId = [this](const StubHT& stub) { static const DataFormat& phiT = dataFormats_->format(Variable::phiT, Process::ht); static const DataFormat& zT = dataFormats_->format(Variable::zT, Process::ht); return (phiT.ttBV(stub.phiT()) + zT.ttBV(stub.zT())).val(); }; std::vector<StubHT> stubs; stubs.reserve(stream.size()); for (const tt::FrameStub& frame : stream) stubs.emplace_back(frame, dataFormats_); for (auto it = stubs.begin(); it != stubs.end();) { const auto start = it; const int id = toTrkId(*it); auto different = [id, toTrkId](const StubHT& stub) { return id != toTrkId(stub); }; it = std::find_if(it, stubs.end(), different); std::vector<TTStubRef> ttStubRefs; ttStubRefs.reserve(std::distance(start, it)); std::transform(start, it, std::back_inserter(ttStubRefs), [](const StubHT& stub) { return stub.frame().first; }); tracks.push_back(ttStubRefs); TTBV hitPattern(0, setup_->numLayers()); for (auto iter = start; iter != it; iter++) hitPattern.set(iter->layer().val(layer.width())); const double cot = dataFormats_->format(Variable::zT, Process::ht).floating(start->zT()) / setup_->chosenRofZ(); hisNumLayers_->Fill(hitPattern.count()); profNumLayers_->Fill(abs(sinh(cot)), hitPattern.count()); for (int layer : hitPattern.ids()) hisLayers_->Fill(layer); } } // void AnalyzerHT::associate(const std::vector<std::vector<TTStubRef>>& tracks, const tt::StubAssociation* ass, std::set<TPPtr>& tps, int& sum) const { for (const std::vector<TTStubRef>& ttStubRefs : tracks) { const std::vector<TPPtr>& tpPtrs = ass->associate(ttStubRefs); if (tpPtrs.empty()) continue; sum++; std::copy(tpPtrs.begin(), tpPtrs.end(), std::inserter(tps, tps.begin())); } } } // namespace trackerTFP DEFINE_FWK_MODULE(trackerTFP::AnalyzerHT);