STEER/ESD/AliESDtrack.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 //           Implementation of the ESD track class
  17 //   ESD = Event Summary Data
  18 //   This is the class to deal with during the phisics analysis of data
  19 //      Origin: Iouri Belikov, CERN
  20 //      e-mail: Jouri.Belikov@cern.ch
  21 //
  22 //
  23 //
  24 //  What do you need to know before starting analysis
  25 //  (by Marian Ivanov: marian.ivanov@cern.ch)
  26 //
  27 //
  28 //   AliESDtrack:
  29 //   1.  What is the AliESDtrack
  30 //   2.  What informations do we store
  31 //   3.  How to use the information for analysis
  32 //
  33 //
  34 //   1.AliESDtrack is the container of the information about the track/particle
  35 //     reconstructed during Barrel Tracking.
  36 //     Content:
  37 //        a.) Track parameters and covariance - AliExternalTrackParam  - snapshots along trajectory at differnt tracking steps
  38 //            current, fIp, fTPCinner, fCp, fOp, fHMPIDp, + friendTrack (fITSout, fTPCout, fTRDin)
  39 //        b.) Flags - per detector status bits
  40 //        c.) Track fit quality -chi2, number of clusters
  41 //        d.) Detector PID information
  42 //        d.) Different detector specific information (e.g number of tracklets in TRD, TOF cluster descriptor ...)
  43 //
  44 //
  45 //     The track information is propagated from one tracking detector to
  46 //     other using the functionality of AliESDtrack (AliExternalTrackParam - current parameters)
  47 //
  48 //     Barrel tracking uses Kalman filtering technique (no global fit model is used). Kalman provides optimal local
  49 //     track parameters estimate at given position under certian assumptions.
  50 //     Following approximations were used:
  51 //           a.) gaussian Multiple scattering
  52 //           b.) gaussian energy loss
  53 //           c.) gaussian error of the space point residuals
  54 //
  55 //     Kalman filter take into account following effects which are
  56 //     difficult to handle using global fit:
  57 //        a.) Multiple scattering
  58 //        b.) Energy loss
  59 //        c.) Non homogenous magnetic field
  60 //
  61 //     In general case, following barrel detectors are contributing to
  62 //     the Kalman track information:
  63 //         a. TPC
  64 //         b. ITS
  65 //         c. TRD
  66 //
  67 //      Track findind/fitting procedure is done in 3 steps:
  68 //         1. Cluster2Track(in)   - inward sequence TPC->ITS
  69 //         2. PropagateBack(out)   - outward sequence ITS->TPC->TRD -> Outer PID detectors
  70 //         3. RefitInward(refit)     - inward sequence TRD->TPC->ITS
  71 //      After each recosntruction step detector status is updated in the data member fFlags
  72 //      fFlags|=k<DetectorName><step> where step={1:in,2:out,3:refit,}
  73 //      For some of detectors a special flags were implemented. Declaration of list of all flags can be find in $ALICE_ROOT/STEER/STEERBase/AliVTrack.h
  74 //
  75 //
  76 //      The current track parameter is updated after each detector (see bellow).
  77 //      In specical cases a track  snapshots (AliExternalTrackParam) are stored
  78 //
  79 //
  80 //
  81 //      For some type of analysis (+visualization) track local parameters at
  82 //      different position are neccesary. A snapshots during the track
  83 //      propagation are created and stored either in track itself (for analysis purposes) or assiciated friend track (for calibration and debugging purposes)
  84 //      (See AliExternalTrackParam class for desctiption of variables and functionality)
  85 //      Snapshots:
  86 //      a.)  Current parameters (AliESDtrack itself)
  87 //         Description:  as obtained in the last succesfull tracking step
  88 //         Contributors: best case TRD->TPC->ITS after RefitInward
  89 //         Recomended way to get track parameters. It includes all of the information.
  90 //         NOTICE - By default the track parameters are stored at the DCA point to the primary vertex.
  91 //                  Optimal for primary tracks, far from optimal for deeply secondary tracks.
  92 //                  To get optimal track parameters at the secondary vertex, OnTheFly V0s with associated track parameters should be used
  93 //
  94 //      b.) Constrained parameters (fCp)
  95 //         Description:
  96 //            Kalman track  updated with the Primary vertex information with corresponding error (soft  constraint - see http://en.wikipedia.org/wiki/Constraint_(mathematics)#Hard_and_soft_constraints)
  97 //         Function:
  98 //            const AliExternalTrackParam *GetConstrainedParam() const {return fCp;}
  99 //         Contributors: best case TRD->TPC->ITS after RefitInward
 100 //         Recommended usage: Use only for tracks selected as primary (check GetConstrainedChi2())
 101 //         NOTICE - not real constraint only soft constraint
 102 //
 103 //      c.) Inner parameters (fIp) -
 104 //         Description: Track parameters at inner wall of the TPC
 105 //         Function:
 106 //           const AliExternalTrackParam *GetInnerParam() const { return fIp;}
 107 //         Contributors: general case TRD->TPC (during RefitInward)
 108 //         Recomended usage: To provide momenta for the TPC PID and to estimate quality of the track determination for further
 109 //
 110 //      d.)  TPCinnerParam (fTPCinner):
 111 //         Description: TPC only parameters at DCA to the primary vertex (corrected for the material between TPC and vertex)
 112 //         Function:
 113 //                const AliExternalTrackParam *GetTPCInnerParam() const {return fTPCInner;}
 114 //         Contributors:  TPC only from in step  1 (Cluster2Tracks)
 115 //         Recomended usage: Requested for HBT study
 116 //                           (smaller correlations as using also ITS information)
 117 //         NOTICE: Optimal for primary, far from optimal for secondary tracks (similar to track parameters a.)
 118 //                 ! We should always use the c.) fIp in case of the TPC PID analysis, or undo material budget correction!
 119 //
 120 //      e.) Outer parameters - (fOp)
 121 //         Description: track parameters during PropagateBack in the last sucessfull propagation
 122 //                  Reason to generate backup  OuterParameters
 123 //                        a.) Local inclination angle bigger than threshold -
 124 //                            Low momenta tracks
 125 //                        b.) Catastrofic (large relative>~20 %)energy loss in material outside of the TPC
 126 //                        c.) No additional  space points contributing to track
 127 //         Function:
 128 //            const AliExternalTrackParam *GetOuterParam() const { return fOp;}
 129 //         Contributors:  general case - ITS-> TPC -> TRD ( during PropagateBack )
 130 //         Recomended usage:
 131 //              a.) Tracking: Starting parameter for Refit inward
 132 //              b.) Visualization
 133 //              c.) QA
 134 //         NOTICE: Should be not used for the physic analysis
 135 //         Function:
 136 //            const AliExternalTrackParam *GetOuterParam() const { return fOp;}
 137 //
 138 //-----------------------------------------------------------------
 139
 140 #include <TMath.h>
 141 #include <TParticle.h>
 142 #include <TDatabasePDG.h>
 143 #include <TMatrixD.h>
 144
 145 #include "AliESDVertex.h"
 146 #include "AliESDtrack.h"
 147 #include "AliESDEvent.h"
 148 #include "AliKalmanTrack.h"
 149 #include "AliVTrack.h"
 150 #include "AliLog.h"
 151 #include "AliTrackPointArray.h"
 152 #include "TPolyMarker3D.h"
 153 #include "AliTrackerBase.h"
 154 #include "AliTPCdEdxInfo.h"
 155 #include "AliDetectorPID.h"
 156 #include "TTreeStream.h"
 157 #include "TObjArray.h"
 158
 159 ClassImp(AliESDtrack)
 160
 161 void SetPIDValues(Double_t * dest, const Double_t * src, Int_t n) {
 162   // This function copies "n" PID weights from "scr" to "dest"
 163   // and normalizes their sum to 1 thus producing conditional probabilities.
 164   // The negative weights are set to 0.
 165   // In case all the weights are non-positive they are replaced by
 166   // uniform probabilities
 167
 168   if (n<=0) return;
 169
 170   Float_t uniform = 1./(Float_t)n;
 171
 172   Float_t sum = 0;
 173   for (Int_t i=0; i<n; i++)
 174     if (src[i]>=0) {
 175       sum+=src[i];
 176       dest[i] = src[i];
 177     }
 178     else {
 179       dest[i] = 0;
 180     }
 181
 182   if(sum>0)
 183     for (Int_t i=0; i<n; i++) dest[i] /= sum;
 184   else
 185     for (Int_t i=0; i<n; i++) dest[i] = uniform;
 186 }
 187
 188 //_______________________________________________________________________
 189 AliESDtrack::AliESDtrack() :
 190   AliExternalTrackParam(),
 191   fCp(0),
 192   fIp(0),
 193   fTPCInner(0),
 194   fOp(0),
 195   fHMPIDp(0),
 196   fFriendTrack(NULL),
 197   fTPCFitMap(159),//number of padrows
 198   fTPCClusterMap(159),//number of padrows
 199   fTPCSharedMap(159),//number of padrows
 200   fFlags(0),
 201   fID(0),
 202   fLabel(0),
 203   fITSLabel(0),
 204   fTPCLabel(0),
 205   fTRDLabel(0),
 206   fTOFLabel(NULL),
 207   fTOFCalChannel(-1),
 208   fTOFindex(-1),
 209   fHMPIDqn(0),
 210   fHMPIDcluIdx(-1),
 211   fCaloIndex(kEMCALNoMatch),
 212   fR(0),
 213   fITSr(0),
 214   fTPCr(0),
 215   fTRDr(0),
 216   fTOFr(0),
 217   fHMPIDr(0),
 218   fHMPIDtrkTheta(0),
 219   fHMPIDtrkPhi(0),
 220   fHMPIDsignal(0),
 221   fTrackTime(0),
 222   fTrackLength(0),
 223   fdTPC(0),fzTPC(0),
 224   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 225   fCchi2TPC(0),
 226   fD(0),fZ(0),
 227   fCdd(0),fCdz(0),fCzz(0),
 228   fCchi2(0),
 229   fITSchi2(0),
 230   fTPCchi2(0),
 231   fTPCchi2Iter1(0),
 232   fTRDchi2(0),
 233   fTOFchi2(0),
 234   fHMPIDchi2(0),
 235   fGlobalChi2(0),
 236   fITSsignal(0),
 237   fTPCsignal(0),
 238   fTPCsignalTuned(0),
 239   fTPCsignalS(0),
 240   fTPCdEdxInfo(0),
 241   fTRDsignal(0),
 242   fTRDQuality(0),
 243   fTRDBudget(0),
 244   fTOFsignal(99999),
 245   fTOFsignalTuned(99999),
 246   fTOFsignalToT(99999),
 247   fTOFsignalRaw(99999),
 248   fTOFsignalDz(999),
 249   fTOFsignalDx(999),
 250   fTOFdeltaBC(999),
 251   fTOFl0l1(999),
 252   fCaloDx(0),
 253   fCaloDz(0),
 254   fHMPIDtrkX(0),
 255   fHMPIDtrkY(0),
 256   fHMPIDmipX(0),
 257   fHMPIDmipY(0),
 258   fTPCncls(0),
 259   fTPCnclsF(0),
 260   fTPCsignalN(0),
 261   fTPCnclsIter1(0),
 262   fTPCnclsFIter1(0),
 263   fITSncls(0),
 264   fITSClusterMap(0),
 265   fITSSharedMap(0),
 266   fTRDncls(0),
 267   fTRDncls0(0),
 268   fTRDntracklets(0),
 269   fTRDNchamberdEdx(0),
 270   fTRDNclusterdEdx(0),
 271   fTRDnSlices(0),
 272   fTRDslices(0x0),
 273   fVertexID(-2),// -2 means an orphan track
 274   fPIDForTracking(AliPID::kPion),
 275   fESDEvent(0),
 276   fCacheNCrossedRows(-10),
 277   fCacheChi2TPCConstrainedVsGlobal(-10),
 278   fCacheChi2TPCConstrainedVsGlobalVertex(0),
 279   fDetectorPID(0x0),
 280   fTrackPhiOnEMCal(-999),
 281   fTrackEtaOnEMCal(-999),
 282   fTrackPtOnEMCal(-999),
 283   fNtofClusters(0),
 284   fTOFcluster(NULL)
 285 {
 286   //
 287   // The default ESD constructor
 288   //
 289   if (!OnlineMode()) fFriendTrack=new AliESDfriendTrack();
 290
 291   Int_t i;
 292   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
 293
 294   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 295   for (i=0; i<3; i++)   { fV0Indexes[i]=0;}
 296   for (i=0;i<kTRDnPlanes;i++) {
 297     fTRDTimBin[i]=0;
 298   }
 299   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
 300   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 301   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 302   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 303 }
 304
 305 bool AliESDtrack::fgkOnlineMode=false;
 306
 307 //_______________________________________________________________________
 308 AliESDtrack::AliESDtrack(const AliESDtrack& track):
 309   AliExternalTrackParam(track),
 310   fCp(0),
 311   fIp(0),
 312   fTPCInner(0),
 313   fOp(0),
 314   fHMPIDp(0),
 315   fFriendTrack(0),
 316   fTPCFitMap(track.fTPCFitMap),
 317   fTPCClusterMap(track.fTPCClusterMap),
 318   fTPCSharedMap(track.fTPCSharedMap),
 319   fFlags(track.fFlags),
 320   fID(track.fID),
 321   fLabel(track.fLabel),
 322   fITSLabel(track.fITSLabel),
 323   fTPCLabel(track.fTPCLabel),
 324   fTRDLabel(track.fTRDLabel),
 325   fTOFLabel(NULL),
 326   fTOFCalChannel(track.fTOFCalChannel),
 327   fTOFindex(track.fTOFindex),
 328   fHMPIDqn(track.fHMPIDqn),
 329   fHMPIDcluIdx(track.fHMPIDcluIdx),
 330   fCaloIndex(track.fCaloIndex),
 331   fR(0),
 332   fITSr(0),
 333   fTPCr(0),
 334   fTRDr(0),
 335   fTOFr(0),
 336   fHMPIDr(0),
 337   fHMPIDtrkTheta(track.fHMPIDtrkTheta),
 338   fHMPIDtrkPhi(track.fHMPIDtrkPhi),
 339   fHMPIDsignal(track.fHMPIDsignal),
 340   fTrackTime(NULL),
 341   fTrackLength(track.fTrackLength),
 342   fdTPC(track.fdTPC),fzTPC(track.fzTPC),
 343   fCddTPC(track.fCddTPC),fCdzTPC(track.fCdzTPC),fCzzTPC(track.fCzzTPC),
 344   fCchi2TPC(track.fCchi2TPC),
 345   fD(track.fD),fZ(track.fZ),
 346   fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
 347   fCchi2(track.fCchi2),
 348   fITSchi2(track.fITSchi2),
 349   fTPCchi2(track.fTPCchi2),
 350   fTPCchi2Iter1(track.fTPCchi2Iter1),
 351   fTRDchi2(track.fTRDchi2),
 352   fTOFchi2(track.fTOFchi2),
 353   fHMPIDchi2(track.fHMPIDchi2),
 354   fGlobalChi2(track.fGlobalChi2),
 355   fITSsignal(track.fITSsignal),
 356   fTPCsignal(track.fTPCsignal),
 357   fTPCsignalTuned(track.fTPCsignalTuned),
 358   fTPCsignalS(track.fTPCsignalS),
 359   fTPCdEdxInfo(0),
 360   fTRDsignal(track.fTRDsignal),
 361   fTRDQuality(track.fTRDQuality),
 362   fTRDBudget(track.fTRDBudget),
 363   fTOFsignal(track.fTOFsignal),
 364   fTOFsignalTuned(track.fTOFsignalTuned),
 365   fTOFsignalToT(track.fTOFsignalToT),
 366   fTOFsignalRaw(track.fTOFsignalRaw),
 367   fTOFsignalDz(track.fTOFsignalDz),
 368   fTOFsignalDx(track.fTOFsignalDx),
 369   fTOFdeltaBC(track.fTOFdeltaBC),
 370   fTOFl0l1(track.fTOFl0l1),
 371   fCaloDx(track.fCaloDx),
 372   fCaloDz(track.fCaloDz),
 373   fHMPIDtrkX(track.fHMPIDtrkX),
 374   fHMPIDtrkY(track.fHMPIDtrkY),
 375   fHMPIDmipX(track.fHMPIDmipX),
 376   fHMPIDmipY(track.fHMPIDmipY),
 377   fTPCncls(track.fTPCncls),
 378   fTPCnclsF(track.fTPCnclsF),
 379   fTPCsignalN(track.fTPCsignalN),
 380   fTPCnclsIter1(track.fTPCnclsIter1),
 381   fTPCnclsFIter1(track.fTPCnclsIter1),
 382   fITSncls(track.fITSncls),
 383   fITSClusterMap(track.fITSClusterMap),
 384   fITSSharedMap(track.fITSSharedMap),
 385   fTRDncls(track.fTRDncls),
 386   fTRDncls0(track.fTRDncls0),
 387   fTRDntracklets(track.fTRDntracklets),
 388   fTRDNchamberdEdx(track.fTRDNchamberdEdx),
 389   fTRDNclusterdEdx(track.fTRDNclusterdEdx),
 390   fTRDnSlices(track.fTRDnSlices),
 391   fTRDslices(0x0),
 392   fVertexID(track.fVertexID),
 393   fPIDForTracking(AliPID::kPion),
 394   fESDEvent(track.fESDEvent),
 395   fCacheNCrossedRows(track.fCacheNCrossedRows),
 396   fCacheChi2TPCConstrainedVsGlobal(track.fCacheChi2TPCConstrainedVsGlobal),
 397   fCacheChi2TPCConstrainedVsGlobalVertex(track.fCacheChi2TPCConstrainedVsGlobalVertex),
 398   fDetectorPID(0x0),
 399   fTrackPhiOnEMCal(track.fTrackPhiOnEMCal),
 400   fTrackEtaOnEMCal(track.fTrackEtaOnEMCal),
 401   fTrackPtOnEMCal(track.fTrackPtOnEMCal),
 402   fNtofClusters(track.fNtofClusters),
 403   fTOFcluster(NULL)
 404 {
 405   //
 406   //copy constructor
 407   //
 408   for (Int_t i=kNITSchi2Std;i--;) fITSchi2Std[i] = track.fITSchi2Std[i];
 409
 410   if(track.fTrackTime){
 411     fTrackTime = new Double32_t[AliPID::kSPECIESC];
 412     for (Int_t i=0;i<AliPID::kSPECIESC;i++) fTrackTime[i]=track.fTrackTime[i];
 413   }
 414
 415   if (track.fR) {
 416     fR = new Double32_t[AliPID::kSPECIES];
 417     for (Int_t i=AliPID::kSPECIES;i--;)  fR[i]=track.fR[i];
 418   }
 419   if (track.fITSr) {
 420     fITSr = new Double32_t[AliPID::kSPECIES];
 421     for (Int_t i=AliPID::kSPECIES;i--;)  fITSr[i]=track.fITSr[i];
 422   }
 423   //
 424   if (track.fTPCr) {
 425     fTPCr = new Double32_t[AliPID::kSPECIES];
 426     for (Int_t i=AliPID::kSPECIES;i--;) fTPCr[i]=track.fTPCr[i];
 427   }
 428
 429   for (Int_t i=0;i<4;i++) {fITSdEdxSamples[i]=track.fITSdEdxSamples[i];}
 430   for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
 431   for (Int_t i=0; i<3;i++)   { fKinkIndexes[i]=track.fKinkIndexes[i];}
 432   for (Int_t i=0; i<3;i++)   { fV0Indexes[i]=track.fV0Indexes[i];}
 433   //
 434   for (Int_t i=0;i<kTRDnPlanes;i++) {
 435     fTRDTimBin[i]=track.fTRDTimBin[i];
 436   }
 437
 438   if (fTRDnSlices) {
 439     fTRDslices=new Double32_t[fTRDnSlices];
 440     for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=track.fTRDslices[i];
 441   }
 442
 443   if (track.fDetectorPID) fDetectorPID = new AliDetectorPID(*track.fDetectorPID);
 444
 445   if (track.fTRDr) {
 446     fTRDr  = new Double32_t[AliPID::kSPECIES];
 447     for (Int_t i=AliPID::kSPECIES;i--;) fTRDr[i]=track.fTRDr[i];
 448   }
 449
 450   if (track.fTOFr) {
 451     fTOFr = new Double32_t[AliPID::kSPECIES];
 452     for (Int_t i=AliPID::kSPECIES;i--;) fTOFr[i]=track.fTOFr[i];
 453   }
 454
 455   if(track.fTOFLabel){
 456     if(!fTOFLabel) fTOFLabel = new Int_t[3];
 457     for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
 458   }
 459
 460   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
 461   for (Int_t i=0;i<12;i++) fITSModule[i]=track.fITSModule[i];
 462
 463   if (track.fHMPIDr) {
 464     fHMPIDr = new Double32_t[AliPID::kSPECIES];
 465     for (Int_t i=AliPID::kSPECIES;i--;) fHMPIDr[i]=track.fHMPIDr[i];
 466   }
 467
 468   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
 469   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
 470   if (track.fTPCInner) fTPCInner=new AliExternalTrackParam(*track.fTPCInner);
 471   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
 472   if (track.fHMPIDp) fHMPIDp=new AliExternalTrackParam(*track.fHMPIDp);
 473   if (track.fTPCdEdxInfo) fTPCdEdxInfo = new AliTPCdEdxInfo(*track.fTPCdEdxInfo);
 474
 475
 476   if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
 477
 478   if(fNtofClusters > 0){
 479     fTOFcluster = new Int_t[fNtofClusters];
 480         for(Int_t i=0;i < fNtofClusters;i++) fTOFcluster[i] = track.fTOFcluster[i];
 481   }
 482 }
 483
 484 //_______________________________________________________________________
 485 AliESDtrack::AliESDtrack(const AliVTrack *track) :
 486   AliExternalTrackParam(track),
 487   fCp(0),
 488   fIp(0),
 489   fTPCInner(0),
 490   fOp(0),
 491   fHMPIDp(0),
 492   fFriendTrack(0),
 493   fTPCFitMap(159),//number of padrows
 494   fTPCClusterMap(159),//number of padrows
 495   fTPCSharedMap(159),//number of padrows
 496   fFlags(0),
 497   fID(),
 498   fLabel(0),
 499   fITSLabel(0),
 500   fTPCLabel(0),
 501   fTRDLabel(0),
 502   fTOFLabel(NULL),
 503   fTOFCalChannel(-1),
 504   fTOFindex(-1),
 505   fHMPIDqn(0),
 506   fHMPIDcluIdx(-1),
 507   fCaloIndex(kEMCALNoMatch),
 508   fR(0),
 509   fITSr(0),
 510   fTPCr(0),
 511   fTRDr(0),
 512   fTOFr(0),
 513   fHMPIDr(0),
 514   fHMPIDtrkTheta(0),
 515   fHMPIDtrkPhi(0),
 516   fHMPIDsignal(0),
 517   fTrackTime(NULL),
 518   fTrackLength(0),
 519   fdTPC(0),fzTPC(0),
 520   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 521   fCchi2TPC(0),
 522   fD(0),fZ(0),
 523   fCdd(0),fCdz(0),fCzz(0),
 524   fCchi2(0),
 525   fITSchi2(0),
 526   fTPCchi2(0),
 527   fTPCchi2Iter1(0),
 528   fTRDchi2(0),
 529   fTOFchi2(0),
 530   fHMPIDchi2(0),
 531   fGlobalChi2(0),
 532   fITSsignal(0),
 533   fTPCsignal(0),
 534   fTPCsignalTuned(0),
 535   fTPCsignalS(0),
 536   fTPCdEdxInfo(0),
 537   fTRDsignal(0),
 538   fTRDQuality(0),
 539   fTRDBudget(0),
 540   fTOFsignal(99999),
 541   fTOFsignalTuned(99999),
 542   fTOFsignalToT(99999),
 543   fTOFsignalRaw(99999),
 544   fTOFsignalDz(999),
 545   fTOFsignalDx(999),
 546   fTOFdeltaBC(999),
 547   fTOFl0l1(999),
 548   fCaloDx(0),
 549   fCaloDz(0),
 550   fHMPIDtrkX(0),
 551   fHMPIDtrkY(0),
 552   fHMPIDmipX(0),
 553   fHMPIDmipY(0),
 554   fTPCncls(0),
 555   fTPCnclsF(0),
 556   fTPCsignalN(0),
 557   fTPCnclsIter1(0),
 558   fTPCnclsFIter1(0),
 559   fITSncls(0),
 560   fITSClusterMap(0),
 561   fITSSharedMap(0),
 562   fTRDncls(0),
 563   fTRDncls0(0),
 564   fTRDntracklets(0),
 565   fTRDNchamberdEdx(0),
 566   fTRDNclusterdEdx(0),
 567   fTRDnSlices(0),
 568   fTRDslices(0x0),
 569   fVertexID(-2),  // -2 means an orphan track
 570   fPIDForTracking(track->GetPIDForTracking()),
 571   fESDEvent(0),
 572   fCacheNCrossedRows(-10),
 573   fCacheChi2TPCConstrainedVsGlobal(-10),
 574   fCacheChi2TPCConstrainedVsGlobalVertex(0),
 575   fDetectorPID(0x0),
 576   fTrackPhiOnEMCal(-999),
 577   fTrackEtaOnEMCal(-999),
 578   fTrackPtOnEMCal(-999),
 579   fNtofClusters(0),
 580   fTOFcluster(NULL)
 581 {
 582   //
 583   // ESD track from AliVTrack.
 584   // This is not a copy constructor !
 585   //
 586
 587   if (track->InheritsFrom("AliExternalTrackParam")) {
 588      AliError("This is not a copy constructor. Use AliESDtrack(const AliESDtrack &) !");
 589      AliWarning("Calling the default constructor...");
 590      AliESDtrack();
 591      return;
 592   }
 593
 594   // Reset all the arrays
 595   Int_t i;
 596   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
 597
 598   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 599   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
 600   for (i=0;i<kTRDnPlanes;i++) {
 601     fTRDTimBin[i]=0;
 602   }
 603   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
 604   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 605   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 606   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 607
 608
 609   // Set ITS cluster map
 610   fITSClusterMap=track->GetITSClusterMap();
 611   fITSSharedMap=0;
 612
 613   fITSncls=0;
 614   for(i=0; i<6; i++) {
 615     if(HasPointOnITSLayer(i)) fITSncls++;
 616   }
 617
 618   // Set TPC ncls
 619   fTPCncls=track->GetTPCNcls();
 620   fTPCnclsF=track->GetTPCNclsF();
 621   // TPC cluster maps
 622   const TBits* bmap = track->GetTPCClusterMapPtr();
 623   if (bmap) SetTPCClusterMap(*bmap);
 624   bmap = GetTPCFitMapPtr();
 625   if (bmap) SetTPCFitMap(*bmap);
 626   bmap = GetTPCSharedMapPtr();
 627   if (bmap) SetTPCSharedMap(*bmap);
 628   //
 629   // Set the combined PID
 630   const Double_t *pid = track->PID();
 631   if(pid) {
 632     fR = new Double32_t[AliPID::kSPECIES];
 633     for (i=AliPID::kSPECIES; i--;) fR[i]=pid[i];
 634   }
 635   //
 636   // calo matched cluster id
 637   SetEMCALcluster(track->GetEMCALcluster());
 638   // AliESD track label
 639   //
 640   // PID info
 641   fITSsignal = track->GetITSsignal();
 642   double itsdEdx[4];
 643   track->GetITSdEdxSamples(itsdEdx);
 644   SetITSdEdxSamples(itsdEdx);
 645   //
 646   SetTPCsignal(track->GetTPCsignal(),fTPCsignalS,track->GetTPCsignalN()); // No signalS in AODPi
 647   AliTPCdEdxInfo * dEdxInfo = track->GetTPCdEdxInfo();
 648   if (dEdxInfo) SetTPCdEdxInfo(new AliTPCdEdxInfo(*dEdxInfo));
 649   //
 650   SetTRDsignal(track->GetTRDsignal());
 651   int ntrdsl = track->GetNumberOfTRDslices();
 652   if (ntrdsl>0) {
 653     SetNumberOfTRDslices((ntrdsl+2)*kTRDnPlanes);
 654     for (int ipl=kTRDnPlanes;ipl--;){
 655       for (int isl=ntrdsl;isl--;) SetTRDslice(track->GetTRDslice(ipl,isl),ipl,isl);
 656       Double_t sp, p = track->GetTRDmomentum(ipl, &sp);
 657       SetTRDmomentum(p, ipl, &sp);
 658     }
 659   }
 660   //
 661   fTRDncls = track->GetTRDncls();
 662   fTRDntracklets &= 0xff & track->GetTRDntrackletsPID();
 663   fTRDchi2 = track->GetTRDchi2();
 664   //
 665   SetTOFsignal(track->GetTOFsignal());
 666   Double_t expt[AliPID::kSPECIESC];
 667   track->GetIntegratedTimes(expt,AliPID::kSPECIESC);
 668   SetIntegratedTimes(expt);
 669   //
 670   SetTrackPhiEtaPtOnEMCal(track->GetTrackPhiOnEMCal(),track->GetTrackEtaOnEMCal(),track->GetTrackPtOnEMCal());
 671   //
 672   SetLabel(track->GetLabel());
 673   // Set the status
 674   SetStatus(track->GetStatus());
 675   //
 676   // Set the ID
 677   SetID(track->GetID());
 678   //
 679 }
 680
 681 //_______________________________________________________________________
 682 AliESDtrack::AliESDtrack(TParticle * part) :
 683   AliExternalTrackParam(),
 684   fCp(0),
 685   fIp(0),
 686   fTPCInner(0),
 687   fOp(0),
 688   fHMPIDp(0),
 689   fFriendTrack(0),
 690   fTPCFitMap(159),//number of padrows
 691   fTPCClusterMap(159),//number of padrows
 692   fTPCSharedMap(159),//number of padrows
 693   fFlags(0),
 694   fID(0),
 695   fLabel(0),
 696   fITSLabel(0),
 697   fTPCLabel(0),
 698   fTRDLabel(0),
 699   fTOFLabel(NULL),
 700   fTOFCalChannel(-1),
 701   fTOFindex(-1),
 702   fHMPIDqn(0),
 703   fHMPIDcluIdx(-1),
 704   fCaloIndex(kEMCALNoMatch),
 705   fR(0),
 706   fITSr(0),
 707   fTPCr(0),
 708   fTRDr(0),
 709   fTOFr(0),
 710   fHMPIDr(0),
 711   fHMPIDtrkTheta(0),
 712   fHMPIDtrkPhi(0),
 713   fHMPIDsignal(0),
 714   fTrackTime(NULL),
 715   fTrackLength(0),
 716   fdTPC(0),fzTPC(0),
 717   fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
 718   fCchi2TPC(0),
 719   fD(0),fZ(0),
 720   fCdd(0),fCdz(0),fCzz(0),
 721   fCchi2(0),
 722   fITSchi2(0),
 723   fTPCchi2(0),
 724   fTPCchi2Iter1(0),
 725   fTRDchi2(0),
 726   fTOFchi2(0),
 727   fHMPIDchi2(0),
 728   fGlobalChi2(0),
 729   fITSsignal(0),
 730   fTPCsignal(0),
 731   fTPCsignalTuned(0),
 732   fTPCsignalS(0),
 733   fTPCdEdxInfo(0),
 734   fTRDsignal(0),
 735   fTRDQuality(0),
 736   fTRDBudget(0),
 737   fTOFsignal(99999),
 738   fTOFsignalTuned(99999),
 739   fTOFsignalToT(99999),
 740   fTOFsignalRaw(99999),
 741   fTOFsignalDz(999),
 742   fTOFsignalDx(999),
 743   fTOFdeltaBC(999),
 744   fTOFl0l1(999),
 745   fCaloDx(0),
 746   fCaloDz(0),
 747   fHMPIDtrkX(0),
 748   fHMPIDtrkY(0),
 749   fHMPIDmipX(0),
 750   fHMPIDmipY(0),
 751   fTPCncls(0),
 752   fTPCnclsF(0),
 753   fTPCsignalN(0),
 754   fTPCnclsIter1(0),
 755   fTPCnclsFIter1(0),
 756   fITSncls(0),
 757   fITSClusterMap(0),
 758   fITSSharedMap(0),
 759   fTRDncls(0),
 760   fTRDncls0(0),
 761   fTRDntracklets(0),
 762   fTRDNchamberdEdx(0),
 763   fTRDNclusterdEdx(0),
 764   fTRDnSlices(0),
 765   fTRDslices(0x0),
 766   fVertexID(-2),  // -2 means an orphan track
 767   fPIDForTracking(AliPID::kPion),
 768   fESDEvent(0),
 769   fCacheNCrossedRows(-10),
 770   fCacheChi2TPCConstrainedVsGlobal(-10),
 771   fCacheChi2TPCConstrainedVsGlobalVertex(0),
 772   fDetectorPID(0x0),
 773   fTrackPhiOnEMCal(-999),
 774   fTrackEtaOnEMCal(-999),
 775   fTrackPtOnEMCal(-999),
 776   fNtofClusters(0),
 777   fTOFcluster(NULL)
 778 {
 779   //
 780   // ESD track from TParticle
 781   //
 782
 783   // Reset all the arrays
 784   Int_t i;
 785   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
 786
 787   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
 788   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
 789   for (i=0;i<kTRDnPlanes;i++) {
 790     fTRDTimBin[i]=0;
 791   }
 792   for (i=0;i<4;i++) {fITSdEdxSamples[i]=0.;}
 793   for (i=0;i<4;i++) {fTPCPoints[i]=0;}
 794   for (i=0;i<10;i++) {fTOFInfo[i]=0;}
 795   for (i=0;i<12;i++) {fITSModule[i]=-1;}
 796
 797   // Calculate the AliExternalTrackParam content
 798
 799   Double_t xref;
 800   Double_t alpha;
 801   Double_t param[5];
 802   Double_t covar[15];
 803
 804   // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
 805   alpha = part->Phi()*180./TMath::Pi();
 806   if (alpha<0) alpha+= 360.;
 807   if (alpha>360) alpha -= 360.;
 808
 809   Int_t sector = (Int_t)(alpha/20.);
 810   alpha = 10. + 20.*sector;
 811   alpha /= 180;
 812   alpha *= TMath::Pi();
 813
 814   // Covariance matrix: no errors, the parameters are exact
 815   for (i=0; i<15; i++) covar[i]=0.;
 816
 817   // Get the vertex of origin and the momentum
 818   TVector3 ver(part->Vx(),part->Vy(),part->Vz());
 819   TVector3 mom(part->Px(),part->Py(),part->Pz());
 820
 821   // Rotate to the local coordinate system (TPC sector)
 822   ver.RotateZ(-alpha);
 823   mom.RotateZ(-alpha);
 824
 825   // X of the referense plane
 826   xref = ver.X();
 827
 828   Int_t pdgCode = part->GetPdgCode();
 829
 830   Double_t charge =
 831     TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
 832
 833   param[0] = ver.Y();
 834   param[1] = ver.Z();
 835   param[2] = TMath::Sin(mom.Phi());
 836   param[3] = mom.Pz()/mom.Pt();
 837   param[4] = TMath::Sign(1/mom.Pt(),charge);
 838
 839   // Set AliExternalTrackParam
 840   Set(xref, alpha, param, covar);
 841
 842   // Set the PID
 843   Int_t indexPID = 99;
 844   if (pdgCode<0) pdgCode = -pdgCode;
 845   for (i=0;i<AliPID::kSPECIESC;i++) if (pdgCode==AliPID::ParticleCode(i)) {indexPID = i; break;}
 846
 847   if (indexPID < AliPID::kSPECIESC) fPIDForTracking = indexPID;
 848
 849   // AliESD track label
 850   SetLabel(part->GetUniqueID());
 851
 852 }
 853
 854 //_______________________________________________________________________
 855 AliESDtrack::~AliESDtrack(){
 856   //
 857   // This is destructor according Coding Conventrions
 858   //
 859   //printf("Delete track\n");
 860   delete fIp;
 861   delete fTPCInner;
 862   delete fOp;
 863   delete fHMPIDp;
 864   delete fCp;
 865   delete fFriendTrack;
 866   delete fTPCdEdxInfo;
 867   if(fTRDnSlices)
 868     delete[] fTRDslices;
 869
 870   //Reset cached values - needed for TClonesArray in AliESDInputHandler
 871   fCacheNCrossedRows = -10.;
 872   fCacheChi2TPCConstrainedVsGlobal = -10.;
 873   if(fCacheChi2TPCConstrainedVsGlobalVertex) fCacheChi2TPCConstrainedVsGlobalVertex = 0;
 874
 875   if(fTOFcluster)
 876     delete[] fTOFcluster;
 877   fTOFcluster = NULL;
 878   fNtofClusters=0;
 879
 880   delete fDetectorPID;
 881
 882   delete[] fR;
 883   delete[] fITSr;
 884   delete[] fTPCr;
 885   delete[] fTRDr;
 886   delete[] fTOFr;
 887   delete[] fHMPIDr;
 888   //
 889   if(fTrackTime) delete[] fTrackTime;
 890   if(fTOFLabel) delete[] fTOFLabel;
 891 }
 892
 893 //_______________________________________________________________________
 894 AliESDtrack &AliESDtrack::operator=(const AliESDtrack &source)
 895 {
 896   // == operator
 897
 898   if(&source == this) return *this;
 899   AliExternalTrackParam::operator=(source);
 900
 901
 902   if(source.fCp){
 903     // we have the trackparam: assign or copy construct
 904     if(fCp)*fCp = *source.fCp;
 905     else fCp = new AliExternalTrackParam(*source.fCp);
 906   }
 907   else{
 908     // no track param delete the old one
 909     delete fCp;
 910     fCp = 0;
 911   }
 912
 913   if(source.fIp){
 914     // we have the trackparam: assign or copy construct
 915     if(fIp)*fIp = *source.fIp;
 916     else fIp = new AliExternalTrackParam(*source.fIp);
 917   }
 918   else{
 919     // no track param delete the old one
 920     delete fIp;
 921     fIp = 0;
 922   }
 923
 924
 925   if(source.fTPCInner){
 926     // we have the trackparam: assign or copy construct
 927     if(fTPCInner) *fTPCInner = *source.fTPCInner;
 928     else fTPCInner = new AliExternalTrackParam(*source.fTPCInner);
 929   }
 930   else{
 931     // no track param delete the old one
 932     delete fTPCInner;
 933     fTPCInner = 0;
 934   }
 935
 936   if(source.fTPCdEdxInfo) {
 937     if(fTPCdEdxInfo) *fTPCdEdxInfo = *source.fTPCdEdxInfo;
 938     fTPCdEdxInfo = new AliTPCdEdxInfo(*source.fTPCdEdxInfo);
 939   }
 940
 941   if(source.fOp){
 942     // we have the trackparam: assign or copy construct
 943     if(fOp) *fOp = *source.fOp;
 944     else fOp = new AliExternalTrackParam(*source.fOp);
 945   }
 946   else{
 947     // no track param delete the old one
 948     delete fOp;
 949     fOp = 0;
 950   }
 951
 952
 953   if(source.fHMPIDp){
 954     // we have the trackparam: assign or copy construct
 955     if(fHMPIDp) *fHMPIDp = *source.fHMPIDp;
 956     else fHMPIDp = new AliExternalTrackParam(*source.fHMPIDp);
 957   }
 958   else{
 959     // no track param delete the old one
 960     delete fHMPIDp;
 961     fHMPIDp = 0;
 962   }
 963
 964   // copy also the friend track
 965   // use copy constructor
 966   if(source.fFriendTrack){
 967     // we have the trackparam: assign or copy construct
 968     delete fFriendTrack; fFriendTrack=new AliESDfriendTrack(*source.fFriendTrack);
 969   }
 970   else{
 971     // no track param delete the old one
 972     delete fFriendTrack; fFriendTrack= 0;
 973   }
 974
 975   fTPCFitMap = source.fTPCFitMap;
 976   fTPCClusterMap = source.fTPCClusterMap;
 977   fTPCSharedMap  = source.fTPCSharedMap;
 978   // the simple stuff
 979   fFlags    = source.fFlags;
 980   fID       = source.fID;
 981   fLabel    = source.fLabel;
 982   fITSLabel = source.fITSLabel;
 983   for(int i = 0; i< 12;++i){
 984     fITSModule[i] = source.fITSModule[i];
 985   }
 986   fTPCLabel = source.fTPCLabel;
 987   fTRDLabel = source.fTRDLabel;
 988   if(source.fTOFLabel){
 989     if(!fTOFLabel) fTOFLabel = new Int_t[3];
 990     for(int i = 0; i< 3;++i){
 991       fTOFLabel[i] = source.fTOFLabel[i];
 992     }
 993   }
 994   fTOFCalChannel = source.fTOFCalChannel;
 995   fTOFindex      = source.fTOFindex;
 996   fHMPIDqn       = source.fHMPIDqn;
 997   fHMPIDcluIdx   = source.fHMPIDcluIdx;
 998   fCaloIndex    = source.fCaloIndex;
 999   for (int i=kNITSchi2Std;i--;) fITSchi2Std[i] = source.fITSchi2Std[i];
1000   for(int i = 0; i< 3;++i){
1001     fKinkIndexes[i] = source.fKinkIndexes[i];
1002     fV0Indexes[i]   = source.fV0Indexes[i];
1003   }
1004
1005   if (source.fR) {
1006     if (!fR) fR = new Double32_t[AliPID::kSPECIES];
1007     for (Int_t i=AliPID::kSPECIES;i--;)  fR[i]=source.fR[i];
1008   }
1009   else {delete[] fR; fR = 0;}
1010
1011   if (source.fITSr) {
1012     if (!fITSr) fITSr = new Double32_t[AliPID::kSPECIES];
1013     for (Int_t i=AliPID::kSPECIES;i--;)  fITSr[i]=source.fITSr[i];
1014   }
1015   else {delete[] fITSr; fITSr = 0;}
1016   //
1017   if (source.fTPCr) {
1018     if (!fTPCr) fTPCr = new Double32_t[AliPID::kSPECIES];
1019     for (Int_t i=AliPID::kSPECIES;i--;) fTPCr[i]=source.fTPCr[i];
1020   }
1021   else {delete[] fTPCr; fTPCr = 0;}
1022
1023   if (source.fTRDr) {
1024     if (!fTRDr) fTRDr  = new Double32_t[AliPID::kSPECIES];
1025     for (Int_t i=AliPID::kSPECIES;i--;) fTRDr[i]=source.fTRDr[i];
1026   }
1027   else {delete[] fTRDr; fTRDr = 0;}
1028
1029   if (source.fTOFr) {
1030     if (!fTOFr) fTOFr = new Double32_t[AliPID::kSPECIES];
1031     for (Int_t i=AliPID::kSPECIES;i--;) fTOFr[i]=source.fTOFr[i];
1032   }
1033   else {delete[] fTOFr; fTOFr = 0;}
1034
1035   if (source.fHMPIDr) {
1036     if (!fHMPIDr) fHMPIDr = new Double32_t[AliPID::kSPECIES];
1037     for (Int_t i=AliPID::kSPECIES;i--;) fHMPIDr[i]=source.fHMPIDr[i];
1038   }
1039   else {delete[] fHMPIDr; fHMPIDr = 0;}
1040
1041   fPIDForTracking = source.fPIDForTracking;
1042
1043   fHMPIDtrkTheta = source.fHMPIDtrkTheta;
1044   fHMPIDtrkPhi   = source.fHMPIDtrkPhi;
1045   fHMPIDsignal   = source.fHMPIDsignal;
1046
1047
1048   if(fTrackTime){
1049     delete[] fTrackTime;
1050   }
1051   if(source.fTrackTime){
1052     fTrackTime = new Double32_t[AliPID::kSPECIESC];
1053     for(Int_t i=0;i < AliPID::kSPECIESC;i++)
1054       fTrackTime[i] = source.fTrackTime[i];
1055   }
1056
1057   fTrackLength   = source. fTrackLength;
1058   fdTPC  = source.fdTPC;
1059   fzTPC  = source.fzTPC;
1060   fCddTPC = source.fCddTPC;
1061   fCdzTPC = source.fCdzTPC;
1062   fCzzTPC = source.fCzzTPC;
1063   fCchi2TPC = source.fCchi2TPC;
1064
1065   fD  = source.fD;
1066   fZ  = source.fZ;
1067   fCdd = source.fCdd;
1068   fCdz = source.fCdz;
1069   fCzz = source.fCzz;
1070   fCchi2     = source.fCchi2;
1071
1072   fITSchi2   = source.fITSchi2;
1073   fTPCchi2   = source.fTPCchi2;
1074   fTPCchi2Iter1   = source.fTPCchi2Iter1;
1075   fTRDchi2   = source.fTRDchi2;
1076   fTOFchi2   = source.fTOFchi2;
1077   fHMPIDchi2 = source.fHMPIDchi2;
1078
1079   fGlobalChi2 = source.fGlobalChi2;
1080
1081   fITSsignal  = source.fITSsignal;
1082   for (Int_t i=0;i<4;i++) {fITSdEdxSamples[i]=source.fITSdEdxSamples[i];}
1083   fTPCsignal  = source.fTPCsignal;
1084   fTPCsignalTuned  = source.fTPCsignalTuned;
1085   fTPCsignalS = source.fTPCsignalS;
1086   for(int i = 0; i< 4;++i){
1087     fTPCPoints[i] = source.fTPCPoints[i];
1088   }
1089   fTRDsignal = source.fTRDsignal;
1090   fTRDNchamberdEdx = source.fTRDNchamberdEdx;
1091   fTRDNclusterdEdx = source.fTRDNclusterdEdx;
1092
1093   for(int i = 0;i < kTRDnPlanes;++i){
1094     fTRDTimBin[i] = source.fTRDTimBin[i];
1095   }
1096
1097   if(fTRDnSlices)
1098     delete[] fTRDslices;
1099   fTRDslices=0;
1100   fTRDnSlices=source.fTRDnSlices;
1101   if (fTRDnSlices) {
1102     fTRDslices=new Double32_t[fTRDnSlices];
1103     for(int j = 0;j < fTRDnSlices;++j) fTRDslices[j] = source.fTRDslices[j];
1104   }
1105
1106   fTRDQuality =   source.fTRDQuality;
1107   fTRDBudget  =   source.fTRDBudget;
1108   fTOFsignal  =   source.fTOFsignal;
1109   fTOFsignalTuned  = source.fTOFsignalTuned;
1110   fTOFsignalToT = source.fTOFsignalToT;
1111   fTOFsignalRaw = source.fTOFsignalRaw;
1112   fTOFsignalDz  = source.fTOFsignalDz;
1113   fTOFsignalDx  = source.fTOFsignalDx;
1114   fTOFdeltaBC   = source.fTOFdeltaBC;
1115   fTOFl0l1      = source.fTOFl0l1;
1116
1117   for(int i = 0;i<10;++i){
1118     fTOFInfo[i] = source.fTOFInfo[i];
1119   }
1120
1121   fHMPIDtrkX = source.fHMPIDtrkX;
1122   fHMPIDtrkY = source.fHMPIDtrkY;
1123   fHMPIDmipX = source.fHMPIDmipX;
1124   fHMPIDmipY = source.fHMPIDmipY;
1125
1126   fTPCncls    = source.fTPCncls;
1127   fTPCnclsF   = source.fTPCnclsF;
1128   fTPCsignalN = source.fTPCsignalN;
1129   fTPCnclsIter1    = source.fTPCnclsIter1;
1130   fTPCnclsFIter1   = source.fTPCnclsFIter1;
1131
1132   fITSncls = source.fITSncls;
1133   fITSClusterMap = source.fITSClusterMap;
1134   fITSSharedMap = source.fITSSharedMap;
1135   fTRDncls   = source.fTRDncls;
1136   fTRDncls0  = source.fTRDncls0;
1137   fTRDntracklets  = source.fTRDntracklets;
1138   fVertexID = source.fVertexID;
1139   fPIDForTracking = source.fPIDForTracking;
1140
1141   fCacheNCrossedRows = source.fCacheNCrossedRows;
1142   fCacheChi2TPCConstrainedVsGlobal = source.fCacheChi2TPCConstrainedVsGlobal;
1143   fCacheChi2TPCConstrainedVsGlobalVertex = source.fCacheChi2TPCConstrainedVsGlobalVertex;
1144
1145   delete fDetectorPID;
1146   fDetectorPID=0x0;
1147   if (source.fDetectorPID) fDetectorPID = new AliDetectorPID(*source.fDetectorPID);
1148
1149   fTrackPhiOnEMCal= source.fTrackPhiOnEMCal;
1150   fTrackEtaOnEMCal= source.fTrackEtaOnEMCal;
1151   fTrackPtOnEMCal= source.fTrackPtOnEMCal;
1152
1153   if(fTOFcluster){
1154     delete[] fTOFcluster;
1155   }
1156   fNtofClusters = source.fNtofClusters;
1157   if(fNtofClusters > 0){
1158     fTOFcluster = new Int_t[fNtofClusters];
1159         for(Int_t i=0;i < fNtofClusters;i++) fTOFcluster[i] = source.fTOFcluster[i];
1160   }
1161
1162   return *this;
1163 }
1164
1165
1166
1167 void AliESDtrack::Copy(TObject &obj) const {
1168
1169   // this overwrites the virtual TOBject::Copy()
1170   // to allow run time copying without casting
1171   // in AliESDEvent
1172
1173   if(this==&obj)return;
1174   AliESDtrack *robj = dynamic_cast<AliESDtrack*>(&obj);
1175   if(!robj)return; // not an AliESDtrack
1176   *robj = *this;
1177
1178 }
1179
1180
1181
1182 void AliESDtrack::AddCalibObject(TObject * object){
1183   //
1184   // add calib object to the list
1185   //
1186   if (!fFriendTrack) fFriendTrack  = new AliESDfriendTrack;
1187   if (!fFriendTrack) return;
1188   fFriendTrack->AddCalibObject(object);
1189 }
1190
1191 TObject *  AliESDtrack::GetCalibObject(Int_t index){
1192   //
1193   // return calib objct at given position
1194   //
1195   if (!fFriendTrack) return 0;
1196   return fFriendTrack->GetCalibObject(index);
1197 }
1198
1199
1200 Bool_t AliESDtrack::FillTPCOnlyTrack(AliESDtrack &track){
1201
1202   // Fills the information of the TPC-only first reconstruction pass
1203   // into the passed ESDtrack object. For consistency fTPCInner is also filled
1204   // again
1205
1206
1207
1208   // For data produced before r26675
1209   // RelateToVertexTPC was not properly called during reco
1210   // so you'll have to call it again, before FillTPCOnlyTrack
1211   //  Float_t p[2],cov[3];
1212   // track->GetImpactParametersTPC(p,cov);
1213   // if(p[0]==0&&p[1]==0) // <- Default values
1214   //  track->RelateToVertexTPC(esd->GetPrimaryVertexTPC(),esd->GetMagneticField(),kVeryBig);
1215
1216
1217   if(!fTPCInner)return kFALSE;
1218
1219   // fill the TPC track params to the global track parameters
1220   track.Set(fTPCInner->GetX(),fTPCInner->GetAlpha(),fTPCInner->GetParameter(),fTPCInner->GetCovariance());
1221   track.fD = fdTPC;
1222   track.fZ = fzTPC;
1223   track.fCdd = fCddTPC;
1224   track.fCdz = fCdzTPC;
1225   track.fCzz = fCzzTPC;
1226
1227   // copy the inner params
1228   if(track.fIp) *track.fIp = *fIp;
1229   else track.fIp = new AliExternalTrackParam(*fIp);
1230
1231   // copy the TPCinner parameters
1232   if(track.fTPCInner) *track.fTPCInner = *fTPCInner;
1233   else track.fTPCInner = new AliExternalTrackParam(*fTPCInner);
1234   track.fdTPC   = fdTPC;
1235   track.fzTPC   = fzTPC;
1236   track.fCddTPC = fCddTPC;
1237   track.fCdzTPC = fCdzTPC;
1238   track.fCzzTPC = fCzzTPC;
1239   track.fCchi2TPC = fCchi2TPC;
1240
1241   // copy all other TPC specific parameters
1242
1243   // replace label by TPC label
1244   track.fLabel    = fTPCLabel;
1245   track.fTPCLabel = fTPCLabel;
1246
1247   track.fTPCchi2 = fTPCchi2;
1248   track.fTPCchi2Iter1 = fTPCchi2Iter1;
1249   track.fTPCsignal = fTPCsignal;
1250   track.fTPCsignalTuned = fTPCsignalTuned;
1251   track.fTPCsignalS = fTPCsignalS;
1252   for(int i = 0;i<4;++i)track.fTPCPoints[i] = fTPCPoints[i];
1253
1254   track.fTPCncls    = fTPCncls;
1255   track.fTPCnclsF   = fTPCnclsF;
1256   track.fTPCsignalN =  fTPCsignalN;
1257   track.fTPCnclsIter1    = fTPCnclsIter1;
1258   track.fTPCnclsFIter1   = fTPCnclsFIter1;
1259
1260   // PID
1261   if (fTPCr) {
1262     if (!track.fTPCr) track.fTPCr = new Double32_t[AliPID::kSPECIES];
1263     for(int i=AliPID::kSPECIES;i--;) track.fTPCr[i] = fTPCr[i];
1264   }
1265   //
1266   if (fR) {
1267     if (!track.fR) track.fR = new Double32_t[AliPID::kSPECIES];
1268     for(int i=AliPID::kSPECIES;i--;) track.fR[i] = fR[i];
1269   }
1270
1271   track.fTPCFitMap = fTPCFitMap;
1272   track.fTPCClusterMap = fTPCClusterMap;
1273   track.fTPCSharedMap = fTPCSharedMap;
1274
1275
1276   // reset the flags
1277   track.fFlags = kTPCin;
1278   track.fID    = fID;
1279
1280   track.fFlags |= fFlags & kTPCpid; //copy the TPCpid status flag
1281
1282   for (Int_t i=0;i<3;i++) track.fKinkIndexes[i] = fKinkIndexes[i];
1283
1284   return kTRUE;
1285
1286 }
1287
1288 //_______________________________________________________________________
1289 void AliESDtrack::MakeMiniESDtrack(){
1290   // Resets everything except
1291   // fFlags: Reconstruction status flags
1292   // fLabel: Track label
1293   // fID:  Unique ID of the track
1294   // Impact parameter information
1295   // fR[AliPID::kSPECIES]: combined "detector response probability"
1296   // Running track parameters in the base class (AliExternalTrackParam)
1297
1298   fTrackLength = 0;
1299
1300   if(fTrackTime)
1301     for (Int_t i=0;i<AliPID::kSPECIESC;i++) fTrackTime[i] = 0;
1302
1303   // Reset track parameters constrained to the primary vertex
1304   delete fCp;fCp = 0;
1305
1306   // Reset track parameters at the inner wall of TPC
1307   delete fIp;fIp = 0;
1308   delete fTPCInner;fTPCInner=0;
1309   // Reset track parameters at the inner wall of the TRD
1310   delete fOp;fOp = 0;
1311   // Reset track parameters at the HMPID
1312   delete fHMPIDp;fHMPIDp = 0;
1313
1314
1315   // Reset ITS track related information
1316   fITSchi2 = 0;
1317   fITSncls = 0;
1318   fITSClusterMap=0;
1319   fITSSharedMap=0;
1320   fITSsignal = 0;
1321   for (Int_t i=0;i<4;i++) fITSdEdxSamples[i] = 0.;
1322   if (fITSr) for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0;
1323   fITSLabel = 0;
1324
1325   // Reset TPC related track information
1326   fTPCchi2 = 0;
1327   fTPCchi2Iter1 = 0;
1328   fTPCncls = 0;
1329   fTPCnclsF = 0;
1330   fTPCnclsIter1 = 0;
1331   fTPCnclsFIter1 = 0;
1332   fTPCFitMap = 0;
1333   fTPCClusterMap = 0;
1334   fTPCSharedMap = 0;
1335   fTPCsignal= 0;
1336   fTPCsignalTuned= 0;
1337   fTPCsignalS= 0;
1338   fTPCsignalN= 0;
1339   if (fTPCr) for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i] = 0;
1340   fTPCLabel=0;
1341   for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
1342   for (Int_t i=0; i<3;i++)   fKinkIndexes[i] = 0;
1343   for (Int_t i=0; i<3;i++)   fV0Indexes[i] = 0;
1344
1345   // Reset TRD related track information
1346   fTRDchi2 = 0;
1347   fTRDncls = 0;
1348   fTRDncls0 = 0;
1349   fTRDsignal = 0;
1350   fTRDNchamberdEdx = 0;
1351   fTRDNclusterdEdx = 0;
1352
1353   for (Int_t i=0;i<kTRDnPlanes;i++) {
1354     fTRDTimBin[i]  = 0;
1355   }
1356   if (fTRDr) for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0;
1357   fTRDLabel = 0;
1358   fTRDQuality  = 0;
1359   fTRDntracklets = 0;
1360   if(fTRDnSlices)
1361     delete[] fTRDslices;
1362   fTRDslices=0x0;
1363   fTRDnSlices=0;
1364   fTRDBudget  = 0;
1365
1366   // Reset TOF related track information
1367   fTOFchi2 = 0;
1368   fTOFindex = -1;
1369   fTOFsignal = 99999;
1370   fTOFCalChannel = -1;
1371   fTOFsignalToT = 99999;
1372   fTOFsignalRaw = 99999;
1373   fTOFsignalDz = 999;
1374   fTOFsignalDx = 999;
1375   fTOFdeltaBC = 999;
1376   fTOFl0l1 = 999;
1377   if (fTOFr) for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
1378   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
1379
1380   // Reset HMPID related track information
1381   fHMPIDchi2 = 0;
1382   fHMPIDqn = 0;
1383   fHMPIDcluIdx = -1;
1384   fHMPIDsignal = 0;
1385   if (fHMPIDr) for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
1386   fHMPIDtrkTheta = 0;
1387   fHMPIDtrkPhi = 0;
1388   fHMPIDtrkX = 0;
1389   fHMPIDtrkY = 0;
1390   fHMPIDmipX = 0;
1391   fHMPIDmipY = 0;
1392   fCaloIndex = kEMCALNoMatch;
1393
1394   // reset global track chi2
1395   fGlobalChi2 = 0;
1396
1397   fVertexID = -2; // an orphan track
1398   fPIDForTracking = AliPID::kPion;
1399   //
1400   delete fFriendTrack; fFriendTrack = 0;
1401 }
1402
1403 //_______________________________________________________________________
1404 Int_t AliESDtrack::GetPID(Bool_t tpcOnly) const
1405 {
1406   // Returns the particle most probable id. For backward compatibility first the prob. arrays
1407   // will be checked, but normally the GetPIDForTracking will be returned
1408   Int_t i;
1409   const Double32_t *prob = 0;
1410   if (tpcOnly) { // check if TPCpid is valid
1411     if (!fTPCr) return GetPIDForTracking();
1412     prob = fTPCr;
1413     for (i=0; i<AliPID::kSPECIES-1; i++) if (prob[i] != prob[i+1]) break;
1414     if (i == AliPID::kSPECIES-1) prob = 0; // not valid, try with combined pid
1415   }
1416   if (!prob) { // either requested TPCpid is not valid or comb.pid is requested
1417     if (!fR) return GetPIDForTracking();
1418     prob = fR;
1419     for (i=0; i<AliPID::kSPECIES-1; i++) if (prob[i] != prob[i+1]) break;
1420     if (i == AliPID::kSPECIES-1) return GetPIDForTracking();  // If all the probabilities are equal, return the pion mass
1421   }
1422   //
1423   Float_t max=0.;
1424   Int_t k=-1;
1425   for (i=0; i<AliPID::kSPECIES; i++) if (prob[i]>max) {k=i; max=prob[i];}
1426   //
1427   if (k==0) { // dE/dx "crossing points" in the TPC
1428     /*
1429     Double_t p=GetP();
1430     if ((p>0.38)&&(p<0.48))
1431       if (prob[0]<prob[3]*10.) return AliPID::kKaon;
1432     if ((p>0.75)&&(p<0.85))
1433       if (prob[0]<prob[4]*10.) return AliPID::kProton;
1434     */
1435     return AliPID::kElectron;
1436   }
1437   if (k==1) return AliPID::kMuon;
1438   if (k==2||k==-1) return AliPID::kPion;
1439   if (k==3) return AliPID::kKaon;
1440   if (k==4) return AliPID::kProton;
1441   //  AliWarning("Undefined PID !");
1442   return GetPIDForTracking();
1443 }
1444
1445 //_______________________________________________________________________
1446 Int_t AliESDtrack::GetTOFBunchCrossing(Double_t b, Bool_t pidTPConly) const
1447 {
1448   // Returns the number of bunch crossings after trigger (assuming 25ns spacing)
1449   const double kSpacing = 25e3; // min interbanch spacing
1450   const double kShift = 0;
1451   Int_t bcid = kTOFBCNA; // defualt one
1452   if (!IsOn(kTOFout)/* || !IsOn(kESDpid)*/) return bcid; // no info
1453   //
1454   double tdif = GetTOFsignal();
1455   if (IsOn(kTIME)) { // integrated time info is there
1456     int pid = GetPID(pidTPConly);
1457     Double_t times[AliPID::kSPECIESC];
1458     // old esd has only AliPID::kSPECIES times
1459     GetIntegratedTimes(times,pid>=AliPID::kSPECIES ? AliPID::kSPECIESC : AliPID::kSPECIES);
1460     tdif -= times[pid];
1461   }
1462   else { // assume integrated time info from TOF radius and momentum
1463     const double kRTOF = 385.;
1464     const double kCSpeed = 3.e-2; // cm/ps
1465     double p = GetP();
1466     if (p<0.01) return bcid;
1467     double m = GetMass(pidTPConly);
1468     double curv = GetC(b);
1469     double path = TMath::Abs(curv)>kAlmost0 ? // account for curvature
1470       2./curv*TMath::ASin(kRTOF*curv/2.)*TMath::Sqrt(1.+GetTgl()*GetTgl()) : kRTOF;
1471     tdif -= path/kCSpeed*TMath::Sqrt(1.+m*m/(p*p));
1472   }
1473   bcid = TMath::Nint((tdif - kShift)/kSpacing);
1474   return bcid;
1475 }
1476
1477 //______________________________________________________________________________
1478 Double_t AliESDtrack::M() const
1479 {
1480   // Returns the assumed mass
1481   // (the pion mass, if the particle can't be identified properly).
1482   static Bool_t printerr=kTRUE;
1483   if (printerr) {
1484      AliWarning("WARNING !!! ... THIS WILL BE PRINTED JUST ONCE !!!");
1485      printerr = kFALSE;
1486      AliWarning("This is the ESD mass. Use it with care !");
1487   }
1488   return GetMass();
1489 }
1490
1491 //______________________________________________________________________________
1492 Double_t AliESDtrack::E() const
1493 {
1494   // Returns the energy of the particle given its assumed mass.
1495   // Assumes the pion mass if the particle can't be identified properly.
1496
1497   Double_t m = M();
1498   Double_t p = P();
1499   return TMath::Sqrt(p*p + m*m);
1500 }
1501
1502 //______________________________________________________________________________
1503 Double_t AliESDtrack::Y() const
1504 {
1505   // Returns the rapidity of a particle given its assumed mass.
1506   // Assumes the pion mass if the particle can't be identified properly.
1507
1508   Double_t e = E();
1509   Double_t pz = Pz();
1510   if (e != TMath::Abs(pz)) { // energy was not equal to pz
1511     return 0.5*TMath::Log((e+pz)/(e-pz));
1512   } else { // energy was equal to pz
1513     return -999.;
1514   }
1515 }
1516
1517 //_______________________________________________________________________
1518 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
1519   //
1520   // This function updates track's running parameters
1521   //
1522   Bool_t rc=kTRUE;
1523
1524   SetStatus(flags);
1525   fLabel=t->GetLabel();
1526
1527   if (t->IsStartedTimeIntegral()) {
1528     SetStatus(kTIME);
1529     Double_t times[AliPID::kSPECIESC];
1530     t->GetIntegratedTimes(times);
1531     SetIntegratedTimes(times);
1532     SetIntegratedLength(t->GetIntegratedLength());
1533   }
1534
1535   Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1536   if (fFriendTrack) {
1537   if (flags==kITSout) fFriendTrack->SetITSOut(*t);
1538   if (flags==kTPCout) fFriendTrack->SetTPCOut(*t);
1539   if (flags==kTRDrefit) fFriendTrack->SetTRDIn(*t);
1540   }
1541
1542   switch (flags) {
1543
1544   case kITSin:
1545     fITSchi2Std[0] = t->GetChi2();
1546     //
1547   case kITSout:
1548     fITSchi2Std[1] = t->GetChi2();
1549   case kITSrefit:
1550     {
1551     fITSchi2Std[2] = t->GetChi2();
1552     fITSClusterMap=0;
1553     fITSncls=t->GetNumberOfClusters();
1554     if (fFriendTrack) {
1555     Int_t* indexITS = new Int_t[AliESDfriendTrack::kMaxITScluster];
1556     for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
1557         indexITS[i]=t->GetClusterIndex(i);
1558
1559         if (i<fITSncls) {
1560           Int_t l=(indexITS[i] & 0xf0000000) >> 28;
1561            SETBIT(fITSClusterMap,l);
1562         }
1563     }
1564     fFriendTrack->SetITSIndices(indexITS,AliESDfriendTrack::kMaxITScluster);
1565     delete [] indexITS;
1566     }
1567
1568     fITSchi2=t->GetChi2();
1569     fITSsignal=t->GetPIDsignal();
1570     fITSLabel = t->GetLabel();
1571     // keep in fOp the parameters outside ITS for ITS stand-alone tracks
1572     if (flags==kITSout) {
1573       if (!fOp) fOp=new AliExternalTrackParam(*t);
1574       else
1575         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1576     }
1577     }
1578     break;
1579
1580   case kTPCin: case kTPCrefit:
1581     {
1582     fTPCLabel = t->GetLabel();
1583     if (flags==kTPCin)  {
1584       fTPCInner=new AliExternalTrackParam(*t);
1585       fTPCnclsIter1=t->GetNumberOfClusters();
1586       fTPCchi2Iter1=t->GetChi2();
1587     }
1588     if (!fIp) fIp=new AliExternalTrackParam(*t);
1589     else
1590       fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1591     }
1592     // Intentionally no break statement; need to set general TPC variables as well
1593   case kTPCout:
1594     {
1595     if (flags & kTPCout){
1596       if (!fOp) fOp=new AliExternalTrackParam(*t);
1597       else
1598         fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1599     }
1600     fTPCncls=t->GetNumberOfClusters();
1601     fTPCchi2=t->GetChi2();
1602
1603     if (fFriendTrack) {  // Copy cluster indices
1604       Int_t* indexTPC = new Int_t[AliESDfriendTrack::kMaxTPCcluster];
1605       for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)
1606         indexTPC[i]=t->GetClusterIndex(i);
1607       fFriendTrack->SetTPCIndices(indexTPC,AliESDfriendTrack::kMaxTPCcluster);
1608       delete [] indexTPC;
1609     }
1610     fTPCsignal=t->GetPIDsignal();
1611     }
1612     break;
1613
1614   case kTRDin: case kTRDrefit:
1615     break;
1616   case kTRDout:
1617     {
1618     fTRDLabel = t->GetLabel();
1619     fTRDchi2  = t->GetChi2();
1620     fTRDncls  = t->GetNumberOfClusters();
1621     if (fFriendTrack) {
1622       Int_t* indexTRD = new Int_t[AliESDfriendTrack::kMaxTRDcluster];
1623       for (Int_t i=0;i<AliESDfriendTrack::kMaxTRDcluster;i++) indexTRD[i]=-2;
1624       for (Int_t i=0;i<6;i++) indexTRD[i]=t->GetTrackletIndex(i);
1625       fFriendTrack->SetTRDIndices(indexTRD,AliESDfriendTrack::kMaxTRDcluster);
1626       delete [] indexTRD;
1627     }
1628
1629     //commented out by Xianguo
1630     //fTRDsignal=t->GetPIDsignal();
1631     }
1632     break;
1633   case kTRDbackup:
1634     if (!fOp) fOp=new AliExternalTrackParam(*t);
1635     else
1636       fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1637     fTRDncls0 = t->GetNumberOfClusters();
1638     break;
1639   case kTOFin:
1640     break;
1641   case kTOFout:
1642     break;
1643   case kTRDStop:
1644     break;
1645   case kHMPIDout:
1646   if (!fHMPIDp) fHMPIDp=new AliExternalTrackParam(*t);
1647     else
1648       fHMPIDp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1649     break;
1650   default:
1651     AliError("Wrong flag !");
1652     return kFALSE;
1653   }
1654
1655   return rc;
1656 }
1657
1658 //_______________________________________________________________________
1659 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
1660   //---------------------------------------------------------------------
1661   // This function returns external representation of the track parameters
1662   //---------------------------------------------------------------------
1663   x=GetX();
1664   for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
1665 }
1666
1667 //_______________________________________________________________________
1668 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
1669   //---------------------------------------------------------------------
1670   // This function returns external representation of the cov. matrix
1671   //---------------------------------------------------------------------
1672   for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
1673 }
1674
1675 //_______________________________________________________________________
1676 Bool_t AliESDtrack::GetConstrainedExternalParameters
1677                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1678   //---------------------------------------------------------------------
1679   // This function returns the constrained external track parameters
1680   //---------------------------------------------------------------------
1681   if (!fCp) return kFALSE;
1682   alpha=fCp->GetAlpha();
1683   x=fCp->GetX();
1684   for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
1685   return kTRUE;
1686 }
1687
1688 //_______________________________________________________________________
1689 Bool_t
1690 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
1691   //---------------------------------------------------------------------
1692   // This function returns the constrained external cov. matrix
1693   //---------------------------------------------------------------------
1694   if (!fCp) return kFALSE;
1695   for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
1696   return kTRUE;
1697 }
1698
1699 Bool_t
1700 AliESDtrack::GetInnerExternalParameters
1701                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1702   //---------------------------------------------------------------------
1703   // This function returns external representation of the track parameters
1704   // at the inner layer of TPC
1705   //---------------------------------------------------------------------
1706   if (!fIp) return kFALSE;
1707   alpha=fIp->GetAlpha();
1708   x=fIp->GetX();
1709   for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
1710   return kTRUE;
1711 }
1712
1713 Bool_t
1714 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
1715  //---------------------------------------------------------------------
1716  // This function returns external representation of the cov. matrix
1717  // at the inner layer of TPC
1718  //---------------------------------------------------------------------
1719   if (!fIp) return kFALSE;
1720   for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
1721   return kTRUE;
1722 }
1723
1724 void
1725 AliESDtrack::SetOuterParam(const AliExternalTrackParam *p, ULong_t flags) {
1726   //
1727   // This is a direct setter for the outer track parameters
1728   //
1729   SetStatus(flags);
1730   if (fOp) delete fOp;
1731   fOp=new AliExternalTrackParam(*p);
1732 }
1733
1734 void
1735 AliESDtrack::SetOuterHmpParam(const AliExternalTrackParam *p, ULong_t flags) {
1736   //
1737   // This is a direct setter for the outer track parameters
1738   //
1739   SetStatus(flags);
1740   if (fHMPIDp) delete fHMPIDp;
1741   fHMPIDp=new AliExternalTrackParam(*p);
1742 }
1743
1744 Bool_t
1745 AliESDtrack::GetOuterExternalParameters
1746                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1747   //---------------------------------------------------------------------
1748   // This function returns external representation of the track parameters
1749   // at the inner layer of TRD
1750   //---------------------------------------------------------------------
1751   if (!fOp) return kFALSE;
1752   alpha=fOp->GetAlpha();
1753   x=fOp->GetX();
1754   for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
1755   return kTRUE;
1756 }
1757
1758 Bool_t
1759 AliESDtrack::GetOuterHmpExternalParameters
1760                  (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1761   //---------------------------------------------------------------------
1762   // This function returns external representation of the track parameters
1763   // at the inner layer of TRD
1764   //---------------------------------------------------------------------
1765   if (!fHMPIDp) return kFALSE;
1766   alpha=fHMPIDp->GetAlpha();
1767   x=fHMPIDp->GetX();
1768   for (Int_t i=0; i<5; i++) p[i]=fHMPIDp->GetParameter()[i];
1769   return kTRUE;
1770 }
1771
1772 Bool_t
1773 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
1774  //---------------------------------------------------------------------
1775  // This function returns external representation of the cov. matrix
1776  // at the inner layer of TRD
1777  //---------------------------------------------------------------------
1778   if (!fOp) return kFALSE;
1779   for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
1780   return kTRUE;
1781 }
1782
1783 Bool_t
1784 AliESDtrack::GetOuterHmpExternalCovariance(Double_t cov[15]) const {
1785  //---------------------------------------------------------------------
1786  // This function returns external representation of the cov. matrix
1787  // at the inner layer of TRD
1788  //---------------------------------------------------------------------
1789   if (!fHMPIDp) return kFALSE;
1790   for (Int_t i=0; i<15; i++) cov[i]=fHMPIDp->GetCovariance()[i];
1791   return kTRUE;
1792 }
1793
1794 Int_t AliESDtrack::GetNcls(Int_t idet) const
1795 {
1796   // Get number of clusters by subdetector index
1797   //
1798   Int_t ncls = 0;
1799   switch(idet){
1800   case 0:
1801     ncls = fITSncls;
1802     break;
1803   case 1:
1804     ncls = fTPCncls;
1805     break;
1806   case 2:
1807     ncls = fTRDncls;
1808     break;
1809   case 3:
1810     if (fTOFindex != -1)
1811       ncls = 1;
1812     break;
1813   case 4: //PHOS
1814     break;
1815   case 5: //HMPID
1816     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1817       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1818         ncls = 1;
1819       }
1820     }
1821     break;
1822   default:
1823     break;
1824   }
1825   return ncls;
1826 }
1827
1828 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
1829 {
1830   // Get cluster index array by subdetector index
1831   //
1832   Int_t ncls = 0;
1833   switch(idet){
1834   case 0:
1835     ncls = GetITSclusters(idx);
1836     break;
1837   case 1:
1838     ncls = GetTPCclusters(idx);
1839     break;
1840   case 2:
1841     ncls = GetTRDclusters(idx);
1842     break;
1843   case 3:
1844     if (fTOFindex != -1) {
1845       idx[0] = fTOFindex;
1846       ncls = 1;
1847     }
1848     break;
1849   case 4: //PHOS
1850     break;
1851   case 5:
1852     if ((fHMPIDcluIdx >= 0) && (fHMPIDcluIdx < 7000000)) {
1853       if ((fHMPIDcluIdx%1000000 != 9999) && (fHMPIDcluIdx%1000000 != 99999)) {
1854         idx[0] = GetHMPIDcluIdx();
1855         ncls = 1;
1856       }
1857     }
1858     break;
1859   case 6: //EMCAL
1860     break;
1861   default:
1862     break;
1863   }
1864   return ncls;
1865 }
1866
1867 //_______________________________________________________________________
1868 void AliESDtrack::GetIntegratedTimes(Double_t *times, Int_t nspec) const
1869 {
1870   // get integrated time for requested N species
1871   if (nspec<1) return;
1872   if(fNtofClusters>0 && GetESDEvent()){
1873     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
1874     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
1875     //
1876     for(int i=tofcl->GetNMatchableTracks();i--;){
1877       if(tofcl->GetTrackIndex(i) == GetID()) {
1878         for (int j=nspec; j--;) times[j]=tofcl->GetIntegratedTime(j,i);
1879         return;
1880       }
1881     }
1882   }
1883   else if(fNtofClusters>0)
1884     AliInfo("No AliESDEvent available here!\n");
1885
1886   // Returns the array with integrated times for each particle hypothesis
1887   if(fTrackTime)
1888     for (int i=nspec; i--;) times[i]=fTrackTime[i];
1889   else
1890     for (int i=AliPID::kSPECIESC; i--;) times[i]=0.0;
1891   //
1892 }
1893 //_______________________________________________________________________
1894 Double_t AliESDtrack::GetIntegratedLength() const{
1895   Int_t index = -1;
1896   if(fNtofClusters>0 && GetESDEvent()){
1897     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
1898     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
1899
1900     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
1901       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
1902     }
1903
1904     if(fNtofClusters>0 && index > -1)
1905       return tofcl->GetLength(index);
1906   }
1907   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
1908
1909   return fTrackLength;
1910 }
1911
1912 //_______________________________________________________________________
1913 void AliESDtrack::SetIntegratedTimes(const Double_t *times)
1914 {
1915   // Sets the array with integrated times for each particle hypotesis
1916   if(!fTrackTime) fTrackTime = new Double32_t[AliPID::kSPECIESC];
1917   for (int i=AliPID::kSPECIESC; i--;) fTrackTime[i]=times[i];
1918 }
1919
1920 //_______________________________________________________________________
1921 void AliESDtrack::SetITSpid(const Double_t *p)
1922 {
1923   // Sets values for the probability of each particle type (in ITS)
1924   if (!fITSr) fITSr = new Double32_t[AliPID::kSPECIESC];
1925   SetPIDValues(fITSr,p,AliPID::kSPECIES);
1926   SetStatus(AliESDtrack::kITSpid);
1927 }
1928
1929 //_______________________________________________________________________
1930 void AliESDtrack::GetITSpid(Double_t *p) const {
1931   // Gets the probability of each particle type (in ITS)
1932   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i] = fITSr ? fITSr[i] : 0;
1933 }
1934
1935 //_______________________________________________________________________
1936 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
1937   //---------------------------------------------------------------------
1938   // This function returns indices of the assgined ITS clusters
1939   //---------------------------------------------------------------------
1940   if (idx && fFriendTrack) {
1941     Int_t *index=fFriendTrack->GetITSindices();
1942     for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) {
1943       if ( (i>=fITSncls) && (i<6) ) idx[i]=-1;
1944       else {
1945         if (index) {
1946           idx[i]=index[i];
1947         }
1948         else idx[i]= -2;
1949       }
1950     }
1951   }
1952   return fITSncls;
1953 }
1954
1955 //_______________________________________________________________________
1956 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
1957                                          Float_t &xloc,Float_t &zloc) const {
1958   //----------------------------------------------------------------------
1959   // This function encodes in the module number also the status of cluster association
1960   // "status" can have the following values:
1961   // 1 "found" (cluster is associated),
1962   // 2 "dead" (module is dead from OCDB),
1963   // 3 "skipped" (module or layer forced to be skipped),
1964   // 4 "outinz" (track out of z acceptance),
1965   // 5 "nocls" (no clusters in the road),
1966   // 6 "norefit" (cluster rejected during refit),
1967   // 7 "deadzspd" (holes in z in SPD)
1968   // Also given are the coordinates of the crossing point of track and module
1969   // (in the local module ref. system)
1970   // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
1971   //----------------------------------------------------------------------
1972
1973   if(fITSModule[ilayer]==-1) {
1974     idet = -1;
1975     status=0;
1976     xloc=-99.; zloc=-99.;
1977     return kFALSE;
1978   }
1979
1980   Int_t module = fITSModule[ilayer];
1981
1982   idet = Int_t(module/1000000);
1983
1984   module -= idet*1000000;
1985
1986   status = Int_t(module/100000);
1987
1988   module -= status*100000;
1989
1990   Int_t signs = Int_t(module/10000);
1991
1992   module-=signs*10000;
1993
1994   Int_t xInt = Int_t(module/100);
1995   module -= xInt*100;
1996
1997   Int_t zInt = module;
1998
1999   if(signs==1) { xInt*=1; zInt*=1; }
2000   if(signs==2) { xInt*=1; zInt*=-1; }
2001   if(signs==3) { xInt*=-1; zInt*=1; }
2002   if(signs==4) { xInt*=-1; zInt*=-1; }
2003
2004   xloc = 0.1*(Float_t)xInt;
2005   zloc = 0.1*(Float_t)zInt;
2006
2007   if(status==4) idet = -1;
2008
2009   return kTRUE;
2010 }
2011
2012 //_______________________________________________________________________
2013 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
2014   //---------------------------------------------------------------------
2015   // This function returns indices of the assgined ITS clusters
2016   //---------------------------------------------------------------------
2017   if (idx && fFriendTrack) {
2018     Int_t *index=fFriendTrack->GetTPCindices();
2019
2020     if (index){
2021       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
2022     }
2023     else {
2024       for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=-2;
2025     }
2026   }
2027   return fTPCncls;
2028 }
2029
2030 //_______________________________________________________________________
2031 Float_t AliESDtrack::GetTPCCrossedRows() const
2032 {
2033   // This function calls GetTPCClusterInfo with some default parameters which are used in the track selection and caches the outcome
2034   // because GetTPCClusterInfo is quite time-consuming
2035
2036   if (fCacheNCrossedRows > -1)
2037     return fCacheNCrossedRows;
2038
2039   fCacheNCrossedRows = GetTPCClusterInfo(2, 1);
2040   return fCacheNCrossedRows;
2041 }
2042
2043 //_______________________________________________________________________
2044 Float_t AliESDtrack::GetTPCClusterInfo(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
2045 {
2046   //
2047   // TPC cluster information
2048   // type 0: get fraction of found/findable clusters with neighbourhood definition
2049   //      1: findable clusters with neighbourhood definition
2050   //      2: found clusters
2051   // bitType:
2052   //      0 - all cluster used
2053   //      1 - clusters  used for the kalman update
2054   // definition of findable clusters:
2055   //            a cluster is defined as findable if there is another cluster
2056   //           within +- nNeighbours pad rows. The idea is to overcome threshold
2057   //           effects with a very simple algorithm.
2058   //
2059
2060
2061   Int_t found=0;
2062   Int_t findable=0;
2063   Int_t last=-nNeighbours;
2064   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
2065
2066   Int_t upperBound=clusterMap.GetNbits();
2067   if (upperBound>row1) upperBound=row1;
2068   for (Int_t i=row0; i<upperBound; ++i){
2069     //look to current row
2070     if (clusterMap[i]) {
2071       last=i;
2072       ++found;
2073       ++findable;
2074       continue;
2075     }
2076     //look to nNeighbours before
2077     if ((i-last)<=nNeighbours) {
2078       ++findable;
2079       continue;
2080     }
2081     //look to nNeighbours after
2082     for (Int_t j=i+1; j<i+1+nNeighbours; ++j){
2083       if (clusterMap[j]){
2084         ++findable;
2085         break;
2086       }
2087     }
2088   }
2089   if (type==2) return found;
2090   if (type==1) return findable;
2091
2092   if (type==0){
2093     Float_t fraction=0;
2094     if (findable>0)
2095       fraction=(Float_t)found/(Float_t)findable;
2096     else
2097       fraction=0;
2098     return fraction;
2099   }
2100   return 0;  // undefined type - default value
2101 }
2102
2103 //_______________________________________________________________________
2104 Float_t AliESDtrack::GetTPCClusterDensity(Int_t nNeighbours/*=3*/, Int_t type/*=0*/, Int_t row0, Int_t row1, Int_t bitType ) const
2105 {
2106   //
2107   // TPC cluster density -  only rows where signal before and after given row are used
2108   //                     -  slower function
2109   // type 0: get fraction of found/findable clusters with neighbourhood definition
2110   //      1: findable clusters with neighbourhood definition
2111   //      2: found clusters
2112   // bitType:
2113   //      0 - all cluster used
2114   //      1 - clusters  used for the kalman update
2115   // definition of findable clusters:
2116   //            a cluster is defined as findable if there is another cluster
2117   //           within +- nNeighbours pad rows. The idea is to overcome threshold
2118   //           effects with a very simple algorithm.
2119   //
2120   Int_t found=0;
2121   Int_t findable=0;
2122   //  Int_t last=-nNeighbours;
2123   const TBits & clusterMap = (bitType%2==0) ? fTPCClusterMap : fTPCFitMap;
2124   Int_t upperBound=clusterMap.GetNbits();
2125   if (upperBound>row1) upperBound=row1;
2126   for (Int_t i=row0; i<upperBound; ++i){
2127     Bool_t isUp=kFALSE;
2128     Bool_t isDown=kFALSE;
2129     for (Int_t idelta=1; idelta<=nNeighbours; idelta++){
2130       if (i-idelta>=0 && clusterMap[i-idelta]) isDown=kTRUE;
2131       if (i+idelta<upperBound && clusterMap[i+idelta]) isUp=kTRUE;
2132     }
2133     if (isUp&&isDown){
2134       ++findable;
2135       if (clusterMap[i]) ++found;
2136     }
2137   }
2138   if (type==2) return found;
2139   if (type==1) return findable;
2140
2141   if (type==0){
2142     Float_t fraction=0;
2143     if (findable>0)
2144       fraction=(Float_t)found/(Float_t)findable;
2145     else
2146       fraction=0;
2147     return fraction;
2148   }
2149   return 0;  // undefined type - default value
2150 }
2151
2152
2153
2154
2155 //_______________________________________________________________________
2156 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
2157   //
2158   // GetDensity of the clusters on given region between row0 and row1
2159   // Dead zone effect takin into acoount
2160   //
2161   if (!fFriendTrack) return 0.0;
2162   Int_t good  = 0;
2163   Int_t found = 0;
2164   //
2165   Int_t *index=fFriendTrack->GetTPCindices();
2166   for (Int_t i=row0;i<=row1;i++){
2167     Int_t idx = index[i];
2168     if (idx!=-1)  good++;             // track outside of dead zone
2169     if (idx>0)    found++;
2170   }
2171   Float_t density=0.5;
2172   if (good>TMath::Max((row1-row0)*0.5,0.0)) density = Float_t(found)/Float_t(good);
2173   return density;
2174 }
2175
2176 //_______________________________________________________________________
2177 void AliESDtrack::SetTPCpid(const Double_t *p) {
2178   // Sets values for the probability of each particle type (in TPC)
2179   if (!fTPCr) fTPCr = new Double32_t[AliPID::kSPECIES];
2180   SetPIDValues(fTPCr,p,AliPID::kSPECIES);
2181   SetStatus(AliESDtrack::kTPCpid);
2182 }
2183
2184 //_______________________________________________________________________
2185 void AliESDtrack::GetTPCpid(Double_t *p) const {
2186   // Gets the probability of each particle type (in TPC)
2187   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i] = fTPCr ? fTPCr[i] : 0;
2188 }
2189
2190 //_______________________________________________________________________
2191 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
2192   //---------------------------------------------------------------------
2193   // This function returns indices of the assgined TRD clusters
2194   //---------------------------------------------------------------------
2195   if (idx && fFriendTrack) {
2196     Int_t *index=fFriendTrack->GetTRDindices();
2197
2198     if (index) {
2199       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
2200     }
2201     else {
2202       for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=-2;
2203     }
2204   }
2205   return fTRDncls;
2206 }
2207
2208 //_______________________________________________________________________
2209 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
2210 //
2211 // This function returns the number of TRD tracklets used in tracking
2212 // and it fills the indices of these tracklets in the array "idx" as they
2213 // are registered in the TRD track list.
2214 //
2215 // Caution :
2216 //   1. The idx array has to be allocated with a size >= AliESDtrack::kTRDnPlanes
2217 //   2. The idx array store not only the index but also the layer of the tracklet.
2218 //      Therefore tracks with TRD gaps contain default values for indices [-1]
2219
2220   if (!fFriendTrack) return 0;
2221   if (!idx) return GetTRDntracklets();
2222   Int_t *index=fFriendTrack->GetTRDindices();
2223   Int_t n = 0;
2224   for (Int_t i=0; i<kTRDnPlanes; i++){
2225     if (index){
2226       if(index[i]>=0) n++;
2227       idx[i]=index[i];
2228     }
2229     else idx[i] = -2;
2230   }
2231   return n;
2232 }
2233
2234 //_______________________________________________________________________
2235 void AliESDtrack::SetTRDpid(const Double_t *p) {
2236   // Sets values for the probability of each particle type (in TRD)
2237   if (!fTRDr) fTRDr = new Double32_t[AliPID::kSPECIES];
2238   SetPIDValues(fTRDr,p,AliPID::kSPECIES);
2239   SetStatus(AliESDtrack::kTRDpid);
2240 }
2241
2242 //_______________________________________________________________________
2243 void AliESDtrack::GetTRDpid(Double_t *p) const {
2244   // Gets the probability of each particle type (in TRD)
2245   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i] = fTRDr ? fTRDr[i]:0;
2246 }
2247
2248 //_______________________________________________________________________
2249 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
2250 {
2251   // Sets the probability of particle type iSpecies to p (in TRD)
2252   if (!fTRDr) {
2253     fTRDr = new Double32_t[AliPID::kSPECIES];
2254     for (Int_t i=AliPID::kSPECIES; i--;) fTRDr[i] = 0;
2255   }
2256   fTRDr[iSpecies] = p;
2257 }
2258
2259 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
2260 {
2261   // Returns the probability of particle type iSpecies (in TRD)
2262   return fTRDr ? fTRDr[iSpecies] : 0;
2263 }
2264
2265 //____________________________________________________
2266 Int_t AliESDtrack::GetNumberOfTRDslices() const
2267 {
2268   // built in backward compatibility
2269   Int_t idx = fTRDnSlices - (kTRDnPlanes<<1);
2270   return idx<18 ? fTRDnSlices/kTRDnPlanes : idx/kTRDnPlanes;
2271 }
2272
2273 //____________________________________________________
2274 Double_t AliESDtrack::GetTRDmomentum(Int_t plane, Double_t *sp) const
2275 {
2276 //Returns momentum estimation and optional its error (sp)
2277 // in TRD layer "plane".
2278
2279   if (!fTRDnSlices) {
2280     AliDebug(2, "No TRD info allocated for this track.");
2281     return -1.;
2282   }
2283   if ((plane<0) || (plane>=kTRDnPlanes)) {
2284     AliWarning(Form("Request for TRD plane[%d] outside range.", plane));
2285     return -1.;
2286   }
2287
2288   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
2289   // Protection for backward compatibility
2290   if(idx<(GetNumberOfTRDslices()*kTRDnPlanes)) return -1.;
2291
2292   if(sp) (*sp) = fTRDslices[idx+kTRDnPlanes];
2293   return fTRDslices[idx];
2294 }
2295
2296 //____________________________________________________
2297 Double_t  AliESDtrack::GetTRDslice(Int_t plane, Int_t slice) const {
2298   //Gets the charge from the slice of the plane
2299
2300   if(!fTRDslices) {
2301     //AliError("No TRD slices allocated for this track !");
2302     return -1.;
2303   }
2304   if ((plane<0) || (plane>=kTRDnPlanes)) {
2305     AliError("Info for TRD plane not available !");
2306     return -1.;
2307   }
2308   Int_t ns=GetNumberOfTRDslices();
2309   if ((slice<-1) || (slice>=ns)) {
2310     //AliError("Wrong TRD slice !");
2311     return -1.;
2312   }
2313
2314   if(slice>=0) return fTRDslices[plane*ns + slice];
2315
2316   // return average of the dEdx measurements
2317   Double_t q=0.; Double32_t *s = &fTRDslices[plane*ns];
2318   for (Int_t i=0; i<ns; i++, s++) if((*s)>0.) q+=(*s);
2319   return q/ns;
2320 }
2321
2322 //____________________________________________________
2323 void  AliESDtrack::SetNumberOfTRDslices(Int_t n) {
2324   //Sets the number of slices used for PID
2325   if (fTRDnSlices) return;
2326
2327   fTRDnSlices=n;
2328   fTRDslices=new Double32_t[fTRDnSlices];
2329
2330   // set-up correctly the allocated memory
2331   memset(fTRDslices, 0, n*sizeof(Double32_t));
2332   for (Int_t i=GetNumberOfTRDslices(); i--;) fTRDslices[i]=-1.;
2333 }
2334
2335 //____________________________________________________
2336 void  AliESDtrack::SetTRDslice(Double_t q, Int_t plane, Int_t slice) {
2337   //Sets the charge q in the slice of the plane
2338   if(!fTRDslices) {
2339     AliError("No TRD slices allocated for this track !");
2340     return;
2341   }
2342   if ((plane<0) || (plane>=kTRDnPlanes)) {
2343     AliError("Info for TRD plane not allocated !");
2344     return;
2345   }
2346   Int_t ns=GetNumberOfTRDslices();
2347   if ((slice<0) || (slice>=ns)) {
2348     AliError(Form("Wrong TRD slice %d/%d, NSlices=%d",plane,slice,ns));
2349     return;
2350   }
2351   Int_t n=plane*ns + slice;
2352   fTRDslices[n]=q;
2353 }
2354
2355
2356 //____________________________________________________
2357 void AliESDtrack::SetTRDmomentum(Double_t p, Int_t plane, Double_t *sp)
2358 {
2359   if(!fTRDslices) {
2360     AliError("No TRD slices allocated for this track !");
2361     return;
2362   }
2363   if ((plane<0) || (plane>=kTRDnPlanes)) {
2364     AliError("Info for TRD plane not allocated !");
2365     return;
2366   }
2367
2368   Int_t idx = fTRDnSlices-(kTRDnPlanes<<1)+plane;
2369   // Protection for backward compatibility
2370   if(idx<GetNumberOfTRDslices()*kTRDnPlanes) return;
2371
2372   if(sp) fTRDslices[idx+kTRDnPlanes] = (*sp);
2373   fTRDslices[idx] = p;
2374 }
2375
2376
2377 //_______________________________________________________________________
2378 void AliESDtrack::SetTOFpid(const Double_t *p) {
2379   // Sets the probability of each particle type (in TOF)
2380   if (!fTOFr) fTOFr = new Double32_t[AliPID::kSPECIES];
2381   SetPIDValues(fTOFr,p,AliPID::kSPECIES);
2382   SetStatus(AliESDtrack::kTOFpid);
2383 }
2384
2385 //_______________________________________________________________________
2386 void AliESDtrack::SetTOFLabel(const Int_t *p) {
2387
2388   if(fNtofClusters>0){
2389     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
2390     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
2391     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
2392
2393     if(hit) hit->SetTOFLabel(p);
2394   }
2395   else{
2396     // Sets  (in TOF)
2397     if(!fTOFLabel) fTOFLabel = new Int_t[3];
2398     for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
2399   }
2400 }
2401
2402 //_______________________________________________________________________
2403 void AliESDtrack::GetTOFpid(Double_t *p) const {
2404   // Gets probabilities of each particle type (in TOF)
2405   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i] = fTOFr ? fTOFr[i]:0;
2406 }
2407
2408 //_______________________________________________________________________
2409 void AliESDtrack::GetTOFLabel(Int_t *p) const {
2410   // Gets (in TOF)
2411   if(fNtofClusters>0){
2412     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
2413     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
2414
2415     for (Int_t i=0; i<3; i++) p[i]=tofcl->GetLabel(i);
2416   }
2417   else{
2418     if(fTOFLabel) for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
2419     else for (int i=3;i--;) p[i] = -1;
2420   }
2421 }
2422
2423 //_______________________________________________________________________
2424 void AliESDtrack::GetTOFInfo(Float_t *info) const {
2425   // Gets (in TOF)
2426   for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
2427 }
2428
2429 //_______________________________________________________________________
2430 void AliESDtrack::SetTOFInfo(Float_t*info) {
2431   // Gets (in TOF)
2432   for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
2433 }
2434
2435
2436
2437 //_______________________________________________________________________
2438 void AliESDtrack::SetHMPIDpid(const Double_t *p) {
2439   // Sets the probability of each particle type (in HMPID)
2440   if (!fHMPIDr) fHMPIDr = new Double32_t[AliPID::kSPECIES];
2441   SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
2442   SetStatus(AliESDtrack::kHMPIDpid);
2443 }
2444
2445 //_______________________________________________________________________
2446 void  AliESDtrack::SetTPCdEdxInfo(AliTPCdEdxInfo * dEdxInfo){
2447   if(fTPCdEdxInfo) delete fTPCdEdxInfo;
2448   fTPCdEdxInfo = dEdxInfo;
2449 }
2450
2451 //_______________________________________________________________________
2452 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
2453   // Gets probabilities of each particle type (in HMPID)
2454   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i] = fHMPIDr ? fHMPIDr[i]:0;
2455 }
2456
2457
2458
2459 //_______________________________________________________________________
2460 void AliESDtrack::SetESDpid(const Double_t *p) {
2461   // Sets the probability of each particle type for the ESD track
2462   if (!fR) fR = new Double32_t[AliPID::kSPECIES];
2463   SetPIDValues(fR,p,AliPID::kSPECIES);
2464   SetStatus(AliESDtrack::kESDpid);
2465 }
2466
2467 //_______________________________________________________________________
2468 void AliESDtrack::GetESDpid(Double_t *p) const {
2469   // Gets probability of each particle type for the ESD track
2470   for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i] = fR ? fR[i]:0;
2471 }
2472
2473 //_______________________________________________________________________
2474 Bool_t AliESDtrack::RelateToVertexTPC(const AliESDVertex *vtx,
2475 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2476   //
2477   // Try to relate the TPC-only track parameters to the vertex "vtx",
2478   // if the (rough) transverse impact parameter is not bigger then "maxd".
2479   //            Magnetic field is "b" (kG).
2480   //
2481   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2482   // b) The impact parameters and their covariance matrix are calculated.
2483   // c) An attempt to constrain the TPC-only params to the vertex is done.
2484   //    The constrained params are returned via "cParam".
2485   //
2486   // In the case of success, the returned value is kTRUE
2487   // otherwise, it's kFALSE)
2488   //
2489
2490   if (!fTPCInner) return kFALSE;
2491   if (!vtx) return kFALSE;
2492
2493   Double_t dz[2],cov[3];
2494   if (!fTPCInner->PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2495
2496   fdTPC = dz[0];
2497   fzTPC = dz[1];
2498   fCddTPC = cov[0];
2499   fCdzTPC = cov[1];
2500   fCzzTPC = cov[2];
2501
2502   Double_t covar[6]; vtx->GetCovMatrix(covar);
2503   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2504   Double_t c[3]={covar[2],0.,covar[5]};
2505
2506   Double_t chi2=GetPredictedChi2(p,c);
2507   if (chi2>kVeryBig) return kFALSE;
2508
2509   fCchi2TPC=chi2;
2510
2511   if (!cParam) return kTRUE;
2512
2513   *cParam = *fTPCInner;
2514   if (!cParam->Update(p,c)) return kFALSE;
2515
2516   return kTRUE;
2517 }
2518
2519 //_______________________________________________________________________
2520 Bool_t AliESDtrack::RelateToVertexTPCBxByBz(const AliESDVertex *vtx,
2521 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2522   //
2523   // Try to relate the TPC-only track parameters to the vertex "vtx",
2524   // if the (rough) transverse impact parameter is not bigger then "maxd".
2525   //
2526   // All three components of the magnetic field ,"b[3]" (kG),
2527   // are taken into account.
2528   //
2529   // a) The TPC-only paramters are extapolated to the DCA to the vertex.
2530   // b) The impact parameters and their covariance matrix are calculated.
2531   // c) An attempt to constrain the TPC-only params to the vertex is done.
2532   //    The constrained params are returned via "cParam".
2533   //
2534   // In the case of success, the returned value is kTRUE
2535   // otherwise, it's kFALSE)
2536   //
2537
2538   if (!fTPCInner) return kFALSE;
2539   if (!vtx) return kFALSE;
2540
2541   Double_t dz[2],cov[3];
2542   if (!fTPCInner->PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2543
2544   fdTPC = dz[0];
2545   fzTPC = dz[1];
2546   fCddTPC = cov[0];
2547   fCdzTPC = cov[1];
2548   fCzzTPC = cov[2];
2549
2550   Double_t covar[6]; vtx->GetCovMatrix(covar);
2551   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2552   Double_t c[3]={covar[2],0.,covar[5]};
2553
2554   Double_t chi2=GetPredictedChi2(p,c);
2555   if (chi2>kVeryBig) return kFALSE;
2556
2557   fCchi2TPC=chi2;
2558
2559   if (!cParam) return kTRUE;
2560
2561   *cParam = *fTPCInner;
2562   if (!cParam->Update(p,c)) return kFALSE;
2563
2564   return kTRUE;
2565 }
2566
2567 //_______________________________________________________________________
2568 Bool_t AliESDtrack::RelateToVertex(const AliESDVertex *vtx,
2569 Double_t b, Double_t maxd, AliExternalTrackParam *cParam) {
2570   //
2571   // Try to relate this track to the vertex "vtx",
2572   // if the (rough) transverse impact parameter is not bigger then "maxd".
2573   //            Magnetic field is "b" (kG).
2574   //
2575   // a) The track gets extapolated to the DCA to the vertex.
2576   // b) The impact parameters and their covariance matrix are calculated.
2577   // c) An attempt to constrain this track to the vertex is done.
2578   //    The constrained params are returned via "cParam".
2579   //
2580   // In the case of success, the returned value is kTRUE
2581   // (otherwise, it's kFALSE)
2582   //
2583
2584   if (!vtx) return kFALSE;
2585
2586   Double_t dz[2],cov[3];
2587   if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
2588
2589   fD = dz[0];
2590   fZ = dz[1];
2591   fCdd = cov[0];
2592   fCdz = cov[1];
2593   fCzz = cov[2];
2594
2595   Double_t covar[6]; vtx->GetCovMatrix(covar);
2596   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2597   Double_t c[3]={covar[2],0.,covar[5]};
2598
2599   Double_t chi2=GetPredictedChi2(p,c);
2600   if (chi2>kVeryBig) return kFALSE;
2601
2602   fCchi2=chi2;
2603
2604
2605   //--- Could now these lines be removed ? ---
2606   delete fCp;
2607   fCp=new AliExternalTrackParam(*this);
2608
2609   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2610   //----------------------------------------
2611
2612   fVertexID = vtx->GetID();
2613
2614   if (!cParam) return kTRUE;
2615
2616   *cParam = *this;
2617   if (!cParam->Update(p,c)) return kFALSE;
2618
2619   return kTRUE;
2620 }
2621
2622 //_______________________________________________________________________
2623 Bool_t AliESDtrack::RelateToVertexBxByBz(const AliESDVertex *vtx,
2624 Double_t b[3], Double_t maxd, AliExternalTrackParam *cParam) {
2625   //
2626   // Try to relate this track to the vertex "vtx",
2627   // if the (rough) transverse impact parameter is not bigger then "maxd".
2628   //            Magnetic field is "b" (kG).
2629   //
2630   // a) The track gets extapolated to the DCA to the vertex.
2631   // b) The impact parameters and their covariance matrix are calculated.
2632   // c) An attempt to constrain this track to the vertex is done.
2633   //    The constrained params are returned via "cParam".
2634   //
2635   // In the case of success, the returned value is kTRUE
2636   // (otherwise, it's kFALSE)
2637   //
2638
2639   if (!vtx) return kFALSE;
2640
2641   Double_t dz[2],cov[3];
2642   if (!PropagateToDCABxByBz(vtx, b, maxd, dz, cov)) return kFALSE;
2643
2644   fD = dz[0];
2645   fZ = dz[1];
2646   fCdd = cov[0];
2647   fCdz = cov[1];
2648   fCzz = cov[2];
2649
2650   Double_t covar[6]; vtx->GetCovMatrix(covar);
2651   Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
2652   Double_t c[3]={covar[2],0.,covar[5]};
2653
2654   Double_t chi2=GetPredictedChi2(p,c);
2655   if (chi2>kVeryBig) return kFALSE;
2656
2657   fCchi2=chi2;
2658
2659
2660   //--- Could now these lines be removed ? ---
2661   delete fCp;
2662   fCp=new AliExternalTrackParam(*this);
2663
2664   if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
2665   //----------------------------------------
2666
2667   fVertexID = vtx->GetID();
2668
2669   if (!cParam) return kTRUE;
2670
2671   *cParam = *this;
2672   if (!cParam->Update(p,c)) return kFALSE;
2673
2674   return kTRUE;
2675 }
2676
2677 //_______________________________________________________________________
2678 void AliESDtrack::Print(Option_t *) const {
2679   // Prints info on the track
2680   AliExternalTrackParam::Print();
2681   printf("ESD track info\n") ;
2682   Double_t p[AliPID::kSPECIES] ;
2683   Int_t index = 0 ;
2684   if( IsOn(kITSpid) ){
2685     printf("From ITS: ") ;
2686     GetITSpid(p) ;
2687     for(index = 0 ; index < AliPID::kSPECIES; index++)
2688       printf("%f, ", p[index]) ;
2689     printf("\n           signal = %f\n", GetITSsignal()) ;
2690   }
2691   if( IsOn(kTPCpid) ){
2692     printf("From TPC: ") ;
2693     GetTPCpid(p) ;
2694     for(index = 0 ; index < AliPID::kSPECIES; index++)
2695       printf("%f, ", p[index]) ;
2696     printf("\n           signal = %f\n", GetTPCsignal()) ;
2697   }
2698   if( IsOn(kTRDpid) ){
2699     printf("From TRD: ") ;
2700     GetTRDpid(p) ;
2701     for(index = 0 ; index < AliPID::kSPECIES; index++)
2702       printf("%f, ", p[index]) ;
2703       printf("\n           signal = %f\n", GetTRDsignal()) ;
2704       printf("\n           NchamberdEdx = %d\n", GetTRDNchamberdEdx()) ;
2705       printf("\n           NclusterdEdx = %d\n", GetTRDNclusterdEdx()) ;
2706   }
2707   if( IsOn(kTOFpid) ){
2708     printf("From TOF: ") ;
2709     GetTOFpid(p) ;
2710     for(index = 0 ; index < AliPID::kSPECIES; index++)
2711       printf("%f, ", p[index]) ;
2712     printf("\n           signal = %f\n", GetTOFsignal()) ;
2713   }
2714   if( IsOn(kHMPIDpid) ){
2715     printf("From HMPID: ") ;
2716     GetHMPIDpid(p) ;
2717     for(index = 0 ; index < AliPID::kSPECIES; index++)
2718       printf("%f, ", p[index]) ;
2719     printf("\n           signal = %f\n", GetHMPIDsignal()) ;
2720   }
2721 }
2722
2723
2724 //
2725 // Draw functionality
2726 // Origin: Marian Ivanov, Marian.Ivanov@cern.ch
2727 //
2728 void AliESDtrack::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
2729   //
2730   // Fill points in the polymarker
2731   //
2732   TObjArray arrayRef;
2733   arrayRef.AddLast(new AliExternalTrackParam(*this));
2734   if (fIp) arrayRef.AddLast(new AliExternalTrackParam(*fIp));
2735   if (fOp) arrayRef.AddLast(new AliExternalTrackParam(*fOp));
2736   if (fHMPIDp) arrayRef.AddLast(new AliExternalTrackParam(*fHMPIDp));
2737   //
2738   Double_t mpos[3]={0,0,0};
2739   Int_t entries=arrayRef.GetEntries();
2740   for (Int_t i=0;i<entries;i++){
2741     Double_t pos[3];
2742     ((AliExternalTrackParam*)arrayRef.At(i))->GetXYZ(pos);
2743     mpos[0]+=pos[0]/entries;
2744     mpos[1]+=pos[1]/entries;
2745     mpos[2]+=pos[2]/entries;
2746   }
2747   // Rotate to the mean position
2748   //
2749   Float_t fi= TMath::ATan2(mpos[1],mpos[0]);
2750   for (Int_t i=0;i<entries;i++){
2751     Bool_t res = ((AliExternalTrackParam*)arrayRef.At(i))->Rotate(fi);
2752     if (!res) delete arrayRef.RemoveAt(i);
2753   }
2754   Int_t counter=0;
2755   for (Double_t r=minR; r<maxR; r+=stepR){
2756     Double_t sweight=0;
2757     Double_t mlpos[3]={0,0,0};
2758     for (Int_t i=0;i<entries;i++){
2759       Double_t point[3]={0,0,0};
2760       AliExternalTrackParam *param = ((AliExternalTrackParam*)arrayRef.At(i));
2761       if (!param) continue;
2762       if (param->GetXYZAt(r,magF,point)){
2763         Double_t weight = 1./(10.+(r-param->GetX())*(r-param->GetX()));
2764         sweight+=weight;
2765         mlpos[0]+=point[0]*weight;
2766         mlpos[1]+=point[1]*weight;
2767         mlpos[2]+=point[2]*weight;
2768       }
2769     }
2770     if (sweight>0){
2771       mlpos[0]/=sweight;
2772       mlpos[1]/=sweight;
2773       mlpos[2]/=sweight;
2774       pol->SetPoint(counter,mlpos[0],mlpos[1], mlpos[2]);
2775       //      printf("xyz\t%f\t%f\t%f\n",mlpos[0], mlpos[1],mlpos[2]);
2776       counter++;
2777     }
2778   }
2779 }
2780
2781 //_______________________________________________________________________
2782 void AliESDtrack::SetITSdEdxSamples(const Double_t s[4]) {
2783   //
2784   // Store the dE/dx samples measured by the two SSD and two SDD layers.
2785   // These samples are corrected for the track segment length.
2786   //
2787   for (Int_t i=0; i<4; i++) fITSdEdxSamples[i]=s[i];
2788 }
2789
2790 //_______________________________________________________________________
2791 void AliESDtrack::GetITSdEdxSamples(Double_t s[4]) const {
2792   //
2793   // Get the dE/dx samples measured by the two SSD and two SDD layers.
2794   // These samples are corrected for the track segment length.
2795   //
2796   for (Int_t i=0; i<4; i++) s[i]=fITSdEdxSamples[i];
2797 }
2798
2799
2800 UShort_t   AliESDtrack::GetTPCnclsS(Int_t i0,Int_t i1) const{
2801   //
2802   // get number of shared TPC clusters
2803   //
2804   return  fTPCSharedMap.CountBits(i0)-fTPCSharedMap.CountBits(i1);
2805 }
2806
2807 UShort_t   AliESDtrack::GetTPCncls(Int_t i0,Int_t i1) const{
2808   //
2809   // get number of TPC clusters
2810   //
2811   return  fTPCClusterMap.CountBits(i0)-fTPCClusterMap.CountBits(i1);
2812 }
2813
2814 //____________________________________________________________________
2815 Double_t AliESDtrack::GetChi2TPCConstrainedVsGlobal(const AliESDVertex* vtx) const
2816 {
2817   // Calculates the chi2 between the TPC track (TPCinner) constrained to the primary vertex and the global track
2818   //
2819   // Returns -1 in case the calculation failed
2820   //
2821   // Value is cached as a non-persistent member.
2822   //
2823   // Code adapted from original code by GSI group (Jacek, Marian, Michael)
2824
2825   // cache, ignoring that a different vertex might be passed
2826   if (fCacheChi2TPCConstrainedVsGlobalVertex == vtx)
2827     return fCacheChi2TPCConstrainedVsGlobal;
2828
2829   fCacheChi2TPCConstrainedVsGlobal = -1;
2830   fCacheChi2TPCConstrainedVsGlobalVertex = vtx;
2831
2832   Double_t x[3];
2833   GetXYZ(x);
2834   Double_t b[3];
2835   AliTrackerBase::GetBxByBz(x,b);
2836
2837   if (!fTPCInner)  {
2838     AliWarning("Could not get TPC Inner Param.");
2839     return fCacheChi2TPCConstrainedVsGlobal;
2840   }
2841
2842   // clone for constraining
2843   AliExternalTrackParam* tpcInnerC = new AliExternalTrackParam(*fTPCInner);
2844   if (!tpcInnerC) {
2845     AliWarning("Clone of TPCInnerParam failed.");
2846     return fCacheChi2TPCConstrainedVsGlobal;
2847   }
2848
2849   // transform to the track reference frame
2850   Bool_t isOK = tpcInnerC->Rotate(GetAlpha());
2851   isOK &= tpcInnerC->PropagateTo(GetX(), b[2]);
2852   if (!isOK) {
2853     delete tpcInnerC;
2854     tpcInnerC = 0;
2855     AliWarning("Rotation/Propagation of track failed.") ;
2856     return fCacheChi2TPCConstrainedVsGlobal;
2857   }
2858
2859   // constrain TPCinner
2860   isOK = tpcInnerC->ConstrainToVertex(vtx, b);
2861
2862   // transform to the track reference frame
2863   isOK &= tpcInnerC->Rotate(GetAlpha());
2864   isOK &= tpcInnerC->PropagateTo(GetX(), b[2]);
2865
2866   if (!isOK) {
2867     AliWarning("ConstrainTPCInner failed.") ;
2868     delete tpcInnerC;
2869     tpcInnerC = 0;
2870     return fCacheChi2TPCConstrainedVsGlobal;
2871   }
2872
2873   // calculate chi2 between vi and vj vectors
2874   // with covi and covj covariance matrices
2875   // chi2ij = (vi-vj)^(T)*(covi+covj)^(-1)*(vi-vj)
2876   TMatrixD deltaT(5,1);
2877   TMatrixD delta(1,5);
2878   TMatrixD covarM(5,5);
2879
2880   for (Int_t ipar=0; ipar<5; ipar++) {
2881     deltaT(ipar,0) = tpcInnerC->GetParameter()[ipar] - GetParameter()[ipar];
2882     delta(0,ipar) = tpcInnerC->GetParameter()[ipar] - GetParameter()[ipar];
2883
2884     for (Int_t jpar=0; jpar<5; jpar++) {
2885       Int_t index = GetIndex(ipar,jpar);
2886       covarM(ipar,jpar) = GetCovariance()[index]+tpcInnerC->GetCovariance()[index];
2887     }
2888   }
2889   // chi2 distance TPC constrained and TPC+ITS
2890   TMatrixD covarMInv = covarM.Invert();
2891   TMatrixD mat2 = covarMInv*deltaT;
2892   TMatrixD chi2 = delta*mat2;
2893
2894   delete tpcInnerC;
2895   tpcInnerC = 0;
2896
2897   fCacheChi2TPCConstrainedVsGlobal = chi2(0,0);
2898   return fCacheChi2TPCConstrainedVsGlobal;
2899 }
2900
2901 void AliESDtrack::SetDetectorPID(const AliDetectorPID *pid)
2902 {
2903   //
2904   // Set the detector PID
2905   //
2906   if (fDetectorPID) delete fDetectorPID;
2907   fDetectorPID=pid;
2908
2909 }
2910
2911 Double_t AliESDtrack::GetLengthInActiveZone( Int_t mode, Double_t deltaY, Double_t deltaZ, Double_t bz, Double_t exbPhi , TTreeSRedirector * pcstream) const {
2912   //
2913   // Input parameters:
2914   //   mode  - type of external track parameters
2915   //   deltaY - user defined "dead region" in cm
2916   //   deltaZ - user defined "active region" in cm (250 cm drift lenght - 14 cm L1 delay
2917   //   bz     - magnetic field
2918   //   exbPhi - optional rotation due to the ExB effect
2919   // return value:
2920   //   the length of the track in cm in "active volume" of the TPC
2921   //
2922   if (mode==0) return GetLengthInActiveZone(this, deltaY,deltaZ,bz, exbPhi,pcstream);
2923   if (mode==1) return GetLengthInActiveZone(fIp, deltaY,deltaZ,bz, exbPhi,pcstream);
2924   if (mode==2) return GetLengthInActiveZone(fOp, deltaY,deltaZ,bz, exbPhi,pcstream);
2925   return 0;
2926 }
2927
2928 Double_t AliESDtrack::GetLengthInActiveZone(const AliExternalTrackParam  *paramT, Double_t deltaY, Double_t deltaZ, Double_t bz, Double_t exbPhi , TTreeSRedirector * pcstream) {
2929   //
2930   // Numerical code to calculate the length of the track in active region of the TPC
2931   // ( can be speed up if somebody wants to invest time - analysical version shoult be possible)
2932   //
2933   // Input parameters:
2934   //   paramT - external track parameters
2935   //   deltaY - user defined "dead region" in cm
2936   //   deltaZ - user defined "active region" in cm (250 cm drift lenght - 14 cm L1 delay
2937   //   bz     - magnetic field
2938   //   exbPhi - optional rotation due to the ExB effect
2939   // return value:
2940   //   the length of the track in cm in "active volume" of the TPC
2941   //
2942   const Double_t rIn=85;
2943   const Double_t rOut=245;
2944   Double_t xyz[3], pxyz[3];
2945   if (paramT->GetXYZAt(rIn,bz,xyz)){
2946     paramT->GetPxPyPzAt(rIn,bz,pxyz);
2947   }else{
2948     paramT->GetXYZ(xyz);
2949     paramT->GetPxPyPz(pxyz);
2950   }
2951   //
2952   Double_t dca   = -paramT->GetD(0,0,bz);  // get impact parameter distance to point (0,0)
2953   Double_t radius= TMath::Abs(1/paramT->GetC(bz));  //
2954   Double_t sign  = paramT->GetSign();
2955   Double_t R0    = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);   // radius at current point
2956   Double_t phiR0 = TMath::ATan2(xyz[1],xyz[0]);                // angle of given point
2957   Double_t dPhiR0= -TMath::ASin((dca*dca-2*dca*radius*sign+R0*R0)/(2*R0*(dca-radius*sign)));
2958   Double_t phi0  = phiR0-(dPhiR0);  // global phi offset to be added
2959   //
2960   //
2961   AliExternalTrackParam paramR=(*paramT);
2962   Double_t length=0;
2963   for (Double_t R=rIn; R<=rOut; R++){
2964     Double_t sinPhi=(dca*dca-2*dca*radius*sign+R*R)/(2*R*(dca-radius*sign));
2965     if (TMath::Abs(sinPhi)>=1) continue;
2966     Double_t dphi     = -TMath::ASin(sinPhi);
2967     Double_t phi      = phi0+dphi;                           // global phi
2968     Int_t    sector   = TMath::Nint(9*phi/(TMath::Pi()));
2969     Double_t dPhiEdge = phi-(sector*TMath::Pi()/9)+exbPhi;   // distance to sector boundary in rphi
2970     Double_t dX   = R*TMath::Cos(phi)-xyz[0];
2971     Double_t dY   = R*TMath::Sin(phi)-xyz[1];
2972     Double_t deltaPhi = 2*TMath::ASin(0.5*TMath::Sqrt(dX*dX+dY*dY)/radius);
2973     Double_t z = xyz[2]+deltaPhi*radius*paramT->GetTgl();
2974     if (TMath::Abs(dPhiEdge*R)>deltaY && TMath::Abs(z)<deltaZ){
2975       length++;
2976     }
2977     //    Double_t deltaZ= dphi*radius;
2978     if (pcstream){
2979       //should we keep debug possibility ?
2980       AliExternalTrackParam paramTcopy=(*paramT);
2981       paramR.Rotate(phi);
2982       paramR.PropagateTo(R,bz);
2983       (*pcstream)<<"debugEdge"<<
2984         "R="<<R<<                   // radius
2985         "dphiEdge="<<dPhiEdge<<     // distance to edge
2986         "phi0="<<phi0<<             // phi0 -phi at the track initial position
2987         "phi="<<phi<<               //
2988         "z="<<z<<
2989         "pT.="<<&paramTcopy<<
2990         "pR.="<<&paramR<<
2991         "\n";
2992     }
2993   }
2994   return length;
2995 }
2996
2997 Double_t AliESDtrack::GetMassForTracking() const
2998 {
2999   int pid = fPIDForTracking;
3000   if (pid<AliPID::kPion) pid = AliPID::kPion;
3001   double m = AliPID::ParticleMass(pid);
3002   return (fPIDForTracking==AliPID::kHe3 || fPIDForTracking==AliPID::kAlpha) ? -m : m;
3003 }
3004
3005
3006 void    AliESDtrack::SetTOFclusterArray(Int_t /*ncluster*/,Int_t */*TOFcluster*/){
3007   AliInfo("Method has to be implemented!");
3008 //   fNtofClusters=ncluster;
3009 //   if(TOFcluster == fTOFcluster) return;
3010 //   if(fTOFcluster){ // reset previous content
3011 //     delete[] fTOFcluster;
3012 //     fTOFcluster = NULL;
3013 //     fNtofClusters=0;
3014 //   }
3015
3016 //   if(ncluster){ // set new content
3017 //     fTOFcluster = new Int_t[fNtofClusters];
3018 //     for(Int_t i=0;i < fNtofClusters;i++) fTOFcluster[i] = TOFcluster[i];
3019 //   }
3020 //   else
3021 //     fTOFcluster = 0;
3022 }
3023
3024 //____________________________________________
3025 void  AliESDtrack::SuppressTOFMatches()
3026 {
3027   // remove reference to this track from TOF clusters
3028   if (!fNtofClusters || !GetESDEvent()) return;
3029   TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3030   for (;fNtofClusters--;) {
3031     AliESDTOFCluster* clTOF = (AliESDTOFCluster*)tofclArray->At(fTOFcluster[fNtofClusters]);
3032     clTOF->SuppressMatchedTrack(GetID());
3033     if (!clTOF->GetNMatchableTracks()) { // remove this cluster
3034       int last = tofclArray->GetEntriesFast()-1;
3035       AliESDTOFCluster* clTOFL = (AliESDTOFCluster*)tofclArray->At(last);
3036       if (last != fTOFcluster[fNtofClusters]) {
3037         *clTOF = *clTOFL; // move last cluster to the place of eliminated one
3038         // fix the references on this cluster
3039         clTOF->FixSelfReferences(last,fTOFcluster[fNtofClusters]);
3040       }
3041       tofclArray->RemoveAt(last);
3042     }
3043   }
3044 }
3045
3046 //____________________________________________
3047 void  AliESDtrack::ReplaceTOFTrackID(int oldID, int newID)
3048 {
3049   // replace the ID in TOF clusters references to this track
3050   if (!fNtofClusters || !GetESDEvent()) return;
3051   TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3052   for (int it=fNtofClusters;it--;) {
3053     AliESDTOFCluster* clTOF = (AliESDTOFCluster*)tofclArray->At(fTOFcluster[it]);
3054     clTOF->ReplaceMatchedTrackID(oldID,newID);
3055   }
3056 }
3057
3058 //____________________________________________
3059 void  AliESDtrack::ReplaceTOFClusterID(int oldID, int newID)
3060 {
3061   // replace the referenc on TOF cluster oldID by newID
3062   if (!fNtofClusters || !GetESDEvent()) return;
3063   for (int it=fNtofClusters;it--;) {
3064     if (fTOFcluster[it] == oldID) {
3065       fTOFcluster[it] = newID;
3066       return;
3067     }
3068   }
3069 }
3070
3071 //____________________________________________
3072 void  AliESDtrack::ReplaceTOFMatchID(int oldID, int newID)
3073 {
3074   // replace in the ESDTOFCluster associated with this track the id of the corresponding
3075   // ESDTOFMatch from oldID to newID
3076   if (!fNtofClusters || !GetESDEvent()) return;
3077   TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3078   for (int it=fNtofClusters;it--;) {
3079     AliESDTOFCluster* clTOF = (AliESDTOFCluster*)tofclArray->At(fTOFcluster[it]);
3080     clTOF->ReplaceMatchID(oldID,newID);
3081   }
3082 }
3083
3084 //____________________________________________
3085 void AliESDtrack::AddTOFcluster(Int_t icl)
3086 {
3087   fNtofClusters++;
3088
3089   Int_t *old = fTOFcluster;
3090   fTOFcluster = new Int_t[fNtofClusters];
3091
3092   for(Int_t i=0;i < fNtofClusters-1;i++) fTOFcluster[i] = old[i];
3093   fTOFcluster[fNtofClusters-1] = icl;
3094
3095   if(fNtofClusters-1)  delete[] old; // delete previous content
3096
3097 }
3098
3099 //____________________________________________
3100 void AliESDtrack::SetTOFsignal(Double_t tof)
3101 {
3102   if(fNtofClusters>0 && GetESDEvent()){
3103     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3104     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3105     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
3106     if(hit) hit->SetTime(tof);
3107   }
3108   else{
3109     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3110     fTOFsignal=tof;
3111   }
3112 }
3113 //____________________________________________
3114 void AliESDtrack::SetTOFCalChannel(Int_t index){
3115   if(fNtofClusters>0 && GetESDEvent()){
3116     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3117     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3118     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
3119     if(hit) hit->SetTOFchannel(index);
3120   }
3121   else{
3122     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3123     fTOFCalChannel=index;
3124   }
3125 }
3126 //____________________________________________
3127 void AliESDtrack::SetTOFsignalToT(Double_t ToT){
3128   if(fNtofClusters>0 && GetESDEvent()){
3129     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3130     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3131     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
3132     if(hit) hit->SetTOT(ToT);
3133   }
3134   else{
3135     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3136     fTOFsignalToT=ToT;
3137   }
3138 }
3139 //____________________________________________
3140 void AliESDtrack::SetTOFsignalRaw(Double_t tof){
3141   if(fNtofClusters>0 && GetESDEvent()){
3142     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3143     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3144     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
3145     if(hit) hit->SetTimeRaw(tof);
3146   }
3147   else{
3148     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3149     fTOFsignalRaw=tof;
3150   }
3151 }
3152 //____________________________________________
3153 void AliESDtrack::SetTOFsignalDz(Double_t dz){
3154   Int_t index = -1;
3155   AliESDTOFCluster *tofcl;
3156
3157   if(fNtofClusters>0 && GetESDEvent()){
3158     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3159     tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3160
3161     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
3162       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
3163     }
3164
3165   }
3166   if(index > -1){
3167     AliESDTOFMatch* match = tofcl->GetTOFMatch(index);
3168     if(match){
3169       match->SetDz(dz);
3170     }
3171   }
3172   else{
3173     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3174     fTOFsignalDz=dz;
3175   }
3176
3177
3178 }
3179 //____________________________________________
3180 void AliESDtrack::SetTOFsignalDx(Double_t dx){
3181   Int_t index = -1;
3182   AliESDTOFCluster *tofcl;
3183
3184   if(fNtofClusters>0 && GetESDEvent()){
3185     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3186     tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3187
3188     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
3189       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
3190     }
3191
3192   }
3193   if(index > -1){
3194     AliESDTOFMatch* match = tofcl->GetTOFMatch(index);
3195     if(match){
3196       match->SetDx(dx);
3197     }
3198   }
3199   else{
3200     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3201     fTOFsignalDx=dx;
3202   }
3203 }
3204 //____________________________________________
3205 void AliESDtrack::SetTOFDeltaBC(Short_t deltaBC){
3206   if(fNtofClusters>0 && GetESDEvent()){
3207     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3208     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3209     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
3210     if(hit) hit->SetDeltaBC(deltaBC);
3211   }
3212   else{
3213     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3214     fTOFdeltaBC=deltaBC;
3215   }
3216 }
3217 //____________________________________________
3218 void AliESDtrack::SetTOFL0L1(Short_t l0l1){
3219   if(fNtofClusters>0 && GetESDEvent()){
3220     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3221     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3222     AliESDTOFHit* hit = tofcl->GetTOFHit(0);
3223     if(hit) hit->SetL0L1Latency(l0l1);
3224   }
3225   else{
3226     if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3227     fTOFl0l1=l0l1;
3228   }
3229 }
3230 //____________________________________________
3231 Double_t AliESDtrack::GetTOFsignal() const
3232 {
3233   if(fNtofClusters>0 && GetESDEvent()){
3234     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3235     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3236
3237     return tofcl->GetTime();
3238   }
3239   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3240
3241   return fTOFsignal;
3242 }
3243
3244 //____________________________________________
3245 Double_t AliESDtrack::GetTOFsignalToT() const
3246 {
3247   if(fNtofClusters>0 && GetESDEvent()){
3248     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3249     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3250
3251     return tofcl->GetTOT();
3252   }
3253   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3254
3255   return fTOFsignalToT;
3256 }
3257
3258 //____________________________________________
3259 Double_t AliESDtrack::GetTOFsignalRaw() const
3260 {
3261   if(fNtofClusters>0 && GetESDEvent()){
3262     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3263     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3264
3265     return tofcl->GetTimeRaw();
3266   }
3267   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3268
3269   return fTOFsignalRaw;
3270 }
3271
3272 //____________________________________________
3273 Double_t AliESDtrack::GetTOFsignalDz() const
3274 {
3275
3276   AliESDTOFCluster *tofcl;
3277
3278   Int_t index = -1;
3279   if(fNtofClusters>0 && GetESDEvent()){
3280     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3281     tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3282
3283     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
3284       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
3285     }
3286   }
3287   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3288
3289   if(fNtofClusters>0 && index > -1){
3290     return tofcl->GetDz(index);
3291   }
3292   return fTOFsignalDz;
3293 }
3294
3295 //____________________________________________
3296 Double_t AliESDtrack::GetTOFsignalDx() const
3297 {
3298   AliESDTOFCluster *tofcl;
3299
3300   Int_t index = -1;
3301   if(fNtofClusters>0 && GetESDEvent()){
3302     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3303     tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3304     for(Int_t i=0;i < tofcl->GetNMatchableTracks();i++){
3305       if(tofcl->GetTrackIndex(i) == GetID()) index = i;
3306     }
3307   }
3308   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3309   if(fNtofClusters>0 && index > -1){
3310     return tofcl->GetDx(index);
3311   }
3312   return fTOFsignalDx;
3313 }
3314
3315 //____________________________________________
3316 Short_t  AliESDtrack::GetTOFDeltaBC() const
3317 {
3318   if(fNtofClusters>0 && GetESDEvent()){
3319     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3320     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3321     return tofcl->GetDeltaBC();
3322   }
3323   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3324
3325   return fTOFdeltaBC;
3326 }
3327
3328 //____________________________________________
3329 Short_t  AliESDtrack::GetTOFL0L1() const
3330 {
3331   if(fNtofClusters>0 && GetESDEvent()){
3332     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3333     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3334
3335     return tofcl->GetL0L1Latency();
3336   }
3337   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3338
3339   return fTOFl0l1;
3340 }
3341
3342 //____________________________________________
3343 Int_t   AliESDtrack::GetTOFCalChannel() const
3344 {
3345   if(fNtofClusters>0 && GetESDEvent()){
3346     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3347     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3348
3349     return tofcl->GetTOFchannel();
3350   }
3351   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3352
3353   return fTOFCalChannel;
3354 }
3355
3356 //____________________________________________
3357 Int_t   AliESDtrack::GetTOFcluster() const
3358 {
3359   if(fNtofClusters>0 && GetESDEvent()){
3360     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3361     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3362
3363     return tofcl->GetClusterIndex();
3364   }
3365   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3366
3367   return fTOFindex;
3368 }
3369
3370 //____________________________________________
3371 Int_t   AliESDtrack::GetTOFclusterN() const
3372 {
3373   return fNtofClusters;
3374 }
3375
3376 //____________________________________________
3377 Bool_t  AliESDtrack::IsTOFHitAlreadyMatched() const{
3378   if(fNtofClusters>0 && GetESDEvent()){
3379     TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3380     AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[0]);
3381
3382     if (tofcl->GetNMatchableTracks() > 1)
3383       return kTRUE;
3384   }
3385   else if(fNtofClusters>0) AliInfo("No AliESDEvent available here!\n");
3386
3387   return kFALSE;
3388 }
3389
3390 //____________________________________________
3391 void AliESDtrack::ReMapTOFcluster(Int_t ncl,Int_t *mapping){
3392   for(Int_t i=0;i<fNtofClusters;i++){
3393     if(fTOFcluster[i]<ncl && fTOFcluster[i]>-1)
3394       fTOFcluster[i] = mapping[fTOFcluster[i]];
3395     else
3396       AliInfo(Form("TOF cluster re-mapping in AliESDtrack: out of range (%i > %i)\n",fTOFcluster[i],ncl));
3397   }
3398 }
3399
3400 //____________________________________________
3401 void AliESDtrack::SortTOFcluster(){
3402   TClonesArray *tofclArray = GetESDEvent()->GetESDTOFClusters();
3403
3404   for(Int_t i=0;i<fNtofClusters-1;i++){
3405     for(Int_t j=i+1;j<fNtofClusters;j++){
3406       AliESDTOFCluster *tofcl = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[i]);
3407       Int_t index1 = -1;
3408       for(Int_t it=0;it < tofcl->GetNMatchableTracks();it++){
3409          if(tofcl->GetTrackIndex(it) == GetID()) index1 = it;
3410       }
3411       Double_t timedist1 = 10000;
3412       for(Int_t isp=0; isp< AliPID::kSPECIESC;isp++){
3413         Double_t timec = TMath::Abs(tofcl->GetTime() - tofcl->GetIntegratedTime(isp));
3414         if(timec < timedist1) timedist1 = timec;
3415       }
3416       timedist1 *= 0.03; // in cm
3417       Double_t radius1 = tofcl->GetDx(index1)*tofcl->GetDx(index1) + tofcl->GetDz(index1)*tofcl->GetDz(index1) + timedist1*timedist1;
3418
3419       AliESDTOFCluster *tofcl2 = (AliESDTOFCluster *) tofclArray->At(fTOFcluster[j]);
3420       Int_t index2 = -1;
3421       for(Int_t it=0;it < tofcl2->GetNMatchableTracks();it++){
3422          if(tofcl2->GetTrackIndex(it) == GetID()) index2 = it;
3423       }
3424       if(index1 == -1 || index2 == -1){
3425       }
3426       Double_t timedist2 = 10000;
3427       for(Int_t isp=0; isp< AliPID::kSPECIESC;isp++){
3428         Double_t timec = TMath::Abs(tofcl2->GetTime() - tofcl2->GetIntegratedTime(isp));
3429         if(timec < timedist2) timedist2 = timec;
3430       }
3431       timedist2 *= 0.03; // in cm
3432       Double_t radius2 = tofcl2->GetDx(index2)*tofcl2->GetDx(index2) + tofcl2->GetDz(index2)*tofcl2->GetDz(index2) + timedist2*timedist2;
3433
3434       if(radius2 < radius1){
3435         Int_t change = fTOFcluster[i];
3436         fTOFcluster[i] = fTOFcluster[j];
3437         fTOFcluster[j] = change;
3438       }
3439     }
3440   }
3441 }
3442
3443 //____________________________________________
3444 const AliTOFHeader* AliESDtrack::GetTOFHeader() const {
3445   return fESDEvent->GetTOFHeader();
3446 }
3447
3448 //___________________________________________
3449 void AliESDtrack::SetID(Short_t id)
3450 {
3451   // set track ID taking care about dependencies
3452   if (fNtofClusters) ReplaceTOFTrackID(fID,id);
3453   fID=id;
3454 }
3455
3456
3457 Double_t  AliESDtrack::GetdEdxInfo(Int_t regionID, Int_t calibID, Int_t qID, Int_t valueID){
3458   //
3459   // Interface to get the calibrated dEdx information
3460   // For details of arguments and return values see
3461   //     AliTPCdEdxInfo::GetdEdxInfo(Int_t regionID, Int_t calibID, Int_t valueID)
3462   //
3463   if (!fTPCdEdxInfo) return 0;
3464   if (!fIp) return 0;
3465   return fTPCdEdxInfo->GetdEdxInfo(fIp, regionID, calibID, qID, valueID);
3466 }
3467
3468
3469 Double_t AliESDtrack::GetdEdxInfoTRD(Int_t method, Double_t p0, Double_t p1, Double_t p2){
3470   //
3471   // Methods
3472   // mean values:
3473   //     0.)   linear
3474   //     1.)   logarithmic
3475   //     2.)   1/sqrt
3476   //     3.)   power()
3477   // time COG:
3478   //     4.)   linear
3479   //     5.)   logarithmic
3480   //     6.)   square
3481   Int_t nSlicesPerLayer=GetNumberOfTRDslices();
3482   Int_t nSlicesAll=GetNumberOfTRDslices()*kTRDnPlanes;
3483
3484   if (method<=3){
3485     Double_t sumAmp=0;
3486     Int_t    sumW=0;
3487     for (Int_t ibin=0; ibin<nSlicesAll; ibin++){
3488       if (fTRDslices[ibin]<=0) continue;
3489       sumW++;
3490       if (method==0) sumAmp+=fTRDslices[ibin];
3491       if (method==1) sumAmp+=TMath::Log(TMath::Abs(fTRDslices[ibin])+p0);
3492       if (method==2) sumAmp+=1/TMath::Sqrt(TMath::Abs(fTRDslices[ibin])+p0);
3493       if (method==3) sumAmp+=TMath::Power(TMath::Abs(fTRDslices[ibin])+p0,p1);
3494     }
3495     if (sumW==0) return 0;
3496     Double_t dEdx=sumAmp/sumW;
3497     if (method==1) dEdx= TMath::Exp(dEdx);
3498     if (method==2) dEdx= 1/(dEdx*dEdx);
3499     if (method==3) dEdx= TMath::Power(dEdx,1/p1);
3500     return dEdx;
3501   }
3502   if (method>3){
3503     Double_t sumWT=0;
3504     Double_t sumW=0;
3505     for (Int_t ibin=0; ibin<nSlicesAll; ibin++){
3506       if (fTRDslices[ibin]<=0) continue;
3507       Double_t time=(ibin%nSlicesPerLayer);
3508       Double_t weight=fTRDslices[ibin];
3509       if (method==5) weight=TMath::Log((weight+p0)/p0);
3510       if (method==6) weight=TMath::Power(weight+p0,p1);
3511       sumWT+=time*weight;
3512       sumW+=weight;
3513     }
3514     if (sumW<=0) return 0;
3515     Double_t meanTime=sumWT/sumW;
3516     return meanTime;
3517   }
3518   return 0;
3519 }