Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:34:23

0001 // system include files
0002 #include <atomic>
0003 #include <memory>
0004 #include <cmath>
0005 #include <iostream>
0006 #include <sstream>
0007 #include <fstream>
0008 
0009 // user include files
0010 #include "FWCore/Framework/interface/Frameworkfwd.h"
0011 #include "FWCore/Framework/interface/stream/EDFilter.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/Framework/interface/Run.h"
0015 #include "FWCore/Framework/interface/LuminosityBlock.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0019 #include "FWCore/Common/interface/TriggerNames.h"
0020 
0021 #include "DataFormats/Common/interface/Handle.h"
0022 #include "DataFormats/Common/interface/Ref.h"
0023 #include "DataFormats/Common/interface/TriggerResults.h"
0024 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0025 #include "DataFormats/MuonReco/interface/Muon.h"
0026 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0027 
0028 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0029 #include "Calibration/IsolatedParticles/interface/CaloPropagateTrack.h"
0030 
0031 #include "MagneticField/Engine/interface/MagneticField.h"
0032 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0033 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0034 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0035 
0036 //#define EDM_ML_DEBUG
0037 //
0038 // class declaration
0039 //
0040 
0041 namespace alCaHEMuonFilter {
0042   struct Counters {
0043     Counters() : nAll_(0), nGood_(0), nFinal_(0) {}
0044     mutable std::atomic<unsigned int> nAll_, nGood_, nFinal_;
0045   };
0046 }  // namespace alCaHEMuonFilter
0047 
0048 class AlCaHEMuonFilter : public edm::stream::EDFilter<edm::GlobalCache<alCaHEMuonFilter::Counters> > {
0049 public:
0050   explicit AlCaHEMuonFilter(edm::ParameterSet const&, const alCaHEMuonFilter::Counters* count);
0051   ~AlCaHEMuonFilter() override;
0052 
0053   static std::unique_ptr<alCaHEMuonFilter::Counters> initializeGlobalCache(edm::ParameterSet const&) {
0054     return std::make_unique<alCaHEMuonFilter::Counters>();
0055   }
0056 
0057   bool filter(edm::Event&, edm::EventSetup const&) override;
0058   void endStream() override;
0059   static void globalEndJob(const alCaHEMuonFilter::Counters* counters);
0060   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0061 
0062 private:
0063   void beginRun(edm::Run const&, edm::EventSetup const&) override;
0064   void endRun(edm::Run const&, edm::EventSetup const&) override;
0065 
0066   // ----------member data ---------------------------
0067   HLTConfigProvider hltConfig_;
0068   unsigned int nRun_, nAll_, nGood_, nFinal_;
0069   const std::vector<std::string> trigNames_;
0070   const std::string processName_;
0071   const edm::InputTag triggerResults_, labelMuon_;
0072   const bool pfCut_;
0073   const double trackIsoCut_, caloIsoCut_, pfIsoCut_, muonptCut_, muonetaCut_;
0074   const int preScale_;
0075   const edm::EDGetTokenT<edm::TriggerResults> tok_trigRes_;
0076   const edm::EDGetTokenT<reco::MuonCollection> tok_Muon_;
0077   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> tok_geom_;
0078   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> tok_magField_;
0079 };
0080 
0081 //
0082 // constants, enums and typedefs
0083 //
0084 
0085 //
0086 // static data member definitions
0087 //
0088 
0089 //
0090 // constructors and destructor
0091 //
0092 AlCaHEMuonFilter::AlCaHEMuonFilter(edm::ParameterSet const& iConfig, const alCaHEMuonFilter::Counters* count)
0093     : nRun_(0),
0094       nAll_(0),
0095       nGood_(0),
0096       nFinal_(0),
0097       trigNames_(iConfig.getParameter<std::vector<std::string> >("triggers")),
0098       processName_(iConfig.getParameter<std::string>("processName")),
0099       triggerResults_(iConfig.getParameter<edm::InputTag>("triggerResultLabel")),
0100       labelMuon_(iConfig.getParameter<edm::InputTag>("muonLabel")),
0101       pfCut_(iConfig.getParameter<bool>("pfCut")),
0102       trackIsoCut_(iConfig.getParameter<double>("trackIsolationCut")),
0103       caloIsoCut_(iConfig.getParameter<double>("caloIsolationCut")),
0104       pfIsoCut_(iConfig.getParameter<double>("pfIsolationCut")),
0105       muonptCut_(iConfig.getParameter<double>("muonPtCut")),
0106       muonetaCut_(iConfig.getParameter<double>("muonEtaCut")),
0107       preScale_(iConfig.getParameter<int>("preScale") > 0 ? iConfig.getParameter<int>("preScale") : 1),
0108       tok_trigRes_(consumes<edm::TriggerResults>(triggerResults_)),
0109       tok_Muon_(consumes<reco::MuonCollection>(labelMuon_)),
0110       tok_geom_(esConsumes<CaloGeometry, CaloGeometryRecord>()),
0111       tok_magField_(esConsumes<MagneticField, IdealMagneticFieldRecord>()) {
0112   // load the parameters and define tokens for access
0113 
0114   edm::LogVerbatim("HEMuon") << "Parameters read from config file \n"
0115                              << "Process " << processName_ << "  Prescale " << preScale_ << "  Isolation Cuts "
0116                              << trackIsoCut_ << ":" << caloIsoCut_ << "\n";
0117   for (unsigned int k = 0; k < trigNames_.size(); ++k)
0118     edm::LogVerbatim("HEMuon") << "Trigger[" << k << "] " << trigNames_[k] << "\n";
0119 }  // AlCaHEMuonFilter::AlCaHEMuonFilter  constructor
0120 
0121 AlCaHEMuonFilter::~AlCaHEMuonFilter() {}
0122 
0123 //
0124 // member functions
0125 //
0126 
0127 // ------------ method called on each new Event  ------------
0128 bool AlCaHEMuonFilter::filter(edm::Event& iEvent, edm::EventSetup const& iSetup) {
0129   bool accept(false);
0130   ++nAll_;
0131 #ifdef EDM_ML_DEBUG
0132   edm::LogVerbatim("HEMuon") << "AlCaHEMuonFilter::Run " << iEvent.id().run() << " Event " << iEvent.id().event()
0133                              << " Luminosity " << iEvent.luminosityBlock() << " Bunch " << iEvent.bunchCrossing()
0134                              << std::endl;
0135 #endif
0136   //Step1: Find if the event passes one of the chosen triggers
0137   /////////////////////////////TriggerResults
0138   edm::Handle<edm::TriggerResults> triggerResults = iEvent.getHandle(tok_trigRes_);
0139   if (triggerResults.isValid()) {
0140     bool ok(false);
0141     std::vector<std::string> modules;
0142     const edm::TriggerNames& triggerNames = iEvent.triggerNames(*triggerResults);
0143     const std::vector<std::string>& triggerNames_ = triggerNames.triggerNames();
0144     for (unsigned int iHLT = 0; iHLT < triggerResults->size(); iHLT++) {
0145       int hlt = triggerResults->accept(iHLT);
0146       for (unsigned int i = 0; i < trigNames_.size(); ++i) {
0147         if (triggerNames_[iHLT].find(trigNames_[i]) != std::string::npos) {
0148           if (hlt > 0) {
0149             ok = true;
0150           }
0151 #ifdef EDM_ML_DEBUG
0152           edm::LogVerbatim("HEMuon") << "AlCaHEMuonFilter::Trigger " << triggerNames_[iHLT] << " Flag " << hlt << ":"
0153                                      << ok << std::endl;
0154 #endif
0155         }
0156       }
0157     }
0158     if (ok) {
0159       //Step2: Get geometry/B-field information
0160       //Get magnetic field
0161       const MagneticField* bField = &(iSetup.getData(tok_magField_));
0162       const CaloGeometry* geo = &(iSetup.getData(tok_geom_));
0163 
0164       // Relevant blocks from iEvent
0165       edm::Handle<reco::MuonCollection> muonHandle = iEvent.getHandle(tok_Muon_);
0166 #ifdef EDM_ML_DEBUG
0167       edm::LogVerbatim("HEMuon") << "AlCaHEMuonFilter::Muon Handle " << muonHandle.isValid();
0168 #endif
0169       if (muonHandle.isValid()) {
0170         for (reco::MuonCollection::const_iterator RecMuon = muonHandle->begin(); RecMuon != muonHandle->end();
0171              ++RecMuon) {
0172 #ifdef EDM_ML_DEBUG
0173           edm::LogVerbatim("HEMuon") << "AlCaHEMuonFilter::Muon:Track " << RecMuon->track().isNonnull()
0174                                      << " innerTrack " << RecMuon->innerTrack().isNonnull() << " outerTrack "
0175                                      << RecMuon->outerTrack().isNonnull() << " globalTrack "
0176                                      << RecMuon->globalTrack().isNonnull() << std::endl;
0177 #endif
0178           if ((RecMuon->track().isNonnull()) && (RecMuon->innerTrack().isNonnull()) &&
0179               (RecMuon->outerTrack().isNonnull()) && (RecMuon->globalTrack().isNonnull())) {
0180             bool ptCut = RecMuon->pt() < muonptCut_;
0181             bool etaCut = fabs(RecMuon->eta()) < muonetaCut_;
0182             if (ptCut || etaCut)
0183               continue;
0184 
0185             const reco::Track* pTrack = (RecMuon->innerTrack()).get();
0186             spr::propagatedTrackID trackID = spr::propagateCALO(pTrack, geo, bField, false);
0187 #ifdef EDM_ML_DEBUG
0188             edm::LogVerbatim("HEMuon") << "AlCaHEMuonFilter::Propagate: ECAL " << trackID.okECAL << " to HCAL "
0189                                        << trackID.okHCAL << std::endl;
0190 #endif
0191             double trackIso = RecMuon->isolationR03().sumPt;
0192             double caloIso = RecMuon->isolationR03().emEt + RecMuon->isolationR03().hadEt;
0193             double isolR04 =
0194                 ((RecMuon->pfIsolationR04().sumChargedHadronPt +
0195                   std::max(0.,
0196                            RecMuon->pfIsolationR04().sumNeutralHadronEt + RecMuon->pfIsolationR04().sumPhotonEt -
0197                                (0.5 * RecMuon->pfIsolationR04().sumPUPt))) /
0198                  RecMuon->pt());
0199             bool isoCut = (pfCut_) ? (isolR04 < pfIsoCut_) : ((trackIso < trackIsoCut_) && (caloIso < caloIsoCut_));
0200             if ((trackID.okECAL) && (trackID.okHCAL) && isoCut) {
0201               accept = true;
0202               break;
0203             }
0204           }
0205         }
0206       }
0207     }
0208   }
0209   // Step 4:  Return the acceptance flag
0210   if (accept) {
0211     ++nGood_;
0212     if (((nGood_ - 1) % preScale_) != 0) {
0213       accept = false;
0214     } else {
0215       ++nFinal_;
0216     }
0217   }
0218   return accept;
0219 
0220 }  // AlCaHEMuonFilter::filter
0221 
0222 // ------------ method called once each job just after ending the event loop  ------------
0223 void AlCaHEMuonFilter::endStream() {
0224   globalCache()->nAll_ += nAll_;
0225   globalCache()->nGood_ += nGood_;
0226   globalCache()->nFinal_ += nFinal_;
0227 }
0228 
0229 void AlCaHEMuonFilter::globalEndJob(const alCaHEMuonFilter::Counters* count) {
0230   edm::LogVerbatim("HEMuon") << "Selects " << count->nFinal_ << " out of " << count->nGood_ << " good events out of "
0231                              << count->nAll_ << " total # of events\n";
0232 }
0233 
0234 // ------------ method called when starting to processes a run  ------------
0235 void AlCaHEMuonFilter::beginRun(edm::Run const& iRun, edm::EventSetup const& iSetup) {
0236   bool changed(false);
0237   bool flag = hltConfig_.init(iRun, iSetup, processName_, changed);
0238   edm::LogVerbatim("HEMuon") << "Run[" << nRun_ << "] " << iRun.run() << " hltconfig.init " << flag << std::endl;
0239 }
0240 
0241 // ------------ method called when ending the processing of a run  ------------
0242 void AlCaHEMuonFilter::endRun(edm::Run const& iRun, edm::EventSetup const&) {
0243   edm::LogVerbatim("HEMuon") << "endRun[" << nRun_ << "] " << iRun.run() << "\n";
0244   nRun_++;
0245 }
0246 
0247 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0248 void AlCaHEMuonFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0249   //The following says we do not know what parameters are allowed so do no validation
0250   // Please change this to state exactly what you do use, even if it is no parameters
0251   edm::ParameterSetDescription desc;
0252   std::vector<std::string> triggers = {"HLT_IsoMu", "HLT_Mu"};
0253   desc.add<std::string>("processName", "HLT");
0254   desc.add<edm::InputTag>("triggerResultLabel", edm::InputTag("TriggerResults", "", "HLT"));
0255   desc.add<edm::InputTag>("muonLabel", edm::InputTag("muons"));
0256   desc.add<std::vector<std::string> >("triggers", triggers);
0257   desc.add<double>("muonPtCut", 5.0);
0258   desc.add<double>("muonEtaCut", 1.305);
0259   desc.add<bool>("pfCut", true);
0260   desc.add<double>("pfIsolationCut", 0.15);
0261   desc.add<double>("trackIsolationCut", 3.0);
0262   desc.add<double>("caloIsolationCut", 5.0);
0263   desc.add<int>("preScale", 1);
0264   descriptions.add("alcaHEMuonFilter", desc);
0265 }
0266 
0267 //define this as a plug-in
0268 #include "FWCore/Framework/interface/MakerMacros.h"
0269 DEFINE_FWK_MODULE(AlCaHEMuonFilter);