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       digitsB->Expand();
1549       for (iDict = 0; iDict < kNDict; iDict++) {
1550         dictionaryA[iDict] = (AliTRDdataArrayI *) fSDigitsManager->GetDictionary(iDet,iDict);
1551         dictionaryB[iDict] = (AliTRDdataArrayI *) mergeSDigitsManager->GetDictionary(iDet,iDict);
1552         dictionaryA[iDict]->Expand();
1553         dictionaryB[iDict]->Expand();
1554       }
1555
1556       // Merge only detectors that contain a signal
1557       Bool_t doMerge = kTRUE;
1558       if (fMergeSignalOnly) {
1559         if (digitsA->GetOverThreshold(0) == 0) {
1560           doMerge = kFALSE;
1561         }
1562       }
1563
1564       if (doMerge) {
1565
1566         AliDebug(1,Form("Merge detector %d of input no.%d",iDet,iMerge+1));
1567
1568         for (Int_t iRow  = 0; iRow  <  nRowMax;   iRow++ ) {
1569           for (Int_t iCol  = 0; iCol  <  nColMax;   iCol++ ) {
1570             for (Int_t iTime = 0; iTime < nTimeTotal; iTime++) {
1571
1572               // Add the amplitudes of the summable digits
1573               Float_t ampA = digitsA->GetDataUnchecked(iRow,iCol,iTime);
1574               Float_t ampB = digitsB->GetDataUnchecked(iRow,iCol,iTime);
1575               ampA += ampB;
1576               digitsA->SetDataUnchecked(iRow,iCol,iTime,ampA);
1577
1578               // Add the mask to the track id if defined.
1579               for (iDict = 0; iDict < kNDict; iDict++) {
1580                 Int_t trackB = dictionaryB[iDict]->GetDataUnchecked(iRow,iCol,iTime);
1581                 if ((fMasks) && (trackB > 0)) {
1582                   for (jDict = 0; jDict < kNDict; jDict++) {
1583                     Int_t trackA = dictionaryA[iDict]->GetDataUnchecked(iRow,iCol,iTime);
1584                     if (trackA == 0) {
1585                       trackA = trackB + fMasks[iMerge];
1586                       dictionaryA[iDict]->SetDataUnchecked(iRow,iCol,iTime,trackA);
1587                     } // if:  track A == 0
1588                   } // for: jDict
1589                 } // if:  fMasks and trackB > 0
1590               } // for: iDict
1591
1592             } // for: iTime
1593           } // for: iCol
1594         } // for: iRow
1595
1596       } // if:  doMerge
1597
1598       mergeSDigitsManager->RemoveDigits(iDet);
1599       mergeSDigitsManager->RemoveDictionaries(iDet);
1600
1601       if (fCompress) {
1602         digitsA->Compress(1,0);
1603         for (iDict = 0; iDict < kNDict; iDict++) {
1604           dictionaryA[iDict]->Compress(1,0);
1605         }
1606       }
1607
1608     } // for: detectors
1609
1610     // The next set of s-digits
1611     mergeSDigitsManager = (AliTRDdigitsManager *)
1612                           fSDigitsManagerList->After(mergeSDigitsManager);
1613
1614   } // while: mergeDigitsManagers
1615
1616   return kTRUE;
1617
1618 }
1619
1620 //_____________________________________________________________________________
1621 Bool_t AliTRDdigitizer::ConvertSDigits()
1622 {
1623   //
1624   // Converts s-digits to normal digits
1625   //
1626
1627   AliTRDdataArrayF *digitsIn = 0x0;
1628
1629   if (!fSDigitsManager->HasSDigits()) {
1630     AliError("No s-digits in digits manager");
1631     return kFALSE;
1632   }
1633
1634   // Loop through the detectors
1635   for (Int_t det = 0; det < AliTRDgeometry::Ndet(); det++) {
1636
1637     // Get the merged s-digits (signals)
1638     digitsIn = (AliTRDdataArrayF *) fSDigitsManager->GetDigits(det);
1639     if (!digitsIn->HasData()) {
1640       AliDebug(2,Form("No digits for det=%d",det));
1641       continue;
1642     }
1643
1644     // Convert the merged sdigits to digits
1645     if (!Signal2ADC(det,digitsIn)) {
1646       continue;
1647     }
1648
1649     // Copy the dictionary information to the output array
1650     if (!CopyDictionary(det)) {
1651       continue;
1652     }
1653
1654     // Delete
1655     fSDigitsManager->RemoveDigits(det);
1656     fSDigitsManager->RemoveDictionaries(det);
1657
1658     // Compress the arrays
1659     CompressOutputArrays(det);
1660
1661   } // for: detector numbers
1662
1663   return kTRUE;
1664
1665 }
1666
1667 //_____________________________________________________________________________
1668 Bool_t AliTRDdigitizer::CopyDictionary(Int_t det)
1669 {
1670   //
1671   // Copies the dictionary information from the s-digits arrays
1672   // to the output arrays
1673   //
1674
1675   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
1676   if (!calibration) {
1677     AliFatal("Could not get calibration object");
1678     return kFALSE;
1679   }
1680
1681   AliDebug(1,Form("Start copying dictionaries for detector=%d",det));
1682
1683   const Int_t kNDict = AliTRDdigitsManager::kNDict;
1684   AliTRDdataArrayI *dictionaryIn[kNDict];
1685   AliTRDdataArrayI *dictionaryOut[kNDict];
1686
1687   Int_t nRowMax    = fGeo->GetPadPlane(det)->GetNrows();
1688   Int_t nColMax    = fGeo->GetPadPlane(det)->GetNcols();
1689   Int_t nTimeTotal = calibration->GetNumberOfTimeBins();
1690
1691   Int_t row  = 0;
1692   Int_t col  = 0;
1693   Int_t time = 0;
1694   Int_t dict = 0;
1695
1696   for (dict = 0; dict < kNDict; dict++) {
1697
1698     dictionaryIn[dict]  = (AliTRDdataArrayI *) fSDigitsManager->GetDictionary(det,dict);
1699     dictionaryIn[dict]->Expand();
1700     dictionaryOut[dict] = (AliTRDdataArrayI *) fDigitsManager->GetDictionary(det,dict);
1701     dictionaryOut[dict]->Allocate(nRowMax,nColMax,nTimeTotal);
1702
1703     for (row = 0; row < nRowMax; row++) {
1704       for (col = 0; col < nColMax; col++) {
1705         for (time = 0; time < nTimeTotal; time++) {
1706
1707           Int_t track = dictionaryIn[dict]->GetDataUnchecked(row,col,time);
1708           dictionaryOut[dict]->SetDataUnchecked(row,col,time,track);
1709
1710         } // for: time
1711       } // for: col
1712     } // for: row
1713
1714   } // for: dictionaries
1715
1716   return kTRUE;
1717
1718 }
1719
1720 //_____________________________________________________________________________
1721 void AliTRDdigitizer::CompressOutputArrays(Int_t det)
1722 {
1723   //
1724   // Compress the output arrays
1725   //
1726
1727   const Int_t kNDict = AliTRDdigitsManager::kNDict;
1728
1729   AliTRDdataArrayS *digits     = 0x0;
1730   AliTRDdataArrayI *dictionary = 0x0;
1731
1732   if (fCompress) {
1733
1734     digits = (AliTRDdataArrayS *) fDigitsManager->GetDigits(det);
1735     digits->Compress(1,0);
1736     for (Int_t dict = 0; dict < kNDict; dict++) {
1737       dictionary = (AliTRDdataArrayI *) fDigitsManager->GetDictionary(det,dict);
1738       dictionary->Compress(1,0);
1739
1740     }
1741
1742   }
1743
1744 }
1745
1746 //_____________________________________________________________________________
1747 Bool_t AliTRDdigitizer::WriteDigits() const
1748 {
1749   //
1750   // Writes out the TRD-digits and the dictionaries
1751   //
1752
1753   // Write parameters
1754   fRunLoader->CdGAFile();
1755
1756   // Store the digits and the dictionary in the tree
1757   return fDigitsManager->WriteDigits();
1758
1759 }
1760
1761 //_____________________________________________________________________________
1762 void AliTRDdigitizer::InitOutput(Int_t iEvent)
1763 {
1764   //
1765   // Initializes the output branches
1766   //
1767
1768   fEvent = iEvent;
1769
1770   if (!fRunLoader) {
1771     AliError("Run Loader is NULL");
1772     return;
1773   }
1774
1775   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
1776   if (!loader) {
1777     AliError("Can not get TRD loader from Run Loader");
1778     return;
1779   }
1780
1781   TTree *tree = 0;
1782
1783   if (fSDigits) {
1784     // If we produce SDigits
1785     tree = loader->TreeS();
1786     if (!tree) {
1787       loader->MakeTree("S");
1788       tree = loader->TreeS();
1789     }
1790   }
1791   else {
1792     // If we produce Digits
1793     tree = loader->TreeD();
1794     if (!tree) {
1795       loader->MakeTree("D");
1796       tree = loader->TreeD();
1797     }
1798   }
1799
1800   fDigitsManager->SetEvent(iEvent);
1801   fDigitsManager->MakeBranch(tree);
1802
1803 }
1804
1805 //_____________________________________________________________________________
1806 Double_t AliTRDdigitizer::TimeStruct(Float_t vdrift, Double_t dist, Double_t z)
1807 {
1808   //
1809   // Applies the time structure of the drift cells (by C.Lippmann).
1810   // The drift time of electrons to the anode wires depends on the
1811   // distance to the wire (z) and on the position in the drift region.
1812   //
1813   // input :
1814   // dist = radial distance from (cathode) pad plane [cm]
1815   // z    = distance from anode wire (parallel to cathode planes) [cm]
1816   //
1817   // output :
1818   // tdrift = the drift time of an electron at the given position
1819   //
1820   // We interpolate between the drift time values at the two drift
1821   // velocities fVDlo and fVDhi, being smaller and larger than
1822   // fDriftVelocity. We use the two stored drift time maps fTimeStruct1
1823   // and fTimeStruct2, calculated for the two mentioned drift velocities.
1824   //
1825
1826   SampleTimeStruct(vdrift);
1827
1828   // Indices:
1829   Int_t r1 = (Int_t)(10 * dist);
1830   if (r1 <  0) r1 =  0;
1831   if (r1 > 37) r1 = 37;
1832   Int_t r2 = r1 + 1;
1833   if (r2 <  0) r2 =  0;
1834   if (r2 > 37) r2 = 37;
1835   const Int_t kz1 = ((Int_t)(100 * z / 2.5));
1836   const Int_t kz2 = kz1 + 1;
1837
1838   if ((r1  <  0) ||
1839       (r1  > 37) ||
1840       (kz1 <  0) ||
1841       (kz1 > 10)) {
1842     AliWarning(Form("Indices out of range: dist=%.2f, z=%.2f, r1=%d, kz1=%d"
1843                    ,dist,z,r1,kz1));
1844   }
1845
1846   const Float_t ky111 = fTimeStruct1[r1+38*kz1];
1847   const Float_t ky221 = ((r2 <= 37) && (kz2 <= 10))
1848                       ? fTimeStruct1[r2+38*kz2]
1849                       : fTimeStruct1[37+38*10];
1850   const Float_t ky121 = (kz2 <= 10)
1851                       ? fTimeStruct1[r1+38*kz2]
1852                       : fTimeStruct1[r1+38*10];
1853   const Float_t ky211 = (r2 <= 37)
1854                       ? fTimeStruct1[r2+38*kz1]
1855                       : fTimeStruct1[37+38*kz1];
1856
1857   // 2D Interpolation, lower drift time map
1858   const Float_t ky11  = (ky211-ky111)*10*dist + ky111 - (ky211-ky111)*r1;
1859   const Float_t ky21  = (ky221-ky121)*10*dist + ky121 - (ky221-ky121)*r1;
1860
1861   const Float_t ky112 = fTimeStruct2[r1+38*kz1];
1862   const Float_t ky222 = ((r2 <= 37) && (kz2 <= 10))
1863                       ? fTimeStruct2[r2+38*kz2]
1864                       : fTimeStruct2[37+38*10];
1865   const Float_t ky122 = (kz2 <= 10)
1866                       ? fTimeStruct2[r1+38*kz2]
1867                       : fTimeStruct2[r1+38*10];
1868   const Float_t ky212 = (r2 <= 37)
1869                       ? fTimeStruct2[r2+38*kz1]
1870                       : fTimeStruct2[37+38*kz1];
1871
1872   // 2D Interpolation, larger drift time map
1873   const Float_t ky12  = (ky212-ky112)*10*dist + ky112 - (ky212-ky112)*r1;
1874   const Float_t ky22  = (ky222-ky122)*10*dist + ky122 - (ky222-ky122)*r1;
1875
1876   // Dist now is the drift distance to the anode wires (negative if electrons are
1877   // between anode wire plane and cathode pad plane)
1878   dist -= AliTRDgeometry::AmThick() / 2.0;
1879
1880   // Get the drift times for the drift velocities fVDlo and fVDhi
1881   const Float_t ktdrift1 = ((TMath::Abs(dist) > 0.005) || (z > 0.005))
1882                          ? (ky21 - ky11) * 100 * z / 2.5 + ky11 - (ky21 - ky11) * kz1
1883                          : 0.0;
1884   const Float_t ktdrift2 = ((TMath::Abs(dist) > 0.005) || (z > 0.005))
1885                          ? (ky22 - ky12) * 100 * z / 2.5 + ky12 - (ky22 - ky12) * kz1
1886                          : 0.0;
1887
1888   // 1D Interpolation between the values at fVDlo and fVDhi
1889   Float_t a = (ktdrift2 - ktdrift1)
1890             / (fVDhi - fVDlo);
1891   Float_t b = ktdrift2 - a * fVDhi;
1892
1893   return a * vdrift + b;
1894
1895 }
1896
1897 //_____________________________________________________________________________
1898 void AliTRDdigitizer::SampleTimeStruct(Float_t vdrift)
1899 {
1900   //
1901   // Samples the timing structure of a drift cell
1902   // Drift Time data calculated with Garfield (by C.Lippmann)
1903   //
1904
1905   // Nothing to do
1906   if (vdrift == fTimeLastVdrift) {
1907     return;
1908   }
1909   fTimeLastVdrift = vdrift;
1910
1911   // Drift time maps are saved for some drift velocity values (in drift region):
1912   Float_t  fVDsmp[8];
1913   fVDsmp[0] = 1.032;
1914   fVDsmp[1] = 1.158;
1915   fVDsmp[2] = 1.299;
1916   fVDsmp[3] = 1.450;
1917   fVDsmp[4] = 1.610;
1918   fVDsmp[5] = 1.783;
1919   fVDsmp[6] = 1.959;
1920   fVDsmp[7] = 2.134;
1921
1922   if      (vdrift < fVDsmp[0]) {
1923     AliWarning(Form("Drift Velocity too small (%.3f<%.3f)",vdrift,fVDsmp[0]));
1924     vdrift = fVDsmp[0];
1925   }
1926   else if (vdrift > fVDsmp[7]) {
1927     AliWarning(Form("Drift Velocity too large (%.3f>%.3f)",vdrift,fVDsmp[6]));
1928     vdrift = fVDsmp[7];
1929   }
1930
1931   const Int_t ktimebin  = 38;
1932   const Int_t kZbin     = 11;
1933
1934   // Garfield simulation at UD = -1500V; vd = 1.032cm/microsec, <driftfield> = 525V/cm
1935   Float_t time1500[ktimebin][kZbin] =
1936       {{0.09170, 0.09205, 0.09306, 0.09475, 0.09716, 0.10035,
1937         0.10445, 0.10993, 0.11838, 0.13986, 0.37858},
1938        {0.06588, 0.06626, 0.06739, 0.06926, 0.07186, 0.07524,
1939         0.07951, 0.08515, 0.09381, 0.11601, 0.35673},
1940        {0.03946, 0.04003, 0.04171, 0.04435, 0.04780, 0.05193,
1941         0.05680, 0.06306, 0.07290, 0.09873, 0.34748},
1942        {0.01151, 0.01283, 0.01718, 0.02282, 0.02880, 0.03479,
1943         0.04098, 0.04910, 0.06413, 0.10567, 0.36897},
1944        {0.01116, 0.01290, 0.01721, 0.02299, 0.02863, 0.03447,
1945         0.04074, 0.04984, 0.06839, 0.11625, 0.37277},
1946        {0.03919, 0.03974, 0.04131, 0.04380, 0.04703, 0.05102,
1947         0.05602, 0.06309, 0.07651, 0.10938, 0.36838},
1948        {0.06493, 0.06560, 0.06640, 0.06802, 0.07051, 0.07392,
1949         0.07853, 0.08510, 0.09690, 0.12621, 0.38058},
1950        {0.09174, 0.09186, 0.09225, 0.09303, 0.09477, 0.00000,
1951         0.11205, 0.11952, 0.13461, 0.16984, 0.43017},
1952        {0.14356, 0.14494, 0.14959, 0.16002, 0.18328, 0.27981,
1953         0.22785, 0.21240, 0.21948, 0.25965, 0.52392},
1954        {0.23120, 0.23366, 0.24046, 0.25422, 0.28071, 0.36914,
1955         0.32999, 0.31208, 0.31772, 0.35804, 0.62249},
1956        {0.32686, 0.32916, 0.33646, 0.35053, 0.37710, 0.46292,
1957         0.42773, 0.40948, 0.41497, 0.45527, 0.71955},
1958        {0.42353, 0.42583, 0.43317, 0.44727, 0.47380, 0.55884,
1959         0.52479, 0.50650, 0.51194, 0.55225, 0.81658},
1960        {0.52038, 0.52271, 0.53000, 0.54415, 0.57064, 0.65545,
1961         0.62172, 0.60341, 0.60885, 0.64915, 0.91339},
1962        {0.61724, 0.61953, 0.62694, 0.64098, 0.66756, 0.75226,
1963         0.71862, 0.70030, 0.70575, 0.74604, 1.01035},
1964        {0.71460, 0.71678, 0.72376, 0.73786, 0.76447, 0.84913,
1965         0.81551, 0.79720, 0.80264, 0.84292, 1.10723},
1966        {0.81101, 0.81334, 0.82066, 0.83475, 0.86127, 0.94599,
1967         0.91240, 0.89408, 0.89952, 0.93981, 1.20409},
1968        {0.90788, 0.91023, 0.91752, 0.93163, 0.95815, 1.04293,
1969         1.00929, 0.99096, 0.99640, 1.03669, 1.30106},
1970        {1.00477, 1.00707, 1.01449, 1.02852, 1.05504, 1.13976,
1971         1.10617, 1.08784, 1.09329, 1.13358, 1.39796},
1972        {1.10166, 1.10397, 1.11130, 1.12541, 1.15257, 1.23672,
1973         1.20307, 1.18472, 1.19018, 1.23046, 1.49486},
1974        {1.19854, 1.20084, 1.20818, 1.22235, 1.24885, 1.33355,
1975         1.29992, 1.28156, 1.28707, 1.32735, 1.59177},
1976        {1.29544, 1.29780, 1.30507, 1.31917, 1.34575, 1.43073,
1977         1.39681, 1.37851, 1.38396, 1.42377, 1.68854},
1978        {1.39236, 1.39462, 1.40205, 1.41607, 1.44259, 1.52745,
1979         1.49368, 1.47541, 1.48083, 1.52112, 1.78546},
1980        {1.49314, 1.49149, 1.49885, 1.51297, 1.53949, 1.62420,
1981         1.59016, 1.57231, 1.57772, 1.61800, 1.88048},
1982        {1.58610, 1.58839, 1.59572, 1.60983, 1.63635, 1.72109,
1983         1.68651, 1.66921, 1.67463, 1.71489, 1.97918},
1984        {1.68400, 1.68529, 1.69261, 1.70671, 1.73331, 1.81830,
1985         1.78341, 1.76605, 1.77150, 1.81179, 2.07608},
1986        {1.77991, 1.78215, 1.78952, 1.80385, 1.83014, 1.91486,
1987         1.88128, 1.86215, 1.86837, 1.90865, 2.17304},
1988        {1.87674, 1.87904, 1.88647, 1.90052, 1.92712, 2.01173,
1989         1.97812, 1.95905, 1.96527, 2.00710, 2.26979},
1990        {1.97369, 1.97594, 1.98326, 1.99869, 2.02442, 2.10865,
1991         2.07501, 2.05666, 2.06214, 2.10243, 2.36669},
1992        {2.07052, 2.07281, 2.08016, 2.09425, 2.12132, 2.20555,
1993         2.17182, 2.15341, 2.15904, 2.19933, 2.46363},
1994        {2.16742, 2.16971, 2.17707, 2.19114, 2.21766, 2.30240,
1995         2.26877, 2.25015, 2.25573, 2.29586, 2.56060},
1996        {2.26423, 2.26659, 2.27396, 2.28803, 2.31456, 2.40828,
1997         2.36567, 2.34705, 2.35282, 2.39765, 2.65744},
1998        {2.36153, 2.36349, 2.37330, 2.38501, 2.41159, 2.49940,
1999         2.46257, 2.44420, 2.44843, 2.48987, 2.75431},
2000        {2.46558, 2.46035, 2.46822, 2.48181, 2.50849, 2.59630,
2001         2.55947, 2.54112, 2.54513, 2.58677, 2.85094},
2002        {2.56248, 2.55723, 2.56486, 2.57871, 2.60520, 2.68998,
2003         2.65626, 2.63790, 2.64316, 2.68360, 2.94813},
2004        {2.65178, 2.65441, 2.66153, 2.67556, 2.70210, 2.78687,
2005         2.75319, 2.73463, 2.74032, 2.78060, 3.04503},
2006        {2.74868, 2.75131, 2.75870, 2.77245, 2.79385, 2.88700,
2007         2.85009, 2.83177, 2.83723, 2.87750, 3.14193},
2008        {2.84574, 2.84789, 2.85560, 2.86935, 2.89075, 2.98060,
2009         2.94576, 2.92868, 2.93356, 2.97436, 3.23868},
2010        {2.94239, 2.94469, 2.95221, 2.96625, 2.99345, 3.07747,
2011         3.04266, 3.02545, 3.03051, 3.07118, 3.33555}};
2012
2013   // Garfield simulation at UD = -1600V; vd = 1.158cm/microsec, <driftfield> = 558V/cm
2014   Float_t time1600[ktimebin][kZbin] =
2015       {{0.09169, 0.09203, 0.09305, 0.09473, 0.09714, 0.10032,
2016         0.10441, 0.10990, 0.11835, 0.13986, 0.37845},
2017        {0.06589, 0.06626, 0.06738, 0.06924, 0.07184, 0.07521,
2018         0.07947, 0.08512, 0.09379, 0.11603, 0.35648},
2019        {0.03947, 0.04003, 0.04171, 0.04434, 0.04778, 0.05190,
2020         0.05678, 0.06304, 0.07292, 0.09876, 0.34759},
2021        {0.01111, 0.01281, 0.01718, 0.02281, 0.02879, 0.03477,
2022         0.04097, 0.04910, 0.06415, 0.10573, 0.36896},
2023        {0.01120, 0.01311, 0.01721, 0.02279, 0.02862, 0.03446,
2024         0.04074, 0.04981, 0.06825, 0.11595, 0.37255},
2025        {0.03919, 0.03980, 0.04132, 0.04380, 0.04704, 0.05102,
2026         0.05602, 0.06302, 0.07633, 0.10896, 0.36743},
2027        {0.06531, 0.06528, 0.06631, 0.06805, 0.07053, 0.07392,
2028         0.07853, 0.08505, 0.09669, 0.12578, 0.37967},
2029        {0.09157, 0.09171, 0.09216, 0.09301, 0.09475, 0.00000,
2030         0.11152, 0.11879, 0.13352, 0.16802, 0.42750},
2031        {0.13977, 0.14095, 0.14509, 0.15433, 0.17534, 0.26406,
2032         0.21660, 0.20345, 0.21113, 0.25067, 0.51434},
2033        {0.21816, 0.22041, 0.22631, 0.23850, 0.26208, 0.34340,
2034         0.30755, 0.29237, 0.29878, 0.33863, 0.60258},
2035        {0.30344, 0.30547, 0.31241, 0.32444, 0.34809, 0.42696,
2036         0.39464, 0.37919, 0.38546, 0.42530, 0.68926},
2037        {0.38969, 0.39164, 0.39810, 0.41059, 0.43441, 0.51246,
2038         0.48112, 0.46562, 0.47191, 0.51172, 0.77558},
2039        {0.47592, 0.47799, 0.48442, 0.49689, 0.52061, 0.59855,
2040         0.56752, 0.55201, 0.55826, 0.59808, 0.86202},
2041        {0.56226, 0.56428, 0.57074, 0.58324, 0.60696, 0.68483,
2042         0.65388, 0.63837, 0.64461, 0.68445, 0.94830},
2043        {0.64881, 0.65063, 0.65709, 0.66958, 0.69331, 0.77117,
2044         0.74023, 0.72472, 0.73098, 0.77079, 1.03486},
2045        {0.73506, 0.73698, 0.74344, 0.75596, 0.77964, 0.85751,
2046         0.82658, 0.81107, 0.81731, 0.85712, 1.12106},
2047        {0.82132, 0.82333, 0.82979, 0.84228, 0.86608, 0.94386,
2048         0.91293, 0.89742, 0.90367, 0.94335, 1.20737},
2049        {0.90767, 0.90968, 0.91614, 0.92863, 0.95236, 1.03021,
2050         0.99928, 0.98377, 0.99001, 1.02984, 1.29371},
2051        {0.99410, 0.99602, 1.00257, 1.01498, 1.03869, 1.11720,
2052         1.08563, 1.07011, 1.07637, 1.11621, 1.37873},
2053        {1.08036, 1.08240, 1.08884, 1.10138, 1.12504, 1.20301,
2054         1.17198, 1.15647, 1.16272, 1.20255, 1.46651},
2055        {1.16670, 1.16872, 1.17525, 1.18783, 1.21139, 1.28934,
2056         1.25833, 1.24281, 1.24909, 1.28889, 1.55275},
2057        {1.25307, 1.25510, 1.26153, 1.27404, 1.29773, 1.37584,
2058         1.34469, 1.32916, 1.33536, 1.37524, 1.63915},
2059        {1.33942, 1.34146, 1.34788, 1.36040, 1.38410, 1.46438,
2060         1.43105, 1.41537, 1.42176, 1.46158, 1.72538},
2061        {1.42579, 1.42782, 1.43458, 1.44674, 1.47043, 1.55085,
2062         1.51675, 1.50168, 1.50810, 1.54793, 1.81174},
2063        {1.51207, 1.51454, 1.52060, 1.53307, 1.55684, 1.63478,
2064         1.60336, 1.58820, 1.59446, 1.63414, 1.89814},
2065        {1.59856, 1.60047, 1.60693, 1.61942, 1.64317, 1.72257,
2066         1.69008, 1.67454, 1.68080, 1.72063, 1.98433},
2067        {1.68481, 1.68682, 1.69330, 1.70584, 1.72949, 1.80752,
2068         1.77643, 1.76089, 1.76716, 1.80692, 2.07069},
2069        {1.77117, 1.77319, 1.77969, 1.79260, 1.81583, 1.89376,
2070         1.86226, 1.84720, 1.85355, 1.89256, 2.15343},
2071        {1.85748, 1.85967, 1.86605, 1.87848, 1.90222, 1.98010,
2072         1.94913, 1.93271, 1.93981, 1.97968, 2.24355},
2073        {1.94386, 1.94587, 1.95233, 1.96484, 1.98854, 2.06646,
2074         2.03542, 2.01755, 2.02617, 2.06604, 2.32993},
2075        {2.03022, 2.03230, 2.03868, 2.05134, 2.07488, 2.15367,
2076         2.12178, 2.10391, 2.11252, 2.15432, 2.41623},
2077        {2.11656, 2.11857, 2.12505, 2.13772, 2.16147, 2.23919,
2078         2.20817, 2.19265, 2.20744, 2.23872, 2.49996},
2079        {2.20291, 2.20611, 2.21137, 2.22387, 2.24758, 2.32563,
2080         2.29450, 2.27901, 2.28525, 2.32507, 2.58897},
2081        {2.28922, 2.29172, 2.29774, 2.31345, 2.33400, 2.41287,
2082         2.38086, 2.36535, 2.37160, 2.40869, 2.67113},
2083        {2.37572, 2.37764, 2.38410, 2.39803, 2.42046, 2.49817,
2084         2.46721, 2.45171, 2.45794, 2.49505, 2.76061},
2085        {2.46190, 2.46396, 2.47043, 2.48340, 2.50665, 2.58453,
2086         2.55357, 2.53728, 2.54430, 2.58407, 2.84816},
2087        {2.54833, 2.55032, 2.55679, 2.56976, 2.59312, 2.67103,
2088         2.63993, 2.62364, 2.63062, 2.67040, 2.93444},
2089        {2.63456, 2.63660, 2.64304, 2.65555, 2.67938, 2.75739,
2090         2.72629, 2.71064, 2.71688, 2.75671, 3.01886}};
2091
2092   // Garfield simulation at UD = -1700V; vd = 1.299cm/microsec, <driftfield> = 590V/cm
2093   Float_t time1700[ktimebin][kZbin] =
2094       {{0.09167, 0.09201, 0.09302, 0.09471, 0.09712, 0.10029,
2095         0.10438, 0.10986, 0.11832, 0.13986, 0.37824},
2096        {0.06591, 0.06626, 0.06736, 0.06923, 0.07183, 0.07519,
2097         0.07944, 0.08511, 0.09378, 0.11603, 0.35625},
2098        {0.03946, 0.04003, 0.04170, 0.04433, 0.04777, 0.05189,
2099         0.05676, 0.06301, 0.07291, 0.09880, 0.34724},
2100        {0.01110, 0.01281, 0.01718, 0.02280, 0.02879, 0.03476,
2101         0.04096, 0.04910, 0.06417, 0.10582, 0.36861},
2102        {0.01115, 0.01294, 0.01721, 0.02276, 0.02862, 0.03447,
2103         0.04074, 0.04980, 0.06812, 0.11565, 0.37231},
2104        {0.03920, 0.03974, 0.04133, 0.04381, 0.04706, 0.05105,
2105         0.05603, 0.06299, 0.07618, 0.10860, 0.36646},
2106        {0.06498, 0.06529, 0.06634, 0.06808, 0.07055, 0.07395,
2107         0.07852, 0.08500, 0.09650, 0.12532, 0.37850},
2108        {0.09143, 0.09159, 0.09207, 0.09297, 0.09473, 0.00000,
2109         0.11102, 0.11809, 0.13245, 0.16627, 0.42496},
2110        {0.13646, 0.13750, 0.14112, 0.14926, 0.16806, 0.24960,
2111         0.20627, 0.19536, 0.20366, 0.24256, 0.50557},
2112        {0.20678, 0.20848, 0.21384, 0.22450, 0.24552, 0.32018,
2113         0.28740, 0.27477, 0.28196, 0.32128, 0.58475},
2114        {0.28287, 0.28461, 0.29020, 0.30108, 0.32224, 0.39467,
2115         0.36500, 0.35217, 0.35926, 0.39860, 0.66194},
2116        {0.35972, 0.36145, 0.36713, 0.37797, 0.39912, 0.47091,
2117         0.44212, 0.42925, 0.43632, 0.47563, 0.73892},
2118        {0.43667, 0.43841, 0.44413, 0.45494, 0.47607, 0.54780,
2119         0.51912, 0.50627, 0.51334, 0.55254, 0.81595},
2120        {0.51365, 0.51540, 0.52101, 0.53193, 0.55305, 0.62463,
2121         0.59617, 0.58328, 0.59035, 0.62965, 0.89303},
2122        {0.59064, 0.59240, 0.59801, 0.60893, 0.63009, 0.70169,
2123         0.67317, 0.66028, 0.66735, 0.70666, 0.96995},
2124        {0.66765, 0.66939, 0.67501, 0.68592, 0.70724, 0.77863,
2125         0.75016, 0.73728, 0.74435, 0.78366, 1.04696},
2126        {0.74464, 0.74636, 0.75200, 0.76293, 0.78405, 0.85561,
2127         0.82716, 0.81427, 0.82133, 0.86064, 1.12396},
2128        {0.82165, 0.82340, 0.82902, 0.83991, 0.86104, 0.93266,
2129         0.90414, 0.89128, 0.89834, 0.93763, 1.20100},
2130        {0.89863, 0.90042, 0.90659, 0.91705, 0.93805, 1.00960,
2131         0.98115, 0.96825, 0.97533, 1.01462, 1.27801},
2132        {0.97563, 0.97740, 0.98310, 0.99391, 1.01504, 1.08659,
2133         1.05814, 1.04526, 1.05233, 1.09163, 1.35503},
2134        {1.05276, 1.05451, 1.06002, 1.07090, 1.09099, 1.16357,
2135         1.13516, 1.12225, 1.12933, 1.16863, 1.43195},
2136        {1.12977, 1.13138, 1.13700, 1.14792, 1.16797, 1.24061,
2137         1.21212, 1.19926, 1.20626, 1.24554, 1.50900},
2138        {1.20664, 1.20839, 1.21400, 1.22490, 1.24606, 1.31772,
2139         1.28914, 1.27382, 1.28329, 1.32262, 1.58550},
2140        {1.28367, 1.28541, 1.29099, 1.30189, 1.32312, 1.39460,
2141         1.36612, 1.34924, 1.36030, 1.39961, 1.66310},
2142        {1.36064, 1.36249, 1.36799, 1.37896, 1.40004, 1.48030,
2143         1.44314, 1.43032, 1.43731, 1.47659, 1.73442},
2144        {1.43762, 1.43937, 1.44497, 1.45618, 1.47704, 1.54932,
2145         1.52012, 1.50725, 1.51430, 1.55357, 1.81708},
2146        {1.51462, 1.51937, 1.52203, 1.53316, 1.55403, 1.62572,
2147         1.59713, 1.58424, 1.59128, 1.63061, 1.89406},
2148        {1.59162, 1.59338, 1.59947, 1.60989, 1.63103, 1.70270,
2149         1.67411, 1.66124, 1.66799, 1.70759, 1.97103},
2150        {1.66874, 1.67037, 1.67597, 1.68687, 1.70814, 1.77969,
2151         1.75112, 1.73806, 1.74530, 1.78457, 2.04794},
2152        {1.74693, 1.74749, 1.75297, 1.76476, 1.78500, 1.85667,
2153         1.82811, 1.81504, 1.82101, 1.86161, 2.12492},
2154        {1.82260, 1.82437, 1.82995, 1.84174, 1.86202, 1.93372,
2155         1.90509, 1.89202, 1.89930, 1.93859, 2.20189},
2156        {1.89964, 1.90135, 1.90693, 1.91789, 1.93952, 2.01080,
2157         1.98207, 1.96921, 1.97628, 2.01384, 2.27887},
2158        {1.97662, 1.97917, 1.98611, 1.99487, 2.01601, 2.08778,
2159         2.05846, 2.04623, 2.05330, 2.09244, 2.35585},
2160        {2.05359, 2.05615, 2.06309, 2.07187, 2.09867, 2.16459,
2161         2.13610, 2.12322, 2.13029, 2.16942, 2.43199},
2162        {2.13063, 2.13233, 2.13795, 2.14886, 2.17008, 2.24199,
2163         2.21310, 2.20020, 2.20727, 2.24659, 2.50983},
2164        {2.20761, 2.20931, 2.21955, 2.22624, 2.24708, 2.32147,
2165         2.29009, 2.27725, 2.28276, 2.32359, 2.58680},
2166        {2.28459, 2.29108, 2.29202, 2.30286, 2.32007, 2.39559,
2167         2.36683, 2.35422, 2.36119, 2.40058, 2.66081},
2168        {2.36153, 2.36806, 2.36889, 2.37985, 2.40092, 2.47828,
2169         2.44381, 2.43099, 2.43819, 2.47750, 2.73779}};
2170
2171   // Garfield simulation at UD = -1800V; vd = 1.450cm/microsec, <driftfield> = 623V/cm
2172   Float_t time1800[ktimebin][kZbin] =
2173       {{0.09166, 0.09199, 0.09300, 0.09470, 0.09709, 0.10026,
2174         0.10434, 0.10983, 0.11831, 0.13987, 0.37802},
2175        {0.06585, 0.06623, 0.06735, 0.06921, 0.07180, 0.07520,
2176         0.07941, 0.08507, 0.09376, 0.11604, 0.35624},
2177        {0.03945, 0.04004, 0.04169, 0.04432, 0.04776, 0.05187,
2178         0.05674, 0.06300, 0.07290, 0.09884, 0.34704},
2179        {0.01108, 0.01287, 0.01717, 0.02280, 0.02880, 0.03476,
2180         0.04095, 0.04909, 0.06419, 0.10589, 0.36843},
2181        {0.01115, 0.01287, 0.01720, 0.02276, 0.02862, 0.03448,
2182         0.04073, 0.04973, 0.06799, 0.11535, 0.37224},
2183        {0.03918, 0.03975, 0.04134, 0.04382, 0.04707, 0.05105,
2184         0.05603, 0.06296, 0.07605, 0.10822, 0.36560},
2185        {0.06498, 0.06532, 0.06635, 0.06809, 0.07058, 0.07395,
2186         0.07855, 0.08495, 0.09632, 0.12488, 0.37730},
2187        {0.09130, 0.09160, 0.09199, 0.09300, 0.09472, 0.00000,
2188         0.11059, 0.11747, 0.13146, 0.16462, 0.42233},
2189        {0.13364, 0.13449, 0.13767, 0.14481, 0.16147, 0.23635,
2190         0.19706, 0.18812, 0.19704, 0.23520, 0.49749},
2191        {0.19698, 0.19844, 0.20311, 0.21236, 0.23082, 0.29920,
2192         0.26936, 0.25927, 0.26732, 0.30601, 0.56871},
2193        {0.26518, 0.26670, 0.27160, 0.28099, 0.29955, 0.36597,
2194         0.33885, 0.32858, 0.33653, 0.37524, 0.63801},
2195        {0.33441, 0.33553, 0.34040, 0.34987, 0.36841, 0.43428,
2196         0.40797, 0.39763, 0.40556, 0.44425, 0.70698},
2197        {0.40296, 0.40447, 0.40934, 0.41881, 0.43737, 0.50306,
2198         0.47695, 0.46662, 0.47455, 0.51329, 0.77600},
2199        {0.47296, 0.47344, 0.47830, 0.48779, 0.50632, 0.57200,
2200         0.54593, 0.53559, 0.54351, 0.58222, 0.84489},
2201        {0.54089, 0.54264, 0.54727, 0.55673, 0.57529, 0.64094,
2202         0.61490, 0.60457, 0.61249, 0.65118, 0.91394},
2203        {0.60987, 0.61138, 0.61624, 0.62573, 0.64428, 0.70989,
2204         0.68397, 0.67354, 0.68147, 0.72015, 0.98291},
2205        {0.67883, 0.68035, 0.68521, 0.69469, 0.71324, 0.77896,
2206         0.75287, 0.74251, 0.75043, 0.78912, 1.04458},
2207        {0.74780, 0.74932, 0.75421, 0.76367, 0.78221, 0.84785,
2208         0.82185, 0.81148, 0.81941, 0.85811, 1.12085},
2209        {0.81690, 0.81830, 0.82316, 0.83263, 0.85120, 0.91683,
2210         0.89077, 0.88045, 0.88837, 0.92707, 1.18976},
2211        {0.88574, 0.88726, 0.89228, 0.90198, 0.92017, 0.98578,
2212         0.95974, 0.94947, 0.95734, 0.99604, 1.25873},
2213        {0.95493, 0.95624, 0.96110, 0.97058, 0.98913, 1.05481,
2214         1.02873, 1.01839, 1.02631, 1.06503, 1.32772},
2215        {1.02392, 1.02524, 1.03008, 1.03955, 1.05810, 1.12378,
2216         1.09757, 1.08605, 1.09530, 1.13399, 1.39669},
2217        {1.09270, 1.09418, 1.09911, 1.10854, 1.12714, 1.19281,
2218         1.16502, 1.15633, 1.16427, 1.20271, 1.46574},
2219        {1.16168, 1.16323, 1.16801, 1.17772, 1.19604, 1.26190,
2220         1.23399, 1.22531, 1.23323, 1.27194, 1.53475},
2221        {1.23061, 1.23214, 1.23698, 1.24669, 1.26503, 1.33073,
2222         1.30461, 1.29428, 1.30220, 1.34091, 1.60372},
2223        {1.29960, 1.30110, 1.30596, 1.31544, 1.33398, 1.39962,
2224         1.37228, 1.36323, 1.37121, 1.40988, 1.67273},
2225        {1.36851, 1.37007, 1.37512, 1.38441, 1.40297, 1.46865,
2226         1.44256, 1.43222, 1.44017, 1.47878, 1.74155},
2227        {1.43752, 1.43904, 1.44773, 1.45338, 1.47220, 1.53759,
2228         1.51136, 1.50119, 1.50914, 1.54775, 1.81050},
2229        {1.50646, 1.50802, 1.51288, 1.52237, 1.54097, 1.60697,
2230         1.58049, 1.57018, 1.57811, 1.61678, 1.87947},
2231        {1.57545, 1.57720, 1.58185, 1.59134, 1.60996, 1.67787,
2232         1.64929, 1.63914, 1.64707, 1.68570, 1.94851},
2233        {1.64442, 1.64617, 1.65081, 1.66035, 1.67893, 1.74684,
2234         1.71826, 1.70745, 1.71604, 1.75310, 2.01748},
2235        {1.71337, 1.71513, 1.71978, 1.72932, 1.74645, 1.81346,
2236         1.78739, 1.77642, 1.78501, 1.82151, 2.08644},
2237        {1.78238, 1.78410, 1.78876, 1.79824, 1.81678, 1.88332,
2238         1.85639, 1.84262, 1.85397, 1.89270, 2.15533},
2239        {1.85135, 1.85306, 1.85778, 1.86728, 1.88580, 1.95615,
2240         1.92536, 1.91171, 1.92283, 1.96165, 2.22428},
2241        {1.92774, 1.92184, 1.92672, 1.93618, 1.95477, 2.02048,
2242         1.99427, 1.98068, 1.99192, 2.03062, 2.29338},
2243        {1.98929, 1.99081, 1.99567, 2.00515, 2.02373, 2.08987,
2244         2.06332, 2.05249, 2.05939, 2.09928, 2.36227},
2245        {2.05829, 2.05978, 2.06464, 2.07414, 2.09272, 2.15850,
2246         2.12928, 2.12194, 2.12987, 2.16825, 2.43083},
2247        {2.12726, 2.12869, 2.13360, 2.14425, 2.16160, 2.22872,
2248         2.20118, 2.19078, 2.19876, 2.23416, 2.50007}};
2249
2250   // Garfield simulation at UD = -1900V; vd = 1.610cm/microsec, <driftfield> = 655V/cm
2251   Float_t time1900[ktimebin][kZbin] =
2252       {{0.09166, 0.09198, 0.09298, 0.09467, 0.09707, 0.10023,
2253         0.10431, 0.10980, 0.11828, 0.13988, 0.37789},
2254        {0.06584, 0.06622, 0.06735, 0.06920, 0.07179, 0.07514,
2255         0.07938, 0.08505, 0.09374, 0.11606, 0.35599},
2256        {0.03945, 0.04002, 0.04169, 0.04432, 0.04775, 0.05185,
2257         0.05672, 0.06298, 0.07290, 0.09888, 0.34695},
2258        {0.01109, 0.01281, 0.01717, 0.02279, 0.02878, 0.03476,
2259         0.04094, 0.04909, 0.06421, 0.10597, 0.36823},
2260        {0.01115, 0.01287, 0.01720, 0.02294, 0.02862, 0.03448,
2261         0.04074, 0.04973, 0.06783, 0.11506, 0.37206},
2262        {0.03940, 0.03975, 0.04135, 0.04386, 0.04708, 0.05106,
2263         0.05604, 0.06293, 0.07592, 0.10787, 0.36484},
2264        {0.06500, 0.06534, 0.06638, 0.06811, 0.07060, 0.07413,
2265         0.07852, 0.08491, 0.09614, 0.12446, 0.37626},
2266        {0.09119, 0.09140, 0.09194, 0.09293, 0.09471, 0.00000,
2267         0.11013, 0.11685, 0.13050, 0.16302, 0.41991},
2268        {0.13111, 0.13190, 0.13466, 0.14091, 0.15554, 0.22409,
2269         0.18846, 0.18167, 0.19113, 0.22854, 0.48995},
2270        {0.18849, 0.18975, 0.19380, 0.20185, 0.21797, 0.28050,
2271         0.25368, 0.24575, 0.25446, 0.29249, 0.55442},
2272        {0.24995, 0.25125, 0.25563, 0.26366, 0.27986, 0.34065,
2273         0.31605, 0.30815, 0.31680, 0.35482, 0.61697},
2274        {0.31187, 0.31324, 0.31745, 0.32580, 0.34205, 0.40217,
2275         0.37825, 0.37031, 0.37897, 0.41696, 0.67890},
2276        {0.37401, 0.37531, 0.37955, 0.38777, 0.40395, 0.46408,
2277         0.44037, 0.43242, 0.44108, 0.47906, 0.74122},
2278        {0.43610, 0.43741, 0.44161, 0.44986, 0.46604, 0.52614,
2279         0.50248, 0.49452, 0.50316, 0.54116, 0.80326},
2280        {0.49820, 0.49988, 0.50372, 0.51196, 0.52814, 0.58822,
2281         0.56459, 0.55661, 0.56527, 0.60326, 0.86526},
2282        {0.56032, 0.56161, 0.56582, 0.57408, 0.59024, 0.65032,
2283         0.62670, 0.61872, 0.62737, 0.66537, 0.92749},
2284        {0.62240, 0.62371, 0.62792, 0.63624, 0.65236, 0.71241,
2285         0.68881, 0.68081, 0.68947, 0.72750, 0.98941},
2286        {0.68449, 0.68581, 0.69002, 0.69828, 0.71444, 0.77452,
2287         0.75089, 0.74295, 0.75157, 0.78957, 1.05157},
2288        {0.74660, 0.74790, 0.75212, 0.76036, 0.77654, 0.83748,
2289         0.81299, 0.80501, 0.81193, 0.85168, 1.11375},
2290        {0.80870, 0.81017, 0.81423, 0.82246, 0.83867, 0.89908,
2291         0.87509, 0.86660, 0.87577, 0.91376, 1.17586},
2292        {0.87080, 0.87233, 0.87632, 0.88458, 0.90074, 0.96083,
2293         0.93718, 0.92922, 0.93787, 0.97588, 1.23794},
2294        {0.93291, 0.93422, 0.93844, 0.94667, 0.96293, 1.02295,
2295         0.99929, 0.99127, 0.99997, 1.03797, 1.29995},
2296        {0.99500, 0.99645, 1.00308, 1.00877, 1.02497, 1.08504,
2297         1.06140, 1.05343, 1.06203, 1.10006, 1.36216},
2298        {1.05712, 1.05926, 1.06262, 1.07092, 1.08706, 1.14754,
2299         1.12350, 1.11550, 1.12417, 1.16218, 1.42427},
2300        {1.11921, 1.12059, 1.12473, 1.13297, 1.14916, 1.21140,
2301         1.18560, 1.17284, 1.18625, 1.22414, 1.48629},
2302        {1.18140, 1.18262, 1.18690, 1.19508, 1.21125, 1.27139,
2303         1.24164, 1.23495, 1.24838, 1.28634, 1.54852},
2304        {1.24340, 1.24473, 1.24901, 1.25732, 1.27336, 1.33358,
2305         1.30793, 1.30179, 1.31047, 1.34848, 1.61066},
2306        {1.30551, 1.30684, 1.31104, 1.32056, 1.33553, 1.39609,
2307         1.37004, 1.36392, 1.37045, 1.41057, 1.67259},
2308        {1.36755, 1.36892, 1.37315, 1.39148, 1.39755, 1.45820,
2309         1.43215, 1.42602, 1.43467, 1.47268, 1.73477},
2310        {1.42966, 1.43101, 1.43549, 1.45359, 1.45976, 1.52031,
2311         1.49601, 1.48811, 1.49677, 1.53477, 1.79691},
2312        {1.49180, 1.49321, 1.49760, 1.51570, 1.52175, 1.58185,
2313         1.55771, 1.55023, 1.55888, 1.59672, 1.85501},
2314        {1.55391, 1.55527, 1.55943, 1.57782, 1.58391, 1.64395,
2315         1.62008, 1.61233, 1.62085, 1.65883, 1.92091},
2316        {1.61599, 1.61732, 1.62154, 1.63993, 1.64612, 1.70608,
2317         1.68237, 1.67108, 1.68301, 1.72110, 1.97931},
2318        {1.67808, 1.67948, 1.68364, 1.70204, 1.70823, 1.76858,
2319         1.74404, 1.73539, 1.74512, 1.78321, 2.04522},
2320        {1.74019, 1.74152, 1.74573, 1.76415, 1.77015, 1.83040,
2321         1.80615, 1.79366, 1.80723, 1.84509, 2.10742},
2322        {1.80235, 1.80362, 1.80783, 1.82626, 1.83227, 1.89246,
2323         1.86795, 1.85405, 1.86938, 1.90720, 2.16953},
2324        {1.86442, 1.86572, 1.86994, 1.88837, 1.89438, 1.95445,
2325         1.93006, 1.92283, 1.93148, 1.96931, 2.23147},
2326        {1.92700, 1.92783, 1.93197, 1.95049, 1.95649, 2.01668,
2327         1.99217, 1.98486, 1.99352, 2.03143, 2.29358}};
2328
2329   // Garfield simulation at UD = -2000V; vd = 1.783cm/microsec, <driftfield> = 688V/cm
2330   Float_t time2000[ktimebin][kZbin] =
2331       {{0.09176, 0.09196, 0.09296, 0.09465, 0.09704, 0.10020,
2332         0.10427, 0.10977, 0.11825, 0.13988, 0.37774},
2333        {0.06583, 0.06620, 0.06735, 0.06918, 0.07177, 0.07513,
2334         0.07936, 0.08503, 0.09372, 0.11606, 0.35586},
2335        {0.03944, 0.04001, 0.04170, 0.04431, 0.04774, 0.05184,
2336         0.05670, 0.06296, 0.07291, 0.09893, 0.34680},
2337        {0.01108, 0.01281, 0.01719, 0.02279, 0.02879, 0.03474,
2338         0.04093, 0.04908, 0.06422, 0.10605, 0.36800},
2339        {0.01114, 0.01287, 0.01720, 0.02276, 0.02863, 0.03449,
2340         0.04073, 0.04970, 0.06774, 0.11478, 0.37179},
2341        {0.03925, 0.03977, 0.04135, 0.04386, 0.04711, 0.05108,
2342         0.05604, 0.06290, 0.07580, 0.10748, 0.36386},
2343        {0.06501, 0.06536, 0.06640, 0.06814, 0.07062, 0.07398,
2344         0.07852, 0.08487, 0.09598, 0.12405, 0.37519},
2345        {0.09109, 0.09127, 0.09188, 0.09292, 0.09472, 0.00000,
2346         0.10964, 0.11630, 0.12960, 0.16150, 0.41765},
2347        {0.12898, 0.12968, 0.13209, 0.13749, 0.15034, 0.21286,
2348         0.18088, 0.17590, 0.18591, 0.22254, 0.48315},
2349        {0.18122, 0.18227, 0.18574, 0.19263, 0.20674, 0.26376,
2350         0.23960, 0.23375, 0.24316, 0.28047, 0.54179},
2351        {0.23674, 0.23784, 0.24142, 0.24847, 0.26264, 0.31810,
2352         0.29602, 0.29008, 0.29944, 0.33675, 0.59795},
2353        {0.29279, 0.29382, 0.29742, 0.30448, 0.31865, 0.37364,
2354         0.35215, 0.34629, 0.35555, 0.39286, 0.65411},
2355        {0.34875, 0.34987, 0.35346, 0.36054, 0.37472, 0.42956,
2356         0.40825, 0.40229, 0.41167, 0.44894, 0.71033},
2357        {0.40484, 0.40594, 0.40954, 0.41660, 0.43077, 0.48560,
2358         0.46433, 0.45840, 0.46772, 0.50500, 0.76632},
2359        {0.46090, 0.46201, 0.46560, 0.47267, 0.48684, 0.54167,
2360         0.52041, 0.51449, 0.52382, 0.56108, 0.82227},
2361        {0.51698, 0.51809, 0.52173, 0.52874, 0.54291, 0.59776,
2362         0.57646, 0.57052, 0.57986, 0.61717, 0.87836},
2363        {0.57306, 0.57418, 0.57782, 0.58513, 0.59899, 0.65380,
2364         0.63255, 0.62661, 0.63594, 0.67325, 0.93460},
2365        {0.62912, 0.63024, 0.63383, 0.64103, 0.65506, 0.70988,
2366         0.68484, 0.68267, 0.69202, 0.72878, 0.99046},
2367        {0.68521, 0.68633, 0.68990, 0.69699, 0.71115, 0.76595,
2368         0.74468, 0.73872, 0.74814, 0.78538, 1.04674},
2369        {0.74127, 0.74239, 0.74605, 0.75303, 0.77022, 0.82204,
2370         0.80078, 0.79484, 0.80416, 0.84147, 1.10261},
2371        {0.79736, 0.79846, 0.80206, 0.80947, 0.82330, 0.87813,
2372         0.85688, 0.85087, 0.86023, 0.89752, 1.15874},
2373        {0.85342, 0.85454, 0.85815, 0.86519, 0.87936, 0.93417,
2374         0.91293, 0.90428, 0.91631, 0.95360, 1.20760},
2375        {0.90949, 0.91061, 0.91423, 0.92128, 0.93544, 0.99026,
2376         0.96807, 0.96305, 0.97239, 1.00967, 1.27078},
2377        {0.96556, 0.96669, 0.97111, 0.97734, 0.99151, 1.04664,
2378         1.02508, 1.01879, 1.02846, 1.06167, 1.32695},
2379        {1.02167, 1.02279, 1.02656, 1.03341, 1.04759, 1.10242,
2380         1.08115, 1.07003, 1.08453, 1.12184, 1.38304},
2381        {1.07776, 1.07883, 1.08242, 1.08950, 1.10384, 1.16422,
2382         1.13725, 1.13133, 1.14061, 1.17793, 1.43910},
2383        {1.13379, 1.13492, 1.13864, 1.14567, 1.15973, 1.21455,
2384         1.19323, 1.18734, 1.19668, 1.23401, 1.49528},
2385        {1.18988, 1.19098, 1.19457, 1.20164, 1.21582, 1.27064,
2386         1.24937, 1.24044, 1.25275, 1.29004, 1.55137},
2387        {1.24592, 1.24706, 1.25087, 1.25773, 1.27188, 1.32670,
2388         1.30544, 1.29953, 1.30883, 1.34613, 1.60743},
2389        {1.30202, 1.30313, 1.30673, 1.31381, 1.32797, 1.38278,
2390         1.36151, 1.35167, 1.36490, 1.40221, 1.66306},
2391        {1.35809, 1.35921, 1.36282, 1.36986, 1.38403, 1.43888,
2392         1.41760, 1.41174, 1.42083, 1.45830, 1.71915},
2393        {1.41419, 1.41528, 1.41890, 1.42595, 1.44011, 1.49496,
2394         1.47368, 1.46769, 1.47706, 1.51436, 1.77523},
2395        {1.47131, 1.47141, 1.47494, 1.48850, 1.49620, 1.55137,
2396         1.52977, 1.51820, 1.53315, 1.57042, 1.83158},
2397        {1.52635, 1.52750, 1.53103, 1.53814, 1.55228, 1.60736,
2398         1.58503, 1.57986, 1.58920, 1.62649, 1.88767},
2399        {1.58418, 1.58355, 1.58711, 1.59526, 1.60833, 1.66316,
2400         1.63345, 1.63261, 1.64556, 1.68204, 1.94359},
2401        {1.64027, 1.63958, 1.64489, 1.65024, 1.66443, 1.71925,
2402         1.69794, 1.69201, 1.70143, 1.73865, 1.99968},
2403        {1.69450, 1.69566, 1.69940, 1.70697, 1.71841, 1.77819,
2404         1.75396, 1.74814, 1.75743, 1.79083, 2.05427},
2405        {1.75054, 1.75221, 1.75527, 1.76306, 1.77662, 1.83428,
2406         1.81006, 1.81173, 1.81345, 1.85076, 2.10289}};
2407
2408   // Garfield simulation at UD = -2100V; vd = 1.959cm/microsec, <driftfield> = 720V/cm
2409   Float_t time2100[ktimebin][kZbin] =
2410       {{0.09160, 0.09194, 0.09294, 0.09462, 0.09701, 0.10017,
2411         0.10424, 0.10974, 0.11823, 0.13988, 0.37762},
2412        {0.06585, 0.06619, 0.06731, 0.06916, 0.07174, 0.07509,
2413         0.07933, 0.08500, 0.09370, 0.11609, 0.35565},
2414        {0.03960, 0.04001, 0.04171, 0.04430, 0.04774, 0.05182,
2415         0.05668, 0.06294, 0.07291, 0.09896, 0.34676},
2416        {0.01109, 0.01280, 0.01716, 0.02279, 0.02876, 0.03474,
2417         0.04096, 0.04908, 0.06424, 0.10612, 0.36790},
2418        {0.01114, 0.01285, 0.01719, 0.02287, 0.02863, 0.03449,
2419         0.04073, 0.04964, 0.06759, 0.11446, 0.37162},
2420        {0.03922, 0.03977, 0.04146, 0.04386, 0.04711, 0.05109,
2421         0.05605, 0.06287, 0.07575, 0.10713, 0.36298},
2422        {0.06504, 0.06538, 0.06641, 0.06818, 0.07064, 0.07426,
2423         0.07852, 0.08483, 0.09581, 0.12363, 0.37424},
2424        {0.09103, 0.09129, 0.09186, 0.09291, 0.09476, 0.00000,
2425         0.10923, 0.11578, 0.12873, 0.16005, 0.41525},
2426        {0.12723, 0.12777, 0.12988, 0.13458, 0.14579, 0.20264,
2427         0.17421, 0.17078, 0.18132, 0.21708, 0.47699},
2428        {0.17508, 0.17601, 0.17897, 0.18487, 0.19698, 0.24881,
2429         0.22737, 0.22337, 0.23348, 0.27000, 0.53032},
2430        {0.22571, 0.22663, 0.22969, 0.23570, 0.24787, 0.29826,
2431         0.27871, 0.27462, 0.28471, 0.32122, 0.58166},
2432        {0.27664, 0.27759, 0.28067, 0.28669, 0.29891, 0.34898,
2433         0.32982, 0.32570, 0.33576, 0.37229, 0.63268},
2434        {0.32766, 0.32862, 0.33170, 0.33778, 0.34988, 0.39973,
2435         0.38088, 0.37675, 0.38680, 0.42333, 0.68159},
2436        {0.37872, 0.37966, 0.38275, 0.38875, 0.40093, 0.45073,
2437         0.43192, 0.42780, 0.43786, 0.47438, 0.73480},
2438        {0.42974, 0.43070, 0.43378, 0.43982, 0.45196, 0.50177,
2439         0.48297, 0.47884, 0.48889, 0.52544, 0.78581},
2440        {0.48081, 0.48175, 0.48482, 0.49084, 0.50302, 0.55290,
2441         0.53398, 0.52988, 0.53994, 0.57647, 0.83687},
2442        {0.53645, 0.53295, 0.53586, 0.54188, 0.55408, 0.60398,
2443         0.58504, 0.58092, 0.59100, 0.62768, 0.88773},
2444        {0.58345, 0.58409, 0.58690, 0.59292, 0.60510, 0.65562,
2445         0.63609, 0.63197, 0.64203, 0.67856, 0.93907},
2446        {0.63397, 0.63490, 0.63795, 0.64403, 0.65613, 0.70612,
2447         0.68714, 0.68301, 0.69294, 0.72955, 0.99000},
2448        {0.68496, 0.68592, 0.68899, 0.69504, 0.70733, 0.75708,
2449         0.73818, 0.73405, 0.74412, 0.78064, 1.04100},
2450        {0.73600, 0.73696, 0.74003, 0.74624, 0.75828, 0.80805,
2451         0.78904, 0.78512, 0.79517, 0.83152, 1.09205},
2452        {0.78709, 0.78801, 0.79108, 0.79709, 0.80931, 0.85906,
2453         0.84027, 0.83614, 0.84621, 0.88269, 1.14058},
2454        {0.83808, 0.83905, 0.84215, 0.84816, 0.86031, 0.91011,
2455         0.89139, 0.88718, 0.89725, 0.93377, 1.19413},
2456        {0.88916, 0.89010, 0.89320, 0.89920, 0.91136, 0.96117,
2457         0.94235, 0.93822, 0.94828, 0.98480, 1.24538},
2458        {0.94036, 0.94113, 0.94422, 0.95023, 0.96241, 1.01220,
2459         0.99310, 0.98927, 0.99933, 1.03585, 1.29629},
2460        {0.99139, 0.99220, 0.99525, 1.00127, 1.01344, 1.06324,
2461         1.04451, 1.04033, 1.04836, 1.08690, 1.34727},
2462        {1.04261, 1.04325, 1.04628, 1.05232, 1.06448, 1.12090,
2463         1.09546, 1.09136, 1.10142, 1.13795, 1.39831},
2464        {1.09331, 1.09429, 1.09742, 1.10336, 1.11557, 1.16547,
2465         1.14658, 1.13642, 1.15247, 1.18898, 1.44936},
2466        {1.14436, 1.14539, 1.14847, 1.15443, 1.16662, 1.21794,
2467         1.19763, 1.19329, 1.20349, 1.23956, 1.50043},
2468        {1.19533, 1.19651, 1.19943, 1.20548, 1.21666, 1.26753,
2469         1.24862, 1.24434, 1.25455, 1.29106, 1.55142},
2470        {1.24638, 1.24756, 1.25046, 1.25648, 1.26764, 1.31858,
2471         1.29967, 1.29538, 1.30499, 1.34211, 1.60250},
2472        {1.29747, 1.29847, 1.30175, 1.30753, 1.31869, 1.36969,
2473         1.35069, 1.34656, 1.35663, 1.39316, 1.64644},
2474        {1.35537, 1.34952, 1.35255, 1.35869, 1.36973, 1.41387,
2475         1.40173, 1.39761, 1.40768, 1.44396, 1.70238},
2476        {1.39956, 1.40056, 1.40380, 1.40961, 1.42178, 1.46492,
2477         1.45278, 1.45423, 1.45872, 1.49522, 1.75557},
2478        {1.45080, 1.45159, 1.45463, 1.46109, 1.47287, 1.52263,
2479         1.50382, 1.50050, 1.50977, 1.54502, 1.80670},
2480        {1.50165, 1.50264, 1.50570, 1.51214, 1.52233, 1.57370,
2481         1.55484, 1.55155, 1.56080, 1.59731, 1.85778},
2482        {1.55269, 1.55364, 1.55675, 1.56274, 1.57491, 1.62598,
2483         1.60590, 1.60259, 1.61185, 1.64836, 1.90883},
2484        {1.60368, 1.60469, 1.60779, 1.61373, 1.62596, 1.67738,
2485         1.65651, 1.65249, 1.66290, 1.69936, 1.95959}};
2486
2487   // Garfield simulation at UD = -2200V; vd = 2.134cm/microsec, <driftfield> = 753V/cm
2488   Float_t time2200[ktimebin][kZbin] =
2489       {{0.09162, 0.09194, 0.09292, 0.09460, 0.09702, 0.10014,
2490         0.10421, 0.10971, 0.11820, 0.13990, 0.37745},
2491        {0.06581, 0.06618, 0.06730, 0.06915, 0.07173, 0.07507,
2492         0.07931, 0.08497, 0.09368, 0.11609, 0.35560},
2493        {0.03947, 0.04001, 0.04167, 0.04429, 0.04772, 0.05183,
2494         0.05667, 0.06293, 0.07292, 0.09900, 0.34673},
2495        {0.01111, 0.01280, 0.01716, 0.02279, 0.02876, 0.03473,
2496         0.04091, 0.04907, 0.06426, 0.10620, 0.36766},
2497        {0.01113, 0.01285, 0.01719, 0.02276, 0.02863, 0.03452,
2498         0.04076, 0.04960, 0.06745, 0.11419, 0.37139},
2499        {0.03923, 0.03978, 0.04137, 0.04387, 0.04713, 0.05110,
2500         0.05605, 0.06284, 0.07551, 0.10677, 0.36210},
2501        {0.06505, 0.06540, 0.06644, 0.06820, 0.07069, 0.07401,
2502         0.07852, 0.08479, 0.09565, 0.12325, 0.37313},
2503        {0.09107, 0.09127, 0.09181, 0.09291, 0.09474, 0.00000,
2504         0.10883, 0.11528, 0.12789, 0.15865, 0.41313},
2505        {0.12559, 0.12622, 0.12800, 0.13206, 0.14166, 0.19331,
2506         0.16832, 0.16632, 0.17724, 0.21218, 0.47098},
2507        {0.16992, 0.17070, 0.17325, 0.17831, 0.18871, 0.23557,
2508         0.21690, 0.21451, 0.22514, 0.26082, 0.52034},
2509        {0.21640, 0.21722, 0.21987, 0.22500, 0.23540, 0.28097,
2510         0.26399, 0.26154, 0.27214, 0.30784, 0.56734},
2511        {0.26318, 0.26400, 0.26679, 0.27181, 0.28220, 0.32739,
2512         0.31090, 0.30842, 0.31902, 0.35474, 0.61415},
2513        {0.31001, 0.31085, 0.31348, 0.31866, 0.32903, 0.37412,
2514         0.35777, 0.35546, 0.36588, 0.40159, 0.66103},
2515        {0.35687, 0.35769, 0.36033, 0.36556, 0.37588, 0.42094,
2516         0.40523, 0.40214, 0.41273, 0.44841, 0.70785},
2517        {0.40372, 0.40457, 0.40723, 0.41234, 0.42273, 0.46778,
2518         0.45148, 0.44903, 0.45961, 0.49526, 0.75486},
2519        {0.45062, 0.45139, 0.45404, 0.45966, 0.46958, 0.51470,
2520         0.49833, 0.49584, 0.50644, 0.54211, 0.80160},
2521        {0.49742, 0.49825, 0.50088, 0.50605, 0.51644, 0.56148,
2522         0.54518, 0.54270, 0.55330, 0.58897, 0.84854},
2523        {0.54427, 0.54510, 0.54774, 0.55290, 0.56329, 0.60846,
2524         0.59203, 0.58955, 0.60014, 0.63578, 0.89528},
2525        {0.59119, 0.59199, 0.59471, 0.59975, 0.61014, 0.65533,
2526         0.63889, 0.63636, 0.64699, 0.68269, 0.94197},
2527        {0.63866, 0.63880, 0.64145, 0.64664, 0.65701, 0.70639,
2528         0.68574, 0.68325, 0.69385, 0.72949, 0.98900},
2529        {0.68483, 0.68566, 0.68831, 0.69347, 0.70386, 0.74890,
2530         0.73260, 0.73010, 0.74069, 0.77638, 1.03320},
2531        {0.73168, 0.73255, 0.73515, 0.74031, 0.75072, 0.79576,
2532         0.77117, 0.77501, 0.78755, 0.82318, 1.08006},
2533        {0.77854, 0.78310, 0.78200, 0.79525, 0.79756, 0.84281,
2534         0.81803, 0.82393, 0.83441, 0.87008, 1.12692},
2535        {0.82541, 0.82642, 0.82916, 0.83528, 0.84442, 0.89086,
2536         0.87315, 0.87079, 0.88125, 0.91694, 1.17648},
2537        {0.87226, 0.87308, 0.87602, 0.88086, 0.89128, 0.93772,
2538         0.92001, 0.91751, 0.92811, 0.95587, 1.22328},
2539        {0.91921, 0.91994, 0.92256, 0.92772, 0.94713, 0.98566,
2540         0.96690, 0.96436, 0.97495, 1.01064, 1.26882},
2541        {0.96790, 0.96679, 0.96941, 0.97463, 0.99399, 1.03001,
2542         1.01376, 1.01112, 1.02181, 1.05749, 1.31568},
2543        {1.01278, 1.01390, 1.01674, 1.02147, 1.03182, 1.07820,
2544         1.06056, 1.05798, 1.06867, 1.10433, 1.36390},
2545        {1.05964, 1.06076, 1.06331, 1.06833, 1.07870, 1.13297,
2546         1.10742, 1.10520, 1.11543, 1.15120, 1.41069},
2547        {1.10664, 1.10762, 1.10997, 1.11519, 1.12556, 1.17531,
2548         1.15427, 1.14620, 1.16229, 1.19805, 1.45758},
2549        {1.15352, 1.15421, 1.15683, 1.16218, 1.17242, 1.21910,
2550         1.20035, 1.19863, 1.20579, 1.24473, 1.50412},
2551        {1.20019, 1.20115, 1.20369, 1.20892, 1.21928, 1.26913,
2552         1.24721, 1.24549, 1.25605, 1.29159, 1.54920},
2553        {1.24707, 1.24846, 1.25052, 1.25602, 1.26608, 1.31558,
2554         1.29448, 1.29232, 1.30293, 1.33675, 1.59798},
2555        {1.29391, 1.29475, 1.29738, 1.30255, 1.31294, 1.36244,
2556         1.34167, 1.33918, 1.34979, 1.38229, 1.64496},
2557        {1.34078, 1.34304, 1.34424, 1.35565, 1.35980, 1.40930,
2558         1.38853, 1.38229, 1.39664, 1.42863, 1.69162},
2559        {1.38762, 1.38847, 1.39110, 1.39627, 1.40666, 1.45183,
2560         1.43539, 1.43289, 1.44348, 1.47549, 1.73876},
2561        {1.43524, 1.43533, 1.43796, 1.44310, 1.45371, 1.49305,
2562         1.48224, 1.47941, 1.49034, 1.52601, 1.78552},
2563        {1.48122, 1.48219, 1.48482, 1.48991, 1.50030, 1.53991,
2564         1.52898, 1.52653, 1.53653, 1.57282, 1.82386}};
2565
2566   if (fTimeStruct1) {
2567     delete [] fTimeStruct1;
2568   }
2569   fTimeStruct1 = new Float_t[ktimebin*kZbin];
2570
2571   if (fTimeStruct2) {
2572     delete [] fTimeStruct2;
2573   }
2574   fTimeStruct2 = new Float_t[ktimebin*kZbin];
2575
2576   for (Int_t ctrt = 0; ctrt < ktimebin; ctrt++) {
2577     for (Int_t ctrz = 0; ctrz <    kZbin; ctrz++) {
2578       if      (vdrift > fVDsmp[6]) {
2579         fTimeStruct1[ctrt+ctrz*ktimebin] = time2100[ctrt][ctrz];
2580         fTimeStruct2[ctrt+ctrz*ktimebin] = time2200[ctrt][ctrz];
2581         fVDlo = fVDsmp[6];
2582         fVDhi = fVDsmp[7];
2583       }
2584       else if (vdrift > fVDsmp[5]) {
2585         fTimeStruct1[ctrt+ctrz*ktimebin] = time2000[ctrt][ctrz];
2586         fTimeStruct2[ctrt+ctrz*ktimebin] = time2100[ctrt][ctrz];
2587         fVDlo = fVDsmp[5];
2588         fVDhi = fVDsmp[6];
2589       }
2590       else if (vdrift > fVDsmp[4]) {
2591         fTimeStruct1[ctrt+ctrz*ktimebin] = time1900[ctrt][ctrz];
2592         fTimeStruct2[ctrt+ctrz*ktimebin] = time2000[ctrt][ctrz];
2593         fVDlo = fVDsmp[4];
2594         fVDhi = fVDsmp[5];
2595       }
2596       else if (vdrift > fVDsmp[3]) {
2597         fTimeStruct1[ctrt+ctrz*ktimebin] = time1800[ctrt][ctrz];
2598         fTimeStruct2[ctrt+ctrz*ktimebin] = time1900[ctrt][ctrz];
2599         fVDlo = fVDsmp[3];
2600         fVDhi = fVDsmp[4];
2601       }
2602       else if (vdrift > fVDsmp[2]) {
2603         fTimeStruct1[ctrt+ctrz*ktimebin] = time1700[ctrt][ctrz];
2604         fTimeStruct2[ctrt+ctrz*ktimebin] = time1800[ctrt][ctrz];
2605         fVDlo = fVDsmp[2];
2606         fVDhi = fVDsmp[3];
2607       }
2608       else if (vdrift > fVDsmp[1]) {
2609         fTimeStruct1[ctrt+ctrz*ktimebin] = time1600[ctrt][ctrz];
2610         fTimeStruct2[ctrt+ctrz*ktimebin] = time1700[ctrt][ctrz];
2611         fVDlo = fVDsmp[1];
2612         fVDhi = fVDsmp[2];
2613       }
2614       else if (vdrift > (fVDsmp[0] - 1.0e-5)) {
2615         fTimeStruct1[ctrt+ctrz*ktimebin] = time1500[ctrt][ctrz];
2616         fTimeStruct2[ctrt+ctrz*ktimebin] = time1600[ctrt][ctrz];
2617         fVDlo = fVDsmp[0];
2618         fVDhi = fVDsmp[1];
2619       }
2620     }
2621   }
2622
2623 }
2624
2625 //_____________________________________________________________________________
2626 void AliTRDdigitizer::RecalcDiffusion(Float_t vdrift)
2627 {
2628   //
2629   // Recalculates the diffusion parameters
2630   //
2631   // The B=0 case is not really included here.
2632   // Should be revisited!
2633   //
2634
2635   if (vdrift == fDiffLastVdrift) {
2636     return;
2637   }
2638
2639   AliTRDSimParam    *simParam    = AliTRDSimParam::Instance();
2640   if (!simParam) {
2641     AliFatal("Could not get simulation parameters");
2642     return;
2643   }
2644
2645   AliTRDCommonParam *commonParam = AliTRDCommonParam::Instance();
2646   if (!commonParam) {
2647     AliFatal("Could not get common parameters");
2648     return;
2649   }
2650
2651   AliTRDcalibDB     *calibration = AliTRDcalibDB::Instance();
2652   if (!calibration) {
2653     AliFatal("Could not get calibration object");
2654     return;
2655   }
2656
2657   // The magnetic field strength
2658   Double_t x[3] = { 0.0, 0.0, 0.0 };
2659   Double_t b[3];
2660   gAlice->Field(x,b);         // b[] is in kilo Gauss
2661   Float_t field = b[2] * 0.1; // Tesla
2662
2663   fDiffLastVdrift = vdrift;
2664
2665   // DiffusionL
2666   const Int_t kNbL = 5;
2667   Float_t p0L[kNbL] = {  0.007440,  0.007493,  0.007513,  0.007672,  0.007831 };
2668   Float_t p1L[kNbL] = {  0.019252,  0.018912,  0.018636,  0.018012,  0.017343 };
2669   Float_t p2L[kNbL] = { -0.005042, -0.004926, -0.004867, -0.004650, -0.004424 };
2670   Float_t p3L[kNbL] = {  0.000195,  0.000189,  0.000195,  0.000182,  0.000169 };
2671
2672   Int_t ibL = ((Int_t) (10 * (field - 0.15)));
2673   ibL       = TMath::Max(   0,ibL);
2674   ibL       = TMath::Min(kNbL,ibL);
2675
2676   fDiffusionL = p0L[ibL]
2677               + p1L[ibL] * vdrift
2678               + p2L[ibL] * vdrift*vdrift
2679               + p3L[ibL] * vdrift*vdrift*vdrift;
2680
2681   // DiffusionT
2682   const Int_t kNbT = 5;
2683   Float_t p0T[kNbT] = {  0.009550,  0.009599,  0.009674,  0.009757,  0.009850 };
2684   Float_t p1T[kNbT] = {  0.006667,  0.006539,  0.006359,  0.006153,  0.005925 };
2685   Float_t p2T[kNbT] = { -0.000853, -0.000798, -0.000721, -0.000635, -0.000541 };
2686   Float_t p3T[kNbT] = {  0.000131,  0.000122,  0.000111,  0.000098,  0.000085 };
2687
2688   Int_t ibT= ((Int_t) (10 * (field - 0.15)));
2689   ibT      = TMath::Max(   0,ibT);
2690   ibT      = TMath::Min(kNbT,ibT);
2691
2692   fDiffusionT = p0T[ibT]
2693               + p1T[ibT] * vdrift
2694               + p2T[ibT] * vdrift*vdrift
2695               + p3T[ibT] * vdrift*vdrift*vdrift;
2696
2697   // OmegaTau
2698   fOmegaTau = calibration->GetOmegaTau(vdrift,field);
2699
2700   // Lorentzfactor
2701   if (commonParam->ExBOn()) {
2702     fLorentzFactor = 1.0 / (1.0 + fOmegaTau*fOmegaTau);
2703   }
2704   else {
2705     fLorentzFactor = 1.0;
2706   }
2707
2708 }
2709
2710 //_____________________________________________________________________________
2711 Float_t AliTRDdigitizer::GetDiffusionL(Float_t vdrift)
2712 {
2713   //
2714   // Returns the longitudinal diffusion coefficient for a given drift
2715   // velocity <vd> and a B-field <b> for Xe/CO2 (15%).
2716   // The values are according to a GARFIELD simulation.
2717   //
2718
2719   RecalcDiffusion(vdrift);
2720
2721   return fDiffusionL;
2722
2723 }
2724
2725 //_____________________________________________________________________________
2726 Float_t AliTRDdigitizer::GetDiffusionT(Float_t vdrift)
2727 {
2728   //
2729   // Returns the transverse diffusion coefficient for a given drift
2730   // velocity <vd> and a B-field <b> for Xe/CO2 (15%).
2731   // The values are according to a GARFIELD simulation.
2732   //
2733
2734   RecalcDiffusion(vdrift);
2735
2736   return fDiffusionT;
2737
2738 }
2739
2740 //_____________________________________________________________________________
2741 Int_t AliTRDdigitizer::Diffusion(Float_t vdrift, Double_t absdriftlength
2742                                , Double_t &lRow, Double_t &lCol, Double_t &lTime)
2743 {
2744   //
2745   // Applies the diffusion smearing to the position of a single electron.
2746   // Depends on absolute drift length.
2747   //
2748
2749   RecalcDiffusion(vdrift);
2750
2751   Float_t driftSqrt = TMath::Sqrt(absdriftlength);
2752   Float_t sigmaT    = driftSqrt * fDiffusionT;
2753   Float_t sigmaL    = driftSqrt * fDiffusionL;
2754   lRow  = gRandom->Gaus(lRow ,sigmaT);
2755   lCol  = gRandom->Gaus(lCol ,sigmaT * GetLorentzFactor(vdrift));
2756   lTime = gRandom->Gaus(lTime,sigmaL * GetLorentzFactor(vdrift));
2757
2758   return 1;
2759
2760 }
2761
2762 //_____________________________________________________________________________
2763 Float_t AliTRDdigitizer::GetLorentzFactor(Float_t vd)
2764 {
2765   //
2766   // Returns the recalculated Lorentz factor
2767   //
2768
2769   RecalcDiffusion(vd);
2770
2771   return fLorentzFactor;
2772
2773 }
2774
2775 //_____________________________________________________________________________
2776 Int_t AliTRDdigitizer::ExB(Float_t vdrift, Double_t driftlength, Double_t &lCol)
2777 {
2778   //
2779   // Applies E x B effects to the position of a single electron.
2780   // Depends on signed drift length.
2781   //
2782
2783   RecalcDiffusion(vdrift);
2784
2785   lCol = lCol + fOmegaTau * driftlength;
2786
2787   return 1;
2788
2789 }