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