Avoid message E-TClonesArray::At: during digitization due to try to access non existi...
[u/mrichter/AliRoot.git] / EMCAL / AliCaloCalibPedestal.cxx
CommitLineData
a235e2bc 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
59ed788d 15//* $Id$ */
a235e2bc 16
17//________________________________________________________________________
18//
19// A help class for monitoring and calibration tools: MOOD, AMORE etc.,
20// It can be created and used a la (ctor):
21/*
22 //Create the object for making the histograms
23 fPedestals = new AliCaloCalibPedestal( fDetType );
24 // AliCaloCalibPedestal knows how many modules we have for PHOS or EMCAL
25 fNumModules = fPedestals->GetModules();
26*/
27// fed an event:
28// fPedestals->ProcessEvent(fCaloRawStream);
29// asked to draw histograms:
30// fPedestals->GetDeadMap(i)->Draw("col");
31// or
32// fPedestals->GetPeakProfileHighGainRatio((i < fNumModules) ? i : fVisibleModule)->Draw("colz");
33// etc.
34// The pseudo-code examples above were from the first implementation in MOOD (summer 2007).
35//________________________________________________________________________
36
59ed788d 37//#include "TCanvas.h"
356c3e0c 38#include "TH1.h"
59ed788d 39#include "TF1.h"
356c3e0c 40#include "TFile.h"
356c3e0c 41#include <fstream>
b07ee441 42#include <sstream>
356c3e0c 43#include <iostream>
c490d8e5 44#include <stdexcept>
3dba9483 45#include <cmath>
356c3e0c 46
30aa89b0 47#include "AliRawReader.h"
32cd4c24 48#include "AliCaloRawStreamV3.h"
a235e2bc 49
356c3e0c 50//The include file
51#include "AliCaloCalibPedestal.h"
52
53ClassImp(AliCaloCalibPedestal)
54
55using namespace std;
56
57// ctor; initialize everything in order to avoid compiler warnings
58AliCaloCalibPedestal::AliCaloCalibPedestal(kDetType detectorType) :
59 TObject(),
60 fPedestalLowGain(),
61 fPedestalHighGain(),
29f94584 62 fPedestalLEDRefLowGain(),
63 fPedestalLEDRefHighGain(),
356c3e0c 64 fPeakMinusPedLowGain(),
65 fPeakMinusPedHighGain(),
59ed788d 66 fPeakMinusPedHighGainHisto(),
356c3e0c 67 fPedestalLowGainDiff(),
68 fPedestalHighGainDiff(),
29f94584 69 fPedestalLEDRefLowGainDiff(),
70 fPedestalLEDRefHighGainDiff(),
356c3e0c 71 fPeakMinusPedLowGainDiff(),
72 fPeakMinusPedHighGainDiff(),
73 fPedestalLowGainRatio(),
74 fPedestalHighGainRatio(),
29f94584 75 fPedestalLEDRefLowGainRatio(),
76 fPedestalLEDRefHighGainRatio(),
356c3e0c 77 fPeakMinusPedLowGainRatio(),
78 fPeakMinusPedHighGainRatio(),
79 fDeadMap(),
419341ea 80 fNEvents(0),
81 fNChanFills(0),
356c3e0c 82 fDeadTowers(0),
83 fNewDeadTowers(0),
84 fResurrectedTowers(0),
85 fReference(0),
86 fDetType(kNone),
87 fColumns(0),
88 fRows(0),
29f94584 89 fLEDRefs(0),
356c3e0c 90 fModules(0),
18db89b7 91 fRowMin(0),
92 fRowMax(0),
93 fRowMultiplier(0),
f4fc542c 94 fCaloString(),
95 fMapping(NULL),
3dba9483 96 fRunNumber(-1),
97 fSelectPedestalSamples(kTRUE),
98 fFirstPedestalSample(0),
59ed788d 99 fLastPedestalSample(15),
100 fDeadThreshold(5),
101 fWarningThreshold(50),
102 fWarningFraction(0.002),
103 fHotSigma(5)
356c3e0c 104{
105 //Default constructor. First we set the detector-type related constants.
106 if (detectorType == kPhos) {
107 fColumns = fgkPhosCols;
108 fRows = fgkPhosRows;
29f94584 109 fLEDRefs = fgkPhosLEDRefs;
356c3e0c 110 fModules = fgkPhosModules;
f4fc542c 111 fCaloString = "PHOS";
18db89b7 112 fRowMin = -1*fRows;
113 fRowMax = 0;
114 fRowMultiplier = -1;
356c3e0c 115 }
116 else {
117 //We'll just trust the enum to keep everything in line, so that if detectorType
118 //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
119 //case, if someone intentionally gives another number
bd3cd9c2 120 fColumns = AliEMCALGeoParams::fgkEMCALCols;
121 fRows = AliEMCALGeoParams::fgkEMCALRows;
29f94584 122 fLEDRefs = AliEMCALGeoParams::fgkEMCALLEDRefs;
bd3cd9c2 123 fModules = AliEMCALGeoParams::fgkEMCALModules;
f4fc542c 124 fCaloString = "EMCAL";
18db89b7 125 fRowMin = 0;
126 fRowMax = fRows;
127 fRowMultiplier = 1;
356c3e0c 128 }
129 fDetType = detectorType;
23ca956b 130
131 // ValidateProfiles(); // not to be done in ctor; info from Axel N.
132}
133
134//_____________________________________________________________________
135void AliCaloCalibPedestal::ValidateProfiles()
136{
137 //Make sure the basic histos exist
138 if (!fPedestalLowGain.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
139 //the same time
140
356c3e0c 141 //Then, loop for the requested number of modules
142 TString title, name;
143 for (int i = 0; i < fModules; i++) {
144 //Pedestals, low gain
145 name = "hPedlowgain";
146 name += i;
147 title = "Pedestals, low gain, module ";
148 title += i;
149 fPedestalLowGain.Add(new TProfile2D(name, title,
150 fColumns, 0.0, fColumns,
18db89b7 151 fRows, fRowMin, fRowMax,"s"));
356c3e0c 152
153 //Pedestals, high gain
154 name = "hPedhighgain";
155 name += i;
156 title = "Pedestals, high gain, module ";
157 title += i;
158 fPedestalHighGain.Add(new TProfile2D(name, title,
159 fColumns, 0.0, fColumns,
18db89b7 160 fRows, fRowMin, fRowMax,"s"));
29f94584 161
162 //LED Ref/Mon pedestals, low gain
163 name = "hPedestalLEDReflowgain";
18db89b7 164 name += i;
29f94584 165 title = "Pedestal LEDRef, low gain, module ";
18db89b7 166 title += i;
29f94584 167 fPedestalLEDRefLowGain.Add(new TProfile(name, title,
168 fLEDRefs, 0.0, fLEDRefs, "s"));
169
170 //LED Ref/Mon pedestals, high gain
171 name = "hPedestalLEDRefhighgain";
18db89b7 172 name += i;
29f94584 173 title = "Pedestal LEDRef, high gain, module ";
18db89b7 174 title += i;
29f94584 175 fPedestalLEDRefHighGain.Add(new TProfile(name, title,
176 fLEDRefs, 0.0, fLEDRefs, "s"));
356c3e0c 177
178 //Peak-Pedestals, low gain
179 name = "hPeakMinusPedlowgain";
180 name += i;
181 title = "Peak-Pedestal, low gain, module ";
182 title += i;
183 fPeakMinusPedLowGain.Add(new TProfile2D(name, title,
184 fColumns, 0.0, fColumns,
18db89b7 185 fRows, fRowMin, fRowMax,"s"));
356c3e0c 186
187 //Peak-Pedestals, high gain
188 name = "hPeakMinusPedhighgain";
189 name += i;
190 title = "Peak-Pedestal, high gain, module ";
191 title += i;
192 fPeakMinusPedHighGain.Add(new TProfile2D(name, title,
193 fColumns, 0.0, fColumns,
18db89b7 194 fRows, fRowMin, fRowMax,"s"));
59ed788d 195
196 //Peak-Pedestals, high gain - TH2F histo
197 name = "hPeakMinusPedhighgainHisto";
198 name += i;
199 title = "Peak-Pedestal, high gain, module ";
200 title += i;
201 fPeakMinusPedHighGainHisto.Add(new TH2F(name, title,
202 fColumns*fRows, 0.0, fColumns*fRows,
203 100, 0, 1000));
204
356c3e0c 205 name = "hDeadMap";
206 name += i;
207 title = "Dead map, module ";
208 title += i;
209 fDeadMap.Add(new TH2D(name, title, fColumns, 0.0, fColumns,
18db89b7 210 fRows, fRowMin, fRowMax));
356c3e0c 211
212 }//end for nModules create the histograms
23ca956b 213
214 CompressAndSetOwner();
215}
216
217//_____________________________________________________________________
218void AliCaloCalibPedestal::CompressAndSetOwner()
219{
356c3e0c 220 //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
221 fPedestalLowGain.Compress();
222 fPedestalHighGain.Compress();
29f94584 223 fPedestalLEDRefLowGain.Compress();
224 fPedestalLEDRefHighGain.Compress();
356c3e0c 225 fPeakMinusPedLowGain.Compress();
226 fPeakMinusPedHighGain.Compress();
59ed788d 227 fPeakMinusPedHighGainHisto.Compress();
356c3e0c 228 fDeadMap.Compress();
29f94584 229
61e4e2e3 230 // set owner ship for everyone
231 fPedestalLowGain.SetOwner(kTRUE);
232 fPedestalHighGain.SetOwner(kTRUE);
233 fPedestalLEDRefLowGain.SetOwner(kTRUE);
234 fPedestalLEDRefHighGain.SetOwner(kTRUE);
235 fPeakMinusPedLowGain.SetOwner(kTRUE);
236 fPeakMinusPedHighGain.SetOwner(kTRUE);
237 fPeakMinusPedHighGainHisto.SetOwner(kTRUE);
238 fPedestalLowGainDiff.SetOwner(kTRUE);
239 fPedestalHighGainDiff.SetOwner(kTRUE);
240 fPedestalLEDRefLowGainDiff.SetOwner(kTRUE);
241 fPedestalLEDRefHighGainDiff.SetOwner(kTRUE);
242 fPeakMinusPedLowGainDiff.SetOwner(kTRUE);
243 fPeakMinusPedHighGainDiff.SetOwner(kTRUE);
244 fPedestalLowGainRatio.SetOwner(kTRUE);
245 fPedestalHighGainRatio.SetOwner(kTRUE);
246 fPedestalLEDRefLowGainRatio.SetOwner(kTRUE);
247 fPedestalLEDRefHighGainRatio.SetOwner(kTRUE);
248 fPeakMinusPedLowGainRatio.SetOwner(kTRUE);
249 fPeakMinusPedHighGainRatio.SetOwner(kTRUE);
250 fDeadMap.SetOwner(kTRUE);
356c3e0c 251}
252
253// dtor
254//_____________________________________________________________________
255AliCaloCalibPedestal::~AliCaloCalibPedestal()
256{
65bec413 257 //dtor
65bec413 258
f62c044a 259 if (fReference) delete fReference;//Delete the reference object, if it has been loaded
65bec413 260
65bec413 261 // delete also TObjArray's
f62c044a 262 fPedestalLowGain.Delete();
263 fPedestalHighGain.Delete();
264 fPedestalLEDRefLowGain.Delete();
265 fPedestalLEDRefHighGain.Delete();
266 fPeakMinusPedLowGain.Delete();
267 fPeakMinusPedHighGain.Delete();
268 fPeakMinusPedHighGainHisto.Delete();
269 fPedestalLowGainDiff.Delete();
270 fPedestalHighGainDiff.Delete();
271 fPedestalLEDRefLowGainDiff.Delete();
272 fPedestalLEDRefHighGainDiff.Delete();
273 fPeakMinusPedLowGainDiff.Delete();
274 fPeakMinusPedHighGainDiff.Delete();
275 fPedestalLowGainRatio.Delete();
276 fPedestalHighGainRatio.Delete();
277 fPedestalLEDRefLowGainRatio.Delete();
278 fPedestalLEDRefHighGainRatio.Delete();
279 fPeakMinusPedLowGainRatio.Delete();
280 fPeakMinusPedHighGainRatio.Delete();
281 fDeadMap.Delete();
65bec413 282
356c3e0c 283}
284
285// copy ctor
286//_____________________________________________________________________
23ca956b 287AliCaloCalibPedestal::AliCaloCalibPedestal(AliCaloCalibPedestal &ped) :
356c3e0c 288 TObject(ped),
289 fPedestalLowGain(),
290 fPedestalHighGain(),
29f94584 291 fPedestalLEDRefLowGain(),
292 fPedestalLEDRefHighGain(),
356c3e0c 293 fPeakMinusPedLowGain(),
294 fPeakMinusPedHighGain(),
59ed788d 295 fPeakMinusPedHighGainHisto(),
356c3e0c 296 fPedestalLowGainDiff(),
297 fPedestalHighGainDiff(),
29f94584 298 fPedestalLEDRefLowGainDiff(),
299 fPedestalLEDRefHighGainDiff(),
356c3e0c 300 fPeakMinusPedLowGainDiff(),
301 fPeakMinusPedHighGainDiff(),
302 fPedestalLowGainRatio(),
303 fPedestalHighGainRatio(),
29f94584 304 fPedestalLEDRefLowGainRatio(),
305 fPedestalLEDRefHighGainRatio(),
356c3e0c 306 fPeakMinusPedLowGainRatio(),
307 fPeakMinusPedHighGainRatio(),
308 fDeadMap(),
419341ea 309 fNEvents(ped.GetNEvents()),
310 fNChanFills(ped.GetNChanFills()),
356c3e0c 311 fDeadTowers(ped.GetDeadTowerCount()),
312 fNewDeadTowers(ped.GetDeadTowerNew()),
313 fResurrectedTowers(ped.GetDeadTowerResurrected()),
314 fReference( 0 ), //! note that we do not try to copy the reference info here
315 fDetType(ped.GetDetectorType()),
316 fColumns(ped.GetColumns()),
317 fRows(ped.GetRows()),
29f94584 318 fLEDRefs(ped.GetLEDRefs()),
356c3e0c 319 fModules(ped.GetModules()),
18db89b7 320 fRowMin(ped.GetRowMin()),
321 fRowMax(ped.GetRowMax()),
322 fRowMultiplier(ped.GetRowMultiplier()),
f4fc542c 323 fCaloString(ped.GetCaloString()),
324 fMapping(NULL), //! note that we are not copying the map info
3dba9483 325 fRunNumber(ped.GetRunNumber()),
326 fSelectPedestalSamples(ped.GetSelectPedestalSamples()),
327 fFirstPedestalSample(ped.GetFirstPedestalSample()),
59ed788d 328 fLastPedestalSample(ped.GetLastPedestalSample()),
329 fDeadThreshold(ped.GetDeadThreshold()),
330 fWarningThreshold(ped.GetWarningThreshold()),
331 fWarningFraction(ped.GetWarningFraction()),
332 fHotSigma(ped.GetHotSigma())
356c3e0c 333{
334 // Then the ObjArray ones; we add the histograms rather than trying TObjArray = assignment
335 //DS: this has not really been tested yet..
336 for (int i = 0; i < fModules; i++) {
337 fPedestalLowGain.Add( ped.GetPedProfileLowGain(i) );
338 fPedestalHighGain.Add( ped.GetPedProfileHighGain(i) );
29f94584 339 fPedestalLEDRefLowGain.Add( ped.GetPedLEDRefProfileLowGain(i) );
340 fPedestalLEDRefHighGain.Add( ped.GetPedLEDRefProfileHighGain(i) );
356c3e0c 341 fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
342 fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
59ed788d 343 fPeakMinusPedHighGainHisto.Add( ped.GetPeakHighGainHisto(i) );
356c3e0c 344
345 fDeadMap.Add( ped.GetDeadMap(i) );
346 }//end for nModules
347
23ca956b 348 CompressAndSetOwner();
356c3e0c 349}
350
351// assignment operator; use copy ctor to make life easy..
352//_____________________________________________________________________
23ca956b 353AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (AliCaloCalibPedestal &source)
356c3e0c 354{
a235e2bc 355 // assignment operator; use copy ctor
356c3e0c 356 if (&source == this) return *this;
357
358 new (this) AliCaloCalibPedestal(source);
359 return *this;
360}
361
362//_____________________________________________________________________
363void AliCaloCalibPedestal::Reset()
364{
23ca956b 365 ValidateProfiles(); // make sure histos/profiles exist
a235e2bc 366 // Reset all arrays/histograms
356c3e0c 367 for (int i = 0; i < fModules; i++) {
368 GetPedProfileLowGain(i)->Reset();
369 GetPedProfileHighGain(i)->Reset();
29f94584 370 GetPedLEDRefProfileLowGain(i)->Reset();
371 GetPedLEDRefProfileHighGain(i)->Reset();
356c3e0c 372 GetPeakProfileLowGain(i)->Reset();
373 GetPeakProfileHighGain(i)->Reset();
59ed788d 374 GetPeakHighGainHisto(i)->Reset();
356c3e0c 375 GetDeadMap(i)->Reset();
376
377 if (!fPedestalLowGainDiff.IsEmpty()) {
378 //This means that the comparison profiles have been created.
379
380 GetPedProfileLowGainDiff(i)->Reset();
381 GetPedProfileHighGainDiff(i)->Reset();
29f94584 382 GetPedLEDRefProfileLowGainDiff(i)->Reset();
383 GetPedLEDRefProfileHighGainDiff(i)->Reset();
356c3e0c 384 GetPeakProfileLowGainDiff(i)->Reset();
385 GetPeakProfileHighGainDiff(i)->Reset();
386
387 GetPedProfileLowGainRatio(i)->Reset();
388 GetPedProfileHighGainRatio(i)->Reset();
29f94584 389 GetPedLEDRefProfileLowGainRatio(i)->Reset();
390 GetPedLEDRefProfileHighGainRatio(i)->Reset();
356c3e0c 391 GetPeakProfileLowGainRatio(i)->Reset();
392 GetPeakProfileHighGainRatio(i)->Reset();
393 }
394 }
419341ea 395 fNEvents = 0;
396 fNChanFills = 0;
356c3e0c 397 fDeadTowers = 0;
398 fNewDeadTowers = 0;
399 fResurrectedTowers = 0;
400
401 //To think about: should fReference be deleted too?... let's not do it this time, at least...
402}
403
b07ee441 404// Parameter/cut handling
405//_____________________________________________________________________
406void AliCaloCalibPedestal::SetParametersFromFile(const char *parameterFile)
40164976 407{
408 // Note: this method is a bit more complicated than it really has to be
409 // - allowing for multiple entries per line, arbitrary order of the
410 // different variables etc. But I wanted to try and do this in as
411 // correct a C++ way as I could (as an exercise).
412
b07ee441 413 static const string delimitor("::");
414
415 // open, check input file
416 ifstream in( parameterFile );
417 if( !in ) {
418 printf("in AliCaloCalibPedestal::SetParametersFromFile - Using default/run_time parameters.\n");
419 return;
420 }
421
b07ee441 422
423 // read in
424 char readline[1024];
425 while ((in.rdstate() & ios::failbit) == 0 ) {
426
427 // Read into the raw char array and then construct a string
428 // to do the searching
429 in.getline(readline, 1024);
430 istringstream s(readline);
431
432 while ( ( s.rdstate() & ios::failbit ) == 0 ) {
433
ab962f7b 434 string keyValue;
435 s >> keyValue;
b07ee441 436
437 // check stream status
b63770f3 438 if( ( s.rdstate() & ios::failbit ) == ios::failbit) break;
b07ee441 439
440 // skip rest of line if comments found
ab962f7b 441 if( keyValue.substr( 0, 2 ) == "//" ) break;
b07ee441 442
ab962f7b 443 // look for "::" in keyValue pair
444 size_t position = keyValue.find( delimitor );
b07ee441 445 if( position == string::npos ) {
ab962f7b 446 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
b07ee441 447 }
448
ab962f7b 449 // split keyValue pair
450 string key( keyValue.substr( 0, position ) );
451 string value( keyValue.substr( position+delimitor.size(),
452 keyValue.size()-delimitor.size() ) );
b07ee441 453
454 // check value does not contain a new delimitor
455 if( value.find( delimitor ) != string::npos ) {
ab962f7b 456 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
b07ee441 457 }
458
459 // debug: check key value pair
460 // printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
461
462 // if the key matches with something we expect, we assign the new value
463 istringstream iss(value);
464 // the comparison strings defined at the beginning of this method
59ed788d 465 if ( (key == "fFirstPedestalSample") || (key == "fLastPedestalSample") || (key == "fDeadThreshold") || (key == "fWarningThreshold") || (key == "fWarningFraction") || (key == "fHotSigma") ) {
b07ee441 466 printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
467
468 if (key == "fFirstPedestalSample") {
469 iss >> fFirstPedestalSample;
470 }
471 else if (key == "fLastPedestalSample") {
472 iss >> fLastPedestalSample;
473 }
59ed788d 474 else if (key == "fDeadThreshold") {
475 iss >> fDeadThreshold;
476 }
477 else if (key == "fWarningThreshold") {
478 iss >> fWarningThreshold;
479 }
480 else if (key == "fWarningFraction") {
481 iss >> fWarningFraction;
482 }
483 else if (key == "fHotSigma") {
484 iss >> fHotSigma;
485 }
486
b07ee441 487 } // some match
488
489 }
490 }
491
492 in.close();
493 return;
494
495}
496
497//_____________________________________________________________________
498void AliCaloCalibPedestal::WriteParametersToFile(const char *parameterFile)
499{
40164976 500 //Write parameters in file.
501
b07ee441 502 static const string delimitor("::");
503 ofstream out( parameterFile );
504 out << "// " << parameterFile << endl;
505 out << "fFirstPedestalSample" << "::" << fFirstPedestalSample << endl;
506 out << "fLastPedestalSample" << "::" << fLastPedestalSample << endl;
59ed788d 507 out << "fDeadThreshold" << "::" << fDeadThreshold << endl;
508 out << "fWarningThreshold" << "::" << fWarningThreshold << endl;
509 out << "fWarningFraction" << "::" << fWarningFraction << endl;
510 out << "fHotSigma" << "::" << fHotSigma << endl;
b07ee441 511
512 out.close();
513 return;
514}
515
356c3e0c 516//_____________________________________________________________________
23ca956b 517Bool_t AliCaloCalibPedestal::AddInfo(AliCaloCalibPedestal *ped)
419341ea 518{
23ca956b 519 ValidateProfiles(); // make sure histos/profiles exist
419341ea 520 // just do this for the basic histograms/profiles that get filled in ProcessEvent
521 // may not have data for all modules, but let's just Add everything..
522 for (int i = 0; i < fModules; i++) {
2ce7ee7f 523 GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
524 GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
fcba0b23 525 GetPedLEDRefProfileLowGain(i)->Add( ped->GetPedLEDRefProfileLowGain(i) );
526 GetPedLEDRefProfileHighGain(i)->Add( ped->GetPedLEDRefProfileHighGain(i) );
2ce7ee7f 527 GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
528 GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
59ed788d 529 GetPeakHighGainHisto(i)->Add( ped->GetPeakHighGainHisto(i) );
530
419341ea 531 }//end for nModules
532
d50632ab 533 // We should also copy other pieces of info: counters and parameters
534 // (not number of columns and rows etc which should be the same)
535 // note that I just assign them here rather than Add them, but we
536 // normally just Add (e.g. in Preprocessor) one object so this should be fine.
537 fNEvents = ped->GetNEvents();
538 fNChanFills = ped->GetNChanFills();
539 fDeadTowers = ped->GetDeadTowerCount();
540 fNewDeadTowers = ped->GetDeadTowerNew();
541 fResurrectedTowers = ped->GetDeadTowerResurrected();
542 fRunNumber = ped->GetRunNumber();
543 fSelectPedestalSamples = ped->GetSelectPedestalSamples();
544 fFirstPedestalSample = ped->GetFirstPedestalSample();
545 fLastPedestalSample = ped->GetLastPedestalSample();
546 fDeadThreshold = ped->GetDeadThreshold();
547 fWarningThreshold = ped->GetWarningThreshold();
548 fWarningFraction = ped->GetWarningFraction();
549 fHotSigma = ped->GetHotSigma();
550
419341ea 551 // DeadMap; Diff profiles etc would need to be redone after this operation
552
553 return kTRUE;//We succesfully added info from the supplied object
554}
555
556//_____________________________________________________________________
f4fc542c 557Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
558{
559 // if fMapping is NULL the rawstream will crate its own mapping
32cd4c24 560 AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
30aa89b0 561 if (fDetType == kEmCal) {
562 rawReader->Select("EMCAL", 0, AliEMCALGeoParams::fgkLastAltroDDL) ; //select EMCAL DDL range
563 }
f4fc542c 564 return ProcessEvent(&rawStream);
565}
566
567//_____________________________________________________________________
32cd4c24 568Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStreamV3 *in)
356c3e0c 569{
a235e2bc 570 // Method to process=analyze one event in the data stream
356c3e0c 571 if (!in) return kFALSE; //Return right away if there's a null pointer
65bec413 572 in->Reset(); // just in case the next customer forgets to check if the stream was reset..
573
419341ea 574 fNEvents++; // one more event
23ca956b 575
576 if (fNEvents==1) ValidateProfiles(); // 1st event, make sure histos/profiles exist
3dba9483 577
578 // indices for the reading
579 int sample = 0;
3dba9483 580 int time = 0;
32cd4c24 581 int i = 0; // sample counter
582 int startBin = 0;
3dba9483 583
32cd4c24 584 // start loop over input stream
585 while (in->NextDDL()) {
586 while (in->NextChannel()) {
3dba9483 587
32cd4c24 588 // counters
bd3cd9c2 589 int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
56613d93 590 int nsamples = 0;
591
29f94584 592 // pedestal samples
593 int nPed = 0;
594 vector<int> pedSamples;
595
32cd4c24 596 while (in->NextBunch()) {
597 const UShort_t *sig = in->GetSignals();
598 startBin = in->GetStartTimeBin();
56613d93 599 nsamples += in->GetBunchLength();
32cd4c24 600 for (i = 0; i < in->GetBunchLength(); i++) {
601 sample = sig[i];
602 time = startBin--;
603
604 // check if it's a min or max value
605 if (sample < min) min = sample;
606 if (sample > max) max = sample;
607
608 // should we add it for the pedestal calculation?
609 if ( (fFirstPedestalSample<=time && time<=fLastPedestalSample) || // sample time in range
610 !fSelectPedestalSamples ) { // or we don't restrict the sample range.. - then we'll take all
29f94584 611 pedSamples.push_back( sig[i] );
612 nPed++;
32cd4c24 613 }
614
615 } // loop over samples in bunch
616 } // loop over bunches
3dba9483 617
56613d93 618 if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
619
32cd4c24 620 // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
356c3e0c 621 int arrayPos = in->GetModule(); //The modules are numbered starting from 0
622 if (arrayPos >= fModules) {
623 //TODO: return an error message, if appopriate (perhaps if debug>0?)
624 return kFALSE;
32cd4c24 625 }
356c3e0c 626 //Debug
627 if (arrayPos < 0 || arrayPos >= fModules) {
628 printf("Oh no: arrayPos = %i.\n", arrayPos);
629 }
32cd4c24 630
419341ea 631 fNChanFills++; // one more channel found, and profile to be filled
356c3e0c 632 //NOTE: coordinates are (column, row) for the profiles
29f94584 633 if ( in->IsLowGain() ) {
356c3e0c 634 //fill the low gain histograms
18db89b7 635 ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
29f94584 636 if (nPed>0) { // only fill pedestal info in case it could be calculated
637 for ( i=0; i<nPed; i++) {
638 ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), pedSamples[i]);
639 }
3dba9483 640 }
356c3e0c 641 }
29f94584 642 else if ( in->IsHighGain() ) {
3dba9483 643 //fill the high gain ones
18db89b7 644 ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
29f94584 645 if (nPed>0) { // only fill pedestal info in case it could be calculated
646 for ( i=0; i<nPed; i++) {
647 ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), pedSamples[i]);
648 }
3dba9483 649 }
59ed788d 650 // for warning checks
651 int idx = in->GetRow() + fRows * in->GetColumn();
652 ((TH2F*)fPeakMinusPedHighGainHisto[arrayPos])->Fill(idx, max - min);
29f94584 653 }
654 else if ( in->IsLEDMonData() ) {
655 // for LED Mon data, the mapping class holds the gain info in the Row variable
656 // and the Strip number in the Column..
657 int gain = in->GetRow();
658 int stripId = in->GetColumn();
659 if (nPed>0 && stripId<fLEDRefs) {
660 if (gain == 0) {
661 for ( i=0; i<nPed; i++) {
662 ((TProfile*)fPedestalLEDRefLowGain[arrayPos])->Fill(stripId, pedSamples[i]);
663 }
664 }
665 else {
666 for ( i=0; i<nPed; i++) {
667 ((TProfile*)fPedestalLEDRefHighGain[arrayPos])->Fill(stripId, pedSamples[i]);
668 }
669 }
670 }
671 }
419341ea 672
56613d93 673 } // nsamples>0 check, some data found for this channel; not only trailer/header
32cd4c24 674 }// end while over channel
675 }//end while over DDL's, of input stream
676
32cd4c24 677
356c3e0c 678 return kTRUE;
679}
680
681//_____________________________________________________________________
682Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
683{
684 //Saves all the histograms (or profiles, to be accurate) to the designated file
23ca956b 685 ValidateProfiles(); // make sure histos/profiles exist
356c3e0c 686 TFile destFile(fileName, "recreate");
687
688 if (destFile.IsZombie()) {
689 return kFALSE;
690 }
691
692 destFile.cd();
693
694 for (int i = 0; i < fModules; i++) {
695 if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
696 fPeakMinusPedLowGain[i]->Write();
697 }
698 if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
699 fPeakMinusPedHighGain[i]->Write();
700 }
701 if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
702 fPedestalLowGain[i]->Write();
703 }
704 if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
705 fPedestalHighGain[i]->Write();
706 }
29f94584 707 if( ((TProfile *)fPedestalLEDRefLowGain[i])->GetEntries() || saveEmptyHistos) {
708 fPedestalLEDRefLowGain[i]->Write();
18db89b7 709 }
29f94584 710 if( ((TProfile *)fPedestalLEDRefHighGain[i])->GetEntries() || saveEmptyHistos) {
711 fPedestalLEDRefHighGain[i]->Write();
18db89b7 712 }
59ed788d 713 if( ((TH2F *)fPeakMinusPedHighGainHisto[i])->GetEntries() || saveEmptyHistos) {
714 fPeakMinusPedHighGainHisto[i]->Write();
715 }
716
356c3e0c 717 }
718
719 destFile.Close();
720
721 return kTRUE;
722}
723
724//_____________________________________________________________________
725Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
726{
727
728 //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
729 TH1::AddDirectory(kFALSE);
730
731 TFile *sourceFile = new TFile(fileName);
732 if (sourceFile->IsZombie()) {
733 return kFALSE;//We couldn't load the reference
734 }
735
736 if (fReference) delete fReference;//Delete the reference object, if it already exists
737 fReference = 0;
738
739 fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
740
741 if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
742 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
743 fReference = 0;
744 return kFALSE;
745 }
746
747 delete sourceFile;
29f94584 748
356c3e0c 749 //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
750 //if we are called by someone who has set it to false...
751 TH1::AddDirectory(kTRUE);
752
753 return kTRUE;//We succesfully loaded the object
754}
755
35eabcf4 756
757//_____________________________________________________________________
758Bool_t AliCaloCalibPedestal::SetReference(AliCaloCalibPedestal *ref)
759{
760 if (fReference) delete fReference;//Delete the reference object, if it already exists
761 fReference = 0;
762
763 fReference = ref;
764
765 if (!fReference || (fReference->GetDetectorType() != fDetType)) {
766 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
767 fReference = 0;
768 return kFALSE;
769 }
770
771 return kTRUE;//We succesfully loaded the object
772}
773
356c3e0c 774//_____________________________________________________________________
775void AliCaloCalibPedestal::ValidateComparisonProfiles()
776{
a235e2bc 777 //Make sure the comparison histos exist
356c3e0c 778 if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
779 //the same time
780
781
782 //Then, loop for the requested number of modules
783 TString title, name;
784 for (int i = 0; i < fModules; i++) {
785 //Pedestals, low gain
786 name = "hPedlowgainDiff";
787 name += i;
788 title = "Pedestals difference, low gain, module ";
789 title += i;
790 fPedestalLowGainDiff.Add(new TProfile2D(name, title,
791 fColumns, 0.0, fColumns,
18db89b7 792 fRows, fRowMin, fRowMax,"s"));
356c3e0c 793
794 //Pedestals, high gain
795 name = "hPedhighgainDiff";
796 name += i;
797 title = "Pedestals difference, high gain, module ";
798 title += i;
799 fPedestalHighGainDiff.Add(new TProfile2D(name, title,
800 fColumns, 0.0, fColumns,
18db89b7 801 fRows, fRowMin, fRowMax,"s"));
802
29f94584 803 //LED Ref/Mon pedestals, low gain
804 name = "hPedestalLEDReflowgainDiff";
805 name += i;
806 title = "Pedestal difference LEDRef, low gain, module ";
807 title += i;
808 fPedestalLEDRefLowGainDiff.Add(new TProfile(name, title,
809 fLEDRefs, 0.0, fLEDRefs, "s"));
810
811 //LED Ref/Mon pedestals, high gain
812 name = "hPedestalLEDRefhighgainDiff";
813 name += i;
814 title = "Pedestal difference LEDRef, high gain, module ";
815 title += i;
816 fPedestalLEDRefHighGainDiff.Add(new TProfile(name, title,
817 fLEDRefs, 0.0, fLEDRefs, "s"));
818
356c3e0c 819 //Peak-Pedestals, high gain
820 name = "hPeakMinusPedhighgainDiff";
821 name += i;
822 title = "Peak-Pedestal difference, high gain, module ";
823 title += i;
824 fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
825 fColumns, 0.0, fColumns,
18db89b7 826 fRows, fRowMin, fRowMax,"s"));
29f94584 827
828 //Peak-Pedestals, low gain
829 name = "hPeakMinusPedlowgainDiff";
830 name += i;
831 title = "Peak-Pedestal difference, low gain, module ";
832 title += i;
833 fPeakMinusPedLowGainDiff.Add(new TProfile2D(name, title,
834 fColumns, 0.0, fColumns,
835 fRows, fRowMin, fRowMax,"s"));
356c3e0c 836
837 //Pedestals, low gain
838 name = "hPedlowgainRatio";
839 name += i;
840 title = "Pedestals ratio, low gain, module ";
841 title += i;
842 fPedestalLowGainRatio.Add(new TProfile2D(name, title,
843 fColumns, 0.0, fColumns,
18db89b7 844 fRows, fRowMin, fRowMax,"s"));
356c3e0c 845
846 //Pedestals, high gain
847 name = "hPedhighgainRatio";
848 name += i;
849 title = "Pedestals ratio, high gain, module ";
850 title += i;
851 fPedestalHighGainRatio.Add(new TProfile2D(name, title,
852 fColumns, 0.0, fColumns,
18db89b7 853 fRows, fRowMin, fRowMax,"s"));
29f94584 854
855 //LED Ref/Mon pedestals, low gain
35eabcf4 856 name = "hPedestalLEDReflowgainRatio";
29f94584 857 name += i;
858 title = "Pedestal ratio LEDRef, low gain, module ";
859 title += i;
860 fPedestalLEDRefLowGainRatio.Add(new TProfile(name, title,
861 fLEDRefs, 0.0, fLEDRefs, "s"));
862
863 //LED Ref/Mon pedestals, high gain
864 name = "hPedestalLEDRefhighgainRatio";
865 name += i;
866 title = "Pedestal ratio LEDRef, high gain, module ";
867 title += i;
868 fPedestalLEDRefHighGainRatio.Add(new TProfile(name, title,
869 fLEDRefs, 0.0, fLEDRefs, "s"));
356c3e0c 870
871 //Peak-Pedestals, low gain
872 name = "hPeakMinusPedlowgainRatio";
873 name += i;
874 title = "Peak-Pedestal ratio, low gain, module ";
875 title += i;
876 fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
877 fColumns, 0.0, fColumns,
18db89b7 878 fRows, fRowMin, fRowMax,"s"));
356c3e0c 879
880 //Peak-Pedestals, high gain
881 name = "hPeakMinusPedhighgainRatio";
882 name += i;
883 title = "Peak-Pedestal ratio, high gain, module ";
884 title += i;
885 fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
886 fColumns, 0.0, fColumns,
18db89b7 887 fRows, fRowMin, fRowMax,"s"));
356c3e0c 888
889 }//end for nModules create the histograms
890}
891
892//_____________________________________________________________________
893void AliCaloCalibPedestal::ComputeDiffAndRatio()
894{
a235e2bc 895 // calculate differences and ratios relative to a reference
23ca956b 896 ValidateProfiles(); // make sure histos/profiles exist
356c3e0c 897 ValidateComparisonProfiles();//Make sure the comparison histos exist
898
899 if (!fReference) {
900 return;//Return if the reference object isn't loaded
901 }
902
29f94584 903 int bin = 0;
904 double diff = 0;
905 double ratio = 1;
356c3e0c 906 for (int i = 0; i < fModules; i++) {
356c3e0c 907 //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
908 //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
29f94584 909 for (int j = 0; j < fColumns; j++) {
910 for (int k = 0; k < fRows; k++) {
911 bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
912
913 if (fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) > 0) {
914 diff = GetPeakProfileHighGain(i)->GetBinContent(bin) - fReference->GetPeakProfileHighGain(i)->GetBinContent(bin);
915 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(bin, diff);
916 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
917 ratio = GetPeakProfileHighGain(i)->GetBinContent(bin) / fReference->GetPeakProfileHighGain(i)->GetBinContent(bin);
918 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->SetBinContent(bin, ratio);
919 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->SetBinEntries(bin, 1);
920 }
921
922 if (fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) > 0) {
923 diff = GetPeakProfileLowGain(i)->GetBinContent(bin) - fReference->GetPeakProfileLowGain(i)->GetBinContent(bin);
924 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(bin, diff);
925 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
926 ratio = GetPeakProfileLowGain(i)->GetBinContent(bin) / fReference->GetPeakProfileLowGain(i)->GetBinContent(bin);
927 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->SetBinContent(bin, ratio);
928 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->SetBinEntries(bin, 1);
929 }
930
931 if (fReference->GetPedProfileHighGain(i)->GetBinContent(bin) > 0) {
932 diff = GetPedProfileHighGain(i)->GetBinContent(bin) - fReference->GetPedProfileHighGain(i)->GetBinContent(bin);
933 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(bin, diff);
934 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
935 ratio = GetPedProfileHighGain(i)->GetBinContent(bin) / fReference->GetPedProfileHighGain(i)->GetBinContent(bin);
936 ((TProfile2D*)fPedestalHighGainRatio[i])->SetBinContent(bin, ratio);
937 ((TProfile2D*)fPedestalHighGainRatio[i])->SetBinEntries(bin, 1);
938 }
939
940 if (fReference->GetPedProfileLowGain(i)->GetBinContent(bin) > 0) {
941 diff = GetPedProfileLowGain(i)->GetBinContent(bin) - fReference->GetPedProfileLowGain(i)->GetBinContent(bin);
942 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(bin, diff);
943 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
944 ratio = GetPedProfileLowGain(i)->GetBinContent(bin) / fReference->GetPedProfileLowGain(i)->GetBinContent(bin);
945 ((TProfile2D*)fPedestalLowGainRatio[i])->SetBinContent(bin, ratio);
946 ((TProfile2D*)fPedestalLowGainRatio[i])->SetBinEntries(bin, 1);
947 }
948
356c3e0c 949 } // rows
950 } // columns
29f94584 951
952 // same for LED Ref/Mon channels
953 for (int j = 0; j <= fLEDRefs; j++) {
954 bin = j+1;//Note that we assume here that all histos have the same structure...
955
956 if (fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) > 0) {
957 diff = GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) - fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin);
958 ((TProfile*)fPedestalLEDRefHighGainDiff[i])->SetBinContent(bin, diff);
959 ((TProfile*)fPedestalLEDRefHighGainDiff[i])->SetBinEntries(bin, 1);
960 ratio = GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) / fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin);
961 ((TProfile*)fPedestalLEDRefHighGainRatio[i])->SetBinContent(bin, ratio);
962 ((TProfile*)fPedestalLEDRefHighGainRatio[i])->SetBinEntries(bin, 1);
963 }
964
965 if (fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) > 0) {
966 diff = GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) - fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin);
967 ((TProfile*)fPedestalLEDRefLowGainDiff[i])->SetBinContent(bin, diff);
968 ((TProfile*)fPedestalLEDRefLowGainDiff[i])->SetBinEntries(bin, 1);
969 ratio = GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) / fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin);
970 ((TProfile*)fPedestalLEDRefLowGainRatio[i])->SetBinContent(bin, ratio);
971 ((TProfile*)fPedestalLEDRefLowGainRatio[i])->SetBinEntries(bin, 1);
972 }
973
974 }
975
356c3e0c 976 } // modules
977
978}
979
980//_____________________________________________________________________
59ed788d 981void AliCaloCalibPedestal::ComputeHotAndWarningTowers(const char * hotMapFile)
982{ // look for hot/noisy towers
23ca956b 983 ValidateProfiles(); // make sure histos/profiles exist
59ed788d 984 ofstream * fout = 0;
985 char name[512];//Quite a long temp buffer, just in case the filename includes a path
986
987 if (hotMapFile) {
988 snprintf(name, 512, "%s.txt", hotMapFile);
989 fout = new ofstream(name);
990 if (!fout->is_open()) {
991 delete fout;
992 fout = 0;//Set the pointer to empty if the file was not opened
993 }
994 }
995
996 for(int i = 0; i < fModules; i++){
997
998 //first we compute the peak-pedestal distribution for each supermodule...
999 if( GetPeakHighGainHisto(i)->GetEntries() > 0 ) {
1000 double min = GetPeakProfileHighGain(i)->GetBinContent(GetPeakProfileHighGain(i)->GetMinimumBin());
1001 double max = GetPeakProfileHighGain(i)->GetBinContent(GetPeakProfileHighGain(i)->GetMaximumBin());
1002 TH1D *hPeakFit = new TH1D(Form("hFit_%d", i), Form("hFit_%d", i), (int)((max-min)*10), min-1, max+1);
1003
1004 for (int j = 1; j <= fColumns; j++) {
1005 for (int k = 1; k <= fRows; k++) {
1006 hPeakFit->Fill(GetPeakProfileHighGain(i)->GetBinContent(j, k));
1007 }
1008 }
1009
c490d8e5 1010 //...and then we try to fit it
1011 double mean = hPeakFit->GetMean();
1012 double sigma = hPeakFit->GetRMS();
1013 try {
1014 hPeakFit->Fit("gaus", "OQ", "", mean - 3*sigma, mean + 3*sigma);
1015 mean = hPeakFit->GetFunction("gaus")->GetParameter(1);
1016 sigma = hPeakFit->GetFunction("gaus")->GetParameter(2);
1017 }
1018 catch (const std::exception & e) {
1019 printf("AliCaloCalibPedestal: TH1D PeakFit exception %s", e.what());
1020 }
59ed788d 1021 //hPeakFit->Draw();
59ed788d 1022
1023 delete hPeakFit;
1024
1025 //Then we look for warm/hot towers
1026 TH2F * hPeak2D = GetPeakHighGainHisto(i);
1027 hPeak2D->GetYaxis()->SetRangeUser( fWarningThreshold, hPeak2D->GetYaxis()->GetBinUpEdge(hPeak2D->GetNbinsY()) );
1028
1029 int idx = 0 ;
1030 int warnCounter = 0;
1031 for (int j = 1; j <= fColumns; j++) {
1032 for (int k = 1; k <= fRows; k++) {
1033 //we start looking for warm/warning towers...
1034 // histogram x-axis index
1035 idx = k-1 + fRows*(j-1); // this is what is used in the Fill call
1036 hPeak2D->GetXaxis()->SetRangeUser(idx, idx);
1037 warnCounter = (int) hPeak2D->Integral();
1038 if(warnCounter > fNEvents * fWarningFraction) {
1039 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kWarning);
1040 /* printf("mod %d col %d row %d warnCounter %d - status %d\n",
1041 i, j-1, k-1, warnCounter, (int) (kWarning)); */
1042 }
1043 //...then we look for hot ones (towers whose values are greater than mean + X*sigma)
1044 if(GetPeakProfileHighGain(i)->GetBinContent(j, k) > mean + fHotSigma*sigma ) {
1045 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kHot);
1046 /* printf("mod %d col %d row %d binc %d - status %d\n",
1047 i, j-1, k-1, (int)(GetPeakProfileHighGain(i)->GetBinContent(j, k)), (int) (kHot)); */
1048 }
1049
1050 //Write the status to the hot/warm map file, if the file is open.
1051 // module - column - row - status (1=dead, 2= warm/warning , 3 = hot, see .h file enum)
1052 if (fout && ((TH2D*)fDeadMap[i])->GetBinContent(j, k) > 1) {
1053
1054 (*fout) << i << " "
1055 << (j - 1) << " "
1056 << (k - 1) << " "
1057 << ((TH2D*)fDeadMap[i])->GetBinContent(j, k) << endl; }
1058
1059 }
1060 }
1061
1062 }
1063 }
1064 return;
1065}
1066
1067//_____________________________________________________________________
1068void AliCaloCalibPedestal::ComputeDeadTowers(const char * deadMapFile)
a235e2bc 1069{
23ca956b 1070 ValidateProfiles(); // make sure histos/profiles exist
a235e2bc 1071 //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
356c3e0c 1072 int countTot = 0;
1073 int countNew = 0;
1074 int countRes = 0;
1075 ofstream * fout = 0;
1076 ofstream * diff = 0;
1077 char name[512];//Quite a long temp buffer, just in case the filename includes a path
1078
1079 if (deadMapFile) {
1080 snprintf(name, 512, "%s.txt", deadMapFile);
1081 fout = new ofstream(name);
1082 snprintf(name, 512, "%sdiff.txt", deadMapFile);
1083 diff = new ofstream(name);
1084 if (!fout->is_open()) {
1085 delete fout;
1086 fout = 0;//Set the pointer to empty if the file was not opened
1087 }
1088 if (!diff->is_open()) {
1089 delete diff;
b48d1356 1090 diff = 0;//Set the pointer to empty if the file was not opened
356c3e0c 1091 }
1092 }
1093
1094 for (int i = 0; i < fModules; i++) {
1095 if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
1096 for (int j = 1; j <= fColumns; j++) {
1097 for (int k = 1; k <= fRows; k++) {
1098
59ed788d 1099 if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {//It's dead
356c3e0c 1100 countTot++;//One more dead total
1101 if (fout) {
1102 (*fout) << i << " "
59ed788d 1103 << (j - 1) << " "
1104 << (k - 1) << " "
356c3e0c 1105 << "1" << " "
1106 << "0" << endl;//Write the status to the deadmap file, if the file is open.
1107 }
1108
59ed788d 1109 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= fDeadThreshold) {
356c3e0c 1110 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
1111 countNew++;//This tower wasn't dead before!
1112 if (diff) {
1113 ( *diff) << i << " "
356c3e0c 1114 << (j - 1) << " "
59ed788d 1115 << (k - 1) << " "
356c3e0c 1116 << "1" << " "
1117 << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
1118 }
1119 }
1120 else {
1121 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
1122 }
1123 }
1124 else { //It's ALIVE!!
1125 //Don't bother with writing the live ones.
59ed788d 1126 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {
356c3e0c 1127 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
1128 countRes++; //This tower was dead before => it's a miracle! :P
1129 if (diff) {
1130 (*diff) << i << " "
356c3e0c 1131 << (j - 1) << " "
59ed788d 1132 << (k - 1) << " "
356c3e0c 1133 << "1" << " "
1134 << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
1135 }
1136 }
1137 else {
1138 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
1139 }
1140 }
1141
1142 }//end for k/rows
1143 }//end for j/columns
1144 }//end if GetEntries >= 0
1145
1146 }//end for modules
1147
1148 if (fout) {
1149 fout->close();
1150 delete fout;
1151 }
1152
1153 fDeadTowers = countTot;
1154 fNewDeadTowers = countNew;
1155 fResurrectedTowers = countRes;
1156}
1157
40164976 1158//_____________________________________________________________________
1159Bool_t AliCaloCalibPedestal::IsBadChannel(int imod, int icol, int irow) const
1160{
23ca956b 1161 // Check if status info histo exists
1162 if (!fDeadMap[imod]) {
1163 return kFALSE;
1164 }
1165
59ed788d 1166 //Check if channel is dead or hot.
1167 Int_t status = (Int_t) ( ((TH2D*)fDeadMap[imod])->GetBinContent(icol,irow) );
1168 if(status == kAlive)
1169 return kFALSE;
1170 else
1171 return kTRUE;
1172
40164976 1173}
1174
1175//_____________________________________________________________________
1176void AliCaloCalibPedestal::SetChannelStatus(int imod, int icol, int irow, int status)
1177{
23ca956b 1178 ValidateProfiles(); // make sure histos/profiles exist
59ed788d 1179 //Set status of channel dead, hot, alive ...
1180 ((TH2D*)fDeadMap[imod])->SetBinContent(icol, irow, status);
40164976 1181}