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