File indexing completed on 2024-04-06 12:18:22
0001
0002
0003
0004
0005 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007
0008 #include "HLTMessages.h"
0009 #include "HLTBitAnalyzer.h"
0010
0011 typedef std::pair<const char*, const edm::InputTag*> MissingCollectionInfo;
0012
0013 template <class T>
0014 static inline bool getCollection(const edm::Event& event,
0015 std::vector<MissingCollectionInfo>& missing,
0016 edm::Handle<T>& handle,
0017 const edm::InputTag& name,
0018 const edm::EDGetTokenT<T> token,
0019 const char* description) {
0020 event.getByToken(token, handle);
0021 bool valid = handle.isValid();
0022 if (not valid) {
0023 missing.push_back(std::make_pair(description, &name));
0024 handle.clear();
0025
0026 }
0027 return valid;
0028 }
0029
0030
0031
0032 HLTBitAnalyzer::HLTBitAnalyzer(edm::ParameterSet const& conf) : hlt_analysis_(conf, consumesCollector(), *this) {
0033
0034
0035
0036 std::cout << " Beginning HLTBitAnalyzer Analysis " << std::endl;
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054 mctruth_ = conf.getParameter<edm::InputTag>("mctruth");
0055 genEventInfo_ = conf.getParameter<edm::InputTag>("genEventInfo");
0056 VertexTagOffline0_ = conf.getParameter<edm::InputTag>("OfflinePrimaryVertices0");
0057 simhits_ = conf.getParameter<edm::InputTag>("simhits");
0058
0059 hltresults_ = conf.getParameter<edm::InputTag>("hltresults");
0060 l1results_ = conf.getParameter<edm::InputTag>("l1results");
0061
0062
0063
0064
0065
0066
0067
0068 pileupInfo_ = edm::InputTag("addPileupInfo");
0069
0070 hltresultsToken_ = consumes<edm::TriggerResults>(hltresults_);
0071 genEventInfoToken_ = consumes<GenEventInfoProduct>(genEventInfo_);
0072
0073
0074
0075
0076
0077
0078
0079
0080
0081
0082
0083
0084
0085 l1resultsToken_ = consumes<GlobalAlgBlkBxCollection>(l1results_);
0086
0087
0088
0089
0090
0091
0092
0093 mctruthToken_ = consumes<reco::CandidateView>(mctruth_);
0094 VertexTagOffline0Token_ = consumes<reco::VertexCollection>(VertexTagOffline0_);
0095 simtracksToken_ = consumes<std::vector<SimTrack> >(simhits_);
0096 simverticesToken_ = consumes<std::vector<SimVertex> >(simhits_);
0097 pileupInfoToken_ = consumes<std::vector<PileupSummaryInfo> >(pileupInfo_);
0098
0099 _UseTFileService = conf.getUntrackedParameter<bool>("UseTFileService", false);
0100
0101 m_file = nullptr;
0102 errCnt = 0;
0103
0104
0105 edm::ParameterSet runParameters = conf.getParameter<edm::ParameterSet>("RunParameters");
0106 _HistName = runParameters.getUntrackedParameter<std::string>("HistogramFile", "test.root");
0107 _isData = runParameters.getUntrackedParameter<bool>("isData", true);
0108
0109
0110 if (_UseTFileService) {
0111 edm::Service<TFileService> fs;
0112 HltTree = fs->make<TTree>("HltTree", "");
0113 } else {
0114 m_file = new TFile(_HistName.c_str(), "RECREATE");
0115 if (m_file)
0116 m_file->cd();
0117 HltTree = new TTree("HltTree", "");
0118 }
0119
0120
0121 hlt_analysis_.setup(conf, HltTree);
0122 if (!_isData) {
0123 mct_analysis_.setup(conf, HltTree);
0124 }
0125 vrt_analysisOffline0_.setup(conf, HltTree, "Offline0");
0126 evt_header_.setup(consumesCollector(), HltTree);
0127 }
0128
0129
0130 void HLTBitAnalyzer::analyze(edm::Event const& iEvent, edm::EventSetup const& iSetup) {
0131 edm::Handle<edm::TriggerResults> hltresults;
0132 edm::Handle<GlobalAlgBlkBxCollection> l1results;
0133
0134
0135
0136
0137
0138
0139
0140
0141
0142
0143
0144 edm::Handle<reco::CandidateView> mctruth;
0145 edm::Handle<GenEventInfoProduct> genEventInfo;
0146 edm::Handle<std::vector<SimTrack> > simTracks;
0147 edm::Handle<std::vector<SimVertex> > simVertices;
0148 edm::Handle<reco::VertexCollection> recoVertexsOffline0;
0149 edm::Handle<std::vector<PileupSummaryInfo> > pupInfo;
0150
0151
0152 std::vector<MissingCollectionInfo> missing;
0153
0154 getCollection(iEvent, missing, hltresults, hltresults_, hltresultsToken_, kHltresults);
0155
0156
0157
0158
0159
0160
0161
0162
0163
0164
0165
0166 getCollection(iEvent, missing, l1results, l1results_, l1resultsToken_, kL1GtRR);
0167
0168
0169
0170
0171
0172
0173 getCollection(iEvent, missing, mctruth, mctruth_, mctruthToken_, kMctruth);
0174 getCollection(iEvent, missing, simTracks, simhits_, simtracksToken_, kSimhit);
0175 getCollection(iEvent, missing, simVertices, simhits_, simverticesToken_, kSimhit);
0176 getCollection(iEvent, missing, genEventInfo, genEventInfo_, genEventInfoToken_, kGenEventInfo);
0177 getCollection(iEvent, missing, pupInfo, pileupInfo_, pileupInfoToken_, kPileupInfo);
0178
0179 getCollection(
0180 iEvent, missing, recoVertexsOffline0, VertexTagOffline0_, VertexTagOffline0Token_, kRecoVerticesOffline0);
0181
0182 if (!_isData) {
0183 ptHat = -1.;
0184 if (genEventInfo.isValid()) {
0185 ptHat = genEventInfo->qScale();
0186 }
0187
0188 weight = genEventInfo->weight();
0189 }
0190
0191 if (not missing.empty() and (errCnt < errMax())) {
0192 errCnt++;
0193 std::stringstream out;
0194 out << "OpenHLT analyser - missing collections:";
0195 for (auto const& entry : missing)
0196 out << "\n\t" << entry.first << ": " << entry.second->encode();
0197 edm::LogPrint("OpenHLT") << out.str() << std::endl;
0198 if (errCnt == errMax())
0199 edm::LogWarning("OpenHLT") << "Maximum error count reached -- No more messages will be printed.";
0200 }
0201
0202
0203 hlt_analysis_.analyze(hltresults,
0204
0205
0206
0207
0208
0209
0210
0211
0212
0213
0214
0215 l1results,
0216
0217
0218
0219
0220
0221 iSetup,
0222 iEvent,
0223 HltTree);
0224
0225 evt_header_.analyze(iEvent, HltTree);
0226
0227 if (!_isData) {
0228 mct_analysis_.analyze(mctruth, ptHat, weight, simTracks, simVertices, pupInfo, HltTree);
0229 }
0230 vrt_analysisOffline0_.analyze(recoVertexsOffline0, HltTree);
0231
0232
0233
0234 if (m_file) {
0235 m_file->cd();
0236 }
0237
0238 HltTree->Fill();
0239 }
0240
0241
0242 void HLTBitAnalyzer::beginRun(edm::Run const& run, edm::EventSetup const& es) { hlt_analysis_.beginRun(run, es); }
0243
0244
0245 void HLTBitAnalyzer::endJob() {
0246 if (!_UseTFileService) {
0247 if (m_file)
0248 m_file->cd();
0249
0250 HltTree->Write();
0251 delete HltTree;
0252 HltTree = nullptr;
0253
0254 if (m_file) {
0255 m_file->Write();
0256 delete m_file;
0257 m_file = nullptr;
0258 }
0259 }
0260 }
0261
0262
0263 #include "FWCore/Framework/interface/MakerMacros.h"
0264 DEFINE_FWK_MODULE(HLTBitAnalyzer);