Added method to read raw data and convert to digits format

diff --git a/EMCAL/AliEMCALGetter.cxx b/EMCAL/AliEMCALGetter.cxx
index 402a535715b1c45ae9fd36cb404d123dc3457408..bcfcae36208b8db3f016b04f6e7a780c7b23ea6a 100644 (file)
--- a/EMCAL/AliEMCALGetter.cxx
+++ b/EMCAL/AliEMCALGetter.cxx
@@ -43,7 +43,8 @@
 #include <TFile.h>
 #include <TROOT.h>
 #include <TSystem.h>
-
+#include <TH1D.h>
+#include <TF1.h>
 
 // --- Standard library ---
 
@@ -56,6 +57,8 @@
 #include "AliMC.h"
 #include "AliRunLoader.h"
 #include "AliStack.h"  
+#include "AliEMCALRawStream.h"
+#include "AliRawReaderFile.h"
 
 ClassImp(AliEMCALGetter)
   
@@ -272,7 +275,10 @@ void AliEMCALGetter::Event(Int_t event, const char* opt)
 
   if( strstr(opt,"P") )
     ReadTreeP() ;
- 
+
+  if( strstr(opt,"W") )
+    ReadRaw(event) ;
+  
 }
 
 
@@ -442,6 +448,85 @@ void AliEMCALGetter::ReadPrimaries()
 }
 
 //____________________________________________________________________________ 
+Int_t AliEMCALGetter::ReadRaw(Int_t event)
+{
+  // reads the raw format data, converts it into digits format and store digits in Digits()
+  // container.
+  
+  AliRawReaderFile rawReader(event) ; 
+  AliEMCALRawStream in(&rawReader);
+  
+  const Int_t kHighGainIdOffset = EMCALGeometry()->GetNTowers()
+    * EMCALGeometry()->GetNPhi()
+    * EMCALGeometry()->GetNZ() 
+    * 2 ;
+  
+  Bool_t first = kTRUE ;
+  
+  TH1D hLowGain("hLowGain", "Low Gain", 1000, 0., EMCAL()->GetRawFormatTimeMax()) ; 
+  TH1D hHighGain("hHighGain", "High Gain", 1000, 0., EMCAL()->GetRawFormatTimeMax()) ; 
+  
+  // fit half the gaussian decay rather than AliEMCAL::RawResponseFunction because thiswould give a floating point
+  // exception during error calculation ... To solve... 
+  TF1 * gauss = new TF1("gauss", "gaus", 
+                       EMCAL()->GetRawFormatTimePeak(), 
+                       EMCAL()->GetRawFormatTimeMax() ) ;   
+  
+  Int_t id = -1;
+  Bool_t hgflag = kFALSE ; 
+  
+  TClonesArray * digits = Digits() ;
+  digits->Clear() ; 
+  Int_t idigit = 0 ; 
+  
+  while ( in.Next() ) { // EMCAL entries loop 
+    
+    if ( in.IsNewId() ) {
+      if (!first) {
+       hLowGain.Fit(gauss, "QRON") ; 
+       Int_t ampL =  static_cast<Int_t>(gauss->Eval(gauss->GetParameter(2)) + 0.5) ; 
+       Double_t timeL  = EMCAL()->GetRawFormatTimePeak() - gauss->GetParameter(2) ;
+       if (timeL < 0 ) // happens with noise 
+         timeL =  EMCAL()->GetRawFormatTimePeak() ;  
+       if (hgflag) {
+         hHighGain.Fit(gauss, "QRON") ; 
+         Int_t ampH =  static_cast<Int_t>(gauss->Eval(gauss->GetParameter(2)) + 0.5) ; 
+         Double_t timeH  = EMCAL()->GetRawFormatTimePeak() - gauss->GetParameter(2) ; 
+         if (timeH < 0 ) // happens with noise 
+           timeH =  EMCAL()->GetRawFormatTimePeak() ;  
+         new((*digits)[idigit]) AliEMCALDigit( -1, -1, id+kHighGainIdOffset, ampH, timeH) ;
+         idigit++ ; 
+       }
+       else {
+         new((*digits)[idigit]) AliEMCALDigit( -1, -1, id, ampL, timeL) ;
+         idigit++ ; 
+       }
+      }
+      first = kFALSE ; 
+      hLowGain.Reset() ; 
+      hHighGain.Reset() ; 
+      id = in.GetId() ; 
+      if (id > 9999 ) { // fixme 
+       hgflag = kTRUE ; 
+      } else 
+       hgflag = kFALSE ; 
+    }
+    if (hgflag)
+      hHighGain.Fill(
+                    in.GetTime() * EMCAL()->GetRawFormatTimeMax() / EMCAL()->GetRawFormatTimeBins(), 
+                    in. GetSignal() * EMCAL()->GetRawFormatHighGainFactor() ) ;
+    else 
+      hLowGain.Fill(
+                   in.GetTime() * EMCAL()->GetRawFormatTimeMax() / EMCAL()->GetRawFormatTimeBins(), 
+                   in. GetSignal() ) ;
+  } // EMCAL entries loop
+  
+  delete gauss ; 
+  
+  return Digits()->GetEntriesFast() ; 
+}
+  
+  //____________________________________________________________________________ 
 Int_t AliEMCALGetter::ReadTreeD()
 {
   // Read the Digits

diff --git a/EMCAL/AliEMCALGetter.h b/EMCAL/AliEMCALGetter.h
index fdf6a5944530e6eea68bec61ce6e3a0cef0a2911..dddd0ae5eebe6c1b5698401f5bd47a3a36bfe9bc 100644 (file)
--- a/EMCAL/AliEMCALGetter.h
+++ b/EMCAL/AliEMCALGetter.h
@@ -172,6 +172,8 @@ class AliEMCALGetter : public TObject {
   Int_t                 WritePID(Option_t* opt="") const {
     return  EmcalLoader()->WritePID(opt) ; }
 
+  //========== Raw ===========
+  Int_t ReadRaw(Int_t event) ; 
 
   void SetDebug(Int_t level) {fgDebug = level;} // Set debug level 
     void PostClusterizer(AliEMCALClusterizer * clu) 

diff --git a/PHOS/AliPHOSGetter.cxx b/PHOS/AliPHOSGetter.cxx
index 6e9bfb620c12cfa3c5a2de2400ab3b7f80a3ab15..c8f9b6365e7ed30b47daac1c193af3eda36e47c0 100644 (file)
--- a/PHOS/AliPHOSGetter.cxx
+++ b/PHOS/AliPHOSGetter.cxx
@@ -44,7 +44,8 @@
 #include <TROOT.h>
 #include <TSystem.h>
 #include <TParticle.h>
-
+#include <TH1D.h>
+#include <TF1.h>
 
 // --- Standard library ---
 
@@ -58,6 +59,8 @@
 #include "AliPHOSLoader.h"
 #include "AliRunLoader.h"
 #include "AliStack.h"  
+#include "AliPHOSRawStream.h"
+#include "AliRawReaderFile.h"
 
 ClassImp(AliPHOSGetter)
   
@@ -303,6 +306,9 @@ void AliPHOSGetter::Event(Int_t event, const char* opt)
 
   if( strstr(opt,"E") )
     ReadTreeE(event) ;
+
+  if( strstr(opt,"W") )
+    ReadRaw(event) ;
  
 
 //   if( strstr(opt,"Q") )
@@ -513,6 +519,90 @@ Bool_t AliPHOSGetter::OpenESDFile()
   return rv ; 
 }
 
+//____________________________________________________________________________ 
+Int_t AliPHOSGetter::ReadRaw(Int_t event)
+{
+  // reads the raw format data, converts it into digits format and store digits in Digits()
+  // container.
+
+  AliRawReaderFile rawReader(event) ; 
+  AliPHOSRawStream in(&rawReader);
+  
+  const Int_t kHighGainIdOffset = PHOSGeometry()->GetNModules()
+    * PHOSGeometry()->GetNPhi()
+    * PHOSGeometry()->GetNZ() 
+    * 2 ;
+
+  Bool_t first = kTRUE ;
+  
+  TH1D hLowGain("hLowGain", "Low Gain", 1000, 0., PHOS()->GetRawFormatTimeMax()) ; 
+  TH1D hHighGain("hHighGain", "High Gain", 1000, 0., PHOS()->GetRawFormatTimeMax()) ; 
+  
+  // fit half the gaussian decay rather than AliPHOS::RawResponseFunction because thiswould give a floating point
+  // exception during error calculation ... To solve... 
+  TF1 * gauss = new TF1("gauss", "gaus", 
+                        PHOS()->GetRawFormatTimePeak(), 
+                        PHOS()->GetRawFormatTimeMax() ) ;   
+
+  Int_t relId[4], id ;
+  Bool_t hgflag = kFALSE ; 
+ 
+  TClonesArray * digits = Digits() ;
+  digits->Clear() ; 
+  Int_t idigit = 0 ; 
+  
+  while ( in.Next() ) { // PHOS entries loop 
+        
+    if ( (in.IsNewRow() || in.IsNewColumn() || in.IsNewModule()) ) {
+      if (!first) {
+       hLowGain.Fit(gauss, "QRON") ; 
+       Int_t ampL =  static_cast<Int_t>(gauss->Eval(gauss->GetParameter(2)) + 0.5) ; 
+       Double_t timeL  = PHOS()->GetRawFormatTimePeak() - gauss->GetParameter(2) ;
+       if (timeL < 0 ) // happens with noise 
+         timeL =  PHOS()->GetRawFormatTimePeak() ;  
+       if (hgflag) {
+         hHighGain.Fit(gauss, "QRON") ; 
+         Int_t ampH =  static_cast<Int_t>(gauss->Eval(gauss->GetParameter(2)) + 0.5) ; 
+         Double_t timeH  = PHOS()->GetRawFormatTimePeak() - gauss->GetParameter(2) ; 
+         if (timeH < 0 ) // happens with noise 
+           timeH =  PHOS()->GetRawFormatTimePeak() ;  
+         new((*digits)[idigit]) AliPHOSDigit( -1, id+kHighGainIdOffset, ampH, timeH) ;
+         idigit++ ; 
+       }
+       else {
+         new((*digits)[idigit]) AliPHOSDigit( -1, id, ampL, timeL) ;
+         idigit++ ; 
+       }
+      }
+      first = kFALSE ; 
+      hLowGain.Reset() ; 
+      hHighGain.Reset() ; 
+      relId[0] = in.GetModule() ;
+      if ( relId[0] >= PHOS()->GetRawFormatHighGainOffset() ) { 
+       relId[0] -=  PHOS()->GetRawFormatHighGainOffset() ; 
+       hgflag = kTRUE ; 
+      } else 
+         hgflag = kFALSE ; 
+      relId[1] = 0 ; 
+      relId[2] = in.GetRow() ; 
+      relId[3] = in.GetColumn() ; 
+      PHOSGeometry()->RelToAbsNumbering(relId, id) ;   
+    }
+    if (hgflag)
+      hHighGain.Fill(
+                    in.GetTime() * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 
+                    in. GetSignal() * PHOS()->GetRawFormatHighGainFactor() ) ;
+    else 
+      hLowGain.Fill(
+                   in.GetTime() * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 
+                   in. GetSignal() ) ;
+  } // PHOS entries loop
+  
+  delete gauss ; 
+  
+  return Digits()->GetEntriesFast() ; 
+}
+
 //____________________________________________________________________________ 
 Int_t AliPHOSGetter::ReadTreeD()
 {

diff --git a/PHOS/AliPHOSGetter.h b/PHOS/AliPHOSGetter.h
index 07240ec9022eb23e0ef35f460dffbe54290c4924..420db469029c78b3e77b623a030d418bd1e77aac 100644 (file)
--- a/PHOS/AliPHOSGetter.h
+++ b/PHOS/AliPHOSGetter.h
@@ -167,7 +167,10 @@ public:
   Int_t                 WriteRecParticles(Option_t* opt="") const { return PhosLoader()->WriteRecParticles(opt) ; }
   Int_t                 WritePID(Option_t* opt="") const {
     return  PhosLoader()->WritePID(opt) ; }
-  
+
+  //========== Raw ===========
+  Int_t ReadRaw(Int_t event) ; 
+
   void SetDebug(Int_t level) {fgDebug = level;} // Set debug level 
   void PostClusterizer(AliPHOSClusterizer * clu) 
     const{PhosLoader()->PostClusterizer(clu) ; }