1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 //________________________________________________________________________
19 // A help class for monitoring and calibration tools: MOOD, AMORE etc.,
20 // It can be created and used a la (ctor):
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();
28 // fPedestals->ProcessEvent(fCaloRawStream);
29 // asked to draw histograms:
30 // fPedestals->GetDeadMap(i)->Draw("col");
32 // fPedestals->GetPeakProfileHighGainRatio((i < fNumModules) ? i : fVisibleModule)->Draw("colz");
34 // The pseudo-code examples above were from the first implementation in MOOD (summer 2007).
35 //________________________________________________________________________
44 #include "AliCaloRawStreamV3.h"
47 #include "AliCaloCalibPedestal.h"
49 ClassImp(AliCaloCalibPedestal)
53 // ctor; initialize everything in order to avoid compiler warnings
54 AliCaloCalibPedestal::AliCaloCalibPedestal(kDetType detectorType) :
58 fPedestalRMSLowGain(),
59 fPedestalRMSHighGain(),
60 fPeakMinusPedLowGain(),
61 fPeakMinusPedHighGain(),
62 fPedestalLowGainDiff(),
63 fPedestalHighGainDiff(),
64 fPeakMinusPedLowGainDiff(),
65 fPeakMinusPedHighGainDiff(),
66 fPedestalLowGainRatio(),
67 fPedestalHighGainRatio(),
68 fPeakMinusPedLowGainRatio(),
69 fPeakMinusPedHighGainRatio(),
75 fResurrectedTowers(0),
87 fSelectPedestalSamples(kTRUE),
88 fFirstPedestalSample(0),
89 fLastPedestalSample(15)
91 //Default constructor. First we set the detector-type related constants.
92 if (detectorType == kPhos) {
93 fColumns = fgkPhosCols;
95 fModules = fgkPhosModules;
102 //We'll just trust the enum to keep everything in line, so that if detectorType
103 //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
104 //case, if someone intentionally gives another number
105 fColumns = AliEMCALGeoParams::fgkEMCALCols;
106 fRows = AliEMCALGeoParams::fgkEMCALRows;
107 fModules = AliEMCALGeoParams::fgkEMCALModules;
108 fCaloString = "EMCAL";
113 fDetType = detectorType;
115 //Then, loop for the requested number of modules
117 for (int i = 0; i < fModules; i++) {
118 //Pedestals, low gain
119 name = "hPedlowgain";
121 title = "Pedestals, low gain, module ";
123 fPedestalLowGain.Add(new TProfile2D(name, title,
124 fColumns, 0.0, fColumns,
125 fRows, fRowMin, fRowMax,"s"));
127 //Pedestals, high gain
128 name = "hPedhighgain";
130 title = "Pedestals, high gain, module ";
132 fPedestalHighGain.Add(new TProfile2D(name, title,
133 fColumns, 0.0, fColumns,
134 fRows, fRowMin, fRowMax,"s"));
135 //All Samples, low gain
136 name = "hPedestalRMSlowgain";
138 title = "Pedestal RMS, low gain, module ";
140 fPedestalRMSLowGain.Add(new TProfile2D(name, title,
141 fColumns, 0.0, fColumns,
142 fRows, fRowMin, fRowMax,"s"));
144 //All Samples, high gain
145 name = "hPedestalRMShighgain";
147 title = "Pedestal RMS, high gain, module ";
149 fPedestalRMSHighGain.Add(new TProfile2D(name, title,
150 fColumns, 0.0, fColumns,
151 fRows, fRowMin, fRowMax,"s"));
153 //Peak-Pedestals, low gain
154 name = "hPeakMinusPedlowgain";
156 title = "Peak-Pedestal, low gain, module ";
158 fPeakMinusPedLowGain.Add(new TProfile2D(name, title,
159 fColumns, 0.0, fColumns,
160 fRows, fRowMin, fRowMax,"s"));
162 //Peak-Pedestals, high gain
163 name = "hPeakMinusPedhighgain";
165 title = "Peak-Pedestal, high gain, module ";
167 fPeakMinusPedHighGain.Add(new TProfile2D(name, title,
168 fColumns, 0.0, fColumns,
169 fRows, fRowMin, fRowMax,"s"));
173 title = "Dead map, module ";
175 fDeadMap.Add(new TH2D(name, title, fColumns, 0.0, fColumns,
176 fRows, fRowMin, fRowMax));
178 }//end for nModules create the histograms
180 //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
181 fPedestalLowGain.Compress();
182 fPedestalHighGain.Compress();
183 fPedestalRMSLowGain.Compress();
184 fPedestalRMSHighGain.Compress();
185 fPeakMinusPedLowGain.Compress();
186 fPeakMinusPedHighGain.Compress();
188 //Make them the owners of the profiles, so we don't need to care about deleting them
189 //fPedestalLowGain.SetOwner();
190 //fPedestalHighGain.SetOwner();
191 //fPeakMinusPedLowGain.SetOwner();
192 //fPeakMinusPedHighGain.SetOwner();
197 //_____________________________________________________________________
198 AliCaloCalibPedestal::~AliCaloCalibPedestal()
200 if (fReference) delete fReference;//Delete the reference object, if it has been loaded
201 //TObjArray will delete the histos/profiles when it is deleted.
205 //_____________________________________________________________________
206 AliCaloCalibPedestal::AliCaloCalibPedestal(const AliCaloCalibPedestal &ped) :
210 fPedestalRMSLowGain(),
211 fPedestalRMSHighGain(),
212 fPeakMinusPedLowGain(),
213 fPeakMinusPedHighGain(),
214 fPedestalLowGainDiff(),
215 fPedestalHighGainDiff(),
216 fPeakMinusPedLowGainDiff(),
217 fPeakMinusPedHighGainDiff(),
218 fPedestalLowGainRatio(),
219 fPedestalHighGainRatio(),
220 fPeakMinusPedLowGainRatio(),
221 fPeakMinusPedHighGainRatio(),
223 fNEvents(ped.GetNEvents()),
224 fNChanFills(ped.GetNChanFills()),
225 fDeadTowers(ped.GetDeadTowerCount()),
226 fNewDeadTowers(ped.GetDeadTowerNew()),
227 fResurrectedTowers(ped.GetDeadTowerResurrected()),
228 fReference( 0 ), //! note that we do not try to copy the reference info here
229 fDetType(ped.GetDetectorType()),
230 fColumns(ped.GetColumns()),
231 fRows(ped.GetRows()),
232 fModules(ped.GetModules()),
233 fRowMin(ped.GetRowMin()),
234 fRowMax(ped.GetRowMax()),
235 fRowMultiplier(ped.GetRowMultiplier()),
236 fCaloString(ped.GetCaloString()),
237 fMapping(NULL), //! note that we are not copying the map info
238 fRunNumber(ped.GetRunNumber()),
239 fSelectPedestalSamples(ped.GetSelectPedestalSamples()),
240 fFirstPedestalSample(ped.GetFirstPedestalSample()),
241 fLastPedestalSample(ped.GetLastPedestalSample())
243 // Then the ObjArray ones; we add the histograms rather than trying TObjArray = assignment
244 //DS: this has not really been tested yet..
245 for (int i = 0; i < fModules; i++) {
246 fPedestalLowGain.Add( ped.GetPedProfileLowGain(i) );
247 fPedestalHighGain.Add( ped.GetPedProfileHighGain(i) );
248 fPedestalRMSLowGain.Add( ped.GetPedRMSProfileLowGain(i) );
249 fPedestalRMSHighGain.Add( ped.GetPedRMSProfileHighGain(i) );
250 fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
251 fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
253 fDeadMap.Add( ped.GetDeadMap(i) );
256 //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
257 fPedestalLowGain.Compress();
258 fPedestalHighGain.Compress();
259 fPedestalRMSLowGain.Compress();
260 fPedestalRMSHighGain.Compress();
261 fPeakMinusPedLowGain.Compress();
262 fPeakMinusPedHighGain.Compress();
266 // assignment operator; use copy ctor to make life easy..
267 //_____________________________________________________________________
268 AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (const AliCaloCalibPedestal &source)
270 // assignment operator; use copy ctor
271 if (&source == this) return *this;
273 new (this) AliCaloCalibPedestal(source);
277 //_____________________________________________________________________
278 void AliCaloCalibPedestal::Reset()
280 // Reset all arrays/histograms
281 for (int i = 0; i < fModules; i++) {
282 GetPedProfileLowGain(i)->Reset();
283 GetPedProfileHighGain(i)->Reset();
284 GetPeakProfileLowGain(i)->Reset();
285 GetPeakProfileHighGain(i)->Reset();
286 GetDeadMap(i)->Reset();
288 if (!fPedestalLowGainDiff.IsEmpty()) {
289 //This means that the comparison profiles have been created.
291 GetPedProfileLowGainDiff(i)->Reset();
292 GetPedProfileHighGainDiff(i)->Reset();
293 GetPeakProfileLowGainDiff(i)->Reset();
294 GetPeakProfileHighGainDiff(i)->Reset();
296 GetPedProfileLowGainRatio(i)->Reset();
297 GetPedProfileHighGainRatio(i)->Reset();
298 GetPeakProfileLowGainRatio(i)->Reset();
299 GetPeakProfileHighGainRatio(i)->Reset();
306 fResurrectedTowers = 0;
308 //To think about: should fReference be deleted too?... let's not do it this time, at least...
311 // Parameter/cut handling
312 //_____________________________________________________________________
313 void AliCaloCalibPedestal::SetParametersFromFile(const char *parameterFile)
315 // Note: this method is a bit more complicated than it really has to be
316 // - allowing for multiple entries per line, arbitrary order of the
317 // different variables etc. But I wanted to try and do this in as
318 // correct a C++ way as I could (as an exercise).
320 static const string delimitor("::");
322 // open, check input file
323 ifstream in( parameterFile );
325 printf("in AliCaloCalibPedestal::SetParametersFromFile - Using default/run_time parameters.\n");
332 while ((in.rdstate() & ios::failbit) == 0 ) {
334 // Read into the raw char array and then construct a string
335 // to do the searching
336 in.getline(readline, 1024);
337 istringstream s(readline);
339 while ( ( s.rdstate() & ios::failbit ) == 0 ) {
344 // check stream status
345 if( s.rdstate() & ios::failbit ) break;
347 // skip rest of line if comments found
348 if( keyValue.substr( 0, 2 ) == "//" ) break;
350 // look for "::" in keyValue pair
351 size_t position = keyValue.find( delimitor );
352 if( position == string::npos ) {
353 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
356 // split keyValue pair
357 string key( keyValue.substr( 0, position ) );
358 string value( keyValue.substr( position+delimitor.size(),
359 keyValue.size()-delimitor.size() ) );
361 // check value does not contain a new delimitor
362 if( value.find( delimitor ) != string::npos ) {
363 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
366 // debug: check key value pair
367 // printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
369 // if the key matches with something we expect, we assign the new value
370 istringstream iss(value);
371 // the comparison strings defined at the beginning of this method
372 if ( (key == "fFirstPedestalSample") || (key == "fLastPedestalSample") ) {
373 printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
375 if (key == "fFirstPedestalSample") {
376 iss >> fFirstPedestalSample;
378 else if (key == "fLastPedestalSample") {
379 iss >> fLastPedestalSample;
391 //_____________________________________________________________________
392 void AliCaloCalibPedestal::WriteParametersToFile(const char *parameterFile)
394 //Write parameters in file.
396 static const string delimitor("::");
397 ofstream out( parameterFile );
398 out << "// " << parameterFile << endl;
399 out << "fFirstPedestalSample" << "::" << fFirstPedestalSample << endl;
400 out << "fLastPedestalSample" << "::" << fLastPedestalSample << endl;
406 //_____________________________________________________________________
407 Bool_t AliCaloCalibPedestal::AddInfo(const AliCaloCalibPedestal *ped)
409 // just do this for the basic histograms/profiles that get filled in ProcessEvent
410 // may not have data for all modules, but let's just Add everything..
411 for (int i = 0; i < fModules; i++) {
412 GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
413 GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
414 GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
415 GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
418 // DeadMap; Diff profiles etc would need to be redone after this operation
420 return kTRUE;//We succesfully added info from the supplied object
423 //_____________________________________________________________________
424 Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
426 // if fMapping is NULL the rawstream will crate its own mapping
427 AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
428 return ProcessEvent(&rawStream);
431 //_____________________________________________________________________
432 Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStreamV3 *in)
434 // Method to process=analyze one event in the data stream
435 if (!in) return kFALSE; //Return right away if there's a null pointer
436 fNEvents++; // one more event
438 // indices for the reading
442 int i = 0; // sample counter
445 // start loop over input stream
446 while (in->NextDDL()) {
447 while (in->NextChannel()) {
450 int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
453 // for the pedestal calculation
454 int sampleSum = 0; // sum of samples
455 int squaredSampleSum = 0; // sum of samples squared
456 int nSum = 0; // number of samples in sum
458 double mean = 0, squaredMean = 0, rms = 0;
460 while (in->NextBunch()) {
461 const UShort_t *sig = in->GetSignals();
462 startBin = in->GetStartTimeBin();
463 nsamples += in->GetBunchLength();
464 for (i = 0; i < in->GetBunchLength(); i++) {
468 // check if it's a min or max value
469 if (sample < min) min = sample;
470 if (sample > max) max = sample;
472 // should we add it for the pedestal calculation?
473 if ( (fFirstPedestalSample<=time && time<=fLastPedestalSample) || // sample time in range
474 !fSelectPedestalSamples ) { // or we don't restrict the sample range.. - then we'll take all
476 squaredSampleSum += sample*sample;
480 } // loop over samples in bunch
481 } // loop over bunches
483 if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
485 // calculate pedesstal estimate: mean of possibly selected samples
487 mean = sampleSum / (1.0 * nSum);
488 squaredMean = squaredSampleSum / (1.0 * nSum);
489 // The variance (rms squared) is equal to the mean of the squares minus the square of the mean..
490 rms = sqrt(squaredMean - mean*mean);
498 // we're done with the calc. for this channel; let's prepare to fill histo
499 gain = -1; // init to not valid value
500 if ( in->IsLowGain() ) {
503 else if ( in->IsHighGain() ) {
507 // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
508 int arrayPos = in->GetModule(); //The modules are numbered starting from 0
509 if (arrayPos >= fModules) {
510 //TODO: return an error message, if appopriate (perhaps if debug>0?)
514 if (arrayPos < 0 || arrayPos >= fModules) {
515 printf("Oh no: arrayPos = %i.\n", arrayPos);
518 fNChanFills++; // one more channel found, and profile to be filled
519 //NOTE: coordinates are (column, row) for the profiles
521 //fill the low gain histograms
522 ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
523 if (nSum>0) { // only fill pedestal info in case it could be calculated
524 ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), mean);
525 ((TProfile2D*)fPedestalRMSLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), rms);
528 else if (gain == 1) {
529 //fill the high gain ones
530 ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
531 if (nSum>0) { // only fill pedestal info in case it could be calculated
532 ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), mean);
533 ((TProfile2D*)fPedestalRMSHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), rms);
537 } // nsamples>0 check, some data found for this channel; not only trailer/header
538 }// end while over channel
539 }//end while over DDL's, of input stream
541 in->Reset(); // just in case the next customer forgets to check if the stream was reset..
546 //_____________________________________________________________________
547 Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
549 //Saves all the histograms (or profiles, to be accurate) to the designated file
551 TFile destFile(fileName, "recreate");
553 if (destFile.IsZombie()) {
559 for (int i = 0; i < fModules; i++) {
560 if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
561 fPeakMinusPedLowGain[i]->Write();
563 if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
564 fPeakMinusPedHighGain[i]->Write();
566 if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
567 fPedestalLowGain[i]->Write();
569 if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
570 fPedestalHighGain[i]->Write();
571 Printf("save %d", i);
573 if( ((TProfile2D *)fPedestalRMSLowGain[i])->GetEntries() || saveEmptyHistos) {
574 fPedestalRMSLowGain[i]->Write();
576 if( ((TProfile2D *)fPedestalRMSHighGain[i])->GetEntries() || saveEmptyHistos) {
577 fPedestalRMSHighGain[i]->Write();
586 //_____________________________________________________________________
587 Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
590 //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
591 TH1::AddDirectory(kFALSE);
593 TFile *sourceFile = new TFile(fileName);
594 if (sourceFile->IsZombie()) {
595 return kFALSE;//We couldn't load the reference
598 if (fReference) delete fReference;//Delete the reference object, if it already exists
601 fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
603 if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
604 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
611 //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
612 //if we are called by someone who has set it to false...
613 TH1::AddDirectory(kTRUE);
615 return kTRUE;//We succesfully loaded the object
618 //_____________________________________________________________________
619 void AliCaloCalibPedestal::ValidateComparisonProfiles()
621 //Make sure the comparison histos exist
622 if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
626 //Then, loop for the requested number of modules
628 for (int i = 0; i < fModules; i++) {
629 //Pedestals, low gain
630 name = "hPedlowgainDiff";
632 title = "Pedestals difference, low gain, module ";
634 fPedestalLowGainDiff.Add(new TProfile2D(name, title,
635 fColumns, 0.0, fColumns,
636 fRows, fRowMin, fRowMax,"s"));
638 //Pedestals, high gain
639 name = "hPedhighgainDiff";
641 title = "Pedestals difference, high gain, module ";
643 fPedestalHighGainDiff.Add(new TProfile2D(name, title,
644 fColumns, 0.0, fColumns,
645 fRows, fRowMin, fRowMax,"s"));
647 //Peak-Pedestals, high gain
648 name = "hPeakMinusPedhighgainDiff";
650 title = "Peak-Pedestal difference, high gain, module ";
652 fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
653 fColumns, 0.0, fColumns,
654 fRows, fRowMin, fRowMax,"s"));
656 //Pedestals, low gain
657 name = "hPedlowgainRatio";
659 title = "Pedestals ratio, low gain, module ";
661 fPedestalLowGainRatio.Add(new TProfile2D(name, title,
662 fColumns, 0.0, fColumns,
663 fRows, fRowMin, fRowMax,"s"));
665 //Pedestals, high gain
666 name = "hPedhighgainRatio";
668 title = "Pedestals ratio, high gain, module ";
670 fPedestalHighGainRatio.Add(new TProfile2D(name, title,
671 fColumns, 0.0, fColumns,
672 fRows, fRowMin, fRowMax,"s"));
674 //Peak-Pedestals, low gain
675 name = "hPeakMinusPedlowgainRatio";
677 title = "Peak-Pedestal ratio, low gain, module ";
679 fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
680 fColumns, 0.0, fColumns,
681 fRows, fRowMin, fRowMax,"s"));
683 //Peak-Pedestals, high gain
684 name = "hPeakMinusPedhighgainRatio";
686 title = "Peak-Pedestal ratio, high gain, module ";
688 fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
689 fColumns, 0.0, fColumns,
690 fRows, fRowMin, fRowMax,"s"));
692 }//end for nModules create the histograms
695 //_____________________________________________________________________
696 void AliCaloCalibPedestal::ComputeDiffAndRatio()
698 // calculate differences and ratios relative to a reference
699 ValidateComparisonProfiles();//Make sure the comparison histos exist
702 return;//Return if the reference object isn't loaded
705 for (int i = 0; i < fModules; i++) {
706 //Compute the ratio of the histograms
708 ((TProfile2D*)fPedestalLowGainRatio[i])->Divide(GetPedProfileLowGain(i), fReference->GetPedProfileLowGain(i), 1.0, 1.0);
709 ((TProfile2D*)fPedestalHighGainRatio[i])->Divide(GetPedProfileHighGain(i), fReference->GetPedProfileHighGain(i), 1.0, 1.0);
710 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->Divide(GetPeakProfileLowGain(i), fReference->GetPeakProfileLowGain(i), 1.0, 1.0);
711 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->Divide(GetPeakProfileHighGain(i), fReference->GetPeakProfileHighGain(i), 1.0, 1.0);
713 //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
714 //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
715 for (int j = 0; j <= fColumns; j++) {
716 for (int k = 0; k <= fRows; k++) {
717 int bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
718 double diff = fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) - GetPeakProfileHighGain(i)->GetBinContent(bin);
719 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
720 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
722 diff = fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) - GetPeakProfileLowGain(i)->GetBinContent(bin);
723 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
724 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
726 diff = fReference->GetPedProfileHighGain(i)->GetBinContent(bin) - GetPedProfileHighGain(i)->GetBinContent(bin);
727 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
728 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
730 diff = fReference->GetPedProfileLowGain(i)->GetBinContent(bin) - GetPedProfileLowGain(i)->GetBinContent(bin);
731 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
732 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
741 //_____________________________________________________________________
742 void AliCaloCalibPedestal::ComputeDeadTowers(int threshold, const char * deadMapFile)
744 //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
750 char name[512];//Quite a long temp buffer, just in case the filename includes a path
753 snprintf(name, 512, "%s.txt", deadMapFile);
754 fout = new ofstream(name);
755 snprintf(name, 512, "%sdiff.txt", deadMapFile);
756 diff = new ofstream(name);
757 if (!fout->is_open()) {
759 fout = 0;//Set the pointer to empty if the file was not opened
761 if (!diff->is_open()) {
763 fout = 0;//Set the pointer to empty if the file was not opened
767 for (int i = 0; i < fModules; i++) {
768 if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
769 for (int j = 1; j <= fColumns; j++) {
770 for (int k = 1; k <= fRows; k++) {
772 if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {//It's dead
773 countTot++;//One more dead total
776 << (fRows - k) << " "
779 << "0" << endl;//Write the status to the deadmap file, if the file is open.
782 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= threshold) {
783 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
784 countNew++;//This tower wasn't dead before!
787 << (fRows - k) << " "
790 << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
794 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
797 else { //It's ALIVE!!
798 //Don't bother with writing the live ones.
800 // (*fout) << i << " "
801 // << (fRows - k) << " "
804 // << "1" << endl;//Write the status to the deadmap file, if the file is open.
805 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {
806 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
807 countRes++; //This tower was dead before => it's a miracle! :P
810 << (fRows - k) << " "
813 << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
817 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
823 }//end if GetEntries >= 0
832 fDeadTowers = countTot;
833 fNewDeadTowers = countNew;
834 fResurrectedTowers = countRes;
837 //_____________________________________________________________________
838 Bool_t AliCaloCalibPedestal::IsBadChannel(int imod, int icol, int irow) const
840 //Check if channel is dead or hot.
842 Int_t status = ((TH2D*)fDeadMap[imod])->GetBinContent(icol,irow);
850 //_____________________________________________________________________
851 void AliCaloCalibPedestal::SetChannelStatus(int imod, int icol, int irow, int status)
853 //Set status of channel dead, hot, alive ...
855 ((TH2D*)fDeadMap[imod])->SetBinContent(icol,irow, status);