Follow rename of method in AliESDtrack.
[u/mrichter/AliRoot.git] / EVE / alice-macros / esd_tracks.C
CommitLineData
5a5a1232 1// $Id$
2
d31ac39d 3Reve::Track* esd_make_track(Reve::TrackRnrStyle* rnrStyle,
4 Int_t index,
5 AliESDtrack* at,
86f12f78 6 AliExternalTrackParam* tp=0)
7{
a8600b56 8 // Helper function
86f12f78 9 Double_t pbuf[3], vbuf[3];
10 Reve::RecTrack rt;
11
12 if(tp == 0) tp = at;
13
d31ac39d 14 rt.label = at->GetLabel();
15 rt.index = index;
86f12f78 16 rt.status = (Int_t) at->GetStatus();
17 rt.sign = tp->GetSign();
18 tp->GetXYZ(vbuf);
19 rt.V.Set(vbuf);
20 tp->GetPxPyPz(pbuf);
21 rt.P.Set(pbuf);
22 Double_t ep = at->GetP(), mc = at->GetMass();
23 rt.beta = ep/TMath::Sqrt(ep*ep + mc*mc);
24
25 Reve::Track* track = new Reve::Track(&rt, rnrStyle);
7be1e8cd 26 //PH The line below is replaced waiting for a fix in Root
27 //PH which permits to use variable siza arguments in CINT
28 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
29 //PH track->SetName(Form("ESDTrack %d", rt.label));
30 //PH track->SetTitle(Form("pT=%.3f, pZ=%.3f; V=(%.3f, %.3f, %.3f)",
31 //PH rt.sign*TMath::Hypot(rt.P.x, rt.P.y), rt.P.z,
32 //PH rt.V.x, rt.V.y, rt.V.z));
33 char form[1000];
54ece6ba 34 sprintf(form,"Track %d", rt.index);
7be1e8cd 35 track->SetName(form);
fc9514f5 36 track->SetStdTitle();
86f12f78 37 return track;
38}
39
a8600b56 40Bool_t gkFixFailedITSExtr = kFALSE;
86f12f78 41
5a5a1232 42Reve::TrackList* esd_tracks(Double_t min_pt=0.1, Double_t max_pt=100)
43{
22aefef8 44 AliESDEvent* esd = Alieve::Event::AssertESD();
5a5a1232 45
46 Double_t minptsq = min_pt*min_pt;
47 Double_t maxptsq = max_pt*max_pt;
48 Double_t ptsq;
49
5a5a1232 50 Reve::TrackList* cont = new Reve::TrackList("ESD Tracks");
51 cont->SetMainColor(Color_t(6));
86f12f78 52 Reve::TrackRnrStyle* rnrStyle = cont->GetRnrStyle();
53 rnrStyle->SetMagField( esd->GetMagneticField() );
5a5a1232 54
5b96ea20 55 gReve->AddRenderElement(cont);
5a5a1232 56
86f12f78 57 Int_t count = 0;
58 Double_t pbuf[3];
d31ac39d 59 for (Int_t n=0; n<esd->GetNumberOfTracks(); n++)
60 {
86f12f78 61 AliESDtrack* at = esd->GetTrack(n);
5a5a1232 62
63 // Here would be sweet to have TObjectFormula.
64 at->GetPxPyPz(pbuf);
65 ptsq = pbuf[0]*pbuf[0] + pbuf[1]*pbuf[1];
66 if(ptsq < minptsq || ptsq > maxptsq)
67 continue;
68
69 ++count;
70
a8600b56 71 // If gkFixFailedITSExtr is TRUE (FALSE by default) and
72 // if ITS refit failed, take track parameters at inner TPC radius.
86f12f78 73 AliExternalTrackParam* tp = at;
74 if (gkFixFailedITSExtr && !at->IsOn(AliESDtrack::kITSrefit)) {
75 tp = at->GetInnerParam();
76 }
5a5a1232 77
d31ac39d 78 Reve::Track* track = esd_make_track(rnrStyle, n, at, tp);
32e219c2 79 track->SetAttLineAttMarker(cont);
80 gReve->AddRenderElement(track, cont);
5a5a1232 81 }
82
7be1e8cd 83 //PH The line below is replaced waiting for a fix in Root
84 //PH which permits to use variable siza arguments in CINT
85 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
86 //PH const Text_t* tooltip = Form("pT ~ (%.2lf, %.2lf), N=%d", min_pt, max_pt, count);
87 char tooltip[1000];
88 sprintf(tooltip,"pT ~ (%.2lf, %.2lf), N=%d", min_pt, max_pt, count);
5a5a1232 89 cont->SetTitle(tooltip); // Not broadcasted automatically ...
5b96ea20 90 cont->UpdateItems();
5a5a1232 91
92 cont->MakeTracks();
32e219c2 93
5b96ea20 94 gReve->Redraw3D();
5a5a1232 95
96 return cont;
97}
98
a8600b56 99/**************************************************************************/
100// esd_tracks_from_array()
101/**************************************************************************/
102
22aefef8 103Reve::TrackList* esd_tracks_from_array(TCollection* col, AliESDEvent* esd=0)
5a5a1232 104{
105 // Retrieves AliESDTrack's from collection.
106 // See example usage with AliAnalysisTrackCuts in the next function.
107
86f12f78 108 if(esd == 0) esd = Alieve::Event::AssertESD();
86f12f78 109
5a5a1232 110 Reve::TrackList* cont = new Reve::TrackList("ESD Tracks");
111 cont->SetMainColor(Color_t(6));
86f12f78 112 Reve::TrackRnrStyle* rnrStyle = cont->GetRnrStyle();
113 rnrStyle->SetMagField( esd->GetMagneticField() );
5a5a1232 114
5b96ea20 115 gReve->AddRenderElement(cont);
5a5a1232 116
5a5a1232 117 Int_t count = 0;
118 TIter next(col);
119 TObject *obj;
d31ac39d 120 while((obj = next()) != 0)
121 {
5a5a1232 122 if(obj->IsA()->InheritsFrom("AliESDtrack") == kFALSE) {
123 Warning("Object '%s', '%s' is not an AliESDtrack.",
124 obj->GetName(), obj->GetTitle());
125 continue;
126 }
127
128 ++count;
86f12f78 129 AliESDtrack* at = (AliESDtrack*) obj;
5a5a1232 130
d31ac39d 131 Reve::Track* track = esd_make_track(rnrStyle, count, at);
32e219c2 132 track->SetAttLineAttMarker(cont);
133 gReve->AddRenderElement(track, cont);
5a5a1232 134 }
135
7be1e8cd 136 //PH The line below is replaced waiting for a fix in Root
137 //PH which permits to use variable siza arguments in CINT
138 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
139 //PH const Text_t* tooltip = Form("N=%d", count);
140 const tooltip[1000];
141 sprintf(tooltip,"N=%d", count);
5a5a1232 142 cont->SetTitle(tooltip); // Not broadcasted automatically ...
5b96ea20 143 cont->UpdateItems();
5a5a1232 144
145 cont->MakeTracks();
32e219c2 146
5b96ea20 147 gReve->Redraw3D();
5a5a1232 148
149 return cont;
150}
151
152void esd_tracks_alianalcuts_demo()
153{
22aefef8 154 AliESDEvent* esd = Alieve::Event::AssertESD();
86f12f78 155 gSystem->Load("libANALYSIS");
5a5a1232 156
157 AliAnalysisTrackCuts atc;
158 atc.SetPtRange(0.1, 5);
159 atc.SetRapRange(-1, 1);
160
86f12f78 161 esd_tracks_from_array(atc.GetAcceptedParticles(esd), esd);
5a5a1232 162}
a8600b56 163
164/**************************************************************************/
165// esd_tracks_vertex_cut
166/**************************************************************************/
167
168Float_t get_sigma_to_vertex(AliESDtrack* esdTrack)
169{
170 // Taken from: PWG0/esdTrackCuts/AliESDtrackCuts.cxx
171 // Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack)
172
173 Float_t b[2];
174 Float_t bRes[2];
175 Float_t bCov[3];
176 esdTrack->GetImpactParameters(b,bCov);
177 if (bCov[0]<=0 || bCov[2]<=0) {
d31ac39d 178 printf("Estimated b resolution lower or equal zero!\n");
a8600b56 179 bCov[0]=0; bCov[2]=0;
180 }
181 bRes[0] = TMath::Sqrt(bCov[0]);
182 bRes[1] = TMath::Sqrt(bCov[2]);
183
184 // -----------------------------------
185 // How to get to a n-sigma cut?
186 //
187 // The accumulated statistics from 0 to d is
188 //
189 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
190 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
191 //
192 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
193 // Can this be expressed in a different way?
194
195 if (bRes[0] == 0 || bRes[1] ==0)
196 return -1;
197
198 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
199
200 // stupid rounding problem screws up everything:
201 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
202 if (TMath::Exp(-d * d / 2) < 1e-10)
203 return 1000;
204
205 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
206 return d;
207}
208
209Reve::RenderElementList* esd_tracks_vertex_cut()
210{
f206464b 211 // Import ESD tracks, separate them into five containers according to
212 // primary-vertex cut and ITS refit status.
a8600b56 213
22aefef8 214 AliESDEvent* esd = Alieve::Event::AssertESD();
a8600b56 215
216 Reve::RenderElementList* cont = new Reve::RenderElementList("ESD Tracks", 0, kTRUE);
217 gReve->AddRenderElement(cont);
f206464b 218 Reve::TrackList *tl[5];
219 Int_t tc[5];
a8600b56 220 Int_t count = 0;
221
222 tl[0] = new Reve::TrackList("Sigma < 3");
a8600b56 223 tc[0] = 0;
224 tl[0]->GetRnrStyle()->SetMagField( esd->GetMagneticField() );
225 tl[0]->SetMainColor(Color_t(3));
32e219c2 226 gReve->AddRenderElement(tl[0], cont);
a8600b56 227
228 tl[1] = new Reve::TrackList("3 < Sigma < 5");
229 tc[1] = 0;
230 tl[1]->GetRnrStyle()->SetMagField( esd->GetMagneticField() );
231 tl[1]->SetMainColor(Color_t(7));
32e219c2 232 gReve->AddRenderElement(tl[1], cont);
a8600b56 233
234 tl[2] = new Reve::TrackList("5 < Sigma");
235 tc[2] = 0;
236 tl[2]->GetRnrStyle()->SetMagField( esd->GetMagneticField() );
237 tl[2]->SetMainColor(Color_t(46));
32e219c2 238 gReve->AddRenderElement(tl[2], cont);
a8600b56 239
846b5c55 240 tl[3] = new Reve::TrackList("no ITS refit; Sigma < 5");
a8600b56 241 tc[3] = 0;
242 tl[3]->GetRnrStyle()->SetMagField( esd->GetMagneticField() );
243 tl[3]->SetMainColor(Color_t(41));
32e219c2 244 gReve->AddRenderElement(tl[3], cont);
a8600b56 245
846b5c55 246 tl[4] = new Reve::TrackList("no ITS refit; Sigma > 5");
f206464b 247 tc[4] = 0;
248 tl[4]->GetRnrStyle()->SetMagField( esd->GetMagneticField() );
249 tl[4]->SetMainColor(Color_t(48));
32e219c2 250 gReve->AddRenderElement(tl[4], cont);
f206464b 251
d31ac39d 252 for (Int_t n=0; n<esd->GetNumberOfTracks(); n++)
253 {
a8600b56 254 AliESDtrack* at = esd->GetTrack(n);
255
256 Float_t s = get_sigma_to_vertex(at);
257 Int_t ti;
258 if (s < 3) ti = 0;
259 else if (s <= 5) ti = 1;
260 else ti = 2;
261
262 AliExternalTrackParam* tp = at;
f206464b 263 // If ITS refit failed, optionally take track parameters at inner
264 // TPC radius and put track in a special container.
a8600b56 265 // This ignores state of gkFixFailedITSExtr (used in esd_tracks()).
266 // Use BOTH functions to compare results.
267 if (!at->IsOn(AliESDtrack::kITSrefit)) {
f206464b 268 // tp = at->GetInnerParam();
269 ti = (ti == 2) ? 4 : 3;
a8600b56 270 }
271
272 Reve::TrackList* tlist = tl[ti];
273 ++tc[ti];
274 ++count;
275
d31ac39d 276 Reve::Track* track = esd_make_track(tlist->GetRnrStyle(), n, at, tp);
32e219c2 277 track->SetAttLineAttMarker(tlist);
d31ac39d 278
7be1e8cd 279 //PH The line below is replaced waiting for a fix in Root
280 //PH which permits to use variable siza arguments in CINT
281 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
282 //PH track->SetName(Form("track %d, sigma=%5.3f", at->GetLabel(), s));
283 char form[1000];
8c62e475 284 sprintf(form,"Track idx=%d, sigma=%5.3f", at->GetID(), s);
7be1e8cd 285 track->SetName(form);
32e219c2 286 gReve->AddRenderElement(track, tlist);
a8600b56 287 }
288
f206464b 289 for (Int_t ti=0; ti<5; ++ti) {
a8600b56 290 Reve::TrackList* tlist = tl[ti];
7be1e8cd 291 //PH The line below is replaced waiting for a fix in Root
292 //PH which permits to use variable siza arguments in CINT
293 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
294 //PH const Text_t* tooltip = Form("N tracks=%d", tc[ti]);
846b5c55 295 //MT Modified somewhat.
296 char buff[1000];
297 sprintf(buff, "%s [%d]", tlist->GetName(), tlist->GetNChildren());
298 tlist->SetName(buff);
299 sprintf(buff, "N tracks=%d", tc[ti]);
300 tlist->SetTitle(buff); // Not broadcasted automatically ...
a8600b56 301 tlist->UpdateItems();
302
303 tlist->MakeTracks();
a8600b56 304 }
7be1e8cd 305 //PH The line below is replaced waiting for a fix in Root
306 //PH which permits to use variable siza arguments in CINT
307 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
308 //PH cont->SetTitle(Form("N all tracks = %d", count));
309 char form[1000];
310 sprintf(form,"N all tracks = %d", count);
311 cont->SetTitle(form);
a8600b56 312 cont->UpdateItems();
313 gReve->Redraw3D();
314
315 return cont;
316}