]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOSGetter.cxx
Restoring raw data fit from version of 29-Aug-2004
[u/mrichter/AliRoot.git] / PHOS / AliPHOSGetter.cxx
index 6e9bfb620c12cfa3c5a2de2400ab3b7f80a3ab15..719daf3c1cfb724c0d487106bc80fe5cfb613c2c 100644 (file)
 #include <TROOT.h>
 #include <TSystem.h>
 #include <TParticle.h>
-
+#include <TF1.h>
+#include <TGraph.h>
+//#include <TCanvas.h>
+//#include <TFrame.h>
 
 // --- Standard library ---
 
@@ -58,6 +61,9 @@
 #include "AliPHOSLoader.h"
 #include "AliRunLoader.h"
 #include "AliStack.h"  
+#include "AliPHOSRawStream.h"
+#include "AliRawReaderFile.h"
+#include "AliLog.h"
 
 ClassImp(AliPHOSGetter)
   
@@ -82,6 +88,7 @@ AliPHOSGetter::AliPHOSGetter(const char* headerFile, const char* version, Option
     if (rl->GetAliRun() == 0x0) {
       rl->LoadgAlice();
       gAlice = rl->GetAliRun(); // should be removed
+      rl->LoadHeader();
     }
   }
   fgPhosLoader = dynamic_cast<AliPHOSLoader*>(rl->GetLoader("PHOSLoader"));
@@ -90,7 +97,6 @@ AliPHOSGetter::AliPHOSGetter(const char* headerFile, const char* version, Option
   else 
     fgPhosLoader->SetTitle(version);
   
-  
   // initialize data members
   SetDebug(0) ; 
   fBTE = 0 ; 
@@ -102,23 +108,31 @@ AliPHOSGetter::AliPHOSGetter(const char* headerFile, const char* version, Option
   fESDFile = 0 ; 
   fESD = 0 ; 
   fESDTree = 0 ; 
+  fRawDigits = kFALSE ;
 }
 
 //____________________________________________________________________________ 
 AliPHOSGetter::~AliPHOSGetter()
 {
   // dtor
-  delete fgPhosLoader ;
-  fgPhosLoader = 0 ;
-  delete fBTE ; 
-  fBTE = 0 ; 
-  fPrimaries->Delete() ; 
-  delete fPrimaries ; 
-  fgObjGetter = 0;
+  if(fgPhosLoader){
+    delete fgPhosLoader ;
+    fgPhosLoader = 0 ;
+  }
+  if(fBTE){
+    delete fBTE ; 
+    fBTE = 0 ;
+  } 
+  if(fPrimaries){
+    fPrimaries->Delete() ; 
+    delete fPrimaries ;
+  } 
   if (fESD) 
     delete fESD ; 
   if (fESDTree) 
-    delete fESDTree ; 
+    delete fESDTree ;
+  fgObjGetter = 0;
 }
 
 //____________________________________________________________________________ 
@@ -254,59 +268,75 @@ void AliPHOSGetter::Event(Int_t event, const char* opt)
 {
   // Reads the content of all Tree's S, D and R
 
-  if ( event >= MaxEvent() ) {
-    Error("Event", "%d not found in TreeE !", event) ; 
-    return ; 
-  }
+//   if ( event >= MaxEvent() ) {
+//     Error("Event", "%d not found in TreeE !", event) ; 
+//     return ; 
+//   }
 
   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 
-  // checks if we are dealing with test-beam data
-  TBranch * btb = rl->TreeE()->GetBranch("AliPHOSBeamTestEvent") ;
-  if(btb){
-    if(!fBTE)
-      fBTE = new AliPHOSBeamTestEvent() ;
-    btb->SetAddress(&fBTE) ;
-    btb->GetEntry(event) ;
-  }
-  else{
-    if(fBTE){
-      delete fBTE ;
-      fBTE = 0 ;
-    }
-  }
+//   // checks if we are dealing with test-beam data
+//   TBranch * btb = rl->TreeE()->GetBranch("AliPHOSBeamTestEvent") ;
+//   if(btb){
+//     if(!fBTE)
+//       fBTE = new AliPHOSBeamTestEvent() ;
+//     btb->SetAddress(&fBTE) ;
+//     btb->GetEntry(event) ;
+//   }
+//   else{
+//     if(fBTE){
+//       delete fBTE ;
+//       fBTE = 0 ;
+//     }
+//   }
 
   // Loads the type of object(s) requested
   
   rl->GetEvent(event) ;
 
-  if( strstr(opt,"X") || (strcmp(opt,"")==0) )
+  if(strstr(opt,"X") || (strcmp(opt,"")==0)){
     ReadPrimaries() ;
-
-  if(strstr(opt,"H") )
+  }
+  
+  if(strstr(opt,"H")  ){
     ReadTreeH();
-
-  if(strstr(opt,"S") )
+  }
+  
+  if(strstr(opt,"S")  ){
     ReadTreeS() ;
-
-  if( strstr(opt,"D") )
+  }
+  
+  if(strstr(opt,"D") ){
     ReadTreeD() ;
-
-  if( strstr(opt,"R") )
+  }
+  
+  if(strstr(opt,"R") ){
     ReadTreeR() ;
+  }
 
-  if( strstr(opt,"T") )
+  if( strstr(opt,"T") ){
     ReadTreeT() ;
+  }    
 
-  if( strstr(opt,"P") )
+  if( strstr(opt,"P") ){
     ReadTreeP() ;
+  }    
 
-  if( strstr(opt,"E") )
+  if( strstr(opt,"E") ){
     ReadTreeE(event) ;
+  }
+
+}
 
-//   if( strstr(opt,"Q") )
-//     ReadTreeQA() ;
+
+//____________________________________________________________________________ 
+void AliPHOSGetter::Event(AliRawReader *rawReader, const char* opt) 
+{
+  // Reads the raw event from rawReader
+
+  if( strstr(opt,"W")  ){
+    ReadRaw(rawReader) ;
+  }    
  
 }
 
@@ -340,8 +370,6 @@ AliPHOSGetter * AliPHOSGetter::Instance(const char* alirunFileName, const char*
   // Creates and returns the pointer of the unique instance
   // Must be called only when the environment has changed
   
-  //::Info("Instance","alirunFileName=%s version=%s openingOption=%s",alirunFileName,version,openingOption);
-  
   if(!fgObjGetter){ // first time the getter is called 
     fgObjGetter = new AliPHOSGetter(alirunFileName, version, openingOption) ;
   }
@@ -413,7 +441,16 @@ TParticle * AliPHOSGetter::Primary(Int_t index) const
 AliPHOS * AliPHOSGetter:: PHOS() const  
 {
   // returns the PHOS object 
-  AliPHOS * phos = dynamic_cast<AliPHOS*>(PhosLoader()->GetModulesFolder()->FindObject("PHOS")) ;  
+  AliPHOSLoader *    loader = 0;
+  static AliPHOSLoader * oldloader = 0;
+  static AliPHOS * phos = 0;
+
+  loader = PhosLoader();
+
+  if ( loader != oldloader) {
+    phos = dynamic_cast<AliPHOS*>(loader->GetModulesFolder()->FindObject("PHOS")) ;
+    oldloader = loader;
+  }
   if (!phos) 
     if (fgDebug)
       Warning("PHOS", "PHOS module not found in module folders: %s", PhosLoader()->GetModulesFolder()->GetName() ) ; 
@@ -513,18 +550,189 @@ Bool_t AliPHOSGetter::OpenESDFile()
   return rv ; 
 }
 
+//____________________________________________________________________________ 
+void AliPHOSGetter::FitRaw(Bool_t lowGainFlag, TGraph * gLowGain, TGraph * gHighGain, TF1* signalF, Int_t & amp, Double_t & time)
+{
+  // Fits the raw signal time distribution 
+
+  const Int_t kNoiseThreshold = 0 ;
+  Double_t timezero1 = 0., timezero2 = 0., timemax = 0. ;
+  Double_t signal = 0., signalmax = 0. ;       
+  Double_t energy = time = 0. ; 
+
+  if (lowGainFlag) {
+    timezero1 = timezero2 = signalmax = timemax = 0. ;
+    signalF->FixParameter(0, PHOS()->GetRawFormatLowCharge()) ; 
+    signalF->FixParameter(1, PHOS()->GetRawFormatLowGain()) ; 
+    Int_t index ; 
+    for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
+      gLowGain->GetPoint(index, time, signal) ; 
+      if (signal > kNoiseThreshold && timezero1 == 0.) 
+       timezero1 = time ;
+      if (signal <= kNoiseThreshold && timezero1 > 0. && timezero2 == 0.)
+       timezero2 = time ; 
+      if (signal > signalmax) {
+       signalmax = signal ; 
+       timemax   = time ; 
+      }
+    }
+    signalmax /= PHOS()->RawResponseFunctionMax(PHOS()->GetRawFormatLowCharge(), 
+                                               PHOS()->GetRawFormatLowGain()) ;
+    if ( timezero1 + PHOS()->GetRawFormatTimePeak() < PHOS()->GetRawFormatTimeMax() * 0.4 ) { // else its noise 
+      signalF->SetParameter(2, signalmax) ; 
+      signalF->SetParameter(3, timezero1) ;                
+      gLowGain->Fit(signalF, "QRON", "", 0., timezero2); //, "QRON") ; 
+      energy = signalF->GetParameter(2) ; 
+      time   = signalF->GetMaximumX() - PHOS()->GetRawFormatTimePeak() - PHOS()->GetRawFormatTimeTrigger() ;
+    }
+  } else {
+    timezero1 = timezero2 = signalmax = timemax = 0. ;
+    signalF->FixParameter(0, PHOS()->GetRawFormatHighCharge()) ; 
+    signalF->FixParameter(1, PHOS()->GetRawFormatHighGain()) ; 
+    Int_t index ; 
+    for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
+      gHighGain->GetPoint(index, time, signal) ;               
+      if (signal > kNoiseThreshold && timezero1 == 0.) 
+       timezero1 = time ;
+      if (signal <= kNoiseThreshold && timezero1 > 0. && timezero2 == 0.)
+       timezero2 = time ; 
+      if (signal > signalmax) {
+       signalmax = signal ;   
+       timemax   = time ; 
+      }
+    }
+    signalmax /= PHOS()->RawResponseFunctionMax(PHOS()->GetRawFormatHighCharge(), 
+                                               PHOS()->GetRawFormatHighGain()) ;;
+    if ( timezero1 + PHOS()->GetRawFormatTimePeak() < PHOS()->GetRawFormatTimeMax() * 0.4 ) { // else its noise  
+      signalF->SetParameter(2, signalmax) ; 
+      signalF->SetParameter(3, timezero1) ;               
+      gHighGain->Fit(signalF, "QRON", "", 0., timezero2) ; 
+      energy = signalF->GetParameter(2) ; 
+      time   = signalF->GetMaximumX() - PHOS()->GetRawFormatTimePeak() - PHOS()->GetRawFormatTimeTrigger() ;
+    }
+  }
+  
+  if (time == 0. && energy == 0.) 
+    amp = 0 ; 
+  else {
+  AliPHOSDigitizer * digitizer = Digitizer() ; 
+  amp = static_cast<Int_t>( (energy - digitizer->GetEMCpedestal()) / digitizer->GetEMCchannel() + 0.5 ) ; 
+  }
+}
+
+//____________________________________________________________________________ 
+Int_t AliPHOSGetter::ReadRaw(AliRawReader *rawReader)
+{
+  // reads the raw format data, converts it into digits format and store digits in Digits()
+  // container.
+
+  AliPHOSRawStream in(rawReader);
+  Bool_t first = kTRUE ;
+  
+  TF1 * signalF = new TF1("signal", AliPHOS::RawResponseFunction, 0, PHOS()->GetRawFormatTimeMax(), 4);
+  signalF->SetParNames("Charge", "Gain", "Amplitude", "TimeZero") ; 
+
+  Int_t relId[4], id ;
+  Bool_t lowGainFlag = kFALSE ; 
+  TClonesArray * digits = Digits() ;
+  digits->Clear() ; 
+  Int_t idigit = 0 ; 
+  Int_t amp = 0 ; 
+  Double_t time = 0. ; 
+
+  TGraph * gLowGain = new TGraph(PHOS()->GetRawFormatTimeBins()) ; 
+  TGraph * gHighGain= new TGraph(PHOS()->GetRawFormatTimeBins()) ;  
+
+  while ( in.Next() ) { // PHOS entries loop 
+    if ( (in.IsNewRow() || in.IsNewColumn() || in.IsNewModule()) ) {
+      if (!first) {
+       FitRaw(lowGainFlag, gLowGain, gHighGain, signalF, amp, time) ; 
+       if (amp > 0) {
+         new((*digits)[idigit]) AliPHOSDigit( -1, id, amp, time) ;     
+         idigit++ ; 
+       }
+       Int_t index ; 
+       for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
+         gLowGain->SetPoint(index, index * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 0) ;  
+         gHighGain->SetPoint(index, index * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 0) ; 
+       } 
+      }
+      first = kFALSE ; 
+      relId[0] = in.GetModule() ;
+      if ( relId[0] >= PHOS()->GetRawFormatLowGainOffset() ) { 
+       relId[0] -=  PHOS()->GetRawFormatLowGainOffset() ; 
+       lowGainFlag = kTRUE ; 
+      } else 
+       lowGainFlag = kFALSE ; 
+      relId[1] = 0 ; 
+      relId[2] = in.GetRow() ; 
+      relId[3] = in.GetColumn() ; 
+      PHOSGeometry()->RelToAbsNumbering(relId, id) ;   
+    }
+    if (lowGainFlag)
+      gLowGain->SetPoint(in.GetTime(), 
+                    in.GetTime()* PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 
+                    in.GetSignal()) ;
+    else 
+      gHighGain->SetPoint(in.GetTime(), 
+                        in.GetTime() * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 
+                        in.GetSignal() ) ;
+
+  } // PHOS entries loop
+
+  FitRaw(lowGainFlag, gLowGain, gHighGain, signalF, amp, time) ; 
+  if (amp > 0 ) {
+    new((*digits)[idigit]) AliPHOSDigit( -1, id, amp, time) ;
+    idigit++ ; 
+  }
+  digits->Sort() ; 
+
+  delete signalF ; 
+  delete gLowGain, gHighGain ; 
+  
+  return digits->GetEntriesFast() ; 
+}
+
+//   TClonesArray * digits = Digits() ;
+//   digits->Clear() ; 
+//   Int_t idigit = 0 ; 
+
+//   while ( in.Next() ) { // PHOS entries loop 
+//     Int_t amp = in.GetSignal() ; 
+//     Double_t time = in.GetTime() ; 
+//     Int_t relId[4], id ;
+
+//     relId[0] = in.GetModule() ;
+//     if ( relId[0] >= PHOS()->GetRawFormatLowGainOffset() ) { 
+//       relId[0] -=  PHOS()->GetRawFormatLowGainOffset() ; 
+//     }
+//     relId[1] = 0 ; 
+//     relId[2] = in.GetRow() ; 
+//     relId[3] = in.GetColumn() ; 
+//     PHOSGeometry()->RelToAbsNumbering(relId, id) ;  
+    
+//     if (amp > 0 && id >=0 ) {
+//       new((*digits)[idigit]) AliPHOSDigit( -1, id, amp, time) ;     
+//       idigit++ ; 
+//     }
+    
+//   } // PHOS entries loop
+  
+//   digits->Sort() ; 
+  
+//   return digits->GetEntriesFast() ; 
+// }
 //____________________________________________________________________________ 
 Int_t AliPHOSGetter::ReadTreeD()
 {
   // Read the Digits
   
   PhosLoader()->CleanDigits() ;    
-  // gets TreeD from the root file (PHOS.Digits.root)
-  // if ( !IsLoaded("D") ) {
-    PhosLoader()->LoadDigits("UPDATE") ;
-    PhosLoader()->LoadDigitizer("UPDATE") ;
-    //  SetLoaded("D") ; 
-    //} 
+  PhosLoader()->LoadDigits("UPDATE") ;
+  PhosLoader()->LoadDigitizer("UPDATE") ;
   return Digits()->GetEntries() ; 
 }