]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliEMCALDigitizer.cxx
changing to AMANDA protocol version 2
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALDigitizer.cxx
index 17e1c66a5730b30d8caaf9a2c424ffaa81b512f5..12ceca89ec337dab28ac2f57188b5b12174bec86 100644 (file)
 //_________________________________________________________________________________
 
 // --- ROOT system ---
-#include "TROOT.h"
-#include "TTree.h"
-#include "TSystem.h"
-#include "TBenchmark.h"
-#include "TList.h"
-#include "TH1.h"
-#include "TBrowser.h"
-#include "TObjectTable.h"
+#include <TROOT.h>
+#include <TTree.h>
+#include <TSystem.h>
+#include <TBenchmark.h>
+#include <TList.h>
+#include <TH1.h>
+#include <TBrowser.h>
+#include <TObjectTable.h>
+#include <TRandom.h>
+#include <cassert>
 
 // --- AliRoot header files ---
+#include "AliLog.h"
 #include "AliRun.h"
 #include "AliRunDigitizer.h"
 #include "AliRunLoader.h"
+#include "AliCDBManager.h"
+#include "AliCDBEntry.h"
 #include "AliEMCALDigit.h"
 #include "AliEMCAL.h"
 #include "AliEMCALLoader.h"
@@ -85,63 +90,128 @@ ClassImp(AliEMCALDigitizer)
 
 
 //____________________________________________________________________________ 
-  AliEMCALDigitizer::AliEMCALDigitizer():AliDigitizer("",""),
-                                      fInput(0),
-                                      fInputFileNames(0x0),
-                                      fEventNames(0x0)
+AliEMCALDigitizer::AliEMCALDigitizer()
+  : AliDigitizer("",""),
+    fDefaultInit(kTRUE),
+    fDigitsInRun(0),
+    fInit(0),
+    fInput(0),
+    fInputFileNames(0x0),
+    fEventNames(0x0),
+    fDigitThreshold(0),
+    fMeanPhotonElectron(0),
+    fPedestal(0),
+    fSlope(0),
+    fPinNoise(0),
+    fTimeResolution(0),
+    fTimeThreshold(0),    
+    fTimeSignalLength(0),
+    fADCchannelEC(0),
+    fADCpedestalEC(0),
+    fNADCEC(0),
+    fEventFolderName(""),
+    fFirstEvent(0),
+    fLastEvent(0),
+    fControlHists(0),
+    fHists(0),fCalibData(0x0)
 {
   // ctor
   InitParameters() ; 
-  fDefaultInit = kTRUE ; 
   fManager = 0 ;                     // We work in the standalong mode
-  fEventFolderName = "" ; 
 }
 
 //____________________________________________________________________________ 
-AliEMCALDigitizer::AliEMCALDigitizer(TString alirunFileName, TString eventFolderName):
-  AliDigitizer("EMCAL"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
-  fInputFileNames(0), fEventNames(0), fEventFolderName(eventFolderName)
+AliEMCALDigitizer::AliEMCALDigitizer(TString alirunFileName, TString eventFolderName)
+  : AliDigitizer("EMCAL"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
+    fDefaultInit(kFALSE),
+    fDigitsInRun(0),
+    fInit(0),
+    fInput(0),
+    fInputFileNames(0), 
+    fEventNames(0), 
+    fDigitThreshold(0),
+    fMeanPhotonElectron(0),
+    fPedestal(0),
+    fSlope(0),
+    fPinNoise(0),
+    fTimeResolution(0),
+    fTimeThreshold(0),
+    fTimeSignalLength(0),
+    fADCchannelEC(0),
+    fADCpedestalEC(0),
+    fNADCEC(0),
+    fEventFolderName(eventFolderName),
+    fFirstEvent(0),
+    fLastEvent(0),
+    fControlHists(0),
+    fHists(0)
 {
   // ctor
-
   InitParameters() ; 
   Init() ;
-  fDefaultInit = kFALSE ; 
   fManager = 0 ;                     // We work in the standalong mode
 }
 
 //____________________________________________________________________________ 
-AliEMCALDigitizer::AliEMCALDigitizer(const AliEMCALDigitizer & d) : AliDigitizer(d)
+AliEMCALDigitizer::AliEMCALDigitizer(const AliEMCALDigitizer & d) 
+  : AliDigitizer(d.GetName(),d.GetTitle()),
+    fDefaultInit(d.fDefaultInit),
+    fDigitsInRun(d.fDigitsInRun),
+    fInit(d.fInit),
+    fInput(d.fInput),
+    fInputFileNames(d.fInputFileNames),
+    fEventNames(d.fEventNames),
+    fDigitThreshold(d.fDigitThreshold),
+    fMeanPhotonElectron(d.fMeanPhotonElectron),
+    fPedestal(d.fPedestal),
+    fSlope(d.fSlope),
+    fPinNoise(d.fPinNoise),
+    fTimeResolution(d.fTimeResolution),
+    fTimeThreshold(d.fTimeThreshold),
+    fTimeSignalLength(d.fTimeSignalLength),
+    fADCchannelEC(d.fADCchannelEC),
+    fADCpedestalEC(d.fADCpedestalEC),
+    fNADCEC(d.fNADCEC),
+    fEventFolderName(d.fEventFolderName),
+    fFirstEvent(d.fFirstEvent),
+    fLastEvent(d.fLastEvent),
+    fControlHists(d.fControlHists),
+    fHists(d.fHists),fCalibData(d.fCalibData)
 {
   // copyy ctor 
-
-  SetName(d.GetName()) ; 
-  SetTitle(d.GetTitle()) ; 
-  fDigitThreshold = d.fDigitThreshold ; 
-  fMeanPhotonElectron = d.fMeanPhotonElectron ; 
-  fPedestal           = d.fPedestal ; 
-  fSlope              = d.fSlope ; 
-  fPinNoise           = d.fPinNoise ; 
-  fTimeResolution     = d.fTimeResolution ; 
-  fTimeThreshold      = d.fTimeThreshold ; 
-  fTimeSignalLength   = d.fTimeSignalLength ; 
-  fADCchannelEC       = d.fADCchannelEC ; 
-  fADCpedestalEC      = d.fADCpedestalEC ; 
-  fNADCEC             = d.fNADCEC ;
-  fEventFolderName    = d.fEventFolderName;
  }
 
 //____________________________________________________________________________ 
-AliEMCALDigitizer::AliEMCALDigitizer(AliRunDigitizer * rd):
- AliDigitizer(rd,"EMCAL"+AliConfig::Instance()->GetDigitizerTaskName()),
- fEventFolderName(0)
+AliEMCALDigitizer::AliEMCALDigitizer(AliRunDigitizer * rd)
+  : AliDigitizer(rd,"EMCAL"+AliConfig::Instance()->GetDigitizerTaskName()),
+    fDefaultInit(kFALSE),
+    fDigitsInRun(0),
+    fInit(0),
+    fInput(0),
+    fInputFileNames(0),
+    fEventNames(0),
+    fDigitThreshold(0.),
+    fMeanPhotonElectron(0),
+    fPedestal(0),
+    fSlope(0.),
+    fPinNoise(0),
+    fTimeResolution(0.),
+    fTimeThreshold(0),
+    fTimeSignalLength(0),
+    fADCchannelEC(0),
+    fADCpedestalEC(0),
+    fNADCEC(0),
+    fEventFolderName(0),
+    fFirstEvent(0),
+    fLastEvent(0),
+    fControlHists(0),
+    fHists(0)
 {
   // ctor Init() is called by RunDigitizer
   fManager = rd ; 
   fEventFolderName = fManager->GetInputFolderName(0) ;
   SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
   InitParameters() ; 
-  fDefaultInit = kFALSE ; 
 }
 
 //____________________________________________________________________________ 
@@ -185,16 +255,19 @@ void AliEMCALDigitizer::Digitize(Int_t event)
   rl->GetEvent(readEvent);
 
   TClonesArray * digits = emcalLoader->Digits() ; 
-  digits->Clear() ;
+  digits->Delete() ;  
 
   // Load Geometry
-  // const AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
-  rl->LoadgAlice(); 
-  AliRun * gAlice = rl->GetAliRun(); 
-  AliEMCAL * emcal  = (AliEMCAL*)gAlice->GetDetector("EMCAL");
-  AliEMCALGeometry * geom = emcal->GetGeometry();
+  AliEMCALGeometry *geom = 0;
+  if (rl->GetAliRun()) {
+    AliEMCAL * emcal  = (AliEMCAL*)rl->GetAliRun()->GetDetector("EMCAL");
+    geom = emcal->GetGeometry();
+  }
+  else 
+    AliFatal("Could not get AliRun from runLoader");
 
   if(isTrd1Geom < 0) { 
+    AliInfo(Form(" get Geometry %s : %s ", geom->GetName(),geom->GetTitle()));
     TString ng(geom->GetName());
     isTrd1Geom = 0;
     if(ng.Contains("SHISH") &&  ng.Contains("TRD1")) isTrd1Geom = 1;
@@ -262,12 +335,12 @@ void AliEMCALDigitizer::Digitize(Int_t event)
   TClonesArray * ticks = new TClonesArray("AliEMCALTick",1000) ;
 
   //Put Noise contribution
-  for(absID = 1; absID <= nEMC; absID++){
+  for(absID = 0; absID < nEMC; absID++){ // Nov 30, 2006 by PAI; was from 1 to nEMC
     Float_t amp = 0 ;
     // amplitude set to zero, noise will be added later
-    new((*digits)[absID-1]) AliEMCALDigit( -1, -1, absID, 0, TimeOfNoise() ) ;
+    new((*digits)[absID]) AliEMCALDigit( -1, -1, absID, 0, TimeOfNoise() ); // absID-1->absID
     //look if we have to add signal?
-    digit = dynamic_cast<AliEMCALDigit *>(digits->At(absID-1)) ;
+    digit = dynamic_cast<AliEMCALDigit *>(digits->At(absID)); // absID-1->absID
     
     if(absID==nextSig){
       //Add SDigits from all inputs    
@@ -354,7 +427,6 @@ void AliEMCALDigitizer::Digitize(Int_t event)
   digits->Compress() ;  
   
   Int_t ndigits = digits->GetEntriesFast() ; 
-  digits->Expand(ndigits) ;
   
   //Set indexes in list of digits and fill hists.
   AliEMCALHistoUtilities::FillH1(fHists, 0, Double_t(ndigits));
@@ -368,6 +440,11 @@ void AliEMCALDigitizer::Digitize(Int_t event)
     AliEMCALHistoUtilities::FillH1(fHists, 2, double(digit->GetAmp()));
     AliEMCALHistoUtilities::FillH1(fHists, 3, double(energy));
     AliEMCALHistoUtilities::FillH1(fHists, 4, double(digit->GetId()));
+    //    if(digit->GetId() == nEMC) {
+    //  printf(" i %i \n", i );
+    //  digit->Dump();
+    //  assert(0);
+    //}
   }
   AliEMCALHistoUtilities::FillH1(fHists, 1, esum);
 }
@@ -376,38 +453,29 @@ void AliEMCALDigitizer::Digitize(Int_t event)
 Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy, Int_t AbsId)
 { 
   // Returns digitized value of the energy in a cell absId
-  // Loader
-  AliRunLoader *rl = AliRunLoader::GetRunLoader();
-  AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>
-    (rl->GetDetectorLoader("EMCAL"));
-  
-  // Load EMCAL Geometry
-  rl->LoadgAlice(); 
-  AliRun * gAlice = rl->GetAliRun(); 
-  AliEMCAL * emcal  = (AliEMCAL*)gAlice->GetDetector("EMCAL");
-  AliEMCALGeometry * geom = emcal->GetGeometry();
+
+  // Load Geometry
+  const AliEMCALGeometry * geom = AliEMCALGeometry::GetInstance();
 
   if (geom==0)
-    AliFatal("Did not get geometry from EMCALLoader") ;
+    AliFatal("Did not get geometry from EMCALLoader");
 
   Int_t iSupMod = -1;
-  Int_t nTower  = -1;
+  Int_t nModule  = -1;
   Int_t nIphi   = -1;
   Int_t nIeta   = -1;
   Int_t iphi    = -1;
   Int_t ieta    = -1;
   Int_t channel = -999; 
 
-  Bool_t bCell = geom->GetCellIndex(AbsId, iSupMod, nTower, nIphi, nIeta) ;
+  Bool_t bCell = geom->GetCellIndex(AbsId, iSupMod, nModule, nIphi, nIeta) ;
   if(!bCell)
-    Error("DigitizeEnergy","Wrong cell id number") ;
-  geom->GetCellPhiEtaIndexInSModule(iSupMod,nTower,nIphi, nIeta,iphi,ieta);
+    Error("DigitizeEnergy","Wrong cell id number : AbsId %i ", AbsId) ;
+  geom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
   
-  if(emcalLoader->CalibData()) {
-    fADCpedestalEC = emcalLoader->CalibData()
-      ->GetADCpedestal(iSupMod,ieta,iphi);
-    fADCchannelEC = emcalLoader->CalibData()
-      ->GetADCchannel(iSupMod,ieta,iphi);
+  if(fCalibData) {
+    fADCpedestalEC = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
+    fADCchannelEC = fCalibData->GetADCchannel(iSupMod,ieta,iphi);
   }
   
   channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalEC)/fADCchannelEC ))  ;
@@ -433,6 +501,7 @@ void AliEMCALDigitizer::Exec(Option_t *option)
   } 
 
   if (strstr(option,"print")) {
+
     Print();
     return ; 
   }
@@ -446,8 +515,9 @@ void AliEMCALDigitizer::Exec(Option_t *option)
   // Post Digitizer to the white board
   emcalLoader->PostDigitizer(this) ;
   
-  if (fLastEvent == -1) 
+  if (fLastEvent == -1)  {
     fLastEvent = rl->GetNumberOfEvents() - 1 ;
+  }
   else if (fManager) 
     fLastEvent = fFirstEvent ; // what is this ??
 
@@ -475,8 +545,8 @@ void AliEMCALDigitizer::Exec(Option_t *option)
 
   if(strstr(option,"tim")){
     gBenchmark->Stop("EMCALDigitizer");
-    printf("Exec: took %f seconds for Digitizing %f seconds per event", 
-        gBenchmark->GetCpuTime("EMCALDigitizer"), gBenchmark->GetCpuTime("EMCALDigitizer")/nEvents ) ;
+    AliInfo(Form("Exec: took %f seconds for Digitizing %f seconds per event", 
+        gBenchmark->GetCpuTime("EMCALDigitizer"), gBenchmark->GetCpuTime("EMCALDigitizer")/nEvents )) ;
   } 
 }
 
@@ -537,26 +607,29 @@ Bool_t AliEMCALDigitizer::Init()
   emcalLoader->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
 
   //PH  Print();
-  
+  //Calibration instance
+  fCalibData = emcalLoader->CalibData();
   return fInit ;    
 }
 
 //____________________________________________________________________________ 
 void AliEMCALDigitizer::InitParameters()
 { 
-  //parameter initialization for digitizer
-  // Tune parameters - 24-nov-04
+  // Parameter initialization for digitizer
+  // Tune parameters - 24-nov-04; Apr 29, 2007
 
-  fMeanPhotonElectron = 3300 ; // electrons per GeV 
-  fPinNoise           = 0.004; 
+  fMeanPhotonElectron = 3300 // electrons per GeV 
+  fPinNoise           = 0.010; // pin noise in GEV from analysis test beam data 
   if (fPinNoise == 0. ) 
     Warning("InitParameters", "No noise added\n") ; 
   fDigitThreshold     = fPinNoise * 3; // 3 * sigma
   fTimeResolution     = 0.3e-9 ; // 300 psc
   fTimeSignalLength   = 1.0e-9 ;
 
-  fADCchannelEC    = 0.00305; // 200./65536 - width of one ADC channel in GeV
-  fADCpedestalEC   = 0.009 ;  // GeV
+  // These defaults are normally not used. 
+  // Values are read from calibration database instead
+  fADCchannelEC    = 0.0153; // Update 24 Apr 2007: 250./16/1024 - width of one ADC channel in GeV
+  fADCpedestalEC   = 0.0 ;  // GeV
   fNADCEC          = (Int_t) TMath::Power(2,16) ;  // number of channels in Tower ADC - 65536
 
   fTimeThreshold      = 0.001*10000000 ; // Means 1 MeV in terms of SDigits amplitude ??
@@ -602,6 +675,7 @@ void AliEMCALDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
       return ;
     }
     // need to increase the arrays
+    // MvL: This code only works when fInput == 1, I think.
     TString tempo = fInputFileNames[fInput-1] ; 
     delete [] fInputFileNames ; 
     fInputFileNames = new TString[fInput+1] ; 
@@ -747,8 +821,10 @@ void AliEMCALDigitizer::WriteDigits()
   // -- create Digits branch
   Int_t bufferSize = 32000 ;    
   TBranch * digitsBranch = 0;
-  if ((digitsBranch = treeD->GetBranch("EMCAL")))
+  if ((digitsBranch = treeD->GetBranch("EMCAL"))) {
     digitsBranch->SetAddress(&digits);
+    AliWarning("Digits Branch already exists. Not all digits will be visible");
+  }
   else
     treeD->Branch("EMCAL","TClonesArray",&digits,bufferSize);
   //digitsBranch->SetTitle(fEventFolderName);