]>
Commit | Line | Data |
---|---|---|
30a892e3 | 1 | ////////////////////////////////////////////////////////////////////// |
2 | // | |
3 | // $Id$ | |
4 | // | |
5 | // Author: Emanuele Simili | |
6 | // | |
7 | ////////////////////////////////////////////////////////////////////// | |
8 | // | |
9 | //_____________________________________________________________ | |
10 | // | |
11 | // Description: | |
12 | // AliFlowEvent is the basic class to perform flow study in ALICE. | |
13 | // The AliFlowEvent object contains data mebers wich summarize the ESD | |
14 | // informations that are most useful for flow study, together with a | |
15 | // collection of tracks (AliFlowTrackCollection) and reconstructed v0s. | |
16 | // The class also implements member functions to calculate Multiplicity, | |
17 | // Mean P and Pt, Q (the event plane vector), Psi (the event plane angle), | |
18 | // normalized q (Q/sqrt(Mult) and Q/sqrt(sum of weights^2)), moreover, | |
19 | // functions to separate random or eta sub-events, to fill the bayesian vector | |
20 | // of particle abundance (for bayesian p.id calculation), to plug-in weights | |
21 | // (to make the tracks' distribution isotropic and to enhance the resolution), | |
22 | // and to select the tracks that are used in the event plane calculation (for | |
23 | // each harmonic, selection, and subevent. | |
24 | // AliFlowEvent supports the standard analysis as well as the Cumulant | |
25 | // Analysis: method to extract the generating functions for the old and new | |
26 | // cumulant methods is also there. | |
27 | // The Flow Event structure (AliFlowEvent, AliFlowTrack, AliFlowV0 and | |
28 | // AliFlowSelection) is independent from AliRoot, i.e. once the FlowEvents | |
29 | // are filled from the AliESD or the MC KineTree (see AliFlowInterface: here | |
30 | // AliRoot is needed), the flow analysis itself can be performed just under | |
31 | // ROOT (see AliFlowAnalysisMaker). | |
32 | // | |
33 | // AliFlowEvent is adapted from the original class StFlowEvent, succesfully | |
34 | // employed to study flow in the STAR experiment at RICH. | |
35 | // Original Authors: Raimond Snellings & Art Poskanzer | |
36 | // | |
37 | ||
38 | #include "AliFlowEvent.h" | |
39 | #include "AliFlowConstants.h" | |
40 | #include "TVector.h" | |
41 | #include "TVector2.h" | |
42 | #include "TVector3.h" | |
43 | #include <iostream> | |
44 | using namespace std; //required for resolving the 'cout' symbol | |
45 | ||
46 | // - Flags & Sets | |
47 | Bool_t AliFlowEvent::fPtWgt = kFALSE ; // gives pT-based weights | |
48 | Bool_t AliFlowEvent::fEtaWgt = kFALSE ; // gives eta-based weights | |
49 | Bool_t AliFlowEvent::fOnePhiWgt = kTRUE ; // kTRUE --> ENABLEs SINGLE WEIGHT ISTOGRAM , kFALSE --> ENABLEs 3 WEIGHT ISTOGRAMS | |
50 | Bool_t AliFlowEvent::fNoWgt = kFALSE ; // No Weight is used | |
51 | // - Eta Sub-Events (later used to calculate the resolution) | |
52 | Bool_t AliFlowEvent::fEtaSubs = kFALSE ; // makes eta subevents | |
53 | ||
54 | ClassImp(AliFlowEvent) ; | |
55 | //----------------------------------------------------------- | |
56 | AliFlowEvent::AliFlowEvent(Int_t lenght) | |
57 | { | |
58 | // Default constructor: initializes the ObjArray of FlowTracks and FlowV0s, | |
59 | // cleans the internal variables, sets all the weights to 1, sets default flags. | |
60 | ||
61 | fEventID = 0 ; | |
62 | fRunID = 0 ; | |
63 | fOrigMult = 0 ; | |
64 | fL0TriggerWord = 0 ; | |
65 | fCentrality = -1 ; | |
66 | fZDCpart = 0 ; | |
67 | for(int zz=0;zz<3;zz++) { fZDCenergy[zz] = 0. ; } | |
68 | ||
69 | // Make a new track collection | |
70 | fTrackCollection = new TObjArray(lenght) ; | |
71 | fV0Collection = new TObjArray(0) ; | |
72 | ||
73 | // Set Weights Arrays to 1 (default) | |
74 | for(int nS=0;nS<Flow::nSels;nS++) | |
75 | { | |
76 | for(int nH=0;nH<Flow::nHars;nH++) | |
77 | { | |
78 | for(int nP=0;nP<Flow::nPhiBins;nP++) | |
79 | { | |
80 | // enable this with: SetOnePhiWgt() | |
81 | fPhiWgt[nS][nH][nP] = 1. ; // cout << nS << nH << nP << " val " << fPhiWgt[nS][nH][nP] << endl ; | |
82 | // enable these with: SetFirstLastPhiWgt() | |
83 | fPhiWgtPlus[nS][nH][nP] = 1. ; // cout << nS << nH << nP << " val " << fPhiWgtPlus[nS][nH][nP] << endl ; | |
84 | fPhiWgtMinus[nS][nH][nP] = 1. ; // cout << nS << nH << nP << " val " << fPhiWgtMinus[nS][nH][nP] << endl ; | |
85 | fPhiWgtCross[nS][nH][nP] = 1. ; // cout << nS << nH << nP << " val " << fPhiWgtCross[nS][nH][nP] << endl ; | |
86 | } | |
87 | } | |
88 | } | |
89 | //for(int nH=0;nH<Flow::nHars;nH++) { fExtPsi[nH] = 0. ; fExtRes[nH] = 0. ; } | |
90 | ||
91 | // The Expected particles abundance is taken directly from Flow::fBayesian[] (see Bayesian P.Id.) | |
92 | ||
93 | fDone = kFALSE ; | |
94 | } | |
95 | //----------------------------------------------------------- | |
96 | AliFlowEvent::~AliFlowEvent() | |
97 | { | |
98 | // Default distructor: deletes the ObjArrays. | |
99 | ||
100 | fTrackCollection->Delete() ; delete fTrackCollection ; | |
101 | fV0Collection->Delete() ; delete fV0Collection ; | |
102 | } | |
103 | //------------------------------------------------------------- | |
104 | ||
105 | //------------------------------------------------------------- | |
106 | Double_t AliFlowEvent::PhiWeightRaw(Int_t selN, Int_t harN, AliFlowTrack* pFlowTrack) const | |
107 | { | |
108 | // Weight for making the event plane isotropic in the lab. | |
109 | ||
110 | float phi = pFlowTrack->Phi() ; if(phi < 0.) { phi += 2*TMath::Pi() ; } | |
111 | Double_t eta = (Double_t)pFlowTrack->Eta() ; | |
112 | int n = (int)((phi/(2*TMath::Pi()))*Flow::nPhiBins); | |
113 | ||
114 | Double_t phiWgt = 1. ; | |
115 | if(OnePhiWgt()) | |
116 | { | |
117 | phiWgt = (Double_t)fPhiWgt[selN][harN][n]; //cout << "Selection " << selN << " ; Harmonic " << harN << " ; PhiBin " << n << " - Wgt = " << phiWgt << endl ; | |
118 | } | |
119 | else if(FirstLastPhiWgt()) | |
120 | { | |
121 | float zFirstPoint = pFlowTrack->ZFirstPoint(); // float zLastPoint = pFlowTrack->ZLastPoint(); | |
122 | ||
123 | if (zFirstPoint > 0. && eta > 0.) { phiWgt = (Double_t)fPhiWgtPlus[selN][harN][n] ; } | |
124 | else if(zFirstPoint < 0. && eta < 0.) { phiWgt = (Double_t)fPhiWgtMinus[selN][harN][n] ; } | |
125 | else { phiWgt = (Double_t)fPhiWgtCross[selN][harN][n] ; } | |
126 | } | |
127 | ||
128 | // ... weight can be implemented as well for ITS , TRD , PHOS , etc. | |
129 | // | |
130 | // if(...) | |
131 | // { | |
132 | // n = (int)((phi/2*TMath::Pi())*Flow::nPhiBins...) ; | |
133 | // ... | |
134 | // } | |
135 | ||
136 | return phiWgt ; | |
137 | } | |
138 | //------------------------------------------------------------- | |
139 | Double_t AliFlowEvent::Weight(Int_t selN, Int_t harN, AliFlowTrack* pFlowTrack) const | |
140 | { | |
141 | // Weight for enhancing the resolution (eta gives sign +/- for Odd Harmonics) | |
142 | ||
143 | bool oddHar = (harN+1) % 2 ; | |
144 | Double_t phiWgt = 1. ; | |
145 | if(PtWgt()) | |
146 | { | |
147 | Double_t pt = (Double_t)pFlowTrack->Pt(); | |
148 | if(pt < Flow::fPtWgtSaturation) { phiWgt *= pt ; } | |
149 | else { phiWgt *= Flow::fPtWgtSaturation ; } // pt weighting going constant | |
150 | } | |
151 | Double_t eta = (Double_t)pFlowTrack->Eta(); | |
152 | Double_t etaAbs = TMath::Abs(eta); | |
153 | if(EtaWgt() && oddHar) { phiWgt *= etaAbs ; } | |
154 | if(oddHar && eta < 0.) { phiWgt *= -1. ; } | |
155 | ||
156 | return phiWgt ; | |
157 | } | |
158 | //------------------------------------------------------------- | |
159 | Double_t AliFlowEvent::PhiWeight(Int_t selN, Int_t harN, AliFlowTrack* pFlowTrack) const | |
160 | { | |
161 | // Weight for making the event plane isotropic in the lab and enhancing the resolution | |
162 | // (it simply rerurns PhiWeightRaw() * Weight()). If fNoWgt = kTRUE, returns +/-1 , | |
163 | // basing on Sign(eta), for odd harmonics . | |
164 | ||
165 | if(NoWgt()) // no weights (but +/- according to eta) | |
166 | { | |
167 | bool oddHar = (harN+1) % 2 ; | |
168 | if(oddHar) { return TMath::Sign((Double_t)1.,(Double_t)pFlowTrack->Eta()) ; } | |
169 | else { return (Double_t)1. ; } | |
170 | } | |
171 | Double_t phiWgtRaw = PhiWeightRaw(selN, harN, pFlowTrack); | |
172 | Double_t weight = Weight(selN, harN, pFlowTrack); | |
173 | if(Flow::fDebug) { cout << "[PhiWeight]: phiWgtRaw = " << phiWgtRaw << " , weight = " << weight << " , eta = " << pFlowTrack->Eta() << endl ; } | |
174 | ||
175 | return phiWgtRaw * weight; | |
176 | } | |
177 | //------------------------------------------------------------- | |
178 | Int_t AliFlowEvent::Mult(AliFlowSelection* pFlowSelect) | |
179 | { | |
180 | // Multiplicity of tracks in the specified Selection | |
181 | ||
182 | if(fDone) | |
183 | { | |
184 | int sub = pFlowSelect->Sub() ; | |
185 | if(sub<0) { return fMult[pFlowSelect->Sel()][pFlowSelect->Har()] ; } | |
186 | else { return fMultSub[sub][pFlowSelect->Sel()][pFlowSelect->Har()] ; } | |
187 | } | |
188 | // - | |
189 | Int_t mult = 0; | |
190 | Int_t itr ; | |
191 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
192 | { | |
193 | AliFlowTrack* pFlowTrack ; | |
194 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
195 | if(pFlowSelect->Select(pFlowTrack)) { mult++ ; } | |
196 | } | |
197 | return mult; | |
198 | } | |
199 | //------------------------------------------------------------- | |
200 | Float_t AliFlowEvent::MeanPt(AliFlowSelection* pFlowSelect) | |
201 | { | |
202 | // Mean pt of tracks in the specified Selection | |
203 | ||
204 | Double_t meanPt = 0. ; | |
205 | Float_t sumPt = 0. ; | |
206 | UInt_t mult = 0 ; | |
207 | Int_t itr ; | |
208 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
209 | { | |
210 | AliFlowTrack* pFlowTrack ; | |
211 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
212 | if (pFlowSelect->Select(pFlowTrack)) | |
213 | { | |
214 | sumPt += pFlowTrack->Pt(); | |
215 | mult++; | |
216 | } | |
217 | } | |
218 | if(mult) { meanPt = sumPt/(float)mult ; } | |
219 | ||
220 | return meanPt ; | |
221 | } | |
222 | //------------------------------------------------------------- | |
223 | TVector2 AliFlowEvent::Q(AliFlowSelection* pFlowSelect) | |
224 | { | |
225 | // Event plane vector for the specified Selection | |
226 | ||
227 | if(fDone) | |
228 | { | |
229 | int sub = pFlowSelect->Sub() ; | |
230 | if(sub<0) { return fQ[pFlowSelect->Sel()][pFlowSelect->Har()] ; } | |
231 | else { return fQSub[sub][pFlowSelect->Sel()][pFlowSelect->Har()] ; } | |
232 | } | |
233 | // - | |
234 | TVector2 mQ ; | |
235 | Double_t mQx=0. , mQy=0. ; | |
236 | int selN = pFlowSelect->Sel() ; | |
237 | int harN = pFlowSelect->Har() ; | |
238 | double order = (double)(harN + 1) ; | |
239 | ||
240 | Int_t itr ; | |
241 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
242 | { | |
243 | AliFlowTrack* pFlowTrack ; | |
244 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
245 | if(pFlowSelect->Select(pFlowTrack)) | |
246 | { | |
247 | double phiWgt = PhiWeight(selN, harN, pFlowTrack); | |
248 | float phi = pFlowTrack->Phi(); | |
249 | mQx += phiWgt * cos(phi * order) ; | |
250 | mQy += phiWgt * sin(phi * order) ; | |
251 | if(Flow::fDebug) { cout << itr << " phi = " << phi << " , wgt = " << phiWgt << endl ; } | |
252 | } | |
253 | } | |
254 | mQ.Set(mQx, mQy); | |
255 | ||
256 | return mQ; | |
257 | } | |
258 | //------------------------------------------------------------- | |
259 | Float_t AliFlowEvent::Psi(AliFlowSelection* pFlowSelect) | |
260 | { | |
261 | // Event plane angle for the specified Selection | |
262 | ||
263 | int harN = pFlowSelect->Har() ; | |
264 | float order = (float)(harN + 1) ; | |
265 | Float_t psi = 0. ; | |
266 | ||
267 | TVector2 mQ = Q(pFlowSelect); | |
268 | if(mQ.Mod()) // if vector is not 0 | |
269 | { | |
270 | psi= mQ.Phi() / order ; | |
271 | if (psi < 0.) { psi += 2*TMath::Pi() / order ; } | |
272 | } | |
273 | return psi; | |
274 | } | |
275 | //------------------------------------------------------------- | |
276 | TVector2 AliFlowEvent::NormQ(AliFlowSelection* pFlowSelect) | |
277 | { | |
278 | // Normalized Q = Q/sqrt(sum of weights^2) for the specified Selection | |
279 | ||
280 | if(fDone) | |
281 | { | |
282 | TVector2 mQ = fQ[pFlowSelect->Sel()][pFlowSelect->Har()] ; | |
283 | double SumOfWeightSqr = fSumOfWeightSqr[pFlowSelect->Sel()][pFlowSelect->Har()] ; | |
284 | if(SumOfWeightSqr) { mQ /= TMath::Sqrt(SumOfWeightSqr) ; } | |
285 | else { mQ.Set(0.,0.) ; } | |
286 | return mQ ; | |
287 | } | |
288 | // - | |
289 | TVector2 mQ ; | |
290 | Double_t mQx=0. , mQy=0. ; | |
291 | double SumOfWeightSqr = 0 ; | |
292 | int selN = pFlowSelect->Sel() ; | |
293 | int harN = pFlowSelect->Har() ; | |
294 | double order = (double)(harN + 1) ; | |
295 | Int_t itr ; | |
296 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
297 | { | |
298 | AliFlowTrack* pFlowTrack ; | |
299 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
300 | if (pFlowSelect->Select(pFlowTrack)) | |
301 | { | |
302 | double phiWgt = PhiWeight(selN, harN, pFlowTrack); | |
303 | SumOfWeightSqr += phiWgt*phiWgt; | |
304 | ||
305 | float phi = pFlowTrack->Phi(); | |
306 | mQx += phiWgt * cos(phi * order); | |
307 | mQy += phiWgt * sin(phi * order); | |
308 | } | |
309 | } | |
310 | if(SumOfWeightSqr) { mQ.Set(mQx/TMath::Sqrt(SumOfWeightSqr), mQy/TMath::Sqrt(SumOfWeightSqr)); } | |
311 | else { mQ.Set(0.,0.); } | |
312 | ||
313 | return mQ; | |
314 | } | |
315 | //------------------------------------------------------------- | |
316 | Float_t AliFlowEvent::q(AliFlowSelection* pFlowSelect) | |
317 | { | |
318 | // Magnitude of normalized Q vector (without pt or eta weighting) for the specified Selection | |
319 | ||
320 | TVector2 mQ ; | |
321 | Double_t mQx = 0. , mQy = 0. ; | |
322 | int selN = pFlowSelect->Sel() ; | |
323 | int harN = pFlowSelect->Har() ; | |
324 | double order = (double)(harN + 1) ; | |
325 | double SumOfWeightSqr = 0 ; | |
326 | ||
327 | Int_t itr ; | |
328 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
329 | { | |
330 | AliFlowTrack* pFlowTrack ; | |
331 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
332 | if(pFlowSelect->Select(pFlowTrack)) | |
333 | { | |
334 | double phiWgt = PhiWeightRaw(selN, harN, pFlowTrack); // Raw | |
335 | SumOfWeightSqr += phiWgt*phiWgt; | |
336 | float phi = pFlowTrack->Phi(); | |
337 | mQx += phiWgt * cos(phi * order); | |
338 | mQy += phiWgt * sin(phi * order); | |
339 | } | |
340 | } | |
341 | if(SumOfWeightSqr) { mQ.Set(mQx/TMath::Sqrt(SumOfWeightSqr), mQy/TMath::Sqrt(SumOfWeightSqr)); } | |
342 | else { mQ.Set(0.,0.); } | |
343 | ||
344 | return mQ.Mod(); | |
345 | } | |
346 | //----------------------------------------------------------------------- | |
347 | Double_t AliFlowEvent::G_New(AliFlowSelection* pFlowSelect, Double_t Zx, Double_t Zy) | |
348 | { | |
349 | // Generating function for the new cumulant method. Eq. 3 in the Practical Guide | |
350 | ||
351 | int selN = pFlowSelect->Sel(); | |
352 | int harN = pFlowSelect->Har(); | |
353 | double order = (double)(harN + 1); | |
354 | ||
355 | double mult = (double)Mult(pFlowSelect); | |
356 | Double_t theG = 1.; | |
357 | ||
358 | Int_t itr ; | |
359 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
360 | { | |
361 | AliFlowTrack* pFlowTrack ; | |
362 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
363 | if (pFlowSelect->Select(pFlowTrack)) | |
364 | { | |
365 | double phiWgt = PhiWeight(selN, harN, pFlowTrack); | |
366 | float phi = pFlowTrack->Phi(); | |
367 | theG *= (1. + (phiWgt/mult) * (2.* Zx * cos(phi * order) + 2.* Zy * sin(phi * order) ) ); | |
368 | } | |
369 | } | |
370 | return theG; | |
371 | } | |
372 | //----------------------------------------------------------------------- | |
373 | Double_t AliFlowEvent::G_Old(AliFlowSelection* pFlowSelect, Double_t Zx, Double_t Zy) | |
374 | { | |
375 | // Generating function for the old cumulant method (if expanded in Taylor | |
376 | // series, one recovers G_New() in new new cumulant method) | |
377 | ||
378 | TVector2 normQ = NormQ(pFlowSelect); | |
379 | ||
380 | return exp(2*Zx*normQ.X() + 2*Zy*normQ.Y()); | |
381 | } | |
382 | //----------------------------------------------------------------------- | |
383 | Double_t AliFlowEvent::SumWeightSquare(AliFlowSelection* pFlowSelect) | |
384 | { | |
385 | // Return sum of weights^2 for the specified Selection (used for normalization) | |
386 | ||
387 | if(fDone) | |
388 | { | |
389 | return fSumOfWeightSqr[pFlowSelect->Sel()][pFlowSelect->Har()] ; | |
390 | } | |
391 | // - | |
392 | int selN = pFlowSelect->Sel(); | |
393 | int harN = pFlowSelect->Har(); | |
394 | Double_t SumOfWeightSqr = 0; | |
395 | ||
396 | Int_t itr ; | |
397 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
398 | { | |
399 | AliFlowTrack* pFlowTrack ; | |
400 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
401 | if (pFlowSelect->Select(pFlowTrack)) | |
402 | { | |
403 | double phiWgt = PhiWeight(selN, harN, pFlowTrack); | |
404 | SumOfWeightSqr += phiWgt*phiWgt; | |
405 | } | |
406 | } | |
407 | if(SumOfWeightSqr < 0.) { return Mult(pFlowSelect) ; } | |
408 | ||
409 | return SumOfWeightSqr; | |
410 | } | |
411 | //------------------------------------------------------------- | |
412 | Double_t AliFlowEvent::WgtMult_q4(AliFlowSelection* pFlowSelect) | |
413 | { | |
414 | // Used only for the old cumulant method, for getting q4 when weight is on. | |
415 | // Replace multiplicity in Eq.(74b) by this quantity when weight is on. | |
416 | // This is derived based on (A4) in the old cumulant paper. | |
417 | ||
418 | int selN = pFlowSelect->Sel(); | |
419 | int harN = pFlowSelect->Har(); | |
420 | double theMult = 0.; | |
421 | double theMeanWj4 = 0.; | |
422 | double theMeanWj2 = 0.; | |
423 | double theSumOfWgtSqr = 0; | |
424 | double phiWgtSq; | |
425 | ||
426 | Int_t itr ; | |
427 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
428 | { | |
429 | AliFlowTrack* pFlowTrack ; | |
430 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
431 | if (pFlowSelect->Select(pFlowTrack)) | |
432 | { | |
433 | double phiWgt = PhiWeight(selN, harN, pFlowTrack); | |
434 | phiWgtSq = phiWgt*phiWgt; | |
435 | theSumOfWgtSqr += phiWgtSq; | |
436 | theMeanWj4 += phiWgtSq*phiWgtSq; | |
437 | theMult += 1.; | |
438 | } | |
439 | } | |
440 | if (theMult <= 0.) return theMult; | |
441 | ||
442 | theMeanWj4 /= theMult; | |
443 | theMeanWj2 = theSumOfWgtSqr / theMult; | |
444 | ||
445 | return (theSumOfWgtSqr*theSumOfWgtSqr)/(theMult*(-theMeanWj4+2*theMeanWj2*theMeanWj2)); | |
446 | } | |
447 | //------------------------------------------------------------- | |
448 | Double_t AliFlowEvent::WgtMult_q6(AliFlowSelection* pFlowSelect) | |
449 | { | |
450 | // Used only for the old cumulant method. For getting q6 when weight is on. | |
451 | // Replace multiplicity in Eq.(74c) by this quantity when weight is on. | |
452 | // This is derived based on (A4) in the old cumulant paper. | |
453 | ||
454 | int selN = pFlowSelect->Sel(); | |
455 | int harN = pFlowSelect->Har(); | |
456 | double theMult = 0.; | |
457 | double theMeanWj6 = 0.; | |
458 | double theMeanWj4 = 0.; | |
459 | double theMeanWj2 = 0.; | |
460 | double theSumOfWgtSqr = 0; | |
461 | double phiWgtSq; | |
462 | ||
463 | Int_t itr ; | |
464 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
465 | { | |
466 | AliFlowTrack* pFlowTrack ; | |
467 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
468 | if (pFlowSelect->Select(pFlowTrack)) | |
469 | { | |
470 | double phiWgt = PhiWeight(selN, harN, pFlowTrack); | |
471 | phiWgtSq = phiWgt*phiWgt; | |
472 | theSumOfWgtSqr += phiWgtSq; | |
473 | theMeanWj4 += phiWgtSq*phiWgtSq; | |
474 | theMeanWj6 += phiWgtSq*phiWgtSq*phiWgtSq; | |
475 | theMult += 1.; | |
476 | } | |
477 | } | |
478 | if (theMult <= 0.) return theMult*theMult; | |
479 | ||
480 | theMeanWj6 /= theMult; | |
481 | theMeanWj4 /= theMult; | |
482 | theMeanWj2 = theSumOfWgtSqr / theMult; | |
483 | ||
484 | return 4.*(theSumOfWgtSqr*theSumOfWgtSqr*theSumOfWgtSqr)/(theMult*(theMeanWj6-9.*theMeanWj2*theMeanWj4+12.*theMeanWj2*theMeanWj2*theMeanWj2)); | |
485 | } | |
486 | //------------------------------------------------------------- | |
487 | void AliFlowEvent::SetSelections(AliFlowSelection* pFlowSelect) | |
488 | { | |
489 | // Sets the selection of tracks used in the Reaction Plane calculation | |
490 | // for the specific Harmonic and Selection - this does not cut trow away | |
491 | // tracks from the event, neither exclude them from flow study. See also | |
492 | // the AliFlowSelection class. | |
493 | // Strategy of Selection: | |
494 | // For the specific Harmonic and Selection, IF cuts | |
495 | // are defined (such that low<high) and IF the track satisfies them, THEN | |
496 | // the respective track's flag (bool AliFlowTrack::mSelection[har][sel]) | |
497 | // is set kTRUE so that the track, from now on, will be included in the | |
498 | // R.P. determination for that selection. | |
499 | // If NO cuts are defined -> ALL Flags are setted kTRUE (all particle | |
500 | // used for all the Reaction Planes -> no difference in Psi[har][sel]). | |
501 | // ------------------------------------------------------------------- | |
502 | // The first selection (all harmonics) is set kTRUE : no conditions. | |
503 | ||
504 | Int_t itr ; | |
505 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
506 | { | |
507 | AliFlowTrack* pFlowTrack ; | |
508 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
509 | pFlowTrack->ResetSelection() ; // this re-sets all the mSelection flags to 0 | |
510 | ||
511 | // * this sets all the selection n.[0] flag kTRUE (all harmonics) * | |
512 | for(int harN=0;harN<Flow::nHars;harN++) { pFlowTrack->SetSelect(harN,0) ; } | |
513 | ||
514 | // Track need to be Constrainable | |
515 | if(pFlowSelect->ConstrainCut() && !pFlowTrack->IsConstrainable()) continue ; | |
516 | ||
517 | // PID - gets rid of the track with AliFlowTrack::Pid() != AliFlowSelection::Pid() (if there) | |
518 | if(pFlowSelect->Pid()[0] != '\0') | |
519 | { | |
520 | if(strstr(pFlowSelect->Pid(), "h")!=0) | |
521 | { | |
522 | int charge = pFlowTrack->Charge() ; | |
523 | if(strcmp("h+",pFlowSelect->Pid())==0 && charge != 1) continue; | |
524 | if(strcmp("h-",pFlowSelect->Pid())==0 && charge != -1) continue; | |
525 | } | |
526 | else | |
527 | { | |
528 | Char_t pid[10]; | |
529 | strcpy(pid, pFlowTrack->Pid()); | |
530 | if(strstr(pid, pFlowSelect->Pid())==0) continue; | |
531 | } | |
532 | } | |
533 | double eta = (double)(pFlowTrack->Eta()); | |
534 | float Pt = pFlowTrack->Pt(); | |
535 | float gDca = pFlowTrack->Dca() ; | |
536 | ||
537 | // Global DCA - gets rid of the track with DCA outside the range | |
538 | if((pFlowSelect->DcaGlobalCutHi()>pFlowSelect->DcaGlobalCutLo()) && (gDca<pFlowSelect->DcaGlobalCutLo() || gDca>pFlowSelect->DcaGlobalCutHi())) continue ; | |
539 | ||
540 | // Pt & Eta - this is done differently for different Harmonic & Selection | |
541 | for(int selN = 1; selN < Flow::nSels; selN++) // not even consider the 0th selection (no cut applied there) | |
542 | { | |
543 | // min. TPC hits required | |
544 | if(pFlowSelect->NhitsCut(selN) && (pFlowTrack->FitPtsTPC()<pFlowSelect->NhitsCut(selN))) continue ; | |
545 | ||
546 | for(int harN = 0; harN < Flow::nHars; harN++) | |
547 | { | |
548 | // Eta - gets rid of the track with Eta outside the range | |
549 | if((pFlowSelect->EtaCutHi(harN%2,selN)>pFlowSelect->EtaCutLo(harN%2,selN)) && (TMath::Abs(eta)<pFlowSelect->EtaCutLo(harN%2,selN) || TMath::Abs(eta)>pFlowSelect->EtaCutHi(harN%2,selN))) continue ; | |
550 | // Pt - gets rid of the track with Pt outside the range | |
551 | if((pFlowSelect->PtCutHi(harN%2,selN)>pFlowSelect->PtCutLo(harN%2,selN)) && (Pt<pFlowSelect->PtCutLo(harN%2,selN) || Pt>pFlowSelect->PtCutHi(harN%2,selN))) continue ; | |
552 | ||
553 | pFlowTrack->SetSelect(harN, selN) ; // if cuts defined (low<high) && track is in the range -> Set [har][sel] Flag ON | |
554 | ||
555 | if(Flow::fDebug) | |
556 | { | |
557 | cout << " harN " << harN%2 << " , selN " << selN << " - si" << endl ; | |
558 | if(pFlowSelect->Pid()[0] != '\0') { cout << " track: pid " << pFlowTrack->Pid() << " = "<< pFlowSelect->Pid() << endl ; } | |
559 | cout << " track: dca " << pFlowSelect->DcaGlobalCutLo() << " < " << gDca << " < " << pFlowSelect->DcaGlobalCutHi() << endl ; | |
560 | cout << " track: eta " << pFlowSelect->EtaCutLo(harN,selN) << " < |" << eta << "| < " << pFlowSelect->EtaCutHi(harN,selN) << endl ; | |
561 | cout << " track: pT " << pFlowSelect->PtCutLo(harN,selN) << " < " << Pt << " < " << pFlowSelect->PtCutHi(harN,selN) << endl ; | |
562 | pFlowTrack->PrintSelection() ; | |
563 | } | |
564 | } | |
565 | } | |
566 | } | |
567 | } | |
568 | //------------------------------------------------------------- | |
569 | void AliFlowEvent::SetPids() | |
570 | { | |
571 | // Re-sets the tracks P.id. (using the current Flow::fBayesian[] array) | |
572 | ||
573 | const Int_t code[] = {11,13,211,321,2212,10010020} ; | |
574 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
575 | { | |
576 | AliFlowTrack* pFlowTrack ; | |
577 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
578 | TVector bayPid = pFlowTrack->PidProbs() ; | |
579 | Int_t maxN = 2 ; // if No id. -> then is a Pi | |
580 | Float_t pid_max = bayPid[2] ; // (if all equal, Pi probability get's the advantage to be the first) | |
581 | for(Int_t nP=0;nP<Flow::nPid;nP++) | |
582 | { | |
583 | if(bayPid[nP]>pid_max) { maxN = nP ; pid_max = bayPid[nP] ; } | |
584 | } | |
585 | Int_t pdg_code = TMath::Sign(code[maxN],pFlowTrack->Charge()) ; | |
586 | pFlowTrack->SetMostLikelihoodPID(pdg_code); | |
587 | } | |
588 | } | |
589 | //------------------------------------------------------------- | |
590 | void AliFlowEvent::MakeSubEvents() | |
591 | { | |
592 | // Make random or eta sub-events | |
593 | ||
594 | if(EtaSubs()) { MakeEtaSubEvents() ; } | |
595 | else { MakeRndSubEvents() ; } | |
596 | } | |
597 | //------------------------------------------------------------- | |
598 | void AliFlowEvent::MakeRndSubEvents() | |
599 | { | |
600 | // Make random subevents | |
601 | ||
602 | int eventMult[Flow::nHars][Flow::nSels] = {{0}}; | |
603 | int harN, selN, subN = 0; | |
604 | ||
605 | // loop to count the total number of tracks for each selection | |
606 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
607 | { | |
608 | AliFlowTrack* pFlowTrack ; | |
609 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
610 | for (selN = 0; selN < Flow::nSels; selN++) | |
611 | { | |
612 | for (harN = 0; harN < Flow::nHars; harN++) | |
613 | { | |
614 | if(pFlowTrack->Select(harN, selN)) { eventMult[harN][selN]++ ; } | |
615 | } | |
616 | } | |
617 | } | |
618 | // loop to set the SubEvent member | |
619 | for (selN = 0; selN < Flow::nSels; selN++) | |
620 | { | |
621 | for (harN = 0; harN < Flow::nHars; harN++) | |
622 | { | |
623 | int subEventMult = eventMult[harN][selN] / Flow::nSubs; | |
624 | if (subEventMult) | |
625 | { | |
626 | subN = 0; | |
627 | int countN = 0; | |
628 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
629 | { | |
630 | AliFlowTrack* pFlowTrack ; | |
631 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
632 | if(pFlowTrack->Select(harN, selN)) | |
633 | { | |
634 | pFlowTrack->SetSubevent(harN, selN, subN); | |
635 | countN++; | |
636 | if (countN % subEventMult == 0.) { subN++ ; } | |
637 | } | |
638 | } | |
639 | } | |
640 | } | |
641 | } | |
642 | return ; | |
643 | } | |
644 | //------------------------------------------------------------- | |
645 | void AliFlowEvent::MakeEtaSubEvents() | |
646 | { | |
647 | // Make subevents for positive and negative eta | |
648 | // (when done, fEtaSubs flag setted to kTRUE). | |
649 | ||
650 | int harN, selN = 0; | |
651 | // loop to set the SubEvent member | |
652 | for (selN = 0; selN < Flow::nSels; selN++) | |
653 | { | |
654 | for (harN = 0; harN < Flow::nHars; harN++) | |
655 | { | |
656 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
657 | { | |
658 | AliFlowTrack* pFlowTrack ; | |
659 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
660 | if(pFlowTrack->Select(harN, selN)) | |
661 | { | |
662 | float eta = pFlowTrack->Eta(); | |
663 | if (eta > 0.) { pFlowTrack->SetSubevent(harN, selN, 0) ; } | |
664 | else { pFlowTrack->SetSubevent(harN, selN, 1) ; } | |
665 | } | |
666 | } | |
667 | } | |
668 | } | |
669 | } | |
670 | //------------------------------------------------------------- | |
671 | void AliFlowEvent::RandomShuffle() | |
672 | { | |
673 | // Randomly re-shuffles the tracks in the array; if a track is not | |
674 | // primary, the reference carried by the reconstructed mother (in | |
675 | // the v0 array) is changed accordingly. | |
676 | ||
677 | Int_t tot = 0 ; | |
678 | UInt_t imax = TrackCollection()->GetEntries() ; | |
679 | TRandom* rnd = new TRandom(0) ; // SetSeed(0) ; | |
680 | TObjArray* newTrackCollection = new TObjArray(imax) ; | |
681 | ||
682 | // re-arranges the ObjArray (TrackCollection()) | |
683 | for(UInt_t itr=0;itr<imax;itr++) | |
684 | { | |
685 | AliFlowTrack* pFlowTrack ; | |
686 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
687 | ||
688 | UInt_t rndNumber = rnd->Integer(imax) ; | |
689 | Bool_t put = kFALSE ; | |
690 | while(!put) | |
691 | { | |
692 | if(!newTrackCollection->At(rndNumber)) | |
693 | { | |
694 | newTrackCollection->AddAt(pFlowTrack, rndNumber) ; | |
695 | put = kTRUE ; tot++ ; if(Flow::fDebug) { cout << " " << itr << " --> " << rndNumber << endl ; } | |
696 | } | |
697 | else | |
698 | { | |
699 | rndNumber++ ; if(rndNumber>=imax) { rndNumber -= imax ; } | |
700 | } | |
701 | } | |
702 | } | |
703 | if(Flow::fDebug) { cout << "* RandomShuffle() : " << tot << "/" << imax << " flow tracks have been shuffled " << endl ; } | |
704 | fTrackCollection = newTrackCollection ; | |
705 | } | |
706 | //----------------------------------------------------------------------- | |
707 | UInt_t AliFlowEvent::Centrality() | |
708 | { | |
709 | // Returns the Centrality Class as stored | |
710 | ||
711 | if(fCentrality < 0) { SetCentrality() ; } | |
712 | return fCentrality ; | |
713 | } | |
714 | //----------------------------------------------------------------------- | |
715 | void AliFlowEvent::SetCentrality(Int_t cent) | |
716 | { | |
717 | // Set the Centrality Classes to "cent" . | |
718 | ||
719 | fCentrality = cent ; | |
720 | } | |
721 | //----------------------------------------------------------------------- | |
722 | void AliFlowEvent::SetCentrality() | |
723 | { | |
724 | // Sets the Centrality Classes basing on Multiplicity at mid rapidity, | |
725 | // ... (ZDC information can be added) . | |
726 | ||
727 | Float_t* cent ; | |
728 | Int_t tracks = MultEta() ; | |
729 | ||
730 | if(RunID() == -1) | |
731 | { | |
732 | cent = Flow::fCentNorm ; | |
733 | //if centrality classes are not defined, does it now (with CentNorm & MaxMult) | |
734 | if(cent[Flow::nCents-1] <= 1) | |
735 | { | |
736 | for(Int_t ic=0;ic<Flow::nCents;ic++) | |
737 | { | |
738 | cent[ic] *= Flow::fMaxMult ; | |
739 | if(Flow::fDebug) { cout << "Centrality[" << ic << "] = " << cent[ic] << " . " << endl ; } | |
740 | } | |
741 | } | |
742 | } | |
743 | else if((RunID() != -1) && (CenterOfMassEnergy() == 5500.)) | |
744 | { | |
745 | cent = (Float_t*)Flow::fCent0 ; | |
746 | } | |
747 | else // other definition of centrality are possible... | |
748 | { | |
749 | cent = (Float_t*)Flow::fCent0 ; | |
750 | } | |
751 | if (tracks < cent[0]) { fCentrality = 0; } | |
752 | else if (tracks < cent[1]) { fCentrality = 1; } | |
753 | else if (tracks < cent[2]) { fCentrality = 2; } | |
754 | else if (tracks < cent[3]) { fCentrality = 3; } | |
755 | else if (tracks < cent[4]) { fCentrality = 4; } | |
756 | else if (tracks < cent[5]) { fCentrality = 5; } | |
757 | else if (tracks < cent[6]) { fCentrality = 6; } | |
758 | else if (tracks < cent[7]) { fCentrality = 7; } | |
759 | else if (tracks < cent[8]) { fCentrality = 8; } | |
760 | else { fCentrality = 9; } | |
761 | ||
762 | if(Flow::fDebug) { cout << " * Centrality Class : " << fCentrality << " . " << endl ; } | |
763 | } | |
764 | //----------------------------------------------------------------------- | |
765 | void AliFlowEvent::Bayesian(Double_t bayes[Flow::nPid]) | |
766 | { | |
767 | // Returns bayesian array of particles' abundances (from Flow::) | |
768 | ||
769 | for(Int_t i=0;i<Flow::nPid;i++) { bayes[i] = Flow::fBayesian[i] ; } | |
770 | } | |
771 | //----------------------------------------------------------------------- | |
772 | TVector AliFlowEvent::Bayesian() | |
773 | { | |
774 | TVector bayes(Flow::nPid) ; | |
775 | for(Int_t i=0;i<Flow::nPid;i++) { bayes[i] = Flow::fBayesian[i] ; } | |
776 | return bayes ; | |
777 | } | |
778 | //----------------------------------------------------------------------- | |
779 | void AliFlowEvent::PrintFlagList() const | |
780 | { | |
781 | // Prints the list of selection cuts ( called in AliFlowInterface::Finish() ) | |
782 | ||
783 | cout << "#######################################################" << endl; | |
784 | cout << "# Weighting and Striping:" << endl; | |
785 | if(PtWgt()) | |
786 | { | |
787 | cout << "# PtWgt = kTRUE " << endl ; // (also for output of PhiWgt file?) | |
788 | cout << "# PtWgt Saturation = " << Flow::fPtWgtSaturation << endl; | |
789 | } | |
790 | else | |
791 | { | |
792 | cout << "# PtWgt = kFALSE" << endl; | |
793 | } | |
794 | if(EtaWgt()) | |
795 | { | |
796 | cout << "# EtaWgt = kTRUE " << endl ; // (also for output of PhiWgt file for odd harmonics?) | |
797 | } | |
798 | else | |
799 | { | |
800 | cout << "# EtaWgt = kFALSE" << endl; | |
801 | } | |
802 | cout << "#######################################################" << endl; | |
803 | } | |
804 | //----------------------------------------------------------------------- | |
805 | void AliFlowEvent::SetEventID(const Int_t& id) | |
806 | { | |
807 | // Sets Event ID and the Event name (name = evtNumber_runId) | |
808 | ||
809 | fEventID = id ; | |
810 | TString name = "" ; | |
811 | name += fEventID ; | |
812 | if(fRunID) { name += "_" ; name += fRunID ; } | |
813 | SetName(name) ; // from TNamed::SetName | |
814 | } | |
815 | //----------------------------------------------------------------------- | |
816 | Int_t AliFlowEvent::MultEta() | |
817 | { | |
818 | // Returns the multiplicity in the interval |eta|<(Flow::fEetaMid), used | |
819 | // for centrality measurement (see centrality classes in fCentrality) . | |
820 | ||
821 | Int_t goodtracks = 0 ; | |
822 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
823 | { | |
824 | AliFlowTrack* pFlowTrack ; | |
825 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
826 | if((pFlowTrack->Charge()) && (TMath::Abs(pFlowTrack->Eta())<Flow::fEtaMid)) { goodtracks++ ; } | |
827 | } | |
828 | return goodtracks ; | |
829 | } | |
830 | //----------------------------------------------------------------------- | |
831 | Int_t AliFlowEvent::UncorrNegMult(Float_t eta) const | |
832 | { | |
833 | // Negative multiplicity in the interval (-eta..eta) | |
834 | // (default is Flow::fEetaGood = 0.9) | |
835 | ||
836 | Int_t negMult = 0 ; | |
837 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
838 | { | |
839 | AliFlowTrack* pFlowTrack ; | |
840 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
841 | if((pFlowTrack->Charge()<0) && (TMath::Abs(pFlowTrack->Eta())<TMath::Abs(eta))) { negMult++ ; } | |
842 | delete pFlowTrack ; | |
843 | } | |
844 | return negMult; | |
845 | } | |
846 | //----------------------------------------------------------------------- | |
847 | Int_t AliFlowEvent::UncorrPosMult(Float_t eta) const | |
848 | { | |
849 | // Positive multiplicity in the interval (-eta..eta) | |
850 | // (default is Flow::fEetaGood = 0.9) | |
851 | ||
852 | Int_t posMult ; | |
853 | for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) | |
854 | { | |
855 | AliFlowTrack* pFlowTrack ; | |
856 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
857 | if((pFlowTrack->Charge()>0) && (TMath::Abs(pFlowTrack->Eta())<TMath::Abs(eta))) { posMult++ ; } | |
858 | delete pFlowTrack ; | |
859 | } | |
860 | return posMult; | |
861 | } | |
862 | //----------------------------------------------------------------------- | |
863 | void AliFlowEvent::VertexPos(Float_t vtx[3]) const | |
864 | { | |
865 | for(Int_t ii=0;ii<3;ii++) { vtx[ii] = fVertexPos[ii] ; } | |
866 | } | |
867 | //----------------------------------------------------------------------- | |
868 | TVector3 AliFlowEvent::VertexPos() const | |
869 | { | |
870 | Float_t vertex[3] ; | |
871 | VertexPos(vertex) ; | |
872 | return TVector3(vertex) ; | |
873 | } | |
874 | //----------------------------------------------------------------------- | |
875 | void AliFlowEvent::SetVertexPos(Float_t v1,Float_t v2,Float_t v3) | |
876 | { | |
877 | fVertexPos[0] = v1 ; fVertexPos[1] = v2 ; fVertexPos[2] = v3 ; | |
878 | } | |
879 | //----------------------------------------------------------------------- | |
880 | void AliFlowEvent::MakeAll() | |
881 | { | |
882 | // calculates all quantities in 1 shoot ... | |
883 | // ... | |
884 | ||
885 | Double_t mQx[Flow::nSels][Flow::nHars] ; | |
886 | Double_t mQy[Flow::nSels][Flow::nHars] ; | |
887 | Double_t mQxSub[Flow::nSubs][Flow::nSels][Flow::nHars] ; | |
888 | Double_t mQySub[Flow::nSubs][Flow::nSels][Flow::nHars] ; | |
889 | // - | |
890 | int selN, harN, subN ; | |
891 | for(selN=0;selN<Flow::nSels;selN++) | |
892 | { | |
893 | for(harN=0;harN<Flow::nHars;harN++) | |
894 | { | |
895 | mQx[selN][harN] = 0. ; | |
896 | mQy[selN][harN] = 0. ; | |
897 | fMult[selN][harN] = 0 ; | |
898 | fSumOfWeightSqr[selN][harN] = 0. ; | |
899 | for(subN=0;subN<Flow::nSubs;subN++) | |
900 | { | |
901 | mQxSub[subN][selN][harN] = 0. ; | |
902 | mQySub[subN][selN][harN] = 0. ; | |
903 | fMultSub[subN][selN][harN] = 0 ; | |
904 | } | |
905 | } | |
906 | } | |
907 | ||
908 | double order = 0. ; | |
909 | double phiWgt = 0. ; | |
910 | float phi = 0. ; | |
911 | // - | |
912 | int itr ; | |
913 | for(itr=0;itr<TrackCollection()->GetEntries();itr++) | |
914 | { | |
915 | AliFlowTrack* pFlowTrack ; | |
916 | pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ; | |
917 | phi = pFlowTrack->Phi(); | |
918 | for(selN=0;selN<Flow::nSels;selN++) | |
919 | { | |
920 | for(harN=0;harN<Flow::nHars;harN++) | |
921 | { | |
922 | order = (double)(harN+1) ; | |
923 | if(pFlowTrack->Select(harN,selN)) | |
924 | { | |
925 | phiWgt = PhiWeight(selN,harN,pFlowTrack) ; | |
926 | fSumOfWeightSqr[selN][harN] += phiWgt*phiWgt ; | |
927 | mQx[selN][harN] += phiWgt * cos(phi * order) ; | |
928 | mQy[selN][harN] += phiWgt * sin(phi * order) ; | |
929 | fMult[selN][harN]++ ; | |
930 | for(subN=0;subN<Flow::nSubs;subN++) | |
931 | { | |
932 | if(pFlowTrack->Select(harN,selN,subN)) | |
933 | { | |
934 | mQxSub[subN][selN][harN] += phiWgt * cos(phi * order) ; | |
935 | mQySub[subN][selN][harN] += phiWgt * sin(phi * order) ; | |
936 | fMultSub[subN][selN][harN]++ ; | |
937 | } | |
938 | } // sub | |
939 | } // if | |
940 | } // har | |
941 | } // sel | |
942 | } //itr | |
943 | ||
944 | for(selN=0;selN<Flow::nSels;selN++) | |
945 | { | |
946 | for(harN=0;harN<Flow::nHars;harN++) | |
947 | { | |
948 | fQ[selN][harN].Set(mQx[selN][harN],mQy[selN][harN]) ; | |
949 | for(subN=0;subN<Flow::nSubs;subN++) { fQSub[subN][selN][harN].Set(mQxSub[subN][selN][harN],mQySub[subN][selN][harN]) ; } | |
950 | } | |
951 | } | |
952 | ||
953 | fDone = kTRUE ; | |
954 | } | |
955 | // //----------------------------------------------------------------------- | |
956 | // Float_t AliFlowEvent::ExtPsi(Int_t harN) const | |
957 | // { | |
958 | // // external R.P. angle (check input source...) | |
959 | // | |
960 | // if(harN<Flow::nHars) { return fExtPsi[harN] ; } | |
961 | // else | |
962 | // { | |
963 | // cout << "AliFlowEvent::ExtPsi(" << harN << ") : harmonic " << harN+1 << " is not there !" << endl ; | |
964 | // return 0. ; | |
965 | // } | |
966 | // } | |
967 | // //----------------------------------------------------------------------- | |
968 | // Float_t AliFlowEvent::ExtRes(Int_t harN) const | |
969 | // { | |
970 | // // external R.P. resolution (check input source...) | |
971 | // | |
972 | // if(harN<Flow::nHars) { return fExtRes[harN] ; } | |
973 | // else | |
974 | // { | |
975 | // cout << "AliFlowEvent::ExtRes(" << harN << ") : harmonic " << harN+1 << " is not there !" << endl ; | |
976 | // return 0. ; | |
977 | // } | |
978 | // } | |
979 | // //----------------------------------------------------------------------- | |
980 | // void AliFlowEvent::SetExtPsi(Int_t harN,Float_t psi) | |
981 | // { | |
982 | // if(harN<Flow::nHars) { fExtPsi[harN] = psi ; } | |
983 | // } | |
984 | // //----------------------------------------------------------------------- | |
985 | // void AliFlowEvent::SetExtRes(Int_t harN,Float_t res) | |
986 | // { | |
987 | // if(harN<Flow::nHars) { fExtRes[harN] = res ; } | |
988 | // } | |
989 | // //----------------------------------------------------------------------- | |
990 | ||
991 | ||
992 | ||
993 | //----------------------------------------------------------------------- | |
994 | // - those will go into the .h as inline functions ..................... | |
995 | //----------------------------------------------------------------------- | |
996 | TObjArray* AliFlowEvent::TrackCollection() const { return fTrackCollection; } | |
997 | TObjArray* AliFlowEvent::V0Collection() const { return fV0Collection; } | |
998 | //----------------------------------------------------------------------- | |
999 | Int_t AliFlowEvent::EventID() const { return fEventID; } | |
1000 | Int_t AliFlowEvent::RunID() const { return fRunID; } | |
1001 | Double_t AliFlowEvent::CenterOfMassEnergy() const { return Flow::fCenterOfMassEnergy ; } | |
1002 | Double_t AliFlowEvent::MagneticField() const { return Flow::fMagneticField ; } | |
1003 | Short_t AliFlowEvent::BeamMassNumberEast() const { return Flow::fBeamMassNumberEast ; } | |
1004 | Short_t AliFlowEvent::BeamMassNumberWest() const { return Flow::fBeamMassNumberWest ; } | |
1005 | UInt_t AliFlowEvent::OrigMult() const { return fOrigMult; } | |
1006 | Long_t AliFlowEvent::L0TriggerWord() const { return fL0TriggerWord; } | |
1007 | Int_t AliFlowEvent::V0Mult() const { return V0Collection()->GetEntries() ; } | |
1008 | Int_t AliFlowEvent::FlowEventMult() const { return TrackCollection()->GetEntries() ; } | |
1009 | Int_t AliFlowEvent::ZDCpart() const { return fZDCpart; } | |
1010 | Float_t AliFlowEvent::ZDCenergy(Int_t npem) const { return fZDCenergy[npem]; } | |
1011 | Float_t AliFlowEvent::PtWgtSaturation() const { return Flow::fPtWgtSaturation; } | |
1012 | Bool_t AliFlowEvent::PtWgt() const { return fPtWgt; } | |
1013 | Bool_t AliFlowEvent::EtaWgt() const { return fEtaWgt; } | |
1014 | Bool_t AliFlowEvent::FirstLastPhiWgt() const { return !fOnePhiWgt ; } | |
1015 | Bool_t AliFlowEvent::OnePhiWgt() const { return fOnePhiWgt ; } | |
1016 | Bool_t AliFlowEvent::NoWgt() const { return fNoWgt; } | |
1017 | Bool_t AliFlowEvent::EtaSubs() const { return fEtaSubs ; } | |
1018 | //----------------------------------------------------------------------- | |
1019 | void AliFlowEvent::SetEtaSubs(Bool_t etasub) { fEtaSubs = etasub ; } | |
1020 | void AliFlowEvent::SetRunID(const Int_t& id) { fRunID = id; } | |
1021 | void AliFlowEvent::SetMagneticField(const Double_t& mf) { Flow::fMagneticField = mf; } | |
1022 | void AliFlowEvent::SetCenterOfMassEnergy(const Double_t& cms) { Flow::fCenterOfMassEnergy = cms; } | |
1023 | void AliFlowEvent::SetBeamMassNumberEast(const Short_t& bme) { Flow::fBeamMassNumberEast = bme; } | |
1024 | void AliFlowEvent::SetBeamMassNumberWest(const Short_t& bmw) { Flow::fBeamMassNumberWest = bmw; } | |
1025 | void AliFlowEvent::SetOrigMult(const UInt_t& tracks) { fOrigMult = tracks; } | |
1026 | void AliFlowEvent::SetL0TriggerWord(const Long_t& trigger) { fL0TriggerWord = trigger; } | |
1027 | void AliFlowEvent::SetZDCpart(Int_t zdcp) { fZDCpart = zdcp ; } | |
1028 | void AliFlowEvent::SetZDCenergy(Float_t n, Float_t p, Float_t em) { fZDCenergy[0] = n ; fZDCenergy[1] = p ; fZDCenergy[2] = em ; } | |
1029 | //----------------------------------------------------------------------- | |
1030 | void AliFlowEvent::SetBayesian(Double_t bayes[Flow::nPid]) | |
1031 | { | |
1032 | for(Int_t i=0;i<Flow::nPid;i++) { Flow::fBayesian[i] = bayes[i] ; } | |
1033 | } | |
1034 | //----------------------------------------------------------------------- | |
1035 | void AliFlowEvent::SetNoWgt(Bool_t nowgt) | |
1036 | { | |
1037 | // Sets no phi weightening, but Wgt = 1*sign(Eta) for odd harmonics . | |
1038 | ||
1039 | fNoWgt = nowgt ; // cout << " Wgt = +1 (positive Eta) or -1 (negative Eta) . " << endl ; | |
1040 | } | |
1041 | //----------------------------------------------------------------------- | |
1042 | void AliFlowEvent::SetOnePhiWgt() { fOnePhiWgt = kTRUE ; } | |
1043 | void AliFlowEvent::SetFirstLastPhiWgt() { fOnePhiWgt = kFALSE ; } | |
1044 | void AliFlowEvent::SetPtWgt(Bool_t PtWgt) { fPtWgt = PtWgt; } | |
1045 | void AliFlowEvent::SetEtaWgt(Bool_t EtaWgt) { fEtaWgt = EtaWgt; } | |
1046 | //----------------------------------------------------------------------- | |
1047 | #ifndef __CINT__ | |
1048 | void AliFlowEvent::SetPhiWeight(const Flow::PhiWgt_t& pPhiWgt) { memcpy (fPhiWgt, pPhiWgt, sizeof(Flow::PhiWgt_t)); } | |
1049 | void AliFlowEvent::SetPhiWeightPlus(const Flow::PhiWgt_t& pPhiWgtPlus) { memcpy (fPhiWgtPlus, pPhiWgtPlus, sizeof(Flow::PhiWgt_t)); } | |
1050 | void AliFlowEvent::SetPhiWeightMinus(const Flow::PhiWgt_t& pPhiWgtMinus) { memcpy (fPhiWgtMinus, pPhiWgtMinus, sizeof(Flow::PhiWgt_t)); } | |
1051 | void AliFlowEvent::SetPhiWeightCross(const Flow::PhiWgt_t& pPhiWgtCross) { memcpy (fPhiWgtCross, pPhiWgtCross, sizeof(Flow::PhiWgt_t)); } | |
1052 | #endif | |
1053 | //----------------------------------------------------------------------- | |
1054 | ////////////////////////////////////////////////////////////////////// |