CMS 3D CMS Logo

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