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 //________________________________________________________________________
43 #include "AliCaloRawStream.h"
46 #include "AliCaloCalibPedestal.h"
48 ClassImp(AliCaloCalibPedestal)
52 // ctor; initialize everything in order to avoid compiler warnings
53 AliCaloCalibPedestal::AliCaloCalibPedestal(kDetType detectorType) :
57 fPedestalRMSLowGain(),
58 fPedestalRMSHighGain(),
59 fPeakMinusPedLowGain(),
60 fPeakMinusPedHighGain(),
61 fPedestalLowGainDiff(),
62 fPedestalHighGainDiff(),
63 fPeakMinusPedLowGainDiff(),
64 fPeakMinusPedHighGainDiff(),
65 fPedestalLowGainRatio(),
66 fPedestalHighGainRatio(),
67 fPeakMinusPedLowGainRatio(),
68 fPeakMinusPedHighGainRatio(),
74 fResurrectedTowers(0),
86 fSelectPedestalSamples(kTRUE),
87 fFirstPedestalSample(0),
88 fLastPedestalSample(15)
90 //Default constructor. First we set the detector-type related constants.
91 if (detectorType == kPhos) {
92 fColumns = fgkPhosCols;
94 fModules = fgkPhosModules;
101 //We'll just trust the enum to keep everything in line, so that if detectorType
102 //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
103 //case, if someone intentionally gives another number
104 fColumns = fgkEmCalCols;
105 fRows = fgkEmCalRows;
106 fModules = fgkEmCalModules;
107 fCaloString = "EMCAL";
112 fDetType = detectorType;
114 //Then, loop for the requested number of modules
116 for (int i = 0; i < fModules; i++) {
117 //Pedestals, low gain
118 name = "hPedlowgain";
120 title = "Pedestals, low gain, module ";
122 fPedestalLowGain.Add(new TProfile2D(name, title,
123 fColumns, 0.0, fColumns,
124 fRows, fRowMin, fRowMax,"s"));
126 //Pedestals, high gain
127 name = "hPedhighgain";
129 title = "Pedestals, high gain, module ";
131 fPedestalHighGain.Add(new TProfile2D(name, title,
132 fColumns, 0.0, fColumns,
133 fRows, fRowMin, fRowMax,"s"));
134 //All Samples, low gain
135 name = "hPedestalRMSlowgain";
137 title = "Pedestal RMS, low gain, module ";
139 fPedestalRMSLowGain.Add(new TProfile2D(name, title,
140 fColumns, 0.0, fColumns,
141 fRows, fRowMin, fRowMax,"s"));
143 //All Samples, high gain
144 name = "hPedestalRMShighgain";
146 title = "Pedestal RMS, high gain, module ";
148 fPedestalRMSHighGain.Add(new TProfile2D(name, title,
149 fColumns, 0.0, fColumns,
150 fRows, fRowMin, fRowMax,"s"));
152 //Peak-Pedestals, low gain
153 name = "hPeakMinusPedlowgain";
155 title = "Peak-Pedestal, low gain, module ";
157 fPeakMinusPedLowGain.Add(new TProfile2D(name, title,
158 fColumns, 0.0, fColumns,
159 fRows, fRowMin, fRowMax,"s"));
161 //Peak-Pedestals, high gain
162 name = "hPeakMinusPedhighgain";
164 title = "Peak-Pedestal, high gain, module ";
166 fPeakMinusPedHighGain.Add(new TProfile2D(name, title,
167 fColumns, 0.0, fColumns,
168 fRows, fRowMin, fRowMax,"s"));
172 title = "Dead map, module ";
174 fDeadMap.Add(new TH2D(name, title, fColumns, 0.0, fColumns,
175 fRows, fRowMin, fRowMax));
177 }//end for nModules create the histograms
179 //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
180 fPedestalLowGain.Compress();
181 fPedestalHighGain.Compress();
182 fPedestalRMSLowGain.Compress();
183 fPedestalRMSHighGain.Compress();
184 fPeakMinusPedLowGain.Compress();
185 fPeakMinusPedHighGain.Compress();
187 //Make them the owners of the profiles, so we don't need to care about deleting them
188 //fPedestalLowGain.SetOwner();
189 //fPedestalHighGain.SetOwner();
190 //fPeakMinusPedLowGain.SetOwner();
191 //fPeakMinusPedHighGain.SetOwner();
196 //_____________________________________________________________________
197 AliCaloCalibPedestal::~AliCaloCalibPedestal()
199 if (fReference) delete fReference;//Delete the reference object, if it has been loaded
200 //TObjArray will delete the histos/profiles when it is deleted.
204 //_____________________________________________________________________
205 AliCaloCalibPedestal::AliCaloCalibPedestal(const AliCaloCalibPedestal &ped) :
209 fPedestalRMSLowGain(),
210 fPedestalRMSHighGain(),
211 fPeakMinusPedLowGain(),
212 fPeakMinusPedHighGain(),
213 fPedestalLowGainDiff(),
214 fPedestalHighGainDiff(),
215 fPeakMinusPedLowGainDiff(),
216 fPeakMinusPedHighGainDiff(),
217 fPedestalLowGainRatio(),
218 fPedestalHighGainRatio(),
219 fPeakMinusPedLowGainRatio(),
220 fPeakMinusPedHighGainRatio(),
222 fNEvents(ped.GetNEvents()),
223 fNChanFills(ped.GetNChanFills()),
224 fDeadTowers(ped.GetDeadTowerCount()),
225 fNewDeadTowers(ped.GetDeadTowerNew()),
226 fResurrectedTowers(ped.GetDeadTowerResurrected()),
227 fReference( 0 ), //! note that we do not try to copy the reference info here
228 fDetType(ped.GetDetectorType()),
229 fColumns(ped.GetColumns()),
230 fRows(ped.GetRows()),
231 fModules(ped.GetModules()),
232 fRowMin(ped.GetRowMin()),
233 fRowMax(ped.GetRowMax()),
234 fRowMultiplier(ped.GetRowMultiplier()),
235 fCaloString(ped.GetCaloString()),
236 fMapping(NULL), //! note that we are not copying the map info
237 fRunNumber(ped.GetRunNumber()),
238 fSelectPedestalSamples(ped.GetSelectPedestalSamples()),
239 fFirstPedestalSample(ped.GetFirstPedestalSample()),
240 fLastPedestalSample(ped.GetLastPedestalSample())
242 // Then the ObjArray ones; we add the histograms rather than trying TObjArray = assignment
243 //DS: this has not really been tested yet..
244 for (int i = 0; i < fModules; i++) {
245 fPedestalLowGain.Add( ped.GetPedProfileLowGain(i) );
246 fPedestalHighGain.Add( ped.GetPedProfileHighGain(i) );
247 fPedestalRMSLowGain.Add( ped.GetPedRMSProfileLowGain(i) );
248 fPedestalRMSHighGain.Add( ped.GetPedRMSProfileHighGain(i) );
249 fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
250 fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
252 fDeadMap.Add( ped.GetDeadMap(i) );
255 //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
256 fPedestalLowGain.Compress();
257 fPedestalHighGain.Compress();
258 fPedestalRMSLowGain.Compress();
259 fPedestalRMSHighGain.Compress();
260 fPeakMinusPedLowGain.Compress();
261 fPeakMinusPedHighGain.Compress();
265 // assignment operator; use copy ctor to make life easy..
266 //_____________________________________________________________________
267 AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (const AliCaloCalibPedestal &source)
269 // assignment operator; use copy ctor
270 if (&source == this) return *this;
272 new (this) AliCaloCalibPedestal(source);
276 //_____________________________________________________________________
277 void AliCaloCalibPedestal::Reset()
279 // Reset all arrays/histograms
280 for (int i = 0; i < fModules; i++) {
281 GetPedProfileLowGain(i)->Reset();
282 GetPedProfileHighGain(i)->Reset();
283 GetPeakProfileLowGain(i)->Reset();
284 GetPeakProfileHighGain(i)->Reset();
285 GetDeadMap(i)->Reset();
287 if (!fPedestalLowGainDiff.IsEmpty()) {
288 //This means that the comparison profiles have been created.
290 GetPedProfileLowGainDiff(i)->Reset();
291 GetPedProfileHighGainDiff(i)->Reset();
292 GetPeakProfileLowGainDiff(i)->Reset();
293 GetPeakProfileHighGainDiff(i)->Reset();
295 GetPedProfileLowGainRatio(i)->Reset();
296 GetPedProfileHighGainRatio(i)->Reset();
297 GetPeakProfileLowGainRatio(i)->Reset();
298 GetPeakProfileHighGainRatio(i)->Reset();
305 fResurrectedTowers = 0;
307 //To think about: should fReference be deleted too?... let's not do it this time, at least...
310 //_____________________________________________________________________
311 Bool_t AliCaloCalibPedestal::AddInfo(const AliCaloCalibPedestal *ped)
313 // just do this for the basic histograms/profiles that get filled in ProcessEvent
314 // may not have data for all modules, but let's just Add everything..
315 for (int i = 0; i < fModules; i++) {
316 GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
317 GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
318 GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
319 GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
322 // DeadMap; Diff profiles etc would need to be redone after this operation
324 return kTRUE;//We succesfully added info from the supplied object
327 //_____________________________________________________________________
328 Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
330 // if fMapping is NULL the rawstream will crate its own mapping
331 AliCaloRawStream rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
332 return ProcessEvent(&rawStream);
335 //_____________________________________________________________________
336 Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStream *in)
338 // Method to process=analyze one event in the data stream
339 if (!in) return kFALSE; //Return right away if there's a null pointer
340 fNEvents++; // one more event
342 // indices for the reading
347 int i = 0; // the sample number in current event.
348 int max = fgkSampleMin, min = fgkSampleMax;//Use these for picking the pedestal
350 // for the pedestal calculation
351 int sampleSum = 0; // sum of samples
352 int squaredSampleSum = 0; // sum of samples squared
353 int nSum = 0; // number of samples in sum
355 double mean = 0, squaredMean = 0, rms = 0;
358 sample = in->GetSignal(); //Get the adc signal
359 if (sample < min) min = sample;
360 if (sample > max) max = sample;
363 // should we add it for the pedestal calculation?
364 time = in->GetTime();
365 if ( (fFirstPedestalSample<=time && time<=fLastPedestalSample) || // sample time in range
366 !fSelectPedestalSamples ) { // or we don't restrict the sample range.. - then we'll take all
368 squaredSampleSum += sample*sample;
372 if ( i >= in->GetTimeLength()) {
373 // calculate pedesstal estimate: mean of possibly selected samples
375 mean = sampleSum / (1.0 * nSum);
376 squaredMean = squaredSampleSum / (1.0 * nSum);
377 // The variance (rms squared) is equal to the mean of the squares minus the square of the mean..
378 rms = sqrt(squaredMean - mean*mean);
386 //If we're here then we're done with this tower
387 gain = -1; // init to not valid value
388 if ( in->IsLowGain() ) {
391 else if ( in->IsHighGain() ) {
395 int arrayPos = in->GetModule(); //The modules are numbered starting from 0
396 if (arrayPos >= fModules) {
397 //TODO: return an error message, if appopriate (perhaps if debug>0?)
402 if (arrayPos < 0 || arrayPos >= fModules) {
403 printf("Oh no: arrayPos = %i.\n", arrayPos);
406 fNChanFills++; // one more channel found, and profile to be filled
407 //NOTE: coordinates are (column, row) for the profiles
409 //fill the low gain histograms
410 ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
411 if (nSum>0) { // only fill pedestal info in case it could be calculated
412 ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), mean);
413 ((TProfile2D*)fPedestalRMSLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), rms);
416 else if (gain == 1) {
417 //fill the high gain ones
418 ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
419 if (nSum>0) { // only fill pedestal info in case it could be calculated
420 ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), mean);
421 ((TProfile2D*)fPedestalRMSHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), rms);
426 max = fgkSampleMin; min = fgkSampleMax;
428 // also pedestal calc counters
429 sampleSum = 0; // sum of samples
430 squaredSampleSum = 0; // sum of samples squared
431 nSum = 0; // number of samples in sum
433 }//End if end of tower
435 }//end while, of stream
440 //_____________________________________________________________________
441 Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
443 //Saves all the histograms (or profiles, to be accurate) to the designated file
445 TFile destFile(fileName, "recreate");
447 if (destFile.IsZombie()) {
453 for (int i = 0; i < fModules; i++) {
454 if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
455 fPeakMinusPedLowGain[i]->Write();
457 if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
458 fPeakMinusPedHighGain[i]->Write();
460 if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
461 fPedestalLowGain[i]->Write();
463 if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
464 fPedestalHighGain[i]->Write();
465 Printf("save %d", i);
467 if( ((TProfile2D *)fPedestalRMSLowGain[i])->GetEntries() || saveEmptyHistos) {
468 fPedestalRMSLowGain[i]->Write();
470 if( ((TProfile2D *)fPedestalRMSHighGain[i])->GetEntries() || saveEmptyHistos) {
471 fPedestalRMSHighGain[i]->Write();
480 //_____________________________________________________________________
481 Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
484 //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
485 TH1::AddDirectory(kFALSE);
487 TFile *sourceFile = new TFile(fileName);
488 if (sourceFile->IsZombie()) {
489 return kFALSE;//We couldn't load the reference
492 if (fReference) delete fReference;//Delete the reference object, if it already exists
495 fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
497 if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
498 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
505 //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
506 //if we are called by someone who has set it to false...
507 TH1::AddDirectory(kTRUE);
509 return kTRUE;//We succesfully loaded the object
512 //_____________________________________________________________________
513 void AliCaloCalibPedestal::ValidateComparisonProfiles()
515 //Make sure the comparison histos exist
516 if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
520 //Then, loop for the requested number of modules
522 for (int i = 0; i < fModules; i++) {
523 //Pedestals, low gain
524 name = "hPedlowgainDiff";
526 title = "Pedestals difference, low gain, module ";
528 fPedestalLowGainDiff.Add(new TProfile2D(name, title,
529 fColumns, 0.0, fColumns,
530 fRows, fRowMin, fRowMax,"s"));
532 //Pedestals, high gain
533 name = "hPedhighgainDiff";
535 title = "Pedestals difference, high gain, module ";
537 fPedestalHighGainDiff.Add(new TProfile2D(name, title,
538 fColumns, 0.0, fColumns,
539 fRows, fRowMin, fRowMax,"s"));
541 //Peak-Pedestals, high gain
542 name = "hPeakMinusPedhighgainDiff";
544 title = "Peak-Pedestal difference, high gain, module ";
546 fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
547 fColumns, 0.0, fColumns,
548 fRows, fRowMin, fRowMax,"s"));
550 //Pedestals, low gain
551 name = "hPedlowgainRatio";
553 title = "Pedestals ratio, low gain, module ";
555 fPedestalLowGainRatio.Add(new TProfile2D(name, title,
556 fColumns, 0.0, fColumns,
557 fRows, fRowMin, fRowMax,"s"));
559 //Pedestals, high gain
560 name = "hPedhighgainRatio";
562 title = "Pedestals ratio, high gain, module ";
564 fPedestalHighGainRatio.Add(new TProfile2D(name, title,
565 fColumns, 0.0, fColumns,
566 fRows, fRowMin, fRowMax,"s"));
568 //Peak-Pedestals, low gain
569 name = "hPeakMinusPedlowgainRatio";
571 title = "Peak-Pedestal ratio, low gain, module ";
573 fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
574 fColumns, 0.0, fColumns,
575 fRows, fRowMin, fRowMax,"s"));
577 //Peak-Pedestals, high gain
578 name = "hPeakMinusPedhighgainRatio";
580 title = "Peak-Pedestal ratio, high gain, module ";
582 fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
583 fColumns, 0.0, fColumns,
584 fRows, fRowMin, fRowMax,"s"));
586 }//end for nModules create the histograms
589 //_____________________________________________________________________
590 void AliCaloCalibPedestal::ComputeDiffAndRatio()
592 // calculate differences and ratios relative to a reference
593 ValidateComparisonProfiles();//Make sure the comparison histos exist
596 return;//Return if the reference object isn't loaded
599 for (int i = 0; i < fModules; i++) {
600 //Compute the ratio of the histograms
602 ((TProfile2D*)fPedestalLowGainRatio[i])->Divide(GetPedProfileLowGain(i), fReference->GetPedProfileLowGain(i), 1.0, 1.0);
603 ((TProfile2D*)fPedestalHighGainRatio[i])->Divide(GetPedProfileHighGain(i), fReference->GetPedProfileHighGain(i), 1.0, 1.0);
604 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->Divide(GetPeakProfileLowGain(i), fReference->GetPeakProfileLowGain(i), 1.0, 1.0);
605 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->Divide(GetPeakProfileHighGain(i), fReference->GetPeakProfileHighGain(i), 1.0, 1.0);
607 //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
608 //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
609 for (int j = 0; j <= fColumns; j++) {
610 for (int k = 0; k <= fRows; k++) {
611 int bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
612 double diff = fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) - GetPeakProfileHighGain(i)->GetBinContent(bin);
613 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
614 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
616 diff = fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) - GetPeakProfileLowGain(i)->GetBinContent(bin);
617 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
618 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
620 diff = fReference->GetPedProfileHighGain(i)->GetBinContent(bin) - GetPedProfileHighGain(i)->GetBinContent(bin);
621 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(j+1, k+1, diff);
622 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
624 diff = fReference->GetPedProfileLowGain(i)->GetBinContent(bin) - GetPedProfileLowGain(i)->GetBinContent(bin);
625 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(j+1, k+1, diff);
626 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
635 //_____________________________________________________________________
636 void AliCaloCalibPedestal::ComputeDeadTowers(int threshold, const char * deadMapFile)
638 //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
644 char name[512];//Quite a long temp buffer, just in case the filename includes a path
647 snprintf(name, 512, "%s.txt", deadMapFile);
648 fout = new ofstream(name);
649 snprintf(name, 512, "%sdiff.txt", deadMapFile);
650 diff = new ofstream(name);
651 if (!fout->is_open()) {
653 fout = 0;//Set the pointer to empty if the file was not opened
655 if (!diff->is_open()) {
657 fout = 0;//Set the pointer to empty if the file was not opened
661 for (int i = 0; i < fModules; i++) {
662 if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
663 for (int j = 1; j <= fColumns; j++) {
664 for (int k = 1; k <= fRows; k++) {
666 if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {//It's dead
667 countTot++;//One more dead total
670 << (fRows - k) << " "
673 << "0" << endl;//Write the status to the deadmap file, if the file is open.
676 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= threshold) {
677 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
678 countNew++;//This tower wasn't dead before!
681 << (fRows - k) << " "
684 << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
688 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
691 else { //It's ALIVE!!
692 //Don't bother with writing the live ones.
694 // (*fout) << i << " "
695 // << (fRows - k) << " "
698 // << "1" << endl;//Write the status to the deadmap file, if the file is open.
699 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < threshold) {
700 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
701 countRes++; //This tower was dead before => it's a miracle! :P
704 << (fRows - k) << " "
707 << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
711 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
717 }//end if GetEntries >= 0
726 fDeadTowers = countTot;
727 fNewDeadTowers = countNew;
728 fResurrectedTowers = countRes;