Tracking in non-uniform nmagnetic field (Yu.Belikov)
[u/mrichter/AliRoot.git] / ANALYSIS / AliTrackPoints.cxx
CommitLineData
c7ffd78f 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
78d7c6d3 18#include "AliTrackPoints.h"
c7ffd78f 19#include "AliLog.h"
20
a5556ea5 21//_________________________________
22////////////////////////////////////////////////////////////
23// //
24// class AliTrackPoints //
25// //
26// used by Anti-Merging cut //
27// contains set of poits the lay on track trajectory //
28// according to reconstructed track parameters - //
29// NOT CLUSTERS POSITIONS!!! //
30// Anti-Merging cut is applied only on tracks coming from //
31// different events (that are use to fill deniminators) //
32// //
33////////////////////////////////////////////////////////////
34
35#include <TClonesArray.h>
36#include <TFile.h>
37#include <TMath.h>
38
39#include "AliESDtrack.h"
a5556ea5 40#include "AliTPCtrack.h"
41#include "AliTrackReference.h"
78d7c6d3 42#include "AliITStrackV2.h"
a5556ea5 43
44ClassImp(AliTrackPoints)
45
46Int_t AliTrackPoints::fgDebug = 0;
78d7c6d3 47
a5556ea5 48AliTrackPoints::AliTrackPoints():
49 fN(0),
50 fX(0x0),
51 fY(0x0),
52 fZ(0x0)
53{
54 //constructor
55}
56/***************************************************************/
57
0f0355ab 58AliTrackPoints::AliTrackPoints(AliTrackPoints::ETypes type, AliESDtrack* track, Float_t mf):
78d7c6d3 59 fN(0),
60 fX(0x0),
61 fY(0x0),
62 fZ(0x0)
63{
0f0355ab 64 //constructor
65 //tupe - what kind of track points should be calculated
66 //mf - magnetic field in [kG] = [T]*10.0
78d7c6d3 67 switch (type)
68 {
69 case kITS:
78d7c6d3 70 fN = 6;
71 fX = new Float_t[fN];
72 fY = new Float_t[fN];
73 fZ = new Float_t[fN];
74 MakeITSPoints(track);
75 break;
76
0f0355ab 77 case kITSInnerFromVertexOuterFromTPC:
78 fN = 6;
79 fX = new Float_t[fN];
80 fY = new Float_t[fN];
81 fZ = new Float_t[fN];
82 MakeITSPointsInnerFromVertexOuterFromTPC(track,mf);
83 break;
84
78d7c6d3 85 default:
86 Info("AliTrackPoints","Not recognized type");
87 }
88
89}
90/***************************************************************/
91
a5556ea5 92AliTrackPoints::AliTrackPoints(Int_t n, AliESDtrack* track, Float_t mf, Float_t dr, Float_t r0):
93 fN(n),
94 fX(new Float_t[fN]),
95 fY(new Float_t[fN]),
96 fZ(new Float_t[fN])
97{
98 //constructor
0f0355ab 99 //mf - magnetic field in kG - needed to calculate curvature out of Pt
a5556ea5 100 //r0 - starting radius
101 //dr - calculate points every dr cm, default every 30cm
102 if (track == 0x0)
103 {
104 Error("AliTrackPoints","ESD track is null");
105 fN = 0;
106 delete [] fX;
107 delete [] fY;
108 delete [] fZ;
109 fX = fY = fZ = 0x0;
110 return;
111 }
112
113 if ( ((track->GetStatus() & AliESDtrack::kTPCrefit) == kFALSE)&&
114 ((track->GetStatus() & AliESDtrack::kTPCin) == kFALSE) )
115 {
116 //could happend: its stand alone tracking
c7ffd78f 117 AliDebug(3,"This ESD track does not contain TPC information");
a5556ea5 118
119 fN = 0;
120 delete [] fX;
121 delete [] fY;
122 delete [] fZ;
123 fX = fY = fZ = 0x0;
124
125 return;
126 }
127
128 Double_t x;
129 Double_t par[5];
130 track->GetInnerExternalParameters(x,par); //get properties of the track
131 if (par[4] == 0)
132 {
133 Error("AliTrackPoints","This ESD track seem not to contain TPC information (curv is 0)");
134 return;
135 }
136
137 if (mf == 0.0)
138 {
139 Error("AliTrackPoints","Zero Magnetic field passed as parameter.");
140 return;
141 }
142
143 Double_t alpha = track->GetInnerAlpha();
144 Double_t cc = 1000./0.299792458/mf;//conversion constant
145 Double_t c=par[4]/cc;
146
147 MakePoints(dr,r0,x,par,c,alpha);
148
149}
150/***************************************************************/
151
152AliTrackPoints::AliTrackPoints(Int_t n, AliTPCtrack* track, Float_t dr, Float_t r0):
153 fN(n),
154 fX(new Float_t[fN]),
155 fY(new Float_t[fN]),
156 fZ(new Float_t[fN])
157{
158 //constructor
159 //r0 starting radius
160 //dr - calculate points every dr cm, default every 30cm
161 if (track == 0x0)
162 {
163 Error("AliTrackPoints","TPC track is null");
164 fN = 0;
165 delete [] fX;
166 delete [] fY;
167 delete [] fZ;
168 fX = fY = fZ = 0x0;
169 return;
170 }
171 track->PropagateTo(r0);
172
173 //* This formation is now fixed in the following way: *
174 //* external param0: local Y-coordinate of a track (cm) *
175 //* external param1: local Z-coordinate of a track (cm) *
176 //* external param2: local sine of the track momentum azimuth angle *
177 //* external param3: tangent of the track momentum dip angle *
178 //* external param4: 1/pt (1/(GeV/c)) *
179
180 Double_t x = 0;
181 Double_t par[5];
182 track->GetExternalParameters(x,par); //get properties of the track
183
184 Double_t alpha = track->GetAlpha();
185 Double_t c=track->GetC();
186 MakePoints(dr,r0,x,par,c,alpha);
187}
78d7c6d3 188/***************************************************************/
189
190AliTrackPoints::~AliTrackPoints()
191{
192 //destructor
193 delete [] fX;
194 delete [] fY;
195 delete [] fZ;
196}
197/***************************************************************/
a5556ea5 198
199void AliTrackPoints::MakePoints( Float_t dr, Float_t r0, Double_t x, Double_t* par, Double_t c, Double_t alpha)
200{
201 //Calculates points starting at radius r0
202 //spacing every dr (in radial direction)
203 // according to track parameters
204 // x - position in sector local reference frame. x i parallel to R and sector is symmetric with respect to x
205 // par - track external parameters; array with 5 elements; look at AliTPCtrack.h or AliESDtrack.h for their meaning
206 // c - track curvature
207 // alpha - sector's rotation angle (phi) == angle needed for local to global transformation
208
209 Double_t y = par[0];
210 Double_t z0 = par[1];
211
212 Double_t phi0local = TMath::ATan2(y,x);
213 Double_t phi0global = phi0local + alpha;
214
215 if (phi0local<0) phi0local+=2*TMath::Pi();
216 if (phi0local>=2.*TMath::Pi()) phi0local-=2*TMath::Pi();
217
218 if (phi0global<0) phi0global+=2*TMath::Pi();
219 if (phi0global>=2.*TMath::Pi()) phi0global-=2*TMath::Pi();
220
221 Double_t r = TMath::Hypot(x,y);
222
223
c7ffd78f 224 AliDebug(9,Form("Radius0 %f, Real Radius %f",r0,r));
a5556ea5 225
c7ffd78f 226 AliDebug(5,Form("Phi Global at first padraw %f, Phi locat %f",phi0global,phi0local));
a5556ea5 227
228 Double_t eta = x*c - par[2] ;//par[2] = fX*C - eta; eta==fP2 ; C==fP4
229
230 //this calculattions are assuming new (current) model
231 Double_t tmp = par[2];
232 tmp = 1. - tmp*tmp;
233 tmp = c*y + TMath::Sqrt(tmp);
234 Double_t dca=(TMath::Hypot(eta,tmp) - 1. )/TMath::Abs(c);
235
236 //Here is old model Cold=Cnew/2.
237 Double_t dcasq = dca*dca;
238 Double_t c2 = c/2.;
239 Double_t cst1 = (1.+c2*dca)*dca;//first constant
240 Double_t cst2 = 1. + 2.*c2*dca;//second constant
241
242 Double_t factorPhi0 = TMath::ASin((c2*r + cst1/r)/cst2);
243 Double_t factorZ0 = TMath::ASin(c2*TMath::Sqrt((r*r-dcasq)/cst2))*par[3]/c2;
244
245 for(Int_t i = 0; i<fN; i++)
246 {
247 Double_t rc = r0 + i*dr;
248 Double_t ftmp = (c2*rc + cst1/rc)/cst2;
249 if (ftmp > 1.0)
250 {
c7ffd78f 251 AliDebug(1,Form("ASin argument > 1 %f:",ftmp));
a5556ea5 252 ftmp=1.0;
253 }
254 else if (ftmp < -1.0)
255 {
c7ffd78f 256 AliDebug(1,Form("ASin argument < -1 %f:",ftmp));
a5556ea5 257 ftmp=-1.0;
258 }
259
260 Double_t factorPhi = TMath::ASin( ftmp );//factor phi od rc
261 Double_t phi = phi0global + factorPhi - factorPhi0;
262
263 ftmp = (rc*rc-dcasq)/cst2;
264 if (ftmp < 0.0)
265 {
c7ffd78f 266 AliDebug(1,Form("Sqrt argument < 0: %f",ftmp));
a5556ea5 267 ftmp=0.0;
268 }
269
270 ftmp = c2*TMath::Sqrt(ftmp);
271 if (ftmp > 1.0)
272 {
c7ffd78f 273 AliDebug(1,Form("ASin argument > 1: %f",ftmp));
a5556ea5 274 ftmp=1.0;
275 }
276 else if (ftmp < -1.0)
277 {
c7ffd78f 278 AliDebug(2,Form("ASin argument < -1: %f",ftmp));
a5556ea5 279 ftmp=-1.0;
280 }
281 Double_t factorZ = TMath::ASin(ftmp)*par[3]/c2;
282 fZ[i] = z0 + factorZ - factorZ0;
283 fX[i] = rc*TMath::Cos(phi);
284 fY[i] = rc*TMath::Sin(phi);
285
c7ffd78f 286 AliDebug(3,Form("AliTrackPoints","X %f Y %f Z %f R asked %f R obtained %f",
287 fX[i],fY[i],fZ[i],rc,TMath::Hypot(fX[i],fY[i])));
a5556ea5 288 }
289}
290/***************************************************************/
291
78d7c6d3 292void AliTrackPoints::MakeITSPoints(AliESDtrack* track)
a5556ea5 293{
78d7c6d3 294//Calculates points in ITS
295// z=R*Pz/Pt
296 AliITStrackV2 itstrack(*track,kTRUE);
297 Double_t x,y,z;
298 static const Double_t r[6] = {4.0, 7.0, 14.9, 23.8, 39.1, 43.6};
299 for (Int_t i = 0; i < 6; i++)
300 {
301 itstrack.GetGlobalXYZat(r[i],x,y,z);
302 fX[i] = x;
303 fY[i] = y;
304 fZ[i] = z;
305// Info("MakeITSPoints","X %f Y %f Z %f R asked %f R obtained %f",
306// fX[i],fY[i],fZ[i],r[i],TMath::Hypot(fX[i],fY[i]));
307 }
308
a5556ea5 309}
a5556ea5 310
78d7c6d3 311/***************************************************************/
0f0355ab 312void AliTrackPoints::MakeITSPointsInnerFromVertexOuterFromTPC(AliESDtrack* track, Float_t mf)
313{
314//makes trackpoints for ITS
315//for 3 inner layers calculates out of the vector at vertex
316//for 3 outer ---------------//------------------ at inner TPC
317
318 static const Double_t r[6] = {4.0, 7.0, 14.9, 23.8, 39.1, 43.6};
319 AliITStrackV2 itstrack(*track,kTRUE);
320 Double_t x,y,z;
321 for (Int_t i = 0; i < 3; i++)
322 {
323 itstrack.GetGlobalXYZat(r[i],x,y,z);
324 fX[i] = x;
325 fY[i] = y;
326 fZ[i] = z;
c7ffd78f 327 AliDebug(3,Form("X %f Y %f Z %f R asked %f R obtained %f",
328 fX[i],fY[i],fZ[i],r[i],TMath::Hypot(fX[i],fY[i])));
0f0355ab 329 }
330
331 for (Int_t i = 3; i < 6; i++)
332 {
333 Float_t ax,ay,az;
334 AliTrackPoints tmptp(1,track,mf,0,r[i]);
335 tmptp.PositionAt(0,ax,ay,az);
336 fX[i] = ax;
337 fY[i] = ay;
338 fZ[i] = az;
c7ffd78f 339 AliDebug(3,Form("X %f Y %f Z %f R asked %f R obtained %f",
340 fX[i],fY[i],fZ[i],r[i],TMath::Hypot(fX[i],fY[i])));
0f0355ab 341 }
342
343}
344
345
346/***************************************************************/
347
a5556ea5 348void AliTrackPoints::PositionAt(Int_t n, Float_t &x,Float_t &y,Float_t &z)
349{
350 //returns position at point n
5b49024d 351 if ((n<0) || (n>=fN))
a5556ea5 352 {
353 Error("PositionAt","Point %d out of range",n);
354 return;
355 }
356
357 x = fX[n];
358 y = fY[n];
359 z = fZ[n];
c7ffd78f 360 AliDebug(2,Form("n %d; X %f; Y %f; Z %f",n,x,y,z));
361
a5556ea5 362}
363/***************************************************************/
364
78d7c6d3 365void AliTrackPoints::Move(Float_t x, Float_t y, Float_t z)
366{
367//Moves all points about vector
368 for (Int_t i = 0; i<fN; i++)
369 {
370 fX[i]+=x;
371 fY[i]+=y;
372 fZ[i]+=z;
373 }
374}
375/***************************************************************/
376
a5556ea5 377Double_t AliTrackPoints::AvarageDistance(const AliTrackPoints& tr)
378{
379 //returns the aritmethic avarage distance between two tracks
380// Info("AvarageDistance","Entered");
381 if ( (fN <= 0) || (tr.fN <=0) )
382 {
c7ffd78f 383 AliDebug(1,"One of tracks is empty");
a5556ea5 384 return -1;
385 }
386
387 if (fN != tr.fN)
388 {
389 Warning("AvarageDistance","Number of points is not equal");
390 return -1;
391 }
392
393 Double_t sum = 0;
394 for (Int_t i = 0; i<fN; i++)
395 {
c7ffd78f 396 AliDebug(10,Form("radii: %f %f",TMath::Hypot(fX[i],fY[i]),TMath::Hypot(tr.fX[i],tr.fY[i])));
a5556ea5 397// Float_t r1sq = fX[i]*fX[i]+fY[i]*fY[i];
398// Float_t r2sq = tr.fX[i]*tr.fX[i]+tr.fY[i]*tr.fY[i];
a5556ea5 399
400
401 Double_t dx = fX[i]-tr.fX[i];
402 Double_t dy = fY[i]-tr.fY[i];
403 Double_t dz = fZ[i]-tr.fZ[i];
404 sum+=TMath::Sqrt(dx*dx + dy*dy + dz*dz);
405
c7ffd78f 406 AliDebug(2,Form("Diff: x ,y z: %f , %f, %f",dx,dy,dz));
407 AliDebug(2,Form("xxyyzz %f %f %f %f %f %f",
408 fX[i],tr.fX[i],fY[i],tr.fY[i],fZ[i],tr.fZ[i]));
a5556ea5 409 }
410
411 Double_t retval = sum/((Double_t)fN);
c7ffd78f 412 AliDebug(1,Form("Avarage distance is %f.",retval));
413
a5556ea5 414 return retval;
415}
416/***************************************************************/
5b49024d 417
418void AliTrackPoints::Print(Option_t* /*option*/) const
419{
420
421 Info("Print","There is %d points",fN);
422 for(Int_t i = 0; i < fN; i++)
423 {
424 Info("Print","%d: %f %f %f",i,fX[i],fY[i],fZ[i]);
425 }
426
427}
428
429/***************************************************************/
a5556ea5 430/***************************************************************/
431/***************************************************************/
432/***************************************************************/
433/***************************************************************/
434/***************************************************************/
435
436#include "AliRun.h"
437#include "AliESD.h"
438#include "AliRunLoader.h"
439#include "AliTPCtrack.h"
440#include "TTree.h"
441#include "TBranch.h"
442#include "TH2D.h"
443#include "TCanvas.h"
444#include "AliMagF.h"
445
446
447
448
78d7c6d3 449void AliTrackPoints::Testesd(Int_t entr,const char* fname )
a5556ea5 450{
a5556ea5 451 delete gAlice;
452 gAlice = 0x0;
453 AliRunLoader* rl = AliRunLoader::Open();
454 rl->LoadgAlice();
455
456 Float_t mf = rl->GetAliRun()->Field()->SolenoidField();
457
458
459 TFile* fFile = TFile::Open(fname);
460
461 if (fFile == 0x0)
462 {
463 printf("testesd: There is no suche a ESD file\n");
464 return;
465 }
466 AliESD* esd = dynamic_cast<AliESD*>(fFile->Get("0"));
467 AliESDtrack *t = esd->GetTrack(entr);
468 if (t == 0x0)
469 {
470 ::Error("testesd","Can not get track %d",entr);
471 return;
472 }
473
474
475 Int_t N = 170;
476 AliTrackPoints* tp = new AliTrackPoints(N,t,mf,1.);
477
478 Float_t xmin = -250;
479 Float_t xmax = 250;
480
481 Float_t ymin = -250;
482 Float_t ymax = 250;
483
484 Float_t zmin = -250;
485 Float_t zmax = 250;
486
487 TH2D* hxy = new TH2D("hxy","hxy",1000,xmin,xmax,1000,ymin,ymax);
488 TH2D* hxyt = new TH2D("hxyt","hxyt",1000,xmin,xmax,1000,ymin,ymax);
489 TH2D* hxyTR = new TH2D("hxyTR","hxyTR",1000,xmin,xmax,1000,ymin,ymax);
490
491 TH2D* hxz = new TH2D("hxz","hxz",1000,xmin,xmax,1000,zmin,zmax);
492 TH2D* hxzt = new TH2D("hxzt","hxzt",1000,xmin,xmax,1000,zmin,zmax);
493 TH2D* hxzTR = new TH2D("hxzTR","hxzTR",1000,xmin,xmax,1000,zmin,zmax);
494
495 hxyt->SetDirectory(0x0);
496 hxy->SetDirectory(0x0);
497 hxyTR->SetDirectory(0x0);
498
499 hxzt->SetDirectory(0x0);
500 hxz->SetDirectory(0x0);
501 hxzTR->SetDirectory(0x0);
502
503 Float_t x,y,z;
504
505 for (Int_t i = 0;i<N;i++)
506 {
507 Double_t r = 84.1+i;
508 tp->PositionAt(i,x,y,z);
509 hxy->Fill(x,y);
510 hxz->Fill(x,z);
511 printf("Rdemanded %f\n",r);
512 printf("tpx %f tpy %f tpz %f Rt =%f\n", x,y,z,TMath::Hypot(x,y));
513
514 }
515
516 rl->LoadTrackRefs();
517 TTree* treeTR = rl->TreeTR();
518 TBranch* b = treeTR->GetBranch("TPC");
519
520 TClonesArray* trackrefs = new TClonesArray("AliTrackReference", 100);
521 AliTrackReference* tref;
522 b->SetAddress(&trackrefs);
523
524 Int_t tlab = TMath::Abs(t->GetLabel());
525
526 Int_t netr = (Int_t)treeTR->GetEntries();
527 printf("Found %d entries in TR tree\n",netr);
528
529 for (Int_t e = 0; e < netr; e++)
530 {
531 treeTR->GetEntry(e);
532 tref = (AliTrackReference*)trackrefs->At(0);
533 if (tref == 0x0) continue;
534 if (tref->GetTrack() != tlab) continue;
535
536 printf("Found %d entries in TR array\n",trackrefs->GetEntries());
537
538 for (Int_t i = 0; i < trackrefs->GetEntries(); i++)
539 {
540 tref = (AliTrackReference*)trackrefs->At(i);
541 if (tref->GetTrack() != tlab) continue;
542 x = tref->X();
543 y = tref->Y();
544 z = tref->Z();
545 printf("Track Ref: x %f y %f z %f\n",tref->X(),tref->Y(),tref->Z());
546
547 hxzTR->Fill(x,z);
548 hxyTR->Fill(x,y);
549 for (Int_t j = 1; j < 10; j++)
550 {
551 hxyTR->Fill(x, y+j*0.1);
552 hxyTR->Fill(x, y-j*0.1);
553 hxyTR->Fill(x+j*0.1,y);
554 hxyTR->Fill(x-j*0.1,y);
555
556 hxzTR->Fill(x,z-j*0.1);
557 hxzTR->Fill(x,z+j*0.1);
558 hxzTR->Fill(x-j*0.1,z);
559 hxzTR->Fill(x+j*0.1,z);
560 }
561 }
562 break;
563 }
564 hxy->Draw("");
565// hxzt->Draw("same");
566 hxyTR->Draw("same");
567
568 delete rl;
569}
570
571/***************************************************************/
572/***************************************************************/
573/***************************************************************/
574
78d7c6d3 575void AliTrackPoints::Testtpc(Int_t entr)
a5556ea5 576{
a5556ea5 577 delete gAlice;
578 gAlice = 0x0;
579 AliRunLoader* rl = AliRunLoader::Open();
580 AliLoader* l = rl->GetLoader("TPCLoader");
581 rl->LoadgAlice();
c84a5e9e 582 AliKalmanTrack::SetFieldMap(rl->GetAliRun()->Field());
a5556ea5 583 l->LoadTracks();
584 AliTPCtrack* t = new AliTPCtrack();
585 TBranch* b=l->TreeT()->GetBranch("tracks");
586 b->SetAddress(&t);
587 l->TreeT()->GetEntry(entr);
588 Int_t N = 160;
589 AliTrackPoints* tp = new AliTrackPoints(N,t,1.);
590
591 Float_t xmin = -250;
592 Float_t xmax = 250;
593
594 Float_t ymin = -250;
595 Float_t ymax = 250;
596
597 Float_t zmin = -250;
598 Float_t zmax = 250;
599
600 TH2D* hxy = new TH2D("hxy","hxy",1000,xmin,xmax,1000,ymin,ymax);
601 TH2D* hxyt = new TH2D("hxyt","hxyt",1000,xmin,xmax,1000,ymin,ymax);
602 TH2D* hxyTR = new TH2D("hxyTR","hxyTR",1000,xmin,xmax,1000,ymin,ymax);
603
604 TH2D* hxz = new TH2D("hxz","hxz",1000,xmin,xmax,1000,zmin,zmax);
605 TH2D* hxzt = new TH2D("hxzt","hxzt",1000,xmin,xmax,1000,zmin,zmax);
606 TH2D* hxzTR = new TH2D("hxzTR","hxzTR",1000,xmin,xmax,1000,zmin,zmax);
607
608 hxyt->SetDirectory(0x0);
609 hxy->SetDirectory(0x0);
610 hxyTR->SetDirectory(0x0);
611
612 hxzt->SetDirectory(0x0);
613 hxz->SetDirectory(0x0);
614 hxzTR->SetDirectory(0x0);
615
616 Float_t x,y,z;
617
618 for (Int_t i = 0;i<N;i++)
619 {
620 Double_t r = 84.1+i;
621 tp->PositionAt(i,x,y,z);
622 hxy->Fill(x,y);
623 hxz->Fill(x,z);
624 printf("Rdemanded %f\n",r);
625 printf("tpx %f tpy %f tpz %f Rt =%f\n", x,y,z,TMath::Hypot(x,y));
5b49024d 626
a5556ea5 627 //BUT they are local!!!!
628 t->PropagateTo(r);
629// Double_t phi = t->Phi();
630 Double_t rl = TMath::Hypot(t->GetX(),t->GetY());//real radius
631
632 Double_t alpha = t->GetAlpha();
633 Double_t salpha = TMath::Sin(alpha);
634 Double_t calpha = TMath::Cos(alpha);
635 x = t->GetX()*calpha - t->GetY()*salpha;
636 y = t->GetX()*salpha + t->GetY()*calpha;
637 z = t->GetZ();
638
639 printf("tx %f ty %f tz %f Rt = %f R from XY %f\n",x,y,z,TMath::Hypot(x,y),rl);
640
641 printf("tpz - tz %f\n",z-t->GetZ());
642 printf("\n");
643 hxyt->Fill(x,y);
644 hxzt->Fill(x,z);
645
646 }
647
648 rl->LoadTrackRefs();
649 TTree* treeTR = rl->TreeTR();
650 b = treeTR->GetBranch("TPC");
651
652 TClonesArray* trackrefs = new TClonesArray("AliTrackReference", 100);
653 AliTrackReference* tref;
654 b->SetAddress(&trackrefs);
655
656 Int_t tlab = TMath::Abs(t->GetLabel());
657
658 Int_t netr = (Int_t)treeTR->GetEntries();
659 printf("Found %d entries in TR tree\n",netr);
660
661 for (Int_t e = 0; e < netr; e++)
662 {
663 treeTR->GetEntry(e);
664 tref = (AliTrackReference*)trackrefs->At(0);
665 if (tref == 0x0) continue;
666 if (tref->GetTrack() != tlab) continue;
667
668 printf("Found %d entries in TR array\n",trackrefs->GetEntries());
669
670 for (Int_t i = 0; i < trackrefs->GetEntries(); i++)
671 {
672 tref = (AliTrackReference*)trackrefs->At(i);
673 if (tref->GetTrack() != tlab) continue;
674 x = tref->X();
675 y = tref->Y();
676 z = tref->Z();
677 printf("Track Ref: x %f y %f z %f\n",tref->X(),tref->Y(),tref->Z());
678
679 hxzTR->Fill(x,z);
680 hxyTR->Fill(x,y);
681 for (Int_t j = 1; j < 10; j++)
682 {
683 hxyTR->Fill(x, y+j*0.1);
684 hxyTR->Fill(x, y-j*0.1);
685 hxyTR->Fill(x+j*0.1,y);
686 hxyTR->Fill(x-j*0.1,y);
687
688 hxzTR->Fill(x,z-j*0.1);
689 hxzTR->Fill(x,z+j*0.1);
690 hxzTR->Fill(x-j*0.1,z);
691 hxzTR->Fill(x+j*0.1,z);
692 }
693 }
694 break;
695 }
696 hxz->Draw("");
697// hxzt->Draw("same");
698 hxzTR->Draw("same");
699
700 delete rl;
701}
78d7c6d3 702