PWG2/FLOW/AliFlowTasks/AliFlowEvent.cxx

   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 *****************************************************************/
  21
  22 #include "Riostream.h"
  23 #include "TList.h"
  24 #include "TH2F.h"
  25 #include "AliMCEvent.h"
  26 #include "AliMCParticle.h"
  27 #include "AliCFManager.h"
  28 #include "AliESDtrack.h"
  29 #include "AliESDPmdTrack.h"
  30 #include "AliESDEvent.h"
  31 #include "AliAODEvent.h"
  32 #include "AliGenCocktailEventHeader.h"
  33 #include "AliGenEposEventHeader.h"
  34 #include "AliGenHijingEventHeader.h"
  35 #include "AliGenGeVSimEventHeader.h"
  36 #include "AliCollisionGeometry.h"
  37 #include "AliMultiplicity.h"
  38 #include "AliFlowTrackCuts.h"
  39 #include "AliFlowEventSimple.h"
  40 #include "AliFlowTrack.h"
  41 #include "AliFlowEvent.h"
  42 #include "AliLog.h"
  43
  44 ClassImp(AliFlowEvent)
  45
  46 //-----------------------------------------------------------------------
  47
  48 AliFlowEvent::AliFlowEvent():
  49   AliFlowEventSimple()
  50 {
  51   //ctor
  52   cout << "AliFlowEvent: Default constructor to be used only by root for io" << endl;
  53 }
  54
  55 //-----------------------------------------------------------------------
  56 AliFlowEvent::AliFlowEvent(const AliFlowEvent& event):
  57   AliFlowEventSimple(event)
  58 {
  59   //cpy ctor
  60 }
  61
  62 //-----------------------------------------------------------------------
  63 AliFlowEvent& AliFlowEvent::operator=(const AliFlowEvent& event)
  64 {
  65   //assignment operator
  66   AliFlowEventSimple::operator=(event);
  67   return *this;
  68 }
  69
  70 //-----------------------------------------------------------------------
  71 AliFlowTrack* AliFlowEvent::GetTrack(Int_t i)
  72 {
  73   //get track i from collection
  74   if (i>=fNumberOfTracks) return NULL;
  75   AliFlowTrack* pTrack = static_cast<AliFlowTrack*>(fTrackCollection->At(i)) ;
  76   return pTrack;
  77 }
  78
  79 //-----------------------------------------------------------------------
  80 void AliFlowEvent::SetMCReactionPlaneAngle(const AliMCEvent* mcEvent)
  81 {
  82   //sets the event plane angle from the proper header in the MC
  83
  84   //COCKTAIL with HIJING
  85   if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Cocktail Header"))   //returns 0 if matches
  86   {
  87     AliGenCocktailEventHeader *headerC = dynamic_cast<AliGenCocktailEventHeader *> (mcEvent-> GenEventHeader());
  88     if (headerC)
  89     {
  90       TList *lhd = headerC->GetHeaders();
  91       if (lhd)
  92       {
  93         AliGenHijingEventHeader *hdh = dynamic_cast<AliGenHijingEventHeader *> (lhd->At(0));
  94         if (hdh) AliFlowEventSimple::SetMCReactionPlaneAngle( hdh->ReactionPlaneAngle() );
  95       }
  96     }
  97   }
  98   //THERMINATOR
  99   else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Therminator"))   //returns 0 if matches
 100   {
 101     AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(mcEvent->GenEventHeader());
 102     if (headerH) AliFlowEventSimple::SetMCReactionPlaneAngle( headerH->ReactionPlaneAngle() );
 103   }
 104   //GEVSIM
 105   else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"GeVSim header"))   //returns 0 if matches
 106   {
 107     AliGenGeVSimEventHeader* headerG = dynamic_cast<AliGenGeVSimEventHeader*>(mcEvent->GenEventHeader());
 108     if (headerG) AliFlowEventSimple::SetMCReactionPlaneAngle( headerG->GetEventPlane() );
 109   }
 110   //HIJING
 111   else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Hijing"))   //returns 0 if matches
 112   {
 113     AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(mcEvent->GenEventHeader());
 114     if (headerH) AliFlowEventSimple::SetMCReactionPlaneAngle( headerH->ReactionPlaneAngle() );
 115   }
 116   //EPOS
 117   else if (!strcmp(mcEvent->GenEventHeader()->GetName(),"EPOS"))
 118   {
 119     AliGenEposEventHeader* headerE = dynamic_cast<AliGenEposEventHeader*>(mcEvent->GenEventHeader());
 120     if (headerE) AliFlowEventSimple::SetMCReactionPlaneAngle( headerE->ReactionPlaneAngle() );
 121   }
 122   //Hydjet
 123   else
 124   {
 125     AliCollisionGeometry* header = dynamic_cast<AliCollisionGeometry*>(mcEvent->GenEventHeader());
 126     if (header) AliFlowEventSimple::SetMCReactionPlaneAngle( header->ReactionPlaneAngle() );
 127   }
 128 }
 129
 130 //-----------------------------------------------------------------------
 131 AliFlowEvent::AliFlowEvent( const AliMCEvent* anInput,
 132                             const AliCFManager* rpCFManager,
 133                             const AliCFManager* poiCFManager):
 134   AliFlowEventSimple(20)
 135 {
 136   //Fills the event from the MC kinematic information
 137
 138   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 139
 140   //loop over tracks
 141   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 142   {
 143     //get input particle
 144     AliMCParticle* pParticle = dynamic_cast<AliMCParticle*>(anInput->GetTrack(itrkN));
 145     if (!pParticle) continue;
 146
 147     //check if pParticle passes the cuts
 148     Bool_t rpOK = kTRUE;
 149     Bool_t poiOK = kTRUE;
 150     if (rpCFManager && poiCFManager)
 151     {
 152       rpOK = rpCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pParticle);
 153       poiOK = poiCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pParticle);
 154     }
 155     if (!(rpOK||poiOK)) continue;
 156
 157     AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 158     pTrack->SetSource(AliFlowTrack::kFromMC);
 159
 160     if (rpOK && rpCFManager)
 161     {
 162       pTrack->SetForRPSelection(kTRUE);
 163       fNumberOfRPs++;
 164     }
 165     if (poiOK && poiCFManager)
 166     {
 167       pTrack->SetForPOISelection(kTRUE);
 168     }
 169
 170     AddTrack(pTrack) ;
 171   }//for all tracks
 172   SetMCReactionPlaneAngle(anInput);
 173 }
 174
 175 //-----------------------------------------------------------------------
 176 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
 177                             const AliCFManager* rpCFManager,
 178                             const AliCFManager* poiCFManager ):
 179   AliFlowEventSimple(20)
 180 {
 181   //Fills the event from the ESD
 182
 183   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 184
 185   //loop over tracks
 186   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 187   {
 188     AliESDtrack* pParticle = anInput->GetTrack(itrkN);   //get input particle
 189
 190     //check if pParticle passes the cuts
 191     Bool_t rpOK = kTRUE;
 192     Bool_t poiOK = kTRUE;
 193     if (rpCFManager && poiCFManager)
 194     {
 195       rpOK = ( rpCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 196                rpCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 197       poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 198                 poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 199     }
 200     if (!(rpOK || poiOK)) continue;
 201
 202     //make new AliFLowTrack
 203     AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 204     pTrack->SetSource(AliFlowTrack::kFromESD);
 205
 206     //marking the particles used for int. flow:
 207     if(rpOK && rpCFManager)
 208     {
 209       pTrack->SetForRPSelection(kTRUE);
 210       fNumberOfRPs++;
 211     }
 212     //marking the particles used for diff. flow:
 213     if(poiOK && poiCFManager)
 214     {
 215       pTrack->SetForPOISelection(kTRUE);
 216     }
 217
 218     AddTrack(pTrack);
 219   }//end of while (itrkN < iNumberOfInputTracks)
 220 }
 221
 222 //-----------------------------------------------------------------------
 223 AliFlowEvent::AliFlowEvent( const AliAODEvent* anInput,
 224                             const AliCFManager* rpCFManager,
 225                             const AliCFManager* poiCFManager):
 226   AliFlowEventSimple(20)
 227 {
 228   //Fills the event from the AOD
 229   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 230
 231   //loop over tracks
 232   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 233   {
 234     AliAODTrack* pParticle = anInput->GetTrack(itrkN);   //get input particle
 235
 236     //check if pParticle passes the cuts
 237     Bool_t rpOK = kTRUE;
 238     Bool_t poiOK = kTRUE;
 239     if (rpCFManager && poiCFManager)
 240     {
 241       rpOK = ( rpCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 242                rpCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 243       poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 244                 poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 245     }
 246     if (!(rpOK || poiOK)) continue;
 247
 248     //make new AliFlowTrack
 249     AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 250     pTrack->SetSource(AliFlowTrack::kFromAOD);
 251
 252     if (rpOK /* && rpCFManager */ ) // to be fixed - with CF managers uncommented only empty events (NULL in header files)
 253     {
 254       pTrack->SetForRPSelection(kTRUE);
 255       fNumberOfRPs++;
 256     }
 257     if (poiOK /* && poiCFManager*/ )
 258     {
 259       pTrack->SetForPOISelection(kTRUE);
 260     }
 261     AddTrack(pTrack);
 262   }
 263
 264   //  if (iSelParticlesRP >= fMinMult && iSelParticlesRP <= fMaxMult)
 265   //  {
 266   //    if ( (++fCount % 100) == 0)
 267   //    {
 268   //      if (!fMCReactionPlaneAngle == 0) cout<<" MC Reaction Plane Angle = "<<  fMCReactionPlaneAngle << endl;
 269   //      else cout<<" MC Reaction Plane Angle = unknown "<< endl;
 270   //      cout<<" iGoodTracks = "<<iGoodTracks<<endl;
 271   //      cout<<" # of RP selected tracks = "<<iSelParticlesRP<<endl;
 272   //      cout<<" # of POI selected tracks = "<<iSelParticlesPOI<<endl;
 273   //      cout << "# " << fCount << " events processed" << endl;
 274   //    }
 275   //    return pEvent;
 276   //  }
 277   //  else
 278   //  {
 279   //    cout<<"Not enough tracks in the FlowEventSimple"<<endl;
 280   //    return 0;
 281   //  }
 282   //}
 283   //else
 284   //{
 285   //  cout<<"Event does not pass multiplicity cuts"<<endl;
 286   //  return 0;
 287   //}
 288
 289 }
 290
 291 //-----------------------------------------------------------------------
 292 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
 293                             const AliMCEvent* anInputMc,
 294                             KineSource anOption,
 295                             const AliCFManager* rpCFManager,
 296                             const AliCFManager* poiCFManager ):
 297   AliFlowEventSimple(20)
 298 {
 299   //fills the event with tracks from the ESD and kinematics from the MC info via the track label
 300   if (anOption==kNoKine)
 301   {
 302     AliFatal("WRONG OPTION IN AliFlowEventMaker::FillTracks(AliESDEvent* anInput, AliMCEvent* anInputMc, KineSource anOption)");
 303     exit(1);
 304   }
 305
 306   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 307
 308   Int_t iNumberOfInputTracksMC = anInputMc->GetNumberOfTracks() ;
 309   if (iNumberOfInputTracksMC==-1)
 310   {
 311     AliError("Skipping Event -- No MC information available for this event");
 312     return;
 313   }
 314
 315   //loop over ESD tracks
 316   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 317   {
 318     AliESDtrack* pParticle = anInput->GetTrack(itrkN);   //get input particle
 319     //get Label
 320     Int_t iLabel = pParticle->GetLabel();
 321     //match to mc particle
 322     AliMCParticle* pMcParticle = (AliMCParticle*) anInputMc->GetTrack(TMath::Abs(iLabel));
 323
 324     //check
 325     if (TMath::Abs(pParticle->GetLabel())!=pMcParticle->Label())
 326       AliWarning(Form("pParticle->GetLabel()!=pMcParticle->Label(), %i, %i", pParticle->GetLabel(), pMcParticle->Label()));
 327
 328     //check if pParticle passes the cuts
 329     Bool_t rpOK = kTRUE;
 330     Bool_t poiOK = kTRUE;
 331     if (rpCFManager && poiCFManager)
 332     {
 333       if(anOption == kESDkine)
 334       {
 335         if (rpCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle,"mcGenCuts1") &&
 336             rpCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle))
 337           rpOK=kTRUE;
 338         if (poiCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle,"mcGenCuts2") &&
 339             poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle))
 340           poiOK=kTRUE;
 341       }
 342       else if (anOption == kMCkine)
 343       {
 344         if (rpCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle))
 345           rpOK=kTRUE;
 346         if (poiCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,pMcParticle))
 347           poiOK=kTRUE;
 348       }
 349     }
 350
 351     if (!(rpOK || poiOK)) continue;
 352
 353     //make new AliFlowTrack
 354     AliFlowTrack* pTrack = NULL;
 355     if(anOption == kESDkine)   //take the PID from the MC & the kinematics from the ESD
 356     {
 357       pTrack = new AliFlowTrack(pParticle);
 358     }
 359     else if (anOption == kMCkine)   //take the PID and kinematics from the MC
 360     {
 361       pTrack = new AliFlowTrack(pMcParticle);
 362     }
 363
 364     if (rpOK && rpCFManager)
 365     {
 366       fNumberOfRPs++;
 367       pTrack->SetForRPSelection();
 368     }
 369     if (poiOK && poiCFManager) pTrack->SetForPOISelection();
 370
 371     AddTrack(pTrack);
 372   }
 373   SetMCReactionPlaneAngle(anInputMc);
 374 }
 375
 376 //-----------------------------------------------------------------------
 377 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
 378                             const AliMultiplicity* anInputTracklets,
 379                             const AliCFManager* poiCFManager ):
 380   AliFlowEventSimple(20)
 381 {
 382
 383   //Select the particles of interest from the ESD
 384   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 385
 386   //loop over tracks
 387   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 388     {
 389       AliESDtrack* pParticle = anInput->GetTrack(itrkN);   //get input particle
 390
 391       //check if pParticle passes the cuts
 392       Bool_t poiOK = kTRUE;
 393       if (poiCFManager)
 394         {
 395           poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 396                     poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 397         }
 398       if (!poiOK) continue;
 399
 400       //make new AliFLowTrack
 401       AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 402
 403       //marking the particles used for the particle of interest (POI) selection:
 404       if(poiOK && poiCFManager)
 405         {
 406           pTrack->SetForPOISelection(kTRUE);
 407           pTrack->SetSource(AliFlowTrack::kFromESD);
 408         }
 409
 410       AddTrack(pTrack);
 411     }//end of while (itrkN < iNumberOfInputTracks)
 412
 413   //Select the reference particles from the SPD tracklets
 414   anInputTracklets = anInput->GetMultiplicity();
 415   Int_t multSPD = anInputTracklets->GetNumberOfTracklets();
 416
 417   //loop over tracklets
 418   for (Int_t itracklet=0; itracklet<multSPD; ++itracklet) {
 419     Float_t thetaTr= anInputTracklets->GetTheta(itracklet);
 420     Float_t phiTr= anInputTracklets->GetPhi(itracklet);
 421     // calculate eta
 422     Float_t etaTr = -TMath::Log(TMath::Tan(thetaTr/2.));
 423
 424     //make new AliFLowTrackSimple
 425     AliFlowTrack* pTrack = new AliFlowTrack();
 426     pTrack->SetPt(0.0);
 427     pTrack->SetEta(etaTr);
 428     pTrack->SetPhi(phiTr);
 429     //marking the particles used for the reference particle (RP) selection:
 430     fNumberOfRPs++;
 431     pTrack->SetForRPSelection(kTRUE);
 432     pTrack->SetSource(AliFlowTrack::kFromTracklet);
 433
 434     //Add the track to the flowevent
 435     AddTrack(pTrack);
 436   }
 437
 438 }
 439
 440 //-----------------------------------------------------------------------
 441 AliFlowEvent::AliFlowEvent( const AliESDEvent* esd,
 442                             const AliCFManager* poiCFManager,
 443                             Bool_t hybrid):
 444   AliFlowEventSimple(20)
 445 {
 446
 447   //Select the particles of interest from the ESD
 448   Int_t iNumberOfInputTracks = esd->GetNumberOfTracks() ;
 449
 450   //Double_t gPt = 0.0, gP = 0.0;
 451   Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};  //The impact parameters and their covariance.
 452   Double_t dca3D = 0.0;
 453
 454   AliESDtrack trackTPC;
 455
 456   //loop over tracks
 457   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 458     {
 459
 460       if (!esd->GetTrack(itrkN)) continue;
 461
 462       Bool_t useTPC = kFALSE;
 463
 464       AliESDtrack* pParticle = esd->GetTrack(itrkN);   //get input particle
 465
 466       //check if pParticle passes the cuts
 467       Bool_t poiOK = kTRUE;
 468
 469       if (poiCFManager)
 470       {
 471         poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 472                   poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 473       }
 474
 475       if (!(poiOK)) continue;
 476
 477       AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)pParticle->GetTPCInnerParam();
 478
 479       if (tpcTrack)
 480       {
 481
 482 //      gPt = tpcTrack->Pt();
 483 //      gP = tpcTrack->P();
 484
 485         useTPC = kTRUE;
 486
 487         const AliESDVertex *vertexSPD = esd->GetPrimaryVertexSPD();
 488         const AliESDVertex *vertexTPC = esd->GetPrimaryVertexTPC();
 489
 490         if(hybrid)
 491           tpcTrack->PropagateToDCA(vertexSPD,esd->GetMagneticField(),100.,dca,cov);
 492         else
 493           tpcTrack->PropagateToDCA(vertexTPC,esd->GetMagneticField(),100.,dca,cov);
 494
 495         dca3D = TMath::Sqrt(TMath::Power(dca[0],2)+TMath::Power(dca[1],2));
 496
 497       }
 498
 499       //make new AliFLowTrack
 500       AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 501
 502       pTrack->SetSource(AliFlowTrack::kFromESD);
 503
 504       //marking the particles used for diff. flow:
 505       if(poiOK && poiCFManager)
 506       {
 507         pTrack->SetForPOISelection(kTRUE);
 508       }
 509
 510       if(useTPC)
 511       {
 512         pTrack->SetForRPSelection(kTRUE);
 513         fNumberOfRPs++;
 514       }
 515
 516       AddTrack(pTrack);
 517
 518     }//end of while (itrkN < iNumberOfInputTracks)
 519
 520 }
 521
 522 //-----------------------------------------------------------------------
 523 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
 524                             const TH2F* anInputFMDhist,
 525                             const AliCFManager* poiCFManager ):
 526   AliFlowEventSimple(20)
 527 {
 528
 529   //Select the particles of interest from the ESD
 530   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 531
 532   //loop over tracks
 533   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 534     {
 535       AliESDtrack* pParticle = anInput->GetTrack(itrkN);   //get input particle
 536
 537       //check if pParticle passes the cuts
 538       Bool_t poiOK = kTRUE;
 539       if (poiCFManager)
 540         {
 541           poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 542                     poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 543         }
 544       if (!poiOK) continue;
 545
 546       //make new AliFLowTrack
 547       AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 548
 549       //marking the particles used for the particle of interest (POI) selection:
 550       if(poiOK && poiCFManager)
 551         {
 552           pTrack->SetForPOISelection(kTRUE);
 553           pTrack->SetSource(AliFlowTrack::kFromESD);
 554         }
 555
 556       AddTrack(pTrack);
 557     }//end of while (itrkN < iNumberOfInputTracks)
 558
 559   //Select the reference particles from the FMD hits
 560   //loop over FMD histogram
 561   Int_t iBinsEta = anInputFMDhist->GetNbinsX();
 562   Int_t iBinsPhi = anInputFMDhist->GetNbinsY();
 563
 564   for (Int_t iEta = 1; iEta <= iBinsEta; iEta++){
 565     Double_t etaFMD = anInputFMDhist->GetXaxis()->GetBinCenter(iEta);
 566     for (Int_t iPhi = 1; iPhi <= iBinsPhi; iPhi++){
 567       Double_t phiFMD = anInputFMDhist->GetYaxis()->GetBinCenter(iPhi);
 568       Double_t weightFMD = anInputFMDhist->GetBinContent(iEta,iPhi);
 569
 570       if (weightFMD > 0.0) { //do not add empty bins
 571         //make new AliFLowTrackSimple
 572         AliFlowTrack* pTrack = new AliFlowTrack();
 573         pTrack->SetPt(0.0);
 574         pTrack->SetEta(etaFMD);
 575         pTrack->SetPhi(phiFMD);
 576         pTrack->SetWeight(weightFMD);
 577         //marking the particles used for the reference particle (RP) selection:
 578         pTrack->TagRP();
 579         fNumberOfRPs++;
 580         pTrack->SetSource(AliFlowTrack::kFromFMD);
 581
 582         //Add the track to the flowevent
 583         AddTrack(pTrack);
 584
 585       }
 586     }
 587   }
 588
 589 }
 590
 591 //-----------------------------------------------------------------------
 592 AliFlowEvent::AliFlowEvent( AliFlowTrackCuts* rpCuts,
 593                             AliFlowTrackCuts* poiCuts ):
 594   AliFlowEventSimple(20)
 595 {
 596   //Fills the event from a vevent: AliESDEvent,AliAODEvent,AliMCEvent
 597   //the input data needs to be attached to the cuts
 598   //we have two cases, if we're cutting the same collection of tracks
 599   //(same param type) then we can have tracks that are both rp and poi
 600   //in the other case we want to have two exclusive sets of rps and pois
 601   //e.g. one tracklets, the other PMD or global - USER IS RESPOSIBLE
 602   //FOR MAKING SURE THEY DONT OVERLAP OR ELSE THE SAME PARTICLE WILL BE
 603   //TAKEN TWICE
 604
 605   AliFlowTrackCuts::trackParameterType sourceRP = rpCuts->GetParamType();
 606   AliFlowTrackCuts::trackParameterType sourcePOI = poiCuts->GetParamType();
 607
 608   if (!rpCuts || !poiCuts) return;
 609   if (sourceRP==sourcePOI)
 610   {
 611     //loop over tracks
 612     for (Int_t i=0; i<rpCuts->GetNumberOfInputObjects(); i++)
 613     {
 614       //get input object (particle)
 615       TObject* particle = rpCuts->GetInputObject(i);
 616
 617       Bool_t rp = rpCuts->IsSelected(particle,i);
 618       Bool_t poi = poiCuts->IsSelected(particle,i);
 619
 620       if (!(rp||poi)) continue;
 621
 622       //make new AliFLowTrack
 623       AliFlowTrack* pTrack = NULL;
 624       if (rp)
 625       {
 626         pTrack = rpCuts->MakeFlowTrack();
 627         if (!pTrack) continue;
 628         pTrack->TagRP(); fNumberOfRPs++;
 629         if (poi) pTrack->TagPOI();
 630       }
 631       else
 632       if (poi)
 633       {
 634         pTrack = poiCuts->MakeFlowTrack();
 635         if (!pTrack) continue;
 636         pTrack->TagPOI();
 637       }
 638
 639       AddTrack(pTrack);
 640     }//end of while (i < numberOfTracks)
 641   }
 642   else if (sourceRP!=sourcePOI)
 643   {
 644     //here we have two different sources of particles, so we fill
 645     //them independently
 646     AliFlowTrack* pTrack = NULL;
 647     //RP
 648     for (Int_t i=0; i<rpCuts->GetNumberOfInputObjects(); i++)
 649     {
 650       TObject* particle = rpCuts->GetInputObject(i);
 651       Bool_t rp = rpCuts->IsSelected(particle,i);
 652       if (!rp) continue;
 653       pTrack = rpCuts->MakeFlowTrack();
 654       if (!pTrack) continue;
 655       pTrack->TagRP(); fNumberOfRPs++;
 656       AddTrack(pTrack);
 657     }
 658     //POI
 659     for (Int_t i=0; i<poiCuts->GetNumberOfInputObjects(); i++)
 660     {
 661       TObject* particle = poiCuts->GetInputObject(i);
 662       Bool_t poi = poiCuts->IsSelected(particle,i);
 663       if (!poi) continue;
 664       pTrack = poiCuts->MakeFlowTrack();
 665       if (!pTrack) continue;
 666       pTrack->TagPOI();
 667       AddTrack(pTrack);
 668     }
 669   }
 670 }
 671
 672 //-------------------------------------------------------------------//
 673 //---- Including PMD tracks as RP --------------------------//
 674
 675 AliFlowEvent::AliFlowEvent( const AliESDEvent* anInput,
 676                             const AliESDPmdTrack *pmdtracks,
 677                             const AliCFManager* poiCFManager ):
 678   AliFlowEventSimple(20)
 679 {
 680   Float_t GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos);
 681   Float_t GetPmdPhi(Float_t xPos, Float_t yPos);
 682   //Select the particles of interest from the ESD
 683   Int_t iNumberOfInputTracks = anInput->GetNumberOfTracks() ;
 684
 685   //loop over tracks
 686   for (Int_t itrkN=0; itrkN<iNumberOfInputTracks; itrkN++)
 687     {
 688       AliESDtrack* pParticle = anInput->GetTrack(itrkN);   //get input particle
 689       //check if pParticle passes the cuts
 690       Bool_t poiOK = kTRUE;
 691       if (poiCFManager)
 692         {
 693           poiOK = ( poiCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,pParticle) &&
 694                     poiCFManager->CheckParticleCuts(AliCFManager::kPartSelCuts,pParticle));
 695         }
 696       if (!poiOK) continue;
 697
 698       //make new AliFLowTrack
 699       AliFlowTrack* pTrack = new AliFlowTrack(pParticle);
 700
 701       //marking the particles used for the particle of interest (POI) selection:
 702       if(poiOK && poiCFManager)
 703         {
 704           pTrack->SetForPOISelection(kTRUE);
 705           pTrack->SetSource(AliFlowTrack::kFromESD);
 706         }
 707
 708       AddTrack(pTrack);
 709     }//end of while (itrkN < iNumberOfInputTracks)
 710
 711   //Select the reference particles from the PMD tracks
 712   Int_t npmdcl = anInput->GetNumberOfPmdTracks();
 713   printf("======There are %d PMD tracks in this event\n-------",npmdcl);
 714   //loop over clusters
 715   for(Int_t iclust=0; iclust < npmdcl; iclust++){
 716     //AliESDPmdTrack *pmdtr = anInput->GetPmdTrack(iclust);
 717     pmdtracks = anInput->GetPmdTrack(iclust);
 718     Int_t   det   = pmdtracks->GetDetector();
 719     //Int_t   smn   = pmdtracks->GetSmn();
 720     Float_t clsX  = pmdtracks->GetClusterX();
 721     Float_t clsY  = pmdtracks->GetClusterY();
 722     Float_t clsZ  = pmdtracks->GetClusterZ();
 723     Float_t ncell = pmdtracks->GetClusterCells();
 724     Float_t adc   = pmdtracks->GetClusterADC();
 725     //Float_t pid   = pmdtracks->GetClusterPID();
 726     Float_t etacls = GetPmdEta(clsX,clsY,clsZ);
 727     Float_t phicls = GetPmdPhi(clsX,clsY);
 728     //make new AliFLowTrackSimple
 729     AliFlowTrack* pTrack = new AliFlowTrack();
 730     //if(det == 0){ //selecting preshower plane only
 731     if(det == 0 && adc > 270 && ncell > 1){ //selecting preshower plane only
 732       //pTrack->SetPt(adc);//cluster adc
 733       pTrack->SetPt(0.0);
 734       pTrack->SetEta(etacls);
 735       pTrack->SetPhi(phicls);
 736       //marking the particles used for the reference particle (RP) selection:
 737       fNumberOfRPs++;
 738       pTrack->SetForRPSelection(kTRUE);
 739       pTrack->SetSource(AliFlowTrack::kFromPMD);
 740       //Add the track to the flowevent
 741       AddTrack(pTrack);
 742     }//if det
 743   }
 744 }
 745 //----------------------------------------------------------------------------//
 746 Float_t GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos)
 747 {
 748   Float_t rpxpy, theta, eta;
 749   rpxpy  = TMath::Sqrt(xPos*xPos + yPos*yPos);
 750   theta  = TMath::ATan2(rpxpy,zPos);
 751   eta    = -TMath::Log(TMath::Tan(0.5*theta));
 752   return eta;
 753 }
 754 //--------------------------------------------------------------------------//
 755 Float_t GetPmdPhi(Float_t xPos, Float_t yPos)
 756 {
 757   Float_t pybypx, phi = 0., phi1;
 758   if(xPos==0)
 759     {
 760       if(yPos>0) phi = 90.;
 761       if(yPos<0) phi = 270.;
 762     }
 763   if(xPos != 0)
 764     {
 765       pybypx = yPos/xPos;
 766       if(pybypx < 0) pybypx = - pybypx;
 767       phi1 = TMath::ATan(pybypx)*180./3.14159;
 768
 769       if(xPos > 0 && yPos > 0) phi = phi1;        // 1st Quadrant
 770       if(xPos < 0 && yPos > 0) phi = 180 - phi1;  // 2nd Quadrant
 771       if(xPos < 0 && yPos < 0) phi = 180 + phi1;  // 3rd Quadrant
 772       if(xPos > 0 && yPos < 0) phi = 360 - phi1;  // 4th Quadrant
 773
 774     }
 775   phi = phi*3.14159/180.;
 776   return   phi;
 777 }
 778 //---------------------------------------------------------------//
 779
 780