CMS 3D CMS Logo

EcalSelectiveReadoutValidation.cc
Go to the documentation of this file.
1 /*
2  * \file EcalSelectiveReadoutValidation.cc
3  *
4  *
5  */
6 
8 
10 
17 
25 
27 
34 
35 using namespace cms;
36 using namespace edm;
37 using namespace std;
38 
39 const double EcalSelectiveReadoutValidation::rad2deg = 45./atan(1.);
40 
42  //EE- DCCs:
43  // 1 2 3 4 5 6 7 8 9
44  34, 32, 33, 33, 32, 34, 33, 34, 33,
45  //EB- DCCs:
46  // 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27
47  68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68,
48  //EB+ DCCs:
49  // 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45
50  68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68, 68,
51  //EE+ DCCs:
52  // 46 47 48 49 50 51 52 53 54
53  32, 33, 33, 32, 34, 33, 34, 33, 34
54 };
55 
57  collNotFoundWarn_(ps.getUntrackedParameter<bool>("warnIfCollectionNotFound", true)),
58  ebDigis_(ps.getParameter<edm::InputTag>("EbDigiCollection"), false,
59  collNotFoundWarn_),
60  eeDigis_(ps.getParameter<edm::InputTag>("EeDigiCollection"), false,
61  collNotFoundWarn_),
62  ebNoZsDigis_(ps.getParameter<edm::InputTag>("EbUnsuppressedDigiCollection"),
63  false, false/*collNotFoundWarn_*/),
64  eeNoZsDigis_(ps.getParameter<edm::InputTag>("EeUnsuppressedDigiCollection"),
65  false, false/*collNotFoundWarn_*/),
66  ebSrFlags_(ps.getParameter<edm::InputTag>("EbSrFlagCollection"), false,
67  collNotFoundWarn_),
68  eeSrFlags_(ps.getParameter<edm::InputTag>("EeSrFlagCollection"), false,
69  collNotFoundWarn_),
70  ebComputedSrFlags_(ps.getParameter<edm::InputTag>("EbSrFlagFromTTCollection"), false,
71  false/*collNotFoundWarn_*/),
72  eeComputedSrFlags_(ps.getParameter<edm::InputTag>("EeSrFlagFromTTCollection"), false,
73  false/*collNotFoundWarn_*/),
74  ebSimHits_(ps.getParameter<edm::InputTag>("EbSimHitCollection"), false,
75  false/*collNotFoundWarn_*/),
76  eeSimHits_(ps.getParameter<edm::InputTag>("EeSimHitCollection"), false,
77  false/*collNotFoundWarn_*/),
78  tps_(ps.getParameter<edm::InputTag>("TrigPrimCollection"), false,
79  collNotFoundWarn_),
80  ebRecHits_(ps.getParameter<edm::InputTag>("EbRecHitCollection"), false,
81  false/*collNotFoundWarn_*/),
82  eeRecHits_(ps.getParameter<edm::InputTag>("EeRecHitCollection"), false,
83  false/*collNotFoundWarn_*/),
84  fedRaw_(ps.getParameter<edm::InputTag>("FEDRawCollection"), false,
85  false/*collNotFoundWarn_*/),
86  tmax(0),
87  tmin(numeric_limits<int64_t>::max()),
88  l1aOfTmin(0),
89  l1aOfTmax(0),
90  triggerTowerMap_(0),
91  localReco_(ps.getParameter<bool>("LocalReco")),
92  weights_(ps.getParameter<vector<double> >("weights")),
93  tpInGeV_(ps.getParameter<bool>("tpInGeV")),
94  firstFIRSample_(ps.getParameter<int>("ecalDccZs1stSample")),
95  useEventRate_(ps.getParameter<bool>("useEventRate")),
96  logErrForDccs_(nDccs_, false),
97  ievt_(0),
98  allHists_(false),
99  histDir_(ps.getParameter<string>("histDir")),
100  withEeSimHit_(false),
101  withEbSimHit_(false){
102 
104  ebDigis_.setToken(collector);
105  eeDigis_.setToken(collector);
106  ebNoZsDigis_.setToken(collector);
107  eeNoZsDigis_.setToken(collector);
108  ebSrFlags_.setToken(collector);
109  eeSrFlags_.setToken(collector);
110  ebComputedSrFlags_.setToken(collector);
111  eeComputedSrFlags_.setToken(collector);
112  ebSimHits_.setToken(collector);
113  eeSimHits_.setToken(collector);
114  tps_.setToken(collector);
115  ebRecHits_.setToken(collector);
116  eeRecHits_.setToken(collector);
117  fedRaw_.setToken(collector);
118 
119  double ebZsThr = ps.getParameter<double>("ebZsThrADCCount");
120  double eeZsThr = ps.getParameter<double>("eeZsThrADCCount");
121 
122  ebZsThr_ = lround(ebZsThr*4);
123  eeZsThr_ = lround(eeZsThr*4);
124 
125  //File to log SRP algorithem inconsistency
127  = ps.getUntrackedParameter<string>("srpAlgoErrorLogFile","");
129 
130  //File to log SRP decision application inconsistency
132  = ps.getUntrackedParameter<string>("srApplicationErrorLogFile","");
134 
135  //FIR ZS weights
136  configFirWeights(ps.getParameter<vector<double> >("dccWeights"));
137 
138  // DQM ROOT output
139  outputFile_ = ps.getUntrackedParameter<string>("outputFile", "");
140 
141  if(outputFile_.size() != 0){
142  LogInfo("OutputInfo") << " Ecal Digi Task histograms will be saved to '"
143  << outputFile_.c_str() << "'";
144  } else{
145  LogInfo("OutputInfo") << " Ecal Digi Task histograms will NOT be saved";
146  }
147 
148  // verbosity switch
149  verbose_ = ps.getUntrackedParameter<bool>("verbose", false);
150 
151  // get hold of back-end interface
152 
153  vector<string>
154  hists(ps.getUntrackedParameter<vector<string> >("histograms",
155  vector<string>(1, "all")));
156 
157  for(vector<string>::iterator it = hists.begin();
158  it!=hists.end(); ++it) histList_.insert(*it);
159  if(histList_.find("all") != histList_.end()) allHists_ = true;
160 
161 }
162 
163 
165  const int32_t bx = event.bunchCrossing();
166  if(bx<1 || bx > 3564) return;
167 
168  int64_t t = event.bunchCrossing() + (event.orbitNumber()-1)*3564;
169 
170  if(t<tmin){
171  tmin = t;
172  l1aOfTmin = event.id().event();
173  }
174 
175  if(t>tmax){
176  tmax = t;
177  l1aOfTmax = event.id().event();
178  }
179 }
180 
182  LogDebug("EcalSrValid") << __FILE__ << ":" << __LINE__ << ": "
183  << "Tmax = " << tmax << " x 25ns; Tmin = " << tmin
184  << " x 25ns; L1A(Tmax) = " << l1aOfTmax << "; L1A(Tmin) = "
185  << l1aOfTmin << "\n";
186  return (double)(l1aOfTmax - l1aOfTmin) / ((tmax-tmin) * 25e-9);
187 }
188 
190  EventSetup const & es){
191 
192  updateL1aRate(event);
193 
194  //retrieves event products:
195  readAllCollections(event);
196 
197  withEeSimHit_ = (eeSimHits_->size()!=0);
198  withEbSimHit_ = (ebSimHits_->size()!=0);
199 
200  if(ievt_<10){
201  edm::LogInfo("EcalSrValid") << "Size of TP collection: " << tps_->size() << std::endl
202  << "Size of EB SRF collection read from data: "
203  << ebSrFlags_->size() << std::endl
204  << "Size of EB SRF collection computed from data TTFs: "
205  << ebComputedSrFlags_->size() << std::endl
206  << "Size of EE SRF collection read from data: "
207  << eeSrFlags_->size() << std::endl
208  << "Size of EE SRF collection computed from data TTFs: "
209  << eeComputedSrFlags_->size() << std::endl;
210  }
211 
212  if(ievt_==0){
213  selectFedsForLog(); //note: must be called after readAllCollection
214  }
215 
216  //computes Et sum trigger tower crystals:
218 
219  //Data Volume
220  analyzeDataVolume(event, es);
221 
222  //EB digis
223  //must be called after analyzeDataVolume because it uses
224  //isRuComplete_ array that this method fills
225  analyzeEB(event, es);
226 
227  //EE digis
228  //must be called after analyzeDataVolume because it uses
229  //isRuComplete_ array that this method fills
230  analyzeEE(event, es);
231 
235 
239 
243 
244  //TP
245  analyzeTP(event, es);
246 
247  if(ebComputedSrFlags_->size()){
249  }
250  if(eeComputedSrFlags_->size()){
252  }
253  nDroppedFRO_ = 0;
254  nIncompleteFRO_ = 0;
255  nCompleteZS_ = 0;
261  ++ievt_;
262 }
263 
264 
266  const edm::EventSetup& es){
267  bool eventError = false;
268  nEeZsErrors_ = 0;
269 
270  for(int iZ0=0; iZ0<nEndcaps; ++iZ0){
271  for(int iX0=0; iX0<nEeX; ++iX0){
272  for(int iY0=0; iY0<nEeY; ++iY0){
273  eeEnergies[iZ0][iX0][iY0].noZsRecE = -numeric_limits<double>::max();
274  eeEnergies[iZ0][iX0][iY0].recE = -numeric_limits<double>::max();
275  eeEnergies[iZ0][iX0][iY0].simE = 0; //must be set to zero.
276  eeEnergies[iZ0][iX0][iY0].simHit = 0;
277  eeEnergies[iZ0][iX0][iY0].gain12 = false;
278  }
279  }
280  }
281 
282  // gets the endcap geometry:
283  edm::ESHandle<CaloGeometry> geoHandle;
284  es.get<CaloGeometryRecord>().get(geoHandle);
285  const CaloSubdetectorGeometry *geometry_p
286  = (*geoHandle).getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
287  CaloSubdetectorGeometry const& geometry = *geometry_p;
288 
289  //EE unsupressed digis:
290  for (unsigned int digis=0; digis<eeNoZsDigis_->size(); ++digis){
291 
292  EEDataFrame frame = (*eeNoZsDigis_)[digis];
293  int iX0 = iXY2cIndex(frame.id().ix());
294  int iY0 = iXY2cIndex(frame.id().iy());
295  int iZ0 = frame.id().zside()>0?1:0;
296 
297  if(iX0<0 || iX0>=nEeX){
298  edm::LogError("EcalSrValid") << "iX0 (= " << iX0 << ") is out of range ("
299  << "[0," << nEeX -1 << "]\n";
300  }
301  if(iY0<0 || iY0>=nEeY){
302  edm::LogError("EcalSrValid") << "iY0 (= " << iY0 << ") is out of range ("
303  << "[0," << nEeY -1 << "]\n";
304  }
305  // cout << "EE no ZS energy computation..." ;
306  eeEnergies[iZ0][iX0][iY0].noZsRecE = frame2Energy(frame);
307 
308  eeEnergies[iZ0][iX0][iY0].gain12 = true;
309  for(int i = 0; i< frame.size(); ++i){
310  const int gain12Code = 0x1;
311  if(frame[i].gainId()!=gain12Code) eeEnergies[iZ0][iX0][iY0].gain12 = false;
312  }
313 
314  const GlobalPoint xtalPos
315  = geometry.getGeometry(frame.id())->getPosition();
316 
317  eeEnergies[iZ0][iX0][iY0].phi = rad2deg*((double)xtalPos.phi());
318  eeEnergies[iZ0][iX0][iY0].eta = xtalPos.eta();
319  }
320 
321  //EE rec hits:
322  if(!localReco_){
324  = eeRecHits_->begin();
325  it != eeRecHits_->end(); ++it){
326  const RecHit& hit = *it;
327  int iX0 = iXY2cIndex(static_cast<const EEDetId&>(hit.id()).ix());
328  int iY0 = iXY2cIndex(static_cast<const EEDetId&>(hit.id()).iy());
329  int iZ0 = static_cast<const EEDetId&>(hit.id()).zside()>0?1:0;
330 
331  if(iX0<0 || iX0>=nEeX){
332  LogError("EcalSrValid") << "iX0 (= " << iX0 << ") is out of range ("
333  << "[0," << nEeX -1 << "]\n";
334  }
335  if(iY0<0 || iY0>=nEeY){
336  LogError("EcalSrValid") << "iY0 (= " << iY0 << ") is out of range ("
337  << "[0," << nEeY -1 << "]\n";
338  }
339  // cout << "EE no ZS energy computation..." ;
340  eeEnergies[iZ0][iX0][iY0].recE = hit.energy();
341  }
342  }
343 
344  //EE sim hits:
345  for(vector<PCaloHit>::const_iterator it = eeSimHits_->begin();
346  it != eeSimHits_->end(); ++it){
347  const PCaloHit& simHit = *it;
348  EEDetId detId(simHit.id());
349  int iX = detId.ix();
350  int iX0 =iXY2cIndex(iX);
351  int iY = detId.iy();
352  int iY0 = iXY2cIndex(iY);
353  int iZ0 = detId.zside()>0?1:0;
354  eeEnergies[iZ0][iX0][iY0].simE += simHit.energy();
355  ++eeEnergies[iZ0][iX0][iY0].simHit;
356  }
357 
358  bool EEcrystalShot[nEeX][nEeY][2];
359  pair<int,int> EExtalCoor[nEeX][nEeY][2];
360 
361  for(int iEeZ=0; iEeZ<2; ++iEeZ){
362  for(int iEeX=0; iEeX<nEeX; ++iEeX){
363  for(int iEeY=0; iEeY<nEeY; ++iEeY){
364  EEcrystalShot[iEeX][iEeY][iEeZ] = false;
365  EExtalCoor[iEeX][iEeY][iEeZ] = make_pair(0,0);
366  }
367  }
368  }
369 
370  //EE suppressed digis
372  it != eeDigis_->end(); ++it){
373  const EEDataFrame& frame = *it;
374  int iX0 = iXY2cIndex(static_cast<const EEDetId&>(frame.id()).ix());
375  int iY0 = iXY2cIndex(static_cast<const EEDetId&>(frame.id()).iy());
376  int iZ0 = static_cast<const EEDetId&>(frame.id()).zside()>0?1:0;
377  if(iX0<0 || iX0>=nEeX){
378  LogError("EcalSrValid") << "iX0 (= " << iX0 << ") is out of range ("
379  << "[0," << nEeX -1 << "]\n";
380  }
381  if(iY0<0 || iY0>=nEeY){
382  LogError("EcalSrValid") << "iY0 (= " << iY0 << ") is out of range ("
383  << "[0," << nEeY -1 << "]\n";
384  }
385 
386  if(!EEcrystalShot[iX0][iY0][iZ0]){
387  EEcrystalShot[iX0][iY0][iZ0] = true;
388  EExtalCoor[iX0][iY0][iZ0] = make_pair(xtalGraphX(frame.id()), xtalGraphY(frame.id()));
389  } else{
390  cout << "Error: several digi for same crystal!";
391  abort();
392  }
393 
394  if(localReco_){
395  eeEnergies[iZ0][iX0][iY0].recE = frame2Energy(frame);
396  }
397 
398  eeEnergies[iZ0][iX0][iY0].gain12 = true;
399  for(int i = 0; i< frame.size(); ++i){
400  const int gain12Code = 0x1;
401  if(frame[i].gainId()!=gain12Code){
402  eeEnergies[iZ0][iX0][iY0].gain12 = false;
403  }
404  }
405 
407  = eeSrFlags_->find(readOutUnitOf(frame.id()));
408 
409  bool highInterest = false;
410 
411 
412  if(srf==eeSrFlags_->end()) continue;
413 
414  if(srf!=eeSrFlags_->end()){
415  highInterest = ((srf->value() & ~EcalSrFlag::SRF_FORCED_MASK)
417  }
418 
419  if(highInterest){
421  } else{
422  int v = dccZsFIR(frame, firWeights_, firstFIRSample_, 0);
423  fill(meEeLiZsFir_, v);
424  if(v < eeZsThr_){
425  eventError = true;
426  ++nEeZsErrors_;
427  pair<int,int> ru = dccCh(frame.id());
428  if(isRuComplete_[ru.first][ru.second-1]) ++nEeZsErrorsType1_;
429  if(nEeZsErrors_ < 3){
430  srApplicationErrorLog_ << event.id() << ", "
431  << "RU " << frame.id() << ", "
432  << "DCC " << ru.first
433  << " Ch : " << ru.second << ": "
434  << "LI channel under ZS threshold.\n";
435  }
436  if(nEeZsErrors_==3){
437  srApplicationErrorLog_ << event.id() << ": "
438  << "more ZS errors for this event...\n";
439  }
440  }
441  }
442  } //next ZS digi.
443 
444  for(int iEeZ=0; iEeZ<2; ++iEeZ){
445  for(int iEeX=0; iEeX<nEeX; ++iEeX){
446  for(int iEeY=0; iEeY<nEeY; ++iEeY){
447  fill(meChOcc_, EExtalCoor[iEeX][iEeY][iEeZ].first, EExtalCoor[iEeX][iEeY][iEeZ].second,EEcrystalShot[iEeX][iEeY][iEeZ] ? 1 : 0);
448  }
449  }
450  }
451 
452  for(int iZ0=0; iZ0<nEndcaps; ++iZ0){
453  for(int iX0=0; iX0<nEeX; ++iX0){
454  for(int iY0=0; iY0<nEeY; ++iY0){
455  double recE = eeEnergies[iZ0][iX0][iY0].recE;
456  if(recE==-numeric_limits<double>::max()) continue; //not a crystal or ZS
457  fill(meEeRecE_, eeEnergies[iZ0][iX0][iY0].recE);
458 
459  fill(meEeEMean_, ievt_+1,
460  eeEnergies[iZ0][iX0][iY0].recE);
461 
462  if(withEeSimHit_){
463  if(!eeEnergies[iZ0][iX0][iY0].simHit){//noise only crystal channel
464  fill(meEeNoise_, eeEnergies[iZ0][iX0][iY0].noZsRecE);
465  } else{
466  fill(meEeSimE_, eeEnergies[iZ0][iX0][iY0].simE);
467  fill(meEeRecEHitXtal_, eeEnergies[iZ0][iX0][iY0].recE);
468  }
469  fill(meEeRecVsSimE_, eeEnergies[iZ0][iX0][iY0].simE,
470  eeEnergies[iZ0][iX0][iY0].recE);
471  fill(meEeNoZsRecVsSimE_, eeEnergies[iZ0][iX0][iY0].simE,
472  eeEnergies[iZ0][iX0][iY0].noZsRecE);
473  }
474  }
475  }
476  }
477 
478  int EEZs1RuCount[2][20][20];
479  int EEFullRuCount[2][20][20];
480  int EEForcedRuCount[2][20][20];
481  for(int iZ(0); iZ < 2; iZ++){
482  for(int iX(0);iX < 20; iX++){
483  for(int iY(0);iY < 20; iY++){
484  EEZs1RuCount[iZ][iX][iY] = 0;
485  EEFullRuCount[iZ][iX][iY] = 0;
486  EEForcedRuCount[iZ][iX][iY] = 0;
487  }
488  }
489  }
490 
491  nEeFROCnt_ = 0;
492  char eeSrfMark[2][20][20];
493  bzero(eeSrfMark, sizeof(eeSrfMark));
494  //Filling RU histo
496  it != eeSrFlags_->end(); ++it){
497  const EESrFlag& srf = *it;
498  // srf.id() is EcalScDetId; 1 <= ix <= 20 1 <= iy <= 20
499  int iX = srf.id().ix();
500  int iY = srf.id().iy();
501  int zside = srf.id().zside(); //-1 for EE-, +1 for EE+
502  if(iX<1 || iY > 100) throw cms::Exception("EcalSelectiveReadoutValidation")
503  << "Found an endcap SRF with an invalid det ID: " << srf.id() << ".\n";
504  ++eeSrfMark[zside>0?1:0][iX-1][iY-1];
505  if(eeSrfMark[zside>0?1:0][iX-1][iY-1] > 1) throw cms::Exception("EcalSelectiveReadoutValidation")
506  << "Duplicate SRF for supercrystal " << srf.id() << ".\n";
507  int flag = srf.value() & ~EcalSrFlag::SRF_FORCED_MASK;
508  if(flag == EcalSrFlag::SRF_ZS1){
509  EEZs1RuCount[zside>0?1:0][iX-1][iY-1] += 1;
510  }
511 
512  if(flag == EcalSrFlag::SRF_FULL){
513  EEFullRuCount[zside>0?1:0][iX-1][iY-1] += 1;
514  ++nEeFROCnt_;
515  }
516 
518  EEForcedRuCount[zside>0?1:0][iX-1][iY-1] += 1;
519  }
520  }
521  for(int iZ(0); iZ < 2; iZ++){
522  int xOffset(iZ == 0 ? -40 : 20);
523  for(int iX(0);iX < 20; iX++){
524  for(int iY(0);iY < 20; iY++){
525  int GraphX = (iX + 1)+ xOffset;
526  int GraphY = (iY+1);
527  fill(meZs1Ru_, GraphX, GraphY, EEZs1RuCount[iZ][iX][iY]);
528  fill(meFullRoRu_, GraphX, GraphY, EEFullRuCount[iZ][iX][iY]);
529  fill(meForcedRu_, GraphX, GraphY, EEForcedRuCount[iZ][iX][iY]);
530  }
531  }
532  }
533 
534  if(eventError) srApplicationErrorLog_ << event.id()
535  << ": " << nEeZsErrors_
536  << " ZS-flagged EE channels under "
537  "the ZS threshold, whose " << nEeZsErrorsType1_
538  << " in a complete RU.\n";
539 } //end of analyzeEE
540 
541 void
543  const edm::EventSetup& es){
544 
545  bool eventError = false;
546  nEbZsErrors_ = 0;
547  vector<pair<int,int> > xtalEtaPhi;
548 
549  xtalEtaPhi.reserve(nEbPhi*nEbEta);
550  for(int iEta0=0; iEta0<nEbEta; ++iEta0){
551  for(int iPhi0=0; iPhi0<nEbPhi; ++iPhi0){
553  ebEnergies[iEta0][iPhi0].recE = -numeric_limits<double>::max();
554  ebEnergies[iEta0][iPhi0].simE = 0; //must be zero.
555  ebEnergies[iEta0][iPhi0].simHit = 0;
556  ebEnergies[iEta0][iPhi0].gain12 = false;
557  xtalEtaPhi.push_back(pair<int,int>(iEta0, iPhi0));
558  }
559  }
560 
561  // get the barrel geometry:
562  edm::ESHandle<CaloGeometry> geoHandle;
563 
564  es.get<CaloGeometryRecord>().get(geoHandle);
565  const CaloSubdetectorGeometry *geometry_p
566  = (*geoHandle).getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
567  CaloSubdetectorGeometry const& geometry = *geometry_p;
568 
569  //EB unsuppressed digis:
571  it != ebNoZsDigis_->end(); ++it){
572  const EBDataFrame& frame = *it;
573  int iEta0 = iEta2cIndex(static_cast<const EBDetId&>(frame.id()).ieta());
574  int iPhi0 = iPhi2cIndex(static_cast<const EBDetId&>(frame.id()).iphi());
575  if(iEta0<0 || iEta0>=nEbEta){
576  stringstream s;
577  s << "EcalSelectiveReadoutValidation: "
578  << "iEta0 (= " << iEta0 << ") is out of range ("
579  << "[0," << nEbEta -1 << "]\n";
580  throw cms::Exception(s.str());
581  }
582  if(iPhi0<0 || iPhi0>=nEbPhi){
583  stringstream s;
584  s << "EcalSelectiveReadoutValidation: "
585  << "iPhi0 (= " << iPhi0 << ") is out of range ("
586  << "[0," << nEbPhi -1 << "]\n";
587  throw cms::Exception(s.str());
588  }
589 
590  ebEnergies[iEta0][iPhi0].noZsRecE = frame2Energy(frame);
591  ebEnergies[iEta0][iPhi0].gain12 = true;
592  for(int i = 0; i< frame.size(); ++i){
593  const int gain12Code = 0x1;
594  if(frame[i].gainId()!=gain12Code) ebEnergies[iEta0][iPhi0].gain12 = false;
595  }
596 
597  const GlobalPoint xtalPos
598  = geometry.getGeometry(frame.id())->getPosition();
599 
600  ebEnergies[iEta0][iPhi0].phi = rad2deg*((double)xtalPos.phi());
601  ebEnergies[iEta0][iPhi0].eta = xtalPos.eta();
602  } //next non-zs digi
603 
604  //EB sim hits
605  for(vector<PCaloHit>::const_iterator it = ebSimHits_->begin();
606  it != ebSimHits_->end(); ++it){
607  const PCaloHit& simHit = *it;
608  EBDetId detId(simHit.id());
609  int iEta = detId.ieta();
610  int iEta0 =iEta2cIndex(iEta);
611  int iPhi = detId.iphi();
612  int iPhi0 = iPhi2cIndex(iPhi);
613  ebEnergies[iEta0][iPhi0].simE += simHit.energy();
614  ++ebEnergies[iEta0][iPhi0].simHit;
615  }
616 
617  bool crystalShot[nEbEta][nEbPhi];
618  pair<int,int> EBxtalCoor[nEbEta][nEbPhi];
619 
620  for(int iEta0=0; iEta0<nEbEta; ++iEta0){
621  for(int iPhi0=0; iPhi0<nEbPhi; ++iPhi0){
622  crystalShot[iEta0][iPhi0] = false;
623  EBxtalCoor[iEta0][iPhi0] = make_pair(0,0);
624  }
625  }
626 
627  int nEbDigi = 0;
628 
630  it != ebDigis_->end(); ++it){
631  ++nEbDigi;
632  const EBDataFrame& frame = *it;
633  int iEta = static_cast<const EBDetId&>(frame.id()).ieta();
634  int iPhi = static_cast<const EBDetId&>(frame.id()).iphi();
635  int iEta0 = iEta2cIndex(iEta);
636  int iPhi0 = iPhi2cIndex(iPhi);
637  if(iEta0<0 || iEta0>=nEbEta){
638  throw (cms::Exception("EcalSelectiveReadoutValidation")
639  << "iEta0 (= " << iEta0 << ") is out of range ("
640  << "[0," << nEbEta -1 << "]");
641  }
642  if(iPhi0<0 || iPhi0>=nEbPhi){
643  throw (cms::Exception("EcalSelectiveReadoutValidation")
644  << "iPhi0 (= " << iPhi0 << ") is out of range ("
645  << "[0," << nEbPhi -1 << "]");
646  }
647  assert(iEta0>=0 && iEta0<nEbEta);
648  assert(iPhi0>=0 && iPhi0<nEbPhi);
649  if(!crystalShot[iEta0][iPhi0]){
650  crystalShot[iEta0][iPhi0] = true;
651  EBxtalCoor[iEta0][iPhi0] = make_pair(xtalGraphX(frame.id()), xtalGraphY(frame.id()));
652  } else{
653  cout << "Error: several digi for same crystal!";
654  abort();
655  }
656  if(localReco_){
657  ebEnergies[iEta0][iPhi0].recE = frame2Energy(frame);
658  }
659 
660  ebEnergies[iEta0][iPhi0].gain12 = true;
661  for(int i = 0; i< frame.size(); ++i){
662  const int gain12Code = 0x1;
663  if(frame[i].gainId()!=gain12Code){
664  ebEnergies[iEta0][iPhi0].gain12 = false;
665  }
666  }
667 
669  = ebSrFlags_->find(readOutUnitOf(frame.id()));
670 
671  bool highInterest = false;
672 
673  // if(srf == ebSrFlags_->end()){
674  // throw cms::Exception("EcalSelectiveReadoutValidation")
675  // << __FILE__ << ":" << __LINE__ << ": SR flag not found";
676  //}
677 
678  if(srf != ebSrFlags_->end()){
679  highInterest = ((srf->value() & ~EcalSrFlag::SRF_FORCED_MASK)
681  }
682 
683 
684  if(highInterest){
686  } else{
687  int v = dccZsFIR(frame, firWeights_, firstFIRSample_, 0);
688  fill(meEbLiZsFir_, v);
689  if(v < ebZsThr_){
690  eventError = true;
691  ++nEbZsErrors_;
692  pair<int,int> ru = dccCh(frame.id());
693  if(isRuComplete_[ru.first][ru.second-1]) ++nEbZsErrorsType1_;
694  if(nEbZsErrors_ < 3){
695  srApplicationErrorLog_ << event.id() << ", "
696  << "RU " << frame.id() << ", "
697  << "DCC " << ru.first
698  << " Ch : " << ru.second << ": "
699  << "LI channel under ZS threshold.\n";
700  }
701  if(nEbZsErrors_==3){
702  srApplicationErrorLog_ << event.id() << ": "
703  << "more ZS errors for this event...\n";
704  }
705  }
706  }
707  } //next EB digi
708 
709  for(int iEta0=0; iEta0<nEbEta; ++iEta0){
710  for(int iPhi0=0; iPhi0<nEbPhi; ++iPhi0)
711  fill(meChOcc_, EBxtalCoor[iEta0][iPhi0].first, EBxtalCoor[iEta0][iPhi0].second,crystalShot[iEta0][iPhi0] ? 1. : 0.);
712  }
713 
714  if(!localReco_){
716  = ebRecHits_->begin();
717  it != ebRecHits_->end(); ++it){
718  ++nEbDigi;
719  const RecHit& hit = *it;
720  int iEta = static_cast<const EBDetId&>(hit.id()).ieta();
721  int iPhi = static_cast<const EBDetId&>(hit.id()).iphi();
722  int iEta0 = iEta2cIndex(iEta);
723  int iPhi0 = iPhi2cIndex(iPhi);
724  if(iEta0<0 || iEta0>=nEbEta){
725  LogError("EcalSrValid") << "iEta0 (= " << iEta0 << ") is out of range ("
726  << "[0," << nEbEta -1 << "]\n";
727  }
728  if(iPhi0<0 || iPhi0>=nEbPhi){
729  LogError("EcalSrValid") << "iPhi0 (= " << iPhi0 << ") is out of range ("
730  << "[0," << nEbPhi -1 << "]\n";
731  }
732  ebEnergies[iEta0][iPhi0].recE = hit.energy();
733  }
734  }
735 
736  for(unsigned int i=0; i<xtalEtaPhi.size(); ++i){
737  int iEta0 = xtalEtaPhi[i].first;
738  int iPhi0= xtalEtaPhi[i].second;
739  energiesEb_t& energies = ebEnergies[iEta0][iPhi0];
740 
741  double recE = energies.recE;
742  if(recE!=-numeric_limits<double>::max()){//not zero suppressed
743  fill(meEbRecE_, ebEnergies[iEta0][iPhi0].recE);
744  fill(meEbEMean_, ievt_+1, recE);
745  } //not zero suppressed
746 
747  if(withEbSimHit_){
748  if(!energies.simHit){//noise only crystal channel
749  fill(meEbNoise_, energies.noZsRecE);
750  } else{
751  fill(meEbSimE_, energies.simE);
752  fill(meEbRecEHitXtal_, energies.recE);
753  }
754  fill(meEbRecVsSimE_, energies.simE, energies.recE);
755  fill(meEbNoZsRecVsSimE_, energies.simE, energies.noZsRecE);
756  }
757  }
758 
759  int EBZs1RuCount[2][17][72];
760  int EBFullRuCount[2][17][72];
761  int EBForcedRuCount[2][17][72];
762  std::pair<int, int> EBtowerCoor[2][17][72];
763  for(int iZ(0); iZ < 2; iZ++){
764  for(int iEta(0);iEta < 17; iEta++){
765  for(int iPhi(0);iPhi < 72; iPhi++){
766  EBZs1RuCount[iZ][iEta][iPhi] = 0;
767  EBFullRuCount[iZ][iEta][iPhi] = 0;
768  EBForcedRuCount[iZ][iEta][iPhi] = 0;
769  }
770  }
771  }
772 
773  //SRF
774  nEbFROCnt_ = 0;
775  char ebSrfMark[2][17][72];
776  bzero(ebSrfMark, sizeof(ebSrfMark));
777  // int idbg = 0;
779  it != ebSrFlags_->end(); ++it){
780  const EBSrFlag& srf = *it;
781  int iEtaAbs = srf.id().ietaAbs();
782  int iPhi = srf.id().iphi();
783  int iZ = srf.id().zside();
784 
785  // cout << "--> " << ++idbg << iEtaAbs << " " << iPhi << " " << iZ
786  // << " " << srf.id() << "\n";
787 
788  if(iEtaAbs < 1 || iEtaAbs > 17
789  || iPhi < 1 || iPhi > 72) throw cms::Exception("EcalSelectiveReadoutValidation")
790  << "Found a barrel SRF with an invalid det ID: " << srf.id() << ".\n";
791  ++ebSrfMark[iZ>0?1:0][iEtaAbs-1][iPhi-1];
792  if(ebSrfMark[iZ>0?1:0][iEtaAbs-1][iPhi-1] > 1) throw cms::Exception("EcalSelectiveReadoutValidation")
793  << "Duplicate SRF for RU " << srf.id() << ".\n";
794 
795  EBtowerCoor[iZ > 0 ? 1 : 0][iEtaAbs - 1][iPhi - 1] = std::pair<int, int>(srf.id().ieta(), srf.id().iphi());
796 
797  int flag = srf.value() & ~EcalSrFlag::SRF_FORCED_MASK;
798  if(flag == EcalSrFlag::SRF_ZS1){
799  EBZs1RuCount[iZ>0?1:0][iEtaAbs-1][iPhi-1] +=1;
800  }
801  if(flag == EcalSrFlag::SRF_FULL){
802  EBFullRuCount[iZ>0?1:0][iEtaAbs-1][iPhi-1] +=1;
803  ++nEbFROCnt_;
804  }
806  EBForcedRuCount[iZ>0?1:0][iEtaAbs-1][iPhi-1] +=1;
807  }
808  }
809  for(int iZ(0); iZ < 2; iZ++){
810  for(int iEta(0);iEta < 17; iEta++){
811  for(int iPhi(0);iPhi < 72; iPhi++){
812  float x(EBtowerCoor[iZ][iEta][iPhi].first);
813  float y(EBtowerCoor[iZ][iEta][iPhi].second);
814  fill(meZs1Ru_, x, y, EBZs1RuCount[iZ][iEta][iPhi]);
815  fill(meFullRoRu_, x, y, EBFullRuCount[iZ][iEta][iPhi]);
816  fill(meForcedRu_, x, y, EBForcedRuCount[iZ][iEta][iPhi]);
817  }
818  }
819  }
820 
821  if(eventError) srApplicationErrorLog_ << event.id()
822  << ": " << nEbZsErrors_
823  << " ZS-flagged EB channels under "
824  "the ZS threshold, whose " << nEbZsErrorsType1_
825  << " in a complete RU.\n";
826 }
827 
829 }
830 
832  // endcap mapping
834  es.get<IdealGeometryRecord>().get(hTriggerTowerMap);
835  triggerTowerMap_ = hTriggerTowerMap.product();
836 
837  //electronics map
839  es.get< EcalMappingRcd >().get(ecalmapping);
840  elecMap_ = ecalmapping.product();
841 
842  initAsciiFile();
843 }
844 
847 }
848 
850 
851  meL1aRate_ = bookFloat(ibooker, "l1aRate_");
852 
853  meDccVol_ = bookProfile(ibooker, "hDccVol", //"EcalDccEventSizeComputed",
854  "ECAL DCC event fragment size;Dcc id; "
855  "<Event size> (kB)", nDccs_, .5, .5+nDccs_);
856 
857  meDccLiVol_ = bookProfile(ibooker, "hDccLiVol",
858  "LI channel payload per DCC;Dcc id; "
859  "<Event size> (kB)", nDccs_, .5, .5+nDccs_);
860 
861  meDccHiVol_ = bookProfile(ibooker, "hDccHiVol",
862  "HI channel payload per DCC;Dcc id; "
863  "<Event size> (kB)", nDccs_, .5, .5+nDccs_);
864 
865  meDccVolFromData_ = bookProfile(ibooker, "hDccVolFromData", //"EcalDccEventSize",
866  "ECAL DCC event fragment size;Dcc id; "
867  "<Event size> (kB)", nDccs_, .5, .5+nDccs_);
868 
869  meVolBLI_ = book1D(ibooker, "hVolBLI",// "EBLowInterestPayload",
870  "ECAL Barrel low interest crystal data payload;"
871  "Event size (kB);Nevts",
872  100, 0., 200.);
873 
874  meVolELI_ = book1D(ibooker, "hVolELI", //"EELowInterestPayload",
875  "Endcap low interest crystal data payload;"
876  "Event size (kB);Nevts",
877  100, 0., 200.);
878 
879  meVolLI_ = book1D(ibooker, "hVolLI", //"EcalLowInterestPayload",
880  "ECAL low interest crystal data payload;"
881  "Event size (kB);Nevts",
882  100, 0., 200.);
883 
884  meVolBHI_ = book1D(ibooker, "hVolBHI", //"EBHighInterestPayload",
885  "Barrel high interest crystal data payload;"
886  "Event size (kB);Nevts",
887  100, 0., 200.);
888 
889  meVolEHI_ = book1D(ibooker, "hVolEHI", //"EEHighInterestPayload",
890  "Endcap high interest crystal data payload;"
891  "Event size (kB);Nevts",
892  100, 0., 200.);
893 
894  meVolHI_ = book1D(ibooker, "hVolHI", //"EcalHighInterestPayload",
895  "ECAL high interest crystal data payload;"
896  "Event size (kB);Nevts",
897  100, 0., 200.);
898 
899  meVolB_ = book1D(ibooker, "hVolB", //"EBEventSize",
900  "Barrel data volume;Event size (kB);Nevts",
901  100, 0., 200.);
902 
903  meVolE_ = book1D(ibooker, "hVolE", //"EEEventSize",
904  "Endcap data volume;Event size (kB);Nevts",
905  100, 0., 200.);
906 
907  meVol_ = book1D(ibooker, "hVol", //"EcalEventSize",
908  "ECAL data volume;Event size (kB);Nevts",
909  100, 0., 200.);
910 
911  meChOcc_ = bookProfile2D(ibooker, "h2ChOcc", //"EcalChannelOccupancy",
912  "ECAL crystal channel occupancy after zero suppression;"
913  "iX -200 / iEta / iX + 100;"
914  "iY / iPhi (starting from -10^{o}!);"
915  "Event count rate",
916  401, -200.5, 200.5,
917  360, .5, 360.5);
918 
919  //TP
920  string tpUnit;
921  if(tpInGeV_) tpUnit = string("GeV"); else tpUnit = string("TP hw unit");
922  string title;
923  title = string("Trigger primitive TT E_{T};E_{T} ")
924  + tpUnit + string(";Event Count");
925  meTp_ = bookProfile(ibooker, "hTp", //"EcalTriggerPrimitiveEt",
926  title.c_str(),
927  (tpInGeV_?100:40), 0., (tpInGeV_?10.:40.));
928 
929  meTtf_ = bookProfile(ibooker, "hTtf", //"EcalTriggerTowerFlag",
930  "Trigger primitive TT flag;Flag number;Event count",
931  8, -.5, 7.5);
932 
933  title = string("Trigger tower flag vs TP;E_{T}(TT) (")
934  + tpUnit + string(");Flag number");
935  meTtfVsTp_ = book2D(ibooker, "h2TtfVsTp",
936  title.c_str(),
937  100, 0., (tpInGeV_?10.:40.),
938  8, -.5, 7.5);
939 
940  meTtfVsEtSum_ = book2D(ibooker, "h2TtfVsEtSum",
941  "Trigger tower flag vs #sumE_{T};"
942  "E_{T}(TT) (GeV);"
943  "TTF",
944  100, 0., 10.,
945  8, -.5, 7.5);
946  title = string("Trigger primitive Et (TP) vs #sumE_{T};"
947  "E_{T} (sum) (GeV);"
948  "E_{T} (TP) (") + tpUnit + string (")");
949 
950  meTpVsEtSum_ = book2D(ibooker, "h2TpVsEtSum",
951  title.c_str(),
952  100, 0., 10.,
953  100, 0., (tpInGeV_?10.:40.));
954 
955  title = string("Trigger primitive E_{T};"
956  "iEta;"
957  "iPhi;"
958  "E_{T} (TP) (") + tpUnit + string (")");
959  meTpMap_ = bookProfile2D(ibooker, "h2Tp",
960  title.c_str(),
961  57, -28.5, 28.5,
962  72, .5, 72.5);
963 
964  //SRF
965  meFullRoRu_ = book2D(ibooker, "h2FRORu", //"EcalFullReadoutSRFlagMap",
966  "Full Read-out readout unit;"
967  "iX - 40 / iEta / iX + 20;"
968  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
969  "Event count",
970  80, -39.5, 40.5,
971  72, .5, 72.5);
972 
973  meFullRoCnt_ = book1D(ibooker, "hFROCnt",
974  "Number of Full-readout-flagged readout units;"
975  "FRO RU count;Event count",
976  300, -.5, 299.5);
977 
978  meEbFullRoCnt_ = book1D(ibooker, "hEbFROCnt",
979  "Number of EB Full-readout-flagged readout units;"
980  "FRO RU count;Event count",
981  200, -.5, 199.5);
982 
983  meEeFullRoCnt_ = book1D(ibooker, "hEeFROCnt",
984  "Number of EE Full-readout-flagged readout units;"
985  "FRO RU count;Event count",
986  200, -.5, 199.5);
987 
988  meZs1Ru_ = book2D(ibooker, "h2Zs1Ru", //"EbZeroSupp1SRFlagMap",
989  "Readout unit with ZS-thr-1 flag;"
990  "iX - 40 / iEta / iX + 20;"
991  "iY0 / iPhi0 (iPhi = 1 at phi = 0 rad);"
992  "Event count",
993  80, -39.5, 40.5,
994  72, .5, 72.5);
995 
996  meForcedRu_ = book2D(ibooker, "h2ForcedRu", //"EcalReadoutUnitForcedBitMap",
997  "ECAL readout unit with forced bit of SR flag on;"
998  "iX - 40 / iEta / iX + 20;"
999  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1000  "Event count",
1001  80, -39.5, 40.5,
1002  72, .5, 72.5);
1003 
1004  meLiTtf_ = bookProfile2D(ibooker, "h2LiTtf", //"EcalLowInterestTriggerTowerFlagMap",
1005  "Low interest trigger tower flags;"
1006  "iEta;"
1007  "iPhi;"
1008  "Event count",
1009  57, -28.5, 28.5,
1010  72, .5, 72.5);
1011 
1012  meMiTtf_ = bookProfile2D(ibooker, "h2MiTtf", //"EcalMidInterestTriggerTowerFlagMap",
1013  "Mid interest trigger tower flags;"
1014  "iEta;"
1015  "iPhi;"
1016  "Event count",
1017  57, -28.5, 28.5,
1018  72, .5, 72.5);
1019 
1020  meHiTtf_ = bookProfile2D(ibooker, "h2HiTtf", //"EcalHighInterestTriggerTowerFlagMap",
1021  "High interest trigger tower flags;"
1022  "iEta;"
1023  "iPhi;"
1024  "Event count",
1025  57, -28.5, 28.5,
1026  72, .5, 72.5);
1027 
1028  meForcedTtf_ = book2D(ibooker, "h2ForcedTtf", //"EcalTtfForcedBitMap",
1029  "Trigger tower flags with forced bit set;"
1030  "iEta;"
1031  "iPhi;"
1032  "Event count",
1033  57, -28.5, 28.5,
1034  72, .5, 72.5);
1035 
1036 
1037  const float ebMinNoise = -1.;
1038  const float ebMaxNoise = 1.;
1039 
1040  const float eeMinNoise = -1.;
1041  const float eeMaxNoise = 1.;
1042 
1043  const float ebMinE = ebMinNoise;
1044  const float ebMaxE = ebMaxNoise;
1045 
1046  const float eeMinE = eeMinNoise;
1047  const float eeMaxE = eeMaxNoise;
1048 
1049  const int evtMax = 500;
1050 
1051  meEbRecE_ = book1D(ibooker, "hEbRecE",
1052  "Crystal reconstructed energy;E (GeV);Event count",
1053  100, ebMinE, ebMaxE);
1054 
1055  meEbEMean_ = bookProfile(ibooker, "hEbEMean",
1056  "EE <E_hit>;event #;<E_hit> (GeV)",
1057  evtMax, .5, evtMax + .5);
1058 
1059  meEbNoise_ = book1D(ibooker, "hEbNoise",
1060  "Crystal noise "
1061  "(rec E of crystal without deposited energy)"
1062  ";Rec E (GeV);Event count",
1063  100, ebMinNoise, ebMaxNoise);
1064 
1065  meEbLiZsFir_ = book1D(ibooker, "zsEbLiFIRemu",
1066  "Emulated ouput of ZS FIR filter for EB "
1067  "low interest crystals;"
1068  "ADC count*4;"
1069  "Event count",
1070  60, -30, 30);
1071 
1072  meEbHiZsFir_ = book1D(ibooker, "zsEbHiFIRemu",
1073  "Emulated ouput of ZS FIR filter for EB "
1074  "high interest crystals;"
1075  "ADC count*4;"
1076  "Event count",
1077  60, -30, 30);
1078 
1079  //TODO: Fill this histogram...
1080 // meEbIncompleteRUZsFir_ = book1D(ibooker, "zsEbIncompleteRUFIRemu",
1081 // "Emulated ouput of ZS FIR filter for EB "
1082 // "incomplete FRO-flagged RU;"
1083 // "ADC count*4;"
1084 // "Event count",
1085 // 60, -30, 30);
1086 
1087  meEbSimE_ = book1D(ibooker, "hEbSimE", "EB hit crystal simulated energy",
1088  100, ebMinE, ebMaxE);
1089 
1090  meEbRecEHitXtal_ = book1D(ibooker, "hEbRecEHitXtal",
1091  "EB rec energy of hit crystals",
1092  100, ebMinE, ebMaxE);
1093 
1094  meEbRecVsSimE_ = book2D(ibooker, "hEbRecVsSimE",
1095  "Crystal simulated vs reconstructed energy;"
1096  "Esim (GeV);Erec GeV);Event count",
1097  100, ebMinE, ebMaxE,
1098  100, ebMinE, ebMaxE);
1099 
1100  meEbNoZsRecVsSimE_ = book2D(ibooker, "hEbNoZsRecVsSimE",
1101  "Crystal no-zs simulated vs reconstructed "
1102  "energy;"
1103  "Esim (GeV);Erec GeV);Event count",
1104  100, ebMinE, ebMaxE,
1105  100, ebMinE, ebMaxE);
1106 
1107  meEeRecE_ = book1D(ibooker, "hEeRecE",
1108  "EE crystal reconstructed energy;E (GeV);"
1109  "Event count",
1110  100, eeMinE, eeMaxE);
1111 
1112  meEeEMean_ = bookProfile(ibooker, "hEeEMean",
1113  "<E_{EE hit}>;event;<E_{hit}> (GeV)",
1114  evtMax, .5, evtMax + .5);
1115 
1116 
1117  meEeNoise_ = book1D(ibooker, "hEeNoise",
1118  "EE crystal noise "
1119  "(rec E of crystal without deposited energy);"
1120  "E (GeV);Event count",
1121  200, eeMinNoise, eeMaxNoise);
1122 
1123  meEeLiZsFir_ = book1D(ibooker, "zsEeLiFIRemu",
1124  "Emulated ouput of ZS FIR filter for EE "
1125  "low interest crystals;"
1126  "ADC count*4;"
1127  "Event count",
1128  60, -30, 30);
1129 
1130  meEeHiZsFir_ = book1D(ibooker, "zsEeHiFIRemu",
1131  "Emulated ouput of ZS FIR filter for EE "
1132  "high interest crystals;"
1133  "ADC count*4;"
1134  "Event count",
1135  60, -30, 30);
1136 
1137  //TODO: Fill this histogram...
1138 // meEeIncompleteRUZsFir_ = book1D(ibooker, "zsEeIncompleteRUFIRemu",
1139 // "Emulated ouput of ZS FIR filter for EE "
1140 // "incomplete FRO-flagged RU;"
1141 // "ADC count*4;"
1142 // "Event count",
1143 // 60, -30, 30);
1144 
1145 
1146  meEeSimE_ = book1D(ibooker, "hEeSimE", "EE hit crystal simulated energy",
1147  100, eeMinE, eeMaxE);
1148 
1149  meEeRecEHitXtal_ = book1D(ibooker, "hEeRecEHitXtal",
1150  "EE rec energy of hit crystals",
1151  100, eeMinE, eeMaxE);
1152 
1153  meEeRecVsSimE_ = book2D(ibooker, "hEeRecVsSimE",
1154  "EE crystal simulated vs reconstructed energy;"
1155  "Esim (GeV);Erec GeV);Event count",
1156  100, eeMinE, eeMaxE,
1157  100, eeMinE, eeMaxE);
1158 
1159  meEeNoZsRecVsSimE_ = book2D(ibooker, "hEeNoZsRecVsSimE",
1160  "EE crystal no-zs simulated vs "
1161  "reconstructed "
1162  "energy;Esim (GeV);Erec GeV);Event count",
1163  100, eeMinE, eeMaxE,
1164  100, eeMinE, eeMaxE);
1165 
1166  meSRFlagsConsistency_ = book2D(ibooker, "hSRAlgoErrorMap",
1167  "TTFlags and SR Flags mismatch;"
1168  "iX - 40 / iEta / iX + 20;"
1169  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1170  "Event count",
1171  80, -39.5, 40.5,
1172  72, .5, 72.5);
1173 
1174  //Readout Units histos (interest/Ncrystals)
1175  meIncompleteFROMap_ = book2D(ibooker, "hIncompleteFROMap",
1176  "Incomplete full-readout-flagged readout units;"
1177  "iX - 40 / iEta / iX + 20;"
1178  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1179  "Event count",
1180  80, -39.5, 40.5,
1181  72, .5, 72.5);
1182 
1183  meIncompleteFROCnt_ = book1D(ibooker, "hIncompleteFROCnt",
1184  "Number of incomplete full-readout-flagged "
1185  "readout units;"
1186  "Number of RUs;Event count;",
1187  200, -.5, 199.5);
1188 
1190  = bookProfile2D(ibooker, "hIncompleteFRORateMap",
1191  "Incomplete full-readout-flagged readout units;"
1192  "iX - 40 / iEta / iX + 20;"
1193  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1194  "Incomplete error rate",
1195  80, -39.5, 40.5,
1196  72, .5, 72.5);
1197 
1198 
1199  meDroppedFROMap_ = book2D(ibooker, "hDroppedFROMap",
1200  "Dropped full-readout-flagged readout units;"
1201  "iX - 40 / iEta / iX + 20;"
1202  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1203  "Event count",
1204  80, -39.5, 40.5,
1205  72, .5, 72.5);
1206 
1207  meDroppedFROCnt_ = book1D(ibooker, "hDroppedFROCnt",
1208  "Number of dropped full-readout-flagged "
1209  "RU count;RU count;Event count",
1210  200, -.5, 199.5);
1211 
1212  meCompleteZSCnt_ = book1D(ibooker, "hCompleteZsCnt",
1213  "Number of zero-suppressed-flagged RU fully "
1214  "readout;"
1215  "RU count;Event count",
1216  200, -.5, 199.5);
1217 
1218  stringstream buf;
1219  buf << "Number of LI EB channels below the " << ebZsThr_/4. << " ADC count ZS threshold;"
1220  "Channel count;Event count",
1221  meEbZsErrCnt_ = book1D(ibooker, "hEbZsErrCnt",
1222  buf.str().c_str(),
1223  200, -.5, 199.5);
1224 
1225  buf.str("");
1226  buf << "Number of LI EE channels below the " << eeZsThr_/4. << " ADC count ZS theshold;"
1227  "Channel count;Event count",
1228  meEeZsErrCnt_ = book1D(ibooker, "hEeZsErrCnt",
1229  buf.str().c_str(),
1230  200, -.5, 199.5);
1231 
1232  meZsErrCnt_ = book1D(ibooker, "hZsErrCnt",
1233  "Number of LI channels below the ZS threshold;"
1234  "Channel count;Event count",
1235  200, -.5, 199.5);
1236 
1237  meEbZsErrType1Cnt_ = book1D(ibooker, "hEbZsErrType1Cnt",
1238  "Number of EB channels below the ZS "
1239  "threshold in a LI but fully readout RU;"
1240  "Channel count;Event count;",
1241  200, -.5, 199.5);
1242 
1243  meEeZsErrType1Cnt_ = book1D(ibooker, "hEeZsErrType1Cnt",
1244  "Number EE channels below the ZS threshold"
1245  " in a LI but fully readout RU;"
1246  "Channel count;Event count",
1247  200, -.5, 199.5);
1248 
1249  meZsErrType1Cnt_ = book1D(ibooker, "hZsErrType1Cnt",
1250  "Number of LI channels below the ZS threshold "
1251  "in a LI but fully readout RU;"
1252  "Channel count;Event count",
1253  200, -.5, 199.5);
1254 
1255 
1257  = bookProfile2D(ibooker, "hDroppedFRORateMap",
1258  "Dropped full-readout-flagged readout units"
1259  "iX - 40 / iEta / iX + 20;"
1260  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1261  "Dropping rate",
1262  80, -39.5, 40.5,
1263  72, .5, 72.5);
1264 
1265  meCompleteZSMap_ = book2D(ibooker, "hCompleteZSMap",
1266  "Complete zero-suppressed-flagged readout units;"
1267  "iX - 40 / iEta / iX + 20;"
1268  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1269  "Event count",
1270  80, -39.5, 40.5,
1271  72, .5, 72.5);
1272 
1274  = bookProfile2D(ibooker, "hCompleteZSRate",
1275  "Complete zero-suppressed-flagged readout units;"
1276  "iX - 40 / iEta / iX + 20;"
1277  "iY / iPhi (iPhi = 1 at phi = 0 rad);"
1278  "Completeness rate",
1279  80, -39.5, 40.5,
1280  72, .5, 72.5);
1281 
1282  //print list of available histograms (must be called after
1283  //the bookXX methods):
1285 
1286  //check the histList parameter:
1287  stringstream s;
1288  for(set<string>::iterator it = histList_.begin();
1289  it!=histList_.end();
1290  ++it){
1291  if(*it!=string("all")
1292  && availableHistList_.find(*it)==availableHistList_.end()){
1293  s << (s.str().size()==0?"":", ") << *it;
1294  }
1295  }
1296  if(s.str().size()!=0){
1297  LogWarning("Configuration")
1298  << "Parameter 'histList' contains some unknown histogram(s). "
1299  "Check spelling. Following name were not found: "
1300  << s.str();
1301  }
1302 }
1303 
1304 void
1306  edm::EventSetup const & es){
1307  int TTFlagCount[8];
1308  int LiTTFlagCount[nTtEta][nTtPhi];
1309  int MiTTFlagCount[nTtEta][nTtPhi];
1310  int HiTTFlagCount[nTtEta][nTtPhi];
1311  for(int iTTFlag(0);iTTFlag <8; iTTFlag++){
1312  TTFlagCount[iTTFlag]=0;
1313  }
1314  for(int iTtEta(0); iTtEta < nTtEta; iTtEta++){
1315  for(int iTtPhi(0); iTtPhi < nTtPhi; iTtPhi++){
1316  LiTTFlagCount[iTtEta][iTtPhi] = 0;
1317  MiTTFlagCount[iTtEta][iTtPhi] = 0;
1318  HiTTFlagCount[iTtEta][iTtPhi] = 0;
1319  }
1320  }
1321  int tpEtCount[100];
1322  for(int iEt(0);iEt < 100; iEt++){
1323  tpEtCount[iEt] = 0;
1324  }
1325 
1327  es.get<EcalTPGPhysicsConstRcd>().get( physHandle );
1328  const EcalTPGPhysicsConstMap & physMap = physHandle.product()->getMap() ;
1329 
1330  edm::ESHandle<EcalTPGLutGroup> lutGrpHandle;
1331  es.get<EcalTPGLutGroupRcd>().get( lutGrpHandle );
1332  const EcalTPGGroups::EcalTPGGroupsMap & lutGrpMap = lutGrpHandle.product()->getMap() ;
1333 
1334  edm::ESHandle<EcalTPGLutIdMap> lutMapHandle;
1335  es.get<EcalTPGLutIdMapRcd>().get( lutMapHandle );
1336  const EcalTPGLutIdMap::EcalTPGLutMap & lutMap = lutMapHandle.product()->getMap() ;
1337 
1338  EcalTPGPhysicsConstMapIterator ebItr(physMap.find(DetId(DetId::Ecal,EcalBarrel).rawId()));
1339  double lsb10bitsEB(ebItr == physMap.end() ? 0. : ebItr->second.EtSat / 1024.);
1340  EcalTPGPhysicsConstMapIterator eeItr(physMap.find(DetId(DetId::Ecal,EcalEndcap).rawId()));
1341  double lsb10bitsEE(eeItr == physMap.end() ? 0. : eeItr->second.EtSat / 1024.);
1342 
1343  for(EcalTrigPrimDigiCollection::const_iterator it = tps_->begin(); it != tps_->end(); ++it){
1344 
1345  double tpEt;
1346  if(tpInGeV_){
1347  EcalTrigTowerDetId const& towerId(it->id());
1348  unsigned int ADC = it->compressedEt();
1349 
1350  double lsb10bits(0.);
1351  if(towerId.subDet() == EcalBarrel) lsb10bits = lsb10bitsEB;
1352  else if(towerId.subDet() == EcalEndcap) lsb10bits = lsb10bitsEE;
1353 
1354  int tpg10bits = 0 ;
1355  EcalTPGGroups::EcalTPGGroupsMapItr itgrp = lutGrpMap.find(towerId.rawId()) ;
1356  uint32_t lutGrp = 999 ;
1357  if (itgrp != lutGrpMap.end()) lutGrp = itgrp->second ;
1358 
1359  EcalTPGLutIdMap::EcalTPGLutMapItr itLut = lutMap.find(lutGrp) ;
1360  if (itLut != lutMap.end()) {
1361  const unsigned int * lut = (itLut->second).getLut() ;
1362  for (unsigned int i=0 ; i<1024 ; i++)
1363  if (ADC == (0xff & lut[i])) {
1364  tpg10bits = i ;
1365  break ;
1366  }
1367  }
1368 
1369  tpEt = lsb10bits * tpg10bits;
1370  }
1371  else{
1372  tpEt = it->compressedEt();
1373  }
1374  int iEta = it->id().ieta();
1375  int iEta0 = iTtEta2cIndex(iEta);
1376  int iPhi = it->id().iphi();
1377  int iPhi0 = iTtPhi2cIndex(iPhi);
1378  double etSum = ttEtSums[iEta0][iPhi0];
1379 
1380  int iE = meTp_->getTProfile()->FindFixBin(tpEt);
1381  if ((iE >= 0) && (iE < 100)) {
1382  ++tpEtCount[iE];
1383  } else {
1384  // FindFixBin might return an overflow bin (outside tpEtCount range).
1385  // To prevent a memory overflow / segfault, these values are ignored.
1386  //std::cout << "EcalSelectiveReadoutValidation: Invalid iE value: " << iE << std::endl;
1387  }
1388 
1389  fill(meTpVsEtSum_, etSum, tpEt);
1390  ++TTFlagCount[it->ttFlag()];
1391  if((it->ttFlag() & 0x3) == 0){
1392  LiTTFlagCount[iEta0][iPhi0] += 1;
1393  }
1394  else if((it->ttFlag() & 0x3) == 1){
1395  MiTTFlagCount[iEta0][iPhi0] += 1;
1396  }
1397  else if((it->ttFlag() & 0x3) == 3){
1398  HiTTFlagCount[iEta0][iPhi0] += 1;
1399  }
1400  if((it->ttFlag() & 0x4)){
1401  fill(meForcedTtf_, iEta, iPhi);
1402  }
1403 
1404  fill(meTtfVsTp_, tpEt, it->ttFlag());
1405  fill(meTtfVsEtSum_, etSum, it->ttFlag());
1406  fill(meTpMap_, iEta, iPhi, tpEt, 1.);
1407  }
1408 
1409  for(int ittflag(0); ittflag < 8; ittflag++){
1410  fill(meTtf_, ittflag, TTFlagCount[ittflag]);
1411  }
1412  for(int iTtEta(0); iTtEta < nTtEta; iTtEta++){
1413  for(int iTtPhi(0); iTtPhi < nTtPhi; iTtPhi++){
1414  fill(meLiTtf_, cIndex2iTtEta(iTtEta), cIndex2iTtPhi(iTtPhi), LiTTFlagCount[iTtEta][iTtPhi]);
1415  fill(meMiTtf_, cIndex2iTtEta(iTtEta), cIndex2iTtPhi(iTtPhi), MiTTFlagCount[iTtEta][iTtPhi]);
1416  fill(meHiTtf_, cIndex2iTtEta(iTtEta), cIndex2iTtPhi(iTtPhi), HiTTFlagCount[iTtEta][iTtPhi]);
1417  }
1418  }
1419  if(tpInGeV_){
1420  for(int iE(0);iE<100; iE++){
1421  fill(meTp_, iE, tpEtCount[iE]);
1422  }
1423  } else{
1424  for(int iE(0); iE<40; iE++){
1425  fill(meTp_, iE, tpEtCount[iE]);
1426  }
1427  }
1428 }
1429 
1431  const EventSetup& es){
1432 
1433  anaDigiInit();
1434 
1435 
1436  //Complete RU, i.e. RU actually fully readout
1437  for(int iDcc = minDccId_; iDcc <= maxDccId_; ++iDcc){
1438  for(int iCh = 1; iCh < nDccRus_[iDcc-minDccId_]; ++iCh){
1439  isRuComplete_[iDcc-minDccId_][iCh-1]
1440  = (nPerRu_[iDcc-minDccId_][iCh-1]==getCrystalCount(iDcc, iCh));
1441  }
1442  }
1443 
1444 
1445  //Barrel
1446  for (unsigned int digis=0; digis<ebDigis_->size(); ++digis){
1447  EBDataFrame ebdf = (*ebDigis_)[digis];
1448  anaDigi(ebdf, *ebSrFlags_);
1449  }
1450 
1451  // Endcap
1452  for (unsigned int digis=0; digis<eeDigis_->size(); ++digis){
1453  EEDataFrame eedf = (*eeDigis_)[digis];
1454  anaDigi(eedf, *eeSrFlags_);
1455  }
1456 
1457  //histos
1458  for(unsigned iDcc0 = 0; iDcc0 < nDccs_; ++iDcc0){
1459  fill(meDccVol_, iDcc0+1, getDccEventSize(iDcc0, nPerDcc_[iDcc0])/kByte_);
1460  fill(meDccLiVol_, iDcc0+1,
1461  getDccSrDependentPayload(iDcc0, nLiRuPerDcc_[iDcc0],
1462  nLiPerDcc_[iDcc0])/kByte_);
1463  fill(meDccHiVol_, iDcc0+1,
1464  getDccSrDependentPayload(iDcc0, nHiRuPerDcc_[iDcc0],
1465  nHiPerDcc_[iDcc0])/kByte_);
1466  const FEDRawDataCollection& raw = *fedRaw_;
1467  fill(meDccVolFromData_, iDcc0+1,
1468  ((double)raw.FEDData(601+iDcc0).size())/kByte_);
1469  }
1470 
1471 
1472  //low interesest channels:
1473  double a = nEbLI_*getBytesPerCrystal()/kByte_; //getEbEventSize(nEbLI_)/kByte_;
1474  fill(meVolBLI_, a);
1475  double b = nEeLI_*getBytesPerCrystal()/kByte_; //getEeEventSize(nEeLI_)/kByte_;
1476  fill(meVolELI_, b);
1477  fill(meVolLI_, a+b);
1478 
1479  //high interest chanels:
1480  a = nEbHI_*getBytesPerCrystal()/kByte_; //getEbEventSize(nEbHI_)/kByte_;
1481  fill(meVolBHI_, a);
1482  b = nEeHI_*getBytesPerCrystal()/kByte_; //getEeEventSize(nEeHI_)/kByte_;
1483  fill(meVolEHI_, b);
1484  fill(meVolHI_, a+b);
1485 
1486  //any-interest channels:
1487  a = getEbEventSize(nEb_)/kByte_;
1488  fill(meVolB_, a);
1489  b = getEeEventSize(nEe_)/kByte_;
1490  fill(meVolE_, b);
1491  fill(meVol_, a+b);
1492 }
1493 
1494 
1495 template<class T, class U>
1497  const U& srFlagColl){
1498  const DetId& xtalId = frame.id();
1499  typedef typename U::key_type RuDetId;
1500  const RuDetId& ruId = readOutUnitOf(frame.id());
1501  typename U::const_iterator srf = srFlagColl.find(ruId);
1502 
1503  bool highInterest = false;
1504  int flag = 0;
1505 
1506  if(srf != srFlagColl.end()){
1507  flag = srf->value() & ~EcalSrFlag::SRF_FORCED_MASK;
1508 
1509  highInterest = (flag == EcalSrFlag::SRF_FULL);
1510 
1511  }
1512 
1513  bool barrel = (xtalId.subdetId()==EcalBarrel);
1514 
1515  pair<int,int> ch = dccCh(xtalId);
1516 
1517  if(barrel){
1518  ++nEb_;
1519  if(highInterest){
1520  ++nEbHI_;
1521  } else{//low interest
1522  ++nEbLI_;
1523  }
1524  int iEta0 = iEta2cIndex(static_cast<const EBDetId&>(xtalId).ieta());
1525  int iPhi0 = iPhi2cIndex(static_cast<const EBDetId&>(xtalId).iphi());
1526  if(!ebRuActive_[iEta0/ebTtEdge][iPhi0/ebTtEdge]){
1527  ++nRuPerDcc_[ch.first-minDccId_];
1528  if(highInterest){
1529  ++nHiRuPerDcc_[ch.first-minDccId_];
1530  } else{
1531  ++nLiRuPerDcc_[ch.first-minDccId_];
1532  }
1533 
1534  ebRuActive_[iEta0/ebTtEdge][iPhi0/ebTtEdge] = true;
1535  }
1536  } else{//endcap
1537  ++nEe_;
1538  if(highInterest){
1539  ++nEeHI_;
1540  } else{//low interest
1541  ++nEeLI_;
1542  }
1543  int iX0 = iXY2cIndex(static_cast<const EEDetId&>(frame.id()).ix());
1544  int iY0 = iXY2cIndex(static_cast<const EEDetId&>(frame.id()).iy());
1545  int iZ0 = static_cast<const EEDetId&>(frame.id()).zside()>0?1:0;
1546 
1547  if(!eeRuActive_[iZ0][iX0/scEdge][iY0/scEdge]){
1548  ++nRuPerDcc_[ch.first-minDccId_];
1549  if(highInterest){
1550  ++nHiRuPerDcc_[ch.first-minDccId_];
1551  } else{
1552  ++nLiRuPerDcc_[ch.first-minDccId_];
1553  }
1554 
1555  eeRuActive_[iZ0][iX0/scEdge][iY0/scEdge] = true;
1556  }
1557  }
1558 
1559  if(ch.second < 1 || ch.second > 68){
1560  throw cms::Exception("EcalSelectiveReadoutValidation")
1561  << "Error in DCC channel retrieval for crystal with detId "
1562  << xtalId.rawId() << "DCC channel out of allowed range [1..68]\n";
1563  }
1564  ++nPerDcc_[ch.first-minDccId_];
1565  ++nPerRu_[ch.first-minDccId_][ch.second-1];
1566  if(highInterest){
1567  ++nHiPerDcc_[ch.first-minDccId_];
1568  } else{//low interest channel
1569  ++nLiPerDcc_[ch.first-minDccId_];
1570  }
1571 }
1572 
1574  nEb_ = 0;
1575  nEe_ = 0;
1576  nEeLI_ = 0;
1577  nEeHI_ = 0;
1578  nEbLI_ = 0;
1579  nEbHI_ = 0;
1580  bzero(nPerDcc_, sizeof(nPerDcc_));
1581  bzero(nLiPerDcc_, sizeof(nLiPerDcc_));
1582  bzero(nHiPerDcc_, sizeof(nHiPerDcc_));
1583  bzero(nRuPerDcc_, sizeof(nRuPerDcc_));
1584  bzero(ebRuActive_, sizeof(ebRuActive_));
1585  bzero(eeRuActive_, sizeof(eeRuActive_));
1586  bzero(nPerRu_, sizeof(nPerRu_));
1587  bzero(nLiRuPerDcc_, sizeof(nLiRuPerDcc_));
1588  bzero(nHiRuPerDcc_, sizeof(nHiRuPerDcc_));
1589 }
1590 
1592  static std::atomic<bool> firstCall {true};
1593  bool expected = true;
1594  if( firstCall.compare_exchange_strong(expected, false) ) {
1595  stringstream buf;
1596  buf << "Weights:";
1597  for(unsigned i=0; i<weights_.size();++i){
1598  buf << "\t" << weights_[i];
1599  }
1600  edm::LogInfo("EcalSrValid") << buf.str() << "\n";
1601  firstCall = false;
1602  }
1603  double adc2GeV = 0.;
1604 
1605  if(typeid(EBDataFrame)==typeid(frame)){//barrel APD
1606  adc2GeV = .035;
1607  } else if(typeid(EEDataFrame)==typeid(frame)){//endcap VPT
1608  adc2GeV = 0.06;
1609  } else{
1610  assert(false);
1611  }
1612 
1613  double acc = 0;
1614 
1615  const int n = min(frame.size(), (int)weights_.size());
1616 
1617  double gainInv[] = {12., 1., 6., 12.};
1618 
1619  for(int i=0; i < n; ++i){
1620  acc += weights_[i]*frame[i].adc()*gainInv[frame[i].gainId()]*adc2GeV;
1621  }
1622  return acc;
1623 }
1624 
1626  return nRuPerDcc_[iDcc0];
1627 }
1628 
1629 pair<int,int> EcalSelectiveReadoutValidation::dccCh(const DetId& detId) const{
1630  if(detId.det()!=DetId::Ecal){
1631  throw cms::Exception("InvalidParameter")
1632  << "Wrong type of DetId passed to the "
1633  "EcalSelectiveReadoutValidation::dccCh(const DetId&). "
1634  "An ECAL DetId was expected.\n";
1635  }
1636 
1637  DetId xtalId;
1638  switch(detId.subdetId()){
1639  case EcalTriggerTower: //Trigger tower
1640  {
1641  const EcalTrigTowerDetId tt = detId;
1642  //pick up one crystal of the trigger tower: they are however all readout by
1643  //the same DCC channel in the barrel.
1644  //Arithmetic is easier on the "c" indices:
1645  const int iTtPhi0 = iTtPhi2cIndex(tt.iphi());
1646  const int iTtEta0 = iTtEta2cIndex(tt.ieta());
1647  const int oneXtalPhi0 = iTtPhi0 * 5;
1648  const int oneXtalEta0 = (iTtEta0 - nOneEeTtEta) * 5;
1649 
1650  xtalId = EBDetId(cIndex2iEta(oneXtalEta0),
1651  cIndex2iPhi(oneXtalPhi0));
1652  }
1653  break;
1654  case EcalEndcap:
1655  if(detId.rawId() & 0x8000){ //Supercrystal
1656  return elecMap_->getDCCandSC(EcalScDetId(detId));
1657  } else { //EE crystal
1658  xtalId = detId;
1659  }
1660  break;
1661  case EcalBarrel: //EB crystal
1662  xtalId = detId;
1663  break;
1664  default:
1665  throw cms::Exception("InvalidParameter")
1666  << "Wrong type of DetId passed to the method "
1667  "EcalSelectiveReadoutValidation::dccCh(const DetId&). "
1668  "A valid EcalTriggerTower, EcalBarrel or EcalEndcap DetId was expected. "
1669  "detid = " << xtalId.rawId() << ".\n";
1670  }
1671 
1672  const EcalElectronicsId& EcalElecId = elecMap_->getElectronicsId(xtalId);
1673 
1674  pair<int,int> result;
1675  result.first = EcalElecId.dccId();
1676 
1677  if(result.first < minDccId_ || result.second > maxDccId_){
1678  throw cms::Exception("OutOfRange")
1679  << "Got an invalid DCC ID, DCCID = " << result.first
1680  << " for DetId 0x" << hex << detId.rawId()
1681  << " and 0x" << xtalId.rawId() << dec << "\n";
1682  }
1683 
1684  result.second = EcalElecId.towerId();
1685 
1686  if(result.second < 1 || result.second > 68){
1687  throw cms::Exception("OutOfRange")
1688  << "Got an invalid DCC channel ID, DCC_CH = " << result.second
1689  << " for DetId 0x" << hex << detId.rawId()
1690  << " and 0x" << xtalId.rawId() << dec << "\n";
1691  }
1692 
1693  return result;
1694 }
1695 
1698  return triggerTowerMap_->towerOf(xtalId);
1699 }
1700 
1703  const EcalElectronicsId& EcalElecId = elecMap_->getElectronicsId(xtalId);
1704  int iDCC= EcalElecId.dccId();
1705  int iDccChan = EcalElecId.towerId();
1706  const bool ignoreSingle = true;
1707  const vector<EcalScDetId> id = elecMap_->getEcalScDetId(iDCC, iDccChan, ignoreSingle);
1708  return id.size()>0?id[0]:EcalScDetId();
1709 }
1710 
1711 void
1713  const EBDigiCollection& ebDigis,
1714  const EEDigiCollection& eeDigis){
1715  //ecal geometry:
1716  static const CaloSubdetectorGeometry* eeGeometry = 0;
1717  static const CaloSubdetectorGeometry* ebGeometry = 0;
1718  if(eeGeometry==0 || ebGeometry==0){
1719  edm::ESHandle<CaloGeometry> geoHandle;
1720  es.get<CaloGeometryRecord>().get(geoHandle);
1721  eeGeometry
1722  = (*geoHandle).getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
1723  ebGeometry
1724  = (*geoHandle).getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
1725  }
1726 
1727  //init etSum array:
1728  for(int iEta0 = 0; iEta0 < nTtEta; ++iEta0){
1729  for(int iPhi0 = 0; iPhi0 < nTtPhi; ++iPhi0){
1730  ttEtSums[iEta0][iPhi0] = 0.;
1731  }
1732  }
1733 
1735  it != ebDigis_->end(); ++it){
1736  const EBDataFrame& frame = *it;
1737  const EcalTrigTowerDetId& ttId = triggerTowerMap_->towerOf(frame.id());
1738 
1739  const int iTtEta0 = iTtEta2cIndex(ttId.ieta());
1740  const int iTtPhi0 = iTtPhi2cIndex(ttId.iphi());
1741  double theta = ebGeometry->getGeometry(frame.id())->getPosition().theta();
1742  double e = frame2EnergyForTp(frame);
1743  if((frame2EnergyForTp(frame,-1) < e) && (frame2EnergyForTp(frame, 1) < e)){
1744  ttEtSums[iTtEta0][iTtPhi0] += e*sin(theta);
1745  }
1746  }
1747 
1748  for(EEDigiCollection::const_iterator it = eeDigis.begin();
1749  it != eeDigis.end(); ++it){
1750  const EEDataFrame& frame = *it;
1751  const EcalTrigTowerDetId& ttId = triggerTowerMap_->towerOf(frame.id());
1752  const int iTtEta0 = iTtEta2cIndex(ttId.ieta());
1753  const int iTtPhi0 = iTtPhi2cIndex(ttId.iphi());
1754 
1755  double theta = eeGeometry->getGeometry(frame.id())->getPosition().theta();
1756  double e = frame2EnergyForTp(frame);
1757  if((frame2EnergyForTp(frame,-1) < e) && (frame2EnergyForTp(frame, 1) < e)){
1758  ttEtSums[iTtEta0][iTtPhi0] += e*sin(theta);
1759  }
1760  }
1761 
1762  //dealing with pseudo-TT in two inner EE eta-ring:
1763  int innerTTEtas[] = {0, 1, 54, 55};
1764  for(unsigned iRing = 0; iRing < sizeof(innerTTEtas)/sizeof(innerTTEtas[0]);
1765  ++iRing){
1766  int iTtEta0 = innerTTEtas[iRing];
1767  //this detector eta-section is divided in only 36 phi bins
1768  //For this eta regions,
1769  //current tower eta numbering scheme is inconsistent. For geometry
1770  //version 133:
1771  //- TT are numbered from 0 to 72 for 36 bins
1772  //- some TT have an even index, some an odd index
1773  //For geometry version 125, there are 72 phi bins.
1774  //The code below should handle both geometry definition.
1775  //If there are 72 input trigger primitives for each inner eta-ring,
1776  //then the average of the trigger primitive of the two pseudo-TT of
1777  //a pair (nEta, nEta+1) is taken as Et of both pseudo TTs.
1778  //If there are only 36 input TTs for each inner eta ring, then half
1779  //of the present primitive of a pseudo TT pair is used as Et of both
1780  //pseudo TTs.
1781 
1782  for(unsigned iTtPhi0 = 0; iTtPhi0 < nTtPhi-1; iTtPhi0 += 2){
1783  double et = .5*(ttEtSums[iTtEta0][iTtPhi0]
1784  +ttEtSums[iTtEta0][iTtPhi0+1]);
1785  //divides the TT into 2 phi bins in order to match with 72 phi-bins SRP
1786  //scheme or average the Et on the two pseudo TTs if the TT is already
1787  //divided into two trigger primitives.
1788  ttEtSums[iTtEta0][iTtPhi0] = et;
1789  ttEtSums[iTtEta0][iTtPhi0+1] = et;
1790  }
1791  }
1792 }
1793 
1794 template<class T>
1796  int offset) const{
1797  //we have to start by 0 in order to handle offset=-1
1798  //(however Fenix FIR has AFAK only 5 taps)
1799  double weights[] = {0., -1/3., -1/3., -1/3., 0., 1.};
1800 
1801  double adc2GeV = 0.;
1802  if(typeid(frame) == typeid(EBDataFrame)){
1803  adc2GeV = 0.035;
1804  } else if(typeid(frame) == typeid(EEDataFrame)){
1805  adc2GeV = 0.060;
1806  } else{ //T is an invalid type!
1807  //TODO: replace message by a cms exception
1808  throw cms::Exception("Severe Error")
1809  << __FILE__ << ":" << __LINE__ << ": "
1810  << "this is a bug. Please report it.\n";
1811  }
1812 
1813  double acc = 0;
1814 
1815  const int n = min<int>(frame.size(), sizeof(weights)/sizeof(weights[0]));
1816 
1817  double gainInv[] = {12., 1., 6., 12};
1818 
1819  for(int i=offset; i < n; ++i){
1820  int iframe = i + offset;
1821  if(iframe>=0 && iframe<frame.size()){
1822  acc += weights[i]*frame[iframe].adc()
1823  *gainInv[frame[iframe].gainId()]*adc2GeV;
1824  }
1825  }
1826  //cout << "\n";
1827  return acc;
1828 }
1829 
1831  if(!registerHist(name, "")) return 0; //this histo is disabled
1832  MonitorElement* result = ibook.bookFloat(name);
1833  if(result==0){
1834  throw cms::Exception("DQM")
1835  << "Failed to book integer DQM monitor element" << name;
1836  }
1837  return result;
1838 }
1839 
1840 
1842  if(!registerHist(name, title)) return 0; //this histo is disabled
1843  MonitorElement* result = ibook.book1D(name, title, nbins, xmin, xmax);
1844  if(result==0){
1845  throw cms::Exception("Histo")
1846  << "Failed to book histogram " << name;
1847  }
1848  return result;
1849 }
1850 
1851 MonitorElement* EcalSelectiveReadoutValidation::book2D(DQMStore::IBooker& ibook, const std::string& name, const std::string& title, int nxbins, double xmin, double xmax, int nybins, double ymin, double ymax){
1852  if(!registerHist(name, title)) return 0; //this histo is disabled
1853  MonitorElement* result = ibook.book2D(name, title, nxbins, xmin, xmax,
1854  nybins, ymin, ymax);
1855  if(result==0){
1856  throw cms::Exception("Histo")
1857  << "Failed to book histogram " << name;
1858  }
1859  return result;
1860 }
1861 
1863  if(!registerHist(name, title)) return 0; //this histo is disabled
1864  MonitorElement* result = ibook.bookProfile(name, title, nbins, xmin, xmax,
1865  0, 0, 0);
1866  if(result==0){
1867  throw cms::Exception("Histo")
1868  << "Failed to book histogram " << name;
1869  }
1870  return result;
1871 }
1872 
1873 MonitorElement* EcalSelectiveReadoutValidation::bookProfile2D(DQMStore::IBooker& ibook, const std::string& name, const std::string& title, int nbinx, double xmin, double xmax, int nbiny, double ymin, double ymax, const char* option){
1874  if(!registerHist(name, title)) return 0; //this histo is disabled
1876  = ibook.bookProfile2D(name,
1877  title,
1878  nbinx, xmin, xmax,
1879  nbiny, ymin, ymax,
1880  0, 0, 0,
1881  option);
1882  if(result==0){
1883  throw cms::Exception("Histo")
1884  << "Failed to book histogram " << name;
1885  }
1886  return result;
1887 }
1888 
1890  const std::string& title){
1891  availableHistList_.insert(pair<string, string>(name, title));
1892  return allHists_ || histList_.find(name)!=histList_.end();
1893 }
1894 
1896  ebRecHits_.read(event);
1897  eeRecHits_.read(event);
1898  ebDigis_.read(event);
1899  eeDigis_.read(event);
1900  ebNoZsDigis_.read(event);
1901  eeNoZsDigis_.read(event);
1902  ebSrFlags_.read(event);
1903  eeSrFlags_.read(event);
1904  ebComputedSrFlags_.read(event);
1905  eeComputedSrFlags_.read(event);
1906  ebSimHits_.read(event);
1907  eeSimHits_.read(event);
1908  tps_.read(event);
1909  fedRaw_.read(event);
1910 }
1911 
1913  LogInfo log("HistoList");
1914  log << "Avalailable histograms (DQM monitor elements): \n";
1915  for(map<string, string>::iterator it = availableHistList_.begin();
1916  it != availableHistList_.end();
1917  ++it){
1918  log << it->first << ": " << it->second << "\n";
1919  }
1920  log << "\nTo include an histogram add its name in the vstring parameter "
1921  "'histograms' of the EcalSelectiveReadoutValidation module\n";
1922 }
1923 
1924 double EcalSelectiveReadoutValidation::getEbEventSize(double nReadXtals) const{
1925  double ruHeaderPayload = 0.;
1926  const int firstEbDcc0 = nEeDccs/2;
1927  for(int iDcc0 = firstEbDcc0; iDcc0 < firstEbDcc0 + nEbDccs; ++iDcc0){
1928  ruHeaderPayload += getRuCount(iDcc0)*8.;
1929  }
1930 
1931  return getDccOverhead(EB)*nEbDccs + nReadXtals*getBytesPerCrystal()
1932  + ruHeaderPayload;
1933 }
1934 
1935 double EcalSelectiveReadoutValidation::getEeEventSize(double nReadXtals) const{
1936  double ruHeaderPayload = 0.;
1937  const unsigned firstEbDcc0 = nEeDccs/2;
1938  for(unsigned iDcc0 = 0; iDcc0 < nDccs_; ++iDcc0){
1939  //skip barrel:
1940  if(iDcc0== firstEbDcc0) iDcc0 += nEbDccs;
1941  ruHeaderPayload += getRuCount(iDcc0)*8.;
1942  }
1943  return getDccOverhead(EE)*nEeDccs + nReadXtals*getBytesPerCrystal()
1944  + ruHeaderPayload;
1945 }
1946 
1947 //This implementation assumes that int is coded on at least 28-bits,
1948 //which in pratice should be always true.
1949 int
1951  const std::vector<int>& firWeights,
1952  int firstFIRSample,
1953  bool* saturated){
1954  const int nFIRTaps = 6;
1955  //FIR filter weights:
1956  const vector<int>& w = firWeights;
1957 
1958  //accumulator used to compute weighted sum of samples
1959  int acc = 0;
1960  bool gain12saturated = false;
1961  const int gain12 = 0x01;
1962  const int lastFIRSample = firstFIRSample + nFIRTaps - 1;
1963  //LogDebug("DccFir") << "DCC FIR operation: ";
1964  int iWeight = 0;
1965  for(int iSample=firstFIRSample-1;
1966  iSample<lastFIRSample; ++iSample, ++iWeight){
1967  if(iSample>=0 && iSample < frame.size()){
1968  EcalMGPASample sample(frame[iSample]);
1969  if(sample.gainId()!=gain12) gain12saturated = true;
1970  LogTrace("DccFir") << (iSample>=firstFIRSample?"+":"") << sample.adc()
1971  << "*(" << w[iWeight] << ")";
1972  acc+=sample.adc()*w[iWeight];
1973  } else{
1974  edm::LogWarning("DccFir") << __FILE__ << ":" << __LINE__ <<
1975  ": Not enough samples in data frame or 'ecalDccZs1stSample' module "
1976  "parameter is not valid...";
1977  }
1978  }
1979  LogTrace("DccFir") << "\n";
1980  //discards the 8 LSBs
1981  //(shift operator cannot be used on negative numbers because
1982  // the result depends on compilator implementation)
1983  acc = (acc>=0)?(acc >> 8):-(-acc >> 8);
1984  //ZS passed if weighted sum acc above ZS threshold or if
1985  //one sample has a lower gain than gain 12 (that is gain 12 output
1986  //is saturated)
1987 
1988  LogTrace("DccFir") << "acc: " << acc << "\n"
1989  << "saturated: " << (gain12saturated?"yes":"no") << "\n";
1990 
1991  if(saturated){
1992  *saturated = gain12saturated;
1993  }
1994 
1995  return gain12saturated?numeric_limits<int>::max():acc;
1996 }
1997 
1998 std::vector<int>
2000  normalizedWeights){
2001  const int nFIRTaps = 6;
2002  vector<int> firWeights(nFIRTaps, 0); //default weight: 0;
2003  const static int maxWeight = 0xEFF; //weights coded on 11+1 signed bits
2004  for(unsigned i=0; i < min((size_t)nFIRTaps,normalizedWeights.size()); ++i){
2005  firWeights[i] = lround(normalizedWeights[i] * (1<<10));
2006  if(abs(firWeights[i])>maxWeight){//overflow
2007  firWeights[i] = firWeights[i]<0?-maxWeight:maxWeight;
2008  }
2009  }
2010  return firWeights;
2011 }
2012 
2013 void
2014 EcalSelectiveReadoutValidation::configFirWeights(const vector<double>& weightsForZsFIR){
2015  bool notNormalized = false;
2016  bool notInt = false;
2017  for(unsigned i=0; i < weightsForZsFIR.size(); ++i){
2018  if(weightsForZsFIR[i] > 1.) notNormalized = true;
2019  if((int)weightsForZsFIR[i]!=weightsForZsFIR[i]) notInt = true;
2020  }
2021  if(notInt && notNormalized){
2022  throw cms::Exception("InvalidParameter")
2023  << "weigtsForZsFIR paramater values are not valid: they "
2024  << "must either be integer and uses the hardware representation "
2025  << "of the weights or less or equal than 1 and used the normalized "
2026  << "representation.";
2027  }
2028  LogInfo log("DccFir");
2029  if(notNormalized){
2030  firWeights_ = vector<int>(weightsForZsFIR.size());
2031  for(unsigned i = 0; i< weightsForZsFIR.size(); ++i){
2032  firWeights_[i] = (int)weightsForZsFIR[i];
2033  }
2034  } else{
2035  firWeights_ = getFIRWeights(weightsForZsFIR);
2036  }
2037 
2038  log << "Input weights for FIR: ";
2039  for(unsigned i = 0; i < weightsForZsFIR.size(); ++i){
2040  log << weightsForZsFIR[i] << "\t";
2041  }
2042 
2043  double s2 = 0.;
2044  log << "\nActual FIR weights: ";
2045  for(unsigned i = 0; i < firWeights_.size(); ++i){
2046  log << firWeights_[i] << "\t";
2047  s2 += firWeights_[i]*firWeights_[i];
2048  }
2049 
2050  s2 = sqrt(s2);
2051  log << "\nNormalized FIR weights after hw representation rounding: ";
2052  for(unsigned i = 0; i < firWeights_.size(); ++i){
2053  log << firWeights_[i] / (double)(1<<10) << "\t";
2054  }
2055 
2056  log <<"\nFirst FIR sample: " << firstFIRSample_;
2057 }
2058 
2060  if(logSrpAlgoErrors_){
2062  if(!srpAlgoErrorLog_.good()){
2063  throw cms::Exception("Output")
2064  << "Failed to open the log file '"
2066  << "' for SRP algorithm result check.\n";
2067  }
2068  }
2069 
2072  if(!srApplicationErrorLog_.good()){
2073  throw cms::Exception("Output")
2074  << "Failed to open the log file '"
2076  << "' for Selective Readout decision application check.\n";
2077  }
2078  }
2079 }
2080 
2081 //Compares two SR flag sorted collections . Both collections
2082 //are sorted by their key (the detid) and following algorithm is based on
2083 //this feature.
2084 template<class T> //T must be either an EBSrFlagCollection or an EESrFlagCollection
2085 void EcalSelectiveReadoutValidation::compareSrfColl(const edm::Event& event, T& srfFromData, T& computedSrf){
2086  typedef typename T::const_iterator SrFlagCollectionConstIt;
2087  typedef typename T::key_type MyRuDetIdType;
2088  SrFlagCollectionConstIt itSrfFromData = srfFromData.begin();
2089  SrFlagCollectionConstIt itComputedSr = computedSrf.begin();
2090 
2091  while(itSrfFromData != srfFromData.end()
2092  || itComputedSr != computedSrf.end()){
2093  MyRuDetIdType inconsistentRu = 0;
2094  bool inconsistent = false;
2095  if(itComputedSr == computedSrf.end() ||
2096  (itSrfFromData != srfFromData.end()
2097  && itSrfFromData->id() < itComputedSr->id())){
2098  //computedSrf is missig a detid found in srfFromData
2099  pair<int, int> ch = dccCh(itSrfFromData->id());
2100  srpAlgoErrorLog_ << event.id() << ": " << itSrfFromData->id()
2101  << ", DCC " << ch.first << " ch " << ch.second
2102  << " found in data (SRF:" << itSrfFromData->flagName()
2103  << ") but not in the set of SRFs computed from the data TTF.\n";
2104  inconsistentRu = itSrfFromData->id();
2105  inconsistent = true;
2106  ++itSrfFromData;
2107  } else if(itSrfFromData==srfFromData.end() ||
2108  (itComputedSr != computedSrf.end()
2109  && itComputedSr->id() < itSrfFromData->id())){
2110  //ebSrFlags is missing a detid found in computedSrf
2111  pair<int, int> ch = dccCh(itComputedSr->id());
2112  if(logErrForDccs_[ch.first-minDccId_]){
2113  srpAlgoErrorLog_ << event.id() << ": " << itComputedSr->id()
2114  << ", DCC " << ch.first << " ch " << ch.second
2115  << " not found in data. Computed SRF: "
2116  << itComputedSr->flagName() << ".\n";
2117  inconsistentRu = itComputedSr->id();
2118  inconsistent = true;
2119  }
2120  ++itComputedSr;
2121  } else{
2122  //*itSrfFromData and *itComputedSr has same detid
2123  if(itComputedSr->value()!=itSrfFromData->value()){
2124  pair<int, int> ch = dccCh(itSrfFromData->id());
2125  srpAlgoErrorLog_ << event.id() << ", "
2126  << itSrfFromData->id()
2127  << ", DCC " << ch.first << " ch " << ch.second
2128  << ", SRF inconsistency: "
2129  << "from data: " << itSrfFromData->flagName()
2130  << ", computed from TTF: "
2131  << itComputedSr->flagName()
2132  << "\n";
2133  inconsistentRu = itComputedSr->id();
2134  inconsistent = true;
2135  }
2136  if(itComputedSr != computedSrf.end()) ++itComputedSr;
2137  if(itSrfFromData != srfFromData.end()) ++itSrfFromData;
2138  }
2139 
2140  if(inconsistent) fill(meSRFlagsConsistency_, ruGraphX(inconsistentRu),
2141  ruGraphY(inconsistentRu));
2142  }
2143 }
2144 
2145 
2147  return elecMap_->getDCCandSC(detId).first;
2148 }
2149 
2151  if(detId.ietaAbs()>17){
2152  throw cms::Exception("InvalidArgument")
2153  << "Argument of EcalSelectiveReadoutValidation::dccId(const EcalTrigTowerDetId&) "
2154  << "must be a barrel trigger tower Id\n";
2155  }
2156  return dccCh(detId).first;
2157 }
2158 
2159 
2160 
2162  logErrForDccs_ = vector<bool>(nDccs_, false);
2163 
2165  it != ebSrFlags_->end();
2166  ++it){
2167  int iDcc = dccId(it->id()) - minDccId_;
2168 
2169  logErrForDccs_.at(iDcc) = true;
2170  }
2171 
2173  it != eeSrFlags_->end();
2174  ++it){
2175  int iDcc = dccId(it->id()) - minDccId_;
2176 
2177  logErrForDccs_.at(iDcc) = true;
2178  }
2179 
2180  stringstream buf;
2181  buf << "List of DCCs found in the first processed event: ";
2182  bool first = true;
2183  for(unsigned iDcc = 0; iDcc < nDccs_; ++iDcc){
2184  if(logErrForDccs_[iDcc]){
2185  buf << (first?"":", ") << (iDcc + minDccId_);
2186  first = false;
2187  }
2188  }
2189  buf << "\nOnly DCCs from this list will be considered for error logging\n";
2190  srpAlgoErrorLog_ << buf.str();
2191  srApplicationErrorLog_<< buf.str();
2192  LogInfo("EcalSrValid") << buf.str();
2193 }
2194 
2195 
2196 template<class T>
2198  T& srfs){
2199  typedef typename T::const_iterator SrFlagCollectionConstIt;
2200  typedef typename T::key_type MyRuDetIdType;
2201 
2202  for(SrFlagCollectionConstIt itSrf = srfs.begin();
2203  itSrf != srfs.end(); ++itSrf){
2204  int flag = itSrf->value() & ~EcalSrFlag::SRF_FORCED_MASK;
2205  pair<int,int> ru = dccCh(itSrf->id());
2206 
2207  if(flag == EcalSrFlag::SRF_FULL){
2208  if(nPerRu_[ru.first-minDccId_][ru.second-1]==getCrystalCount(ru.first, ru.second)){ //no error
2209  fill(meIncompleteFRORateMap_, ruGraphX(itSrf->id()),
2210  ruGraphY(itSrf->id()), 0);
2212  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 0);
2213  } else if(nPerRu_[ru.first-minDccId_][ru.second-1]==0) {//tower dropped!
2215  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 0);
2217  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 1);
2218  fill(meDroppedFROMap_, ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 1);
2219  ++nDroppedFRO_;
2220  srApplicationErrorLog_ << event.id() << ": Flag of RU "
2221  << itSrf->id() << " (DCC " << ru.first
2222  << " ch " << ru.second << ") is 'Full readout' "
2223  << "while none of its channel was read out\n";
2224  } else{ //tower partially read out
2226  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 1);
2228  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 0);
2230  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 1);
2231  ++nIncompleteFRO_;
2232  srApplicationErrorLog_ << event.id() << ": Flag of RU"
2233  << itSrf->id() << " (DCC " << ru.first
2234  << " ch " << ru.second << ") is 'Full readout' "
2235  << "while only "
2236  << nPerRu_[ru.first-minDccId_][ru.second-1]
2237  << " / " << getCrystalCount(ru.first, ru.second)
2238  << " channels were read out.\n";
2239  }
2240  }
2241 
2242  if(flag == EcalSrFlag::SRF_ZS1 || flag == EcalSrFlag::SRF_ZS2){
2243  if(nPerRu_[ru.first-minDccId_][ru.second-1]
2244  ==getCrystalCount(ru.first, ru.second)){
2245  //ZS readout unit whose every channel was read
2246 
2247  fill(meCompleteZSMap_, ruGraphX(itSrf->id()), ruGraphY(itSrf->id()));
2249  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 1);
2250 
2251  ++nCompleteZS_;
2252  } else{
2254  ruGraphX(itSrf->id()), ruGraphY(itSrf->id()), 0);
2255  }
2256  }
2257  }
2258 }
2259 
2261  if(iDcc < minDccId_ || iDcc > maxDccId_){ //invalid DCC
2262  return 0;
2263  } else if (10 <= iDcc && iDcc <= 45) {//EB
2264  return 25;
2265  } else { //EE
2266  int iDccPhi;
2267  if(iDcc < 10) iDccPhi = iDcc;
2268  else iDccPhi = iDcc - 45;
2269  switch(iDccPhi*100+iDccCh){
2270  case 110:
2271  case 232:
2272  case 312:
2273  case 412:
2274  case 532:
2275  case 610:
2276  case 830:
2277  case 806:
2278  //inner partials at 12, 3, and 9 o'clock
2279  return 20;
2280  case 134:
2281  case 634:
2282  case 827:
2283  case 803:
2284  return 10;
2285  case 330:
2286  case 430:
2287  return 20;
2288  case 203:
2289  case 503:
2290  case 721:
2291  case 921:
2292  return 21;
2293  default:
2294  return 25;
2295  }
2296  }
2297 }
#define LogDebug(id)
void configFirWeights(const std::vector< double > &weightsForZsFIR)
void endRun(const edm::Run &r, const edm::EventSetup &c) override
double getDccSrDependentPayload(int iDcc0, double nReadRus, double nReadXtals) const
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
const EcalTPGGroupsMap & getMap() const
Definition: EcalTPGGroups.h:28
CollHandle< EBDigiCollection > ebNoZsDigis_
int value() const
Definition: EcalSrFlag.h:44
std::ofstream srApplicationErrorLog_
Output ascii file for unconsistency between Xtals and RU Flags.
int ix() const
Definition: EEDetId.h:76
double getDccOverhead(subdet_t subdet) const
static const int nEeX
EE crystal grid size along X.
key_type id() const
Definition: EBDataFrame.h:31
boost::transform_iterator< IterHelp, boost::counting_iterator< int > > const_iterator
int xtalGraphX(const EEDetId &id) const
double getDccEventSize(int iDcc0, double nReadXtals) const
int nIncompleteFRO_
Counter of FRO-flagged RU only partial data.
const double w
Definition: UKUtility.cc:23
void dqmBeginRun(edm::Run const &, edm::EventSetup const &) override
int gainId(sample_type sample)
get the gainId (2 bits)
MonitorElement * bookProfile(Args &&...args)
Definition: DQMStore.h:157
static const int SRF_ZS2
Definition: EcalSrFlag.h:21
static const int nEbPhi
number of crystals along Phi in EB
MonitorElement * bookFloat(DQMStore::IBooker &, const std::string &name)
std::map< uint32_t, uint32_t >::const_iterator EcalTPGGroupsMapItr
Definition: EcalTPGGroups.h:23
Ecal readout channel identification [32:20] Unused (so far) [19:13] DCC id [12:6] tower [5:3] strip [...
static const int SRF_FORCED_MASK
Definition: EcalSrFlag.h:29
static const double rad2deg
Conversion factor from radian to degree.
energiesEe_t eeEnergies[nEndcaps][nEeX][nEeY]
double energy() const
Definition: PCaloHit.h:29
int simHit
energy reconstructed from zero-suppressed digi
CollHandle< EcalTrigPrimDigiCollection > tps_
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
void analyzeEB(const edm::Event &event, const edm::EventSetup &es)
std::vector< EcalRecHit >::const_iterator const_iterator
Geom::Theta< T > theta() const
int simHit
energy reconstructed from zero-suppressed digi
int nEeZsErrors_
Counter of EE ZS errors (LI channel below ZS threshold)
CollHandle< FEDRawDataCollection > fedRaw_
int towerId() const
get the tower id
void analyzeEE(const edm::Event &event, const edm::EventSetup &es)
const_iterator begin() const
static const int nEeDccs
number of DCCs for EE
std::pair< int, int > getDCCandSC(EcalScDetId id) const
int gainId() const
get the gainId (2 bits)
unsigned ttId(DetId const &)
size_t size() const
Lenght of the data buffer in bytes.
Definition: FEDRawData.h:47
int size() const
Definition: EcalDataFrame.h:26
bool registerHist(const std::string &name, const std::string &title)
EcalTrigTowerDetId towerOf(const DetId &id) const
Get the tower id for this det id (or null if not known)
int ieta() const
get the tower ieta
void setTtEtSums(const edm::EventSetup &es, const EBDigiCollection &ebDigis, const EEDigiCollection &eeDigis)
int zside() const
get the z-side of the tower (1/-1)
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
void Fill(long long x)
static const int nEbEta
number of crystals along Eta in EB
uint32_t rawId() const
get the raw id
Definition: DetId.h:43
U second(std::pair< T, U > const &p)
double frame2EnergyForTp(const T &frame, int offset=0) const
EcalElectronicsId getElectronicsId(const DetId &id) const
Get the electronics id for this det id.
CollHandle< EESrFlagCollection > eeComputedSrFlags_
static const int nTtPhi
Number of Trigger Towers along Phi.
const EcalTPGLutMap & getMap() const
int xtalGraphY(const EEDetId &id) const
CollHandle< EEDigiCollection > eeDigis_
const FEDRawData & FEDData(int fedid) const
retrieve data for fed
void checkSrApplication(const edm::Event &event, T &srfs)
std::string outputFile_
Output file for histograms.
CollHandle< RecHitCollection > eeRecHits_
int ix() const
Definition: EcalScDetId.h:71
int ietaAbs() const
get the absolute value of the tower ieta
MonitorElement * book2D(DQMStore::IBooker &, const std::string &name, const std::string &title, int nxbins, double xmin, double xmax, int nybins, double ymin, double ymax)
MonitorElement * bookProfile2D(Args &&...args)
Definition: DQMStore.h:163
MonitorElement * book1D(DQMStore::IBooker &, const std::string &name, const std::string &title, int nbins, double xmin, double xmax)
void analyzeDataVolume(const edm::Event &e, const edm::EventSetup &es)
ConsumesCollector consumesCollector()
Use a ConsumesCollector to gather consumes information from helper functions.
T sqrt(T t)
Definition: SSEVec.h:18
void analyze(edm::Event const &e, edm::EventSetup const &c) override
Analyzes the event.
void anaDigi(const T &frame, const U &srFlagColl)
CollHandle< EBSrFlagCollection > ebSrFlags_
int zside() const
Definition: EEDetId.h:70
void updateL1aRate(const edm::Event &event)
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:115
static const int SRF_FULL
Definition: EcalSrFlag.h:24
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int nCompleteZS_
Counter of ZS-flagged RU fully read out.
float energy() const
Definition: EcalRecHit.h:68
CollHandle< std::vector< PCaloHit > > eeSimHits_
static int dccZsFIR(const EcalDataFrame &frame, const std::vector< int > &firWeights, int firstFIRSample, bool *saturated=0)
int iy() const
Definition: EEDetId.h:82
int nEeFROCnt_
Counter of EE FRO-flagged RUs.
int ieta() const
get the crystal ieta
Definition: EBDetId.h:51
unsigned int id() const
Definition: PCaloHit.h:43
T min(T a, T b)
Definition: MathUtil.h:58
void bookHistograms(DQMStore::IBooker &i, edm::Run const &, edm::EventSetup const &) override
MonitorElement * bookProfile(DQMStore::IBooker &, const std::string &name, const std::string &title, int nbins, double xmin, double xmax)
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
int dccId() const
get the DCC (Ecal Local DCC value not global one) id
unsigned towerId(DetId const &)
#define LogTrace(id)
int iy() const
Definition: EcalScDetId.h:77
void read(const edm::Event &event)
Definition: CollHandle.h:62
static const int SRF_ZS1
Definition: EcalSrFlag.h:18
key_type id() const
Definition: EEDataFrame.h:28
std::vector< EcalScDetId > getEcalScDetId(int DCCid, int DCC_Channel, bool ignoreSingleCrystal=true) const
CollHandle< EESrFlagCollection > eeSrFlags_
static const int nDccRus_[nDccs_]
number of RUs for each DCC
static const int nOneEeTtEta
Number of Trigger Towers in an endcap along Eta.
int nEbFROCnt_
Counter of EB FRO-flagged RUs.
static const int nEbDccs
number of DCCs for EB
CollHandle< std::vector< PCaloHit > > ebSimHits_
int nEbZsErrors_
Counter of EB ZS errors (LI channel below ZS threshold)
static const int scEdge
Number of crystals along a supercrystal edge.
double recE
energy reconstructed from unsuppressed digi
static const double tmax[3]
Definition: DetId.h:18
static std::vector< int > getFIRWeights(const std::vector< double > &normalizedWeights)
int iphi() const
get the tower iphi
static const int ebTtEdge
Number of crystals along an EB TT.
EcalSelectiveReadoutValidation(const edm::ParameterSet &ps)
Constructor.
const std::map< uint32_t, Item > & getMap() const
DetId id() const
get the id
Definition: EcalRecHit.h:77
double getEeEventSize(double nReadXtals) const
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:133
static const int kByte_
number of bytes in 1 kByte:
const EcalScDetId & id() const
Definition: EESrFlag.h:37
bool ebRuActive_[nEbEta/ebTtEdge][nEbPhi/ebTtEdge]
const T & get() const
Definition: EventSetup.h:56
double b
Definition: hdecay.h:120
int zside() const
Definition: EcalScDetId.h:65
et
define resolution functions of each parameter
EcalTrigTowerDetId readOutUnitOf(const EBDetId &xtalId) const
CollHandle< EEDigiCollection > eeNoZsDigis_
double recE
energy reconstructed from unsuppressed digi
const EcalTrigTowerDetId & id() const
Definition: EBSrFlag.h:36
int dccId(const EcalScDetId &detId) const
std::map< std::string, std::string > availableHistList_
const_iterator end() const
T eta() const
Definition: PV3DBase.h:76
std::map< uint32_t, EcalTPGLut > EcalTPGLutMap
int ruGraphX(const EcalScDetId &id) const
HLT enums.
static const int nTtEta
Number of Trigger Towers along Eta.
double a
Definition: hdecay.h:121
TProfile * getTProfile(void) const
int nDroppedFRO_
Counter of FRO-flagged RU dropped from data.
void setToken(edm::ConsumesCollector &collector)
Definition: CollHandle.h:49
void fill(MonitorElement *me, float x)
std::map< uint32_t, EcalTPGPhysicsConst::Item > EcalTPGPhysicsConstMap
MonitorElement * bookFloat(Args &&...args)
Definition: DQMStore.h:109
int ruGraphY(const EcalScDetId &id) const
void compareSrfColl(const edm::Event &event, T &srfFromData, T &computedSrf)
std::pair< int, int > dccCh(const DetId &xtalId) const
double frame2Energy(const EcalDataFrame &frame) const
std::ofstream srpAlgoErrorLog_
Output ascii file for unconsistency on SR flags.
static const int nEeY
EE crystal grid size along Y.
CollHandle< EBDigiCollection > ebDigis_
std::map< uint32_t, EcalTPGPhysicsConst::Item >::const_iterator EcalTPGPhysicsConstMapIterator
double getEbEventSize(double nReadXtals) const
Detector det() const
get the detector field from this detid
Definition: DetId.h:35
long double T
MonitorElement * bookProfile2D(DQMStore::IBooker &, const std::string &name, const std::string &title, int nbinx, double xmin, double xmax, int nbiny, double ymin, double ymax, const char *option="")
static const unsigned nDccs_
Total number of DCCs.
T const * product() const
Definition: ESHandle.h:86
static const int nEndcaps
number of endcaps
std::map< uint32_t, uint32_t > EcalTPGGroupsMap
Definition: EcalTPGGroups.h:22
std::map< uint32_t, EcalTPGLut >::const_iterator EcalTPGLutMapItr
const EcalTrigTowerConstituentsMap * triggerTowerMap_
void analyzeTP(const edm::Event &event, const edm::EventSetup &es)
const EcalElectronicsMapping * elecMap_
CollHandle< RecHitCollection > ebRecHits_
Definition: event.py:1
Definition: Run.h:42
int adc() const
get the ADC sample (12 bits)
bool eeRuActive_[nEndcaps][nEeX/scEdge][nEeY/scEdge]
CollHandle< EBSrFlagCollection > ebComputedSrFlags_