LCOV - code coverage report
Current view: top level - TRD/TRDsim - AliTRDdigitizer.cxx (source / functions) Hit Total Coverage
Test: coverage.info Lines: 544 906 60.0 %
Date: 2016-06-14 17:26:59 Functions: 28 39 71.8 %

          Line data    Source code
       1             : /**************************************************************************
       2             :  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
       3             :  *                                                                        *
       4             :  * Author: The ALICE Off-line Project.                                    *
       5             :  * Contributors are mentioned in the code where appropriate.              *
       6             :  *                                                                        *
       7             :  * Permission to use, copy, modify and distribute this software and its   *
       8             :  * documentation strictly for non-commercial purposes is hereby granted   *
       9             :  * without fee, provided that the above copyright notice appears in all   *
      10             :  * copies and that both the copyright notice and this permission notice   *
      11             :  * appear in the supporting documentation. The authors make no claims     *
      12             :  * about the suitability of this software for any purpose. It is          *
      13             :  * provided "as is" without express or implied warranty.                  *
      14             :  **************************************************************************/
      15             : 
      16             : /* $Id$ */
      17             : 
      18             : ////////////////////////////////////////////////////////////////////////////
      19             : //                                                                        //
      20             : //  Creates and handles digits from TRD hits                              //
      21             : //                                                                        //
      22             : //  Authors: C. Blume (blume@ikf.uni-frankfurt.de)                        //
      23             : //           C. Lippmann                                                  //
      24             : //           B. Vulpescu                                                  //
      25             : //                                                                        //
      26             : //  The following effects are included:                                   //
      27             : //      - Diffusion                                                       //
      28             : //      - ExB effects                                                     //
      29             : //      - Gas gain including fluctuations                                 //
      30             : //      - Pad-response (simple Gaussian approximation)                    //
      31             : //      - Time-response                                                   //
      32             : //      - Electronics noise                                               //
      33             : //      - Electronics gain                                                //
      34             : //      - Digitization                                                    //
      35             : //      - Zero suppression                                                //
      36             : //                                                                        //
      37             : ////////////////////////////////////////////////////////////////////////////
      38             : 
      39             : #include <TGeoManager.h>
      40             : #include <TList.h>
      41             : #include <TMath.h>
      42             : #include <TRandom.h>
      43             : #include <TTree.h>
      44             : 
      45             : #include "AliRun.h"
      46             : #include "AliMC.h"
      47             : #include "AliRunLoader.h"
      48             : #include "AliLoader.h"
      49             : #include "AliConfig.h"
      50             : #include "AliDigitizationInput.h"
      51             : #include "AliRunLoader.h"
      52             : #include "AliLoader.h"
      53             : #include "AliLog.h"
      54             : 
      55             : #include "AliTRD.h"
      56             : #include "AliTRDhit.h"
      57             : #include "AliTRDdigitizer.h"
      58             : #include "AliTRDarrayDictionary.h"
      59             : #include "AliTRDarrayADC.h"
      60             : #include "AliTRDarraySignal.h"
      61             : #include "AliTRDdigitsManager.h"
      62             : #include "AliTRDgeometry.h"
      63             : #include "AliTRDpadPlane.h"
      64             : #include "AliTRDcalibDB.h"
      65             : #include "AliTRDSimParam.h"
      66             : #include "AliTRDCommonParam.h"
      67             : #include "AliTRDfeeParam.h"
      68             : #include "AliTRDmcmSim.h"
      69             : #include "AliTRDdigitsParam.h"
      70             : 
      71             : #include "AliTRDCalROC.h"
      72             : #include "AliTRDCalDet.h"
      73             : #include "AliTRDCalOnlineGainTableROC.h"
      74             : 
      75          12 : ClassImp(AliTRDdigitizer)
      76             : 
      77             : //_____________________________________________________________________________
      78             : AliTRDdigitizer::AliTRDdigitizer()
      79           0 :   :AliDigitizer()
      80           0 :   ,fRunLoader(0)
      81           0 :   ,fDigitsManager(0)
      82           0 :   ,fSDigitsManager(0)
      83           0 :   ,fSDigitsManagerList(0)
      84           0 :   ,fTRD(0)
      85           0 :   ,fGeo(0)
      86           0 :   ,fMcmSim(new AliTRDmcmSim)
      87           0 :   ,fEvent(0)
      88           0 :   ,fMasks(0)
      89           0 :   ,fCompress(kTRUE)
      90           0 :   ,fSDigits(kFALSE)
      91           0 :   ,fMergeSignalOnly(kFALSE)
      92           0 : {
      93             :   //
      94             :   // AliTRDdigitizer default constructor
      95             :   //
      96             :   
      97           0 : }
      98             : 
      99             : //_____________________________________________________________________________
     100             : AliTRDdigitizer::AliTRDdigitizer(const Text_t *name, const Text_t *title)              
     101           1 :   :AliDigitizer(name,title)
     102           1 :   ,fRunLoader(0)
     103           1 :   ,fDigitsManager(0)
     104           1 :   ,fSDigitsManager(0)
     105           1 :   ,fSDigitsManagerList(0)
     106           1 :   ,fTRD(0)
     107           1 :   ,fGeo(0)
     108           3 :   ,fMcmSim(new AliTRDmcmSim)
     109           1 :   ,fEvent(0)
     110           1 :   ,fMasks(0)
     111           1 :   ,fCompress(kTRUE)
     112           1 :   ,fSDigits(kFALSE)
     113           1 :   ,fMergeSignalOnly(kFALSE)
     114           5 : {
     115             :   //
     116             :   // AliTRDdigitizer constructor
     117             :   //
     118             : 
     119           2 : }
     120             : 
     121             : //_____________________________________________________________________________
     122             : AliTRDdigitizer::AliTRDdigitizer(AliDigitizationInput* digInput
     123             :                                , const Text_t *name, const Text_t *title)
     124           0 :   :AliDigitizer(digInput,name,title)
     125           0 :   ,fRunLoader(0)
     126           0 :   ,fDigitsManager(0)
     127           0 :   ,fSDigitsManager(0)
     128           0 :   ,fSDigitsManagerList(0)
     129           0 :   ,fTRD(0)
     130           0 :   ,fGeo(0)
     131           0 :   ,fMcmSim(new AliTRDmcmSim)
     132           0 :   ,fEvent(0)
     133           0 :   ,fMasks(0)
     134           0 :   ,fCompress(kTRUE)
     135           0 :   ,fSDigits(kFALSE)
     136           0 :   ,fMergeSignalOnly(kFALSE)
     137           0 : {
     138             :   //
     139             :   // AliTRDdigitizer constructor
     140             :   //
     141             : 
     142           0 : }
     143             : 
     144             : //_____________________________________________________________________________
     145             : AliTRDdigitizer::AliTRDdigitizer(AliDigitizationInput* digInput)
     146           1 :   :AliDigitizer(digInput,"AliTRDdigitizer","TRD digitizer")
     147           1 :   ,fRunLoader(0)
     148           1 :   ,fDigitsManager(0)
     149           1 :   ,fSDigitsManager(0)
     150           1 :   ,fSDigitsManagerList(0)
     151           1 :   ,fTRD(0)
     152           1 :   ,fGeo(0)
     153           3 :   ,fMcmSim(new AliTRDmcmSim)
     154           1 :   ,fEvent(0)
     155           1 :   ,fMasks(0)
     156           1 :   ,fCompress(kTRUE)
     157           1 :   ,fSDigits(kFALSE)
     158           1 :   ,fMergeSignalOnly(kFALSE)
     159           5 : {
     160             :   //
     161             :   // AliTRDdigitizer constructor
     162             :   //
     163             : 
     164           2 : }
     165             : 
     166             : //_____________________________________________________________________________
     167             : AliTRDdigitizer::AliTRDdigitizer(const AliTRDdigitizer &d)
     168           0 :   :AliDigitizer(d)
     169           0 :   ,fRunLoader(0)
     170           0 :   ,fDigitsManager(0)
     171           0 :   ,fSDigitsManager(0)
     172           0 :   ,fSDigitsManagerList(0)
     173           0 :   ,fTRD(0)
     174           0 :   ,fGeo(0)
     175           0 :   ,fMcmSim(new AliTRDmcmSim)
     176           0 :   ,fEvent(0)
     177           0 :   ,fMasks(0)
     178           0 :   ,fCompress(d.fCompress)
     179           0 :   ,fSDigits(d.fSDigits)
     180           0 :   ,fMergeSignalOnly(d.fMergeSignalOnly)
     181           0 : {
     182             :   //
     183             :   // AliTRDdigitizer copy constructor
     184             :   //
     185             : 
     186           0 : }
     187             : 
     188             : //_____________________________________________________________________________
     189             : AliTRDdigitizer::~AliTRDdigitizer()
     190          10 : {
     191             :   //
     192             :   // AliTRDdigitizer destructor
     193             :   //
     194             : 
     195           4 :   delete fDigitsManager;
     196           2 :   fDigitsManager      = 0;
     197             : 
     198             :   // s-digitsmanager will be deleted via list
     199           2 :   fSDigitsManager     = 0;
     200           2 :   if (fSDigitsManagerList) {
     201           2 :     fSDigitsManagerList->Delete();
     202           4 :     delete fSDigitsManagerList;
     203             :   }
     204           2 :   fSDigitsManagerList = 0;
     205             : 
     206           3 :   delete [] fMasks;
     207           2 :   fMasks       = 0;
     208             : 
     209           4 :   delete fMcmSim;
     210           2 :   fMcmSim = 0;
     211             : 
     212           4 :   delete fGeo;
     213           2 :   fGeo = 0;
     214             : 
     215           5 : }
     216             : 
     217             : //_____________________________________________________________________________
     218             : AliTRDdigitizer &AliTRDdigitizer::operator=(const AliTRDdigitizer &d)
     219             : {
     220             :   //
     221             :   // Assignment operator
     222             :   //
     223             : 
     224           0 :   if (this != &d) {
     225           0 :     ((AliTRDdigitizer &) d).Copy(*this);
     226           0 :   }
     227             : 
     228           0 :   return *this;
     229             : 
     230             : }
     231             : 
     232             : //_____________________________________________________________________________
     233             : void AliTRDdigitizer::Copy(TObject &d) const
     234             : {
     235             :   //
     236             :   // Copy function
     237             :   //
     238             : 
     239           0 :   ((AliTRDdigitizer &) d).fRunLoader          = 0;
     240           0 :   ((AliTRDdigitizer &) d).fDigitsManager      = 0;
     241           0 :   ((AliTRDdigitizer &) d).fSDigitsManager     = 0;
     242           0 :   ((AliTRDdigitizer &) d).fSDigitsManagerList = 0;
     243           0 :   ((AliTRDdigitizer &) d).fTRD                = 0;
     244           0 :   ((AliTRDdigitizer &) d).fGeo                = 0;
     245           0 :   ((AliTRDdigitizer &) d).fEvent              = 0;
     246           0 :   ((AliTRDdigitizer &) d).fMasks              = 0;
     247           0 :   ((AliTRDdigitizer &) d).fCompress           = fCompress;
     248           0 :   ((AliTRDdigitizer &) d).fSDigits            = fSDigits;
     249           0 :   ((AliTRDdigitizer &) d).fMergeSignalOnly    = fMergeSignalOnly;
     250             : 
     251           0 : }
     252             : 
     253             : //_____________________________________________________________________________
     254             : void AliTRDdigitizer::Digitize(const Option_t* option)
     255             : {
     256             :   //
     257             :   // Executes the merging
     258             :   //
     259             : 
     260             :   Int_t iInput;
     261             : 
     262             :   AliTRDdigitsManager *sdigitsManager;
     263             : 
     264           8 :   TString optionString = option;
     265           8 :   if (optionString.Contains("deb")) {
     266           0 :     AliLog::SetClassDebugLevel("AliTRDdigitizer",1);
     267           0 :     AliInfo("Called with debug option");
     268             :   }
     269             : 
     270             :   // The AliRoot file is already connected by the manager
     271             :   AliRunLoader *inrl = 0x0;
     272             :   
     273           4 :   if (gAlice) {
     274          20 :     AliDebug(1,"AliRun object found on file.");
     275             :   }
     276             :   else {
     277           0 :     inrl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
     278           0 :     inrl->LoadgAlice();
     279           0 :     gAlice = inrl->GetAliRun();
     280           0 :     if (!gAlice) {
     281           0 :       AliError("Could not find AliRun object.");
     282           0 :       return;
     283             :     }
     284             :   }
     285             :                                                                            
     286           4 :   Int_t nInput = fDigInput->GetNinputs();
     287           8 :   fMasks       = new Int_t[nInput];
     288          16 :   for (iInput = 0; iInput < nInput; iInput++) {
     289           4 :     fMasks[iInput] = fDigInput->GetMask(iInput);
     290             :   }
     291             : 
     292             :   //
     293             :   // Initialization
     294             :   //
     295             : 
     296           8 :   AliRunLoader *orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
     297             : 
     298           8 :   if (InitDetector()) {
     299             : 
     300           4 :     AliLoader *ogime = orl->GetLoader("TRDLoader");
     301             : 
     302             :     TTree *tree = 0;
     303           4 :     if (fSDigits) { 
     304             :       // If we produce SDigits
     305           0 :       tree = ogime->TreeS();
     306           0 :       if (!tree) {
     307           0 :         ogime->MakeTree("S");
     308           0 :         tree = ogime->TreeS();
     309           0 :       }
     310             :     }
     311             :     else {
     312             :       // If we produce Digits
     313           4 :       tree = ogime->TreeD();
     314           4 :       if (!tree) {
     315           4 :         ogime->MakeTree("D");
     316           4 :         tree = ogime->TreeD();
     317           4 :       }
     318             :     }
     319             : 
     320           4 :     MakeBranch(tree);
     321             : 
     322           4 :   }
     323             :  
     324          24 :   for (iInput = 0; iInput < nInput; iInput++) {
     325             : 
     326          28 :     AliDebug(1,Form("Add input stream %d",iInput));
     327             : 
     328             :     // Check if the input tree exists
     329          12 :     inrl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(iInput));
     330           4 :     AliLoader *gime = inrl->GetLoader("TRDLoader");
     331             : 
     332           4 :     TTree *treees = gime->TreeS();
     333           4 :     if (treees == 0x0) {
     334           2 :       if (gime->LoadSDigits()) {
     335           0 :         AliError(Form("Error Occured while loading S. Digits for input %d.",iInput));
     336           0 :         return;
     337             :       }
     338           1 :       treees = gime->TreeS();
     339           1 :     }
     340             :     
     341           4 :     if (treees == 0x0) {
     342           0 :       AliError(Form("Input stream %d does not exist",iInput));
     343           0 :       return;
     344             :     } 
     345             : 
     346             :     // Read the s-digits via digits manager
     347           8 :     sdigitsManager = new AliTRDdigitsManager();
     348           4 :     sdigitsManager->SetSDigits(kTRUE);
     349             :     
     350          12 :     AliRunLoader *rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(iInput));
     351           4 :     AliLoader *gimme = rl->GetLoader("TRDLoader");
     352           8 :     if (!gimme->TreeS()) 
     353             :       {
     354           0 :         gimme->LoadSDigits();
     355             :       }
     356             : 
     357           8 :     sdigitsManager->ReadDigits(gimme->TreeS());
     358             :    
     359             :     // Add the s-digits to the input list 
     360           4 :     AddSDigitsManager(sdigitsManager);
     361             : 
     362           4 :   }
     363             : 
     364             :   // Convert the s-digits to normal digits
     365          20 :   AliDebug(1,"Do the conversion");
     366           4 :   SDigits2Digits();
     367             : 
     368             :   // Store the digits
     369          20 :   AliDebug(1,"Write the digits");
     370           4 :   fDigitsManager->WriteDigits();
     371             : 
     372             :   // Write parameters
     373           4 :   orl->CdGAFile();
     374             : 
     375             :   // Clean up
     376           4 :   DeleteSDigitsManager();
     377             : 
     378          20 :   AliDebug(1,"Done");
     379             : 
     380           8 : }
     381             : 
     382             : //_____________________________________________________________________________
     383             : Bool_t AliTRDdigitizer::Open(const Char_t *file, Int_t nEvent)
     384             : {
     385             :   //
     386             :   // Opens a ROOT-file with TRD-hits and reads in the hit-tree
     387             :   //
     388             :   // Connect the AliRoot file containing Geometry, Kine, and Hits
     389             :   //  
     390             : 
     391           0 :   TString evfoldname = AliConfig::GetDefaultEventFolderName();
     392             : 
     393           0 :   fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
     394           0 :   if (!fRunLoader) {
     395           0 :     fRunLoader = AliRunLoader::Open(file,evfoldname,"UPDATE");
     396           0 :   }  
     397           0 :   if (!fRunLoader) {
     398           0 :     AliError(Form("Can not open session for file %s.",file));
     399           0 :     return kFALSE;
     400             :   }
     401             :    
     402           0 :   if (!fRunLoader->GetAliRun()) {
     403           0 :     fRunLoader->LoadgAlice();
     404             :   }
     405           0 :   gAlice = fRunLoader->GetAliRun();
     406             :   
     407           0 :   if (gAlice) {
     408           0 :     AliDebug(1,"AliRun object found on file.");
     409             :   }
     410             :   else {
     411           0 :     AliError("Could not find AliRun object.");
     412           0 :     return kFALSE;
     413             :   }
     414             : 
     415           0 :   fEvent = nEvent;
     416             : 
     417           0 :   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
     418           0 :   if (!loader) {
     419           0 :     AliError("Can not get TRD loader from Run Loader");
     420           0 :     return kFALSE;
     421             :   }
     422             :   
     423           0 :   if (InitDetector()) {
     424             :     TTree *tree = 0;
     425           0 :     if (fSDigits) { 
     426             :       // If we produce SDigits
     427           0 :       tree = loader->TreeS();
     428           0 :       if (!tree) {
     429           0 :         loader->MakeTree("S");
     430           0 :         tree = loader->TreeS();
     431           0 :       }
     432             :     }
     433             :     else {
     434             :       // If we produce Digits
     435           0 :       tree = loader->TreeD();
     436           0 :       if (!tree) {
     437           0 :         loader->MakeTree("D");
     438           0 :         tree = loader->TreeD();
     439           0 :       }
     440             :     }
     441           0 :     return MakeBranch(tree);
     442             :   }
     443             :   else {
     444           0 :     return kFALSE;
     445             :   }
     446             : 
     447           0 : }
     448             : 
     449             : //_____________________________________________________________________________
     450             : Bool_t AliTRDdigitizer::Open(AliRunLoader * const runLoader, Int_t nEvent)
     451             : {
     452             :   //
     453             :   // Opens a ROOT-file with TRD-hits and reads in the hit-tree
     454             :   //
     455             :   // Connect the AliRoot file containing Geometry, Kine, and Hits
     456             :   //  
     457             : 
     458           8 :   fRunLoader = runLoader;
     459           4 :   if (!fRunLoader) {
     460           0 :     AliError("RunLoader does not exist");
     461           0 :     return kFALSE;
     462             :   }
     463             :    
     464           4 :   if (!fRunLoader->GetAliRun()) {
     465           0 :     fRunLoader->LoadgAlice();
     466           0 :   }
     467           4 :   gAlice = fRunLoader->GetAliRun();
     468             :   
     469           4 :   if (gAlice) {
     470          12 :     AliDebug(1,"AliRun object found on file.");
     471             :   }
     472             :   else {
     473           0 :     AliError("Could not find AliRun object.");
     474           0 :     return kFALSE;
     475             :   }
     476             : 
     477           4 :   fEvent = nEvent;
     478             : 
     479           4 :   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
     480           4 :   if (!loader) {
     481           0 :     AliError("Can not get TRD loader from Run Loader");
     482           0 :     return kFALSE;
     483             :   }
     484             :   
     485           4 :   if (InitDetector()) {
     486             :     TTree *tree = 0;
     487           4 :     if (fSDigits) { 
     488             :       // If we produce SDigits
     489           4 :       tree = loader->TreeS();
     490           4 :       if (!tree) {
     491           4 :         loader->MakeTree("S");
     492           4 :         tree = loader->TreeS();
     493           4 :       }
     494             :     }
     495             :     else {
     496             :       // If we produce Digits
     497           0 :       tree = loader->TreeD();
     498           0 :       if (!tree) {
     499           0 :         loader->MakeTree("D");
     500           0 :         tree = loader->TreeD();
     501           0 :       }
     502             :     }
     503           4 :     return MakeBranch(tree);
     504             :   }
     505             :   else {
     506           0 :     return kFALSE;
     507             :   }
     508             : 
     509           4 : }
     510             : 
     511             : //_____________________________________________________________________________
     512             : Bool_t AliTRDdigitizer::InitDetector()
     513             : {
     514             :   //
     515             :   // Sets the pointer to the TRD detector and the geometry
     516             :   //
     517             : 
     518             :   // Get the pointer to the detector class and check for version 1
     519          18 :   fTRD = (AliTRD *) gAlice->GetDetector("TRD");
     520           9 :   if (!fTRD) {
     521           0 :     AliFatal("No TRD module found");
     522           0 :     exit(1);
     523             :   }
     524           9 :   if (fTRD->IsVersion() != 1) {
     525           0 :     AliFatal("TRD must be version 1 (slow simulator)");
     526           0 :     exit(1);
     527             :   }
     528             : 
     529             :   // Get the geometry
     530          18 :   fGeo = new AliTRDgeometry();
     531             : 
     532             :   // Create a digits manager
     533           9 :   if (fDigitsManager) {
     534          14 :     delete fDigitsManager;
     535             :   }
     536          18 :   fDigitsManager = new AliTRDdigitsManager();
     537           9 :   fDigitsManager->SetSDigits(fSDigits);
     538           9 :   fDigitsManager->CreateArrays();
     539           9 :   fDigitsManager->SetEvent(fEvent);
     540             : 
     541             :   // The list for the input s-digits manager to be merged
     542           9 :   if (fSDigitsManagerList) {
     543           7 :     fSDigitsManagerList->Delete();
     544           7 :   } 
     545             :   else {
     546           4 :     fSDigitsManagerList = new TList();
     547             :   }
     548             : 
     549           9 :   return kTRUE;
     550             : 
     551           0 : }
     552             : 
     553             : //_____________________________________________________________________________
     554             : Bool_t AliTRDdigitizer::MakeBranch(TTree *tree) const
     555             : {
     556             :   // 
     557             :   // Create the branches for the digits array
     558             :   //
     559             : 
     560          16 :   return fDigitsManager->MakeBranch(tree);
     561             : 
     562             : }
     563             : 
     564             : //_____________________________________________________________________________
     565             : void AliTRDdigitizer::AddSDigitsManager(AliTRDdigitsManager *man)
     566             : {
     567             :   //
     568             :   // Add a digits manager for s-digits to the input list.
     569             :   //
     570             : 
     571           8 :   fSDigitsManagerList->Add(man);
     572             : 
     573           4 : }
     574             : 
     575             : //_____________________________________________________________________________
     576             : void AliTRDdigitizer::DeleteSDigitsManager()
     577             : {
     578             :   //
     579             :   // Removes digits manager from the input list.
     580             :   //
     581             : 
     582           8 :   fSDigitsManagerList->Delete();
     583             : 
     584           4 : }
     585             : 
     586             : //_____________________________________________________________________________
     587             : Bool_t AliTRDdigitizer::MakeDigits()
     588             : {
     589             :   //
     590             :   // Creates digits.
     591             :   //
     592             : 
     593          16 :   AliDebug(1,"Start creating digits");
     594             : 
     595           4 :   if (!fGeo) {
     596           0 :     AliError("No geometry defined");
     597           0 :     return kFALSE;
     598             :   }
     599             : 
     600           4 :   AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
     601           4 :   if (!calibration) {
     602           0 :     AliFatal("Could not get calibration object");
     603           0 :     return kFALSE;
     604             :   }
     605             : 
     606           4 :   const Int_t kNdet  = AliTRDgeometry::Ndet();
     607             : 
     608           4 :   Float_t **hits = new Float_t*[kNdet];
     609           4 :   Int_t    *nhit = new Int_t[kNdet];
     610           4 :   memset(nhit,0,kNdet*sizeof(Int_t));
     611             : 
     612             :   AliTRDarraySignal *signals = 0x0;
     613             : 
     614             :   // Check the number of time bins from simParam against OCDB,
     615             :   // if OCDB value is not supposed to be used.
     616             :   // As default, the value from OCDB is taken
     617           4 :   if (AliTRDSimParam::Instance()->GetNTBoverwriteOCDB()) {
     618           0 :     if (calibration->GetNumberOfTimeBinsDCS() != AliTRDSimParam::Instance()->GetNTimeBins()) {
     619           0 :       AliWarning(Form("Number of time bins is different to OCDB value [SIM=%d, OCDB=%d]"
     620             :                      ,AliTRDSimParam::Instance()->GetNTimeBins()
     621             :                      ,calibration->GetNumberOfTimeBinsDCS()));
     622           0 :     }
     623             :     // Save the values for the raw data headers
     624           0 :     fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDSimParam::Instance()->GetNTimeBins());
     625           0 :   }
     626             :   else {
     627             :     // Get the OCDB values
     628           4 :     Int_t nTB = calibration->GetNumberOfTimeBinsDCS();
     629           4 :     if (nTB < 0) { // Currently -1 gets returned for "undefined" and "mixed",
     630             :                    // one might go back to -1 undefined and -2 mixed?
     631           0 :       AliError("No useful DCS information available for this run! Using standard values.");
     632             :       // // We fall back to the standard OCDB object, 
     633             :       // // cache the current run number..
     634             :       // Long64_t run = calibration->GetRun();
     635             :       // calibration->SetRun(0);
     636             :       // nTB = calibration->GetNumberOfTimeBinsDCS();
     637             :       // // ..to set it again
     638             :       // calibration->SetRun(run);
     639             :       // // If there's no standard OCDB object, we can still fail
     640             :       // if (nTB < 0) {
     641             :       //        AliFatal("No standard object found in the OCDB!");
     642             :       // }
     643           0 :       nTB = AliTRDSimParam::Instance()->GetNTimeBins();
     644           0 :     }
     645             :     // Save the values for the raw data headers
     646           4 :     fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(nTB);
     647             :   }
     648             : 
     649             :   // Save the values for the raw data headers
     650           4 :   fDigitsManager->GetDigitsParam()->SetADCbaselineAll(AliTRDSimParam::Instance()->GetADCbaseline());
     651             :  
     652             :   // Sort all hits according to detector number
     653           4 :   if (!SortHits(hits,nhit)) {
     654           0 :     AliError("Sorting hits failed");
     655           0 :     delete [] hits;
     656           0 :     delete [] nhit;
     657           0 :     return kFALSE;
     658             :   }
     659             : 
     660             :   // Loop through all detectors
     661        4332 :   for (Int_t det = 0; det < kNdet; det++) {
     662             : 
     663             :     // Detectors that are switched off, not installed, etc.
     664        2996 :     if ((!calibration->IsChamberNoData(det))    &&
     665        2160 :         ( fGeo->ChamberInGeometry(det))         &&
     666         836 :         (nhit[det] > 0)) {
     667             : 
     668         188 :       signals = new AliTRDarraySignal();
     669             : 
     670             :       // Convert the hits of the current detector to detector signals
     671         188 :       if (!ConvertHits(det,hits[det],nhit[det],signals)) {
     672           0 :         AliError(Form("Conversion of hits failed for detector=%d",det));
     673           0 :         delete [] hits;
     674           0 :         delete [] nhit;
     675           0 :         delete signals;
     676             :         signals = 0x0;
     677           0 :         return kFALSE;
     678             :       }
     679             : 
     680             :       // Convert the detector signals to digits or s-digits
     681         188 :       if (!ConvertSignals(det,signals)) {
     682           0 :         AliError(Form("Conversion of signals failed for detector=%d",det));
     683           0 :         delete [] hits;
     684           0 :         delete [] nhit;
     685           0 :         delete signals;
     686             :         signals = 0x0;
     687           0 :         return kFALSE;
     688             :       }
     689             : 
     690             :       // Delete the signals array
     691         376 :       delete signals;
     692             :       signals = 0x0;
     693             : 
     694         188 :     } // if: detector status
     695             : 
     696        2348 :     delete [] hits[det];
     697             : 
     698             :   } // for: detector
     699             : 
     700           4 :   if (!fSDigits) {
     701           0 :     if (AliDataLoader *trklLoader 
     702           0 :           = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets")) {
     703           0 :       if (trklLoader->Tree())
     704           0 :         trklLoader->WriteData("OVERWRITE");
     705             :     }
     706           0 :   }
     707             :     
     708           8 :   delete [] hits;
     709           8 :   delete [] nhit;
     710             : 
     711           4 :   return kTRUE;
     712             : 
     713           4 : }
     714             : 
     715             : //_____________________________________________________________________________
     716             : Bool_t AliTRDdigitizer::SortHits(Float_t **hits, Int_t *nhit)
     717             : {
     718             :   //
     719             :   // Read all the hits and sorts them according to detector number
     720             :   // in the output array <hits>.
     721             :   //
     722             : 
     723          16 :   AliDebug(1,"Start sorting hits");
     724             : 
     725           4 :   const Int_t kNdet = AliTRDgeometry::Ndet();
     726             :   // Size of the hit vector
     727             :   const Int_t kNhit = 6;
     728             : 
     729             :   Float_t *xyz      = 0;
     730             :   Int_t    nhitTrk  = 0;
     731             : 
     732           4 :   Int_t   *lhit     = new Int_t[kNdet];
     733           4 :   memset(lhit,0,kNdet*sizeof(Int_t));
     734             : 
     735        4328 :   for (Int_t det = 0; det < kNdet; det++) {
     736        2160 :     hits[det] = 0x0;
     737             :   }
     738             : 
     739           4 :   AliLoader *gimme = fRunLoader->GetLoader("TRDLoader");
     740           4 :   if (!gimme->TreeH()) {
     741           0 :     gimme->LoadHits();
     742           0 :   }
     743           4 :   TTree *hitTree = gimme->TreeH();
     744           4 :   if (hitTree == 0x0) {
     745           0 :     AliError("Can not get TreeH");
     746           0 :     delete [] lhit;
     747           0 :     return kFALSE;
     748             :   }
     749           4 :   fTRD->SetTreeAddress();
     750             : 
     751             :   // Get the number of entries in the hit tree
     752             :   // (Number of primary particles creating a hit somewhere)
     753           4 :   Int_t nTrk = (Int_t) hitTree->GetEntries();
     754          12 :   AliDebug(1,Form("Found %d tracks",nTrk));
     755             : 
     756             :   // Loop through all the tracks in the tree
     757         232 :   for (Int_t iTrk = 0; iTrk < nTrk; iTrk++) {
     758             : 
     759         112 :     gAlice->GetMCApp()->ResetHits();
     760         112 :     hitTree->GetEvent(iTrk);
     761             : 
     762         112 :     if (!fTRD->Hits()) {
     763           0 :       AliError(Form("No hits array for track = %d",iTrk));
     764           0 :       continue;
     765             :     }
     766             : 
     767             :     // Number of hits for this track
     768         112 :     nhitTrk = fTRD->Hits()->GetEntriesFast();
     769             : 
     770             :     Int_t hitCnt = 0;
     771             :     // Loop through the TRD hits
     772         112 :     AliTRDhit *hit = (AliTRDhit *) fTRD->FirstHit(-1);
     773       46674 :     while (hit) {
     774             : 
     775       23225 :       hitCnt++;
     776             : 
     777             :       // Don't analyze test hits
     778       23225 :       if (((Int_t) hit->GetCharge()) != 0) {
     779             : 
     780       22224 :         Int_t   trk  = hit->Track();
     781       22224 :         Int_t   det  = hit->GetDetector();
     782       22224 :         Int_t   q    = hit->GetCharge();
     783       22224 :         Float_t x    = hit->X();
     784       22224 :         Float_t y    = hit->Y();
     785       22224 :         Float_t z    = hit->Z();
     786       22224 :         Float_t time = hit->GetTime();
     787             : 
     788       22224 :         if (nhit[det] == lhit[det]) {
     789             :           // Inititialization of new detector
     790         189 :           xyz        = new Float_t[kNhit*(nhitTrk+lhit[det])];
     791         189 :           if (hits[det]) {
     792           1 :             memcpy(xyz,hits[det],sizeof(Float_t)*kNhit*lhit[det]);
     793           2 :             delete [] hits[det];
     794             :           }
     795         189 :           lhit[det] += nhitTrk;
     796         189 :           hits[det]  = xyz;
     797         189 :         }
     798             :         else {
     799       22035 :           xyz        = hits[det];
     800             :         }
     801       22224 :         xyz[nhit[det]*kNhit+0] = x;
     802       22224 :         xyz[nhit[det]*kNhit+1] = y;
     803       22224 :         xyz[nhit[det]*kNhit+2] = z;
     804       22224 :         xyz[nhit[det]*kNhit+3] = q;
     805       22224 :         xyz[nhit[det]*kNhit+4] = trk;
     806       22224 :         xyz[nhit[det]*kNhit+5] = time;
     807       22224 :         nhit[det]++;
     808             : 
     809       22224 :       } // if: charge != 0
     810             : 
     811       23225 :       hit = (AliTRDhit *) fTRD->NextHit();   
     812             : 
     813             :     } // for: hits of one track
     814             : 
     815         112 :   } // for: tracks
     816             : 
     817           8 :   delete [] lhit;
     818             : 
     819             :   return kTRUE;
     820             : 
     821           4 : }
     822             : 
     823             : //_____________________________________________________________________________
     824             : Bool_t AliTRDdigitizer::ConvertHits(Int_t det
     825             :                                   , const Float_t * const hits
     826             :                                   , Int_t nhit
     827             :                                   , AliTRDarraySignal *signals)
     828             : {
     829             :   //
     830             :   // Converts the detectorwise sorted hits to detector signals
     831             :   //
     832             : 
     833         752 :   AliDebug(1,Form("Start converting hits for detector=%d (nhits=%d)",det,nhit));
     834             : 
     835             :   // Number of pads included in the pad response
     836             :   const Int_t kNpad      = 3;
     837             :   // Number of track dictionary arrays
     838             :   const Int_t kNdict     = AliTRDdigitsManager::kNDict;
     839             :   // Size of the hit vector
     840             :   const Int_t kNhit      = 6;
     841             : 
     842             :   // Width of the amplification region
     843         188 :   const Float_t kAmWidth = AliTRDgeometry::AmThick();
     844             :   // Width of the drift region
     845         188 :   const Float_t kDrWidth = AliTRDgeometry::DrThick();
     846             :   // Drift + amplification region 
     847         188 :   const Float_t kDrMin   =          - 0.5 * kAmWidth;
     848         188 :   const Float_t kDrMax   = kDrWidth + 0.5 * kAmWidth;
     849             :   
     850             :   Int_t    iPad          = 0;
     851             :   Int_t    dict          = 0;
     852             :   Int_t    timeBinTRFend = 1;
     853             : 
     854         188 :   Double_t pos[3];
     855         188 :   Double_t loc[3];
     856         188 :   Double_t padSignal[kNpad];
     857         188 :   Double_t signalOld[kNpad];
     858             : 
     859         188 :   AliTRDarrayDictionary *dictionary[kNdict];
     860             : 
     861         188 :   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
     862         188 :   AliTRDCommonParam *commonParam = AliTRDCommonParam::Instance();
     863         188 :   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
     864             : 
     865         188 :   if (!commonParam) {
     866           0 :     AliFatal("Could not get common parameterss");
     867           0 :     return kFALSE;
     868             :   }
     869         188 :   if (!simParam) {
     870           0 :     AliFatal("Could not get simulation parameters");
     871           0 :     return kFALSE;
     872             :   }  
     873         188 :   if (!calibration) {
     874           0 :     AliFatal("Could not get calibration object");  
     875           0 :     return kFALSE;
     876             :   }
     877             : 
     878             :   // Get the detector wise calibration objects
     879             :   AliTRDCalROC       *calVdriftROC      = 0;
     880             :   Float_t             calVdriftDetValue = 0.0;
     881         188 :   const AliTRDCalDet *calVdriftDet      = calibration->GetVdriftDet();  
     882             :   AliTRDCalROC       *calT0ROC          = 0;
     883             :   Float_t             calT0DetValue     = 0.0;
     884         188 :   const AliTRDCalDet *calT0Det          = calibration->GetT0Det();  
     885             :   Double_t            calExBDetValue    = 0.0;
     886         188 :   const AliTRDCalDet *calExBDet         = calibration->GetExBDet();
     887             : 
     888         188 :   if (simParam->TRFOn()) {
     889         376 :     timeBinTRFend = ((Int_t) (simParam->GetTRFhi() 
     890         376 :                   * commonParam->GetSamplingFrequency())) - 1;
     891         188 :   }
     892             : 
     893         188 :   Int_t   nTimeTotal   = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
     894         188 :   Float_t samplingRate = commonParam->GetSamplingFrequency();
     895         188 :   Float_t elAttachProp = simParam->GetElAttachProp() / 100.0; 
     896             : 
     897         188 :   AliTRDpadPlane *padPlane = fGeo->GetPadPlane(det);
     898         188 :   Int_t   layer   = fGeo->GetLayer(det);         //update
     899         188 :   Float_t row0    = padPlane->GetRow0ROC();
     900         188 :   Int_t   nRowMax = padPlane->GetNrows();
     901         188 :   Int_t   nColMax = padPlane->GetNcols();
     902             : 
     903             :   // Create a new array for the signals
     904         188 :   signals->Allocate(nRowMax,nColMax,nTimeTotal);
     905             : 
     906             :   // Create a new array for the dictionary
     907        1504 :   for (dict = 0; dict < kNdict; dict++) {       
     908         564 :     dictionary[dict] = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(det,dict);
     909         564 :     dictionary[dict]->Allocate(nRowMax,nColMax,nTimeTotal);
     910             :   }      
     911             : 
     912             :   // Loop through the hits in this detector
     913       44824 :   for (Int_t hit = 0; hit < nhit; hit++) {
     914             : 
     915       22224 :     pos[0]          = hits[hit*kNhit+0];
     916       22224 :     pos[1]          = hits[hit*kNhit+1];
     917       22224 :     pos[2]          = hits[hit*kNhit+2];
     918       22224 :     Float_t q       = hits[hit*kNhit+3];
     919       22224 :     Float_t hittime = hits[hit*kNhit+5];
     920       22224 :     Int_t   track   = ((Int_t) hits[hit*kNhit+4]);
     921             : 
     922             :     Int_t   inDrift = 1;
     923             : 
     924             :     // Find the current volume with the geo manager
     925       22224 :     gGeoManager->SetCurrentPoint(pos);
     926       22224 :     gGeoManager->FindNode();      
     927       22224 :     if (strstr(gGeoManager->GetPath(),"/UK")) {
     928             :       inDrift = 0;
     929        4200 :     }
     930             : 
     931             :     // Get the calibration objects
     932       22224 :     calVdriftROC      = calibration->GetVdriftROC(det);
     933       22224 :     calVdriftDetValue = calVdriftDet->GetValue(det);
     934       22224 :     calT0ROC          = calibration->GetT0ROC(det);
     935       22224 :     calT0DetValue     = calT0Det->GetValue(det);
     936       22224 :     calExBDetValue    = calExBDet->GetValue(det);
     937             : 
     938             :     // Go to the local coordinate system:
     939             :     // loc[0] - col  direction in amplification or driftvolume
     940             :     // loc[1] - row  direction in amplification or driftvolume
     941             :     // loc[2] - time direction in amplification or driftvolume
     942       22224 :     gGeoManager->MasterToLocal(pos,loc);
     943       22224 :     if (inDrift) {
     944             :       // Relative to middle of amplification region
     945       18024 :       loc[2] = loc[2] - kDrWidth/2.0 - kAmWidth/2.0;
     946       18024 :     } 
     947             : 
     948             :     // The driftlength [cm] (w/o diffusion yet !).
     949             :     // It is negative if the hit is between pad plane and anode wires.
     950       22224 :     Double_t driftlength = -1.0 * loc[2];
     951             : 
     952             :     // Stupid patch to take care of TR photons that are absorbed
     953             :     // outside the chamber volume. A real fix would actually need
     954             :     // a more clever implementation of the TR hit generation
     955       22224 :     if (q < 0.0) {
     956          38 :       if ((loc[1] < padPlane->GetRowEndROC()) ||
     957          19 :           (loc[1] > padPlane->GetRow0ROC())) {
     958           0 :         continue;
     959             :       }
     960          38 :       if ((driftlength < kDrMin) ||
     961          19 :           (driftlength > kDrMax)) {
     962           0 :         continue;
     963             :       }
     964             :     }
     965             : 
     966             :     // Get row and col of unsmeared electron to retrieve drift velocity
     967             :     // The pad row (z-direction)
     968       22224 :     Int_t    rowE         = padPlane->GetPadRowNumberROC(loc[1]);
     969       22224 :     if (rowE < 0) {
     970           0 :       continue;
     971             :     }
     972       22224 :     Double_t rowOffset    = padPlane->GetPadRowOffsetROC(rowE,loc[1]);
     973             :     // The pad column (rphi-direction)
     974       22224 :     Double_t offsetTilt   = padPlane->GetTiltOffset(rowOffset);
     975       22224 :     Int_t    colE         = padPlane->GetPadColNumber(loc[0]+offsetTilt);
     976       22224 :     if (colE < 0) {
     977           0 :       continue;   
     978             :     }
     979             :     Double_t colOffset    = 0.0;
     980             : 
     981             :     // Normalized drift length
     982       22224 :     Float_t  driftvelocity  = calVdriftDetValue * calVdriftROC->GetValue(colE,rowE);
     983       22224 :     Double_t absdriftlength = TMath::Abs(driftlength);
     984       22224 :     if (commonParam->ExBOn()) {
     985       22224 :       absdriftlength /= TMath::Sqrt(1.0 / (1.0 + calExBDetValue*calExBDetValue));
     986       22224 :     }
     987             : 
     988             :     // Loop over all electrons of this hit
     989             :     // TR photons produce hits with negative charge
     990       22224 :     Int_t nEl = ((Int_t) TMath::Abs(q));
     991     1194272 :     for (Int_t iEl = 0; iEl < nEl; iEl++) {
     992             : 
     993             :       // Now the real local coordinate system of the ROC
     994             :       // column direction: locC
     995             :       // row direction:    locR 
     996             :       // time direction:   locT
     997             :       // locR and locC are identical to the coordinates of the corresponding
     998             :       // volumina of the drift or amplification region.
     999             :       // locT is defined relative to the wire plane (i.e. middle of amplification
    1000             :       // region), meaning locT = 0, and is negative for hits coming from the
    1001             :       // drift region. 
    1002      574912 :       Double_t locC = loc[0];
    1003      574912 :       Double_t locR = loc[1];
    1004      574912 :       Double_t locT = loc[2];
    1005             : 
    1006             :       // Electron attachment
    1007      574912 :       if (simParam->ElAttachOn()) {
    1008           0 :         if (gRandom->Rndm() < (absdriftlength * elAttachProp)) {
    1009           0 :           continue;
    1010             :         }
    1011             :       }
    1012             :           
    1013             :       // Apply the diffusion smearing
    1014      574912 :       if (simParam->DiffusionOn()) {
    1015      574912 :         if (!(Diffusion(driftvelocity,absdriftlength,calExBDetValue,locR,locC,locT))) {
    1016           0 :           continue;
    1017             :         }
    1018             :       }
    1019             : 
    1020             :       // Apply E x B effects (depends on drift direction)
    1021      574912 :       if (commonParam->ExBOn()) {
    1022      574912 :         locC = locC + calExBDetValue * driftlength;
    1023      574912 :       }
    1024             : 
    1025             :       // The electron position after diffusion and ExB in pad coordinates.
    1026             :       // The pad row (z-direction)
    1027      574912 :       rowE       = padPlane->GetPadRowNumberROC(locR);
    1028      574920 :       if (rowE < 0) continue;
    1029      574904 :       rowOffset  = padPlane->GetPadRowOffsetROC(rowE,locR);
    1030             : 
    1031             :       // The pad column (rphi-direction)
    1032      574904 :       offsetTilt = padPlane->GetTiltOffset(rowOffset);
    1033      574904 :       colE       = padPlane->GetPadColNumber(locC+offsetTilt);
    1034      574904 :       if (colE < 0) continue;         
    1035      574904 :       colOffset  = padPlane->GetPadColOffset(colE,locC+offsetTilt);
    1036             : 
    1037             :       // Also re-retrieve drift velocity because col and row may have changed
    1038      574904 :       driftvelocity = calVdriftDetValue * calVdriftROC->GetValue(colE,rowE);
    1039      574904 :       Float_t t0    = calT0DetValue     + calT0ROC->GetValue(colE,rowE);
    1040             : 
    1041             :       // Convert the position to drift time [mus], using either constant drift velocity or
    1042             :       // time structure of drift cells (non-isochronity, GARFIELD calculation).
    1043             :       // Also add absolute time of hits to take pile-up events into account properly
    1044             :       Double_t drifttime;
    1045      574904 :       if (simParam->TimeStructOn()) {
    1046             :         // Get z-position with respect to anode wire
    1047      574904 :         Double_t zz  =  row0 - locR + padPlane->GetAnodeWireOffset();
    1048      574904 :         zz -= ((Int_t)(2 * zz)) / 2.0;
    1049      574904 :         if (zz > 0.25) {
    1050      288352 :           zz  = 0.5 - zz;
    1051      288352 :         }
    1052             :         // Use drift time map (GARFIELD)
    1053      574904 :         drifttime = commonParam->TimeStruct(driftvelocity,0.5*kAmWidth-1.0*locT,zz)
    1054      574904 :                   + hittime;
    1055      574904 :       } 
    1056             :       else {
    1057             :         // Use constant drift velocity
    1058           0 :         drifttime = TMath::Abs(locT) / driftvelocity
    1059           0 :                   + hittime;
    1060             :       }
    1061             : 
    1062             :       // Apply the gas gain including fluctuations
    1063             :       Double_t ggRndm = 0.0;
    1064      574904 :       do {
    1065      574904 :         ggRndm = gRandom->Rndm();
    1066      574904 :       } while (ggRndm <= 0);
    1067      574904 :       Double_t signal = -(simParam->GetGasGain()) * TMath::Log(ggRndm);
    1068             : 
    1069             :       // Apply the pad response 
    1070      574904 :       if (simParam->PRFOn()) {
    1071             :         // The distance of the electron to the center of the pad 
    1072             :         // in units of pad width
    1073      574904 :         Double_t dist = (colOffset - 0.5*padPlane->GetColSize(colE))
    1074      574904 :                       / padPlane->GetColSize(colE);
    1075             :         // This is a fixed parametrization, i.e. not dependent on
    1076             :         // calibration values !
    1077      574904 :         if (!(calibration->PadResponse(signal,dist,layer,padSignal))) continue;
    1078      574904 :       }
    1079             :       else {
    1080           0 :         padSignal[0] = 0.0;
    1081           0 :         padSignal[1] = signal;
    1082           0 :         padSignal[2] = 0.0;
    1083             :       }
    1084             : 
    1085             :       // The time bin (always positive), with t0 distortion
    1086      574904 :       Double_t timeBinIdeal = drifttime * samplingRate + t0;
    1087             :       // Protection 
    1088      574904 :       if (TMath::Abs(timeBinIdeal) > 2*nTimeTotal) {
    1089             :         timeBinIdeal = 2 * nTimeTotal;
    1090           0 :       }
    1091      574904 :       Int_t    timeBinTruncated = ((Int_t) timeBinIdeal);
    1092             :       // The distance of the position to the middle of the timebin
    1093      574904 :       Double_t timeOffset       = ((Float_t) timeBinTruncated 
    1094      574904 :                                 + 0.5 - timeBinIdeal) / samplingRate;
    1095             :           
    1096             :       // Sample the time response inside the drift region
    1097             :       // + additional time bins before and after.
    1098             :       // The sampling is done always in the middle of the time bin
    1099    21584996 :       for (Int_t iTimeBin = TMath::Max(timeBinTruncated,0)
    1100    10792498 :           ;iTimeBin < TMath::Min(timeBinTruncated+timeBinTRFend,nTimeTotal)
    1101    10217594 :           ;iTimeBin++) {
    1102             : 
    1103             :         // Apply the time response
    1104             :         Double_t timeResponse = 1.0;
    1105             :         Double_t crossTalk    = 0.0;
    1106    10217594 :         Double_t time         = (iTimeBin - timeBinTruncated) / samplingRate + timeOffset;
    1107             : 
    1108    10217594 :         if (simParam->TRFOn()) {
    1109    10217594 :           timeResponse = simParam->TimeResponse(time);
    1110    10217594 :         }
    1111    10217594 :         if (simParam->CTOn()) {
    1112    10217594 :           crossTalk    = simParam->CrossTalk(time);
    1113    10217594 :         }
    1114             : 
    1115    10217594 :         signalOld[0] = 0.0;
    1116    10217594 :         signalOld[1] = 0.0;
    1117    10217594 :         signalOld[2] = 0.0;
    1118             : 
    1119    81592430 :         for (iPad = 0; iPad < kNpad; iPad++) {
    1120             : 
    1121    30652782 :           Int_t colPos = colE + iPad - 1;
    1122    30690682 :           if (colPos <        0) continue;
    1123    30689043 :           if (colPos >= nColMax) break;
    1124             : 
    1125             :           // Add the signals
    1126    30540721 :           signalOld[iPad] = signals->GetData(rowE,colPos,iTimeBin);
    1127             : 
    1128    61081442 :           if (colPos != colE) {
    1129             :             // Cross talk added to non-central pads
    1130    50863848 :             signalOld[iPad] += padSignal[iPad] 
    1131    20323127 :                              * (timeResponse + crossTalk);
    1132    20323127 :           } 
    1133             :           else {
    1134             :             // W/o cross talk at central pad
    1135    10217594 :             signalOld[iPad] += padSignal[iPad] 
    1136    10217594 :                              * timeResponse;
    1137             :           }
    1138             : 
    1139    30540721 :           signals->SetData(rowE,colPos,iTimeBin,signalOld[iPad]);
    1140             : 
    1141             :           // Store the track index in the dictionary
    1142             :           // Note: We store index+1 in order to allow the array to be compressed
    1143             :           // Note2: Taking out the +1 in track
    1144    30540721 :           if (signalOld[iPad] > 0.0) { 
    1145    83070874 :             for (dict = 0; dict < kNdict; dict++) {
    1146    40709089 :               Int_t oldTrack = dictionary[dict]->GetData(rowE,colPos,iTimeBin); 
    1147    70351597 :               if (oldTrack == track) break;
    1148    11066581 :               if (oldTrack ==  -1 ) {
    1149       71865 :                 dictionary[dict]->SetData(rowE,colPos,iTimeBin,track);
    1150       71865 :                 break;
    1151             :               }
    1152    10994716 :             }
    1153             :           }
    1154             : 
    1155    30540721 :         } // Loop: pads
    1156             : 
    1157             :       } // Loop: time bins
    1158             : 
    1159     1149816 :     } // Loop: electrons of a single hit
    1160             : 
    1161       22224 :   } // Loop: hits
    1162             : 
    1163         564 :   AliDebug(2,Form("Finished analyzing %d hits",nhit));
    1164             : 
    1165             :   return kTRUE;
    1166             : 
    1167         188 : }
    1168             : 
    1169             : //_____________________________________________________________________________
    1170             : Bool_t AliTRDdigitizer::ConvertSignals(Int_t det, AliTRDarraySignal *signals)
    1171             : {
    1172             :   //
    1173             :   // Convert signals to digits
    1174             :   //
    1175             : 
    1176         752 :   AliDebug(1,Form("Start converting the signals for detector %d",det));
    1177             : 
    1178         188 :   if (fSDigits) {
    1179             :     // Convert the signal array to s-digits
    1180         188 :     if (!Signal2SDigits(det,signals)) {
    1181           0 :       return kFALSE;
    1182             :     }
    1183             :   }
    1184             :   else {
    1185             :     // Convert the signal array to digits
    1186           0 :     if (!Signal2ADC(det,signals)) {
    1187           0 :       return kFALSE;
    1188             :     }
    1189             :     // Run digital processing for digits
    1190           0 :     RunDigitalProcessing(det);
    1191             :   }
    1192             : 
    1193             :   // Compress the arrays
    1194         188 :   CompressOutputArrays(det);
    1195             : 
    1196         188 :   return kTRUE;
    1197             : 
    1198         188 : }
    1199             : 
    1200             : //_____________________________________________________________________________
    1201             : Bool_t AliTRDdigitizer::Signal2ADC(Int_t det, AliTRDarraySignal *signals)
    1202             : {
    1203             :   //
    1204             :   // Converts the sampled electron signals to ADC values for a given chamber
    1205             :   //
    1206             : 
    1207         752 :   AliDebug(1,Form("Start converting signals to ADC values for detector=%d",det));
    1208             : 
    1209         188 :   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
    1210         188 :   if (!calibration) {
    1211           0 :     AliFatal("Could not get calibration object");
    1212           0 :     return kFALSE;
    1213             :   }
    1214             : 
    1215         188 :   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
    1216         188 :   if (!simParam) {
    1217           0 :     AliFatal("Could not get simulation parameters");
    1218           0 :     return kFALSE;
    1219             :   }  
    1220             : 
    1221             :   // Converts number of electrons to fC
    1222             :   const Double_t kEl2fC = 1.602e-19 * 1.0e15; 
    1223             : 
    1224             :   // Coupling factor
    1225         376 :   Double_t coupling     = simParam->GetPadCoupling() 
    1226         188 :                         * simParam->GetTimeCoupling();
    1227             :   // Electronics conversion factor
    1228             :   Double_t convert      = kEl2fC 
    1229         188 :                         * simParam->GetChipGain();
    1230             :   // ADC conversion factor
    1231         376 :   Double_t adcConvert   = simParam->GetADCoutRange()
    1232         188 :                         / simParam->GetADCinRange();
    1233             :   // The electronics baseline in mV
    1234         188 :   Double_t baseline     = simParam->GetADCbaseline() 
    1235         188 :                         / adcConvert;
    1236             :   // The electronics baseline in electrons
    1237             :   Double_t baselineEl   = baseline
    1238         188 :                         / convert;
    1239             : 
    1240             :   Int_t row  = 0;
    1241             :   Int_t col  = 0;
    1242             :   Int_t time = 0;
    1243             : 
    1244         188 :   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
    1245         188 :   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
    1246         188 :   Int_t nTimeTotal = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
    1247         188 :   if (fSDigitsManager->GetDigitsParam()->GetNTimeBins(det)) {
    1248         188 :     nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(det);
    1249             :   }
    1250             :   else {
    1251           0 :     AliFatal("Could not get number of time bins");
    1252           0 :     return kFALSE;
    1253             :   }
    1254             : 
    1255             :   // The gain factor calibration objects
    1256         188 :   const AliTRDCalDet *calGainFactorDet      = calibration->GetGainFactorDet();  
    1257             :   AliTRDCalROC       *calGainFactorROC      = 0x0;
    1258             :   Float_t             calGainFactorDetValue = 0.0;
    1259             : 
    1260             :   AliTRDarrayADC     *digits                = 0x0;
    1261             : 
    1262         188 :   if (!signals) {
    1263           0 :     AliError(Form("Signals array for detector %d does not exist\n",det));
    1264           0 :     return kFALSE;
    1265             :   }
    1266         188 :   if (signals->HasData()) {
    1267             :     // Expand the container if neccessary
    1268         188 :     signals->Expand();   
    1269         188 :   }
    1270             :   else {
    1271             :     // Create missing containers
    1272           0 :     signals->Allocate(nRowMax,nColMax,nTimeTotal);      
    1273             :   }
    1274             : 
    1275             :   // Get the container for the digits of this detector
    1276         188 :   if (fDigitsManager->HasSDigits()) {
    1277           0 :     AliError("Digits manager has s-digits");
    1278           0 :     return kFALSE;
    1279             :   }
    1280             : 
    1281         188 :   digits = (AliTRDarrayADC *) fDigitsManager->GetDigits(det);
    1282             :   // Allocate memory space for the digits buffer
    1283         188 :   if (!digits->HasData()) {
    1284         188 :     digits->Allocate(nRowMax,nColMax,nTimeTotal);
    1285         188 :   }
    1286             : 
    1287             :   // Get the calibration objects
    1288         188 :   calGainFactorROC      = calibration->GetGainFactorROC(det);
    1289         188 :   calGainFactorDetValue = calGainFactorDet->GetValue(det);
    1290             : 
    1291             :   // Create the digits for this chamber
    1292        5576 :   for (row  = 0; row  <  nRowMax; row++ ) {
    1293      754000 :     for (col  = 0; col  <  nColMax; col++ ) {
    1294             : 
    1295             :       // halfchamber masking
    1296      374400 :       Int_t iMcm            = (Int_t)(col/18); // current group of 18 col pads
    1297      374400 :       Int_t halfchamberside = (iMcm>3 ? 1 : 0); // 0=Aside, 1=Bside
    1298             :       // Halfchambers that are switched off, masked by calibration
    1299      374400 :       if (calibration->IsHalfChamberNoData(det, halfchamberside)) 
    1300           0 :         continue;
    1301             : 
    1302             :       // Check whether pad is masked
    1303             :       // Bridged pads are not considered yet!!!
    1304      748800 :       if (calibration->IsPadMasked(det,col,row) || 
    1305      374400 :           calibration->IsPadNotConnected(det,col,row)) {
    1306           0 :         continue;
    1307             :       }
    1308             : 
    1309             :       // The gain factors
    1310             :       Float_t padgain    = calGainFactorDetValue 
    1311      374400 :                          * calGainFactorROC->GetValue(col,row);
    1312      374400 :       if (padgain <= 0) {
    1313           0 :         AliError(Form("Not a valid gain %f, %d %d %d",padgain,det,col,row));
    1314           0 :       }
    1315             : 
    1316    20966400 :       for (time = 0; time < nTimeTotal; time++) {
    1317             : 
    1318             :         // Get the signal amplitude
    1319    10108800 :         Float_t signalAmp = signals->GetData(row,col,time);
    1320             :         // Pad and time coupling
    1321    10108800 :         signalAmp *= coupling;
    1322             :         // Gain factors
    1323    10108800 :         signalAmp *= padgain;
    1324             : 
    1325             :         // Add the noise, starting from minus ADC baseline in electrons
    1326    20217600 :         signalAmp  = TMath::Max((Double_t) gRandom->Gaus(signalAmp,simParam->GetNoise())
    1327    10108800 :                                ,-baselineEl);
    1328             : 
    1329             :         // Convert to mV
    1330    10108800 :         signalAmp *= convert;
    1331             :         // Add ADC baseline in mV
    1332    10108800 :         signalAmp += baseline;
    1333             : 
    1334             :         // Convert to ADC counts. Set the overflow-bit fADCoutRange if the
    1335             :         // signal is larger than fADCinRange
    1336             :         Short_t adc  = 0;
    1337    10108800 :         if (signalAmp >= simParam->GetADCinRange()) {
    1338          14 :           adc = ((Short_t) simParam->GetADCoutRange());
    1339          14 :         }
    1340             :         else {
    1341    10108786 :           adc = TMath::Nint(signalAmp * adcConvert);
    1342             :         }
    1343             : 
    1344             :         // Saving all digits
    1345    10108800 :         digits->SetData(row,col,time,adc);
    1346             : 
    1347             :       } // for: time
    1348             : 
    1349      374400 :     } // for: col
    1350             :   } // for: row
    1351             : 
    1352         188 :   return kTRUE;
    1353             : 
    1354         188 : }
    1355             : 
    1356             : //_____________________________________________________________________________
    1357             : Bool_t AliTRDdigitizer::Signal2SDigits(Int_t det, AliTRDarraySignal *signals)
    1358             : {
    1359             :   //
    1360             :   // Converts the sampled electron signals to s-digits
    1361             :   //
    1362             : 
    1363         752 :   AliDebug(1,Form("Start converting signals to s-digits for detector=%d",det));
    1364             : 
    1365         188 :   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
    1366         188 :   if (!calibration) {
    1367           0 :     AliFatal("Could not get calibration object");
    1368           0 :     return kFALSE;
    1369             :   }
    1370             : 
    1371             :   Int_t row  = 0;
    1372             :   Int_t col  = 0;
    1373             :   Int_t time = 0;
    1374             : 
    1375         188 :   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
    1376         188 :   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
    1377         188 :   Int_t nTimeTotal = fDigitsManager->GetDigitsParam()->GetNTimeBins(det);
    1378             : 
    1379             :   // Get the container for the digits of this detector
    1380         188 :   if (!fDigitsManager->HasSDigits()) {
    1381           0 :     AliError("Digits manager has no s-digits");
    1382           0 :     return kFALSE;
    1383             :   }
    1384             : 
    1385         188 :   AliTRDarraySignal *digits = (AliTRDarraySignal *) fDigitsManager->GetSDigits(det);
    1386             :   // Allocate memory space for the digits buffer
    1387         188 :   if (!digits->HasData()) {
    1388         188 :     digits->Allocate(nRowMax,nColMax,nTimeTotal);
    1389         188 :   }
    1390             : 
    1391             :   // Create the sdigits for this chamber
    1392        5576 :   for (row  = 0; row  <  nRowMax; row++ ) {
    1393      754000 :     for (col  = 0; col  <  nColMax; col++ ) {
    1394             : 
    1395             :       // halfchamber masking
    1396      374400 :       Int_t iMcm            = (Int_t)(col/18); // current group of 18 col pads
    1397      374400 :       Int_t halfchamberside = (iMcm>3 ? 1 : 0); // 0=Aside, 1=Bside
    1398             :       // Halfchambers that are switched off, masked by calibration
    1399      374400 :       if (calibration->IsHalfChamberNoData(det, halfchamberside))
    1400           0 :         continue;
    1401             : 
    1402    20966400 :       for (time = 0; time < nTimeTotal; time++) {
    1403    10108800 :         digits->SetData(row,col,time,signals->GetData(row,col,time));
    1404             :       } // for: time
    1405      374400 :     } // for: col
    1406             :   } // for: row
    1407             : 
    1408             :   return kTRUE;
    1409             : 
    1410         188 : }
    1411             : 
    1412             : //_____________________________________________________________________________
    1413             : Bool_t AliTRDdigitizer::Digits2SDigits(AliTRDdigitsManager * const manDig
    1414             :                                      , AliTRDdigitsManager * const manSDig)
    1415             : {
    1416             :   //
    1417             :   // Converts digits into s-digits. Needed for embedding into real data.
    1418             :   //
    1419             : 
    1420           0 :   AliDebug(1,"Start converting digits to s-digits");
    1421             : 
    1422           0 :   if (!fGeo) {
    1423           0 :     fGeo = new AliTRDgeometry();
    1424           0 :   }
    1425             : 
    1426           0 :   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
    1427           0 :   if (!calibration) {
    1428           0 :     AliFatal("Could not get calibration object");
    1429           0 :     return kFALSE;
    1430             :   }
    1431             : 
    1432           0 :   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
    1433           0 :   if (!simParam) {
    1434           0 :     AliFatal("Could not get simulation parameters");
    1435           0 :     return kFALSE;
    1436             :   }  
    1437             : 
    1438             :   // Converts number of electrons to fC
    1439             :   const Double_t kEl2fC = 1.602e-19 * 1.0e15; 
    1440             : 
    1441             :   // Coupling factor
    1442           0 :   Double_t coupling     = simParam->GetPadCoupling() 
    1443           0 :                         * simParam->GetTimeCoupling();
    1444             :   // Electronics conversion factor
    1445             :   Double_t convert      = kEl2fC 
    1446           0 :                         * simParam->GetChipGain();
    1447             :   // ADC conversion factor
    1448           0 :   Double_t adcConvert   = simParam->GetADCoutRange()
    1449           0 :                         / simParam->GetADCinRange();
    1450             :   // The electronics baseline in mV
    1451           0 :   Double_t baseline     = simParam->GetADCbaseline() 
    1452           0 :                         / adcConvert;
    1453             :   // The electronics baseline in electrons
    1454             :   //Double_t baselineEl   = baseline
    1455             :   //                      / convert;
    1456             : 
    1457             :   // The gainfactor calibration objects
    1458             :   // Not used since these digits are supposed to be from real raw data
    1459             :   //const AliTRDCalDet *calGainFactorDet      = calibration->GetGainFactorDet();  
    1460             :   //AliTRDCalROC       *calGainFactorROC      = 0;
    1461             :   //Float_t             calGainFactorDetValue = 0.0;
    1462             : 
    1463             :   Int_t row  = 0;
    1464             :   Int_t col  = 0;
    1465             :   Int_t time = 0;
    1466             : 
    1467           0 :   for (Int_t det = 0; det < AliTRDgeometry::Ndet(); det++) {
    1468             : 
    1469           0 :     Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
    1470           0 :     Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
    1471           0 :     Int_t nTimeTotal = manDig->GetDigitsParam()->GetNTimeBins(det);
    1472             : 
    1473             :     // Get the calibration objects
    1474             :     //calGainFactorROC      = calibration->GetGainFactorROC(det);
    1475             :     //calGainFactorDetValue = calGainFactorDet->GetValue(det);
    1476             : 
    1477             :     // Get the digits
    1478           0 :     AliTRDarrayADC *digits = (AliTRDarrayADC *) manDig->GetDigits(det);
    1479             : 
    1480           0 :     if (!manSDig->HasSDigits()) {
    1481           0 :       AliError("SDigits manager has no s-digits");
    1482           0 :       return kFALSE;
    1483             :     }
    1484             :     // Get the s-digits
    1485           0 :     AliTRDarraySignal     *sdigits = (AliTRDarraySignal *)     manSDig->GetSDigits(det);
    1486           0 :     AliTRDarrayDictionary *tracks0 = (AliTRDarrayDictionary *) manSDig->GetDictionary(det,0);
    1487           0 :     AliTRDarrayDictionary *tracks1 = (AliTRDarrayDictionary *) manSDig->GetDictionary(det,1);
    1488           0 :     AliTRDarrayDictionary *tracks2 = (AliTRDarrayDictionary *) manSDig->GetDictionary(det,2);
    1489             :     // Allocate memory space for the digits buffer
    1490           0 :     sdigits->Allocate(nRowMax,nColMax,nTimeTotal);
    1491           0 :     tracks0->Allocate(nRowMax,nColMax,nTimeTotal);
    1492           0 :     tracks1->Allocate(nRowMax,nColMax,nTimeTotal);
    1493           0 :     tracks2->Allocate(nRowMax,nColMax,nTimeTotal);
    1494             : 
    1495             :     // Keep the digits param
    1496           0 :     manSDig->GetDigitsParam()->SetNTimeBinsAll(manDig->GetDigitsParam()->GetNTimeBins(0));
    1497           0 :     manSDig->GetDigitsParam()->SetADCbaselineAll(manDig->GetDigitsParam()->GetADCbaseline(0));
    1498             : 
    1499           0 :     if (digits->HasData()) {
    1500             : 
    1501           0 :       digits->Expand();
    1502             : 
    1503             :       // Create the sdigits for this chamber
    1504           0 :       for (row  = 0; row  <  nRowMax; row++ ) {
    1505           0 :         for (col  = 0; col  <  nColMax; col++ ) {
    1506             : 
    1507             :           // The gain factors
    1508             :           //Float_t padgain = calGainFactorDetValue 
    1509             :           //                * calGainFactorROC->GetValue(col,row);
    1510             : 
    1511           0 :           for (time = 0; time < nTimeTotal; time++) {
    1512             : 
    1513           0 :             Short_t  adcVal = digits->GetData(row,col,time);
    1514           0 :             Double_t signal = (Double_t) adcVal;
    1515             :             // ADC -> signal in mV
    1516           0 :             signal /= adcConvert;
    1517             :             // Subtract baseline in mV
    1518           0 :             signal -= baseline;
    1519             :             // Signal in mV -> signal in #electrons
    1520           0 :             signal /= convert;
    1521             :             // Gain factor
    1522             :             //signal /= padgain; // Not needed for real data
    1523             :             // Pad and time coupling
    1524           0 :             signal /= coupling;
    1525             : 
    1526           0 :             sdigits->SetData(row,col,time,signal);
    1527           0 :             tracks0->SetData(row,col,time,0);
    1528           0 :             tracks1->SetData(row,col,time,0);
    1529           0 :             tracks2->SetData(row,col,time,0);
    1530             : 
    1531             :           } // for: time
    1532             : 
    1533             :         } // for: col
    1534             :       } // for: row
    1535             :   
    1536             :     } // if: has data
    1537             : 
    1538           0 :     sdigits->Compress(0);
    1539           0 :     tracks0->Compress();
    1540           0 :     tracks1->Compress();
    1541           0 :     tracks2->Compress();
    1542             : 
    1543             :     // No compress just remove
    1544           0 :     manDig->RemoveDigits(det);
    1545           0 :     manDig->RemoveDictionaries(det);      
    1546             : 
    1547           0 :   } // for: det
    1548             : 
    1549           0 :   return kTRUE;
    1550             : 
    1551           0 : }
    1552             : 
    1553             : //_____________________________________________________________________________
    1554             : Bool_t AliTRDdigitizer::SDigits2Digits()
    1555             : {
    1556             :   //
    1557             :   // Merges the input s-digits and converts them to normal digits
    1558             :   //
    1559             : 
    1560           8 :   if (!MergeSDigits()) {
    1561           0 :     return kFALSE;
    1562             :   }
    1563             : 
    1564           4 :   return ConvertSDigits();
    1565             : 
    1566           4 : }
    1567             : 
    1568             : //_____________________________________________________________________________
    1569             : Bool_t AliTRDdigitizer::MergeSDigits()
    1570             : {
    1571             :   //
    1572             :   // Merges the input s-digits:
    1573             :   //   - The amplitude of the different inputs are summed up.
    1574             :   //   - Of the track IDs from the input dictionaries only one is
    1575             :   //     kept for each input. This works for maximal 3 different merged inputs.
    1576             :   //
    1577             : 
    1578             :   // Number of track dictionary arrays
    1579             :   const Int_t kNDict = AliTRDdigitsManager::kNDict;
    1580             : 
    1581           8 :   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
    1582           4 :   if (!simParam) {
    1583           0 :     AliFatal("Could not get simulation parameters");
    1584           0 :     return kFALSE;
    1585             :   }
    1586             :   
    1587           4 :   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
    1588           4 :   if (!calibration) {
    1589           0 :     AliFatal("Could not get calibration object");
    1590           0 :     return kFALSE;
    1591             :   }
    1592             :   
    1593             :   Int_t iDict = 0;
    1594             :   Int_t jDict = 0;
    1595             : 
    1596             :   AliTRDarraySignal     *digitsA;
    1597             :   AliTRDarraySignal     *digitsB;
    1598           4 :   AliTRDarrayDictionary *dictionaryA[kNDict];
    1599           4 :   AliTRDarrayDictionary *dictionaryB[kNDict];
    1600             : 
    1601             :   AliTRDdigitsManager   *mergeSDigitsManager = 0x0;
    1602             :   // Get the first s-digits
    1603           4 :   fSDigitsManager = (AliTRDdigitsManager *) fSDigitsManagerList->First();
    1604           4 :   if (!fSDigitsManager) { 
    1605           0 :     AliError("No SDigits manager");
    1606           0 :     return kFALSE;
    1607             :   }
    1608             : 
    1609             :   // Loop through the other sets of s-digits
    1610           4 :   mergeSDigitsManager = (AliTRDdigitsManager *) fSDigitsManagerList->After(fSDigitsManager);
    1611             : 
    1612           8 :   if (mergeSDigitsManager) {
    1613           4 :     AliDebug(1,Form("Merge %d input files.",fSDigitsManagerList->GetSize()));
    1614             :   }
    1615             :   else {
    1616          12 :     AliDebug(1,"Only one input file.");
    1617             :   }
    1618             :   
    1619             :   Int_t iMerge = 0;
    1620             : 
    1621           8 :   while (mergeSDigitsManager) {
    1622             : 
    1623           0 :     iMerge++;
    1624             : 
    1625             :     // Loop through the detectors
    1626           0 :     for (Int_t iDet = 0; iDet < AliTRDgeometry::Ndet(); iDet++) {
    1627             : 
    1628           0 :       Int_t nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(iDet);
    1629           0 :       if (mergeSDigitsManager->GetDigitsParam()->GetNTimeBins(iDet) != nTimeTotal) {
    1630           0 :         AliError(Form("Mismatch in the number of time bins [%d,%d] in detector %d"
    1631             :                      ,nTimeTotal
    1632             :                      ,mergeSDigitsManager->GetDigitsParam()->GetNTimeBins(iDet)
    1633             :                      ,iDet));
    1634           0 :         return kFALSE;
    1635             :       }
    1636             : 
    1637           0 :       Int_t nRowMax = fGeo->GetPadPlane(iDet)->GetNrows();
    1638           0 :       Int_t nColMax = fGeo->GetPadPlane(iDet)->GetNcols();
    1639             :           
    1640             :       // Loop through the pixels of one detector and add the signals
    1641           0 :       digitsA = (AliTRDarraySignal *) fSDigitsManager->GetSDigits(iDet);    
    1642           0 :       digitsB = (AliTRDarraySignal *) mergeSDigitsManager->GetSDigits(iDet); 
    1643           0 :       digitsA->Expand();  
    1644           0 :       if (!digitsA->HasData()) continue;
    1645           0 :       digitsB->Expand();    
    1646           0 :       if (!digitsB->HasData()) continue;
    1647             :           
    1648           0 :       for (iDict = 0; iDict < kNDict; iDict++) {
    1649           0 :         dictionaryA[iDict] = (AliTRDarrayDictionary *) fSDigitsManager->GetDictionary(iDet,iDict);
    1650           0 :         dictionaryB[iDict] = (AliTRDarrayDictionary *) mergeSDigitsManager->GetDictionary(iDet,iDict);
    1651           0 :         dictionaryA[iDict]->Expand();  
    1652           0 :         dictionaryB[iDict]->Expand();
    1653             :       }
    1654             : 
    1655             :       // Merge only detectors that contain a signal
    1656             :       Bool_t doMerge = kTRUE;
    1657           0 :       if (fMergeSignalOnly) {
    1658           0 :         if (digitsA->GetOverThreshold(0) == 0) {                             
    1659             :           doMerge = kFALSE;
    1660           0 :         }
    1661             :       }
    1662             :           
    1663           0 :       if (doMerge) {
    1664             :               
    1665           0 :         AliDebug(1,Form("Merge detector %d of input no.%d",iDet,iMerge+1));
    1666             :               
    1667           0 :         for (Int_t iRow  = 0; iRow  <  nRowMax;   iRow++ ) {
    1668           0 :           for (Int_t iCol  = 0; iCol  <  nColMax;   iCol++ ) {
    1669           0 :             for (Int_t iTime = 0; iTime < nTimeTotal; iTime++) {
    1670             :                 
    1671             :               // Add the amplitudes of the summable digits 
    1672           0 :               Float_t ampA = digitsA->GetData(iRow,iCol,iTime);
    1673           0 :               Float_t ampB = digitsB->GetData(iRow,iCol,iTime);
    1674           0 :               ampA += ampB;
    1675           0 :               digitsA->SetData(iRow,iCol,iTime,ampA);
    1676             : 
    1677             :               // Add the mask to the track id if defined.
    1678           0 :               for (iDict = 0; iDict < kNDict; iDict++) {
    1679           0 :                 Int_t trackB = dictionaryB[iDict]->GetData(iRow,iCol,iTime);
    1680           0 :                 if ((fMasks) && (trackB > 0))  {
    1681           0 :                   for (jDict = 0; jDict < kNDict; jDict++) { 
    1682           0 :                     Int_t trackA = dictionaryA[iDict]->GetData(iRow,iCol,iTime); 
    1683           0 :                     if (trackA == 0) {
    1684           0 :                       trackA = trackB + fMasks[iMerge];
    1685           0 :                       dictionaryA[iDict]->SetData(iRow,iCol,iTime,trackA);  
    1686           0 :                     } // if:  track A == 0
    1687             :                   } // for: jDict
    1688             :                 } // if:  fMasks and trackB > 0
    1689             :               } // for: iDict
    1690             : 
    1691             :             } // for: iTime
    1692             :           } // for: iCol
    1693             :         } // for: iRow
    1694             : 
    1695           0 :       } // if:  doMerge
    1696             : 
    1697           0 :       mergeSDigitsManager->RemoveDigits(iDet);
    1698           0 :       mergeSDigitsManager->RemoveDictionaries(iDet);
    1699             :   
    1700           0 :       if (fCompress) {
    1701           0 :         digitsA->Compress(0); 
    1702           0 :         for (iDict = 0; iDict < kNDict; iDict++) {                                     
    1703           0 :           dictionaryA[iDict]->Compress();
    1704             :         }
    1705             :       }
    1706             :       
    1707           0 :     } // for: detectors    
    1708             :       
    1709             :     // The next set of s-digits
    1710           0 :     mergeSDigitsManager = (AliTRDdigitsManager *) fSDigitsManagerList->After(mergeSDigitsManager);
    1711             :     
    1712             :   } // while: mergeDigitsManagers
    1713             :   
    1714           4 :   return kTRUE;
    1715             : 
    1716           8 : }
    1717             : 
    1718             : //_____________________________________________________________________________
    1719             : Bool_t AliTRDdigitizer::ConvertSDigits()
    1720             : {
    1721             :   //
    1722             :   // Converts s-digits to normal digits
    1723             :   //
    1724             : 
    1725             :   AliTRDarraySignal *digitsIn = 0x0;
    1726             : 
    1727           8 :   if (!fSDigitsManager->HasSDigits()) {
    1728           0 :     AliError("No s-digits in digits manager");
    1729           0 :     return kFALSE;
    1730             :   }
    1731             : 
    1732             :   // Loop through the detectors
    1733        4328 :   for (Int_t det = 0; det < AliTRDgeometry::Ndet(); det++) {
    1734             : 
    1735             :     // Get the merged s-digits (signals)
    1736        2160 :     digitsIn = (AliTRDarraySignal *) fSDigitsManager->GetSDigits(det);
    1737        2160 :     if (!digitsIn->HasData()) {
    1738        5916 :       AliDebug(2,Form("No digits for det=%d",det));
    1739             :       continue;
    1740             :     }
    1741             :     
    1742             :     // Convert the merged sdigits to digits
    1743         188 :     if (!Signal2ADC(det,digitsIn)) {
    1744             :       continue;
    1745             :     }
    1746             : 
    1747             :     // Copy the dictionary information to the output array
    1748         188 :     if (!CopyDictionary(det)) {
    1749             :       continue;
    1750             :     }
    1751             : 
    1752             :     // Delete 
    1753         188 :     fSDigitsManager->RemoveDigits(det);
    1754         188 :     fSDigitsManager->RemoveDictionaries(det);
    1755             : 
    1756             :     // Run digital processing
    1757         188 :     RunDigitalProcessing(det);
    1758             : 
    1759             :     // Compress the arrays
    1760         188 :     CompressOutputArrays(det);
    1761             : 
    1762         188 :   } // for: detector numbers
    1763             : 
    1764           4 :   if (AliDataLoader *trklLoader = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets")) {
    1765           4 :     if (trklLoader->Tree())
    1766           4 :       trklLoader->WriteData("OVERWRITE");
    1767             :   }
    1768             : 
    1769             :   // Save the values for the raw data headers
    1770           8 :   if (AliTRDSimParam::Instance()->GetNTBoverwriteOCDB()) {
    1771           4 :     fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDSimParam::Instance()->GetNTimeBins());
    1772           0 :   }
    1773             :   else {
    1774           4 :     fDigitsManager->GetDigitsParam()->SetNTimeBinsAll(AliTRDcalibDB::Instance()->GetNumberOfTimeBinsDCS());
    1775             :   }
    1776           4 :   fDigitsManager->GetDigitsParam()->SetADCbaselineAll(AliTRDSimParam::Instance()->GetADCbaseline());
    1777             : 
    1778           4 :   return kTRUE;
    1779             : 
    1780           4 : }
    1781             : 
    1782             : //_____________________________________________________________________________
    1783             : Bool_t AliTRDdigitizer::CopyDictionary(Int_t det)
    1784             : {
    1785             :   //
    1786             :   // Copies the dictionary information from the s-digits arrays
    1787             :   // to the output arrays
    1788             :   //
    1789             : 
    1790         376 :   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
    1791         188 :   if (!calibration) {
    1792           0 :     AliFatal("Could not get calibration object");
    1793           0 :     return kFALSE;
    1794             :   }
    1795             : 
    1796         564 :   AliDebug(1,Form("Start copying dictionaries for detector=%d",det));
    1797             : 
    1798             :   const Int_t kNDict = AliTRDdigitsManager::kNDict;
    1799         188 :   AliTRDarrayDictionary *dictionaryIn[kNDict];
    1800         188 :   AliTRDarrayDictionary *dictionaryOut[kNDict];
    1801             : 
    1802         188 :   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
    1803         188 :   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
    1804         188 :   Int_t nTimeTotal = fSDigitsManager->GetDigitsParam()->GetNTimeBins(det);
    1805             : 
    1806             :   Int_t row  = 0;
    1807             :   Int_t col  = 0;
    1808             :   Int_t time = 0;
    1809             :   Int_t dict = 0;
    1810             : 
    1811        1504 :   for (dict = 0; dict < kNDict; dict++) {
    1812             : 
    1813         564 :     dictionaryIn[dict]  = (AliTRDarrayDictionary *) fSDigitsManager->GetDictionary(det,dict);
    1814         564 :     dictionaryIn[dict]->Expand();
    1815         564 :     dictionaryOut[dict] = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(det,dict);
    1816         564 :     dictionaryOut[dict]->Allocate(nRowMax,nColMax,nTimeTotal);
    1817             : 
    1818       16728 :     for (row = 0; row < nRowMax; row++) {
    1819     2262000 :       for (col = 0; col < nColMax; col++) {
    1820    62899200 :         for (time = 0; time < nTimeTotal; time++) {
    1821    30326400 :           Int_t track = dictionaryIn[dict]->GetData(row,col,time);
    1822    30326400 :           dictionaryOut[dict]->SetData(row,col,time,track);
    1823             :         } // for: time
    1824             :       } // for: col
    1825             :     } // for: row
    1826             :     
    1827             :   } // for: dictionaries
    1828             :   
    1829             :   return kTRUE;
    1830             : 
    1831         376 : }
    1832             : 
    1833             : //_____________________________________________________________________________
    1834             : void AliTRDdigitizer::CompressOutputArrays(Int_t det)
    1835             : {
    1836             :   //
    1837             :   // Compress the output arrays
    1838             :   //
    1839             : 
    1840             :   const Int_t kNDict = AliTRDdigitsManager::kNDict;
    1841             :   AliTRDarrayDictionary *dictionary = 0x0;
    1842             : 
    1843         752 :   if (fCompress) {
    1844             : 
    1845         376 :     if (!fSDigits) {
    1846             :       AliTRDarrayADC *digits = 0x0;  
    1847         188 :       digits = (AliTRDarrayADC *) fDigitsManager->GetDigits(det);
    1848         188 :       digits->Compress();
    1849         188 :     }
    1850             : 
    1851         376 :     if (fSDigits) {
    1852             :       AliTRDarraySignal *digits = 0x0; 
    1853         188 :       digits = (AliTRDarraySignal *) fDigitsManager->GetSDigits(det);
    1854         188 :       digits->Compress(0);
    1855         188 :     }
    1856             : 
    1857        3008 :     for (Int_t dict = 0; dict < kNDict; dict++) {
    1858        1128 :       dictionary = (AliTRDarrayDictionary *) fDigitsManager->GetDictionary(det,dict);
    1859        1128 :       dictionary->Compress();
    1860             :     }
    1861             : 
    1862         376 :   }
    1863             : 
    1864         376 : }
    1865             : 
    1866             : //_____________________________________________________________________________
    1867             : Bool_t AliTRDdigitizer::WriteDigits() const
    1868             : {
    1869             :   //
    1870             :   // Writes out the TRD-digits and the dictionaries
    1871             :   //
    1872             : 
    1873             :   // Write parameters
    1874           8 :   fRunLoader->CdGAFile();
    1875             : 
    1876             :   // Store the digits and the dictionary in the tree
    1877           4 :   return fDigitsManager->WriteDigits();
    1878             : 
    1879             : }
    1880             : 
    1881             : //_____________________________________________________________________________
    1882             : void AliTRDdigitizer::InitOutput(Int_t iEvent)
    1883             : {
    1884             :   //
    1885             :   // Initializes the output branches
    1886             :   //
    1887             : 
    1888           0 :   fEvent = iEvent;
    1889             :    
    1890           0 :   if (!fRunLoader) {
    1891           0 :     AliError("Run Loader is NULL");
    1892           0 :     return;  
    1893             :   }
    1894             : 
    1895           0 :   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
    1896           0 :   if (!loader) {
    1897           0 :     AliError("Can not get TRD loader from Run Loader");
    1898           0 :     return;
    1899             :   }
    1900             : 
    1901             :   TTree *tree = 0;
    1902             :   
    1903           0 :   if (fSDigits) { 
    1904             :     // If we produce SDigits
    1905           0 :     tree = loader->TreeS();
    1906           0 :     if (!tree) {
    1907           0 :       loader->MakeTree("S");
    1908           0 :       tree = loader->TreeS();
    1909           0 :     }
    1910             :   }
    1911             :   else {
    1912             :     // If we produce Digits
    1913           0 :     tree = loader->TreeD();
    1914           0 :     if (!tree) {
    1915           0 :       loader->MakeTree("D");
    1916           0 :       tree = loader->TreeD();
    1917           0 :     }
    1918             :   }
    1919           0 :   fDigitsManager->SetEvent(iEvent);
    1920           0 :   fDigitsManager->MakeBranch(tree);
    1921             : 
    1922           0 : }
    1923             :   
    1924             : //_____________________________________________________________________________
    1925             : Int_t AliTRDdigitizer::Diffusion(Float_t vdrift, Double_t absdriftlength
    1926             :                                , Double_t exbvalue
    1927             :                                , Double_t &lRow, Double_t &lCol, Double_t &lTime)
    1928             : {
    1929             :   //
    1930             :   // Applies the diffusion smearing to the position of a single electron.
    1931             :   // Depends on absolute drift length.
    1932             :   //
    1933             :   
    1934     1149824 :   Float_t diffL = 0.0;
    1935      574912 :   Float_t diffT = 0.0;
    1936             : 
    1937      574912 :   if (AliTRDCommonParam::Instance()->GetDiffCoeff(diffL,diffT,vdrift)) {
    1938             : 
    1939      574912 :     Float_t driftSqrt = TMath::Sqrt(absdriftlength);
    1940      574912 :     Float_t sigmaT    = driftSqrt * diffT;
    1941      574912 :     Float_t sigmaL    = driftSqrt * diffL;
    1942      574912 :     lRow  = gRandom->Gaus(lRow ,sigmaT);
    1943     1149824 :     if (AliTRDCommonParam::Instance()->ExBOn()) {
    1944     1149824 :       lCol  = gRandom->Gaus(lCol ,sigmaT * 1.0 / (1.0 + exbvalue*exbvalue));
    1945      574912 :       lTime = gRandom->Gaus(lTime,sigmaL * 1.0 / (1.0 + exbvalue*exbvalue));
    1946      574912 :     }
    1947             :     else {
    1948           0 :       lCol  = gRandom->Gaus(lCol ,sigmaT);
    1949           0 :       lTime = gRandom->Gaus(lTime,sigmaL);
    1950             :     }
    1951             : 
    1952             :     return 1;
    1953             : 
    1954             :   }
    1955             :   else {
    1956             : 
    1957           0 :     return 0;
    1958             : 
    1959             :   }
    1960             : 
    1961      574912 : }
    1962             :   
    1963             : //_____________________________________________________________________________
    1964             : void AliTRDdigitizer::RunDigitalProcessing(Int_t det)
    1965             : {
    1966             :   //
    1967             :   // Run the digital processing in the TRAP
    1968             :   //
    1969             : 
    1970         376 :   AliTRDfeeParam *feeParam = AliTRDfeeParam::Instance();
    1971             : 
    1972         188 :   AliTRDarrayADC *digits = fDigitsManager->GetDigits(det);
    1973         188 :   if (!digits)
    1974           0 :     return;
    1975             : 
    1976             :   //Call the methods in the mcm class using the temporary array as input  
    1977             :   // process the data in the same order as in hardware
    1978        1128 :   for (Int_t side = 0; side <= 1; side++) {
    1979        3352 :     for(Int_t rob = side; rob < digits->GetNrow() / 2; rob += 2) {
    1980       44200 :       for(Int_t mcm = 0; mcm < 16; mcm++) {
    1981       20800 :         fMcmSim->Init(det, rob, mcm);
    1982       20800 :         fMcmSim->SetDataByPad(digits, fDigitsManager);
    1983       20800 :         fMcmSim->Filter();
    1984       20800 :         if (feeParam->GetTracklet()) {
    1985       20800 :           fMcmSim->Tracklet();
    1986       20800 :           fMcmSim->StoreTracklets();
    1987       20800 :         }
    1988       20800 :         fMcmSim->ZSMapping();
    1989       20800 :         fMcmSim->WriteData(digits);
    1990             :       }
    1991             :     }
    1992             :   }
    1993         376 : }
    1994             : 

Generated by: LCOV version 1.11