1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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
24 // What do you need to know before starting analysis
25 // (by Marian Ivanov: marian.ivanov@cern.ch)
29 // 1. What is the AliESDtrack
30 // 2. What informations do we store
31 // 3. How to use the information for analysis
34 // 1.AliESDtrack is the container of the information about the track/particle
35 // reconstructed during Barrel Tracking.
36 // The track information is propagated from one tracking detector to
37 // other using the functionality of AliESDtrack - Current parameters.
39 // No global fit model is used.
40 // Barrel tracking use Kalman filtering technique, it gives optimal local
41 // track parameters at given point under certian assumptions.
43 // Kalman filter take into account additional effect which are
44 // difficult to handle using global fit.
46 // a.) Multiple scattering
48 // c.) Non homogenous magnetic field
50 // In general case, following barrel detectors are contributing to
51 // the Kalman track information:
56 // In general 3 reconstruction itteration are performed:
57 // 1. Find tracks - sequence TPC->ITS
58 // 2. PropagateBack - sequence ITS->TPC->TRD -> Outer PID detectors
59 // 3. Refit invward - sequence TRD->TPC->ITS
60 // The current tracks are updated after each detector (see bellow).
61 // In specical cases a track sanpshots are stored.
64 // For some type of analysis (+visualization) track local parameters at
65 // different position are neccesary. A snapshots during the track
66 // propagation are created.
67 // (See AliExternalTrackParam class for desctiption of variables and
70 // a. Current parameters - class itself (AliExternalTrackParam)
71 // Contributors: general case TRD->TPC->ITS
72 // Preferable usage: Decission - primary or secondary track
73 // NOTICE - By default the track parameters are stored at the DCA point
74 // to the primary vertex. optimal for primary tracks,
75 // far from optimal for secondary tracks.
76 // b. Constrained parameters - Kalman information updated with
77 // the Primary vertex information
78 // Contributors: general case TRD->TPC->ITS
79 // Preferable usage: Use only for tracks selected as primary
80 // NOTICE - not real constrain - taken as additional measurement
81 // with corresponding error
83 // const AliExternalTrackParam *GetConstrainedParam() const {return fCp;}
84 // c. Inner parameters - Track parameters at inner wall of the TPC
85 // Contributors: general case TRD->TPC
87 // const AliExternalTrackParam *GetInnerParam() const { return fIp;}
89 // d. TPCinnerparam - contributors - TPC only
91 // Preferable usage: Requested for HBT study
92 // (smaller correlations as using also ITS information)
93 // NOTICE - the track parameters are propagated to the DCA to
95 // Optimal for primary, far from optimal for secondary tracks
97 // const AliExternalTrackParam *GetTPCInnerParam() const {return fTPCInner;}
99 // e. Outer parameters -
100 // Contributors- general case - ITS-> TPC -> TRD
101 // The last point - Outer parameters radius is determined
102 // e.a) Local inclination angle bigger than threshold -
103 // Low momenta tracks
104 // e.a) Catastrofic energy losss in material
105 // e.b) Not further improvement (no space points)
107 // a.) Tracking: Starting parameter for Refit inward
110 // NOTICE: Should be not used for the physic analysis
112 // const AliExternalTrackParam *GetOuterParam() const { return fOp;}
114 //-----------------------------------------------------------------
117 #include <TParticle.h>
119 #include "AliESDVertex.h"
120 #include "AliESDtrack.h"
121 #include "AliKalmanTrack.h"
123 #include "AliTrackPointArray.h"
124 #include "TPolyMarker3D.h"
126 ClassImp(AliESDtrack)
128 void SetPIDValues(Double_t * dest, const Double_t * src, Int_t n) {
129 // This function copies "n" PID weights from "scr" to "dest"
130 // and normalizes their sum to 1 thus producing conditional probabilities.
131 // The negative weights are set to 0.
132 // In case all the weights are non-positive they are replaced by
133 // uniform probabilities
137 Float_t uniform = 1./(Float_t)n;
140 for (Int_t i=0; i<n; i++)
150 for (Int_t i=0; i<n; i++) dest[i] /= sum;
152 for (Int_t i=0; i<n; i++) dest[i] = uniform;
155 //_______________________________________________________________________
156 AliESDtrack::AliESDtrack() :
157 AliExternalTrackParam(),
162 fFriendTrack(new AliESDfriendTrack()),
163 fTPCClusterMap(159),//number of padrows
164 fTPCSharedMap(159),//number of padrows
175 fEMCALindex(kEMCALNoMatch),
181 fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
183 fCdd(0),fCdz(0),fCzz(0),
217 // The default ESD constructor
220 for (i=0; i<AliPID::kSPECIES; i++) {
230 for (i=0; i<3; i++) { fKinkIndexes[i]=0;}
231 for (i=0; i<3; i++) { fV0Indexes[i]=0;}
232 for (i=0;i<kTRDnPlanes;i++) {
235 for (i=0;i<4;i++) {fTPCPoints[i]=0;}
236 for (i=0;i<3;i++) {fTOFLabel[i]=0;}
237 for (i=0;i<10;i++) {fTOFInfo[i]=0;}
238 for (i=0;i<12;i++) {fITSModule[i]=-1;}
241 //_______________________________________________________________________
242 AliESDtrack::AliESDtrack(const AliESDtrack& track):
243 AliExternalTrackParam(track),
249 fTPCClusterMap(track.fTPCClusterMap),
250 fTPCSharedMap(track.fTPCSharedMap),
251 fFlags(track.fFlags),
253 fLabel(track.fLabel),
254 fITSLabel(track.fITSLabel),
255 fTPCLabel(track.fTPCLabel),
256 fTRDLabel(track.fTRDLabel),
257 fTOFCalChannel(track.fTOFCalChannel),
258 fTOFindex(track.fTOFindex),
259 fHMPIDqn(track.fHMPIDqn),
260 fHMPIDcluIdx(track.fHMPIDcluIdx),
261 fEMCALindex(track.fEMCALindex),
262 fHMPIDtrkTheta(track.fHMPIDtrkTheta),
263 fHMPIDtrkPhi(track.fHMPIDtrkPhi),
264 fHMPIDsignal(track.fHMPIDsignal),
265 fTrackLength(track.fTrackLength),
266 fdTPC(track.fdTPC),fzTPC(track.fzTPC),
267 fCddTPC(track.fCddTPC),fCdzTPC(track.fCdzTPC),fCzzTPC(track.fCzzTPC),
268 fD(track.fD),fZ(track.fZ),
269 fCdd(track.fCdd),fCdz(track.fCdz),fCzz(track.fCzz),
270 fCchi2(track.fCchi2),
271 fITSchi2(track.fITSchi2),
272 fTPCchi2(track.fTPCchi2),
273 fTRDchi2(track.fTRDchi2),
274 fTOFchi2(track.fTOFchi2),
275 fHMPIDchi2(track.fHMPIDchi2),
276 fITSsignal(track.fITSsignal),
277 fTPCsignal(track.fTPCsignal),
278 fTPCsignalS(track.fTPCsignalS),
279 fTRDsignal(track.fTRDsignal),
280 fTRDQuality(track.fTRDQuality),
281 fTRDBudget(track.fTRDBudget),
282 fTOFsignal(track.fTOFsignal),
283 fTOFsignalToT(track.fTOFsignalToT),
284 fTOFsignalRaw(track.fTOFsignalRaw),
285 fTOFsignalDz(track.fTOFsignalDz),
286 fHMPIDtrkX(track.fHMPIDtrkX),
287 fHMPIDtrkY(track.fHMPIDtrkY),
288 fHMPIDmipX(track.fHMPIDmipX),
289 fHMPIDmipY(track.fHMPIDmipY),
290 fTPCncls(track.fTPCncls),
291 fTPCnclsF(track.fTPCnclsF),
292 fTPCsignalN(track.fTPCsignalN),
293 fITSncls(track.fITSncls),
294 fITSClusterMap(track.fITSClusterMap),
295 fTRDncls(track.fTRDncls),
296 fTRDncls0(track.fTRDncls0),
297 fTRDpidQuality(track.fTRDpidQuality),
298 fTRDnSlices(track.fTRDnSlices),
304 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i]=track.fTrackTime[i];
305 for (Int_t i=0;i<AliPID::kSPECIES;i++) fR[i]=track.fR[i];
307 for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=track.fITSr[i];
309 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=track.fTPCr[i];
310 for (Int_t i=0;i<4;i++) {fTPCPoints[i]=track.fTPCPoints[i];}
311 for (Int_t i=0; i<3;i++) { fKinkIndexes[i]=track.fKinkIndexes[i];}
312 for (Int_t i=0; i<3;i++) { fV0Indexes[i]=track.fV0Indexes[i];}
314 for (Int_t i=0;i<kTRDnPlanes;i++) {
315 fTRDTimBin[i]=track.fTRDTimBin[i];
319 fTRDslices=new Double32_t[fTRDnSlices];
320 for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=track.fTRDslices[i];
323 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i]=track.fTRDr[i];
324 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
325 for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
326 for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
327 for (Int_t i=0;i<12;i++) fITSModule[i]=track.fITSModule[i];
328 for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i]=track.fHMPIDr[i];
330 if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
331 if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
332 if (track.fTPCInner) fTPCInner=new AliExternalTrackParam(*track.fTPCInner);
333 if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
335 if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
338 //_______________________________________________________________________
339 AliESDtrack::AliESDtrack(TParticle * part) :
340 AliExternalTrackParam(),
346 fTPCClusterMap(159),//number of padrows
347 fTPCSharedMap(159),//number of padrows
358 fEMCALindex(kEMCALNoMatch),
364 fCddTPC(0),fCdzTPC(0),fCzzTPC(0),
366 fCdd(0),fCdz(0),fCzz(0),
399 // ESD track from TParticle
402 // Reset all the arrays
404 for (i=0; i<AliPID::kSPECIES; i++) {
414 for (i=0; i<3; i++) { fKinkIndexes[i]=0;}
415 for (i=0; i<3; i++) { fV0Indexes[i]=-1;}
416 for (i=0;i<kTRDnPlanes;i++) {
419 for (i=0;i<4;i++) {fTPCPoints[i]=0;}
420 for (i=0;i<3;i++) {fTOFLabel[i]=0;}
421 for (i=0;i<10;i++) {fTOFInfo[i]=0;}
422 for (i=0;i<12;i++) {fITSModule[i]=-1;}
424 // Calculate the AliExternalTrackParam content
431 // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
432 alpha = part->Phi()*180./TMath::Pi();
433 if (alpha<0) alpha+= 360.;
434 if (alpha>360) alpha -= 360.;
436 Int_t sector = (Int_t)(alpha/20.);
437 alpha = 10. + 20.*sector;
439 alpha *= TMath::Pi();
441 // Covariance matrix: no errors, the parameters are exact
442 for (i=0; i<15; i++) covar[i]=0.;
444 // Get the vertex of origin and the momentum
445 TVector3 ver(part->Vx(),part->Vy(),part->Vz());
446 TVector3 mom(part->Px(),part->Py(),part->Pz());
448 // Rotate to the local coordinate system (TPC sector)
452 // X of the referense plane
455 Int_t pdgCode = part->GetPdgCode();
458 TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
462 param[2] = TMath::Sin(mom.Phi());
463 param[3] = mom.Pz()/mom.Pt();
464 param[4] = TMath::Sign(1/mom.Pt(),charge);
466 // Set AliExternalTrackParam
467 Set(xref, alpha, param, covar);
472 switch (TMath::Abs(pdgCode)) {
498 // If the particle is not e,mu,pi,K or p the PID probabilities are set to 0
499 if (indexPID < AliPID::kSPECIES) {
505 fHMPIDr[indexPID]=1.;
508 // AliESD track label
509 SetLabel(part->GetUniqueID());
513 //_______________________________________________________________________
514 AliESDtrack::~AliESDtrack(){
516 // This is destructor according Coding Conventrions
518 //printf("Delete track\n");
528 AliESDtrack &AliESDtrack::operator=(const AliESDtrack &source){
531 if(&source == this) return *this;
532 AliExternalTrackParam::operator=(source);
536 // we have the trackparam: assign or copy construct
537 if(fCp)*fCp = *source.fCp;
538 else fCp = new AliExternalTrackParam(*source.fCp);
541 // no track param delete the old one
547 // we have the trackparam: assign or copy construct
548 if(fIp)*fIp = *source.fIp;
549 else fIp = new AliExternalTrackParam(*source.fIp);
552 // no track param delete the old one
558 if(source.fTPCInner){
559 // we have the trackparam: assign or copy construct
560 if(fTPCInner) *fTPCInner = *source.fTPCInner;
561 else fTPCInner = new AliExternalTrackParam(*source.fTPCInner);
564 // no track param delete the old one
565 if(fTPCInner)delete fTPCInner;
571 // we have the trackparam: assign or copy construct
572 if(fOp) *fOp = *source.fOp;
573 else fOp = new AliExternalTrackParam(*source.fOp);
576 // no track param delete the old one
581 // copy also the friend track
582 // use copy constructor
583 if(source.fFriendTrack){
584 // we have the trackparam: assign or copy construct
585 delete fFriendTrack; fFriendTrack=new AliESDfriendTrack(*source.fFriendTrack);
588 // no track param delete the old one
589 delete fFriendTrack; fFriendTrack= 0;
592 fTPCClusterMap = source.fTPCClusterMap;
593 fTPCSharedMap = source.fTPCSharedMap;
595 fFlags = source.fFlags;
597 fLabel = source.fLabel;
598 fITSLabel = source.fITSLabel;
599 for(int i = 0; i< 12;++i){
600 fITSModule[i] = source.fITSModule[i];
602 fTPCLabel = source.fTPCLabel;
603 fTRDLabel = source.fTRDLabel;
604 for(int i = 0; i< 3;++i){
605 fTOFLabel[i] = source.fTOFLabel[i];
607 fTOFCalChannel = source.fTOFCalChannel;
608 fTOFindex = source.fTOFindex;
609 fHMPIDqn = source.fHMPIDqn;
610 fHMPIDcluIdx = source.fHMPIDcluIdx;
611 fEMCALindex = source.fEMCALindex;
613 for(int i = 0; i< 3;++i){
614 fKinkIndexes[i] = source.fKinkIndexes[i];
615 fV0Indexes[i] = source.fV0Indexes[i];
618 for(int i = 0; i< AliPID::kSPECIES;++i){
619 fR[i] = source.fR[i];
620 fITSr[i] = source.fITSr[i];
621 fTPCr[i] = source.fTPCr[i];
622 fTRDr[i] = source.fTRDr[i];
623 fTOFr[i] = source.fTOFr[i];
624 fHMPIDr[i] = source.fHMPIDr[i];
625 fTrackTime[i] = source.fTrackTime[i];
628 fHMPIDtrkTheta = source.fHMPIDtrkTheta;
629 fHMPIDtrkPhi = source.fHMPIDtrkPhi;
630 fHMPIDsignal = source.fHMPIDsignal;
633 fTrackLength = source. fTrackLength;
634 fdTPC = source.fdTPC;
635 fzTPC = source.fzTPC;
636 fCddTPC = source.fCddTPC;
637 fCdzTPC = source.fCdzTPC;
638 fCzzTPC = source.fCzzTPC;
645 fCchi2 = source.fCchi2;
646 fITSchi2 = source.fITSchi2;
647 fTPCchi2 = source.fTPCchi2;
648 fTRDchi2 = source.fTRDchi2;
649 fTOFchi2 = source.fTOFchi2;
650 fHMPIDchi2 = source.fHMPIDchi2;
653 fITSsignal = source.fITSsignal;
654 fTPCsignal = source.fTPCsignal;
655 fTPCsignalS = source.fTPCsignalS;
656 for(int i = 0; i< 4;++i){
657 fTPCPoints[i] = source.fTPCPoints[i];
659 fTRDsignal = source.fTRDsignal;
661 for(int i = 0;i < kTRDnPlanes;++i){
662 fTRDTimBin[i] = source.fTRDTimBin[i];
668 fTRDnSlices=source.fTRDnSlices;
670 fTRDslices=new Double32_t[fTRDnSlices];
671 for(int j = 0;j < fTRDnSlices;++j) fTRDslices[j] = source.fTRDslices[j];
674 fTRDQuality = source.fTRDQuality;
675 fTRDBudget = source.fTRDBudget;
676 fTOFsignal = source.fTOFsignal;
677 fTOFsignalToT = source.fTOFsignalToT;
678 fTOFsignalRaw = source.fTOFsignalRaw;
679 fTOFsignalDz = source.fTOFsignalDz;
681 for(int i = 0;i<10;++i){
682 fTOFInfo[i] = source.fTOFInfo[i];
685 fHMPIDtrkX = source.fHMPIDtrkX;
686 fHMPIDtrkY = source.fHMPIDtrkY;
687 fHMPIDmipX = source.fHMPIDmipX;
688 fHMPIDmipY = source.fHMPIDmipY;
690 fTPCncls = source.fTPCncls;
691 fTPCnclsF = source.fTPCnclsF;
692 fTPCsignalN = source.fTPCsignalN;
694 fITSncls = source.fITSncls;
695 fITSClusterMap = source.fITSClusterMap;
696 fTRDncls = source.fTRDncls;
697 fTRDncls0 = source.fTRDncls0;
698 fTRDpidQuality = source.fTRDpidQuality;
704 void AliESDtrack::Copy(TObject &obj) const {
706 // this overwrites the virtual TOBject::Copy()
707 // to allow run time copying without casting
710 if(this==&obj)return;
711 AliESDtrack *robj = dynamic_cast<AliESDtrack*>(&obj);
712 if(!robj)return; // not an AliESDtrack
719 void AliESDtrack::AddCalibObject(TObject * object){
721 // add calib object to the list
723 if (!fFriendTrack) fFriendTrack = new AliESDfriendTrack;
724 fFriendTrack->AddCalibObject(object);
727 TObject * AliESDtrack::GetCalibObject(Int_t index){
729 // return calib objct at given position
731 if (!fFriendTrack) return 0;
732 return fFriendTrack->GetCalibObject(index);
736 const Bool_t AliESDtrack::FillTPCOnlyTrack(AliESDtrack &track){
738 // Fills the information of the TPC-only first reconstruction pass
739 // into the passed ESDtrack object. For consistency fTPCInner is also filled
742 if(!fTPCInner)return kFALSE;
744 // fill the TPC track params to the global track parameters
745 track.Set(fTPCInner->GetX(),fTPCInner->GetAlpha(),fTPCInner->GetParameter(),fTPCInner->GetCovariance());
748 track.fCdd = fCddTPC;
749 track.fCdz = fCdzTPC;
750 track.fCzz = fCzzTPC;
752 // copy the TPCinner parameters
753 if(track.fTPCInner) *track.fTPCInner = *fTPCInner;
754 else track.fTPCInner = new AliExternalTrackParam(*fTPCInner);
757 track.fCddTPC = fCddTPC;
758 track.fCdzTPC = fCdzTPC;
759 track.fCzzTPC = fCzzTPC;
762 // copy all other TPC specific parameters
764 // replace label by TPC label
765 track.fLabel = fTPCLabel;
766 track.fTPCLabel = fTPCLabel;
768 track.fTPCchi2 = fTPCchi2;
769 track.fTPCsignal = fTPCsignal;
770 track.fTPCsignalS = fTPCsignalS;
771 for(int i = 0;i<4;++i)track.fTPCPoints[i] = fTPCPoints[i];
773 track.fTPCncls = fTPCncls;
774 track.fTPCnclsF = fTPCnclsF;
775 track.fTPCsignalN = fTPCsignalN;
778 for(int i=0;i<AliPID::kSPECIES;++i){
779 track.fTPCr[i] = fTPCr[i];
780 // combined PID is TPC only!
781 track.fR[i] = fTPCr[i];
783 track.fTPCClusterMap = fTPCClusterMap;
784 track.fTPCSharedMap = fTPCSharedMap;
788 track.fFlags = kTPCin;
792 for (Int_t i=0;i<3;i++) track.fKinkIndexes[i] = fKinkIndexes[i];
798 //_______________________________________________________________________
799 void AliESDtrack::MakeMiniESDtrack(){
800 // Resets everything except
801 // fFlags: Reconstruction status flags
802 // fLabel: Track label
803 // fID: Unique ID of the track
804 // Impact parameter information
805 // fR[AliPID::kSPECIES]: combined "detector response probability"
806 // Running track parameters in the base class (AliExternalTrackParam)
810 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i] = 0;
812 // Reset track parameters constrained to the primary vertex
816 // Reset track parameters at the inner wall of TPC
818 delete fTPCInner;fTPCInner=0;
819 // Reset track parameters at the inner wall of the TRD
823 // Reset ITS track related information
828 for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0;
831 // Reset TPC related track information
840 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTPCr[i]=0;
842 for (Int_t i=0;i<4;i++) fTPCPoints[i] = 0;
843 for (Int_t i=0; i<3;i++) fKinkIndexes[i] = 0;
844 for (Int_t i=0; i<3;i++) fV0Indexes[i] = 0;
846 // Reset TRD related track information
851 for (Int_t i=0;i<kTRDnPlanes;i++) {
854 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTRDr[i] = 0;
864 // Reset TOF related track information
872 for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i] = 0;
873 for (Int_t i=0;i<3;i++) fTOFLabel[i] = 0;
874 for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
876 // Reset HMPID related track information
881 for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
888 fEMCALindex = kEMCALNoMatch;
890 delete fFriendTrack; fFriendTrack = 0;
892 //_______________________________________________________________________
893 Double_t AliESDtrack::GetMass() const {
894 // Returns the mass of the most probable particle type
897 for (Int_t i=0; i<AliPID::kSPECIES; i++) {
898 if (fR[i]>max) {k=i; max=fR[i];}
900 if (k==0) { // dE/dx "crossing points" in the TPC
902 if ((p>0.38)&&(p<0.48))
903 if (fR[0]<fR[3]*10.) return AliPID::ParticleMass(AliPID::kKaon);
904 if ((p>0.75)&&(p<0.85))
905 if (fR[0]<fR[4]*10.) return AliPID::ParticleMass(AliPID::kProton);
908 if (k==1) return AliPID::ParticleMass(AliPID::kMuon);
909 if (k==2||k==-1) return AliPID::ParticleMass(AliPID::kPion);
910 if (k==3) return AliPID::ParticleMass(AliPID::kKaon);
911 if (k==4) return AliPID::ParticleMass(AliPID::kProton);
912 AliWarning("Undefined mass !");
913 return AliPID::ParticleMass(AliPID::kPion);
916 //______________________________________________________________________________
917 Double_t AliESDtrack::E() const
919 // Returns the energy of the particle given its assumed mass.
920 // Assumes the pion mass if the particle can't be identified properly.
924 return TMath::Sqrt(p*p + m*m);
927 //______________________________________________________________________________
928 Double_t AliESDtrack::Y() const
930 // Returns the rapidity of a particle given its assumed mass.
931 // Assumes the pion mass if the particle can't be identified properly.
935 if (e != TMath::Abs(pz)) { // energy was not equal to pz
936 return 0.5*TMath::Log((e+pz)/(e-pz));
937 } else { // energy was equal to pz
942 //_______________________________________________________________________
943 Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
945 // This function updates track's running parameters
951 fLabel=t->GetLabel();
953 if (t->IsStartedTimeIntegral()) {
955 Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
956 SetIntegratedLength(t->GetIntegratedLength());
959 Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
963 case kITSin: case kITSout: case kITSrefit:
965 fITSncls=t->GetNumberOfClusters();
966 index=fFriendTrack->GetITSindices();
967 for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
968 index[i]=t->GetClusterIndex(i);
970 Int_t l=(index[i] & 0xf0000000) >> 28;
971 SETBIT(fITSClusterMap,l);
974 fITSchi2=t->GetChi2();
975 fITSsignal=t->GetPIDsignal();
976 fITSLabel = t->GetLabel();
977 // keep in fOp the parameters outside ITS for ITS stand-alone tracks
978 if (flags==kITSout) {
979 if (!fOp) fOp=new AliExternalTrackParam(*t);
981 fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
985 case kTPCin: case kTPCrefit:
986 fTPCLabel = t->GetLabel();
987 if (flags==kTPCin) fTPCInner=new AliExternalTrackParam(*t);
988 if (!fIp) fIp=new AliExternalTrackParam(*t);
990 fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
992 index=fFriendTrack->GetTPCindices();
993 if (flags & kTPCout){
994 if (!fOp) fOp=new AliExternalTrackParam(*t);
996 fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
998 fTPCncls=t->GetNumberOfClusters();
999 fTPCchi2=t->GetChi2();
1001 {//prevrow must be declared in separate namespace, otherwise compiler cries:
1002 //"jump to case label crosses initialization of `Int_t prevrow'"
1004 // for (Int_t i=0;i<fTPCncls;i++)
1005 for (Int_t i=0;i<AliESDfriendTrack::kMaxTPCcluster;i++)
1007 index[i]=t->GetClusterIndex(i);
1008 Int_t idx = index[i];
1010 if (idx<0) continue;
1012 // Piotr's Cluster Map for HBT
1013 // ### please change accordingly if cluster array is changing
1014 // to "New TPC Tracking" style (with gaps in array)
1015 Int_t sect = (idx&0xff000000)>>24;
1016 Int_t row = (idx&0x00ff0000)>>16;
1017 if (sect > 18) row +=63; //if it is outer sector, add number of inner sectors
1019 fTPCClusterMap.SetBitNumber(row,kTRUE);
1021 //Fill the gap between previous row and this row with 0 bits
1022 //In case ### pleas change it as well - just set bit 0 in case there
1023 //is no associated clusters for current "i"
1026 prevrow = row;//if previous bit was not assigned yet == this is the first one
1029 { //we don't know the order (inner to outer or reverse)
1030 //just to be save in case it is going to change
1043 for (Int_t j = n+1; j < m; j++)
1045 fTPCClusterMap.SetBitNumber(j,kFALSE);
1049 // End Of Piotr's Cluster Map for HBT
1052 fTPCsignal=t->GetPIDsignal();
1055 case kTRDout: case kTRDin: case kTRDrefit:
1056 index = fFriendTrack->GetTRDindices();
1057 fTRDLabel = t->GetLabel();
1058 fTRDchi2 = t->GetChi2();
1059 fTRDncls = t->GetNumberOfClusters();
1060 for (Int_t i=0;i<6;i++) index[i]=t->GetTrackletIndex(i);
1062 fTRDsignal=t->GetPIDsignal();
1065 if (!fOp) fOp=new AliExternalTrackParam(*t);
1067 fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
1068 fTRDncls0 = t->GetNumberOfClusters();
1077 AliError("Wrong flag !");
1084 //_______________________________________________________________________
1085 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
1086 //---------------------------------------------------------------------
1087 // This function returns external representation of the track parameters
1088 //---------------------------------------------------------------------
1090 for (Int_t i=0; i<5; i++) p[i]=GetParameter()[i];
1093 //_______________________________________________________________________
1094 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
1095 //---------------------------------------------------------------------
1096 // This function returns external representation of the cov. matrix
1097 //---------------------------------------------------------------------
1098 for (Int_t i=0; i<15; i++) cov[i]=AliExternalTrackParam::GetCovariance()[i];
1101 //_______________________________________________________________________
1102 Bool_t AliESDtrack::GetConstrainedExternalParameters
1103 (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1104 //---------------------------------------------------------------------
1105 // This function returns the constrained external track parameters
1106 //---------------------------------------------------------------------
1107 if (!fCp) return kFALSE;
1108 alpha=fCp->GetAlpha();
1110 for (Int_t i=0; i<5; i++) p[i]=fCp->GetParameter()[i];
1114 //_______________________________________________________________________
1116 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
1117 //---------------------------------------------------------------------
1118 // This function returns the constrained external cov. matrix
1119 //---------------------------------------------------------------------
1120 if (!fCp) return kFALSE;
1121 for (Int_t i=0; i<15; i++) c[i]=fCp->GetCovariance()[i];
1126 AliESDtrack::GetInnerExternalParameters
1127 (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1128 //---------------------------------------------------------------------
1129 // This function returns external representation of the track parameters
1130 // at the inner layer of TPC
1131 //---------------------------------------------------------------------
1132 if (!fIp) return kFALSE;
1133 alpha=fIp->GetAlpha();
1135 for (Int_t i=0; i<5; i++) p[i]=fIp->GetParameter()[i];
1140 AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const {
1141 //---------------------------------------------------------------------
1142 // This function returns external representation of the cov. matrix
1143 // at the inner layer of TPC
1144 //---------------------------------------------------------------------
1145 if (!fIp) return kFALSE;
1146 for (Int_t i=0; i<15; i++) cov[i]=fIp->GetCovariance()[i];
1151 AliESDtrack::SetOuterParam(const AliExternalTrackParam *p, ULong_t flags) {
1153 // This is a direct setter for the outer track parameters
1156 if (fOp) delete fOp;
1157 fOp=new AliExternalTrackParam(*p);
1161 AliESDtrack::GetOuterExternalParameters
1162 (Double_t &alpha, Double_t &x, Double_t p[5]) const {
1163 //---------------------------------------------------------------------
1164 // This function returns external representation of the track parameters
1165 // at the inner layer of TRD
1166 //---------------------------------------------------------------------
1167 if (!fOp) return kFALSE;
1168 alpha=fOp->GetAlpha();
1170 for (Int_t i=0; i<5; i++) p[i]=fOp->GetParameter()[i];
1175 AliESDtrack::GetOuterExternalCovariance(Double_t cov[15]) const {
1176 //---------------------------------------------------------------------
1177 // This function returns external representation of the cov. matrix
1178 // at the inner layer of TRD
1179 //---------------------------------------------------------------------
1180 if (!fOp) return kFALSE;
1181 for (Int_t i=0; i<15; i++) cov[i]=fOp->GetCovariance()[i];
1185 Int_t AliESDtrack::GetNcls(Int_t idet) const
1187 // Get number of clusters by subdetector index
1201 if (fTOFindex != -1)
1210 Int_t AliESDtrack::GetClusters(Int_t idet, Int_t *idx) const
1212 // Get cluster index array by subdetector index
1217 ncls = GetITSclusters(idx);
1220 ncls = GetTPCclusters(idx);
1223 ncls = GetTRDclusters(idx);
1226 if (fTOFindex != -1) {
1234 if (fHMPIDcluIdx != 0) {
1235 idx[0] = GetHMPIDcluIdx();
1247 //_______________________________________________________________________
1248 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
1249 // Returns the array with integrated times for each particle hypothesis
1250 for (Int_t i=0; i<AliPID::kSPECIES; i++) times[i]=fTrackTime[i];
1253 //_______________________________________________________________________
1254 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
1255 // Sets the array with integrated times for each particle hypotesis
1256 for (Int_t i=0; i<AliPID::kSPECIES; i++) fTrackTime[i]=times[i];
1259 //_______________________________________________________________________
1260 void AliESDtrack::SetITSpid(const Double_t *p) {
1261 // Sets values for the probability of each particle type (in ITS)
1262 SetPIDValues(fITSr,p,AliPID::kSPECIES);
1263 SetStatus(AliESDtrack::kITSpid);
1266 //_______________________________________________________________________
1267 void AliESDtrack::GetITSpid(Double_t *p) const {
1268 // Gets the probability of each particle type (in ITS)
1269 for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fITSr[i];
1272 //_______________________________________________________________________
1273 Char_t AliESDtrack::GetITSclusters(Int_t *idx) const {
1274 //---------------------------------------------------------------------
1275 // This function returns indices of the assgined ITS clusters
1276 //---------------------------------------------------------------------
1278 Int_t *index=fFriendTrack->GetITSindices();
1279 for (Int_t i=0; i<AliESDfriendTrack::kMaxITScluster; i++) {
1280 if ( (i>=fITSncls) && (i<6) ) idx[i]=-1;
1281 else idx[i]=index[i];
1287 //_______________________________________________________________________
1288 Bool_t AliESDtrack::GetITSModuleIndexInfo(Int_t ilayer,Int_t &idet,Int_t &status,
1289 Float_t &xloc,Float_t &zloc) const {
1290 //----------------------------------------------------------------------
1291 // This function encodes in the module number also the status of cluster association
1292 // "status" can have the following values:
1293 // 1 "found" (cluster is associated),
1294 // 2 "dead" (module is dead from OCDB),
1295 // 3 "skipped" (module or layer forced to be skipped),
1296 // 4 "outinz" (track out of z acceptance),
1297 // 5 "nocls" (no clusters in the road),
1298 // 6 "norefit" (cluster rejected during refit),
1299 // 7 "deadzspd" (holes in z in SPD)
1300 // Also given are the coordinates of the crossing point of track and module
1301 // (in the local module ref. system)
1302 // WARNING: THIS METHOD HAS TO BE SYNCHRONIZED WITH AliITStrackV2::GetModuleIndexInfo()!
1303 //----------------------------------------------------------------------
1305 if(fITSModule[ilayer]==-1) {
1306 AliError("fModule was not set !");
1309 xloc=-99.; zloc=-99.;
1313 Int_t module = fITSModule[ilayer];
1315 idet = Int_t(module/1000000);
1317 module -= idet*1000000;
1319 status = Int_t(module/100000);
1321 module -= status*100000;
1323 Int_t signs = Int_t(module/10000);
1325 module-=signs*10000;
1327 Int_t xInt = Int_t(module/100);
1330 Int_t zInt = module;
1332 if(signs==1) { xInt*=1; zInt*=1; }
1333 if(signs==2) { xInt*=1; zInt*=-1; }
1334 if(signs==3) { xInt*=-1; zInt*=1; }
1335 if(signs==4) { xInt*=-1; zInt*=-1; }
1337 xloc = 0.1*(Float_t)xInt;
1338 zloc = 0.1*(Float_t)zInt;
1340 if(status==4) idet = -1;
1345 //_______________________________________________________________________
1346 UShort_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
1347 //---------------------------------------------------------------------
1348 // This function returns indices of the assgined ITS clusters
1349 //---------------------------------------------------------------------
1351 Int_t *index=fFriendTrack->GetTPCindices();
1352 for (Int_t i=0; i<AliESDfriendTrack::kMaxTPCcluster; i++) idx[i]=index[i];
1357 Double_t AliESDtrack::GetTPCdensity(Int_t row0, Int_t row1) const{
1359 // GetDensity of the clusters on given region between row0 and row1
1360 // Dead zone effect takin into acoount
1365 Int_t *index=fFriendTrack->GetTPCindices();
1366 for (Int_t i=row0;i<=row1;i++){
1367 Int_t idx = index[i];
1368 if (idx!=-1) good++; // track outside of dead zone
1371 Float_t density=0.5;
1372 if (good>(row1-row0)*0.5) density = Float_t(found)/Float_t(good);
1376 //_______________________________________________________________________
1377 void AliESDtrack::SetTPCpid(const Double_t *p) {
1378 // Sets values for the probability of each particle type (in TPC)
1379 SetPIDValues(fTPCr,p,AliPID::kSPECIES);
1380 SetStatus(AliESDtrack::kTPCpid);
1383 //_______________________________________________________________________
1384 void AliESDtrack::GetTPCpid(Double_t *p) const {
1385 // Gets the probability of each particle type (in TPC)
1386 for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTPCr[i];
1389 //_______________________________________________________________________
1390 UChar_t AliESDtrack::GetTRDclusters(Int_t *idx) const {
1391 //---------------------------------------------------------------------
1392 // This function returns indices of the assgined TRD clusters
1393 //---------------------------------------------------------------------
1395 Int_t *index=fFriendTrack->GetTRDindices();
1396 for (Int_t i=0; i<AliESDfriendTrack::kMaxTRDcluster; i++) idx[i]=index[i];
1401 //_______________________________________________________________________
1402 UChar_t AliESDtrack::GetTRDtracklets(Int_t *idx) const {
1403 //---------------------------------------------------------------------
1404 // This function returns indices of the assigned TRD tracklets
1405 //---------------------------------------------------------------------
1407 Int_t *index=fFriendTrack->GetTRDindices();
1408 for (Int_t i=0; i<6/*AliESDfriendTrack::kMaxTRDcluster*/; i++) idx[i]=index[i];
1413 //_______________________________________________________________________
1414 void AliESDtrack::SetTRDpid(const Double_t *p) {
1415 // Sets values for the probability of each particle type (in TRD)
1416 SetPIDValues(fTRDr,p,AliPID::kSPECIES);
1417 SetStatus(AliESDtrack::kTRDpid);
1420 //_______________________________________________________________________
1421 void AliESDtrack::GetTRDpid(Double_t *p) const {
1422 // Gets the probability of each particle type (in TRD)
1423 for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTRDr[i];
1426 //_______________________________________________________________________
1427 void AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
1429 // Sets the probability of particle type iSpecies to p (in TRD)
1430 fTRDr[iSpecies] = p;
1433 Double_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
1435 // Returns the probability of particle type iSpecies (in TRD)
1436 return fTRDr[iSpecies];
1439 void AliESDtrack::SetNumberOfTRDslices(Int_t n) {
1440 //Sets the number of slices used for PID
1441 if (fTRDnSlices != 0) return;
1442 fTRDnSlices=kTRDnPlanes*n;
1443 fTRDslices=new Double32_t[fTRDnSlices];
1444 for (Int_t i=0; i<fTRDnSlices; i++) fTRDslices[i]=-1.;
1447 void AliESDtrack::SetTRDslice(Double_t q, Int_t plane, Int_t slice) {
1448 //Sets the charge q in the slice of the plane
1449 Int_t ns=GetNumberOfTRDslices();
1451 AliError("No TRD slices allocated for this track !");
1455 if ((plane<0) || (plane>=kTRDnPlanes)) {
1456 AliError("Wrong TRD plane !");
1459 if ((slice<0) || (slice>=ns)) {
1460 AliError("Wrong TRD slice !");
1463 Int_t n=plane*ns + slice;
1467 Double_t AliESDtrack::GetTRDslice(Int_t plane, Int_t slice) const {
1468 //Gets the charge from the slice of the plane
1469 Int_t ns=GetNumberOfTRDslices();
1471 //AliError("No TRD slices allocated for this track !");
1475 if ((plane<0) || (plane>=kTRDnPlanes)) {
1476 AliError("Wrong TRD plane !");
1479 if ((slice<-1) || (slice>=ns)) {
1480 //AliError("Wrong TRD slice !");
1486 for (Int_t i=0; i<ns; i++) q+=fTRDslices[plane*ns + i];
1490 return fTRDslices[plane*ns + slice];
1494 //_______________________________________________________________________
1495 void AliESDtrack::SetTOFpid(const Double_t *p) {
1496 // Sets the probability of each particle type (in TOF)
1497 SetPIDValues(fTOFr,p,AliPID::kSPECIES);
1498 SetStatus(AliESDtrack::kTOFpid);
1501 //_______________________________________________________________________
1502 void AliESDtrack::SetTOFLabel(const Int_t *p) {
1504 for (Int_t i=0; i<3; i++) fTOFLabel[i]=p[i];
1507 //_______________________________________________________________________
1508 void AliESDtrack::GetTOFpid(Double_t *p) const {
1509 // Gets probabilities of each particle type (in TOF)
1510 for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fTOFr[i];
1513 //_______________________________________________________________________
1514 void AliESDtrack::GetTOFLabel(Int_t *p) const {
1516 for (Int_t i=0; i<3; i++) p[i]=fTOFLabel[i];
1519 //_______________________________________________________________________
1520 void AliESDtrack::GetTOFInfo(Float_t *info) const {
1522 for (Int_t i=0; i<10; i++) info[i]=fTOFInfo[i];
1525 //_______________________________________________________________________
1526 void AliESDtrack::SetTOFInfo(Float_t*info) {
1528 for (Int_t i=0; i<10; i++) fTOFInfo[i]=info[i];
1533 //_______________________________________________________________________
1534 void AliESDtrack::SetHMPIDpid(const Double_t *p) {
1535 // Sets the probability of each particle type (in HMPID)
1536 SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
1537 SetStatus(AliESDtrack::kHMPIDpid);
1540 //_______________________________________________________________________
1541 void AliESDtrack::GetHMPIDpid(Double_t *p) const {
1542 // Gets probabilities of each particle type (in HMPID)
1543 for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
1548 //_______________________________________________________________________
1549 void AliESDtrack::SetESDpid(const Double_t *p) {
1550 // Sets the probability of each particle type for the ESD track
1551 SetPIDValues(fR,p,AliPID::kSPECIES);
1552 SetStatus(AliESDtrack::kESDpid);
1555 //_______________________________________________________________________
1556 void AliESDtrack::GetESDpid(Double_t *p) const {
1557 // Gets probability of each particle type for the ESD track
1558 for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fR[i];
1561 //_______________________________________________________________________
1562 Bool_t AliESDtrack::RelateToVertexTPC
1563 (const AliESDVertex *vtx, Double_t b, Double_t maxd) {
1565 // Try to relate the TPC-only track paramters to the vertex "vtx",
1566 // if the (rough) transverse impact parameter is not bigger then "maxd".
1567 // Magnetic field is "b" (kG).
1569 // a) The TPC-only paramters are extapolated to the DCA to the vertex.
1570 // b) The impact parameters and their covariance matrix are calculated.
1573 if (!fTPCInner) return kFALSE;
1574 if (!vtx) return kFALSE;
1576 Double_t dz[2],cov[3];
1577 if (!fTPCInner->PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
1588 //_______________________________________________________________________
1589 Bool_t AliESDtrack::RelateToVertex
1590 (const AliESDVertex *vtx, Double_t b, Double_t maxd) {
1592 // Try to relate this track to the vertex "vtx",
1593 // if the (rough) transverse impact parameter is not bigger then "maxd".
1594 // Magnetic field is "b" (kG).
1596 // a) The track gets extapolated to the DCA to the vertex.
1597 // b) The impact parameters and their covariance matrix are calculated.
1598 // c) An attempt to constrain this track to the vertex is done.
1600 // In the case of success, the returned value is kTRUE
1601 // (otherwise, it's kFALSE)
1604 if (!vtx) return kFALSE;
1606 Double_t dz[2],cov[3];
1607 if (!PropagateToDCA(vtx, b, maxd, dz, cov)) return kFALSE;
1615 Double_t covar[6]; vtx->GetCovMatrix(covar);
1616 Double_t p[2]={GetParameter()[0]-dz[0],GetParameter()[1]-dz[1]};
1617 Double_t c[3]={covar[2],0.,covar[5]};
1619 Double_t chi2=GetPredictedChi2(p,c);
1620 if (chi2>77.) return kFALSE;
1623 fCp=new AliExternalTrackParam(*this);
1625 if (!fCp->Update(p,c)) {delete fCp; fCp=0; return kFALSE;}
1631 //_______________________________________________________________________
1632 void AliESDtrack::Print(Option_t *) const {
1633 // Prints info on the track
1634 AliExternalTrackParam::Print();
1635 printf("ESD track info\n") ;
1636 Double_t p[AliPID::kSPECIESN] ;
1638 if( IsOn(kITSpid) ){
1639 printf("From ITS: ") ;
1641 for(index = 0 ; index < AliPID::kSPECIES; index++)
1642 printf("%f, ", p[index]) ;
1643 printf("\n signal = %f\n", GetITSsignal()) ;
1645 if( IsOn(kTPCpid) ){
1646 printf("From TPC: ") ;
1648 for(index = 0 ; index < AliPID::kSPECIES; index++)
1649 printf("%f, ", p[index]) ;
1650 printf("\n signal = %f\n", GetTPCsignal()) ;
1652 if( IsOn(kTRDpid) ){
1653 printf("From TRD: ") ;
1655 for(index = 0 ; index < AliPID::kSPECIES; index++)
1656 printf("%f, ", p[index]) ;
1657 printf("\n signal = %f\n", GetTRDsignal()) ;
1659 if( IsOn(kTOFpid) ){
1660 printf("From TOF: ") ;
1662 for(index = 0 ; index < AliPID::kSPECIES; index++)
1663 printf("%f, ", p[index]) ;
1664 printf("\n signal = %f\n", GetTOFsignal()) ;
1666 if( IsOn(kHMPIDpid) ){
1667 printf("From HMPID: ") ;
1669 for(index = 0 ; index < AliPID::kSPECIES; index++)
1670 printf("%f, ", p[index]) ;
1671 printf("\n signal = %f\n", GetHMPIDsignal()) ;
1677 // Draw functionality
1678 // Origin: Marian Ivanov, Marian.Ivanov@cern.ch
1680 void AliESDtrack::FillPolymarker(TPolyMarker3D *pol, Float_t magF, Float_t minR, Float_t maxR, Float_t stepR){
1682 // Fill points in the polymarker
1685 arrayRef.AddLast(new AliExternalTrackParam(*this));
1686 if (fIp) arrayRef.AddLast(new AliExternalTrackParam(*fIp));
1687 if (fOp) arrayRef.AddLast(new AliExternalTrackParam(*fOp));
1689 Double_t mpos[3]={0,0,0};
1690 Int_t entries=arrayRef.GetEntries();
1691 for (Int_t i=0;i<entries;i++){
1693 ((AliExternalTrackParam*)arrayRef.At(i))->GetXYZ(pos);
1694 mpos[0]+=pos[0]/entries;
1695 mpos[1]+=pos[1]/entries;
1696 mpos[2]+=pos[2]/entries;
1698 // Rotate to the mean position
1700 Float_t fi= TMath::ATan2(mpos[1],mpos[0]);
1701 for (Int_t i=0;i<entries;i++){
1702 Bool_t res = ((AliExternalTrackParam*)arrayRef.At(i))->Rotate(fi);
1703 if (!res) delete arrayRef.RemoveAt(i);
1706 for (Double_t r=minR; r<maxR; r+=stepR){
1708 Double_t mlpos[3]={0,0,0};
1709 for (Int_t i=0;i<entries;i++){
1710 Double_t point[3]={0,0,0};
1711 AliExternalTrackParam *param = ((AliExternalTrackParam*)arrayRef.At(i));
1712 if (!param) continue;
1713 if (param->GetXYZAt(r,magF,point)){
1714 Double_t weight = 1./(10.+(r-param->GetX())*(r-param->GetX()));
1716 mlpos[0]+=point[0]*weight;
1717 mlpos[1]+=point[1]*weight;
1718 mlpos[2]+=point[2]*weight;
1725 pol->SetPoint(counter,mlpos[0],mlpos[1], mlpos[2]);
1726 printf("xyz\t%f\t%f\t%f\n",mlpos[0], mlpos[1],mlpos[2]);