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