AliPhysics  ff0b22e (ff0b22e)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliFlowEvent.cxx
Go to the documentation of this file.
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  AliFlowEvent: Event container for flow analysis
18 
19  origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
20  mods: Redmer A. Bertens (rbertens@cern.ch)
21 *****************************************************************/
22 
23 #include "Riostream.h"
24 #include "TFile.h"
25 #include "TList.h"
26 #include "TH1.h"
27 #include "TH1F.h"
28 #include "TH2F.h"
29 #include "TH3F.h"
30 #include "TArrayD.h"
31 #include "TProfile.h"
32 #include "AliMCEvent.h"
33 #include "AliMCParticle.h"
34 #include "AliCFManager.h"
35 #include "AliESDtrack.h"
36 #include "AliESDPmdTrack.h"
37 #include "AliESDEvent.h"
38 #include "AliAODEvent.h"
39 #include "AliOADBContainer.h"
40 #include "AliGenCocktailEventHeader.h"
41 #include "AliGenEposEventHeader.h"
42 #include "AliGenHijingEventHeader.h"
43 #include "AliGenGeVSimEventHeader.h"
44 #include "AliCollisionGeometry.h"
45 #include "AliMultiplicity.h"
46 #include "AliMultSelection.h"
47 #include "AliFlowTrackCuts.h"
48 #include "AliFlowEventSimple.h"
49 #include "AliFlowTrack.h"
50 #include "AliFlowVector.h"
51 #include "AliFlowEvent.h"
52 #include "AliLog.h"
53 
54 using std::endl;
55 using std::cout;
56 ClassImp(AliFlowEvent)
57 
58 //-----------------------------------------------------------------------
59 AliFlowEvent::AliFlowEvent():
60  AliFlowEventSimple(), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
61 {
62  // constructor
63  for(Int_t i(0); i < 9; i++) {
64  for(Int_t j(0); j < 2; j++) {
65  for(Int_t k(0); k < 2; k++) {
66  fMeanQ[i][j][k] = 0.;
67  fWidthQ[i][j][k] = 0.;
68  fMeanQv3[i][j][k] = 0.;
69  fWidthQv3[i][j][k] = 0.;
70  }
71  }
72  }
73  for(Int_t i(0); i < 5; i++) {
74  fQxavsV0[i] = 0x0;
75  fQyavsV0[i] = 0x0;
76  fQxcvsV0[i] = 0x0;
77  fQycvsV0[i] = 0x0;
78  }
79  //ctor
80  cout << "AliFlowEvent: Default constructor to be used only by root for io" << endl;
81 }
82 
83 //-----------------------------------------------------------------------
84 AliFlowEvent::AliFlowEvent(Int_t n):
85  AliFlowEventSimple(n), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
86 {
87  // constructor
88  for(Int_t i(0); i < 9; i++) {
89  for(Int_t j(0); j < 2; j++) {
90  for(Int_t k(0); k < 2; k++) {
91  fMeanQ[i][j][k] = 0.;
92  fWidthQ[i][j][k] = 0.;
93  fMeanQv3[i][j][k] = 0.;
94  fWidthQv3[i][j][k] = 0.;
95  }
96  }
97  }
98  for(Int_t i(0); i < 5; i++) {
99  fQxavsV0[i] = 0x0;
100  fQyavsV0[i] = 0x0;
101  fQxcvsV0[i] = 0x0;
102  fQycvsV0[i] = 0x0;
103  }
104 }
105 
106 //-----------------------------------------------------------------------
107 AliFlowEvent::AliFlowEvent(const AliFlowEvent& event):
108  AliFlowEventSimple(event), fApplyRecentering(event.fApplyRecentering), fDivSigma(event.fDivSigma), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
109 {
110  // copy constructor
111  for(Int_t i(0); i < 9; i++) {
112  for(Int_t j(0); j < 2; j++) {
113  for(Int_t k(0); k < 2; k++) {
114  fMeanQ[i][j][k] = 0.;
115  fWidthQ[i][j][k] = 0.;
116  fMeanQv3[i][j][k] = 0.;
117  fWidthQv3[i][j][k] = 0.;
118  }
119  }
120  }
121  for(Int_t i(0); i < 5; i++) {
122  fQxavsV0[i] = 0x0;
123  fQyavsV0[i] = 0x0;
124  fQxcvsV0[i] = 0x0;
125  fQycvsV0[i] = 0x0;
126  }
127 }
128 
129 //-----------------------------------------------------------------------
130 AliFlowEvent& AliFlowEvent::operator=(const AliFlowEvent& event)
131 {
132  //assignment operator
133  if (&event==this) return *this; // check self-assignment
134 
135  fApplyRecentering = event.fApplyRecentering;
136  fCachedRun = event.fCachedRun;
137  fVZEROcentralityBin = event.fVZEROcentralityBin;
138  fDivSigma = event.fDivSigma;
139  fEvent = 0x0; // should never be copied
140  fChi2A = 0x0; // do not clone these; if 0x0 they will be retrieved from the rp cuts object
141  fChi2C = 0x0;
142  fChi3A = 0x0;
143  fChi3C = 0x0;
144  for(Int_t i(0); i < 9; i++) {
145  for(Int_t j(0); j < 2; j++) {
146  for(Int_t k(0); k < 2; k++) {
147  fMeanQ[i][j][k] = event.fMeanQ[i][j][k];
148  fWidthQ[i][j][k] = event.fWidthQ[i][j][k];
149  fMeanQv3[i][j][k] = event.fMeanQv3[i][j][k];
150  fWidthQv3[i][j][k] = event.fWidthQv3[i][j][k];
151  }
152  }
153  }
154  for(Int_t i(0); i < 5; i++) {
155  fQxavsV0[i] = 0x0;
156  fQyavsV0[i] = 0x0;
157  fQxcvsV0[i] = 0x0;
158  fQycvsV0[i] = 0x0;
159  }
160 
161  AliFlowEventSimple::operator=(event);
162  return *this;
163 }
164 //-----------------------------------------------------------------------
165 AliFlowTrack* AliFlowEvent::GetTrack(Int_t i)
166 {
167  //get track i from collection
168  if (i>=fNumberOfTracks) return NULL;
169  AliFlowTrack* pTrack = static_cast<AliFlowTrack*>(fTrackCollection->At(i)) ;
170  return pTrack;
171 }
172 
173 //-----------------------------------------------------------------------
174 void AliFlowEvent::SetMCReactionPlaneAngle(const AliMCEvent* mcEvent)
175 {
176  //sets the event plane angle from the proper header in the MC
177 
178  //COCKTAIL with HIJING
179  if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Cocktail Header")) //returns 0 if matches
180  {
181  AliGenCocktailEventHeader *headerC = dynamic_cast<AliGenCocktailEventHeader *> (mcEvent-> GenEventHeader());
182  if (headerC)
183  {
184  TList *lhd = headerC->GetHeaders();
185  if (lhd)
186  {
187  AliGenHijingEventHeader *hdh = dynamic_cast<AliGenHijingEventHeader *> (lhd->At(0));
188  if (hdh) AliFlowEventSimple::SetMCReactionPlaneAngle( hdh->ReactionPlaneAngle() );
189  }
190  }
191  }
192  //THERMINATOR
193  else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Therminator")) //returns 0 if matches
194  {
195  AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(mcEvent->GenEventHeader());
196  if (headerH) AliFlowEventSimple::SetMCReactionPlaneAngle( headerH->ReactionPlaneAngle() );
197  }
198  //GEVSIM
199  else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"GeVSim header")) //returns 0 if matches
200  {
201  AliGenGeVSimEventHeader* headerG = dynamic_cast<AliGenGeVSimEventHeader*>(mcEvent->GenEventHeader());
202  if (headerG) AliFlowEventSimple::SetMCReactionPlaneAngle( headerG->GetEventPlane() );
203  }
204  //HIJING
205  else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Hijing")) //returns 0 if matches
206  {
207  AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(mcEvent->GenEventHeader());
208  if (headerH) AliFlowEventSimple::SetMCReactionPlaneAngle( headerH->ReactionPlaneAngle() );
209  }
210  //AMPT
211  else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Ampt")) //returns 0 if matches
212  {
213  AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(mcEvent->GenEventHeader());
214  if (headerH) AliFlowEventSimple::SetMCReactionPlaneAngle( headerH->ReactionPlaneAngle() );
215  }
216  //EPOS
217  else if (!strcmp(mcEvent->GenEventHeader()->GetName(),"EPOS"))
218  {
219  AliGenEposEventHeader* headerE = dynamic_cast<AliGenEposEventHeader*>(mcEvent->GenEventHeader());
220  if (headerE) AliFlowEventSimple::SetMCReactionPlaneAngle( headerE->ReactionPlaneAngle() );
221  }
222  //Hydjet
223  else
224  {
225  AliCollisionGeometry* header = dynamic_cast<AliCollisionGeometry*>(mcEvent->GenEventHeader());
226  if (header) AliFlowEventSimple::SetMCReactionPlaneAngle( header->ReactionPlaneAngle() );
227  }
228 }
229 
230 //-----------------------------------------------------------------------
231 AliFlowEvent::AliFlowEvent( const AliMCEvent* anInput,
232  const AliCFManager* rpCFManager,
233  const AliCFManager* poiCFManager):
234  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
235 {
236  // constructor
237  for(Int_t i(0); i < 9; i++) {
238  for(Int_t j(0); j < 2; j++) {
239  for(Int_t k(0); k < 2; k++) {
240  fMeanQ[i][j][k] = 0.;
241  fWidthQ[i][j][k] = 0.;
242  fMeanQv3[i][j][k] = 0.;
243  fWidthQv3[i][j][k] = 0.;
244  }
245  }
246  }
247  for(Int_t i(0); i < 5; i++) {
248  fQxavsV0[i] = 0x0;
249  fQyavsV0[i] = 0x0;
250  fQxcvsV0[i] = 0x0;
251  fQycvsV0[i] = 0x0;
252  }
253 
254  //Fills the event from the MC kinematic information
255  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
256 
257  //loop over tracks
258  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
259  {
260  //get input particle
261  AliMCParticle* pParticle = dynamic_cast<AliMCParticle*>(anInput->GetTrack(itrkN));
262  if (!pParticle) continue;
263 
264  //check if pParticle passes the cuts
265  Bool_t rpOK = kTRUE;
266  Bool_t poiOK = kTRUE;
267  if (rpCFManager && poiCFManager)
268  {
269  rpOK = rpCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pParticle);
270  poiOK = poiCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pParticle);
271  }
272  if (!(rpOK||poiOK)) continue;
273 
274  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
275  pTrack->SetSource(AliFlowTrack::kFromMC);
276 
277  if (rpOK && rpCFManager)
278  {
279  pTrack->SetForRPSelection(kTRUE);
280  IncrementNumberOfPOIs(0);
281  }
282  if (poiOK && poiCFManager)
283  {
284  pTrack->SetForPOISelection(kTRUE);
285  IncrementNumberOfPOIs(1);
286  }
287 
288  AddTrack(pTrack) ;
289  }//for all tracks
290  SetMCReactionPlaneAngle(anInput);
291 }
292 
293 //-----------------------------------------------------------------------
294 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
295  const AliCFManager* rpCFManager,
296  const AliCFManager* poiCFManager ):
297  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
298 {
299  // constructor
300  for(Int_t i(0); i < 9; i++) {
301  for(Int_t j(0); j < 2; j++) {
302  for(Int_t k(0); k < 2; k++) {
303  fMeanQ[i][j][k] = 0.;
304  fWidthQ[i][j][k] = 0.;
305  fMeanQv3[i][j][k] = 0.;
306  fWidthQv3[i][j][k] = 0.;
307  }
308  }
309  }
310  for(Int_t i(0); i < 5; i++) {
311  fQxavsV0[i] = 0x0;
312  fQyavsV0[i] = 0x0;
313  fQxcvsV0[i] = 0x0;
314  fQycvsV0[i] = 0x0;
315  }
316  //set run number
317  if(anInput->GetRunNumber()) fRun = anInput->GetRunNumber();
318 
319  //Fills the event from the ESD
320 
321  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
322 
323  //loop over tracks
324  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
325  {
326  AliESDtrack* pParticle = anInput->GetTrack(itrkN); //get input particle
327 
328  //check if pParticle passes the cuts
329  Bool_t rpOK = kTRUE;
330  Bool_t poiOK = kTRUE;
331  if (rpCFManager && poiCFManager)
332  {
333  rpOK = ( rpCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
334  rpCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
335  poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
336  poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
337  }
338  if (!(rpOK || poiOK)) continue;
339 
340  //make new AliFLowTrack
341  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
342  pTrack->SetSource(AliFlowTrack::kFromESD);
343 
344  //marking the particles used for int. flow:
345  if(rpOK && rpCFManager)
346  {
347  pTrack->SetForRPSelection(kTRUE);
348  IncrementNumberOfPOIs(0);
349  }
350  //marking the particles used for diff. flow:
351  if(poiOK && poiCFManager)
352  {
353  pTrack->SetForPOISelection(kTRUE);
354  IncrementNumberOfPOIs(1);
355  }
356 
357  AddTrack(pTrack);
358  }//end of while (itrkN < iNumberOfInputTracks)
359 }
360 
361 //-----------------------------------------------------------------------
362 AliFlowEvent::AliFlowEvent( const AliAODEvent* anInput,
363  const AliCFManager* rpCFManager,
364  const AliCFManager* poiCFManager):
365  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
366 {
367  // constructor
368  for(Int_t i(0); i < 9; i++) {
369  for(Int_t j(0); j < 2; j++) {
370  for(Int_t k(0); k < 2; k++) {
371  fMeanQ[i][j][k] = 0.;
372  fWidthQ[i][j][k] = 0.;
373  fMeanQv3[i][j][k] = 0.;
374  fWidthQv3[i][j][k] = 0.;
375  }
376  }
377  }
378  for(Int_t i(0); i < 5; i++) {
379  fQxavsV0[i] = 0x0;
380  fQyavsV0[i] = 0x0;
381  fQxcvsV0[i] = 0x0;
382  fQycvsV0[i] = 0x0;
383  }
384 
385  //set run number
386  if(anInput->GetRunNumber()) fRun = anInput->GetRunNumber();
387 
388  //Fills the event from the AOD
389  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
390 
391  //loop over tracks
392  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
393  {
394  AliAODTrack* pParticle = dynamic_cast<AliAODTrack*>(anInput->GetTrack(itrkN));
395  if(!pParticle) AliFatal("Not a standard AOD"); //get input particle
396 
397  //check if pParticle passes the cuts
398  Bool_t rpOK = kTRUE;
399  Bool_t poiOK = kTRUE;
400  if (rpCFManager && poiCFManager)
401  {
402  rpOK = ( rpCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
403  rpCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
404  poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
405  poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
406  }
407  if (!(rpOK || poiOK)) continue;
408 
409  //make new AliFlowTrack
410  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
411  pTrack->SetSource(AliFlowTrack::kFromAOD);
412 
413  if (rpOK /* && rpCFManager */ ) // to be fixed - with CF managers uncommented only empty events (NULL in header files)
414  {
415  pTrack->SetForRPSelection(kTRUE);
416  IncrementNumberOfPOIs(0);
417  }
418  if (poiOK /* && poiCFManager*/ )
419  {
420  pTrack->SetForPOISelection(kTRUE);
421  IncrementNumberOfPOIs(1);
422  }
423  AddTrack(pTrack);
424  }
425 
426  // if (iSelParticlesRP >= fMinMult && iSelParticlesRP <= fMaxMult)
427  // {
428  // if ( (++fCount % 100) == 0)
429  // {
430  // if (!fMCReactionPlaneAngle == 0) cout<<" MC Reaction Plane Angle = "<< fMCReactionPlaneAngle << endl;
431  // else cout<<" MC Reaction Plane Angle = unknown "<< endl;
432  // cout<<" iGoodTracks = "<<iGoodTracks<<endl;
433  // cout<<" # of RP selected tracks = "<<iSelParticlesRP<<endl;
434  // cout<<" # of POI selected tracks = "<<iSelParticlesPOI<<endl;
435  // cout << "# " << fCount << " events processed" << endl;
436  // }
437  // return pEvent;
438  // }
439  // else
440  // {
441  // cout<<"Not enough tracks in the FlowEventSimple"<<endl;
442  // return 0;
443  // }
444  //}
445  //else
446  //{
447  // cout<<"Event does not pass multiplicity cuts"<<endl;
448  // return 0;
449  //}
450 
451 }
452 
453 //-----------------------------------------------------------------------
454 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
455  const AliMCEvent* anInputMc,
456  KineSource anOption,
457  const AliCFManager* rpCFManager,
458  const AliCFManager* poiCFManager ):
459  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
460 {
461  // constructor
462  for(Int_t i(0); i < 9; i++) {
463  for(Int_t j(0); j < 2; j++) {
464  for(Int_t k(0); k < 2; k++) {
465  fMeanQ[i][j][k] = 0.;
466  fWidthQ[i][j][k] = 0.;
467  fMeanQv3[i][j][k] = 0.;
468  fWidthQv3[i][j][k] = 0.;
469  }
470  }
471  }
472  for(Int_t i(0); i < 5; i++) {
473  fQxavsV0[i] = 0x0;
474  fQyavsV0[i] = 0x0;
475  fQxcvsV0[i] = 0x0;
476  fQycvsV0[i] = 0x0;
477  }
478  //fills the event with tracks from the ESD and kinematics from the MC info via the track label
479  if (anOption==kNoKine)
480  {
481  AliFatal("WRONG OPTION IN AliFlowEventMaker::FillTracks(AliESDEvent* anInput, AliMCEvent* anInputMc, KineSource anOption)");
482  exit(1);
483  }
484 
485  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
486 
487  Int_t iNumberOfInputTracksMC = anInputMc->GetNumberOfTracks() ;
488  if (iNumberOfInputTracksMC==-1)
489  {
490  AliError("Skipping Event -- No MC information available for this event");
491  return;
492  }
493 
494  //loop over ESD tracks
495  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
496  {
497  AliESDtrack* pParticle = anInput->GetTrack(itrkN); //get input particle
498  //get Label
499  Int_t iLabel = pParticle->GetLabel();
500  //match to mc particle
501  AliMCParticle* pMcParticle = (AliMCParticle*) anInputMc->GetTrack(TMath::Abs(iLabel));
502 
503  //check
504  if (TMath::Abs(pParticle->GetLabel())!=pMcParticle->Label())
505  AliWarning(Form("pParticle->GetLabel()!=pMcParticle->Label(), %i, %i", pParticle->GetLabel(), pMcParticle->Label()));
506 
507  //check if pParticle passes the cuts
508  Bool_t rpOK = kFALSE;
509  Bool_t poiOK = kFALSE;
510  if (rpCFManager && poiCFManager)
511  {
512  if(anOption == kESDkine)
513  {
514  if (rpCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle,"mcGenCuts1") &&
515  rpCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle))
516  rpOK=kTRUE;
517  if (poiCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle,"mcGenCuts2") &&
518  poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle))
519  poiOK=kTRUE;
520  }
521  else if (anOption == kMCkine)
522  {
523  if (rpCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle))
524  rpOK=kTRUE;
525  if (poiCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle))
526  poiOK=kTRUE;
527  }
528  }
529 
530  if (!(rpOK || poiOK)) continue;
531 
532  //make new AliFlowTrack
533  AliFlowTrack* pTrack = NULL;
534  if(anOption == kESDkine) //take the PID from the MC & the kinematics from the ESD
535  {
536  pTrack = new AliFlowTrack(pParticle);
537  }
538  else if (anOption == kMCkine) //take the PID and kinematics from the MC
539  {
540  pTrack = new AliFlowTrack(pMcParticle);
541  }
542 
543  if (rpOK && rpCFManager)
544  {
545  IncrementNumberOfPOIs(0);
546  pTrack->SetForRPSelection();
547  }
548  if (poiOK && poiCFManager)
549  {
550  IncrementNumberOfPOIs(1);
551  pTrack->SetForPOISelection();
552  }
553 
554  AddTrack(pTrack);
555  }
556  SetMCReactionPlaneAngle(anInputMc);
557 }
558 
559 //-----------------------------------------------------------------------
560 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
561  const AliMultiplicity* anInputTracklets,
562  const AliCFManager* poiCFManager ):
563  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
564 {
565  // constructor
566  for(Int_t i(0); i < 9; i++) {
567  for(Int_t j(0); j < 2; j++) {
568  for(Int_t k(0); k < 2; k++) {
569  fMeanQ[i][j][k] = 0.;
570  fWidthQ[i][j][k] = 0.;
571  fMeanQv3[i][j][k] = 0.;
572  fWidthQv3[i][j][k] = 0.;
573  }
574  }
575  }
576  for(Int_t i(0); i < 5; i++) {
577  fQxavsV0[i] = 0x0;
578  fQyavsV0[i] = 0x0;
579  fQxcvsV0[i] = 0x0;
580  fQycvsV0[i] = 0x0;
581  }
582 
583  //Select the particles of interest from the ESD
584  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
585 
586  //loop over tracks
587  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
588  {
589  AliESDtrack* pParticle = anInput->GetTrack(itrkN); //get input particle
590 
591  //check if pParticle passes the cuts
592  Bool_t poiOK = kTRUE;
593  if (poiCFManager)
594  {
595  poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
596  poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
597  }
598  if (!poiOK) continue;
599 
600  //make new AliFLowTrack
601  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
602 
603  //marking the particles used for the particle of interest (POI) selection:
604  if(poiOK && poiCFManager)
605  {
606  IncrementNumberOfPOIs(1);
607  pTrack->SetForPOISelection(kTRUE);
608  pTrack->SetSource(AliFlowTrack::kFromESD);
609  }
610 
611  AddTrack(pTrack);
612  }//end of while (itrkN < iNumberOfInputTracks)
613 
614  //Select the reference particles from the SPD tracklets
615  anInputTracklets = anInput->GetMultiplicity();
616  Int_t multSPD = anInputTracklets->GetNumberOfTracklets();
617 
618  //loop over tracklets
619  for (Int_t itracklet=0; itracklet<multSPD; ++itracklet) {
620  Float_t thetaTr= anInputTracklets->GetTheta(itracklet);
621  Float_t phiTr= anInputTracklets->GetPhi(itracklet);
622  // calculate eta
623  Float_t etaTr = -TMath::Log(TMath::Tan(thetaTr/2.));
624 
625  //make new AliFLowTrackSimple
626  AliFlowTrack* pTrack = new AliFlowTrack();
627  pTrack->SetPt(0.0);
628  pTrack->SetEta(etaTr);
629  pTrack->SetPhi(phiTr);
630  //marking the particles used for the reference particle (RP) selection:
631  IncrementNumberOfPOIs(0);
632  pTrack->SetForRPSelection(kTRUE);
633  pTrack->SetSource(AliFlowTrack::kFromTracklet);
634 
635  //Add the track to the flowevent
636  AddTrack(pTrack);
637  }
638 
639 }
640 
641 //-----------------------------------------------------------------------
642 AliFlowEvent::AliFlowEvent( const AliESDEvent* esd,
643  const AliCFManager* poiCFManager,
644  Bool_t hybrid):
645  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
646 {
647  // constructor
648  for(Int_t i(0); i < 9; i++) {
649  for(Int_t j(0); j < 2; j++) {
650  for(Int_t k(0); k < 2; k++) {
651  fMeanQ[i][j][k] = 0.;
652  fWidthQ[i][j][k] = 0.;
653  fMeanQv3[i][j][k] = 0.;
654  fWidthQv3[i][j][k] = 0.;
655  }
656  }
657  }
658  for(Int_t i(0); i < 5; i++) {
659  fQxavsV0[i] = 0x0;
660  fQyavsV0[i] = 0x0;
661  fQxcvsV0[i] = 0x0;
662  fQycvsV0[i] = 0x0;
663  }
664  //Select the particles of interest from the ESD
665  Int_t iNumberOfInputTracks = esd->GetNumberOfTracks() ;
666 
667  //Double_t gPt = 0.0, gP = 0.0;
668  Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0}; //The impact parameters and their covariance.
669 // Double_t dca3D = 0.0; FIXME unused variable
670 
671  AliESDtrack trackTPC;
672 
673  //loop over tracks
674  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
675  {
676 
677  if (!esd->GetTrack(itrkN)) continue;
678 
679  Bool_t useTPC = kFALSE;
680 
681  AliESDtrack* pParticle = esd->GetTrack(itrkN); //get input particle
682 
683  //check if pParticle passes the cuts
684  Bool_t poiOK = kTRUE;
685 
686  if (poiCFManager)
687  {
688  poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
689  poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
690  }
691 
692  if (!(poiOK)) continue;
693 
694  AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)pParticle->GetTPCInnerParam();
695 
696  if (tpcTrack)
697  {
698 
699 // gPt = tpcTrack->Pt();
700 // gP = tpcTrack->P();
701 
702  useTPC = kTRUE;
703 
704  const AliESDVertex *vertexSPD = esd->GetPrimaryVertexSPD();
705  const AliESDVertex *vertexTPC = esd->GetPrimaryVertexTPC();
706 
707  AliExternalTrackParam copy(*tpcTrack);
708  if(hybrid)
709  copy.PropagateToDCA(vertexSPD,esd->GetMagneticField(),100.,dca,cov);
710  else
711  copy.PropagateToDCA(vertexTPC,esd->GetMagneticField(),100.,dca,cov);
712 
713 // dca3D = TMath::Sqrt(TMath::Power(dca[0],2)+TMath::Power(dca[1],2)); FIXME unused variable
714 
715  }
716 
717  //make new AliFLowTrack
718  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
719 
720  pTrack->SetSource(AliFlowTrack::kFromESD);
721 
722  //marking the particles used for diff. flow:
723  if(poiOK && poiCFManager)
724  {
725  pTrack->SetForPOISelection(kTRUE);
726  IncrementNumberOfPOIs(1);
727  }
728 
729  if(useTPC)
730  {
731  pTrack->SetForRPSelection(kTRUE);
732  IncrementNumberOfPOIs(0);
733  }
734 
735  AddTrack(pTrack);
736 
737  }//end of while (itrkN < iNumberOfInputTracks)
738 
739 }
740 
741 //-----------------------------------------------------------------------
742 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
743  const TH2F* anInputFMDhist,
744  const AliCFManager* poiCFManager ):
745  AliFlowEventSimple(20), fApplyRecentering(-1), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
746 {
747  // constructor
748  for(Int_t i(0); i < 9; i++) {
749  for(Int_t j(0); j < 2; j++) {
750  for(Int_t k(0); k < 2; k++) {
751  fMeanQ[i][j][k] = 0.;
752  fWidthQ[i][j][k] = 0.;
753  fMeanQv3[i][j][k] = 0.;
754  fWidthQv3[i][j][k] = 0.;
755  }
756  }
757  }
758  for(Int_t i(0); i < 5; i++) {
759  fQxavsV0[i] = 0x0;
760  fQyavsV0[i] = 0x0;
761  fQxcvsV0[i] = 0x0;
762  fQycvsV0[i] = 0x0;
763  }
764 
765  //Select the particles of interest from the ESD
766  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
767 
768  //loop over tracks
769  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
770  {
771  AliESDtrack* pParticle = anInput->GetTrack(itrkN); //get input particle
772 
773  //check if pParticle passes the cuts
774  Bool_t poiOK = kTRUE;
775  if (poiCFManager)
776  {
777  poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
778  poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
779  }
780  if (!poiOK) continue;
781 
782  //make new AliFLowTrack
783  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
784 
785  //marking the particles used for the particle of interest (POI) selection:
786  if(poiOK && poiCFManager)
787  {
788  IncrementNumberOfPOIs(1);
789  pTrack->SetForPOISelection(kTRUE);
790  pTrack->SetSource(AliFlowTrack::kFromESD);
791  }
792 
793  AddTrack(pTrack);
794  }//end of while (itrkN < iNumberOfInputTracks)
795 
796  //Select the reference particles from the FMD hits
797  //loop over FMD histogram
798  Int_t iBinsEta = anInputFMDhist->GetNbinsX();
799  Int_t iBinsPhi = anInputFMDhist->GetNbinsY();
800 
801  for (Int_t iEta = 1; iEta <= iBinsEta; iEta++){
802  Double_t etaFMD = anInputFMDhist->GetXaxis()->GetBinCenter(iEta);
803  for (Int_t iPhi = 1; iPhi <= iBinsPhi; iPhi++){
804  Double_t phiFMD = anInputFMDhist->GetYaxis()->GetBinCenter(iPhi);
805  Double_t weightFMD = anInputFMDhist->GetBinContent(iEta,iPhi);
806 
807  if (weightFMD > 0.0) { //do not add empty bins
808  //make new AliFLowTrackSimple
809  AliFlowTrack* pTrack = new AliFlowTrack();
810  pTrack->SetPt(0.0);
811  pTrack->SetEta(etaFMD);
812  pTrack->SetPhi(phiFMD);
813  pTrack->SetWeight(weightFMD);
814  //marking the particles used for the reference particle (RP) selection:
815  pTrack->TagRP();
816  IncrementNumberOfPOIs(0);
817  pTrack->SetSource(AliFlowTrack::kFromFMD);
818 
819  //Add the track to the flowevent
820  AddTrack(pTrack);
821 
822  }
823  }
824  }
825 
826 }
827 
828 //-----------------------------------------------------------------------
829 void AliFlowEvent::FindDaughters(Bool_t keepDaughtersInRPselection)
830 {
831  //each flow track holds it's esd track index as well as its daughters esd index.
832  //fill the array of daughters for every track with the pointers to flow tracks
833  //to associate the mothers with daughters directly
834  for (Int_t iTrack=0; iTrack<fMothersCollection->GetEntriesFast(); iTrack++)
835  {
836  AliFlowTrack* mother = static_cast<AliFlowTrack*>(fMothersCollection->At(iTrack));
837  if (!mother) continue;
838  if (mother->GetNDaughters()<1) continue;
839  for (Int_t iDaughterCandidate=0; iDaughterCandidate<fNumberOfTracks; iDaughterCandidate++)
840  {
841  AliFlowTrack* daughterCandidate = static_cast<AliFlowTrack*>(fTrackCollection->At(iDaughterCandidate));
842  Int_t esdIndexDaughterCandidate = daughterCandidate->GetID();
843  for (Int_t iDaughter=0; iDaughter<mother->GetNDaughters(); iDaughter++)
844  {
845  Int_t esdIndexDaughter = mother->GetIDDaughter(iDaughter);
846  if (esdIndexDaughter==esdIndexDaughterCandidate)
847  {
848  mother->SetDaughter(iDaughter,daughterCandidate);
849  daughterCandidate->SetForRPSelection(keepDaughtersInRPselection);
850  }
851  }
852  }
853  }
854 }
855 
856 //-----------------------------------------------------------------------
857 void AliFlowEvent::Fill( AliFlowTrackCuts* rpCuts,
858  AliFlowTrackCuts* poiCuts )
859 {
860  //Fills the event from a vevent: AliESDEvent,AliAODEvent,AliMCEvent
861  //the input data needs to be attached to the cuts
862  //we have two cases, if we're cutting the same collection of tracks
863  //(same param type) then we can have tracks that are both rp and poi
864  //in the other case we want to have two exclusive sets of rps and pois
865  //e.g. one tracklets, the other PMD or global - USER IS RESPOSIBLE
866  //FOR MAKING SURE THEY DONT OVERLAP OR ELSE THE SAME PARTICLE WILL BE
867  //TAKEN TWICE
868 
869  ClearFast();
870  if (!rpCuts || !poiCuts) return;
871  AliFlowTrackCuts::trackParameterType sourceRP = rpCuts->GetParamType();
872  AliFlowTrackCuts::trackParameterType sourcePOI = poiCuts->GetParamType();
873  AliFlowTrack* pTrack=NULL;
874 
875  //set run
876  if(rpCuts->GetRun()) fRun = rpCuts->GetRun();
877 
878  // if the source for rp's or poi's is the VZERO detector, get the calibration
879  // and set the calibration parameters
880  if (sourceRP == AliFlowTrackCuts::kBetaVZERO) {
881  SetBetaVZEROCalibrationForTrackCuts(rpCuts);
882  fDivSigma = rpCuts->GetDivSigma();
883  if(!rpCuts->GetApplyRecentering()) {
884  // if the user does not want to recenter, switch the flag
885  fApplyRecentering = -1;
886  }
887  // note: this flag is used in the overloaded implementation of Get2Qsub()
888  // and tells the function to use as Qsub vectors the recentered Q-vectors
889  // from the VZERO oadb file or from the event header
890  }
891  if (sourcePOI == AliFlowTrackCuts::kBetaVZERO) {
892  // probably no-one will choose vzero tracks as poi's ...
893  SetBetaVZEROCalibrationForTrackCuts(poiCuts);
894  fDivSigma = poiCuts->GetDivSigma();
895  }
896 
897  if (sourceRP == AliFlowTrackCuts::kDeltaVZERO) {
898  SetDeltaVZEROCalibrationForTrackCuts(rpCuts);
899  fDivSigma = rpCuts->GetDivSigma();
900  if(!rpCuts->GetApplyRecentering()) {
901  // if the user does not want to recenter, switch the flag
902  fApplyRecentering = -1;
903  }
904  }
905  if (sourcePOI == AliFlowTrackCuts::kDeltaVZERO) {
906  // probably no-one will choose vzero tracks as poi's ...
907  SetDeltaVZEROCalibrationForTrackCuts(poiCuts);
908  fDivSigma = poiCuts->GetDivSigma();
909  }
910 
911  if (sourceRP == AliFlowTrackCuts::kVZERO) {
912  SetVZEROCalibrationForTrackCuts(rpCuts);
913  if(!rpCuts->GetApplyRecentering()) {
914  // if the user does not want to recenter, switch the flag
915  fApplyRecentering = -1;
916  }
917  // note: this flag is used in the overloaded implementation of Get2Qsub()
918  // and tells the function to use as Qsub vectors the recentered Q-vectors
919  // from the VZERO oadb file or from the event header
920  }
921  if (sourcePOI == AliFlowTrackCuts::kVZERO) {
922  // probably no-one will choose vzero tracks as poi's ...
923  SetVZEROCalibrationForTrackCuts(poiCuts);
924  }
925 
926  if (sourceRP==sourcePOI)
927  {
928  //loop over tracks
929  Int_t numberOfInputObjects = rpCuts->GetNumberOfInputObjects();
930  for (Int_t i=0; i<numberOfInputObjects; i++)
931  {
932  //get input object (particle)
933  TObject* particle = rpCuts->GetInputObject(i);
934 
935  Bool_t rp = rpCuts->IsSelected(particle,i);
936  Bool_t poi = poiCuts->IsSelected(particle,i);
937 
938  if (!(rp||poi)) continue;
939 
940  //make new AliFlowTrack
941  if (rp)
942  {
943  pTrack = rpCuts->FillFlowTrack(fTrackCollection,fNumberOfTracks);
944  if (!pTrack) continue;
945  pTrack->Tag(0); IncrementNumberOfPOIs(0);
946  if (poi) {pTrack->Tag(1); IncrementNumberOfPOIs(1);}
947  if (pTrack->GetNDaughters()>0) fMothersCollection->Add(pTrack);
948  }
949  else if (poi)
950  {
951  pTrack = poiCuts->FillFlowTrack(fTrackCollection,fNumberOfTracks);
952  if (!pTrack) continue;
953  pTrack->Tag(1); IncrementNumberOfPOIs(1);
954  if (pTrack->GetNDaughters()>0) fMothersCollection->Add(pTrack);
955  }
956  fNumberOfTracks++;
957  }//end of while (i < numberOfTracks)
958  }
959  else if (sourceRP!=sourcePOI)
960  {
961  //here we have two different sources of particles, so we fill
962  //them independently
963  //POI
964  for (Int_t i=0; i<poiCuts->GetNumberOfInputObjects(); i++)
965  {
966  TObject* particle = poiCuts->GetInputObject(i);
967  Bool_t poi = poiCuts->IsSelected(particle,i);
968  if (!poi) continue;
969  pTrack = poiCuts->FillFlowTrack(fTrackCollection,fNumberOfTracks);
970  if (!pTrack) continue;
971  pTrack->Tag(1);
972  IncrementNumberOfPOIs(1);
973  fNumberOfTracks++;
974  if (pTrack->GetNDaughters()>0) fMothersCollection->Add(pTrack);
975  }
976  //RP
977  Int_t numberOfInputObjects = rpCuts->GetNumberOfInputObjects();
978  for (Int_t i=0; i<numberOfInputObjects; i++)
979  {
980  TObject* particle = rpCuts->GetInputObject(i);
981  Bool_t rp = rpCuts->IsSelected(particle,i);
982  if (!rp) continue;
983  pTrack = rpCuts->FillFlowTrack(fTrackCollection,fNumberOfTracks);
984  if (!pTrack) continue;
985  pTrack->Tag(0);
986  IncrementNumberOfPOIs(0);
987  fNumberOfTracks++;
988  if (pTrack->GetNDaughters()>0) fMothersCollection->Add(pTrack);
989  }
990  }
991 }
992 
993 //-----------------------------------------------------------------------
994 void AliFlowEvent::InsertTrack(AliFlowTrack *track) {
995  // adds a flow track at the end of the container
996  AliFlowTrack *pTrack = ReuseTrack( fNumberOfTracks++ );
997  *pTrack = *track;
998  if (track->GetNDaughters()>0)
999  {
1000  fMothersCollection->Add(pTrack);
1001  }
1002  return;
1003 }
1004 
1005 //-----------------------------------------------------------------------
1006 AliFlowTrack* AliFlowEvent::ReuseTrack(Int_t i)
1007 {
1008  //try to reuse an existing track, if empty, make new one
1009  AliFlowTrack* pTrack = static_cast<AliFlowTrack*>(fTrackCollection->At(i));
1010  if (pTrack)
1011  {
1012  pTrack->Clear();
1013  }
1014  else
1015  {
1016  pTrack = new AliFlowTrack();
1017  fTrackCollection->AddAtAndExpand(pTrack,i);
1018  }
1019  return pTrack;
1020 }
1021 
1022 //-----------------------------------------------------------------------
1023 AliFlowEvent::AliFlowEvent( AliFlowTrackCuts* rpCuts,
1024  AliFlowTrackCuts* poiCuts ):
1025  AliFlowEventSimple(20), fApplyRecentering(kFALSE), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
1026 
1027 {
1028  // constructor
1029  for(Int_t i(0); i < 9; i++) {
1030  for(Int_t j(0); j < 2; j++) {
1031  for(Int_t k(0); k < 2; k++) {
1032  fMeanQ[i][j][k] = 0.;
1033  fWidthQ[i][j][k] = 0.;
1034  fMeanQv3[i][j][k] = 0.;
1035  fWidthQv3[i][j][k] = 0.;
1036  }
1037  }
1038  }
1039  for(Int_t i(0); i < 5; i++) {
1040  fQxavsV0[i] = 0x0;
1041  fQyavsV0[i] = 0x0;
1042  fQxcvsV0[i] = 0x0;
1043  fQycvsV0[i] = 0x0;
1044  }
1045  //Fills the event from a vevent: AliESDEvent,AliAODEvent,AliMCEvent
1046  //the input data needs to be attached to the cuts
1047  Fill(rpCuts,poiCuts);
1048 }
1049 
1050 //-------------------------------------------------------------------//
1051 //---- Including PMD tracks as RP --------------------------//
1052 
1053 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
1054  const AliESDPmdTrack *pmdtracks,
1055  const AliCFManager* poiCFManager ):
1056  AliFlowEventSimple(20), fApplyRecentering(kFALSE), fDivSigma(kTRUE), fCachedRun(-1), fVZEROcentralityBin(-1), fEvent(0x0), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0)
1057 
1058 {
1059  // constructor
1060  for(Int_t i(0); i < 9; i++) {
1061  for(Int_t j(0); j < 2; j++) {
1062  for(Int_t k(0); k < 2; k++) {
1063  fMeanQ[i][j][k] = 0.;
1064  fWidthQ[i][j][k] = 0.;
1065  fMeanQv3[i][j][k] = 0.;
1066  fWidthQv3[i][j][k] = 0.;
1067  }
1068  }
1069  }
1070  for(Int_t i(0); i < 5; i++) {
1071  fQxavsV0[i] = 0x0;
1072  fQyavsV0[i] = 0x0;
1073  fQxcvsV0[i] = 0x0;
1074  fQycvsV0[i] = 0x0;
1075  }
1076 
1077  Float_t GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos);
1078  Float_t GetPmdPhi(Float_t xPos, Float_t yPos);
1079  //Select the particles of interest from the ESD
1080  Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
1081 
1082  //loop over tracks
1083  for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
1084  {
1085  AliESDtrack* pParticle = anInput->GetTrack(itrkN); //get input particle
1086  //check if pParticle passes the cuts
1087  Bool_t poiOK = kTRUE;
1088  if (poiCFManager)
1089  {
1090  poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
1091  poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
1092  }
1093  if (!poiOK) continue;
1094 
1095  //make new AliFLowTrack
1096  AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
1097 
1098  //marking the particles used for the particle of interest (POI) selection:
1099  if(poiOK && poiCFManager)
1100  {
1101  IncrementNumberOfPOIs(1);
1102  pTrack->SetForPOISelection(kTRUE);
1103  pTrack->SetSource(AliFlowTrack::kFromESD);
1104  }
1105 
1106  AddTrack(pTrack);
1107  }//end of while (itrkN < iNumberOfInputTracks)
1108 
1109  //Select the reference particles from the PMD tracks
1110  Int_t npmdcl = anInput->GetNumberOfPmdTracks();
1111  printf("======There are %d PMD tracks in this event\n-------",npmdcl);
1112  //loop over clusters
1113  for(Int_t iclust=0; iclust < npmdcl; iclust++){
1114  //AliESDPmdTrack *pmdtr = anInput->GetPmdTrack(iclust);
1115  pmdtracks = anInput->GetPmdTrack(iclust);
1116  Int_t det = pmdtracks->GetDetector();
1117  //Int_t smn = pmdtracks->GetSmn();
1118  Float_t clsX = pmdtracks->GetClusterX();
1119  Float_t clsY = pmdtracks->GetClusterY();
1120  Float_t clsZ = pmdtracks->GetClusterZ();
1121  Float_t ncell = pmdtracks->GetClusterCells();
1122  Float_t adc = pmdtracks->GetClusterADC();
1123  //Float_t pid = pmdtracks->GetClusterPID();
1124  Float_t etacls = GetPmdEta(clsX,clsY,clsZ);
1125  Float_t phicls = GetPmdPhi(clsX,clsY);
1126  //make new AliFLowTrackSimple
1127  AliFlowTrack* pTrack = new AliFlowTrack();
1128  //if(det == 0){ //selecting preshower plane only
1129  if(det == 0 && adc > 270 && ncell > 1){ //selecting preshower plane only
1130  //pTrack->SetPt(adc);//cluster adc
1131  pTrack->SetPt(0.0);
1132  pTrack->SetEta(etacls);
1133  pTrack->SetPhi(phicls);
1134  //marking the particles used for the reference particle (RP) selection:
1135  IncrementNumberOfPOIs(0);
1136  pTrack->SetForRPSelection(kTRUE);
1137  pTrack->SetSource(AliFlowTrack::kFromPMD);
1138  //Add the track to the flowevent
1139  AddTrack(pTrack);
1140  }//if det
1141  }
1142 }
1143 //----------------------------------------------------------------------------//
1144 Float_t GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos)
1145 {
1146  Float_t rpxpy, theta, eta;
1147  rpxpy = TMath::Sqrt(xPos*xPos + yPos*yPos);
1148  theta = TMath::ATan2(rpxpy,zPos);
1149  eta = -TMath::Log(TMath::Tan(0.5*theta));
1150  return eta;
1151 }
1152 //--------------------------------------------------------------------------//
1153 Float_t GetPmdPhi(Float_t xPos, Float_t yPos)
1154 {
1155  Float_t pybypx, phi = 0., phi1;
1156  if(xPos==0)
1157  {
1158  if(yPos>0) phi = 90.;
1159  if(yPos<0) phi = 270.;
1160  }
1161  if(xPos != 0)
1162  {
1163  pybypx = yPos/xPos;
1164  if(pybypx < 0) pybypx = - pybypx;
1165  phi1 = TMath::ATan(pybypx)*180./3.14159;
1166 
1167  if(xPos > 0 && yPos > 0) phi = phi1; // 1st Quadrant
1168  if(xPos < 0 && yPos > 0) phi = 180 - phi1; // 2nd Quadrant
1169  if(xPos < 0 && yPos < 0) phi = 180 + phi1; // 3rd Quadrant
1170  if(xPos > 0 && yPos < 0) phi = 360 - phi1; // 4th Quadrant
1171 
1172  }
1173  phi = phi*3.14159/180.;
1174  return phi;
1175 }
1176 //---------------------------------------------------------------//
1177 AliFlowVector AliFlowEvent::GetQ(Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
1178 {
1179  // start with the uncalibrated Q-vector. if we're not looking at vzero data, this will be returned
1180  // so this is completely backward compatible
1181  AliFlowVector vQ = AliFlowEventSimple::GetQ(n, weightsList, usePhiWeights, usePhiWeights, useEtaWeights);
1182  // see if we want to re-center 2010 style
1183  if(fApplyRecentering == 2010) {
1184  AliFlowVector Qarray[2];
1185  AliFlowEvent::Get2Qsub(Qarray, n, weightsList, usePhiWeights, usePtWeights, useEtaWeights);
1186  AliFlowVector vA = Qarray[0];
1187  AliFlowVector vB = Qarray[1];
1188 
1189  // now we have the calibrated sub-event q-vectors, these will be merged using a resolution derived weight
1190  Double_t chiA(1.), chiC(1.), dQX(0.), dQY(0.);
1191  if(n==2) {
1192  chiA = fChi2A->At(fVZEROcentralityBin);
1193  chiC = fChi2C->At(fVZEROcentralityBin);
1194  } else if (n==3) {
1195  chiA = fChi3A->At(fVZEROcentralityBin);
1196  chiC = fChi3C->At(fVZEROcentralityBin);
1197  }
1198 
1199  // combine the vzera and vzeroc signal, note that vzeroa is stored in vB and vzeroc in vA
1200  dQX = chiA*chiA*vB.X()+chiC*chiC*vA.X();
1201  dQY = chiA*chiA*vB.Y()+chiC*chiC*vA.Y();
1202  vQ.Set(dQX, dQY);
1203 
1204  } else if (fApplyRecentering == 2011) { // 11h style recentering
1205  Double_t dQX = 0.;
1206  Double_t dQY = 0.;
1207 
1208  if(fEvent && fEvent->GetEventplane()) {
1209  fEvent->GetEventplane()->CalculateVZEROEventPlane(fEvent, 10, n, dQX, dQY);
1210  // set the new q-vectors (which in this case means REPLACING)
1211  vQ.Set(dQX, dQY);
1212  } // if for some reason the info from the event header is not available, the AliFlowTrackCuts object
1213  // has provided the equalized multiplicity info so this should still be relatively safe
1214  }
1215 
1216  // return the Q-vector
1217  return vQ;
1218 }
1219  //---------------------------------------------------------------//
1220 void AliFlowEvent::Get2Qsub(AliFlowVector* Qarray, Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
1221 {
1222  // get q vectors for the subevents. if no recentering is necessary, get the guy from the flow event simple
1223  AliFlowEventSimple::Get2Qsub(Qarray, n, weightsList, usePhiWeights, usePtWeights, useEtaWeights);
1224  AliFlowVector vA = Qarray[0];
1225  AliFlowVector vB = Qarray[1];
1226 
1227  // else get the recentering from the cached info
1228  if (fApplyRecentering == 2010) // 10h style recentering, implemented for n=2 and n=3
1229  {
1230  // first retrieve the q-vectors from the AliFlowEventSimple:: routine
1231  // extract the information form the current flow vectors
1232  Double_t Qxc(vA.X()); // IMPORTANT: user is responsible for the sign of eta
1233  Double_t Qyc(vA.Y()); // vzeroC has negative pseudorapidity and is taken as subevent A
1234  Double_t Qxa(vB.X()); // vzeroA has positive pseudorapidity and is taken as subevent B
1235  Double_t Qya(vB.Y());
1236  // init some values for the corrections
1237 
1238  // default values for vector a (VZEROA)
1239  Double_t Qxamean(0);
1240  Double_t Qxarms(1);
1241  Double_t Qyamean(0);
1242  Double_t Qyarms(1);
1243  // default values for vector b (VZEROC)
1244  Double_t Qxcmean(0);
1245  Double_t Qxcrms(1);
1246  Double_t Qycmean(0);
1247  Double_t Qycrms(1);
1248  // note that defaults are chosen such that for n!=2||n!=3 re-centering is a null-operation
1249 
1250  if( n == 2) { // second order symmetry
1251  Qxamean = fMeanQ[fVZEROcentralityBin][1][0];
1252  Qxarms = fWidthQ[fVZEROcentralityBin][1][0];
1253  Qyamean = fMeanQ[fVZEROcentralityBin][1][1];
1254  Qyarms = fWidthQ[fVZEROcentralityBin][1][1];
1255 
1256  Qxcmean = fMeanQ[fVZEROcentralityBin][0][0];
1257  Qxcrms = fWidthQ[fVZEROcentralityBin][0][0];
1258  Qycmean = fMeanQ[fVZEROcentralityBin][0][1];
1259  Qycrms = fWidthQ[fVZEROcentralityBin][0][1];
1260  } else if (n == 3) { // third order symmetry
1261  Qxamean = fMeanQv3[fVZEROcentralityBin][1][0];
1262  Qxarms = fWidthQv3[fVZEROcentralityBin][1][0];
1263  Qyamean = fMeanQv3[fVZEROcentralityBin][1][1];
1264  Qyarms = fWidthQv3[fVZEROcentralityBin][1][1];
1265 
1266  Qxcmean = fMeanQv3[fVZEROcentralityBin][0][0];
1267  Qxcrms = fWidthQv3[fVZEROcentralityBin][0][0];
1268  Qycmean = fMeanQv3[fVZEROcentralityBin][0][1];
1269  Qycrms = fWidthQv3[fVZEROcentralityBin][0][1];
1270  }
1271  // do the correction
1272  Double_t QxaCor = (Qxa - Qxamean)/Qxarms;
1273  Double_t QyaCor = (Qya - Qyamean)/Qyarms;
1274  Double_t QxcCor = (Qxc - Qxcmean)/Qxcrms;
1275  Double_t QycCor = (Qyc - Qycmean)/Qycrms;
1276  // update the vector
1277  vA.Set(QxcCor, QycCor);
1278  vB.Set(QxaCor, QyaCor);
1279  } else if (fApplyRecentering == 2011) { // 11h style recentering
1280 
1281  Double_t QxaCor = 0.;
1282  Double_t QyaCor = 0.;
1283  Double_t QxcCor = 0.;
1284  Double_t QycCor = 0.;
1285 
1286  if(fEvent && fEvent->GetEventplane()) {
1287  fEvent->GetEventplane()->CalculateVZEROEventPlane(fEvent, 8, n, QxaCor, QyaCor);
1288  fEvent->GetEventplane()->CalculateVZEROEventPlane(fEvent, 9, n, QxcCor, QycCor);
1289 
1290  // set the new q-vectors (which in this case means REPLACING)
1291  vA.Set(QxcCor, QycCor);
1292  vB.Set(QxaCor, QyaCor);
1293  } // if for some reason the info from the event header is not available, the AliFlowTrackCuts object
1294  // has provided the equalized multiplicity info so this should still be relatively safe
1295  } else if (fApplyRecentering == 999) {
1296  // experimental VZERO recentering
1297  // first retrieve the q-vectors from the AliFlowEventSimple:: routine
1298  // extract the information form the current flow vectors
1299  Double_t Qxc(vA.X()); // IMPORTANT: user is responsible for the sign of eta
1300  Double_t Qyc(vA.Y()); // vzeroC has negative pseudorapidity and is taken as subevent A
1301  Double_t Qxa(vB.X()); // vzeroA has positive pseudorapidity and is taken as subevent B
1302  Double_t Qya(vB.Y());
1303  // init some values for the corrections
1304 
1305  Double_t Cen = fEvent->GetCentrality()->GetCentralityPercentile("V0M");
1306  // default values for vector a (VZEROA)
1307  Double_t Qxamean(fQxavsV0[n-1]->GetBinContent(fQxavsV0[n-1]->FindBin(Cen)));
1308  Double_t Qxarms(fQxavsV0[n-1]->GetBinError(fQxavsV0[n-1]->FindBin(Cen)));
1309  Double_t Qyamean(fQyavsV0[n-1]->GetBinContent(fQyavsV0[n-1]->FindBin(Cen)));
1310  Double_t Qyarms(fQyavsV0[n-1]->GetBinError(fQyavsV0[n-1]->FindBin(Cen)));
1311  // default values for vector b (VZEROC)
1312  Double_t Qxcmean(fQxcvsV0[n-1]->GetBinContent(fQxcvsV0[n-1]->FindBin(Cen)));
1313  Double_t Qxcrms(fQxcvsV0[n-1]->GetBinError(fQxcvsV0[n-1]->FindBin(Cen)));
1314  Double_t Qycmean(fQycvsV0[n-1]->GetBinContent(fQycvsV0[n-1]->FindBin(Cen)));
1315  Double_t Qycrms(fQycvsV0[n-1]->GetBinError(fQycvsV0[n-1]->FindBin(Cen)));
1316 
1317  // just a precaution ...
1318  if(n > 5) {
1319  Qxamean = 0;
1320  Qxarms = 1;
1321  Qyamean = 0;
1322  Qyarms = 1;
1323  Qxcmean = 0;
1324  Qxcrms = 1;
1325  Qycmean = 0; // this effectively disables the calibration
1326  Qycrms = 1;
1327  }
1328 
1329  Double_t QxaR = Qxa - Qxamean;
1330  Double_t QyaR = Qya - Qyamean;
1331  Double_t QxcR = Qxc - Qxcmean;
1332  Double_t QycR = Qyc - Qycmean;
1333  if(fDivSigma && Qxarms>0. && Qyarms>0. && Qxcrms>0. && Qycrms>0.) {
1334  QxaR /= Qxarms;
1335  QyaR /= Qyarms;
1336  QxcR /= Qxcrms;
1337  QycR /= Qycrms;
1338  }
1339  // update the vector
1340  vA.Set(QxcR, QycR);
1341  vB.Set(QxaR, QyaR);
1342 
1343  } else if (fApplyRecentering == 666) {
1344  // experimental VZERO recentering for 11h Full TPC Flow data, harmonics 1,2,3
1345  // first retrieve the q-vectors from the AliFlowEventSimple:: routine
1346  // extract the information form the current flow vectors
1347  if(n < 4) {
1348  Double_t Qxc(vA.X()); // IMPORTANT: user is responsible for the sign of eta
1349  Double_t Qyc(vA.Y()); // vzeroC has negative pseudorapidity and is taken as subevent A
1350  Double_t Qxa(vB.X()); // vzeroA has positive pseudorapidity and is taken as subevent B
1351  Double_t Qya(vB.Y());
1352  Double_t MultC = vA.GetMult();
1353  Double_t MultA = vB.GetMult();
1354  // init some values for the corrections
1355  Double_t Cen = - 999;
1356  //Added by Bernhard Hohlweger - bhohlweg@cern.ch
1357  if(fEvent->GetRunNumber() < 209122){
1358  //For Run1 Data the Old Centrality Percentile Method is available whereas for Run2 a new method was implemented
1359  //Cut was done for the first run of the LHC15a period
1360  Cen = fEvent->GetCentrality()->GetCentralityPercentile("V0M");
1361  }else{
1362  AliMultSelection *MultSelection = 0x0;
1363  MultSelection = (AliMultSelection * ) fEvent->FindListObject("MultSelection");
1364  if( !MultSelection) {
1365  //If you get this warning (and EventCentrality -999) please check that the AliMultSelectionTask actually ran (before your task)
1366  AliWarning("AliMultSelection not found, did you Run AliMultSelectionTask? \n");
1367  }else{
1368  Cen = MultSelection->GetMultiplicityPercentile("V0M");
1369  }
1370  }
1371  //14062016 Bernhard Hohlweger
1372  // default values for vector a (VZEROA)
1373  Double_t Qxamean(fQxavsV0[n-1]->GetBinContent(fQxavsV0[n-1]->FindBin(Cen)));
1374  Double_t Qxarms(fQxavsV0[n-1]->GetBinError(fQxavsV0[n-1]->FindBin(Cen)));
1375  Double_t Qyamean(fQyavsV0[n-1]->GetBinContent(fQyavsV0[n-1]->FindBin(Cen)));
1376  Double_t Qyarms(fQyavsV0[n-1]->GetBinError(fQyavsV0[n-1]->FindBin(Cen)));
1377  // default values for vector b (VZEROC)
1378  Double_t Qxcmean(fQxcvsV0[n-1]->GetBinContent(fQxcvsV0[n-1]->FindBin(Cen)));
1379  Double_t Qxcrms(fQxcvsV0[n-1]->GetBinError(fQxcvsV0[n-1]->FindBin(Cen)));
1380  Double_t Qycmean(fQycvsV0[n-1]->GetBinContent(fQycvsV0[n-1]->FindBin(Cen)));
1381  Double_t Qycrms(fQycvsV0[n-1]->GetBinError(fQycvsV0[n-1]->FindBin(Cen)));
1382 
1383  Double_t QxaR = Qxa - Qxamean*MultA;
1384  Double_t QyaR = Qya - Qyamean*MultA;
1385  Double_t QxcR = Qxc - Qxcmean*MultC;
1386  Double_t QycR = Qyc - Qycmean*MultC;
1387  if(fDivSigma && Qxarms>0. && Qyarms>0. && Qxcrms>0. && Qycrms>0.) {
1388  QxaR /= Qxarms;
1389  QyaR /= Qyarms;
1390  QxcR /= Qxcrms;
1391  QycR /= Qycrms;
1392  }
1393  // update the vector
1394  vA.Set(QxcR, QycR);
1395  vB.Set(QxaR, QyaR);
1396 
1397  } else {
1398  cout << " WARNING: recentering not possible for harmonic " << n << " (delta calibration) " << endl;
1399  } // end of if(n < 4)
1400  }
1401 
1402  Qarray[0] = vA;
1403  Qarray[1] = vB;
1404 }
1405 //_____________________________________________________________________________
1406 void AliFlowEvent::SetVZEROCalibrationForTrackCuts(AliFlowTrackCuts* cuts) {
1407  // open calibration info, copied from AliAnalyisTaskVnV0.cxx
1408  fEvent = cuts->GetEvent();
1409  if(!fEvent) return; // coverity. we need to know the event to get the runnumber and centrlaity
1410  // get the vzero centrality percentile (cc dependent calibration)
1411  Float_t v0Centr(fEvent->GetCentrality()->GetCentralityPercentile("V0M"));
1412  if(v0Centr < 5) fVZEROcentralityBin = 0;
1413  else if(v0Centr < 10) fVZEROcentralityBin = 1;
1414  else if(v0Centr < 20) fVZEROcentralityBin = 2;
1415  else if(v0Centr < 30) fVZEROcentralityBin = 3;
1416  else if(v0Centr < 40) fVZEROcentralityBin = 4;
1417  else if(v0Centr < 50) fVZEROcentralityBin = 5;
1418  else if(v0Centr < 60) fVZEROcentralityBin = 6;
1419  else if(v0Centr < 70) fVZEROcentralityBin = 7;
1420  else fVZEROcentralityBin = 8;
1421 
1422  // if this event is from the same run as the previous event
1423  // we can use the cached calibration values, no need to re-open the
1424  // aodb file, else cache the new run
1425  Int_t run(fEvent->GetRunNumber());
1426  if(fCachedRun == run) return;
1427  else fCachedRun = run;
1428 
1429  // relevant for 2010 data: check if the proper chi weights for merging vzero a and vzero c ep are present
1430  // if not, use sane defaults. centrality binning is equal to that given in the fVZEROcentralityBin snippet
1431  //
1432  // chi values can be calculated using the static helper function
1433  // AliAnalysisTaskJetV2::CalculateEventPlaneChi(Double_t res) where res is the event plane
1434  // resolution in a given centrality bin
1435  //
1436  // the resolutions that were used for these defaults are
1437  // Double_t R2VZEROA[] = {.35, .40, .48, .50, .48, .45, .38, .26, .16};
1438  // Double_t R2VZEROC[] = {.45, .60, .70, .73, .68, .60, .40, .36, .17};
1439  //
1440  // Double_t R3VZEROA[] = {.22, .23, .22, .19, .15, .12, .08, .00, .00};
1441  // Double_t R3VZEROC[] = {.30, .30, .28, .25, .22, .17, .11, .00, .00};
1442  // this might need a bit of updating as they were read 'by-eye' from a performance plot ..
1443  Double_t chiC2[] = {0.771423, 1.10236, 1.38116, 1.48077, 1.31964, 1.10236, 0.674622, 0.600403, 0.273865};
1444  Double_t chiA2[] = {0.582214, 0.674622, 0.832214, 0.873962, 0.832214, 0.771423, 0.637146, 0.424255, 0.257385};
1445  Double_t chiC3[] = {0.493347, 0.493347, 0.458557, 0.407166, 0.356628, 0.273865, 0.176208, 6.10352e-05, 6.10352e-05};
1446  Double_t chiA3[] = {0.356628, 0.373474, 0.356628, 0.306702, 0.24115, 0.192322, 0.127869, 6.10352e-05, 6.10352e-05};
1447 
1448  // this may seem redundant but in this way the cuts object is owner of the arrays
1449  // even if they're created here (so we won't get into trouble with dtor, assigmnet and copying)
1450  if(!cuts->GetChi2A()) cuts->SetChi2A(new TArrayD(9, chiA2));
1451  if(!cuts->GetChi2C()) cuts->SetChi2C(new TArrayD(9, chiC2));
1452  if(!cuts->GetChi3A()) cuts->SetChi3A(new TArrayD(9, chiA3));
1453  if(!cuts->GetChi3C()) cuts->SetChi3C(new TArrayD(9, chiC3));
1454 
1455  if(!fChi2A) fChi2A = cuts->GetChi2A();
1456  if(!fChi2C) fChi2C = cuts->GetChi2C();
1457  if(!fChi3A) fChi3A = cuts->GetChi3A();
1458  if(!fChi3C) fChi3C = cuts->GetChi3C();
1459 
1460  TFile *foadb = TFile::Open("$ALICE_PHYSICS/OADB/PWGCF/VZERO/VZEROcalibEP.root");
1461  if(!foadb){
1462  printf("OADB file $ALICE_PHYSICS/OADB/PWGCF/VZERO/VZEROcalibEP.root cannot be opened, CALIBRATION FAILED !");
1463  return;
1464  }
1465 
1466  AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr");
1467  if(!cont){
1468  printf("OADB object hMultV0BefCorr is not available in the file\n");
1469  return;
1470  }
1471  if(!(cont->GetObject(run))){
1472  // if the multiplicity correction cannot be found for the specified run,
1473  // loop over the 11h runs to see if it's 11h data
1474  Int_t runs11h[] = {170593, 170572, 170556, 170552, 170546, 170390, 170389, 170388, 170387, 170315, 170313, 170312, 170311, 170309, 170308, 170306, 170270, 170269, 170268, 170267, 170264, 170230, 170228, 170208, 170207, 170205, 170204, 170203, 170195, 170193, 170163, 170162, 170159, 170155, 170152, 170091, 170089, 170088, 170085, 170084, 170083, 170081, 170040, 170038, 170036, 170027, 169981, 169975, 169969, 169965, 169961, 169956, 169926, 169924, 169923, 169922, 169919, 169918, 169914, 169859, 169858, 169855, 169846, 169838, 169837, 169835, 169683, 169628, 169591, 169590, 169588, 169587, 169586, 169584, 169557, 169555, 169554, 169553, 169550, 169515, 169512, 169506, 169504, 169498, 169475, 169420, 169419, 169418, 169417, 169415, 169411, 169238, 169236, 169167, 169160, 169156, 169148, 169145, 169144, 169143, 169138, 169099, 169094, 169091, 169045, 169044, 169040, 169035, 168992, 168988, 168984, 168826, 168777, 168514, 168512, 168511, 168467, 168464, 168461, 168460, 168458, 168362, 168361, 168356, 168342, 168341, 168325, 168322, 168318, 168311, 168310, 168213, 168212, 168208, 168207, 168206, 168205, 168204, 168203, 168181, 168177, 168175, 168173, 168172, 168171, 168115, 168108, 168107, 168105, 168104, 168103, 168076, 168069, 168068, 168066, 167988, 167987, 167986, 167985, 167921, 167920, 167915, 167909, 167903, 167902, 167818, 167814, 167813, 167808, 167807, 167806, 167713, 167712, 167711, 167706, 167693};
1475  for(Int_t r(0); r < 176; r++) {
1476  if(run == runs11h[r]) {
1477  printf(" > run has been identified as 11h < \n");
1478  if(cuts->GetVZEROgainEqualizationPerRing()) {
1479  // enable or disable rings through the weights, weight 1. is enabled, 0. is disabled
1480  // start with the vzero c rings (segments 0 through 31)
1481  (cuts->GetUseVZERORing(0)) ? cuts->SetVZEROCpol(0, 1.) : cuts->SetVZEROCpol(0, 0.);
1482  (cuts->GetUseVZERORing(1)) ? cuts->SetVZEROCpol(1, 1.) : cuts->SetVZEROCpol(1, 0.);
1483  (cuts->GetUseVZERORing(2)) ? cuts->SetVZEROCpol(2, 1.) : cuts->SetVZEROCpol(2, 0.);
1484  (cuts->GetUseVZERORing(3)) ? cuts->SetVZEROCpol(3, 1.) : cuts->SetVZEROCpol(3, 0.);
1485  // same for vzero a
1486  (cuts->GetUseVZERORing(4)) ? cuts->SetVZEROApol(0, 1.) : cuts->SetVZEROApol(0, 0.);
1487  (cuts->GetUseVZERORing(5)) ? cuts->SetVZEROApol(1, 1.) : cuts->SetVZEROApol(1, 0.);
1488  (cuts->GetUseVZERORing(6)) ? cuts->SetVZEROApol(2, 1.) : cuts->SetVZEROApol(2, 0.);
1489  (cuts->GetUseVZERORing(7)) ? cuts->SetVZEROApol(3, 1.) : cuts->SetVZEROApol(3, 0.);
1490  } else {
1491  // else enable all rings, which is also default
1492  for(Int_t i(0); i < 4; i++) cuts->SetVZEROCpol(i, 1.);
1493  for(Int_t i(0); i < 4; i++) cuts->SetVZEROApol(i, 1.);
1494  }
1495  // pass a NULL pointer to the track cuts object, the NULL pointer will identify 11h runs
1496  cuts->SetVZEROgainEqualisation(NULL);
1497  fApplyRecentering = 2011;
1498  return; // the rest of the steps are not necessary
1499  }
1500  }
1501  // the run has not been identified as lhc11h data, so we assume a template calibration
1502  printf("OADB object hMultVZEROBefCorr is not available for run %i (used default run 137366)\n",run);
1503  run = 137366;
1504  }
1505  printf(" > run has been identified as 10h < \n");
1506  // step 1) get the proper multiplicity weights from the vzero signal
1507  TProfile* fMultVZERO = ((TH2F *) cont->GetObject(run))->ProfileX();
1508 
1509  TF1 *fpol0 = new TF1("fpol0","pol0");
1510  if(cuts->GetVZEROgainEqualizationPerRing()) {
1511  // do the calibration per ring
1512  // start with the vzero c rings (segments 0 through 31)
1513  fMultVZERO->Fit(fpol0, "N0", "", 0, 8);
1514  (cuts->GetUseVZERORing(0)) ? cuts->SetVZEROCpol(0, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(0, 0.);
1515  fMultVZERO->Fit(fpol0, "N0", "", 8, 16);
1516  (cuts->GetUseVZERORing(1)) ? cuts->SetVZEROCpol(1, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(1, 0.);
1517  fMultVZERO->Fit(fpol0, "N0", "", 16, 24);
1518  (cuts->GetUseVZERORing(2)) ? cuts->SetVZEROCpol(2, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(2, 0.);
1519  fMultVZERO->Fit(fpol0, "N0", "", 24, 32);
1520  (cuts->GetUseVZERORing(3)) ? cuts->SetVZEROCpol(3, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(3, 0.);
1521  // same thing for vero A
1522  fMultVZERO->Fit(fpol0, "N0", "", 32, 40);
1523  (cuts->GetUseVZERORing(4)) ? cuts->SetVZEROApol(0, fpol0->GetParameter(0)) : cuts->SetVZEROApol(0, 0.);
1524  fMultVZERO->Fit(fpol0, "N0", "", 40, 48);
1525  (cuts->GetUseVZERORing(5)) ? cuts->SetVZEROApol(1, fpol0->GetParameter(0)) : cuts->SetVZEROApol(1, 0.);
1526  fMultVZERO->Fit(fpol0, "N0", "", 48, 56);
1527  (cuts->GetUseVZERORing(6)) ? cuts->SetVZEROApol(2, fpol0->GetParameter(0)) : cuts->SetVZEROApol(2, 0.);
1528  fMultVZERO->Fit(fpol0, "N0", "", 56, 64);
1529  (cuts->GetUseVZERORing(7)) ? cuts->SetVZEROApol(3, fpol0->GetParameter(0)) : cuts->SetVZEROApol(3, 0.);
1530  } else {
1531  // do the calibration in one go. the calibration will still be
1532  // stored per ring, but each ring has the same weight now
1533  fMultVZERO->Fit(fpol0,"N0","",0,31);
1534  for(Int_t i(0); i < 4; i++) cuts->SetVZEROCpol(i, fpol0->GetParameter(0));
1535  fMultVZERO->Fit(fpol0,"N0","",32,64);
1536  for(Int_t i(0); i < 4; i++) cuts->SetVZEROApol(i, fpol0->GetParameter(0));
1537  }
1538  // the parameters to weigh the vzero track cuts have been extracted now,
1539  // so we can pass them to the current track cuts obect
1540  cuts->SetVZEROgainEqualisation(fMultVZERO); // passed as a TH1
1541 
1542  // step 2) reweight the q-vectors that will be called by flow methods which use
1543  // subevents
1544  // underlying assumption is that subevent a uses VZEROA
1545  // and subevent b uses VZEROC
1546  for(Int_t iside=0;iside<2;iside++){
1547  for(Int_t icoord=0;icoord<2;icoord++){
1548  for(Int_t i=0;i < 9;i++){
1549  char namecont[100];
1550  if(iside==0 && icoord==0)
1551  snprintf(namecont,100,"hQxc2_%i",i);
1552  else if(iside==1 && icoord==0)
1553  snprintf(namecont,100,"hQxa2_%i",i);
1554  else if(iside==0 && icoord==1)
1555  snprintf(namecont,100,"hQyc2_%i",i);
1556  else if(iside==1 && icoord==1)
1557  snprintf(namecont,100,"hQya2_%i",i);
1558 
1559  cont = (AliOADBContainer*) foadb->Get(namecont);
1560  if(!cont){
1561  printf("OADB object %s is not available in the file\n",namecont);
1562  return;
1563  }
1564 
1565  if(!(cont->GetObject(run))){
1566  printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
1567  run = 137366;
1568  }
1569 
1570  // after grabbing all the info, set the CORRECTION TERMS to
1571  // the 2nd and 3rd order qsub-vectors
1572  // we do this here for all centralities, so that subsequent events
1573  // can grab the correction from these cached values
1574  fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
1575  fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
1576 
1577  //for v3
1578  if(iside==0 && icoord==0)
1579  snprintf(namecont,100,"hQxc3_%i",i);
1580  else if(iside==1 && icoord==0)
1581  snprintf(namecont,100,"hQxa3_%i",i);
1582  else if(iside==0 && icoord==1)
1583  snprintf(namecont,100,"hQyc3_%i",i);
1584  else if(iside==1 && icoord==1)
1585  snprintf(namecont,100,"hQya3_%i",i);
1586 
1587  cont = (AliOADBContainer*) foadb->Get(namecont);
1588  if(!cont){
1589  printf("OADB object %s is not available in the file\n",namecont);
1590  return;
1591  }
1592 
1593  if(!(cont->GetObject(run))){
1594  printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
1595  run = 137366;
1596  }
1597  fMeanQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
1598  fWidthQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
1599 
1600  }
1601  }
1602  }
1603  // set the recentering style (might be switched back to -1 if recentering is disabeled)
1604  fApplyRecentering = 2010;
1605 }
1606 //-----------------------------------------------------------------------
1607 void AliFlowEvent::SetBetaVZEROCalibrationForTrackCuts(AliFlowTrackCuts* cuts) {
1608  // implementation of beta vzero calibration
1609 
1610  fEvent = cuts->GetEvent();
1611  if(!fEvent) return; // coverity. we need to know the event to get the runnumber and centrlaity
1612  // get the vzero centrality percentile (cc dependent calibration)
1613 
1614 
1615  // if this event is from the same run as the previous event
1616  // we can use the cached calibration values, no need to re-open the
1617  // aodb file, else cache the new run
1618  Int_t run(fEvent->GetRunNumber());
1619  if(fCachedRun == run) return;
1620  else fCachedRun = run;
1621 
1622  // check if the proper chi weights for merging vzero a and vzero c ep are present
1623  // if not, use sane defaults. centrality binning is equal to that given in the fVZEROcentralityBin snippet
1624  //
1625  // chi values can be calculated using the static helper function
1626  // AliAnalysisTaskJetV2::CalculateEventPlaneChi(Double_t res) where res is the event plane
1627  // resolution in a given centrality bin
1628  //
1629  // the resolutions that were used for these defaults are
1630  // Double_t R2VZEROA[] = {.35, .40, .48, .50, .48, .45, .38, .26, .16};
1631  // Double_t R2VZEROC[] = {.45, .60, .70, .73, .68, .60, .40, .36, .17};
1632  //
1633  // Double_t R3VZEROA[] = {.22, .23, .22, .19, .15, .12, .08, .00, .00};
1634  // Double_t R3VZEROC[] = {.30, .30, .28, .25, .22, .17, .11, .00, .00};
1635  Double_t chiC2[] = {0.771423, 1.10236, 1.38116, 1.48077, 1.31964, 1.10236, 0.674622, 0.600403, 0.273865};
1636  Double_t chiA2[] = {0.582214, 0.674622, 0.832214, 0.873962, 0.832214, 0.771423, 0.637146, 0.424255, 0.257385};
1637  Double_t chiC3[] = {0.493347, 0.493347, 0.458557, 0.407166, 0.356628, 0.273865, 0.176208, 6.10352e-05, 6.10352e-05};
1638  Double_t chiA3[] = {0.356628, 0.373474, 0.356628, 0.306702, 0.24115, 0.192322, 0.127869, 6.10352e-05, 6.10352e-05};
1639 
1640  // this may seem redundant but in this way the cuts object is owner of the arrays
1641  // even if they're created here (so we won't get into trouble with dtor, assigmnet and copying)
1642  if(!cuts->GetChi2A()) cuts->SetChi2A(new TArrayD(9, chiA2));
1643  if(!cuts->GetChi2C()) cuts->SetChi2C(new TArrayD(9, chiC2));
1644  if(!cuts->GetChi3A()) cuts->SetChi3A(new TArrayD(9, chiA3));
1645  if(!cuts->GetChi3C()) cuts->SetChi3C(new TArrayD(9, chiC3));
1646 
1647  if(!fChi2A) fChi2A = cuts->GetChi2A();
1648  if(!fChi2C) fChi2C = cuts->GetChi2C();
1649  if(!fChi3A) fChi3A = cuts->GetChi3A();
1650  if(!fChi3C) fChi3C = cuts->GetChi3C();
1651 
1652 
1653  TFile *foadb = TFile::Open("$ALICE_PHYSICS/PWGCF/FLOW/database/calibV0_filtered.root");
1654  if(!foadb){
1655  printf("OADB file $ALICE_PHYSICS/PWGCF/FLOW/database/calibV0_filtered.root cannot be opened, CALIBRATION FAILED !");
1656  return;
1657  }
1658 
1659  // first get the mutiplicity of the vzero channels before gain equalization
1660  AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr_filtered");
1661  if(!cont){
1662  printf("OADB object hMultV0BefCorr is not available in the file\n");
1663  return;
1664  }
1665  if(!(cont->GetObject(run))) {
1666  // if the multiplicity correction cannot be found for the specified run,
1667  // loop over the 11h runs to see if it's 11h data
1668  Int_t runs11h[] = {170593, 170572, 170556, 170552, 170546, 170390, 170389, 170388, 170387, 170315, 170313, 170312, 170311, 170309, 170308, 170306, 170270, 170269, 170268, 170267, 170264, 170230, 170228, 170208, 170207, 170205, 170204, 170203, 170195, 170193, 170163, 170162, 170159, 170155, 170152, 170091, 170089, 170088, 170085, 170084, 170083, 170081, 170040, 170038, 170036, 170027, 169981, 169975, 169969, 169965, 169961, 169956, 169926, 169924, 169923, 169922, 169919, 169918, 169914, 169859, 169858, 169855, 169846, 169838, 169837, 169835, 169683, 169628, 169591, 169590, 169588, 169587, 169586, 169584, 169557, 169555, 169554, 169553, 169550, 169515, 169512, 169506, 169504, 169498, 169475, 169420, 169419, 169418, 169417, 169415, 169411, 169238, 169236, 169167, 169160, 169156, 169148, 169145, 169144, 169143, 169138, 169099, 169094, 169091, 169045, 169044, 169040, 169035, 168992, 168988, 168984, 168826, 168777, 168514, 168512, 168511, 168467, 168464, 168461, 168460, 168458, 168362, 168361, 168356, 168342, 168341, 168325, 168322, 168318, 168311, 168310, 168213, 168212, 168208, 168207, 168206, 168205, 168204, 168203, 168181, 168177, 168175, 168173, 168172, 168171, 168115, 168108, 168107, 168105, 168104, 168103, 168076, 168069, 168068, 168066, 167988, 167987, 167986, 167985, 167921, 167920, 167915, 167909, 167903, 167902, 167818, 167814, 167813, 167808, 167807, 167806, 167713, 167712, 167711, 167706, 167693};
1669  for(Int_t r(0); r < 176; r++) {
1670  if(run == runs11h[r]) {
1671  printf(" > run has been identified as 11h < \n");
1672  if(cuts->GetVZEROgainEqualizationPerRing()) {
1673  // enable or disable rings through the weights, weight 1. is enabled, 0. is disabled
1674  // start with the vzero c rings (segments 0 through 31)
1675  (cuts->GetUseVZERORing(0)) ? cuts->SetVZEROCpol(0, 1.) : cuts->SetVZEROCpol(0, 0.);
1676  (cuts->GetUseVZERORing(1)) ? cuts->SetVZEROCpol(1, 1.) : cuts->SetVZEROCpol(1, 0.);
1677  (cuts->GetUseVZERORing(2)) ? cuts->SetVZEROCpol(2, 1.) : cuts->SetVZEROCpol(2, 0.);
1678  (cuts->GetUseVZERORing(3)) ? cuts->SetVZEROCpol(3, 1.) : cuts->SetVZEROCpol(3, 0.);
1679  // same for vzero a
1680  (cuts->GetUseVZERORing(4)) ? cuts->SetVZEROApol(0, 1.) : cuts->SetVZEROApol(0, 0.);
1681  (cuts->GetUseVZERORing(5)) ? cuts->SetVZEROApol(1, 1.) : cuts->SetVZEROApol(1, 0.);
1682  (cuts->GetUseVZERORing(6)) ? cuts->SetVZEROApol(2, 1.) : cuts->SetVZEROApol(2, 0.);
1683  (cuts->GetUseVZERORing(7)) ? cuts->SetVZEROApol(3, 1.) : cuts->SetVZEROApol(3, 0.);
1684  } else {
1685  // else enable all rings, which is also default
1686  for(Int_t i(0); i < 4; i++) cuts->SetVZEROCpol(i, 1.);
1687  for(Int_t i(0); i < 4; i++) cuts->SetVZEROApol(i, 1.);
1688  }
1689  // pass a NULL pointer to the track cuts object, the NULL pointer will identify 11h runs
1690  cuts->SetVZEROgainEqualisation(NULL);
1691  fApplyRecentering = 2011;
1692  return; // the rest of the steps are not necessary
1693  }
1694  }
1695  // the run has not been identified as lhc11h data, so we assume a template calibration
1696  printf("OADB object hMultVZEROBefCorr is not available for run %i (used default run 138275)\n",run);
1697  run = 138275;
1698  }
1699  printf(" > run has been identified as 10h < \n");
1700  // step 0) get the profile which contains average multiplicity per VZERO channel
1701  TProfile* fMultVZERO = static_cast<TProfile*>(cont->GetObject(run));
1702 
1703  TF1 *fpol0 = new TF1("fpol0","pol0");
1704  // step 1) extract the proper weights from the profile. Q-vector recentering relies on
1705  // ring-by-ring gain equalization, so the terms are extracted ring-by-ring here
1706 
1707  // start with the vzero c rings (segments 0 through 31)
1708  fMultVZERO->Fit(fpol0, "N0", "", 0, 8);
1709  (cuts->GetUseVZERORing(0)) ? cuts->SetVZEROCpol(0, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(0, 0.);
1710  fMultVZERO->Fit(fpol0, "N0", "", 8, 16);
1711  (cuts->GetUseVZERORing(1)) ? cuts->SetVZEROCpol(1, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(1, 0.);
1712  fMultVZERO->Fit(fpol0, "N0", "", 16, 24);
1713  (cuts->GetUseVZERORing(2)) ? cuts->SetVZEROCpol(2, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(2, 0.);
1714  fMultVZERO->Fit(fpol0, "N0", "", 24, 32);
1715  (cuts->GetUseVZERORing(3)) ? cuts->SetVZEROCpol(3, fpol0->GetParameter(0)) : cuts->SetVZEROCpol(3, 0.);
1716  // same thing for vero A
1717  fMultVZERO->Fit(fpol0, "N0", "", 32, 40);
1718  (cuts->GetUseVZERORing(4)) ? cuts->SetVZEROApol(0, fpol0->GetParameter(0)) : cuts->SetVZEROApol(0, 0.);
1719  fMultVZERO->Fit(fpol0, "N0", "", 40, 48);
1720  (cuts->GetUseVZERORing(5)) ? cuts->SetVZEROApol(1, fpol0->GetParameter(0)) : cuts->SetVZEROApol(1, 0.);
1721  fMultVZERO->Fit(fpol0, "N0", "", 48, 56);
1722  (cuts->GetUseVZERORing(6)) ? cuts->SetVZEROApol(2, fpol0->GetParameter(0)) : cuts->SetVZEROApol(2, 0.);
1723  fMultVZERO->Fit(fpol0, "N0", "", 56, 64);
1724  (cuts->GetUseVZERORing(7)) ? cuts->SetVZEROApol(3, fpol0->GetParameter(0)) : cuts->SetVZEROApol(3, 0.);
1725 
1726  // the parameters to weigh the vzero track cuts have been extracted now,
1727  // so we can pass them to the current track cuts obect
1728  cuts->SetVZEROgainEqualisation(fMultVZERO); // passed as a TH1
1729 
1730  // step 2) extract the calibration histograms from the database and
1731  // pass them to the cuts object
1732  //
1733  // first index of the oadb array is the harmonic n, the second index is either qax, qay, qcx, qcy
1734  AliOADBContainer* h[5][4];
1735  for(Int_t i(0); i < 5; i++) {
1736  h[i][0] = (AliOADBContainer*)foadb->Get(Form("hQxa%i_filtered", i+1));
1737  if(h[i][0]) fQxavsV0[i] = static_cast<TH1F*>(h[i][0]->GetObject(run));
1738  h[i][1] = (AliOADBContainer*)foadb->Get(Form("hQya%i_filtered", i+1));
1739  if(h[i][1]) fQyavsV0[i] = static_cast<TH1F*>(h[i][1]->GetObject(run));
1740  h[i][2] = (AliOADBContainer*)foadb->Get(Form("hQxc%i_filtered", i+1));
1741  if(h[i][2]) fQxcvsV0[i] = static_cast<TH1F*>(h[i][2]->GetObject(run));
1742  h[i][3] = (AliOADBContainer*)foadb->Get(Form("hQyc%i_filtered", i+1));
1743  if(h[i][3]) fQycvsV0[i] = static_cast<TH1F*>(h[i][3]->GetObject(run));
1744  }
1745 
1746  // set the recentering style (might be switched back to -1 if recentering is disabeled)
1747  // FIXME as an ugly hack, for now I mark this as 999 to denote the experimental nature
1748  // of this and use it transparently without disrupting the existing calbiration
1749  fApplyRecentering = 999;
1750 }
1751 
1752 //-----------------------------------------------------------------------------
1753 
1754 void AliFlowEvent::SetDeltaVZEROCalibrationForTrackCuts(AliFlowTrackCuts* cuts) {
1755  // implementation of delta vzero calibration 2011
1756 
1757  fEvent = cuts->GetEvent();
1758  if(!fEvent) return; // coverity. we need to know the event to get the runnumber and centrlaity
1759  // get the vzero centrality percentile (cc dependent calibration)
1760 
1761 
1762  // if this event is from the same run as the previous event
1763  // we can use the cached calibration values, no need to re-open the
1764  // aodb file, else cache the new run
1765  Int_t run(fEvent->GetRunNumber());
1766  if(fCachedRun == run) return;
1767  else fCachedRun = run;
1768 
1769  // check if the proper chi weights for merging vzero a and vzero c ep are present
1770  // if not, use sane defaults. centrality binning is equal to that given in the fVZEROcentralityBin snippet
1771  //
1772  // chi values can be calculated using the static helper function
1773  // AliAnalysisTaskJetV2::CalculateEventPlaneChi(Double_t res) where res is the event plane
1774  // resolution in a given centrality bin
1775  //
1776  // the resolutions that were used for these defaults are
1777  // Double_t R2VZEROA[] = {.35, .40, .48, .50, .48, .45, .38, .26, .16};
1778  // Double_t R2VZEROC[] = {.45, .60, .70, .73, .68, .60, .40, .36, .17};
1779  //
1780  // Double_t R3VZEROA[] = {.22, .23, .22, .19, .15, .12, .08, .00, .00};
1781  // Double_t R3VZEROC[] = {.30, .30, .28, .25, .22, .17, .11, .00, .00};
1782  Double_t chiC2[] = {0.771423, 1.10236, 1.38116, 1.48077, 1.31964, 1.10236, 0.674622, 0.600403, 0.273865};
1783  Double_t chiA2[] = {0.582214, 0.674622, 0.832214, 0.873962, 0.832214, 0.771423, 0.637146, 0.424255, 0.257385};
1784  Double_t chiC3[] = {0.493347, 0.493347, 0.458557, 0.407166, 0.356628, 0.273865, 0.176208, 6.10352e-05, 6.10352e-05};
1785  Double_t chiA3[] = {0.356628, 0.373474, 0.356628, 0.306702, 0.24115, 0.192322, 0.127869, 6.10352e-05, 6.10352e-05};
1786 
1787  // this may seem redundant but in this way the cuts object is owner of the arrays
1788  // even if they're created here (so we won't get into trouble with dtor, assigmnet and copying)
1789  if(!cuts->GetChi2A()) cuts->SetChi2A(new TArrayD(9, chiA2));
1790  if(!cuts->GetChi2C()) cuts->SetChi2C(new TArrayD(9, chiC2));
1791  if(!cuts->GetChi3A()) cuts->SetChi3A(new TArrayD(9, chiA3));
1792  if(!cuts->GetChi3C()) cuts->SetChi3C(new TArrayD(9, chiC3));
1793 
1794  if(!fChi2A) fChi2A = cuts->GetChi2A();
1795  if(!fChi2C) fChi2C = cuts->GetChi2C();
1796  if(!fChi3A) fChi3A = cuts->GetChi3A();
1797  if(!fChi3C) fChi3C = cuts->GetChi3C();
1798 
1799  // get the calibration file for the gain equalization
1800  TFile *foadb = TFile::Open("alien:///alice/cern.ch/user/j/jmargutt/gainVZERO.LHC11h.root");
1801  if(!foadb){
1802  printf("file alien:///alice/cern.ch/user/j/jmargutt/gainVZERO.LHC11h.root cannot be opened, CALIBRATION FAILED !");
1803  return;
1804  }
1805  TH3F* Weights = dynamic_cast<TH3F*>(foadb->FindObjectAny("LHC11h")->FindObject("gHistVZEROChannelGainEqualizationMap"));
1806  if(!Weights){
1807  printf("gHistVZEROChannelGainEqualizationMap is not available in the file\n");
1808  return;
1809  }
1810 
1811  if(cuts->GetVZEROgainEqualizationPerRing()) {
1812  // enable or disable rings through the weights, weight 1. is enabled, 0. is disabled
1813  // start with the vzero c rings (segments 0 through 31)
1814  (cuts->GetUseVZERORing(0)) ? cuts->SetVZEROCpol(0, 1.) : cuts->SetVZEROCpol(0, 0.);
1815  (cuts->GetUseVZERORing(1)) ? cuts->SetVZEROCpol(1, 1.) : cuts->SetVZEROCpol(1, 0.);
1816  (cuts->GetUseVZERORing(2)) ? cuts->SetVZEROCpol(2, 1.) : cuts->SetVZEROCpol(2, 0.);
1817  (cuts->GetUseVZERORing(3)) ? cuts->SetVZEROCpol(3, 1.) : cuts->SetVZEROCpol(3, 0.);
1818  // same for vzero a
1819  (cuts->GetUseVZERORing(4)) ? cuts->SetVZEROApol(0, 1.) : cuts->SetVZEROApol(0, 0.);
1820  (cuts->GetUseVZERORing(5)) ? cuts->SetVZEROApol(1, 1.) : cuts->SetVZEROApol(1, 0.);
1821  (cuts->GetUseVZERORing(6)) ? cuts->SetVZEROApol(2, 1.) : cuts->SetVZEROApol(2, 0.);
1822  (cuts->GetUseVZERORing(7)) ? cuts->SetVZEROApol(3, 1.) : cuts->SetVZEROApol(3, 0.);
1823  } else {
1824  // else enable all rings, which is also default
1825  for(Int_t i(0); i < 4; i++) cuts->SetVZEROCpol(i, 1.);
1826  for(Int_t i(0); i < 4; i++) cuts->SetVZEROApol(i, 1.);
1827  }
1828 
1829  // loop over the 11h runs to see if it's 11h data (Full TPC Flow)
1830  Int_t runs11h[] = {167915, 168115, 168460, 169035, 169238, 169859, 170228, 167920, 168310, 168464, 169091, 169411, 169923, 170230, 167985, 168311, 168467, 169094, 169415, 170027, 170268, 167987, 168322, 168511, 169138, 169417, 170081, 170269, 167988, 168325, 168512, 169144, 169835, 170155, 170270, 168069, 168341, 168514, 169145, 169837, 170159, 170306, 168076, 168342, 168777, 169148, 169838, 170163, 170308, 168105, 168361, 168826, 169156, 169846, 170193, 170309, 168107, 168362, 168988, 169160, 169855, 170203, 168108, 168458, 168992, 169167, 169858, 170204};
1831  Int_t RunBin = -1;
1832  for(Int_t r(0); r < 68; r++) {
1833  if(run == runs11h[r]) {
1834  RunBin = r+1;
1835  }
1836  }
1837  if(RunBin>0) {
1838  printf(" > run has been identified as 11h Full TPC Flow < \n");
1839  } else {
1840  printf(" > run has NOT been identified as 11h Full TPC Flow, use default for 11h < \n");
1841  // pass a NULL pointer to the track cuts object, the NULL pointer will identify 11h runs
1842  cuts->SetVZEROgainEqualisation(NULL);
1843  fApplyRecentering = 2011;
1844  return; // the rest of the steps are not necessary
1845  }
1846 
1847  // step 0) get the TH2D which contains the correction factor per VZERO channel per centrality bin
1848  // and pass it to the current track cuts obect
1849  Weights->GetXaxis()->SetRange(RunBin,RunBin);
1850  TH2D* fMultVZEROCen = dynamic_cast<TH2D*>(Weights->Project3D("zy"));
1851  cuts->SetVZEROgainEqualisationCen(fMultVZEROCen); // passed as a TH2
1852 
1853  // get the calibration file for the q-vector recentering
1854  TFile *fqvec = TFile::Open("alien:///alice/cern.ch/user/j/jmargutt/recenteringVZERO.LHC11h.root");
1855  if(!fqvec){
1856  printf("file alien:///alice/cern.ch/user/j/jmargutt/recenteringVZERO.LHC11h.root cannot be opened, CALIBRATION FAILED !");
1857  return;
1858  }
1859 
1860  // first index of the oadb array is the harmonic n, the second index is either qax, qay, qcx, qcy
1861  TH2D* h[3][4];
1862  for(Int_t i(0); i < 3; i++) {
1863  h[i][0] = (TH2D*)(fqvec->FindObjectAny("LHC11h")->FindObject(Form("gHistVZEROAQ%ixRecenteringMap",i+1)));
1864  if(h[i][0]) fQxavsV0[i] = static_cast<TH1D*>(h[i][0]->ProjectionY(Form("fQxavsV0[%i]",i),RunBin,RunBin));
1865  h[i][1] = (TH2D*)(fqvec->FindObjectAny("LHC11h")->FindObject(Form("gHistVZEROAQ%iyRecenteringMap",i+1)));
1866  if(h[i][1]) fQyavsV0[i] = static_cast<TH1D*>(h[i][1]->ProjectionY(Form("fQyavsV0[%i]",i),RunBin,RunBin));
1867  h[i][2] = (TH2D*)(fqvec->FindObjectAny("LHC11h")->FindObject(Form("gHistVZEROCQ%ixRecenteringMap",i+1)));
1868  if(h[i][2]) fQxcvsV0[i] = static_cast<TH1D*>(h[i][2]->ProjectionY(Form("fQxcvsV0[%i]",i),RunBin,RunBin));
1869  h[i][3] = (TH2D*)(fqvec->FindObjectAny("LHC11h")->FindObject(Form("gHistVZEROCQ%iyRecenteringMap",i+1)));
1870  if(h[i][3]) fQycvsV0[i] = static_cast<TH1D*>(h[i][3]->ProjectionY(Form("fQycvsV0[%i]",i),RunBin,RunBin));
1871  }
1872 
1873  // set the recentering style (might be switched back to -1 if recentering is disabeled)
1874  // FIXME as an ugly hack, for now I mark this as 666 to denote the experimental nature
1875  // of this and use it transparently without disrupting the existing calbiration
1876  fApplyRecentering = 666;
1877 }
1878 
1879 //-----------------------------------------------------------------------------
1880 
1881 void AliFlowEvent::ClearFast()
1882 {
1883  //clear the event without releasing any memory
1884  //note that cached run number of recentering settigns are not clear
1885  //(see AliFlowEvent::ClearCachedRun() )
1886  AliFlowEventSimple::ClearFast();
1887 }
1888 
1889 //_____________________________________________________________________________
1890 
1891 void AliFlowEvent::ClearCachedRun()
1892 {
1893  //clear the cached run (not in clear fast as cache needs to be persistent in most cases )
1894  fCachedRun=0;
1895  fApplyRecentering=0;
1896 }
AliFlowEvent::GetQ
virtual AliFlowVector GetQ(Int_t n=2, TList *weightsList=NULL, Bool_t usePhiWeights=kFALSE, Bool_t usePtWeights=kFALSE, Bool_t useEtaWeights=kFALSE)
Definition: AliFlowEvent.cxx:1177
AliFlowTrack.h
GetPmdEta
Float_t GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos)
Definition: AliFlowEvent.cxx:1144
AliFlowEvent.h
AliFlowEvent::fChi2A
TArrayD * fChi2A
current event
Definition: AliFlowEvent.h:117
AliFlowEvent::ClearCachedRun
virtual void ClearCachedRun()
Definition: AliFlowEvent.cxx:1891
AliFlowEvent::FindDaughters
void FindDaughters(Bool_t keepDaughtersInRPselection=kFALSE)
Definition: AliFlowEvent.cxx:829
AliFlowEvent::SetVZEROCalibrationForTrackCuts
void SetVZEROCalibrationForTrackCuts(AliFlowTrackCuts *cuts)
Definition: AliFlowEvent.cxx:1406
AliFlowTrackCuts::SetChi3C
void SetChi3C(TArrayD *Chi3C)
Definition: AliFlowTrackCuts.h:328
AliFlowTrackSimple::SetDaughter
virtual void SetDaughter(Int_t, AliFlowTrackSimple *)
Definition: AliFlowTrackSimple.h:131
AliFlowEventSimple::GetQ
virtual AliFlowVector GetQ(Int_t n=2, TList *weightsList=NULL, Bool_t usePhiWeights=kFALSE, Bool_t usePtWeights=kFALSE, Bool_t useEtaWeights=kFALSE)
Definition: AliFlowEventSimple.cxx:453
AliFlowTrackSimple::SetEta
void SetEta(Double_t eta)
Definition: AliFlowTrackSimple.h:174
AliFlowTrack::Clear
virtual void Clear(Option_t *o="")
Definition: AliFlowTrack.h:42
AliFlowTrackCuts::GetInputObject
TObject * GetInputObject(Int_t i)
Definition: AliFlowTrackCuts.cxx:2673
AliFlowTrackSimple::GetNDaughters
virtual Int_t GetNDaughters() const
Definition: AliFlowTrackSimple.h:128
AliFlowTrackCuts.h
AliFlowTrackCuts::SetChi3A
void SetChi3A(TArrayD *Chi3A)
Definition: AliFlowTrackCuts.h:326
AliFlowEvent::fChi3C
TArrayD * fChi3C
chi vs cent for vzero A ep_3
Definition: AliFlowEvent.h:120
AliFlowEvent::fChi2C
TArrayD * fChi2C
chi vs cent for vzero A ep_2
Definition: AliFlowEvent.h:118
AliFlowVector::GetMult
Double_t GetMult() const
Definition: AliFlowVector.h:46
AliFlowEvent::SetMCReactionPlaneAngle
void SetMCReactionPlaneAngle(const AliMCEvent *mcEvent)
Definition: AliFlowEvent.cxx:174
AliFlowEventSimple::operator=
AliFlowEventSimple & operator=(const AliFlowEventSimple &anEvent)
Definition: AliFlowEventSimple.cxx:228
AliFlowTrackCuts::GetRun
Int_t GetRun() const
Definition: AliFlowTrackCuts.h:339
AliFlowEvent::operator=
AliFlowEvent & operator=(const AliFlowEvent &event)
Definition: AliFlowEvent.cxx:130
AliFlowTrackCuts::GetChi2A
TArrayD * GetChi2A()
Definition: AliFlowTrackCuts.h:330
AliFlowTrackCuts::SetChi2A
void SetChi2A(TArrayD *Chi2A)
Definition: AliFlowTrackCuts.h:325
AliFlowTrackSimple::GetIDDaughter
virtual Int_t GetIDDaughter(Int_t) const
Definition: AliFlowTrackSimple.h:130
AliFlowEventSimple::AddTrack
void AddTrack(AliFlowTrackSimple *track)
Definition: AliFlowEventSimple.cxx:412
AliFlowTrackCuts::IsSelected
virtual Bool_t IsSelected(TObject *obj, Int_t id=-666)
Definition: AliFlowTrackCuts.cxx:730
AliFlowEvent::fEvent
AliVEvent * fEvent
recentering
Definition: AliFlowEvent.h:116
AliFlowTrackCuts::SetVZEROgainEqualisationCen
void SetVZEROgainEqualisationCen(TH2 *g)
Definition: AliFlowTrackCuts.h:307
ClassImp
ClassImp(AliFlowEvent) AliFlowEvent
Definition: AliFlowEvent.cxx:56
AliFlowTrackCuts::GetChi2C
TArrayD * GetChi2C()
Definition: AliFlowTrackCuts.h:332
AliFlowEventSimple.h
AliFlowTrackCuts::GetParamType
trackParameterType GetParamType() const
Definition: AliFlowTrackCuts.h:226
AliFlowTrackSimple::TagRP
void TagRP(Bool_t b=kTRUE)
Definition: AliFlowTrackSimple.h:53
AliFlowEvent::fDivSigma
Bool_t fDivSigma
Definition: AliFlowEvent.h:103
AliFlowEventSimple::IncrementNumberOfPOIs
void IncrementNumberOfPOIs(Int_t poiType=1)
Definition: AliFlowEventSimple.cxx:220
AliFlowEvent::GetTrack
AliFlowTrack * GetTrack(Int_t i)
Definition: AliFlowEvent.cxx:165
AliFlowTrackCuts::GetVZEROgainEqualizationPerRing
Bool_t GetVZEROgainEqualizationPerRing() const
Definition: AliFlowTrackCuts.h:314
AliFlowTrackSimple::SetForRPSelection
void SetForRPSelection(Bool_t b=kTRUE)
Definition: AliFlowTrackSimple.h:188
AliFlowTrackCuts::GetApplyRecentering
Bool_t GetApplyRecentering() const
Definition: AliFlowTrackCuts.h:312
AliFlowTrackCuts::GetNumberOfInputObjects
Int_t GetNumberOfInputObjects() const
Definition: AliFlowTrackCuts.cxx:2623
AliFlowEvent::SetDeltaVZEROCalibrationForTrackCuts
void SetDeltaVZEROCalibrationForTrackCuts(AliFlowTrackCuts *cuts)
Definition: AliFlowEvent.cxx:1754
AliFlowEvent::Fill
void Fill(AliFlowTrackCuts *rpCuts, AliFlowTrackCuts *poiCuts)
Definition: AliFlowEvent.cxx:857
AliFlowTrackCuts::FillFlowTrack
Bool_t FillFlowTrack(AliFlowTrack *track) const
Definition: AliFlowTrackCuts.cxx:2131
AliFlowEvent::fQxcvsV0
TH1 * fQxcvsV0[5]
recentering
Definition: AliFlowEvent.h:113
AliFlowTrackCuts::SetVZEROCpol
void SetVZEROCpol(Int_t ring, Float_t f)
Definition: AliFlowTrackCuts.h:309
AliFlowTrackCuts::GetChi3A
TArrayD * GetChi3A()
Definition: AliFlowTrackCuts.h:331
AliFlowTrackCuts::GetDivSigma
Bool_t GetDivSigma() const
Definition: AliFlowTrackCuts.h:316
AliFlowEvent::fWidthQ
Float_t fWidthQ[9][2][2]
recentering
Definition: AliFlowEvent.h:107
AliFlowTrackCuts::SetVZEROgainEqualisation
void SetVZEROgainEqualisation(TH1 *g)
Definition: AliFlowTrackCuts.h:306
GetPmdPhi
Float_t GetPmdPhi(Float_t xPos, Float_t yPos)
Definition: AliFlowEvent.cxx:1153
AliFlowTrackSimple::Tag
void Tag(Int_t n, Bool_t b=kTRUE)
Definition: AliFlowTrackSimple.h:55
AliFlowTrackCuts::GetEvent
AliVEvent * GetEvent() const
Definition: AliFlowTrackCuts.h:245
AliFlowEventSimple::fTrackCollection
TObjArray * fTrackCollection
Definition: AliFlowEventSimple.h:151
AliFlowTrack::SetSource
void SetSource(trackSource s)
Definition: AliFlowTrack.h:37
AliFlowEvent::fMeanQv3
Float_t fMeanQv3[9][2][2]
recentering
Definition: AliFlowEvent.h:108
AliFlowEventSimple::Get2Qsub
virtual void Get2Qsub(AliFlowVector *Qarray, Int_t n=2, TList *weightsList=NULL, Bool_t usePhiWeights=kFALSE, Bool_t usePtWeights=kFALSE, Bool_t useEtaWeights=kFALSE)
Definition: AliFlowEventSimple.cxx:570
AliFlowTrackSimple::SetPhi
void SetPhi(Double_t phi)
Definition: AliFlowTrackSimple.h:178
AliFlowEvent::fVZEROcentralityBin
Int_t fVZEROcentralityBin
cached calibration info for vzero
Definition: AliFlowEvent.h:105
AliFlowTrackSimple::SetPt
void SetPt(Double_t pt)
Definition: AliFlowTrackSimple.h:176
AliFlowTrackCuts::GetUseVZERORing
Bool_t GetUseVZERORing(Int_t i) const
Definition: AliFlowTrackCuts.h:324
AliFlowTrackCuts::GetChi3C
TArrayD * GetChi3C()
Definition: AliFlowTrackCuts.h:333
AliFlowTrackSimple::SetWeight
void SetWeight(Double_t weight)
Definition: AliFlowTrackSimple.h:180
AliFlowEvent::SetBetaVZEROCalibrationForTrackCuts
void SetBetaVZEROCalibrationForTrackCuts(AliFlowTrackCuts *cuts)
Definition: AliFlowEvent.cxx:1607
AliFlowTrackSimple::SetForPOISelection
void SetForPOISelection(Bool_t b=kTRUE)
Definition: AliFlowTrackSimple.h:190
AliFlowEventSimple::fMothersCollection
TObjArray * fMothersCollection
Definition: AliFlowEventSimple.h:174
AliFlowEventSimple::SetMCReactionPlaneAngle
void SetMCReactionPlaneAngle(Double_t fPhiRP)
Definition: AliFlowEventSimple.h:81
AliFlowEvent::fMeanQ
Float_t fMeanQ[9][2][2]
centrality bin for the current event
Definition: AliFlowEvent.h:106
AliFlowVector.h
AliFlowEvent::fQycvsV0
TH1 * fQycvsV0[5]
recentering
Definition: AliFlowEvent.h:114
AliFlowEvent::ReuseTrack
AliFlowTrack * ReuseTrack(Int_t i)
Definition: AliFlowEvent.cxx:1006
AliFlowEvent::Get2Qsub
virtual void Get2Qsub(AliFlowVector *Qarray, Int_t n=2, TList *weightsList=0x0, Bool_t usePhiWeights=0x0, Bool_t usePtWeights=0x0, Bool_t useEtaWeights=0x0)
Definition: AliFlowEvent.cxx:1220
AliFlowEvent::fQyavsV0
TH1 * fQyavsV0[5]
recentering
Definition: AliFlowEvent.h:112
AliFlowEvent::fCachedRun
Int_t fCachedRun
Definition: AliFlowEvent.h:104
AliFlowEvent::InsertTrack
void InsertTrack(AliFlowTrack *)
Definition: AliFlowEvent.cxx:994
AliFlowEvent::ClearFast
virtual void ClearFast()
Definition: AliFlowEvent.cxx:1881
AliFlowTrackCuts::SetChi2C
void SetChi2C(TArrayD *Chi2C)
Definition: AliFlowTrackCuts.h:327
AliFlowEvent::fQxavsV0
TH1 * fQxavsV0[5]
recentering
Definition: AliFlowEvent.h:111
AliFlowEvent::fChi3A
TArrayD * fChi3A
chi vs cent for vzero C ep_2
Definition: AliFlowEvent.h:119
AliFlowEvent::fWidthQv3
Float_t fWidthQv3[9][2][2]
recentering
Definition: AliFlowEvent.h:109
AliFlowTrackCuts::SetVZEROApol
void SetVZEROApol(Int_t ring, Float_t f)
Definition: AliFlowTrackCuts.h:308
AliFlowEvent::fApplyRecentering
Int_t fApplyRecentering
Definition: AliFlowEvent.h:102