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