]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliCaloCalibSignal.cxx
add getter for the array with digits and initialization inside, just for analysis...
[u/mrichter/AliRoot.git] / EMCAL / AliCaloCalibSignal.cxx
index ab0a1b0d3676e134726c536df5e6453f0d6f15bc..4f2d2020496a48efb90f920a99b92866f8faf50a 100644 (file)
 */
 // fed an event:
 //  fSignals->ProcessEvent(fCaloRawStream,fRawEventHeaderBase);
-// asked to draw graphs or profiles:
-//  fSignals->GetGraphAmpVsTimeHighGain(module,column,row)->Draw("ap");
-// or
-//  fSignals->GetProfAmpVsTimeHighGain(module,column,row)->Draw();
+// get some info:
+//  fSignals->GetXXX..()
 // etc.
 //________________________________________________________________________
+#include <string>
+#include <sstream>
+#include <fstream>
 
+#include "TProfile.h"
 #include "TFile.h"
 
-#include "AliRawEventHeaderBase.h"
-#include "AliCaloRawStream.h"
+#include "AliRawReader.h"
+#include "AliCaloRawStreamV3.h"
 
 //The include file
 #include "AliCaloCalibSignal.h"
@@ -45,6 +47,13 @@ ClassImp(AliCaloCalibSignal)
 
 using namespace std;
 
+// variables for TTree filling; not sure if they should be static or not
+static int fChannelNum = 0; // for regular towers
+static int fRefNum = 0; // for LED
+static double fAmp = 0;
+static double fAvgAmp = 0;
+static double fRMS = 0;
+
 // ctor; initialize everything in order to avoid compiler warnings
 // put some reasonable defaults
 AliCaloCalibSignal::AliCaloCalibSignal(kDetType detectorType) :  
@@ -52,365 +61,655 @@ AliCaloCalibSignal::AliCaloCalibSignal(kDetType detectorType) :
   fDetType(kNone),
   fColumns(0),
   fRows(0),
+  fLEDRefs(0),
   fModules(0),
+  fCaloString(),
+  fMapping(NULL),
   fRunNumber(-1),
   fStartTime(0),
-  fAmpCut(50),
-  fReqFractionAboveAmpCutVal(0.8),
+  fAmpCut(40), // min. 40 ADC counts as default
+  fReqFractionAboveAmpCutVal(0.6), // 60% in a strip, per default
   fReqFractionAboveAmp(kTRUE),
+  fAmpCutLEDRef(100), // min. 100 ADC counts as default
+  fReqLEDRefAboveAmpCutVal(kTRUE),
   fHour(0),
   fLatestHour(0),
   fUseAverage(kTRUE),
   fSecInAverage(1800), 
+  fDownscale(10), 
   fNEvents(0),
-  fNAcceptedEvents(0)
+  fNAcceptedEvents(0),
+  fTreeAmpVsTime(NULL),
+  fTreeAvgAmpVsTime(NULL),
+  fTreeLEDAmpVsTime(NULL),
+  fTreeLEDAvgAmpVsTime(NULL)
 {
   //Default constructor. First we set the detector-type related constants.
   if (detectorType == kPhos) {
     fColumns = fgkPhosCols;
     fRows = fgkPhosRows;
+    fLEDRefs = fgkPhosLEDRefs;
     fModules = fgkPhosModules;
+    fCaloString = "PHOS";
   } 
   else {
     //We'll just trust the enum to keep everything in line, so that if detectorType
     //isn't kPhos then it is kEmCal. Note, however, that this is not necessarily the
     //case, if someone intentionally gives another number
-    fColumns = fgkEmCalCols;
-    fRows = fgkEmCalRows;
-    fModules = fgkEmCalModules;
+    fColumns = AliEMCALGeoParams::fgkEMCALCols;
+    fRows = AliEMCALGeoParams::fgkEMCALRows;
+    fLEDRefs = AliEMCALGeoParams::fgkEMCALLEDRefs;
+    fModules = AliEMCALGeoParams::fgkEMCALModules;
+    fCaloString = "EMCAL";
   }
 
   fDetType = detectorType;
 
-  // Set the number of points for each Amp vs. Time graph to 0
-  memset(fNHighGain, 0, sizeof(fNHighGain));
-  memset(fNLowGain, 0, sizeof(fNLowGain));
-
-  CreateGraphs(); // set up the TGraphs
-
-  // init TProfiles to NULL=0 also
-  memset(fProfAmpVsTimeHighGain, 0, sizeof(fProfAmpVsTimeHighGain));
-  memset(fProfAmpVsTimeLowGain, 0, sizeof(fProfAmpVsTimeLowGain));
+  ResetInfo(); // trees and counters
 } 
 
 // dtor
 //_____________________________________________________________________
 AliCaloCalibSignal::~AliCaloCalibSignal()
 {
-  ClearObjects();
+  DeleteTrees();
 }
 
 //_____________________________________________________________________
-void AliCaloCalibSignal::ClearObjects()
+void AliCaloCalibSignal::DeleteTrees()
 {
-  // delete what was created in the ctor (TGraphs), and possible later (TProfiles)
-  for (int i=0; i<fgkMaxTowers; i++) {
-    if ( fGraphAmpVsTimeHighGain[i] ) { delete fGraphAmpVsTimeHighGain[i]; }
-    if ( fGraphAmpVsTimeLowGain[i] ) { delete fGraphAmpVsTimeLowGain[i]; }
-    if ( fProfAmpVsTimeHighGain[i] ) { delete fProfAmpVsTimeHighGain[i]; }
-    if ( fProfAmpVsTimeLowGain[i] ) { delete fProfAmpVsTimeLowGain[i]; }
-  }
-  // set pointers
-  memset(fGraphAmpVsTimeHighGain, 0, sizeof(fGraphAmpVsTimeHighGain));
-  memset(fGraphAmpVsTimeLowGain, 0, sizeof(fGraphAmpVsTimeLowGain));
-  memset(fProfAmpVsTimeHighGain, 0, sizeof(fProfAmpVsTimeHighGain));
-  memset(fProfAmpVsTimeLowGain, 0, sizeof(fProfAmpVsTimeLowGain));
+  // delete what was created in the ctor (TTrees)
+  if (fTreeAmpVsTime) delete fTreeAmpVsTime;
+  if (fTreeAvgAmpVsTime) delete fTreeAvgAmpVsTime;
+  if (fTreeLEDAmpVsTime) delete fTreeLEDAmpVsTime;
+  if (fTreeLEDAvgAmpVsTime) delete fTreeLEDAvgAmpVsTime;
+  // and reset pointers
+  fTreeAmpVsTime = NULL;
+  fTreeAvgAmpVsTime = NULL;
+  fTreeLEDAmpVsTime = NULL;
+  fTreeLEDAvgAmpVsTime = NULL;
 
   return;
 }
 
 // copy ctor
 //_____________________________________________________________________
-AliCaloCalibSignal::AliCaloCalibSignal(const AliCaloCalibSignal &sig) :
-  TObject(sig),
-  fDetType(sig.GetDetectorType()),
-  fColumns(sig.GetColumns()),
-  fRows(sig.GetRows()),
-  fModules(sig.GetModules()),
-  fRunNumber(sig.GetRunNumber()),
-  fStartTime(sig.GetStartTime()),
-  fAmpCut(sig.GetAmpCut()),
-  fReqFractionAboveAmpCutVal(sig.GetReqFractionAboveAmpCutVal()),
-  fReqFractionAboveAmp(sig.GetReqFractionAboveAmp()),
-  fHour(sig.GetHour()),
-  fLatestHour(sig.GetLatestHour()),
-  fUseAverage(sig.GetUseAverage()),
-  fSecInAverage(sig.GetSecInAverage()),
-  fNEvents(sig.GetNEvents()),
-  fNAcceptedEvents(sig.GetNAcceptedEvents())
-{
-  // also the TGraph contents
-  AddInfo(&sig);
-}
-
+//AliCaloCalibSignal::AliCaloCalibSignal(const AliCaloCalibSignal &sig) :
+//  TObject(sig),
+//  fDetType(sig.GetDetectorType()),
+//  fColumns(sig.GetColumns()),
+//  fRows(sig.GetRows()),
+//  fLEDRefs(sig.GetLEDRefs()),
+//  fModules(sig.GetModules()),
+//  fCaloString(sig.GetCaloString()),
+//  fMapping(), //! note that we are not copying the map info
+//  fRunNumber(sig.GetRunNumber()),
+//  fStartTime(sig.GetStartTime()),
+//  fAmpCut(sig.GetAmpCut()),
+//  fReqFractionAboveAmpCutVal(sig.GetReqFractionAboveAmpCutVal()),
+//  fReqFractionAboveAmp(sig.GetReqFractionAboveAmp()),
+//  fAmpCutLEDRef(sig.GetAmpCutLEDRef()),
+//  fReqLEDRefAboveAmpCutVal(sig.GetReqLEDRefAboveAmpCutVal()),
+//  fHour(sig.GetHour()),
+//  fLatestHour(sig.GetLatestHour()),
+//  fUseAverage(sig.GetUseAverage()),
+//  fSecInAverage(sig.GetSecInAverage()),
+//  fDownscale(sig.GetDownscale()),
+//  fNEvents(sig.GetNEvents()),
+//  fNAcceptedEvents(sig.GetNAcceptedEvents()),
+//  fTreeAmpVsTime(),
+//  fTreeAvgAmpVsTime(),
+//  fTreeLEDAmpVsTime(),
+//  fTreeLEDAvgAmpVsTime()
+//{
+//  // also the TTree contents
+//  AddInfo(&sig);
+//  for (Int_t i = 0; i<fgkMaxTowers; i++) {
+//      fNHighGain[i] = sig.fNHighGain[i];
+//      fNLowGain[i]  = sig.fNLowGain[i]; 
+//  }
+//  for (Int_t i = 0; i<(2*fgkMaxRefs); i++) {
+//    fNRef[i] = sig.fNRef[i]; 
+//  }
+//  
+//  
+//}
+//
 // assignment operator; use copy ctor to make life easy..
 //_____________________________________________________________________
-AliCaloCalibSignal& AliCaloCalibSignal::operator = (const AliCaloCalibSignal &source)
-{
-  // assignment operator; use copy ctor
-  if (&source == this) return *this;
-
-  new (this) AliCaloCalibSignal(source);
-  return *this;
-}
+//AliCaloCalibSignal& AliCaloCalibSignal::operator = (const AliCaloCalibSignal &source)
+//{
+//  // assignment operator; use copy ctor
+//  if (&source == this) return *this;
+//
+//  new (this) AliCaloCalibSignal(source);
+//  return *this;
+//}
 
 //_____________________________________________________________________
-void AliCaloCalibSignal::CreateGraphs()
+void AliCaloCalibSignal::CreateTrees()
 {
-  //Then, loop for the requested number of modules
-  TString title, name;
-  for (int i = 0; i < fModules; i++) {
-    
-    // Amplitude vs. Time graph for each channel
-    for(int ic=0;ic < fColumns;ic++){
-      for(int ir=0;ir < fRows;ir++){
-       
-       int id = GetTowerNum(i, ic, ir);
-       
-       // high gain graph
-       name = "fGraphAmpVsTimeHighGain_"; name += i;
-       name += "_"; name += ic;
-       name += "_"; name += ir;
-       title = "Amp vs. Time High Gain Mod "; title += i;
-       title += " Col "; title += ic;
-       title += " Row "; title += ir;
-       
-       fGraphAmpVsTimeHighGain[id] = new TGraph();
-       fGraphAmpVsTimeHighGain[id]->SetName(name);
-       fGraphAmpVsTimeHighGain[id]->SetTitle(title);
-       fGraphAmpVsTimeHighGain[id]->SetMarkerStyle(20);
-       
-       // Low Gain
-       name = "fGraphAmpVsTimeLowGain_"; name += i;
-       name += "_"; name += ic;
-       name += "_"; name += ir;
-       title = "Amp vs. Time Low Gain Mod "; title += i;
-       title += " Col "; title += ic;
-       title += " Row "; title += ir;
-       
-       fGraphAmpVsTimeLowGain[id] = new TGraph();
-       fGraphAmpVsTimeLowGain[id]->SetName(name);
-       fGraphAmpVsTimeLowGain[id]->SetTitle(title);
-       fGraphAmpVsTimeLowGain[id]->SetMarkerStyle(20);
-       
-      }
-    }
+  // initialize trees
+  // first, regular version
+  fTreeAmpVsTime = new TTree("fTreeAmpVsTime","Amplitude vs. Time Tree Variables");
+
+  fTreeAmpVsTime->Branch("fChannelNum", &fChannelNum, "fChannelNum/I");
+  fTreeAmpVsTime->Branch("fHour", &fHour, "fHour/D");
+  fTreeAmpVsTime->Branch("fAmp", &fAmp, "fAmp/D");
+
+  // then, average version
+  fTreeAvgAmpVsTime = new TTree("fTreeAvgAmpVsTime","Average Amplitude vs. Time Tree Variables");
+
+  fTreeAvgAmpVsTime->Branch("fChannelNum", &fChannelNum, "fChannelNum/I");
+  fTreeAvgAmpVsTime->Branch("fHour", &fHour, "fHour/D");
+  fTreeAvgAmpVsTime->Branch("fAvgAmp", &fAvgAmp, "fAvgAmp/D");
+  fTreeAvgAmpVsTime->Branch("fRMS", &fRMS, "fRMS/D");
+
+  // then same for LED..
+  fTreeLEDAmpVsTime = new TTree("fTreeLEDAmpVsTime","LED Amplitude vs. Time Tree Variables");
+  fTreeLEDAmpVsTime->Branch("fRefNum", &fRefNum, "fRefNum/I");
+  fTreeLEDAmpVsTime->Branch("fHour", &fHour, "fHour/D");
+  fTreeLEDAmpVsTime->Branch("fAmp", &fAmp, "fAmp/D");
+
+  fTreeLEDAvgAmpVsTime = new TTree("fTreeLEDAvgAmpVsTime","Average LED Amplitude vs. Time Tree Variables");
+  fTreeLEDAvgAmpVsTime->Branch("fRefNum", &fRefNum, "fRefNum/I");
+  fTreeLEDAvgAmpVsTime->Branch("fHour", &fHour, "fHour/D");
+  fTreeLEDAvgAmpVsTime->Branch("fAvgAmp", &fAvgAmp, "fAvgAmp/D");
+  fTreeLEDAvgAmpVsTime->Branch("fRMS", &fRMS, "fRMS/D");
 
-  }//end for nModules 
+  return;
 }
 
 //_____________________________________________________________________
-void AliCaloCalibSignal::Reset()
-{
+void AliCaloCalibSignal::ResetInfo()
+{ // reset trees and counters
   Zero(); // set all counters to 0
-  ClearObjects(); // delete previous TGraphs and TProfiles
-  CreateGraphs(); // and create some new ones
+  DeleteTrees(); // delete previous stuff
+  CreateTrees(); // and create some new ones
   return;
 }
 
 //_____________________________________________________________________
 void AliCaloCalibSignal::Zero()
 {
-  // set all counters to 0; not cuts etc.though
+  // set all counters to 0; not cuts etc. though
   fHour = 0;
   fLatestHour = 0;
   fNEvents = 0;
   fNAcceptedEvents = 0;
+
+  // Set the number of points for each tower: Amp vs. Time
+  memset(fNHighGain, 0, sizeof(fNHighGain));
+  memset(fNLowGain, 0, sizeof(fNLowGain));
+  // and LED reference
+  memset(fNRef, 0, sizeof(fNRef));
+
   return;
 }
 
 //_____________________________________________________________________
-Bool_t AliCaloCalibSignal::CheckFractionAboveAmp(int *AmpVal, int nTotChan)
-{
-  int nAbove = 0;
+Bool_t AliCaloCalibSignal::CheckFractionAboveAmp(const int *iAmpVal, 
+                                                int resultArray[]) const
+{ // check fraction of towers, per column, that are above amplitude cut
+  Bool_t returnCode = false;
     
-  int TowerNum = 0;
+  int iTowerNum = 0;
+  double fraction = 0;
   for (int i = 0; i<fModules; i++) {
     for (int j = 0; j<fColumns; j++) {
+      int nAbove = 0;
       for (int k = 0; k<fRows; k++) {
-       TowerNum = GetTowerNum(i,j,k);
-       if (AmpVal[TowerNum] > fAmpCut) { 
+       iTowerNum = GetTowerNum(i,j,k);
+       if (iAmpVal[iTowerNum] > fAmpCut) { 
          nAbove++;
        }
       }
+      resultArray[i*fColumns +j] = 0; // init. to denied
+      if (nAbove > 0) {
+       fraction = (1.0*nAbove) / fRows;
+       /*
+       printf("DS mod %d col %d nAbove %d fraction %3.2f\n",
+              i, j, nAbove, fraction);
+       */
+       if (fraction > fReqFractionAboveAmpCutVal) {
+         resultArray[i*fColumns + j] = nAbove;
+         returnCode = true;
+       }  
+      }
     }
-  }
+  } // modules loop
   
-  double fraction = (1.0*nAbove) / nTotChan;
+  return returnCode;
+}
+
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::CheckLEDRefAboveAmp(const int *iAmpVal, 
+                                              int resultArray[]) const
+{ // check which LEDRef/Mon strips are above amplitude cut
+  Bool_t returnCode = false;
+    
+  int iRefNum = 0;
+  int gain = 1; // look at high gain; this should be rather saturated usually..
+  for (int i = 0; i<fModules; i++) {
+    for (int j = 0; j<fLEDRefs; j++) {
+      iRefNum = GetRefNum(i, j, gain);
+      if (iAmpVal[iRefNum] > fAmpCutLEDRef) { 
+       resultArray[i*fLEDRefs +j] = 1; // enough signal
+       returnCode = true;
+      }
+      else {
+       resultArray[i*fLEDRefs +j] = 0; // not enough signal
+      }
+      
+      /*
+      printf("DS mod %d LEDRef %d ampVal %d\n",
+            i, j, iAmpVal[iRefNum]);
+      */
+    } // LEDRefs
+  } // modules loop
   
-  if (fraction > fReqFractionAboveAmpCutVal) {  
-    return true;
-  }
-  else return false;
+  return returnCode;
 }
 
+// Parameter/cut handling
 //_____________________________________________________________________
-Bool_t AliCaloCalibSignal::AddInfo(const AliCaloCalibSignal *sig)
-{
-  // just do this for the basic graphs/profiles that get filled in ProcessEvent
-  // may not have data for all channels, but let's just Add everything..
-  for (int i = 0; i < fModules; i++) {
-    for (int j = 0; j < fColumns; j++) {
-      for (int k = 0; k < fRows; k++) {
+void AliCaloCalibSignal::SetParametersFromFile(const char *parameterFile)
+{ // set parameters from file
+  static const string delimitor("::");
        
-       int id = GetTowerNum(i,j,k);
-       if(fUseAverage){
-         GetProfAmpVsTimeHighGain(id)->Add(sig->GetProfAmpVsTimeHighGain(id));
-         GetProfAmpVsTimeLowGain(id)->Add(sig->GetProfAmpVsTimeLowGain(id));
+  // open, check input file
+  ifstream in( parameterFile );
+  if( !in ) {
+    printf("in AliCaloCalibSignal::SetParametersFromFile - Using default/run_time parameters.\n");
+    return;
+  } 
+
+  // Note: this method is a bit more complicated than it really has to be
+  // - allowing for multiple entries per line, arbitrary order of the
+  // different variables etc. But I wanted to try and do this in as
+  // correct a C++ way as I could (as an exercise).
+
+  // read in
+  char readline[1024];
+  while ((in.rdstate() & ios::failbit) == 0 ) {
+    
+    // Read into the raw char array and then construct a string
+    // to do the searching
+    in.getline(readline, 1024);
+    istringstream s(readline);         
+               
+    while ( ( s.rdstate() & ios::failbit ) == 0 ) {
+                       
+      string keyValue; 
+      s >> keyValue;
+      
+      // check stream status
+      if( ( s.rdstate() & ios::failbit ) == ios::failbit ) break;
+                       
+      // skip rest of line if comments found
+      if( keyValue.substr( 0, 2 ) == "//" ) break;
+                       
+      // look for "::" in keyValue pair
+      size_t position = keyValue.find( delimitor );
+      if( position == string::npos ) {
+       printf("wrong format for key::value pair: %s\n", keyValue.c_str());
+      }
+                               
+      // split keyValue pair
+      string key( keyValue.substr( 0, position ) );
+      string value( keyValue.substr( position+delimitor.size(), 
+                                     keyValue.size()-delimitor.size() ) );
+                       
+      // check value does not contain a new delimitor
+      if( value.find( delimitor ) != string::npos ) {
+       printf("wrong format for key::value pair: %s\n", keyValue.c_str());
+      }
+      
+      // debug: check key value pair
+      // printf("AliCaloCalibSignal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
+
+      // if the key matches with something we expect, we assign the new value
+      if ( (key == "fAmpCut") || (key == "fReqFractionAboveAmpCutVal") ||
+          (key == "fAmpCutLEDRef") || (key == "fSecInAverage") || 
+          (key == "fDownscale") ) {
+       istringstream iss(value);
+       printf("AliCaloCalibSignal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
+
+       if (key == "fAmpCut") { 
+         iss >> fAmpCut; 
+       }
+       else if (key == "fReqFractionAboveAmpCutVal") { 
+         iss >> fReqFractionAboveAmpCutVal; 
+       }
+       else if (key == "fAmpCutLEDRef") { 
+         iss >> fAmpCutLEDRef; 
+       }
+       else if (key == "fSecInAverage") { 
+         iss >> fSecInAverage; 
        }
-       else{     
-         //DS
-//       sig->GetGraphAmpVsTimeHighGain(i,j,k);
-//       sig->GetGraphAmpVsTimeLowGain(i,j,k);
+       else if (key == "fDownscale") { 
+         iss >> fDownscale; 
        }
+      } // some match found/expected
+
+    }          
+  }
+
+  in.close();
+  return;      
+}
+
+//_____________________________________________________________________
+void AliCaloCalibSignal::WriteParametersToFile(const char *parameterFile)
+{ // write parameters to file
+  static const string delimitor("::");
+  ofstream out( parameterFile );
+  out << "// " << parameterFile << endl;
+  out << "fAmpCut" << "::" << fAmpCut << endl;
+  out << "fReqFractionAboveAmpCutVal" << "::" << fReqFractionAboveAmpCutVal << endl;
+  out << "fAmpCutLEDRef" << "::" << fAmpCutLEDRef << endl;
+  out << "fSecInAverage" << "::" << fSecInAverage << endl;
+  out << "fDownscale" << "::" << fDownscale << endl;
+
+  out.close();
+  return;
+}
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::AddInfo(const AliCaloCalibSignal *sig)
+{ 
+  // 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.
+
+  // add info from sig's TTrees to ours..
+  TTree *sigAmp = sig->GetTreeAmpVsTime();
+  TTree *sigAvgAmp = sig->GetTreeAvgAmpVsTime();
+
+  // we could try some merging via TList or what also as a more elegant approach
+  // but I wanted with the stupid/simple and hopefully safe approach of looping
+  // over what we want to add..
+
+  // associate variables for sigAmp and sigAvgAmp:
+  sigAmp->SetBranchAddress("fChannelNum",&fChannelNum);
+  sigAmp->SetBranchAddress("fHour",&fHour);
+  sigAmp->SetBranchAddress("fAmp",&fAmp);
+
+  // loop over the trees.. note that since we use the same variables we should not need
+  // to do any assignments between the getting and filling
+  for (int i=0; i<sigAmp->GetEntries(); i++) {
+    sigAmp->GetEntry(i);
+    fTreeAmpVsTime->Fill();
+  }
+
+  sigAvgAmp->SetBranchAddress("fChannelNum",&fChannelNum);
+  sigAvgAmp->SetBranchAddress("fHour",&fHour);
+  sigAvgAmp->SetBranchAddress("fAvgAmp",&fAvgAmp);
+  sigAvgAmp->SetBranchAddress("fRMS",&fRMS);
+
+  for (int i=0; i<sigAvgAmp->GetEntries(); i++) {
+    sigAvgAmp->GetEntry(i);
+    fTreeAvgAmpVsTime->Fill();
+  }
+
+  // also LED.. 
+  TTree *sigLEDAmp = sig->GetTreeLEDAmpVsTime();
+  TTree *sigLEDAvgAmp = sig->GetTreeLEDAvgAmpVsTime();
+
+  // associate variables for sigAmp and sigAvgAmp:
+  sigLEDAmp->SetBranchAddress("fRefNum",&fRefNum);
+  sigLEDAmp->SetBranchAddress("fHour",&fHour);
+  sigLEDAmp->SetBranchAddress("fAmp",&fAmp);
+
+  // loop over the trees.. note that since we use the same variables we should not need
+  // to do any assignments between the getting and filling
+  for (int i=0; i<sigLEDAmp->GetEntries(); i++) {
+    sigLEDAmp->GetEntry(i);
+    fTreeLEDAmpVsTime->Fill();
+  }
 
-      }//end for nModules 
-    }//end for nColumns
-  }//end for nRows
+  sigLEDAvgAmp->SetBranchAddress("fRefNum",&fRefNum);
+  sigLEDAvgAmp->SetBranchAddress("fHour",&fHour);
+  sigLEDAvgAmp->SetBranchAddress("fAvgAmp",&fAvgAmp);
+  sigLEDAvgAmp->SetBranchAddress("fRMS",&fRMS);
 
-  return kTRUE;//We succesfully added info from the supplied object
+  for (int i=0; i<sigLEDAvgAmp->GetEntries(); i++) {
+    sigLEDAvgAmp->GetEntry(i);
+    fTreeLEDAvgAmpVsTime->Fill();
+  }
+
+  // We should also copy other pieces of info: counters and parameters 
+  // (not number of columns and rows etc which should be the same)
+  // note that I just assign them here rather than Add them, but we
+  // normally just Add (e.g. in Preprocessor) one object so this should be fine.
+  fRunNumber = sig->GetRunNumber();
+  fStartTime = sig->GetStartTime();
+  fAmpCut = sig->GetAmpCut();
+  fReqFractionAboveAmpCutVal = sig->GetReqFractionAboveAmpCutVal();
+  fReqFractionAboveAmp = sig->GetReqFractionAboveAmp();
+  fAmpCutLEDRef = sig->GetAmpCutLEDRef();
+  fReqLEDRefAboveAmpCutVal = sig->GetReqLEDRefAboveAmpCutVal();
+  fHour = sig->GetHour();
+  fLatestHour = sig->GetLatestHour();
+  fUseAverage = sig->GetUseAverage();
+  fSecInAverage = sig->GetSecInAverage();
+  fDownscale = sig->GetDownscale();
+  fNEvents = sig->GetNEvents();
+  fNAcceptedEvents = sig->GetNAcceptedEvents();
+
+  return kTRUE;//We hopefully succesfully added info from the supplied object
 }
 
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::ProcessEvent(AliRawReader *rawReader)
+{
+  // if fMapping is NULL the rawstream will crate its own mapping
+  AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);  
+  if (fDetType == kEmCal) {
+    rawReader->Select("EMCAL", 0, AliEMCALGeoParams::fgkLastAltroDDL) ; //select EMCAL DDL range 
+  }
+  return ProcessEvent( &rawStream, rawReader->GetTimestamp() );
+}
 
 //_____________________________________________________________________
-Bool_t AliCaloCalibSignal::ProcessEvent(AliCaloRawStream *in, AliRawEventHeaderBase *aliHeader)
+Bool_t AliCaloCalibSignal::ProcessEvent(AliCaloRawStreamV3 *in, UInt_t Timestamp)
 { 
   // Method to process=analyze one event in the data stream
   if (!in) return kFALSE; //Return right away if there's a null pointer
   
   fNEvents++; // one more event
 
-  // PHOS has more towers than EMCAL, so use PHOS numbers to set array sizes
-  int AmpValHighGain[fgkMaxTowers];
-  int AmpValLowGain[fgkMaxTowers];
+  if ( (fNEvents%fDownscale)!=0 ) return kFALSE; // mechanism to skip some of the input events, if we want
+
+  // use maximum numbers to set array sizes
+  int iAmpValHighGain[fgkMaxTowers];
+  int iAmpValLowGain[fgkMaxTowers];
+  memset(iAmpValHighGain, 0, sizeof(iAmpValHighGain));
+  memset(iAmpValLowGain, 0, sizeof(iAmpValLowGain));
 
-  memset(AmpValHighGain, 0, sizeof(AmpValHighGain));
-  memset(AmpValLowGain, 0, sizeof(AmpValLowGain));
+  // also for LED reference
+  int iLEDAmpVal[fgkMaxRefs * 2]; // factor 2 is for the two gain values
+  memset(iLEDAmpVal, 0, sizeof(iLEDAmpVal));
 
-  int sample, i = 0; //The sample temp, and the sample number in current event.
-  int max = fgkSampleMin, min = fgkSampleMax;//Use these for picking the signal
-  int gain = 0;
+  int sample = 0; // temporary value
+  int gain = 0; // high or low gain
   
-  // Number of Low and High gain channels for this event:
+  // Number of Low and High gain, and LED Ref, channels for this event:
   int nLowChan = 0; 
   int nHighChan = 0; 
+  int nLEDRefChan = 0;
 
-  int TowerNum = 0; // array index for TGraphs etc.
+  int iTowerNum = 0; // array index for regular towers
+  int iRefNum = 0; // array index for LED references
 
   // loop first to get the fraction of channels with amplitudes above cut
-  while (in->Next()) {
-    sample = in->GetSignal(); //Get the adc signal
-    if (sample < min) min = sample;
-    if (sample > max) max = sample;
-    i++;
-    if ( i >= in->GetTimeLength()) {
+
+  while (in->NextDDL()) {
+    while (in->NextChannel()) {
+
+      // counters
+      int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
+      int nsamples = 0;
+
+      while (in->NextBunch()) {
+       const UShort_t *sig = in->GetSignals();
+       nsamples += in->GetBunchLength();
+       for (Int_t i = 0; i < in->GetBunchLength(); i++) {
+         sample = sig[i];
+
+         // check if it's a min or max value
+         if (sample < min) min = sample;
+         if (sample > max) max = sample;
+
+       } // loop over samples in bunch
+      } // loop over bunches
+
+      if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
+
+      gain = -1; // init to not valid value
       //If we're here then we're done with this tower
-      gain = 1 - in->IsLowGain();
-      
+      if ( in->IsLowGain() ) {
+       gain = 0;
+      }
+      else if ( in->IsHighGain() ) {
+       gain = 1;
+      }
+      else if ( in->IsLEDMonData() ) {
+       gain = in->GetRow(); // gain coded in (in RCU/Altro mapping) as Row info for LED refs..
+      }
+
+      // it should be enough to check the SuperModule info for each DDL really, but let's keep it here for now
       int arrayPos = in->GetModule(); //The modules are numbered starting from 0
-      if (arrayPos >= fModules) {
-       //TODO: return an error message, if appopriate (perhaps if debug>0?)
-       return kFALSE;
-      } 
-      
       //Debug
       if (arrayPos < 0 || arrayPos >= fModules) {
-       printf("Oh no: arrayPos = %i.\n", arrayPos); 
+       printf("AliCaloCalibSignal::ProcessEvent = Oh no: arrayPos = %i.\n", arrayPos); 
+       return kFALSE;
       }
 
-      // get tower number for AmpVal array
-      TowerNum = GetTowerNum(arrayPos, in->GetColumn(), in->GetRow()); 
-
-      if (gain == 0) {
-       // fill amplitude into the array           
-        AmpValLowGain[TowerNum] = max - min;
-       nLowChan++;
-      } 
-      else if (gain==1) {//fill the high gain ones
-       // fill amplitude into the array
-       AmpValHighGain[TowerNum] = max - min;
-       nHighChan++;
-      }//end if gain
-
-      
-      max = fgkSampleMin; min = fgkSampleMax;
-      i = 0;
-      
-    }//End if end of tower
+      if ( in->IsHighGain() || in->IsLowGain() ) { // regular tower
+       // get tower number for AmpVal array
+       iTowerNum = GetTowerNum(arrayPos, in->GetColumn(), in->GetRow()); 
+
+       if (gain == 0) {
+         // fill amplitude into the array         
+         iAmpValLowGain[iTowerNum] = max - min;
+         nLowChan++;
+       } 
+       else if (gain==1) {//fill the high gain ones
+         // fill amplitude into the array
+         iAmpValHighGain[iTowerNum] = max - min;
+         nHighChan++;
+       }//end if gain
+      } // regular tower
+      else if ( in->IsLEDMonData() ) { // LED ref.; 
+       // strip # is coded is 'column' in the channel maps 
+       iRefNum = GetRefNum(arrayPos, in->GetColumn(), gain); 
+       iLEDAmpVal[iRefNum] = max - min;
+       nLEDRefChan++;
+      } // end of LED ref
+
+      } // nsamples>0 check, some data found for this channel; not only trailer/header      
+    } // end while over channel 
    
-  }//end while, of stream
+  }//end while over DDL's, of input stream
   
-  // now check if it was a led event, only use high gain (that should be sufficient)
+  in->Reset(); // just in case the next customer forgets to check if the stream was reset..
+
+  // now check if it was an LED event, using the LED Reference info per strip
+
+  // by default all columns are accepted (init check to > 0)
+  int checkResultArray[AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALCols];
+  for (int ia=0; ia<(AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALCols); ia++) { 
+    checkResultArray[ia] = 1; 
+  }
   if (fReqFractionAboveAmp) {
-    bool ok = CheckFractionAboveAmp(AmpValHighGain, nHighChan);
+    bool ok = false;
+    if (nHighChan > 0) { 
+      ok = CheckFractionAboveAmp(iAmpValHighGain, checkResultArray); 
+    }
+    if (!ok) return false; // skip event
+  }
+
+  // by default all columns are accepted (init check to > 0)
+  int checkResultArrayLEDRef[AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALLEDRefs];
+  for (int ia=0; ia<(AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALLEDRefs); ia++) { 
+    checkResultArrayLEDRef[ia] = 1; 
+  }
+  if (fReqLEDRefAboveAmpCutVal) {
+    bool ok = false;
+    if (nLEDRefChan > 0) { 
+      ok = CheckLEDRefAboveAmp(iLEDAmpVal, checkResultArrayLEDRef); 
+    }
     if (!ok) return false; // skip event
   }
 
   fNAcceptedEvents++; // one more event accepted
 
   if (fStartTime == 0) { // if start-timestamp wasn't set,we'll pick it up from the first event we encounter
-    fStartTime = aliHeader->Get("Timestamp");
+    fStartTime = Timestamp;
   }
 
-  fHour = (aliHeader->Get("Timestamp")-fStartTime)/(double)fgkNumSecInHr;
+  fHour = (Timestamp - fStartTime)/(double)fgkNumSecInHr;
   if (fLatestHour < fHour) {
     fLatestHour = fHour; 
   }
   
-  // it is a led event, now fill graphs (maybe profiles later)
-  for(int i=0;i<fModules;i++){
-    for(int j=0;j<fColumns;j++){
-      for(int k=0;k<fRows;k++){
+  // it is a led event, now fill TTree
+  // We also do the activity check for LEDRefs/Strips, but need to translate between column
+  // and strip indices for that; based on these relations: 
+  // iStrip = AliEMCALGeoParams::GetStripModule(iSM, iCol);
+  // iStrip = (iSM%2==0) ? iCol/2 : AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - iCol/2;
+  // which leads to
+  // iColFirst = (iSM%2==0) ? iStrip*2 : (AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - iStrip)*2;
+
+  for(int i=0; i<fModules; i++){
+    for(int j=0; j<fColumns; j++) {
+      int iStrip = (i%2==0) ? j/2 : AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - j/2;
+      if (checkResultArray[i*fColumns + j]>0  && checkResultArrayLEDRef[i*fLEDRefs + iStrip]>0) { // column passed check 
+      for(int k=0; k<fRows; k++){
        
-       TowerNum = GetTowerNum(i, j, k); 
+       iTowerNum = GetTowerNum(i, j, k); 
 
-       if(AmpValHighGain[TowerNum]) {
-         fGraphAmpVsTimeHighGain[TowerNum]->SetPoint(fNHighGain[TowerNum]++,fHour,AmpValHighGain[TowerNum]);
+       if(iAmpValHighGain[iTowerNum]) {
+         fAmp = iAmpValHighGain[iTowerNum];
+         fChannelNum = GetChannelNum(i,j,k,1);
+         fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValHighGain[iTowerNum]);
+         fNHighGain[iTowerNum]++;
        }
-       if(AmpValLowGain[TowerNum]) {
-         fGraphAmpVsTimeLowGain[TowerNum]->SetPoint(fNLowGain[TowerNum]++,fHour,AmpValLowGain[TowerNum]);
+       if(iAmpValLowGain[iTowerNum]) {
+         fAmp = iAmpValLowGain[iTowerNum];
+         fChannelNum = GetChannelNum(i,j,k,0);
+         fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValLowGain[iTowerNum]);
+         fNLowGain[iTowerNum]++;
        }
-      }
+      } // rows
+      } // column passed check, and LED Ref for strip passed check (if any)
+    } // columns
+
+    // also LED refs
+    for(int j=0; j<fLEDRefs; j++){
+      int iColFirst = (i%2==0) ? j*2 : (AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - j)*2; //CHECKME!!!
+      if ( ((checkResultArray[i*fColumns + iColFirst]>0) || (checkResultArray[i*fColumns + iColFirst + 1]>0)) && // at least one column in strip passed check 
+          (checkResultArrayLEDRef[i*fLEDRefs + j]>0) ) { // and LED Ref passed checks
+       for (gain=0; gain<2; gain++) {
+         fRefNum = GetRefNum(i, j, gain); 
+         if (iLEDAmpVal[fRefNum]) {
+           fAmp = iLEDAmpVal[fRefNum];
+           fTreeLEDAmpVsTime->Fill();//fRefNum,fHour,fAmp);
+           fNRef[fRefNum]++;
+         }
+       } // gain
+      } // at least one column in strip passed check, and LED Ref passed check (if any) 
     }
-  }
+
+  } // modules
   
   return kTRUE;
 }
 
 //_____________________________________________________________________
-void AliCaloCalibSignal::CreateProfile(int imod, int ic, int ir, int towerId, int gain,
-                                      int nbins, double min, double max)
-{ //! create/setup a TProfile
-  TString title, name;   
-  if (gain == 0) { 
-    name = "fProfAmpVsTimeLowGain_";   
-    title = "Amp vs. Time Low Gain Mod "; 
-  } 
-  else if (gain == 1) { 
-    name = "fProfAmpVsTimeHighGain_"; 
-    title = "Amp vs. Time High Gain Mod "; 
-  } 
-  name += imod;
-  name += "_"; name += ic;
-  name += "_"; name += ir;
-  title += imod;
-  title += " Col "; title += ic;
-  title += " Row "; title += ir;
-           
-  // use "s" option for RMS
-  if (gain==0) { 
-    fProfAmpVsTimeLowGain[towerId] = new TProfile(name,title, nbins, min, max,"s");
-  }
-  else if (gain==1) {
-    fProfAmpVsTimeHighGain[towerId] = new TProfile(name,title, nbins, min, max,"s");
-  }
-
-  return;
-}
-//_____________________________________________________________________
-Bool_t AliCaloCalibSignal::Save(TString fileName, Bool_t saveEmptyGraphs)
+Bool_t AliCaloCalibSignal::Save(TString fileName)
 {
-  //Saves all the histograms (or profiles, to be accurate) to the designated file
+  //Saves all the TTrees to the designated file
   
   TFile destFile(fileName, "recreate");
   
@@ -420,101 +719,187 @@ Bool_t AliCaloCalibSignal::Save(TString fileName, Bool_t saveEmptyGraphs)
   
   destFile.cd();
 
-  // setup variables for the TProfile plot
+  // save the trees
+  fTreeAmpVsTime->Write();
+  fTreeLEDAmpVsTime->Write();
+  if (fUseAverage) { 
+    Analyze(); // get the latest and greatest averages
+    fTreeAvgAmpVsTime->Write();
+    fTreeLEDAvgAmpVsTime->Write();
+  }
+
+  destFile.Close();
+  
+  return kTRUE;
+}
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::Analyze()
+{
+  // Fill the tree holding the average values
+  if (!fUseAverage) { return kFALSE; }
+
+  // Reset the average TTree if Analyze has already been called earlier,
+  // meaning that the TTree could have been partially filled
+  if (fTreeAvgAmpVsTime->GetEntries() > 0) {
+    fTreeAvgAmpVsTime->Reset();
+  }
+
+  //0: setup variables for the TProfile plots that we'll use to do the averages
   int numProfBins = 0;
   double timeMin = 0;
   double timeMax = 0;
-  if (fUseAverage) {
-    if (fSecInAverage > 0) {
-      numProfBins = (int)( (fLatestHour*fgkNumSecInHr)/fSecInAverage + 1 ); // round-off
-    }
-    numProfBins += 2; // add extra buffer : first and last
-    double binSize = 1.0*fSecInAverage / fgkNumSecInHr;
-    timeMin = - binSize;
-    timeMax = timeMin + numProfBins*binSize;
+  if (fSecInAverage > 0) {
+    numProfBins = (int)( (fLatestHour*fgkNumSecInHr)/fSecInAverage + 1 ); // round-off
   }
+  numProfBins += 2; // add extra buffer : first and last
+  double binSize = 1.0*fSecInAverage / fgkNumSecInHr;
+  timeMin = - binSize;
+  timeMax = timeMin + numProfBins*binSize;
+
+  //1: set up TProfiles for the towers that had data
+  TProfile * profile[fgkMaxTowers*2]; // *2 is since we include both high and low gains
+  memset(profile, 0, sizeof(profile));
+  const Int_t buffersize = 200;
+  char name[buffersize]; // for profile id and title
+  int iTowerNum = 0;
 
-  int numGraphPoints= 0;
-  int TowerNum = 0;    
-  for (int i = 0; i < fModules; i++) {
-    
-    for(int ic=0;ic < fColumns;ic++){
-      for(int ir=0;ir < fRows;ir++){
-
-       TowerNum = GetTowerNum(i, ic, ir); 
-
-       // 1st: high gain
-       numGraphPoints= fGraphAmpVsTimeHighGain[TowerNum]->GetN();
-       if( numGraphPoints>0 || saveEmptyGraphs) {
-         
-         // average the graphs points over time if requested and put them in a profile plot
-         if(fUseAverage && numGraphPoints>0) {
-           
-           // get the values
-           double *graphX = fGraphAmpVsTimeHighGain[TowerNum]->GetX();
-           double *graphY = fGraphAmpVsTimeHighGain[TowerNum]->GetY();
-
-           // create the TProfile: 1 is for High gain              
-           CreateProfile(i, ic, ir, TowerNum, 1,
-                         numProfBins, timeMin, timeMax);
-
-           // loop over graph points and fill profile
-           for(int ip=0; ip < numGraphPoints; ip++){
-             fProfAmpVsTimeHighGain[TowerNum]->Fill(graphX[ip],graphY[ip]);
-           }
-           
-           fProfAmpVsTimeHighGain[TowerNum]->GetXaxis()->SetTitle("Hours");
-           fProfAmpVsTimeHighGain[TowerNum]->GetYaxis()->SetTitle("MaxAmplitude - Pedestal");
-           fProfAmpVsTimeHighGain[TowerNum]->Write();
+  for (int i = 0; i<fModules; i++) {
+    for (int ic=0; ic<fColumns; ic++){
+      for (int ir=0; ir<fRows; ir++) {
+
+       iTowerNum = GetTowerNum(i, ic, ir);
+       // high gain
+       if (fNHighGain[iTowerNum] > 0) {
+         fChannelNum = GetChannelNum(i, ic, ir, 1); 
+         snprintf(name,buffersize,"profileChan%d", fChannelNum);
+         profile[fChannelNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
+       }
 
-         }
-          else{
-            //otherwise, just save the graphs and forget the profiling
-            fGraphAmpVsTimeHighGain[TowerNum]->GetXaxis()->SetTitle("Hours");
-            fGraphAmpVsTimeHighGain[TowerNum]->GetYaxis()->SetTitle("MaxAmplitude - Pedestal");
-            fGraphAmpVsTimeHighGain[TowerNum]->Write();
-          }
-         
-       } // low gain graph info should be saved in one form or another
-       
-       // 2nd: now go to the low gain case
-       numGraphPoints= fGraphAmpVsTimeLowGain[TowerNum]->GetN();
-       if( numGraphPoints>0 || saveEmptyGraphs) {
-         
-         // average the graphs points over time if requested and put them in a profile plot
-         if(fUseAverage && numGraphPoints>0) {
-           
-           double *graphX = fGraphAmpVsTimeLowGain[TowerNum]->GetX();
-           double *graphY = fGraphAmpVsTimeLowGain[TowerNum]->GetY();
-           
-           // create the TProfile: 0 is for Low gain       
-           CreateProfile(i, ic, ir, TowerNum, 0,
-                         numProfBins, timeMin, timeMax);
-
-           // loop over graph points and fill profile
-           for(int ip=0; ip < numGraphPoints; ip++){
-             fProfAmpVsTimeLowGain[TowerNum]->Fill(graphX[ip],graphY[ip]);
-           }
-           
-           fProfAmpVsTimeLowGain[TowerNum]->GetXaxis()->SetTitle("Hours");
-           fProfAmpVsTimeLowGain[TowerNum]->GetYaxis()->SetTitle("MaxAmplitude - Pedestal");
-           fProfAmpVsTimeLowGain[TowerNum]->Write();
+       // same for low gain
+       if (fNLowGain[iTowerNum] > 0) {
+         fChannelNum = GetChannelNum(i, ic, ir, 0); 
+         snprintf(name,buffersize,"profileChan%d", fChannelNum);
+         profile[fChannelNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
+       }
 
-         }
-         else{
-            //otherwise, just save the graphs and forget the profiling
-           fGraphAmpVsTimeLowGain[TowerNum]->GetXaxis()->SetTitle("Hours");
-           fGraphAmpVsTimeLowGain[TowerNum]->GetYaxis()->SetTitle("MaxAmplitude - Pedestal");
-           fGraphAmpVsTimeLowGain[TowerNum]->Write();
-         }
-         
-       } // low gain graph info should be saved in one form or another
+      } // rows
+    } // columns
+  } // modules
+
+  //2: fill profiles by looping over tree
+  // Set addresses for tree-readback also
+  fTreeAmpVsTime->SetBranchAddress("fChannelNum", &fChannelNum);
+  fTreeAmpVsTime->SetBranchAddress("fHour", &fHour);
+  fTreeAmpVsTime->SetBranchAddress("fAmp", &fAmp);
+
+  for (int ient=0; ient<fTreeAmpVsTime->GetEntries(); ient++) {
+    fTreeAmpVsTime->GetEntry(ient);
+    if (profile[fChannelNum]) { 
+      // profile should always have been created above, for active channels
+      profile[fChannelNum]->Fill(fHour, fAmp);
+    }
+  }
 
-      } // fRows
-    } // fColumns
+  // re-associating the branch addresses here seems to be needed for OK 'average' storage          
+  fTreeAvgAmpVsTime->SetBranchAddress("fChannelNum", &fChannelNum);
+  fTreeAvgAmpVsTime->SetBranchAddress("fHour", &fHour);
+  fTreeAvgAmpVsTime->SetBranchAddress("fAvgAmp", &fAvgAmp);
+  fTreeAvgAmpVsTime->SetBranchAddress("fRMS", &fRMS);
+
+  //3: fill avg tree by looping over the profiles
+  for (fChannelNum = 0; fChannelNum<(fgkMaxTowers*2); fChannelNum++) {
+    if (profile[fChannelNum]) { // profile was created
+      if (profile[fChannelNum]->GetEntries() > 0) { // profile had some entries
+       for(int it=0; it<numProfBins; it++) {
+         if (profile[fChannelNum]->GetBinEntries(it+1) > 0) {
+           fAvgAmp = profile[fChannelNum]->GetBinContent(it+1);
+           fHour = profile[fChannelNum]->GetBinCenter(it+1);
+           fRMS = profile[fChannelNum]->GetBinError(it+1);
+           fTreeAvgAmpVsTime->Fill();
+         } // some entries for this bin
+       } // loop over bins
+      } // some entries for this profile
+    } // profile exists  
+  } // loop over all possible channels
+
+
+  // and finally, go through same exercise for LED also.. 
+
+  //1: set up TProfiles for the towers that had data
+  TProfile * profileLED[fgkMaxRefs*2]; // *2 is since we include both high and low gains
+  memset(profileLED, 0, sizeof(profileLED));
+
+  for (int i = 0; i<fModules; i++) {
+    for(int j=0; j<fLEDRefs; j++){
+      for (int gain=0; gain<2; gain++) {
+       fRefNum = GetRefNum(i, j, gain);
+       if (fNRef[fRefNum] > 0) { 
+         snprintf(name, buffersize, "profileLEDRef%d", fRefNum);
+         profileLED[fRefNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
+       } 
+      }// gain
+    } 
+  } // modules
+
+  //2: fill profiles by looping over tree
+  // Set addresses for tree-readback also
+  fTreeLEDAmpVsTime->SetBranchAddress("fRefNum", &fRefNum);
+  fTreeLEDAmpVsTime->SetBranchAddress("fHour", &fHour);
+  fTreeLEDAmpVsTime->SetBranchAddress("fAmp", &fAmp);
+
+  for (int ient=0; ient<fTreeLEDAmpVsTime->GetEntries(); ient++) {
+    fTreeLEDAmpVsTime->GetEntry(ient);
+    if (profileLED[fRefNum]) { 
+      // profile should always have been created above, for active channels
+      profileLED[fRefNum]->Fill(fHour, fAmp);
+    }
+  }
+
+  // re-associating the branch addresses here seems to be needed for OK 'average' storage          
+  fTreeLEDAvgAmpVsTime->SetBranchAddress("fRefNum", &fRefNum);
+  fTreeLEDAvgAmpVsTime->SetBranchAddress("fHour", &fHour);
+  fTreeLEDAvgAmpVsTime->SetBranchAddress("fAvgAmp", &fAvgAmp);
+  fTreeLEDAvgAmpVsTime->SetBranchAddress("fRMS", &fRMS);
+
+  //3: fill avg tree by looping over the profiles
+  for (fRefNum = 0; fRefNum<(fgkMaxRefs*2); fRefNum++) {
+    if (profileLED[fRefNum]) { // profile was created
+      if (profileLED[fRefNum]->GetEntries() > 0) { // profile had some entries
+       for(int it=0; it<numProfBins; it++) {
+         if (profileLED[fRefNum]->GetBinEntries(it+1) > 0) {
+           fAvgAmp = profileLED[fRefNum]->GetBinContent(it+1);
+           fHour = profileLED[fRefNum]->GetBinCenter(it+1);
+           fRMS = profileLED[fRefNum]->GetBinError(it+1);
+           fTreeLEDAvgAmpVsTime->Fill();
+         } // some entries for this bin
+       } // loop over bins
+      } // some entries for this profile
+    } // profile exists  
+  } // loop over all possible channels
+
+  // OK, we're done..
 
-  } // fModules
-  destFile.Close();
-  
   return kTRUE;
 }
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::DecodeChannelNum(const int chanId, 
+                                           int *imod, int *icol, int *irow, int *igain) const  
+{ // return the module, column, row, and gain for a given channel number
+  *igain = chanId/(fModules*fColumns*fRows);
+  *imod = (chanId/(fColumns*fRows)) % fModules;
+  *icol = (chanId/fRows) % fColumns;
+  *irow = chanId % fRows;
+  return kTRUE;
+} 
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::DecodeRefNum(const int refId, 
+                                       int *imod, int *istripMod, int *igain) const 
+{ // return the module, stripModule, and gain for a given reference number
+  *igain = refId/(fModules*fLEDRefs);
+  *imod = (refId/(fLEDRefs)) % fModules;
+  *istripMod = refId % fLEDRefs;
+  return kTRUE;
+}