CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Analyzer_minbias.cc
Go to the documentation of this file.
1 // system include files
2 #include <memory>
3 #include <string>
4 #include <iostream>
5 
6 // user include files
7 #include "Geometry/Records/interface/IdealGeometryRecord.h"
8 #include "Calibration/HcalCalibAlgos/interface/Analyzer_minbias.h"
9 #include "DataFormats/Provenance/interface/Provenance.h"
10 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
11 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
12 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
13 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
14 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
15 #include "DataFormats/L1GlobalTrigger/interface/L1GtfeWord.h"
16 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
17 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerRecord.h"
18 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMapRecord.h"
19 #include "CondFormats/L1TObjects/interface/L1GtTriggerMenu.h"
20 #include "CondFormats/DataRecord/interface/L1GtTriggerMenuRcd.h"
21 
22 #include "TFile.h"
23 #include "TH1.h"
24 #include "TH2.h"
25 #include <fstream>
26 #include <sstream>
27 
28 #include "CondFormats/HcalObjects/interface/HcalRespCorrs.h"
29 #include "CondFormats/DataRecord/interface/HcalRespCorrsRcd.h"
30 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMap.h"
31 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMapRecord.h"
32 using namespace std;
33 using namespace reco;
34 //
35 // constructors and destructor
36 //
37 namespace cms{
38 Analyzer_minbias::Analyzer_minbias(const edm::ParameterSet& iConfig)
39 {
40  // get name of output file with histogramms
41  fOutputFileName = iConfig.getUntrackedParameter<string>("HistOutFile");
42  // get names of modules, producing object collections
43 
44  hbherecoMB = iConfig.getParameter<edm::InputTag>("hbheInputMB");
45  horecoMB = iConfig.getParameter<edm::InputTag>("hoInputMB");
46  hfrecoMB = iConfig.getParameter<edm::InputTag>("hfInputMB");
47 
48  hbherecoNoise = iConfig.getParameter<edm::InputTag>("hbheInputNoise");
49  horecoNoise = iConfig.getParameter<edm::InputTag>("hoInputNoise");
50  hfrecoNoise = iConfig.getParameter<edm::InputTag>("hfInputNoise");
51 
52  theRecalib = iConfig.getParameter<bool>("Recalib");
53 
54 //
55 //
56  for(int i=0; i<73; i++)
57  {
58  for(int j=0; j<43; j++)
59  {
60  noise_min[i][j] = 0.;
61  noise_pl[i][j] = 0.;
62  }
63  }
64 //
65 //
66 
67 }
68 
69 Analyzer_minbias::~Analyzer_minbias()
70 {
71 
72  // do anything here that needs to be done at desctruction time
73  // (e.g. close files, deallocate resources etc.)
74 
75 }
76 
77 void Analyzer_minbias::beginRun( const edm::Run& r, const edm::EventSetup& iSetup)
78 {
79  nevent_run = 0;
80 }
81 void Analyzer_minbias::endRun( const edm::Run& r, const edm::EventSetup& iSetup)
82 {
83  std::cout<<" Runnumber "<<r.run()<<" Nevents "<<nevent_run<<std::endl;
84 }
85 
86 void Analyzer_minbias::beginJob()
87 {
88 
89  hOutputFile = new TFile( fOutputFileName.c_str(), "RECREATE" ) ;
90 
91  myTree = new TTree("RecJet","RecJet Tree");
92  myTree->Branch("mydet", &mydet, "mydet/I");
93  myTree->Branch("mysubd", &mysubd, "mysubd/I");
94  myTree->Branch("depth", &depth, "depth/I");
95  myTree->Branch("ieta", &ieta, "ieta/I");
96  myTree->Branch("iphi", &iphi, "iphi/I");
97  myTree->Branch("eta", &eta, "eta/F");
98  myTree->Branch("phi", &phi, "phi/F");
99 
100  myTree->Branch("mom0_MB", &mom0_MB, "mom0_MB/F");
101  myTree->Branch("mom1_MB", &mom1_MB, "mom1_MB/F");
102  myTree->Branch("mom2_MB", &mom2_MB, "mom2_MB/F");
103  myTree->Branch("mom4_MB", &mom4_MB, "mom4_MB/F");
104 
105  myTree->Branch("mom0_Noise", &mom0_Noise, "mom0_Noise/F");
106  myTree->Branch("mom1_Noise", &mom1_Noise, "mom1_Noise/F");
107  myTree->Branch("mom2_Noise", &mom2_Noise, "mom2_Noise/F");
108  myTree->Branch("mom4_Noise", &mom4_Noise, "mom4_Noise/F");
109 
110  myTree->Branch("mom0_Diff", &mom0_Diff, "mom0_Diff/F");
111  myTree->Branch("mom1_Diff", &mom1_Diff, "mom1_Diff/F");
112  myTree->Branch("mom2_Diff", &mom2_Diff, "mom2_Diff/F");
113 
114  myTree->Branch("occup", &occup, "occup/F");
115 
116  std::cout<<" Before ordering Histos "<<std::endl;
117 
118  char str0[15];
119  char str1[15];
120 
121  char str10[15];
122  char str11[15];
123 
124  int k=0;
125  nevent = 0;
126 // Size of collections
127 
128  hHBHEsize_vs_run = new TH2F("hHBHEsize_vs_run","hHBHEsize_vs_run",500,111500.,112000.,6101,-100.5,6000.5);
129  hHFsize_vs_run = new TH2F("hHFsize_vs_run","hHFsize_vs_run",500,111500.,112000.,6101,-100.5,6000.5);
130 
131  for(int i=1;i<73;i++){
132  for(int j=1;j<43;j++){
133 
134  meannoise_pl[i][j] = 0.;
135  meannoise_min[i][j] = 0.;
136 
137 // for(int l=1;l<5;l++){
138  k = i*1000+j;
139  sprintf(str0,"mpl%d",k);
140  sprintf(str1,"mmin%d",k);
141 
142  sprintf(str10,"vpl%d",k);
143  sprintf(str11,"vmin%d",k);
144 // cout<<" "<<i<<" "<<j<<endl;
145  if( j < 30 )
146  {
147 // first order moment
148  hCalo1[i][j] = new TH1F(str0, "h0", 320, -10., 10.);
149  hCalo2[i][j] = new TH1F(str1, "h1", 320, -10., 10.);
150 
151 // second order moment
152  hCalo1mom2[i][j] = new TH1F(str10, "h10", 320, 0., 20.);
153  hCalo2mom2[i][j] = new TH1F(str11, "h11", 320, 0., 20.);
154  }
155  else
156  {
157 // HF
158 // first order moment
159 // cout<<" "<<i<<" "<<j<<" "<<k<<endl;
160  if(j < 40)
161  {
162  hCalo1[i][j] = new TH1F(str0, "h0", 320, -10., 10.);
163  hCalo2[i][j] = new TH1F(str1, "h1", 320, -10., 10.);
164 //
165 // second order moment
166  hCalo1mom2[i][j] = new TH1F(str10, "h10", 320, 0., 40.);
167  hCalo2mom2[i][j] = new TH1F(str11, "h11", 320, 0., 40.);
168  }
169  else
170  {
171  hCalo1[i][j] = new TH1F(str0,"h0" , 320, -10., 10.);
172  hCalo2[i][j] = new TH1F(str1, "h1", 320, -10., 10.);
173 
174 // second order moment
175  hCalo1mom2[i][j] = new TH1F(str10, "h10", 320, 0., 120.);
176  hCalo2mom2[i][j] = new TH1F(str11, "h11", 320, 0., 120.);
177 
178  }
179  } // HE/HF boundary
180 // } // l
181  } // j
182  } // i
183 
184 
185  hbheNoiseE = new TH1F("hbheNoiseE","hbheNoiseE", 320, -10., 10.);
186  hfNoiseE = new TH1F("hfNoiseE","hfNoiseE", 320, -10., 10.);
187  hbheSignalE = new TH1F("hbheSignalE","hbheSignalE", 320, -10., 10.);
188  hfSignalE = new TH1F("hfSignalE","hfSignalE", 320, -10., 10.);
189 
190 
191  std::cout<<" After ordering Histos "<<std::endl;
192 
193  std::string ccc = "noise_0.dat";
194 
195  myout_hcal = new ofstream(ccc.c_str());
196  if(!myout_hcal) cout << " Output file not open!!! "<<endl;
197 
198 //
199  for (int i=0; i<5;i++)
200  {
201  for (int j=0; j<5;j++)
202  {
203  for (int k=0; k<73;k++)
204  {
205  for (int l=0; l<43;l++)
206  {
207  theMBFillDetMapPl0[i][j][k][l] = 0.;
208  theMBFillDetMapPl1[i][j][k][l] = 0.;
209  theMBFillDetMapPl2[i][j][k][l] = 0.;
210  theMBFillDetMapPl4[i][j][k][l] = 0.;
211 
212  theMBFillDetMapMin0[i][j][k][l] = 0.;
213  theMBFillDetMapMin1[i][j][k][l] = 0.;
214  theMBFillDetMapMin2[i][j][k][l] = 0.;
215  theMBFillDetMapMin4[i][j][k][l] = 0.;
216 
217 
218  theNSFillDetMapPl0[i][j][k][l] = 0.;
219  theNSFillDetMapPl1[i][j][k][l] = 0.;
220  theNSFillDetMapPl2[i][j][k][l] = 0.;
221  theNSFillDetMapPl4[i][j][k][l] = 0.;
222 
223  theNSFillDetMapMin0[i][j][k][l] = 0.;
224  theNSFillDetMapMin1[i][j][k][l] = 0.;
225  theNSFillDetMapMin2[i][j][k][l] = 0.;
226  theNSFillDetMapMin4[i][j][k][l] = 0.;
227 
228  theDFFillDetMapPl0[i][j][k][l] = 0.;
229  theDFFillDetMapPl1[i][j][k][l] = 0.;
230  theDFFillDetMapPl2[i][j][k][l] = 0.;
231  theDFFillDetMapMin0[i][j][k][l] = 0.;
232  theDFFillDetMapMin1[i][j][k][l] = 0.;
233  theDFFillDetMapMin2[i][j][k][l] = 0.;
234  }
235  }
236  }
237  }
238 
239  return ;
240 }
241 //
242 // EndJob
243 //
244 void Analyzer_minbias::endJob()
245 {
246  int ii=0;
247 
248  for (int i=1; i<5;i++)
249  {
250  for (int j=1; j<5;j++)
251  {
252  for (int k=1; k<73;k++)
253  {
254  for (int l=1; l<43;l++)
255  {
256  if(theMBFillDetMapPl0[i][j][k][l] > 0)
257  {
258  mom0_MB = theMBFillDetMapPl0[i][j][k][l];
259  mom1_MB = theMBFillDetMapPl1[i][j][k][l];
260  mom2_MB = theMBFillDetMapPl2[i][j][k][l];
261  mom4_MB = theMBFillDetMapPl4[i][j][k][l];
262  mom0_Noise = theNSFillDetMapPl0[i][j][k][l];
263  mom1_Noise = theNSFillDetMapPl1[i][j][k][l];
264  mom2_Noise = theNSFillDetMapPl2[i][j][k][l];
265  mom4_Noise = theNSFillDetMapPl4[i][j][k][l];
266  mom0_Diff = theDFFillDetMapPl0[i][j][k][l];
267  mom1_Diff = theDFFillDetMapPl1[i][j][k][l];
268  mom2_Diff = theDFFillDetMapPl2[i][j][k][l];
269 
270  mysubd = i;
271  depth = j;
272  ieta = l;
273  iphi = k;
274  cout<<" Result Plus= "<<mysubd<<" "<<ieta<<" "<<iphi<<" mom0 "<<mom0_MB<<" mom1 "<<mom1_MB<<" mom2 "<<mom2_MB<<endl;
275  myTree->Fill();
276  ii++;
277  } // Pl > 0
278 
279 
280  if(theMBFillDetMapMin0[i][j][k][l] > 0)
281  {
282  mom0_MB = theMBFillDetMapMin0[i][j][k][l];
283  mom1_MB = theMBFillDetMapMin1[i][j][k][l];
284  mom2_MB = theMBFillDetMapMin2[i][j][k][l];
285  mom4_MB = theMBFillDetMapMin4[i][j][k][l];
286  mom0_Noise = theNSFillDetMapMin0[i][j][k][l];
287  mom1_Noise = theNSFillDetMapMin1[i][j][k][l];
288  mom2_Noise = theNSFillDetMapMin2[i][j][k][l];
289  mom4_Noise = theNSFillDetMapMin4[i][j][k][l];
290  mom0_Diff = theDFFillDetMapMin0[i][j][k][l];
291  mom1_Diff = theDFFillDetMapMin1[i][j][k][l];
292  mom2_Diff = theDFFillDetMapMin2[i][j][k][l];
293 
294 
295  mysubd = i;
296  depth = j;
297  ieta = -1*l;
298  iphi = k;
299  cout<<" Result Minus= "<<mysubd<<" "<<ieta<<" "<<iphi<<" mom0 "<<mom0_MB<<" mom1 "<<mom1_MB<<" mom2 "<<mom2_MB<<endl;
300  myTree->Fill();
301  ii++;
302 
303  } // Min>0
304  } // ieta
305  } // iphi
306  } // depth
307  } //subd
308 
309 
310 
311  cout<<" Number of cells "<<ii<<endl;
312 
313  hOutputFile->Write();
314 
315  hOutputFile->cd();
316 
317  myTree->Write();
318 
319  hHBHEsize_vs_run->Write() ;
320  hHFsize_vs_run->Write() ;
321 
322  for(int i=1;i<73;i++){
323  for(int j=1;j<43;j++){
324  hCalo1[i][j]->Write();
325  hCalo2[i][j]->Write();
326  hCalo1mom2[i][j]->Write();
327  hCalo2mom2[i][j]->Write();
328  }
329  }
330 
331  hbheNoiseE->Write() ;
332  hfNoiseE->Write() ;
333  hbheSignalE->Write() ;
334  hfSignalE->Write() ;
335 
336 
337  hOutputFile->Close() ;
338 
339  cout<<" File is closed "<<endl;
340 
341  return ;
342 }
343 
344 
345 //
346 // member functions
347 //
348 
349 // ------------ method called to produce the data ------------
350 void
351 Analyzer_minbias::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
352 {
353 
354  std::cout<<" Start Analyzer_minbias::analyze "<<nevent<<std::endl;
355  nevent++;
356  nevent_run++;
357  using namespace edm;
358 
359  float rnnum = (float)iEvent.run();
360 /*
361  std::vector<Provenance const*> theProvenance;
362  iEvent.getAllProvenance(theProvenance);
363 
364  for( std::vector<Provenance const*>::const_iterator ip = theProvenance.begin();
365  ip != theProvenance.end(); ip++)
366  {
367  cout<<" Print all module/label names "<<(**ip).moduleName()<<" "<<(**ip).moduleLabel()<<
368  " "<<(**ip).productInstanceName()<<endl;
369  }
370 */
371 
372  edm::ESHandle<L1GtTriggerMenu> menuRcd;
373  iSetup.get<L1GtTriggerMenuRcd>().get(menuRcd) ;
374  const L1GtTriggerMenu* menu = menuRcd.product();
375  const AlgorithmMap& bitMap = menu->gtAlgorithmMap();
376 
377  edm::Handle<L1GlobalTriggerReadoutRecord> gtRecord;
378  iEvent.getByLabel("gtDigisAlCaMB", gtRecord);
379 
380  if (!gtRecord.isValid()) {
381 
382 // LogDebug("L1GlobalTriggerRecordProducer")
383 // << "\n\n Error: no L1GlobalTriggerReadoutRecord found with input tag "
384 // << m_l1GtReadoutRecord
385 // << "\n Returning empty L1GlobalTriggerRecord.\n\n"
386 // << std::endl;
387  cout<<" No L1 trigger record "<<endl;
388  } else {
389 
390  const DecisionWord dWord = gtRecord->decisionWord();
391 
392  for (CItAlgo itAlgo = bitMap.begin(); itAlgo != bitMap.end(); itAlgo++)
393  {
394  bool decision=menu->gtAlgorithmResult(itAlgo->first,dWord);
395  if(decision == 1) std::cout<<" Trigger "<<itAlgo->first<<" "<<decision<<std::endl;
396  }
397 
398  }
399 
400  const HcalRespCorrs* myRecalib=0;
401  if( theRecalib ) {
402 // Radek:
403  edm::ESHandle <HcalRespCorrs> recalibCorrs;
404  iSetup.get<HcalRespCorrsRcd>().get("recalibrate",recalibCorrs);
405  myRecalib = recalibCorrs.product();
406 // end
407  } // theRecalib
408 
409 // Noise part for HB HE
410 
411  double tmpNSFillDetMapPl1[5][5][73][43];
412  double tmpNSFillDetMapMin1[5][5][73][43];
413 
414  for (int i=0; i<5;i++)
415  {
416  for (int j=0; j<5;j++)
417  {
418  for (int k=0; k<73;k++)
419  {
420  for (int l=0; l<43;l++)
421  {
422  tmpNSFillDetMapPl1[i][j][k][l] = 0.;
423  tmpNSFillDetMapMin1[i][j][k][l] = 0.;
424  }
425  }
426  }
427  }
428  edm::Handle<HBHERecHitCollection> hbheNormal;
429  iEvent.getByLabel("hbhereco", hbheNormal);
430  if(!hbheNormal.isValid()){
431  cout<<" hbheNormal failed "<<endl;
432  } else {
433  cout<<" The size of the normal collection "<<hbheNormal->size()<<endl;
434  }
435 
436 
437  edm::Handle<HBHERecHitCollection> hbheNS;
438  iEvent.getByLabel(hbherecoNoise, hbheNS);
439 
440 
441  if(!hbheNS.isValid()){
442  LogDebug("") << "HcalCalibAlgos: Error! can't get hbhe product!" << std::endl;
443  cout<<" No HBHE MS "<<endl;
444  return ;
445  }
446 
447 
448  const HBHERecHitCollection HithbheNS = *(hbheNS.product());
449  cout<<" HBHE NS size of collection "<<HithbheNS.size()<<endl;
450  hHBHEsize_vs_run->Fill(rnnum,(float)HithbheNS.size());
451 
452  if(HithbheNS.size()!= 5184) {
453  cout<<" HBHE problem "<<rnnum<<" "<<HithbheNS.size()<<endl;
454  return;
455  }
456  edm::Handle<HBHERecHitCollection> hbheMB;
457  iEvent.getByLabel(hbherecoMB, hbheMB);
458 
459  if(!hbheMB.isValid()){
460  LogDebug("") << "HcalCalibAlgos: Error! can't get hbhe product!" << std::endl;
461  cout<<" No HBHE MB"<<endl;
462  return ;
463  }
464 
465  const HBHERecHitCollection HithbheMB = *(hbheMB.product());
466  cout<<" HBHE MB size of collection "<<HithbheMB.size()<<endl;
467  if(HithbheMB.size()!= 5184) {
468  cout<<" HBHE problem "<<rnnum<<" "<<HithbheMB.size()<<endl;
469  return;
470  }
471 
472  edm::Handle<HFRecHitCollection> hfNS;
473  iEvent.getByLabel(hfrecoNoise, hfNS);
474 
475  if(!hfNS.isValid()){
476  LogDebug("") << "HcalCalibAlgos: Error! can't get hbhe product!" << std::endl;
477  cout<<" No HF NS "<<endl;
478  return ;
479  }
480 
481  const HFRecHitCollection HithfNS = *(hfNS.product());
482  cout<<" HFE NS size of collection "<<HithfNS.size()<<endl;
483  hHFsize_vs_run->Fill(rnnum,(float)HithfNS.size());
484  if(HithfNS.size()!= 1728) {
485  cout<<" HF problem "<<rnnum<<" "<<HithfNS.size()<<endl;
486  return;
487  }
488 
489  edm::Handle<HFRecHitCollection> hfMB;
490  iEvent.getByLabel(hfrecoMB, hfMB);
491 
492  if(!hfMB.isValid()){
493  LogDebug("") << "HcalCalibAlgos: Error! can't get hbhe product!" << std::endl;
494  cout<<" No HBHE MB"<<endl;
495  return ;
496  }
497 
498  const HFRecHitCollection HithfMB = *(hfMB.product());
499  cout<<" HF MB size of collection "<<HithfMB.size()<<endl;
500  if(HithfMB.size()!= 1728) {
501  cout<<" HF problem "<<rnnum<<" "<<HithfMB.size()<<endl;
502  return;
503  }
504 
505 
506 
507  for(HBHERecHitCollection::const_iterator hbheItr=HithbheNS.begin(); hbheItr!=HithbheNS.end(); hbheItr++)
508  {
509 // Recalibration of energy
510  float icalconst=1.;
511  DetId mydetid = hbheItr->id().rawId();
512  if( theRecalib ) icalconst=myRecalib->getValues(mydetid)->getValue();
513 
514  HBHERecHit aHit(hbheItr->id(),hbheItr->energy()*icalconst,hbheItr->time());
515 
516  double energyhit = aHit.energy();
517 
518  DetId id = (*hbheItr).detid();
519  HcalDetId hid=HcalDetId(id);
520 
521 
522 
523  int mysu = ((hid).rawId()>>25)&0x7;
524  if( hid.ieta() > 0 ) {
525  theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()]+ 1.;
526  theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+energyhit;
527  theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,2);
528  theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,4);
529 
530  tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
531 
532 
533  } else {
534  theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+ 1.;
535  theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+energyhit;
536  theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,2);
537  theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,4);
538 
539 
540  tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
541 
542  }
543 
544  if(hid.depth() == 1) {
545  hbheNoiseE->Fill(energyhit);
546 
547  if(energyhit<-2.) std::cout<<" Run "<<rnnum<<" ieta,iphi "<<hid.ieta()<<" "<<hid.iphi()<<energyhit<<std::endl;
548 
549  // if( hid.ieta() > 0 ) {
550  // hCalo1[hid.iphi()][hid.ieta()]->Fill(energyhit-noise_pl[hid.iphi()][hid.ieta()]);
551  // hCalo1mom2[hid.iphi()][hid.ieta()]->Fill(pow(energyhit,2));
552  // } else {
553  // hCalo2[hid.iphi()][abs(hid.ieta())]->Fill(energyhit-noise_min[hid.iphi()][abs(hid.ieta())]);
554  // hCalo2mom2[hid.iphi()][abs(hid.ieta())]->Fill(pow(energyhit,2));
555  // } // eta><0
556 
557  } // depth=1
558 
559 
560  } // HBHE_NS
561 
562 
563 // Signal part for HB HE
564 
565 
566  for(HBHERecHitCollection::const_iterator hbheItr=HithbheMB.begin(); hbheItr!=HithbheMB.end(); hbheItr++)
567  {
568 // Recalibration of energy
569  float icalconst=1.;
570  DetId mydetid = hbheItr->id().rawId();
571  if( theRecalib ) icalconst=myRecalib->getValues(mydetid)->getValue();
572 
573  HBHERecHit aHit(hbheItr->id(),hbheItr->energy()*icalconst,hbheItr->time());
574 
575  double energyhit = aHit.energy();
576 
577  DetId id = (*hbheItr).detid();
578  HcalDetId hid=HcalDetId(id);
579 
580  int mysu = ((hid).rawId()>>25)&0x7;
581  if( hid.ieta() > 0 ) {
582  theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()]+ 1.;
583  theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+energyhit;
584  theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,2);
585  theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,4);
586  float mydiff = energyhit - tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
587 
588 
589  theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+mydiff;
590  theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(mydiff,2);
591  } else {
592  theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+ 1.;
593  theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+energyhit;
594  theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,2);
595  theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,4);
596 
597 
598  float mydiff = energyhit - tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
599  theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+mydiff;
600  theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(mydiff,2);
601  }
602 
603 
604  if(hid.depth() == 1) {
605 
606  hbheSignalE->Fill(energyhit);
607 
608  if( hid.ieta() > 0 ) {
609  hCalo1[hid.iphi()][hid.ieta()]->Fill(energyhit);
610  hCalo1mom2[hid.iphi()][hid.ieta()]->Fill(pow(energyhit,2));
611  } else {
612  hCalo2[hid.iphi()][abs(hid.ieta())]->Fill(energyhit);
613  hCalo2mom2[hid.iphi()][abs(hid.ieta())]->Fill(pow(energyhit,2));
614  } // eta><0
615 
616  } // depth=1
617 
618 
619  } // HBHE_MB
620 
621 // HF
622 
623  for(HFRecHitCollection::const_iterator hbheItr=HithfNS.begin(); hbheItr!=HithfNS.end(); hbheItr++)
624  {
625 // Recalibration of energy
626  float icalconst=1.;
627  DetId mydetid = hbheItr->id().rawId();
628  if( theRecalib ) icalconst=myRecalib->getValues(mydetid)->getValue();
629 
630  HFRecHit aHit(hbheItr->id(),hbheItr->energy()*icalconst,hbheItr->time());
631 
632  double energyhit = aHit.energy();
633 //
634 // Remove PMT hits
635 //
636  if(fabs(energyhit) > 40. ) continue;
637 
638  DetId id = (*hbheItr).detid();
639  HcalDetId hid=HcalDetId(id);
640 
641  int mysu = ((hid).rawId()>>25)&0x7;
642  if( hid.ieta() > 0 ) {
643  theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()]+ 1.;
644  theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+energyhit;
645  theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,2);
646  theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,4);
647 
648  tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
649 
650 
651  } else {
652  theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+ 1.;
653  theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+energyhit;
654  theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,2);
655  theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,4);
656 
657 
658  tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
659 
660  }
661 
662  if(hid.depth() == 1) {
663  hfNoiseE->Fill(energyhit);
664 
665  //if( hid.ieta() > 0 ) {
666  // hCalo1[hid.iphi()][hid.ieta()]->Fill(energyhit-noise_pl[hid.iphi()][hid.ieta()]);
667  // hCalo1mom2[hid.iphi()][hid.ieta()]->Fill(pow(energyhit,2));
668  //} else {
669  // hCalo2[hid.iphi()][abs(hid.ieta())]->Fill(energyhit-noise_min[hid.iphi()][abs(hid.ieta())]);
670  // hCalo2mom2[hid.iphi()][abs(hid.ieta())]->Fill(pow(energyhit,2));
671  //} // eta><0
672 
673  } // depth=1
674 
675  } // HBHE_NS
676 
677 
678 // Signal part for HB HE
679 
680  for(HFRecHitCollection::const_iterator hbheItr=HithfMB.begin(); hbheItr!=HithfMB.end(); hbheItr++)
681  {
682 // Recalibration of energy
683  float icalconst=1.;
684  DetId mydetid = hbheItr->id().rawId();
685  if( theRecalib ) icalconst=myRecalib->getValues(mydetid)->getValue();
686 
687  HFRecHit aHit(hbheItr->id(),hbheItr->energy()*icalconst,hbheItr->time());
688 
689  double energyhit = aHit.energy();
690 //
691 // Remove PMT hits
692 //
693  if(fabs(energyhit) > 40. ) continue;
694 
695  DetId id = (*hbheItr).detid();
696  HcalDetId hid=HcalDetId(id);
697 
698  int mysu = ((hid).rawId()>>25)&0x7;
699  if( hid.ieta() > 0 ) {
700  theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()]+ 1.;
701  theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+energyhit;
702  theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,2);
703  theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow(energyhit,4);
704 
705 
706  theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+energyhit-tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
707  theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
708  theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow((energyhit-tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]),2);
709  } else {
710  theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+ 1.;
711  theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+energyhit;
712  theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,2);
713  theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] = theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())]+pow(energyhit,4);
714 
715  theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()]+energyhit-tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
716  theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
717  theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()]+pow((energyhit-tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()]),2);
718  }
719 
720 
721  if(hid.depth() == 1) {
722  hfSignalE->Fill(energyhit);
723 
724  if( hid.ieta() > 0 )
725  {
726  hCalo1[hid.iphi()][hid.ieta()]->Fill(energyhit);
727  hCalo1mom2[hid.iphi()][hid.ieta()]->Fill(pow(energyhit,2));
728  } else
729  {
730  hCalo2[hid.iphi()][abs(hid.ieta())]->Fill(energyhit);
731  hCalo2mom2[hid.iphi()][abs(hid.ieta())]->Fill(pow(energyhit,2));
732  } // eta><0
733 
734  } // depth=1
735 
736  } // HF_MB
737 
738  std::cout<<" Event is finished "<<std::endl;
739 }
740 }
741 //define this as a plug-in
742 //DEFINE_ANOTHER_FWK_MODULE(Analyzer_minbias)
743 
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:501
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
i
int i
Definition: DBlmapReader.cc:9
hcalLocalRecoNZS_cff.horecoMB
tuple horecoMB
Definition: hcalLocalRecoNZS_cff.py:10
edm::RunBase::run
RunNumber_t run() const
Definition: RunBase.h:42
AlgorithmMap
std::map< std::string, L1GtAlgorithm > AlgorithmMap
map containing the algorithms
Definition: L1GtTriggerMenuFwd.h:31
Analyzer_minbias.h
edm::SortedCollection::const_iterator
std::vector< T >::const_iterator const_iterator
Definition: SortedCollection.h:80
abs
#define abs(x)
Definition: mlp_lapack.h:159
eta
T eta() const
Definition: Basic3DVectorLD.h:172
dt_dqm_sourceclient_common_cff.reco
tuple reco
Definition: dt_dqm_sourceclient_common_cff.py:105
bk::beginJob
void beginJob()
Definition: Breakpoints.cc:15
CItAlgo
AlgorithmMap::const_iterator CItAlgo
iterators through map containing the algorithms
Definition: L1GtTriggerMenuFwd.h:38
DetId::rawId
uint32_t rawId() const
get the raw id
Definition: DetId.h:45
nevent
int nevent
Definition: AMPTWrapper.h:74
HcalDetId::depth
int depth() const
get the tower depth
Definition: HcalDetId.h:42
iEvent
int iEvent
Definition: GenABIO.cc:243
HcalObjRepresent::Fill
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
Definition: HcalObjRepresent.h:611
DecisionWord
std::vector< bool > DecisionWord
typedefs
Definition: L1GlobalTriggerReadoutSetupFwd.h:36
CaloRecHit::energy
float energy() const
Definition: CaloRecHit.h:19
HcalDetId::ieta
int ieta() const
get the cell ieta
Definition: HcalDetId.h:38
edm::Event::run
RunNumber_t run() const
Definition: Event.h:67
hcalLocalRecoNZS_cff.hfrecoMB
tuple hfrecoMB
Definition: hcalLocalRecoNZS_cff.py:7
j
int j
Definition: DBlmapReader.cc:9
ALCARECOHcalCalMinBias_cff.hbherecoNoise
tuple hbherecoNoise
Definition: ALCARECOHcalCalMinBias_cff.py:10
ALCARECOHcalCalMinBias_cff.hfrecoNoise
tuple hfrecoNoise
Definition: ALCARECOHcalCalMinBias_cff.py:12
edm::Event::getByLabel
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
gen::k
int k[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:81
HcalDetId::iphi
int iphi() const
get the cell iphi
Definition: HcalDetId.h:40
DetId
Definition: DetId.h:20
edm::EventSetup::get
const T & get() const
Definition: EventSetup.h:55
edm::ESHandle::product
T const * product() const
Definition: ESHandle.h:62
L1GtTriggerMenu::gtAlgorithmMap
const AlgorithmMap & gtAlgorithmMap() const
get / set the algorithm map (by name)
Definition: L1GtTriggerMenu.h:239
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
ALCARECOHcalCalMinBias_cff.horecoNoise
tuple horecoNoise
Definition: ALCARECOHcalCalMinBias_cff.py:16
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
hcalLocalRecoNZS_cff.hbherecoMB
tuple hbherecoMB
Definition: hcalLocalRecoNZS_cff.py:4
edm::Run
Definition: Run.h:33
phi
Definition: DDAxes.h:10
Generated for CMSSW Reference Manual by doxygen 1.8.5