TRD/AliTRDdigitizer.cxx

   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 //  Authors: C. Blume (blume@ikf.uni-frankfurt.de)                           //
  22 //           C. Lippmann                                                     //
  23 //           B. Vulpescu                                                     //
  24 //                                                                           //
  25 //  The following effects are included:                                      //
  26 //      - Diffusion                                                          //
  27 //      - ExB effects                                                        //
  28 //      - Gas gain including fluctuations                                    //
  29 //      - Pad-response (simple Gaussian approximation)                       //
  30 //      - Time-response                                                      //
  31 //      - Electronics noise                                                  //
  32 //      - Electronics gain                                                   //
  33 //      - Digitization                                                       //
  34 //      - ADC threshold                                                      //
  35 //  The corresponding parameter can be adjusted via the various              //
  36 //  Set-functions. If these parameters are not explicitly set, default       //
  37 //  values are used (see Init-function).                                     //
  38 //                                                                           //
  39 ///////////////////////////////////////////////////////////////////////////////
  40
  41 #include <stdlib.h>
  42
  43 #include <TMath.h>
  44 #include <TVector.h>
  45 #include <TRandom.h>
  46 #include <TROOT.h>
  47 #include <TTree.h>
  48 #include <TFile.h>
  49 #include <TF1.h>
  50 #include <TList.h>
  51 #include <TTask.h>
  52 #include <TGeoManager.h>
  53
  54 #include "AliRun.h"
  55 #include "AliRunLoader.h"
  56 #include "AliLoader.h"
  57 #include "AliConfig.h"
  58 #include "AliMagF.h"
  59 #include "AliRunDigitizer.h"
  60 #include "AliRunLoader.h"
  61 #include "AliLoader.h"
  62 #include "AliLog.h"
  63
  64 #include "AliTRD.h"
  65 #include "AliTRDhit.h"
  66 #include "AliTRDdigitizer.h"
  67 #include "AliTRDdataArrayI.h"
  68 #include "AliTRDdataArrayF.h"
  69 #include "AliTRDsegmentArray.h"
  70 #include "AliTRDdigitsManager.h"
  71 #include "AliTRDgeometry.h"
  72 #include "AliTRDpadPlane.h"
  73 #include "AliTRDcalibDB.h"
  74 #include "AliTRDSimParam.h"
  75 #include "AliTRDCommonParam.h"
  76
  77 #include "Cal/AliTRDCalROC.h"
  78 #include "Cal/AliTRDCalDet.h"
  79
  80 ClassImp(AliTRDdigitizer)
  81
  82 //_____________________________________________________________________________
  83 AliTRDdigitizer::AliTRDdigitizer()
  84   :AliDigitizer()
  85   ,fRunLoader(0)
  86   ,fDigitsManager(0)
  87   ,fSDigitsManager(0)
  88   ,fSDigitsManagerList(0)
  89   ,fTRD(0)
  90   ,fGeo(0)
  91   ,fEvent(0)
  92   ,fMasks(0)
  93   ,fCompress(kTRUE)
  94   ,fSDigits(kFALSE)
  95   ,fSDigitsScale(0.0)
  96   ,fMergeSignalOnly(kFALSE)
  97   ,fDiffLastVdrift(0)
  98   ,fDiffusionT(0.0)
  99   ,fDiffusionL(0.0)
 100   ,fOmegaTau(0.0)
 101   ,fLorentzFactor(0.0)
 102   ,fTimeLastVdrift(0)
 103   ,fTimeStruct1(0)
 104   ,fTimeStruct2(0)
 105   ,fVDlo(0)
 106   ,fVDhi(0)
 107 {
 108   //
 109   // AliTRDdigitizer default constructor
 110   //
 111
 112   Init();
 113
 114 }
 115
 116 //_____________________________________________________________________________
 117 AliTRDdigitizer::AliTRDdigitizer(const Text_t *name, const Text_t *title)
 118   :AliDigitizer(name,title)
 119   ,fRunLoader(0)
 120   ,fDigitsManager(0)
 121   ,fSDigitsManager(0)
 122   ,fSDigitsManagerList(0)
 123   ,fTRD(0)
 124   ,fGeo(0)
 125   ,fEvent(0)
 126   ,fMasks(0)
 127   ,fCompress(kTRUE)
 128   ,fSDigits(kFALSE)
 129   ,fSDigitsScale(0.0)
 130   ,fMergeSignalOnly(kFALSE)
 131   ,fDiffLastVdrift(0)
 132   ,fDiffusionT(0.0)
 133   ,fDiffusionL(0.0)
 134   ,fOmegaTau(0.0)
 135   ,fLorentzFactor(0.0)
 136   ,fTimeLastVdrift(0)
 137   ,fTimeStruct1(0)
 138   ,fTimeStruct2(0)
 139   ,fVDlo(0)
 140   ,fVDhi(0)
 141 {
 142   //
 143   // AliTRDdigitizer constructor
 144   //
 145
 146   Init();
 147
 148 }
 149
 150 //_____________________________________________________________________________
 151 AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager
 152                                , const Text_t *name, const Text_t *title)
 153   :AliDigitizer(manager,name,title)
 154   ,fRunLoader(0)
 155   ,fDigitsManager(0)
 156   ,fSDigitsManager(0)
 157   ,fSDigitsManagerList(0)
 158   ,fTRD(0)
 159   ,fGeo(0)
 160   ,fEvent(0)
 161   ,fMasks(0)
 162   ,fCompress(kTRUE)
 163   ,fSDigits(kFALSE)
 164   ,fSDigitsScale(0.0)
 165   ,fMergeSignalOnly(kFALSE)
 166   ,fDiffLastVdrift(0)
 167   ,fDiffusionT(0.0)
 168   ,fDiffusionL(0.0)
 169   ,fOmegaTau(0.0)
 170   ,fLorentzFactor(0.0)
 171   ,fTimeLastVdrift(0)
 172   ,fTimeStruct1(0)
 173   ,fTimeStruct2(0)
 174   ,fVDlo(0)
 175   ,fVDhi(0)
 176 {
 177   //
 178   // AliTRDdigitizer constructor
 179   //
 180
 181   Init();
 182
 183 }
 184
 185 //_____________________________________________________________________________
 186 AliTRDdigitizer::AliTRDdigitizer(AliRunDigitizer *manager)
 187   :AliDigitizer(manager,"AliTRDdigitizer","TRD digitizer")
 188   ,fRunLoader(0)
 189   ,fDigitsManager(0)
 190   ,fSDigitsManager(0)
 191   ,fSDigitsManagerList(0)
 192   ,fTRD(0)
 193   ,fGeo(0)
 194   ,fEvent(0)
 195   ,fMasks(0)
 196   ,fCompress(kTRUE)
 197   ,fSDigits(kFALSE)
 198   ,fSDigitsScale(0.0)
 199   ,fMergeSignalOnly(kFALSE)
 200   ,fDiffLastVdrift(0)
 201   ,fDiffusionT(0.0)
 202   ,fDiffusionL(0.0)
 203   ,fOmegaTau(0.0)
 204   ,fLorentzFactor(0.0)
 205   ,fTimeLastVdrift(0)
 206   ,fTimeStruct1(0)
 207   ,fTimeStruct2(0)
 208   ,fVDlo(0)
 209   ,fVDhi(0)
 210 {
 211   //
 212   // AliTRDdigitizer constructor
 213   //
 214
 215   Init();
 216
 217 }
 218
 219 //_____________________________________________________________________________
 220 Bool_t AliTRDdigitizer::Init()
 221 {
 222   //
 223   // Initialize the digitizer with default values
 224   //
 225
 226   fRunLoader          = 0;
 227   fDigitsManager      = 0;
 228   fSDigitsManager     = 0;
 229   fSDigitsManagerList = 0;
 230   fTRD                = 0;
 231   fGeo                = 0;
 232
 233   fEvent              = 0;
 234   fMasks              = 0;
 235   fCompress           = kTRUE;
 236   fSDigits            = kFALSE;
 237   fSDigitsScale       = 100.0;
 238   fMergeSignalOnly    = kFALSE;
 239
 240   fTimeLastVdrift     = -1;
 241   fTimeStruct1        =  0;
 242   fTimeStruct2        =  0;
 243   fVDlo               =  0;
 244   fVDhi               =  0;
 245
 246   fDiffLastVdrift     = -1;
 247   fDiffusionT         =  0.0;
 248   fDiffusionL         =  0.0;
 249   fLorentzFactor      =  0.0;
 250
 251   return AliDigitizer::Init();
 252
 253 }
 254
 255 //_____________________________________________________________________________
 256 AliTRDdigitizer::AliTRDdigitizer(const AliTRDdigitizer &d)
 257   :AliDigitizer(d)
 258   ,fRunLoader(0)
 259   ,fDigitsManager(0)
 260   ,fSDigitsManager(0)
 261   ,fSDigitsManagerList(0)
 262   ,fTRD(0)
 263   ,fGeo(0)
 264   ,fEvent(0)
 265   ,fMasks(0)
 266   ,fCompress(d.fCompress)
 267   ,fSDigits(d.fSDigits)
 268   ,fSDigitsScale(d.fSDigitsScale)
 269   ,fMergeSignalOnly(d.fMergeSignalOnly)
 270   ,fDiffLastVdrift(-1)
 271   ,fDiffusionT(0.0)
 272   ,fDiffusionL(0.0)
 273   ,fOmegaTau(0.0)
 274   ,fLorentzFactor(0.0)
 275   ,fTimeLastVdrift(-1)
 276   ,fTimeStruct1(0)
 277   ,fTimeStruct2(0)
 278   ,fVDlo(0)
 279   ,fVDhi(0)
 280 {
 281   //
 282   // AliTRDdigitizer copy constructor
 283   //
 284
 285   // Do not copy timestructs, just invalidate lastvdrift.
 286   // Next time they are requested, they get recalculated
 287   if (((AliTRDdigitizer &) d).fTimeStruct1) {
 288     delete [] ((AliTRDdigitizer &) d).fTimeStruct1;
 289     ((AliTRDdigitizer &) d).fTimeStruct1 = 0;
 290   }
 291   if (((AliTRDdigitizer &) d).fTimeStruct2) {
 292     delete [] ((AliTRDdigitizer &) d).fTimeStruct2;
 293     ((AliTRDdigitizer &) d).fTimeStruct2 = 0;
 294   }
 295
 296 }
 297
 298 //_____________________________________________________________________________
 299 AliTRDdigitizer::~AliTRDdigitizer()
 300 {
 301   //
 302   // AliTRDdigitizer destructor
 303   //
 304
 305   if (fDigitsManager) {
 306     delete fDigitsManager;
 307     fDigitsManager      = 0;
 308   }
 309
 310   if (fDigitsManager) {
 311     delete fSDigitsManager;
 312     fSDigitsManager     = 0;
 313   }
 314
 315   if (fSDigitsManagerList) {
 316     fSDigitsManagerList->Delete();
 317     delete fSDigitsManagerList;
 318     fSDigitsManagerList = 0;
 319   }
 320
 321   if (fMasks) {
 322     delete [] fMasks;
 323     fMasks       = 0;
 324   }
 325
 326   if (fTimeStruct1) {
 327     delete [] fTimeStruct1;
 328     fTimeStruct1 = 0;
 329   }
 330
 331   if (fTimeStruct2) {
 332     delete [] fTimeStruct2;
 333     fTimeStruct2 = 0;
 334   }
 335
 336   if (fGeo) {
 337     delete fGeo;
 338     fGeo = 0;
 339   }
 340
 341 }
 342
 343 //_____________________________________________________________________________
 344 AliTRDdigitizer &AliTRDdigitizer::operator=(const AliTRDdigitizer &d)
 345 {
 346   //
 347   // Assignment operator
 348   //
 349
 350   if (this != &d) ((AliTRDdigitizer &) d).Copy(*this);
 351
 352   return *this;
 353
 354 }
 355
 356 //_____________________________________________________________________________
 357 void AliTRDdigitizer::Copy(TObject &d) const
 358 {
 359   //
 360   // Copy function
 361   //
 362
 363   ((AliTRDdigitizer &) d).fRunLoader          = 0;
 364   ((AliTRDdigitizer &) d).fDigitsManager      = 0;
 365   ((AliTRDdigitizer &) d).fSDigitsManager     = 0;
 366   ((AliTRDdigitizer &) d).fSDigitsManagerList = 0;
 367   ((AliTRDdigitizer &) d).fTRD                = 0;
 368   ((AliTRDdigitizer &) d).fGeo                = 0;
 369   ((AliTRDdigitizer &) d).fEvent              = 0;
 370   ((AliTRDdigitizer &) d).fMasks              = 0;
 371   ((AliTRDdigitizer &) d).fCompress           = fCompress;
 372   ((AliTRDdigitizer &) d).fSDigits            = fSDigits;
 373   ((AliTRDdigitizer &) d).fSDigitsScale       = fSDigitsScale;
 374   ((AliTRDdigitizer &) d).fMergeSignalOnly    = fMergeSignalOnly;
 375
 376   AliTRDdigitizer& target = (AliTRDdigitizer &) d;
 377
 378   // Do not copy timestructs, just invalidate lastvdrift.
 379   // Next time they are requested, they get recalculated
 380   if (target.fTimeStruct1) {
 381     delete [] target.fTimeStruct1;
 382     target.fTimeStruct1 = 0;
 383   }
 384   if (target.fTimeStruct2) {
 385     delete [] target.fTimeStruct2;
 386     target.fTimeStruct2 = 0;
 387   }
 388   target.fTimeLastVdrift = -1;
 389
 390 }
 391
 392 //_____________________________________________________________________________
 393 void AliTRDdigitizer::Exec(Option_t *option)
 394 {
 395   //
 396   // Executes the merging
 397   //
 398
 399   Int_t iInput;
 400
 401   AliTRDdigitsManager *sdigitsManager;
 402
 403   TString optionString = option;
 404   if (optionString.Contains("deb")) {
 405     AliLog::SetClassDebugLevel("AliTRDdigitizer",1);
 406     AliInfo("Called with debug option");
 407   }
 408
 409   // The AliRoot file is already connected by the manager
 410   AliRunLoader *inrl;
 411
 412   if (gAlice) {
 413     AliDebug(1,"AliRun object found on file.");
 414   }
 415   else {
 416     inrl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 417     inrl->LoadgAlice();
 418     gAlice = inrl->GetAliRun();
 419     if (!gAlice) {
 420       AliError("Could not find AliRun object.")
 421       return;
 422     }
 423   }
 424
 425   Int_t nInput = fManager->GetNinputs();
 426   fMasks       = new Int_t[nInput];
 427   for (iInput = 0; iInput < nInput; iInput++) {
 428     fMasks[iInput] = fManager->GetMask(iInput);
 429   }
 430
 431   //
 432   // Initialization
 433   //
 434
 435   AliRunLoader *orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 436
 437   if (InitDetector()) {
 438
 439     AliLoader *ogime = orl->GetLoader("TRDLoader");
 440
 441     TTree *tree = 0;
 442     if (fSDigits) {
 443       // If we produce SDigits
 444       tree = ogime->TreeS();
 445       if (!tree) {
 446         ogime->MakeTree("S");
 447         tree = ogime->TreeS();
 448       }
 449     }
 450     else {
 451       // If we produce Digits
 452       tree = ogime->TreeD();
 453       if (!tree) {
 454         ogime->MakeTree("D");
 455         tree = ogime->TreeD();
 456       }
 457     }
 458
 459     MakeBranch(tree);
 460
 461   }
 462
 463   for (iInput = 0; iInput < nInput; iInput++) {
 464
 465     AliDebug(1,Form("Add input stream %d",iInput));
 466
 467     // Check if the input tree exists
 468     inrl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput));
 469     AliLoader *gime = inrl->GetLoader("TRDLoader");
 470
 471     TTree *treees = gime->TreeS();
 472     if (treees == 0x0) {
 473       if (gime->LoadSDigits()) {
 474         AliError(Form("Error Occured while loading S. Digits for input %d.",iInput));
 475         return;
 476       }
 477       treees = gime->TreeS();
 478     }
 479
 480     if (treees == 0x0) {
 481       AliError(Form("Input stream %d does not exist",iInput));
 482       return;
 483     }
 484
 485     // Read the s-digits via digits manager
 486     sdigitsManager = new AliTRDdigitsManager();
 487     sdigitsManager->SetSDigits(kTRUE);
 488
 489     AliRunLoader *rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput));
 490     AliLoader *gimme = rl->GetLoader("TRDLoader");
 491     if (!gimme->TreeS()) {
 492       gimme->LoadSDigits();
 493     }
 494     sdigitsManager->ReadDigits(gimme->TreeS());
 495
 496     // Add the s-digits to the input list
 497     AddSDigitsManager(sdigitsManager);
 498
 499   }
 500
 501   // Convert the s-digits to normal digits
 502   AliDebug(1,"Do the conversion");
 503   SDigits2Digits();
 504
 505   // Store the digits
 506   AliDebug(1,"Write the digits");
 507   fDigitsManager->WriteDigits();
 508
 509   // Write parameters
 510   orl->CdGAFile();
 511
 512   AliDebug(1,"Done");
 513
 514   DeleteSDigitsManager();
 515
 516 }
 517
 518 //_____________________________________________________________________________
 519 Bool_t AliTRDdigitizer::Open(const Char_t *file, Int_t nEvent)
 520 {
 521   //
 522   // Opens a ROOT-file with TRD-hits and reads in the hit-tree
 523   //
 524   // Connect the AliRoot file containing Geometry, Kine, and Hits
 525   //
 526
 527   TString evfoldname = AliConfig::GetDefaultEventFolderName();
 528
 529   fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
 530   if (!fRunLoader) {
 531     fRunLoader = AliRunLoader::Open(file,evfoldname,"UPDATE");
 532   }
 533   if (!fRunLoader) {
 534     AliError(Form("Can not open session for file %s.",file));
 535     return kFALSE;
 536   }
 537
 538   if (!fRunLoader->GetAliRun()) {
 539     fRunLoader->LoadgAlice();
 540   }
 541   gAlice = fRunLoader->GetAliRun();
 542
 543   if (gAlice) {
 544     AliDebug(1,"AliRun object found on file.");
 545   }
 546   else {
 547     AliError("Could not find AliRun object.");
 548     return kFALSE;
 549   }
 550
 551   fEvent = nEvent;
 552
 553   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
 554   if (!loader) {
 555     AliError("Can not get TRD loader from Run Loader");
 556     return kFALSE;
 557   }
 558
 559   if (InitDetector()) {
 560     TTree *tree = 0;
 561     if (fSDigits) {
 562       // If we produce SDigits
 563       tree = loader->TreeS();
 564       if (!tree) {
 565         loader->MakeTree("S");
 566         tree = loader->TreeS();
 567       }
 568     }
 569     else {
 570       // If we produce Digits
 571       tree = loader->TreeD();
 572       if (!tree) {
 573         loader->MakeTree("D");
 574         tree = loader->TreeD();
 575       }
 576     }
 577     return MakeBranch(tree);
 578   }
 579   else {
 580     return kFALSE;
 581   }
 582
 583 }
 584
 585 //_____________________________________________________________________________
 586 Bool_t AliTRDdigitizer::Open(AliRunLoader *runLoader, Int_t nEvent)
 587 {
 588   //
 589   // Opens a ROOT-file with TRD-hits and reads in the hit-tree
 590   //
 591   // Connect the AliRoot file containing Geometry, Kine, and Hits
 592   //
 593
 594   fRunLoader = runLoader;
 595   if (!fRunLoader) {
 596     AliError("RunLoader does not exist");
 597     return kFALSE;
 598   }
 599
 600   if (!fRunLoader->GetAliRun()) {
 601     fRunLoader->LoadgAlice();
 602   }
 603   gAlice = fRunLoader->GetAliRun();
 604
 605   if (gAlice) {
 606     AliDebug(1,"AliRun object found on file.");
 607   }
 608   else {
 609     AliError("Could not find AliRun object.");
 610     return kFALSE;
 611   }
 612
 613   fEvent = nEvent;
 614
 615   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
 616   if (!loader) {
 617     AliError("Can not get TRD loader from Run Loader");
 618     return kFALSE;
 619   }
 620
 621   if (InitDetector()) {
 622     TTree *tree = 0;
 623     if (fSDigits) {
 624       // If we produce SDigits
 625       tree = loader->TreeS();
 626       if (!tree) {
 627         loader->MakeTree("S");
 628         tree = loader->TreeS();
 629       }
 630     }
 631     else {
 632       // If we produce Digits
 633       tree = loader->TreeD();
 634       if (!tree) {
 635         loader->MakeTree("D");
 636         tree = loader->TreeD();
 637       }
 638     }
 639     return MakeBranch(tree);
 640   }
 641   else {
 642     return kFALSE;
 643   }
 644
 645 }
 646
 647 //_____________________________________________________________________________
 648 Bool_t AliTRDdigitizer::InitDetector()
 649 {
 650   //
 651   // Sets the pointer to the TRD detector and the geometry
 652   //
 653
 654   // Get the pointer to the detector class and check for version 1
 655   fTRD = (AliTRD *) gAlice->GetDetector("TRD");
 656   if (!fTRD) {
 657     AliFatal("No TRD module found");
 658     exit(1);
 659   }
 660   if (fTRD->IsVersion() != 1) {
 661     AliFatal("TRD must be version 1 (slow simulator)");
 662     exit(1);
 663   }
 664
 665   // Get the geometry
 666   fGeo = new AliTRDgeometry();
 667
 668   // Create a digits manager
 669   delete fDigitsManager;
 670   fDigitsManager = new AliTRDdigitsManager();
 671   fDigitsManager->SetSDigits(fSDigits);
 672   fDigitsManager->CreateArrays();
 673   fDigitsManager->SetEvent(fEvent);
 674
 675   // The list for the input s-digits manager to be merged
 676   if (fSDigitsManagerList) {
 677     fSDigitsManagerList->Delete();
 678   }
 679   else {
 680     fSDigitsManagerList = new TList();
 681   }
 682
 683   return kTRUE;
 684
 685 }
 686
 687 //_____________________________________________________________________________
 688 Bool_t AliTRDdigitizer::MakeBranch(TTree *tree) const
 689 {
 690   //
 691   // Create the branches for the digits array
 692   //
 693
 694   return fDigitsManager->MakeBranch(tree);
 695
 696 }
 697
 698 //_____________________________________________________________________________
 699 Bool_t AliTRDdigitizer::MakeDigits()
 700 {
 701   //
 702   // Creates digits.
 703   //
 704
 705   ///////////////////////////////////////////////////////////////
 706   // Parameter
 707   ///////////////////////////////////////////////////////////////
 708
 709   // Converts number of electrons to fC
 710   const Double_t kEl2fC  = 1.602e-19 * 1.0e15;
 711
 712   ///////////////////////////////////////////////////////////////
 713
 714   // Number of pads included in the pad response
 715   const Int_t kNpad      = 3;
 716
 717   // Number of track dictionary arrays
 718   const Int_t kNDict     = AliTRDdigitsManager::kNDict;
 719
 720   // Width of the amplification region
 721   const Float_t kAmWidth = AliTRDgeometry::AmThick();
 722   // Width of the drift region
 723   const Float_t kDrWidth = AliTRDgeometry::DrThick();
 724   // Drift + amplification region
 725   const Float_t kDrMin   =          - 0.5 * kAmWidth;
 726   const Float_t kDrMax   = kDrWidth + 0.5 * kAmWidth;
 727
 728   Int_t    iRow;
 729   Int_t    iCol;
 730   Int_t    iTime;
 731   Int_t    iPad;
 732   Int_t    iDict  = 0;
 733   Int_t    nBytes = 0;
 734
 735   Int_t    totalSizeDigits = 0;
 736   Int_t    totalSizeDict0  = 0;
 737   Int_t    totalSizeDict1  = 0;
 738   Int_t    totalSizeDict2  = 0;
 739
 740   Int_t    timeBinTRFend   = 1;
 741
 742   Double_t pos[3];
 743   Double_t loc[3];
 744   Double_t padSignal[kNpad];
 745   Double_t signalOld[kNpad];
 746
 747   AliTRDdataArrayF *signals  = 0;
 748   AliTRDdataArrayI *digits   = 0;
 749   AliTRDdataArrayI *dictionary[kNDict];
 750
 751   AliTRDpadPlane   *padPlane = 0;
 752
 753   AliTRDCalROC     *calVdriftROC          = 0;
 754   Float_t           calVdriftDetValue     = 0.0;
 755   AliTRDCalROC     *calT0ROC              = 0;
 756   Float_t           calT0DetValue         = 0.0;
 757   AliTRDCalROC     *calGainFactorROC      = 0;
 758   Float_t           calGainFactorDetValue = 0.0;
 759
 760   if (!gGeoManager) {
 761     AliFatal("No geometry manager!");
 762   }
 763
 764   if (!fGeo) {
 765     AliError("No geometry defined");
 766     return kFALSE;
 767   }
 768   fGeo->ReadGeoMatrices();
 769
 770   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
 771   if (!simParam) {
 772     AliFatal("Could not get simulation parameters");
 773     return kFALSE;
 774   }
 775
 776   AliTRDCommonParam *commonParam = AliTRDCommonParam::Instance();
 777   if (!commonParam) {
 778     AliFatal("Could not get common parameterss");
 779     return kFALSE;
 780   }
 781
 782   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
 783   if (!calibration) {
 784     AliFatal("Could not get calibration object");
 785     return kFALSE;
 786   }
 787
 788   AliDebug(1,"Start creating digits");
 789
 790   // Create a container for the amplitudes
 791   AliTRDsegmentArray *signalsArray = new AliTRDsegmentArray("AliTRDdataArrayF"
 792                                                            ,AliTRDgeometry::Ndet());
 793
 794   AliLoader *gimme = fRunLoader->GetLoader("TRDLoader");
 795   if (!gimme->TreeH()) {
 796     gimme->LoadHits();
 797   }
 798   TTree *hitTree = gimme->TreeH();
 799   if (hitTree == 0x0) {
 800     AliError("Can not get TreeH");
 801     return kFALSE;
 802   }
 803   fTRD->SetTreeAddress();
 804
 805   // Get the number of entries in the hit tree
 806   // (Number of primary particles creating a hit somewhere)
 807   Int_t nTrack = (Int_t) hitTree->GetEntries();
 808   AliDebug(1,Form("Found %d primary particles",nTrack));
 809   AliDebug(1,Form("Sampling = %.0fMHz"        ,commonParam->GetSamplingFrequency()));
 810   AliDebug(1,Form("Gain     = %d"             ,((Int_t) simParam->GetGasGain())));
 811   AliDebug(1,Form("Noise    = %d"             ,((Int_t) simParam->GetNoise())));
 812   if (simParam->TimeStructOn()) {
 813     AliDebug(1,"Time Structure of drift cells implemented.");
 814   }
 815   else {
 816     AliDebug(1,"Constant drift velocity in drift cells.");
 817   }
 818   if (simParam->TRFOn()) {
 819     timeBinTRFend = ((Int_t) (simParam->GetTRFhi()
 820                   * commonParam->GetSamplingFrequency())) - 1;
 821     AliDebug(1,Form("Sample the TRF up to bin %d",timeBinTRFend));
 822   }
 823
 824   // Get the detector wise calibration objects
 825   const AliTRDCalDet *calVdriftDet     = calibration->GetVdriftDet();
 826   const AliTRDCalDet *calT0Det         = calibration->GetT0Det();
 827   const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
 828
 829   Int_t   detectorOld  = -1;
 830   Int_t   countHits    =  0;
 831
 832   Int_t   nTimeTotal   = calibration->GetNumberOfTimeBins();
 833   Float_t samplingRate = commonParam->GetSamplingFrequency();
 834   Float_t elAttachProp = simParam->GetElAttachProp() / 100.0;
 835
 836   // Loop through all entries in the tree
 837   for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
 838
 839     gAlice->ResetHits();
 840     nBytes += hitTree->GetEvent(iTrack);
 841
 842     // Loop through the TRD hits
 843     Int_t iHit = 0;
 844     AliTRDhit *hit = (AliTRDhit *) fTRD->FirstHit(-1);
 845     while (hit) {
 846
 847       countHits++;
 848       iHit++;
 849
 850       Float_t q        = hit->GetCharge();
 851       // Don't analyze test hits
 852       if (((Int_t) q) == 0) {
 853         hit = (AliTRDhit *) fTRD->NextHit();
 854         continue;
 855       }
 856
 857       pos[0]           = hit->X();
 858       pos[1]           = hit->Y();
 859       pos[2]           = hit->Z();
 860       Int_t   track    = hit->Track();
 861       Int_t   detector = hit->GetDetector();
 862       Float_t hittime  = hit->GetTime();
 863       Int_t   plane    = fGeo->GetPlane(detector);
 864       Int_t   chamber  = fGeo->GetChamber(detector);
 865       padPlane         = fGeo->GetPadPlane(plane,chamber);
 866       Float_t row0     = padPlane->GetRow0ROC();
 867       Int_t   nRowMax  = padPlane->GetNrows();
 868       Int_t   nColMax  = padPlane->GetNcols();
 869       Int_t   inDrift  = 1;
 870
 871       // Find the current volume with the geo manager
 872       gGeoManager->SetCurrentPoint(pos);
 873       gGeoManager->FindNode();
 874       if (strstr(gGeoManager->GetPath(),"/UK")) {
 875         inDrift = 0;
 876       }
 877
 878       if (detector != detectorOld) {
 879
 880         // Compress the old one if enabled
 881         if ((fCompress) && (detectorOld > -1)) {
 882           signals->Compress(1,0);
 883           for (iDict = 0; iDict < kNDict; iDict++) {
 884             dictionary[iDict]->Compress(1,0);
 885           }
 886         }
 887         // Get the new container
 888         signals = (AliTRDdataArrayF *) signalsArray->At(detector);
 889         if (signals->GetNtime() == 0) {
 890           // Allocate a new one if not yet existing
 891           signals->Allocate(nRowMax,nColMax,nTimeTotal);
 892         }
 893         else {
 894           // Expand an existing one
 895           if (fCompress) {
 896             signals->Expand();
 897           }
 898         }
 899         // The same for the dictionary
 900         for (iDict = 0; iDict < kNDict; iDict++) {
 901           dictionary[iDict] = fDigitsManager->GetDictionary(detector,iDict);
 902           if (dictionary[iDict]->GetNtime() == 0) {
 903             dictionary[iDict]->Allocate(nRowMax,nColMax,nTimeTotal);
 904           }
 905           else {
 906             if (fCompress) dictionary[iDict]->Expand();
 907           }
 908         }
 909
 910         // Get the calibration objects
 911         calVdriftROC      = calibration->GetVdriftROC(detector);
 912         calVdriftDetValue = calVdriftDet->GetValue(detector);
 913         calT0ROC          = calibration->GetT0ROC(detector);
 914         calT0DetValue     = calT0Det->GetValue(detector);
 915
 916         detectorOld = detector;
 917
 918       }
 919
 920       // Go to the local coordinate system:
 921       // loc[0] - col  direction in amplification or driftvolume
 922       // loc[1] - row  direction in amplification or driftvolume
 923       // loc[2] - time direction in amplification or driftvolume
 924       gGeoManager->MasterToLocal(pos,loc);
 925       if (inDrift) {
 926         // Relative to middle of amplification region
 927         loc[2] = loc[2] - kDrWidth/2.0 - kAmWidth/2.0;
 928       }
 929
 930       // The driftlength [cm] (w/o diffusion yet !).
 931       // It is negative if the hit is between pad plane and anode wires.
 932       Double_t driftlength = -1.0 * loc[2];
 933
 934       // Stupid patch to take care of TR photons that are absorbed
 935       // outside the chamber volume. A real fix would actually need
 936       // a more clever implementation of the TR hit generation
 937       if (q < 0.0) {
 938         if ((loc[1] < padPlane->GetRowEndROC()) ||
 939             (loc[1] > padPlane->GetRow0ROC())) {
 940           hit = (AliTRDhit *) fTRD->NextHit();
 941           continue;
 942         }
 943         if ((driftlength < kDrMin) ||
 944             (driftlength > kDrMax)) {
 945           hit = (AliTRDhit *) fTRD->NextHit();
 946           continue;
 947         }
 948       }
 949
 950       // Get row and col of unsmeared electron to retrieve drift velocity
 951       // The pad row (z-direction)
 952       Int_t    rowE         = padPlane->GetPadRowNumberROC(loc[1]);
 953       if (rowE < 0) {
 954         hit = (AliTRDhit *) fTRD->NextHit();
 955         continue;
 956       }
 957       Double_t rowOffset    = padPlane->GetPadRowOffsetROC(rowE,loc[1]);
 958
 959       // The pad column (rphi-direction)
 960       Double_t offsetTilt   = padPlane->GetTiltOffset(rowOffset);
 961       Int_t    colE         = padPlane->GetPadColNumber(loc[0]+offsetTilt);
 962       if (colE < 0) {
 963         hit = (AliTRDhit *) fTRD->NextHit();
 964         continue;
 965       }
 966       Double_t colOffset    = padPlane->GetPadColOffset(colE,loc[0]+offsetTilt);
 967
 968       Float_t driftvelocity = calVdriftDetValue * calVdriftROC->GetValue(colE,rowE);
 969
 970       // Normalized drift length
 971       Double_t absdriftlength = TMath::Abs(driftlength);
 972       if (commonParam->ExBOn()) {
 973         absdriftlength /= TMath::Sqrt(GetLorentzFactor(driftvelocity));
 974       }
 975
 976       // Loop over all electrons of this hit
 977       // TR photons produce hits with negative charge
 978       Int_t nEl = ((Int_t) TMath::Abs(q));
 979       for (Int_t iEl = 0; iEl < nEl; iEl++) {
 980
 981         // Now the real local coordinate system of the ROC
 982         // column direction: locC
 983         // row direction:    locR
 984         // time direction:   locT
 985         // locR and locC are identical to the coordinates of the corresponding
 986         // volumina of the drift or amplification region.
 987         // locT is defined relative to the wire plane (i.e. middle of amplification
 988         // region), meaming locT = 0, and is negative for hits coming from the
 989         // drift region.
 990         Double_t locC = loc[0];
 991         Double_t locR = loc[1];
 992         Double_t locT = loc[2];
 993
 994         // Electron attachment
 995         if (simParam->ElAttachOn()) {
 996           if (gRandom->Rndm() < (absdriftlength * elAttachProp)) {
 997             continue;
 998           }
 999         }
1000
1001         // Apply the diffusion smearing
1002         if (simParam->DiffusionOn()) {
1003           if (!(Diffusion(driftvelocity,absdriftlength,locR,locC,locT))) {
1004             continue;
1005           }
1006         }
1007
1008         // Apply E x B effects (depends on drift direction)
1009         if (commonParam->ExBOn()) {
1010           if (!(ExB(driftvelocity,driftlength,locC))) {
1011             continue;
1012           }
1013         }
1014
1015         // The electron position after diffusion and ExB in pad coordinates.
1016         // The pad row (z-direction)
1017         rowE       = padPlane->GetPadRowNumberROC(locR);
1018         if (rowE < 0) continue;
1019         rowOffset  = padPlane->GetPadRowOffsetROC(rowE,locR);
1020
1021         // The pad column (rphi-direction)
1022         offsetTilt = padPlane->GetTiltOffset(rowOffset);
1023         colE       = padPlane->GetPadColNumber(locC+offsetTilt);
1024         if (colE < 0) continue;
1025         colOffset  = padPlane->GetPadColOffset(colE,locC+offsetTilt);
1026
1027         // Also re-retrieve drift velocity because col and row may have changed
1028         driftvelocity = calVdriftDetValue * calVdriftROC->GetValue(colE,rowE);
1029         Float_t t0    = calT0DetValue     + calT0ROC->GetValue(colE,rowE);
1030
1031         // Convert the position to drift time [mus], using either constant drift velocity or
1032         // time structure of drift cells (non-isochronity, GARFIELD calculation).
1033         // Also add absolute time of hits to take pile-up events into account properly
1034         Double_t drifttime;
1035         if (simParam->TimeStructOn()) {
1036           // Get z-position with respect to anode wire
1037           Double_t zz  =  row0 - locR + simParam->GetAnodeWireOffset();
1038           zz -= ((Int_t)(2 * zz)) / 2.0;
1039           if (zz > 0.25) {
1040             zz  = 0.5 - zz;
1041           }
1042           // Use drift time map (GARFIELD)
1043           drifttime = TimeStruct(driftvelocity,0.5*kAmWidth-1.0*locT,zz)
1044                     + hittime;
1045         }
1046         else {
1047           // Use constant drift velocity
1048           drifttime = -1.0 * locT / driftvelocity
1049                     + hittime;
1050         }
1051
1052         // Apply the gas gain including fluctuations
1053         Double_t ggRndm = 0.0;
1054         do {
1055           ggRndm = gRandom->Rndm();
1056         } while (ggRndm <= 0);
1057         Int_t signal = (Int_t) (-(simParam->GetGasGain()) * TMath::Log(ggRndm));
1058
1059         // Apply the pad response
1060         if (simParam->PRFOn()) {
1061           // The distance of the electron to the center of the pad
1062           // in units of pad width
1063           Double_t dist = (colOffset - 0.5*padPlane->GetColSize(colE))
1064                         / padPlane->GetColSize(colE);
1065           // This is still the fixed parametrization, i.e. not dependent on
1066           // calibration values !!!!
1067           if (!(calibration->PadResponse(signal,dist,plane,padSignal))) continue;
1068         }
1069         else {
1070           padSignal[0] = 0.0;
1071           padSignal[1] = signal;
1072           padSignal[2] = 0.0;
1073         }
1074
1075         // The time bin (always positive), with t0 distortion
1076         Double_t timeBinIdeal = drifttime * samplingRate + t0;
1077         // Protection
1078         if (TMath::Abs(timeBinIdeal) > 2*nTimeTotal) {
1079           timeBinIdeal = 2 * nTimeTotal;
1080         }
1081         Int_t    timeBinTruncated = (Int_t) timeBinIdeal;
1082         // The distance of the position to the middle of the timebin
1083         Double_t timeOffset       = ((Float_t) timeBinTruncated
1084                                   + 0.5 - timeBinIdeal) / samplingRate;
1085
1086         // Sample the time response inside the drift region
1087         // + additional time bins before and after.
1088         // The sampling is done always in the middle of the time bin
1089         for (Int_t iTimeBin = TMath::Max(timeBinTruncated,0)
1090             ;iTimeBin < TMath::Min(timeBinTruncated+timeBinTRFend,nTimeTotal)
1091             ;iTimeBin++) {
1092
1093           // Apply the time response
1094           Double_t timeResponse = 1.0;
1095           Double_t crossTalk    = 0.0;
1096           Double_t time         = (iTimeBin - timeBinTruncated) / samplingRate + timeOffset;
1097           if (simParam->TRFOn()) {
1098             timeResponse = simParam->TimeResponse(time);
1099           }
1100           if (simParam->CTOn()) {
1101             crossTalk    = simParam->CrossTalk(time);
1102           }
1103
1104           signalOld[0] = 0.0;
1105           signalOld[1] = 0.0;
1106           signalOld[2] = 0.0;
1107
1108           for (iPad = 0; iPad < kNpad; iPad++) {
1109
1110             Int_t colPos = colE + iPad - 1;
1111             if (colPos <        0) continue;
1112             if (colPos >= nColMax) break;
1113
1114             // Add the signals
1115             Int_t iCurrentTimeBin = iTimeBin;
1116             signalOld[iPad]  = signals->GetDataUnchecked(rowE,colPos,iCurrentTimeBin);
1117             if (colPos != colE) {
1118               signalOld[iPad] += padSignal[iPad] * (timeResponse + crossTalk);
1119             }
1120             else {
1121               signalOld[iPad] += padSignal[iPad] * timeResponse;
1122             }
1123             signals->SetDataUnchecked(rowE,colPos,iCurrentTimeBin,signalOld[iPad]);
1124
1125             // Store the track index in the dictionary
1126             // Note: We store index+1 in order to allow the array to be compressed
1127             if (signalOld[iPad] > 0) {
1128               for (iDict = 0; iDict < kNDict; iDict++) {
1129                 Int_t oldTrack = dictionary[iDict]->GetDataUnchecked(rowE
1130                                                                     ,colPos
1131                                                                     ,iCurrentTimeBin);
1132                 if (oldTrack == track+1) break;
1133                 if (oldTrack ==       0) {
1134                   dictionary[iDict]->SetDataUnchecked(rowE,colPos,iCurrentTimeBin,track+1);
1135                   break;
1136                 }
1137               }
1138             }
1139
1140           } // Loop: pads
1141
1142         } // Loop: time bins
1143
1144       } // Loop: electrons of a single hit
1145
1146       hit = (AliTRDhit *) fTRD->NextHit();
1147
1148     } // Loop: hits of one primary track
1149
1150   } // Loop: primary tracks
1151
1152   AliDebug(1,Form("Finished analyzing %d hits",countHits));
1153
1154   //____________________________________________________________________
1155   //
1156   // Create (s)digits
1157   //____________________________________________________________________
1158   //
1159
1160   // The coupling factor
1161   Double_t coupling = simParam->GetPadCoupling()
1162                     * simParam->GetTimeCoupling();
1163
1164   // The conversion factor
1165   Double_t convert  = kEl2fC
1166                     * simParam->GetChipGain();
1167
1168   // Loop through all chambers to finalize the digits
1169   Int_t iDetBeg = 0;
1170   Int_t iDetEnd = AliTRDgeometry::Ndet();
1171   for (Int_t iDet = iDetBeg; iDet < iDetEnd; iDet++) {
1172
1173     Int_t plane      = fGeo->GetPlane(iDet);
1174     Int_t sector     = fGeo->GetSector(iDet);
1175     Int_t chamber    = fGeo->GetChamber(iDet);
1176     Int_t nRowMax    = fGeo->GetRowMax(plane,chamber,sector);
1177     Int_t nColMax    = fGeo->GetColMax(plane);
1178
1179     Double_t *inADC  = new Double_t[nTimeTotal];
1180     Double_t *outADC = new Double_t[nTimeTotal];
1181
1182     // Add a container for the digits of this detector
1183     digits = fDigitsManager->GetDigits(iDet);
1184     // Allocate memory space for the digits buffer
1185     if (digits->GetNtime() == 0) {
1186       digits->Allocate(nRowMax,nColMax,nTimeTotal);
1187     }
1188
1189     // Get the signal container
1190     signals = (AliTRDdataArrayF *) signalsArray->At(iDet);
1191     if (signals->GetNtime() == 0) {
1192       // Create missing containers
1193       signals->Allocate(nRowMax,nColMax,nTimeTotal);
1194     }
1195     else {
1196       // Expand the container if neccessary
1197       if (fCompress) signals->Expand();
1198     }
1199     // Create the missing dictionary containers
1200     for (iDict = 0; iDict < kNDict; iDict++) {
1201       dictionary[iDict] = fDigitsManager->GetDictionary(iDet,iDict);
1202       if (dictionary[iDict]->GetNtime() == 0) {
1203         dictionary[iDict]->Allocate(nRowMax,nColMax,nTimeTotal);
1204       }
1205     }
1206
1207     Int_t nDigits = 0;
1208
1209     // Don't create noise in detectors that are switched off / not installed, etc.
1210     if (( calibration->IsChamberInstalled(iDet)) &&
1211         (!calibration->IsChamberMasked(iDet))    &&
1212         ( fGeo->GetSMstatus(sector))) {
1213
1214       // Get the calibration objects
1215       calGainFactorROC      = calibration->GetGainFactorROC(iDet);
1216       calGainFactorDetValue = calGainFactorDet->GetValue(iDet);
1217
1218       // Create the digits for this chamber
1219       for (iRow  = 0; iRow  <  nRowMax;   iRow++ ) {
1220         for (iCol  = 0; iCol  <  nColMax;   iCol++ ) {
1221
1222           // Check whether pad is masked
1223           // Bridged pads are not considered yet!!!
1224           if (calibration->IsPadMasked(iDet,iCol,iRow)) {
1225             continue;
1226           }
1227
1228           // Create summable digits
1229           if (fSDigits) {
1230
1231             for (iTime = 0; iTime < nTimeTotal; iTime++) {
1232               Float_t signalAmp = signals->GetDataUnchecked(iRow,iCol,iTime);
1233               signalAmp *= fSDigitsScale;
1234               signalAmp  = TMath::Min(signalAmp,(Float_t) 1.0e9);
1235               Int_t adc  = (Int_t) signalAmp;
1236               if (adc > 0) nDigits++;
1237               digits->SetDataUnchecked(iRow,iCol,iTime,adc);
1238             }
1239
1240           }
1241           // Create normal digits
1242           else {
1243
1244             Float_t padgain = calGainFactorDetValue
1245                             * calGainFactorROC->GetValue(iCol,iRow);
1246             if (padgain <= 0) {
1247               AliError(Form("Not a valid gain %f, %d %d %d",padgain,iDet,iCol,iRow));
1248             }
1249
1250             for (iTime = 0; iTime < nTimeTotal; iTime++) {
1251
1252               Float_t signalAmp = signals->GetDataUnchecked(iRow,iCol,iTime);
1253
1254               // Pad and time coupling
1255               signalAmp *= coupling;
1256               // Gain factors
1257               signalAmp *= padgain;
1258
1259               // Add the noise, starting from minus ADC baseline in electrons
1260               Double_t baselineEl = simParam->GetADCbaseline()
1261                                   * (simParam->GetADCinRange() / simParam->GetADCoutRange())
1262                                   / convert;
1263               signalAmp  = TMath::Max((Double_t) gRandom->Gaus(signalAmp,simParam->GetNoise())
1264                                      ,-baselineEl);
1265               // Convert to mV
1266               signalAmp *= convert;
1267               // Add ADC baseline in mV
1268               signalAmp += simParam->GetADCbaseline()
1269                          * (simParam->GetADCinRange() / simParam->GetADCoutRange());
1270               // Convert to ADC counts. Set the overflow-bit fADCoutRange if the
1271               // signal is larger than fADCinRange
1272               Int_t adc  = 0;
1273               if (signalAmp >= simParam->GetADCinRange()) {
1274                 adc = ((Int_t) simParam->GetADCoutRange());
1275               }
1276               else {
1277                 adc = TMath::Nint(signalAmp * (simParam->GetADCoutRange()
1278                                             / simParam->GetADCinRange()));
1279               }
1280
1281               inADC[iTime]  = adc;
1282               outADC[iTime] = adc;
1283
1284             }
1285
1286             for (iTime = 0; iTime < nTimeTotal; iTime++) {
1287               // Store the amplitude of the digit if above threshold
1288               if (outADC[iTime] > (simParam->GetADCbaseline() + simParam->GetADCthreshold())) {
1289                 nDigits++;
1290                 digits->SetDataUnchecked(iRow,iCol,iTime,((Int_t) outADC[iTime]));
1291               }
1292             }
1293
1294           }
1295
1296         }
1297       }
1298
1299     }
1300
1301     // Compress the arrays
1302     digits->Compress(1,0);
1303     for (iDict = 0; iDict < kNDict; iDict++) {
1304       dictionary[iDict]->Compress(1,0);
1305     }
1306
1307     totalSizeDigits += digits->GetSize();
1308     totalSizeDict0  += dictionary[0]->GetSize();
1309     totalSizeDict1  += dictionary[1]->GetSize();
1310     totalSizeDict2  += dictionary[2]->GetSize();
1311
1312     if (nDigits > 0) {
1313       Float_t nPixel = nRowMax * nColMax * nTimeTotal;
1314       AliDebug(1,Form("Found %d digits in detector %d (%3.0f)."
1315                      ,nDigits,iDet
1316                      ,100.0 * ((Float_t) nDigits) / nPixel));
1317     }
1318
1319     if (fCompress) {
1320       signals->Compress(1,0);
1321     }
1322
1323     delete [] inADC;
1324     delete [] outADC;
1325
1326   }
1327
1328   if (signalsArray) {
1329     delete signalsArray;
1330     signalsArray = 0;
1331   }
1332
1333   AliDebug(1,Form("Total number of analyzed hits = %d",countHits));
1334   AliDebug(1,Form("Total digits data size = %d, %d, %d, %d",totalSizeDigits
1335                                                            ,totalSizeDict0
1336                                                            ,totalSizeDict1
1337                                                            ,totalSizeDict2));
1338
1339   return kTRUE;
1340
1341 }
1342
1343 //_____________________________________________________________________________
1344 void AliTRDdigitizer::AddSDigitsManager(AliTRDdigitsManager *man)
1345 {
1346   //
1347   // Add a digits manager for s-digits to the input list.
1348   //
1349
1350   fSDigitsManagerList->Add(man);
1351
1352 }
1353
1354 //_____________________________________________________________________________
1355 void AliTRDdigitizer::DeleteSDigitsManager()
1356 {
1357   //
1358   // Removes digits manager from the input list.
1359   //
1360
1361   fSDigitsManagerList->Delete();
1362
1363 }
1364
1365 //_____________________________________________________________________________
1366 Bool_t AliTRDdigitizer::ConvertSDigits()
1367 {
1368   //
1369   // Converts s-digits to normal digits
1370   //
1371
1372   // Number of track dictionary arrays
1373   const Int_t    kNDict = AliTRDdigitsManager::kNDict;
1374
1375   // Converts number of electrons to fC
1376   const Double_t kEl2fC = 1.602e-19 * 1.0e15;
1377
1378   Int_t iDict = 0;
1379   Int_t iRow;
1380   Int_t iCol;
1381   Int_t iTime;
1382
1383   AliTRDCalROC *calGainFactorROC      = 0;
1384   Float_t       calGainFactorDetValue = 0.0;
1385
1386   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
1387   if (!simParam) {
1388     AliFatal("Could not get simulation parameters");
1389     return kFALSE;
1390   }
1391
1392   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
1393   if (!calibration) {
1394     AliFatal("Could not get calibration object");
1395     return kFALSE;
1396   }
1397
1398   Double_t sDigitsScale = 1.0 / GetSDigitsScale();
1399   Double_t noise        = simParam->GetNoise();
1400   Double_t padCoupling  = simParam->GetPadCoupling();
1401   Double_t timeCoupling = simParam->GetTimeCoupling();
1402   Double_t chipGain     = simParam->GetChipGain();
1403   Double_t coupling     = padCoupling * timeCoupling;
1404   Double_t convert      = kEl2fC * chipGain;
1405   Double_t adcInRange   = simParam->GetADCinRange();
1406   Double_t adcOutRange  = simParam->GetADCoutRange();
1407   Int_t    adcThreshold = simParam->GetADCthreshold();
1408   Int_t    adcBaseline  = simParam->GetADCbaseline();
1409   Int_t    nTimeTotal   = calibration->GetNumberOfTimeBins();
1410
1411   AliTRDdataArrayI *digitsIn;
1412   AliTRDdataArrayI *digitsOut;
1413   AliTRDdataArrayI *dictionaryIn[kNDict];
1414   AliTRDdataArrayI *dictionaryOut[kNDict];
1415
1416   // Get the detector wise calibration objects
1417   const AliTRDCalDet *calGainFactorDet = calibration->GetGainFactorDet();
1418
1419   // Loop through the detectors
1420   for (Int_t iDet = 0; iDet < AliTRDgeometry::Ndet(); iDet++) {
1421
1422     Int_t plane      = fGeo->GetPlane(iDet);
1423     Int_t sector     = fGeo->GetSector(iDet);
1424     Int_t chamber    = fGeo->GetChamber(iDet);
1425     Int_t nRowMax    = fGeo->GetRowMax(plane,chamber,sector);
1426     Int_t nColMax    = fGeo->GetColMax(plane);
1427
1428     Double_t *inADC  = new Double_t[nTimeTotal];
1429     Double_t *outADC = new Double_t[nTimeTotal];
1430
1431     digitsIn  = fSDigitsManager->GetDigits(iDet);
1432     digitsIn->Expand();
1433     digitsOut = fDigitsManager->GetDigits(iDet);
1434     digitsOut->Allocate(nRowMax,nColMax,nTimeTotal);
1435     for (iDict = 0; iDict < kNDict; iDict++) {
1436       dictionaryIn[iDict]  = fSDigitsManager->GetDictionary(iDet,iDict);
1437       dictionaryIn[iDict]->Expand();
1438       dictionaryOut[iDict] = fDigitsManager->GetDictionary(iDet,iDict);
1439       dictionaryOut[iDict]->Allocate(nRowMax,nColMax,nTimeTotal);
1440     }
1441
1442     // Don't create noise in detectors that are switched off / not installed, etc.
1443     if (( calibration->IsChamberInstalled(iDet)) &&
1444         (!calibration->IsChamberMasked(iDet))    &&
1445         ( fGeo->GetSMstatus(sector))) {
1446
1447       // Get the calibration objects
1448       calGainFactorROC      = calibration->GetGainFactorROC(iDet);
1449       calGainFactorDetValue = calGainFactorDet->GetValue(iDet);
1450
1451       for (iRow  = 0; iRow  <  nRowMax;   iRow++ ) {
1452         for (iCol  = 0; iCol  <  nColMax;   iCol++ ) {
1453
1454           // Check whether pad is masked
1455           // Bridged pads are not considered yet!!!
1456           if (calibration->IsPadMasked(iDet,iCol,iRow)) {
1457             continue;
1458           }
1459
1460           Float_t padgain = calGainFactorDetValue
1461                           * calGainFactorROC->GetValue(iCol,iRow);
1462           if (padgain <= 0) {
1463             AliError(Form("Not a valid gain %f, %d %d %d",padgain,iDet,iCol,iRow));
1464           }
1465
1466           for (iTime = 0; iTime < nTimeTotal; iTime++) {
1467
1468             // Scale s-digits to normal digits
1469             Double_t signal = (Double_t) digitsIn->GetDataUnchecked(iRow,iCol,iTime);
1470             signal         *= sDigitsScale;
1471             // Pad and time coupling
1472             signal *= coupling;
1473             // Gain factors
1474             signal *= padgain;
1475             // Add the noise, starting from minus ADC baseline in electrons
1476             Double_t baselineEl = adcBaseline * (adcInRange / adcOutRange) / convert;
1477             signal  = TMath::Max((Double_t) gRandom->Gaus(signal,noise),-baselineEl);
1478             // Convert to mV
1479             signal *= convert;
1480             // add ADC baseline in mV
1481             signal += adcBaseline * (adcInRange / adcOutRange);
1482             // Convert to ADC counts. Set the overflow-bit adcOutRange if the
1483             // signal is larger than adcInRange
1484             Int_t adc  = 0;
1485             if (signal >= adcInRange) {
1486               adc = ((Int_t) adcOutRange);
1487             }
1488             else {
1489               adc = TMath::Nint(signal * (adcOutRange / adcInRange));
1490             }
1491             inADC[iTime]  = adc;
1492             outADC[iTime] = adc;
1493
1494           }
1495
1496           for (iTime = 0; iTime < nTimeTotal; iTime++) {
1497             // Store the amplitude of the digit if above threshold
1498             if (outADC[iTime] > (adcBaseline + adcThreshold)) {
1499               digitsOut->SetDataUnchecked(iRow,iCol,iTime,((Int_t) outADC[iTime]));
1500               // Copy the dictionary
1501               for (iDict = 0; iDict < kNDict; iDict++) {
1502                 Int_t track = dictionaryIn[iDict]->GetDataUnchecked(iRow,iCol,iTime);
1503                 dictionaryOut[iDict]->SetDataUnchecked(iRow,iCol,iTime,track);
1504               }
1505             }
1506           }
1507
1508         }
1509       }
1510
1511     }
1512
1513     if (fCompress) {
1514       digitsIn->Compress(1,0);
1515       digitsOut->Compress(1,0);
1516       for (iDict = 0; iDict < kNDict; iDict++) {
1517         dictionaryIn[iDict]->Compress(1,0);
1518         dictionaryOut[iDict]->Compress(1,0);
1519       }
1520     }
1521
1522     delete [] inADC;
1523     delete [] outADC;
1524
1525   }
1526
1527   return kTRUE;
1528
1529 }
1530
1531 //_____________________________________________________________________________
1532 Bool_t AliTRDdigitizer::MergeSDigits()
1533 {
1534   //
1535   // Merges the input s-digits:
1536   //   - The amplitude of the different inputs are summed up.
1537   //   - Of the track IDs from the input dictionaries only one is
1538   //     kept for each input. This works for maximal 3 different merged inputs.
1539   //
1540
1541   // Number of track dictionary arrays
1542   const Int_t kNDict = AliTRDdigitsManager::kNDict;
1543
1544   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
1545   if (!simParam) {
1546     AliFatal("Could not get simulation parameters");
1547     return kFALSE;
1548   }
1549
1550   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
1551   if (!calibration) {
1552     AliFatal("Could not get calibration object");
1553     return kFALSE;
1554   }
1555
1556   Int_t iDict = 0;
1557   Int_t jDict = 0;
1558
1559   AliTRDdataArrayI *digitsA;
1560   AliTRDdataArrayI *digitsB;
1561   AliTRDdataArrayI *dictionaryA[kNDict];
1562   AliTRDdataArrayI *dictionaryB[kNDict];
1563
1564   // Get the first s-digits
1565   fSDigitsManager = (AliTRDdigitsManager *) fSDigitsManagerList->First();
1566   if (!fSDigitsManager) {
1567     AliError("No SDigits manager");
1568     return kFALSE;
1569   }
1570
1571   // Loop through the other sets of s-digits
1572   AliTRDdigitsManager *mergeSDigitsManager;
1573   mergeSDigitsManager = (AliTRDdigitsManager *)
1574                         fSDigitsManagerList->After(fSDigitsManager);
1575
1576   if (mergeSDigitsManager) {
1577     AliDebug(1,Form("Merge %d input files.",fSDigitsManagerList->GetSize()));
1578   }
1579   else {
1580     AliDebug(1,"Only one input file.");
1581   }
1582
1583   Int_t nTimeTotal = calibration->GetNumberOfTimeBins();
1584
1585   Int_t iMerge = 0;
1586   while (mergeSDigitsManager) {
1587
1588     iMerge++;
1589
1590     // Loop through the detectors
1591     for (Int_t iDet = 0; iDet < AliTRDgeometry::Ndet(); iDet++) {
1592
1593       Int_t plane   = fGeo->GetPlane(iDet);
1594       Int_t sector  = fGeo->GetSector(iDet);
1595       Int_t chamber = fGeo->GetChamber(iDet);
1596       Int_t nRowMax = fGeo->GetRowMax(plane,chamber,sector);
1597       Int_t nColMax = fGeo->GetColMax(plane);
1598
1599       // Loop through the pixels of one detector and add the signals
1600       digitsA = fSDigitsManager->GetDigits(iDet);
1601       digitsB = mergeSDigitsManager->GetDigits(iDet);
1602       digitsA->Expand();
1603       digitsB->Expand();
1604       for (iDict = 0; iDict < kNDict; iDict++) {
1605         dictionaryA[iDict] = fSDigitsManager->GetDictionary(iDet,iDict);
1606         dictionaryB[iDict] = mergeSDigitsManager->GetDictionary(iDet,iDict);
1607         dictionaryA[iDict]->Expand();
1608         dictionaryB[iDict]->Expand();
1609       }
1610
1611       // Merge only detectors that contain a signal
1612       Bool_t doMerge = kTRUE;
1613       if (fMergeSignalOnly) {
1614         if (digitsA->GetOverThreshold(0) == 0) {
1615           doMerge = kFALSE;
1616         }
1617       }
1618
1619       if (doMerge) {
1620
1621         AliDebug(1,Form("Merge detector %d of input no.%d",iDet,iMerge+1));
1622
1623         for (Int_t iRow  = 0; iRow  <  nRowMax;   iRow++ ) {
1624           for (Int_t iCol  = 0; iCol  <  nColMax;   iCol++ ) {
1625             for (Int_t iTime = 0; iTime < nTimeTotal; iTime++) {
1626
1627               // Add the amplitudes of the summable digits
1628               Int_t ampA = digitsA->GetDataUnchecked(iRow,iCol,iTime);
1629               Int_t ampB = digitsB->GetDataUnchecked(iRow,iCol,iTime);
1630               ampA += ampB;
1631               digitsA->SetDataUnchecked(iRow,iCol,iTime,ampA);
1632
1633               // Add the mask to the track id if defined.
1634               for (iDict = 0; iDict < kNDict; iDict++) {
1635                 Int_t trackB = dictionaryB[iDict]->GetDataUnchecked(iRow,iCol,iTime);
1636                 if ((fMasks) && (trackB > 0)) {
1637                   for (jDict = 0; jDict < kNDict; jDict++) {
1638                     Int_t trackA = dictionaryA[iDict]->GetDataUnchecked(iRow,iCol,iTime);
1639                     if (trackA == 0) {
1640                       trackA = trackB + fMasks[iMerge];
1641                       dictionaryA[iDict]->SetDataUnchecked(iRow,iCol,iTime,trackA);
1642                     }
1643                   }
1644                 }
1645               }
1646
1647             }
1648           }
1649         }
1650
1651       }
1652
1653       if (fCompress) {
1654         digitsA->Compress(1,0);
1655         digitsB->Compress(1,0);
1656         for (iDict = 0; iDict < kNDict; iDict++) {
1657           dictionaryA[iDict]->Compress(1,0);
1658           dictionaryB[iDict]->Compress(1,0);
1659         }
1660       }
1661
1662     }
1663
1664     // The next set of s-digits
1665     mergeSDigitsManager = (AliTRDdigitsManager *)
1666                           fSDigitsManagerList->After(mergeSDigitsManager);
1667
1668   }
1669
1670   return kTRUE;
1671
1672 }
1673
1674 //_____________________________________________________________________________
1675 Bool_t AliTRDdigitizer::SDigits2Digits()
1676 {
1677   //
1678   // Merges the input s-digits and converts them to normal digits
1679   //
1680
1681   if (!MergeSDigits()) {
1682     return kFALSE;
1683   }
1684
1685   return ConvertSDigits();
1686
1687 }
1688
1689 //_____________________________________________________________________________
1690 Bool_t AliTRDdigitizer::WriteDigits() const
1691 {
1692   //
1693   // Writes out the TRD-digits and the dictionaries
1694   //
1695
1696   // Write parameters
1697   fRunLoader->CdGAFile();
1698
1699   // Store the digits and the dictionary in the tree
1700   return fDigitsManager->WriteDigits();
1701
1702 }
1703
1704 //_____________________________________________________________________________
1705 void AliTRDdigitizer::InitOutput(Int_t iEvent)
1706 {
1707   //
1708   // Initializes the output branches
1709   //
1710
1711   fEvent = iEvent;
1712
1713   if (!fRunLoader) {
1714     AliError("Run Loader is NULL");
1715     return;
1716   }
1717
1718   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
1719   if (!loader) {
1720     AliError("Can not get TRD loader from Run Loader");
1721     return;
1722   }
1723
1724   TTree *tree = 0;
1725
1726   if (fSDigits) {
1727     // If we produce SDigits
1728     tree = loader->TreeS();
1729     if (!tree) {
1730       loader->MakeTree("S");
1731       tree = loader->TreeS();
1732     }
1733   }
1734   else {
1735     // If we produce Digits
1736     tree = loader->TreeD();
1737     if (!tree) {
1738       loader->MakeTree("D");
1739       tree = loader->TreeD();
1740     }
1741   }
1742
1743   fDigitsManager->SetEvent(iEvent);
1744   fDigitsManager->MakeBranch(tree);
1745
1746 }
1747
1748 //_____________________________________________________________________________
1749 Double_t AliTRDdigitizer::TimeStruct(Float_t vdrift, Double_t dist, Double_t z)
1750 {
1751   //
1752   // Applies the time structure of the drift cells (by C.Lippmann).
1753   // The drift time of electrons to the anode wires depends on the
1754   // distance to the wire (z) and on the position in the drift region.
1755   //
1756   // input :
1757   // dist = radial distance from (cathode) pad plane [cm]
1758   // z    = distance from anode wire (parallel to cathode planes) [cm]
1759   //
1760   // output :
1761   // tdrift = the drift time of an electron at the given position
1762   //
1763   // We interpolate between the drift time values at the two drift
1764   // velocities fVDlo and fVDhi, being smaller and larger than
1765   // fDriftVelocity. We use the two stored drift time maps fTimeStruct1
1766   // and fTimeStruct2, calculated for the two mentioned drift velocities.
1767   //
1768
1769   SampleTimeStruct(vdrift);
1770
1771   // Indices:
1772   Int_t r1 = (Int_t)(10 * dist);
1773   if (r1 <  0) r1 =  0;
1774   if (r1 > 37) r1 = 37;
1775   Int_t r2 = r1 + 1;
1776   if (r2 <  0) r2 =  0;
1777   if (r2 > 37) r2 = 37;
1778   const Int_t kz1 = ((Int_t)(100 * z / 2.5));
1779   const Int_t kz2 = kz1 + 1;
1780
1781   if ((r1  <  0) ||
1782       (r1  > 37) ||
1783       (kz1 <  0) ||
1784       (kz1 > 10)) {
1785     AliWarning(Form("Indices out of range: dist=%.2f, z=%.2f, r1=%d, kz1=%d"
1786                    ,dist,z,r1,kz1));
1787   }
1788
1789   const Float_t ky111 = fTimeStruct1[r1+38*kz1];
1790   const Float_t ky221 = ((r2 <= 37) && (kz2 <= 10))
1791                       ? fTimeStruct1[r2+38*kz2]
1792                       : fTimeStruct1[37+38*10];
1793   const Float_t ky121 = (kz2 <= 10)
1794                       ? fTimeStruct1[r1+38*kz2]
1795                       : fTimeStruct1[r1+38*10];
1796   const Float_t ky211 = (r2 <= 37)
1797                       ? fTimeStruct1[r2+38*kz1]
1798                       : fTimeStruct1[37+38*kz1];
1799
1800   // 2D Interpolation, lower drift time map
1801   const Float_t ky11  = (ky211-ky111)*10*dist + ky111 - (ky211-ky111)*r1;
1802   const Float_t ky21  = (ky221-ky121)*10*dist + ky121 - (ky221-ky121)*r1;
1803
1804   const Float_t ky112 = fTimeStruct2[r1+38*kz1];
1805   const Float_t ky222 = ((r2 <= 37) && (kz2 <= 10))
1806                       ? fTimeStruct2[r2+38*kz2]
1807                       : fTimeStruct2[37+38*10];
1808   const Float_t ky122 = (kz2 <= 10)
1809                       ? fTimeStruct2[r1+38*kz2]
1810                       : fTimeStruct2[r1+38*10];
1811   const Float_t ky212 = (r2 <= 37)
1812                       ? fTimeStruct2[r2+38*kz1]
1813                       : fTimeStruct2[37+38*kz1];
1814
1815   // 2D Interpolation, larger drift time map
1816   const Float_t ky12  = (ky212-ky112)*10*dist + ky112 - (ky212-ky112)*r1;
1817   const Float_t ky22  = (ky222-ky122)*10*dist + ky122 - (ky222-ky122)*r1;
1818
1819   // Dist now is the drift distance to the anode wires (negative if electrons are
1820   // between anode wire plane and cathode pad plane)
1821   dist -= AliTRDgeometry::AmThick() / 2.0;
1822
1823   // Get the drift times for the drift velocities fVDlo and fVDhi
1824   const Float_t ktdrift1 = ((TMath::Abs(dist) > 0.005) || (z > 0.005))
1825                          ? (ky21 - ky11) * 100 * z / 2.5 + ky11 - (ky21 - ky11) * kz1
1826                          : 0.0;
1827   const Float_t ktdrift2 = ((TMath::Abs(dist) > 0.005) || (z > 0.005))
1828                          ? (ky22 - ky12) * 100 * z / 2.5 + ky12 - (ky22 - ky12) * kz1
1829                          : 0.0;
1830
1831   // 1D Interpolation between the values at fVDlo and fVDhi
1832   Float_t a = (ktdrift2 - ktdrift1)
1833             / (fVDhi - fVDlo);
1834   Float_t b = ktdrift2 - a * fVDhi;
1835
1836   return a * vdrift + b;
1837
1838 }
1839
1840 //_____________________________________________________________________________
1841 void AliTRDdigitizer::SampleTimeStruct(Float_t vdrift)
1842 {
1843   //
1844   // Samples the timing structure of a drift cell
1845   // Drift Time data calculated with Garfield (by C.Lippmann)
1846   //
1847
1848   // Nothing to do
1849   if (vdrift == fTimeLastVdrift) {
1850     return;
1851   }
1852   fTimeLastVdrift = vdrift;
1853
1854   // Drift time maps are saved for some drift velocity values (in drift region):
1855   Float_t  fVDsmp[8];
1856   fVDsmp[0] = 1.032;
1857   fVDsmp[1] = 1.158;
1858   fVDsmp[2] = 1.299;
1859   fVDsmp[3] = 1.450;
1860   fVDsmp[4] = 1.610;
1861   fVDsmp[5] = 1.783;
1862   fVDsmp[6] = 1.959;
1863   fVDsmp[7] = 2.134;
1864
1865   if      (vdrift < fVDsmp[0]) {
1866     AliWarning(Form("Drift Velocity too small (%.3f<%.3f)",vdrift,fVDsmp[0]));
1867     vdrift = fVDsmp[0];
1868   }
1869   else if (vdrift > fVDsmp[7]) {
1870     AliWarning(Form("Drift Velocity too large (%.3f>%.3f)",vdrift,fVDsmp[6]));
1871     vdrift = fVDsmp[7];
1872   }
1873
1874   const Int_t ktimebin  = 38;
1875   const Int_t kZbin     = 11;
1876
1877   // Garfield simulation at UD = -1500V; vd = 1.032cm/microsec, <driftfield> = 525V/cm
1878   Float_t time1500[ktimebin][kZbin] =
1879       {{0.09170, 0.09205, 0.09306, 0.09475, 0.09716, 0.10035,
1880         0.10445, 0.10993, 0.11838, 0.13986, 0.37858},
1881        {0.06588, 0.06626, 0.06739, 0.06926, 0.07186, 0.07524,
1882         0.07951, 0.08515, 0.09381, 0.11601, 0.35673},
1883        {0.03946, 0.04003, 0.04171, 0.04435, 0.04780, 0.05193,
1884         0.05680, 0.06306, 0.07290, 0.09873, 0.34748},
1885        {0.01151, 0.01283, 0.01718, 0.02282, 0.02880, 0.03479,
1886         0.04098, 0.04910, 0.06413, 0.10567, 0.36897},
1887        {0.01116, 0.01290, 0.01721, 0.02299, 0.02863, 0.03447,
1888         0.04074, 0.04984, 0.06839, 0.11625, 0.37277},
1889        {0.03919, 0.03974, 0.04131, 0.04380, 0.04703, 0.05102,
1890         0.05602, 0.06309, 0.07651, 0.10938, 0.36838},
1891        {0.06493, 0.06560, 0.06640, 0.06802, 0.07051, 0.07392,
1892         0.07853, 0.08510, 0.09690, 0.12621, 0.38058},
1893        {0.09174, 0.09186, 0.09225, 0.09303, 0.09477, 0.00000,
1894         0.11205, 0.11952, 0.13461, 0.16984, 0.43017},
1895        {0.14356, 0.14494, 0.14959, 0.16002, 0.18328, 0.27981,
1896         0.22785, 0.21240, 0.21948, 0.25965, 0.52392},
1897        {0.23120, 0.23366, 0.24046, 0.25422, 0.28071, 0.36914,
1898         0.32999, 0.31208, 0.31772, 0.35804, 0.62249},
1899        {0.32686, 0.32916, 0.33646, 0.35053, 0.37710, 0.46292,
1900         0.42773, 0.40948, 0.41497, 0.45527, 0.71955},
1901        {0.42353, 0.42583, 0.43317, 0.44727, 0.47380, 0.55884,
1902         0.52479, 0.50650, 0.51194, 0.55225, 0.81658},
1903        {0.52038, 0.52271, 0.53000, 0.54415, 0.57064, 0.65545,
1904         0.62172, 0.60341, 0.60885, 0.64915, 0.91339},
1905        {0.61724, 0.61953, 0.62694, 0.64098, 0.66756, 0.75226,
1906         0.71862, 0.70030, 0.70575, 0.74604, 1.01035},
1907        {0.71460, 0.71678, 0.72376, 0.73786, 0.76447, 0.84913,
1908         0.81551, 0.79720, 0.80264, 0.84292, 1.10723},
1909        {0.81101, 0.81334, 0.82066, 0.83475, 0.86127, 0.94599,
1910         0.91240, 0.89408, 0.89952, 0.93981, 1.20409},
1911        {0.90788, 0.91023, 0.91752, 0.93163, 0.95815, 1.04293,
1912         1.00929, 0.99096, 0.99640, 1.03669, 1.30106},
1913        {1.00477, 1.00707, 1.01449, 1.02852, 1.05504, 1.13976,
1914         1.10617, 1.08784, 1.09329, 1.13358, 1.39796},
1915        {1.10166, 1.10397, 1.11130, 1.12541, 1.15257, 1.23672,
1916         1.20307, 1.18472, 1.19018, 1.23046, 1.49486},
1917        {1.19854, 1.20084, 1.20818, 1.22235, 1.24885, 1.33355,
1918         1.29992, 1.28156, 1.28707, 1.32735, 1.59177},
1919        {1.29544, 1.29780, 1.30507, 1.31917, 1.34575, 1.43073,
1920         1.39681, 1.37851, 1.38396, 1.42377, 1.68854},
1921        {1.39236, 1.39462, 1.40205, 1.41607, 1.44259, 1.52745,
1922         1.49368, 1.47541, 1.48083, 1.52112, 1.78546},
1923        {1.49314, 1.49149, 1.49885, 1.51297, 1.53949, 1.62420,
1924         1.59016, 1.57231, 1.57772, 1.61800, 1.88048},
1925        {1.58610, 1.58839, 1.59572, 1.60983, 1.63635, 1.72109,
1926         1.68651, 1.66921, 1.67463, 1.71489, 1.97918},
1927        {1.68400, 1.68529, 1.69261, 1.70671, 1.73331, 1.81830,
1928         1.78341, 1.76605, 1.77150, 1.81179, 2.07608},
1929        {1.77991, 1.78215, 1.78952, 1.80385, 1.83014, 1.91486,
1930         1.88128, 1.86215, 1.86837, 1.90865, 2.17304},
1931        {1.87674, 1.87904, 1.88647, 1.90052, 1.92712, 2.01173,
1932         1.97812, 1.95905, 1.96527, 2.00710, 2.26979},
1933        {1.97369, 1.97594, 1.98326, 1.99869, 2.02442, 2.10865,
1934         2.07501, 2.05666, 2.06214, 2.10243, 2.36669},
1935        {2.07052, 2.07281, 2.08016, 2.09425, 2.12132, 2.20555,
1936         2.17182, 2.15341, 2.15904, 2.19933, 2.46363},
1937        {2.16742, 2.16971, 2.17707, 2.19114, 2.21766, 2.30240,
1938         2.26877, 2.25015, 2.25573, 2.29586, 2.56060},
1939        {2.26423, 2.26659, 2.27396, 2.28803, 2.31456, 2.40828,
1940         2.36567, 2.34705, 2.35282, 2.39765, 2.65744},
1941        {2.36153, 2.36349, 2.37330, 2.38501, 2.41159, 2.49940,
1942         2.46257, 2.44420, 2.44843, 2.48987, 2.75431},
1943        {2.46558, 2.46035, 2.46822, 2.48181, 2.50849, 2.59630,
1944         2.55947, 2.54112, 2.54513, 2.58677, 2.85094},
1945        {2.56248, 2.55723, 2.56486, 2.57871, 2.60520, 2.68998,
1946         2.65626, 2.63790, 2.64316, 2.68360, 2.94813},
1947        {2.65178, 2.65441, 2.66153, 2.67556, 2.70210, 2.78687,
1948         2.75319, 2.73463, 2.74032, 2.78060, 3.04503},
1949        {2.74868, 2.75131, 2.75870, 2.77245, 2.79385, 2.88700,
1950         2.85009, 2.83177, 2.83723, 2.87750, 3.14193},
1951        {2.84574, 2.84789, 2.85560, 2.86935, 2.89075, 2.98060,
1952         2.94576, 2.92868, 2.93356, 2.97436, 3.23868},
1953        {2.94239, 2.94469, 2.95221, 2.96625, 2.99345, 3.07747,
1954         3.04266, 3.02545, 3.03051, 3.07118, 3.33555}};
1955
1956   // Garfield simulation at UD = -1600V; vd = 1.158cm/microsec, <driftfield> = 558V/cm
1957   Float_t time1600[ktimebin][kZbin] =
1958       {{0.09169, 0.09203, 0.09305, 0.09473, 0.09714, 0.10032,
1959         0.10441, 0.10990, 0.11835, 0.13986, 0.37845},
1960        {0.06589, 0.06626, 0.06738, 0.06924, 0.07184, 0.07521,
1961         0.07947, 0.08512, 0.09379, 0.11603, 0.35648},
1962        {0.03947, 0.04003, 0.04171, 0.04434, 0.04778, 0.05190,
1963         0.05678, 0.06304, 0.07292, 0.09876, 0.34759},
1964        {0.01111, 0.01281, 0.01718, 0.02281, 0.02879, 0.03477,
1965         0.04097, 0.04910, 0.06415, 0.10573, 0.36896},
1966        {0.01120, 0.01311, 0.01721, 0.02279, 0.02862, 0.03446,
1967         0.04074, 0.04981, 0.06825, 0.11595, 0.37255},
1968        {0.03919, 0.03980, 0.04132, 0.04380, 0.04704, 0.05102,
1969         0.05602, 0.06302, 0.07633, 0.10896, 0.36743},
1970        {0.06531, 0.06528, 0.06631, 0.06805, 0.07053, 0.07392,
1971         0.07853, 0.08505, 0.09669, 0.12578, 0.37967},
1972        {0.09157, 0.09171, 0.09216, 0.09301, 0.09475, 0.00000,
1973         0.11152, 0.11879, 0.13352, 0.16802, 0.42750},
1974        {0.13977, 0.14095, 0.14509, 0.15433, 0.17534, 0.26406,
1975         0.21660, 0.20345, 0.21113, 0.25067, 0.51434},
1976        {0.21816, 0.22041, 0.22631, 0.23850, 0.26208, 0.34340,
1977         0.30755, 0.29237, 0.29878, 0.33863, 0.60258},
1978        {0.30344, 0.30547, 0.31241, 0.32444, 0.34809, 0.42696,
1979         0.39464, 0.37919, 0.38546, 0.42530, 0.68926},
1980        {0.38969, 0.39164, 0.39810, 0.41059, 0.43441, 0.51246,
1981         0.48112, 0.46562, 0.47191, 0.51172, 0.77558},
1982        {0.47592, 0.47799, 0.48442, 0.49689, 0.52061, 0.59855,
1983         0.56752, 0.55201, 0.55826, 0.59808, 0.86202},
1984        {0.56226, 0.56428, 0.57074, 0.58324, 0.60696, 0.68483,
1985         0.65388, 0.63837, 0.64461, 0.68445, 0.94830},
1986        {0.64881, 0.65063, 0.65709, 0.66958, 0.69331, 0.77117,
1987         0.74023, 0.72472, 0.73098, 0.77079, 1.03486},
1988        {0.73506, 0.73698, 0.74344, 0.75596, 0.77964, 0.85751,
1989         0.82658, 0.81107, 0.81731, 0.85712, 1.12106},
1990        {0.82132, 0.82333, 0.82979, 0.84228, 0.86608, 0.94386,
1991         0.91293, 0.89742, 0.90367, 0.94335, 1.20737},
1992        {0.90767, 0.90968, 0.91614, 0.92863, 0.95236, 1.03021,
1993         0.99928, 0.98377, 0.99001, 1.02984, 1.29371},
1994        {0.99410, 0.99602, 1.00257, 1.01498, 1.03869, 1.11720,
1995         1.08563, 1.07011, 1.07637, 1.11621, 1.37873},
1996        {1.08036, 1.08240, 1.08884, 1.10138, 1.12504, 1.20301,
1997         1.17198, 1.15647, 1.16272, 1.20255, 1.46651},
1998        {1.16670, 1.16872, 1.17525, 1.18783, 1.21139, 1.28934,
1999         1.25833, 1.24281, 1.24909, 1.28889, 1.55275},
2000        {1.25307, 1.25510, 1.26153, 1.27404, 1.29773, 1.37584,
2001         1.34469, 1.32916, 1.33536, 1.37524, 1.63915},
2002        {1.33942, 1.34146, 1.34788, 1.36040, 1.38410, 1.46438,
2003         1.43105, 1.41537, 1.42176, 1.46158, 1.72538},
2004        {1.42579, 1.42782, 1.43458, 1.44674, 1.47043, 1.55085,
2005         1.51675, 1.50168, 1.50810, 1.54793, 1.81174},
2006        {1.51207, 1.51454, 1.52060, 1.53307, 1.55684, 1.63478,
2007         1.60336, 1.58820, 1.59446, 1.63414, 1.89814},
2008        {1.59856, 1.60047, 1.60693, 1.61942, 1.64317, 1.72257,
2009         1.69008, 1.67454, 1.68080, 1.72063, 1.98433},
2010        {1.68481, 1.68682, 1.69330, 1.70584, 1.72949, 1.80752,
2011         1.77643, 1.76089, 1.76716, 1.80692, 2.07069},
2012        {1.77117, 1.77319, 1.77969, 1.79260, 1.81583, 1.89376,
2013         1.86226, 1.84720, 1.85355, 1.89256, 2.15343},
2014        {1.85748, 1.85967, 1.86605, 1.87848, 1.90222, 1.98010,
2015         1.94913, 1.93271, 1.93981, 1.97968, 2.24355},
2016        {1.94386, 1.94587, 1.95233, 1.96484, 1.98854, 2.06646,
2017         2.03542, 2.01755, 2.02617, 2.06604, 2.32993},
2018        {2.03022, 2.03230, 2.03868, 2.05134, 2.07488, 2.15367,
2019         2.12178, 2.10391, 2.11252, 2.15432, 2.41623},
2020        {2.11656, 2.11857, 2.12505, 2.13772, 2.16147, 2.23919,
2021         2.20817, 2.19265, 2.20744, 2.23872, 2.49996},
2022        {2.20291, 2.20611, 2.21137, 2.22387, 2.24758, 2.32563,
2023         2.29450, 2.27901, 2.28525, 2.32507, 2.58897},
2024        {2.28922, 2.29172, 2.29774, 2.31345, 2.33400, 2.41287,
2025         2.38086, 2.36535, 2.37160, 2.40869, 2.67113},
2026        {2.37572, 2.37764, 2.38410, 2.39803, 2.42046, 2.49817,
2027         2.46721, 2.45171, 2.45794, 2.49505, 2.76061},
2028        {2.46190, 2.46396, 2.47043, 2.48340, 2.50665, 2.58453,
2029         2.55357, 2.53728, 2.54430, 2.58407, 2.84816},
2030        {2.54833, 2.55032, 2.55679, 2.56976, 2.59312, 2.67103,
2031         2.63993, 2.62364, 2.63062, 2.67040, 2.93444},
2032        {2.63456, 2.63660, 2.64304, 2.65555, 2.67938, 2.75739,
2033         2.72629, 2.71064, 2.71688, 2.75671, 3.01886}};
2034
2035   // Garfield simulation at UD = -1700V; vd = 1.299cm/microsec, <driftfield> = 590V/cm
2036   Float_t time1700[ktimebin][kZbin] =
2037       {{0.09167, 0.09201, 0.09302, 0.09471, 0.09712, 0.10029,
2038         0.10438, 0.10986, 0.11832, 0.13986, 0.37824},
2039        {0.06591, 0.06626, 0.06736, 0.06923, 0.07183, 0.07519,
2040         0.07944, 0.08511, 0.09378, 0.11603, 0.35625},
2041        {0.03946, 0.04003, 0.04170, 0.04433, 0.04777, 0.05189,
2042         0.05676, 0.06301, 0.07291, 0.09880, 0.34724},
2043        {0.01110, 0.01281, 0.01718, 0.02280, 0.02879, 0.03476,
2044         0.04096, 0.04910, 0.06417, 0.10582, 0.36861},
2045        {0.01115, 0.01294, 0.01721, 0.02276, 0.02862, 0.03447,
2046         0.04074, 0.04980, 0.06812, 0.11565, 0.37231},
2047        {0.03920, 0.03974, 0.04133, 0.04381, 0.04706, 0.05105,
2048         0.05603, 0.06299, 0.07618, 0.10860, 0.36646},
2049        {0.06498, 0.06529, 0.06634, 0.06808, 0.07055, 0.07395,
2050         0.07852, 0.08500, 0.09650, 0.12532, 0.37850},
2051        {0.09143, 0.09159, 0.09207, 0.09297, 0.09473, 0.00000,
2052         0.11102, 0.11809, 0.13245, 0.16627, 0.42496},
2053        {0.13646, 0.13750, 0.14112, 0.14926, 0.16806, 0.24960,
2054         0.20627, 0.19536, 0.20366, 0.24256, 0.50557},
2055        {0.20678, 0.20848, 0.21384, 0.22450, 0.24552, 0.32018,
2056         0.28740, 0.27477, 0.28196, 0.32128, 0.58475},
2057        {0.28287, 0.28461, 0.29020, 0.30108, 0.32224, 0.39467,
2058         0.36500, 0.35217, 0.35926, 0.39860, 0.66194},
2059        {0.35972, 0.36145, 0.36713, 0.37797, 0.39912, 0.47091,
2060         0.44212, 0.42925, 0.43632, 0.47563, 0.73892},
2061        {0.43667, 0.43841, 0.44413, 0.45494, 0.47607, 0.54780,
2062         0.51912, 0.50627, 0.51334, 0.55254, 0.81595},
2063        {0.51365, 0.51540, 0.52101, 0.53193, 0.55305, 0.62463,
2064         0.59617, 0.58328, 0.59035, 0.62965, 0.89303},
2065        {0.59064, 0.59240, 0.59801, 0.60893, 0.63009, 0.70169,
2066         0.67317, 0.66028, 0.66735, 0.70666, 0.96995},
2067        {0.66765, 0.66939, 0.67501, 0.68592, 0.70724, 0.77863,
2068         0.75016, 0.73728, 0.74435, 0.78366, 1.04696},
2069        {0.74464, 0.74636, 0.75200, 0.76293, 0.78405, 0.85561,
2070         0.82716, 0.81427, 0.82133, 0.86064, 1.12396},
2071        {0.82165, 0.82340, 0.82902, 0.83991, 0.86104, 0.93266,
2072         0.90414, 0.89128, 0.89834, 0.93763, 1.20100},
2073        {0.89863, 0.90042, 0.90659, 0.91705, 0.93805, 1.00960,
2074         0.98115, 0.96825, 0.97533, 1.01462, 1.27801},
2075        {0.97563, 0.97740, 0.98310, 0.99391, 1.01504, 1.08659,
2076         1.05814, 1.04526, 1.05233, 1.09163, 1.35503},
2077        {1.05276, 1.05451, 1.06002, 1.07090, 1.09099, 1.16357,
2078         1.13516, 1.12225, 1.12933, 1.16863, 1.43195},
2079        {1.12977, 1.13138, 1.13700, 1.14792, 1.16797, 1.24061,
2080         1.21212, 1.19926, 1.20626, 1.24554, 1.50900},
2081        {1.20664, 1.20839, 1.21400, 1.22490, 1.24606, 1.31772,
2082         1.28914, 1.27382, 1.28329, 1.32262, 1.58550},
2083        {1.28367, 1.28541, 1.29099, 1.30189, 1.32312, 1.39460,
2084         1.36612, 1.34924, 1.36030, 1.39961, 1.66310},
2085        {1.36064, 1.36249, 1.36799, 1.37896, 1.40004, 1.48030,
2086         1.44314, 1.43032, 1.43731, 1.47659, 1.73442},
2087        {1.43762, 1.43937, 1.44497, 1.45618, 1.47704, 1.54932,
2088         1.52012, 1.50725, 1.51430, 1.55357, 1.81708},
2089        {1.51462, 1.51937, 1.52203, 1.53316, 1.55403, 1.62572,
2090         1.59713, 1.58424, 1.59128, 1.63061, 1.89406},
2091        {1.59162, 1.59338, 1.59947, 1.60989, 1.63103, 1.70270,
2092         1.67411, 1.66124, 1.66799, 1.70759, 1.97103},
2093        {1.66874, 1.67037, 1.67597, 1.68687, 1.70814, 1.77969,
2094         1.75112, 1.73806, 1.74530, 1.78457, 2.04794},
2095        {1.74693, 1.74749, 1.75297, 1.76476, 1.78500, 1.85667,
2096         1.82811, 1.81504, 1.82101, 1.86161, 2.12492},
2097        {1.82260, 1.82437, 1.82995, 1.84174, 1.86202, 1.93372,
2098         1.90509, 1.89202, 1.89930, 1.93859, 2.20189},
2099        {1.89964, 1.90135, 1.90693, 1.91789, 1.93952, 2.01080,
2100         1.98207, 1.96921, 1.97628, 2.01384, 2.27887},
2101        {1.97662, 1.97917, 1.98611, 1.99487, 2.01601, 2.08778,
2102         2.05846, 2.04623, 2.05330, 2.09244, 2.35585},
2103        {2.05359, 2.05615, 2.06309, 2.07187, 2.09867, 2.16459,
2104         2.13610, 2.12322, 2.13029, 2.16942, 2.43199},
2105        {2.13063, 2.13233, 2.13795, 2.14886, 2.17008, 2.24199,
2106         2.21310, 2.20020, 2.20727, 2.24659, 2.50983},
2107        {2.20761, 2.20931, 2.21955, 2.22624, 2.24708, 2.32147,
2108         2.29009, 2.27725, 2.28276, 2.32359, 2.58680},
2109        {2.28459, 2.29108, 2.29202, 2.30286, 2.32007, 2.39559,
2110         2.36683, 2.35422, 2.36119, 2.40058, 2.66081},
2111        {2.36153, 2.36806, 2.36889, 2.37985, 2.40092, 2.47828,
2112         2.44381, 2.43099, 2.43819, 2.47750, 2.73779}};
2113
2114   // Garfield simulation at UD = -1800V; vd = 1.450cm/microsec, <driftfield> = 623V/cm
2115   Float_t time1800[ktimebin][kZbin] =
2116       {{0.09166, 0.09199, 0.09300, 0.09470, 0.09709, 0.10026,
2117         0.10434, 0.10983, 0.11831, 0.13987, 0.37802},
2118        {0.06585, 0.06623, 0.06735, 0.06921, 0.07180, 0.07520,
2119         0.07941, 0.08507, 0.09376, 0.11604, 0.35624},
2120        {0.03945, 0.04004, 0.04169, 0.04432, 0.04776, 0.05187,
2121         0.05674, 0.06300, 0.07290, 0.09884, 0.34704},
2122        {0.01108, 0.01287, 0.01717, 0.02280, 0.02880, 0.03476,
2123         0.04095, 0.04909, 0.06419, 0.10589, 0.36843},
2124        {0.01115, 0.01287, 0.01720, 0.02276, 0.02862, 0.03448,
2125         0.04073, 0.04973, 0.06799, 0.11535, 0.37224},
2126        {0.03918, 0.03975, 0.04134, 0.04382, 0.04707, 0.05105,
2127         0.05603, 0.06296, 0.07605, 0.10822, 0.36560},
2128        {0.06498, 0.06532, 0.06635, 0.06809, 0.07058, 0.07395,
2129         0.07855, 0.08495, 0.09632, 0.12488, 0.37730},
2130        {0.09130, 0.09160, 0.09199, 0.09300, 0.09472, 0.00000,
2131         0.11059, 0.11747, 0.13146, 0.16462, 0.42233},
2132        {0.13364, 0.13449, 0.13767, 0.14481, 0.16147, 0.23635,
2133         0.19706, 0.18812, 0.19704, 0.23520, 0.49749},
2134        {0.19698, 0.19844, 0.20311, 0.21236, 0.23082, 0.29920,
2135         0.26936, 0.25927, 0.26732, 0.30601, 0.56871},
2136        {0.26518, 0.26670, 0.27160, 0.28099, 0.29955, 0.36597,
2137         0.33885, 0.32858, 0.33653, 0.37524, 0.63801},
2138        {0.33441, 0.33553, 0.34040, 0.34987, 0.36841, 0.43428,
2139         0.40797, 0.39763, 0.40556, 0.44425, 0.70698},
2140        {0.40296, 0.40447, 0.40934, 0.41881, 0.43737, 0.50306,
2141         0.47695, 0.46662, 0.47455, 0.51329, 0.77600},
2142        {0.47296, 0.47344, 0.47830, 0.48779, 0.50632, 0.57200,
2143         0.54593, 0.53559, 0.54351, 0.58222, 0.84489},
2144        {0.54089, 0.54264, 0.54727, 0.55673, 0.57529, 0.64094,
2145         0.61490, 0.60457, 0.61249, 0.65118, 0.91394},
2146        {0.60987, 0.61138, 0.61624, 0.62573, 0.64428, 0.70989,
2147         0.68397, 0.67354, 0.68147, 0.72015, 0.98291},
2148        {0.67883, 0.68035, 0.68521, 0.69469, 0.71324, 0.77896,
2149         0.75287, 0.74251, 0.75043, 0.78912, 1.04458},
2150        {0.74780, 0.74932, 0.75421, 0.76367, 0.78221, 0.84785,
2151         0.82185, 0.81148, 0.81941, 0.85811, 1.12085},
2152        {0.81690, 0.81830, 0.82316, 0.83263, 0.85120, 0.91683,
2153         0.89077, 0.88045, 0.88837, 0.92707, 1.18976},
2154        {0.88574, 0.88726, 0.89228, 0.90198, 0.92017, 0.98578,
2155         0.95974, 0.94947, 0.95734, 0.99604, 1.25873},
2156        {0.95493, 0.95624, 0.96110, 0.97058, 0.98913, 1.05481,
2157         1.02873, 1.01839, 1.02631, 1.06503, 1.32772},
2158        {1.02392, 1.02524, 1.03008, 1.03955, 1.05810, 1.12378,
2159         1.09757, 1.08605, 1.09530, 1.13399, 1.39669},
2160        {1.09270, 1.09418, 1.09911, 1.10854, 1.12714, 1.19281,
2161         1.16502, 1.15633, 1.16427, 1.20271, 1.46574},
2162        {1.16168, 1.16323, 1.16801, 1.17772, 1.19604, 1.26190,
2163         1.23399, 1.22531, 1.23323, 1.27194, 1.53475},
2164        {1.23061, 1.23214, 1.23698, 1.24669, 1.26503, 1.33073,
2165         1.30461, 1.29428, 1.30220, 1.34091, 1.60372},
2166        {1.29960, 1.30110, 1.30596, 1.31544, 1.33398, 1.39962,
2167         1.37228, 1.36323, 1.37121, 1.40988, 1.67273},
2168        {1.36851, 1.37007, 1.37512, 1.38441, 1.40297, 1.46865,
2169         1.44256, 1.43222, 1.44017, 1.47878, 1.74155},
2170        {1.43752, 1.43904, 1.44773, 1.45338, 1.47220, 1.53759,
2171         1.51136, 1.50119, 1.50914, 1.54775, 1.81050},
2172        {1.50646, 1.50802, 1.51288, 1.52237, 1.54097, 1.60697,
2173         1.58049, 1.57018, 1.57811, 1.61678, 1.87947},
2174        {1.57545, 1.57720, 1.58185, 1.59134, 1.60996, 1.67787,
2175         1.64929, 1.63914, 1.64707, 1.68570, 1.94851},
2176        {1.64442, 1.64617, 1.65081, 1.66035, 1.67893, 1.74684,
2177         1.71826, 1.70745, 1.71604, 1.75310, 2.01748},
2178        {1.71337, 1.71513, 1.71978, 1.72932, 1.74645, 1.81346,
2179         1.78739, 1.77642, 1.78501, 1.82151, 2.08644},
2180        {1.78238, 1.78410, 1.78876, 1.79824, 1.81678, 1.88332,
2181         1.85639, 1.84262, 1.85397, 1.89270, 2.15533},
2182        {1.85135, 1.85306, 1.85778, 1.86728, 1.88580, 1.95615,
2183         1.92536, 1.91171, 1.92283, 1.96165, 2.22428},
2184        {1.92774, 1.92184, 1.92672, 1.93618, 1.95477, 2.02048,
2185         1.99427, 1.98068, 1.99192, 2.03062, 2.29338},
2186        {1.98929, 1.99081, 1.99567, 2.00515, 2.02373, 2.08987,
2187         2.06332, 2.05249, 2.05939, 2.09928, 2.36227},
2188        {2.05829, 2.05978, 2.06464, 2.07414, 2.09272, 2.15850,
2189         2.12928, 2.12194, 2.12987, 2.16825, 2.43083},
2190        {2.12726, 2.12869, 2.13360, 2.14425, 2.16160, 2.22872,
2191         2.20118, 2.19078, 2.19876, 2.23416, 2.50007}};
2192
2193   // Garfield simulation at UD = -1900V; vd = 1.610cm/microsec, <driftfield> = 655V/cm
2194   Float_t time1900[ktimebin][kZbin] =
2195       {{0.09166, 0.09198, 0.09298, 0.09467, 0.09707, 0.10023,
2196         0.10431, 0.10980, 0.11828, 0.13988, 0.37789},
2197        {0.06584, 0.06622, 0.06735, 0.06920, 0.07179, 0.07514,
2198         0.07938, 0.08505, 0.09374, 0.11606, 0.35599},
2199        {0.03945, 0.04002, 0.04169, 0.04432, 0.04775, 0.05185,
2200         0.05672, 0.06298, 0.07290, 0.09888, 0.34695},
2201        {0.01109, 0.01281, 0.01717, 0.02279, 0.02878, 0.03476,
2202         0.04094, 0.04909, 0.06421, 0.10597, 0.36823},
2203        {0.01115, 0.01287, 0.01720, 0.02294, 0.02862, 0.03448,
2204         0.04074, 0.04973, 0.06783, 0.11506, 0.37206},
2205        {0.03940, 0.03975, 0.04135, 0.04386, 0.04708, 0.05106,
2206         0.05604, 0.06293, 0.07592, 0.10787, 0.36484},
2207        {0.06500, 0.06534, 0.06638, 0.06811, 0.07060, 0.07413,
2208         0.07852, 0.08491, 0.09614, 0.12446, 0.37626},
2209        {0.09119, 0.09140, 0.09194, 0.09293, 0.09471, 0.00000,
2210         0.11013, 0.11685, 0.13050, 0.16302, 0.41991},
2211        {0.13111, 0.13190, 0.13466, 0.14091, 0.15554, 0.22409,
2212         0.18846, 0.18167, 0.19113, 0.22854, 0.48995},
2213        {0.18849, 0.18975, 0.19380, 0.20185, 0.21797, 0.28050,
2214         0.25368, 0.24575, 0.25446, 0.29249, 0.55442},
2215        {0.24995, 0.25125, 0.25563, 0.26366, 0.27986, 0.34065,
2216         0.31605, 0.30815, 0.31680, 0.35482, 0.61697},
2217        {0.31187, 0.31324, 0.31745, 0.32580, 0.34205, 0.40217,
2218         0.37825, 0.37031, 0.37897, 0.41696, 0.67890},
2219        {0.37401, 0.37531, 0.37955, 0.38777, 0.40395, 0.46408,
2220         0.44037, 0.43242, 0.44108, 0.47906, 0.74122},
2221        {0.43610, 0.43741, 0.44161, 0.44986, 0.46604, 0.52614,
2222         0.50248, 0.49452, 0.50316, 0.54116, 0.80326},
2223        {0.49820, 0.49988, 0.50372, 0.51196, 0.52814, 0.58822,
2224         0.56459, 0.55661, 0.56527, 0.60326, 0.86526},
2225        {0.56032, 0.56161, 0.56582, 0.57408, 0.59024, 0.65032,
2226         0.62670, 0.61872, 0.62737, 0.66537, 0.92749},
2227        {0.62240, 0.62371, 0.62792, 0.63624, 0.65236, 0.71241,
2228         0.68881, 0.68081, 0.68947, 0.72750, 0.98941},
2229        {0.68449, 0.68581, 0.69002, 0.69828, 0.71444, 0.77452,
2230         0.75089, 0.74295, 0.75157, 0.78957, 1.05157},
2231        {0.74660, 0.74790, 0.75212, 0.76036, 0.77654, 0.83748,
2232         0.81299, 0.80501, 0.81193, 0.85168, 1.11375},
2233        {0.80870, 0.81017, 0.81423, 0.82246, 0.83867, 0.89908,
2234         0.87509, 0.86660, 0.87577, 0.91376, 1.17586},
2235        {0.87080, 0.87233, 0.87632, 0.88458, 0.90074, 0.96083,
2236         0.93718, 0.92922, 0.93787, 0.97588, 1.23794},
2237        {0.93291, 0.93422, 0.93844, 0.94667, 0.96293, 1.02295,
2238         0.99929, 0.99127, 0.99997, 1.03797, 1.29995},
2239        {0.99500, 0.99645, 1.00308, 1.00877, 1.02497, 1.08504,
2240         1.06140, 1.05343, 1.06203, 1.10006, 1.36216},
2241        {1.05712, 1.05926, 1.06262, 1.07092, 1.08706, 1.14754,
2242         1.12350, 1.11550, 1.12417, 1.16218, 1.42427},
2243        {1.11921, 1.12059, 1.12473, 1.13297, 1.14916, 1.21140,
2244         1.18560, 1.17284, 1.18625, 1.22414, 1.48629},
2245        {1.18140, 1.18262, 1.18690, 1.19508, 1.21125, 1.27139,
2246         1.24164, 1.23495, 1.24838, 1.28634, 1.54852},
2247        {1.24340, 1.24473, 1.24901, 1.25732, 1.27336, 1.33358,
2248         1.30793, 1.30179, 1.31047, 1.34848, 1.61066},
2249        {1.30551, 1.30684, 1.31104, 1.32056, 1.33553, 1.39609,
2250         1.37004, 1.36392, 1.37045, 1.41057, 1.67259},
2251        {1.36755, 1.36892, 1.37315, 1.39148, 1.39755, 1.45820,
2252         1.43215, 1.42602, 1.43467, 1.47268, 1.73477},
2253        {1.42966, 1.43101, 1.43549, 1.45359, 1.45976, 1.52031,
2254         1.49601, 1.48811, 1.49677, 1.53477, 1.79691},
2255        {1.49180, 1.49321, 1.49760, 1.51570, 1.52175, 1.58185,
2256         1.55771, 1.55023, 1.55888, 1.59672, 1.85501},
2257        {1.55391, 1.55527, 1.55943, 1.57782, 1.58391, 1.64395,
2258         1.62008, 1.61233, 1.62085, 1.65883, 1.92091},
2259        {1.61599, 1.61732, 1.62154, 1.63993, 1.64612, 1.70608,
2260         1.68237, 1.67108, 1.68301, 1.72110, 1.97931},
2261        {1.67808, 1.67948, 1.68364, 1.70204, 1.70823, 1.76858,
2262         1.74404, 1.73539, 1.74512, 1.78321, 2.04522},
2263        {1.74019, 1.74152, 1.74573, 1.76415, 1.77015, 1.83040,
2264         1.80615, 1.79366, 1.80723, 1.84509, 2.10742},
2265        {1.80235, 1.80362, 1.80783, 1.82626, 1.83227, 1.89246,
2266         1.86795, 1.85405, 1.86938, 1.90720, 2.16953},
2267        {1.86442, 1.86572, 1.86994, 1.88837, 1.89438, 1.95445,
2268         1.93006, 1.92283, 1.93148, 1.96931, 2.23147},
2269        {1.92700, 1.92783, 1.93197, 1.95049, 1.95649, 2.01668,
2270         1.99217, 1.98486, 1.99352, 2.03143, 2.29358}};
2271
2272   // Garfield simulation at UD = -2000V; vd = 1.783cm/microsec, <driftfield> = 688V/cm
2273   Float_t time2000[ktimebin][kZbin] =
2274       {{0.09176, 0.09196, 0.09296, 0.09465, 0.09704, 0.10020,
2275         0.10427, 0.10977, 0.11825, 0.13988, 0.37774},
2276        {0.06583, 0.06620, 0.06735, 0.06918, 0.07177, 0.07513,
2277         0.07936, 0.08503, 0.09372, 0.11606, 0.35586},
2278        {0.03944, 0.04001, 0.04170, 0.04431, 0.04774, 0.05184,
2279         0.05670, 0.06296, 0.07291, 0.09893, 0.34680},
2280        {0.01108, 0.01281, 0.01719, 0.02279, 0.02879, 0.03474,
2281         0.04093, 0.04908, 0.06422, 0.10605, 0.36800},
2282        {0.01114, 0.01287, 0.01720, 0.02276, 0.02863, 0.03449,
2283         0.04073, 0.04970, 0.06774, 0.11478, 0.37179},
2284        {0.03925, 0.03977, 0.04135, 0.04386, 0.04711, 0.05108,
2285         0.05604, 0.06290, 0.07580, 0.10748, 0.36386},
2286        {0.06501, 0.06536, 0.06640, 0.06814, 0.07062, 0.07398,
2287         0.07852, 0.08487, 0.09598, 0.12405, 0.37519},
2288        {0.09109, 0.09127, 0.09188, 0.09292, 0.09472, 0.00000,
2289         0.10964, 0.11630, 0.12960, 0.16150, 0.41765},
2290        {0.12898, 0.12968, 0.13209, 0.13749, 0.15034, 0.21286,
2291         0.18088, 0.17590, 0.18591, 0.22254, 0.48315},
2292        {0.18122, 0.18227, 0.18574, 0.19263, 0.20674, 0.26376,
2293         0.23960, 0.23375, 0.24316, 0.28047, 0.54179},
2294        {0.23674, 0.23784, 0.24142, 0.24847, 0.26264, 0.31810,
2295         0.29602, 0.29008, 0.29944, 0.33675, 0.59795},
2296        {0.29279, 0.29382, 0.29742, 0.30448, 0.31865, 0.37364,
2297         0.35215, 0.34629, 0.35555, 0.39286, 0.65411},
2298        {0.34875, 0.34987, 0.35346, 0.36054, 0.37472, 0.42956,
2299         0.40825, 0.40229, 0.41167, 0.44894, 0.71033},
2300        {0.40484, 0.40594, 0.40954, 0.41660, 0.43077, 0.48560,
2301         0.46433, 0.45840, 0.46772, 0.50500, 0.76632},
2302        {0.46090, 0.46201, 0.46560, 0.47267, 0.48684, 0.54167,
2303         0.52041, 0.51449, 0.52382, 0.56108, 0.82227},
2304        {0.51698, 0.51809, 0.52173, 0.52874, 0.54291, 0.59776,
2305         0.57646, 0.57052, 0.57986, 0.61717, 0.87836},
2306        {0.57306, 0.57418, 0.57782, 0.58513, 0.59899, 0.65380,
2307         0.63255, 0.62661, 0.63594, 0.67325, 0.93460},
2308        {0.62912, 0.63024, 0.63383, 0.64103, 0.65506, 0.70988,
2309         0.68484, 0.68267, 0.69202, 0.72878, 0.99046},
2310        {0.68521, 0.68633, 0.68990, 0.69699, 0.71115, 0.76595,
2311         0.74468, 0.73872, 0.74814, 0.78538, 1.04674},
2312        {0.74127, 0.74239, 0.74605, 0.75303, 0.77022, 0.82204,
2313         0.80078, 0.79484, 0.80416, 0.84147, 1.10261},
2314        {0.79736, 0.79846, 0.80206, 0.80947, 0.82330, 0.87813,
2315         0.85688, 0.85087, 0.86023, 0.89752, 1.15874},
2316        {0.85342, 0.85454, 0.85815, 0.86519, 0.87936, 0.93417,
2317         0.91293, 0.90428, 0.91631, 0.95360, 1.20760},
2318        {0.90949, 0.91061, 0.91423, 0.92128, 0.93544, 0.99026,
2319         0.96807, 0.96305, 0.97239, 1.00967, 1.27078},
2320        {0.96556, 0.96669, 0.97111, 0.97734, 0.99151, 1.04664,
2321         1.02508, 1.01879, 1.02846, 1.06167, 1.32695},
2322        {1.02167, 1.02279, 1.02656, 1.03341, 1.04759, 1.10242,
2323         1.08115, 1.07003, 1.08453, 1.12184, 1.38304},
2324        {1.07776, 1.07883, 1.08242, 1.08950, 1.10384, 1.16422,
2325         1.13725, 1.13133, 1.14061, 1.17793, 1.43910},
2326        {1.13379, 1.13492, 1.13864, 1.14567, 1.15973, 1.21455,
2327         1.19323, 1.18734, 1.19668, 1.23401, 1.49528},
2328        {1.18988, 1.19098, 1.19457, 1.20164, 1.21582, 1.27064,
2329         1.24937, 1.24044, 1.25275, 1.29004, 1.55137},
2330        {1.24592, 1.24706, 1.25087, 1.25773, 1.27188, 1.32670,
2331         1.30544, 1.29953, 1.30883, 1.34613, 1.60743},
2332        {1.30202, 1.30313, 1.30673, 1.31381, 1.32797, 1.38278,
2333         1.36151, 1.35167, 1.36490, 1.40221, 1.66306},
2334        {1.35809, 1.35921, 1.36282, 1.36986, 1.38403, 1.43888,
2335         1.41760, 1.41174, 1.42083, 1.45830, 1.71915},
2336        {1.41419, 1.41528, 1.41890, 1.42595, 1.44011, 1.49496,
2337         1.47368, 1.46769, 1.47706, 1.51436, 1.77523},
2338        {1.47131, 1.47141, 1.47494, 1.48850, 1.49620, 1.55137,
2339         1.52977, 1.51820, 1.53315, 1.57042, 1.83158},
2340        {1.52635, 1.52750, 1.53103, 1.53814, 1.55228, 1.60736,
2341         1.58503, 1.57986, 1.58920, 1.62649, 1.88767},
2342        {1.58418, 1.58355, 1.58711, 1.59526, 1.60833, 1.66316,
2343         1.63345, 1.63261, 1.64556, 1.68204, 1.94359},
2344        {1.64027, 1.63958, 1.64489, 1.65024, 1.66443, 1.71925,
2345         1.69794, 1.69201, 1.70143, 1.73865, 1.99968},
2346        {1.69450, 1.69566, 1.69940, 1.70697, 1.71841, 1.77819,
2347         1.75396, 1.74814, 1.75743, 1.79083, 2.05427},
2348        {1.75054, 1.75221, 1.75527, 1.76306, 1.77662, 1.83428,
2349         1.81006, 1.81173, 1.81345, 1.85076, 2.10289}};
2350
2351   // Garfield simulation at UD = -2100V; vd = 1.959cm/microsec, <driftfield> = 720V/cm
2352   Float_t time2100[ktimebin][kZbin] =
2353       {{0.09160, 0.09194, 0.09294, 0.09462, 0.09701, 0.10017,
2354         0.10424, 0.10974, 0.11823, 0.13988, 0.37762},
2355        {0.06585, 0.06619, 0.06731, 0.06916, 0.07174, 0.07509,
2356         0.07933, 0.08500, 0.09370, 0.11609, 0.35565},
2357        {0.03960, 0.04001, 0.04171, 0.04430, 0.04774, 0.05182,
2358         0.05668, 0.06294, 0.07291, 0.09896, 0.34676},
2359        {0.01109, 0.01280, 0.01716, 0.02279, 0.02876, 0.03474,
2360         0.04096, 0.04908, 0.06424, 0.10612, 0.36790},
2361        {0.01114, 0.01285, 0.01719, 0.02287, 0.02863, 0.03449,
2362         0.04073, 0.04964, 0.06759, 0.11446, 0.37162},
2363        {0.03922, 0.03977, 0.04146, 0.04386, 0.04711, 0.05109,
2364         0.05605, 0.06287, 0.07575, 0.10713, 0.36298},
2365        {0.06504, 0.06538, 0.06641, 0.06818, 0.07064, 0.07426,
2366         0.07852, 0.08483, 0.09581, 0.12363, 0.37424},
2367        {0.09103, 0.09129, 0.09186, 0.09291, 0.09476, 0.00000,
2368         0.10923, 0.11578, 0.12873, 0.16005, 0.41525},
2369        {0.12723, 0.12777, 0.12988, 0.13458, 0.14579, 0.20264,
2370         0.17421, 0.17078, 0.18132, 0.21708, 0.47699},
2371        {0.17508, 0.17601, 0.17897, 0.18487, 0.19698, 0.24881,
2372         0.22737, 0.22337, 0.23348, 0.27000, 0.53032},
2373        {0.22571, 0.22663, 0.22969, 0.23570, 0.24787, 0.29826,
2374         0.27871, 0.27462, 0.28471, 0.32122, 0.58166},
2375        {0.27664, 0.27759, 0.28067, 0.28669, 0.29891, 0.34898,
2376         0.32982, 0.32570, 0.33576, 0.37229, 0.63268},
2377        {0.32766, 0.32862, 0.33170, 0.33778, 0.34988, 0.39973,
2378         0.38088, 0.37675, 0.38680, 0.42333, 0.68159},
2379        {0.37872, 0.37966, 0.38275, 0.38875, 0.40093, 0.45073,
2380         0.43192, 0.42780, 0.43786, 0.47438, 0.73480},
2381        {0.42974, 0.43070, 0.43378, 0.43982, 0.45196, 0.50177,
2382         0.48297, 0.47884, 0.48889, 0.52544, 0.78581},
2383        {0.48081, 0.48175, 0.48482, 0.49084, 0.50302, 0.55290,
2384         0.53398, 0.52988, 0.53994, 0.57647, 0.83687},
2385        {0.53645, 0.53295, 0.53586, 0.54188, 0.55408, 0.60398,
2386         0.58504, 0.58092, 0.59100, 0.62768, 0.88773},
2387        {0.58345, 0.58409, 0.58690, 0.59292, 0.60510, 0.65562,
2388         0.63609, 0.63197, 0.64203, 0.67856, 0.93907},
2389        {0.63397, 0.63490, 0.63795, 0.64403, 0.65613, 0.70612,
2390         0.68714, 0.68301, 0.69294, 0.72955, 0.99000},
2391        {0.68496, 0.68592, 0.68899, 0.69504, 0.70733, 0.75708,
2392         0.73818, 0.73405, 0.74412, 0.78064, 1.04100},
2393        {0.73600, 0.73696, 0.74003, 0.74624, 0.75828, 0.80805,
2394         0.78904, 0.78512, 0.79517, 0.83152, 1.09205},
2395        {0.78709, 0.78801, 0.79108, 0.79709, 0.80931, 0.85906,
2396         0.84027, 0.83614, 0.84621, 0.88269, 1.14058},
2397        {0.83808, 0.83905, 0.84215, 0.84816, 0.86031, 0.91011,
2398         0.89139, 0.88718, 0.89725, 0.93377, 1.19413},
2399        {0.88916, 0.89010, 0.89320, 0.89920, 0.91136, 0.96117,
2400         0.94235, 0.93822, 0.94828, 0.98480, 1.24538},
2401        {0.94036, 0.94113, 0.94422, 0.95023, 0.96241, 1.01220,
2402         0.99310, 0.98927, 0.99933, 1.03585, 1.29629},
2403        {0.99139, 0.99220, 0.99525, 1.00127, 1.01344, 1.06324,
2404         1.04451, 1.04033, 1.04836, 1.08690, 1.34727},
2405        {1.04261, 1.04325, 1.04628, 1.05232, 1.06448, 1.12090,
2406         1.09546, 1.09136, 1.10142, 1.13795, 1.39831},
2407        {1.09331, 1.09429, 1.09742, 1.10336, 1.11557, 1.16547,
2408         1.14658, 1.13642, 1.15247, 1.18898, 1.44936},
2409        {1.14436, 1.14539, 1.14847, 1.15443, 1.16662, 1.21794,
2410         1.19763, 1.19329, 1.20349, 1.23956, 1.50043},
2411        {1.19533, 1.19651, 1.19943, 1.20548, 1.21666, 1.26753,
2412         1.24862, 1.24434, 1.25455, 1.29106, 1.55142},
2413        {1.24638, 1.24756, 1.25046, 1.25648, 1.26764, 1.31858,
2414         1.29967, 1.29538, 1.30499, 1.34211, 1.60250},
2415        {1.29747, 1.29847, 1.30175, 1.30753, 1.31869, 1.36969,
2416         1.35069, 1.34656, 1.35663, 1.39316, 1.64644},
2417        {1.35537, 1.34952, 1.35255, 1.35869, 1.36973, 1.41387,
2418         1.40173, 1.39761, 1.40768, 1.44396, 1.70238},
2419        {1.39956, 1.40056, 1.40380, 1.40961, 1.42178, 1.46492,
2420         1.45278, 1.45423, 1.45872, 1.49522, 1.75557},
2421        {1.45080, 1.45159, 1.45463, 1.46109, 1.47287, 1.52263,
2422         1.50382, 1.50050, 1.50977, 1.54502, 1.80670},
2423        {1.50165, 1.50264, 1.50570, 1.51214, 1.52233, 1.57370,
2424         1.55484, 1.55155, 1.56080, 1.59731, 1.85778},
2425        {1.55269, 1.55364, 1.55675, 1.56274, 1.57491, 1.62598,
2426         1.60590, 1.60259, 1.61185, 1.64836, 1.90883},
2427        {1.60368, 1.60469, 1.60779, 1.61373, 1.62596, 1.67738,
2428         1.65651, 1.65249, 1.66290, 1.69936, 1.95959}};
2429
2430   // Garfield simulation at UD = -2200V; vd = 2.134cm/microsec, <driftfield> = 753V/cm
2431   Float_t time2200[ktimebin][kZbin] =
2432       {{0.09162, 0.09194, 0.09292, 0.09460, 0.09702, 0.10014,
2433         0.10421, 0.10971, 0.11820, 0.13990, 0.37745},
2434        {0.06581, 0.06618, 0.06730, 0.06915, 0.07173, 0.07507,
2435         0.07931, 0.08497, 0.09368, 0.11609, 0.35560},
2436        {0.03947, 0.04001, 0.04167, 0.04429, 0.04772, 0.05183,
2437         0.05667, 0.06293, 0.07292, 0.09900, 0.34673},
2438        {0.01111, 0.01280, 0.01716, 0.02279, 0.02876, 0.03473,
2439         0.04091, 0.04907, 0.06426, 0.10620, 0.36766},
2440        {0.01113, 0.01285, 0.01719, 0.02276, 0.02863, 0.03452,
2441         0.04076, 0.04960, 0.06745, 0.11419, 0.37139},
2442        {0.03923, 0.03978, 0.04137, 0.04387, 0.04713, 0.05110,
2443         0.05605, 0.06284, 0.07551, 0.10677, 0.36210},
2444        {0.06505, 0.06540, 0.06644, 0.06820, 0.07069, 0.07401,
2445         0.07852, 0.08479, 0.09565, 0.12325, 0.37313},
2446        {0.09107, 0.09127, 0.09181, 0.09291, 0.09474, 0.00000,
2447         0.10883, 0.11528, 0.12789, 0.15865, 0.41313},
2448        {0.12559, 0.12622, 0.12800, 0.13206, 0.14166, 0.19331,
2449         0.16832, 0.16632, 0.17724, 0.21218, 0.47098},
2450        {0.16992, 0.17070, 0.17325, 0.17831, 0.18871, 0.23557,
2451         0.21690, 0.21451, 0.22514, 0.26082, 0.52034},
2452        {0.21640, 0.21722, 0.21987, 0.22500, 0.23540, 0.28097,
2453         0.26399, 0.26154, 0.27214, 0.30784, 0.56734},
2454        {0.26318, 0.26400, 0.26679, 0.27181, 0.28220, 0.32739,
2455         0.31090, 0.30842, 0.31902, 0.35474, 0.61415},
2456        {0.31001, 0.31085, 0.31348, 0.31866, 0.32903, 0.37412,
2457         0.35777, 0.35546, 0.36588, 0.40159, 0.66103},
2458        {0.35687, 0.35769, 0.36033, 0.36556, 0.37588, 0.42094,
2459         0.40523, 0.40214, 0.41273, 0.44841, 0.70785},
2460        {0.40372, 0.40457, 0.40723, 0.41234, 0.42273, 0.46778,
2461         0.45148, 0.44903, 0.45961, 0.49526, 0.75486},
2462        {0.45062, 0.45139, 0.45404, 0.45966, 0.46958, 0.51470,
2463         0.49833, 0.49584, 0.50644, 0.54211, 0.80160},
2464        {0.49742, 0.49825, 0.50088, 0.50605, 0.51644, 0.56148,
2465         0.54518, 0.54270, 0.55330, 0.58897, 0.84854},
2466        {0.54427, 0.54510, 0.54774, 0.55290, 0.56329, 0.60846,
2467         0.59203, 0.58955, 0.60014, 0.63578, 0.89528},
2468        {0.59119, 0.59199, 0.59471, 0.59975, 0.61014, 0.65533,
2469         0.63889, 0.63636, 0.64699, 0.68269, 0.94197},
2470        {0.63866, 0.63880, 0.64145, 0.64664, 0.65701, 0.70639,
2471         0.68574, 0.68325, 0.69385, 0.72949, 0.98900},
2472        {0.68483, 0.68566, 0.68831, 0.69347, 0.70386, 0.74890,
2473         0.73260, 0.73010, 0.74069, 0.77638, 1.03320},
2474        {0.73168, 0.73255, 0.73515, 0.74031, 0.75072, 0.79576,
2475         0.77117, 0.77501, 0.78755, 0.82318, 1.08006},
2476        {0.77854, 0.78310, 0.78200, 0.79525, 0.79756, 0.84281,
2477         0.81803, 0.82393, 0.83441, 0.87008, 1.12692},
2478        {0.82541, 0.82642, 0.82916, 0.83528, 0.84442, 0.89086,
2479         0.87315, 0.87079, 0.88125, 0.91694, 1.17648},
2480        {0.87226, 0.87308, 0.87602, 0.88086, 0.89128, 0.93772,
2481         0.92001, 0.91751, 0.92811, 0.95587, 1.22328},
2482        {0.91921, 0.91994, 0.92256, 0.92772, 0.94713, 0.98566,
2483         0.96690, 0.96436, 0.97495, 1.01064, 1.26882},
2484        {0.96790, 0.96679, 0.96941, 0.97463, 0.99399, 1.03001,
2485         1.01376, 1.01112, 1.02181, 1.05749, 1.31568},
2486        {1.01278, 1.01390, 1.01674, 1.02147, 1.03182, 1.07820,
2487         1.06056, 1.05798, 1.06867, 1.10433, 1.36390},
2488        {1.05964, 1.06076, 1.06331, 1.06833, 1.07870, 1.13297,
2489         1.10742, 1.10520, 1.11543, 1.15120, 1.41069},
2490        {1.10664, 1.10762, 1.10997, 1.11519, 1.12556, 1.17531,
2491         1.15427, 1.14620, 1.16229, 1.19805, 1.45758},
2492        {1.15352, 1.15421, 1.15683, 1.16218, 1.17242, 1.21910,
2493         1.20035, 1.19863, 1.20579, 1.24473, 1.50412},
2494        {1.20019, 1.20115, 1.20369, 1.20892, 1.21928, 1.26913,
2495         1.24721, 1.24549, 1.25605, 1.29159, 1.54920},
2496        {1.24707, 1.24846, 1.25052, 1.25602, 1.26608, 1.31558,
2497         1.29448, 1.29232, 1.30293, 1.33675, 1.59798},
2498        {1.29391, 1.29475, 1.29738, 1.30255, 1.31294, 1.36244,
2499         1.34167, 1.33918, 1.34979, 1.38229, 1.64496},
2500        {1.34078, 1.34304, 1.34424, 1.35565, 1.35980, 1.40930,
2501         1.38853, 1.38229, 1.39664, 1.42863, 1.69162},
2502        {1.38762, 1.38847, 1.39110, 1.39627, 1.40666, 1.45183,
2503         1.43539, 1.43289, 1.44348, 1.47549, 1.73876},
2504        {1.43524, 1.43533, 1.43796, 1.44310, 1.45371, 1.49305,
2505         1.48224, 1.47941, 1.49034, 1.52601, 1.78552},
2506        {1.48122, 1.48219, 1.48482, 1.48991, 1.50030, 1.53991,
2507         1.52898, 1.52653, 1.53653, 1.57282, 1.82386}};
2508
2509   if (fTimeStruct1) {
2510     delete [] fTimeStruct1;
2511   }
2512   fTimeStruct1 = new Float_t[ktimebin*kZbin];
2513
2514   if (fTimeStruct2) {
2515     delete [] fTimeStruct2;
2516   }
2517   fTimeStruct2 = new Float_t[ktimebin*kZbin];
2518
2519   for (Int_t ctrt = 0; ctrt < ktimebin; ctrt++) {
2520     for (Int_t ctrz = 0; ctrz <    kZbin; ctrz++) {
2521       if      (vdrift > fVDsmp[6]) {
2522         fTimeStruct1[ctrt+ctrz*ktimebin] = time2100[ctrt][ctrz];
2523         fTimeStruct2[ctrt+ctrz*ktimebin] = time2200[ctrt][ctrz];
2524         fVDlo = fVDsmp[6];
2525         fVDhi = fVDsmp[7];
2526       }
2527       else if (vdrift > fVDsmp[5]) {
2528         fTimeStruct1[ctrt+ctrz*ktimebin] = time2000[ctrt][ctrz];
2529         fTimeStruct2[ctrt+ctrz*ktimebin] = time2100[ctrt][ctrz];
2530         fVDlo = fVDsmp[5];
2531         fVDhi = fVDsmp[6];
2532       }
2533       else if (vdrift > fVDsmp[4]) {
2534         fTimeStruct1[ctrt+ctrz*ktimebin] = time1900[ctrt][ctrz];
2535         fTimeStruct2[ctrt+ctrz*ktimebin] = time2000[ctrt][ctrz];
2536         fVDlo = fVDsmp[4];
2537         fVDhi = fVDsmp[5];
2538       }
2539       else if (vdrift > fVDsmp[3]) {
2540         fTimeStruct1[ctrt+ctrz*ktimebin] = time1800[ctrt][ctrz];
2541         fTimeStruct2[ctrt+ctrz*ktimebin] = time1900[ctrt][ctrz];
2542         fVDlo = fVDsmp[3];
2543         fVDhi = fVDsmp[4];
2544       }
2545       else if (vdrift > fVDsmp[2]) {
2546         fTimeStruct1[ctrt+ctrz*ktimebin] = time1700[ctrt][ctrz];
2547         fTimeStruct2[ctrt+ctrz*ktimebin] = time1800[ctrt][ctrz];
2548         fVDlo = fVDsmp[2];
2549         fVDhi = fVDsmp[3];
2550       }
2551       else if (vdrift > fVDsmp[1]) {
2552         fTimeStruct1[ctrt+ctrz*ktimebin] = time1600[ctrt][ctrz];
2553         fTimeStruct2[ctrt+ctrz*ktimebin] = time1700[ctrt][ctrz];
2554         fVDlo = fVDsmp[1];
2555         fVDhi = fVDsmp[2];
2556       }
2557       else if (vdrift > (fVDsmp[0] - 1.0e-5)) {
2558         fTimeStruct1[ctrt+ctrz*ktimebin] = time1500[ctrt][ctrz];
2559         fTimeStruct2[ctrt+ctrz*ktimebin] = time1600[ctrt][ctrz];
2560         fVDlo = fVDsmp[0];
2561         fVDhi = fVDsmp[1];
2562       }
2563     }
2564   }
2565
2566 }
2567
2568 //_____________________________________________________________________________
2569 void AliTRDdigitizer::RecalcDiffusion(Float_t vdrift)
2570 {
2571   //
2572   // Recalculates the diffusion parameters
2573   //
2574   // The B=0 case is not really included here.
2575   // Should be revisited!
2576   //
2577
2578   if (vdrift == fDiffLastVdrift) {
2579     return;
2580   }
2581
2582   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
2583   if (!simParam) {
2584     AliFatal("Could not get simulation parameters");
2585     return;
2586   }
2587
2588   AliTRDCommonParam *commonParam = AliTRDCommonParam::Instance();
2589   if (!commonParam) {
2590     AliFatal("Could not get common parameters");
2591     return;
2592   }
2593
2594   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
2595   if (!calibration) {
2596     AliFatal("Could not get calibration object");
2597     return;
2598   }
2599
2600   // The magnetic field strength
2601   Double_t x[3] = { 0.0, 0.0, 0.0 };
2602   Double_t b[3];
2603   gAlice->Field(x,b);         // b[] is in kilo Gauss
2604   Float_t field = b[2] * 0.1; // Tesla
2605
2606   fDiffLastVdrift = vdrift;
2607
2608   // DiffusionL
2609   const Int_t kNbL = 5;
2610   Float_t p0L[kNbL] = {  0.007440,  0.007493,  0.007513,  0.007672,  0.007831 };
2611   Float_t p1L[kNbL] = {  0.019252,  0.018912,  0.018636,  0.018012,  0.017343 };
2612   Float_t p2L[kNbL] = { -0.005042, -0.004926, -0.004867, -0.004650, -0.004424 };
2613   Float_t p3L[kNbL] = {  0.000195,  0.000189,  0.000195,  0.000182,  0.000169 };
2614
2615   Int_t ibL = ((Int_t) (10 * (field - 0.15)));
2616   ibL       = TMath::Max(   0,ibL);
2617   ibL       = TMath::Min(kNbL,ibL);
2618
2619   fDiffusionL = p0L[ibL]
2620               + p1L[ibL] * vdrift
2621               + p2L[ibL] * vdrift*vdrift
2622               + p3L[ibL] * vdrift*vdrift*vdrift;
2623
2624   // DiffusionT
2625   const Int_t kNbT = 5;
2626   Float_t p0T[kNbT] = {  0.009550,  0.009599,  0.009674,  0.009757,  0.009850 };
2627   Float_t p1T[kNbT] = {  0.006667,  0.006539,  0.006359,  0.006153,  0.005925 };
2628   Float_t p2T[kNbT] = { -0.000853, -0.000798, -0.000721, -0.000635, -0.000541 };
2629   Float_t p3T[kNbT] = {  0.000131,  0.000122,  0.000111,  0.000098,  0.000085 };
2630
2631   Int_t ibT= ((Int_t) (10 * (field - 0.15)));
2632   ibT      = TMath::Max(   0,ibT);
2633   ibT      = TMath::Min(kNbT,ibT);
2634
2635   fDiffusionT = p0T[ibT]
2636               + p1T[ibT] * vdrift
2637               + p2T[ibT] * vdrift*vdrift
2638               + p3T[ibT] * vdrift*vdrift*vdrift;
2639
2640   // OmegaTau
2641   fOmegaTau = calibration->GetOmegaTau(vdrift,field);
2642
2643   // Lorentzfactor
2644   if (commonParam->ExBOn()) {
2645     fLorentzFactor = 1.0 / (1.0 + fOmegaTau*fOmegaTau);
2646   }
2647   else {
2648     fLorentzFactor = 1.0;
2649   }
2650
2651 }
2652
2653 //_____________________________________________________________________________
2654 Float_t AliTRDdigitizer::GetDiffusionL(Float_t vdrift)
2655 {
2656   //
2657   // Returns the longitudinal diffusion coefficient for a given drift
2658   // velocity <vd> and a B-field <b> for Xe/CO2 (15%).
2659   // The values are according to a GARFIELD simulation.
2660   //
2661
2662   RecalcDiffusion(vdrift);
2663
2664   return fDiffusionL;
2665
2666 }
2667
2668 //_____________________________________________________________________________
2669 Float_t AliTRDdigitizer::GetDiffusionT(Float_t vdrift)
2670 {
2671   //
2672   // Returns the transverse diffusion coefficient for a given drift
2673   // velocity <vd> and a B-field <b> for Xe/CO2 (15%).
2674   // The values are according to a GARFIELD simulation.
2675   //
2676
2677   RecalcDiffusion(vdrift);
2678
2679   return fDiffusionT;
2680
2681 }
2682
2683 //_____________________________________________________________________________
2684 Int_t AliTRDdigitizer::Diffusion(Float_t vdrift, Double_t absdriftlength
2685                                , Double_t &lRow, Double_t &lCol, Double_t &lTime)
2686 {
2687   //
2688   // Applies the diffusion smearing to the position of a single electron.
2689   // Depends on absolute drift length.
2690   //
2691
2692   RecalcDiffusion(vdrift);
2693
2694   Float_t driftSqrt = TMath::Sqrt(absdriftlength);
2695   Float_t sigmaT    = driftSqrt * fDiffusionT;
2696   Float_t sigmaL    = driftSqrt * fDiffusionL;
2697   lRow  = gRandom->Gaus(lRow ,sigmaT);
2698   lCol  = gRandom->Gaus(lCol ,sigmaT * GetLorentzFactor(vdrift));
2699   lTime = gRandom->Gaus(lTime,sigmaL * GetLorentzFactor(vdrift));
2700
2701   return 1;
2702
2703 }
2704
2705 //_____________________________________________________________________________
2706 Float_t AliTRDdigitizer::GetLorentzFactor(Float_t vd)
2707 {
2708   //
2709   // Returns the recalculated Lorentz factor
2710   //
2711
2712   RecalcDiffusion(vd);
2713
2714   return fLorentzFactor;
2715
2716 }
2717
2718 //_____________________________________________________________________________
2719 Int_t AliTRDdigitizer::ExB(Float_t vdrift, Double_t driftlength, Double_t &lCol)
2720 {
2721   //
2722   // Applies E x B effects to the position of a single electron.
2723   // Depends on signed drift length.
2724   //
2725
2726   RecalcDiffusion(vdrift);
2727
2728   lCol = lCol + fOmegaTau * driftlength;
2729
2730   return 1;
2731
2732 }