SensorThickness was defined twice. Set inner chip thickness to 250mum to bypass bug...
[u/mrichter/AliRoot.git] / EMCAL / AliCaloCalibPedestal.cxx
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  **************************************************************************/
15 //* $Id$ */
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
37 //#include "TCanvas.h"
38 #include "TH1.h"
39 #include "TF1.h"
40 #include "TFile.h"
41 #include <fstream>
42 #include <sstream>
43 #include <iostream>
44 #include <stdexcept>
45 #include <cmath>
46
47 #include "AliRawReader.h"
48 #include "AliCaloRawStreamV3.h"
49
50 //The include file
51 #include "AliCaloCalibPedestal.h"
52
53 ClassImp(AliCaloCalibPedestal)
54
55 using namespace std;
56
57 // ctor; initialize everything in order to avoid compiler warnings
58 AliCaloCalibPedestal::AliCaloCalibPedestal(kDetType detectorType) :  
59   TObject(),
60   fPedestalLowGain(),
61   fPedestalHighGain(),
62   fPedestalLEDRefLowGain(),
63   fPedestalLEDRefHighGain(),
64   fPeakMinusPedLowGain(),
65   fPeakMinusPedHighGain(),
66   fPeakMinusPedHighGainHisto(),
67   fPedestalLowGainDiff(),
68   fPedestalHighGainDiff(),
69   fPedestalLEDRefLowGainDiff(),
70   fPedestalLEDRefHighGainDiff(),
71   fPeakMinusPedLowGainDiff(),
72   fPeakMinusPedHighGainDiff(),
73   fPedestalLowGainRatio(),
74   fPedestalHighGainRatio(),
75   fPedestalLEDRefLowGainRatio(),
76   fPedestalLEDRefHighGainRatio(),
77   fPeakMinusPedLowGainRatio(),
78   fPeakMinusPedHighGainRatio(),
79   fDeadMap(),
80   fNEvents(0),
81   fNChanFills(0),
82   fDeadTowers(0),
83   fNewDeadTowers(0),
84   fResurrectedTowers(0),
85   fReference(0),
86   fDetType(kNone),
87   fColumns(0),
88   fRows(0),
89   fLEDRefs(0),
90   fModules(0),
91   fRowMin(0),
92   fRowMax(0),
93   fRowMultiplier(0),
94   fCaloString(),
95   fMapping(NULL),
96   fRunNumber(-1),
97   fSelectPedestalSamples(kTRUE), 
98   fFirstPedestalSample(0),
99   fLastPedestalSample(15),
100   fDeadThreshold(5),
101   fWarningThreshold(50),
102   fWarningFraction(0.002),
103   fHotSigma(5)
104 {
105   //Default constructor. First we set the detector-type related constants.
106   if (detectorType == kPhos) {
107     fColumns = fgkPhosCols;
108     fRows = fgkPhosRows;
109     fLEDRefs = fgkPhosLEDRefs;
110     fModules = fgkPhosModules;
111     fCaloString = "PHOS";
112     fRowMin = -1*fRows;
113     fRowMax = 0;
114     fRowMultiplier = -1;
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
120     fColumns = AliEMCALGeoParams::fgkEMCALCols;
121     fRows = AliEMCALGeoParams::fgkEMCALRows;
122     fLEDRefs = AliEMCALGeoParams::fgkEMCALLEDRefs;
123     fModules = AliEMCALGeoParams::fgkEMCALModules;
124     fCaloString = "EMCAL";
125     fRowMin = 0;
126     fRowMax = fRows;
127     fRowMultiplier = 1;
128   } 
129   fDetType = detectorType;
130
131   // ValidateProfiles(); // not to be done in ctor; info from Axel N. 
132 }
133
134 //_____________________________________________________________________
135 void 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
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, 
151                                         fRows, fRowMin, fRowMax,"s"));
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, 
160                                          fRows, fRowMin, fRowMax,"s"));
161
162     //LED Ref/Mon pedestals, low gain
163     name = "hPedestalLEDReflowgain";
164     name += i;
165     title = "Pedestal LEDRef, low gain, module ";
166     title += i; 
167     fPedestalLEDRefLowGain.Add(new TProfile(name, title,
168                                             fLEDRefs, 0.0, fLEDRefs, "s"));
169     
170     //LED Ref/Mon pedestals, high gain
171     name = "hPedestalLEDRefhighgain";
172     name += i;
173     title = "Pedestal LEDRef, high gain, module ";
174     title += i; 
175     fPedestalLEDRefHighGain.Add(new TProfile(name, title,
176                                              fLEDRefs, 0.0, fLEDRefs, "s"));
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, 
185                                             fRows, fRowMin, fRowMax,"s"));
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, 
194                                              fRows, fRowMin, fRowMax,"s"));
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  
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, 
210                           fRows, fRowMin, fRowMax));
211   
212   }//end for nModules create the histograms
213
214   CompressAndSetOwner();
215 }
216
217 //_____________________________________________________________________
218 void AliCaloCalibPedestal::CompressAndSetOwner()
219
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();
223   fPedestalLEDRefLowGain.Compress();
224   fPedestalLEDRefHighGain.Compress();
225   fPeakMinusPedLowGain.Compress();
226   fPeakMinusPedHighGain.Compress();
227   fPeakMinusPedHighGainHisto.Compress();
228   fDeadMap.Compress();
229
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);
251 }
252
253 // dtor
254 //_____________________________________________________________________
255 AliCaloCalibPedestal::~AliCaloCalibPedestal()
256 {
257   //dtor
258   
259   if (fReference) delete fReference;//Delete the reference object, if it has been loaded
260   
261   // delete also TObjArray's 
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();
282   
283 }
284
285 // copy ctor
286 //_____________________________________________________________________
287 AliCaloCalibPedestal::AliCaloCalibPedestal(AliCaloCalibPedestal &ped) :
288   TObject(ped),
289   fPedestalLowGain(),
290   fPedestalHighGain(),
291   fPedestalLEDRefLowGain(),
292   fPedestalLEDRefHighGain(),
293   fPeakMinusPedLowGain(),
294   fPeakMinusPedHighGain(),
295   fPeakMinusPedHighGainHisto(),
296   fPedestalLowGainDiff(),
297   fPedestalHighGainDiff(),
298   fPedestalLEDRefLowGainDiff(),
299   fPedestalLEDRefHighGainDiff(),
300   fPeakMinusPedLowGainDiff(),
301   fPeakMinusPedHighGainDiff(),
302   fPedestalLowGainRatio(),
303   fPedestalHighGainRatio(),
304   fPedestalLEDRefLowGainRatio(),
305   fPedestalLEDRefHighGainRatio(),
306   fPeakMinusPedLowGainRatio(),
307   fPeakMinusPedHighGainRatio(),
308   fDeadMap(),
309   fNEvents(ped.GetNEvents()),
310   fNChanFills(ped.GetNChanFills()),
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()),
318   fLEDRefs(ped.GetLEDRefs()),
319   fModules(ped.GetModules()),
320   fRowMin(ped.GetRowMin()),
321   fRowMax(ped.GetRowMax()),
322   fRowMultiplier(ped.GetRowMultiplier()),
323   fCaloString(ped.GetCaloString()),
324   fMapping(NULL), //! note that we are not copying the map info
325   fRunNumber(ped.GetRunNumber()),
326   fSelectPedestalSamples(ped.GetSelectPedestalSamples()),
327   fFirstPedestalSample(ped.GetFirstPedestalSample()),
328   fLastPedestalSample(ped.GetLastPedestalSample()),
329   fDeadThreshold(ped.GetDeadThreshold()),
330   fWarningThreshold(ped.GetWarningThreshold()),
331   fWarningFraction(ped.GetWarningFraction()),
332   fHotSigma(ped.GetHotSigma())
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) );
339     fPedestalLEDRefLowGain.Add( ped.GetPedLEDRefProfileLowGain(i) );
340     fPedestalLEDRefHighGain.Add( ped.GetPedLEDRefProfileHighGain(i) );
341     fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
342     fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
343     fPeakMinusPedHighGainHisto.Add( ped.GetPeakHighGainHisto(i) );
344
345     fDeadMap.Add( ped.GetDeadMap(i) );  
346   }//end for nModules 
347  
348   CompressAndSetOwner();
349 }
350
351 // assignment operator; use copy ctor to make life easy..
352 //_____________________________________________________________________
353 AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (AliCaloCalibPedestal &source)
354 {
355   // assignment operator; use copy ctor
356   if (&source == this) return *this;
357
358   new (this) AliCaloCalibPedestal(source);
359   return *this;
360 }
361
362 //_____________________________________________________________________
363 void AliCaloCalibPedestal::Reset()
364 {   // Reset all arrays/histograms
365   ValidateProfiles(); // make sure histos/profiles exist
366   for (int i = 0; i < fModules; i++) {
367     GetPedProfileLowGain(i)->Reset();
368     GetPedProfileHighGain(i)->Reset();
369     GetPedLEDRefProfileLowGain(i)->Reset();
370     GetPedLEDRefProfileHighGain(i)->Reset();
371     GetPeakProfileLowGain(i)->Reset();
372     GetPeakProfileHighGain(i)->Reset();
373     GetPeakHighGainHisto(i)->Reset();
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();
381       GetPedLEDRefProfileLowGainDiff(i)->Reset();
382       GetPedLEDRefProfileHighGainDiff(i)->Reset();
383       GetPeakProfileLowGainDiff(i)->Reset();
384       GetPeakProfileHighGainDiff(i)->Reset();
385       
386       GetPedProfileLowGainRatio(i)->Reset();
387       GetPedProfileHighGainRatio(i)->Reset();
388       GetPedLEDRefProfileLowGainRatio(i)->Reset();
389       GetPedLEDRefProfileHighGainRatio(i)->Reset();
390       GetPeakProfileLowGainRatio(i)->Reset();
391       GetPeakProfileHighGainRatio(i)->Reset();
392     }
393   }
394   fNEvents = 0;
395   fNChanFills = 0;
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
403 // Parameter/cut handling
404 //_____________________________________________________________________
405 void AliCaloCalibPedestal::SetParametersFromFile(const char *parameterFile)
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
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
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                         
433       string keyValue; 
434       s >> keyValue;
435       
436       // check stream status
437       if( ( s.rdstate() & ios::failbit ) == ios::failbit) break;
438                         
439       // skip rest of line if comments found
440       if( keyValue.substr( 0, 2 ) == "//" ) break;
441                         
442       // look for "::" in keyValue pair
443       size_t position = keyValue.find( delimitor );
444       if( position == string::npos ) {
445         printf("wrong format for key::value pair: %s\n", keyValue.c_str());
446       }
447                                 
448       // split keyValue pair
449       string key( keyValue.substr( 0, position ) );
450       string value( keyValue.substr( position+delimitor.size(), 
451                                       keyValue.size()-delimitor.size() ) );
452                         
453       // check value does not contain a new delimitor
454       if( value.find( delimitor ) != string::npos ) {
455         printf("wrong format for key::value pair: %s\n", keyValue.c_str());
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
464       if ( (key == "fFirstPedestalSample") || (key == "fLastPedestalSample") || (key == "fDeadThreshold") || (key == "fWarningThreshold") || (key == "fWarningFraction") || (key == "fHotSigma") ) {
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         }
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
486       } // some match
487
488     }           
489   }
490
491   in.close();
492   return;
493         
494 }
495
496 //_____________________________________________________________________
497 void AliCaloCalibPedestal::WriteParametersToFile(const char *parameterFile)
498 {
499   //Write parameters in file.
500         
501   static const string delimitor("::");
502   ofstream out( parameterFile );
503   out << "// " << parameterFile << endl;
504   out << "fFirstPedestalSample" << "::" << fFirstPedestalSample << endl;
505   out << "fLastPedestalSample" << "::" << fLastPedestalSample << endl;
506   out << "fDeadThreshold" << "::" << fDeadThreshold << endl;
507   out << "fWarningThreshold" << "::" << fWarningThreshold << endl;
508   out << "fWarningFraction" << "::" << fWarningFraction << endl;
509   out << "fHotSigma" << "::" << fHotSigma << endl;
510
511   out.close();
512   return;
513 }
514
515 //_____________________________________________________________________
516 Bool_t AliCaloCalibPedestal::AddInfo(AliCaloCalibPedestal *ped)
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..
520   ValidateProfiles(); // make sure histos/profiles exist
521
522   for (int i = 0; i < fModules; i++) {
523     GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
524     GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
525     GetPedLEDRefProfileLowGain(i)->Add( ped->GetPedLEDRefProfileLowGain(i) );
526     GetPedLEDRefProfileHighGain(i)->Add( ped->GetPedLEDRefProfileHighGain(i) );
527     GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
528     GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
529     GetPeakHighGainHisto(i)->Add( ped->GetPeakHighGainHisto(i) );
530   }//end for nModules 
531
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
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 //_____________________________________________________________________
556 Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
557
558   // if fMapping is NULL the rawstream will crate its own mapping
559   AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
560   if (fDetType == kEmCal) {
561     rawReader->Select("EMCAL", 0, AliEMCALGeoParams::fgkLastAltroDDL) ; //select EMCAL DDL range 
562   }
563   return ProcessEvent(&rawStream);
564 }
565
566 //_____________________________________________________________________
567 Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStreamV3 *in)
568
569   // Method to process=analyze one event in the data stream
570   if (!in) return kFALSE; //Return right away if there's a null pointer
571   in->Reset(); // just in case the next customer forgets to check if the stream was reset..
572
573   fNEvents++; // one more event
574
575   if (fNEvents==1) ValidateProfiles(); // 1st event, make sure histos/profiles exist
576   
577   // indices for the reading
578   int sample = 0;
579   int time = 0;
580   int i = 0; // sample counter
581   int startBin = 0;
582
583   // start loop over input stream 
584   while (in->NextDDL()) {
585     while (in->NextChannel()) {
586
587       // counters
588       int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
589       int nsamples = 0;
590
591       // pedestal samples
592       int nPed = 0;
593       vector<int> pedSamples; 
594
595       while (in->NextBunch()) {
596         const UShort_t *sig = in->GetSignals();
597         startBin = in->GetStartTimeBin();
598         nsamples += in->GetBunchLength();
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 
610             pedSamples.push_back( sig[i] );
611             nPed++;
612           }
613           
614         } // loop over samples in bunch
615       } // loop over bunches
616
617       if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
618
619       // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
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;
624       }     
625       //Debug
626       if (arrayPos < 0 || arrayPos >= fModules) {
627         printf("Oh no: arrayPos = %i.\n", arrayPos); 
628       }
629       
630       fNChanFills++; // one more channel found, and profile to be filled
631       //NOTE: coordinates are (column, row) for the profiles
632       if ( in->IsLowGain() ) {
633         //fill the low gain histograms
634         ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
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           }
639         }
640       } 
641       else if ( in->IsHighGain() ) {    
642         //fill the high gain ones
643         ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
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           }       
648         }
649         // for warning checks
650         int idx = in->GetRow() + fRows * in->GetColumn();
651         ((TH2F*)fPeakMinusPedHighGainHisto[arrayPos])->Fill(idx, max - min);
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       }
671
672       } // nsamples>0 check, some data found for this channel; not only trailer/header
673     }// end while over channel   
674   }//end while over DDL's, of input stream 
675
676  
677   return kTRUE;
678 }
679
680 //_____________________________________________________________________
681 Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
682 {
683   //Saves all the histograms (or profiles, to be accurate) to the designated file
684   ValidateProfiles(); // make sure histos/profiles exist
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     }
706     if( ((TProfile *)fPedestalLEDRefLowGain[i])->GetEntries() || saveEmptyHistos) {
707       fPedestalLEDRefLowGain[i]->Write();
708     }
709     if( ((TProfile *)fPedestalLEDRefHighGain[i])->GetEntries() || saveEmptyHistos) {
710       fPedestalLEDRefHighGain[i]->Write();
711     }
712     if( ((TH2F *)fPeakMinusPedHighGainHisto[i])->GetEntries() || saveEmptyHistos) { 
713       fPeakMinusPedHighGainHisto[i]->Write();
714     }
715
716   } 
717   
718   destFile.Close();
719   
720   return kTRUE;
721 }
722
723 //_____________________________________________________________________
724 Bool_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;
747
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
755
756 //_____________________________________________________________________
757 Bool_t AliCaloCalibPedestal::SetReference(AliCaloCalibPedestal *ref)
758 { // set reference object
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
773 //_____________________________________________________________________
774 void AliCaloCalibPedestal::ValidateComparisonProfiles()
775 {
776   //Make sure the comparison histos exist
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, 
791                                             fRows, fRowMin, fRowMax,"s"));
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, 
800                                              fRows, fRowMin, fRowMax,"s"));
801
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
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, 
825                                                  fRows, fRowMin, fRowMax,"s"));
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"));
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, 
843                                              fRows, fRowMin, fRowMax,"s"));
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, 
852                                               fRows, fRowMin, fRowMax,"s"));
853
854     //LED Ref/Mon pedestals, low gain
855     name = "hPedestalLEDReflowgainRatio";
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"));
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, 
877                                                  fRows, fRowMin, fRowMax,"s"));
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, 
886                                                   fRows, fRowMin, fRowMax,"s"));
887     
888   }//end for nModules create the histograms
889 }
890
891 //_____________________________________________________________________
892 void AliCaloCalibPedestal::ComputeDiffAndRatio()
893 { // calculate differences and ratios relative to a reference
894   ValidateProfiles(); // make sure histos/profiles exist
895   ValidateComparisonProfiles();//Make sure the comparison histos exist
896  
897   if (!fReference) {
898     return;//Return if the reference object isn't loaded
899   }
900
901   int bin = 0;
902   double diff = 0;
903   double ratio = 1;
904   for (int i = 0; i < fModules; i++) {
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:
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
947       } // rows
948     } // columns
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
974   } // modules
975  
976 }
977
978 //_____________________________________________________________________
979 void AliCaloCalibPedestal::ComputeHotAndWarningTowers(const char * hotMapFile)
980 { // look for hot/noisy towers
981   ValidateProfiles(); // make sure histos/profiles exist
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
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       }      
1019       //hPeakFit->Draw();
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 //_____________________________________________________________________
1066 void AliCaloCalibPedestal::ComputeDeadTowers(const char * deadMapFile)
1067 {
1068   ValidateProfiles(); // make sure histos/profiles exist
1069   //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
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;
1088       diff = 0;//Set the pointer to empty if the file was not opened
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
1097           if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {//It's dead
1098             countTot++;//One more dead total
1099             if (fout) {
1100               (*fout) << i << " " 
1101                       << (j - 1) << " " 
1102                       << (k - 1) << " " 
1103                       << "1" << " " 
1104                       << "0" << endl;//Write the status to the deadmap file, if the file is open.
1105             }
1106             
1107             if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= fDeadThreshold) {
1108               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased); 
1109               countNew++;//This tower wasn't dead before!
1110               if (diff) {
1111                 ( *diff) << i << " " 
1112                          << (j - 1) << " " 
1113                          << (k - 1) << " " 
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.
1124             if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {
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 << " " 
1129                         << (j - 1) << " " 
1130                         << (k - 1) << " " 
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
1156 //_____________________________________________________________________
1157 Bool_t AliCaloCalibPedestal::IsBadChannel(int imod, int icol, int irow) const
1158 {
1159   // Check if status info histo exists
1160   if (!fDeadMap[imod]) { 
1161     return kFALSE;
1162   }
1163
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   
1171 }
1172
1173 //_____________________________________________________________________
1174 void AliCaloCalibPedestal::SetChannelStatus(int imod, int icol, int irow, int status)
1175 {
1176   ValidateProfiles(); // make sure histos/profiles exist
1177   //Set status of channel dead, hot, alive ...  
1178   ((TH2D*)fDeadMap[imod])->SetBinContent(icol, irow, status);   
1179 }