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 //________________________________________________________________________
37 //#include "TCanvas.h"
47 #include "AliRawReader.h"
48 #include "AliCaloRawStreamV3.h"
51 #include "AliCaloCalibPedestal.h"
53 ClassImp(AliCaloCalibPedestal)
57 // ctor; initialize everything in order to avoid compiler warnings
58 AliCaloCalibPedestal::AliCaloCalibPedestal(kDetType detectorType) :
62 fPedestalLEDRefLowGain(),
63 fPedestalLEDRefHighGain(),
64 fPeakMinusPedLowGain(),
65 fPeakMinusPedHighGain(),
66 fPeakMinusPedHighGainHisto(),
67 fPedestalLowGainDiff(),
68 fPedestalHighGainDiff(),
69 fPedestalLEDRefLowGainDiff(),
70 fPedestalLEDRefHighGainDiff(),
71 fPeakMinusPedLowGainDiff(),
72 fPeakMinusPedHighGainDiff(),
73 fPedestalLowGainRatio(),
74 fPedestalHighGainRatio(),
75 fPedestalLEDRefLowGainRatio(),
76 fPedestalLEDRefHighGainRatio(),
77 fPeakMinusPedLowGainRatio(),
78 fPeakMinusPedHighGainRatio(),
84 fResurrectedTowers(0),
97 fSelectPedestalSamples(kTRUE),
98 fFirstPedestalSample(0),
99 fLastPedestalSample(15),
101 fWarningThreshold(50),
102 fWarningFraction(0.002),
105 //Default constructor. First we set the detector-type related constants.
106 if (detectorType == kPhos) {
107 fColumns = fgkPhosCols;
109 fLEDRefs = fgkPhosLEDRefs;
110 fModules = fgkPhosModules;
111 fCaloString = "PHOS";
117 //We'll just trust the enum to keep everything in line, so that if detectorType
118 //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
119 //case, if someone intentionally gives another number
120 fColumns = AliEMCALGeoParams::fgkEMCALCols;
121 fRows = AliEMCALGeoParams::fgkEMCALRows;
122 fLEDRefs = AliEMCALGeoParams::fgkEMCALLEDRefs;
123 fModules = AliEMCALGeoParams::fgkEMCALModules;
124 fCaloString = "EMCAL";
129 fDetType = detectorType;
131 // ValidateProfiles(); // not to be done in ctor; info from Axel N.
134 //_____________________________________________________________________
135 void AliCaloCalibPedestal::ValidateProfiles()
137 //Make sure the basic histos exist
138 if (!fPedestalLowGain.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
141 //Then, loop for the requested number of modules
143 for (int i = 0; i < fModules; i++) {
144 //Pedestals, low gain
145 name = "hPedlowgain";
147 title = "Pedestals, low gain, module ";
149 fPedestalLowGain.Add(new TProfile2D(name, title,
150 fColumns, 0.0, fColumns,
151 fRows, fRowMin, fRowMax,"s"));
153 //Pedestals, high gain
154 name = "hPedhighgain";
156 title = "Pedestals, high gain, module ";
158 fPedestalHighGain.Add(new TProfile2D(name, title,
159 fColumns, 0.0, fColumns,
160 fRows, fRowMin, fRowMax,"s"));
162 //LED Ref/Mon pedestals, low gain
163 name = "hPedestalLEDReflowgain";
165 title = "Pedestal LEDRef, low gain, module ";
167 fPedestalLEDRefLowGain.Add(new TProfile(name, title,
168 fLEDRefs, 0.0, fLEDRefs, "s"));
170 //LED Ref/Mon pedestals, high gain
171 name = "hPedestalLEDRefhighgain";
173 title = "Pedestal LEDRef, high gain, module ";
175 fPedestalLEDRefHighGain.Add(new TProfile(name, title,
176 fLEDRefs, 0.0, fLEDRefs, "s"));
178 //Peak-Pedestals, low gain
179 name = "hPeakMinusPedlowgain";
181 title = "Peak-Pedestal, low gain, module ";
183 fPeakMinusPedLowGain.Add(new TProfile2D(name, title,
184 fColumns, 0.0, fColumns,
185 fRows, fRowMin, fRowMax,"s"));
187 //Peak-Pedestals, high gain
188 name = "hPeakMinusPedhighgain";
190 title = "Peak-Pedestal, high gain, module ";
192 fPeakMinusPedHighGain.Add(new TProfile2D(name, title,
193 fColumns, 0.0, fColumns,
194 fRows, fRowMin, fRowMax,"s"));
196 //Peak-Pedestals, high gain - TH2F histo
197 name = "hPeakMinusPedhighgainHisto";
199 title = "Peak-Pedestal, high gain, module ";
201 fPeakMinusPedHighGainHisto.Add(new TH2F(name, title,
202 fColumns*fRows, 0.0, fColumns*fRows,
207 title = "Dead map, module ";
209 fDeadMap.Add(new TH2D(name, title, fColumns, 0.0, fColumns,
210 fRows, fRowMin, fRowMax));
212 }//end for nModules create the histograms
214 CompressAndSetOwner();
217 //_____________________________________________________________________
218 void AliCaloCalibPedestal::CompressAndSetOwner()
220 //Compress the arrays, in order to remove the empty objects (a 16 slot array is created by default)
221 fPedestalLowGain.Compress();
222 fPedestalHighGain.Compress();
223 fPedestalLEDRefLowGain.Compress();
224 fPedestalLEDRefHighGain.Compress();
225 fPeakMinusPedLowGain.Compress();
226 fPeakMinusPedHighGain.Compress();
227 fPeakMinusPedHighGainHisto.Compress();
230 // set owner ship for everyone
231 fPedestalLowGain.SetOwner(kTRUE);
232 fPedestalHighGain.SetOwner(kTRUE);
233 fPedestalLEDRefLowGain.SetOwner(kTRUE);
234 fPedestalLEDRefHighGain.SetOwner(kTRUE);
235 fPeakMinusPedLowGain.SetOwner(kTRUE);
236 fPeakMinusPedHighGain.SetOwner(kTRUE);
237 fPeakMinusPedHighGainHisto.SetOwner(kTRUE);
238 fPedestalLowGainDiff.SetOwner(kTRUE);
239 fPedestalHighGainDiff.SetOwner(kTRUE);
240 fPedestalLEDRefLowGainDiff.SetOwner(kTRUE);
241 fPedestalLEDRefHighGainDiff.SetOwner(kTRUE);
242 fPeakMinusPedLowGainDiff.SetOwner(kTRUE);
243 fPeakMinusPedHighGainDiff.SetOwner(kTRUE);
244 fPedestalLowGainRatio.SetOwner(kTRUE);
245 fPedestalHighGainRatio.SetOwner(kTRUE);
246 fPedestalLEDRefLowGainRatio.SetOwner(kTRUE);
247 fPedestalLEDRefHighGainRatio.SetOwner(kTRUE);
248 fPeakMinusPedLowGainRatio.SetOwner(kTRUE);
249 fPeakMinusPedHighGainRatio.SetOwner(kTRUE);
250 fDeadMap.SetOwner(kTRUE);
254 //_____________________________________________________________________
255 AliCaloCalibPedestal::~AliCaloCalibPedestal()
259 if (fReference) delete fReference;//Delete the reference object, if it has been loaded
261 // delete also TObjArray's
262 fPedestalLowGain.Delete();
263 fPedestalHighGain.Delete();
264 fPedestalLEDRefLowGain.Delete();
265 fPedestalLEDRefHighGain.Delete();
266 fPeakMinusPedLowGain.Delete();
267 fPeakMinusPedHighGain.Delete();
268 fPeakMinusPedHighGainHisto.Delete();
269 fPedestalLowGainDiff.Delete();
270 fPedestalHighGainDiff.Delete();
271 fPedestalLEDRefLowGainDiff.Delete();
272 fPedestalLEDRefHighGainDiff.Delete();
273 fPeakMinusPedLowGainDiff.Delete();
274 fPeakMinusPedHighGainDiff.Delete();
275 fPedestalLowGainRatio.Delete();
276 fPedestalHighGainRatio.Delete();
277 fPedestalLEDRefLowGainRatio.Delete();
278 fPedestalLEDRefHighGainRatio.Delete();
279 fPeakMinusPedLowGainRatio.Delete();
280 fPeakMinusPedHighGainRatio.Delete();
286 //_____________________________________________________________________
287 AliCaloCalibPedestal::AliCaloCalibPedestal(AliCaloCalibPedestal &ped) :
291 fPedestalLEDRefLowGain(),
292 fPedestalLEDRefHighGain(),
293 fPeakMinusPedLowGain(),
294 fPeakMinusPedHighGain(),
295 fPeakMinusPedHighGainHisto(),
296 fPedestalLowGainDiff(),
297 fPedestalHighGainDiff(),
298 fPedestalLEDRefLowGainDiff(),
299 fPedestalLEDRefHighGainDiff(),
300 fPeakMinusPedLowGainDiff(),
301 fPeakMinusPedHighGainDiff(),
302 fPedestalLowGainRatio(),
303 fPedestalHighGainRatio(),
304 fPedestalLEDRefLowGainRatio(),
305 fPedestalLEDRefHighGainRatio(),
306 fPeakMinusPedLowGainRatio(),
307 fPeakMinusPedHighGainRatio(),
309 fNEvents(ped.GetNEvents()),
310 fNChanFills(ped.GetNChanFills()),
311 fDeadTowers(ped.GetDeadTowerCount()),
312 fNewDeadTowers(ped.GetDeadTowerNew()),
313 fResurrectedTowers(ped.GetDeadTowerResurrected()),
314 fReference( 0 ), //! note that we do not try to copy the reference info here
315 fDetType(ped.GetDetectorType()),
316 fColumns(ped.GetColumns()),
317 fRows(ped.GetRows()),
318 fLEDRefs(ped.GetLEDRefs()),
319 fModules(ped.GetModules()),
320 fRowMin(ped.GetRowMin()),
321 fRowMax(ped.GetRowMax()),
322 fRowMultiplier(ped.GetRowMultiplier()),
323 fCaloString(ped.GetCaloString()),
324 fMapping(NULL), //! note that we are not copying the map info
325 fRunNumber(ped.GetRunNumber()),
326 fSelectPedestalSamples(ped.GetSelectPedestalSamples()),
327 fFirstPedestalSample(ped.GetFirstPedestalSample()),
328 fLastPedestalSample(ped.GetLastPedestalSample()),
329 fDeadThreshold(ped.GetDeadThreshold()),
330 fWarningThreshold(ped.GetWarningThreshold()),
331 fWarningFraction(ped.GetWarningFraction()),
332 fHotSigma(ped.GetHotSigma())
334 // Then the ObjArray ones; we add the histograms rather than trying TObjArray = assignment
335 //DS: this has not really been tested yet..
336 for (int i = 0; i < fModules; i++) {
337 fPedestalLowGain.Add( ped.GetPedProfileLowGain(i) );
338 fPedestalHighGain.Add( ped.GetPedProfileHighGain(i) );
339 fPedestalLEDRefLowGain.Add( ped.GetPedLEDRefProfileLowGain(i) );
340 fPedestalLEDRefHighGain.Add( ped.GetPedLEDRefProfileHighGain(i) );
341 fPeakMinusPedLowGain.Add( ped.GetPeakProfileLowGain(i) );
342 fPeakMinusPedHighGain.Add( ped.GetPeakProfileHighGain(i) );
343 fPeakMinusPedHighGainHisto.Add( ped.GetPeakHighGainHisto(i) );
345 fDeadMap.Add( ped.GetDeadMap(i) );
348 CompressAndSetOwner();
351 // assignment operator; use copy ctor to make life easy..
352 //_____________________________________________________________________
353 AliCaloCalibPedestal& AliCaloCalibPedestal::operator = (AliCaloCalibPedestal &source)
355 // assignment operator; use copy ctor
356 if (&source == this) return *this;
358 new (this) AliCaloCalibPedestal(source);
362 //_____________________________________________________________________
363 void AliCaloCalibPedestal::Reset()
365 ValidateProfiles(); // make sure histos/profiles exist
366 // Reset all arrays/histograms
367 for (int i = 0; i < fModules; i++) {
368 GetPedProfileLowGain(i)->Reset();
369 GetPedProfileHighGain(i)->Reset();
370 GetPedLEDRefProfileLowGain(i)->Reset();
371 GetPedLEDRefProfileHighGain(i)->Reset();
372 GetPeakProfileLowGain(i)->Reset();
373 GetPeakProfileHighGain(i)->Reset();
374 GetPeakHighGainHisto(i)->Reset();
375 GetDeadMap(i)->Reset();
377 if (!fPedestalLowGainDiff.IsEmpty()) {
378 //This means that the comparison profiles have been created.
380 GetPedProfileLowGainDiff(i)->Reset();
381 GetPedProfileHighGainDiff(i)->Reset();
382 GetPedLEDRefProfileLowGainDiff(i)->Reset();
383 GetPedLEDRefProfileHighGainDiff(i)->Reset();
384 GetPeakProfileLowGainDiff(i)->Reset();
385 GetPeakProfileHighGainDiff(i)->Reset();
387 GetPedProfileLowGainRatio(i)->Reset();
388 GetPedProfileHighGainRatio(i)->Reset();
389 GetPedLEDRefProfileLowGainRatio(i)->Reset();
390 GetPedLEDRefProfileHighGainRatio(i)->Reset();
391 GetPeakProfileLowGainRatio(i)->Reset();
392 GetPeakProfileHighGainRatio(i)->Reset();
399 fResurrectedTowers = 0;
401 //To think about: should fReference be deleted too?... let's not do it this time, at least...
404 // Parameter/cut handling
405 //_____________________________________________________________________
406 void AliCaloCalibPedestal::SetParametersFromFile(const char *parameterFile)
408 // Note: this method is a bit more complicated than it really has to be
409 // - allowing for multiple entries per line, arbitrary order of the
410 // different variables etc. But I wanted to try and do this in as
411 // correct a C++ way as I could (as an exercise).
413 static const string delimitor("::");
415 // open, check input file
416 ifstream in( parameterFile );
418 printf("in AliCaloCalibPedestal::SetParametersFromFile - Using default/run_time parameters.\n");
425 while ((in.rdstate() & ios::failbit) == 0 ) {
427 // Read into the raw char array and then construct a string
428 // to do the searching
429 in.getline(readline, 1024);
430 istringstream s(readline);
432 while ( ( s.rdstate() & ios::failbit ) == 0 ) {
437 // check stream status
438 if( ( s.rdstate() & ios::failbit ) == ios::failbit) break;
440 // skip rest of line if comments found
441 if( keyValue.substr( 0, 2 ) == "//" ) break;
443 // look for "::" in keyValue pair
444 size_t position = keyValue.find( delimitor );
445 if( position == string::npos ) {
446 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
449 // split keyValue pair
450 string key( keyValue.substr( 0, position ) );
451 string value( keyValue.substr( position+delimitor.size(),
452 keyValue.size()-delimitor.size() ) );
454 // check value does not contain a new delimitor
455 if( value.find( delimitor ) != string::npos ) {
456 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
459 // debug: check key value pair
460 // printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
462 // if the key matches with something we expect, we assign the new value
463 istringstream iss(value);
464 // the comparison strings defined at the beginning of this method
465 if ( (key == "fFirstPedestalSample") || (key == "fLastPedestalSample") || (key == "fDeadThreshold") || (key == "fWarningThreshold") || (key == "fWarningFraction") || (key == "fHotSigma") ) {
466 printf("AliCaloCalibPedestal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
468 if (key == "fFirstPedestalSample") {
469 iss >> fFirstPedestalSample;
471 else if (key == "fLastPedestalSample") {
472 iss >> fLastPedestalSample;
474 else if (key == "fDeadThreshold") {
475 iss >> fDeadThreshold;
477 else if (key == "fWarningThreshold") {
478 iss >> fWarningThreshold;
480 else if (key == "fWarningFraction") {
481 iss >> fWarningFraction;
483 else if (key == "fHotSigma") {
497 //_____________________________________________________________________
498 void AliCaloCalibPedestal::WriteParametersToFile(const char *parameterFile)
500 //Write parameters in file.
502 static const string delimitor("::");
503 ofstream out( parameterFile );
504 out << "// " << parameterFile << endl;
505 out << "fFirstPedestalSample" << "::" << fFirstPedestalSample << endl;
506 out << "fLastPedestalSample" << "::" << fLastPedestalSample << endl;
507 out << "fDeadThreshold" << "::" << fDeadThreshold << endl;
508 out << "fWarningThreshold" << "::" << fWarningThreshold << endl;
509 out << "fWarningFraction" << "::" << fWarningFraction << endl;
510 out << "fHotSigma" << "::" << fHotSigma << endl;
516 //_____________________________________________________________________
517 Bool_t AliCaloCalibPedestal::AddInfo(AliCaloCalibPedestal *ped)
519 ValidateProfiles(); // make sure histos/profiles exist
520 // just do this for the basic histograms/profiles that get filled in ProcessEvent
521 // may not have data for all modules, but let's just Add everything..
522 for (int i = 0; i < fModules; i++) {
523 GetPedProfileLowGain(i)->Add( ped->GetPedProfileLowGain(i) );
524 GetPedProfileHighGain(i)->Add( ped->GetPedProfileHighGain(i) );
525 GetPedLEDRefProfileLowGain(i)->Add( ped->GetPedLEDRefProfileLowGain(i) );
526 GetPedLEDRefProfileHighGain(i)->Add( ped->GetPedLEDRefProfileHighGain(i) );
527 GetPeakProfileLowGain(i)->Add( ped->GetPeakProfileLowGain(i) );
528 GetPeakProfileHighGain(i)->Add( ped->GetPeakProfileHighGain(i) );
529 GetPeakHighGainHisto(i)->Add( ped->GetPeakHighGainHisto(i) );
533 // We should also copy other pieces of info: counters and parameters
534 // (not number of columns and rows etc which should be the same)
535 // note that I just assign them here rather than Add them, but we
536 // normally just Add (e.g. in Preprocessor) one object so this should be fine.
537 fNEvents = ped->GetNEvents();
538 fNChanFills = ped->GetNChanFills();
539 fDeadTowers = ped->GetDeadTowerCount();
540 fNewDeadTowers = ped->GetDeadTowerNew();
541 fResurrectedTowers = ped->GetDeadTowerResurrected();
542 fRunNumber = ped->GetRunNumber();
543 fSelectPedestalSamples = ped->GetSelectPedestalSamples();
544 fFirstPedestalSample = ped->GetFirstPedestalSample();
545 fLastPedestalSample = ped->GetLastPedestalSample();
546 fDeadThreshold = ped->GetDeadThreshold();
547 fWarningThreshold = ped->GetWarningThreshold();
548 fWarningFraction = ped->GetWarningFraction();
549 fHotSigma = ped->GetHotSigma();
551 // DeadMap; Diff profiles etc would need to be redone after this operation
553 return kTRUE;//We succesfully added info from the supplied object
556 //_____________________________________________________________________
557 Bool_t AliCaloCalibPedestal::ProcessEvent(AliRawReader *rawReader)
559 // if fMapping is NULL the rawstream will crate its own mapping
560 AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
561 if (fDetType == kEmCal) {
562 rawReader->Select("EMCAL", 0, AliEMCALGeoParams::fgkLastAltroDDL) ; //select EMCAL DDL range
564 return ProcessEvent(&rawStream);
567 //_____________________________________________________________________
568 Bool_t AliCaloCalibPedestal::ProcessEvent(AliCaloRawStreamV3 *in)
570 // Method to process=analyze one event in the data stream
571 if (!in) return kFALSE; //Return right away if there's a null pointer
572 in->Reset(); // just in case the next customer forgets to check if the stream was reset..
574 fNEvents++; // one more event
576 if (fNEvents==1) ValidateProfiles(); // 1st event, make sure histos/profiles exist
578 // indices for the reading
581 int i = 0; // sample counter
584 // start loop over input stream
585 while (in->NextDDL()) {
586 while (in->NextChannel()) {
589 int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
594 vector<int> pedSamples;
596 while (in->NextBunch()) {
597 const UShort_t *sig = in->GetSignals();
598 startBin = in->GetStartTimeBin();
599 nsamples += in->GetBunchLength();
600 for (i = 0; i < in->GetBunchLength(); i++) {
604 // check if it's a min or max value
605 if (sample < min) min = sample;
606 if (sample > max) max = sample;
608 // should we add it for the pedestal calculation?
609 if ( (fFirstPedestalSample<=time && time<=fLastPedestalSample) || // sample time in range
610 !fSelectPedestalSamples ) { // or we don't restrict the sample range.. - then we'll take all
611 pedSamples.push_back( sig[i] );
615 } // loop over samples in bunch
616 } // loop over bunches
618 if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
620 // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
621 int arrayPos = in->GetModule(); //The modules are numbered starting from 0
622 if (arrayPos >= fModules) {
623 //TODO: return an error message, if appopriate (perhaps if debug>0?)
627 if (arrayPos < 0 || arrayPos >= fModules) {
628 printf("Oh no: arrayPos = %i.\n", arrayPos);
631 fNChanFills++; // one more channel found, and profile to be filled
632 //NOTE: coordinates are (column, row) for the profiles
633 if ( in->IsLowGain() ) {
634 //fill the low gain histograms
635 ((TProfile2D*)fPeakMinusPedLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
636 if (nPed>0) { // only fill pedestal info in case it could be calculated
637 for ( i=0; i<nPed; i++) {
638 ((TProfile2D*)fPedestalLowGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), pedSamples[i]);
642 else if ( in->IsHighGain() ) {
643 //fill the high gain ones
644 ((TProfile2D*)fPeakMinusPedHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), max - min);
645 if (nPed>0) { // only fill pedestal info in case it could be calculated
646 for ( i=0; i<nPed; i++) {
647 ((TProfile2D*)fPedestalHighGain[arrayPos])->Fill(in->GetColumn(), fRowMultiplier*in->GetRow(), pedSamples[i]);
650 // for warning checks
651 int idx = in->GetRow() + fRows * in->GetColumn();
652 ((TH2F*)fPeakMinusPedHighGainHisto[arrayPos])->Fill(idx, max - min);
654 else if ( in->IsLEDMonData() ) {
655 // for LED Mon data, the mapping class holds the gain info in the Row variable
656 // and the Strip number in the Column..
657 int gain = in->GetRow();
658 int stripId = in->GetColumn();
659 if (nPed>0 && stripId<fLEDRefs) {
661 for ( i=0; i<nPed; i++) {
662 ((TProfile*)fPedestalLEDRefLowGain[arrayPos])->Fill(stripId, pedSamples[i]);
666 for ( i=0; i<nPed; i++) {
667 ((TProfile*)fPedestalLEDRefHighGain[arrayPos])->Fill(stripId, pedSamples[i]);
673 } // nsamples>0 check, some data found for this channel; not only trailer/header
674 }// end while over channel
675 }//end while over DDL's, of input stream
681 //_____________________________________________________________________
682 Bool_t AliCaloCalibPedestal::SaveHistograms(TString fileName, Bool_t saveEmptyHistos)
684 //Saves all the histograms (or profiles, to be accurate) to the designated file
685 ValidateProfiles(); // make sure histos/profiles exist
686 TFile destFile(fileName, "recreate");
688 if (destFile.IsZombie()) {
694 for (int i = 0; i < fModules; i++) {
695 if( ((TProfile2D *)fPeakMinusPedLowGain[i])->GetEntries() || saveEmptyHistos) {
696 fPeakMinusPedLowGain[i]->Write();
698 if( ((TProfile2D *)fPeakMinusPedHighGain[i])->GetEntries() || saveEmptyHistos) {
699 fPeakMinusPedHighGain[i]->Write();
701 if( ((TProfile2D *)fPedestalLowGain[i])->GetEntries() || saveEmptyHistos) {
702 fPedestalLowGain[i]->Write();
704 if( ((TProfile2D *)fPedestalHighGain[i])->GetEntries() || saveEmptyHistos) {
705 fPedestalHighGain[i]->Write();
707 if( ((TProfile *)fPedestalLEDRefLowGain[i])->GetEntries() || saveEmptyHistos) {
708 fPedestalLEDRefLowGain[i]->Write();
710 if( ((TProfile *)fPedestalLEDRefHighGain[i])->GetEntries() || saveEmptyHistos) {
711 fPedestalLEDRefHighGain[i]->Write();
713 if( ((TH2F *)fPeakMinusPedHighGainHisto[i])->GetEntries() || saveEmptyHistos) {
714 fPeakMinusPedHighGainHisto[i]->Write();
724 //_____________________________________________________________________
725 Bool_t AliCaloCalibPedestal::LoadReferenceCalib(TString fileName, TString objectName)
728 //Make sure that the histograms created when loading the object are not destroyed as the file object is destroyed
729 TH1::AddDirectory(kFALSE);
731 TFile *sourceFile = new TFile(fileName);
732 if (sourceFile->IsZombie()) {
733 return kFALSE;//We couldn't load the reference
736 if (fReference) delete fReference;//Delete the reference object, if it already exists
739 fReference = (AliCaloCalibPedestal*)sourceFile->Get(objectName);
741 if (!fReference || !(fReference->InheritsFrom(AliCaloCalibPedestal::Class())) || (fReference->GetDetectorType() != fDetType)) {
742 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
749 //Reset the histogram ownership behaviour. NOTE: a better workaround would be good, since this may accidentally set AddDirectory to true, even
750 //if we are called by someone who has set it to false...
751 TH1::AddDirectory(kTRUE);
753 return kTRUE;//We succesfully loaded the object
757 //_____________________________________________________________________
758 Bool_t AliCaloCalibPedestal::SetReference(AliCaloCalibPedestal *ref)
760 if (fReference) delete fReference;//Delete the reference object, if it already exists
765 if (!fReference || (fReference->GetDetectorType() != fDetType)) {
766 if (fReference) delete fReference;//Delete the object, in case we had an object of the wrong type
771 return kTRUE;//We succesfully loaded the object
774 //_____________________________________________________________________
775 void AliCaloCalibPedestal::ValidateComparisonProfiles()
777 //Make sure the comparison histos exist
778 if (!fPedestalLowGainDiff.IsEmpty()) return; //The profiles already exist. We just check one, because they're all created at
782 //Then, loop for the requested number of modules
784 for (int i = 0; i < fModules; i++) {
785 //Pedestals, low gain
786 name = "hPedlowgainDiff";
788 title = "Pedestals difference, low gain, module ";
790 fPedestalLowGainDiff.Add(new TProfile2D(name, title,
791 fColumns, 0.0, fColumns,
792 fRows, fRowMin, fRowMax,"s"));
794 //Pedestals, high gain
795 name = "hPedhighgainDiff";
797 title = "Pedestals difference, high gain, module ";
799 fPedestalHighGainDiff.Add(new TProfile2D(name, title,
800 fColumns, 0.0, fColumns,
801 fRows, fRowMin, fRowMax,"s"));
803 //LED Ref/Mon pedestals, low gain
804 name = "hPedestalLEDReflowgainDiff";
806 title = "Pedestal difference LEDRef, low gain, module ";
808 fPedestalLEDRefLowGainDiff.Add(new TProfile(name, title,
809 fLEDRefs, 0.0, fLEDRefs, "s"));
811 //LED Ref/Mon pedestals, high gain
812 name = "hPedestalLEDRefhighgainDiff";
814 title = "Pedestal difference LEDRef, high gain, module ";
816 fPedestalLEDRefHighGainDiff.Add(new TProfile(name, title,
817 fLEDRefs, 0.0, fLEDRefs, "s"));
819 //Peak-Pedestals, high gain
820 name = "hPeakMinusPedhighgainDiff";
822 title = "Peak-Pedestal difference, high gain, module ";
824 fPeakMinusPedHighGainDiff.Add(new TProfile2D(name, title,
825 fColumns, 0.0, fColumns,
826 fRows, fRowMin, fRowMax,"s"));
828 //Peak-Pedestals, low gain
829 name = "hPeakMinusPedlowgainDiff";
831 title = "Peak-Pedestal difference, low gain, module ";
833 fPeakMinusPedLowGainDiff.Add(new TProfile2D(name, title,
834 fColumns, 0.0, fColumns,
835 fRows, fRowMin, fRowMax,"s"));
837 //Pedestals, low gain
838 name = "hPedlowgainRatio";
840 title = "Pedestals ratio, low gain, module ";
842 fPedestalLowGainRatio.Add(new TProfile2D(name, title,
843 fColumns, 0.0, fColumns,
844 fRows, fRowMin, fRowMax,"s"));
846 //Pedestals, high gain
847 name = "hPedhighgainRatio";
849 title = "Pedestals ratio, high gain, module ";
851 fPedestalHighGainRatio.Add(new TProfile2D(name, title,
852 fColumns, 0.0, fColumns,
853 fRows, fRowMin, fRowMax,"s"));
855 //LED Ref/Mon pedestals, low gain
856 name = "hPedestalLEDReflowgainRatio";
858 title = "Pedestal ratio LEDRef, low gain, module ";
860 fPedestalLEDRefLowGainRatio.Add(new TProfile(name, title,
861 fLEDRefs, 0.0, fLEDRefs, "s"));
863 //LED Ref/Mon pedestals, high gain
864 name = "hPedestalLEDRefhighgainRatio";
866 title = "Pedestal ratio LEDRef, high gain, module ";
868 fPedestalLEDRefHighGainRatio.Add(new TProfile(name, title,
869 fLEDRefs, 0.0, fLEDRefs, "s"));
871 //Peak-Pedestals, low gain
872 name = "hPeakMinusPedlowgainRatio";
874 title = "Peak-Pedestal ratio, low gain, module ";
876 fPeakMinusPedLowGainRatio.Add(new TProfile2D(name, title,
877 fColumns, 0.0, fColumns,
878 fRows, fRowMin, fRowMax,"s"));
880 //Peak-Pedestals, high gain
881 name = "hPeakMinusPedhighgainRatio";
883 title = "Peak-Pedestal ratio, high gain, module ";
885 fPeakMinusPedHighGainRatio.Add(new TProfile2D(name, title,
886 fColumns, 0.0, fColumns,
887 fRows, fRowMin, fRowMax,"s"));
889 }//end for nModules create the histograms
892 //_____________________________________________________________________
893 void AliCaloCalibPedestal::ComputeDiffAndRatio()
895 // calculate differences and ratios relative to a reference
896 ValidateProfiles(); // make sure histos/profiles exist
897 ValidateComparisonProfiles();//Make sure the comparison histos exist
900 return;//Return if the reference object isn't loaded
906 for (int i = 0; i < fModules; i++) {
907 //For computing the difference, we cannot simply do TProfile2D->Add(), because that subtracts the sum of all entries,
908 //which means that the mean of the new profile will not be the difference of the means. So do it by hand:
909 for (int j = 0; j < fColumns; j++) {
910 for (int k = 0; k < fRows; k++) {
911 bin = ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->GetBin(j+1, k+1);//Note that we assume here that all histos have the same structure...
913 if (fReference->GetPeakProfileHighGain(i)->GetBinContent(bin) > 0) {
914 diff = GetPeakProfileHighGain(i)->GetBinContent(bin) - fReference->GetPeakProfileHighGain(i)->GetBinContent(bin);
915 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinContent(bin, diff);
916 ((TProfile2D*)fPeakMinusPedHighGainDiff[i])->SetBinEntries(bin, 1);
917 ratio = GetPeakProfileHighGain(i)->GetBinContent(bin) / fReference->GetPeakProfileHighGain(i)->GetBinContent(bin);
918 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->SetBinContent(bin, ratio);
919 ((TProfile2D*)fPeakMinusPedHighGainRatio[i])->SetBinEntries(bin, 1);
922 if (fReference->GetPeakProfileLowGain(i)->GetBinContent(bin) > 0) {
923 diff = GetPeakProfileLowGain(i)->GetBinContent(bin) - fReference->GetPeakProfileLowGain(i)->GetBinContent(bin);
924 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinContent(bin, diff);
925 ((TProfile2D*)fPeakMinusPedLowGainDiff[i])->SetBinEntries(bin, 1);
926 ratio = GetPeakProfileLowGain(i)->GetBinContent(bin) / fReference->GetPeakProfileLowGain(i)->GetBinContent(bin);
927 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->SetBinContent(bin, ratio);
928 ((TProfile2D*)fPeakMinusPedLowGainRatio[i])->SetBinEntries(bin, 1);
931 if (fReference->GetPedProfileHighGain(i)->GetBinContent(bin) > 0) {
932 diff = GetPedProfileHighGain(i)->GetBinContent(bin) - fReference->GetPedProfileHighGain(i)->GetBinContent(bin);
933 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinContent(bin, diff);
934 ((TProfile2D*)fPedestalHighGainDiff[i])->SetBinEntries(bin, 1);
935 ratio = GetPedProfileHighGain(i)->GetBinContent(bin) / fReference->GetPedProfileHighGain(i)->GetBinContent(bin);
936 ((TProfile2D*)fPedestalHighGainRatio[i])->SetBinContent(bin, ratio);
937 ((TProfile2D*)fPedestalHighGainRatio[i])->SetBinEntries(bin, 1);
940 if (fReference->GetPedProfileLowGain(i)->GetBinContent(bin) > 0) {
941 diff = GetPedProfileLowGain(i)->GetBinContent(bin) - fReference->GetPedProfileLowGain(i)->GetBinContent(bin);
942 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinContent(bin, diff);
943 ((TProfile2D*)fPedestalLowGainDiff[i])->SetBinEntries(bin, 1);
944 ratio = GetPedProfileLowGain(i)->GetBinContent(bin) / fReference->GetPedProfileLowGain(i)->GetBinContent(bin);
945 ((TProfile2D*)fPedestalLowGainRatio[i])->SetBinContent(bin, ratio);
946 ((TProfile2D*)fPedestalLowGainRatio[i])->SetBinEntries(bin, 1);
952 // same for LED Ref/Mon channels
953 for (int j = 0; j <= fLEDRefs; j++) {
954 bin = j+1;//Note that we assume here that all histos have the same structure...
956 if (fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) > 0) {
957 diff = GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) - fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin);
958 ((TProfile*)fPedestalLEDRefHighGainDiff[i])->SetBinContent(bin, diff);
959 ((TProfile*)fPedestalLEDRefHighGainDiff[i])->SetBinEntries(bin, 1);
960 ratio = GetPedLEDRefProfileHighGain(i)->GetBinContent(bin) / fReference->GetPedLEDRefProfileHighGain(i)->GetBinContent(bin);
961 ((TProfile*)fPedestalLEDRefHighGainRatio[i])->SetBinContent(bin, ratio);
962 ((TProfile*)fPedestalLEDRefHighGainRatio[i])->SetBinEntries(bin, 1);
965 if (fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) > 0) {
966 diff = GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) - fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin);
967 ((TProfile*)fPedestalLEDRefLowGainDiff[i])->SetBinContent(bin, diff);
968 ((TProfile*)fPedestalLEDRefLowGainDiff[i])->SetBinEntries(bin, 1);
969 ratio = GetPedLEDRefProfileLowGain(i)->GetBinContent(bin) / fReference->GetPedLEDRefProfileLowGain(i)->GetBinContent(bin);
970 ((TProfile*)fPedestalLEDRefLowGainRatio[i])->SetBinContent(bin, ratio);
971 ((TProfile*)fPedestalLEDRefLowGainRatio[i])->SetBinEntries(bin, 1);
980 //_____________________________________________________________________
981 void AliCaloCalibPedestal::ComputeHotAndWarningTowers(const char * hotMapFile)
982 { // look for hot/noisy towers
983 ValidateProfiles(); // make sure histos/profiles exist
985 char name[512];//Quite a long temp buffer, just in case the filename includes a path
988 snprintf(name, 512, "%s.txt", hotMapFile);
989 fout = new ofstream(name);
990 if (!fout->is_open()) {
992 fout = 0;//Set the pointer to empty if the file was not opened
996 for(int i = 0; i < fModules; i++){
998 //first we compute the peak-pedestal distribution for each supermodule...
999 if( GetPeakHighGainHisto(i)->GetEntries() > 0 ) {
1000 double min = GetPeakProfileHighGain(i)->GetBinContent(GetPeakProfileHighGain(i)->GetMinimumBin());
1001 double max = GetPeakProfileHighGain(i)->GetBinContent(GetPeakProfileHighGain(i)->GetMaximumBin());
1002 TH1D *hPeakFit = new TH1D(Form("hFit_%d", i), Form("hFit_%d", i), (int)((max-min)*10), min-1, max+1);
1004 for (int j = 1; j <= fColumns; j++) {
1005 for (int k = 1; k <= fRows; k++) {
1006 hPeakFit->Fill(GetPeakProfileHighGain(i)->GetBinContent(j, k));
1010 //...and then we try to fit it
1011 double mean = hPeakFit->GetMean();
1012 double sigma = hPeakFit->GetRMS();
1014 hPeakFit->Fit("gaus", "OQ", "", mean - 3*sigma, mean + 3*sigma);
1015 mean = hPeakFit->GetFunction("gaus")->GetParameter(1);
1016 sigma = hPeakFit->GetFunction("gaus")->GetParameter(2);
1018 catch (const std::exception & e) {
1019 printf("AliCaloCalibPedestal: TH1D PeakFit exception %s", e.what());
1025 //Then we look for warm/hot towers
1026 TH2F * hPeak2D = GetPeakHighGainHisto(i);
1027 hPeak2D->GetYaxis()->SetRangeUser( fWarningThreshold, hPeak2D->GetYaxis()->GetBinUpEdge(hPeak2D->GetNbinsY()) );
1030 int warnCounter = 0;
1031 for (int j = 1; j <= fColumns; j++) {
1032 for (int k = 1; k <= fRows; k++) {
1033 //we start looking for warm/warning towers...
1034 // histogram x-axis index
1035 idx = k-1 + fRows*(j-1); // this is what is used in the Fill call
1036 hPeak2D->GetXaxis()->SetRangeUser(idx, idx);
1037 warnCounter = (int) hPeak2D->Integral();
1038 if(warnCounter > fNEvents * fWarningFraction) {
1039 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kWarning);
1040 /* printf("mod %d col %d row %d warnCounter %d - status %d\n",
1041 i, j-1, k-1, warnCounter, (int) (kWarning)); */
1043 //...then we look for hot ones (towers whose values are greater than mean + X*sigma)
1044 if(GetPeakProfileHighGain(i)->GetBinContent(j, k) > mean + fHotSigma*sigma ) {
1045 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kHot);
1046 /* printf("mod %d col %d row %d binc %d - status %d\n",
1047 i, j-1, k-1, (int)(GetPeakProfileHighGain(i)->GetBinContent(j, k)), (int) (kHot)); */
1050 //Write the status to the hot/warm map file, if the file is open.
1051 // module - column - row - status (1=dead, 2= warm/warning , 3 = hot, see .h file enum)
1052 if (fout && ((TH2D*)fDeadMap[i])->GetBinContent(j, k) > 1) {
1057 << ((TH2D*)fDeadMap[i])->GetBinContent(j, k) << endl; }
1067 //_____________________________________________________________________
1068 void AliCaloCalibPedestal::ComputeDeadTowers(const char * deadMapFile)
1070 ValidateProfiles(); // make sure histos/profiles exist
1071 //Computes the number of dead towers etc etc into memory, after this you can call the GetDead... -functions
1075 ofstream * fout = 0;
1076 ofstream * diff = 0;
1077 char name[512];//Quite a long temp buffer, just in case the filename includes a path
1080 snprintf(name, 512, "%s.txt", deadMapFile);
1081 fout = new ofstream(name);
1082 snprintf(name, 512, "%sdiff.txt", deadMapFile);
1083 diff = new ofstream(name);
1084 if (!fout->is_open()) {
1086 fout = 0;//Set the pointer to empty if the file was not opened
1088 if (!diff->is_open()) {
1090 diff = 0;//Set the pointer to empty if the file was not opened
1094 for (int i = 0; i < fModules; i++) {
1095 if (GetPeakProfileHighGain(i)->GetEntries() > 0) { //don't care about empty histos
1096 for (int j = 1; j <= fColumns; j++) {
1097 for (int k = 1; k <= fRows; k++) {
1099 if (GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {//It's dead
1100 countTot++;//One more dead total
1106 << "0" << endl;//Write the status to the deadmap file, if the file is open.
1109 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) >= fDeadThreshold) {
1110 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kRecentlyDeceased);
1111 countNew++;//This tower wasn't dead before!
1113 ( *diff) << i << " "
1117 << "0" << endl;//Write the status to the deadmap difference file, if the file is open.
1121 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kDead);//This has been dead before. Nothing new
1124 else { //It's ALIVE!!
1125 //Don't bother with writing the live ones.
1126 if (fReference && fReference->GetPeakProfileHighGain(i)->GetBinContent(j, k) < fDeadThreshold) {
1127 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kResurrected);
1128 countRes++; //This tower was dead before => it's a miracle! :P
1134 << "1" << endl;//Write the status to the deadmap difference file, if the file is open.
1138 ((TH2D*)fDeadMap[i])->SetBinContent(j, k, kAlive);
1143 }//end for j/columns
1144 }//end if GetEntries >= 0
1153 fDeadTowers = countTot;
1154 fNewDeadTowers = countNew;
1155 fResurrectedTowers = countRes;
1158 //_____________________________________________________________________
1159 Bool_t AliCaloCalibPedestal::IsBadChannel(int imod, int icol, int irow) const
1161 // Check if status info histo exists
1162 if (!fDeadMap[imod]) {
1166 //Check if channel is dead or hot.
1167 Int_t status = (Int_t) ( ((TH2D*)fDeadMap[imod])->GetBinContent(icol,irow) );
1168 if(status == kAlive)
1175 //_____________________________________________________________________
1176 void AliCaloCalibPedestal::SetChannelStatus(int imod, int icol, int irow, int status)
1178 ValidateProfiles(); // make sure histos/profiles exist
1179 //Set status of channel dead, hot, alive ...
1180 ((TH2D*)fDeadMap[imod])->SetBinContent(icol, irow, status);