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