TPC/AliTPC.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 //  Time Projection Chamber                                                  //
  21 //  This class contains the basic functions for the Time Projection Chamber  //
  22 //  detector. Functions specific to one particular geometry are              //
  23 //  contained in the derived classes                                         //
  24 //                                                                           //
  25 //Begin_Html
  26 /*
  27 <img src="picts/AliTPCClass.gif">
  28 */
  29 //End_Html
  30 //                                                                           //
  31 //                                                                          //
  32 ///////////////////////////////////////////////////////////////////////////////
  33
  34 //
  35
  36 #include <Riostream.h>
  37 #include <stdlib.h>
  38
  39 #include <TF2.h>
  40 #include <TFile.h>
  41 #include <TGeoGlobalMagField.h>
  42 #include <TInterpreter.h>
  43 #include <TMath.h>
  44 #include <TMatrixF.h>
  45 #include <TObjectTable.h>
  46 #include <TParticle.h>
  47 #include <TROOT.h>
  48 #include <TRandom.h>
  49 #include <TStopwatch.h>
  50 #include <TString.h>
  51 #include <TSystem.h>
  52 #include <TTree.h>
  53 #include <TVector.h>
  54 #include <TVirtualMC.h>
  55 #include <TParameter.h>
  56
  57 #include "AliDigits.h"
  58 #include "AliMagF.h"
  59 #include "AliRun.h"
  60 #include "AliRunLoader.h"
  61 #include "AliSimDigits.h"
  62 #include "AliTPC.h"
  63 #include "AliTPC.h"
  64 #include "AliTPCDigitsArray.h"
  65 #include "AliTPCLoader.h"
  66 #include "AliTPCPRF2D.h"
  67 #include "AliTPCParamSR.h"
  68 #include "AliTPCRF1D.h"
  69 #include "AliTPCTrackHitsV2.h"
  70 #include "AliTrackReference.h"
  71 #include "AliMC.h"
  72 #include "AliStack.h"
  73 #include "AliTPCDigitizer.h"
  74 #include "AliTPCBuffer.h"
  75 #include "AliTPCDDLRawData.h"
  76 #include "AliLog.h"
  77 #include "AliTPCcalibDB.h"
  78 #include "AliTPCCalPad.h"
  79 #include "AliTPCCalROC.h"
  80 #include "AliTPCExB.h"
  81 #include "AliRawReader.h"
  82 #include "AliTPCRawStreamV3.h"
  83 #include "TTreeStream.h"
  84
  85 ClassImp(AliTPC)
  86 //_____________________________________________________________________________
  87   AliTPC::AliTPC():AliDetector(),
  88                    fDefaults(0),
  89                    fSens(0),
  90                    fNsectors(0),
  91                    fDigitsArray(0),
  92                    fTPCParam(0),
  93                    fTrackHits(0),
  94                    fHitType(0),
  95                    fDigitsSwitch(0),
  96                    fSide(0),
  97                    fPrimaryIonisation(0),
  98                    fNoiseDepth(0),
  99                    fNoiseTable(0),
 100                    fCurrentNoise(0),
 101                    fActiveSectors(0),
 102                    fGainFactor(1.),
 103                    fDebugStreamer(0),
 104                    fLHCclockPhaseSw(0),
 105                    fIsGEM(0)
 106
 107 {
 108   //
 109   // Default constructor
 110   //
 111   fIshunt   = 0;
 112   for(Int_t i=0;i<4;i++) fCurrentIndex[i]=0;
 113
 114   //  fTrackHitsOld = 0;
 115 #if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
 116   fHitType = 4; // ROOT containers
 117 #else
 118   fHitType = 2; //default CONTAINERS - based on ROOT structure
 119 #endif
 120 }
 121
 122 //_____________________________________________________________________________
 123 AliTPC::AliTPC(const char *name, const char *title)
 124   : AliDetector(name,title),
 125                    fDefaults(0),
 126                    fSens(0),
 127                    fNsectors(0),
 128                    fDigitsArray(0),
 129                    fTPCParam(0),
 130                    fTrackHits(0),
 131                    fHitType(0),
 132                    fDigitsSwitch(0),
 133                    fSide(0),
 134                    fPrimaryIonisation(0),
 135                    fNoiseDepth(0),
 136                    fNoiseTable(0),
 137                    fCurrentNoise(0),
 138                    fActiveSectors(0),
 139                    fGainFactor(1.),
 140                    fDebugStreamer(0),
 141     fLHCclockPhaseSw(0),
 142     fIsGEM(0)
 143
 144 {
 145   //
 146   // Standard constructor
 147   //
 148
 149   //
 150   // Initialise arrays of hits and digits
 151   fHits     = new TClonesArray("AliTPChit",  176);
 152   gAlice->GetMCApp()->AddHitList(fHits);
 153   //
 154   fTrackHits = new AliTPCTrackHitsV2;
 155   fTrackHits->SetHitPrecision(0.002);
 156   fTrackHits->SetStepPrecision(0.003);
 157   fTrackHits->SetMaxDistance(100);
 158
 159   //fTrackHitsOld = new AliTPCTrackHits;  //MI - 13.09.2000
 160   //fTrackHitsOld->SetHitPrecision(0.002);
 161   //fTrackHitsOld->SetStepPrecision(0.003);
 162   //fTrackHitsOld->SetMaxDistance(100);
 163
 164
 165 #if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
 166   fHitType = 4; // ROOT containers
 167 #else
 168   fHitType = 2;
 169 #endif
 170
 171   for(Int_t i=0;i<4;i++) fCurrentIndex[i]=0;
 172
 173   //
 174   fIshunt     =  0;
 175   //
 176   // Initialise color attributes
 177   //PH SetMarkerColor(kYellow);
 178
 179   //
 180   //  Set TPC parameters
 181   //
 182
 183
 184   if (!strcmp(title,"Default")) {
 185     //fTPCParam = new AliTPCParamSR;
 186     fTPCParam = AliTPCcalibDB::Instance()->GetParameters();
 187   } else {
 188     AliWarning("In Config.C you must set non-default parameters.");
 189     fTPCParam=0;
 190   }
 191
 192 }
 193 void AliTPC::CreateDebugStremer(){
 194   //
 195   // Create Debug streamer to check simulation
 196   //
 197   fDebugStreamer = new TTreeSRedirector("TPCSimdebug.root");
 198 }
 199 //_____________________________________________________________________________
 200 AliTPC::~AliTPC()
 201 {
 202   //
 203   // TPC destructor
 204   //
 205
 206   fIshunt   = 0;
 207   delete fHits;
 208   delete fDigits;
 209   //delete fTPCParam;
 210   delete fTrackHits; //MI 15.09.2000
 211   //  delete fTrackHitsOld; //MI 10.12.2001
 212
 213   fDigitsArray = 0x0;
 214   delete [] fNoiseTable;
 215   delete [] fActiveSectors;
 216   if (fDebugStreamer) delete fDebugStreamer;
 217 }
 218
 219 //_____________________________________________________________________________
 220 void AliTPC::AddHit(Int_t track, Int_t *vol, Float_t *hits)
 221 {
 222   //
 223   // Add a hit to the list
 224   //
 225   if (fHitType&1){
 226     TClonesArray &lhits = *fHits;
 227     new(lhits[fNhits++]) AliTPChit(fIshunt,track,vol,hits);
 228   }
 229   if (fHitType>1)
 230     AddHit2(track,vol,hits);
 231 }
 232
 233 //_____________________________________________________________________________
 234 void AliTPC::CreateMaterials()
 235 {
 236   //-----------------------------------------------
 237   // Create Materials for for TPC simulations
 238   //-----------------------------------------------
 239
 240   //-----------------------------------------------------------------
 241   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
 242   //-----------------------------------------------------------------
 243
 244    Int_t iSXFLD=((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Integ();
 245   Float_t sXMGMX=((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Max();
 246
 247   Float_t amat[5]; // atomic numbers
 248   Float_t zmat[5]; // z
 249   Float_t wmat[5]; // proportions
 250
 251   Float_t density;
 252
 253
 254
 255   //***************** Gases *************************
 256
 257
 258   //--------------------------------------------------------------
 259   // gases - air and CO2
 260   //--------------------------------------------------------------
 261
 262   // CO2
 263
 264   amat[0]=12.011;
 265   amat[1]=15.9994;
 266
 267   zmat[0]=6.;
 268   zmat[1]=8.;
 269
 270   wmat[0]=0.2729;
 271   wmat[1]=0.7271;
 272
 273   density=0.001754609;
 274
 275
 276   AliMixture(10,"CO2",amat,zmat,density,2,wmat);
 277   //
 278   // Air
 279   //
 280   amat[0]=15.9994;
 281   amat[1]=14.007;
 282   //
 283   zmat[0]=8.;
 284   zmat[1]=7.;
 285   //
 286   wmat[0]=0.233;
 287   wmat[1]=0.767;
 288   //
 289   density=0.001205;
 290
 291   AliMixture(11,"Air",amat,zmat,density,2,wmat);
 292
 293   //----------------------------------------------------------------
 294   // drift gases
 295   //----------------------------------------------------------------
 296
 297   //
 298   // Drift gases 1 - nonsensitive, 2 - sensitive
 299   // Ne-CO2 (90-10) (volume) values at 20deg and 1 atm.
 300   // rho(Ne) = 0.839 g/cm^3, rho(CO2) = 1.842 g/cm^3
 301
 302
 303   amat[0]= 20.18;
 304   amat[1]=12.011;
 305   amat[2]=15.9994;
 306   // amat[3]=14.007;
 307
 308   zmat[0]= 10.;
 309   zmat[1]=6.;
 310   zmat[2]=8.;
 311   // zmat[3]=7.;
 312
 313   //wmat[0]=0.756992632;
 314   wmat[0]=0.8038965;
 315   //wmat[1]=0.056235789;
 316   wmat[1]= 0.053519;
 317   //wmat[2]=0.128469474;
 318   wmat[2]= 0.1425743;
 319   // wmat[3]=0.058395789;
 320
 321   density=0.0009393;
 322
 323   AliMixture(12,"Ne-CO2-1",amat,zmat,density,3,wmat);
 324   AliMixture(13,"Ne-CO2-2",amat,zmat,density,3,wmat);
 325   AliMixture(35,"Ne-CO2-3",amat,zmat,density,3,wmat);
 326   //----------------------------------------------------------------------
 327   //               solid materials
 328   //----------------------------------------------------------------------
 329
 330
 331   // Kevlar C14H22O2N2
 332
 333   amat[0] = 12.011;
 334   amat[1] = 1.;
 335   amat[2] = 15.999;
 336   amat[3] = 14.006;
 337
 338   zmat[0] = 6.;
 339   zmat[1] = 1.;
 340   zmat[2] = 8.;
 341   zmat[3] = 7.;
 342
 343   wmat[0] = 14.;
 344   wmat[1] = 22.;
 345   wmat[2] = 2.;
 346   wmat[3] = 2.;
 347
 348   density = 1.45;
 349
 350   AliMixture(14,"Kevlar",amat,zmat,density,-4,wmat);
 351
 352   // NOMEX
 353
 354   amat[0] = 12.011;
 355   amat[1] = 1.;
 356   amat[2] = 15.999;
 357   amat[3] = 14.006;
 358
 359   zmat[0] = 6.;
 360   zmat[1] = 1.;
 361   zmat[2] = 8.;
 362   zmat[3] = 7.;
 363
 364   wmat[0] = 14.;
 365   wmat[1] = 22.;
 366   wmat[2] = 2.;
 367   wmat[3] = 2.;
 368
 369   density = 0.029;
 370
 371   AliMixture(15,"NOMEX",amat,zmat,density,-4,wmat);
 372
 373   // Makrolon C16H18O3
 374
 375   amat[0] = 12.011;
 376   amat[1] = 1.;
 377   amat[2] = 15.999;
 378
 379   zmat[0] = 6.;
 380   zmat[1] = 1.;
 381   zmat[2] = 8.;
 382
 383   wmat[0] = 16.;
 384   wmat[1] = 18.;
 385   wmat[2] = 3.;
 386
 387   density = 1.2;
 388
 389   AliMixture(16,"Makrolon",amat,zmat,density,-3,wmat);
 390
 391   // Tedlar C2H3F
 392
 393   amat[0] = 12.011;
 394   amat[1] = 1.;
 395   amat[2] = 18.998;
 396
 397   zmat[0] = 6.;
 398   zmat[1] = 1.;
 399   zmat[2] = 9.;
 400
 401   wmat[0] = 2.;
 402   wmat[1] = 3.;
 403   wmat[2] = 1.;
 404
 405   density = 1.71;
 406
 407   AliMixture(17, "Tedlar",amat,zmat,density,-3,wmat);
 408
 409   // Mylar C5H4O2
 410
 411   amat[0]=12.011;
 412   amat[1]=1.;
 413   amat[2]=15.9994;
 414
 415   zmat[0]=6.;
 416   zmat[1]=1.;
 417   zmat[2]=8.;
 418
 419   wmat[0]=5.;
 420   wmat[1]=4.;
 421   wmat[2]=2.;
 422
 423   density = 1.39;
 424
 425   AliMixture(18, "Mylar",amat,zmat,density,-3,wmat);
 426   // material for "prepregs"
 427   // Epoxy - C14 H20 O3
 428   // Quartz SiO2
 429   // Carbon C
 430   // prepreg1 60% C-fiber, 40% epoxy (vol)
 431   amat[0]=12.011;
 432   amat[1]=1.;
 433   amat[2]=15.994;
 434
 435   zmat[0]=6.;
 436   zmat[1]=1.;
 437   zmat[2]=8.;
 438
 439   wmat[0]=0.923;
 440   wmat[1]=0.023;
 441   wmat[2]=0.054;
 442
 443   density=1.859;
 444
 445   AliMixture(19, "Prepreg1",amat,zmat,density,3,wmat);
 446
 447   //prepreg2 60% glass-fiber, 40% epoxy
 448
 449   amat[0]=12.01;
 450   amat[1]=1.;
 451   amat[2]=15.994;
 452   amat[3]=28.086;
 453
 454   zmat[0]=6.;
 455   zmat[1]=1.;
 456   zmat[2]=8.;
 457   zmat[3]=14.;
 458
 459   wmat[0]=0.194;
 460   wmat[1]=0.023;
 461   wmat[2]=0.443;
 462   wmat[3]=0.34;
 463
 464   density=1.82;
 465
 466   AliMixture(20, "Prepreg2",amat,zmat,density,4,wmat);
 467
 468   //prepreg3 50% glass-fiber, 50% epoxy
 469
 470   amat[0]=12.01;
 471   amat[1]=1.;
 472   amat[2]=15.994;
 473   amat[3]=28.086;
 474
 475   zmat[0]=6.;
 476   zmat[1]=1.;
 477   zmat[2]=8.;
 478   zmat[3]=14.;
 479
 480   wmat[0]=0.257;
 481   wmat[1]=0.03;
 482   wmat[2]=0.412;
 483   wmat[3]=0.3;
 484
 485   density=1.725;
 486
 487   AliMixture(21, "Prepreg3",amat,zmat,density,4,wmat);
 488
 489   // G10 60% SiO2 40% epoxy
 490
 491   amat[0]=12.01;
 492   amat[1]=1.;
 493   amat[2]=15.994;
 494   amat[3]=28.086;
 495
 496   zmat[0]=6.;
 497   zmat[1]=1.;
 498   zmat[2]=8.;
 499   zmat[3]=14.;
 500
 501   wmat[0]=0.194;
 502   wmat[1]=0.023;
 503   wmat[2]=0.443;
 504   wmat[3]=0.340;
 505
 506   density=1.7;
 507
 508   AliMixture(22, "G10",amat,zmat,density,4,wmat);
 509
 510   // Al
 511
 512   amat[0] = 26.98;
 513   zmat[0] = 13.;
 514
 515   density = 2.7;
 516
 517   AliMaterial(23,"Al",amat[0],zmat[0],density,999.,999.);
 518
 519   // Si (for electronics
 520
 521   amat[0] = 28.086;
 522   zmat[0] = 14.;
 523
 524   density = 2.33;
 525
 526   AliMaterial(24,"Si",amat[0],zmat[0],density,999.,999.);
 527
 528   // Cu
 529
 530   amat[0] = 63.546;
 531   zmat[0] = 29.;
 532
 533   density = 8.96;
 534
 535   AliMaterial(25,"Cu",amat[0],zmat[0],density,999.,999.);
 536
 537   // brass
 538
 539   amat[0] = 63.546;
 540   zmat[0] = 29.;
 541   //
 542   amat[1]= 65.409;
 543   zmat[1]= 30.;
 544   //
 545   wmat[0]= 0.6;
 546   wmat[1]= 0.4;
 547
 548   //
 549   density = 8.23;
 550
 551
 552   //
 553   AliMixture(33,"Brass",amat,zmat,density,2,wmat);
 554
 555   // Epoxy - C14 H20 O3
 556
 557   amat[0]=12.011;
 558   amat[1]=1.;
 559   amat[2]=15.9994;
 560
 561   zmat[0]=6.;
 562   zmat[1]=1.;
 563   zmat[2]=8.;
 564
 565   wmat[0]=14.;
 566   wmat[1]=20.;
 567   wmat[2]=3.;
 568
 569   density=1.25;
 570
 571   AliMixture(26,"Epoxy",amat,zmat,density,-3,wmat);
 572   //
 573   // epoxy film - 90% epoxy, 10% glass fiber
 574   //
 575   amat[0]=12.01;
 576   amat[1]=1.;
 577   amat[2]=15.994;
 578   amat[3]=28.086;
 579
 580   zmat[0]=6.;
 581   zmat[1]=1.;
 582   zmat[2]=8.;
 583   zmat[3]=14.;
 584
 585   wmat[0]=0.596;
 586   wmat[1]=0.071;
 587   wmat[2]=0.257;
 588   wmat[3]=0.076;
 589
 590
 591   density=1.345;
 592
 593   AliMixture(34, "Epoxy-film",amat,zmat,density,4,wmat);
 594
 595   // Plexiglas  C5H8O2
 596
 597   amat[0]=12.011;
 598   amat[1]=1.;
 599   amat[2]=15.9994;
 600
 601   zmat[0]=6.;
 602   zmat[1]=1.;
 603   zmat[2]=8.;
 604
 605   wmat[0]=5.;
 606   wmat[1]=8.;
 607   wmat[2]=2.;
 608
 609   density=1.18;
 610
 611   AliMixture(27,"Plexiglas",amat,zmat,density,-3,wmat);
 612
 613   // Carbon
 614
 615   amat[0]=12.011;
 616   zmat[0]=6.;
 617   density= 2.265;
 618
 619   AliMaterial(28,"C",amat[0],zmat[0],density,999.,999.);
 620
 621   // Fe (steel for the inner heat screen)
 622
 623   amat[0]=55.845;
 624
 625   zmat[0]=26.;
 626
 627   density=7.87;
 628
 629   AliMaterial(29,"Fe",amat[0],zmat[0],density,999.,999.);
 630   //
 631   // Peek - (C6H4-O-OC6H4-O-C6H4-CO)n
 632   amat[0]=12.011;
 633   amat[1]=1.;
 634   amat[2]=15.9994;
 635
 636   zmat[0]=6.;
 637   zmat[1]=1.;
 638   zmat[2]=8.;
 639
 640   wmat[0]=19.;
 641   wmat[1]=12.;
 642   wmat[2]=3.;
 643   //
 644   density=1.3;
 645   //
 646   AliMixture(30,"Peek",amat,zmat,density,-3,wmat);
 647   //
 648   //  Ceramics - Al2O3
 649   //
 650   amat[0] = 26.98;
 651   amat[1]= 15.9994;
 652
 653   zmat[0] = 13.;
 654   zmat[1]=8.;
 655
 656   wmat[0]=2.;
 657   wmat[1]=3.;
 658
 659   density = 3.97;
 660
 661   AliMixture(31,"Alumina",amat,zmat,density,-2,wmat);
 662
 663   //
 664   // liquids
 665   //
 666
 667   // water
 668
 669   amat[0]=1.;
 670   amat[1]=15.9994;
 671
 672   zmat[0]=1.;
 673   zmat[1]=8.;
 674
 675   wmat[0]=2.;
 676   wmat[1]=1.;
 677
 678   density=1.;
 679
 680   AliMixture(32,"Water",amat,zmat,density,-2,wmat);
 681
 682
 683   //----------------------------------------------------------
 684   // tracking media for gases
 685   //----------------------------------------------------------
 686
 687   AliMedium(0, "Air", 11, 0, iSXFLD, sXMGMX, 10., 999., .1, .01, .1);
 688   AliMedium(1, "Ne-CO2-1", 12, 0, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
 689   AliMedium(2, "Ne-CO2-2", 13, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
 690   AliMedium(3,"CO2",10,0, iSXFLD, sXMGMX, 10., 999.,.1, .001, .001);
 691   AliMedium(20, "Ne-CO2-3", 35, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
 692   //-----------------------------------------------------------
 693   // tracking media for solids
 694   //-----------------------------------------------------------
 695
 696   AliMedium(4,"Al",23,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 697   AliMedium(5,"Kevlar",14,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 698   AliMedium(6,"Nomex",15,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 699   AliMedium(7,"Makrolon",16,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 700   AliMedium(8,"Mylar",18,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 701   AliMedium(9,"Tedlar",17,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 702   //
 703   AliMedium(10,"Prepreg1",19,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 704   AliMedium(11,"Prepreg2",20,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 705   AliMedium(12,"Prepreg3",21,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 706   AliMedium(13,"Epoxy",26,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 707
 708   AliMedium(14,"Cu",25,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 709   AliMedium(15,"Si",24,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 710   AliMedium(16,"G10",22,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 711   AliMedium(17,"Plexiglas",27,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 712   AliMedium(18,"Steel",29,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 713   AliMedium(19,"Peek",30,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 714   AliMedium(21,"Alumina",31,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 715   AliMedium(22,"Water",32,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 716   AliMedium(23,"Brass",33,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
 717   AliMedium(24,"Epoxyfm",34,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
 718 }
 719
 720 void AliTPC::GenerNoise(Int_t tablesize)
 721 {
 722   //
 723   //Generate table with noise
 724   //
 725   if (fTPCParam==0) {
 726     // error message
 727     fNoiseDepth=0;
 728     return;
 729   }
 730   if (fNoiseTable)  delete[] fNoiseTable;
 731   fNoiseTable = new Float_t[tablesize];
 732   fNoiseDepth = tablesize;
 733   fCurrentNoise =0; //!index of the noise in  the noise table
 734
 735   Float_t norm = fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac();
 736   for (Int_t i=0;i<tablesize;i++) fNoiseTable[i]= gRandom->Gaus(0,norm);
 737 }
 738
 739 Float_t AliTPC::GetNoise()
 740 {
 741   // get noise from table
 742   //  if ((fCurrentNoise%10)==0)
 743   //  fCurrentNoise= gRandom->Rndm()*fNoiseDepth;
 744   if (fCurrentNoise>=fNoiseDepth) fCurrentNoise=0;
 745   return fNoiseTable[fCurrentNoise++];
 746   //gRandom->Gaus(0, fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac());
 747 }
 748
 749
 750 Bool_t  AliTPC::IsSectorActive(Int_t sec) const
 751 {
 752   //
 753   // check if the sector is active
 754   if (!fActiveSectors) return kTRUE;
 755   else return fActiveSectors[sec];
 756 }
 757
 758 void    AliTPC::SetActiveSectors(Int_t * sectors, Int_t n)
 759 {
 760   // activate interesting sectors
 761   SetTreeAddress();//just for security
 762   if (!fActiveSectors) fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
 763   for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
 764   for (Int_t i=0;i<n;i++)
 765     if ((sectors[i]>=0) && sectors[i]<fTPCParam->GetNSector())  fActiveSectors[sectors[i]]=kTRUE;
 766
 767 }
 768
 769 void    AliTPC::SetActiveSectors(Int_t flag)
 770 {
 771   //
 772   // activate sectors which were hitted by tracks
 773   //loop over tracks
 774   SetTreeAddress();//just for security
 775   if (fHitType==0) return;  // if Clones hit - not short volume ID information
 776   if (!fActiveSectors) fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
 777   if (flag) {
 778     for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kTRUE;
 779     return;
 780   }
 781   for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
 782   //TBranch * branch=0;
 783   if (fLoader->TreeH() == 0x0)
 784    {
 785      AliFatal("Can not find TreeH in folder");
 786      return;
 787    }
 788   //if (fHitType>1) branch = fLoader->TreeH()->GetBranch("TPC2");
 789   if (fHitType>1) fLoader->TreeH()->GetBranch("TPC2");
 790   //else branch = fLoader->TreeH()->GetBranch("TPC");
 791   else fLoader->TreeH()->GetBranch("TPC");
 792   Stat_t ntracks = fLoader->TreeH()->GetEntries();
 793   // loop over all hits
 794   AliDebug(1,Form("Got %d tracks", (Int_t) ntracks));
 795
 796   for(Int_t track=0;track<ntracks;track++) {
 797     ResetHits();
 798     //
 799     if (fTrackHits && fHitType&4) {
 800       TBranch * br1 = fLoader->TreeH()->GetBranch("fVolumes");
 801       TBranch * br2 = fLoader->TreeH()->GetBranch("fNVolumes");
 802       br1->GetEvent(track);
 803       br2->GetEvent(track);
 804       Int_t *volumes = fTrackHits->GetVolumes();
 805       for (Int_t j=0;j<fTrackHits->GetNVolumes(); j++) {
 806         if (volumes[j]>-1 && volumes[j]<fTPCParam->GetNSector()) {
 807           fActiveSectors[volumes[j]]=kTRUE;
 808         }
 809         else {
 810             AliError(Form("Volume %d -> sector number %d is outside (0..%d)",
 811                           j,
 812                           volumes[j],
 813                           fTPCParam->GetNSector()));
 814         }
 815       }
 816     }
 817
 818     //
 819 //     if (fTrackHitsOld && fHitType&2) {
 820 //       TBranch * br = fLoader->TreeH()->GetBranch("fTrackHitsInfo");
 821 //       br->GetEvent(track);
 822 //       AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
 823 //       for (UInt_t j=0;j<ar->GetSize();j++){
 824 //      fActiveSectors[((AliTrackHitsInfo*)ar->At(j))->fVolumeID] =kTRUE;
 825 //       }
 826 //     }
 827   }
 828 }
 829
 830
 831
 832
 833 //_____________________________________________________________________________
 834 void AliTPC::Digits2Raw()
 835 {
 836 // convert digits of the current event to raw data
 837
 838   static const Int_t kThreshold = 0;
 839
 840   fLoader->LoadDigits();
 841   TTree* digits = fLoader->TreeD();
 842   if (!digits) {
 843     AliError("No digits tree");
 844     return;
 845   }
 846
 847   //
 848   AliSimDigits digarr;
 849   AliSimDigits* digrow = &digarr;
 850   digits->GetBranch("Segment")->SetAddress(&digrow);
 851
 852   const char* fileName = "AliTPCDDL.dat";
 853   AliTPCBuffer* buffer  = new AliTPCBuffer(fileName);
 854   //Verbose level
 855   // 0: Silent
 856   // 1: cout messages
 857   // 2: txt files with digits
 858   //BE CAREFUL, verbose level 2 MUST be used only for debugging and
 859   //it is highly suggested to use this mode only for debugging digits files
 860   //reasonably small, because otherwise the size of the txt files can reach
 861   //quickly several MB wasting time and disk space.
 862   buffer->SetVerbose(0);
 863
 864   Int_t nEntries = Int_t(digits->GetEntries());
 865   Int_t previousSector = -1;
 866   Int_t subSector = 0;
 867   for (Int_t i = 0; i < nEntries; i++) {
 868     digits->GetEntry(i);
 869     Int_t sector, row;
 870     fTPCParam->AdjustSectorRow(digarr.GetID(), sector, row);
 871     if(previousSector != sector) {
 872       subSector = 0;
 873       previousSector = sector;
 874     }
 875
 876     if (sector < 36) { //inner sector [0;35]
 877       if (row != 30) {
 878         //the whole row is written into the output file
 879         buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
 880                                sector, subSector, row);
 881       } else {
 882         //only the pads in the range [37;48] are written into the output file
 883         buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 1,
 884                                sector, subSector, row);
 885         subSector = 1;
 886         //only the pads outside the range [37;48] are written into the output file
 887         buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 2,
 888                                sector, subSector, row);
 889       }//end else
 890
 891     } else { //outer sector [36;71]
 892       if (row == 54) subSector = 2;
 893       if ((row != 27) && (row != 76)) {
 894         buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
 895                                sector, subSector, row);
 896       } else if (row == 27) {
 897         //only the pads outside the range [43;46] are written into the output file
 898         buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 2,
 899                                  sector, subSector, row);
 900         subSector = 1;
 901         //only the pads in the range [43;46] are written into the output file
 902         buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 1,
 903                                  sector, subSector, row);
 904       } else if (row == 76) {
 905         //only the pads outside the range [33;88] are written into the output file
 906         buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 2,
 907                                sector, subSector, row);
 908         subSector = 3;
 909         //only the pads in the range [33;88] are written into the output file
 910         buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 1,
 911                                  sector, subSector, row);
 912       }
 913     }//end else
 914   }//end for
 915
 916   delete buffer;
 917   fLoader->UnloadDigits();
 918
 919   AliTPCDDLRawData rawWriter;
 920   rawWriter.SetVerbose(0);
 921
 922   rawWriter.RawData(fileName);
 923   gSystem->Unlink(fileName);
 924
 925 }
 926
 927
 928 //_____________________________________________________________________________
 929 Bool_t AliTPC::Raw2SDigits(AliRawReader* rawReader){
 930   // Converts the TPC raw data into summable digits
 931   // The method is used for merging simulated and
 932   // real data events
 933   if (fLoader->TreeS() == 0x0 ) {
 934     fLoader->MakeTree("S");
 935   }
 936
 937   if(fDefaults == 0) SetDefaults();  // check if the parameters are set
 938
 939   //setup TPCDigitsArray
 940   if(GetDigitsArray()) delete GetDigitsArray();
 941
 942   AliTPCDigitsArray *arr = new AliTPCDigitsArray;
 943   arr->SetClass("AliSimDigits");
 944   arr->Setup(fTPCParam);
 945   arr->MakeTree(fLoader->TreeS());
 946
 947   SetDigitsArray(arr);
 948
 949   // set zero suppression to "0"
 950   fTPCParam->SetZeroSup(0);
 951
 952   // Loop over sectors
 953   const Int_t kmaxTime = fTPCParam->GetMaxTBin();
 954   const Int_t kNIS = fTPCParam->GetNInnerSector();
 955   const Int_t kNOS = fTPCParam->GetNOuterSector();
 956   const Int_t kNS = kNIS + kNOS;
 957
 958   // Setup storage
 959   AliTPCROC * roc = AliTPCROC::Instance();
 960   Int_t nRowsMax = roc->GetNRows(roc->GetNSector()-1);
 961   Int_t nPadsMax = roc->GetNPads(roc->GetNSector()-1,nRowsMax-1);
 962   Short_t** allBins = new Short_t*[nRowsMax];
 963   for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
 964     Int_t maxBin = kmaxTime*nPadsMax;
 965     allBins[iRow] = new Short_t[maxBin];
 966     memset(allBins[iRow],0,sizeof(Short_t)*maxBin);
 967   }
 968
 969   for(Int_t iSector = 0; iSector < kNS; iSector++) {
 970
 971     Int_t nRows = fTPCParam->GetNRow(iSector);
 972     Int_t nDDLs = 0, indexDDL = 0;
 973     if (iSector < kNIS) {
 974       nDDLs = 2;
 975       indexDDL = iSector * 2;
 976     }
 977     else {
 978       nDDLs = 4;
 979       indexDDL = (iSector-kNIS) * 4 + kNIS * 2;
 980     }
 981
 982     // Load the raw data for corresponding DDLs
 983     rawReader->Reset();
 984
 985     AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
 986     AliTPCRawStreamV3 input(rawReader,(AliAltroMapping**)mapping);
 987     rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
 988
 989     // Clean storage
 990     for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
 991       Int_t maxBin = kmaxTime*nPadsMax;
 992       memset(allBins[iRow],0,sizeof(Short_t)*maxBin);
 993     }
 994
 995     // Begin loop over altro data
 996     while (input.NextDDL()) {
 997
 998       if (input.GetSector() != iSector)
 999         AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSector,input.GetSector()));
1000
1001       //loop over pads
1002       while ( input.NextChannel() ) {
1003
1004         Int_t iRow = input.GetRow();
1005         if (iRow < 0 || iRow >= nRows)
1006           AliFatal(Form("Pad-row index (%d) outside the range (%d -> %d) !",
1007                         iRow, 0, nRows -1));
1008         Int_t iPad = input.GetPad();
1009
1010         Int_t maxPad = fTPCParam->GetNPads(iSector,iRow);
1011
1012         if (iPad < 0 || iPad >= maxPad)
1013           AliFatal(Form("Pad index (%d) outside the range (%d -> %d) !",
1014                         iPad, 0, maxPad -1));
1015
1016         //loop over bunches
1017         while ( input.NextBunch() ){
1018           Int_t  startTbin    = (Int_t)input.GetStartTimeBin();
1019           Int_t  bunchlength  = (Int_t)input.GetBunchLength();
1020           const UShort_t *sig = input.GetSignals();
1021           for (Int_t iTime = 0; iTime<bunchlength; iTime++){
1022             Int_t iTimeBin=startTbin-iTime;
1023             if ( iTimeBin < 0 || iTimeBin >= kmaxTime) {
1024               continue;
1025               //AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
1026               //               iTimeBin, 0, kmaxTime -1));
1027             }
1028
1029             Int_t maxBin = kmaxTime*maxPad;
1030             if (((iPad*kmaxTime+iTimeBin) >= maxBin) ||
1031                 ((iPad*kmaxTime+iTimeBin) < 0))
1032               AliFatal(Form("Index outside the allowed range"
1033                             " Sector=%d Row=%d Pad=%d Timebin=%d"
1034                             " (Max.index=%d)",iSector,iRow,iPad,iTimeBin,maxBin));
1035             allBins[iRow][iPad*kmaxTime+iTimeBin] = sig[iTime];
1036           }
1037         }
1038       } // End loop over altro data
1039     }
1040
1041     // Now fill the digits array
1042     if (fDigitsArray->GetTree()==0) {
1043       AliFatal("Tree not set in fDigitsArray");
1044     }
1045
1046     for (Int_t iRow = 0; iRow < nRows; iRow++) {
1047       AliDigits * dig = fDigitsArray->CreateRow(iSector,iRow);
1048       Int_t maxPad = fTPCParam->GetNPads(iSector,iRow);
1049       for(Int_t iPad = 0; iPad < maxPad; iPad++) {
1050         for(Int_t iTimeBin = 0; iTimeBin < kmaxTime; iTimeBin++) {
1051           Short_t q = allBins[iRow][iPad*kmaxTime + iTimeBin];
1052           if (q <= 0) continue;
1053           q *= 16;
1054           dig->SetDigitFast((Short_t)q,iTimeBin,iPad);
1055         }
1056       }
1057       fDigitsArray->StoreRow(iSector,iRow);
1058       Int_t ndig = dig->GetDigitSize();
1059
1060       AliDebug(10,
1061                Form("*** Sector, row, compressed digits %d %d %d ***\n",
1062                     iSector,iRow,ndig));
1063
1064       fDigitsArray->ClearRow(iSector,iRow);
1065
1066     } // end of the sector digitization
1067   }
1068   // get LHC clock phase from the digits tree
1069
1070   TParameter<float> *ph;
1071   Float_t phase;
1072   TTree *digtree = fLoader->TreeD();
1073   //
1074   if(digtree){ // if TreeD exists
1075     ph = (TParameter<float>*)digtree->GetUserInfo()->FindObject("lhcphase0");
1076     phase = ph->GetVal();
1077   }
1078   else{ //TreeD does not exist
1079     phase = 0.;
1080   }
1081     //
1082     // store lhc clock phase in S-digits tree
1083     //
1084     fLoader->TreeS()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",phase));
1085    //
1086    fLoader->WriteSDigits("OVERWRITE");
1087
1088   if(GetDigitsArray()) delete GetDigitsArray();
1089   SetDigitsArray(0x0);
1090
1091   // cleanup storage
1092   for (Int_t iRow = 0; iRow < nRowsMax; iRow++)
1093     delete [] allBins[iRow];
1094   delete [] allBins;
1095
1096   return kTRUE;
1097 }
1098
1099 //______________________________________________________________________
1100 AliDigitizer* AliTPC::CreateDigitizer(AliDigitizationInput* digInput) const
1101 {
1102   return new AliTPCDigitizer(digInput);
1103 }
1104 //__
1105 void AliTPC::SDigits2Digits2(Int_t /*eventnumber*/)
1106 {
1107   //create digits from summable digits
1108   GenerNoise(500000); //create teble with noise
1109
1110   //conect tree with sSDigits
1111   TTree *t = fLoader->TreeS();
1112
1113   if (t == 0x0) {
1114     fLoader->LoadSDigits("READ");
1115     t = fLoader->TreeS();
1116     if (t == 0x0) {
1117       AliError("Can not get input TreeS");
1118       return;
1119     }
1120   }
1121
1122   if (fLoader->TreeD() == 0x0) fLoader->MakeTree("D");
1123
1124   AliSimDigits digarr, *dummy=&digarr;
1125   TBranch* sdb = t->GetBranch("Segment");
1126   if (sdb == 0x0) {
1127     AliError("Can not find branch with segments in TreeS.");
1128     return;
1129   }
1130
1131   sdb->SetAddress(&dummy);
1132
1133   Stat_t nentries = t->GetEntries();
1134
1135   // set zero suppression
1136
1137   fTPCParam->SetZeroSup(2);
1138
1139   // get zero suppression
1140
1141   Int_t zerosup = fTPCParam->GetZeroSup();
1142
1143   //make tree with digits
1144
1145   AliTPCDigitsArray *arr = new AliTPCDigitsArray;
1146   arr->SetClass("AliSimDigits");
1147   arr->Setup(fTPCParam);
1148   arr->MakeTree(fLoader->TreeD());
1149
1150   AliTPCParam * par = fTPCParam;
1151
1152   //Loop over segments of the TPC
1153
1154   for (Int_t n=0; n<nentries; n++) {
1155     t->GetEvent(n);
1156     Int_t sec, row;
1157     if (!par->AdjustSectorRow(digarr.GetID(),sec,row)) {
1158       AliWarning(Form("Invalid segment ID ! %d",digarr.GetID()));
1159       continue;
1160     }
1161     if (!IsSectorActive(sec)) continue;
1162
1163     AliSimDigits * digrow =(AliSimDigits*) arr->CreateRow(sec,row);
1164     Int_t nrows = digrow->GetNRows();
1165     Int_t ncols = digrow->GetNCols();
1166
1167     digrow->ExpandBuffer();
1168     digarr.ExpandBuffer();
1169     digrow->ExpandTrackBuffer();
1170     digarr.ExpandTrackBuffer();
1171
1172
1173     Short_t * pamp0 = digarr.GetDigits();
1174     Int_t   * ptracks0 = digarr.GetTracks();
1175     Short_t * pamp1 = digrow->GetDigits();
1176     Int_t   * ptracks1 = digrow->GetTracks();
1177     Int_t  nelems =nrows*ncols;
1178     Int_t saturation = fTPCParam->GetADCSat() - 1;
1179     //use internal structure of the AliDigits - for speed reason
1180     //if you cahnge implementation
1181     //of the Alidigits - it must be rewriten -
1182     for (Int_t i= 0; i<nelems; i++){
1183       Float_t q = TMath::Nint(Float_t(*pamp0)/16.+GetNoise());
1184       if (q>zerosup){
1185         if (q>saturation) q=saturation;
1186         *pamp1=(Short_t)q;
1187
1188         ptracks1[0]=ptracks0[0];
1189         ptracks1[nelems]=ptracks0[nelems];
1190         ptracks1[2*nelems]=ptracks0[2*nelems];
1191       }
1192       pamp0++;
1193       pamp1++;
1194       ptracks0++;
1195       ptracks1++;
1196     }
1197
1198     arr->StoreRow(sec,row);
1199     arr->ClearRow(sec,row);
1200   }
1201
1202
1203   //write results
1204   fLoader->WriteDigits("OVERWRITE");
1205
1206   delete arr;
1207 }
1208 //__________________________________________________________________
1209 void AliTPC::SetDefaults(){
1210   //
1211   // setting the defaults
1212   //
1213
1214   // Set response functions
1215
1216   //
1217   AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
1218   rl->CdGAFile();
1219   //AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
1220   //gDirectory->Get("75x40_100x60");
1221   AliTPCParamSR *param = (AliTPCParamSR*)AliTPCcalibDB::Instance()->GetParameters();
1222   if(!param){
1223     AliFatal("No TPC parameters found");
1224     return;
1225   }
1226   if (!param->IsGeoRead()){
1227       //
1228       // read transformation matrices for gGeoManager
1229       //
1230       param->ReadGeoMatrices();
1231     }
1232
1233
1234
1235   AliTPCPRF2D    * prfinner   = new AliTPCPRF2D;
1236   AliTPCPRF2D    * prfouter1   = new AliTPCPRF2D;
1237   AliTPCPRF2D    * prfouter2   = new AliTPCPRF2D;
1238
1239
1240   //AliTPCRF1D     * rf    = new AliTPCRF1D(kTRUE);
1241   //rf->SetGauss(param->GetZSigma(),param->GetZWidth(),1.);
1242   //rf->SetOffset(3*param->GetZSigma());
1243   //rf->Update();
1244   //
1245   // Use gamma 4
1246   //
1247   char  strgamma4[1000];
1248   //sprintf(strgamma4,"AliTPCRF1D::Gamma4((x-0.135+%f)*%f,55,160)",3*param->GetZSigma(), 1000000000*param->GetTSample()/param->GetZWidth());
1249
1250   snprintf(strgamma4,1000,"AliTPCRF1D::Gamma4((x-0.135+%f)*%f,55,160)",3*param->GetZSigma(), 1000000000*param->GetTSample()/param->GetZWidth());
1251   TF1 * fgamma4 = new TF1("fgamma4",strgamma4, -1,1);
1252   AliTPCRF1D     * rf    = new AliTPCRF1D(kTRUE,1000);
1253   rf->SetParam(fgamma4,param->GetZWidth(), 1,0.2);
1254   rf->SetOffset(3*param->GetZSigma());
1255   rf->Update();
1256   TDirectory *savedir=gDirectory;
1257
1258   if (fIsGEM==0){
1259     printf ("TPC MWPC readout\n");
1260     TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
1261     if (!f->IsOpen())
1262       AliFatal("Can't open $ALICE_ROOT/TPC/AliTPCprf2d.root !");
1263
1264     TString s;
1265     prfinner->Read("prf_07504_Gati_056068_d02");
1266     //PH Set different names
1267     s=prfinner->GetGRF()->GetName();
1268     s+="in";
1269     prfinner->GetGRF()->SetName(s.Data());
1270
1271     prfouter1->Read("prf_10006_Gati_047051_d03");
1272     s=prfouter1->GetGRF()->GetName();
1273     s+="out1";
1274     prfouter1->GetGRF()->SetName(s.Data());
1275
1276     prfouter2->Read("prf_15006_Gati_047051_d03");
1277     s=prfouter2->GetGRF()->GetName();
1278     s+="out2";
1279     prfouter2->GetGRF()->SetName(s.Data());
1280     f->Close();
1281   }
1282
1283   if (fIsGEM==1){
1284     printf ("TPC GEM readout\n");
1285     TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2dGEM.root");
1286     if (!f->IsOpen())
1287       AliFatal("Can't open $ALICE_ROOT/TPC/AliTPCprf2dGEM.root !");
1288
1289     TString s;
1290     prfinner->Read("prf0");
1291     //PH Set different names
1292     s=prfinner->GetGRF()->GetName();
1293     s+="in";
1294     prfinner->GetGRF()->SetName(s.Data());
1295
1296     prfouter1->Read("prf1");
1297     s=prfouter1->GetGRF()->GetName();
1298     s+="out1";
1299     prfouter1->GetGRF()->SetName(s.Data());
1300
1301     prfouter2->Read("prf2");
1302     s=prfouter2->GetGRF()->GetName();
1303     s+="out2";
1304     prfouter2->GetGRF()->SetName(s.Data());
1305     f->Close();
1306   }
1307   savedir->cd();
1308
1309   param->SetInnerPRF(prfinner);
1310   param->SetOuter1PRF(prfouter1);
1311   param->SetOuter2PRF(prfouter2);
1312   param->SetTimeRF(rf);
1313
1314   // set fTPCParam
1315
1316   SetParam(param);
1317
1318
1319   fDefaults = 1;
1320
1321 }
1322 //__________________________________________________________________
1323 void AliTPC::Hits2Digits()
1324 {
1325   //
1326   // creates digits from hits
1327   //
1328   if (!fTPCParam->IsGeoRead()){
1329     //
1330     // read transformation matrices for gGeoManager
1331     //
1332     fTPCParam->ReadGeoMatrices();
1333   }
1334
1335   fLoader->LoadHits("read");
1336   fLoader->LoadDigits("recreate");
1337   AliRunLoader* runLoader = fLoader->GetRunLoader();
1338
1339   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1340     //PH    runLoader->GetEvent(iEvent);
1341     Hits2Digits(iEvent);
1342   }
1343
1344   fLoader->UnloadHits();
1345   fLoader->UnloadDigits();
1346 }
1347 //__________________________________________________________________
1348 void AliTPC::Hits2Digits(Int_t eventnumber)
1349 {
1350  //----------------------------------------------------
1351  // Loop over all sectors for a single event
1352  //----------------------------------------------------
1353   AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
1354   rl->GetEvent(eventnumber);
1355   SetActiveSectors();
1356   if (fLoader->TreeH() == 0x0) {
1357     if(fLoader->LoadHits()) {
1358       AliError("Can not load hits.");
1359     }
1360   }
1361   SetTreeAddress();
1362
1363   if (fLoader->TreeD() == 0x0 ) {
1364     fLoader->MakeTree("D");
1365     if (fLoader->TreeD() == 0x0 ) {
1366       AliError("Can not get TreeD");
1367       return;
1368     }
1369   }
1370
1371   if(fDefaults == 0) SetDefaults();  // check if the parameters are set
1372   GenerNoise(500000); //create teble with noise
1373
1374   //setup TPCDigitsArray
1375
1376   if(GetDigitsArray()) delete GetDigitsArray();
1377
1378   AliTPCDigitsArray *arr = new AliTPCDigitsArray;
1379   arr->SetClass("AliSimDigits");
1380   arr->Setup(fTPCParam);
1381
1382   arr->MakeTree(fLoader->TreeD());
1383   SetDigitsArray(arr);
1384
1385   fDigitsSwitch=0; // standard digits
1386   // here LHC clock phase
1387   Float_t lhcph = 0.;
1388   switch (fLHCclockPhaseSw){
1389   case 0:
1390     // no phase
1391     lhcph=0.;
1392     break;
1393   case 1:
1394     // random phase
1395     lhcph = (Int_t)(gRandom->Rndm()/0.25);
1396     break;
1397   case 2:
1398     lhcph=0.;
1399     // not implemented yet
1400     break;
1401   }
1402   // adding phase to the TreeD user info
1403   fLoader->TreeD()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",lhcph));
1404   //
1405   for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++)
1406     if (IsSectorActive(isec)) {
1407       AliDebug(1,Form("Hits2Digits: Sector %d is active.",isec));
1408       Hits2DigitsSector(isec);
1409     }
1410     else {
1411       AliDebug(1,Form("Hits2Digits: Sector %d is NOT active.",isec));
1412     }
1413
1414   fLoader->WriteDigits("OVERWRITE");
1415
1416 //this line prevents the crash in the similar one
1417 //on the beginning of this method
1418 //destructor attempts to reset the tree, which is deleted by the loader
1419 //need to be redesign
1420   if(GetDigitsArray()) delete GetDigitsArray();
1421   SetDigitsArray(0x0);
1422
1423 }
1424
1425 //__________________________________________________________________
1426 void AliTPC::Hits2SDigits2(Int_t eventnumber)
1427 {
1428
1429   //-----------------------------------------------------------
1430   //   summable digits - 16 bit "ADC", no noise, no saturation
1431   //-----------------------------------------------------------
1432
1433   //----------------------------------------------------
1434   // Loop over all sectors for a single event
1435   //----------------------------------------------------
1436
1437   AliRunLoader* rl = fLoader->GetRunLoader();
1438
1439   rl->GetEvent(eventnumber);
1440   if (fLoader->TreeH() == 0x0) {
1441     if(fLoader->LoadHits()) {
1442       AliError("Can not load hits.");
1443       return;
1444     }
1445   }
1446   SetTreeAddress();
1447
1448
1449   if (fLoader->TreeS() == 0x0 ) {
1450     fLoader->MakeTree("S");
1451   }
1452
1453   if(fDefaults == 0) SetDefaults();
1454
1455   GenerNoise(500000); //create table with noise
1456   //setup TPCDigitsArray
1457
1458   if(GetDigitsArray()) delete GetDigitsArray();
1459
1460
1461   AliTPCDigitsArray *arr = new AliTPCDigitsArray;
1462   arr->SetClass("AliSimDigits");
1463   arr->Setup(fTPCParam);
1464   arr->MakeTree(fLoader->TreeS());
1465
1466   SetDigitsArray(arr);
1467
1468   fDigitsSwitch=1; // summable digits
1469
1470     // set zero suppression to "0"
1471   // here LHC clock phase
1472   Float_t lhcph = 0.;
1473   switch (fLHCclockPhaseSw){
1474   case 0:
1475     // no phase
1476     lhcph=0.;
1477     break;
1478   case 1:
1479     // random phase
1480     lhcph = (Int_t)(gRandom->Rndm()/0.25);
1481     break;
1482   case 2:
1483     lhcph=0.;
1484     // not implemented yet
1485     break;
1486   }
1487   // adding phase to the TreeS user info
1488
1489   fLoader->TreeS()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",lhcph));
1490
1491   fTPCParam->SetZeroSup(0);
1492
1493   for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++)
1494     if (IsSectorActive(isec)) {
1495       Hits2DigitsSector(isec);
1496     }
1497
1498   fLoader->WriteSDigits("OVERWRITE");
1499
1500 //this line prevents the crash in the similar one
1501 //on the beginning of this method
1502 //destructor attempts to reset the tree, which is deleted by the loader
1503 //need to be redesign
1504   if(GetDigitsArray()) delete GetDigitsArray();
1505   SetDigitsArray(0x0);
1506 }
1507 //__________________________________________________________________
1508
1509 void AliTPC::Hits2SDigits()
1510 {
1511
1512   //-----------------------------------------------------------
1513   //   summable digits - 16 bit "ADC", no noise, no saturation
1514   //-----------------------------------------------------------
1515
1516   if (!fTPCParam->IsGeoRead()){
1517     //
1518     // read transformation matrices for gGeoManager
1519     //
1520     fTPCParam->ReadGeoMatrices();
1521   }
1522
1523   fLoader->LoadHits("read");
1524   fLoader->LoadSDigits("recreate");
1525   AliRunLoader* runLoader = fLoader->GetRunLoader();
1526
1527   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1528     runLoader->GetEvent(iEvent);
1529     SetTreeAddress();
1530     SetActiveSectors();
1531     Hits2SDigits2(iEvent);
1532   }
1533
1534   fLoader->UnloadHits();
1535   fLoader->UnloadSDigits();
1536   if (fDebugStreamer) {
1537     delete fDebugStreamer;
1538     fDebugStreamer=0;
1539   }
1540 }
1541 //_____________________________________________________________________________
1542
1543 void AliTPC::Hits2DigitsSector(Int_t isec)
1544 {
1545   //-------------------------------------------------------------------
1546   // TPC conversion from hits to digits.
1547   //-------------------------------------------------------------------
1548
1549   //-----------------------------------------------------------------
1550   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1551   //-----------------------------------------------------------------
1552
1553   //-------------------------------------------------------
1554   //  Get the access to the track hits
1555   //-------------------------------------------------------
1556
1557   // check if the parameters are set - important if one calls this method
1558   // directly, not from the Hits2Digits
1559
1560   if(fDefaults == 0) SetDefaults();
1561
1562   TTree *tH = fLoader->TreeH(); // pointer to the hits tree
1563   if (tH == 0x0) {
1564     AliFatal("Can not find TreeH in folder");
1565     return;
1566   }
1567
1568   Stat_t ntracks = tH->GetEntries();
1569
1570     Int_t nrows =fTPCParam->GetNRow(isec);
1571
1572     TObjArray **row=new TObjArray* [nrows+2]; // 2 extra rows for cross talk
1573     for(Int_t j=0;j<nrows+2;j++) row[j]=0;
1574
1575     MakeSector(isec,nrows,tH,ntracks,row);
1576
1577     //--------------------------------------------------------
1578     //   Digitize this sector, row by row
1579     //   row[i] is the pointer to the TObjArray of TVectors,
1580     //   each one containing electrons accepted on this
1581     //   row, assigned into tracks
1582     //--------------------------------------------------------
1583
1584     Int_t i;
1585
1586     if (fDigitsArray->GetTree()==0) {
1587       AliFatal("Tree not set in fDigitsArray");
1588     }
1589
1590     for (i=0;i<nrows;i++){
1591
1592       AliDigits * dig = fDigitsArray->CreateRow(isec,i);
1593
1594       DigitizeRow(i,isec,row);
1595
1596       fDigitsArray->StoreRow(isec,i);
1597
1598       Int_t ndig = dig->GetDigitSize();
1599
1600       AliDebug(10,
1601                Form("*** Sector, row, compressed digits %d %d %d ***\n",
1602                     isec,i,ndig));
1603
1604       fDigitsArray->ClearRow(isec,i);
1605
1606
1607     } // end of the sector digitization
1608
1609     for(i=0;i<nrows+2;i++){
1610       row[i]->Delete();
1611       delete row[i];
1612     }
1613
1614     delete [] row; // delete the array of pointers to TObjArray-s
1615
1616
1617 } // end of Hits2DigitsSector
1618
1619
1620 //_____________________________________________________________________________
1621 void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
1622 {
1623   //-----------------------------------------------------------
1624   // Single row digitization, coupling from the neighbouring
1625   // rows taken into account
1626   //-----------------------------------------------------------
1627
1628   //-----------------------------------------------------------------
1629   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1630   // Modified: Marian Ivanov GSI Darmstadt, m.ivanov@gsi.de
1631   //-----------------------------------------------------------------
1632
1633   Float_t zerosup = fTPCParam->GetZeroSup();
1634   AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
1635   AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
1636   AliTPCCalROC * gainROC = gainTPC->GetCalROC(isec);  // pad gains per given sector
1637   AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(isec);  // noise per given sector
1638
1639
1640   fCurrentIndex[1]= isec;
1641
1642
1643   Int_t nofPads = fTPCParam->GetNPads(isec,irow);
1644   Int_t nofTbins = fTPCParam->GetMaxTBin();
1645   Int_t indexRange[4];
1646   //
1647   //  Integrated signal for this row
1648   //  and a single track signal
1649   //
1650
1651   TMatrixF *m1 = new TMatrixF(0,nofPads,0,nofTbins); // integrated
1652   TMatrixF *m2 = new TMatrixF(0,nofPads,0,nofTbins); // single
1653   //
1654   TMatrixF &total  = *m1;
1655
1656   //  Array of pointers to the label-signal list
1657
1658   Int_t nofDigits = nofPads*nofTbins; // number of digits for this row
1659   Float_t  **pList = new Float_t* [nofDigits];
1660
1661   Int_t lp;
1662   Int_t i1;
1663   for(lp=0;lp<nofDigits;lp++)pList[lp]=0; // set all pointers to NULL
1664   //
1665   //calculate signal
1666   //
1667   Int_t row1=irow;
1668   Int_t row2=irow+2;
1669   for (Int_t row= row1;row<=row2;row++){
1670     Int_t nTracks= rows[row]->GetEntries();
1671     for (i1=0;i1<nTracks;i1++){
1672       fCurrentIndex[2]= row;
1673       fCurrentIndex[3]=irow+1;
1674       if (row==irow+1){
1675         m2->Zero();  // clear single track signal matrix
1676         Float_t trackLabel = GetSignal(rows[row],i1,m2,m1,indexRange);
1677         GetList(trackLabel,nofPads,m2,indexRange,pList);
1678       }
1679       else   GetSignal(rows[row],i1,0,m1,indexRange);
1680     }
1681   }
1682
1683   Int_t tracks[3];
1684
1685   AliDigits *dig = fDigitsArray->GetRow(isec,irow);
1686   Int_t gi=-1;
1687   Float_t fzerosup = zerosup+0.5;
1688   for(Int_t it=0;it<nofTbins;it++){
1689     for(Int_t ip=0;ip<nofPads;ip++){
1690       gi++;
1691       Float_t q=total(ip,it);
1692       if(fDigitsSwitch == 0){
1693         Float_t gain = gainROC->GetValue(irow,ip);  // get gain for given - pad-row pad
1694         Float_t noisePad = noiseROC->GetValue(irow,ip);
1695         //
1696         q*=gain;
1697         q+=GetNoise()*noisePad;
1698         if(q <=fzerosup) continue; // do not fill zeros
1699         q = TMath::Nint(q);
1700         if(q >= fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat() - 1;  // saturation
1701
1702       }
1703
1704       else {
1705         if(q <= 0.) continue; // do not fill zeros
1706         if(q>2000.) q=2000.;
1707         q *= 16.;
1708         q = TMath::Nint(q);
1709       }
1710
1711       //
1712       //  "real" signal or electronic noise (list = -1)?
1713       //
1714
1715       for(Int_t j1=0;j1<3;j1++){
1716         tracks[j1] = (pList[gi]) ?(Int_t)(*(pList[gi]+j1)) : -2;
1717       }
1718
1719 //Begin_Html
1720 /*
1721   <A NAME="AliDigits"></A>
1722   using of AliDigits object
1723 */
1724 //End_Html
1725       dig->SetDigitFast((Short_t)q,it,ip);
1726       if (fDigitsArray->IsSimulated()) {
1727         ((AliSimDigits*)dig)->SetTrackIDFast(tracks[0],it,ip,0);
1728         ((AliSimDigits*)dig)->SetTrackIDFast(tracks[1],it,ip,1);
1729         ((AliSimDigits*)dig)->SetTrackIDFast(tracks[2],it,ip,2);
1730       }
1731
1732     } // end of loop over time buckets
1733   }  // end of lop over pads
1734   //
1735   // test
1736   //
1737   //
1738
1739   // glitch filters if normal simulated digits
1740   //
1741   if(!fDigitsSwitch) ((AliSimDigits*)dig)->GlitchFilter();
1742   //
1743   //  This row has been digitized, delete nonused stuff
1744   //
1745
1746   for(lp=0;lp<nofDigits;lp++){
1747     if(pList[lp]) delete [] pList[lp];
1748   }
1749
1750   delete [] pList;
1751
1752   delete m1;
1753   delete m2;
1754
1755 } // end of DigitizeRow
1756
1757 //_____________________________________________________________________________
1758
1759 Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
1760              TMatrixF *m1, TMatrixF *m2,Int_t *indexRange)
1761 {
1762
1763   //---------------------------------------------------------------
1764   //  Calculates 2-D signal (pad,time) for a single track,
1765   //  returns a pointer to the signal matrix and the track label
1766   //  No digitization is performed at this level!!!
1767   //---------------------------------------------------------------
1768
1769   //-----------------------------------------------------------------
1770   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1771   // Modified: Marian Ivanov
1772   //-----------------------------------------------------------------
1773
1774   TVector *tv;
1775
1776   tv = (TVector*)p1->At(ntr); // pointer to a track
1777   TVector &v = *tv;
1778
1779   Float_t label = v(0);
1780   Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3]-1))/2;
1781
1782   Int_t nElectrons = (tv->GetNrows()-1)/5;
1783   indexRange[0]=9999; // min pad
1784   indexRange[1]=-1; // max pad
1785   indexRange[2]=9999; //min time
1786   indexRange[3]=-1; // max time
1787
1788   TMatrixF &signal = *m1;
1789   TMatrixF &total = *m2;
1790   //
1791   // Get LHC clock phase
1792   //
1793   TParameter<float> *ph;
1794   if(fDigitsSwitch){// s-digits
1795     ph = (TParameter<float>*)fLoader->TreeS()->GetUserInfo()->FindObject("lhcphase0");
1796   }
1797   else{ // normal digits
1798     ph = (TParameter<float>*)fLoader->TreeD()->GetUserInfo()->FindObject("lhcphase0");
1799   }
1800   //  Loop over all electrons
1801   //
1802   for(Int_t nel=0; nel<nElectrons; nel++){
1803     Int_t idx=nel*5;
1804     Float_t aval =  v(idx+4);
1805     Float_t eltoadcfac=aval*fTPCParam->GetTotalNormFac();
1806     Float_t xyz[4]={v(idx+1),v(idx+2),v(idx+3),v(idx+5)};
1807     Int_t n = ((AliTPCParamSR*)fTPCParam)->CalcResponseFast(xyz,fCurrentIndex,
1808                                                             fCurrentIndex[3],ph->GetVal());
1809
1810     Int_t *index = fTPCParam->GetResBin(0);
1811     Float_t *weight = & (fTPCParam->GetResWeight(0));
1812
1813     if (n>0) for (Int_t i =0; i<n; i++){
1814       Int_t pad=index[1]+centralPad;  //in digit coordinates central pad has coordinate 0
1815
1816       if (pad>=0){
1817         Int_t time=index[2];
1818         Float_t qweight = *(weight)*eltoadcfac;
1819
1820         if (m1!=0) signal(pad,time)+=qweight;
1821         total(pad,time)+=qweight;
1822         if (indexRange[0]>pad) indexRange[0]=pad;
1823         if (indexRange[1]<pad) indexRange[1]=pad;
1824         if (indexRange[2]>time) indexRange[2]=time;
1825         if (indexRange[3]<time) indexRange[3]=time;
1826
1827         index+=3;
1828         weight++;
1829
1830       }
1831     }
1832   } // end of loop over electrons
1833
1834   return label; // returns track label when finished
1835 }
1836
1837 //_____________________________________________________________________________
1838 void AliTPC::GetList(Float_t label,Int_t np,TMatrixF *m,
1839                      Int_t *indexRange, Float_t **pList)
1840 {
1841   //----------------------------------------------------------------------
1842   //  Updates the list of tracks contributing to digits for a given row
1843   //----------------------------------------------------------------------
1844
1845   //-----------------------------------------------------------------
1846   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1847   //-----------------------------------------------------------------
1848
1849   TMatrixF &signal = *m;
1850
1851   // lop over nonzero digits
1852
1853   for(Int_t it=indexRange[2];it<indexRange[3]+1;it++){
1854     for(Int_t ip=indexRange[0];ip<indexRange[1]+1;ip++){
1855
1856
1857       // accept only the contribution larger than 500 electrons (1/2 s_noise)
1858
1859       if(signal(ip,it)<0.5) continue;
1860
1861       Int_t globalIndex = it*np+ip; // globalIndex starts from 0!
1862
1863       if(!pList[globalIndex]){
1864
1865         //
1866         // Create new list (6 elements - 3 signals and 3 labels),
1867         //
1868
1869         pList[globalIndex] = new Float_t [6];
1870
1871         // set list to -1
1872
1873         *pList[globalIndex] = -1.;
1874         *(pList[globalIndex]+1) = -1.;
1875         *(pList[globalIndex]+2) = -1.;
1876         *(pList[globalIndex]+3) = -1.;
1877         *(pList[globalIndex]+4) = -1.;
1878         *(pList[globalIndex]+5) = -1.;
1879
1880         *pList[globalIndex] = label;
1881         *(pList[globalIndex]+3) = signal(ip,it);
1882       }
1883       else {
1884
1885         // check the signal magnitude
1886
1887         Float_t highest = *(pList[globalIndex]+3);
1888         Float_t middle = *(pList[globalIndex]+4);
1889         Float_t lowest = *(pList[globalIndex]+5);
1890
1891         //
1892         //  compare the new signal with already existing list
1893         //
1894
1895         if(signal(ip,it)<lowest) continue; // neglect this track
1896
1897         //
1898
1899         if (signal(ip,it)>highest){
1900           *(pList[globalIndex]+5) = middle;
1901           *(pList[globalIndex]+4) = highest;
1902           *(pList[globalIndex]+3) = signal(ip,it);
1903
1904           *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
1905           *(pList[globalIndex]+1) = *pList[globalIndex];
1906           *pList[globalIndex] = label;
1907         }
1908         else if (signal(ip,it)>middle){
1909           *(pList[globalIndex]+5) = middle;
1910           *(pList[globalIndex]+4) = signal(ip,it);
1911
1912           *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
1913           *(pList[globalIndex]+1) = label;
1914         }
1915         else{
1916           *(pList[globalIndex]+5) = signal(ip,it);
1917           *(pList[globalIndex]+2) = label;
1918         }
1919       }
1920
1921     } // end of loop over pads
1922   } // end of loop over time bins
1923
1924 }//end of GetList
1925 //___________________________________________________________________
1926 void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
1927                         Stat_t ntracks,TObjArray **row)
1928 {
1929
1930   //-----------------------------------------------------------------
1931   // Prepares the sector digitization, creates the vectors of
1932   // tracks for each row of this sector. The track vector
1933   // contains the track label and the position of electrons.
1934   //-----------------------------------------------------------------
1935
1936   //
1937   // The trasport of the electrons through TPC drift volume
1938   //    Drift (drift velocity + velocity map(not yet implemented)))
1939   //    Application of the random processes (diffusion, gas gain)
1940   //    Systematic effects (ExB effect in drift volume + ROCs)
1941   //
1942   // Algorithm:
1943   // Loop over primary electrons:
1944   //    Creation of the secondary electrons
1945   //    Loop over electrons (primary+ secondaries)
1946   //        Global coordinate frame:
1947   //          1. Skip electrons if attached
1948   //          2. ExB effect in drift volume
1949   //             a.) Simulation   calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
1950   //             b.) Reconstruction -  calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
1951   //          3. Generation of gas gain (Random - Exponential distribution)
1952   //          4. TransportElectron function (diffusion)
1953   //
1954   //        5. Conversion to the local coordinate frame  pad-row, pad, timebin
1955   //        6. Apply Time0 shift - AliTPCCalPad class
1956   //            a.) Plus sign in simulation
1957   //            b.) Minus sign in reconstruction
1958   // end of loop
1959   //
1960   //-----------------------------------------------------------------
1961   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1962   // Origin: Marian Ivanov,  marian.ivanov@cern.ch
1963   //-----------------------------------------------------------------
1964   AliTPCcalibDB* const calib=AliTPCcalibDB::Instance();
1965   if (gAlice){ // Set correctly the magnetic field in the ExB calculation
1966     if (!calib->GetExB()){
1967       AliMagF * field = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField());
1968       if (field) {
1969         calib->SetExBField(field);
1970       }
1971     }
1972   }
1973
1974   Float_t gasgain = fTPCParam->GetGasGain();
1975   gasgain = gasgain/fGainFactor;
1976   Int_t i;
1977   Float_t xyz[5];
1978
1979   AliTPChit *tpcHit; // pointer to a sigle TPC hit
1980   //MI change
1981   TBranch * branch=0;
1982   if (fHitType>1) branch = TH->GetBranch("TPC2");
1983   else branch = TH->GetBranch("TPC");
1984
1985
1986   //----------------------------------------------
1987   // Create TObjArray-s, one for each row,
1988   // each TObjArray will store the TVectors
1989   // of electrons, one TVectors per each track.
1990   //----------------------------------------------
1991
1992   Int_t *nofElectrons = new Int_t [nrows+2]; // electron counter for each row
1993   TVector **tracks = new TVector* [nrows+2]; //pointers to the track vectors
1994
1995   for(i=0; i<nrows+2; i++){
1996     row[i] = new TObjArray;
1997     nofElectrons[i]=0;
1998     tracks[i]=0;
1999   }
2000
2001
2002
2003   //--------------------------------------------------------------------
2004   //  Loop over tracks, the "track" contains the full history
2005   //--------------------------------------------------------------------
2006
2007   Int_t previousTrack,currentTrack;
2008   previousTrack = -1; // nothing to store so far!
2009
2010   for(Int_t track=0;track<ntracks;track++){
2011     Bool_t isInSector=kTRUE;
2012     ResetHits();
2013     isInSector = TrackInVolume(isec,track);
2014     if (!isInSector) continue;
2015     //MI change
2016     branch->GetEntry(track); // get next track
2017
2018     //M.I. changes
2019
2020     tpcHit = (AliTPChit*)FirstHit(-1);
2021
2022     //--------------------------------------------------------------
2023     //  Loop over hits
2024     //--------------------------------------------------------------
2025
2026
2027     while(tpcHit){
2028
2029       Int_t sector=tpcHit->fSector; // sector number
2030       if(sector != isec){
2031         tpcHit = (AliTPChit*) NextHit();
2032         continue;
2033       }
2034
2035       // Remove hits which arrive before the TPC opening gate signal
2036       if(((fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
2037           /fTPCParam->GetDriftV()+tpcHit->Time())<fTPCParam->GetGateDelay()) {
2038         tpcHit = (AliTPChit*) NextHit();
2039         continue;
2040       }
2041
2042       currentTrack = tpcHit->Track(); // track number
2043
2044       if(currentTrack != previousTrack){
2045
2046         // store already filled fTrack
2047
2048         for(i=0;i<nrows+2;i++){
2049           if(previousTrack != -1){
2050             if(nofElectrons[i]>0){
2051               TVector &v = *tracks[i];
2052               v(0) = previousTrack;
2053               tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
2054               row[i]->Add(tracks[i]);
2055             }
2056             else {
2057               delete tracks[i]; // delete empty TVector
2058               tracks[i]=0;
2059             }
2060           }
2061
2062           nofElectrons[i]=0;
2063           tracks[i] = new TVector(601); // TVectors for the next fTrack
2064
2065         } // end of loop over rows
2066
2067         previousTrack=currentTrack; // update track label
2068       }
2069
2070       Int_t qI = (Int_t) (tpcHit->fQ); // energy loss (number of electrons)
2071
2072       //---------------------------------------------------
2073       //  Calculate the electron attachment probability
2074       //---------------------------------------------------
2075
2076
2077       Float_t time = 1.e6*(fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
2078         /fTPCParam->GetDriftV();
2079       // in microseconds!
2080       Float_t attProb = fTPCParam->GetAttCoef()*
2081         fTPCParam->GetOxyCont()*time; //  fraction!
2082
2083       //-----------------------------------------------
2084       //  Loop over electrons
2085       //-----------------------------------------------
2086       Int_t index[3];
2087       index[1]=isec;
2088       for(Int_t nel=0;nel<qI;nel++){
2089         // skip if electron lost due to the attachment
2090         if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
2091
2092         //
2093         // ExB effect
2094         //
2095         Double_t dxyz0[3],dxyz1[3];
2096         dxyz0[0]=tpcHit->X();
2097         dxyz0[1]=tpcHit->Y();
2098         dxyz0[2]=tpcHit->Z();
2099         if (calib->GetExB()){
2100           calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
2101         }else{
2102           AliError("Not valid ExB calibration");
2103           dxyz1[0]=tpcHit->X();
2104           dxyz1[1]=tpcHit->Y();
2105           dxyz1[2]=tpcHit->Z();
2106         }
2107         xyz[0]=dxyz1[0];
2108         xyz[1]=dxyz1[1];
2109         xyz[2]=dxyz1[2];
2110         //
2111         //
2112         //
2113         // protection for the nonphysical avalanche size (10**6 maximum)
2114         //
2115         Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
2116         xyz[3]= (Float_t) (-gasgain*TMath::Log(rn));
2117         index[0]=1;
2118
2119         TransportElectron(xyz,index);
2120         Int_t rowNumber;
2121         Int_t padrow = fTPCParam->GetPadRow(xyz,index);
2122         //
2123         // Add Time0 correction due unisochronity
2124         // xyz[0] - pad row coordinate
2125         // xyz[1] - pad coordinate
2126         // xyz[2] - is in now time bin coordinate system
2127         Float_t correction =0;
2128         if (calib->GetPadTime0()){
2129           if (!calib->GetPadTime0()->GetCalROC(isec)) continue;
2130           Int_t npads = fTPCParam->GetNPads(isec,padrow);
2131           //      Int_t pad  = TMath::Nint(xyz[1]+fTPCParam->GetNPads(isec,TMath::Nint(xyz[0]))*0.5);
2132           // pad numbering from -npads/2 .. npads/2-1
2133           Int_t pad  = TMath::Nint(xyz[1]+npads/2);
2134           if (pad<0) pad=0;
2135           if (pad>=npads) pad=npads-1;
2136           correction = calib->GetPadTime0()->GetCalROC(isec)->GetValue(padrow,pad);
2137           //      printf("%d\t%d\t%d\t%f\n",isec,padrow,pad,correction);
2138           if (fDebugStreamer){
2139             (*fDebugStreamer)<<"Time0"<<
2140               "isec="<<isec<<
2141               "padrow="<<padrow<<
2142               "pad="<<pad<<
2143               "x0="<<xyz[0]<<
2144               "x1="<<xyz[1]<<
2145               "x2="<<xyz[2]<<
2146               "hit.="<<tpcHit<<
2147               "cor="<<correction<<
2148               "\n";
2149           }
2150         }
2151         xyz[2]+=correction;
2152         xyz[2]+=fTPCParam->GetNTBinsL1();    // adding Level 1 time bin offset
2153         //
2154         // Electron track time (for pileup simulation)
2155         xyz[2]+=tpcHit->Time()/fTPCParam->GetTSample(); // adding time of flight
2156         xyz[4] =0;
2157
2158         //
2159         // row 0 - cross talk from the innermost row
2160         // row fNRow+1 cross talk from the outermost row
2161         rowNumber = index[2]+1;
2162         //transform position to local digit coordinates
2163         //relative to nearest pad row
2164         if ((rowNumber<0)||rowNumber>fTPCParam->GetNRow(isec)+1) continue;
2165         /*      Float_t x1,y1;
2166         if (isec <fTPCParam->GetNInnerSector()) {
2167           x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
2168           y1 = fTPCParam->GetYInner(rowNumber);
2169         }
2170         else{
2171           x1=xyz[1]*fTPCParam->GetOuterPadPitchWidth();
2172           y1 = fTPCParam->GetYOuter(rowNumber);
2173         }
2174         // gain inefficiency at the wires edges - linear
2175         x1=TMath::Abs(x1);
2176         y1-=1.;
2177         if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.)); */
2178
2179         nofElectrons[rowNumber]++;
2180         //----------------------------------
2181         // Expand vector if necessary
2182         //----------------------------------
2183         if(nofElectrons[rowNumber]>120){
2184           Int_t range = tracks[rowNumber]->GetNrows();
2185           if((nofElectrons[rowNumber])>(range-1)/5){
2186
2187             tracks[rowNumber]->ResizeTo(range+500); // Add 100 electrons
2188           }
2189         }
2190
2191         TVector &v = *tracks[rowNumber];
2192         Int_t idx = 5*nofElectrons[rowNumber]-4;
2193         Real_t * position = &(((TVector&)v)(idx)); //make code faster
2194         memcpy(position,xyz,5*sizeof(Float_t));
2195
2196       } // end of loop over electrons
2197
2198       tpcHit = (AliTPChit*)NextHit();
2199
2200     } // end of loop over hits
2201   } // end of loop over tracks
2202
2203     //
2204     //   store remaining track (the last one) if not empty
2205     //
2206
2207   for(i=0;i<nrows+2;i++){
2208     if(nofElectrons[i]>0){
2209       TVector &v = *tracks[i];
2210       v(0) = previousTrack;
2211       tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
2212       row[i]->Add(tracks[i]);
2213     }
2214     else{
2215       delete tracks[i];
2216       tracks[i]=0;
2217     }
2218   }
2219
2220   delete [] tracks;
2221   delete [] nofElectrons;
2222
2223 } // end of MakeSector
2224
2225
2226 //_____________________________________________________________________________
2227 void AliTPC::Init()
2228 {
2229   //
2230   // Initialise TPC detector after definition of geometry
2231   //
2232   AliDebug(1,"*********************************************");
2233 }
2234
2235 //_____________________________________________________________________________
2236 void AliTPC::ResetDigits()
2237 {
2238   //
2239   // Reset number of digits and the digits array for this detector
2240   //
2241   fNdigits   = 0;
2242   if (fDigits)   fDigits->Clear();
2243 }
2244
2245
2246
2247 //_____________________________________________________________________________
2248 void AliTPC::SetSens(Int_t sens)
2249 {
2250
2251   //-------------------------------------------------------------
2252   // Activates/deactivates the sensitive strips at the center of
2253   // the pad row -- this is for the space-point resolution calculations
2254   //-------------------------------------------------------------
2255
2256   //-----------------------------------------------------------------
2257   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
2258   //-----------------------------------------------------------------
2259
2260   fSens = sens;
2261 }
2262
2263
2264 void AliTPC::SetSide(Float_t side=0.)
2265 {
2266   // choice of the TPC side
2267
2268   fSide = side;
2269
2270 }
2271 //_____________________________________________________________________________
2272
2273 void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
2274 {
2275   //
2276   // electron transport taking into account:
2277   // 1. diffusion,
2278   // 2.ExB at the wires
2279   // 3. nonisochronity
2280   //
2281   // xyz and index must be already transformed to system 1
2282   //
2283
2284   fTPCParam->Transform1to2(xyz,index);  // mis-alignment applied in this step
2285
2286   //add diffusion
2287   Float_t driftl=xyz[2];
2288   if(driftl<0.01) driftl=0.01;
2289   driftl=TMath::Sqrt(driftl);
2290   Float_t sigT = driftl*(fTPCParam->GetDiffT());
2291   Float_t sigL = driftl*(fTPCParam->GetDiffL());
2292   xyz[0]=gRandom->Gaus(xyz[0],sigT);
2293   xyz[1]=gRandom->Gaus(xyz[1],sigT);
2294   xyz[2]=gRandom->Gaus(xyz[2],sigL);
2295
2296   // ExB
2297
2298   if (fTPCParam->GetMWPCReadout()==kTRUE){
2299     Float_t dx = fTPCParam->Transform2to2NearestWire(xyz,index);
2300     xyz[1]+=dx*(fTPCParam->GetOmegaTau());
2301   }
2302   //add nonisochronity (not implemented yet)
2303
2304
2305 }
2306
2307 ClassImp(AliTPChit)
2308   //______________________________________________________________________
2309   AliTPChit::AliTPChit()
2310             :AliHit(),
2311              fSector(0),
2312              fPadRow(0),
2313              fQ(0),
2314              fTime(0)
2315 {
2316   //
2317   // default
2318   //
2319
2320 }
2321 //_____________________________________________________________________________
2322 AliTPChit::AliTPChit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits)
2323           :AliHit(shunt,track),
2324              fSector(0),
2325              fPadRow(0),
2326              fQ(0),
2327              fTime(0)
2328 {
2329   //
2330   // Creates a TPC hit object
2331   //
2332   fSector     = vol[0];
2333   fPadRow     = vol[1];
2334   fX          = hits[0];
2335   fY          = hits[1];
2336   fZ          = hits[2];
2337   fQ          = hits[3];
2338   fTime       = hits[4];
2339 }
2340
2341 //________________________________________________________________________
2342 // Additional code because of the AliTPCTrackHitsV2
2343
2344 void AliTPC::MakeBranch(Option_t *option)
2345 {
2346   //
2347   // Create a new branch in the current Root Tree
2348   // The branch of fHits is automatically split
2349   // MI change 14.09.2000
2350   AliDebug(1,"");
2351   if (fHitType<2) return;
2352   char branchname[10];
2353   //sprintf(branchname,"%s2",GetName());
2354   snprintf(branchname,10,"%s2",GetName());
2355   //
2356   // Get the pointer to the header
2357   const char *cH = strstr(option,"H");
2358   //
2359   if (fTrackHits   && fLoader->TreeH() && cH && fHitType&4) {
2360     AliDebug(1,"Making branch for Type 4 Hits");
2361     fLoader->TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
2362   }
2363
2364 //   if (fTrackHitsOld   && fLoader->TreeH() && cH && fHitType&2) {
2365 //     AliDebug(1,"Making branch for Type 2 Hits");
2366 //     AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHitsOld,
2367 //                                                    fLoader->TreeH(),fBufferSize,99);
2368 //     fLoader->TreeH()->GetListOfBranches()->Add(branch);
2369 //   }
2370 }
2371
2372 void AliTPC::SetTreeAddress()
2373 {
2374   //Sets tree address for hits
2375   if (fHitType<=1) {
2376     if (fHits == 0x0 ) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
2377     AliDetector::SetTreeAddress();
2378   }
2379   if (fHitType>1) SetTreeAddress2();
2380 }
2381
2382 void AliTPC::SetTreeAddress2()
2383 {
2384   //
2385   // Set branch address for the TrackHits Tree
2386   //
2387   AliDebug(1,"");
2388
2389   TBranch *branch;
2390   char branchname[20];
2391   //sprintf(branchname,"%s2",GetName());
2392   snprintf(branchname,20,"%s2",GetName());
2393   //
2394   // Branch address for hit tree
2395   TTree *treeH = fLoader->TreeH();
2396   if ((treeH)&&(fHitType&4)) {
2397     branch = treeH->GetBranch(branchname);
2398     if (branch) {
2399       branch->SetAddress(&fTrackHits);
2400       AliDebug(1,"fHitType&4 Setting");
2401     }
2402     else
2403       AliDebug(1,"fHitType&4 Failed (can not find branch)");
2404
2405   }
2406  //  if ((treeH)&&(fHitType&2)) {
2407 //     branch = treeH->GetBranch(branchname);
2408 //     if (branch) {
2409 //       branch->SetAddress(&fTrackHitsOld);
2410 //       AliDebug(1,"fHitType&2 Setting");
2411 //     }
2412 //     else
2413 //       AliDebug(1,"fHitType&2 Failed (can not find branch)");
2414 //   }
2415 }
2416
2417 void AliTPC::FinishPrimary()
2418 {
2419   if (fTrackHits &&fHitType&4)      fTrackHits->FlushHitStack();
2420   //  if (fTrackHitsOld && fHitType&2)  fTrackHitsOld->FlushHitStack();
2421 }
2422
2423
2424 void AliTPC::AddHit2(Int_t track, Int_t *vol, Float_t *hits)
2425 {
2426   //
2427   // add hit to the list
2428
2429   Int_t rtrack;
2430   if (fIshunt) {
2431     int primary = gAlice->GetMCApp()->GetPrimary(track);
2432     gAlice->GetMCApp()->Particle(primary)->SetBit(kKeepBit);
2433     rtrack=primary;
2434   } else {
2435     rtrack=track;
2436     gAlice->GetMCApp()->FlagTrack(track);
2437   }
2438   if (fTrackHits && fHitType&4)
2439     fTrackHits->AddHitKartez(vol[0],rtrack, hits[0],
2440                              hits[1],hits[2],(Int_t)hits[3],hits[4]);
2441  //  if (fTrackHitsOld &&fHitType&2 )
2442 //     fTrackHitsOld->AddHitKartez(vol[0],rtrack, hits[0],
2443 //                                 hits[1],hits[2],(Int_t)hits[3]);
2444
2445 }
2446
2447 void AliTPC::ResetHits()
2448 {
2449   if (fHitType&1) AliDetector::ResetHits();
2450   if (fHitType>1) ResetHits2();
2451 }
2452
2453 void AliTPC::ResetHits2()
2454 {
2455   //
2456   //reset hits
2457   if (fTrackHits && fHitType&4) fTrackHits->Clear();
2458   // if (fTrackHitsOld && fHitType&2) fTrackHitsOld->Clear();
2459
2460 }
2461
2462 AliHit* AliTPC::FirstHit(Int_t track)
2463 {
2464   if (fHitType>1) return FirstHit2(track);
2465   return AliDetector::FirstHit(track);
2466 }
2467 AliHit* AliTPC::NextHit()
2468 {
2469   //
2470   // gets next hit
2471   //
2472   if (fHitType>1) return NextHit2();
2473
2474   return AliDetector::NextHit();
2475 }
2476
2477 AliHit* AliTPC::FirstHit2(Int_t track)
2478 {
2479   //
2480   // Initialise the hit iterator
2481   // Return the address of the first hit for track
2482   // If track>=0 the track is read from disk
2483   // while if track<0 the first hit of the current
2484   // track is returned
2485   //
2486   if(track>=0) {
2487     gAlice->GetMCApp()->ResetHits();
2488     fLoader->TreeH()->GetEvent(track);
2489   }
2490   //
2491   if (fTrackHits && fHitType&4) {
2492     fTrackHits->First();
2493     return fTrackHits->GetHit();
2494   }
2495  //  if (fTrackHitsOld && fHitType&2) {
2496 //     fTrackHitsOld->First();
2497 //     return fTrackHitsOld->GetHit();
2498 //   }
2499
2500   else return 0;
2501 }
2502
2503 AliHit* AliTPC::NextHit2()
2504 {
2505   //
2506   //Return the next hit for the current track
2507
2508
2509 //   if (fTrackHitsOld && fHitType&2) {
2510 //     fTrackHitsOld->Next();
2511 //     return fTrackHitsOld->GetHit();
2512 //   }
2513   if (fTrackHits) {
2514     fTrackHits->Next();
2515     return fTrackHits->GetHit();
2516   }
2517   else
2518     return 0;
2519 }
2520
2521 void AliTPC::RemapTrackHitIDs(Int_t *map)
2522 {
2523   //
2524   // remapping
2525   //
2526   if (!fTrackHits) return;
2527
2528 //   if (fTrackHitsOld && fHitType&2){
2529 //     AliObjectArray * arr = fTrackHitsOld->fTrackHitsInfo;
2530 //     for (UInt_t i=0;i<arr->GetSize();i++){
2531 //       AliTrackHitsInfo * info = (AliTrackHitsInfo *)(arr->At(i));
2532 //       info->fTrackID = map[info->fTrackID];
2533 //     }
2534 //   }
2535 //  if (fTrackHitsOld && fHitType&4){
2536   if (fTrackHits && fHitType&4){
2537     TClonesArray * arr = fTrackHits->GetArray();;
2538     for (Int_t i=0;i<arr->GetEntriesFast();i++){
2539       AliTrackHitsParamV2 * info = (AliTrackHitsParamV2 *)(arr->At(i));
2540       info->SetTrackID(map[info->GetTrackID()]);
2541     }
2542   }
2543 }
2544
2545 Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
2546 {
2547   //return bool information - is track in given volume
2548   //load only part of the track information
2549   //return true if current track is in volume
2550   //
2551   //  return kTRUE;
2552  //  if (fTrackHitsOld && fHitType&2) {
2553 //     TBranch * br = fLoader->TreeH()->GetBranch("fTrackHitsInfo");
2554 //     br->GetEvent(track);
2555 //     AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
2556 //     for (UInt_t j=0;j<ar->GetSize();j++){
2557 //       if (  ((AliTrackHitsInfo*)ar->At(j))->fVolumeID==id) return kTRUE;
2558 //     }
2559 //   }
2560
2561   if (fTrackHits && fHitType&4) {
2562     TBranch * br1 = fLoader->TreeH()->GetBranch("fVolumes");
2563     TBranch * br2 = fLoader->TreeH()->GetBranch("fNVolumes");
2564     br2->GetEvent(track);
2565     br1->GetEvent(track);
2566     Int_t *volumes = fTrackHits->GetVolumes();
2567     Int_t nvolumes = fTrackHits->GetNVolumes();
2568     if (!volumes && nvolumes>0) {
2569       AliWarning(Form("Problematic track\t%d\t%d",track,nvolumes));
2570       return kFALSE;
2571     }
2572     for (Int_t j=0;j<nvolumes; j++)
2573       if (volumes[j]==id) return kTRUE;
2574   }
2575
2576   if (fHitType&1) {
2577     TBranch * br = fLoader->TreeH()->GetBranch("fSector");
2578     br->GetEvent(track);
2579     for (Int_t j=0;j<fHits->GetEntriesFast();j++){
2580       if (  ((AliTPChit*)fHits->At(j))->fSector==id) return kTRUE;
2581     }
2582   }
2583   return kFALSE;
2584
2585 }
2586
2587
2588 AliLoader* AliTPC::MakeLoader(const char* topfoldername)
2589 {
2590   //Makes TPC loader
2591   fLoader = new AliTPCLoader(GetName(),topfoldername);
2592   return fLoader;
2593 }
2594
2595 ////////////////////////////////////////////////////////////////////////
2596 AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
2597 //
2598 // load TPC paarmeters from a given file or create new if the object
2599 // is not found there
2600 // 12/05/2003 This method should be moved to the AliTPCLoader
2601 // and one has to decide where to store the TPC parameters
2602 // M.Kowalski
2603   char paramName[50];
2604   //sprintf(paramName,"75x40_100x60_150x60");
2605   snprintf(paramName,50,"75x40_100x60_150x60");
2606   AliTPCParam *paramTPC=(AliTPCParam*)file->Get(paramName);
2607   if (paramTPC) {
2608     AliDebugClass(1,Form("TPC parameters %s found.",paramName));
2609   } else {
2610     AliWarningClass("TPC parameters not found. Create new (they may be incorrect)");
2611     //paramTPC = new AliTPCParamSR;
2612     paramTPC = AliTPCcalibDB::Instance()->GetParameters();
2613     if (!paramTPC->IsGeoRead()){
2614       //
2615       // read transformation matrices for gGeoManager
2616       //
2617       paramTPC->ReadGeoMatrices();
2618     }
2619
2620   }
2621   return paramTPC;
2622
2623 // the older version of parameters can be accessed with this code.
2624 // In some cases, we have old parameters saved in the file but
2625 // digits were created with new parameters, it can be distinguish
2626 // by the name of TPC TreeD. The code here is just for the case
2627 // we would need to compare with old data, uncomment it if needed.
2628 //
2629 //  char paramName[50];
2630 //  sprintf(paramName,"75x40_100x60");
2631 //  AliTPCParam *paramTPC=(AliTPCParam*)in->Get(paramName);
2632 //  if (paramTPC) {
2633 //    cout<<"TPC parameters "<<paramName<<" found."<<endl;
2634 //  } else {
2635 //    sprintf(paramName,"75x40_100x60_150x60");
2636 //    paramTPC=(AliTPCParam*)in->Get(paramName);
2637 //    if (paramTPC) {
2638 //      cout<<"TPC parameters "<<paramName<<" found."<<endl;
2639 //    } else {
2640 //      cerr<<"TPC parameters not found. Create new (they may be incorrect)."
2641 //          <<endl;
2642 //      paramTPC = new AliTPCParamSR;
2643 //    }
2644 //  }
2645 //  return paramTPC;
2646
2647 }
2648
2649