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