]> 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 c4c097afc727af81b641d751a8667dc4a14cdce4..43b9a6342aed118c45c471f1aca9e88fd9b11c71 100644 (file)
 #include "TROOT.h"
 #include "TFolder.h"
 #include "TObjString.h"
+#include "TGeometry.h"
 #include "TBenchmark.h"
 
 // --- Standard library ---
 #include <iomanip.h>
 
 // --- AliRoot header files ---
-
 #include "AliRun.h"
+#include "AliHeader.h"
+#include "AliStream.h"
 #include "AliRunDigitizer.h"
 #include "AliPHOSDigit.h"
 #include "AliPHOS.h"
@@ -85,51 +87,51 @@ ClassImp(AliPHOSDigitizer)
   AliPHOSDigitizer::AliPHOSDigitizer() 
 {
   // ctor
+  InitParameters() ; 
+  fDefaultInit = kTRUE ;
+  fManager = 0 ;                     // We work in the standalong mode
 
-  fPinNoise           = 0.01 ;
-  fEMCDigitThreshold  = 0.01 ;
-  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
-}
+ }
 
 //____________________________________________________________________________ 
-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
-
-
+     
+   fSplitFile = 0 ; 
 }
 
 //____________________________________________________________________________
@@ -149,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();
   
@@ -176,13 +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()) ) ) {
-      cout << "INFO: AliPHOSDigitizer::Digitize -> Adding SDigits " 
-          << GetName() << " from " << folder->GetName() << endl ; 
+      TString fileName(folder->GetName()) ;
+      fileName.ReplaceAll("_","/") ;
+//       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 ; 
@@ -209,10 +213,10 @@ void AliPHOSDigitizer::Digitize(const Int_t event)
     Float_t noise = gRandom->Gaus(0., fPinNoise) ; 
     new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
     //look if we have to add signal?
+    digit = (AliPHOSDigit *) digits->At(absID-1) ;
     if(absID==nextSig){
       //Add SDigits from all inputs 
-      digit = (AliPHOSDigit *) digits->At(absID-1) ;
-
       ticks->Clear() ;
       Int_t contrib = 0 ;
       Float_t a = digit->GetAmp() ;
@@ -235,8 +239,8 @@ void AliPHOSDigitizer::Digitize(const Int_t event)
          if(fManager)
            primaryoffset = fManager->GetMask(i) ; 
          else
-           primaryoffset = i ;
-         curSDigit->ShiftPrimary(i) ;
+           primaryoffset = 10000000*i ;
+         curSDigit->ShiftPrimary(primaryoffset) ;
          
          a = curSDigit->GetAmp() ;
          b = a /fTimeSignalLength ;
@@ -327,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
@@ -346,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)
 {
@@ -367,6 +372,7 @@ Int_t AliPHOSDigitizer::DigitizeEnergy(Float_t energy, Int_t absId)
   }
   return chanel ;
 }
+
 //____________________________________________________________________________
 void AliPHOSDigitizer::Exec(Option_t *option) 
 { 
@@ -374,7 +380,6 @@ void AliPHOSDigitizer::Exec(Option_t *option)
 
   if(strcmp(GetName(), "") == 0 )   
     Init() ;
-  
   if (strstr(option,"print")) {
     Print("");
     return ; 
@@ -392,58 +397,58 @@ void AliPHOSDigitizer::Exec(Option_t *option)
   if(fManager){
     treeD = fManager->GetTreeD() ;
     nevents = 1 ;    // Will process only one event
+    
+    //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 ;
+       
+       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 {
-    gAlice->GetEvent(0) ;
-    nevents = (Int_t) gAlice->TreeE()->GetEntries() ;
-    treeD=gAlice->TreeD() ;
-  }
-  if(treeD == 0 ){
-    cerr << " AliPHOSDigitizer :: Can not find TreeD " << endl ;
-    return ;
+  else { //PHOS standalone
+    if(gime->BranchExists("Digits") ) 
+      return ;
+    nevents=gime->MaxEvent() ;
   }
-
-  //Check, if this branch already exits
-  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 ;
-    
-    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 ; 
-  }   
-
   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") ;
@@ -458,7 +463,7 @@ void AliPHOSDigitizer::Exec(Option_t *option)
     //increment the total number of Digits per run 
     fDigitsInRun += gime->Digits()->GetEntriesFast() ;  
   }
-  
+   
   if(strstr(option,"tim")){
     gBenchmark->Stop("PHOSDigitizer");
     cout << "AliPHOSDigitizer:" << endl ;
@@ -488,41 +493,23 @@ Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks)
   }
   return time ;
 }
+
 //____________________________________________________________________________ 
 Bool_t AliPHOSDigitizer::Init()
 {
-  fPinNoise           = 0.01 ;
-  fEMCDigitThreshold  = 0.01 ;
-  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
-  //  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
@@ -530,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){
@@ -539,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)
 {
@@ -565,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() ;
@@ -574,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
@@ -593,21 +627,33 @@ void AliPHOSDigitizer::MixWith(const char* headerFile)
 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 ;
@@ -680,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 ) ;
-  
 }
 
 //____________________________________________________________________________
@@ -725,54 +748,30 @@ void AliPHOSDigitizer::WriteDigits(Int_t event)
 
   AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
   const TClonesArray * digits = gime->Digits(GetName()) ; 
- TTree * treeD ;
-
-  if(fManager)
 TTree * treeD ;
+  
+  if(fManager) 
     treeD = fManager->GetTreeD() ;
- else
+  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 new branches
-  // -- generate file name if necessary
-  char * file =0;
-  if(gSystem->Getenv("CONFIG_SPLIT_FILE")){ //generating file name
-    file = new char[strlen(gAlice->GetBaseFile())+20] ;
-    sprintf(file,"%s/PHOS.Digits.root",gAlice->GetBaseFile()) ;
   }
-
-  TDirectory *cwd = gDirectory;
+  
   // -- create Digits branch
   Int_t bufferSize = 32000 ;    
   TBranch * digitsBranch = treeD->Branch("PHOS",&digits,bufferSize);
   digitsBranch->SetTitle(GetName());
-  if (file) {
-    digitsBranch->SetFile(file);
-    TIter next( digitsBranch->GetListOfBranches());
-    TBranch * sbr ;
-    while ((sbr=(TBranch*)next())) {
-      sbr->SetFile(file);
-    }   
-    cwd->cd();
-  } 
-    
+  
   // -- 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());
-  if (file) {
-    digitizerBranch->SetFile(file);
-    TIter next( digitizerBranch->GetListOfBranches());
-    TBranch * sbr;
-    while ((sbr=(TBranch*)next())) {
-      sbr->SetFile(file);
-    }   
-    cwd->cd();
-  }
-  
-  digitsBranch->Fill() ;      
-  digitizerBranch->Fill() ;
 
-  treeD->Write(0,kOverwrite) ;  
+  digitsBranch->Fill() ;
+  digitizerBranch->Fill() ; 
+  treeD->AutoSave() ; 
 }