]>
Commit | Line | Data |
---|---|---|
661eeab1 | 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 | ||
17 | /*********************************************************************** | |
18 | * This code defines the reconstructed v0 visualized with EVE | |
19 | * | |
20 | * Ludovic Gaudichet (gaudichet@to.infn.it) | |
21 | ************************************************************************/ | |
22 | ||
23 | #include "Track.h" | |
24 | #include "V0.h" | |
25 | #include "MCHelixLine.hi" | |
26 | ||
27 | #include <TPolyLine3D.h> | |
28 | #include <TPolyMarker3D.h> | |
29 | #include <TColor.h> | |
30 | ||
31 | #include <Reve/RGTopFrame.h> | |
32 | #include <TCanvas.h> | |
33 | #include <TH1F.h> | |
34 | #include <TH2F.h> | |
35 | ||
36 | #include <vector> | |
37 | ||
38 | using namespace Reve; | |
39 | ||
40 | ||
41 | ||
42 | /*********************************************************************** | |
43 | * | |
44 | * V0 class | |
45 | * | |
46 | ************************************************************************/ | |
47 | ||
48 | const Float_t V0::fgkMassPion2 = 0.13956995*0.13956995; | |
49 | const Float_t V0::fgkMassProton2 = 0.93827231*0.93827231; | |
50 | ||
51 | ClassImp(Reve::V0) | |
52 | ||
53 | V0::V0() : | |
54 | RenderElement(), | |
55 | TPolyMarker3D(1), | |
56 | fV_neg(), | |
57 | fP_neg(), | |
58 | fV_pos(), | |
59 | fP_pos(), | |
60 | fV_v0(), | |
61 | fV0_birth(), | |
62 | fBeta_neg(0), | |
63 | fBeta_pos(0), | |
64 | fLabel_neg(0), | |
65 | fLabel_pos(0), | |
66 | fPathMarksNeg(), | |
67 | fPathMarksPos(), | |
68 | fRnrStyle(0), | |
69 | fPolyLineNeg(), | |
70 | fPolyLinePos(), | |
71 | fPolyLineV0(), | |
72 | fESDIndex(-1), | |
73 | fDaughterDCA(999), | |
74 | fCosPointingAngle(999), | |
75 | fDecayLength(0) | |
76 | {} | |
77 | ||
78 | ||
79 | V0::V0(Reve::RecTrack* tNeg, Reve::RecTrack* tPos, | |
80 | Reve::RecV0* v0, TrackRnrStyle* rs) : | |
81 | RenderElement(), | |
82 | TPolyMarker3D(1), | |
83 | fV_neg(v0->V_neg), | |
84 | fP_neg(v0->P_neg ), | |
85 | fV_pos(v0->V_pos ), | |
86 | fP_pos(v0->P_pos), | |
87 | fV_v0(v0->V_ca), | |
88 | fV0_birth(v0->V0_birth), | |
89 | fBeta_neg(tNeg->beta), | |
90 | fBeta_pos(tPos->beta), | |
91 | fLabel_neg(v0->d_label[0]), | |
92 | fLabel_pos(v0->d_label[1]), | |
93 | fPathMarksNeg(), | |
94 | fPathMarksPos(), | |
95 | fRnrStyle(rs), | |
96 | fPolyLineNeg(), | |
97 | fPolyLinePos(), | |
98 | fPolyLineV0(), | |
99 | fESDIndex(-1), | |
100 | fDaughterDCA(999), | |
101 | fCosPointingAngle(999), | |
102 | fDecayLength(0) | |
103 | { | |
104 | fMarkerColor = fRnrStyle->GetColor(); | |
105 | fPolyLineV0.SetLineColor(fMarkerColor); | |
106 | fPolyLinePos.SetLineColor(2); // red | |
107 | fPolyLineNeg.SetLineColor(7); // light blue | |
108 | ||
109 | fMainColorPtr = &fMarkerColor; | |
110 | fMarkerStyle = 20; | |
111 | fMarkerColor = 5; | |
112 | fMarkerSize = 0.3; | |
113 | } | |
114 | ||
115 | ||
116 | V0::~V0() | |
117 | { | |
118 | for (vpPathMark_i i=fPathMarksNeg.begin(); i!=fPathMarksNeg.end(); ++i) | |
119 | delete *i; | |
120 | for (vpPathMark_i i=fPathMarksPos.begin(); i!=fPathMarksPos.end(); ++i) | |
121 | delete *i; | |
122 | } | |
123 | ||
124 | ||
125 | void V0::Reset(TPolyLine3D* polyLine) { | |
126 | //polyLine->SetPolyLine(n_points); | |
127 | polyLine->SetPolyLine(0); | |
128 | } | |
129 | ||
130 | //______________________________________________________________________ | |
131 | void V0::SetDecayLength(Float_t primx, Float_t primy, Float_t primz) { | |
132 | ||
133 | Float_t dx = fV_v0.x-primx; | |
134 | Float_t dy = fV_v0.y-primy; | |
135 | Float_t dz = fV_v0.z-primz; | |
136 | ||
137 | fDecayLength = sqrt(dx*dx+dy*dy+dz*dz); | |
138 | ||
139 | // This is probably wrong but I can only do this for now | |
140 | Float_t distNorm = fDecayLength/GetMomentum(); | |
141 | fV0_birth.x = fV_v0.x - distNorm*GetPx(); | |
142 | fV0_birth.y = fV_v0.y - distNorm*GetPy(); | |
143 | fV0_birth.z = fV_v0.z - distNorm*GetPz(); | |
144 | } | |
145 | ||
146 | ||
147 | ||
148 | //______________________________________________________________________ | |
149 | void V0::MakeTrack(vpPathMark_t& pathMark, Reve::Vector& vtx, Reve::Vector& p, | |
150 | Int_t charge, Float_t beta, TPolyLine3D& polyLine) { | |
151 | ||
152 | TrackRnrStyle& RS((fRnrStyle != 0) ? *fRnrStyle : TrackRnrStyle::fgDefStyle); | |
153 | ||
154 | Float_t px = p.x, py = p.y, pz = p.z; | |
155 | ||
156 | MCVertex mc_v0; | |
157 | mc_v0.x = vtx.x; | |
158 | mc_v0.y = vtx.y; | |
159 | mc_v0.z = vtx.z; | |
160 | mc_v0.t = 0; | |
161 | ||
162 | std::vector<MCVertex> track_points; | |
163 | Bool_t decay = kFALSE; | |
164 | ||
165 | if ((TMath::Abs(vtx.z) > RS.fMaxZ) || (vtx.x*vtx.x + vtx.y*vtx.y > RS.fMaxR*RS.fMaxR)) | |
166 | goto make_polyline; | |
167 | ||
168 | if (TMath::Abs(RS.fMagField) > 1e-5) { | |
169 | ||
170 | // Charged particle in magnetic field | |
171 | ||
172 | Float_t a = RS.fgkB2C * RS.fMagField * charge; | |
173 | ||
174 | MCHelix helix(fRnrStyle, &mc_v0, TMath::C()*beta, &track_points, a); //m->cm | |
175 | helix.Init(TMath::Sqrt(px*px+py*py), pz); | |
176 | ||
177 | if(!pathMark.empty()){ | |
178 | for(std::vector<Reve::PathMark*>::iterator i=pathMark.begin(); | |
179 | i!=pathMark.end(); ++i) { | |
180 | ||
181 | Reve::PathMark* pm = *i; | |
182 | ||
183 | if(RS.fFitDaughters && pm->type == Reve::PathMark::Daughter){ | |
184 | if(TMath::Abs(pm->V.z) > RS.fMaxZ | |
185 | || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR ) | |
186 | goto helix_bounds; | |
187 | ||
188 | //printf("%s fit daughter \n", fName.Data()); | |
189 | helix.LoopToVertex(p.x, p.y, p.z, pm->V.x, pm->V.y, pm->V.z); | |
190 | p.x -= pm->P.x; | |
191 | p.y -= pm->P.y; | |
192 | p.z -= pm->P.z; | |
193 | } | |
194 | if(RS.fFitDecay && pm->type == Reve::PathMark::Decay){ | |
195 | ||
196 | if(TMath::Abs(pm->V.z) > RS.fMaxZ | |
197 | || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR ) | |
198 | goto helix_bounds; | |
199 | helix.LoopToVertex(p.x, p.y, p.z, pm->V.x, pm->V.y, pm->V.z); | |
200 | decay = true; | |
201 | break; | |
202 | } | |
203 | } | |
204 | } | |
205 | helix_bounds: | |
206 | //go to bounds | |
207 | if(!decay || RS.fFitDecay == kFALSE){ | |
208 | helix.LoopToBounds(px,py,pz); | |
209 | // printf("%s loop to bounds \n",fName.Data() ); | |
210 | } | |
211 | ||
212 | } else { | |
213 | ||
214 | // Neutral particle or no field | |
215 | ||
216 | MCLine line(fRnrStyle, &mc_v0, TMath::C()*beta, &track_points); | |
217 | ||
218 | if(!pathMark.empty()) { | |
219 | for(std::vector<Reve::PathMark*>::iterator i=pathMark.begin(); | |
220 | i!=pathMark.end(); ++i) { | |
221 | Reve::PathMark* pm = *i; | |
222 | ||
223 | if(RS.fFitDaughters && pm->type == Reve::PathMark::Daughter){ | |
224 | if(TMath::Abs(pm->V.z) > RS.fMaxZ | |
225 | || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR ) | |
226 | goto line_bounds; | |
227 | line.GotoVertex(pm->V.x, pm->V.y, pm->V.z); | |
228 | p.x -= pm->P.x; | |
229 | p.y -= pm->P.y; | |
230 | p.z -= pm->P.z; | |
231 | } | |
232 | ||
233 | if(RS.fFitDecay && pm->type == Reve::PathMark::Decay){ | |
234 | if(TMath::Abs(pm->V.z) > RS.fMaxZ | |
235 | || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR ) | |
236 | goto line_bounds; | |
237 | line.GotoVertex(pm->V.x, pm->V.y, pm->V.z); | |
238 | decay = true; | |
239 | break; | |
240 | } | |
241 | } | |
242 | } | |
243 | ||
244 | line_bounds: | |
245 | if(!decay || RS.fFitDecay == kFALSE) | |
246 | line.GotoBounds(px,py,pz); | |
247 | ||
248 | } | |
249 | make_polyline: | |
250 | Reset(&polyLine); | |
251 | for(std::vector<MCVertex>::iterator i=track_points.begin(); | |
252 | i!=track_points.end(); ++i) { | |
253 | polyLine.SetNextPoint(i->x,i->y, i->z); | |
254 | } | |
255 | ||
256 | } | |
257 | ||
258 | ||
259 | //______________________________________________________________________ | |
260 | void V0::MakeV0path() { | |
261 | ||
262 | MCVertex mc_v0; | |
263 | mc_v0.x = fV_v0.x; | |
264 | mc_v0.y = fV_v0.y; | |
265 | mc_v0.z = fV_v0.z; | |
266 | mc_v0.t = 0; | |
267 | ||
268 | std::vector<MCVertex> track_points; | |
269 | MCLine line(fRnrStyle, &mc_v0, TMath::C()*0.99, &track_points); | |
270 | ||
271 | line.GotoVertex(fV0_birth.x,fV0_birth.y,fV0_birth.z); | |
272 | ||
273 | Reset(&fPolyLineV0); | |
274 | for(std::vector<MCVertex>::iterator i=track_points.begin(); | |
275 | i!=track_points.end(); ++i) { | |
276 | fPolyLineV0.SetNextPoint(i->x,i->y, i->z); | |
277 | } | |
278 | } | |
279 | ||
280 | ||
281 | //______________________________________________________________________ | |
282 | void V0::MakeV0() | |
283 | { | |
284 | SetNextPoint(fV_v0.x, fV_v0.y, fV_v0.z); | |
285 | MakeTrack(fPathMarksNeg, fV_neg, fP_neg, -1, fBeta_neg, fPolyLineNeg); | |
286 | MakeTrack(fPathMarksPos, fV_pos, fP_pos, 1, fBeta_pos, fPolyLinePos); | |
287 | MakeV0path(); | |
288 | //fN = fPolyLineNeg.GetN(); | |
289 | } | |
290 | ||
291 | ||
292 | //______________________________________________________________________ | |
87738a65 | 293 | Float_t V0::GetAlphaArmenteros() const |
294 | { | |
661eeab1 | 295 | Float_t posXv0 = fP_pos.x*GetPx() + fP_pos.y*GetPy() + fP_pos.z*GetPz(); |
296 | Float_t negXv0 = fP_neg.x*GetPx() + fP_neg.y*GetPy() + fP_neg.z*GetPz(); | |
297 | ||
298 | if (posXv0+negXv0 > 1.e-39) | |
299 | return (posXv0-negXv0)/(posXv0+negXv0); | |
300 | else return -999; | |
301 | } | |
302 | ||
303 | //______________________________________________________________________ | |
87738a65 | 304 | Float_t V0::GetPtArmenteros() const |
305 | { | |
661eeab1 | 306 | Float_t posXv0 = fP_pos.x*GetPx() + fP_pos.y*GetPy() + fP_pos.z*GetPz(); |
307 | Float_t v0mom2 = GetP2(); | |
308 | ||
309 | if (v0mom2 > 1.e-39) | |
310 | return TMath::Sqrt( GetPosP2() - posXv0*posXv0/v0mom2 ) ; | |
311 | else return -999; | |
312 | } | |
313 | ||
314 | ||
315 | ||
661eeab1 | 316 | /*********************************************************************** |
317 | * | |
318 | * V0List class | |
319 | * | |
320 | ************************************************************************/ | |
321 | ||
322 | ClassImp(Reve::V0List) | |
323 | ||
324 | //______________________________________________________________________ | |
325 | V0List::V0List() : | |
5b457dfa | 326 | RenderElement(), |
661eeab1 | 327 | fTitle(), |
328 | fRnrStyle(0), | |
329 | fRnrDaughters(kTRUE), | |
330 | fRnrV0vtx(kTRUE), | |
331 | fRnrV0path(kTRUE), | |
332 | fNegColor(0), | |
333 | fPosColor(0) | |
334 | { | |
335 | for (Int_t i=0; i<fgkNcutVar; i++) | |
336 | fHist[i] = 0; | |
337 | for (Int_t i=0; i<fgkNcutVar2D; i++) | |
338 | fHist2D[i] = 0; | |
339 | } | |
340 | ||
341 | //______________________________________________________________________ | |
342 | V0List::V0List(TrackRnrStyle* rs) : | |
5b457dfa | 343 | RenderElement(), |
661eeab1 | 344 | fTitle(), |
345 | fRnrStyle(rs), | |
346 | fRnrDaughters(kTRUE), | |
347 | fRnrV0vtx(kTRUE), | |
348 | fRnrV0path(kTRUE), | |
349 | fNegColor(0), | |
350 | fPosColor(0) | |
351 | { | |
352 | Init(); | |
353 | } | |
354 | ||
355 | //______________________________________________________________________ | |
356 | V0List::V0List(const Text_t* name, TrackRnrStyle* rs) : | |
5b457dfa | 357 | RenderElement(), |
661eeab1 | 358 | fTitle(), |
359 | fRnrStyle(rs), | |
360 | fRnrDaughters(kTRUE), | |
361 | fRnrV0vtx(kTRUE), | |
362 | fRnrV0path(kTRUE), | |
363 | fNegColor(0), | |
364 | fPosColor(0) | |
365 | { | |
366 | Init(); | |
367 | SetName(name); | |
368 | } | |
369 | ||
370 | //______________________________________________________________________ | |
371 | void V0List::Init() | |
372 | { | |
373 | ||
374 | if (fRnrStyle== 0) fRnrStyle = new TrackRnrStyle; | |
375 | SetMainColorPtr(&fRnrStyle->fColor); | |
376 | ||
377 | fMin[0] = 0; fMax[0] = 10; // pt | |
378 | fMin[1] = 0; fMax[1] = 10; // K0s mass | |
379 | fMin[2] = 0; fMax[2] = 10; // lambda mass | |
380 | fMin[3] = 0; fMax[3] = 10; // anti-lambda mass | |
381 | fMin[4] = 0; fMax[4] = 10; // daughter DCA | |
382 | fMin[5] = 0.8; fMax[5] = 1; // cos of pointing angle | |
383 | ||
384 | fMin[6] = 0; fMax[6] = 200; // radius | |
385 | fMin[7] = -2; fMax[7] = 2; // PseudoRapidity | |
386 | fMin[8] = 0; fMax[8] = 10; // NegPt | |
387 | fMin[9] = -2; fMax[9] = 2; // NegPseudoRapidity | |
388 | fMin[10] = 0; fMax[10] = 10; // PosPt | |
389 | fMin[11] = -2; fMax[11] = 2; // PosPseudoRapidity | |
390 | fMin[12] = 0; fMax[12] = 1e5; // index, would be good to be able to adjust it | |
391 | ||
392 | char *ch = "ptV0"; | |
393 | fHist[0] = new TH1F(ch,ch, 100, fMin[0], fMax[0]); | |
394 | ch = "K0sMass"; | |
395 | fHist[1] = new TH1F(ch,ch, 100, fMin[1], fMax[1]); | |
396 | ch = "LambdaMass"; | |
397 | fHist[2] = new TH1F(ch,ch, 100, fMin[2], fMax[2]); | |
398 | ch = "AntiLambdaMass"; | |
399 | fHist[3] = new TH1F(ch,ch, 100, fMin[3], fMax[3]); | |
400 | ch = "daughterDCA"; | |
401 | fHist[4] = new TH1F(ch,ch, 100, fMin[4], fMax[4]); | |
402 | ch = "cosPointingAngle"; | |
403 | fHist[5] = new TH1F(ch,ch, 100, fMin[5], fMax[5]); | |
404 | ch = "radius"; | |
405 | fHist[6] = new TH1F(ch,ch, 100, fMin[6], fMax[6]); | |
406 | ch = "PseudoRapidity"; | |
407 | fHist[7] = new TH1F(ch,ch, 100, fMin[7], fMax[7]); | |
408 | ||
409 | ch = "NegPt"; | |
410 | fHist[8] = new TH1F(ch,ch, 100, fMin[8], fMax[8]); | |
411 | ch = "NegPseudoRapidity"; | |
412 | fHist[9] = new TH1F(ch,ch, 100, fMin[9], fMax[9]); | |
413 | ch = "PosPt"; | |
414 | fHist[10] = new TH1F(ch,ch, 100, fMin[10], fMax[10]); | |
415 | ch = "PosPseudoRapidity"; | |
416 | fHist[11] = new TH1F(ch,ch, 100, fMin[11], fMax[11]); | |
417 | ch = "v0Index"; | |
418 | fHist[12] = new TH1F(ch,ch, 100, fMin[12], fMax[12]); | |
419 | ||
420 | ||
421 | fMinX[0] = -1.2; | |
422 | fMaxX[0] = 1.2; | |
423 | fMinY[0] = 0; | |
424 | fMaxY[0] = 0.4; | |
425 | ch = "ArmenterosPodolansky"; | |
426 | fHist2D[0] = new TH2F(ch,ch, 70, fMinX[0], fMaxX[0], 70, | |
427 | fMinY[0], fMaxY[0]); | |
428 | ||
429 | for (Int_t i=0; i<fgkNcutVar; i++) { | |
430 | fHist[i]->GetXaxis()->SetLabelSize(0.07); | |
431 | fHist[i]->GetYaxis()->SetLabelSize(0.07); | |
432 | fHist[i]->SetStats(0); | |
433 | } | |
434 | for (Int_t i=0; i<fgkNcutVar2D; i++) { | |
435 | fHist2D[i]->GetXaxis()->SetLabelSize(0.07); | |
436 | fHist2D[i]->GetYaxis()->SetLabelSize(0.07); | |
437 | fHist2D[i]->SetStats(0); | |
438 | } | |
439 | } | |
440 | ||
441 | //______________________________________________________________________ | |
442 | V0List::~V0List() { | |
443 | ||
444 | for (Int_t i=0; i<fgkNcutVar; i++) | |
445 | if (fHist[i]) delete fHist[i]; | |
446 | for (Int_t i=0; i<fgkNcutVar2D; i++) | |
447 | if (fHist2D[i]) delete fHist2D[i]; | |
448 | ||
449 | } | |
450 | ||
451 | //______________________________________________________________________ | |
452 | void V0List::Paint(Option_t* option) { | |
5b457dfa | 453 | if(fRnrSelf) { |
661eeab1 | 454 | |
455 | if(fRnrV0vtx) { | |
e6ac3950 | 456 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
5b457dfa | 457 | if((*i)->GetRnrSelf()) { |
661eeab1 | 458 | ((V0*)(*i))->Paint(option); |
459 | } | |
460 | } | |
461 | } | |
462 | ||
463 | if(fRnrDaughters) { | |
e6ac3950 | 464 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
5b457dfa | 465 | if((*i)->GetRnrSelf()) { |
661eeab1 | 466 | ((V0*)(*i))->PaintDaughters(option); |
467 | } | |
468 | } | |
469 | } | |
470 | ||
471 | if(fRnrV0path) { | |
e6ac3950 | 472 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
5b457dfa | 473 | if((*i)->GetRnrSelf()) { |
661eeab1 | 474 | ((V0*)(*i))->PaintPath(option); |
475 | } | |
476 | } | |
477 | } | |
478 | } | |
479 | } | |
480 | ||
481 | ||
482 | //______________________________________________________________________ | |
483 | ||
484 | void V0List::AddElement(RenderElement* el) { | |
485 | ||
486 | static const Exc_t eH("V0List::AddElement "); | |
487 | if (dynamic_cast<V0*>(el) == 0) | |
488 | throw(eH + "new element not a V0."); | |
5b457dfa | 489 | RenderElement::AddElement(el); |
661eeab1 | 490 | } |
491 | ||
492 | ||
493 | ||
494 | void V0List::SetRnrV0vtx(Bool_t rnr) { | |
495 | fRnrV0vtx = rnr; | |
496 | gReve->Redraw3D(); | |
497 | } | |
498 | ||
499 | void V0List::SetRnrV0path(Bool_t rnr) { | |
500 | fRnrV0path = rnr; | |
501 | gReve->Redraw3D(); | |
502 | } | |
503 | ||
504 | void V0List::SetRnrDaughters(Bool_t rnr) { | |
505 | fRnrDaughters = rnr; | |
506 | gReve->Redraw3D(); | |
507 | } | |
508 | ||
509 | ||
510 | //______________________________________________________________________ | |
511 | ||
512 | void V0List::MakeV0s() { | |
e6ac3950 | 513 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 514 | ((V0*)(*i))->MakeV0(); |
515 | } | |
516 | gReve->Redraw3D(); | |
517 | } | |
518 | ||
519 | ||
520 | void V0List::MakeMarkers() { | |
521 | gReve->Redraw3D(); | |
522 | } | |
523 | ||
524 | ||
525 | //_________________________________________________________________________ | |
526 | void V0List::AdjustHist(Int_t iHist) { | |
527 | ||
8f14e54c | 528 | if ((iHist<0)||(iHist>=fgkNcutVar)) return; |
661eeab1 | 529 | if (! fHist[iHist]) return; |
530 | ||
531 | TString name = fHist[iHist]->GetName(); | |
532 | Int_t nBin = fHist[iHist]->GetXaxis()->GetNbins(); | |
533 | delete fHist[iHist]; | |
534 | fHist[iHist] = new TH1F(name.Data(), name.Data(), nBin, GetMin(iHist), | |
535 | GetMax(iHist)); | |
536 | fHist[iHist]->GetXaxis()->SetLabelSize(0.07); | |
537 | fHist[iHist]->GetYaxis()->SetLabelSize(0.07); | |
538 | fHist[iHist]->SetStats(0); | |
539 | ||
540 | } | |
541 | ||
542 | ||
543 | //______________________________________________________________________ | |
544 | void V0List::UnFill(V0* v0) { | |
545 | ||
546 | Int_t bin = fHist[0]->GetXaxis()->FindBin(v0->GetPt()); | |
547 | fHist[0]->SetBinContent( bin, fHist[0]->GetBinContent(bin)-1 ); | |
548 | ||
549 | bin = fHist[1]->GetXaxis()->FindBin( v0->GetK0mass() ); | |
550 | fHist[1]->SetBinContent( bin, fHist[1]->GetBinContent(bin)-1 ); | |
551 | ||
552 | bin = fHist[2]->GetXaxis()->FindBin( v0->GetLamMass() ); | |
553 | fHist[2]->SetBinContent( bin, fHist[2]->GetBinContent(bin)-1 ); | |
554 | ||
555 | bin = fHist[3]->GetXaxis()->FindBin( v0->GetAntiLamMass() ); | |
556 | fHist[3]->SetBinContent( bin, fHist[3]->GetBinContent(bin)-1 ); | |
557 | ||
558 | bin = fHist[4]->GetXaxis()->FindBin( v0->GetDaughterDCA() ); | |
559 | fHist[4]->SetBinContent( bin, fHist[4]->GetBinContent(bin)-1 ); | |
560 | ||
561 | bin = fHist[5]->GetXaxis()->FindBin( v0->GetCosPointingAngle() ); | |
562 | fHist[5]->SetBinContent( bin, fHist[5]->GetBinContent(bin)-1 ); | |
563 | ||
564 | ||
565 | bin = fHist[6]->GetXaxis()->FindBin( v0->GetRadius() );// radius | |
566 | fHist[6]->SetBinContent( bin, fHist[6]->GetBinContent(bin)-1 ); | |
567 | ||
568 | bin = fHist[7]->GetXaxis()->FindBin( v0->GetPseudoRapidity() ); // PseudoRapidity | |
569 | fHist[7]->SetBinContent( bin, fHist[7]->GetBinContent(bin)-1 ); | |
570 | ||
571 | bin = fHist[8]->GetXaxis()->FindBin( v0->GetNegPt() );// NegPt | |
572 | fHist[8]->SetBinContent( bin, fHist[8]->GetBinContent(bin)-1 ); | |
573 | ||
574 | bin = fHist[9]->GetXaxis()->FindBin( v0->GetNegPseudoRapidity() );// NegPseudoRapidity | |
575 | fHist[9]->SetBinContent( bin, fHist[9]->GetBinContent(bin)-1 ); | |
576 | ||
577 | bin = fHist[10]->GetXaxis()->FindBin( v0->GetPosPt() ); // PosPt | |
578 | fHist[10]->SetBinContent( bin, fHist[10]->GetBinContent(bin)-1 ); | |
579 | ||
580 | bin = fHist[11]->GetXaxis()->FindBin( v0->GetPosPseudoRapidity() ); // PosPseudoRapidity | |
581 | fHist[11]->SetBinContent( bin, fHist[11]->GetBinContent(bin)-1 ); | |
582 | ||
583 | bin = fHist[12]->GetXaxis()->FindBin( v0->GetESDIndex() ); // ESD index | |
584 | fHist[12]->SetBinContent( bin, fHist[12]->GetBinContent(bin)-1 ); | |
585 | ||
586 | //--- | |
587 | bin = fHist2D[0]->GetXaxis()->FindBin( v0->GetAlphaArmenteros() ); | |
588 | Int_t binY = fHist2D[0]->GetYaxis()->FindBin( v0->GetPtArmenteros() ); | |
589 | fHist2D[0]->SetBinContent( bin, binY, fHist2D[0]->GetBinContent(bin,binY)-1 ); | |
590 | } | |
591 | ||
592 | ||
593 | //______________________________________________________________________ | |
594 | void V0List::Filter(V0* v0) { | |
595 | ||
596 | Float_t pt = v0->GetPt(); | |
597 | if ((pt<fMin[0])||(pt>fMax[0])) return; | |
598 | ||
599 | Float_t k0sMass = v0->GetK0mass(); | |
600 | if ( (k0sMass<fMin[1])||(k0sMass>fMax[1]) ) return; | |
601 | ||
602 | Float_t lamMass = v0->GetLamMass(); | |
603 | if ( (lamMass<fMin[2])||(lamMass>fMax[2]) ) return; | |
604 | ||
605 | Float_t alamMass = v0->GetAntiLamMass(); | |
606 | if ( (alamMass<fMin[3])||(alamMass>fMax[3]) ) return; | |
607 | ||
608 | Float_t daughtDCA = v0->GetDaughterDCA(); | |
609 | if ( (daughtDCA<fMin[4])||(daughtDCA>fMax[4]) ) return; | |
610 | ||
611 | Float_t cosPointing = v0->GetCosPointingAngle(); | |
612 | if ( (cosPointing<fMin[5])||(cosPointing>fMax[5]) ) return; | |
613 | ||
614 | ||
615 | Float_t radius = v0->GetRadius(); | |
616 | if ( (radius<fMin[6])||(radius>fMax[6]) ) return; | |
617 | ||
618 | Float_t pseudoRapidity = v0->GetPseudoRapidity(); | |
619 | if ( (pseudoRapidity<fMin[7])||(pseudoRapidity>fMax[7]) ) return; | |
620 | ||
621 | Float_t negPt = v0->GetNegPt(); | |
622 | if ( (negPt<fMin[8])||(negPt>fMax[8]) ) return; | |
623 | ||
624 | Float_t negPseudoRapidity = v0->GetNegPseudoRapidity(); | |
625 | if ( (negPseudoRapidity<fMin[9])||(negPseudoRapidity>fMax[9]) ) return; | |
626 | ||
627 | Float_t posPt = v0->GetPosPt(); | |
628 | if ( (posPt<fMin[10])||(posPt>fMax[10]) ) return; | |
629 | ||
630 | Float_t posPseudoRapidity = v0->GetPosPseudoRapidity(); | |
631 | if ( (posPseudoRapidity<fMin[11])||(posPseudoRapidity>fMax[11]) ) return; | |
632 | ||
633 | Float_t esdIndex = v0->GetESDIndex(); | |
634 | if ( (esdIndex<fMin[12])||(esdIndex>fMax[12]) ) return; | |
635 | ||
636 | Float_t alphaArm = v0->GetAlphaArmenteros(); | |
637 | // if ( (alphaArm<fMinX[0])||(alphaArm>fMaxX[0]) ) return; | |
638 | ||
639 | Float_t ptArm = v0->GetPtArmenteros(); | |
640 | // if ( (ptArm<fMinY[0])||(ptArm>fMaxY[0]) ) return; | |
641 | ||
5b457dfa | 642 | v0->SetRnrSelf(kTRUE); |
661eeab1 | 643 | fHist[0]->Fill(pt); |
644 | fHist[1]->Fill(k0sMass); | |
645 | fHist[2]->Fill(lamMass); | |
646 | fHist[3]->Fill(alamMass); | |
647 | fHist[4]->Fill(daughtDCA); | |
648 | fHist[5]->Fill(cosPointing); | |
649 | fHist[6]->Fill(radius); | |
650 | fHist[7]->Fill(pseudoRapidity); | |
651 | fHist[8]->Fill(negPt); | |
652 | fHist[9]->Fill(negPseudoRapidity); | |
653 | fHist[10]->Fill(posPt); | |
654 | fHist[11]->Fill(posPseudoRapidity); | |
655 | fHist[12]->Fill(esdIndex); | |
656 | fHist2D[0]->Fill(alphaArm, ptArm); | |
657 | } | |
658 | ||
659 | //______________________________________________________________________ | |
660 | void V0List::FilterAll() { | |
661 | ||
662 | for (Int_t i=0; i<fgkNcutVar; i++) | |
663 | fHist[i]->Reset(); | |
664 | ||
665 | for (Int_t i=0; i<fgkNcutVar2D; i++) | |
666 | fHist2D[i]->Reset(); | |
667 | ||
668 | V0* myV0; | |
e6ac3950 | 669 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 670 | |
671 | myV0 = (V0*)(*i); | |
672 | Filter(myV0); | |
673 | } | |
674 | } | |
675 | ||
676 | ||
677 | //______________________________________________________________________ | |
678 | void V0List::GetV0IndexRange(Int_t &imin, Int_t &imax) { | |
679 | ||
680 | Int_t index; | |
681 | V0* myV0; | |
e6ac3950 | 682 | List_i i=fChildren.begin(); |
661eeab1 | 683 | myV0 = (V0*)(*i); |
684 | index = myV0->GetESDIndex(); | |
685 | imin = index; | |
686 | imax = index; | |
687 | ||
688 | for(; i!=fChildren.end(); ++i) { | |
689 | ||
690 | myV0 = (V0*)(*i); | |
691 | index = myV0->GetESDIndex(); | |
692 | if (index<imin) imin = index; | |
693 | if (index>imax) imax = index; | |
694 | } | |
695 | } | |
696 | ||
697 | ||
698 | //______________________________________________________________________ | |
699 | void V0List::PtFilter(Float_t min, Float_t max) { | |
700 | ||
701 | fMin[0] = min; | |
702 | fMax[0] = max; | |
703 | ||
704 | Float_t val; | |
705 | Bool_t wasSelected; | |
706 | Bool_t isSelected; | |
707 | V0* myV0; | |
708 | ||
e6ac3950 | 709 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 710 | |
711 | myV0 = (V0*)(*i); | |
712 | val = myV0->GetPt(); | |
5b457dfa | 713 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 714 | isSelected = ( (val>=min) && (val<=max) ); |
715 | ||
716 | if (wasSelected) { | |
717 | if (! isSelected) { | |
718 | UnFill(myV0); | |
5b457dfa | 719 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 720 | } |
721 | } else { | |
722 | if (isSelected) Filter(myV0); | |
723 | } | |
724 | } | |
725 | } | |
726 | ||
727 | ||
728 | //______________________________________________________________________ | |
729 | void V0List::K0sMFilter(Float_t min, Float_t max) { | |
730 | ||
731 | fMin[1] = min; | |
732 | fMax[1] = max; | |
733 | ||
734 | Float_t val; | |
735 | Bool_t wasSelected; | |
736 | Bool_t isSelected; | |
737 | V0* myV0; | |
738 | ||
e6ac3950 | 739 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 740 | |
741 | myV0 = (V0*)(*i); | |
742 | val = myV0->GetK0mass(); | |
5b457dfa | 743 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 744 | isSelected = ( (val>=min) && (val<=max) ); |
745 | ||
746 | if (wasSelected) { | |
747 | if (! isSelected) { | |
748 | UnFill(myV0); | |
5b457dfa | 749 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 750 | } |
751 | } else { | |
752 | if (isSelected) Filter(myV0); | |
753 | } | |
754 | } | |
755 | } | |
756 | ||
757 | //______________________________________________________________________ | |
758 | void V0List::LamMFilter(Float_t min, Float_t max) { | |
759 | ||
760 | fMin[2] = min; | |
761 | fMax[2] = max; | |
762 | ||
763 | Float_t val; | |
764 | Bool_t wasSelected; | |
765 | Bool_t isSelected; | |
766 | V0* myV0; | |
767 | ||
e6ac3950 | 768 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 769 | |
770 | myV0 = (V0*)(*i); | |
771 | val = myV0->GetLamMass(); | |
5b457dfa | 772 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 773 | isSelected = ( (val>=min) && (val<=max) ); |
774 | ||
775 | if (wasSelected) { | |
776 | if (! isSelected) { | |
777 | UnFill(myV0); | |
5b457dfa | 778 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 779 | } |
780 | } else { | |
781 | if (isSelected) Filter(myV0); | |
782 | } | |
783 | } | |
784 | } | |
785 | ||
786 | ||
787 | ||
788 | //______________________________________________________________________ | |
789 | void V0List::ALamMFilter(Float_t min, Float_t max) { | |
790 | ||
791 | fMin[3] = min; | |
792 | fMax[3] = max; | |
793 | ||
794 | Float_t val; | |
795 | Bool_t wasSelected; | |
796 | Bool_t isSelected; | |
797 | V0* myV0; | |
798 | ||
e6ac3950 | 799 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 800 | |
801 | myV0 = (V0*)(*i); | |
802 | val = myV0->GetAntiLamMass(); | |
5b457dfa | 803 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 804 | isSelected = ( (val>=min) && (val<=max) ); |
805 | ||
806 | if (wasSelected) { | |
807 | if (! isSelected) { | |
808 | UnFill(myV0); | |
5b457dfa | 809 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 810 | } |
811 | } else { | |
812 | if (isSelected) Filter(myV0); | |
813 | } | |
814 | } | |
815 | } | |
816 | ||
817 | //______________________________________________________________________ | |
818 | void V0List::CosPointingFilter(Float_t min, Float_t max) { | |
819 | ||
820 | fMin[5] = min; | |
821 | fMax[5] = max; | |
822 | ||
823 | Float_t val; | |
824 | Bool_t wasSelected; | |
825 | Bool_t isSelected; | |
826 | V0* myV0; | |
827 | ||
e6ac3950 | 828 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 829 | |
830 | myV0 = (V0*)(*i); | |
831 | val = myV0->GetCosPointingAngle(); | |
5b457dfa | 832 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 833 | isSelected = ( (val>=min) && (val<=max) ); |
834 | ||
835 | if (wasSelected) { | |
836 | if (! isSelected) { | |
837 | UnFill(myV0); | |
5b457dfa | 838 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 839 | } |
840 | } else { | |
841 | if (isSelected) Filter(myV0); | |
842 | } | |
843 | } | |
844 | } | |
845 | ||
846 | //______________________________________________________________________ | |
847 | void V0List::DaughterDCAFilter(Float_t min, Float_t max) { | |
848 | ||
849 | fMin[4] = min; | |
850 | fMax[4] = max; | |
851 | ||
852 | Float_t val; | |
853 | Bool_t wasSelected; | |
854 | Bool_t isSelected; | |
855 | V0* myV0; | |
856 | ||
e6ac3950 | 857 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 858 | |
859 | myV0 = (V0*)(*i); | |
860 | val = myV0->GetDaughterDCA(); | |
5b457dfa | 861 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 862 | isSelected = ( (val>=min) && (val<=max) ); |
863 | ||
864 | if (wasSelected) { | |
865 | if (! isSelected) { | |
866 | UnFill(myV0); | |
5b457dfa | 867 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 868 | } |
869 | } else { | |
870 | if (isSelected) Filter(myV0); | |
871 | } | |
872 | } | |
873 | } | |
874 | ||
875 | //______________________________________________________________________ | |
876 | void V0List::RadiusFilter(Float_t min, Float_t max) { | |
877 | ||
878 | fMin[6] = min; | |
879 | fMax[6] = max; | |
880 | ||
881 | Float_t val; | |
882 | Bool_t wasSelected; | |
883 | Bool_t isSelected; | |
884 | V0* myV0; | |
885 | ||
e6ac3950 | 886 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 887 | |
888 | myV0 = (V0*)(*i); | |
889 | val = myV0->GetRadius(); | |
5b457dfa | 890 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 891 | isSelected = ( (val>=min) && (val<=max) ); |
892 | ||
893 | if (wasSelected) { | |
894 | if (! isSelected) { | |
895 | UnFill(myV0); | |
5b457dfa | 896 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 897 | } |
898 | } else { | |
899 | if (isSelected) Filter(myV0); | |
900 | } | |
901 | } | |
902 | } | |
903 | ||
904 | //______________________________________________________________________ | |
905 | void V0List::EtaFilter(Float_t min, Float_t max) { | |
906 | ||
907 | fMin[7] = min; | |
908 | fMax[7] = max; | |
909 | ||
910 | Float_t val; | |
911 | Bool_t wasSelected; | |
912 | Bool_t isSelected; | |
913 | V0* myV0; | |
914 | ||
e6ac3950 | 915 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 916 | |
917 | myV0 = (V0*)(*i); | |
918 | val = myV0->GetPseudoRapidity(); | |
5b457dfa | 919 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 920 | isSelected = ( (val>=min) && (val<=max) ); |
921 | ||
922 | if (wasSelected) { | |
923 | if (! isSelected) { | |
924 | UnFill(myV0); | |
5b457dfa | 925 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 926 | } |
927 | } else { | |
928 | if (isSelected) Filter(myV0); | |
929 | } | |
930 | } | |
931 | } | |
932 | ||
933 | //______________________________________________________________________ | |
934 | void V0List::NegPtFilter(Float_t min, Float_t max) { | |
935 | ||
936 | fMin[8] = min; | |
937 | fMax[8] = max; | |
938 | ||
939 | Float_t val; | |
940 | Bool_t wasSelected; | |
941 | Bool_t isSelected; | |
942 | V0* myV0; | |
943 | ||
e6ac3950 | 944 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 945 | |
946 | myV0 = (V0*)(*i); | |
947 | val = myV0->GetNegPt(); | |
5b457dfa | 948 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 949 | isSelected = ( (val>=min) && (val<=max) ); |
950 | ||
951 | if (wasSelected) { | |
952 | if (! isSelected) { | |
953 | UnFill(myV0); | |
5b457dfa | 954 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 955 | } |
956 | } else { | |
957 | if (isSelected) Filter(myV0); | |
958 | } | |
959 | } | |
960 | } | |
961 | ||
962 | //______________________________________________________________________ | |
963 | void V0List::NegEtaFilter(Float_t min, Float_t max) { | |
964 | ||
965 | fMin[9] = min; | |
966 | fMax[9] = max; | |
967 | ||
968 | Float_t val; | |
969 | Bool_t wasSelected; | |
970 | Bool_t isSelected; | |
971 | V0* myV0; | |
972 | ||
e6ac3950 | 973 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 974 | |
975 | myV0 = (V0*)(*i); | |
976 | val = myV0->GetNegPseudoRapidity(); | |
5b457dfa | 977 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 978 | isSelected = ( (val>=min) && (val<=max) ); |
979 | ||
980 | if (wasSelected) { | |
981 | if (! isSelected) { | |
982 | UnFill(myV0); | |
5b457dfa | 983 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 984 | } |
985 | } else { | |
986 | if (isSelected) Filter(myV0); | |
987 | } | |
988 | } | |
989 | } | |
990 | ||
991 | //______________________________________________________________________ | |
992 | void V0List::PosPtFilter(Float_t min, Float_t max) { | |
993 | ||
994 | fMin[10] = min; | |
995 | fMax[10] = max; | |
996 | ||
997 | Float_t val; | |
998 | Bool_t wasSelected; | |
999 | Bool_t isSelected; | |
1000 | V0* myV0; | |
1001 | ||
e6ac3950 | 1002 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 1003 | |
1004 | myV0 = (V0*)(*i); | |
1005 | val = myV0->GetPosPt(); | |
5b457dfa | 1006 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 1007 | isSelected = ( (val>=min) && (val<=max) ); |
1008 | ||
1009 | if (wasSelected) { | |
1010 | if (! isSelected) { | |
1011 | UnFill(myV0); | |
5b457dfa | 1012 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 1013 | } |
1014 | } else { | |
1015 | if (isSelected) Filter(myV0); | |
1016 | } | |
1017 | } | |
1018 | } | |
1019 | ||
1020 | //______________________________________________________________________ | |
1021 | void V0List::PosEtaFilter(Float_t min, Float_t max) { | |
1022 | ||
1023 | fMin[11] = min; | |
1024 | fMax[11] = max; | |
1025 | ||
1026 | Float_t val; | |
1027 | Bool_t wasSelected; | |
1028 | Bool_t isSelected; | |
1029 | V0* myV0; | |
1030 | ||
e6ac3950 | 1031 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 1032 | |
1033 | myV0 = (V0*)(*i); | |
1034 | val = myV0->GetPosPseudoRapidity(); | |
5b457dfa | 1035 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 1036 | isSelected = ( (val>=min) && (val<=max) ); |
1037 | ||
1038 | if (wasSelected) { | |
1039 | if (! isSelected) { | |
1040 | UnFill(myV0); | |
5b457dfa | 1041 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 1042 | } |
1043 | } else { | |
1044 | if (isSelected) Filter(myV0); | |
1045 | } | |
1046 | } | |
1047 | } | |
1048 | ||
1049 | //______________________________________________________________________ | |
1050 | void V0List::IndexFilter(Float_t min, Float_t max) { | |
1051 | ||
1052 | fMin[12] = min; | |
1053 | fMax[12] = max; | |
1054 | ||
1055 | Float_t val; | |
1056 | Bool_t wasSelected; | |
1057 | Bool_t isSelected; | |
1058 | V0* myV0; | |
1059 | ||
e6ac3950 | 1060 | for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) { |
661eeab1 | 1061 | |
1062 | myV0 = (V0*)(*i); | |
1063 | val = myV0->GetESDIndex(); | |
5b457dfa | 1064 | wasSelected = myV0->GetRnrSelf(); |
661eeab1 | 1065 | isSelected = ( (val>=min) && (val<=max) ); |
1066 | ||
1067 | if (wasSelected) { | |
1068 | if (! isSelected) { | |
1069 | UnFill(myV0); | |
5b457dfa | 1070 | myV0->SetRnrSelf(isSelected); |
661eeab1 | 1071 | } |
1072 | } else { | |
1073 | if (isSelected) Filter(myV0); | |
1074 | } | |
1075 | } | |
1076 | } | |
1077 |