]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDdigitizer.cxx
Not ready for previous changes, revert to stable code.
[u/mrichter/AliRoot.git] / TRD / AliTRDdigitizer.cxx
index 583c240e8bbef236acd0c8729ab898e09948bc79..b1f4dc3e35b6484157c6547fb1b76e969b137b96 100644 (file)
@@ -47,7 +47,7 @@
 #include "AliRunLoader.h"
 #include "AliLoader.h"
 #include "AliConfig.h"
-#include "AliRunDigitizer.h"
+#include "AliDigitizationInput.h"
 #include "AliRunLoader.h"
 #include "AliLoader.h"
 #include "AliLog.h"
@@ -70,6 +70,7 @@
 
 #include "Cal/AliTRDCalROC.h"
 #include "Cal/AliTRDCalDet.h"
+#include "Cal/AliTRDCalOnlineGainTableROC.h"
 
 ClassImp(AliTRDdigitizer)
 
@@ -118,9 +119,9 @@ AliTRDdigitizer::AliTRDdigitizer(const Text_t *name, const Text_t *title)
 }
 
 //_____________________________________________________________________________
-AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager
+AliTRDdigitizer::AliTRDdigitizer(AliDigitizationInput* digInput
                                , const Text_t *name, const Text_t *title)
-  :AliDigitizer(manager,name,title)
+  :AliDigitizer(digInput,name,title)
   ,fRunLoader(0)
   ,fDigitsManager(0)
   ,fSDigitsManager(0)
@@ -141,8 +142,8 @@ AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager
 }
 
 //_____________________________________________________________________________
-AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager)
-  :AliDigitizer(manager,"AliTRDdigitizer","TRD digitizer")
+AliTRDdigitizer::AliTRDdigitizer(AliDigitizationInput* digInput)
+  :AliDigitizer(digInput,"AliTRDdigitizer","TRD digitizer")
   ,fRunLoader(0)
   ,fDigitsManager(0)
   ,fSDigitsManager(0)
@@ -191,36 +192,25 @@ AliTRDdigitizer::~AliTRDdigitizer()
   // AliTRDdigitizer destructor
   //
 
-  if (fDigitsManager) {
-    delete fDigitsManager;
-    fDigitsManager      = 0;
-  }
-
-  if (fSDigitsManager) {
-    // s-digitsmanager will be deleted via list
-    fSDigitsManager     = 0;
-  }
+  delete fDigitsManager;
+  fDigitsManager      = 0;
 
+  // s-digitsmanager will be deleted via list
+  fSDigitsManager     = 0;
   if (fSDigitsManagerList) {
     fSDigitsManagerList->Delete();
     delete fSDigitsManagerList;
-    fSDigitsManagerList = 0;
   }
+  fSDigitsManagerList = 0;
 
-  if (fMasks) {
-    delete [] fMasks;
-    fMasks       = 0;
-  }
+  delete [] fMasks;
+  fMasks       = 0;
 
-  if (fMcmSim) {
-    delete fMcmSim;
-    fMcmSim = 0;
-  }
+  delete fMcmSim;
+  fMcmSim = 0;
 
-  if (fGeo) {
-    delete fGeo;
-    fGeo = 0;
-  }
+  delete fGeo;
+  fGeo = 0;
 
 }
 
@@ -261,7 +251,7 @@ void AliTRDdigitizer::Copy(TObject &d) const
 }
 
 //_____________________________________________________________________________
-void AliTRDdigitizer::Exec(const Option_t * const option)
+void AliTRDdigitizer::Digitize(const Option_t* option)
 {
   //
   // Executes the merging
@@ -284,26 +274,26 @@ void AliTRDdigitizer::Exec(const Option_t * const option)
     AliDebug(1,"AliRun object found on file.");
   }
   else {
-    inrl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
+    inrl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
     inrl->LoadgAlice();
     gAlice = inrl->GetAliRun();
     if (!gAlice) {
-      AliError("Could not find AliRun object.")
+      AliError("Could not find AliRun object.");
       return;
     }
   }
                                                                            
-  Int_t nInput = fManager->GetNinputs();
+  Int_t nInput = fDigInput->GetNinputs();
   fMasks       = new Int_t[nInput];
   for (iInput = 0; iInput < nInput; iInput++) {
-    fMasks[iInput] = fManager->GetMask(iInput);
+    fMasks[iInput] = fDigInput->GetMask(iInput);
   }
 
   //
   // Initialization
   //
 
-  AliRunLoader *orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
+  AliRunLoader *orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
 
   if (InitDetector()) {
 
@@ -336,7 +326,7 @@ void AliTRDdigitizer::Exec(const Option_t * const option)
     AliDebug(1,Form("Add input stream %d",iInput));
 
     // Check if the input tree exists
-    inrl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput));
+    inrl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(iInput));
     AliLoader *gime = inrl->GetLoader("TRDLoader");
 
     TTree *treees = gime->TreeS();
@@ -357,7 +347,7 @@ void AliTRDdigitizer::Exec(const Option_t * const option)
     sdigitsManager = new AliTRDdigitsManager();
     sdigitsManager->SetSDigits(kTRUE);
     
-    AliRunLoader *rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput));
+    AliRunLoader *rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(iInput));
     AliLoader *gimme = rl->GetLoader("TRDLoader");
     if (!gimme->TreeS()) 
       {
@@ -559,6 +549,7 @@ Bool_t AliTRDdigitizer::InitDetector()
   return kTRUE;
 
 }
+
 //_____________________________________________________________________________
 Bool_t AliTRDdigitizer::MakeBranch(TTree *tree) const
 {
@@ -616,19 +607,52 @@ Bool_t AliTRDdigitizer::MakeDigits()
 
   Float_t **hits = new Float_t*[kNdet];
   Int_t    *nhit = new Int_t[kNdet];
+  memset(nhit,0,kNdet*sizeof(Int_t));
 
   AliTRDarraySignal *signals = 0x0;
 
-  // 
-  if (calibration->GetNumberOfTimeBinsDCS() != AliTRDSimParam::Instance()->GetNTimeBins()) {
-    AliWarning(Form("Number of time bins is different to OCDB value [SIM=%d, OCDB=%d]"
-                  ,AliTRDSimParam::Instance()->GetNTimeBins()
-                   ,calibration->GetNumberOfTimeBinsDCS()));
+  // Check the number of time bins from simParam against OCDB,
+  // if OCDB value is not supposed to be used.
+  // As default, the value from OCDB is taken
+  if (AliTRDSimParam::Instance()->GetNTBoverwriteOCDB()) {
+    if (calibration->GetNumberOfTimeBinsDCS() != AliTRDSimParam::Instance()->GetNTimeBins()) {
+      AliWarning(Form("Number of time bins is different to OCDB value [SIM=%d, OCDB=%d]"
+                     ,AliTRDSimParam::Instance()->GetNTimeBins()
+                     ,calibration->GetNumberOfTimeBinsDCS()));
+    }
+    // Save the values for the raw data headers
+    fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDSimParam::Instance()->GetNTimeBins());
+  }
+  else {
+    // Get the OCDB values
+    Int_t NTB=calibration->GetNumberOfTimeBinsDCS();
+    if (NTB < 0) { // Currently -1 gets returned for "undefined" and "mixed",
+                   // one might go back to -1 undefined and -2 mixed?
+      AliError("No useful DCS information available for this run! Using standard values.");
+      // We fall back to the standard OCDB object, 
+      // cache the current run number..
+      Long64_t run = calibration->GetRun();
+      calibration->SetRun(0);
+      NTB=calibration->GetNumberOfTimeBinsDCS();
+      // ..to set it again
+      calibration->SetRun(run);
+      // If there's no standard OCDB object, we can still fail
+      if(NTB < 0){
+       AliFatal("No standard object found in the OCDB!");
+      }
+    }
+    // Save the values for the raw data headers
+    fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(NTB);
   }
+
+  // Save the values for the raw data headers
+  fDigitsManager->GetDigitsParam()->SetADCbaselineAll(AliTRDSimParam::Instance()->GetADCbaseline());
  
   // Sort all hits according to detector number
   if (!SortHits(hits,nhit)) {
     AliError("Sorting hits failed");
+    delete [] hits;
+    delete [] nhit;
     return kFALSE;
   }
 
@@ -636,8 +660,7 @@ Bool_t AliTRDdigitizer::MakeDigits()
   for (Int_t det = 0; det < kNdet; det++) {
 
     // Detectors that are switched off, not installed, etc.
-    if (( calibration->IsChamberInstalled(det)) &&
-       (!calibration->IsChamberMasked(det))    &&
+    if ((!calibration->IsChamberNoData(det))    &&
         ( fGeo->ChamberInGeometry(det))         &&
         (nhit[det] > 0)) {
 
@@ -646,11 +669,20 @@ Bool_t AliTRDdigitizer::MakeDigits()
       // Convert the hits of the current detector to detector signals
       if (!ConvertHits(det,hits[det],nhit[det],signals)) {
        AliError(Form("Conversion of hits failed for detector=%d",det));
+        delete [] hits;
+        delete [] nhit;
+        delete signals;
+        signals = 0x0;
         return kFALSE;
       }
+
       // Convert the detector signals to digits or s-digits
       if (!ConvertSignals(det,signals)) {
        AliError(Form("Conversion of signals failed for detector=%d",det));
+        delete [] hits;
+        delete [] nhit;
+        delete signals;
+        signals = 0x0;
        return kFALSE;
       }
 
@@ -665,7 +697,8 @@ Bool_t AliTRDdigitizer::MakeDigits()
   } // for: detector
 
   if (!fSDigits) {
-    if (AliDataLoader *trklLoader = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets")) {
+    if (AliDataLoader *trklLoader 
+          = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets")) {
       if (trklLoader->Tree())
         trklLoader->WriteData("OVERWRITE");
     }
@@ -674,10 +707,6 @@ Bool_t AliTRDdigitizer::MakeDigits()
   delete [] hits;
   delete [] nhit;
 
-  // Save the values for the raw data headers
-  fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDSimParam::Instance()->GetNTimeBins());
-  fDigitsManager->GetDigitsParam()->SetADCbaselineAll(AliTRDSimParam::Instance()->GetADCbaseline());
-
   return kTRUE;
 
 }
@@ -700,11 +729,10 @@ Bool_t AliTRDdigitizer::SortHits(Float_t **hits, Int_t *nhit)
   Int_t    nhitTrk  = 0;
 
   Int_t   *lhit     = new Int_t[kNdet];
+  memset(lhit,0,kNdet*sizeof(Int_t));
 
   for (Int_t det = 0; det < kNdet; det++) {
-    lhit[det] = 0;
-    nhit[det] = 0;
-    hits[det] = 0;
+    hits[det] = 0x0;
   }
 
   AliLoader *gimme = fRunLoader->GetLoader("TRDLoader");
@@ -714,6 +742,7 @@ Bool_t AliTRDdigitizer::SortHits(Float_t **hits, Int_t *nhit)
   TTree *hitTree = gimme->TreeH();
   if (hitTree == 0x0) {
     AliError("Can not get TreeH");
+    delete [] lhit;
     return kFALSE;
   }
   fTRD->SetTreeAddress();
@@ -852,13 +881,15 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det
   AliTRDCalROC       *calT0ROC          = 0;
   Float_t             calT0DetValue     = 0.0;
   const AliTRDCalDet *calT0Det          = calibration->GetT0Det();  
+  Double_t            calExBDetValue    = 0.0;
+  const AliTRDCalDet *calExBDet         = calibration->GetExBDet();
 
   if (simParam->TRFOn()) {
     timeBinTRFend = ((Int_t) (simParam->GetTRFhi() 
                   * commonParam->GetSamplingFrequency())) - 1;
   }
 
-  Int_t   nTimeTotal   = simParam->GetNTimeBins();
+  Int_t   nTimeTotal   = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
   Float_t samplingRate = commonParam->GetSamplingFrequency();
   Float_t elAttachProp = simParam->GetElAttachProp() / 100.0; 
 
@@ -901,6 +932,7 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det
     calVdriftDetValue = calVdriftDet->GetValue(det);
     calT0ROC          = calibration->GetT0ROC(det);
     calT0DetValue     = calT0Det->GetValue(det);
+    calExBDetValue    = calExBDet->GetValue(det);
 
     // Go to the local coordinate system:
     // loc[0] - col  direction in amplification or driftvolume
@@ -949,7 +981,7 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det
     Float_t  driftvelocity  = calVdriftDetValue * calVdriftROC->GetValue(colE,rowE);
     Double_t absdriftlength = TMath::Abs(driftlength);
     if (commonParam->ExBOn()) {
-      absdriftlength /= TMath::Sqrt(GetLorentzFactor(driftvelocity));
+      absdriftlength /= TMath::Sqrt(1.0 / (1.0 + calExBDetValue*calExBDetValue));
     }
 
     // Loop over all electrons of this hit
@@ -979,16 +1011,14 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det
           
       // Apply the diffusion smearing
       if (simParam->DiffusionOn()) {
-        if (!(Diffusion(driftvelocity,absdriftlength,locR,locC,locT))) {
+        if (!(Diffusion(driftvelocity,absdriftlength,calExBDetValue,locR,locC,locT))) {
           continue;
        }
       }
 
       // Apply E x B effects (depends on drift direction)
-      if (commonParam->ExBOn()) { 
-        if (!(ExB(driftvelocity,driftlength,locC))) {
-          continue;
-       }
+      if (commonParam->ExBOn()) {
+        locC = locC + calExBDetValue * driftlength;
       }
 
       // The electron position after diffusion and ExB in pad coordinates.
@@ -1212,14 +1242,21 @@ Bool_t AliTRDdigitizer::Signal2ADC(Int_t det, AliTRDarraySignal *signals)
 
   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
-  Int_t nTimeTotal = simParam->GetNTimeBins();
+  Int_t nTimeTotal = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
+  if (fSDigitsManager->GetDigitsParam()->GetNTimeBins(det)) {
+    nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(det);
+  }
+  else {
+    AliFatal("Could not get number of time bins");
+    return kFALSE;
+  }
 
-  // The gainfactor calibration objects
+  // The gain factor calibration objects
   const AliTRDCalDet *calGainFactorDet      = calibration->GetGainFactorDet();  
-  AliTRDCalROC       *calGainFactorROC      = 0;
+  AliTRDCalROC       *calGainFactorROC      = 0x0;
   Float_t             calGainFactorDetValue = 0.0;
 
-  AliTRDarrayADC     *digits = 0x0;
+  AliTRDarrayADC     *digits                = 0x0;
 
   if (!signals) {
     AliError(Form("Signals array for detector %d does not exist\n",det));
@@ -1256,13 +1293,14 @@ Bool_t AliTRDdigitizer::Signal2ADC(Int_t det, AliTRDarraySignal *signals)
 
       // Check whether pad is masked
       // Bridged pads are not considered yet!!!
-      if (calibration->IsPadMasked(det,col,row) || calibration->IsPadNotConnected(det,col,row)) {
+      if (calibration->IsPadMasked(det,col,row) || 
+          calibration->IsPadNotConnected(det,col,row)) {
         continue;
       }
 
       // The gain factors
-      Float_t padgain = calGainFactorDetValue 
-                      * calGainFactorROC->GetValue(col,row);
+      Float_t padgain    = calGainFactorDetValue 
+                         * calGainFactorROC->GetValue(col,row);
       if (padgain <= 0) {
         AliError(Form("Not a valid gain %f, %d %d %d",padgain,det,col,row));
       }
@@ -1328,10 +1366,9 @@ Bool_t AliTRDdigitizer::Signal2SDigits(Int_t det, AliTRDarraySignal *signals)
 
   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
-  Int_t nTimeTotal = AliTRDSimParam::Instance()->GetNTimeBins();
+  Int_t nTimeTotal = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
 
   // Get the container for the digits of this detector
-
   if (!fDigitsManager->HasSDigits()) {
     AliError("Digits manager has no s-digits");
     return kFALSE;
@@ -1402,6 +1439,7 @@ Bool_t AliTRDdigitizer::Digits2SDigits(AliTRDdigitsManager * const manDig
   //                      / convert;
 
   // The gainfactor calibration objects
+  // Not used since these digits are supposed to be from real raw data
   //const AliTRDCalDet *calGainFactorDet      = calibration->GetGainFactorDet();  
   //AliTRDCalROC       *calGainFactorROC      = 0;
   //Float_t             calGainFactorDetValue = 0.0;
@@ -1713,7 +1751,12 @@ Bool_t AliTRDdigitizer::ConvertSDigits()
   }
 
   // Save the values for the raw data headers
-  fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDSimParam::Instance()->GetNTimeBins());
+  if (AliTRDSimParam::Instance()->GetNTBoverwriteOCDB()) {
+    fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDSimParam::Instance()->GetNTimeBins());
+  }
+  else {
+    fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDcalibDB::Instance()->GetNumberOfTimeBinsDCS());
+  }
   fDigitsManager->GetDigitsParam()->SetADCbaselineAll(AliTRDSimParam::Instance()->GetADCbaseline());
 
   return kTRUE;
@@ -1742,7 +1785,7 @@ Bool_t AliTRDdigitizer::CopyDictionary(Int_t det)
 
   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
-  Int_t nTimeTotal = AliTRDSimParam::Instance()->GetNTimeBins();
+  Int_t nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(det);
 
   Int_t row  = 0;
   Int_t col  = 0;
@@ -1864,6 +1907,7 @@ void AliTRDdigitizer::InitOutput(Int_t iEvent)
   
 //_____________________________________________________________________________
 Int_t AliTRDdigitizer::Diffusion(Float_t vdrift, Double_t absdriftlength
+                               , Double_t exbvalue
                                , Double_t &lRow, Double_t &lCol, Double_t &lTime)
 {
   //
@@ -1880,8 +1924,14 @@ Int_t AliTRDdigitizer::Diffusion(Float_t vdrift, Double_t absdriftlength
     Float_t sigmaT    = driftSqrt * diffT;
     Float_t sigmaL    = driftSqrt * diffL;
     lRow  = gRandom->Gaus(lRow ,sigmaT);
-    lCol  = gRandom->Gaus(lCol ,sigmaT * GetLorentzFactor(vdrift));
-    lTime = gRandom->Gaus(lTime,sigmaL * GetLorentzFactor(vdrift));
+    if (AliTRDCommonParam::Instance()->ExBOn()) {
+      lCol  = gRandom->Gaus(lCol ,sigmaT * 1.0 / (1.0 + exbvalue*exbvalue));
+      lTime = gRandom->Gaus(lTime,sigmaL * 1.0 / (1.0 + exbvalue*exbvalue));
+    }
+    else {
+      lCol  = gRandom->Gaus(lCol ,sigmaT);
+      lTime = gRandom->Gaus(lTime,sigmaL);
+    }
 
     return 1;
 
@@ -1893,40 +1943,7 @@ Int_t AliTRDdigitizer::Diffusion(Float_t vdrift, Double_t absdriftlength
   }
 
 }
-
-//_____________________________________________________________________________
-Float_t AliTRDdigitizer::GetLorentzFactor(Float_t vd)
-{
-  //
-  // Returns the Lorentz factor
-  //
-
-  Double_t omegaTau      = AliTRDCommonParam::Instance()->GetOmegaTau(vd);
-  Double_t lorentzFactor = 1.0;
-  if (AliTRDCommonParam::Instance()->ExBOn()) {
-    lorentzFactor = 1.0 / (1.0 + omegaTau*omegaTau);
-  }
-
-  return lorentzFactor;
-
-}
   
-//_____________________________________________________________________________
-Int_t AliTRDdigitizer::ExB(Float_t vdrift, Double_t driftlength, Double_t &lCol)
-{
-  //
-  // Applies E x B effects to the position of a single electron.
-  // Depends on signed drift length.
-  //
-
-  lCol = lCol 
-       + AliTRDCommonParam::Instance()->GetOmegaTau(vdrift) 
-       * driftlength;
-
-  return 1;
-
-}
-
 //_____________________________________________________________________________
 void AliTRDdigitizer::RunDigitalProcessing(Int_t det)
 {
@@ -1941,19 +1958,20 @@ void AliTRDdigitizer::RunDigitalProcessing(Int_t det)
     return;
 
   //Call the methods in the mcm class using the temporary array as input  
-  for(Int_t rob = 0; rob < digits->GetNrow() / 2; rob++)
-  {
-    for(Int_t mcm = 0; mcm < 16; mcm++)
-    {
-      fMcmSim->Init(det, rob, mcm); 
-      fMcmSim->SetDataByPad(digits, fDigitsManager);
-      fMcmSim->Filter();
-      if (feeParam->GetTracklet()) {
-        fMcmSim->Tracklet();
-        fMcmSim->StoreTracklets();
+  // process the data in the same order as in hardware
+  for (Int_t side = 0; side <= 1; side++) {
+    for(Int_t rob = side; rob < digits->GetNrow() / 2; rob += 2) {
+      for(Int_t mcm = 0; mcm < 16; mcm++) {
+       fMcmSim->Init(det, rob, mcm); 
+       fMcmSim->SetDataByPad(digits, fDigitsManager);
+       fMcmSim->Filter();
+       if (feeParam->GetTracklet()) {
+         fMcmSim->Tracklet();
+         fMcmSim->StoreTracklets();
+       }
+       fMcmSim->ZSMapping();
+       fMcmSim->WriteData(digits);
       }
-      fMcmSim->ZSMapping();
-      fMcmSim->WriteData(digits);
     }
   }
 }