]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDdigitizer.cxx
Add a protection againts runs with missing DCS information in the OCDB
[u/mrichter/AliRoot.git] / TRD / AliTRDdigitizer.cxx
index dd75f254b69955c6b740c499e782f839aaf664c4..5979125772fed1afcdd78956d705bbcbd4924740 100644 (file)
 //                                                                        //
 ////////////////////////////////////////////////////////////////////////////
 
-#include <TF1.h>
-#include <TFile.h>
-#include <TGeoGlobalMagField.h>
 #include <TGeoManager.h>
 #include <TList.h>
 #include <TMath.h>
-#include <TROOT.h>
 #include <TRandom.h>
-#include <TTask.h>
 #include <TTree.h>
-#include <TVector.h>
 
 #include "AliRun.h"
 #include "AliMC.h"
 #include "AliRunLoader.h"
 #include "AliLoader.h"
 #include "AliConfig.h"
-#include "AliMagF.h"
-#include "AliRunDigitizer.h"
+#include "AliDigitizationInput.h"
 #include "AliRunLoader.h"
 #include "AliLoader.h"
 #include "AliLog.h"
+
 #include "AliTRD.h"
 #include "AliTRDhit.h"
 #include "AliTRDdigitizer.h"
-
 #include "AliTRDarrayDictionary.h"
 #include "AliTRDarrayADC.h"
 #include "AliTRDarraySignal.h"
 #include "AliTRDCommonParam.h"
 #include "AliTRDfeeParam.h"
 #include "AliTRDmcmSim.h"
+#include "AliTRDdigitsParam.h"
+
 #include "Cal/AliTRDCalROC.h"
 #include "Cal/AliTRDCalDet.h"
+#include "Cal/AliTRDCalOnlineGainTableROC.h"
 
 ClassImp(AliTRDdigitizer)
 
@@ -87,6 +83,7 @@ AliTRDdigitizer::AliTRDdigitizer()
   ,fSDigitsManagerList(0)
   ,fTRD(0)
   ,fGeo(0)
+  ,fMcmSim(new AliTRDmcmSim)
   ,fEvent(0)
   ,fMasks(0)
   ,fCompress(kTRUE)
@@ -97,8 +94,6 @@ AliTRDdigitizer::AliTRDdigitizer()
   // AliTRDdigitizer default constructor
   //
   
-  Init();
-
 }
 
 //_____________________________________________________________________________
@@ -110,6 +105,7 @@ AliTRDdigitizer::AliTRDdigitizer(const Text_t *name, const Text_t *title)
   ,fSDigitsManagerList(0)
   ,fTRD(0)
   ,fGeo(0)
+  ,fMcmSim(new AliTRDmcmSim)
   ,fEvent(0)
   ,fMasks(0)
   ,fCompress(kTRUE)
@@ -120,20 +116,19 @@ AliTRDdigitizer::AliTRDdigitizer(const Text_t *name, const Text_t *title)
   // AliTRDdigitizer constructor
   //
 
-  Init();
-
 }
 
 //_____________________________________________________________________________
-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)
   ,fSDigitsManagerList(0)
   ,fTRD(0)
   ,fGeo(0)
+  ,fMcmSim(new AliTRDmcmSim)
   ,fEvent(0)
   ,fMasks(0)
   ,fCompress(kTRUE)
@@ -144,19 +139,18 @@ AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager
   // AliTRDdigitizer constructor
   //
 
-  Init();
-
 }
 
 //_____________________________________________________________________________
-AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager)
-  :AliDigitizer(manager,"AliTRDdigitizer","TRD digitizer")
+AliTRDdigitizer::AliTRDdigitizer(AliDigitizationInput* digInput)
+  :AliDigitizer(digInput,"AliTRDdigitizer","TRD digitizer")
   ,fRunLoader(0)
   ,fDigitsManager(0)
   ,fSDigitsManager(0)
   ,fSDigitsManagerList(0)
   ,fTRD(0)
   ,fGeo(0)
+  ,fMcmSim(new AliTRDmcmSim)
   ,fEvent(0)
   ,fMasks(0)
   ,fCompress(kTRUE)
@@ -167,32 +161,6 @@ AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager)
   // AliTRDdigitizer constructor
   //
 
-  Init();
-
-}
-
-//_____________________________________________________________________________
-Bool_t AliTRDdigitizer::Init()
-{
-  //
-  // Initialize the digitizer with default values
-  //
-
-  fRunLoader          = 0;
-  fDigitsManager      = 0;
-  fSDigitsManager     = 0;
-  fSDigitsManagerList = 0;
-  fTRD                = 0;
-  fGeo                = 0;
-
-  fEvent              = 0;
-  fMasks              = 0;
-  fCompress           = kTRUE;
-  fSDigits            = kFALSE;
-  fMergeSignalOnly    = kFALSE;
-  
-  return AliDigitizer::Init();
-
 }
 
 //_____________________________________________________________________________
@@ -204,6 +172,7 @@ AliTRDdigitizer::AliTRDdigitizer(const AliTRDdigitizer &d)
   ,fSDigitsManagerList(0)
   ,fTRD(0)
   ,fGeo(0)
+  ,fMcmSim(new AliTRDmcmSim)
   ,fEvent(0)
   ,fMasks(0)
   ,fCompress(d.fCompress)
@@ -223,31 +192,25 @@ AliTRDdigitizer::~AliTRDdigitizer()
   // AliTRDdigitizer destructor
   //
 
-  if (fDigitsManager) {
-    delete fDigitsManager;
-    fDigitsManager      = 0;
-  }
-
-  if (fDigitsManager) {              //typo? fSDigitsManager?
-    delete fSDigitsManager;
-    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 (fGeo) {
-    delete fGeo;
-    fGeo = 0;
-  }
+  delete fMcmSim;
+  fMcmSim = 0;
+
+  delete fGeo;
+  fGeo = 0;
 
 }
 
@@ -288,7 +251,7 @@ void AliTRDdigitizer::Copy(TObject &d) const
 }
 
 //_____________________________________________________________________________
-void AliTRDdigitizer::Exec(Option_t *option)
+void AliTRDdigitizer::Digitize(const Option_t* option)
 {
   //
   // Executes the merging
@@ -311,26 +274,26 @@ void AliTRDdigitizer::Exec(Option_t *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()) {
 
@@ -363,7 +326,7 @@ void AliTRDdigitizer::Exec(Option_t *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();
@@ -384,7 +347,7 @@ void AliTRDdigitizer::Exec(Option_t *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()) 
       {
@@ -484,7 +447,7 @@ Bool_t AliTRDdigitizer::Open(const Char_t *file, Int_t nEvent)
 }
 
 //_____________________________________________________________________________
-Bool_t AliTRDdigitizer::Open(AliRunLoader *runLoader, Int_t nEvent)
+Bool_t AliTRDdigitizer::Open(AliRunLoader * const runLoader, Int_t nEvent)
 {
   //
   // Opens a ROOT-file with TRD-hits and reads in the hit-tree
@@ -586,6 +549,7 @@ Bool_t AliTRDdigitizer::InitDetector()
   return kTRUE;
 
 }
+
 //_____________________________________________________________________________
 Bool_t AliTRDdigitizer::MakeBranch(TTree *tree) const
 {
@@ -643,12 +607,38 @@ 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;
+
+  // 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 {
+    // Save the values for the raw data headers
+    if (calibration->GetNumberOfTimeBinsDCS() == -1) {
+      AliFatal("No useful DCS information available for this run!");
+    }
+    fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(calibration->GetNumberOfTimeBinsDCS());
+  }
+
+  // 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;
   }
 
@@ -656,8 +646,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)) {
 
@@ -666,11 +655,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;
       }
 
@@ -684,6 +682,14 @@ Bool_t AliTRDdigitizer::MakeDigits()
 
   } // for: detector
 
+  if (!fSDigits) {
+    if (AliDataLoader *trklLoader 
+          = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets")) {
+      if (trklLoader->Tree())
+        trklLoader->WriteData("OVERWRITE");
+    }
+  }
+    
   delete [] hits;
   delete [] nhit;
 
@@ -709,11 +715,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");
@@ -723,6 +728,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();
@@ -800,7 +806,9 @@ Bool_t AliTRDdigitizer::SortHits(Float_t **hits, Int_t *nhit)
 }
 
 //_____________________________________________________________________________
-Bool_t AliTRDdigitizer::ConvertHits(Int_t det, Float_t *hits, Int_t nhit
+Bool_t AliTRDdigitizer::ConvertHits(Int_t det
+                                  , const Float_t * const hits
+                                  , Int_t nhit
                                   , AliTRDarraySignal *signals)
 {
   //
@@ -859,13 +867,15 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det, Float_t *hits, Int_t nhit
   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   = calibration->GetNumberOfTimeBins();
+  Int_t   nTimeTotal   = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
   Float_t samplingRate = commonParam->GetSamplingFrequency();
   Float_t elAttachProp = simParam->GetElAttachProp() / 100.0; 
 
@@ -908,6 +918,7 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det, Float_t *hits, Int_t nhit
     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
@@ -956,7 +967,7 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det, Float_t *hits, Int_t nhit
     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
@@ -986,16 +997,14 @@ Bool_t AliTRDdigitizer::ConvertHits(Int_t det, Float_t *hits, Int_t nhit
           
       // 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.
@@ -1162,6 +1171,8 @@ Bool_t AliTRDdigitizer::ConvertSignals(Int_t det, AliTRDarraySignal *signals)
     if (!Signal2ADC(det,signals)) {
       return kFALSE;
     }
+    // Run digital processing for digits
+    RunDigitalProcessing(det);
   }
 
   // Compress the arrays
@@ -1217,14 +1228,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 = calibration->GetNumberOfTimeBins();
+  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));
@@ -1261,13 +1279,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)) {
+      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));
       }
@@ -1308,9 +1327,6 @@ Bool_t AliTRDdigitizer::Signal2ADC(Int_t det, AliTRDarraySignal *signals)
     } // for: col
   } // for: row
 
-  // Run the digital processing in the MCM
-  RunDigitalProcessing(digits, det);
-
   return kTRUE;
 
 }
@@ -1336,10 +1352,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 = calibration->GetNumberOfTimeBins();
+  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;
@@ -1364,6 +1379,147 @@ Bool_t AliTRDdigitizer::Signal2SDigits(Int_t det, AliTRDarraySignal *signals)
 
 }
 
+//_____________________________________________________________________________
+Bool_t AliTRDdigitizer::Digits2SDigits(AliTRDdigitsManager * const manDig
+                                     , AliTRDdigitsManager * const manSDig)
+{
+  //
+  // Converts digits into s-digits. Needed for embedding into real data.
+  //
+
+  AliDebug(1,"Start converting digits to s-digits");
+
+  if (!fGeo) {
+    fGeo = new AliTRDgeometry();
+  }
+
+  AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
+  if (!calibration) {
+    AliFatal("Could not get calibration object");
+    return kFALSE;
+  }
+
+  AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
+  if (!simParam) {
+    AliFatal("Could not get simulation parameters");
+    return kFALSE;
+  }  
+
+  // Converts number of electrons to fC
+  const Double_t kEl2fC = 1.602e-19 * 1.0e15; 
+
+  // Coupling factor
+  Double_t coupling     = simParam->GetPadCoupling() 
+                        * simParam->GetTimeCoupling();
+  // Electronics conversion factor
+  Double_t convert      = kEl2fC 
+                        * simParam->GetChipGain();
+  // ADC conversion factor
+  Double_t adcConvert   = simParam->GetADCoutRange()
+                        / simParam->GetADCinRange();
+  // The electronics baseline in mV
+  Double_t baseline     = simParam->GetADCbaseline() 
+                        / adcConvert;
+  // The electronics baseline in electrons
+  //Double_t baselineEl   = baseline
+  //                      / 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;
+
+  Int_t row  = 0;
+  Int_t col  = 0;
+  Int_t time = 0;
+
+  for (Int_t det = 0; det < AliTRDgeometry::Ndet(); det++) {
+
+    Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
+    Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
+    Int_t nTimeTotal = manDig->GetDigitsParam()->GetNTimeBins(det);
+
+    // Get the calibration objects
+    //calGainFactorROC      = calibration->GetGainFactorROC(det);
+    //calGainFactorDetValue = calGainFactorDet->GetValue(det);
+
+    // Get the digits
+    AliTRDarrayADC *digits = (AliTRDarrayADC *) manDig->GetDigits(det);
+
+    if (!manSDig->HasSDigits()) {
+      AliError("SDigits manager has no s-digits");
+      return kFALSE;
+    }
+    // Get the s-digits
+    AliTRDarraySignal     *sdigits = (AliTRDarraySignal *)     manSDig->GetSDigits(det);
+    AliTRDarrayDictionary *tracks0 = (AliTRDarrayDictionary *) manSDig->GetDictionary(det,0);
+    AliTRDarrayDictionary *tracks1 = (AliTRDarrayDictionary *) manSDig->GetDictionary(det,1);
+    AliTRDarrayDictionary *tracks2 = (AliTRDarrayDictionary *) manSDig->GetDictionary(det,2);
+    // Allocate memory space for the digits buffer
+    sdigits->Allocate(nRowMax,nColMax,nTimeTotal);
+    tracks0->Allocate(nRowMax,nColMax,nTimeTotal);
+    tracks1->Allocate(nRowMax,nColMax,nTimeTotal);
+    tracks2->Allocate(nRowMax,nColMax,nTimeTotal);
+
+    // Keep the digits param
+    manSDig->GetDigitsParam()->SetNTimeBinsAll(manDig->GetDigitsParam()->GetNTimeBins(0));
+    manSDig->GetDigitsParam()->SetADCbaselineAll(manDig->GetDigitsParam()->GetADCbaseline(0));
+
+    if (digits->HasData()) {
+
+      digits->Expand();
+
+      // Create the sdigits for this chamber
+      for (row  = 0; row  <  nRowMax; row++ ) {
+        for (col  = 0; col  <  nColMax; col++ ) {
+
+          // The gain factors
+          //Float_t padgain = calGainFactorDetValue 
+          //                * calGainFactorROC->GetValue(col,row);
+
+          for (time = 0; time < nTimeTotal; time++) {
+
+            Short_t  adcVal = digits->GetData(row,col,time);
+            Double_t signal = (Double_t) adcVal;
+            // ADC -> signal in mV
+            signal /= adcConvert;
+            // Subtract baseline in mV
+            signal -= baseline;
+            // Signal in mV -> signal in #electrons
+            signal /= convert;
+            // Gain factor
+            //signal /= padgain; // Not needed for real data
+            // Pad and time coupling
+            signal /= coupling;
+
+            sdigits->SetData(row,col,time,signal);
+            tracks0->SetData(row,col,time,0);
+            tracks1->SetData(row,col,time,0);
+            tracks2->SetData(row,col,time,0);
+
+          } // for: time
+
+        } // for: col
+      } // for: row
+  
+    } // if: has data
+
+    sdigits->Compress(0);
+    tracks0->Compress();
+    tracks1->Compress();
+    tracks2->Compress();
+
+    // No compress just remove
+    manDig->RemoveDigits(det);
+    manDig->RemoveDictionaries(det);      
+
+  } // for: det
+
+  return kTRUE;
+
+}
+
 //_____________________________________________________________________________
 Bool_t AliTRDdigitizer::SDigits2Digits()
 {
@@ -1430,16 +1586,24 @@ Bool_t AliTRDdigitizer::MergeSDigits()
     AliDebug(1,"Only one input file.");
   }
   
-  Int_t nTimeTotal = calibration->GetNumberOfTimeBins();  
   Int_t iMerge = 0;
 
   while (mergeSDigitsManager) {
 
     iMerge++;
-      
+
     // Loop through the detectors
     for (Int_t iDet = 0; iDet < AliTRDgeometry::Ndet(); iDet++) {
 
+      Int_t nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(iDet);
+      if (mergeSDigitsManager->GetDigitsParam()->GetNTimeBins(iDet) != nTimeTotal) {
+        AliError(Form("Mismatch in the number of time bins [%d,%d] in detector %d"
+                     ,nTimeTotal
+                    ,mergeSDigitsManager->GetDigitsParam()->GetNTimeBins(iDet)
+                    ,iDet));
+        return kFALSE;
+      }
+
       Int_t nRowMax = fGeo->GetPadPlane(iDet)->GetNrows();
       Int_t nColMax = fGeo->GetPadPlane(iDet)->GetNcols();
          
@@ -1559,11 +1723,27 @@ Bool_t AliTRDdigitizer::ConvertSDigits()
     fSDigitsManager->RemoveDigits(det);
     fSDigitsManager->RemoveDictionaries(det);
 
+    // Run digital processing
+    RunDigitalProcessing(det);
+
     // Compress the arrays
     CompressOutputArrays(det);
 
   } // for: detector numbers
 
+  if (AliDataLoader *trklLoader = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets")) {
+    if (trklLoader->Tree())
+      trklLoader->WriteData("OVERWRITE");
+  }
+
+  // Save the values for the raw data headers
+  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;
 
@@ -1591,7 +1771,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 = calibration->GetNumberOfTimeBins();
+  Int_t nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(det);
 
   Int_t row  = 0;
   Int_t col  = 0;
@@ -1713,6 +1893,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)
 {
   //
@@ -1729,8 +1910,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;
 
@@ -1742,42 +1929,9 @@ 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(AliTRDarrayADC *digits, Int_t det)
+void AliTRDdigitizer::RunDigitalProcessing(Int_t det)
 {
   //
   // Run the digital processing in the TRAP
@@ -1785,23 +1939,26 @@ void AliTRDdigitizer::RunDigitalProcessing(AliTRDarrayADC *digits, Int_t det)
 
   AliTRDfeeParam *feeParam = AliTRDfeeParam::Instance();
 
-  //Create and initialize the mcm object 
-  AliTRDmcmSim* mcmfast = new AliTRDmcmSim(); 
+  AliTRDarrayADC *digits = fDigitsManager->GetDigits(det);
+  if (!digits)
+    return;
 
   //Call the methods in the mcm class using the temporary array as input  
-  for(Int_t rob = 0; rob < digits->GetNrow() / 4; rob++)
-  {
-    for(Int_t mcm = 0; mcm < 16; mcm++)
-    {
-      mcmfast->Init(det, rob, mcm); 
-      mcmfast->SetData(digits);
-      mcmfast->Filter();
-      if (feeParam->GetTracklet())
-        mcmfast->Tracklet();
-      mcmfast->ZSMapping();
-      mcmfast->WriteData(digits);
+  // 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);
+      }
     }
   }
-
-  delete mcmfast;
 }
+