Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-10-01 22:40:14

0001 // system include files
0002 #include <memory>
0003 #include <string>
0004 // user include files
0005 #include "FWCore/Framework/interface/Frameworkfwd.h"
0006 #include "FWCore/Framework/interface/global/EDProducer.h"
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "FWCore/Framework/interface/EventSetup.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 
0014 #include "DataFormats/DetId/interface/DetId.h"
0015 #include "DataFormats/EgammaCandidates/interface/Photon.h"
0016 #include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
0017 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0018 #include "DataFormats/HcalRecHit/interface/HBHERecHit.h"
0019 #include "DataFormats/HcalRecHit/interface/HFRecHit.h"
0020 #include "DataFormats/HcalRecHit/interface/HORecHit.h"
0021 #include "DataFormats/Common/interface/ValueMap.h"
0022 #include "DataFormats/Common/interface/TriggerResults.h"
0023 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0024 #include "DataFormats/HLTReco/interface/TriggerObject.h"
0025 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0026 #include "DataFormats/METReco/interface/PFMET.h"
0027 #include "DataFormats/METReco/interface/PFMETCollection.h"
0028 #include "DataFormats/VertexReco/interface/Vertex.h"
0029 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0030 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0031 
0032 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0033 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0034 #include <iostream>
0035 
0036 namespace alcaHcalDiJet {
0037   struct Counters {
0038     Counters() : nAll_(0), nSelect_(0) {}
0039     mutable std::atomic<unsigned int> nAll_, nSelect_;
0040   };
0041 }  // namespace alcaHcalDiJet
0042 
0043 //
0044 // class declaration
0045 //
0046 
0047 class AlCaDiJetsProducer : public edm::global::EDProducer<edm::RunCache<alcaHcalDiJet::Counters>> {
0048 public:
0049   explicit AlCaDiJetsProducer(const edm::ParameterSet&);
0050   ~AlCaDiJetsProducer() override = default;
0051 
0052   std::shared_ptr<alcaHcalDiJet::Counters> globalBeginRun(edm::Run const&, edm::EventSetup const&) const override {
0053     return std::make_shared<alcaHcalDiJet::Counters>();
0054   }
0055   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0056   void globalEndRun(edm::Run const& iRun, edm::EventSetup const&) const override;
0057   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0058 
0059 private:
0060   bool select(reco::PFJetCollection) const;
0061 
0062   // ----------member data ---------------------------
0063 
0064   edm::InputTag labelPFJet_, labelHBHE_, labelHF_, labelHO_, labelPFCandidate_, labelVertex_;  //labelTrigger_,
0065   double minPtJet_;
0066 
0067   edm::EDGetTokenT<reco::PFJetCollection> tok_PFJet_;
0068   edm::EDGetTokenT<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>> tok_HBHE_;
0069   edm::EDGetTokenT<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>> tok_HF_;
0070   edm::EDGetTokenT<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>> tok_HO_;
0071   //edm::EDGetTokenT<edm::TriggerResults>                                                   tok_TrigRes_;
0072   edm::EDGetTokenT<reco::PFCandidateCollection> tok_PFCand_;
0073   edm::EDGetTokenT<reco::VertexCollection> tok_Vertex_;
0074 };
0075 
0076 AlCaDiJetsProducer::AlCaDiJetsProducer(const edm::ParameterSet& iConfig) {
0077   // Take input
0078   labelPFJet_ = iConfig.getParameter<edm::InputTag>("PFjetInput");
0079   labelHBHE_ = iConfig.getParameter<edm::InputTag>("HBHEInput");
0080   labelHF_ = iConfig.getParameter<edm::InputTag>("HFInput");
0081   labelHO_ = iConfig.getParameter<edm::InputTag>("HOInput");
0082   //labelTrigger_    = iConfig.getParameter<edm::InputTag>("TriggerResults");
0083   labelPFCandidate_ = iConfig.getParameter<edm::InputTag>("particleFlowInput");
0084   labelVertex_ = iConfig.getParameter<edm::InputTag>("VertexInput");
0085   minPtJet_ = iConfig.getParameter<double>("MinPtJet");
0086 
0087   tok_PFJet_ = consumes<reco::PFJetCollection>(labelPFJet_);
0088   tok_HBHE_ = consumes<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>>(labelHBHE_);
0089   tok_HF_ = consumes<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>>(labelHF_);
0090   tok_HO_ = consumes<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>>(labelHO_);
0091   //tok_TrigRes_= consumes<edm::TriggerResults>(labelTrigger_);
0092   tok_PFCand_ = consumes<reco::PFCandidateCollection>(labelPFCandidate_);
0093   tok_Vertex_ = consumes<reco::VertexCollection>(labelVertex_);
0094 
0095   // register your products
0096   produces<reco::PFJetCollection>(labelPFJet_.encode());
0097   produces<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>>(labelHBHE_.encode());
0098   produces<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>>(labelHF_.encode());
0099   produces<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>>(labelHO_.encode());
0100   //produces<edm::TriggerResults>(labelTrigger_.encode());
0101   produces<reco::PFCandidateCollection>(labelPFCandidate_.encode());
0102   produces<reco::VertexCollection>(labelVertex_.encode());
0103 }
0104 
0105 void AlCaDiJetsProducer::globalEndRun(edm::Run const& iRun, edm::EventSetup const&) const {
0106   edm::LogVerbatim("AlcaDiJets") << "Accepts " << runCache(iRun.index())->nSelect_ << " events from a total of "
0107                                  << runCache(iRun.index())->nAll_ << " events";
0108 }
0109 
0110 bool AlCaDiJetsProducer::select(reco::PFJetCollection jt) const {
0111   if (jt.size() < 2)
0112     return false;
0113   if (((jt.at(0)).pt()) < minPtJet_)
0114     return false;
0115   return true;
0116 }
0117 
0118 // ------------ method called to produce the data  ------------
0119 void AlCaDiJetsProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0120   ++(runCache(iEvent.getRun().index())->nAll_);
0121 
0122   // Access the collections from iEvent
0123   edm::Handle<reco::PFJetCollection> pfjet = iEvent.getHandle(tok_PFJet_);
0124   if (!pfjet.isValid()) {
0125     edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelPFJet_;
0126     return;
0127   }
0128   const reco::PFJetCollection pfjets = *(pfjet.product());
0129 
0130   edm::Handle<reco::PFCandidateCollection> pfc = iEvent.getHandle(tok_PFCand_);
0131   if (!pfc.isValid()) {
0132     edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelPFCandidate_;
0133     return;
0134   }
0135   const reco::PFCandidateCollection pfcand = *(pfc.product());
0136 
0137   edm::Handle<reco::VertexCollection> vt = iEvent.getHandle(tok_Vertex_);
0138   if (!vt.isValid()) {
0139     edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelVertex_;
0140     return;
0141   }
0142   const reco::VertexCollection vtx = *(vt.product());
0143 
0144   auto hbhe = iEvent.getHandle(tok_HBHE_);
0145   if (!hbhe.isValid()) {
0146     edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelHBHE_;
0147     return;
0148   }
0149   const edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>> Hithbhe = *(hbhe.product());
0150 
0151   auto ho = iEvent.getHandle(tok_HO_);
0152   if (!ho.isValid()) {
0153     edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelHO_;
0154     return;
0155   }
0156   const edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>> Hitho = *(ho.product());
0157 
0158   auto hf = iEvent.getHandle(tok_HF_);
0159   if (!hf.isValid()) {
0160     edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelHF_;
0161     return;
0162   }
0163   const edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>> Hithf = *(hf.product());
0164 
0165   // See if this event is useful
0166   bool accept = select(pfjets);
0167   if (accept) {
0168     ++(runCache(iEvent.getRun().index())->nSelect_);
0169 
0170     //Copy from standard place
0171     auto miniPFjetCollection = std::make_unique<reco::PFJetCollection>();
0172     for (reco::PFJetCollection::const_iterator pfjetItr = pfjets.begin(); pfjetItr != pfjets.end(); pfjetItr++) {
0173       miniPFjetCollection->push_back(*pfjetItr);
0174     }
0175 
0176     auto miniPFCandCollection = std::make_unique<reco::PFCandidateCollection>();
0177     for (reco::PFCandidateCollection::const_iterator pfcItr = pfcand.begin(); pfcItr != pfcand.end(); pfcItr++) {
0178       miniPFCandCollection->push_back(*pfcItr);
0179     }
0180 
0181     auto miniVtxCollection = std::make_unique<reco::VertexCollection>();
0182     for (reco::VertexCollection::const_iterator vtxItr = vtx.begin(); vtxItr != vtx.end(); vtxItr++) {
0183       miniVtxCollection->push_back(*vtxItr);
0184     }
0185 
0186     auto miniHBHECollection =
0187         std::make_unique<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>>();
0188     for (edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>::const_iterator hbheItr =
0189              Hithbhe.begin();
0190          hbheItr != Hithbhe.end();
0191          hbheItr++) {
0192       miniHBHECollection->push_back(*hbheItr);
0193     }
0194 
0195     auto miniHOCollection = std::make_unique<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>>();
0196     for (edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>::const_iterator hoItr = Hitho.begin();
0197          hoItr != Hitho.end();
0198          hoItr++) {
0199       miniHOCollection->push_back(*hoItr);
0200     }
0201 
0202     auto miniHFCollection = std::make_unique<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>>();
0203     for (edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>::const_iterator hfItr = Hithf.begin();
0204          hfItr != Hithf.end();
0205          hfItr++) {
0206       miniHFCollection->push_back(*hfItr);
0207     }
0208 
0209     //Put them in the event
0210     iEvent.put(std::move(miniPFjetCollection), labelPFJet_.encode());
0211     iEvent.put(std::move(miniHBHECollection), labelHBHE_.encode());
0212     iEvent.put(std::move(miniHFCollection), labelHF_.encode());
0213     iEvent.put(std::move(miniHOCollection), labelHO_.encode());
0214     //iEvent.put(std::move(miniTriggerCollection),     labelTrigger_.encode());
0215     iEvent.put(std::move(miniPFCandCollection), labelPFCandidate_.encode());
0216     iEvent.put(std::move(miniVtxCollection), labelVertex_.encode());
0217   }
0218   return;
0219 }
0220 
0221 void AlCaDiJetsProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0222   edm::ParameterSetDescription desc;
0223   desc.add<edm::InputTag>("PFjetInput", edm::InputTag("ak4PFJetsCHS")),
0224       desc.add<edm::InputTag>("HBHEInput", edm::InputTag("hbhereco"));
0225   desc.add<edm::InputTag>("HFInput", edm::InputTag("hfreco"));
0226   desc.add<edm::InputTag>("HOInput", edm::InputTag("horeco"));
0227   desc.add<edm::InputTag>("particleFlowInput", edm::InputTag("particleFlow"));
0228   desc.add<edm::InputTag>("VertexInput", edm::InputTag("offlinePrimaryVertices"));
0229   desc.add<double>("MinPtJet", 20.0);
0230   descriptions.add("alcaDiJetsProducer", desc);
0231 }
0232 
0233 DEFINE_FWK_MODULE(AlCaDiJetsProducer);