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