1 //_________________________________
2 ////////////////////////////////////////////////////////////
4 // class AliTrackPoints //
6 // used by Anti-Merging cut //
7 // contains set of poits the lay on track trajectory //
8 // according to reconstructed track parameters - //
9 // NOT CLUSTERS POSITIONS!!! //
10 // Anti-Merging cut is applied only on tracks coming from //
11 // different events (that are use to fill deniminators) //
13 ////////////////////////////////////////////////////////////
15 #include <TClonesArray.h>
19 #include "AliESDtrack.h"
20 #include "AliTrackPoints.h"
21 #include "AliTPCtrack.h"
22 #include "AliTrackReference.h"
24 ClassImp(AliTrackPoints)
26 Int_t AliTrackPoints::fgDebug = 0;
27 AliTrackPoints::AliTrackPoints():
35 /***************************************************************/
37 AliTrackPoints::AliTrackPoints(Int_t n, AliESDtrack* track, Float_t mf, Float_t dr, Float_t r0):
44 //mf - magnetic field in kG - needed to calculated curvature out of Pt
45 //r0 - starting radius
46 //dr - calculate points every dr cm, default every 30cm
49 Error("AliTrackPoints","ESD track is null");
58 if ( ((track->GetStatus() & AliESDtrack::kTPCrefit) == kFALSE)&&
59 ((track->GetStatus() & AliESDtrack::kTPCin) == kFALSE) )
61 //could happend: its stand alone tracking
63 Warning("AliTrackPoints","This ESD track does not contain TPC information");
76 track->GetInnerExternalParameters(x,par); //get properties of the track
79 Error("AliTrackPoints","This ESD track seem not to contain TPC information (curv is 0)");
85 Error("AliTrackPoints","Zero Magnetic field passed as parameter.");
89 Double_t alpha = track->GetInnerAlpha();
90 Double_t cc = 1000./0.299792458/mf;//conversion constant
93 MakePoints(dr,r0,x,par,c,alpha);
96 /***************************************************************/
98 AliTrackPoints::AliTrackPoints(Int_t n, AliTPCtrack* track, Float_t dr, Float_t r0):
106 //dr - calculate points every dr cm, default every 30cm
109 Error("AliTrackPoints","TPC track is null");
117 track->PropagateTo(r0);
119 //* This formation is now fixed in the following way: *
120 //* external param0: local Y-coordinate of a track (cm) *
121 //* external param1: local Z-coordinate of a track (cm) *
122 //* external param2: local sine of the track momentum azimuth angle *
123 //* external param3: tangent of the track momentum dip angle *
124 //* external param4: 1/pt (1/(GeV/c)) *
128 track->GetExternalParameters(x,par); //get properties of the track
130 Double_t alpha = track->GetAlpha();
131 Double_t c=track->GetC();
132 MakePoints(dr,r0,x,par,c,alpha);
135 void AliTrackPoints::MakePoints( Float_t dr, Float_t r0, Double_t x, Double_t* par, Double_t c, Double_t alpha)
137 //Calculates points starting at radius r0
138 //spacing every dr (in radial direction)
139 // according to track parameters
140 // x - position in sector local reference frame. x i parallel to R and sector is symmetric with respect to x
141 // par - track external parameters; array with 5 elements; look at AliTPCtrack.h or AliESDtrack.h for their meaning
142 // c - track curvature
143 // alpha - sector's rotation angle (phi) == angle needed for local to global transformation
146 Double_t z0 = par[1];
148 Double_t phi0local = TMath::ATan2(y,x);
149 Double_t phi0global = phi0local + alpha;
151 if (phi0local<0) phi0local+=2*TMath::Pi();
152 if (phi0local>=2.*TMath::Pi()) phi0local-=2*TMath::Pi();
154 if (phi0global<0) phi0global+=2*TMath::Pi();
155 if (phi0global>=2.*TMath::Pi()) phi0global-=2*TMath::Pi();
157 Double_t r = TMath::Hypot(x,y);
161 Info("AliTrackPoints","Radius0 %f, Real Radius %f",r0,r);
164 Info("AliTrackPoints","Phi Global at first padraw %f, Phi locat %f",phi0global,phi0local);
166 Double_t eta = x*c - par[2] ;//par[2] = fX*C - eta; eta==fP2 ; C==fP4
168 //this calculattions are assuming new (current) model
169 Double_t tmp = par[2];
171 tmp = c*y + TMath::Sqrt(tmp);
172 Double_t dca=(TMath::Hypot(eta,tmp) - 1. )/TMath::Abs(c);
174 //Here is old model Cold=Cnew/2.
175 Double_t dcasq = dca*dca;
177 Double_t cst1 = (1.+c2*dca)*dca;//first constant
178 Double_t cst2 = 1. + 2.*c2*dca;//second constant
180 Double_t factorPhi0 = TMath::ASin((c2*r + cst1/r)/cst2);
181 Double_t factorZ0 = TMath::ASin(c2*TMath::Sqrt((r*r-dcasq)/cst2))*par[3]/c2;
183 for(Int_t i = 0; i<fN; i++)
185 Double_t rc = r0 + i*dr;
186 Double_t ftmp = (c2*rc + cst1/rc)/cst2;
190 Warning("AliTrackPoints","ASin argument > 1 %f:",ftmp);
193 else if (ftmp < -1.0)
196 Warning("AliTrackPoints","ASin argument < -1 %f:",ftmp);
200 Double_t factorPhi = TMath::ASin( ftmp );//factor phi od rc
201 Double_t phi = phi0global + factorPhi - factorPhi0;
203 ftmp = (rc*rc-dcasq)/cst2;
207 Warning("AliTrackPoints","Sqrt argument < 0: %f",ftmp);
211 ftmp = c2*TMath::Sqrt(ftmp);
215 Warning("AliTrackPoints","ASin argument > 1: %f",ftmp);
218 else if (ftmp < -1.0)
221 Warning("AliTrackPoints","ASin argument < -1: %f",ftmp);
224 Double_t factorZ = TMath::ASin(ftmp)*par[3]/c2;
225 fZ[i] = z0 + factorZ - factorZ0;
226 fX[i] = rc*TMath::Cos(phi);
227 fY[i] = rc*TMath::Sin(phi);
229 if ( GetDebug() > 2 )
231 Info("AliTrackPoints","X %f Y %f Z %f R asked %f R obtained %f",
232 fX[i],fY[i],fZ[i],rc,TMath::Hypot(fX[i],fY[i]));
236 /***************************************************************/
238 AliTrackPoints::~AliTrackPoints()
245 /***************************************************************/
247 void AliTrackPoints::PositionAt(Int_t n, Float_t &x,Float_t &y,Float_t &z)
249 //returns position at point n
252 Error("PositionAt","Point %d out of range",n);
259 if ( GetDebug() > 1 )
261 Info("AliTrackPoints","n %d; X %f; Y %f; Z %f",n,x,y,z);
264 /***************************************************************/
266 Double_t AliTrackPoints::AvarageDistance(const AliTrackPoints& tr)
268 //returns the aritmethic avarage distance between two tracks
269 // Info("AvarageDistance","Entered");
270 if ( (fN <= 0) || (tr.fN <=0) )
272 if (GetDebug()) Warning("AvarageDistance","One of tracks is empty");
278 Warning("AvarageDistance","Number of points is not equal");
283 for (Int_t i = 0; i<fN; i++)
287 // Float_t r1sq = fX[i]*fX[i]+fY[i]*fY[i];
288 // Float_t r2sq = tr.fX[i]*tr.fX[i]+tr.fY[i]*tr.fY[i];
289 Float_t r1sq = TMath::Hypot(fX[i],fY[i]);
290 Float_t r2sq = TMath::Hypot(tr.fX[i],tr.fY[i]);
291 Info("AvarageDistance","radii: %f %f",r1sq,r2sq);
295 Double_t dx = fX[i]-tr.fX[i];
296 Double_t dy = fY[i]-tr.fY[i];
297 Double_t dz = fZ[i]-tr.fZ[i];
298 sum+=TMath::Sqrt(dx*dx + dy*dy + dz*dz);
302 Info("AvarageDistance","Diff: x ,y z: %f , %f, %f",dx,dy,dz);
303 Info("AvarageDistance","xxyyzz %f %f %f %f %f %f",
304 fX[i],tr.fX[i],fY[i],tr.fY[i],fZ[i],tr.fZ[i]);
308 Double_t retval = sum/((Double_t)fN);
311 Info("AvarageDistance","Avarage distance is %f.",retval);
315 /***************************************************************/
316 /***************************************************************/
317 /***************************************************************/
318 /***************************************************************/
319 /***************************************************************/
320 /***************************************************************/
324 #include "AliRunLoader.h"
325 #include "AliTPCtrack.h"
335 void AliTrackPoints::TestESD(Int_t entr,const char* fname )
337 //Tests creation of track points for ESD tracks
340 AliRunLoader* rl = AliRunLoader::Open();
343 Float_t mf = rl->GetAliRun()->Field()->SolenoidField();
346 TFile* fFile = TFile::Open(fname);
350 printf("testesd: There is no suche a ESD file\n");
353 AliESD* esd = dynamic_cast<AliESD*>(fFile->Get("0"));
354 AliESDtrack *t = esd->GetTrack(entr);
357 ::Error("testesd","Can not get track %d",entr);
363 AliTrackPoints* tp = new AliTrackPoints(N,t,mf,1.);
374 TH2D* hxy = new TH2D("hxy","hxy",1000,xmin,xmax,1000,ymin,ymax);
375 TH2D* hxyt = new TH2D("hxyt","hxyt",1000,xmin,xmax,1000,ymin,ymax);
376 TH2D* hxyTR = new TH2D("hxyTR","hxyTR",1000,xmin,xmax,1000,ymin,ymax);
378 TH2D* hxz = new TH2D("hxz","hxz",1000,xmin,xmax,1000,zmin,zmax);
379 TH2D* hxzt = new TH2D("hxzt","hxzt",1000,xmin,xmax,1000,zmin,zmax);
380 TH2D* hxzTR = new TH2D("hxzTR","hxzTR",1000,xmin,xmax,1000,zmin,zmax);
382 hxyt->SetDirectory(0x0);
383 hxy->SetDirectory(0x0);
384 hxyTR->SetDirectory(0x0);
386 hxzt->SetDirectory(0x0);
387 hxz->SetDirectory(0x0);
388 hxzTR->SetDirectory(0x0);
392 for (Int_t i = 0;i<N;i++)
395 tp->PositionAt(i,x,y,z);
398 printf("Rdemanded %f\n",r);
399 printf("tpx %f tpy %f tpz %f Rt =%f\n", x,y,z,TMath::Hypot(x,y));
404 TTree* treeTR = rl->TreeTR();
405 TBranch* b = treeTR->GetBranch("TPC");
407 TClonesArray* trackrefs = new TClonesArray("AliTrackReference", 100);
408 AliTrackReference* tref;
409 b->SetAddress(&trackrefs);
411 Int_t tlab = TMath::Abs(t->GetLabel());
413 Int_t netr = (Int_t)treeTR->GetEntries();
414 printf("Found %d entries in TR tree\n",netr);
416 for (Int_t e = 0; e < netr; e++)
419 tref = (AliTrackReference*)trackrefs->At(0);
420 if (tref == 0x0) continue;
421 if (tref->GetTrack() != tlab) continue;
423 printf("Found %d entries in TR array\n",trackrefs->GetEntries());
425 for (Int_t i = 0; i < trackrefs->GetEntries(); i++)
427 tref = (AliTrackReference*)trackrefs->At(i);
428 if (tref->GetTrack() != tlab) continue;
432 printf("Track Ref: x %f y %f z %f\n",tref->X(),tref->Y(),tref->Z());
436 for (Int_t j = 1; j < 10; j++)
438 hxyTR->Fill(x, y+j*0.1);
439 hxyTR->Fill(x, y-j*0.1);
440 hxyTR->Fill(x+j*0.1,y);
441 hxyTR->Fill(x-j*0.1,y);
443 hxzTR->Fill(x,z-j*0.1);
444 hxzTR->Fill(x,z+j*0.1);
445 hxzTR->Fill(x-j*0.1,z);
446 hxzTR->Fill(x+j*0.1,z);
452 // hxzt->Draw("same");
458 /***************************************************************/
459 /***************************************************************/
460 /***************************************************************/
462 void AliTrackPoints::TestTPC(Int_t entr)
464 //Tests creation of track points for TPC tracks
467 AliRunLoader* rl = AliRunLoader::Open();
468 AliLoader* l = rl->GetLoader("TPCLoader");
470 AliKalmanTrack::SetConvConst(100/0.299792458/0.2/rl->GetAliRun()->Field()->Factor());
472 AliTPCtrack* t = new AliTPCtrack();
473 TBranch* b=l->TreeT()->GetBranch("tracks");
475 l->TreeT()->GetEntry(entr);
477 AliTrackPoints* tp = new AliTrackPoints(N,t,1.);
488 TH2D* hxy = new TH2D("hxy","hxy",1000,xmin,xmax,1000,ymin,ymax);
489 TH2D* hxyt = new TH2D("hxyt","hxyt",1000,xmin,xmax,1000,ymin,ymax);
490 TH2D* hxyTR = new TH2D("hxyTR","hxyTR",1000,xmin,xmax,1000,ymin,ymax);
492 TH2D* hxz = new TH2D("hxz","hxz",1000,xmin,xmax,1000,zmin,zmax);
493 TH2D* hxzt = new TH2D("hxzt","hxzt",1000,xmin,xmax,1000,zmin,zmax);
494 TH2D* hxzTR = new TH2D("hxzTR","hxzTR",1000,xmin,xmax,1000,zmin,zmax);
496 hxyt->SetDirectory(0x0);
497 hxy->SetDirectory(0x0);
498 hxyTR->SetDirectory(0x0);
500 hxzt->SetDirectory(0x0);
501 hxz->SetDirectory(0x0);
502 hxzTR->SetDirectory(0x0);
506 for (Int_t i = 0;i<N;i++)
509 tp->PositionAt(i,x,y,z);
512 printf("Rdemanded %f\n",r);
513 printf("tpx %f tpy %f tpz %f Rt =%f\n", x,y,z,TMath::Hypot(x,y));
515 //BUT they are local!!!!
517 // Double_t phi = t->Phi();
518 Double_t rl = TMath::Hypot(t->GetX(),t->GetY());//real radius
520 Double_t alpha = t->GetAlpha();
521 Double_t salpha = TMath::Sin(alpha);
522 Double_t calpha = TMath::Cos(alpha);
523 x = t->GetX()*calpha - t->GetY()*salpha;
524 y = t->GetX()*salpha + t->GetY()*calpha;
527 printf("tx %f ty %f tz %f Rt = %f R from XY %f\n",x,y,z,TMath::Hypot(x,y),rl);
529 printf("tpz - tz %f\n",z-t->GetZ());
537 TTree* treeTR = rl->TreeTR();
538 b = treeTR->GetBranch("TPC");
540 TClonesArray* trackrefs = new TClonesArray("AliTrackReference", 100);
541 AliTrackReference* tref;
542 b->SetAddress(&trackrefs);
544 Int_t tlab = TMath::Abs(t->GetLabel());
546 Int_t netr = (Int_t)treeTR->GetEntries();
547 printf("Found %d entries in TR tree\n",netr);
549 for (Int_t e = 0; e < netr; e++)
552 tref = (AliTrackReference*)trackrefs->At(0);
553 if (tref == 0x0) continue;
554 if (tref->GetTrack() != tlab) continue;
556 printf("Found %d entries in TR array\n",trackrefs->GetEntries());
558 for (Int_t i = 0; i < trackrefs->GetEntries(); i++)
560 tref = (AliTrackReference*)trackrefs->At(i);
561 if (tref->GetTrack() != tlab) continue;
565 printf("Track Ref: x %f y %f z %f\n",tref->X(),tref->Y(),tref->Z());
569 for (Int_t j = 1; j < 10; j++)
571 hxyTR->Fill(x, y+j*0.1);
572 hxyTR->Fill(x, y-j*0.1);
573 hxyTR->Fill(x+j*0.1,y);
574 hxyTR->Fill(x-j*0.1,y);
576 hxzTR->Fill(x,z-j*0.1);
577 hxzTR->Fill(x,z+j*0.1);
578 hxzTR->Fill(x-j*0.1,z);
579 hxzTR->Fill(x+j*0.1,z);
585 // hxzt->Draw("same");