CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
CaloTowerAnalyzer.cc
Go to the documentation of this file.
2 // author: Bobby Scurlock, University of Florida
3 // first version 12/18/2006
4 // modified: Mike Schmitt
5 // date: 02.28.2007
6 // note: code rewrite
7 
11 
18 
26 
29 
34 
36 
39 
43 
44 //#include "Geometry/Vector/interface/GlobalPoint.h"
47 
48 #include <vector>
49 #include <utility>
50 #include <ostream>
51 #include <fstream>
52 #include <string>
53 #include <algorithm>
54 #include <cmath>
55 #include <memory>
56 #include <TLorentzVector.h>
57 
58 using namespace reco;
59 using namespace edm;
60 using namespace std;
61 
63  caloTowersLabel_ = consumes<edm::View<reco::Candidate> >(iConfig.getParameter<edm::InputTag>("CaloTowersLabel"));
64  HLTResultsLabel_ = consumes<edm::TriggerResults>(iConfig.getParameter<edm::InputTag>("HLTResultsLabel"));
65  HBHENoiseFilterResultLabel_ = consumes<bool>(iConfig.getParameter<edm::InputTag>("HBHENoiseFilterResultLabel"));
66 
67  if (iConfig.exists("HLTBitLabels"))
68  HLTBitLabel_ = iConfig.getParameter<std::vector<edm::InputTag> >("HLTBitLabels");
69 
70  debug_ = iConfig.getParameter<bool>("Debug");
71  finebinning_ = iConfig.getUntrackedParameter<bool>("FineBinning");
72  allhist_ = iConfig.getUntrackedParameter<bool>("AllHist");
73  FolderName_ = iConfig.getUntrackedParameter<std::string>("FolderName");
74 
75  hltselection_ = iConfig.getUntrackedParameter<bool>("HLTSelection");
76 }
77 
78 void CaloTowerAnalyzer::dqmbeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) { Nevents = 0; }
79 
81  ibooker.setCurrentFolder(FolderName_);
82 
83  //Store number of events which pass each HLT bit
84  for (unsigned int i = 0; i < HLTBitLabel_.size(); i++) {
85  if (!HLTBitLabel_[i].label().empty()) {
86  hCT_NEvents_HLT.push_back(
87  ibooker.book1D("METTask_CT_" + HLTBitLabel_[i].label(), HLTBitLabel_[i].label(), 2, -0.5, 1.5));
88  }
89  }
90 
91  //--Store number of events used
92  hCT_Nevents = ibooker.book1D("METTask_CT_Nevents", "", 1, 0, 1);
93  //--Data integrated over all events and stored by CaloTower(ieta,iphi)
94  hCT_et_ieta_iphi = ibooker.book2D("METTask_CT_et_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
95  hCT_et_ieta_iphi->setOption("colz");
96  hCT_et_ieta_iphi->setAxisTitle("ieta", 1);
97  hCT_et_ieta_iphi->setAxisTitle("ephi", 2);
98 
99  hCT_emEt_ieta_iphi = ibooker.book2D("METTask_CT_emEt_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
100  hCT_emEt_ieta_iphi->setOption("colz");
101  hCT_emEt_ieta_iphi->setAxisTitle("ieta", 1);
102  hCT_emEt_ieta_iphi->setAxisTitle("ephi", 2);
103  hCT_hadEt_ieta_iphi = ibooker.book2D("METTask_CT_hadEt_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
104  hCT_hadEt_ieta_iphi->setOption("colz");
105  hCT_hadEt_ieta_iphi->setAxisTitle("ieta", 1);
106  hCT_hadEt_ieta_iphi->setAxisTitle("ephi", 2);
107  hCT_outerEt_ieta_iphi = ibooker.book2D("METTask_CT_outerEt_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
108  hCT_outerEt_ieta_iphi->setOption("colz");
109  hCT_outerEt_ieta_iphi->setAxisTitle("ieta", 1);
110  hCT_outerEt_ieta_iphi->setAxisTitle("ephi", 2);
111  hCT_Occ_ieta_iphi = ibooker.book2D("METTask_CT_Occ_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
112  hCT_Occ_ieta_iphi->setOption("colz");
113  hCT_Occ_ieta_iphi->setAxisTitle("ieta", 1);
114  hCT_Occ_ieta_iphi->setAxisTitle("ephi", 2);
115 
116  hCT_Occ_EM_Et_ieta_iphi = ibooker.book2D("METTask_CT_Occ_EM_Et_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
117  hCT_Occ_EM_Et_ieta_iphi->setOption("colz");
118  hCT_Occ_EM_Et_ieta_iphi->setAxisTitle("ieta", 1);
119  hCT_Occ_EM_Et_ieta_iphi->setAxisTitle("ephi", 2);
120 
121  hCT_Occ_HAD_Et_ieta_iphi = ibooker.book2D("METTask_CT_Occ_HAD_Et_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
122  hCT_Occ_HAD_Et_ieta_iphi->setOption("colz");
123  hCT_Occ_HAD_Et_ieta_iphi->setAxisTitle("ieta", 1);
124  hCT_Occ_HAD_Et_ieta_iphi->setAxisTitle("ephi", 2);
125 
126  hCT_Occ_Outer_Et_ieta_iphi = ibooker.book2D("METTask_CT_Occ_Outer_Et_ieta_iphi", "", 83, -41, 42, 72, 1, 73);
127  hCT_Occ_Outer_Et_ieta_iphi->setOption("colz");
128  hCT_Occ_Outer_Et_ieta_iphi->setAxisTitle("ieta", 1);
129  hCT_Occ_Outer_Et_ieta_iphi->setAxisTitle("ephi", 2);
130 
131  //--Data over eta-rings
132 
133  // CaloTower values
134  if (allhist_) {
135  if (finebinning_) {
136  hCT_etvsieta = ibooker.book2D("METTask_CT_etvsieta", "", 83, -41, 42, 10001, 0, 1001);
137  hCT_Minetvsieta = ibooker.book2D("METTask_CT_Minetvsieta", "", 83, -41, 42, 10001, 0, 1001);
138  hCT_Maxetvsieta = ibooker.book2D("METTask_CT_Maxetvsieta", "", 83, -41, 42, 10001, 0, 1001);
139  hCT_emEtvsieta = ibooker.book2D("METTask_CT_emEtvsieta", "", 83, -41, 42, 10001, 0, 1001);
140  hCT_hadEtvsieta = ibooker.book2D("METTask_CT_hadEtvsieta", "", 83, -41, 42, 10001, 0, 1001);
141  hCT_outerEtvsieta = ibooker.book2D("METTask_CT_outerEtvsieta", "", 83, -41, 42, 10001, 0, 1001);
142  // Integrated over phi
143 
144  hCT_Occvsieta = ibooker.book2D("METTask_CT_Occvsieta", "", 83, -41, 42, 84, 0, 84);
145  hCT_SETvsieta = ibooker.book2D("METTask_CT_SETvsieta", "", 83, -41, 42, 20001, 0, 2001);
146  hCT_METvsieta = ibooker.book2D("METTask_CT_METvsieta", "", 83, -41, 42, 20001, 0, 2001);
147  hCT_METPhivsieta = ibooker.book2D("METTask_CT_METPhivsieta", "", 83, -41, 42, 80, -4, 4);
148  hCT_MExvsieta = ibooker.book2D("METTask_CT_MExvsieta", "", 83, -41, 42, 10001, -500, 501);
149  hCT_MEyvsieta = ibooker.book2D("METTask_CT_MEyvsieta", "", 83, -41, 42, 10001, -500, 501);
150  } else {
151  if (allhist_) {
152  hCT_etvsieta = ibooker.book2D("METTask_CT_etvsieta", "", 83, -41, 42, 200, -0.5, 999.5);
153  hCT_Minetvsieta = ibooker.book2D("METTask_CT_Minetvsieta", "", 83, -41, 42, 200, -0.5, 999.5);
154  hCT_Maxetvsieta = ibooker.book2D("METTask_CT_Maxetvsieta", "", 83, -41, 42, 200, -0.5, 999.5);
155  hCT_emEtvsieta = ibooker.book2D("METTask_CT_emEtvsieta", "", 83, -41, 42, 200, -0.5, 999.5);
156  hCT_hadEtvsieta = ibooker.book2D("METTask_CT_hadEtvsieta", "", 83, -41, 42, 200, -0.5, 999.5);
157  hCT_outerEtvsieta = ibooker.book2D("METTask_CT_outerEtvsieta", "", 83, -41, 42, 80, -0.5, 399.5);
158  // Integrated over phi
159  }
160 
161  hCT_Occvsieta = ibooker.book2D("METTask_CT_Occvsieta", "", 83, -41, 42, 73, -0.5, 72.5);
162  hCT_SETvsieta = ibooker.book2D("METTask_CT_SETvsieta", "", 83, -41, 42, 200, -0.5, 1999.5);
163  hCT_METvsieta = ibooker.book2D("METTask_CT_METvsieta", "", 83, -41, 42, 200, -0.5, 1999.5);
164  hCT_METPhivsieta = ibooker.book2D("METTask_CT_METPhivsieta", "", 83, -41, 42, 80, -4, 4);
165  hCT_MExvsieta = ibooker.book2D("METTask_CT_MExvsieta", "", 83, -41, 42, 100, -499.5, 499.5);
166  hCT_MEyvsieta = ibooker.book2D("METTask_CT_MEyvsieta", "", 83, -41, 42, 100, -499.5, 499.5);
167  }
168  } // allhist
169 }
170 
172  // Get HLT Results
173  edm::Handle<edm::TriggerResults> TheHLTResults;
174  iEvent.getByToken(HLTResultsLabel_, TheHLTResults);
175 
176  // **** Get the TriggerResults container
177  //triggerResultsToken_= consumes<edm::TriggerResults>(edm::InputTag(theTriggerResultsLabel));
178  //edm::Handle<edm::TriggerResults> triggerResults;
179  //iEvent.getByToken(triggerResultsToken_, triggerResults);
180 
181  bool EventPasses = true;
182  // Make sure handle is valid
183  if (TheHLTResults.isValid() && hltselection_) {
184  //Get HLT Names
185  const edm::TriggerNames& TheTriggerNames = iEvent.triggerNames(*TheHLTResults);
186 
187  for (unsigned int index = 0; index < HLTBitLabel_.size(); index++) {
188  if (!HLTBitLabel_[index].label().empty()) {
189  //Change the default value since HLT requirement has been issued by the user
190  if (index == 0)
191  EventPasses = false;
192  //Get the HLT bit and check to make sure it is valid
193  unsigned int bit = TheTriggerNames.triggerIndex(HLTBitLabel_[index].label());
194  if (bit < TheHLTResults->size()) {
195  //If any of the HLT names given by the user accept, then the event passes
196  if (TheHLTResults->accept(bit) && !TheHLTResults->error(bit)) {
197  EventPasses = true;
198  hCT_NEvents_HLT[index]->Fill(1);
199  } else
200  hCT_NEvents_HLT[index]->Fill(0);
201  } else {
202  edm::LogInfo("OutputInfo") << "The HLT Trigger Name : " << HLTBitLabel_[index].label()
203  << " is not valid for this trigger table " << std::endl;
204  }
205  }
206  }
207  }
208 
209  if (!EventPasses && hltselection_)
210  return;
211 
212  //----------GREG & CHRIS' idea---///
213  float ETTowerMin = -1; //GeV
214  float METRingMin = -2; // GeV
215 
216  Nevents++;
217  hCT_Nevents->Fill(0);
218 
219  // ==========================================================
220  // Retrieve!
221  // ==========================================================
222 
224  iEvent.getByToken(caloTowersLabel_, towers);
225 
226  if ((!towers.isValid())) {
227  //DD:fix print label
228  //edm::LogInfo("")<<"CaloTowers "<< caloTowersLabel_<<" not found!"<<std::endl;
229  return;
230  }
231 
232  //HBHENoiseFilterResultToken_=consumes<HBHENoiseFilter>(edm::InputTag(HBHENoiseFilterResultLabel_));
233  edm::Handle<bool> HBHENoiseFilterResultHandle;
234  iEvent.getByToken(HBHENoiseFilterResultLabel_, HBHENoiseFilterResultHandle);
235  bool HBHENoiseFilterResult = *HBHENoiseFilterResultHandle;
236  if (!HBHENoiseFilterResultHandle.isValid()) {
237  LogDebug("") << "CaloTowerAnalyzer: Could not find HBHENoiseFilterResult" << std::endl;
238  }
239 
240  bool bHcalNoiseFilter = HBHENoiseFilterResult;
241 
242  if (!bHcalNoiseFilter)
243  return;
244 
245  edm::View<Candidate>::const_iterator towerCand = towers->begin();
246 
247  // ==========================================================
248  // Fill Histograms!
249  // ==========================================================
250 
251  int CTmin_iphi = 99, CTmax_iphi = -99;
252  int CTmin_ieta = 99, CTmax_ieta = -99;
253 
254  TLorentzVector vMET_EtaRing[83];
255  int ActiveRing[83];
256  int NActiveTowers[83];
257  double SET_EtaRing[83];
258  double MinEt_EtaRing[83];
259  double MaxEt_EtaRing[83];
260  for (int i = 0; i < 83; i++) {
261  ActiveRing[i] = 0;
262  NActiveTowers[i] = 0;
263  SET_EtaRing[i] = 0;
264  MinEt_EtaRing[i] = 0;
265  MaxEt_EtaRing[i] = 0;
266  }
267 
268  // for (CaloTowerCollection::const_iterator calotower = towers->begin(); calotower != towers->end(); calotower++)
269  for (; towerCand != towers->end(); towerCand++) {
270  const Candidate* candidate = &(*towerCand);
271  if (candidate) {
272  const CaloTower* calotower = dynamic_cast<const CaloTower*>(candidate);
273  if (calotower) {
274  //math::RhoEtaPhiVector Momentum = calotower->momentum();
275  double Tower_ET = calotower->et();
276  //double Tower_Energy = calotower->energy();
277  // double Tower_Eta = calotower->eta();
278  double Tower_Phi = calotower->phi();
279  //double Tower_EMEnergy = calotower->emEnergy();
280  //double Tower_HadEnergy = calotower->hadEnergy();
281  double Tower_OuterEt = calotower->outerEt();
282  double Tower_EMEt = calotower->emEt();
283  double Tower_HadEt = calotower->hadEt();
284  //int Tower_EMLV1 = calotower->emLvl1();
285  //int Tower_HadLV1 = calotower->hadLv11();
286  int Tower_ieta = calotower->id().ieta();
287  int Tower_iphi = calotower->id().iphi();
288  int EtaRing = 41 + Tower_ieta;
289  ActiveRing[EtaRing] = 1;
290  NActiveTowers[EtaRing]++;
291  SET_EtaRing[EtaRing] += Tower_ET;
292  TLorentzVector v_;
293  v_.SetPtEtaPhiE(Tower_ET, 0, Tower_Phi, Tower_ET);
294  if (Tower_ET > ETTowerMin)
295  vMET_EtaRing[EtaRing] -= v_;
296 
297  // Fill Histograms
298  hCT_Occ_ieta_iphi->Fill(Tower_ieta, Tower_iphi);
299  if (calotower->emEt() > 0 && calotower->emEt() + calotower->hadEt() > 0.3)
300  hCT_Occ_EM_Et_ieta_iphi->Fill(Tower_ieta, Tower_iphi);
301  if (calotower->hadEt() > 0 && calotower->emEt() + calotower->hadEt() > 0.3)
302  hCT_Occ_HAD_Et_ieta_iphi->Fill(Tower_ieta, Tower_iphi);
303  if (calotower->outerEt() > 0 && calotower->emEt() + calotower->hadEt() > 0.3)
304  hCT_Occ_Outer_Et_ieta_iphi->Fill(Tower_ieta, Tower_iphi);
305 
306  hCT_et_ieta_iphi->Fill(Tower_ieta, Tower_iphi, Tower_ET);
307  hCT_emEt_ieta_iphi->Fill(Tower_ieta, Tower_iphi, Tower_EMEt);
308  hCT_hadEt_ieta_iphi->Fill(Tower_ieta, Tower_iphi, Tower_HadEt);
309  hCT_outerEt_ieta_iphi->Fill(Tower_ieta, Tower_iphi, Tower_OuterEt);
310 
311  if (allhist_) {
312  hCT_etvsieta->Fill(Tower_ieta, Tower_ET);
313  hCT_emEtvsieta->Fill(Tower_ieta, Tower_EMEt);
314  hCT_hadEtvsieta->Fill(Tower_ieta, Tower_HadEt);
315  hCT_outerEtvsieta->Fill(Tower_ieta, Tower_OuterEt);
316  }
317 
318  if (Tower_ET > MaxEt_EtaRing[EtaRing])
319  MaxEt_EtaRing[EtaRing] = Tower_ET;
320  if (Tower_ET < MinEt_EtaRing[EtaRing] && Tower_ET > 0)
321  MinEt_EtaRing[EtaRing] = Tower_ET;
322 
323  if (Tower_ieta < CTmin_ieta)
324  CTmin_ieta = Tower_ieta;
325  if (Tower_ieta > CTmax_ieta)
326  CTmax_ieta = Tower_ieta;
327  if (Tower_iphi < CTmin_iphi)
328  CTmin_iphi = Tower_iphi;
329  if (Tower_iphi > CTmax_iphi)
330  CTmax_iphi = Tower_iphi;
331  } //end if (calotower) ..
332  } // end if(candidate) ...
333 
334  } // end loop over towers
335 
336  // Fill eta-ring MET quantities
337  if (allhist_) {
338  for (int iEtaRing = 0; iEtaRing < 83; iEtaRing++) {
339  hCT_Minetvsieta->Fill(iEtaRing - 41, MinEt_EtaRing[iEtaRing]);
340  hCT_Maxetvsieta->Fill(iEtaRing - 41, MaxEt_EtaRing[iEtaRing]);
341 
342  if (ActiveRing[iEtaRing]) {
343  if (vMET_EtaRing[iEtaRing].Pt() > METRingMin) {
344  hCT_METPhivsieta->Fill(iEtaRing - 41, vMET_EtaRing[iEtaRing].Phi());
345  hCT_MExvsieta->Fill(iEtaRing - 41, vMET_EtaRing[iEtaRing].Px());
346  hCT_MEyvsieta->Fill(iEtaRing - 41, vMET_EtaRing[iEtaRing].Py());
347  hCT_METvsieta->Fill(iEtaRing - 41, vMET_EtaRing[iEtaRing].Pt());
348  }
349  hCT_SETvsieta->Fill(iEtaRing - 41, SET_EtaRing[iEtaRing]);
350  hCT_Occvsieta->Fill(iEtaRing - 41, NActiveTowers[iEtaRing]);
351  }
352  } // ietaring
353  } // allhist_
354 
355  edm::LogInfo("OutputInfo") << "CT ieta range: " << CTmin_ieta << " " << CTmax_ieta;
356  edm::LogInfo("OutputInfo") << "CT iphi range: " << CTmin_iphi << " " << CTmax_iphi;
357 }
T getUntrackedParameter(std::string const &, T const &) const
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:32
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
double hadEt() const
Definition: CaloTower.h:136
virtual void setOption(const char *option)
double outerEt() const
Definition: CaloTower.h:137
bool exists(std::string const &parameterName) const
checks if a parameter exists
unsigned int triggerIndex(std::string_view name) const
Definition: TriggerNames.cc:52
edm::TriggerNames const & triggerNames(edm::TriggerResults const &triggerResults) const override
Definition: Event.cc:275
char const * label
virtual void dqmbeginRun(const edm::Run &, const edm::EventSetup &)
CaloTowerAnalyzer(const edm::ParameterSet &)
int iEvent
Definition: GenABIO.cc:224
int iphi() const
get the tower iphi
bool isValid() const
Definition: HandleBase.h:70
Log< level::Info, false > LogInfo
CaloTowerDetId id() const
Definition: CaloTower.h:123
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:177
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
void analyze(const edm::Event &, const edm::EventSetup &) override
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:86
double et(double vtxZ) const
Definition: CaloTower.h:150
int ieta() const
get the tower ieta
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
double phi() const final
momentum azimuthal angle
double emEt() const
Definition: CaloTower.h:135
tuple size
Write out results.
Definition: Run.h:45
#define LogDebug(id)