ATLAS Offline Software
Loading...
Searching...
No Matches
TileDigitsMonTool.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3*/
4
5// ********************************************************************
6//
7// NAME: TileDigitsMonTool.cxx
8// PACKAGE:
9//
10// AUTHOR: Alexander Solodkov
11// Luca.Fiorini@cern.ch
12//
13// ********************************************************************
14
15#include "TileDigitsMonTool.h"
16
23
24#include "TH1S.h"
25#include "TH1F.h"
26#include "TH2F.h"
27#include "TH1D.h"
28#include "TProfile.h"
29#include "TVirtualPad.h"
30#include "TCanvas.h"
31#include "TStyle.h"
32#include "TLine.h"
33#include "TPaveLabel.h"
34#include "TObjArray.h"
35#include "TSpectrum.h"
36#include "TText.h"
37#include "TExec.h"
38#include "TColor.h"
39#include "TPaletteAxis.h"
40#include "TList.h"
41#include "TTree.h"
42
43#include <iostream>
44#include <sstream>
45#include <iomanip>
46#include <cmath>
47#include <algorithm>
48
49
50/*---------------------------------------------------------*/
51TileDigitsMonTool::TileDigitsMonTool(const std::string & type, const std::string & name, const IInterface* parent)
52 : TilePaterMonTool(type, name, parent)
53 , m_tileToolNoiseSample("TileCondToolNoiseSample")
54 , m_cispar(0)
56 , m_allHistsFilled(false)
57 //, hp(-1)
58 //, hb(-1)
59 , m_tileInfo(0)
60 , m_is12bit(false)
61 , m_shiftnbins(0)
62/*---------------------------------------------------------*/
63{
64 declareInterface<IMonitorToolBase>(this);
65 declareInterface<ITileStuckBitsProbsTool>(this);
66
67 declareProperty("bookAllDrawers", m_bookAll = false);
68 declareProperty("book2D", m_book2D = false);
69
70 // run type 1 - phys, 2 - las, 4 - ped, 8 - cis, 9- mono
71 declareProperty("runType", m_runType = 8);
72 declareProperty("bigain", m_bigain = true);
73 declareProperty("NSamples", m_nSamples = 7);
74 declareProperty("TileRawChannelContainerDSP", m_contNameDSP = "TileRawChannelCnt");
75 declareProperty("TileDigitsContainer", m_digitsContainerName = "TileDigitsCnt");
76 declareProperty("FillPedestalDifference", m_fillPedestalDifference = true);
77 declareProperty("TileInfoName", m_infoName = "TileInfo");
78 declareProperty("TileDQstatus", m_DQstatusKey = "TileDQstatus");
79 declareProperty("ZeroLimitHG", m_zeroLimitHG = 0);
80 declareProperty("SaturationLimitHG", m_saturationLimitHG = 1023);
81
82 m_path = "/Tile/Digits"; //ROOT file directory
83}
84
85/*---------------------------------------------------------*/
87/*---------------------------------------------------------*/
88{
89}
90
91/*---------------------------------------------------------*/
93/*---------------------------------------------------------*/
94{
95 ATH_MSG_INFO("in initialize()");
96
98
100 m_allHistsFilled = false;
101 m_data = std::make_unique<Data>();
102 memset(m_data->m_sumPed1, 0, sizeof(m_data->m_sumPed1));
103 memset(m_data->m_sumPed2, 0, sizeof(m_data->m_sumPed2));
104 memset(m_data->m_sumRms1, 0, sizeof(m_data->m_sumRms1));
105 memset(m_data->m_sumRms2, 0, sizeof(m_data->m_sumRms2));
106 memset(m_data->m_meanAmp, 0, sizeof(m_data->m_meanAmp));
107 memset(m_data->m_meanAmp_ij, 0, sizeof(m_data->m_meanAmp_ij));
108 memset(m_data->m_nEvents_i, 0, sizeof(m_data->m_nEvents_i));
109 memset(m_data->m_nEvents_ij, 0, sizeof(m_data->m_nEvents_ij));
110 memset(m_data->m_stuck_probs, 0, sizeof(m_data->m_stuck_probs));
111
112 //For test stuck_bits_maker
113 //hp = 1;
114 //hb = -1;
115
116 // TileInfo
117 CHECK( detStore()->retrieve(m_tileInfo, m_infoName) );
118 m_i_ADCmax = m_tileInfo->ADCmax();
121
122 m_shiftnbins = 796;
123 if (m_i_ADCmax == 4095){
124 m_is12bit = true;
126 }
127
128 CHECK( m_DQstatusKey.initialize() );
129
131}
132
133/*---------------------------------------------------------*/
135/*---------------------------------------------------------*/
136{
137
138 if (msgLvl(MSG::DEBUG)) {
139 msg(MSG::DEBUG) << "in bookHists()" << endmsg;
140 msg(MSG::DEBUG) << "Using base path " << m_path << endmsg;
141 }
142
143 if (m_bookAll) {
144 for (int ros = 1; ros < 5; ++ros) {
145 for (int drawer = 0; drawer < 64; ++drawer) {
146 bookHists(ros, drawer);
147 }
148 }
149 }
150
151 m_data->m_outInHighGain.push_back(book1D("Summary", "OverFlowInHighGain", "Injected charge with overflow in high gain", 10010, -0.5, 1000.5));
152
153 //SetBookStatus(true);
154
155 return StatusCode::SUCCESS;
156}
157
158/*---------------------------------------------------------*/
159void TileDigitsMonTool::bookHists(int ros, int drawer)
160/*---------------------------------------------------------*/
161{
162 const char *part[5] = { "AUX", "LBA", "LBC", "EBA", "EBC" };
163 const char *gain[6] = { "_lo", "_hi", "", " low gain", " high gain", "" };
164
174
175 MsgStream log(msgSvc(), name());
176
177 std::ostringstream sStr;
178 sStr << part[ros] << std::setfill('0') << std::setw(2) << drawer + 1 << std::setfill(' ');
179 std::string moduleName = sStr.str();
180 std::string subDir = moduleName;
181 std::string histName, histTitle;
182
183 ATH_MSG_DEBUG("in bookHists() for module " << moduleName);
184
187 int mingain = (m_bigain) ? 0 : 2;
188 int maxgain = (m_bigain) ? 2 : 3;
189
190 for (int gn = mingain; gn < maxgain; ++gn) {
191 //Book DMU header errors histograms
192 int adc = gn & 1;
193
194 sStr.str("");
195 sStr << moduleName << "_DMU_head_" << gain[gn] << "_err";
196 histName = sStr.str();
197
198 sStr.str("");
199 sStr << moduleName << " DMU Header " << gain[3 + gn] << " errors";
200 histTitle = sStr.str();
201
202 m_data->m_hist2[ros][drawer][adc].push_back(book2F(subDir, histName, histTitle, 16, 0., 16., 8, -0.5, 7.5));
203 }
204
205 for (int ch = 0; ch < 48; ++ch) {
206
207 // int pmt=m_cabling->channel2hole(ros,ch); // 1:48, negative for non-connected channels
208 // int pmt0 = ((pmt > 0)? pmt-1 : pmt+1); // 0:47, negaitve for non-connected channels
209 sStr.str("");
210 sStr << std::setfill('0') << std::setw(2) << ch << std::setfill(' ');
211 std::string sCh = sStr.str(); // 00:47, always positive
212
213 for (int gn = mingain; gn < maxgain; ++gn) {
214
215 int adc = gn & 1;
216
217 // first sample
218 sStr.str("");
219 sStr << moduleName << "_ch_" << sCh << gain[gn] << "_sample0";
220 histName = sStr.str();
221
222 sStr.str("");
223 sStr << moduleName << " CH " << ch << gain[3 + gn] << " sample[0]";
224 histTitle = sStr.str();
225
226 m_data->m_hist1[ros][drawer][ch][adc].push_back(book1S(subDir, histName, histTitle, 151, -0.5, 150.5));
227
228 // RMS over samples in event
229 sStr.str("");
230 sStr << moduleName << "_ch_" << sCh << gain[gn] << "_meanRMS";
231 histName = sStr.str();
232
233 sStr.str("");
234 sStr << moduleName << " CH " << ch << gain[3 + gn] << " mean RMS in event";
235 histTitle = sStr.str();
236
237 m_data->m_hist1[ros][drawer][ch][adc].push_back(book1S(subDir, histName, histTitle, 101, -0.05, 10.05));
238
239 //all samples to find stuck bits
240 sStr.str("");
241 sStr << moduleName << "_ch_" << sCh << gain[gn] << "_samples";
242 histName = sStr.str();
243
244 sStr.str("");
245 sStr << moduleName << " CH " << ch << gain[3 + gn] << " all samples";
246 histTitle = sStr.str();
247
248 m_data->m_hist1[ros][drawer][ch][adc].push_back(book1S(subDir, histName, histTitle, m_i_ADCmax + 1, -0.5, m_i_ADCmax + 0.5));
249
250 //eventnumber % 32 for n/32 gain switch failure
251 sStr.str("");
252 sStr << moduleName << "_ch_" << sCh << gain[gn] << "_evtn_mod32";
253 histName = sStr.str();
254
255 sStr.str("");
256 sStr << moduleName << " CH " << ch << gain[3 + gn] << " event number % 32";
257 histTitle = sStr.str();
258
259 m_data->m_histC[ros][drawer][ch][adc].push_back(book1C(subDir, histName, histTitle, 32, -0.5, 31.5));
260
261 // average profile for a given channel/gain
262 sStr.str("");
263 sStr << moduleName << "_ch_" << sCh << gain[gn] << "_prof";
264 histName = sStr.str();
265
266 sStr.str("");
267 sStr << moduleName << " CH " << ch << gain[3 + gn] << " average profile";
268 histTitle = sStr.str();
269
270 m_data->m_histP[ros][drawer][ch][adc].push_back(bookProfile(subDir, histName, histTitle, m_nSamples, 0.0, m_nSamples * 1.0, -0.5, m_i_ADCmax + 0.5));
271
272 // shifted samples histograms
273
274 double shiftxmin = 0.;
275 double shiftxmax = m_shiftnbins;
276 int shiftnbins = m_shiftnbins;
277
278 sStr.str("");
279 sStr << moduleName << "_ch_" << sCh << gain[gn] << "_samples_shifted";
280 histName = sStr.str();
281
282 sStr.str("");
283 sStr << moduleName << " CH " << ch << gain[3 + gn] << " all samples shifted";
284 histTitle = sStr.str();
285
286 m_data->m_shifted_hist[ros][drawer][ch][adc] = book1S(subDir, histName, histTitle, shiftnbins, shiftxmin, shiftxmax);
287
288 if (ch < 16) { // use first 16 channels to store BCID/CRC errors (one per TileDMU)
289
290 sStr.str("");
291 sStr << std::setfill('0') << std::setw(2) << ch << std::setfill(' ');
292 std::string sDmu = sStr.str();
293
294 sStr.str("");
295 sStr << moduleName << "_dmu_" << sDmu << gain[gn] << "_BCID";
296 histName = sStr.str();
297
298 sStr.str("");
299 sStr << moduleName << " DMU0 " << ch << gain[3 + gn] << " BCID";
300 histTitle = sStr.str();
301
302 m_data->m_hist1[ros][drawer][ch][adc].push_back(book1S(subDir, histName, histTitle, 100, 2050.5, 2150.5));
303
304 sStr.str("");
305 sStr << moduleName << "_dmu_" << sDmu << gain[gn] << "_BCID_err";
306 histName = sStr.str();
307
308 sStr.str("");
309 sStr << moduleName << " DMU0 " << ch << gain[3 + gn] << " BCID errors";
310 histTitle = sStr.str();
311
312 m_data->m_hist_DMUerr[ros][drawer][ch][adc].push_back(book1I(subDir, histName, histTitle, 3, -0.5, 2.5));
313
314 if (adc) continue; // don't book CRC for high gain
315
316 sStr.str("");
317 sStr << moduleName << "_dmu_" << sDmu << "_CRC_err";
318 histName = sStr.str();
319
320 sStr.str("");
321 sStr << moduleName << " DMU0 " << ch << " CRC errors";
322 histTitle = sStr.str();
323
324 m_data->m_hist_DMUerr[ros][drawer][ch][0].push_back(book1I(subDir, histName, histTitle, 5, -0.5, 4.5));
325 }
326 }
327 }
328
329// and now global histograms for all channels in a drawer
330
331 // global CRC
332 sStr.str("");
333 sStr << moduleName << "_CRC_global";
334 histName = sStr.str();
335
336 sStr.str("");
337 sStr << moduleName << " global CRC errors";
338 histTitle = sStr.str();
339
340 m_data->m_hist0[ros][drawer].push_back(book1S("Summary", histName, histTitle, 3, -0.5, 2.5));
341
342}
343
344/*---------------------------------------------------------*/
346/*---------------------------------------------------------*/
347{
348// static counter = 0;
349// bool extra = ((counter%m_updateCounter) == 0);
350// ++counter;
351
352 ATH_MSG_DEBUG("in fillHists()");
353
354 const TileDQstatus* dqStatus = SG::makeHandle (m_DQstatusKey).get();
355
356 // array of 16 CIS parameters
357 m_cispar = dqStatus->cispar();
359
360 // std::cout << "Calib Mode=" << dqStatus->calibMode() << "\n";
361
362 const TileDigitsContainer* digitsContainer;
363 CHECK(evtStore()->retrieve(digitsContainer, m_digitsContainerName));
364
365 for (const TileDigitsCollection* digitsCollection : *digitsContainer) {
366
367 if (digitsCollection->empty()) continue;
368
369 HWIdentifier adc_id = digitsCollection->front()->adc_HWID();
370 int ros = m_tileHWID->ros(adc_id);
371 int drawer = m_tileHWID->drawer(adc_id);
372
373 if (m_data->m_hist1[ros][drawer][0][0].size() == 0) {
374 m_bigain = true;
375 m_nSamples = digitsCollection->front()->NtimeSamples();
376 if (!m_bookAll) bookHists(ros, drawer);
377 }
378
379 uint32_t status = digitsCollection->getFragStatus();
380 if (status != TileFragStatus::ALL_OK) {
381 float bin = 99.;
383 bin = 6.;
384 } else if (status & (TileFragStatus::NO_FRAG | TileFragStatus::NO_ROB)) {
385 bin = 7.;
386 }
387 for (int dmu = 0; dmu < 16; dmu++) {
388 m_data->m_corrup[ros][drawer][0][dmu] = true;
389 m_data->m_corrup[ros][drawer][1][dmu] = true;
390 m_data->m_hist2[ros][drawer][0][0]->Fill(dmu + 0.5, bin, 1.);
391 m_data->m_hist2[ros][drawer][1][0]->Fill(dmu + 0.5, bin, 1.);
392 }
393 continue;
394 }
395
396 std::vector<uint32_t> headerVec = digitsCollection->getFragChipHeaderWords();
397
398 int headsize = headerVec.size();
399 if (headsize > 16) {
400 ATH_MSG_ERROR("length of BCID vector " << headsize << " - greater than 16 !");
401 headsize = 16;
402 }
403
404 for (int dmu = 0; dmu < headsize; dmu++) {
405 m_data->m_corrup[ros][drawer][0][dmu] = DMUheaderCheck(&headerVec, ros, drawer, 0, dmu); //tests and fills
406 m_data->m_corrup[ros][drawer][1][dmu] = DMUheaderCheck(&headerVec, ros, drawer, 1, dmu);
407 }
408
409 for (int dmu = headsize; dmu < 16; dmu++) {
410 m_data->m_corrup[ros][drawer][0][dmu] = false;
411 m_data->m_corrup[ros][drawer][1][dmu] = false;
412 }
413 }
414
415 for (const TileDigitsCollection* digitsCollection : *digitsContainer) {
416
417 if (digitsCollection->empty()) continue;
418
419 HWIdentifier adc_id = digitsCollection->front()->adc_HWID();
420 int ros = m_tileHWID->ros(adc_id);
421 int drawer = m_tileHWID->drawer(adc_id);
422
423 if (m_data->m_hist1[ros][drawer][0][0].size() == 0) {
424 //m_bigain = (dqStatus->calibMode() == 1); // true if bi-gain run
425 // For the time being, we fill both high and low gain plots, as it was requiered by Tomas
426 m_bigain = true;
427 m_nSamples = digitsCollection->front()->NtimeSamples(); // number of samples
428 // bookHists(ros,drawer);
429 }
430
431 std::vector<uint32_t> headerVec = digitsCollection->getFragChipHeaderWords();
432
433 int headsize = headerVec.size();
434 //if (headsize > 16 ) {
435 // log << MSG::ERROR << "legth of BCID vector " << headsize << " - greater than 16 !" << endreq;
436 // headsize = 16;
437 //}
438
439 //coverity[STACK_USE]
440 double mean_tmp[48][2][16];
441 memset(mean_tmp, 0, sizeof(mean_tmp));
442
443 // Below, the charge conversion for 12 bit is just an approximation;
444 // 4095. can be changed later to gain precision if needed.
445 double charge = 0;
446 if (m_is12bit){ // if 12 bit ADCs
447 charge = m_cispar[6] * m_cispar[7] * (2. * 4.096 / 4095.);
448 }
449 else{
450 charge = m_cispar[6] * m_cispar[7] * (2. * 4.096 / 1023.);
451 }
452
453 for (const TileDigits* tileDigits : *digitsCollection) {
454
455 adc_id = tileDigits->adc_HWID();
456 int chan = m_tileHWID->channel(adc_id);
457 int gain = (m_bigain) ? m_tileHWID->adc(adc_id) : 0; // ignore gain in monogain run
458
459 m_data->m_histC[ros][drawer][chan][gain][0]->Fill(m_nEventsTileMon % 32, 1.0);
460
461 std::vector<float> vdigits = tileDigits->samples();
462
463 double meansamp = 0.0;
464 double rmssamp = 0.0;
465 unsigned int dsize = vdigits.size();
466 if (dsize > 16) {
467 ATH_MSG_ERROR("length of digits vector " << dsize << " - greater than 16 !");
468 dsize = 16;
469 }
470 //std::cout << "m_runType= "<< m_runType << "\n";
471 for (unsigned int i = 0; i < dsize; ++i) {
472 double dig = vdigits[i];
473 meansamp += dig;
474 rmssamp += dig * dig;
475 mean_tmp[chan][gain][i] = dig;
476
477 if (!m_data->m_corrup[ros][drawer][gain][chan / 3]) {
478 if (m_runType != CisRun) {
479 m_data->m_histP[ros][drawer][chan][gain][0]->Fill(i + 0.1, dig, 1.0);
480 m_data->m_hist1[ros][drawer][chan][gain][2]->Fill(dig, 1.0);
481 } else if (m_cispar[6] > 0. && (dig > 0 || gain < 1 || charge < 12.)) { //LF If CIS run, Protection to avoid zero amplitudes due to 0 injected charge
482 m_data->m_histP[ros][drawer][chan][gain][0]->Fill(i + 0.1, dig, 1.0);
483 m_data->m_hist1[ros][drawer][chan][gain][2]->Fill(dig, 1.0);
484 }
485
486 if (gain == 1 && dig > m_ADCmaxMinusEps) // overflow in high gain, find charge is it
487 m_data->m_outInHighGain[0]->Fill(charge);
488 }
489 }
490
491 if (dsize > 0 && !m_data->m_corrup[ros][drawer][gain][chan / 3]) {
492 double ped = vdigits[0];
493 m_data->m_hist1[ros][drawer][chan][gain][0]->Fill(ped, 1.0);
494 m_data->m_sumPed1[ros][drawer][chan][gain] += ped;
495 m_data->m_sumPed2[ros][drawer][chan][gain] += ped * ped;
496 //if ( (chan==23)&&(ros==1)&&(drawer==18)&&(gain==0)) {
497 //std::cout << "ped=" << ped << "\tm_sumPed1=" << m_data->m_sumPed1[ros][drawer][chan][gain] << "\n";
498 //std::cout << "ped^2=" << ped*ped << "\tm_sumPed2=" << m_data->m_sumPed2[ros][drawer][chan][gain] << "\n";
499 //}
500 if (dsize > 1) {
501 meansamp /= dsize;
502 rmssamp = rmssamp / dsize - meansamp * meansamp;
503 rmssamp = (rmssamp > 0.0) ? sqrt(rmssamp * dsize / (dsize - 1)) : 0.0;
504 m_data->m_hist1[ros][drawer][chan][gain][1]->Fill(rmssamp, 1.0);
505 m_data->m_sumRms1[ros][drawer][chan][gain] += rmssamp;
506 m_data->m_sumRms2[ros][drawer][chan][gain] += rmssamp * rmssamp;
507 }
508 }
509 } // digits
510
511 //For cor&cov
512 for (int sample = 0; sample < m_nSamples; ++sample) {
513 for (int gain = 0; gain < 2; ++gain) {
514 for (int ch_i = 0; ch_i < 48; ++ch_i)
515 if (!m_data->m_corrup[ros][drawer][gain][ch_i / 3]) {
516 m_data->m_nEvents_i[ros][drawer][gain][ch_i]++;
517 m_data->m_meanAmp[ros][drawer][gain][ch_i] += mean_tmp[ch_i][gain][sample];
518 for (int ch_j = 0; ch_j < 48; ++ch_j)
519 if (!m_data->m_corrup[ros][drawer][gain][ch_j / 3]) {
520 m_data->m_nEvents_ij[ros][drawer][gain][ch_i][ch_j]++;
521 m_data->m_meanAmp_ij[ros][drawer][gain][ch_i][ch_j] += mean_tmp[ch_i][gain][sample] * mean_tmp[ch_j][gain][sample];
522 }
523 }
524 }
525 }
526
527 // BCID
528 uint32_t bcid = digitsCollection->getRODBCID();
529
530 for (int ch = 0; ch < headsize; ++ch) {
531 uint32_t bcid_ch = (headerVec[ch] & 0xFFF);
532 m_data->m_hist1[ros][drawer][ch][0][3]->Fill(bcid_ch, 1.0);
533 if ((bcid_ch == bcid) || (bcid_ch == bcid - 1)) // Conservative condition to be consistent with both run before Feb07 and
534 m_data->m_hist_DMUerr[ros][drawer][ch][0][0]->Fill(1.0, 1.0); // runs after Feb07. Introducing a RunNum variable it could be more strict.
535 else if ((bcid == 0) && ((bcid_ch == 3563) || (bcid_ch == 3564))) // if bcid==0 then bcid_ch should be 3563 (0xDEB)
536 m_data->m_hist_DMUerr[ros][drawer][ch][0][0]->Fill(1.0, 1.0); // but sometimes 3564 (0xDEC) is observed.
537 // if (bcid_ch == bcid) // Now allow only exact bcid: ROD BCID = DMU BCID+1
538 // m_data->m_hist_DMUerr[ros][drawer][ch][0][0]->Fill(1.0,1.0); // Apr 2013
539 // else if (bcid_ch == 0)
540 // m_data->m_hist_DMUerr[ros][drawer][ch][0][0]->Fill(0.0,1.0);
541 // else
542 // m_data->m_hist_DMUerr[ros][drawer][ch][0][0]->Fill(2.0,1.0);
543 }
544
545 //DMUheaderCheck(&headerVec,headsize,ros,drawer,0);
546
547 if ((m_bigain) && ((digitsCollection->getFragChipHeaderWordsHigh()).size() > 0)) {// LF we need to avoid that headsize is set to zero in monogain runs
548
549 headerVec = digitsCollection->getFragChipHeaderWordsHigh();
550 headsize = headerVec.size();
551 if (headsize > 16) {
552 ATH_MSG_ERROR("length of BCIDhigh vector " << headsize << " - greater than 16 !");
553 headsize = 16;
554 }
555 }
556 // if monogain run with m_bigain==1, we fill the BCID plots we the same info
557 for (int ch = 0; ch < headsize; ++ch) {
558 uint32_t bcid_ch = (headerVec[ch] & 0xFFF);
559 m_data->m_hist1[ros][drawer][ch][1][3]->Fill(bcid_ch, 1.0);
560 if ((bcid_ch == bcid) || (bcid_ch == bcid - 1)) // BCID from TileDMU should match BCID from ROD header or
561 m_data->m_hist_DMUerr[ros][drawer][ch][1][0]->Fill(1.0, 1.0); // bcid-1 for runs after Feb07.
562 else if ((bcid == 0) && ((bcid_ch == 3563) || (bcid_ch == 3564))) // if bcid==0 then bcid_ch should be 3563 (0xDEB)
563 m_data->m_hist_DMUerr[ros][drawer][ch][1][0]->Fill(1.0, 1.0); // but sometimes 3564 (0xDEC) is observed.
564 // if (bcid_ch == bcid) // Now allow only exact bcid
565 // m_data->m_hist_DMUerr[ros][drawer][ch][0][0]->Fill(1.0,1.0); // Apr 2013
566 // else if (bcid_ch == 0)
567 // m_data->m_hist_DMUerr[ros][drawer][ch][1][0]->Fill(0.0,1.0);
568 // else
569 // m_data->m_hist_DMUerr[ros][drawer][ch][1][0]->Fill(2.0,1.0);
570 }
571
572 //DMUheaderCheck(&headerVec,headsize,ros,drawer,1);
573
574 if (dqStatus->calibMode() == 1) {
575 // global and DMU CRC check
576 uint32_t crc32 = digitsCollection->getFragCRC();
577 uint32_t CRCmask = digitsCollection->getFragDMUMask();
578
579 CRCcheck(dqStatus, crc32, CRCmask, headsize, ros, drawer);
580 }
581 } //loop over drawers
582 //
583
584 if (dqStatus->calibMode() == 0) {
585 if (RODCRCcalc(dqStatus).isFailure()) ATH_MSG_WARNING("Not possible to check CRC from ROD");
586 }
587
588 return StatusCode::SUCCESS;
589}
590
591/*---------------------------------------------------------*/
593/*---------------------------------------------------------*/
594{
595
596 ATH_MSG_INFO("in finalHists()");
597 if (m_allHistsFilled) return StatusCode::SUCCESS;
598 m_allHistsFilled = true;
599
600 const EventContext &ctx = Gaudi::Hive::currentContext();
601
602 const char *part[5] = { "AUX", "LBA", "LBC", "EBA", "EBC" };
603 const char *gain[6] = { "_lo", "_hi", "", " low gain", " high gain", "" };
604
605 std::string pedestalTitile(" Pedestal[0]");
606 if (m_fillPedestalDifference) pedestalTitile += " - pedestal from DB";
607
608 const char *HistName[8] = { "_ped", "_rms_lfr", "_rms_hfr", "_bits"
609 , pedestalTitile.c_str(), " RMS noise low frequency", " RMS noise high frequency", " Stuck bits and zero amp" };
610
611 const char *ErrName[4] = { "_bcid", "_crc", " BCID errors", " CRC errors" };
612 const char *HistName2[4] = { "_covar", "_corr", " covariance", " correlation" };
613
614 // for bigain run book 2 histograms per channel
615 // for monogain run book just one histogram per channel
616 int mingain = (m_bigain) ? 0 : 2;
617 int maxgain = (m_bigain) ? 2 : 3;
618
619 for (unsigned int ros = 1; ros < TileCalibUtils::MAX_ROS; ++ros) {
620 for (unsigned int drawer = 0; drawer < TileCalibUtils::MAX_DRAWER ; ++drawer) {
621 if (m_data->m_hist1[ros][drawer][0][0].size() != 0) {
622 unsigned int drawerIdx = TileCalibUtils::getDrawerIdx(ros, drawer);
623
624 std::ostringstream sStr;
625 std::string moduleName = TileCalibUtils::getDrawerString(ros, drawer);
626 std::string subDir = "Summary";
627 std::string histName, histTitle;
628
629 for (int gn = mingain; gn < maxgain; ++gn) {
630
631 int adc = gn & 1;
632
633 for (int type = 0; type < 4; ++type) {
634 sStr.str("");
635 sStr << moduleName << gain[gn] << HistName[type];
636 histName = sStr.str();
637 sStr.str("");
638 sStr << moduleName << gain[3 + gn] << HistName[4 + type];
639 histTitle = sStr.str();
640 m_data->m_final_hist1[ros][drawer][adc].push_back(book1F(subDir, histName, histTitle, 48, 0.0, 48.0));
641 }
642
643 // stuck bits and saturations / 0s
644 sStr.str("");
645 sStr << moduleName << gain[gn] << "_stucks";
646 histName = sStr.str();
647 sStr.str("");
648 sStr << moduleName << gain[3 + gn] << " Stuck bits and saturation";
649 histTitle = sStr.str();
650 m_data->m_final_hist_stucks[ros][drawer][adc] = book2C(subDir, histName, histTitle, 48, 0.0, 48.0, 7, 0., 7.);
651 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(1, "SB 0");
652 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(2, "SB 1,2");
653 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(3, "SB 3,4");
654 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(4, "SB 5-9");
655 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(5, "zeros");
656 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(6, "saturation");
657 m_data->m_final_hist_stucks[ros][drawer][adc]->GetYaxis()->SetBinLabel(7, "n % 32");
658
659 for (unsigned int channel = 0; channel < TileCalibUtils::MAX_CHAN; ++channel) {
660
661 // int pmt = abs(m_cabling->channel2hole(ros,channel));
662
663 double Ped = 0.0, PedRMS = 0.0, RMSHi = 0.0;
664
665 int nevents = int(m_data->m_hist1[ros][drawer][channel][adc][0]->GetEntries()); //ROOT GetEntries return a Double_t.
666 bool fillPedAndRms = (nevents > 0);
667
668 if ((nevents != 0) && (nevents != m_nEventsTileMon)) {
681 ATH_MSG_VERBOSE( "Number of entries in histo " << m_data->m_hist1[ros][drawer][channel][adc][0]->GetTitle()
682 << " doesnt match n. of events! Using the first one in mean and RMS calculation" );
683
684 ATH_MSG_VERBOSE( "Number of entries in histo =" << nevents << "\tNumber of events= " << m_nEventsTileMon );
685
686 } else { //all ok
687 nevents = m_nEventsTileMon;
688 }
689
690 if (fillPedAndRms) {
691
692 if (nevents > 0) {
693 Ped = m_data->m_sumPed1[ros][drawer][channel][adc] / nevents;
694 RMSHi = m_data->m_sumRms1[ros][drawer][channel][adc] / nevents;
695
696 if (nevents > 1) {
697 PedRMS = m_data->m_sumPed2[ros][drawer][channel][adc] / nevents - Ped * Ped;
698 PedRMS = (PedRMS > 0.0) ? sqrt(PedRMS * nevents / (nevents - 1)) : 0.0;
699 }
700 }
701
702 // if ( (pmt==24)&&(ros==1)&&(drawer==18)&&(adc==0)) {
703 //std::cout << "Evt = " << m_nEventsTileMon << "\t channel=" << channel << " pmt=" << pmt << "\n";
704 //std::cout << "Ped = " << Ped << "\n";
705 //std::cout << "E(x^2) = " << m_data->m_sumPed2[ros][drawer][channel][adc] << "\n";
706 //std::cout << "PedRMS = " << PedRMS << "\n";
707
708 //}
709
711
712 if (!isDisconnected(ros, drawer, channel)) {
713 m_data->m_final_hist1[ros][drawer][adc][0]->SetBinContent(channel + 1, Ped - m_tileToolNoiseSample->getPed(drawerIdx, channel, adc, TileRawChannelUnit::ADCcounts, ctx));
714 m_data->m_final_hist1[ros][drawer][adc][0]->SetBinError(channel + 1, PedRMS);
715 }
716
717 } else {
718 m_data->m_final_hist1[ros][drawer][adc][0]->SetBinContent(channel + 1, Ped);
719 m_data->m_final_hist1[ros][drawer][adc][0]->SetBinError(channel + 1, PedRMS);
720 }
721
722 m_data->m_final_hist1[ros][drawer][adc][1]->SetBinContent(channel + 1, PedRMS);
723 m_data->m_final_hist1[ros][drawer][adc][2]->SetBinContent(channel + 1, RMSHi);
724
725 }
726
727 // For tests
728 //stuckBits_maker(m_data->m_hist1[ros][drawer][channel][adc][2]);
729
730 TH1S * hist = m_data->m_hist1[ros][drawer][channel][adc][2];
731 double weight = 0.0;
732 for (int i = 2; i < m_i_ADCmax + 2; ++i) { // first bin is excluded
733 weight += hist->GetBinContent(i);
734 }
735 // if we have weight=0 after the loop, then there is no non-zero amplitudes
736 if (weight > 0.0) {
737 if (hist->GetBinContent(1) > 0) {
738 weight = 0.5; // sometimes zeros
739 } else {
740 weight = 1.0; // assume that all OK
741 }
742
743 int stuck_bit = stuckBits(hist, adc);
744 if (stuck_bit > 0) {
745 weight += 1.0; // some stuck bits
746 }
747 }
748 m_data->m_final_hist1[ros][drawer][adc][3]->SetBinContent(channel + 1, weight);
749 stuckBits_Amp2(hist, m_data->m_histC[ros][drawer][channel][adc][0], adc,
750 m_data->m_final_hist_stucks[ros][drawer][adc], channel,
751 m_data->m_stuck_probs[ros][drawer][channel][adc]);
752 } // end of loop over channels
753
754 // BCID
755 sStr.str("");
756 sStr << moduleName << gain[gn] << ErrName[0];
757 histName = sStr.str();
758 sStr.str("");
759 sStr << moduleName << gain[3 + gn] << ErrName[2];
760 histTitle = sStr.str();
761 m_data->m_final_hist1[ros][drawer][adc].push_back(book1F(subDir, histName, histTitle, 16, 0.0, 16.0));
762
763 // CRC only for one gain
764 if (adc == 0) {
765 sStr.str("");
766 sStr << moduleName << ErrName[1];
767 histName = sStr.str();
768 sStr.str("");
769 sStr << moduleName << ErrName[3];
770 histTitle = sStr.str();
771 m_data->m_final_hist1[ros][drawer][adc].push_back(book1F(subDir, histName, histTitle, 16, 0.0, 16.0));
772 }
773
774 for (int ch = 0; ch < 16; ++ch) {
775
776 int dmu = ch + 1;
777
778 // BCID and CRC errors
779 for (int id = 0; id < 2; ++id) {
780 int bin0 = lround(m_data->m_hist_DMUerr[ros][drawer][ch][adc][id]->GetBinContent(0));
781 int bin1 = lround(m_data->m_hist_DMUerr[ros][drawer][ch][adc][id]->GetBinContent(1));
782 int bin2 = lround(m_data->m_hist_DMUerr[ros][drawer][ch][adc][id]->GetBinContent(2));
783 int bin3 = lround(m_data->m_hist_DMUerr[ros][drawer][ch][adc][id]->GetBinContent(3));
784 int bin4 = lround(m_data->m_hist_DMUerr[ros][drawer][ch][adc][id]->GetBinContent(4));
785 int bin5 = lround(m_data->m_hist_DMUerr[ros][drawer][ch][adc][id]->GetBinContent(5));
786 double weight = -1.0;
787
788 if (bin0 + bin3 + bin4 + bin5 > 0) {
789 if ((bin4 == 0) && (bin5 == 0)) {
790 weight = 2.0; // BCID: all other errors; CRC: FE only errors
791 } else if ((bin3 == 0) && (bin5 == 0)) {
792 weight = 3.0; // CRC only: ROD only errors
793 } else { // if there are ROD and FE errors (not necessary at the same time),
794 //weight is 2.5
795 weight = 2.5; // CRC only: ROD and FE errors
796 }
797 } else {
798 if (bin2 == 0) {
799 weight = 0.0; // all zeros
800 } else if (bin1 == 0) {
801 weight = 1.0; // all OK
802 } else {
803 weight = 1.0 - (double) bin1 / m_nEventsTileMon; // percentage of good events
804 if (weight > 0.8) weight = 0.8; // to see clearly even one event with zeros
805 }
806 }
807 m_data->m_final_hist1[ros][drawer][adc][id + 4]->SetBinContent(dmu, weight);
808
809 if (adc) break; // no CRC histogram in high gain
810 }
811 } // end of loop over 16 DMU
812
813 //Correlation and Covariance
814 if (m_runType == PedRun) {
815
816 for (int type = 0; type < 2; ++type) {
817 sStr.str("");
818 sStr << moduleName << gain[gn] << HistName2[type];
819 histName = sStr.str();
820 sStr.str("");
821 sStr << moduleName << gain[3 + gn] << HistName2[2 + type];
822 histTitle = sStr.str();
823 m_data->m_final_hist2[ros][drawer][adc].push_back(book2F(subDir, histName, histTitle, 48, 0.0, 48.0, 48, 0.0, 48.0));
824 }
825
826 if (m_nEventsTileMon * m_nSamples > 0) {
827
828 for (int ch_i = 0; ch_i < 48; ++ch_i) {
829 if (m_data->m_nEvents_i[ros][drawer][adc][ch_i] > 0)
830 m_data->m_meanAmp[ros][drawer][adc][ch_i] /= m_data->m_nEvents_i[ros][drawer][adc][ch_i];
831 for (int ch_j = 0; ch_j < 48; ++ch_j)
832 if (m_data->m_nEvents_ij[ros][drawer][adc][ch_i][ch_j] > 0)
833 m_data->m_meanAmp_ij[ros][drawer][adc][ch_i][ch_j] /= m_data->m_nEvents_ij[ros][drawer][adc][ch_i][ch_j];
834 }
835
836 //coverity[STACK_USE]
837 double covar[48][48];
838 //coverity[STACK_USE]
839 double corr[48][48];
840 double mean_rms[48];
841 double mean_cov_ii = 0.; // mean values of covar in diag terms
842 double mean_cov_ij = 0.; // mean values of covar in off diag terms
843
844 for (int ch_i = 0; ch_i < 48; ++ch_i) {
845 // int pmt_i = abs(m_cabling->channel2hole(ros,ch_i))-1;
846 for (int ch_j = 0; ch_j < 48; ++ch_j) {
847 // int pmt_j = abs(m_cabling->channel2hole(ros,ch_j))-1;
848 covar[ch_i][ch_j] = m_data->m_meanAmp_ij[ros][drawer][adc][ch_i][ch_j] - m_data->m_meanAmp[ros][drawer][adc][ch_i] * m_data->m_meanAmp[ros][drawer][adc][ch_j];
849
850 if (ch_j < ch_i) {
851 mean_cov_ij += covar[ch_i][ch_j]; //LF: we take C_ij with its sign
852 }
853 }
854
855 mean_cov_ii += covar[ch_i][ch_i];
856
857 if (covar[ch_i][ch_i] != 0.0)
858 mean_rms[ch_i] = covar[ch_i][ch_i] / sqrt(TMath::Abs(covar[ch_i][ch_i]));
859 else
860 mean_rms[ch_i] = 0.;
861 }
862
863 if (mean_cov_ii != 0.) {
864 m_data->m_cov_ratio[ros][drawer][adc] = (2. * mean_cov_ij) / (mean_cov_ii * 47.); //(2*cov_ij/(48*47))/(cov_ii/48)
865 } else {
866 m_data->m_cov_ratio[ros][drawer][adc] = 0.;
867 }
868 // std::cout << "m_cov_ratio["<<adc<<"]= " << m_data->m_cov_ratio[ros][drawer][adc] << "\n";
869
870 for (int i = 0; i < 48; ++i) {
871 for (int j = 0; j < 48; ++j) {
872 if ((mean_rms[i] == 0.) || (mean_rms[j] == 0.)) {
873 corr[i][j] = 0.;
874 } else {
875 corr[i][j] = covar[i][j] / mean_rms[i] / mean_rms[j];
876 }
877 m_data->m_final_hist2[ros][drawer][adc][0]->SetBinContent(i + 1, j + 1, covar[i][j]);
878 m_data->m_final_hist2[ros][drawer][adc][1]->SetBinContent(i + 1, j + 1, corr[i][j]);
879 }
880 }
881 } // end if m_nEventsTileMon > 0
882
883 } // end if PedRun
884 } // end of loop over gn
885
886 drawHists(ros, drawer, moduleName);
887
888 //if RunType is Cis, the Kolmogorov test of the sample distribution is performed
889 if (m_runType == CisRun) {
890
891 for (int gn = mingain; gn < maxgain; ++gn) {
892 for (int i = 0; i < 48; ++i) {
893 TH1S *h = (TH1S*) m_data->m_hist1[ros][drawer][i][gn][2]->Clone("temphist");
894 h->SetDirectory(0);
895 shiftHisto(h, ros, drawer, i, gn);
896 delete h;
897 }
898 double shiftxmin = 0.;
899 double shiftxmax = m_shiftnbins;
900 int shiftnbins = m_shiftnbins;
901
902 std::ostringstream s;
903 s << part[ros] << std::setfill('0') << std::setw(2) << drawer + 1 << std::setfill(' ');
904 std::string moduleName = s.str();
905 std::string subDir = moduleName;
906
907 s.str("");
908 s << moduleName << gain[gn] << "_samples_ref";
909 std::string histName = s.str();
910
911 s.str("");
912 s << moduleName << gain[3 + gn] << " all samples reference";
913 std::string histTitle = s.str();
914
915 m_data->m_shifted_hist[ros][drawer][48][gn] = book1S(subDir, histName, histTitle, shiftnbins, shiftxmin, shiftxmax);
916
917 statTestHistos(ros, drawer, gn);
918 } // if runTYpe==CisRun
919 }
920 }
921
922 } //loop over drawer
923 } //loop over ros
924
925
926 return StatusCode::SUCCESS;
927}
928
929/*---------------------------------------------------------*/
930const uint8_t * TileDigitsMonTool::stuckBitProb (int ros, int module, int channel, int gain) const
931/*---------------------------------------------------------*/
932{
933 return m_data->m_stuck_probs[ros][module][channel][gain];
934}
935
936
937/*---------------------------------------------------------*/
938StatusCode TileDigitsMonTool::checkHists(bool /* fromFinalize */)
939/*---------------------------------------------------------*/
940{
941
942 ATH_MSG_INFO("in checkHists()");
943
944 return StatusCode::SUCCESS;
945}
946
947/*---------------------------------------------------------*/
948int TileDigitsMonTool::define_palette(int ncolors, int *colors) {
949 /*---------------------------------------------------------*/
950 Float_t alpha = 1.; // 0- transparent, 1- opaque
951 Int_t ind;
952 // set Deep Sea Blue palette (dark to light)
953 if (ncolors == 51 && colors == 0) {
954 TColor::InitializeColors();
955 const int nRGBs = 5;
956 Double_t stops[nRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
957 Double_t red[nRGBs] = { 0.00, 0.09, 0.18, 0.09, 0.00 };
958 Double_t green[nRGBs] = { 0.01, 0.02, 0.39, 0.68, 0.97 };
959 Double_t blue[nRGBs] = { 0.17, 0.39, 0.62, 0.79, 0.97 };
960 ind = TColor::CreateGradientColorTable(nRGBs, stops, red, green, blue, m_NCont, alpha);
961 }
962
963 // set Grey Scale palette (light to dark)
964 else if (ncolors == 52 && colors == 0) {
965 TColor::InitializeColors();
966 const int nRGBs = 3;
967 Double_t stops[nRGBs] = { 0.00, 0.50, 1.00 };
968 Double_t red[nRGBs] = { 0.90, 0.50, 0.10 };
969 Double_t green[nRGBs] = { 0.90, 0.50, 0.10 };
970 Double_t blue[nRGBs] = { 0.90, 0.50, 0.10 };
971 ind = TColor::CreateGradientColorTable(nRGBs, stops, red, green, blue, m_NCont, alpha);
972 } else
973 return -1;
974 return ind;
975}
976
977/*---------------------------------------------------------*/
978void TileDigitsMonTool::drawHists(int ros, int drawer, const std::string& moduleName)
979/*---------------------------------------------------------*/
980{
981
982 ATH_MSG_VERBOSE("in drawHists()");
983
984 int maxgain = (m_bigain) ? 2 : 1;
985 double ms = (m_bigain) ? 0.75 : 1.0; // marker size
986
987 static bool palette_init_done = false;
988 static int pal51[m_NCont];
989 static int pal52[m_NCont];
990 static TExec *ex_pal51;
991 static TExec *ex_pal52;
992 bool do_plots = m_savePng || m_savePs || m_saveSvg;
993
994 if (palette_init_done == false) {
995 palette_init_done = true;
996 int i;
997 int ind = define_palette(51);
998 if (ind < 0) ind = 0; // fall-back colours
999 for (i = 0; i < m_NCont; i++)
1000 pal51[i] = ind + i;
1001 ind = define_palette(52);
1002 if (ind < 0) ind = 0; // fall-back colours
1003 for (i = 0; i < m_NCont; i++)
1004 pal52[i] = ind + i;
1005 TString s;
1006 s.Form("gStyle->SetPalette(%d,(int*)%ld);", m_NCont, (long) pal51);
1007 ex_pal51 = new TExec("ex_pal51", s.Data());
1008 s.Form("gStyle->SetPalette(%d,(int*)%ld);", m_NCont, (long) pal52);
1009 ex_pal52 = new TExec("ex_pal52", s.Data());
1010 }
1011
1012 TLine * maxline[4];
1013 double line[8] = { 1.2, 2.5, 1.0, 1.0, 48.0, 48.0, 48.0, 16.0 };
1014 double miny[12] = { 20.0, 20.0, 0.0, 0.0, 0.0, 0.0, -0.1, -0.1, -0.1, -0.1, -0.1, -0.1 };
1015 double maxy[12] = { 80.0, 80.0, 1.5, 3.0, 1.5, 3.0, 2.2, 2.2, 2.2, 2.2, 3.2, 2.2 };
1016 // ? , ? , ? , ? , ? , ? , ? , ? , BC hi,BC lo, CRC , ?
1018 miny[0] = miny[1] = -5.0;
1019 maxy[0] = maxy[1] = 5.0;
1020 }
1021
1022 TText *t = new TText();
1023 t->SetTextAlign(32);
1024 t->SetTextSize(0.06);
1025 t->SetTextFont(72);
1026
1027 if (maxgain == 1) {
1028 line[0] = line[1];
1029 for (int ind = 0; ind < 12; ind += 2) {
1030 miny[ind] = miny[ind + 1] = std::min(miny[ind], miny[ind + 1]);
1031 maxy[ind] = maxy[ind + 1] = std::max(maxy[ind], maxy[ind + 1]);
1032 }
1033 }
1034
1035 for (int type = 0; type < 4; ++type) {
1036 maxline[type] = new TLine(0.0, line[type], line[4 + type], line[type]);
1037 maxline[type]->SetLineColor(2 + type / 2);
1038 }
1039
1040 TCanvas * Can = NULL;
1041 if (do_plots) {
1042 Can = new TCanvas("dig_noise", "dig_noise", 402 * maxgain, 588);
1043 Can->Divide(maxgain, 3);
1044 gStyle->SetOptStat(0);
1045 gStyle->SetTitleFontSize(0.1);
1046 }
1047
1048 for (int gain = 0; gain < maxgain; ++gain) {
1049 for (int type = 0; type < 3; ++type) {
1050
1051 if (do_plots) {
1052 TVirtualPad * pad = Can->cd(type * maxgain + gain + 1);
1053 pad->SetTopMargin(0.15);
1054 pad->SetGridx();
1055 }
1056 int ind = type * 2 + gain;
1057
1058 if (m_data->m_final_hist1[ros][drawer][gain][type]->GetMaximum() < 0.9 * maxy[ind]) //if maximum is below reasonable limit use a default scale
1059 m_data->m_final_hist1[ros][drawer][gain][type]->SetMaximum(maxy[ind]);
1060 if (m_data->m_final_hist1[ros][drawer][gain][type]->GetMinimum() > (miny[ind] + 0.1 * TMath::Abs(miny[ind]))) // if minimum is above reasonable limit
1061 m_data->m_final_hist1[ros][drawer][gain][type]->SetMinimum(miny[ind]); // use a default scale
1062
1063 m_data->m_final_hist1[ros][drawer][gain][type]->SetMarkerStyle(21);
1064 m_data->m_final_hist1[ros][drawer][gain][type]->SetMarkerSize(ms);
1065 m_data->m_final_hist1[ros][drawer][gain][type]->SetLabelSize(0.08, "X");
1066 m_data->m_final_hist1[ros][drawer][gain][type]->SetLabelSize(0.08, "Y");
1067 if (do_plots) {
1068 m_data->m_final_hist1[ros][drawer][gain][type]->Draw("P0");
1069 if (type > 0) maxline[gain]->Draw();
1070 }
1071 }
1072 }
1073
1074 if (m_savePng) {
1075 Can->Print(TString(moduleName + "_dig_noise.png"), "png");
1076 }
1077 if (m_savePs) {
1078 Can->Print(TString(moduleName + "_dig_noise.ps"), "ps");
1079 }
1080 if (m_saveSvg) {
1081 Can->Print(TString(moduleName + "_dig_noise.svg"), "svg");
1082 }
1083 if (do_plots) {
1084 delete Can;
1085
1086 Can = new TCanvas("dig_errors", "dig_errors", 402 * maxgain, 588);
1087 Can->Divide(maxgain, 4);
1088 gStyle->SetOptStat(1);
1089 gStyle->SetTitleFontSize(0.1);
1090 }
1091
1092 for (int gain = 0; gain < maxgain; ++gain) {
1093 for (int type = 3; type < 6; ++type) {
1094
1095 TVirtualPad * pad = NULL;
1096 if (do_plots) {
1097 pad = Can->cd((type - 3) * maxgain + gain + 1);
1098 pad->SetTopMargin(0.15);
1099 }
1100 int ind = type * 2 + gain;
1101
1102 if (type == 5 && gain == 1) {
1103 if (do_plots) pad->SetLogy();
1104
1105 m_data->m_hist0[ros][drawer][0]->SetStats(kFALSE);
1106 gStyle->SetStatFontSize(0.1);
1107 m_data->m_hist0[ros][drawer][0]->SetStats(kTRUE);
1108 m_data->m_hist0[ros][drawer][0]->SetLabelSize(0.08, "X");
1109 m_data->m_hist0[ros][drawer][0]->SetLabelSize(0.08, "Y");
1110 if (do_plots) m_data->m_hist0[ros][drawer][0]->Draw(); // global CRC in bottom right corner
1111
1112 } else {
1113 if (do_plots) pad->SetGridx();
1114
1115 if (m_data->m_final_hist1[ros][drawer][gain][type]->GetMaximum() < 0.9 * maxy[ind]) //if maximum is below reasonable limit use a default scale
1116 m_data->m_final_hist1[ros][drawer][gain][type]->SetMaximum(maxy[ind]);
1117 if (m_data->m_final_hist1[ros][drawer][gain][type]->GetMinimum() > (miny[ind] + 0.1 * TMath::Abs(miny[ind]))) // if minimum is above reasonable limit
1118 m_data->m_final_hist1[ros][drawer][gain][type]->SetMinimum(miny[ind]); // use a default scale
1119
1120 m_data->m_final_hist1[ros][drawer][gain][type]->SetMarkerStyle(21);
1121 m_data->m_final_hist1[ros][drawer][gain][type]->SetMarkerSize(ms);
1122 m_data->m_final_hist1[ros][drawer][gain][type]->SetStats(kFALSE);
1123 m_data->m_final_hist1[ros][drawer][gain][type]->SetLabelSize(0.08, "X");
1124 m_data->m_final_hist1[ros][drawer][gain][type]->SetLabelSize(0.08, "Y");
1125 if (type != 3) {
1126 if (do_plots) m_data->m_final_hist1[ros][drawer][gain][type]->Draw("P0");
1127 } else {
1128 m_data->m_final_hist_stucks[ros][drawer][gain]->SetStats(kFALSE);
1129 m_data->m_final_hist_stucks[ros][drawer][gain]->SetMarkerSize(3.);
1130 m_data->m_final_hist_stucks[ros][drawer][gain]->SetMarkerColor(4);
1131 m_data->m_final_hist_stucks[ros][drawer][gain]->SetLabelSize(0.08, "X");
1132 m_data->m_final_hist_stucks[ros][drawer][gain]->SetLabelSize(0.1, "Y");
1133 gStyle->SetNumberContours(m_NCont);
1134 if (do_plots) {
1135 m_data->m_final_hist_stucks[ros][drawer][gain]->Draw("textzcol");
1136 ex_pal52->Draw();
1137 m_data->m_final_hist_stucks[ros][drawer][gain]->Draw("textcolzsame");
1138 pad->Update();
1139 TPaletteAxis *palette = (TPaletteAxis*) m_data->m_final_hist_stucks[ros][drawer][gain]->GetListOfFunctions()->FindObject("palette");
1140 if (palette != NULL) palette->SetLabelSize(0.07);
1141 pad->Modified();
1142 }
1143 }
1144 if (type == 4 || (type == 5 && gain == 0)) {
1145 m_data->m_final_hist1[ros][drawer][gain][type]->GetYaxis()->SetLabelOffset(-0.85); // do not draw default lables
1146 if (do_plots) {
1147 t->DrawText(-0.2, 1., "all OK");
1148 t->DrawText(-0.2, .0, "wrong");
1149 if (type == 4)
1150 t->DrawText(-0.2, 2., "mismatch");
1151 else {
1152 t->DrawText(-0.2, 2., "FE fail");
1153 t->DrawText(-0.2, 2.5, "FE+ROD fail");
1154 t->DrawText(-0.2, 3., "ROD fail");
1155 }
1156 }
1157 }
1158 /* if (type == 3)
1159 maxline[2]->Draw();
1160 else */
1161 if (type != 3) if (do_plots) maxline[3]->Draw();
1162 }
1163 }
1164 //LF Now plot DMU Header Errors histograms
1165 TVirtualPad * pad = NULL;
1166 if (do_plots) {
1167 pad = Can->cd(6 + gain + 1);
1168 pad->SetTopMargin(0.15);
1169 pad->SetGridx();
1170 }
1171 m_data->m_hist2[ros][drawer][gain][0]->SetMarkerSize(3.);
1172 m_data->m_hist2[ros][drawer][gain][0]->SetMarkerColor(2);
1173 m_data->m_hist2[ros][drawer][gain][0]->SetLabelSize(0.11, "Y");
1174 m_data->m_hist2[ros][drawer][gain][0]->SetLabelSize(0.08, "X");
1175 m_data->m_hist2[ros][drawer][gain][0]->SetStats(kFALSE);
1176 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(1, "OK ");
1177 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(2, "Format");
1178 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(3, "Parity");
1179 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(4, "Memory");
1180 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(5, "SingleStr");
1181 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(6, "DbleStr");
1182 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(7, "DummyFrag");
1183 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetBinLabel(8, "NoDataFrag");
1184 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetTickLength(0.01);
1185 m_data->m_hist2[ros][drawer][gain][0]->GetYaxis()->SetLabelOffset(0.001);
1186 if (do_plots) {
1187 gStyle->SetNumberContours(m_NCont);
1188 m_data->m_hist2[ros][drawer][gain][0]->Draw("textzcol");
1189 ex_pal51->Draw();
1190 m_data->m_hist2[ros][drawer][gain][0]->Draw("textzcolsame");
1191 pad->Update();
1192 TPaletteAxis *palette = (TPaletteAxis*) m_data->m_hist2[ros][drawer][gain][0]->GetListOfFunctions()->FindObject("palette");
1193 if (palette != NULL) palette->SetLabelSize(0.07);
1194 pad->Modified();
1195 }
1196
1197 }
1198
1199 if (m_savePng) {
1200 Can->Print(TString(moduleName + "_dig_errors.png"), "png");
1201 }
1202 if (m_savePs) {
1203 Can->Print(TString(moduleName + "_dig_errors.ps"), "ps");
1204 }
1205 if (m_saveSvg) {
1206 Can->Print(TString(moduleName + "_dig_errors.svg"), "svg");
1207 }
1208 if (do_plots) delete Can;
1209
1210 for (int type = 0; type < 4; ++type) {
1211 delete maxline[type];
1212 }
1213
1214 delete t;
1215
1216 if (m_runType == PedRun) {
1217
1218 if (do_plots) {
1219 Can = new TCanvas("correlation", "correlation", 402 * maxgain, 588);
1220 Can->Divide(maxgain, 2);
1221 gStyle->SetOptStat(0);
1222 gStyle->SetPalette(1);
1223 gStyle->SetTitleFontSize(0.1);
1224 }
1225
1226 TPaveLabel covar_label[2];
1227 covar_label[0] = TPaveLabel(37., 49., 54., 56., "");
1228 covar_label[1] = TPaveLabel(37., 49., 54., 56., "");
1229
1230 for (int gain = 0; gain < maxgain; ++gain) {
1231 for (int type = 0; type < 2; ++type) {
1232 if (do_plots) {
1233 TVirtualPad * pad = Can->cd(type * maxgain + gain + 1);
1234 pad->SetTopMargin(0.15);
1235 }
1236 m_data->m_final_hist2[ros][drawer][gain][type]->SetLabelSize(0.06, "X");
1237 m_data->m_final_hist2[ros][drawer][gain][type]->SetLabelSize(0.06, "Y");
1238 if (do_plots) {
1239 m_data->m_final_hist2[ros][drawer][gain][type]->Draw("COLZ");
1240 if (type == 0) {
1241 std::ostringstream label_text;
1242 label_text << "<C_ij>/<C_ii>= " << std::setprecision(4) << m_data->m_cov_ratio[ros][drawer][gain] << std::setprecision(4);
1243 covar_label[gain].SetLabel(label_text.str().c_str());
1244 covar_label[gain].SetTextSize(0.4);
1245 covar_label[gain].Draw();
1246 }
1247 }
1248 }
1249 }
1250
1251 if (m_savePng) {
1252 Can->Print(TString(moduleName + "_correlation.png"), "png");
1253 }
1254 if (m_savePs) {
1255 Can->Print(TString(moduleName + "_correlation.ps"), "ps");
1256 }
1257 if (do_plots) delete Can;
1258 }
1259
1260 gStyle->SetOptStat(1111);
1261}
1262
1263int TileDigitsMonTool::stuckBits(TH1S * hist, int adc) {
1264 // function returns number of stuck bits
1265 int NSB = 0;
1266
1267 // If all bins at given hypothetical stuck_position (0,1) and stuck_bit (0,...,9)
1268 // are empty then this bin is stuck in this position. Do this procedure only for
1269 // bins between the first and the last non-zeroth one.
1270
1271 // Find first/last non-zero bin in the histogram
1272 int MinBin = 1;
1273 for (; MinBin < m_i_ADCmax + 2; ++MinBin) {
1274 if (hist->GetBinContent(MinBin) > 0) {
1275 break;
1276 }
1277 }
1278 if (MinBin == m_i_ADCmax + 2) return 0; // empty histogram, nothing to do
1279
1280 int MaxBin = m_i_ADCmax + 1;
1281 for (; MaxBin > 0; --MaxBin) {
1282 if (hist->GetBinContent(MaxBin) > 0) {
1283 break;
1284 }
1285 }
1286
1287 // bins in hist are counted from 1 to m_i_ADCmax + 1, but actual X value for those bins are 0-m_i_ADCmax
1288
1289 // if there is nothing in first half of histogram
1290 // or there is nothing in second half and there is a sharp edge at 512
1291 // it can be that upper most bit is stuck in 0 or 1
1292 if (MinBin == (m_i_ADCmax + 1) / 2 + 1 || (MaxBin == (m_i_ADCmax + 1) / 2 && hist->GetBinContent(MaxBin) > 3)) {
1293 ++NSB;
1294 // std::cout << "Bit 9 is stuck" << std::endl;
1295 }
1296
1297 --MinBin; // shift to zero bin
1298 ++MaxBin; // shift to zero bin
1299
1300 for (int hyp_pos = 0; hyp_pos < 2; ++hyp_pos) { // bit stuck at value 0 or 1
1301 for (int hyp_bit = 0; hyp_bit < 10; ++hyp_bit) { // which bit is stuck
1302 double sum = 0;
1303 int bin_counter = 0;
1304
1305 int win_length = (1 << hyp_bit); // length of one interval with all zeros if bit is stuck
1306 int BeginBin = (1 - hyp_pos) * win_length; // beginning of first interval to check
1307
1308 while (BeginBin + 1 < MaxBin) {
1309 int EndBin = BeginBin + win_length; // end of interval to check
1310
1311 if (EndBin > MinBin) {
1312 for (int bin = BeginBin; bin < EndBin; ++bin) {
1313 //bin is from 0 to m_i_ADCmax here - must be changed in GetBinConent to (bin+1)
1314 if (MinBin < (bin + 1) && (bin + 1) < MaxBin) {
1315 sum += hist->GetBinContent(bin + 1);
1316 ++bin_counter;
1317 }
1318 }
1319 }
1320 BeginBin = EndBin + win_length; // beginning of next interval
1321 }
1322 if (sum == 0 && bin_counter > 0) {
1323 // std::cout << "Bit "<<hyp_bit<<" is stuck at "<< hyp_pos << std::endl;
1324 ++NSB;
1325 }
1326 }
1327 }
1328
1329 if ((NSB == 0) && ((m_runType == CisRun) || (m_runType == CisRamp))) {
1330
1331 NSB = stuckBits_Amp(hist, adc); // if no stuck bits found, lets try another crude procedure
1332 }
1333
1334 return NSB;
1335}
1336
1337/* not used
1338 int TileDigitsMonTool::stuckBits_maker(TH1S * hist)
1339 {
1340
1341 // If all bins at given hypothetical stuck_position and stuck_bit are empty then
1342 // this bin is stuck in this position. Do this procedure only for
1343 // bins less last bin
1344 hp = 1 - hp;
1345 ++hb;
1346 if (hb > 15) hb = 0;
1347 // Stuck bits maker
1348 hist->SetBinContent(1,0);
1349 if (hb < 9) {
1350 int hyp_pos = hp;
1351 int hyp_bit = hb;
1352 int max_Nsample = (m_i_ADCmax + 1)/(1<<(hyp_bit+1));
1353 for (int bin_sample = 0; bin_sample < max_Nsample; ++bin_sample) {
1354 //beginning of sample in histogram
1355 int BeginBin = (1 - hyp_pos)*(1<<hyp_bit)+bin_sample*(1<<(hyp_bit+1));
1356 //end of sample in histogram
1357 int EndBin = BeginBin + (1<<hyp_bit);
1358 for (int bin = BeginBin; bin < EndBin; ++bin) {
1359 //variable bin from 0 to m_i_ADCmax - must be changed in GetBunConent to (bin+1)
1360 if ((bin+1) < m_i_ADCmax + 2) hist->SetBinContent(bin+1,0);
1361 }
1362 }
1363 }
1364 return -1;
1365 }
1366 */
1367
1368// New function by Tibor Zenis
1369// It checks how many bit flips happen compared to the expected ones
1370// Only works for CIS and Laser because the adc range is not sufficiently excersized
1371// in the other run types.
1372int TileDigitsMonTool::stuckBits_Amp(TH1S * hist, int /*adc*/) {
1373 int i, b, c, cc = 2;
1374 int bc[10], bv[10];
1375 double bs[10];
1376 int f; // flip-flop of non-zero contents of previous bin
1377 int cm, cp; // previous and next bin contents
1378 double prob; // probability
1379 for (b = 0; b < 10; b++) {
1380 bc[b] = 0;
1381 bs[b] = 0.;
1382 bv[b] = 1 << b;
1383 }
1384 cp = hist->GetBinContent(1);
1385 f = !!cp;
1386 for (i = 1; i < m_i_ADCmax + 1; i++) {
1387 cm = cc;
1388 c = cc = cp;
1389 cp = (i < m_i_ADCmax) ? hist->GetBinContent(i + 2) : 0;
1390 prob = 1.;
1391 if (c > 0 && c < 0.125 * (cp + cm - std::sqrt((double) cp) - std::sqrt((double) cm))) {
1392 prob = 1. - (double) c / (0.125 * (cp + cm - std::sqrt((double) cp) - std::sqrt((double) cm)));
1393 c = 0;
1394 }
1395 if (c || f) // count bit flips (bc[b]) in bins with non-zero contents and adjacent bins
1396 for (b = 0; b < 10; b++) {
1397 if ((i & bv[b]) != ((i - 1) & bv[b]))
1398 bc[b]++;
1399 else
1400 break; // higher bit can flip only of the previous one does
1401 }
1402 if (f > 0 && c > 0) continue;
1403 if (f == 0 && c == 0) continue;
1404 if (c || f) // bin content flips (zero - non-zero)
1405 for (b = 0; b < 10; b++) { // count bits flipped at this position
1406 if ((i & bv[b]) != ((i - 1) & bv[b]))
1407 bs[b] += (i & bv[b]) ? (c ? prob : -prob) : (c ? -prob : prob);
1408 else
1409 break;
1410 }
1411 f = !f;
1412 }
1413
1414 int NSB = 0;
1415 for (b = 0; b < 10; b++)
1416 // if (bc[b] > 2 && abs(bs[b]) * 2 > bc[b])
1417 if (((bc[b] > 2) && (std::abs(bs[b]) == bc[b])) || ((std::abs(bs[b]) > 7) && (std::abs(bs[b]) * 3 > bc[b]))) {
1418 NSB++;
1419 }
1420
1421 return NSB;
1422
1423}
1424
1425/* version 2. */
1426int TileDigitsMonTool::stuckBits_Amp2(TH1S * hist, TH1C *modhist, int gain, TH2C *outhist, int ch, uint8_t *stuck_probs) {
1427
1428 if (hist->GetEntries() < 1000) return 0; /* too few events (1000 / N_samples) in histogram, ignore */
1429
1430 int i, b, c, cc = 0; /* bin counter, but position, bin content */
1431 double cont;
1432 int last_non0, first_non0 = -1;
1433 int bc[10] = { 0 }, bv[10] = { 1, 2, 4, 8, 16, 32, 64, 128, 256, 512 }; // bit flips at signal regions, bit values
1434 double bs[10] = { 0. }; // bit flips at signal edges
1435 int bac[10] = { 0 }, ba0[10] = { 0 }, ba1[10] = { 0 }, bas[10] = { 0 }; // for stuck bits between first_non0, last_non0 using old technique
1436 int f; // flip-flop of non-zero contents of previous bin
1437 int cm, cp; // previous and next bin contents
1438 double prob; // probability
1439 int zero_limit, saturation_limit;
1440 zero_limit = gain ? m_zeroLimitHG : 0;
1441 saturation_limit = gain ? m_saturationLimitHG : m_i_ADCmax;
1442 cp = hist->GetBinContent(1);
1443 f = !!cp;
1444 if (f) first_non0 = 0;
1445 for (last_non0 = m_i_ADCmax; last_non0 >= 0; last_non0--)
1446 if (hist->GetBinContent(last_non0 + 1) > 0) break;
1447 if (last_non0 < 0) // empty histogram
1448 return 0;
1449 for (i = 1; i <= last_non0; i++) {
1450 cm = cc;
1451 c = cc = cp;
1452 cp = (i < m_i_ADCmax) ? hist->GetBinContent(i + 2) : 0;
1453 if (first_non0 < 0) {
1454 if (cc > 0)
1455 first_non0 = i;
1456 else
1457 continue;
1458 }
1459 for (b = 0; b < 10; b++) { // count bits flipped at this position
1460 if ((i & bv[b]) != ((i - 1) & bv[b])) {
1461 bac[b]++;
1462 if (bas[b] > 0) {
1463 if ((i - 1) & bv[b])
1464 ba1[b]++;
1465 else
1466 ba0[b]++;
1467 bas[b] = 0;
1468 }
1469 }
1470 if (cc > 0) bas[b]++;
1471 }
1472 prob = 1.;
1473 cont = sqrt((double) (c + cm));
1474 if ((c == 0 || cm == 0) && (c > 0 || cm > 0)) {
1475 prob = erf((cont - sqrt(cont)) / sqrt(2 * cont));
1476 }
1477 /* bin content is lower than usual */
1478 if (c > 0 && c < 0.25 * (cp + cm - std::sqrt((double) cp) - std::sqrt((double) cm))) {
1479 prob = 1. - (double) c / (0.25 * (cp + cm - std::sqrt((double) cp) - std::sqrt((double) cm)));
1480 c = 0;
1481 }
1482 if (c || f) // count bit flips (bc[b]) in bins with non-zero contents or adjacent bins
1483 for (b = 0; b < 10; b++) {
1484 if ((i & bv[b]) != ((i - 1) & bv[b]))
1485 bc[b]++;
1486 else
1487 break; // a higher bit can flip only of the lower one does
1488 }
1489 if (f > 0 && c > 0) continue;
1490 if (f == 0 && c == 0) continue;
1491 if (c || f) // bin content flips (zero - non-zero)
1492 for (b = 0; b < 10; b++) { // count bits flipped at this position
1493 if ((i & bv[b]) != ((i - 1) & bv[b]))
1494 bs[b] += (i & bv[b]) ? (c ? prob : -prob) : (c ? -prob : prob);
1495 else
1496 break;
1497 }
1498 f = !f;
1499 }
1500
1501 for (b = 0; b < 10; b++) // count bits flipped at this position
1502 if ((i & bv[b]) != ((i - 1) & bv[b])) {
1503 bac[b]++;
1504 if (bas[b] > 0) {
1505 if ((i - 1) & bv[b])
1506 ba1[b]++;
1507 else
1508 ba0[b]++;
1509 }
1510 }
1511
1512 int is_stack = 0;
1513 double sb_prob[4] = {0., 0., 0., 0.};
1514 static const int sb_map[10] = {0, 1, 1, 2, 2, 3, 3, 3, 3, 3}; // mapping of SB to histogram rows
1515 for (b = 0; b < 10; b++) {
1516 if ((ba0[b] == 0 || ba1[b] == 0) && bac[b] > 2 && (ba0[b] + ba1[b] >= bac[b] / 2 || ba0[b] + ba1[b] > 2)) {
1517 is_stack = 1;
1518 if (outhist != NULL) {
1519 sb_prob[sb_map[b]] = 1.;
1520 }
1521 if (stuck_probs != NULL)
1522 stuck_probs[b] = ba0[b] == 0 ? 100u : 200u;
1523 continue;
1524 }
1525 double bs1 = std::fabs(bs[b]) - sqrt(std::fabs(bs[b]));
1526 if (bs1 < 0.) bs1 = 0.;
1527 if ((bs1 > 0.5 * bc[b]) || (bc[b] > 7 && bs1 * 3 > bc[b])) is_stack = 1;
1528 if (outhist != NULL && bc[b] > 0) {
1529 // if (sb_prob[sb_map[b]] < 100. * bs1 / bc[b]) sb_prob[sb_map[b]] = 100. * bs1 / bc[b];
1530 sb_prob[sb_map[b]] = 1. - (1. - sb_prob[sb_map[b]]) * (1. - 1. * bs1 / bc[b]); // prod of probs. of not-stuck
1531 }
1532 if (stuck_probs != NULL)
1533 {
1534 stuck_probs[b] = (uint8_t) (100. * bs1 / bc[b]);
1535 if (bs[b] > 0)
1536 {
1537 stuck_probs[b] += 100u;
1538 if (stuck_probs[b] == 100u)
1539 stuck_probs[b] = 0u;
1540 }
1541 }
1542 }
1543 if ((first_non0 >= (m_i_ADCmax + 1) / 2 && first_non0 < m_i_ADCmax)
1544 || (last_non0 == (m_i_ADCmax + 1) / 2 - 1 && hist->GetBinContent(last_non0) > 3)) {
1545 is_stack = 1;
1546
1547 sb_prob[3] = 1.;
1548 if (stuck_probs != NULL)
1549 stuck_probs[9] = first_non0 >= 512 ? 200u : 100u;
1550 }
1551 if (outhist != NULL) {
1552 outhist->Fill((double) ch, 0., 100. * sb_prob[0]);
1553 outhist->Fill((double) ch, 1., 100. * sb_prob[1]);
1554 outhist->Fill((double) ch, 2., 100. * sb_prob[2]);
1555 outhist->Fill((double) ch, 3., 100. * sb_prob[3]);
1556
1557 if (first_non0 >= saturation_limit)
1558 outhist->Fill((double) ch, 5., 100.);
1559 else if (last_non0 >= saturation_limit) {
1560 double frac;
1561 int saturation_entries = 0;
1562 int maxADC = m_i_ADCmax + 1;
1563 for (i = saturation_limit + 1; i <= maxADC; ++i)
1564 saturation_entries += hist->GetBinContent(i);
1565 frac = 100. * (double) saturation_entries / hist->GetEntries();
1566 if (frac > 0. && frac < 1.) frac = 1.;
1567 if (frac > 99. && frac < 100.) frac = 99.;
1568 outhist->Fill((double) ch, 5., frac);
1569 }
1570 if (last_non0 <= zero_limit)
1571 outhist->Fill((double) ch, 4., 100.);
1572 else if (first_non0 <= zero_limit) {
1573 double frac;
1574 int zero_entries = 0;
1575 for (i = 1; i <= zero_limit + 1; ++i)
1576 zero_entries += hist->GetBinContent(i);
1577 frac = 100. * (double) zero_entries / hist->GetEntries();
1578 if (frac > 0. && frac < 1.) frac = 1.;
1579 if (frac > 99. && frac < 100.) frac = 99.;
1580 outhist->Fill((double) ch, 4., frac);
1581 }
1582 double entries, empty_cut, full_cut;
1583 int mod32empty = 0;
1584 int mod32full = 0;
1585 bool enough = false;
1586 entries = modhist -> GetEntries();
1587 empty_cut = entries < 4096. ? entries / (32 * 4) : 32.;
1588 full_cut = entries / 32.;
1589 if (full_cut > 126.)
1590 full_cut = 126.;
1591 if (full_cut - empty_cut < 1.)
1592 full_cut = empty_cut + 1.;
1593
1594 for (i = 1; i <= 32; ++i) {
1595 if (modhist->GetBinContent(i) <= empty_cut)
1596 ++mod32empty;
1597 else if (modhist->GetBinContent(i) > full_cut)
1598 ++mod32full;
1599 if (modhist->GetBinContent(i) > 64)
1600 enough = true;
1601 }
1602 if (mod32empty + mod32full == 32 || enough) //skip in the case of small number of events
1603 if (mod32empty != 0 && mod32full != 0) //some but not all n%32 positions have events
1604 outhist->Fill((double) ch, 6., mod32empty);
1605
1606 outhist->SetMaximum(100.);
1607 }
1608 return is_stack;
1609}
1610
1611
1613 uint32_t crc32, uint32_t crcMask, int headsize, int ros, int drawer) {
1623
1624 // protection against non-exisitng histograms
1625 if (m_data->m_hist0[ros][drawer].size() == 0) return;
1626
1627 //LF CRC check
1628 //std::string ros_name[2] = {"LBA","LBC"};
1629 //std::cout << "Module= " << ros_name[ros-1] << drawer+1 << " CRCmask= 0x" << std::setfill('0') << std::setw(8) << std::hex << crcMask << std::setfill(' ') << std::dec <<"\n";
1630
1631 //END CRC check
1632 // array to adjust order of DMU's of EB* FE CRC status in DMU chip mask word
1633 // this is needed because FE and ROD fill DMU mask word in different order
1634 int chFE_ext[16] = { 0, 1, 2, 3, 4, 5, 6, 7, 12, 13, 8, 9, 10, 11, 14, 15 };
1635 int chFE_sp[16] = { 11, 0, 1, 2, 3, 4, 5, 6, 12, 13, 7, 8, 9, 10, 14, 15 };
1636 // global CRC
1637 // uint32_t crc32 = (*collItr)->getFragCRC();
1638 uint32_t crc0 = crc32 & 0xFFFF;
1639 uint32_t crc1 = crc32 >> 16;
1640 if (crc32 == 0)
1641 m_data->m_hist0[ros][drawer][0]->Fill(0.0, 1.0);
1642 else if (crc0 == crc1) // two half-words should match
1643 m_data->m_hist0[ros][drawer][0]->Fill(1.0, 1.0);
1644 else
1645 m_data->m_hist0[ros][drawer][0]->Fill(2.0, 1.0);
1646
1647 // CRC per TileDMU.For mono gain we have it from ROD and FE.
1648 // In bi gain, it is just a placeholder with global CRC information.
1649
1650 if (headsize < 16) { // How to handle EB? All bits are filled anyway or not?
1651 }
1652
1653 if (dqStatus->calibMode() == 1) {
1654 for (int ch = 0; ch < headsize; ++ch) {
1655 if (crc32 == 0)
1656 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(0.0, 1.0);
1657 else if (crc0 == crc1)
1658 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(1.0, 1.0);
1659 else
1660 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(2.0, 1.0);
1661
1662 }
1663 } else {
1669 if (crc32 == 0) { //std::cout << "Global crc is zero\n";
1670 for (int ch = 0; ch < headsize; ++ch) {
1671 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(0.0, 1.0);
1672 }
1673 } else if (crcMask == 0xFFFFFFFF) {
1674 for (int ch = 0; ch < headsize; ++ch) {
1675 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(1.0, 1.0);
1676 }
1677 } else {
1678 uint32_t fe_crc = crcMask & 0xFFFF;
1679 uint32_t rod_crc = crcMask >> 16;
1680 for (int ch = 0; ch < headsize; ++ch) {
1681 int flag = 0;
1682 if (ros == 1 || ros == 2) // LB* use simple FECRC map
1683 {
1684 if ((fe_crc >> ch & 0x1) == 0x0) flag++;
1685 } else if (ros == 3 || ros == 4) // EB* use different FECRC map
1686 {
1687 if ((ros == 3 && drawer == 14) || (ros == 4 && drawer == 17)) // Special Check for EBA15,EBC18
1688 {
1689 if ((fe_crc >> chFE_sp[ch] & 0x1) == 0x0) flag++;
1690 } else {
1691 if ((fe_crc >> chFE_ext[ch] & 0x1) == 0x0) flag++;
1692 }
1693 }
1694 if ((rod_crc >> ch & 0x1) == 0x0) flag += 2;
1695
1696 //std::cout << "flag for chan "<< ch << "=" << flag << "\n";
1697
1698 switch (flag) {
1699 case 0: //TileDMU is fine
1700 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(1.0, 1.0);
1701 break;
1702 case 1: //fe error only
1703 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(2.0, 1.0);
1704 break;
1705 case 2: // rod errors
1706 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(3.0, 1.0);
1707 break;
1708 default: // fe+rod errors
1709 m_data->m_hist_DMUerr[ros][drawer][ch][0][1]->Fill(4.0, 1.0);
1710 } // end switch case
1711
1712 } // end loop on chips
1713
1714 } // end if on crcMask
1715
1716 } // end if on runtype
1717}
1718
1737//}
1738//New Method to check for header errors. Checks per channel and returns True for data corrupted, or False for data not corrupted.
1739
1740/*-------------------------------------------------------*/
1741bool TileDigitsMonTool::DMUheaderCheck(std::vector<uint32_t>* headerVec, int ros, int drawer, int gain, int dmu)
1742/*-------------------------------------------------------*/
1743{
1744 bool err = false;
1745
1746 if (DMUheaderFormatCheck((*headerVec)[dmu])) {
1747 m_data->m_hist2[ros][drawer][gain][0]->Fill(dmu + 0.5, 1., 1.);
1748 err = true;
1749 return err;
1750 }
1751 if (DMUheaderParityCheck((*headerVec)[dmu])) {
1752 m_data->m_hist2[ros][drawer][gain][0]->Fill(dmu + 0.5, 2., 1.);
1753 err = true;
1754 return err;
1755 }
1756 if (((*headerVec)[dmu] >> 25) & 0x1) {
1757 //Memory Parity Error
1758 m_data->m_hist2[ros][drawer][gain][0]->Fill(dmu + 0.5, 3., 1.);
1759 err = true;
1760 }
1761 if (((*headerVec)[dmu] >> 24) & 0x1) {
1762 //Single Strobe Error
1763 m_data->m_hist2[ros][drawer][gain][0]->Fill(dmu + 0.5, 4., 1.);
1764 err = true;
1765 }
1766 if (((*headerVec)[dmu] >> 23) & 0x1) {
1767 //Double Strobe Error
1768 m_data->m_hist2[ros][drawer][gain][0]->Fill(dmu + 0.5, 5., 1.);
1769 err = true;
1770 }
1771
1772 if (!err) m_data->m_hist2[ros][drawer][gain][0]->Fill(dmu + 0.5, 0., 1.);
1773
1774 int fragId = m_tileHWID->frag(ros, drawer);
1775 if (std::binary_search(m_fragIDsToIgnoreDMUerrors.begin(), m_fragIDsToIgnoreDMUerrors.end(), fragId)) {
1776 err = false;
1777 }
1778
1779 return err;
1780}
1781
1782StatusCode TileDigitsMonTool::RODCRCcalc(const TileDQstatus* dqStatus) {
1793
1794 const TileRawChannelContainer* rawChannelContainer;
1795 CHECK(evtStore()->retrieve(rawChannelContainer, m_contNameDSP));
1796
1797 for (const TileRawChannelCollection* rawChannelCollection : *rawChannelContainer) {
1798
1799 if (rawChannelCollection->empty()) continue;
1800
1801 HWIdentifier hwid = rawChannelCollection->front()->adc_HWID(); //take the first channel in the drawer
1802 int ros = m_tileHWID->ros(hwid); //take the ros and drawer from the first channel
1803 int drawer = m_tileHWID->drawer(hwid);
1804
1805 uint32_t crc32 = rawChannelCollection->getFragGlobalCRC() & 1;
1806
1807 if (crc32 == 0) {
1808 crc32 = 0xFFFFFFFF;
1809 } else { // means OK CRC match
1810 crc32 = 0xFFFF;
1811 } //means NOT OK, CRC mismatch
1812
1813 uint32_t CRCmask = rawChannelCollection->getFragRODChipMask();
1814 CRCmask = CRCmask << 16; // ROD mask is stored in the 16 most significant bits ce
1815 CRCmask += rawChannelCollection->getFragFEChipMask();
1816
1817 CRCcheck(dqStatus, crc32, CRCmask, 16, ros, drawer); //reuse the same funcion used for digits
1818
1819 }
1820
1821 return StatusCode::SUCCESS;
1822}
1823
1824/*---------------------------------------------------------*/
1825void TileDigitsMonTool::shiftHisto(TH1S *hist, int ros, int drawer, int ch, int gain)
1826/*---------------------------------------------------------*/
1827{
1828
1829 hist->GetXaxis()->SetRange(5, 100);
1830
1831 TSpectrum s(8, 1.);
1832
1833 s.Search(hist, 2, "goff");
1834 int shift = int(s.GetPositionX()[0]);
1835 if (shift > 0) {
1836 int factor = m_is12bit ? 4 : 1;
1837 int xmax = std::max(1025 * factor, 846 * factor + shift);
1838 int xmin = std::max(1, shift + 50 * factor);
1839 for (int bin = xmin; bin < xmax; ++bin) {
1840 double c = hist->GetBinContent(bin);
1841 if (c > 0) m_data->m_shifted_hist[ros][drawer][ch][gain]->SetBinContent(bin - shift - 50 * factor, c);
1842 }
1843 }
1844}
1845
1846/*---------------------------------------------------------*/
1847void TileDigitsMonTool::statTestHistos(int ros, int drawer, int gain)
1848/*---------------------------------------------------------*/
1849
1850{
1851 std::vector<TH1S*> refbld;
1852 std::vector<TH1S*> newrefbld;
1853 TH1F *ref = new TH1F("ref", "ref", m_shiftnbins, 0., m_shiftnbins);
1854 ref->SetDirectory(0);
1855 TH1F *ref1 = new TH1F("ref1", "ref1", m_shiftnbins, 0., m_shiftnbins);
1856 ref1->SetDirectory(0);
1857
1858 for (int i = 0; i < 48; i++) {
1859
1860 TH1S *h = m_data->m_shifted_hist[ros][drawer][i][gain];
1861 float integ = h->Integral(200, 600);
1862 if (integ > 0) {
1863 refbld.push_back(h);
1864 }
1865
1866 }
1867
1868 int entr = refbld.size();
1869 const double inv_entr = (entr > 0 ? 1. / static_cast<double>(entr) : 1);
1870 for (int i = 0; i < entr; i++) {
1871 TH1S *h = refbld.at(i);
1872 ref->Add(h, inv_entr);
1873 }
1874
1875 // refbld.push_back(ref);
1876 TH1S *obj = 0;
1877
1878 for (int i = 0; i < entr; i++) {
1879 obj = refbld.at(i);
1880 float kol = obj->KolmogorovTest(ref);
1881 //std::cout<<"The value of the Kolmogorov test for "<< ros << drawer << gain << " pmt:" << i << " is "<< kol << std::endl;
1882 if (kol > 0.5) {
1883 newrefbld.push_back(obj);
1884 }
1885 }
1886
1887 int ent = newrefbld.size();
1888 const double inv_ent = (ent > 0 ? 1. / static_cast<double>(ent) : 1);
1889 for (int i = 0; i < ent; i++) {
1890 obj = newrefbld.at(i);
1891 ref1->Add(obj, inv_ent);
1892 }
1893
1894 // newrefbld->Add(ref1);
1895
1896 if (ent > 6)
1897 m_data->m_shifted_hist[ros][drawer][48][gain]->Add(ref1, 1.);
1898 else
1899 m_data->m_shifted_hist[ros][drawer][48][gain]->Add(ref, 1.);
1900
1901 delete ref;
1902 delete ref1;
1903
1904}
1905
1906
1908
1909 tree->Branch("StuckBitsProb", m_data->m_stuck_probs, "StuckBitsProb[5][64][48][2][10]/b");
1910}
const boost::regex ref(r_ef)
#define endmsg
#define ATH_MSG_ERROR(x)
#define ATH_MSG_INFO(x)
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
#define ATH_MSG_DEBUG(x)
double charge(const T &p)
Definition AtlasPID.h:997
#define CHECK(...)
Evaluate an expression and check for errors.
Handle class for reading from StoreGate.
TGraphErrors * GetEntries(TH2F *histo)
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
ServiceHandle< StoreGateSvc > & evtStore()
const ServiceHandle< StoreGateSvc > & detStore() const
bool msgLvl(const MSG::Level lvl) const
Header file for AthHistogramAlgorithm.
static const unsigned int MAX_ROS
Number of ROSs.
static std::string getDrawerString(unsigned int ros, unsigned int drawer)
Return the drawer name, e.g.
static const unsigned int MAX_DRAWER
Number of drawers in ROS 1-4.
static unsigned int getDrawerIdx(unsigned int ros, unsigned int drawer)
Returns a drawer hash.
static const unsigned int MAX_CHAN
Number of channels in drawer.
Class that holds Data Quality fragment information and provides functions to extract the data quality...
uint32_t calibMode() const
Calibration mode.
const uint32_t * cispar() const
CIS parameters.
bool DMUheaderFormatCheck(uint32_t header)
Function to check that the DMU header format is correct bit_31 of the DMU header must be 1 and bit_17...
int define_palette(int ncolors, int *colors=NULL)
const TileInfo * m_tileInfo
int stuckBits_Amp2(TH1S *hist, TH1C *modhist, int adc, TH2C *outhist=NULL, int ch=0, uint8_t *stuck_probs=NULL)
SG::ReadHandleKey< TileDQstatus > m_DQstatusKey
virtual StatusCode bookHists() override
Calls bookHists( true, true, true ) and initializes lumiBlock and run numbers.
void CRCcheck(const TileDQstatus *dqStatus, uint32_t crc32, uint32_t crcMask, int headsize, int ros, int drawer)
Method to check global CRC and DMU CRC.
void shiftHisto(TH1S *hist, int ros, int drawer, int ch, int gain)
virtual const uint8_t * stuckBitProb(int ros, int module, int channel, int gain) const override
virtual StatusCode initialize() override
void drawHists(int ros, int drawer, const std::string &moduleName)
virtual void saveStuckBitsProbabilities(TTree *tree) override
int stuckBits(TH1S *hist, int adc)
Method to find stuckbits in Read-Out ADC channel.
bool DMUheaderParityCheck(uint32_t header)
Function to check that the DMU header parity is correct Parity of the DMU header should be odd Return...
StatusCode RODCRCcalc(const TileDQstatus *dqStatus)
virtual StatusCode checkHists(bool fromFinalize) override
This implementation does nothing; equivalent functionality may be provided by procHists(....
TileDigitsMonTool(const std::string &type, const std::string &name, const IInterface *parent)
std::unique_ptr< Data > m_data
static const int m_NCont
bool DMUheaderCheck(std::vector< uint32_t > *headerVec, int ros, int drawer, int gain, int dmu)
DMU header format as defined in http://www.sysf.physto.se/~klere/tile-dmu/header.html Bits: 1pllllesd...
ToolHandle< TileCondToolNoiseSample > m_tileToolNoiseSample
tool which provided noise values
int stuckBits_Amp(TH1S *hist, int adc)
A crude method to check Read-Out ADC channel stuckbits.
const uint32_t * m_cispar
void statTestHistos(int ros, int drawer, int gain)
std::string m_digitsContainerName
virtual StatusCode fillHists() override
Calls fillHists( bool, bool, bool ); if an eventBlock,lumiBlock, or run has turned over,...
virtual StatusCode finalHists() override
Calls procHists( true, true, true ).
bool isDisconnected(int ros, int drawer, int ch)
TH2F * book2F(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, int ny, double ymin, double ymax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
std::vector< int > m_fragIDsToIgnoreDMUerrors
TProfile * bookProfile(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
virtual StatusCode initialize() override
TH1S * book1S(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
TH1F * book1F(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
const TileHWID * m_tileHWID
TH1D * book1D(const std::string &nam, const std::string &tit, int nx, double xmin, double xmax)
Implicit version of book1D.
TilePaterMonTool(const std::string &type, const std::string &name, const IInterface *parent)
TH1C * book1C(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
TH2C * book2C(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, int ny, double ymin, double ymax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
TH1I * book1I(const std::string &dir, const std::string &nam, const std::string &tit, int nx, double xmin, double xmax, Interval_t interval=run, MgmtAttr_t attribute=ATTRIB_MANAGED, const std::string &trigChain="", const std::string &mergeAlgo="")
double xmax
Definition listroot.cxx:61
double entries
Definition listroot.cxx:49
double xmin
Definition listroot.cxx:60
SG::ReadCondHandle< T > makeHandle(const SG::ReadCondHandleKey< T > &key, const EventContext &ctx=Gaudi::Hive::currentContext())
MsgStream & msg
Definition testRead.cxx:32
TChain * tree