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