TRD/AliTRDdigitizer.cxx

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