]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOSDigitizer.cxx
Extending to beam test analysis
[u/mrichter/AliRoot.git] / PHOS / AliPHOSDigitizer.cxx
index 4da8bd1df594a379277e5adc445ecdda90f0be68..43b9a6342aed118c45c471f1aca9e88fd9b11c71 100644 (file)
@@ -68,9 +68,9 @@
 #include <iomanip.h>
 
 // --- AliRoot header files ---
-
 #include "AliRun.h"
 #include "AliHeader.h"
+#include "AliStream.h"
 #include "AliRunDigitizer.h"
 #include "AliPHOSDigit.h"
 #include "AliPHOS.h"
@@ -87,56 +87,51 @@ ClassImp(AliPHOSDigitizer)
   AliPHOSDigitizer::AliPHOSDigitizer() 
 {
   // ctor
+  InitParameters() ; 
+  fDefaultInit = kTRUE ;
+  fManager = 0 ;                     // We work in the standalong mode
 
-  fPinNoise           = 0.004 ;
-  fEMCDigitThreshold  = 0.012 ;
-  fCPVNoise           = 0.01;
-  fCPVDigitThreshold  = 0.09 ;
-  fTimeResolution     = 0.5e-9 ;
-  fTimeSignalLength   = 1.0e-9 ;
-  fDigitsInRun  = 0 ; 
-  fADCchanelEmc = 0.0015;        // width of one ADC channel in GeV
-  fADCpedestalEmc = 0.005 ;      //
-  fNADCemc = (Int_t) TMath::Power(2,16) ;  // number of channels in EMC ADC
-
-  fADCchanelCpv = 0.0012 ;          // width of one ADC channel in CPV 'popugais'
-  fADCpedestalCpv = 0.012 ;         // 
-  fNADCcpv = (Int_t) TMath::Power(2,12);      // number of channels in CPV ADC
-
-  fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
-  fManager = 0 ;                        // We work in the standalong mode
-  fSplitFile= 0 ; 
-}
+ }
 
 //____________________________________________________________________________ 
-AliPHOSDigitizer::AliPHOSDigitizer(const char *headerFile,const char * name)
+AliPHOSDigitizer::AliPHOSDigitizer(const char *headerFile, const char * name, const Bool_t toSplit)
 {
   // ctor
-  SetName(name) ;
+
   SetTitle(headerFile) ;
+  SetName(name) ;
   fManager = 0 ;                     // We work in the standalong mode
   fSplitFile= 0 ; 
+  InitParameters() ; 
+  fToSplit = toSplit ;
   Init() ;
-  
+  fDefaultInit = kFALSE ; 
 }
 
 //____________________________________________________________________________ 
 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * ard):AliDigitizer(ard)
 {
   // ctor
-  SetName("");     //Will call init in the digitizing
-  SetTitle("aliroot") ;  
+  SetTitle(ard->GetInputFileName(0,0)) ;
+  InitParameters() ; 
+  fDefaultInit = kFALSE ; 
+  fSplitFile   = 0 ; 
+
+  if (ard->GetOutputFile()) {
+    SetName(ard->GetOutputFile().Data());
+    fToSplit = kTRUE ;
+  } else {
+    SetName("Default") ;
+    fToSplit = kFALSE ;
+  }
 }
 
 //____________________________________________________________________________ 
   AliPHOSDigitizer::~AliPHOSDigitizer()
 {
   // dtor
-
- if (fSplitFile)       
-   if ( fSplitFile->IsOpen() )         
-     fSplitFile->Close() ;
-
+     
+   fSplitFile = 0 ; 
 }
 
 //____________________________________________________________________________
@@ -156,7 +151,6 @@ void AliPHOSDigitizer::Digitize(const Int_t event)
   digits->Clear() ;
 
   const AliPHOSGeometry *geom = gime->PHOSGeometry() ; 
-
   //Making digits with noise, first EMC
   Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
   
@@ -183,15 +177,16 @@ void AliPHOSDigitizer::Digitize(const Int_t event)
   TClonesArray * sdigits = 0 ;
   Int_t input = 0 ;
   TObjArray * sdigArray = new TObjArray(2) ;
-  while ( (folder = (TFolder*)next()) ) 
+  while ( (folder = (TFolder*)next()) ) 
     if ( (sdigits = (TClonesArray*)folder->FindObject(GetName()) ) ) {
       TString fileName(folder->GetName()) ;
       fileName.ReplaceAll("_","/") ;
-      cout << "INFO: AliPHOSDigitizer::Digitize -> Adding SDigits " 
-          << GetName() << " from " << fileName << endl ; 
+//       cout << "INFO: AliPHOSDigitizer::Digitize -> Adding SDigits " 
+//        << GetName() << " from " << fileName << endl ; 
       sdigArray->AddAt(sdigits, input) ;
       input++ ;
     }
+  }
 
   //Find the first crystall with signal
   Int_t nextSig = 200000 ; 
@@ -336,7 +331,7 @@ void AliPHOSDigitizer::Digitize(const Int_t event)
   
   //remove digits below thresholds
   for(i = 0; i < nEMC ; i++){
-    digit = (AliPHOSDigit*) digits->At(i) ;
+    digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
     if(sDigitizer->Calibrate( digit->GetAmp() ) < fEMCDigitThreshold)
       digits->RemoveAt(i) ;
     else
@@ -355,13 +350,14 @@ void AliPHOSDigitizer::Digitize(const Int_t event)
 
   //Set indexes in list of digits and make true digitization of the energy
   for (i = 0 ; i < ndigits ; i++) { 
-    AliPHOSDigit * digit = (AliPHOSDigit *) digits->At(i) ; 
+    digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ; 
     digit->SetIndexInList(i) ;     
     Float_t energy = sDigitizer->Calibrate(digit->GetAmp()) ;
     digit->SetAmp(DigitizeEnergy(energy,digit->GetId()) ) ;
   }
 
 }
+
 //____________________________________________________________________________
 Int_t AliPHOSDigitizer::DigitizeEnergy(Float_t energy, Int_t absId)
 {
@@ -376,6 +372,7 @@ Int_t AliPHOSDigitizer::DigitizeEnergy(Float_t energy, Int_t absId)
   }
   return chanel ;
 }
+
 //____________________________________________________________________________
 void AliPHOSDigitizer::Exec(Option_t *option) 
 { 
@@ -383,7 +380,6 @@ void AliPHOSDigitizer::Exec(Option_t *option)
 
   if(strcmp(GetName(), "") == 0 )   
     Init() ;
-  
   if (strstr(option,"print")) {
     Print("");
     return ; 
@@ -401,57 +397,58 @@ void AliPHOSDigitizer::Exec(Option_t *option)
   if(fManager){
     treeD = fManager->GetTreeD() ;
     nevents = 1 ;    // Will process only one event
-  }
-  else {
-    gAlice->GetEvent(0) ;
-    nevents = (Int_t) gAlice->TreeE()->GetEntries() ;
-    treeD=gAlice->TreeD() ;
-  }
-
-
-  //Check, if this branch already exits
-  if (treeD) { 
-    TObjArray * lob = (TObjArray*)treeD->GetListOfBranches() ;
-    TIter next(lob) ; 
-    TBranch * branch = 0 ;  
-    Bool_t phosfound = kFALSE, digitizerfound = kFALSE ; 
     
-    while ( (branch = (TBranch*)next()) && (!phosfound || !digitizerfound) ) {
-      if ( (strcmp(branch->GetName(), "PHOS")==0) && 
-          (strcmp(branch->GetTitle(), GetName())==0) ) 
-       phosfound = kTRUE ;
+    //Check, if this branch already exits
+    if (treeD) { 
+      TObjArray * lob = (TObjArray*)treeD->GetListOfBranches() ;
+      TIter next(lob) ; 
+      TBranch * branch = 0 ;  
+      Bool_t phosfound = kFALSE, digitizerfound = kFALSE ; 
       
-      else if ( (strcmp(branch->GetName(), "AliPHOSDigitizer")==0) && 
-               (strcmp(branch->GetTitle(), GetName())==0) ) 
-       digitizerfound = kTRUE ; 
-    }
-    
-    if ( phosfound ) {
-      cerr << "WARNING: AliPHOSDigitizer -> Digits branch with name " << GetName() 
-          << " already exits" << endl ;
-      return ; 
-    }   
-    if ( digitizerfound ) {
-      cerr << "WARNING: AliPHOSDigitizer -> Digitizer branch with name " << GetName() 
-          << " already exits" << endl ;
-      return ; 
+      while ( (branch = (TBranch*)next()) && (!phosfound || !digitizerfound) ) {
+       if ( (strcmp(branch->GetName(), "PHOS")==0) && 
+            (strcmp(branch->GetTitle(), GetName())==0) ) 
+         phosfound = kTRUE ;
+       
+       else if ( (strcmp(branch->GetName(), "AliPHOSDigitizer")==0) && 
+                 (strcmp(branch->GetTitle(), GetName())==0) ) 
+         digitizerfound = kTRUE ; 
+      }
+      
+      if ( phosfound ) {
+       cerr << "WARNING: AliPHOSDigitizer -> Digits branch with name " << GetName() 
+            << " already exits" << endl ;
+       return ; 
+      }   
+      if ( digitizerfound ) {
+       cerr << "WARNING: AliPHOSDigitizer -> Digitizer branch with name " << GetName() 
+            << " already exits" << endl ;
+       return ; 
+      }
     }   
   }
-
+  else { //PHOS standalone
+    if(gime->BranchExists("Digits") ) 
+      return ;
+    nevents=gime->MaxEvent() ;
+  }
+  
   Int_t ievent ;
   
   for(ievent = 0; ievent < nevents; ievent++){
   
     if(fManager){
+
       Int_t input ;
       for(input = 0 ; input < fManager->GetNinputs(); input ++){
        TTree * treeS = fManager->GetInputTreeS(input) ;
        if(!treeS){
-         cerr << "AliPHOSDigitizer -> No Input " << endl ;
+         cerr << "AliPHOSDigitizer::Exec -> No Input " << endl ;
          return ;
        }
        gime->ReadTreeS(treeS,input) ;
       }
+
     }
     else
       gime->Event(ievent,"S") ;
@@ -466,11 +463,7 @@ void AliPHOSDigitizer::Exec(Option_t *option)
     //increment the total number of Digits per run 
     fDigitsInRun += gime->Digits()->GetEntriesFast() ;  
   }
-  
-  if (fSplitFile) 
-    if ( fSplitFile->IsOpen() ) 
-      fSplitFile->Close() ; 
+   
   if(strstr(option,"tim")){
     gBenchmark->Stop("PHOSDigitizer");
     cout << "AliPHOSDigitizer:" << endl ;
@@ -500,41 +493,23 @@ Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks)
   }
   return time ;
 }
+
 //____________________________________________________________________________ 
 Bool_t AliPHOSDigitizer::Init()
 {
-  fPinNoise           = 0.004 ;
-  fEMCDigitThreshold  = 0.012 ;
-  fCPVNoise           = 0.01;
-  fCPVDigitThreshold  = 0.09 ;
-  fTimeResolution     = 0.5e-9 ;
-  fTimeSignalLength   = 1.0e-9 ;
-  fDigitsInRun  = 0 ; 
-  fADCchanelEmc = 0.0015;        // width of one ADC channel in GeV
-  fADCpedestalEmc = 0.005 ;      //
-  fNADCemc = (Int_t) TMath::Power(2,16) ;  // number of channels in EMC ADC
+  // Makes all memory allocations
 
-  fADCchanelCpv = 0.0012 ;          // width of one ADC channel in CPV 'popugais'
-  fADCpedestalCpv = 0.012 ;         // 
-  fNADCcpv = (Int_t) TMath::Power(2,12);      // number of channels in CPV ADC
-
-  fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
-    
-  if(fManager)
-    SetTitle("aliroot") ;
-  else if (strcmp(GetTitle(),"")==0) 
-   SetTitle("galice.root") ;
+  if( strcmp(GetTitle(), "") == 0 )
+    SetTitle("galice.root") ;
 
-  if( strcmp(GetName(), "") == 0 )
-    SetName("Default") ;
-  
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance(GetTitle(), GetName()) ; 
+  AliPHOSGetter * gime = AliPHOSGetter::GetInstance(GetTitle(), GetName(), fToSplit) ; 
   if ( gime == 0 ) {
     cerr << "ERROR: AliPHOSDigitizer::Init -> Could not obtain the Getter object !" << endl ; 
     return kFALSE;
   } 
   
   const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
+
   fEmcCrystals = geom->GetNModules() *  geom->GetNCristalsInModule() ;
   
   // Post Digits to the white board
@@ -542,6 +517,29 @@ Bool_t AliPHOSDigitizer::Init()
   
   // Post Digitizer to the white board
   gime->PostDigitizer(this) ;
+
+  fSplitFile = 0 ;
+  if(fToSplit){
+    // construct the name of the file as /path/PHOS.SDigits.root
+    //First - extract full path if necessary
+    TString digitsFileName(GetTitle()) ;
+    Ssiz_t islash = digitsFileName.Last('/') ;
+    if(islash<digitsFileName.Length())
+      digitsFileName.Remove(islash+1,digitsFileName.Length()) ;
+    else
+      digitsFileName="" ;
+    // Next - append the file name 
+    digitsFileName+="PHOS.Digits." ;
+    if((strcmp(GetName(),"Default")!=0)&&(strcmp(GetName(),"")!=0)){
+      digitsFileName+=GetName() ;
+      digitsFileName+="." ;
+    }
+    digitsFileName+="root" ;
+    // Finally - check if the file already opened or open the file
+    fSplitFile = static_cast<TFile*>(gROOT->GetFile(digitsFileName.Data()));   
+    if(!fSplitFile)
+      fSplitFile =  TFile::Open(digitsFileName.Data(),"update") ;
+  }
   
   //Mark that we will use current header file
   if(!fManager){
@@ -551,6 +549,28 @@ Bool_t AliPHOSDigitizer::Init()
   return kTRUE ;
 }
 
+//____________________________________________________________________________ 
+void AliPHOSDigitizer::InitParameters()
+{
+  fPinNoise           = 0.004 ;
+  fEMCDigitThreshold  = 0.012 ;
+  fCPVNoise           = 0.01;
+  fCPVDigitThreshold  = 0.09 ;
+  fTimeResolution     = 0.5e-9 ;
+  fTimeSignalLength   = 1.0e-9 ;
+  fDigitsInRun  = 0 ; 
+  fADCchanelEmc = 0.0015;        // width of one ADC channel in GeV
+  fADCpedestalEmc = 0.005 ;      //
+  fNADCemc = (Int_t) TMath::Power(2,16) ;  // number of channels in EMC ADC
+
+  fADCchanelCpv = 0.0012 ;          // width of one ADC channel in CPV 'popugais'
+  fADCpedestalCpv = 0.012 ;         // 
+  fNADCcpv = (Int_t) TMath::Power(2,12);      // number of channels in CPV ADC
+
+  fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
+    
+}
+
 //__________________________________________________________________
 void AliPHOSDigitizer::MixWith(const char* headerFile)
 {
@@ -577,7 +597,10 @@ void AliPHOSDigitizer::MixWith(const char* headerFile)
   // check if the specified SDigits do not already exist on the White Board:
   // //Folders/RunMC/Event/Data/PHOS/SDigits/headerFile/sdigitsname
 
-  TString path = "Folders/RunMC/Event/Data/PHOS/SDigits" ; 
+  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
+  TString path = gime->SDigitsFolder()->GetName() ; 
+
+  // before it was ???? "Folders/RunMC/Event/Data/PHOS/SDigits" ; 
   path += headerFile ; 
   path += "/" ; 
   path += GetName() ;
@@ -586,7 +609,6 @@ void AliPHOSDigitizer::MixWith(const char* headerFile)
     return;
   }
 
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;
   gime->PostSDigits(GetName(),headerFile) ;
   
   // check if the requested file is already open or exist and if SDigits Branch exist
@@ -601,77 +623,37 @@ void AliPHOSDigitizer::MixWith(const char* headerFile)
   
 }
 
-//__________________________________________________________________
-void AliPHOSDigitizer::SetSplitFile(const TString splitFileName)
-{
-  // Diverts the Digits in a file separate from the hits file
-  
-  // I guess it is not going to work if we do merging
-  if (fManager) {
-    cerr << "ERROR: AliPHOSDigitizer::SetSplitFile -> Not yet available in case of merging activated " << endl ;  
-    return ; 
-  }
-
-  TDirectory * cwd = gDirectory ;
-  if ( !(gAlice->GetTreeDFileName() == splitFileName) ) {
-    if (gAlice->GetTreeDFile() )  
-      gAlice->GetTreeDFile()->Close() ; 
-  }
-
-  fSplitFile = gAlice->InitTreeFile("D",splitFileName.Data());
-  fSplitFile->cd() ; 
-  if ( !fSplitFile->Get("gAlice") ) 
-    gAlice->Write();
-  
-  TTree *treeE  = gAlice->TreeE();
-  if (!treeE) {
-    cerr << "ERROR: AliPHOSDigitizer::SetSplitFile -> No TreeE found "<<endl;
-    abort() ;
-  }      
-  
-  // copy TreeE
-  if ( !fSplitFile->Get("TreeE") ) {
-    AliHeader *header = new AliHeader();
-    treeE->SetBranchAddress("Header", &header);
-    treeE->SetBranchStatus("*",1);
-    TTree *treeENew =  treeE->CloneTree();
-    treeENew->Write();
-  }
-
-  // copy AliceGeom
-  if ( !fSplitFile->Get("AliceGeom") ) {
-    TGeometry *AliceGeom = static_cast<TGeometry*>(cwd->Get("AliceGeom"));
-    if (!AliceGeom) {
-      cerr << "ERROR: AliPHOSDigitizer::SetSplitFile -> AliceGeom was not found in the input file "<<endl;
-      abort() ;
-    }
-    AliceGeom->Write();
-  }
-  
-  gAlice->MakeTree("D",fSplitFile);
-  cwd->cd() ; 
-  cout << "INFO: AliPHOSDigitizer::SetSPlitMode -> Digits will be stored in " << splitFileName.Data() << endl ; 
-}
-
 //__________________________________________________________________
 void AliPHOSDigitizer::Print(Option_t* option)const {
   // Print Digitizer's parameters
   if( strcmp(GetName(), "") != 0 ){
-    
+
+
     cout << "------------------- "<< GetName() << " -------------" << endl ;
-    cout << "Digitizing sDigits from file(s): " <<endl ;
-    
-     TCollection * folderslist = ((TFolder*)gROOT->FindObjectAny("Folders/RunMC/Event/Data/PHOS/SDigits"))->GetListOfFolders() ; 
-    TIter next(folderslist) ; 
-    TFolder * folder = 0 ; 
-    
-    while ( (folder = (TFolder*)next()) ) {
-      if ( folder->FindObject(GetName())  ) 
-       cout << "Adding SDigits " << GetName() << " from " << folder->GetName() << endl ; 
-    }
-    cout << endl ;
-    cout << "Writing digits to " << GetTitle() << endl ;
-    
+
+    const Int_t nStreams = GetNInputStreams() ; 
+    if (nStreams) {
+      Int_t index = 0 ;  
+      for (index = 0 ; index < nStreams ; index++)  
+       cout << "Adding SDigits " << GetName() << " from " <<  fManager->GetInputFileName(index, 0) << endl ; 
+      
+      cout << endl ;
+      cout << "Writing digits to " <<   fManager->GetInputFileName(0, 0) << endl ;   
+    } else { 
+//       AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;  
+//       gime->Folder("sdigits")  ;
+//       cout << "Digitizing sDigits from file(s): " <<endl ;
+//       TCollection * folderslist = gime->Folder("sdigits")->GetListOfFolders() ; 
+//       TIter next(folderslist) ; 
+//       TFolder * folder = 0 ; 
+      
+//       while ( (folder = (TFolder*)next()) ) {
+//     if ( folder->FindObject(GetName())  ) 
+      //    cout << "Adding SDigits " << GetName() << " from " << GetSDigitsFileName() << endl ; 
+//      }
+      cout << endl ;
+      cout << "Writing digits to " << GetTitle() << endl ;
+    }    
     cout << endl ;
     cout << "With following parameters: " << endl ;
     cout << "     Electronics noise in EMC (fPinNoise) = " << fPinNoise << endl ;
@@ -744,35 +726,12 @@ void AliPHOSDigitizer::PrintDigits(Option_t * option){
 
 }
 
-//__________________________________________________________________
-void AliPHOSDigitizer::SetSDigitsBranch(const char* title)
-{
-  // we set title (comment) of the SDigits branch in the first! header file
-  if( strcmp(GetName(), "") == 0 )
-    Init() ;
-
-  AliPHOSGetter::GetInstance()->SDigits()->SetName(title) ; 
-}
 //__________________________________________________________________
 Float_t AliPHOSDigitizer::TimeOfNoise(void)
 {  // Calculates the time signal generated by noise
   //to be rewritten, now returns just big number
   return 1. ;
 
-}
-//____________________________________________________________________________
-void AliPHOSDigitizer::Reset() 
-{ 
-  // sets current event number to the first simulated event
-
-  if( strcmp(GetName(), "") == 0 )
-    Init() ;
-
- //  Int_t inputs ;
-//   for(inputs = 0; inputs < fNinputs ;inputs++)
-//       fIevent->AddAt(-1, inputs ) ;
-  
 }
 
 //____________________________________________________________________________
@@ -789,31 +748,30 @@ void AliPHOSDigitizer::WriteDigits(Int_t event)
 
   AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
   const TClonesArray * digits = gime->Digits(GetName()) ; 
- TTree * treeD ;
-
- if(fManager) 
-   treeD = fManager->GetTreeD() ;
- else {
-   if (!gAlice->TreeD() ) 
-     gAlice->MakeTree("D", fSplitFile);
-   treeD = gAlice->TreeD();
- }
-
-
 TTree * treeD ;
+  
+  if(fManager) 
+    treeD = fManager->GetTreeD() ;
+  else {
+    if((gAlice->TreeD() == 0) || (fSplitFile))  // we should not create TreeD if it is already here
+      gAlice->MakeTree("D", fSplitFile);       // We overwrite TreeD in split file in the case of second reconstruction
+    if(fSplitFile)
+      fSplitFile->cd() ;
+    treeD = gAlice->TreeD();
+  }
+  
   // -- create Digits branch
   Int_t bufferSize = 32000 ;    
   TBranch * digitsBranch = treeD->Branch("PHOS",&digits,bufferSize);
   digitsBranch->SetTitle(GetName());
-    
+  
   // -- Create Digitizer branch
   Int_t splitlevel = 0 ;
-  AliPHOSDigitizer * d = gime->Digitizer(GetName()) ;
+  const AliPHOSDigitizer * d = dynamic_cast<const AliPHOSDigitizer *>(gime->Digitizer()) ;
   TBranch * digitizerBranch = treeD->Branch("AliPHOSDigitizer", "AliPHOSDigitizer", &d,bufferSize,splitlevel); 
   digitizerBranch->SetTitle(GetName());
 
   digitsBranch->Fill() ;
   digitizerBranch->Fill() ; 
   treeD->AutoSave() ; 
 }