CMS 3D CMS Logo

AlCaDiJetsProducer.cc
Go to the documentation of this file.
1 // system include files
2 #include <memory>
3 #include <string>
4 // user include files
13 
31 
34 #include <iostream>
35 
37  struct Counters {
38  Counters() : nAll_(0), nSelect_(0) {}
39  mutable std::atomic<unsigned int> nAll_, nSelect_;
40  };
41 } // namespace alCaHcalDiJetsProducer
42 
43 //
44 // class declaration
45 //
46 
47 class AlCaDiJetsProducer : public edm::global::EDProducer<edm::RunCache<alCaHcalDiJetsProducer::Counters>> {
48 public:
49  explicit AlCaDiJetsProducer(const edm::ParameterSet&);
50  ~AlCaDiJetsProducer() override = default;
51 
52  std::shared_ptr<alCaHcalDiJetsProducer::Counters> globalBeginRun(edm::Run const&,
53  edm::EventSetup const&) const override {
54  return std::make_shared<alCaHcalDiJetsProducer::Counters>();
55  }
56  void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
57  void globalEndRun(edm::Run const& iRun, edm::EventSetup const&) const override;
58  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
59 
60 private:
61  bool select(reco::PFJetCollection) const;
62 
63  // ----------member data ---------------------------
64 
66  double minPtJet_;
67 
72  //edm::EDGetTokenT<edm::TriggerResults> tok_TrigRes_;
75 };
76 
78  // Take input
79  labelPFJet_ = iConfig.getParameter<edm::InputTag>("PFjetInput");
80  labelHBHE_ = iConfig.getParameter<edm::InputTag>("HBHEInput");
81  labelHF_ = iConfig.getParameter<edm::InputTag>("HFInput");
82  labelHO_ = iConfig.getParameter<edm::InputTag>("HOInput");
83  //labelTrigger_ = iConfig.getParameter<edm::InputTag>("TriggerResults");
84  labelPFCandidate_ = iConfig.getParameter<edm::InputTag>("particleFlowInput");
85  labelVertex_ = iConfig.getParameter<edm::InputTag>("VertexInput");
86  minPtJet_ = iConfig.getParameter<double>("MinPtJet");
87 
88  tok_PFJet_ = consumes<reco::PFJetCollection>(labelPFJet_);
89  tok_HBHE_ = consumes<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>>(labelHBHE_);
90  tok_HF_ = consumes<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>>(labelHF_);
91  tok_HO_ = consumes<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>>(labelHO_);
92  //tok_TrigRes_= consumes<edm::TriggerResults>(labelTrigger_);
93  tok_PFCand_ = consumes<reco::PFCandidateCollection>(labelPFCandidate_);
94  tok_Vertex_ = consumes<reco::VertexCollection>(labelVertex_);
95 
96  // register your products
97  produces<reco::PFJetCollection>(labelPFJet_.encode());
98  produces<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>>(labelHBHE_.encode());
99  produces<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>>(labelHF_.encode());
100  produces<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>>(labelHO_.encode());
101  //produces<edm::TriggerResults>(labelTrigger_.encode());
102  produces<reco::PFCandidateCollection>(labelPFCandidate_.encode());
103  produces<reco::VertexCollection>(labelVertex_.encode());
104 }
105 
107  edm::LogVerbatim("AlcaDiJets") << "Accepts " << runCache(iRun.index())->nSelect_ << " events from a total of "
108  << runCache(iRun.index())->nAll_ << " events";
109 }
110 
112  if (jt.size() < 2)
113  return false;
114  if (((jt.at(0)).pt()) < minPtJet_)
115  return false;
116  return true;
117 }
118 
119 // ------------ method called to produce the data ------------
121  ++(runCache(iEvent.getRun().index())->nAll_);
122 
123  // Access the collections from iEvent
125  if (!pfjet.isValid()) {
126  edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelPFJet_;
127  return;
128  }
129  const reco::PFJetCollection pfjets = *(pfjet.product());
130 
132  if (!pfc.isValid()) {
133  edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelPFCandidate_;
134  return;
135  }
136  const reco::PFCandidateCollection pfcand = *(pfc.product());
137 
139  if (!vt.isValid()) {
140  edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelVertex_;
141  return;
142  }
143  const reco::VertexCollection vtx = *(vt.product());
144 
145  auto hbhe = iEvent.getHandle(tok_HBHE_);
146  if (!hbhe.isValid()) {
147  edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelHBHE_;
148  return;
149  }
151 
152  auto ho = iEvent.getHandle(tok_HO_);
153  if (!ho.isValid()) {
154  edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelHO_;
155  return;
156  }
158 
159  auto hf = iEvent.getHandle(tok_HF_);
160  if (!hf.isValid()) {
161  edm::LogWarning("AlCaDiJets") << "AlCaDiJetsProducer: Error! can't get product " << labelHF_;
162  return;
163  }
165 
166  // See if this event is useful
167  bool accept = select(pfjets);
168  if (accept) {
169  ++(runCache(iEvent.getRun().index())->nSelect_);
170 
171  //Copy from standard place
172  auto miniPFjetCollection = std::make_unique<reco::PFJetCollection>();
173  for (reco::PFJetCollection::const_iterator pfjetItr = pfjets.begin(); pfjetItr != pfjets.end(); pfjetItr++) {
174  miniPFjetCollection->push_back(*pfjetItr);
175  }
176 
177  auto miniPFCandCollection = std::make_unique<reco::PFCandidateCollection>();
178  for (reco::PFCandidateCollection::const_iterator pfcItr = pfcand.begin(); pfcItr != pfcand.end(); pfcItr++) {
179  miniPFCandCollection->push_back(*pfcItr);
180  }
181 
182  auto miniVtxCollection = std::make_unique<reco::VertexCollection>();
183  for (reco::VertexCollection::const_iterator vtxItr = vtx.begin(); vtxItr != vtx.end(); vtxItr++) {
184  miniVtxCollection->push_back(*vtxItr);
185  }
186 
187  auto miniHBHECollection =
188  std::make_unique<edm::SortedCollection<HBHERecHit, edm::StrictWeakOrdering<HBHERecHit>>>();
190  Hithbhe.begin();
191  hbheItr != Hithbhe.end();
192  hbheItr++) {
193  miniHBHECollection->push_back(*hbheItr);
194  }
195 
196  auto miniHOCollection = std::make_unique<edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>>();
197  for (edm::SortedCollection<HORecHit, edm::StrictWeakOrdering<HORecHit>>::const_iterator hoItr = Hitho.begin();
198  hoItr != Hitho.end();
199  hoItr++) {
200  miniHOCollection->push_back(*hoItr);
201  }
202 
203  auto miniHFCollection = std::make_unique<edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>>();
204  for (edm::SortedCollection<HFRecHit, edm::StrictWeakOrdering<HFRecHit>>::const_iterator hfItr = Hithf.begin();
205  hfItr != Hithf.end();
206  hfItr++) {
207  miniHFCollection->push_back(*hfItr);
208  }
209 
210  //Put them in the event
211  iEvent.put(std::move(miniPFjetCollection), labelPFJet_.encode());
212  iEvent.put(std::move(miniHBHECollection), labelHBHE_.encode());
213  iEvent.put(std::move(miniHFCollection), labelHF_.encode());
214  iEvent.put(std::move(miniHOCollection), labelHO_.encode());
215  //iEvent.put(std::move(miniTriggerCollection), labelTrigger_.encode());
216  iEvent.put(std::move(miniPFCandCollection), labelPFCandidate_.encode());
217  iEvent.put(std::move(miniVtxCollection), labelVertex_.encode());
218  }
219  return;
220 }
221 
224  desc.add<edm::InputTag>("PFjetInput", edm::InputTag("ak4PFJetsCHS")),
225  desc.add<edm::InputTag>("HBHEInput", edm::InputTag("hbhereco"));
226  desc.add<edm::InputTag>("HFInput", edm::InputTag("hfreco"));
227  desc.add<edm::InputTag>("HOInput", edm::InputTag("horeco"));
228  desc.add<edm::InputTag>("particleFlowInput", edm::InputTag("particleFlow"));
229  desc.add<edm::InputTag>("VertexInput", edm::InputTag("offlinePrimaryVertices"));
230  desc.add<double>("MinPtJet", 20.0);
231  descriptions.add("alcaDiJetsProducer", desc);
232 }
233 
Log< level::Info, true > LogVerbatim
edm::EDGetTokenT< edm::SortedCollection< HBHERecHit, edm::StrictWeakOrdering< HBHERecHit > > > tok_HBHE_
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
edm::InputTag labelPFCandidate_
std::string encode() const
Definition: InputTag.cc:159
~AlCaDiJetsProducer() override=default
edm::EDGetTokenT< edm::SortedCollection< HORecHit, edm::StrictWeakOrdering< HORecHit > > > tok_HO_
T const * product() const
Definition: Handle.h:70
edm::InputTag labelHBHE_
std::atomic< unsigned int > nAll_
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
edm::EDGetTokenT< reco::PFJetCollection > tok_PFJet_
std::atomic< unsigned int > nSelect_
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:31
void globalEndRun(edm::Run const &iRun, edm::EventSetup const &) const override
edm::EDGetTokenT< edm::SortedCollection< HFRecHit, edm::StrictWeakOrdering< HFRecHit > > > tok_HF_
bool select(reco::PFJetCollection) const
int iEvent
Definition: GenABIO.cc:224
edm::InputTag labelPFJet_
edm::EDGetTokenT< reco::VertexCollection > tok_Vertex_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
const_iterator begin() const
edm::EDGetTokenT< reco::PFCandidateCollection > tok_PFCand_
void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override
const_iterator end() const
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
void add(std::string const &label, ParameterSetDescription const &psetDescription)
bool isValid() const
Definition: HandleBase.h:70
std::vector< PFJet > PFJetCollection
collection of PFJet objects
RunIndex index() const
Definition: Run.cc:28
std::shared_ptr< alCaHcalDiJetsProducer::Counters > globalBeginRun(edm::Run const &, edm::EventSetup const &) const override
Log< level::Warning, false > LogWarning
AlCaDiJetsProducer(const edm::ParameterSet &)
edm::InputTag labelVertex_
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
def move(src, dest)
Definition: eostools.py:511
Definition: Run.h:45