new unfolding for clusterization
[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
258 printf("Dtor\n");
61e4e2e3 259
65bec413 260
261 if (fReference){ printf("Ref\n"); delete fReference;}//Delete the reference object, if it has been loaded
262
61e4e2e3 263
65bec413 264 printf("Delete\n");
265
266 // delete also TObjArray's
267 fPedestalLowGain.Delete(); printf("D 1\n");
268 fPedestalHighGain.Delete();printf("D 2\n");
269 fPedestalLEDRefLowGain.Delete();printf("D 3\n");
270 fPedestalLEDRefHighGain.Delete();printf("D 4\n");
271 fPeakMinusPedLowGain.Delete();printf("D 5\n");
272 fPeakMinusPedHighGain.Delete();printf("D 6\n");
273 fPeakMinusPedHighGainHisto.Delete();printf("D 7\n");
274 fPedestalLowGainDiff.Delete();printf("D 8\n");
275 fPedestalHighGainDiff.Delete();printf("D 9\n");
276 fPedestalLEDRefLowGainDiff.Delete();printf("D 10\n");
277 fPedestalLEDRefHighGainDiff.Delete();printf("D 11\n");
278 fPeakMinusPedLowGainDiff.Delete();printf("D 12\n");
279 fPeakMinusPedHighGainDiff.Delete();printf("D 13\n");
280 fPedestalLowGainRatio.Delete();printf("D 14\n");
281 fPedestalHighGainRatio.Delete();printf("D 15\n");
282 fPedestalLEDRefLowGainRatio.Delete();printf("D 16\n");
283 fPedestalLEDRefHighGainRatio.Delete();printf("D 17\n");
284 fPeakMinusPedLowGainRatio.Delete();printf("D 18\n");
285 fPeakMinusPedHighGainRatio.Delete();printf("D 19\n");
286 fDeadMap.Delete();printf("D 20\n");
287
356c3e0c 288}
289
290// copy ctor
291//_____________________________________________________________________
23ca956b 292AliCaloCalibPedestal::AliCaloCalibPedestal(AliCaloCalibPedestal &ped) :
356c3e0c 293 TObject(ped),
294 fPedestalLowGain(),
295 fPedestalHighGain(),
29f94584 296 fPedestalLEDRefLowGain(),
297 fPedestalLEDRefHighGain(),
356c3e0c 298 fPeakMinusPedLowGain(),
299 fPeakMinusPedHighGain(),
59ed788d 300 fPeakMinusPedHighGainHisto(),
356c3e0c 301 fPedestalLowGainDiff(),
302 fPedestalHighGainDiff(),
29f94584 303 fPedestalLEDRefLowGainDiff(),
304 fPedestalLEDRefHighGainDiff(),
356c3e0c 305 fPeakMinusPedLowGainDiff(),
306 fPeakMinusPedHighGainDiff(),
307 fPedestalLowGainRatio(),
308 fPedestalHighGainRatio(),
29f94584 309 fPedestalLEDRefLowGainRatio(),
310 fPedestalLEDRefHighGainRatio(),
356c3e0c 311 fPeakMinusPedLowGainRatio(),
312 fPeakMinusPedHighGainRatio(),
313 fDeadMap(),
419341ea 314 fNEvents(ped.GetNEvents()),
315 fNChanFills(ped.GetNChanFills()),
356c3e0c 316 fDeadTowers(ped.GetDeadTowerCount()),
317 fNewDeadTowers(ped.GetDeadTowerNew()),
318 fResurrectedTowers(ped.GetDeadTowerResurrected()),
319 fReference( 0 ), //! note that we do not try to copy the reference info here
320 fDetType(ped.GetDetectorType()),
321 fColumns(ped.GetColumns()),
322 fRows(ped.GetRows()),
29f94584 323 fLEDRefs(ped.GetLEDRefs()),
356c3e0c 324 fModules(ped.GetModules()),
18db89b7 325 fRowMin(ped.GetRowMin()),
326 fRowMax(ped.GetRowMax()),
327 fRowMultiplier(ped.GetRowMultiplier()),
f4fc542c 328 fCaloString(ped.GetCaloString()),
329 fMapping(NULL), //! note that we are not copying the map info
3dba9483 330 fRunNumber(ped.GetRunNumber()),
331 fSelectPedestalSamples(ped.GetSelectPedestalSamples()),
332 fFirstPedestalSample(ped.GetFirstPedestalSample()),
59ed788d 333 fLastPedestalSample(ped.GetLastPedestalSample()),
334 fDeadThreshold(ped.GetDeadThreshold()),
335 fWarningThreshold(ped.GetWarningThreshold()),
336 fWarningFraction(ped.GetWarningFraction()),
337 fHotSigma(ped.GetHotSigma())
356c3e0c 338{
339 // Then the ObjArray ones; we add the histograms rather than trying TObjArray = assignment
340 //DS: this has not really been tested yet..
341 for (int i = 0; i < fModules; i++) {
342 fPedestalLowGain.Add( ped.GetPedProfileLowGain(i) );
343 fPedestalHighGain.Add( ped.GetPedProfileHighGain(i) );
29f94584 344 fPedestalLEDRefLowGain.Add( ped.GetPedLEDRefProfileLowGain(i) );
345 fPedestalLEDRefHighGain.Add( ped.GetPedLEDRefProfileHighGain(i) );
356c3e0c 346 fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
347 fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
59ed788d 348 fPeakMinusPedHighGainHisto.Add( ped.GetPeakHighGainHisto(i) );
356c3e0c 349
350 fDeadMap.Add( ped.GetDeadMap(i) );
351 }//end for nModules
352
23ca956b 353 CompressAndSetOwner();
356c3e0c 354}
355
356// assignment operator; use copy ctor to make life easy..
357//_____________________________________________________________________
23ca956b 358AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (AliCaloCalibPedestal &source)
356c3e0c 359{
a235e2bc 360 // assignment operator; use copy ctor
356c3e0c 361 if (&source == this) return *this;
362
363 new (this) AliCaloCalibPedestal(source);
364 return *this;
365}
366
367//_____________________________________________________________________
368void AliCaloCalibPedestal::Reset()
369{
23ca956b 370 ValidateProfiles(); // make sure histos/profiles exist
a235e2bc 371 // Reset all arrays/histograms
356c3e0c 372 for (int i = 0; i < fModules; i++) {
373 GetPedProfileLowGain(i)->Reset();
374 GetPedProfileHighGain(i)->Reset();
29f94584 375 GetPedLEDRefProfileLowGain(i)->Reset();
376 GetPedLEDRefProfileHighGain(i)->Reset();
356c3e0c 377 GetPeakProfileLowGain(i)->Reset();
378 GetPeakProfileHighGain(i)->Reset();
59ed788d 379 GetPeakHighGainHisto(i)->Reset();
356c3e0c 380 GetDeadMap(i)->Reset();
381
382 if (!fPedestalLowGainDiff.IsEmpty()) {
383 //This means that the comparison profiles have been created.
384
385 GetPedProfileLowGainDiff(i)->Reset();
386 GetPedProfileHighGainDiff(i)->Reset();
29f94584 387 GetPedLEDRefProfileLowGainDiff(i)->Reset();
388 GetPedLEDRefProfileHighGainDiff(i)->Reset();
356c3e0c 389 GetPeakProfileLowGainDiff(i)->Reset();
390 GetPeakProfileHighGainDiff(i)->Reset();
391
392 GetPedProfileLowGainRatio(i)->Reset();
393 GetPedProfileHighGainRatio(i)->Reset();
29f94584 394 GetPedLEDRefProfileLowGainRatio(i)->Reset();
395 GetPedLEDRefProfileHighGainRatio(i)->Reset();
356c3e0c 396 GetPeakProfileLowGainRatio(i)->Reset();
397 GetPeakProfileHighGainRatio(i)->Reset();
398 }
399 }
419341ea 400 fNEvents = 0;
401 fNChanFills = 0;
356c3e0c 402 fDeadTowers = 0;
403 fNewDeadTowers = 0;
404 fResurrectedTowers = 0;
405
406 //To think about: should fReference be deleted too?... let's not do it this time, at least...
407}
408
b07ee441 409// Parameter/cut handling
410//_____________________________________________________________________
411void AliCaloCalibPedestal::SetParametersFromFile(const char *parameterFile)
40164976 412{
413 // Note: this method is a bit more complicated than it really has to be
414 // - allowing for multiple entries per line, arbitrary order of the
415 // different variables etc. But I wanted to try and do this in as
416 // correct a C++ way as I could (as an exercise).
417
b07ee441 418 static const string delimitor("::");
419
420 // open, check input file
421 ifstream in( parameterFile );
422 if( !in ) {
423 printf("in AliCaloCalibPedestal::SetParametersFromFile - Using default/run_time parameters.\n");
424 return;
425 }
426
b07ee441 427
428 // read in
429 char readline[1024];
430 while ((in.rdstate() & ios::failbit) == 0 ) {
431
432 // Read into the raw char array and then construct a string
433 // to do the searching
434 in.getline(readline, 1024);
435 istringstream s(readline);
436
437 while ( ( s.rdstate() & ios::failbit ) == 0 ) {
438
ab962f7b 439 string keyValue;
440 s >> keyValue;
b07ee441 441
442 // check stream status
b63770f3 443 if( ( s.rdstate() & ios::failbit ) == ios::failbit) break;
b07ee441 444
445 // skip rest of line if comments found
ab962f7b 446 if( keyValue.substr( 0, 2 ) == "//" ) break;
b07ee441 447
ab962f7b 448 // look for "::" in keyValue pair
449 size_t position = keyValue.find( delimitor );
b07ee441 450 if( position == string::npos ) {
ab962f7b 451 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
b07ee441 452 }
453
ab962f7b 454 // split keyValue pair
455 string key( keyValue.substr( 0, position ) );
456 string value( keyValue.substr( position+delimitor.size(),
457 keyValue.size()-delimitor.size() ) );
b07ee441 458
459 // check value does not contain a new delimitor
460 if( value.find( delimitor ) != string::npos ) {
ab962f7b 461 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
b07ee441 462 }
463
464 // debug: check key value pair
465 // printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
466
467 // if the key matches with something we expect, we assign the new value
468 istringstream iss(value);
469 // the comparison strings defined at the beginning of this method
59ed788d 470 if ( (key == "fFirstPedestalSample") || (key == "fLastPedestalSample") || (key == "fDeadThreshold") || (key == "fWarningThreshold") || (key == "fWarningFraction") || (key == "fHotSigma") ) {
b07ee441 471 printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
472
473 if (key == "fFirstPedestalSample") {
474 iss >> fFirstPedestalSample;
475 }
476 else if (key == "fLastPedestalSample") {
477 iss >> fLastPedestalSample;
478 }
59ed788d 479 else if (key == "fDeadThreshold") {
480 iss >> fDeadThreshold;
481 }
482 else if (key == "fWarningThreshold") {
483 iss >> fWarningThreshold;
484 }
485 else if (key == "fWarningFraction") {
486 iss >> fWarningFraction;
487 }
488 else if (key == "fHotSigma") {
489 iss >> fHotSigma;
490 }
491
b07ee441 492 } // some match
493
494 }
495 }
496
497 in.close();
498 return;
499
500}
501
502//_____________________________________________________________________
503void AliCaloCalibPedestal::WriteParametersToFile(const char *parameterFile)
504{
40164976 505 //Write parameters in file.
506
b07ee441 507 static const string delimitor("::");
508 ofstream out( parameterFile );
509 out << "// " << parameterFile << endl;
510 out << "fFirstPedestalSample" << "::" << fFirstPedestalSample << endl;
511 out << "fLastPedestalSample" << "::" << fLastPedestalSample << endl;
59ed788d 512 out << "fDeadThreshold" << "::" << fDeadThreshold << endl;
513 out << "fWarningThreshold" << "::" << fWarningThreshold << endl;
514 out << "fWarningFraction" << "::" << fWarningFraction << endl;
515 out << "fHotSigma" << "::" << fHotSigma << endl;
b07ee441 516
517 out.close();
518 return;
519}
520
356c3e0c 521//_____________________________________________________________________
23ca956b 522Bool_t AliCaloCalibPedestal::AddInfo(AliCaloCalibPedestal *ped)
419341ea 523{
23ca956b 524 ValidateProfiles(); // make sure histos/profiles exist
419341ea 525 // just do this for the basic histograms/profiles that get filled in ProcessEvent
526 // may not have data for all modules, but let's just Add everything..
527 for (int i = 0; i < fModules; i++) {
2ce7ee7f 528 GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
529 GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
fcba0b23 530 GetPedLEDRefProfileLowGain(i)->Add( ped->GetPedLEDRefProfileLowGain(i) );
531 GetPedLEDRefProfileHighGain(i)->Add( ped->GetPedLEDRefProfileHighGain(i) );
2ce7ee7f 532 GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
533 GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
59ed788d 534 GetPeakHighGainHisto(i)->Add( ped->GetPeakHighGainHisto(i) );
535
419341ea 536 }//end for nModules
537
538 // DeadMap; Diff profiles etc would need to be redone after this operation
539
540 return kTRUE;//We succesfully added info from the supplied object
541}
542
543//_____________________________________________________________________
f4fc542c 544Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
545{
546 // if fMapping is NULL the rawstream will crate its own mapping
32cd4c24 547 AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
30aa89b0 548 if (fDetType == kEmCal) {
549 rawReader->Select("EMCAL", 0, AliEMCALGeoParams::fgkLastAltroDDL) ; //select EMCAL DDL range
550 }
f4fc542c 551 return ProcessEvent(&rawStream);
552}
553
554//_____________________________________________________________________
32cd4c24 555Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStreamV3 *in)
356c3e0c 556{
a235e2bc 557 // Method to process=analyze one event in the data stream
356c3e0c 558 if (!in) return kFALSE; //Return right away if there's a null pointer
65bec413 559 in->Reset(); // just in case the next customer forgets to check if the stream was reset..
560
419341ea 561 fNEvents++; // one more event
23ca956b 562
563 if (fNEvents==1) ValidateProfiles(); // 1st event, make sure histos/profiles exist
3dba9483 564
565 // indices for the reading
566 int sample = 0;
3dba9483 567 int time = 0;
32cd4c24 568 int i = 0; // sample counter
569 int startBin = 0;
3dba9483 570
32cd4c24 571 // start loop over input stream
572 while (in->NextDDL()) {
573 while (in->NextChannel()) {
3dba9483 574
32cd4c24 575 // counters
bd3cd9c2 576 int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
56613d93 577 int nsamples = 0;
578
29f94584 579 // pedestal samples
580 int nPed = 0;
581 vector<int> pedSamples;
582
32cd4c24 583 while (in->NextBunch()) {
584 const UShort_t *sig = in->GetSignals();
585 startBin = in->GetStartTimeBin();
56613d93 586 nsamples += in->GetBunchLength();
32cd4c24 587 for (i = 0; i < in->GetBunchLength(); i++) {
588 sample = sig[i];
589 time = startBin--;
590
591 // check if it's a min or max value
592 if (sample < min) min = sample;
593 if (sample > max) max = sample;
594
595 // should we add it for the pedestal calculation?
596 if ( (fFirstPedestalSample<=time && time<=fLastPedestalSample) || // sample time in range
597 !fSelectPedestalSamples ) { // or we don't restrict the sample range.. - then we'll take all
29f94584 598 pedSamples.push_back( sig[i] );
599 nPed++;
32cd4c24 600 }
601
602 } // loop over samples in bunch
603 } // loop over bunches
3dba9483 604
56613d93 605 if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
606
32cd4c24 607 // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
356c3e0c 608 int arrayPos = in->GetModule(); //The modules are numbered starting from 0
609 if (arrayPos >= fModules) {
610 //TODO: return an error message, if appopriate (perhaps if debug>0?)
611 return kFALSE;
32cd4c24 612 }
356c3e0c 613 //Debug
614 if (arrayPos < 0 || arrayPos >= fModules) {
615 printf("Oh no: arrayPos = %i.\n", arrayPos);
616 }
32cd4c24 617
419341ea 618 fNChanFills++; // one more channel found, and profile to be filled
356c3e0c 619 //NOTE: coordinates are (column, row) for the profiles
29f94584 620 if ( in->IsLowGain() ) {
356c3e0c 621 //fill the low gain histograms
18db89b7 622 ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
29f94584 623 if (nPed>0) { // only fill pedestal info in case it could be calculated
624 for ( i=0; i<nPed; i++) {
625 ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), pedSamples[i]);
626 }
3dba9483 627 }
356c3e0c 628 }
29f94584 629 else if ( in->IsHighGain() ) {
3dba9483 630 //fill the high gain ones
18db89b7 631 ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
29f94584 632 if (nPed>0) { // only fill pedestal info in case it could be calculated
633 for ( i=0; i<nPed; i++) {
634 ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), pedSamples[i]);
635 }
3dba9483 636 }
59ed788d 637 // for warning checks
638 int idx = in->GetRow() + fRows * in->GetColumn();
639 ((TH2F*)fPeakMinusPedHighGainHisto[arrayPos])->Fill(idx, max - min);
29f94584 640 }
641 else if ( in->IsLEDMonData() ) {
642 // for LED Mon data, the mapping class holds the gain info in the Row variable
643 // and the Strip number in the Column..
644 int gain = in->GetRow();
645 int stripId = in->GetColumn();
646 if (nPed>0 && stripId<fLEDRefs) {
647 if (gain == 0) {
648 for ( i=0; i<nPed; i++) {
649 ((TProfile*)fPedestalLEDRefLowGain[arrayPos])->Fill(stripId, pedSamples[i]);
650 }
651 }
652 else {
653 for ( i=0; i<nPed; i++) {
654 ((TProfile*)fPedestalLEDRefHighGain[arrayPos])->Fill(stripId, pedSamples[i]);
655 }
656 }
657 }
658 }
419341ea 659
56613d93 660 } // nsamples>0 check, some data found for this channel; not only trailer/header
32cd4c24 661 }// end while over channel
662 }//end while over DDL's, of input stream
663
32cd4c24 664
356c3e0c 665 return kTRUE;
666}
667
668//_____________________________________________________________________
669Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
670{
671 //Saves all the histograms (or profiles, to be accurate) to the designated file
23ca956b 672 ValidateProfiles(); // make sure histos/profiles exist
356c3e0c 673 TFile destFile(fileName, "recreate");
674
675 if (destFile.IsZombie()) {
676 return kFALSE;
677 }
678
679 destFile.cd();
680
681 for (int i = 0; i < fModules; i++) {
682 if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
683 fPeakMinusPedLowGain[i]->Write();
684 }
685 if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
686 fPeakMinusPedHighGain[i]->Write();
687 }
688 if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
689 fPedestalLowGain[i]->Write();
690 }
691 if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
692 fPedestalHighGain[i]->Write();
693 }
29f94584 694 if( ((TProfile *)fPedestalLEDRefLowGain[i])->GetEntries() || saveEmptyHistos) {
695 fPedestalLEDRefLowGain[i]->Write();
18db89b7 696 }
29f94584 697 if( ((TProfile *)fPedestalLEDRefHighGain[i])->GetEntries() || saveEmptyHistos) {
698 fPedestalLEDRefHighGain[i]->Write();
18db89b7 699 }
59ed788d 700 if( ((TH2F *)fPeakMinusPedHighGainHisto[i])->GetEntries() || saveEmptyHistos) {
701 fPeakMinusPedHighGainHisto[i]->Write();
702 }
703
356c3e0c 704 }
705
706 destFile.Close();
707
708 return kTRUE;
709}
710
711//_____________________________________________________________________
712Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
713{
714
715 //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
716 TH1::AddDirectory(kFALSE);
717
718 TFile *sourceFile = new TFile(fileName);
719 if (sourceFile->IsZombie()) {
720 return kFALSE;//We couldn't load the reference
721 }
722
723 if (fReference) delete fReference;//Delete the reference object, if it already exists
724 fReference = 0;
725
726 fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
727
728 if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
729 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
730 fReference = 0;
731 return kFALSE;
732 }
733
734 delete sourceFile;
29f94584 735
356c3e0c 736 //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
737 //if we are called by someone who has set it to false...
738 TH1::AddDirectory(kTRUE);
739
740 return kTRUE;//We succesfully loaded the object
741}
742
35eabcf4 743
744//_____________________________________________________________________
745Bool_t AliCaloCalibPedestal::SetReference(AliCaloCalibPedestal *ref)
746{
747 if (fReference) delete fReference;//Delete the reference object, if it already exists
748 fReference = 0;
749
750 fReference = ref;
751
752 if (!fReference || (fReference->GetDetectorType() != fDetType)) {
753 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
754 fReference = 0;
755 return kFALSE;
756 }
757
758 return kTRUE;//We succesfully loaded the object
759}
760
356c3e0c 761//_____________________________________________________________________
762void AliCaloCalibPedestal::ValidateComparisonProfiles()
763{
a235e2bc 764 //Make sure the comparison histos exist
356c3e0c 765 if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
766 //the same time
767
768
769 //Then, loop for the requested number of modules
770 TString title, name;
771 for (int i = 0; i < fModules; i++) {
772 //Pedestals, low gain
773 name = "hPedlowgainDiff";
774 name += i;
775 title = "Pedestals difference, low gain, module ";
776 title += i;
777 fPedestalLowGainDiff.Add(new TProfile2D(name, title,
778 fColumns, 0.0, fColumns,
18db89b7 779 fRows, fRowMin, fRowMax,"s"));
356c3e0c 780
781 //Pedestals, high gain
782 name = "hPedhighgainDiff";
783 name += i;
784 title = "Pedestals difference, high gain, module ";
785 title += i;
786 fPedestalHighGainDiff.Add(new TProfile2D(name, title,
787 fColumns, 0.0, fColumns,
18db89b7 788 fRows, fRowMin, fRowMax,"s"));
789
29f94584 790 //LED Ref/Mon pedestals, low gain
791 name = "hPedestalLEDReflowgainDiff";
792 name += i;
793 title = "Pedestal difference LEDRef, low gain, module ";
794 title += i;
795 fPedestalLEDRefLowGainDiff.Add(new TProfile(name, title,
796 fLEDRefs, 0.0, fLEDRefs, "s"));
797
798 //LED Ref/Mon pedestals, high gain
799 name = "hPedestalLEDRefhighgainDiff";
800 name += i;
801 title = "Pedestal difference LEDRef, high gain, module ";
802 title += i;
803 fPedestalLEDRefHighGainDiff.Add(new TProfile(name, title,
804 fLEDRefs, 0.0, fLEDRefs, "s"));
805
356c3e0c 806 //Peak-Pedestals, high gain
807 name = "hPeakMinusPedhighgainDiff";
808 name += i;
809 title = "Peak-Pedestal difference, high gain, module ";
810 title += i;
811 fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
812 fColumns, 0.0, fColumns,
18db89b7 813 fRows, fRowMin, fRowMax,"s"));
29f94584 814
815 //Peak-Pedestals, low gain
816 name = "hPeakMinusPedlowgainDiff";
817 name += i;
818 title = "Peak-Pedestal difference, low gain, module ";
819 title += i;
820 fPeakMinusPedLowGainDiff.Add(new TProfile2D(name, title,
821 fColumns, 0.0, fColumns,
822 fRows, fRowMin, fRowMax,"s"));
356c3e0c 823
824 //Pedestals, low gain
825 name = "hPedlowgainRatio";
826 name += i;
827 title = "Pedestals ratio, low gain, module ";
828 title += i;
829 fPedestalLowGainRatio.Add(new TProfile2D(name, title,
830 fColumns, 0.0, fColumns,
18db89b7 831 fRows, fRowMin, fRowMax,"s"));
356c3e0c 832
833 //Pedestals, high gain
834 name = "hPedhighgainRatio";
835 name += i;
836 title = "Pedestals ratio, high gain, module ";
837 title += i;
838 fPedestalHighGainRatio.Add(new TProfile2D(name, title,
839 fColumns, 0.0, fColumns,
18db89b7 840 fRows, fRowMin, fRowMax,"s"));
29f94584 841
842 //LED Ref/Mon pedestals, low gain
35eabcf4 843 name = "hPedestalLEDReflowgainRatio";
29f94584 844 name += i;
845 title = "Pedestal ratio LEDRef, low gain, module ";
846 title += i;
847 fPedestalLEDRefLowGainRatio.Add(new TProfile(name, title,
848 fLEDRefs, 0.0, fLEDRefs, "s"));
849
850 //LED Ref/Mon pedestals, high gain
851 name = "hPedestalLEDRefhighgainRatio";
852 name += i;
853 title = "Pedestal ratio LEDRef, high gain, module ";
854 title += i;
855 fPedestalLEDRefHighGainRatio.Add(new TProfile(name, title,
856 fLEDRefs, 0.0, fLEDRefs, "s"));
356c3e0c 857
858 //Peak-Pedestals, low gain
859 name = "hPeakMinusPedlowgainRatio";
860 name += i;
861 title = "Peak-Pedestal ratio, low gain, module ";
862 title += i;
863 fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
864 fColumns, 0.0, fColumns,
18db89b7 865 fRows, fRowMin, fRowMax,"s"));
356c3e0c 866
867 //Peak-Pedestals, high gain
868 name = "hPeakMinusPedhighgainRatio";
869 name += i;
870 title = "Peak-Pedestal ratio, high gain, module ";
871 title += i;
872 fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
873 fColumns, 0.0, fColumns,
18db89b7 874 fRows, fRowMin, fRowMax,"s"));
356c3e0c 875
876 }//end for nModules create the histograms
877}
878
879//_____________________________________________________________________
880void AliCaloCalibPedestal::ComputeDiffAndRatio()
881{
a235e2bc 882 // calculate differences and ratios relative to a reference
23ca956b 883 ValidateProfiles(); // make sure histos/profiles exist
356c3e0c 884 ValidateComparisonProfiles();//Make sure the comparison histos exist
885
886 if (!fReference) {
887 return;//Return if the reference object isn't loaded
888 }
889
29f94584 890 int bin = 0;
891 double diff = 0;
892 double ratio = 1;
356c3e0c 893 for (int i = 0; i < fModules; i++) {
356c3e0c 894 //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
895 //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
29f94584 896 for (int j = 0; j < fColumns; j++) {
897 for (int k = 0; k < fRows; k++) {
898 bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
899
900 if (fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) > 0) {
901 diff = GetPeakProfileHighGain(i)->GetBinContent(bin) - fReference->GetPeakProfileHighGain(i)->GetBinContent(bin);
902 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(bin, diff);
903 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
904 ratio = GetPeakProfileHighGain(i)->GetBinContent(bin) / fReference->GetPeakProfileHighGain(i)->GetBinContent(bin);
905 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->SetBinContent(bin, ratio);
906 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->SetBinEntries(bin, 1);
907 }
908
909 if (fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) > 0) {
910 diff = GetPeakProfileLowGain(i)->GetBinContent(bin) - fReference->GetPeakProfileLowGain(i)->GetBinContent(bin);
911 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(bin, diff);
912 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
913 ratio = GetPeakProfileLowGain(i)->GetBinContent(bin) / fReference->GetPeakProfileLowGain(i)->GetBinContent(bin);
914 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->SetBinContent(bin, ratio);
915 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->SetBinEntries(bin, 1);
916 }
917
918 if (fReference->GetPedProfileHighGain(i)->GetBinContent(bin) > 0) {
919 diff = GetPedProfileHighGain(i)->GetBinContent(bin) - fReference->GetPedProfileHighGain(i)->GetBinContent(bin);
920 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(bin, diff);
921 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
922 ratio = GetPedProfileHighGain(i)->GetBinContent(bin) / fReference->GetPedProfileHighGain(i)->GetBinContent(bin);
923 ((TProfile2D*)fPedestalHighGainRatio[i])->SetBinContent(bin, ratio);
924 ((TProfile2D*)fPedestalHighGainRatio[i])->SetBinEntries(bin, 1);
925 }
926
927 if (fReference->GetPedProfileLowGain(i)->GetBinContent(bin) > 0) {
928 diff = GetPedProfileLowGain(i)->GetBinContent(bin) - fReference->GetPedProfileLowGain(i)->GetBinContent(bin);
929 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(bin, diff);
930 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
931 ratio = GetPedProfileLowGain(i)->GetBinContent(bin) / fReference->GetPedProfileLowGain(i)->GetBinContent(bin);
932 ((TProfile2D*)fPedestalLowGainRatio[i])->SetBinContent(bin, ratio);
933 ((TProfile2D*)fPedestalLowGainRatio[i])->SetBinEntries(bin, 1);
934 }
935
356c3e0c 936 } // rows
937 } // columns
29f94584 938
939 // same for LED Ref/Mon channels
940 for (int j = 0; j <= fLEDRefs; j++) {
941 bin = j+1;//Note that we assume here that all histos have the same structure...
942
943 if (fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) > 0) {
944 diff = GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) - fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin);
945 ((TProfile*)fPedestalLEDRefHighGainDiff[i])->SetBinContent(bin, diff);
946 ((TProfile*)fPedestalLEDRefHighGainDiff[i])->SetBinEntries(bin, 1);
947 ratio = GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) / fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin);
948 ((TProfile*)fPedestalLEDRefHighGainRatio[i])->SetBinContent(bin, ratio);
949 ((TProfile*)fPedestalLEDRefHighGainRatio[i])->SetBinEntries(bin, 1);
950 }
951
952 if (fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) > 0) {
953 diff = GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) - fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin);
954 ((TProfile*)fPedestalLEDRefLowGainDiff[i])->SetBinContent(bin, diff);
955 ((TProfile*)fPedestalLEDRefLowGainDiff[i])->SetBinEntries(bin, 1);
956 ratio = GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) / fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin);
957 ((TProfile*)fPedestalLEDRefLowGainRatio[i])->SetBinContent(bin, ratio);
958 ((TProfile*)fPedestalLEDRefLowGainRatio[i])->SetBinEntries(bin, 1);
959 }
960
961 }
962
356c3e0c 963 } // modules
964
965}
966
967//_____________________________________________________________________
59ed788d 968void AliCaloCalibPedestal::ComputeHotAndWarningTowers(const char * hotMapFile)
969{ // look for hot/noisy towers
23ca956b 970 ValidateProfiles(); // make sure histos/profiles exist
59ed788d 971 ofstream * fout = 0;
972 char name[512];//Quite a long temp buffer, just in case the filename includes a path
973
974 if (hotMapFile) {
975 snprintf(name, 512, "%s.txt", hotMapFile);
976 fout = new ofstream(name);
977 if (!fout->is_open()) {
978 delete fout;
979 fout = 0;//Set the pointer to empty if the file was not opened
980 }
981 }
982
983 for(int i = 0; i < fModules; i++){
984
985 //first we compute the peak-pedestal distribution for each supermodule...
986 if( GetPeakHighGainHisto(i)->GetEntries() > 0 ) {
987 double min = GetPeakProfileHighGain(i)->GetBinContent(GetPeakProfileHighGain(i)->GetMinimumBin());
988 double max = GetPeakProfileHighGain(i)->GetBinContent(GetPeakProfileHighGain(i)->GetMaximumBin());
989 TH1D *hPeakFit = new TH1D(Form("hFit_%d", i), Form("hFit_%d", i), (int)((max-min)*10), min-1, max+1);
990
991 for (int j = 1; j <= fColumns; j++) {
992 for (int k = 1; k <= fRows; k++) {
993 hPeakFit->Fill(GetPeakProfileHighGain(i)->GetBinContent(j, k));
994 }
995 }
996
c490d8e5 997 //...and then we try to fit it
998 double mean = hPeakFit->GetMean();
999 double sigma = hPeakFit->GetRMS();
1000 try {
1001 hPeakFit->Fit("gaus", "OQ", "", mean - 3*sigma, mean + 3*sigma);
1002 mean = hPeakFit->GetFunction("gaus")->GetParameter(1);
1003 sigma = hPeakFit->GetFunction("gaus")->GetParameter(2);
1004 }
1005 catch (const std::exception & e) {
1006 printf("AliCaloCalibPedestal: TH1D PeakFit exception %s", e.what());
1007 }
59ed788d 1008 //hPeakFit->Draw();
59ed788d 1009
1010 delete hPeakFit;
1011
1012 //Then we look for warm/hot towers
1013 TH2F * hPeak2D = GetPeakHighGainHisto(i);
1014 hPeak2D->GetYaxis()->SetRangeUser( fWarningThreshold, hPeak2D->GetYaxis()->GetBinUpEdge(hPeak2D->GetNbinsY()) );
1015
1016 int idx = 0 ;
1017 int warnCounter = 0;
1018 for (int j = 1; j <= fColumns; j++) {
1019 for (int k = 1; k <= fRows; k++) {
1020 //we start looking for warm/warning towers...
1021 // histogram x-axis index
1022 idx = k-1 + fRows*(j-1); // this is what is used in the Fill call
1023 hPeak2D->GetXaxis()->SetRangeUser(idx, idx);
1024 warnCounter = (int) hPeak2D->Integral();
1025 if(warnCounter > fNEvents * fWarningFraction) {
1026 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kWarning);
1027 /* printf("mod %d col %d row %d warnCounter %d - status %d\n",
1028 i, j-1, k-1, warnCounter, (int) (kWarning)); */
1029 }
1030 //...then we look for hot ones (towers whose values are greater than mean + X*sigma)
1031 if(GetPeakProfileHighGain(i)->GetBinContent(j, k) > mean + fHotSigma*sigma ) {
1032 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kHot);
1033 /* printf("mod %d col %d row %d binc %d - status %d\n",
1034 i, j-1, k-1, (int)(GetPeakProfileHighGain(i)->GetBinContent(j, k)), (int) (kHot)); */
1035 }
1036
1037 //Write the status to the hot/warm map file, if the file is open.
1038 // module - column - row - status (1=dead, 2= warm/warning , 3 = hot, see .h file enum)
1039 if (fout && ((TH2D*)fDeadMap[i])->GetBinContent(j, k) > 1) {
1040
1041 (*fout) << i << " "
1042 << (j - 1) << " "
1043 << (k - 1) << " "
1044 << ((TH2D*)fDeadMap[i])->GetBinContent(j, k) << endl; }
1045
1046 }
1047 }
1048
1049 }
1050 }
1051 return;
1052}
1053
1054//_____________________________________________________________________
1055void AliCaloCalibPedestal::ComputeDeadTowers(const char * deadMapFile)
a235e2bc 1056{
23ca956b 1057 ValidateProfiles(); // make sure histos/profiles exist
a235e2bc 1058 //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
356c3e0c 1059 int countTot = 0;
1060 int countNew = 0;
1061 int countRes = 0;
1062 ofstream * fout = 0;
1063 ofstream * diff = 0;
1064 char name[512];//Quite a long temp buffer, just in case the filename includes a path
1065
1066 if (deadMapFile) {
1067 snprintf(name, 512, "%s.txt", deadMapFile);
1068 fout = new ofstream(name);
1069 snprintf(name, 512, "%sdiff.txt", deadMapFile);
1070 diff = new ofstream(name);
1071 if (!fout->is_open()) {
1072 delete fout;
1073 fout = 0;//Set the pointer to empty if the file was not opened
1074 }
1075 if (!diff->is_open()) {
1076 delete diff;
1077 fout = 0;//Set the pointer to empty if the file was not opened
1078 }
1079 }
1080
1081 for (int i = 0; i < fModules; i++) {
1082 if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
1083 for (int j = 1; j <= fColumns; j++) {
1084 for (int k = 1; k <= fRows; k++) {
1085
59ed788d 1086 if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {//It's dead
356c3e0c 1087 countTot++;//One more dead total
1088 if (fout) {
1089 (*fout) << i << " "
59ed788d 1090 << (j - 1) << " "
1091 << (k - 1) << " "
356c3e0c 1092 << "1" << " "
1093 << "0" << endl;//Write the status to the deadmap file, if the file is open.
1094 }
1095
59ed788d 1096 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= fDeadThreshold) {
356c3e0c 1097 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
1098 countNew++;//This tower wasn't dead before!
1099 if (diff) {
1100 ( *diff) << i << " "
356c3e0c 1101 << (j - 1) << " "
59ed788d 1102 << (k - 1) << " "
356c3e0c 1103 << "1" << " "
1104 << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
1105 }
1106 }
1107 else {
1108 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
1109 }
1110 }
1111 else { //It's ALIVE!!
1112 //Don't bother with writing the live ones.
59ed788d 1113 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {
356c3e0c 1114 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
1115 countRes++; //This tower was dead before => it's a miracle! :P
1116 if (diff) {
1117 (*diff) << i << " "
356c3e0c 1118 << (j - 1) << " "
59ed788d 1119 << (k - 1) << " "
356c3e0c 1120 << "1" << " "
1121 << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
1122 }
1123 }
1124 else {
1125 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
1126 }
1127 }
1128
1129 }//end for k/rows
1130 }//end for j/columns
1131 }//end if GetEntries >= 0
1132
1133 }//end for modules
1134
1135 if (fout) {
1136 fout->close();
1137 delete fout;
1138 }
1139
1140 fDeadTowers = countTot;
1141 fNewDeadTowers = countNew;
1142 fResurrectedTowers = countRes;
1143}
1144
40164976 1145//_____________________________________________________________________
1146Bool_t AliCaloCalibPedestal::IsBadChannel(int imod, int icol, int irow) const
1147{
23ca956b 1148 // Check if status info histo exists
1149 if (!fDeadMap[imod]) {
1150 return kFALSE;
1151 }
1152
59ed788d 1153 //Check if channel is dead or hot.
1154 Int_t status = (Int_t) ( ((TH2D*)fDeadMap[imod])->GetBinContent(icol,irow) );
1155 if(status == kAlive)
1156 return kFALSE;
1157 else
1158 return kTRUE;
1159
40164976 1160}
1161
1162//_____________________________________________________________________
1163void AliCaloCalibPedestal::SetChannelStatus(int imod, int icol, int irow, int status)
1164{
23ca956b 1165 ValidateProfiles(); // make sure histos/profiles exist
59ed788d 1166 //Set status of channel dead, hot, alive ...
1167 ((TH2D*)fDeadMap[imod])->SetBinContent(icol, irow, status);
40164976 1168}