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 "TH1.h"
  38 #include "TFile.h"
  39
  40 #include <fstream>
  41 #include <iostream>
  42
  43 #include "AliCaloRawStream.h"
  44
  45 //The include file
  46 #include "AliCaloCalibPedestal.h"
  47
  48 ClassImp(AliCaloCalibPedestal)
  49
  50 using namespace std;
  51
  52 // ctor; initialize everything in order to avoid compiler warnings
  53 AliCaloCalibPedestal::AliCaloCalibPedestal(kDetType detectorType) :
  54   TObject(),
  55   fPedestalLowGain(),
  56   fPedestalHighGain(),
  57   fSampleLowGain(),
  58   fSampleHighGain(),
  59   fPeakMinusPedLowGain(),
  60   fPeakMinusPedHighGain(),
  61   fPedestalLowGainDiff(),
  62   fPedestalHighGainDiff(),
  63   fPeakMinusPedLowGainDiff(),
  64   fPeakMinusPedHighGainDiff(),
  65   fPedestalLowGainRatio(),
  66   fPedestalHighGainRatio(),
  67   fPeakMinusPedLowGainRatio(),
  68   fPeakMinusPedHighGainRatio(),
  69   fDeadMap(),
  70   fNEvents(0),
  71   fNChanFills(0),
  72   fDeadTowers(0),
  73   fNewDeadTowers(0),
  74   fResurrectedTowers(0),
  75   fReference(0),
  76   fDetType(kNone),
  77   fColumns(0),
  78   fRows(0),
  79   fModules(0),
  80   fRowMin(0),
  81   fRowMax(0),
  82   fRowMultiplier(0),
  83   fCaloString(),
  84   fMapping(NULL),
  85   fRunNumber(-1)
  86 {
  87   //Default constructor. First we set the detector-type related constants.
  88   if (detectorType == kPhos) {
  89     fColumns = fgkPhosCols;
  90     fRows = fgkPhosRows;
  91     fModules = fgkPhosModules;
  92     fCaloString = "PHOS";
  93     fRowMin = -1*fRows;
  94     fRowMax = 0;
  95     fRowMultiplier = -1;
  96   }
  97   else {
  98     //We'll just trust the enum to keep everything in line, so that if detectorType
  99     //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
 100     //case, if someone intentionally gives another number
 101     fColumns = fgkEmCalCols;
 102     fRows = fgkEmCalRows;
 103     fModules = fgkEmCalModules;
 104     fCaloString = "EMCAL";
 105     fRowMin = 0;
 106     fRowMax = fRows;
 107     fRowMultiplier = 1;
 108   }
 109   fDetType = detectorType;
 110
 111   //Then, loop for the requested number of modules
 112   TString title, name;
 113   for (int i = 0; i < fModules; i++) {
 114     //Pedestals, low gain
 115     name = "hPedlowgain";
 116     name += i;
 117     title = "Pedestals, low gain, module ";
 118     title += i;
 119     fPedestalLowGain.Add(new TProfile2D(name, title,
 120                                         fColumns, 0.0, fColumns,
 121                                         fRows, fRowMin, fRowMax,"s"));
 122
 123     //Pedestals, high gain
 124     name = "hPedhighgain";
 125     name += i;
 126     title = "Pedestals, high gain, module ";
 127     title += i;
 128     fPedestalHighGain.Add(new TProfile2D(name, title,
 129                                          fColumns, 0.0, fColumns,
 130                                          fRows, fRowMin, fRowMax,"s"));
 131     //All Samples, low gain
 132     name = "hSamplelowgain";
 133     name += i;
 134     title = "All Samples, low gain, module ";
 135     title += i;
 136     fSampleLowGain.Add(new TProfile2D(name, title,
 137                                         fColumns, 0.0, fColumns,
 138                                         fRows, fRowMin, fRowMax,"s"));
 139
 140     //All Samples, high gain
 141     name = "hSamplehighgain";
 142     name += i;
 143     title = "All Samples, high gain, module ";
 144     title += i;
 145     fSampleHighGain.Add(new TProfile2D(name, title,
 146                                          fColumns, 0.0, fColumns,
 147                                          fRows, fRowMin, fRowMax,"s"));
 148
 149     //Peak-Pedestals, low gain
 150     name = "hPeakMinusPedlowgain";
 151     name += i;
 152     title = "Peak-Pedestal, low gain, module ";
 153     title += i;
 154     fPeakMinusPedLowGain.Add(new TProfile2D(name, title,
 155                                             fColumns, 0.0, fColumns,
 156                                             fRows, fRowMin, fRowMax,"s"));
 157
 158     //Peak-Pedestals, high gain
 159     name = "hPeakMinusPedhighgain";
 160     name += i;
 161     title = "Peak-Pedestal, high gain, module ";
 162     title += i;
 163     fPeakMinusPedHighGain.Add(new TProfile2D(name, title,
 164                                              fColumns, 0.0, fColumns,
 165                                              fRows, fRowMin, fRowMax,"s"));
 166
 167     name = "hDeadMap";
 168     name += i;
 169     title = "Dead map, module ";
 170     title += i;
 171     fDeadMap.Add(new TH2D(name, title, fColumns, 0.0, fColumns,
 172                           fRows, fRowMin, fRowMax));
 173
 174   }//end for nModules create the histograms
 175
 176   //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
 177   fPedestalLowGain.Compress();
 178   fPedestalHighGain.Compress();
 179   fSampleLowGain.Compress();
 180   fSampleHighGain.Compress();
 181   fPeakMinusPedLowGain.Compress();
 182   fPeakMinusPedHighGain.Compress();
 183   fDeadMap.Compress();
 184   //Make them the owners of the profiles, so we don't need to care about deleting them
 185   //fPedestalLowGain.SetOwner();
 186   //fPedestalHighGain.SetOwner();
 187   //fPeakMinusPedLowGain.SetOwner();
 188   //fPeakMinusPedHighGain.SetOwner();
 189
 190 }
 191
 192 // dtor
 193 //_____________________________________________________________________
 194 AliCaloCalibPedestal::~AliCaloCalibPedestal()
 195 {
 196   if (fReference) delete fReference;//Delete the reference object, if it has been loaded
 197   //TObjArray will delete the histos/profiles when it is deleted.
 198 }
 199
 200 // copy ctor
 201 //_____________________________________________________________________
 202 AliCaloCalibPedestal::AliCaloCalibPedestal(const AliCaloCalibPedestal &ped) :
 203   TObject(ped),
 204   fPedestalLowGain(),
 205   fPedestalHighGain(),
 206   fSampleLowGain(),
 207   fSampleHighGain(),
 208   fPeakMinusPedLowGain(),
 209   fPeakMinusPedHighGain(),
 210   fPedestalLowGainDiff(),
 211   fPedestalHighGainDiff(),
 212   fPeakMinusPedLowGainDiff(),
 213   fPeakMinusPedHighGainDiff(),
 214   fPedestalLowGainRatio(),
 215   fPedestalHighGainRatio(),
 216   fPeakMinusPedLowGainRatio(),
 217   fPeakMinusPedHighGainRatio(),
 218   fDeadMap(),
 219   fNEvents(ped.GetNEvents()),
 220   fNChanFills(ped.GetNChanFills()),
 221   fDeadTowers(ped.GetDeadTowerCount()),
 222   fNewDeadTowers(ped.GetDeadTowerNew()),
 223   fResurrectedTowers(ped.GetDeadTowerResurrected()),
 224   fReference( 0 ), //! note that we do not try to copy the reference info here
 225   fDetType(ped.GetDetectorType()),
 226   fColumns(ped.GetColumns()),
 227   fRows(ped.GetRows()),
 228   fModules(ped.GetModules()),
 229   fRowMin(ped.GetRowMin()),
 230   fRowMax(ped.GetRowMax()),
 231   fRowMultiplier(ped.GetRowMultiplier()),
 232   fCaloString(ped.GetCaloString()),
 233   fMapping(NULL), //! note that we are not copying the map info
 234   fRunNumber(ped.GetRunNumber())
 235 {
 236   // Then the ObjArray ones; we add the histograms rather than trying TObjArray = assignment
 237   //DS: this has not really been tested yet..
 238   for (int i = 0; i < fModules; i++) {
 239     fPedestalLowGain.Add( ped.GetPedProfileLowGain(i) );
 240     fPedestalHighGain.Add( ped.GetPedProfileHighGain(i) );
 241     fSampleLowGain.Add( ped.GetSampleProfileLowGain(i) );
 242     fSampleHighGain.Add( ped.GetSampleProfileHighGain(i) );
 243     fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
 244     fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
 245
 246     fDeadMap.Add( ped.GetDeadMap(i) );
 247   }//end for nModules
 248
 249   //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
 250   fPedestalLowGain.Compress();
 251   fPedestalHighGain.Compress();
 252   fSampleLowGain.Compress();
 253   fSampleHighGain.Compress();
 254   fPeakMinusPedLowGain.Compress();
 255   fPeakMinusPedHighGain.Compress();
 256   fDeadMap.Compress();
 257 }
 258
 259 // assignment operator; use copy ctor to make life easy..
 260 //_____________________________________________________________________
 261 AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (const AliCaloCalibPedestal &source)
 262 {
 263   // assignment operator; use copy ctor
 264   if (&source == this) return *this;
 265
 266   new (this) AliCaloCalibPedestal(source);
 267   return *this;
 268 }
 269
 270 //_____________________________________________________________________
 271 void AliCaloCalibPedestal::Reset()
 272 {
 273   // Reset all arrays/histograms
 274   for (int i = 0; i < fModules; i++) {
 275     GetPedProfileLowGain(i)->Reset();
 276     GetPedProfileHighGain(i)->Reset();
 277     GetPeakProfileLowGain(i)->Reset();
 278     GetPeakProfileHighGain(i)->Reset();
 279     GetDeadMap(i)->Reset();
 280
 281     if (!fPedestalLowGainDiff.IsEmpty()) {
 282       //This means that the comparison profiles have been created.
 283
 284       GetPedProfileLowGainDiff(i)->Reset();
 285       GetPedProfileHighGainDiff(i)->Reset();
 286       GetPeakProfileLowGainDiff(i)->Reset();
 287       GetPeakProfileHighGainDiff(i)->Reset();
 288
 289       GetPedProfileLowGainRatio(i)->Reset();
 290       GetPedProfileHighGainRatio(i)->Reset();
 291       GetPeakProfileLowGainRatio(i)->Reset();
 292       GetPeakProfileHighGainRatio(i)->Reset();
 293     }
 294   }
 295   fNEvents = 0;
 296   fNChanFills = 0;
 297   fDeadTowers = 0;
 298   fNewDeadTowers = 0;
 299   fResurrectedTowers = 0;
 300
 301   //To think about: should fReference be deleted too?... let's not do it this time, at least...
 302 }
 303
 304 //_____________________________________________________________________
 305 Bool_t AliCaloCalibPedestal::AddInfo(const AliCaloCalibPedestal *ped)
 306 {
 307   // just do this for the basic histograms/profiles that get filled in ProcessEvent
 308   // may not have data for all modules, but let's just Add everything..
 309   for (int i = 0; i < fModules; i++) {
 310     GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
 311     GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
 312     GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
 313     GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
 314   }//end for nModules
 315
 316   // DeadMap; Diff profiles etc would need to be redone after this operation
 317
 318   return kTRUE;//We succesfully added info from the supplied object
 319 }
 320
 321 //_____________________________________________________________________
 322 Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
 323 {
 324   // if fMapping is NULL the rawstream will crate its own mapping
 325   AliCaloRawStream rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
 326   return ProcessEvent(&rawStream);
 327 }
 328
 329 //_____________________________________________________________________
 330 Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStream *in)
 331 {
 332   // Method to process=analyze one event in the data stream
 333   if (!in) return kFALSE; //Return right away if there's a null pointer
 334
 335   fNEvents++; // one more event
 336   int sample, i = 0; //The sample temp, and the sample number in current event.
 337   int max = fgkSampleMin, min = fgkSampleMax;//Use these for picking the pedestal
 338   int gain = 0;
 339
 340   while (in->Next()) {
 341     sample = in->GetSignal(); //Get the adc signal
 342     if (sample < min) min = sample;
 343     if (sample > max) max = sample;
 344     i++;
 345     if ( i >= in->GetTimeLength()) {
 346       //If we're here then we're done with this tower
 347       gain = 1 - in->IsLowGain();
 348
 349       int arrayPos = in->GetModule(); //The modules are numbered starting from 0
 350       if (arrayPos >= fModules) {
 351         //TODO: return an error message, if appopriate (perhaps if debug>0?)
 352         return kFALSE;
 353       }
 354
 355       //Debug
 356       if (arrayPos < 0 || arrayPos >= fModules) {
 357         printf("Oh no: arrayPos = %i.\n", arrayPos);
 358       }
 359
 360       fNChanFills++; // one more channel found, and profile to be filled
 361       //NOTE: coordinates are (column, row) for the profiles
 362       if (gain == 0) {
 363         //fill the low gain histograms
 364         ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), min);
 365         ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
 366         ((TProfile2D*)fSampleLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), sample);
 367       }
 368       else {//fill the high gain ones
 369         ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), min);
 370         ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
 371         ((TProfile2D*)fSampleHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), sample);
 372       }//end if gain
 373
 374       max = fgkSampleMin; min = fgkSampleMax;
 375       i = 0;
 376
 377     }//End if end of tower
 378
 379   }//end while, of stream
 380
 381   return kTRUE;
 382 }
 383
 384 //_____________________________________________________________________
 385 Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
 386 {
 387   //Saves all the histograms (or profiles, to be accurate) to the designated file
 388
 389   TFile destFile(fileName, "recreate");
 390
 391   if (destFile.IsZombie()) {
 392     return kFALSE;
 393   }
 394
 395   destFile.cd();
 396
 397   for (int i = 0; i < fModules; i++) {
 398     if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
 399       fPeakMinusPedLowGain[i]->Write();
 400     }
 401     if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
 402       fPeakMinusPedHighGain[i]->Write();
 403     }
 404     if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
 405       fPedestalLowGain[i]->Write();
 406     }
 407     if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
 408       fPedestalHighGain[i]->Write();
 409     }
 410     if( ((TProfile2D *)fSampleLowGain[i])->GetEntries() || saveEmptyHistos) {
 411       fSampleLowGain[i]->Write();
 412     }
 413     if( ((TProfile2D *)fSampleHighGain[i])->GetEntries() || saveEmptyHistos) {
 414       fSampleHighGain[i]->Write();
 415     }
 416   }
 417
 418   destFile.Close();
 419
 420   return kTRUE;
 421 }
 422
 423 //_____________________________________________________________________
 424 Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
 425 {
 426
 427   //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
 428   TH1::AddDirectory(kFALSE);
 429
 430   TFile *sourceFile = new TFile(fileName);
 431   if (sourceFile->IsZombie()) {
 432     return kFALSE;//We couldn't load the reference
 433   }
 434
 435   if (fReference) delete fReference;//Delete the reference object, if it already exists
 436   fReference = 0;
 437
 438   fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
 439
 440   if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
 441     if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
 442     fReference = 0;
 443     return kFALSE;
 444   }
 445
 446   delete sourceFile;
 447
 448   //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
 449   //if we are called by someone who has set it to false...
 450   TH1::AddDirectory(kTRUE);
 451
 452   return kTRUE;//We succesfully loaded the object
 453 }
 454
 455 //_____________________________________________________________________
 456 void AliCaloCalibPedestal::ValidateComparisonProfiles()
 457 {
 458   //Make sure the comparison histos exist
 459   if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
 460   //the same time
 461
 462
 463   //Then, loop for the requested number of modules
 464   TString title, name;
 465   for (int i = 0; i < fModules; i++) {
 466     //Pedestals, low gain
 467     name = "hPedlowgainDiff";
 468     name += i;
 469     title = "Pedestals difference, low gain, module ";
 470     title += i;
 471     fPedestalLowGainDiff.Add(new TProfile2D(name, title,
 472                                             fColumns, 0.0, fColumns,
 473                                             fRows, fRowMin, fRowMax,"s"));
 474
 475     //Pedestals, high gain
 476     name = "hPedhighgainDiff";
 477     name += i;
 478     title = "Pedestals difference, high gain, module ";
 479     title += i;
 480     fPedestalHighGainDiff.Add(new TProfile2D(name, title,
 481                                              fColumns, 0.0, fColumns,
 482                                              fRows, fRowMin, fRowMax,"s"));
 483
 484     //Peak-Pedestals, high gain
 485     name = "hPeakMinusPedhighgainDiff";
 486     name += i;
 487     title = "Peak-Pedestal difference, high gain, module ";
 488     title += i;
 489     fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
 490                                                  fColumns, 0.0, fColumns,
 491                                                  fRows, fRowMin, fRowMax,"s"));
 492
 493     //Pedestals, low gain
 494     name = "hPedlowgainRatio";
 495     name += i;
 496     title = "Pedestals ratio, low gain, module ";
 497     title += i;
 498     fPedestalLowGainRatio.Add(new TProfile2D(name, title,
 499                                              fColumns, 0.0, fColumns,
 500                                              fRows, fRowMin, fRowMax,"s"));
 501
 502     //Pedestals, high gain
 503     name = "hPedhighgainRatio";
 504     name += i;
 505     title = "Pedestals ratio, high gain, module ";
 506     title += i;
 507     fPedestalHighGainRatio.Add(new TProfile2D(name, title,
 508                                               fColumns, 0.0, fColumns,
 509                                               fRows, fRowMin, fRowMax,"s"));
 510
 511     //Peak-Pedestals, low gain
 512     name = "hPeakMinusPedlowgainRatio";
 513     name += i;
 514     title = "Peak-Pedestal ratio, low gain, module ";
 515     title += i;
 516     fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
 517                                                  fColumns, 0.0, fColumns,
 518                                                  fRows, fRowMin, fRowMax,"s"));
 519
 520     //Peak-Pedestals, high gain
 521     name = "hPeakMinusPedhighgainRatio";
 522     name += i;
 523     title = "Peak-Pedestal ratio, high gain, module ";
 524     title += i;
 525     fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
 526                                                   fColumns, 0.0, fColumns,
 527                                                   fRows, fRowMin, fRowMax,"s"));
 528
 529   }//end for nModules create the histograms
 530 }
 531
 532 //_____________________________________________________________________
 533 void AliCaloCalibPedestal::ComputeDiffAndRatio()
 534 {
 535   // calculate differences and ratios relative to a reference
 536   ValidateComparisonProfiles();//Make sure the comparison histos exist
 537
 538   if (!fReference) {
 539     return;//Return if the reference object isn't loaded
 540   }
 541
 542   for (int i = 0; i < fModules; i++) {
 543     //Compute the ratio of the histograms
 544
 545     ((TProfile2D*)fPedestalLowGainRatio[i])->Divide(GetPedProfileLowGain(i), fReference->GetPedProfileLowGain(i), 1.0, 1.0);
 546     ((TProfile2D*)fPedestalHighGainRatio[i])->Divide(GetPedProfileHighGain(i), fReference->GetPedProfileHighGain(i), 1.0, 1.0);
 547     ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->Divide(GetPeakProfileLowGain(i), fReference->GetPeakProfileLowGain(i), 1.0, 1.0);
 548     ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->Divide(GetPeakProfileHighGain(i), fReference->GetPeakProfileHighGain(i), 1.0, 1.0);
 549
 550     //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
 551     //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
 552     for (int j = 0; j <= fColumns; j++) {
 553       for (int k = 0; k <= fRows; k++) {
 554         int bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
 555         double diff = fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) - GetPeakProfileHighGain(i)->GetBinContent(bin);
 556         ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
 557         ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
 558
 559         diff = fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) - GetPeakProfileLowGain(i)->GetBinContent(bin);
 560         ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
 561         ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
 562
 563         diff = fReference->GetPedProfileHighGain(i)->GetBinContent(bin) - GetPedProfileHighGain(i)->GetBinContent(bin);
 564         ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
 565         ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
 566
 567         diff = fReference->GetPedProfileLowGain(i)->GetBinContent(bin) - GetPedProfileLowGain(i)->GetBinContent(bin);
 568         ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
 569         ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
 570
 571       } // rows
 572     } // columns
 573
 574   } // modules
 575
 576 }
 577
 578 //_____________________________________________________________________
 579 void AliCaloCalibPedestal::ComputeDeadTowers(int threshold, const char * deadMapFile)
 580 {
 581   //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
 582   int countTot = 0;
 583   int countNew = 0;
 584   int countRes = 0;
 585   ofstream * fout = 0;
 586   ofstream * diff = 0;
 587   char name[512];//Quite a long temp buffer, just in case the filename includes a path
 588
 589   if (deadMapFile) {
 590     snprintf(name, 512, "%s.txt", deadMapFile);
 591     fout = new ofstream(name);
 592     snprintf(name, 512, "%sdiff.txt", deadMapFile);
 593     diff = new ofstream(name);
 594     if (!fout->is_open()) {
 595       delete fout;
 596       fout = 0;//Set the pointer to empty if the file was not opened
 597     }
 598     if (!diff->is_open()) {
 599       delete diff;
 600       fout = 0;//Set the pointer to empty if the file was not opened
 601     }
 602   }
 603
 604   for (int i = 0; i < fModules; i++) {
 605     if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
 606       for (int j = 1; j <= fColumns; j++) {
 607         for (int k = 1; k <= fRows; k++) {
 608
 609           if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {//It's dead
 610             countTot++;//One more dead total
 611             if (fout) {
 612               (*fout) << i << " "
 613                       << (fRows - k) << " "
 614                       << (j-1) << " "
 615                       << "1" << " "
 616                       << "0" << endl;//Write the status to the deadmap file, if the file is open.
 617             }
 618
 619             if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= threshold) {
 620               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
 621               countNew++;//This tower wasn't dead before!
 622               if (diff) {
 623                 ( *diff) << i << " "
 624                          << (fRows - k) << " "
 625                          << (j - 1) << " "
 626                          << "1" << " "
 627                          << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
 628               }
 629             }
 630             else {
 631               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
 632             }
 633           }
 634           else { //It's ALIVE!!
 635             //Don't bother with writing the live ones.
 636             //if (fout)
 637             //  (*fout) << i << " "
 638             //     << (fRows - k) << " "
 639             //     << (j - 1) << " "
 640             //     << "1" << " "
 641             //     << "1" << endl;//Write the status to the deadmap file, if the file is open.
 642             if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {
 643               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
 644               countRes++; //This tower was dead before => it's a miracle! :P
 645               if (diff) {
 646                 (*diff) << i << " "
 647                         << (fRows - k) << " "
 648                         << (j - 1) << " "
 649                         << "1" << " "
 650                         << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
 651               }
 652             }
 653             else {
 654               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
 655             }
 656           }
 657
 658         }//end for k/rows
 659       }//end for j/columns
 660     }//end if GetEntries >= 0
 661
 662   }//end for modules
 663
 664  if (fout) {
 665    fout->close();
 666    delete fout;
 667  }
 668
 669  fDeadTowers = countTot;
 670  fNewDeadTowers = countNew;
 671  fResurrectedTowers = countRes;
 672 }
 673