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; // init to not valid value
 348       if ( in->IsLowGain() ) {
 349         gain = 0;
 350       }
 351       else if ( in->IsHighGain() ) {
 352         gain = 1;
 353       }
 354
 355       int arrayPos = in->GetModule(); //The modules are numbered starting from 0
 356       if (arrayPos >= fModules) {
 357         //TODO: return an error message, if appopriate (perhaps if debug>0?)
 358         return kFALSE;
 359       }
 360
 361       //Debug
 362       if (arrayPos < 0 || arrayPos >= fModules) {
 363         printf("Oh no: arrayPos = %i.\n", arrayPos);
 364       }
 365
 366       fNChanFills++; // one more channel found, and profile to be filled
 367       //NOTE: coordinates are (column, row) for the profiles
 368       if (gain == 0) {
 369         //fill the low gain histograms
 370         ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), min);
 371         ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
 372         ((TProfile2D*)fSampleLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), sample);
 373       }
 374       else if (gain == 1) {//fill the high gain ones
 375         ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), min);
 376         ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
 377         ((TProfile2D*)fSampleHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), sample);
 378       }//end if gain
 379
 380       max = fgkSampleMin; min = fgkSampleMax;
 381       i = 0;
 382
 383     }//End if end of tower
 384
 385   }//end while, of stream
 386
 387   return kTRUE;
 388 }
 389
 390 //_____________________________________________________________________
 391 Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
 392 {
 393   //Saves all the histograms (or profiles, to be accurate) to the designated file
 394
 395   TFile destFile(fileName, "recreate");
 396
 397   if (destFile.IsZombie()) {
 398     return kFALSE;
 399   }
 400
 401   destFile.cd();
 402
 403   for (int i = 0; i < fModules; i++) {
 404     if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
 405       fPeakMinusPedLowGain[i]->Write();
 406     }
 407     if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
 408       fPeakMinusPedHighGain[i]->Write();
 409     }
 410     if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
 411       fPedestalLowGain[i]->Write();
 412     }
 413     if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
 414       fPedestalHighGain[i]->Write();
 415     }
 416     if( ((TProfile2D *)fSampleLowGain[i])->GetEntries() || saveEmptyHistos) {
 417       fSampleLowGain[i]->Write();
 418     }
 419     if( ((TProfile2D *)fSampleHighGain[i])->GetEntries() || saveEmptyHistos) {
 420       fSampleHighGain[i]->Write();
 421     }
 422   }
 423
 424   destFile.Close();
 425
 426   return kTRUE;
 427 }
 428
 429 //_____________________________________________________________________
 430 Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
 431 {
 432
 433   //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
 434   TH1::AddDirectory(kFALSE);
 435
 436   TFile *sourceFile = new TFile(fileName);
 437   if (sourceFile->IsZombie()) {
 438     return kFALSE;//We couldn't load the reference
 439   }
 440
 441   if (fReference) delete fReference;//Delete the reference object, if it already exists
 442   fReference = 0;
 443
 444   fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
 445
 446   if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
 447     if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
 448     fReference = 0;
 449     return kFALSE;
 450   }
 451
 452   delete sourceFile;
 453
 454   //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
 455   //if we are called by someone who has set it to false...
 456   TH1::AddDirectory(kTRUE);
 457
 458   return kTRUE;//We succesfully loaded the object
 459 }
 460
 461 //_____________________________________________________________________
 462 void AliCaloCalibPedestal::ValidateComparisonProfiles()
 463 {
 464   //Make sure the comparison histos exist
 465   if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
 466   //the same time
 467
 468
 469   //Then, loop for the requested number of modules
 470   TString title, name;
 471   for (int i = 0; i < fModules; i++) {
 472     //Pedestals, low gain
 473     name = "hPedlowgainDiff";
 474     name += i;
 475     title = "Pedestals difference, low gain, module ";
 476     title += i;
 477     fPedestalLowGainDiff.Add(new TProfile2D(name, title,
 478                                             fColumns, 0.0, fColumns,
 479                                             fRows, fRowMin, fRowMax,"s"));
 480
 481     //Pedestals, high gain
 482     name = "hPedhighgainDiff";
 483     name += i;
 484     title = "Pedestals difference, high gain, module ";
 485     title += i;
 486     fPedestalHighGainDiff.Add(new TProfile2D(name, title,
 487                                              fColumns, 0.0, fColumns,
 488                                              fRows, fRowMin, fRowMax,"s"));
 489
 490     //Peak-Pedestals, high gain
 491     name = "hPeakMinusPedhighgainDiff";
 492     name += i;
 493     title = "Peak-Pedestal difference, high gain, module ";
 494     title += i;
 495     fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
 496                                                  fColumns, 0.0, fColumns,
 497                                                  fRows, fRowMin, fRowMax,"s"));
 498
 499     //Pedestals, low gain
 500     name = "hPedlowgainRatio";
 501     name += i;
 502     title = "Pedestals ratio, low gain, module ";
 503     title += i;
 504     fPedestalLowGainRatio.Add(new TProfile2D(name, title,
 505                                              fColumns, 0.0, fColumns,
 506                                              fRows, fRowMin, fRowMax,"s"));
 507
 508     //Pedestals, high gain
 509     name = "hPedhighgainRatio";
 510     name += i;
 511     title = "Pedestals ratio, high gain, module ";
 512     title += i;
 513     fPedestalHighGainRatio.Add(new TProfile2D(name, title,
 514                                               fColumns, 0.0, fColumns,
 515                                               fRows, fRowMin, fRowMax,"s"));
 516
 517     //Peak-Pedestals, low gain
 518     name = "hPeakMinusPedlowgainRatio";
 519     name += i;
 520     title = "Peak-Pedestal ratio, low gain, module ";
 521     title += i;
 522     fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
 523                                                  fColumns, 0.0, fColumns,
 524                                                  fRows, fRowMin, fRowMax,"s"));
 525
 526     //Peak-Pedestals, high gain
 527     name = "hPeakMinusPedhighgainRatio";
 528     name += i;
 529     title = "Peak-Pedestal ratio, high gain, module ";
 530     title += i;
 531     fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
 532                                                   fColumns, 0.0, fColumns,
 533                                                   fRows, fRowMin, fRowMax,"s"));
 534
 535   }//end for nModules create the histograms
 536 }
 537
 538 //_____________________________________________________________________
 539 void AliCaloCalibPedestal::ComputeDiffAndRatio()
 540 {
 541   // calculate differences and ratios relative to a reference
 542   ValidateComparisonProfiles();//Make sure the comparison histos exist
 543
 544   if (!fReference) {
 545     return;//Return if the reference object isn't loaded
 546   }
 547
 548   for (int i = 0; i < fModules; i++) {
 549     //Compute the ratio of the histograms
 550
 551     ((TProfile2D*)fPedestalLowGainRatio[i])->Divide(GetPedProfileLowGain(i), fReference->GetPedProfileLowGain(i), 1.0, 1.0);
 552     ((TProfile2D*)fPedestalHighGainRatio[i])->Divide(GetPedProfileHighGain(i), fReference->GetPedProfileHighGain(i), 1.0, 1.0);
 553     ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->Divide(GetPeakProfileLowGain(i), fReference->GetPeakProfileLowGain(i), 1.0, 1.0);
 554     ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->Divide(GetPeakProfileHighGain(i), fReference->GetPeakProfileHighGain(i), 1.0, 1.0);
 555
 556     //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
 557     //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
 558     for (int j = 0; j <= fColumns; j++) {
 559       for (int k = 0; k <= fRows; k++) {
 560         int bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
 561         double diff = fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) - GetPeakProfileHighGain(i)->GetBinContent(bin);
 562         ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
 563         ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
 564
 565         diff = fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) - GetPeakProfileLowGain(i)->GetBinContent(bin);
 566         ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
 567         ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
 568
 569         diff = fReference->GetPedProfileHighGain(i)->GetBinContent(bin) - GetPedProfileHighGain(i)->GetBinContent(bin);
 570         ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
 571         ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
 572
 573         diff = fReference->GetPedProfileLowGain(i)->GetBinContent(bin) - GetPedProfileLowGain(i)->GetBinContent(bin);
 574         ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
 575         ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
 576
 577       } // rows
 578     } // columns
 579
 580   } // modules
 581
 582 }
 583
 584 //_____________________________________________________________________
 585 void AliCaloCalibPedestal::ComputeDeadTowers(int threshold, const char * deadMapFile)
 586 {
 587   //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
 588   int countTot = 0;
 589   int countNew = 0;
 590   int countRes = 0;
 591   ofstream * fout = 0;
 592   ofstream * diff = 0;
 593   char name[512];//Quite a long temp buffer, just in case the filename includes a path
 594
 595   if (deadMapFile) {
 596     snprintf(name, 512, "%s.txt", deadMapFile);
 597     fout = new ofstream(name);
 598     snprintf(name, 512, "%sdiff.txt", deadMapFile);
 599     diff = new ofstream(name);
 600     if (!fout->is_open()) {
 601       delete fout;
 602       fout = 0;//Set the pointer to empty if the file was not opened
 603     }
 604     if (!diff->is_open()) {
 605       delete diff;
 606       fout = 0;//Set the pointer to empty if the file was not opened
 607     }
 608   }
 609
 610   for (int i = 0; i < fModules; i++) {
 611     if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
 612       for (int j = 1; j <= fColumns; j++) {
 613         for (int k = 1; k <= fRows; k++) {
 614
 615           if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {//It's dead
 616             countTot++;//One more dead total
 617             if (fout) {
 618               (*fout) << i << " "
 619                       << (fRows - k) << " "
 620                       << (j-1) << " "
 621                       << "1" << " "
 622                       << "0" << endl;//Write the status to the deadmap file, if the file is open.
 623             }
 624
 625             if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= threshold) {
 626               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
 627               countNew++;//This tower wasn't dead before!
 628               if (diff) {
 629                 ( *diff) << i << " "
 630                          << (fRows - k) << " "
 631                          << (j - 1) << " "
 632                          << "1" << " "
 633                          << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
 634               }
 635             }
 636             else {
 637               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
 638             }
 639           }
 640           else { //It's ALIVE!!
 641             //Don't bother with writing the live ones.
 642             //if (fout)
 643             //  (*fout) << i << " "
 644             //     << (fRows - k) << " "
 645             //     << (j - 1) << " "
 646             //     << "1" << " "
 647             //     << "1" << endl;//Write the status to the deadmap file, if the file is open.
 648             if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {
 649               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
 650               countRes++; //This tower was dead before => it's a miracle! :P
 651               if (diff) {
 652                 (*diff) << i << " "
 653                         << (fRows - k) << " "
 654                         << (j - 1) << " "
 655                         << "1" << " "
 656                         << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
 657               }
 658             }
 659             else {
 660               ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
 661             }
 662           }
 663
 664         }//end for k/rows
 665       }//end for j/columns
 666     }//end if GetEntries >= 0
 667
 668   }//end for modules
 669
 670  if (fout) {
 671    fout->close();
 672    delete fout;
 673  }
 674
 675  fDeadTowers = countTot;
 676  fNewDeadTowers = countNew;
 677  fResurrectedTowers = countRes;
 678 }
 679