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 **************************************************************************/
15 /* $Id: AliCaloCalibSignal.cxx $ */
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 fSignals = new AliCaloCalibSignal( fDetType );
24 // AliCaloCalibSignal knows how many modules we have for PHOS or EMCAL
25 fNumModules = fSignals->GetModules();
28 // fSignals->ProcessEvent(fCaloRawStream,fRawEventHeaderBase);
30 // fSignals->GetXXX..()
32 //________________________________________________________________________
40 #include "AliRawReader.h"
41 #include "AliCaloRawStreamV3.h"
44 #include "AliCaloCalibSignal.h"
46 ClassImp(AliCaloCalibSignal)
50 // variables for TTree filling; not sure if they should be static or not
51 static int fChannelNum = 0; // for regular towers
52 static int fRefNum = 0; // for LED
53 static double fAmp = 0;
54 static double fAvgAmp = 0;
55 static double fRMS = 0;
57 // ctor; initialize everything in order to avoid compiler warnings
58 // put some reasonable defaults
59 AliCaloCalibSignal::AliCaloCalibSignal(kDetType detectorType) :
70 fAmpCut(40), // min. 40 ADC counts as default
71 fReqFractionAboveAmpCutVal(0.6), // 60% in a strip, per default
72 fReqFractionAboveAmp(kTRUE),
73 fAmpCutLEDRef(100), // min. 100 ADC counts as default
74 fReqLEDRefAboveAmpCutVal(kTRUE),
82 fTreeAvgAmpVsTime(NULL),
83 fTreeLEDAmpVsTime(NULL),
84 fTreeLEDAvgAmpVsTime(NULL)
86 //Default constructor. First we set the detector-type related constants.
87 if (detectorType == kPhos) {
88 fColumns = fgkPhosCols;
90 fLEDRefs = fgkPhosLEDRefs;
91 fModules = fgkPhosModules;
95 //We'll just trust the enum to keep everything in line, so that if detectorType
96 //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
97 //case, if someone intentionally gives another number
98 fColumns = AliEMCALGeoParams::fgkEMCALCols;
99 fRows = AliEMCALGeoParams::fgkEMCALRows;
100 fLEDRefs = AliEMCALGeoParams::fgkEMCALLEDRefs;
101 fModules = AliEMCALGeoParams::fgkEMCALModules;
102 fCaloString = "EMCAL";
105 fDetType = detectorType;
107 ResetInfo(); // trees and counters
111 //_____________________________________________________________________
112 AliCaloCalibSignal::~AliCaloCalibSignal()
117 //_____________________________________________________________________
118 void AliCaloCalibSignal::DeleteTrees()
120 // delete what was created in the ctor (TTrees)
121 if (fTreeAmpVsTime) delete fTreeAmpVsTime;
122 if (fTreeAvgAmpVsTime) delete fTreeAvgAmpVsTime;
123 if (fTreeLEDAmpVsTime) delete fTreeLEDAmpVsTime;
124 if (fTreeLEDAvgAmpVsTime) delete fTreeLEDAvgAmpVsTime;
125 // and reset pointers
126 fTreeAmpVsTime = NULL;
127 fTreeAvgAmpVsTime = NULL;
128 fTreeLEDAmpVsTime = NULL;
129 fTreeLEDAvgAmpVsTime = NULL;
135 //_____________________________________________________________________
136 //AliCaloCalibSignal::AliCaloCalibSignal(const AliCaloCalibSignal &sig) :
138 // fDetType(sig.GetDetectorType()),
139 // fColumns(sig.GetColumns()),
140 // fRows(sig.GetRows()),
141 // fLEDRefs(sig.GetLEDRefs()),
142 // fModules(sig.GetModules()),
143 // fCaloString(sig.GetCaloString()),
144 // fMapping(), //! note that we are not copying the map info
145 // fRunNumber(sig.GetRunNumber()),
146 // fStartTime(sig.GetStartTime()),
147 // fAmpCut(sig.GetAmpCut()),
148 // fReqFractionAboveAmpCutVal(sig.GetReqFractionAboveAmpCutVal()),
149 // fReqFractionAboveAmp(sig.GetReqFractionAboveAmp()),
150 // fAmpCutLEDRef(sig.GetAmpCutLEDRef()),
151 // fReqLEDRefAboveAmpCutVal(sig.GetReqLEDRefAboveAmpCutVal()),
152 // fHour(sig.GetHour()),
153 // fLatestHour(sig.GetLatestHour()),
154 // fUseAverage(sig.GetUseAverage()),
155 // fSecInAverage(sig.GetSecInAverage()),
156 // fNEvents(sig.GetNEvents()),
157 // fNAcceptedEvents(sig.GetNAcceptedEvents()),
159 // fTreeAvgAmpVsTime(),
160 // fTreeLEDAmpVsTime(),
161 // fTreeLEDAvgAmpVsTime()
163 // // also the TTree contents
165 // for (Int_t i = 0; i<fgkMaxTowers; i++) {
166 // fNHighGain[i] = sig.fNHighGain[i];
167 // fNLowGain[i] = sig.fNLowGain[i];
169 // for (Int_t i = 0; i<(2*fgkMaxRefs); i++) {
170 // fNRef[i] = sig.fNRef[i];
176 // assignment operator; use copy ctor to make life easy..
177 //_____________________________________________________________________
178 //AliCaloCalibSignal& AliCaloCalibSignal::operator = (const AliCaloCalibSignal &source)
180 // // assignment operator; use copy ctor
181 // if (&source == this) return *this;
183 // new (this) AliCaloCalibSignal(source);
187 //_____________________________________________________________________
188 void AliCaloCalibSignal::CreateTrees()
191 // first, regular version
192 fTreeAmpVsTime = new TTree("fTreeAmpVsTime","Amplitude vs. Time Tree Variables");
194 fTreeAmpVsTime->Branch("fChannelNum", &fChannelNum, "fChannelNum/I");
195 fTreeAmpVsTime->Branch("fHour", &fHour, "fHour/D");
196 fTreeAmpVsTime->Branch("fAmp", &fAmp, "fAmp/D");
198 // then, average version
199 fTreeAvgAmpVsTime = new TTree("fTreeAvgAmpVsTime","Average Amplitude vs. Time Tree Variables");
201 fTreeAvgAmpVsTime->Branch("fChannelNum", &fChannelNum, "fChannelNum/I");
202 fTreeAvgAmpVsTime->Branch("fHour", &fHour, "fHour/D");
203 fTreeAvgAmpVsTime->Branch("fAvgAmp", &fAvgAmp, "fAvgAmp/D");
204 fTreeAvgAmpVsTime->Branch("fRMS", &fRMS, "fRMS/D");
206 // then same for LED..
207 fTreeLEDAmpVsTime = new TTree("fTreeLEDAmpVsTime","LED Amplitude vs. Time Tree Variables");
208 fTreeLEDAmpVsTime->Branch("fRefNum", &fRefNum, "fRefNum/I");
209 fTreeLEDAmpVsTime->Branch("fHour", &fHour, "fHour/D");
210 fTreeLEDAmpVsTime->Branch("fAmp", &fAmp, "fAmp/D");
212 fTreeLEDAvgAmpVsTime = new TTree("fTreeLEDAvgAmpVsTime","Average LED Amplitude vs. Time Tree Variables");
213 fTreeLEDAvgAmpVsTime->Branch("fRefNum", &fRefNum, "fRefNum/I");
214 fTreeLEDAvgAmpVsTime->Branch("fHour", &fHour, "fHour/D");
215 fTreeLEDAvgAmpVsTime->Branch("fAvgAmp", &fAvgAmp, "fAvgAmp/D");
216 fTreeLEDAvgAmpVsTime->Branch("fRMS", &fRMS, "fRMS/D");
221 //_____________________________________________________________________
222 void AliCaloCalibSignal::ResetInfo()
223 { // reset trees and counters
224 Zero(); // set all counters to 0
225 DeleteTrees(); // delete previous stuff
226 CreateTrees(); // and create some new ones
230 //_____________________________________________________________________
231 void AliCaloCalibSignal::Zero()
233 // set all counters to 0; not cuts etc. though
237 fNAcceptedEvents = 0;
239 // Set the number of points for each tower: Amp vs. Time
240 memset(fNHighGain, 0, sizeof(fNHighGain));
241 memset(fNLowGain, 0, sizeof(fNLowGain));
243 memset(fNRef, 0, sizeof(fNRef));
248 //_____________________________________________________________________
249 Bool_t AliCaloCalibSignal::CheckFractionAboveAmp(const int *iAmpVal,
251 { // check fraction of towers, per column, that are above amplitude cut
252 Bool_t returnCode = false;
256 for (int i = 0; i<fModules; i++) {
257 for (int j = 0; j<fColumns; j++) {
259 for (int k = 0; k<fRows; k++) {
260 iTowerNum = GetTowerNum(i,j,k);
261 if (iAmpVal[iTowerNum] > fAmpCut) {
265 resultArray[i*fColumns +j] = 0; // init. to denied
267 fraction = (1.0*nAbove) / fRows;
269 printf("DS mod %d col %d nAbove %d fraction %3.2f\n",
270 i, j, nAbove, fraction);
272 if (fraction > fReqFractionAboveAmpCutVal) {
273 resultArray[i*fColumns + j] = nAbove;
284 //_____________________________________________________________________
285 Bool_t AliCaloCalibSignal::CheckLEDRefAboveAmp(const int *iAmpVal,
287 { // check which LEDRef/Mon strips are above amplitude cut
288 Bool_t returnCode = false;
291 int gain = 1; // look at high gain; this should be rather saturated usually..
292 for (int i = 0; i<fModules; i++) {
293 for (int j = 0; j<fLEDRefs; j++) {
294 iRefNum = GetRefNum(i, j, gain);
295 if (iAmpVal[iRefNum] > fAmpCutLEDRef) {
296 resultArray[i*fLEDRefs +j] = 1; // enough signal
300 resultArray[i*fLEDRefs +j] = 0; // not enough signal
304 printf("DS mod %d LEDRef %d ampVal %d\n",
305 i, j, iAmpVal[iRefNum]);
313 // Parameter/cut handling
314 //_____________________________________________________________________
315 void AliCaloCalibSignal::SetParametersFromFile(const char *parameterFile)
316 { // set parameters from file
317 static const string delimitor("::");
319 // open, check input file
320 ifstream in( parameterFile );
322 printf("in AliCaloCalibSignal::SetParametersFromFile - Using default/run_time parameters.\n");
326 // Note: this method is a bit more complicated than it really has to be
327 // - allowing for multiple entries per line, arbitrary order of the
328 // different variables etc. But I wanted to try and do this in as
329 // correct a C++ way as I could (as an exercise).
333 while ((in.rdstate() & ios::failbit) == 0 ) {
335 // Read into the raw char array and then construct a string
336 // to do the searching
337 in.getline(readline, 1024);
338 istringstream s(readline);
340 while ( ( s.rdstate() & ios::failbit ) == 0 ) {
345 // check stream status
346 if( ( s.rdstate() & ios::failbit ) == ios::failbit ) break;
348 // skip rest of line if comments found
349 if( keyValue.substr( 0, 2 ) == "//" ) break;
351 // look for "::" in keyValue pair
352 size_t position = keyValue.find( delimitor );
353 if( position == string::npos ) {
354 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
357 // split keyValue pair
358 string key( keyValue.substr( 0, position ) );
359 string value( keyValue.substr( position+delimitor.size(),
360 keyValue.size()-delimitor.size() ) );
362 // check value does not contain a new delimitor
363 if( value.find( delimitor ) != string::npos ) {
364 printf("wrong format for key::value pair: %s\n", keyValue.c_str());
367 // debug: check key value pair
368 // printf("AliCaloCalibSignal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
370 // if the key matches with something we expect, we assign the new value
371 if ( (key == "fAmpCut") || (key == "fReqFractionAboveAmpCutVal") ||
372 (key == "fAmpCutLEDRef") || (key == "fSecInAverage") ) {
373 istringstream iss(value);
374 printf("AliCaloCalibSignal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
376 if (key == "fAmpCut") {
379 else if (key == "fReqFractionAboveAmpCutVal") {
380 iss >> fReqFractionAboveAmpCutVal;
382 else if (key == "fAmpCutLEDRef") {
383 iss >> fAmpCutLEDRef;
385 else if (key == "fSecInAverage") {
386 iss >> fSecInAverage;
388 } // some match found/expected
397 //_____________________________________________________________________
398 void AliCaloCalibSignal::WriteParametersToFile(const char *parameterFile)
399 { // write parameters to file
400 static const string delimitor("::");
401 ofstream out( parameterFile );
402 out << "// " << parameterFile << endl;
403 out << "fAmpCut" << "::" << fAmpCut << endl;
404 out << "fReqFractionAboveAmpCutVal" << "::" << fReqFractionAboveAmpCutVal << endl;
405 out << "fAmpCutLEDRef" << "::" << fAmpCutLEDRef << endl;
406 out << "fSecInAverage" << "::" << fSecInAverage << endl;
412 //_____________________________________________________________________
413 Bool_t AliCaloCalibSignal::AddInfo(const AliCaloCalibSignal *sig)
415 // note/FIXME: we are not yet adding correctly the info for fN{HighGain,LowGain,Ref} here - but consider this a feature for now (20080905): we'll do Analyze() unless entries were found for a tower in this original object.
417 // add info from sig's TTrees to ours..
418 TTree *sigAmp = sig->GetTreeAmpVsTime();
419 TTree *sigAvgAmp = sig->GetTreeAvgAmpVsTime();
421 // we could try some merging via TList or what also as a more elegant approach
422 // but I wanted with the stupid/simple and hopefully safe approach of looping
423 // over what we want to add..
425 // associate variables for sigAmp and sigAvgAmp:
426 sigAmp->SetBranchAddress("fChannelNum",&fChannelNum);
427 sigAmp->SetBranchAddress("fHour",&fHour);
428 sigAmp->SetBranchAddress("fAmp",&fAmp);
430 // loop over the trees.. note that since we use the same variables we should not need
431 // to do any assignments between the getting and filling
432 for (int i=0; i<sigAmp->GetEntries(); i++) {
434 fTreeAmpVsTime->Fill();
437 sigAvgAmp->SetBranchAddress("fChannelNum",&fChannelNum);
438 sigAvgAmp->SetBranchAddress("fHour",&fHour);
439 sigAvgAmp->SetBranchAddress("fAvgAmp",&fAvgAmp);
440 sigAvgAmp->SetBranchAddress("fRMS",&fRMS);
442 for (int i=0; i<sigAvgAmp->GetEntries(); i++) {
443 sigAvgAmp->GetEntry(i);
444 fTreeAvgAmpVsTime->Fill();
448 TTree *sigLEDAmp = sig->GetTreeLEDAmpVsTime();
449 TTree *sigLEDAvgAmp = sig->GetTreeLEDAvgAmpVsTime();
451 // associate variables for sigAmp and sigAvgAmp:
452 sigLEDAmp->SetBranchAddress("fRefNum",&fRefNum);
453 sigLEDAmp->SetBranchAddress("fHour",&fHour);
454 sigLEDAmp->SetBranchAddress("fAmp",&fAmp);
456 // loop over the trees.. note that since we use the same variables we should not need
457 // to do any assignments between the getting and filling
458 for (int i=0; i<sigLEDAmp->GetEntries(); i++) {
459 sigLEDAmp->GetEntry(i);
460 fTreeLEDAmpVsTime->Fill();
463 sigLEDAvgAmp->SetBranchAddress("fRefNum",&fRefNum);
464 sigLEDAvgAmp->SetBranchAddress("fHour",&fHour);
465 sigLEDAvgAmp->SetBranchAddress("fAvgAmp",&fAvgAmp);
466 sigLEDAvgAmp->SetBranchAddress("fRMS",&fRMS);
468 for (int i=0; i<sigLEDAvgAmp->GetEntries(); i++) {
469 sigLEDAvgAmp->GetEntry(i);
470 fTreeLEDAvgAmpVsTime->Fill();
473 // We should also copy other pieces of info: counters and parameters
474 // (not number of columns and rows etc which should be the same)
475 // note that I just assign them here rather than Add them, but we
476 // normally just Add (e.g. in Preprocessor) one object so this should be fine.
477 fRunNumber = sig->GetRunNumber();
478 fStartTime = sig->GetStartTime();
479 fAmpCut = sig->GetAmpCut();
480 fReqFractionAboveAmpCutVal = sig->GetReqFractionAboveAmpCutVal();
481 fReqFractionAboveAmp = sig->GetReqFractionAboveAmp();
482 fAmpCutLEDRef = sig->GetAmpCutLEDRef();
483 fReqLEDRefAboveAmpCutVal = sig->GetReqLEDRefAboveAmpCutVal();
484 fHour = sig->GetHour();
485 fLatestHour = sig->GetLatestHour();
486 fUseAverage = sig->GetUseAverage();
487 fSecInAverage = sig->GetSecInAverage();
488 fNEvents = sig->GetNEvents();
489 fNAcceptedEvents = sig->GetNAcceptedEvents();
491 return kTRUE;//We hopefully succesfully added info from the supplied object
494 //_____________________________________________________________________
495 Bool_t AliCaloCalibSignal::ProcessEvent(AliRawReader *rawReader)
497 // if fMapping is NULL the rawstream will crate its own mapping
498 AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
499 if (fDetType == kEmCal) {
500 rawReader->Select("EMCAL", 0, AliEMCALGeoParams::fgkLastAltroDDL) ; //select EMCAL DDL range
502 return ProcessEvent( &rawStream, rawReader->GetTimestamp() );
505 //_____________________________________________________________________
506 Bool_t AliCaloCalibSignal::ProcessEvent(AliCaloRawStreamV3 *in, UInt_t Timestamp)
508 // Method to process=analyze one event in the data stream
509 if (!in) return kFALSE; //Return right away if there's a null pointer
511 fNEvents++; // one more event
513 // use maximum numbers to set array sizes
514 int iAmpValHighGain[fgkMaxTowers];
515 int iAmpValLowGain[fgkMaxTowers];
516 memset(iAmpValHighGain, 0, sizeof(iAmpValHighGain));
517 memset(iAmpValLowGain, 0, sizeof(iAmpValLowGain));
519 // also for LED reference
520 int iLEDAmpVal[fgkMaxRefs * 2]; // factor 2 is for the two gain values
521 memset(iLEDAmpVal, 0, sizeof(iLEDAmpVal));
523 int sample = 0; // temporary value
524 int gain = 0; // high or low gain
526 // Number of Low and High gain, and LED Ref, channels for this event:
531 int iTowerNum = 0; // array index for regular towers
532 int iRefNum = 0; // array index for LED references
534 // loop first to get the fraction of channels with amplitudes above cut
536 while (in->NextDDL()) {
537 while (in->NextChannel()) {
540 int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
543 while (in->NextBunch()) {
544 const UShort_t *sig = in->GetSignals();
545 nsamples += in->GetBunchLength();
546 for (Int_t i = 0; i < in->GetBunchLength(); i++) {
549 // check if it's a min or max value
550 if (sample < min) min = sample;
551 if (sample > max) max = sample;
553 } // loop over samples in bunch
554 } // loop over bunches
556 if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
558 gain = -1; // init to not valid value
559 //If we're here then we're done with this tower
560 if ( in->IsLowGain() ) {
563 else if ( in->IsHighGain() ) {
566 else if ( in->IsLEDMonData() ) {
567 gain = in->GetRow(); // gain coded in (in RCU/Altro mapping) as Row info for LED refs..
570 // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
571 int arrayPos = in->GetModule(); //The modules are numbered starting from 0
573 if (arrayPos < 0 || arrayPos >= fModules) {
574 printf("AliCaloCalibSignal::ProcessEvent = Oh no: arrayPos = %i.\n", arrayPos);
578 if ( in->IsHighGain() || in->IsLowGain() ) { // regular tower
579 // get tower number for AmpVal array
580 iTowerNum = GetTowerNum(arrayPos, in->GetColumn(), in->GetRow());
583 // fill amplitude into the array
584 iAmpValLowGain[iTowerNum] = max - min;
587 else if (gain==1) {//fill the high gain ones
588 // fill amplitude into the array
589 iAmpValHighGain[iTowerNum] = max - min;
593 else if ( in->IsLEDMonData() ) { // LED ref.;
594 // strip # is coded is 'column' in the channel maps
595 iRefNum = GetRefNum(arrayPos, in->GetColumn(), gain);
596 iLEDAmpVal[iRefNum] = max - min;
600 } // nsamples>0 check, some data found for this channel; not only trailer/header
601 } // end while over channel
603 }//end while over DDL's, of input stream
605 in->Reset(); // just in case the next customer forgets to check if the stream was reset..
607 // now check if it was an LED event, using the LED Reference info per strip
609 // by default all columns are accepted (init check to > 0)
610 int checkResultArray[AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALCols];
611 for (int ia=0; ia<(AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALCols); ia++) {
612 checkResultArray[ia] = 1;
614 if (fReqFractionAboveAmp) {
617 ok = CheckFractionAboveAmp(iAmpValHighGain, checkResultArray);
619 if (!ok) return false; // skip event
622 // by default all columns are accepted (init check to > 0)
623 int checkResultArrayLEDRef[AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALLEDRefs];
624 for (int ia=0; ia<(AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALLEDRefs); ia++) {
625 checkResultArrayLEDRef[ia] = 1;
627 if (fReqLEDRefAboveAmpCutVal) {
629 if (nLEDRefChan > 0) {
630 ok = CheckLEDRefAboveAmp(iLEDAmpVal, checkResultArrayLEDRef);
632 if (!ok) return false; // skip event
635 fNAcceptedEvents++; // one more event accepted
637 if (fStartTime == 0) { // if start-timestamp wasn't set,we'll pick it up from the first event we encounter
638 fStartTime = Timestamp;
641 fHour = (Timestamp - fStartTime)/(double)fgkNumSecInHr;
642 if (fLatestHour < fHour) {
646 // it is a led event, now fill TTree
647 // We also do the activity check for LEDRefs/Strips, but need to translate between column
648 // and strip indices for that; based on these relations:
649 // iStrip = AliEMCALGeoParams::GetStripModule(iSM, iCol);
650 // iStrip = (iSM%2==0) ? iCol/2 : AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - iCol/2;
652 // iColFirst = (iSM%2==0) ? iStrip*2 : (AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - iStrip)*2;
654 for(int i=0; i<fModules; i++){
655 for(int j=0; j<fColumns; j++) {
656 int iStrip = (i%2==0) ? j/2 : AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - j/2;
657 if (checkResultArray[i*fColumns + j]>0 && checkResultArrayLEDRef[i*fLEDRefs + iStrip]>0) { // column passed check
658 for(int k=0; k<fRows; k++){
660 iTowerNum = GetTowerNum(i, j, k);
662 if(iAmpValHighGain[iTowerNum]) {
663 fAmp = iAmpValHighGain[iTowerNum];
664 fChannelNum = GetChannelNum(i,j,k,1);
665 fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValHighGain[iTowerNum]);
666 fNHighGain[iTowerNum]++;
668 if(iAmpValLowGain[iTowerNum]) {
669 fAmp = iAmpValLowGain[iTowerNum];
670 fChannelNum = GetChannelNum(i,j,k,0);
671 fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValLowGain[iTowerNum]);
672 fNLowGain[iTowerNum]++;
675 } // column passed check, and LED Ref for strip passed check (if any)
679 for(int j=0; j<fLEDRefs; j++){
680 int iColFirst = (i%2==0) ? j*2 : (AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - j)*2; //CHECKME!!!
681 if ( ((checkResultArray[i*fColumns + iColFirst]>0) || (checkResultArray[i*fColumns + iColFirst + 1]>0)) && // at least one column in strip passed check
682 (checkResultArrayLEDRef[i*fLEDRefs + j]>0) ) { // and LED Ref passed checks
683 for (gain=0; gain<2; gain++) {
684 fRefNum = GetRefNum(i, j, gain);
685 if (iLEDAmpVal[fRefNum]) {
686 fAmp = iLEDAmpVal[fRefNum];
687 fTreeLEDAmpVsTime->Fill();//fRefNum,fHour,fAmp);
691 } // at least one column in strip passed check, and LED Ref passed check (if any)
699 //_____________________________________________________________________
700 Bool_t AliCaloCalibSignal::Save(TString fileName)
702 //Saves all the TTrees to the designated file
704 TFile destFile(fileName, "recreate");
706 if (destFile.IsZombie()) {
713 fTreeAmpVsTime->Write();
714 fTreeLEDAmpVsTime->Write();
716 Analyze(); // get the latest and greatest averages
717 fTreeAvgAmpVsTime->Write();
718 fTreeLEDAvgAmpVsTime->Write();
726 //_____________________________________________________________________
727 Bool_t AliCaloCalibSignal::Analyze()
729 // Fill the tree holding the average values
730 if (!fUseAverage) { return kFALSE; }
732 // Reset the average TTree if Analyze has already been called earlier,
733 // meaning that the TTree could have been partially filled
734 if (fTreeAvgAmpVsTime->GetEntries() > 0) {
735 fTreeAvgAmpVsTime->Reset();
738 //0: setup variables for the TProfile plots that we'll use to do the averages
742 if (fSecInAverage > 0) {
743 numProfBins = (int)( (fLatestHour*fgkNumSecInHr)/fSecInAverage + 1 ); // round-off
745 numProfBins += 2; // add extra buffer : first and last
746 double binSize = 1.0*fSecInAverage / fgkNumSecInHr;
748 timeMax = timeMin + numProfBins*binSize;
750 //1: set up TProfiles for the towers that had data
751 TProfile * profile[fgkMaxTowers*2]; // *2 is since we include both high and low gains
752 memset(profile, 0, sizeof(profile));
753 const Int_t buffersize = 200;
754 char name[buffersize]; // for profile id and title
757 for (int i = 0; i<fModules; i++) {
758 for (int ic=0; ic<fColumns; ic++){
759 for (int ir=0; ir<fRows; ir++) {
761 iTowerNum = GetTowerNum(i, ic, ir);
763 if (fNHighGain[iTowerNum] > 0) {
764 fChannelNum = GetChannelNum(i, ic, ir, 1);
765 snprintf(name,buffersize,"profileChan%d", fChannelNum);
766 profile[fChannelNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
770 if (fNLowGain[iTowerNum] > 0) {
771 fChannelNum = GetChannelNum(i, ic, ir, 0);
772 snprintf(name,buffersize,"profileChan%d", fChannelNum);
773 profile[fChannelNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
780 //2: fill profiles by looping over tree
781 // Set addresses for tree-readback also
782 fTreeAmpVsTime->SetBranchAddress("fChannelNum", &fChannelNum);
783 fTreeAmpVsTime->SetBranchAddress("fHour", &fHour);
784 fTreeAmpVsTime->SetBranchAddress("fAmp", &fAmp);
786 for (int ient=0; ient<fTreeAmpVsTime->GetEntries(); ient++) {
787 fTreeAmpVsTime->GetEntry(ient);
788 if (profile[fChannelNum]) {
789 // profile should always have been created above, for active channels
790 profile[fChannelNum]->Fill(fHour, fAmp);
794 // re-associating the branch addresses here seems to be needed for OK 'average' storage
795 fTreeAvgAmpVsTime->SetBranchAddress("fChannelNum", &fChannelNum);
796 fTreeAvgAmpVsTime->SetBranchAddress("fHour", &fHour);
797 fTreeAvgAmpVsTime->SetBranchAddress("fAvgAmp", &fAvgAmp);
798 fTreeAvgAmpVsTime->SetBranchAddress("fRMS", &fRMS);
800 //3: fill avg tree by looping over the profiles
801 for (fChannelNum = 0; fChannelNum<(fgkMaxTowers*2); fChannelNum++) {
802 if (profile[fChannelNum]) { // profile was created
803 if (profile[fChannelNum]->GetEntries() > 0) { // profile had some entries
804 for(int it=0; it<numProfBins; it++) {
805 if (profile[fChannelNum]->GetBinEntries(it+1) > 0) {
806 fAvgAmp = profile[fChannelNum]->GetBinContent(it+1);
807 fHour = profile[fChannelNum]->GetBinCenter(it+1);
808 fRMS = profile[fChannelNum]->GetBinError(it+1);
809 fTreeAvgAmpVsTime->Fill();
810 } // some entries for this bin
812 } // some entries for this profile
814 } // loop over all possible channels
817 // and finally, go through same exercise for LED also..
819 //1: set up TProfiles for the towers that had data
820 TProfile * profileLED[fgkMaxRefs*2]; // *2 is since we include both high and low gains
821 memset(profileLED, 0, sizeof(profileLED));
823 for (int i = 0; i<fModules; i++) {
824 for(int j=0; j<fLEDRefs; j++){
825 for (int gain=0; gain<2; gain++) {
826 fRefNum = GetRefNum(i, j, gain);
827 if (fNRef[fRefNum] > 0) {
828 snprintf(name, buffersize, "profileLEDRef%d", fRefNum);
829 profileLED[fRefNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
835 //2: fill profiles by looping over tree
836 // Set addresses for tree-readback also
837 fTreeLEDAmpVsTime->SetBranchAddress("fRefNum", &fRefNum);
838 fTreeLEDAmpVsTime->SetBranchAddress("fHour", &fHour);
839 fTreeLEDAmpVsTime->SetBranchAddress("fAmp", &fAmp);
841 for (int ient=0; ient<fTreeLEDAmpVsTime->GetEntries(); ient++) {
842 fTreeLEDAmpVsTime->GetEntry(ient);
843 if (profileLED[fRefNum]) {
844 // profile should always have been created above, for active channels
845 profileLED[fRefNum]->Fill(fHour, fAmp);
849 // re-associating the branch addresses here seems to be needed for OK 'average' storage
850 fTreeLEDAvgAmpVsTime->SetBranchAddress("fRefNum", &fRefNum);
851 fTreeLEDAvgAmpVsTime->SetBranchAddress("fHour", &fHour);
852 fTreeLEDAvgAmpVsTime->SetBranchAddress("fAvgAmp", &fAvgAmp);
853 fTreeLEDAvgAmpVsTime->SetBranchAddress("fRMS", &fRMS);
855 //3: fill avg tree by looping over the profiles
856 for (fRefNum = 0; fRefNum<(fgkMaxRefs*2); fRefNum++) {
857 if (profileLED[fRefNum]) { // profile was created
858 if (profileLED[fRefNum]->GetEntries() > 0) { // profile had some entries
859 for(int it=0; it<numProfBins; it++) {
860 if (profileLED[fRefNum]->GetBinEntries(it+1) > 0) {
861 fAvgAmp = profileLED[fRefNum]->GetBinContent(it+1);
862 fHour = profileLED[fRefNum]->GetBinCenter(it+1);
863 fRMS = profileLED[fRefNum]->GetBinError(it+1);
864 fTreeLEDAvgAmpVsTime->Fill();
865 } // some entries for this bin
867 } // some entries for this profile
869 } // loop over all possible channels
876 //_____________________________________________________________________
877 Bool_t AliCaloCalibSignal::DecodeChannelNum(const int chanId,
878 int *imod, int *icol, int *irow, int *igain) const
879 { // return the module, column, row, and gain for a given channel number
880 *igain = chanId/(fModules*fColumns*fRows);
881 *imod = (chanId/(fColumns*fRows)) % fModules;
882 *icol = (chanId/fRows) % fColumns;
883 *irow = chanId % fRows;
887 //_____________________________________________________________________
888 Bool_t AliCaloCalibSignal::DecodeRefNum(const int refId,
889 int *imod, int *istripMod, int *igain) const
890 { // return the module, stripModule, and gain for a given reference number
891 *igain = refId/(fModules*fLEDRefs);
892 *imod = (refId/(fLEDRefs)) % fModules;
893 *istripMod = refId % fLEDRefs;