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   if (gAlice){ // Set correctly the magnetic field in the ExB calculation
1853     AliMagF * field = gAlice->Field();
1854     if (field) {
1855       calib->SetExBField(field->SolenoidField());
1856     }
1857   }
1858
1859   Float_t gasgain = fTPCParam->GetGasGain();
1860   gasgain = gasgain/fGainFactor;
1861   Int_t i;
1862   Float_t xyz[5];
1863
1864   AliTPChit *tpcHit; // pointer to a sigle TPC hit
1865   //MI change
1866   TBranch * branch=0;
1867   if (fHitType>1) branch = TH->GetBranch("TPC2");
1868   else branch = TH->GetBranch("TPC");
1869
1870
1871   //----------------------------------------------
1872   // Create TObjArray-s, one for each row,
1873   // each TObjArray will store the TVectors
1874   // of electrons, one TVectors per each track.
1875   //----------------------------------------------
1876
1877   Int_t *nofElectrons = new Int_t [nrows+2]; // electron counter for each row
1878   TVector **tracks = new TVector* [nrows+2]; //pointers to the track vectors
1879
1880   for(i=0; i<nrows+2; i++){
1881     row[i] = new TObjArray;
1882     nofElectrons[i]=0;
1883     tracks[i]=0;
1884   }
1885
1886
1887
1888   //--------------------------------------------------------------------
1889   //  Loop over tracks, the "track" contains the full history
1890   //--------------------------------------------------------------------
1891
1892   Int_t previousTrack,currentTrack;
1893   previousTrack = -1; // nothing to store so far!
1894
1895   for(Int_t track=0;track<ntracks;track++){
1896     Bool_t isInSector=kTRUE;
1897     ResetHits();
1898     isInSector = TrackInVolume(isec,track);
1899     if (!isInSector) continue;
1900     //MI change
1901     branch->GetEntry(track); // get next track
1902
1903     //M.I. changes
1904
1905     tpcHit = (AliTPChit*)FirstHit(-1);
1906
1907     //--------------------------------------------------------------
1908     //  Loop over hits
1909     //--------------------------------------------------------------
1910
1911
1912     while(tpcHit){
1913
1914       Int_t sector=tpcHit->fSector; // sector number
1915       if(sector != isec){
1916         tpcHit = (AliTPChit*) NextHit();
1917         continue;
1918       }
1919
1920       // Remove hits which arrive before the TPC opening gate signal
1921       if(((fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
1922           /fTPCParam->GetDriftV()+tpcHit->Time())<fTPCParam->GetGateDelay()) {
1923         tpcHit = (AliTPChit*) NextHit();
1924         continue;
1925       }
1926
1927       currentTrack = tpcHit->Track(); // track number
1928
1929       if(currentTrack != previousTrack){
1930
1931         // store already filled fTrack
1932
1933         for(i=0;i<nrows+2;i++){
1934           if(previousTrack != -1){
1935             if(nofElectrons[i]>0){
1936               TVector &v = *tracks[i];
1937               v(0) = previousTrack;
1938               tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
1939               row[i]->Add(tracks[i]);
1940             }
1941             else {
1942               delete tracks[i]; // delete empty TVector
1943               tracks[i]=0;
1944             }
1945           }
1946
1947           nofElectrons[i]=0;
1948           tracks[i] = new TVector(601); // TVectors for the next fTrack
1949
1950         } // end of loop over rows
1951
1952         previousTrack=currentTrack; // update track label
1953       }
1954
1955       Int_t qI = (Int_t) (tpcHit->fQ); // energy loss (number of electrons)
1956
1957       //---------------------------------------------------
1958       //  Calculate the electron attachment probability
1959       //---------------------------------------------------
1960
1961
1962       Float_t time = 1.e6*(fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
1963         /fTPCParam->GetDriftV();
1964       // in microseconds!
1965       Float_t attProb = fTPCParam->GetAttCoef()*
1966         fTPCParam->GetOxyCont()*time; //  fraction!
1967
1968       //-----------------------------------------------
1969       //  Loop over electrons
1970       //-----------------------------------------------
1971       Int_t index[3];
1972       index[1]=isec;
1973       for(Int_t nel=0;nel<qI;nel++){
1974         // skip if electron lost due to the attachment
1975         if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
1976
1977         //
1978         // ExB effect
1979         //
1980         Double_t dxyz0[3],dxyz1[3];
1981         dxyz0[0]=tpcHit->X();
1982         dxyz0[1]=tpcHit->Y();
1983         dxyz0[2]=tpcHit->Z();
1984         if (calib->GetExB()){
1985           calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
1986         }else{
1987           AliError("Not valid ExB calibration");
1988           dxyz1[0]=tpcHit->X();
1989           dxyz1[1]=tpcHit->Y();
1990           dxyz1[2]=tpcHit->Z();
1991         }
1992         xyz[0]=dxyz1[0];
1993         xyz[1]=dxyz1[1];
1994         xyz[2]=dxyz1[2];
1995         //
1996         //
1997         //
1998         // protection for the nonphysical avalanche size (10**6 maximum)
1999         //
2000         Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
2001         xyz[3]= (Float_t) (-gasgain*TMath::Log(rn));
2002         index[0]=1;
2003
2004         TransportElectron(xyz,index);
2005         Int_t rowNumber;
2006         fTPCParam->GetPadRow(xyz,index);
2007         //
2008         // Add Time0 correction due unisochronity
2009         // xyz[0] - pad row coordinate
2010         // xyz[1] - pad coordinate
2011         // xyz[2] - is in now time bin coordinate system
2012         Float_t correction =0;
2013         if (calib->GetPadTime0()){
2014           if (!calib->GetPadTime0()->GetCalROC(isec)) continue;
2015           Int_t npads = fTPCParam->GetNPads(isec,TMath::Nint(xyz[0]));
2016           //      Int_t pad  = TMath::Nint(xyz[1]+fTPCParam->GetNPads(isec,TMath::Nint(xyz[0]))*0.5);
2017           Int_t pad  = TMath::Nint(xyz[1]);
2018           if (pad<0) pad=0;
2019           if (pad>=npads) pad=npads-1;
2020           correction = calib->GetPadTime0()->GetCalROC(isec)->GetValue(TMath::Nint(xyz[0]),pad);
2021
2022         }
2023         xyz[2]+=correction;
2024         //
2025         // Electron track time (for pileup simulation)
2026         xyz[4] = tpcHit->Time()/fTPCParam->GetTSample();
2027         // row 0 - cross talk from the innermost row
2028         // row fNRow+1 cross talk from the outermost row
2029         rowNumber = index[2]+1;
2030         //transform position to local digit coordinates
2031         //relative to nearest pad row
2032         if ((rowNumber<0)||rowNumber>fTPCParam->GetNRow(isec)+1) continue;
2033         /*      Float_t x1,y1;
2034         if (isec <fTPCParam->GetNInnerSector()) {
2035           x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
2036           y1 = fTPCParam->GetYInner(rowNumber);
2037         }
2038         else{
2039           x1=xyz[1]*fTPCParam->GetOuterPadPitchWidth();
2040           y1 = fTPCParam->GetYOuter(rowNumber);
2041         }
2042         // gain inefficiency at the wires edges - linear
2043         x1=TMath::Abs(x1);
2044         y1-=1.;
2045         if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.)); */
2046
2047         nofElectrons[rowNumber]++;
2048         //----------------------------------
2049         // Expand vector if necessary
2050         //----------------------------------
2051         if(nofElectrons[rowNumber]>120){
2052           Int_t range = tracks[rowNumber]->GetNrows();
2053           if((nofElectrons[rowNumber])>(range-1)/5){
2054
2055             tracks[rowNumber]->ResizeTo(range+500); // Add 100 electrons
2056           }
2057         }
2058
2059         TVector &v = *tracks[rowNumber];
2060         Int_t idx = 5*nofElectrons[rowNumber]-4;
2061         Real_t * position = &(((TVector&)v)(idx)); //make code faster
2062         memcpy(position,xyz,5*sizeof(Float_t));
2063
2064       } // end of loop over electrons
2065
2066       tpcHit = (AliTPChit*)NextHit();
2067
2068     } // end of loop over hits
2069   } // end of loop over tracks
2070
2071     //
2072     //   store remaining track (the last one) if not empty
2073     //
2074
2075   for(i=0;i<nrows+2;i++){
2076     if(nofElectrons[i]>0){
2077       TVector &v = *tracks[i];
2078       v(0) = previousTrack;
2079       tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
2080       row[i]->Add(tracks[i]);
2081     }
2082     else{
2083       delete tracks[i];
2084       tracks[i]=0;
2085     }
2086   }
2087
2088   delete [] tracks;
2089   delete [] nofElectrons;
2090
2091 } // end of MakeSector
2092
2093
2094 //_____________________________________________________________________________
2095 void AliTPC::Init()
2096 {
2097   //
2098   // Initialise TPC detector after definition of geometry
2099   //
2100   AliDebug(1,"*********************************************");
2101 }
2102
2103 //_____________________________________________________________________________
2104 void AliTPC::MakeBranch(Option_t* option)
2105 {
2106   //
2107   // Create Tree branches for the TPC.
2108   //
2109   AliDebug(1,"");
2110   Int_t buffersize = 4000;
2111   char branchname[10];
2112   sprintf(branchname,"%s",GetName());
2113
2114   const char *h = strstr(option,"H");
2115
2116   if ( h && (fHitType<=1) && (fHits == 0x0)) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
2117
2118   AliDetector::MakeBranch(option);
2119
2120   const char *d = strstr(option,"D");
2121
2122   if (fDigits   && fLoader->TreeD() && d) {
2123     MakeBranchInTree(gAlice->TreeD(), branchname, &fDigits, buffersize, 0);
2124   }
2125
2126   if (fHitType>1) MakeBranch2(option,0); // MI change 14.09.2000
2127 }
2128
2129 //_____________________________________________________________________________
2130 void AliTPC::ResetDigits()
2131 {
2132   //
2133   // Reset number of digits and the digits array for this detector
2134   //
2135   fNdigits   = 0;
2136   if (fDigits)   fDigits->Clear();
2137 }
2138
2139
2140
2141 //_____________________________________________________________________________
2142 void AliTPC::SetSens(Int_t sens)
2143 {
2144
2145   //-------------------------------------------------------------
2146   // Activates/deactivates the sensitive strips at the center of
2147   // the pad row -- this is for the space-point resolution calculations
2148   //-------------------------------------------------------------
2149
2150   //-----------------------------------------------------------------
2151   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
2152   //-----------------------------------------------------------------
2153
2154   fSens = sens;
2155 }
2156
2157
2158 void AliTPC::SetSide(Float_t side=0.)
2159 {
2160   // choice of the TPC side
2161
2162   fSide = side;
2163
2164 }
2165 //_____________________________________________________________________________
2166
2167 void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
2168 {
2169   //
2170   // electron transport taking into account:
2171   // 1. diffusion,
2172   // 2.ExB at the wires
2173   // 3. nonisochronity
2174   //
2175   // xyz and index must be already transformed to system 1
2176   //
2177
2178   fTPCParam->Transform1to2(xyz,index);
2179
2180   //add diffusion
2181   Float_t driftl=xyz[2];
2182   if(driftl<0.01) driftl=0.01;
2183   driftl=TMath::Sqrt(driftl);
2184   Float_t sigT = driftl*(fTPCParam->GetDiffT());
2185   Float_t sigL = driftl*(fTPCParam->GetDiffL());
2186   xyz[0]=gRandom->Gaus(xyz[0],sigT);
2187   xyz[1]=gRandom->Gaus(xyz[1],sigT);
2188   xyz[2]=gRandom->Gaus(xyz[2],sigL);
2189
2190   // ExB
2191
2192   if (fTPCParam->GetMWPCReadout()==kTRUE){
2193     Float_t dx = fTPCParam->Transform2to2NearestWire(xyz,index);
2194     xyz[1]+=dx*(fTPCParam->GetOmegaTau());
2195   }
2196   //add nonisochronity (not implemented yet)
2197
2198
2199 }
2200
2201 ClassImp(AliTPChit)
2202   //______________________________________________________________________
2203   AliTPChit::AliTPChit()
2204             :AliHit(),
2205              fSector(0),
2206              fPadRow(0),
2207              fQ(0),
2208              fTime(0)
2209 {
2210   //
2211   // default
2212   //
2213
2214 }
2215 //_____________________________________________________________________________
2216 AliTPChit::AliTPChit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits)
2217           :AliHit(shunt,track),
2218              fSector(0),
2219              fPadRow(0),
2220              fQ(0),
2221              fTime(0)
2222 {
2223   //
2224   // Creates a TPC hit object
2225   //
2226   fSector     = vol[0];
2227   fPadRow     = vol[1];
2228   fX          = hits[0];
2229   fY          = hits[1];
2230   fZ          = hits[2];
2231   fQ          = hits[3];
2232   fTime       = hits[4];
2233 }
2234
2235 //________________________________________________________________________
2236 // Additional code because of the AliTPCTrackHitsV2
2237
2238 void AliTPC::MakeBranch2(Option_t *option,const char */*file*/)
2239 {
2240   //
2241   // Create a new branch in the current Root Tree
2242   // The branch of fHits is automatically split
2243   // MI change 14.09.2000
2244   AliDebug(1,"");
2245   if (fHitType<2) return;
2246   char branchname[10];
2247   sprintf(branchname,"%s2",GetName());
2248   //
2249   // Get the pointer to the header
2250   const char *cH = strstr(option,"H");
2251   //
2252   if (fTrackHits   && TreeH() && cH && fHitType&4) {
2253     AliDebug(1,"Making branch for Type 4 Hits");
2254     TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
2255   }
2256
2257 //   if (fTrackHitsOld   && TreeH() && cH && fHitType&2) {
2258 //     AliDebug(1,"Making branch for Type 2 Hits");
2259 //     AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHitsOld,
2260 //                                                    TreeH(),fBufferSize,99);
2261 //     TreeH()->GetListOfBranches()->Add(branch);
2262 //   }
2263 }
2264
2265 void AliTPC::SetTreeAddress()
2266 {
2267   //Sets tree address for hits
2268   if (fHitType<=1) {
2269     if (fHits == 0x0 ) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
2270     AliDetector::SetTreeAddress();
2271   }
2272   if (fHitType>1) SetTreeAddress2();
2273 }
2274
2275 void AliTPC::SetTreeAddress2()
2276 {
2277   //
2278   // Set branch address for the TrackHits Tree
2279   //
2280   AliDebug(1,"");
2281
2282   TBranch *branch;
2283   char branchname[20];
2284   sprintf(branchname,"%s2",GetName());
2285   //
2286   // Branch address for hit tree
2287   TTree *treeH = TreeH();
2288   if ((treeH)&&(fHitType&4)) {
2289     branch = treeH->GetBranch(branchname);
2290     if (branch) {
2291       branch->SetAddress(&fTrackHits);
2292       AliDebug(1,"fHitType&4 Setting");
2293     }
2294     else
2295       AliDebug(1,"fHitType&4 Failed (can not find branch)");
2296
2297   }
2298  //  if ((treeH)&&(fHitType&2)) {
2299 //     branch = treeH->GetBranch(branchname);
2300 //     if (branch) {
2301 //       branch->SetAddress(&fTrackHitsOld);
2302 //       AliDebug(1,"fHitType&2 Setting");
2303 //     }
2304 //     else
2305 //       AliDebug(1,"fHitType&2 Failed (can not find branch)");
2306 //   }
2307 }
2308
2309 void AliTPC::FinishPrimary()
2310 {
2311   if (fTrackHits &&fHitType&4)      fTrackHits->FlushHitStack();
2312   //  if (fTrackHitsOld && fHitType&2)  fTrackHitsOld->FlushHitStack();
2313 }
2314
2315
2316 void AliTPC::AddHit2(Int_t track, Int_t *vol, Float_t *hits)
2317 {
2318   //
2319   // add hit to the list
2320   Int_t rtrack;
2321   if (fIshunt) {
2322     int primary = gAlice->GetMCApp()->GetPrimary(track);
2323     gAlice->GetMCApp()->Particle(primary)->SetBit(kKeepBit);
2324     rtrack=primary;
2325   } else {
2326     rtrack=track;
2327     gAlice->GetMCApp()->FlagTrack(track);
2328   }
2329   if (fTrackHits && fHitType&4)
2330     fTrackHits->AddHitKartez(vol[0],rtrack, hits[0],
2331                              hits[1],hits[2],(Int_t)hits[3],hits[4]);
2332  //  if (fTrackHitsOld &&fHitType&2 )
2333 //     fTrackHitsOld->AddHitKartez(vol[0],rtrack, hits[0],
2334 //                                 hits[1],hits[2],(Int_t)hits[3]);
2335
2336 }
2337
2338 void AliTPC::ResetHits()
2339 {
2340   if (fHitType&1) AliDetector::ResetHits();
2341   if (fHitType>1) ResetHits2();
2342 }
2343
2344 void AliTPC::ResetHits2()
2345 {
2346   //
2347   //reset hits
2348   if (fTrackHits && fHitType&4) fTrackHits->Clear();
2349   // if (fTrackHitsOld && fHitType&2) fTrackHitsOld->Clear();
2350
2351 }
2352
2353 AliHit* AliTPC::FirstHit(Int_t track)
2354 {
2355   if (fHitType>1) return FirstHit2(track);
2356   return AliDetector::FirstHit(track);
2357 }
2358 AliHit* AliTPC::NextHit()
2359 {
2360   //
2361   // gets next hit
2362   //
2363   if (fHitType>1) return NextHit2();
2364
2365   return AliDetector::NextHit();
2366 }
2367
2368 AliHit* AliTPC::FirstHit2(Int_t track)
2369 {
2370   //
2371   // Initialise the hit iterator
2372   // Return the address of the first hit for track
2373   // If track>=0 the track is read from disk
2374   // while if track<0 the first hit of the current
2375   // track is returned
2376   //
2377   if(track>=0) {
2378     gAlice->ResetHits();
2379     TreeH()->GetEvent(track);
2380   }
2381   //
2382   if (fTrackHits && fHitType&4) {
2383     fTrackHits->First();
2384     return fTrackHits->GetHit();
2385   }
2386  //  if (fTrackHitsOld && fHitType&2) {
2387 //     fTrackHitsOld->First();
2388 //     return fTrackHitsOld->GetHit();
2389 //   }
2390
2391   else return 0;
2392 }
2393
2394 AliHit* AliTPC::NextHit2()
2395 {
2396   //
2397   //Return the next hit for the current track
2398
2399
2400 //   if (fTrackHitsOld && fHitType&2) {
2401 //     fTrackHitsOld->Next();
2402 //     return fTrackHitsOld->GetHit();
2403 //   }
2404   if (fTrackHits) {
2405     fTrackHits->Next();
2406     return fTrackHits->GetHit();
2407   }
2408   else
2409     return 0;
2410 }
2411
2412 void AliTPC::LoadPoints(Int_t)
2413 {
2414   //
2415   Int_t a = 0;
2416
2417   if(fHitType==1) AliDetector::LoadPoints(a);
2418   else LoadPoints2(a);
2419 }
2420
2421
2422 void AliTPC::RemapTrackHitIDs(Int_t *map)
2423 {
2424   //
2425   // remapping
2426   //
2427   if (!fTrackHits) return;
2428
2429 //   if (fTrackHitsOld && fHitType&2){
2430 //     AliObjectArray * arr = fTrackHitsOld->fTrackHitsInfo;
2431 //     for (UInt_t i=0;i<arr->GetSize();i++){
2432 //       AliTrackHitsInfo * info = (AliTrackHitsInfo *)(arr->At(i));
2433 //       info->fTrackID = map[info->fTrackID];
2434 //     }
2435 //   }
2436 //  if (fTrackHitsOld && fHitType&4){
2437   if (fTrackHits && fHitType&4){
2438     TClonesArray * arr = fTrackHits->GetArray();;
2439     for (Int_t i=0;i<arr->GetEntriesFast();i++){
2440       AliTrackHitsParamV2 * info = (AliTrackHitsParamV2 *)(arr->At(i));
2441       info->SetTrackID(map[info->GetTrackID()]);
2442     }
2443   }
2444 }
2445
2446 Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
2447 {
2448   //return bool information - is track in given volume
2449   //load only part of the track information
2450   //return true if current track is in volume
2451   //
2452   //  return kTRUE;
2453  //  if (fTrackHitsOld && fHitType&2) {
2454 //     TBranch * br = TreeH()->GetBranch("fTrackHitsInfo");
2455 //     br->GetEvent(track);
2456 //     AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
2457 //     for (UInt_t j=0;j<ar->GetSize();j++){
2458 //       if (  ((AliTrackHitsInfo*)ar->At(j))->fVolumeID==id) return kTRUE;
2459 //     }
2460 //   }
2461
2462   if (fTrackHits && fHitType&4) {
2463     TBranch * br1 = TreeH()->GetBranch("fVolumes");
2464     TBranch * br2 = TreeH()->GetBranch("fNVolumes");
2465     br2->GetEvent(track);
2466     br1->GetEvent(track);
2467     Int_t *volumes = fTrackHits->GetVolumes();
2468     Int_t nvolumes = fTrackHits->GetNVolumes();
2469     if (!volumes && nvolumes>0) {
2470       AliWarning(Form("Problematic track\t%d\t%d",track,nvolumes));
2471       return kFALSE;
2472     }
2473     for (Int_t j=0;j<nvolumes; j++)
2474       if (volumes[j]==id) return kTRUE;
2475   }
2476
2477   if (fHitType&1) {
2478     TBranch * br = TreeH()->GetBranch("fSector");
2479     br->GetEvent(track);
2480     for (Int_t j=0;j<fHits->GetEntriesFast();j++){
2481       if (  ((AliTPChit*)fHits->At(j))->fSector==id) return kTRUE;
2482     }
2483   }
2484   return kFALSE;
2485
2486 }
2487
2488 //_____________________________________________________________________________
2489 void AliTPC::LoadPoints2(Int_t)
2490 {
2491   //
2492   // Store x, y, z of all hits in memory
2493   //
2494   //  if (fTrackHits == 0 && fTrackHitsOld==0) return;
2495   if (fTrackHits == 0 ) return;
2496   //
2497   Int_t nhits =0;
2498   if (fHitType&4) nhits = fTrackHits->GetEntriesFast();
2499   //  if (fHitType&2) nhits = fTrackHitsOld->GetEntriesFast();
2500
2501   if (nhits == 0) return;
2502   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
2503   if (fPoints == 0) fPoints = new TObjArray(tracks);
2504   AliHit *ahit;
2505   //
2506   Int_t *ntrk=new Int_t[tracks];
2507   Int_t *limi=new Int_t[tracks];
2508   Float_t **coor=new Float_t*[tracks];
2509   for(Int_t i=0;i<tracks;i++) {
2510     ntrk[i]=0;
2511     coor[i]=0;
2512     limi[i]=0;
2513   }
2514   //
2515   AliPoints *points = 0;
2516   Float_t *fp=0;
2517   Int_t trk;
2518   Int_t chunk=nhits/4+1;
2519   //
2520   // Loop over all the hits and store their position
2521   //
2522   ahit = FirstHit2(-1);
2523   while (ahit){
2524     trk=ahit->GetTrack();
2525     if(ntrk[trk]==limi[trk]) {
2526       //
2527       // Initialise a new track
2528       fp=new Float_t[3*(limi[trk]+chunk)];
2529       if(coor[trk]) {
2530         memcpy(fp,coor[trk],sizeof(Float_t)*3*limi[trk]);
2531         delete [] coor[trk];
2532       }
2533       limi[trk]+=chunk;
2534       coor[trk] = fp;
2535     } else {
2536       fp = coor[trk];
2537     }
2538     fp[3*ntrk[trk]  ] = ahit->X();
2539     fp[3*ntrk[trk]+1] = ahit->Y();
2540     fp[3*ntrk[trk]+2] = ahit->Z();
2541     ntrk[trk]++;
2542     ahit = NextHit2();
2543   }
2544
2545
2546
2547   //
2548   for(trk=0; trk<tracks; ++trk) {
2549     if(ntrk[trk]) {
2550       points = new AliPoints();
2551       points->SetMarkerColor(kYellow); //PH kYellow it the default in TPC
2552       points->SetMarkerSize(1);//PH Default size=1
2553       points->SetDetector(this);
2554       points->SetParticle(trk);
2555       points->SetPolyMarker(ntrk[trk],coor[trk],1); // Default style=1
2556       fPoints->AddAt(points,trk);
2557       delete [] coor[trk];
2558       coor[trk]=0;
2559     }
2560   }
2561   delete [] coor;
2562   delete [] ntrk;
2563   delete [] limi;
2564 }
2565
2566
2567 //_____________________________________________________________________________
2568 void AliTPC::LoadPoints3(Int_t)
2569 {
2570   //
2571   // Store x, y, z of all hits in memory
2572   // - only intersection point with pad row
2573   if (fTrackHits == 0) return;
2574   //
2575   Int_t nhits = fTrackHits->GetEntriesFast();
2576   if (nhits == 0) return;
2577   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
2578   if (fPoints == 0) fPoints = new TObjArray(2*tracks);
2579   fPoints->Expand(2*tracks);
2580   AliHit *ahit;
2581   //
2582   Int_t *ntrk=new Int_t[tracks];
2583   Int_t *limi=new Int_t[tracks];
2584   Float_t **coor=new Float_t*[tracks];
2585   for(Int_t i=0;i<tracks;i++) {
2586     ntrk[i]=0;
2587     coor[i]=0;
2588     limi[i]=0;
2589   }
2590   //
2591   AliPoints *points = 0;
2592   Float_t *fp=0;
2593   Int_t trk;
2594   Int_t chunk=nhits/4+1;
2595   //
2596   // Loop over all the hits and store their position
2597   //
2598   ahit = FirstHit2(-1);
2599
2600   Int_t lastrow = -1;
2601   while (ahit){
2602     trk=ahit->GetTrack();
2603     Float_t  x[3]={ahit->X(),ahit->Y(),ahit->Z()};
2604     Int_t    index[3]={1,((AliTPChit*)ahit)->fSector,0};
2605     Int_t    currentrow = fTPCParam->GetPadRow(x,index) ;
2606     if (currentrow!=lastrow){
2607       lastrow = currentrow;
2608       //later calculate intersection point
2609       if(ntrk[trk]==limi[trk]) {
2610         //
2611         // Initialise a new track
2612         fp=new Float_t[3*(limi[trk]+chunk)];
2613         if(coor[trk]) {
2614           memcpy(fp,coor[trk],sizeof(Float_t)*3*limi[trk]);
2615           delete [] coor[trk];
2616         }
2617         limi[trk]+=chunk;
2618         coor[trk] = fp;
2619       } else {
2620         fp = coor[trk];
2621       }
2622       fp[3*ntrk[trk]  ] = ahit->X();
2623       fp[3*ntrk[trk]+1] = ahit->Y();
2624       fp[3*ntrk[trk]+2] = ahit->Z();
2625       ntrk[trk]++;
2626     }
2627     ahit = NextHit2();
2628   }
2629
2630   //
2631   for(trk=0; trk<tracks; ++trk) {
2632     if(ntrk[trk]) {
2633       points = new AliPoints();
2634       points->SetMarkerColor(kMagenta); //PH kYellow + 1 is kMagenta
2635       points->SetMarkerStyle(5);
2636       points->SetMarkerSize(0.2);
2637       points->SetDetector(this);
2638       points->SetParticle(trk);
2639       points->SetPolyMarker(ntrk[trk],coor[trk],30);
2640       fPoints->AddAt(points,tracks+trk);
2641       delete [] coor[trk];
2642       coor[trk]=0;
2643     }
2644   }
2645   delete [] coor;
2646   delete [] ntrk;
2647   delete [] limi;
2648 }
2649
2650
2651
2652 AliLoader* AliTPC::MakeLoader(const char* topfoldername)
2653 {
2654   //Makes TPC loader
2655   fLoader = new AliTPCLoader(GetName(),topfoldername);
2656   return fLoader;
2657 }
2658
2659 ////////////////////////////////////////////////////////////////////////
2660 AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
2661 //
2662 // load TPC paarmeters from a given file or create new if the object
2663 // is not found there
2664 // 12/05/2003 This method should be moved to the AliTPCLoader
2665 // and one has to decide where to store the TPC parameters
2666 // M.Kowalski
2667   char paramName[50];
2668   sprintf(paramName,"75x40_100x60_150x60");
2669   AliTPCParam *paramTPC=(AliTPCParam*)file->Get(paramName);
2670   if (paramTPC) {
2671     AliDebugClass(1,Form("TPC parameters %s found.",paramName));
2672   } else {
2673     AliWarningClass("TPC parameters not found. Create new (they may be incorrect)");
2674     paramTPC = new AliTPCParamSR;
2675   }
2676   return paramTPC;
2677
2678 // the older version of parameters can be accessed with this code.
2679 // In some cases, we have old parameters saved in the file but
2680 // digits were created with new parameters, it can be distinguish
2681 // by the name of TPC TreeD. The code here is just for the case
2682 // we would need to compare with old data, uncomment it if needed.
2683 //
2684 //  char paramName[50];
2685 //  sprintf(paramName,"75x40_100x60");
2686 //  AliTPCParam *paramTPC=(AliTPCParam*)in->Get(paramName);
2687 //  if (paramTPC) {
2688 //    cout<<"TPC parameters "<<paramName<<" found."<<endl;
2689 //  } else {
2690 //    sprintf(paramName,"75x40_100x60_150x60");
2691 //    paramTPC=(AliTPCParam*)in->Get(paramName);
2692 //    if (paramTPC) {
2693 //      cout<<"TPC parameters "<<paramName<<" found."<<endl;
2694 //    } else {
2695 //      cerr<<"TPC parameters not found. Create new (they may be incorrect)."
2696 //          <<endl;
2697 //      paramTPC = new AliTPCParamSR;
2698 //    }
2699 //  }
2700 //  return paramTPC;
2701
2702 }
2703
2704