Improved version, kinetree included
authorsnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 23 Mar 2007 16:32:54 +0000 (16:32 +0000)
committersnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 23 Mar 2007 16:32:54 +0000 (16:32 +0000)
22 files changed:
PWG2/FLOW/AliFlowAnalyser.cxx
PWG2/FLOW/AliFlowAnalyser.h
PWG2/FLOW/AliFlowEvent.cxx
PWG2/FLOW/AliFlowEvent.h
PWG2/FLOW/AliFlowKineMaker.cxx [new file with mode: 0644]
PWG2/FLOW/AliFlowKineMaker.h [new file with mode: 0644]
PWG2/FLOW/AliFlowMaker.cxx
PWG2/FLOW/AliFlowMaker.h
PWG2/FLOW/AliFlowSelection.cxx
PWG2/FLOW/AliFlowSelection.h
PWG2/FLOW/AliFlowTrack.cxx
PWG2/FLOW/AliFlowTrack.h
PWG2/FLOW/AliFlowV0.cxx
PWG2/FLOW/AliFlowV0.h
PWG2/FLOW/AliFlowWeighter.cxx
PWG2/FLOW/AliSelectorFoF.cxx
PWG2/FLOW/macros/compilePWG2
PWG2/FLOW/macros/kOpen.C [new file with mode: 0644]
PWG2/FLOW/macros/testAnal.C
PWG2/FLOW/macros/testFoF.C
PWG2/FLOW/macros/testKiner.C [new file with mode: 0644]
PWG2/FLOW/macros/testMaker.C

index 4a1028f..509f8bb 100644 (file)
@@ -11,7 +11,7 @@
 //         the AliFlowAnalyser provides the methods to perform an Event 
 // plane flow analysis over AliFlowEvents. 
 // - The method Init() generates the histograms.
-// - The method Analyse(AliFlowEvent*) calculates/extracts flow observables,  
+// - The method Analyze(AliFlowEvent*) calculates/extracts flow observables,  
 //  fills some histograms and performs the E.P. analysis.
 // - The method Resolution() calculates the resolution correction and 
 //  applies it to observed v_n values.
@@ -84,7 +84,10 @@ using namespace std; //required for resolving the 'cout' symbol
 
 ClassImp(AliFlowAnalyser) 
 //-----------------------------------------------------------------------
-AliFlowAnalyser::AliFlowAnalyser(const AliFlowSelection* flowSelect)
+AliFlowAnalyser::AliFlowAnalyser(const AliFlowSelection* flowSelect):
+  fHistFile(0x0), fPhiWgtFile(0x0),
+  fFlowEvent(0x0), fFlowTrack(0x0), fFlowV0(0x0), fFlowSelect(0x0), fFlowTracks(0x0), fFlowV0s(0x0), 
+  fPhiWgtHistList(0x0), fVnResHistList(0x0)
 {
  // default constructor (selection given or default selection)
 
@@ -102,18 +105,19 @@ AliFlowAnalyser::AliFlowAnalyser(const AliFlowSelection* flowSelect)
  fPhiMax  = 2*TMath::Pi() ; 
 
  // flags 
- fTrackLoop    = kTRUE ;  // main loop for charged tracks
- fV0loop       = kTRUE ;  // correlation analysis is done also for neutral secundary vertex
- fShuffle      = kFALSE ; // randomly reshuffles tracks
- fV1Ep1Ep2     = kFALSE ; // disabled by default               
- fEtaSub       = kFALSE ; // eta subevents
- fReadPhiWgt   = kFALSE ; // if kTRUE & if flowPhiWgt file is there -> Phi Weights are used                                            
- fBayWgt       = kFALSE ; // Bayesian Weights for P.Id. used   
- fRePid        = kFALSE ; // recalculates the P.Id. (becomes kTRUE if new bayesian weights are plugged in)                     
- fPtWgt        = kFALSE ; // pT as a weight
- fEtaWgt       = kFALSE ; // eta as a weight
- fOnePhiWgt    = kTRUE  ; // one/three phi-wgt histogram(s)
-//  fMaxLabel  = 1000 ;     // for labels histogram
+ fTrackLoop     = kTRUE ;  // main loop for charged tracks
+ fV0loop        = kTRUE ;  // correlation analysis is done also for neutral secundary vertex
+ fShuffle       = kFALSE ; // randomly reshuffles tracks
+ fV1Ep1Ep2      = kFALSE ; // disabled by default               
+ fEtaSub        = kFALSE ; // eta subevents
+ fReadPhiWgt    = kFALSE ; // if kTRUE & if flowPhiWgt file is there -> Phi Weights are used                                           
+ fBayWgt        = kFALSE ; // Bayesian Weights for P.Id. used   
+ fRePid         = kFALSE ; // recalculates the P.Id. (becomes kTRUE if new bayesian weights are plugged in)                    
+ fPtWgt         = kFALSE ; // pT as a weight
+ fEtaWgt        = kFALSE ; // eta as a weight
+ fOnePhiWgt     = kTRUE  ; // one/three phi-wgt histogram(s)
+ fCustomRespFunc = kFALSE ; // custom "detector response function" used for P.Id
+// fMaxLabel    = 1000 ;    // for labels histogram
 
  fPhiWgtHistList = 0 ;
  fVnResHistList = 0 ;
@@ -1926,7 +1930,7 @@ void AliFlowAnalyser::Weightening()
 //-----------------------------------------------------------------------
 // ###
 //-----------------------------------------------------------------------
-Bool_t AliFlowAnalyser::Analyse(AliFlowEvent* flowEvent)         
+Bool_t AliFlowAnalyser::Analyze(AliFlowEvent* flowEvent)         
 {
  // Runs the analysis on the AliFlowEvent (* fFlowEvent). 
  // This method can be inserted in a loop over a collection of 
@@ -1949,6 +1953,8 @@ Bool_t AliFlowAnalyser::Analyse(AliFlowEvent* flowEvent)
   if(fReadPhiWgt) { FillEvtPhiWgt(fFlowEvent) ; }          // phi and bayesian weights are filled previous to the loop (FillWgtArrays(TFile*))
   else                   { fFlowEvent->SetNoWgt() ; }              // phi weights can be used or not , this plugs them into the event
 
+  fFlowEvent->SetCustomRespFunc(fCustomRespFunc) ;         // a custom "detector response function" is used for assigning P.Id
+
   if(fBayWgt)    { FillBayesianWgt(fFlowEvent) ; }         // bayesian weights can be used or not , this plugs them into the event
   if(fRePid)     { fFlowEvent->SetPids() ; }               // re-calculate all p.id. hypotesis with the (new) bayesian array
 
index 69861d2..84fa870 100644 (file)
@@ -49,7 +49,7 @@ public:
   Bool_t   Finish() ;                                          // Saves histograms, Closes stuff
 
  // Analysis of 1 event (can be called from outside)
-  Bool_t   Analyse(AliFlowEvent* flowEvent = 0) ;              // Fills the defaults histograms (init them first!) and performs the calculation for the given event         
+  Bool_t   Analyze(AliFlowEvent* flowEvent = 0) ;              // Fills the defaults histograms (init them first!) and performs the calculation for the given event         
 
  // Resolution corrections to v_n (call it at the end of the evts loop)
   Bool_t   Resolution() ;                                      // Calculates resolution and mean flow values
@@ -69,6 +69,7 @@ public:
   void     SetUseFirstLastPhiWgt(Bool_t flw = kTRUE)           { fOnePhiWgt = !flw ; }         // uses 3 wgt histograms
   void    SetFlowForV0(Bool_t v0 = kTRUE)                      { fV0loop = v0 ; }              // Enables Flow study for v0
   void    SetTrackLoop(Bool_t trkl = kTRUE)                    { fTrackLoop = trkl ; }         // Enables Tracks loop (keep it kTRUE)
+  void    SetCustomRespFunc(Bool_t crf = kTRUE)                { fCustomRespFunc = crf ; }     // Enables to use a custom detector response function
   //void     SetDebugg(Int_t db = 1) ;                                 // set the cout's for debug (default is 1)
 
  // Histograms
@@ -129,6 +130,7 @@ public:
   Bool_t          fReadPhiWgt ;                                //! Phi Weights are applied to Phi distrib. (default is false)
   Bool_t          fBayWgt ;                                    //! Bayesian Weights are applied to P.Id. (default is false) 
   Bool_t          fRePid ;                                     //! Re-Calculates the P.Id. basing on the bayesian wgts (if plugged in)
+  Bool_t          fCustomRespFunc ;                            //! A custom "detector response function" is used for P.Id (default is false -> the combined response function from the ESD will be used)
 
   Bool_t          fPtWgt ;                                     //! flag to use pT as a weight for RP determination
   Bool_t          fEtaWgt ;                                    //! flag to use eta as a weight for RP determination
index f2a3ab1..93abd5e 100644 (file)
@@ -54,14 +54,16 @@ using namespace std; //required for resolving the 'cout' symbol
 // - Flags & Sets
 Bool_t  AliFlowEvent::fgPtWgt         = kFALSE ;  // gives pT-based weights
 Bool_t  AliFlowEvent::fgEtaWgt        = kFALSE ;  // gives eta-based weights
-Bool_t  AliFlowEvent::fgOnePhiWgt       = kTRUE  ;  // kTRUE --> ENABLEs SINGLE WEIGHT ISTOGRAM , kFALSE --> ENABLEs 3 WEIGHT ISTOGRAMS
+Bool_t  AliFlowEvent::fgOnePhiWgt      = kTRUE  ;  // kTRUE --> ENABLEs SINGLE WEIGHT ISTOGRAM , kFALSE --> ENABLEs 3 WEIGHT ISTOGRAMS
 Bool_t  AliFlowEvent::fgNoWgt         = kFALSE ;  // No Weight is used
 // - Eta Sub-Events (later used to calculate the resolution)
-Bool_t  AliFlowEvent::fgEtaSubs               = kFALSE ;  // makes eta subevents
+Bool_t  AliFlowEvent::fgEtaSubs        = kFALSE ;  // makes eta subevents
+Bool_t  AliFlowEvent::fgCustomRespFunc = kFALSE ;  // custom "detector response function" is used for P.Id (see AliFlowTrack)
 
 ClassImp(AliFlowEvent) 
 //-----------------------------------------------------------
-AliFlowEvent::AliFlowEvent(Int_t lenght)  
+AliFlowEvent::AliFlowEvent(Int_t lenght):
+  fTrackCollection(0x0), fV0Collection(0x0)    
 {
  // Default constructor: initializes the ObjArray of FlowTracks and FlowV0s, 
  // cleans the internal variables, sets all the weights to 1, sets default flags.
@@ -570,22 +572,27 @@ void AliFlowEvent::SetSelections(AliFlowSelection* pFlowSelect)
 //-------------------------------------------------------------
 void AliFlowEvent::SetPids()
 {
- // Re-sets the tracks P.id. (using the current AliFlowConstants::fgBayesian[] array)
+ // Re-sets the tracks P.id. (using the current AliFlowConstants::fgBayesian[] array).
+ // To use the Raw P.Id (just the detector response function), set fgBayesian[] = {1,1,1,1,1,1}.
  
  const Int_t kCode[] = {11,13,211,321,2212,10010020} ;
  for(Int_t itr=0;itr<TrackCollection()->GetEntries();itr++) 
  {
   AliFlowTrack* pFlowTrack ;
   pFlowTrack = (AliFlowTrack*)TrackCollection()->At(itr) ;
-  TVector bayPid = pFlowTrack->PidProbs() ;
-  Int_t maxN = 2 ;                // if No id. -> then is a Pi
-  Float_t pidMax = bayPid[2] ;    // (if all equal, Pi probability get's the advantage to be the first)
+
+  Float_t bayPid[AliFlowConstants::kPid] ;
+  if(fgCustomRespFunc)  { pFlowTrack->PidProbsC(bayPid) ; }
+  else                         { pFlowTrack->PidProbs(bayPid)  ; }
+
+  Int_t maxN = 2 ;               // if No id. -> then is a Pi
+  Float_t pidMax = bayPid[2] ;    // (if all equal, Pi probability was the first)
   for(Int_t nP=0;nP<AliFlowConstants::kPid;nP++)
   {
    if(bayPid[nP]>pidMax) { maxN = nP ; pidMax = bayPid[nP] ; }
   }
   Int_t pdgCode = TMath::Sign(kCode[maxN],pFlowTrack->Charge()) ;     
-  pFlowTrack->SetMostLikelihoodPID(pdgCode);                   
+  pFlowTrack->SetMostLikelihoodPID(pdgCode) ;                  
  }
 }
 //-------------------------------------------------------------
index 2909c7a..5e0155a 100644 (file)
@@ -75,6 +75,7 @@ public:
   Bool_t     OnePhiWgt() const       { return fgOnePhiWgt ; }                        // Returns flag for using just one phi weight
   Bool_t     NoWgt() const           { return fgNoWgt; }                             // returns kTRUE if weight are NOT used
   Bool_t     EtaSubs() const         { return fgEtaSubs ; }                          // Returns flag for eta sub-events
+  Bool_t     CustomRespFunc() const   { return fgCustomRespFunc ; }
 
  // Gets
   Int_t      EventID() const               { return fEventID; }                          // Returns ID of the event
@@ -96,10 +97,10 @@ public:
   Float_t    ExtPsi(Int_t harN = 0) const    { if(harN<AliFlowConstants::kHars) { return fExtPsi[harN] ; } else { return 0. ; } } // external RP angle
   Float_t    ExtRes(Int_t harN = 0) const    { if(harN<AliFlowConstants::kHars) { return fExtRes[harN] ; } else { return 0. ; } } // external RP resolution
 
-  Double_t CenterOfMassEnergy() const       { return AliFlowConstants::fgCenterOfMassEnergy ; }   // Returns center of mass energy (5.5 TeV)
-  Double_t MagneticField() const            { return AliFlowConstants::fgMagneticField ; }        // Returns magnetic field value
-  Short_t  BeamMassNumberEast() const       { return AliFlowConstants::fgBeamMassNumberEast ; }   // Returns beam mass (Pb = 208)
-  Short_t  BeamMassNumberWest() const       { return AliFlowConstants::fgBeamMassNumberWest ; }   // Returns beam mass (Pb = 208)
+  Double_t   CenterOfMassEnergy() const             { return AliFlowConstants::fgCenterOfMassEnergy ; }   // Returns center of mass energy (5.5 TeV)
+  Double_t   MagneticField() const          { return AliFlowConstants::fgMagneticField ; }        // Returns magnetic field value
+  Short_t    BeamMassNumberEast() const             { return AliFlowConstants::fgBeamMassNumberEast ; }   // Returns beam mass (Pb = 208)
+  Short_t    BeamMassNumberWest() const             { return AliFlowConstants::fgBeamMassNumberWest ; }   // Returns beam mass (Pb = 208)
 
  // Sets
   void     SetEventID(const Int_t& id)                             { fEventID = id ; }  // Sets Event ID and the Event name (name = evtNumber_runId)
@@ -122,6 +123,7 @@ public:
   static void SetPtWgt(Bool_t ptWgt = kTRUE)                       { fgPtWgt = ptWgt; }
   static void SetEtaWgt(Bool_t etaWgt = kTRUE)                     { fgEtaWgt = etaWgt ; }
   static void SetNoWgt(Bool_t nowgt = kTRUE)                       { fgNoWgt = nowgt ; }  // still for odd harmonics: Wgt = +1 (positive Eta) or -1 (negative Eta)
+  static void SetCustomRespFunc(Bool_t crf = kTRUE)                { fgCustomRespFunc = crf ; }
                
   void    SetBayesian(Double_t bayes[AliFlowConstants::kPid]) const { for(Int_t i=0;i<AliFlowConstants::kPid;i++) { AliFlowConstants::fgBayesian[i] = bayes[i] ; } } // Set the Bayesian vector of particle abundances
   void    SetMagneticField(const Double_t& mf) const               { AliFlowConstants::fgMagneticField = mf; }
@@ -173,14 +175,15 @@ private:
   static Bool_t       fgOnePhiWgt;                               //! flag for phi weights (just one hist)
   static Bool_t       fgNoWgt;                                  //! No Weights (Wgt == 1)
   static Bool_t       fgEtaSubs;                                 //! Flag for making Eta Subevents
+  static Bool_t       fgCustomRespFunc ;                        //! A custom "detector response function" is used for P.Id
 
  // shortcuts (to speed up the execution)
-  Bool_t   fDone ;                                                                                     //! flag setted kTRUE when the loop is done
-  TVector2 fQ[AliFlowConstants::kSels][AliFlowConstants::kHars];                                       //! flow vector
-  UInt_t   fMult[AliFlowConstants::kSels][AliFlowConstants::kHars];                                    //! multiplicity
-  Float_t  fSumOfWeightSqr[AliFlowConstants::kSels][AliFlowConstants::kHars];                                  //! Sqrt(Sum(wgt)) ~ Sqrt(Mult)
-  TVector2 fQSub[AliFlowConstants::kSubs][AliFlowConstants::kSels][AliFlowConstants::kHars];           //! flow vector subs
-  UInt_t   fMultSub[AliFlowConstants::kSubs][AliFlowConstants::kSels][AliFlowConstants::kHars];        //! multiplicity subs
+  Bool_t   fDone ;                                                                             //! flag setted kTRUE when the loop is done
+  TVector2 fQ[AliFlowConstants::kSels][AliFlowConstants::kHars];                               //! flow vector
+  UInt_t   fMult[AliFlowConstants::kSels][AliFlowConstants::kHars];                            //! multiplicity
+  Float_t  fSumOfWeightSqr[AliFlowConstants::kSels][AliFlowConstants::kHars];                          //! Sqrt(Sum(wgt)) ~ Sqrt(Mult)
+  TVector2 fQSub[AliFlowConstants::kSubs][AliFlowConstants::kSels][AliFlowConstants::kHars];    //! flow vector subs
+  UInt_t   fMultSub[AliFlowConstants::kSubs][AliFlowConstants::kSels][AliFlowConstants::kHars]; //! multiplicity subs
 
   ClassDef(AliFlowEvent,2) ;                    // macro for rootcint
 };
diff --git a/PWG2/FLOW/AliFlowKineMaker.cxx b/PWG2/FLOW/AliFlowKineMaker.cxx
new file mode 100644 (file)
index 0000000..1671804
--- /dev/null
@@ -0,0 +1,576 @@
+//////////////////////////////////////////////////////////////////////
+//
+// $Id$
+//
+// Author: Emanuele Simili
+//
+//////////////////////////////////////////////////////////////////////
+//_____________________________________________________________
+//
+// Description: 
+//        AliFlowKineMaker provides the method to create AliFlowEvent(s)
+// creates AliFlowEvent from the KineTree . 
+// TParticle(s) is translated into AliFlowTrack(s), with exact momentum, 
+// P.Id., etc. Very basic track cuts are applyed (like primaries). 
+// The present class can be used in a simple AliRoot macro or in a 
+// more complex enviroment such as AliSelector or AliTask.
+//
+//////////////////////////////////////////////////////////////////////
+
+#ifndef ALIFLOWKINEMAKER_CXX
+#define ALIFLOWKINEMAKER_CXX
+
+// ROOT things
+#include <TROOT.h>
+#include <TFile.h>
+#include <TString.h>
+#include <TMath.h>
+#include <TTree.h>
+#include "TClonesArray.h"
+#include "TParticle.h"
+#include "TParticlePDG.h"
+//#include "TDatabasePDG.h"
+
+// AliRoot things (...not used here, but in the macro)
+//#include "AliRun.h"
+//#include "AliRunLoader.h"
+//#include "AliStack.h"
+
+// Flow things
+#include "AliFlowEvent.h"
+#include "AliFlowTrack.h"
+#include "AliFlowV0.h"
+#include "AliFlowConstants.h"
+#include "AliFlowKineMaker.h"
+
+// ANSI things
+#include <stdlib.h>
+using namespace std; //required for resolving the 'cout' symbol
+
+ClassImp(AliFlowKineMaker) 
+//-----------------------------------------------------------------------
+AliFlowKineMaker::AliFlowKineMaker():
+  fKTree(0x0), fParticle(0x0), fParticlePDG(0x0),
+  fFlowEvent(0x0), fFlowTrack(0x0), fFlowV0(0x0)
+{
+ // default constructor 
+ // resets counters , sets defaults
+ fNewAli = kFALSE ;
+ // flags
+ fLoopParts = kTRUE ;
+ fCounter = 0 ; 
+ // loop variable
+ fRunID = 0 ;       
+ fNumberOfParticles = 0 ;
+ fEventNumber = 0 ; 
+ fPartNumber = 0 ; 
+ fEventNumber = 0 ;
+ fMagField = 0. ;
+ // counters
+ fGoodTracks = 0  ;  
+ fGoodV0s    = 0  ;    
+ fGoodTracksEta = 0  ;
+ fPosiTracks = 0  ;  
+ fNegaTracks = 0  ;  
+ fUnconstrained = 0  ;
+ fCutParts   = 0  ;    
+ for(Int_t bb=0;bb<5;bb++) { fBayesianAll[bb] = 0 ; } ;
+ fSumAll = 0 ; 
+ fCharge = 0 ;
+
+ // trak cuts
+ fPrimary = kTRUE ;
+ fAbsEta = 2.1 ;                       
+ fElow = 0.001 ; fEup = 1000. ;   
+ fLabel[0] = 0 ; fLabel[1] = -1 ;
+
+//  // TGeant3::AddParticlesToPdgDataBase() ---  Stolen From TGeant3.cxx ----(
+//  TDatabasePDG *pdgDB = TDatabasePDG::Instance();
+//  const Int_t kion=10000000;
+//  const Int_t kspe=50000000;
+//  const Double_t kAu2Gev=0.9314943228;
+//  const Double_t khSlash = 1.0545726663e-27;
+//  const Double_t kErg2Gev = 1/1.6021773349e-3;
+//  const Double_t khShGev = khSlash*kErg2Gev;
+//  const Double_t kYear2Sec = 3600*24*365.25;
+//  // Ions
+//  pdgDB->AddParticle("Deuteron","Deuteron",2*kAu2Gev+8.071e-3,kTRUE,0,3,"Ion",kion+10020);
+//  pdgDB->AddParticle("Triton","Triton",3*kAu2Gev+14.931e-3,kFALSE,khShGev/(12.33*kYear2Sec),3,"Ion",kion+10030);
+//  pdgDB->AddParticle("Alpha","Alpha",4*kAu2Gev+2.424e-3,kTRUE,khShGev/(12.33*kYear2Sec),6,"Ion",kion+20040);
+//  pdgDB->AddParticle("HE3","HE3",3*kAu2Gev+14.931e-3,kFALSE,0,6,"Ion",kion+20030);
+//  // Special particles
+//  pdgDB->AddParticle("Cherenkov","Cherenkov",0,kFALSE,0,0,"Special",kspe+50);
+//  pdgDB->AddParticle("FeedbackPhoton","FeedbackPhoton",0,kFALSE,0,0,"Special",kspe+51);
+//  // ----------------------------------------------------------------------)
+}
+//-----------------------------------------------------------------------
+AliFlowKineMaker::~AliFlowKineMaker()
+{
+ // default destructor (no actions) 
+}
+//-----------------------------------------------------------------------
+
+//-----------------------------------------------------------------------
+AliFlowEvent* AliFlowKineMaker::FillFlowEvent(TTree* fKTree)
+{
+ // From the MC KineTree (input) fills the AliFlowEvent (output) . 
+ // It loops on the stored TParticles and calls the methods to fill the 
+ // arrays in the AliFlowEvent (charged -> tracks , neutral -> v0s) . 
+
+ fFlowEvent = new AliFlowEvent() ; if(!fFlowEvent) { return 0 ; }
+ //cout << " -evt- " << fFlowEvent << endl ;
+
+ fRunID = -1 ;
+ fEventNumber = -1 ;
+ fNumberOfParticles = fKTree->GetEntries() ; 
+ //
+ cout << " *evt n. " << fEventNumber << " (run " << fRunID << ")  -  tracks/v0s : " << fNumberOfParticles << endl ;
+
+ // Event id 
+ fFlowEvent->SetRunID(fRunID) ;               
+ fFlowEvent->SetEventID(fEventNumber) ;         
+ fFlowEvent->SetOrigMult((UInt_t)fNumberOfParticles) ;
+
+ // Run information (fixed - ???)
+ fMagField = 4 ; // (?)
+ fFlowEvent->SetMagneticField(fMagField) ;     
+ fFlowEvent->SetCenterOfMassEnergy(AliFlowConstants::fgCenterOfMassEnergy) ;  
+ fFlowEvent->SetBeamMassNumberEast(AliFlowConstants::fgBeamMassNumberEast) ;  
+ fFlowEvent->SetBeamMassNumberWest(AliFlowConstants::fgBeamMassNumberWest) ;  
+
+ // Trigger information (now is: ULon64_t - some trigger mask)
+ fFlowEvent->SetL0TriggerWord(-1); 
+ // Get primary vertex position
+ fVertex[0] = 0. ; fVertex[1] = 0. ; fVertex[2] = 0. ;         // fVertex = // ?! how to get primary vertex !?
+ fFlowEvent->SetVertexPos(fVertex[0],fVertex[1],fVertex[2]) ; 
+
+ // Zero Degree Calorimeter information
+ Int_t zdcp = (Int_t)(TMath::Sqrt(TMath::Sqrt(fNumberOfParticles))) ;
+ Float_t zdce[3] ; zdce[0] = -1 ; zdce[1] = -1 ; zdce[2] = -1 ;
+ fFlowEvent->SetZDCpart(zdcp);                         
+ fFlowEvent->SetZDCenergy(zdce[0],zdce[1],zdce[2]);    
+
+ fKTree->SetBranchAddress("Particles",&fParticle) ;
+
+ // Track (& V0) loop
+ if(fLoopParts)
+ {
+  Int_t badPart = 0 ;   
+  for(fPartNumber=0;fPartNumber<fNumberOfParticles;fPartNumber++) 
+  {
+   fKTree->GetEntry(fPartNumber) ;
+   if(CheckTrack(fParticle))
+   {
+    // fParticlePDG = fParticle->GetPDG() ; 
+    fCharge = (Int_t)((fParticle->GetPDG()->Charge())/3) ; // cout << fCharge << endl ;         
+    // fCharge = (Int_t)(TMath::Sign(1,(fParticle->GetPdgCode()))) ;            
+  
+    if(TMath::Abs(fCharge) > 0)
+    {
+     FillFlowTrack(fParticle) ;   
+     fGoodTracks++ ;                                  
+    }
+    else if(fCharge == 0)
+    {
+     FillFlowV0(fParticle) ;
+     fGoodV0s++ ;                       
+    }
+   }
+   else { badPart++ ; continue ; }
+  }
+  fCutParts += badPart ;
+ } // cout << " -particle number- :  " << fPartNumber << endl ;
+
+ // Evt setting stuff
+ fFlowEvent->SetCentrality();  
+                               
+ fCounter++ ; 
+ return fFlowEvent ;
+}
+//----------------------------------------------------------------------
+AliFlowTrack* AliFlowKineMaker::FillFlowTrack(TParticle* fParticle)
+{
+ // From a charged TParticle (input) fills the AliFlowTrack (output) .
+
+ TString name = "" ; name += fPartNumber ;
+ Int_t idx = fFlowEvent->TrackCollection()->GetEntries() ;
+ fFlowTrack = (AliFlowTrack*)(fFlowEvent->TrackCollection()->New(idx)) ;
+ fFlowTrack->SetName(name.Data()) ;
+ // cout << " -tr- " << name.Data() << "(" << idx << ")"  << endl ;
+
+ // TParticle label (link: KineTree-ESD)
+ Int_t label = TMath::Abs(fPartNumber);
+ fFlowTrack->SetLabel(label) ;                         
+
+ // signed DCA from ESDtrack
+ Float_t x = fParticle->Vx() ; 
+ Float_t y = fParticle->Vy() ; 
+ Float_t z = fParticle->Vz() ; 
+ Float_t xy = TMath::Sqrt(x*x + y*y) ; 
+ fFlowTrack->SetDcaSigned(xy,z) ;                  
+
+ // error on the DCA = 0
+ fFlowTrack->SetDcaError(0.,0.,0.) ;               
+
+ // UnConstrained (global) first
+ Double_t gD[3] ;                              
+ gD[0] = fParticle->Px() ; gD[1] = fParticle->Py() ; gD[2] = fParticle->Pz() ;                 
+ // -
+ Float_t phiGl = (Float_t)Phi(gD) ;  
+ if(phiGl<0) { phiGl += 2*TMath::Pi() ; }
+ fFlowTrack->SetPhiGlobal(phiGl) ;             
+ Float_t ptGl = (Float_t)Pt(gD) ;  if(ptGl<=0) { cout << " !!! ptGlobal = " << ptGl << endl ; }
+ fFlowTrack->SetPtGlobal(ptGl) ;               
+ Float_t etaGl = (Float_t)Eta(gD) ; 
+ fFlowTrack->SetEtaGlobal(etaGl) ;             
+
+ // Constrained (same, if primary)
+ if((fParticle->IsPrimary()) && (Norm(gD)!=0.))  
+ {                                                     
+  fFlowTrack->SetPhi(phiGl) ;                          
+  fFlowTrack->SetPt(ptGl) ;                            
+  fFlowTrack->SetEta(etaGl) ;                          
+
+  // number of constrainable tracks with |eta| < AliFlowConstants::fgEtaGood (0.9)
+  if(TMath::Abs(etaGl) < AliFlowConstants::fgEtaGood)  { fGoodTracksEta++ ; }
+ }
+ else  // in case Constriction impossible for track, fill the UnConstrained (global)
+ {
+  fUnconstrained++ ;   
+  fFlowTrack->SetPhi(0.) ;
+  fFlowTrack->SetPt(0.) ;   
+  fFlowTrack->SetEta(0.) ; 
+ }          
+
+ // positive - negative tracks
+ //Int_t fCharge = (Int_t)(fParticle->GetPDG()->Charge()) ; 
+ fFlowTrack->SetCharge(fCharge) ;              
+ if(fCharge>0)                 { fPosiTracks++ ; }
+ else if(fCharge<0)    { fNegaTracks++ ; }
+ else                  { return 0 ; }
+
+ // Track parametrization (p at, hits, clusters, dE/dx) 
+ Double_t pVecAt[3] ; 
+ for(Int_t gg=0;gg<3;gg++) { pVecAt[gg] = gD[gg] ; }
+ Float_t pAt = (Float_t)Norm(pVecAt) ;
+ Int_t nClus[9] ; Int_t nHits[9] ; 
+ nClus[0] = 1 ;   nHits[0] = 0 ;                           // ITS - pixel
+ nClus[1] = 1 ;   nHits[1] = 0 ;  
+ nClus[2] = 1 ;   nHits[2] = 0 ;                           // ITS - drift
+ nClus[3] = 1 ;   nHits[3] = 0 ;  
+ nClus[4] = 1 ;   nHits[4] = 0 ;                           // ITS - strips
+ nClus[5] = 1 ;   nHits[5] = 0 ;  
+ nClus[6] = 160 ; nHits[6] = 0 ;                           // TPC
+ nClus[7] = 130 ; nHits[7] = 0 ;                           // TRD
+ nClus[8] = 1 ;   nHits[8] = 0 ;                           // TOF
+
+ Int_t pdgcode = 0 ;
+ Float_t dEdx[4] ; for(Int_t de=0;de<4;de++) { dEdx[de] = -1. ; }
+ Float_t detResFun[6] ; for(Int_t de=0;de<6;de++) { detResFun[de] = 0. ; }
+ Float_t zFirst = fVertex[2] ; 
+ Float_t zLast  = 1000. ;
+ Float_t rFirst = TMath::Sqrt((fVertex[0]*fVertex[0]) + (fVertex[1]*fVertex[1])) ; 
+ Float_t rLast  = 1000. ;
+
+// // Geometrical acceptance (calculated assuming straight tracks) 
+//
+//  Float_t rDet[9][2] ; Float_t zDet[9] ;      
+//  rDet[0][0] = 3.9 ;   rDet[0][1] = 3.9 ;   zDet[0] = 14.1/2. ;     // ITS - pixel
+//  rDet[1][0] = 7.6 ;   rDet[1][1] = 7.6 ;   zDet[1] = 14.1/2. ;    
+//  rDet[2][0] = 15.0 ;  rDet[2][1] = 15.0 ;  zDet[2] = 22.2/2. ;     // ITS - drift
+//  rDet[3][0] = 23.9 ;  rDet[3][1] = 23.9 ;  zDet[3] = 29.7/2. ;    
+//  rDet[4][0] = 37.8 ;  rDet[4][1] = 38.4 ;  zDet[4] = 43.1/2. ;     // ITS - strips
+//  rDet[5][0] = 42.8 ;  rDet[5][1] = 43.4 ;  zDet[5] = 48.9/2. ;    
+//  rDet[6][0] = 84.5 ;  rDet[6][1] = 246.6 ; zDet[6] = 500./2. ;     // TPC
+//  rDet[7][0] = 290. ;  rDet[7][1] = 370. ;  zDet[7] = 700./2. ;     // TRD
+//  rDet[8][0] = 370. ;  rDet[8][1] = 399. ;  zDet[8] = 745./2. ;     // TOF
+//
+//  Float_t Atheta = fParticle->Pz()/fParticle->Pt() ; 
+//  Float_t z0 = fParticle->Vz() ; 
+//  Float_t r0 = TMath::Sqrt((fParticle->Vx()*fParticle->Vx())+(fParticle->Vy()*fParticle->Vy())) ;
+//  if((fParticle->Vx()*fParticle->Px()+fParticle->Vy()*fParticle->Py())>0)     { r0 *= 1. ; }   // sign given basing on track direction in respect to position
+//  else if((fParticle->Vx()*fParticle->Px()+fParticle->Vy()*fParticle->Py())<0) { r0 *= -1.; }
+//  else                                                                        { r0  = 0. ; }
+//
+//  // rFirst = rDet[0][0] ; rLast  = rDet[0][0] ;
+//  zFirst = z0 + Atheta * (rDet[0][0] - r0) ; 
+//  zLast  = z0 + Atheta * (rDet[4][1] - r0) ;
+//  Float_t Pin  = 0. ; Float_t Pout = 0. ; Float_t Rout = 0. ; 
+//  for(int dd=0;dd<9;dd++)
+//  {
+//   Pin  = z0 + Atheta * (rDet[dd][0] - r0) ; 
+//   if(Pin<zDet[dd]) 
+//   {
+//    Pout = z0 + Atheta * (rDet[dd][1] - r0) ; 
+//    if(TMath::Abs(Pout<zDet[dd]))                    // track gets in and out inside acceptance -> full hits
+//    { 
+//     nHits[dd] = nClus[dd] ; 
+//     Rout = rDet[dd][1] ; 
+//     rLast = TMath::Abs(Rout) ;              
+//    } 
+//    else                                             // track goes out from one side -> SOME hits (...)
+//    {
+//     Rout = r0 + ((TMath::Sign(zDet[dd],eta)-z0)/Atheta) ; 
+//     rLast = TMath::Abs(Rout) ;      
+//     Float_t proportion = TMath::Abs((rLast-rDet[dd][0])/(rDet[dd][1]-rDet[dd][0])) ; proportion *= nClus[dd] ;      
+//     nHits[dd] = (Int_t)proportion ; 
+//    }                        
+//   }
+//   else                                              // track does not get in -> zero hits
+//   {
+//    nHits[0] = 0 ; rFirst = 0. ; //rLast = 0. ; 
+//   }         
+//  }
+//
+//  if(nHits[7])               { dEdx[0] = 1. ; }  // implement bethe-block for TPC
+//  if(nHits[5] || nHits[6])   { dEdx[1] = 1. ; }  // implement bethe-block for ITS
+//  if(nHits[8])               { dEdx[2] = 1. ; }  // implement transition-radiation for TRD
+//  if(nHits[9])               { dEdx[3] = 1. ; }  // implement time of flight for TOF
+//
+
+ // P.id. (basing on the pdg code from MC -> an exact P.Id (probability=1) is given for [e , mu , pi , K , p , d], others are 0) 
+ pdgcode = fParticle->GetPdgCode() ;                        // cout << FlowDebug << "PDG code = " << pdgcode << endl ; 
+ if(TMath::Abs(pdgcode) == 11)              { detResFun[0] = 1. ; fBayesianAll[0]++ ; }
+ else if(TMath::Abs(pdgcode) == 13)         { detResFun[1] = 1. ; fBayesianAll[1]++ ; }
+ else if(TMath::Abs(pdgcode) == 211)        { detResFun[2] = 1. ; fBayesianAll[2]++ ; } 
+ else if(TMath::Abs(pdgcode) == 321)        { detResFun[3] = 1. ; fBayesianAll[3]++ ; } 
+ else if(TMath::Abs(pdgcode) == 2212)       { detResFun[4] = 1. ; fBayesianAll[4]++ ; }
+ else if(TMath::Abs(pdgcode) == 10010020)    { detResFun[5] = 1. ; fBayesianAll[5]++ ; }
+ else          { for(Int_t de=0;de<6;de++)  { detResFun[de] = 0.2 ; } }
+ fSumAll++ ;
+
+ // Fill the (fake) track parameters (fit , P.Id. ...)
+ fFlowTrack->SetMostLikelihoodPID(pdgcode); 
+
+ fFlowTrack->SetElectronPositronProb(detResFun[0]);              
+ fFlowTrack->SetMuonPlusMinusProb(detResFun[1]);
+ fFlowTrack->SetPionPlusMinusProb(detResFun[2]);
+ fFlowTrack->SetKaonPlusMinusProb(detResFun[3]);
+ fFlowTrack->SetProtonPbarProb(detResFun[4]);
+ fFlowTrack->SetDeuteriumAntiDeuteriumProb(detResFun[5]);    // *!* implement P.Id. for Deuterium
+
+ fFlowTrack->SetZFirstPoint(zFirst) ;  fFlowTrack->SetZLastPoint(zLast) ;           
+ fFlowTrack->SetTrackLength(TMath::Sqrt((zLast-zFirst)*(zLast-zFirst)+(rLast-rFirst)*(rLast-rFirst))) ;
+ fFlowTrack->SetChi2(0.) ;
+
+ // Fill the (fake) detector information (nHits, dE/dx, det, resp. func.) 
+ fFlowTrack->SetFitPtsTPC(nHits[6]) ;               // cout << FlowDebug << "nHits TPC = " << nHits[6] << endl ; 
+ fFlowTrack->SetMaxPtsTPC(nClus[6]) ;               // cout << FlowDebug << "nClus = " << nClus[6] << endl ;  
+ fFlowTrack->SetChi2TPC(nHits[6]/nClus[6]) ;        // cout << FlowDebug << "Chi2 = " << nHits[6]/nClus[6] << endl ;        
+ fFlowTrack->SetDedxTPC(dEdx[0]) ;                  // cout << FlowDebug << "Dedx = " << dEdx << endl ; 
+ fFlowTrack->SetPatTPC(pAt) ;                       // cout << FlowDebug << "p = " << pAt << endl ;                            
+ fFlowTrack->SetRespFunTPC(detResFun) ;             // cout << FlowDebug << "response function = " << detResFun << endl ;  
+ // -
+ Int_t nITShits = 0 ; for(int dd=0;dd<6;dd++) { nITShits += nHits[dd] ; } 
+ Int_t nITSclus = 6 ; 
+ fFlowTrack->SetFitPtsITS(nITShits) ;               // cout << FlowDebug << "nHits ITS = " << nITShits << endl ;                                   
+ fFlowTrack->SetMaxPtsITS(nITSclus) ;               // cout << FlowDebug << "nClus = " << nITSclus << endl ; 
+ fFlowTrack->SetChi2ITS(nITShits/nITSclus) ;        // cout << FlowDebug << "Chi2 = " << nITShits/nITSclus << endl ; 
+ fFlowTrack->SetDedxITS(dEdx[1]) ;                  // cout << FlowDebug << "Dedx = " << dEdx << endl ; 
+ fFlowTrack->SetPatITS(pAt) ;                       // cout << FlowDebug << "p = " << pAt << endl ;                                    
+ fFlowTrack->SetRespFunITS(detResFun) ;             // cout << FlowDebug << "response function = " << detResFun << endl ;  
+ // -
+ fFlowTrack->SetNhitsTRD(nHits[7]) ;                // cout << FlowDebug << "nHits TOF = " << nHits[7] << endl ;                                                   
+ fFlowTrack->SetMaxPtsTRD(nClus[7]) ;               // cout << FlowDebug << "nClus = " << nClus[7] << endl ;  
+ fFlowTrack->SetChi2TRD(nHits[7]/nClus[7]) ;        // cout << FlowDebug << "Chi2 = " << nHits[7]/nClus[7] << endl ; 
+ fFlowTrack->SetSigTRD(dEdx[2]) ;                   // cout << FlowDebug << "Dedx = " << dEdx << endl ; 
+ fFlowTrack->SetPatTRD(pAt) ;                       // cout << FlowDebug << "p = " << pAt << endl ;                                    
+ fFlowTrack->SetRespFunTRD(detResFun) ;             // cout << FlowDebug << "response function = " << detResFun << endl ;  
+ // -
+ fFlowTrack->SetNhitsTOF(nHits[8]) ;                // cout << FlowDebug << "nHits TOF = " << nHits[8] << endl ;                                                    
+ fFlowTrack->SetMaxPtsTOF(nClus[8]) ;               // cout << FlowDebug << "nClus = " << nClus[8] << endl ; 
+ fFlowTrack->SetChi2TOF(nHits[8]/nClus[8]) ;        // cout << FlowDebug << "Chi2 = " << nHits[8]/nClus[8] << endl ; 
+ fFlowTrack->SetTofTOF(dEdx[3]) ;                   // cout << FlowDebug << "Dedx = " << 1. << endl ; 
+ fFlowTrack->SetPatTOF(pAt) ;                       // cout << FlowDebug << "p = " << pAt << endl ;                                    
+ fFlowTrack->SetRespFunTOF(detResFun) ;             // cout << FlowDebug << "response function = " << detResFun << endl ;  
+ return fFlowTrack ; 
+}
+//-----------------------------------------------------------------------
+AliFlowV0* AliFlowKineMaker::FillFlowV0(TParticle* fParticle)
+{
+ // From a neutral TParticle (input) fills the AliFlowV0 (output) .
+
+ TString name = "" ; name += fPartNumber ;
+ Int_t idx = fFlowEvent->V0Collection()->GetEntries() ;
+ fFlowV0 = (AliFlowV0*)(fFlowEvent->V0Collection()->New(idx)) ;
+ fFlowV0->SetName(name.Data()) ;
+ // cout << " -v0- " << name.Data() << "(" << idx << ")"  << endl ;
+
+ // TParticle label (link: KineTree-ESD)
+ Int_t label = TMath::Abs(fPartNumber);
+ fFlowV0->SetLabel(label) ;
+
+ // reconstructed position of the V0 
+ Double_t xyz[3] ;             
+ xyz[0] = fParticle->Vx() ; 
+ xyz[1] = fParticle->Vy() ; 
+ xyz[2] = fParticle->Vz() ; 
+ fFlowV0->SetCrossPoint(xyz[0],xyz[1],xyz[2]) ;
+
+ // V0's impact parameter & error (chi2 , DCA , sigma , pointing angle)  
+ fFlowV0->SetDca((Float_t)Norm(xyz)) ;
+ fFlowV0->SetSigma(0.) ;
+ fFlowV0->SetCosPointingAngle(1.) ;  
+ fFlowV0->SetDaughtersDca(0.) ;
+ fFlowV0->SetChi2(0.) ;
+
+ // reconstructed momentum of the V0
+ Double_t pxyz[3] ;
+ pxyz[0] = fParticle->Px() ; pxyz[1] = fParticle->Py() ; pxyz[2] = fParticle->Pz() ;                   
+ Float_t phi = (Float_t)Phi(pxyz) ; if(phi<0) { phi += 2*TMath::Pi() ; }
+ fFlowV0->SetPhi(phi) ;                
+ Float_t pt = (Float_t)Pt(pxyz) ; 
+ fFlowV0->SetPt(pt) ;          
+ Float_t eta = (Float_t)Eta(pxyz) ; 
+ fFlowV0->SetEta(eta) ;        
+
+ // P.id. 
+ Int_t pdgCode = fParticle->GetPdgCode() ;
+ fFlowV0->SetMostLikelihoodPID(pdgCode);       
+
+ // mass 
+ fFlowV0->SetVmass((Float_t)fParticle->GetMass()) ; 
+ // daughters (should be taken on puprose from the KineTree, and wrote into the flow event)
+ Int_t nDaughters = fParticle->GetNDaughters() ;
+ Int_t d1 = fParticle->GetDaughter(nDaughters-1) ;
+ Int_t d2 = fParticle->GetDaughter(nDaughters-2) ;
+//
+// AliFlowTrack* pos = (AliFlowTrack*)fFlowEvent->TrackCollection()->At(pN) ; 
+// AliFlowTrack* neg = (AliFlowTrack*)fFlowEvent->TrackCollection()->At(nN) ; 
+// fFlowV0->SetDaughters(pos,neg) ;
+//
+ // d1 + d2 ; // warning: statement with no effect :)
+
+ return fFlowV0 ;
+}
+//-----------------------------------------------------------------------
+Bool_t AliFlowKineMaker::CheckTrack(TParticle* fParticle) const
+{
+ // applies track cuts (pE , eta , label)
+
+ Float_t eta   = (Float_t)fParticle->Eta() ;
+ Float_t pE    = (Float_t)fParticle->P() ;
+ Int_t   label = -1 ;   // check how to assign this !
+ Bool_t  prim  = fParticle->IsPrimary() ;
+ if(fAbsEta && (eta > fAbsEta))                                         { return kFALSE ; }
+ if((fElow < fEup) && ((pE<fElow) || (pE>fEup)))                        { return kFALSE ; }
+ if((fLabel[0] < fLabel[1]) && ((label<fLabel[0]) || (label>fLabel[1]))) { return kFALSE ; }
+ if(fPrimary && !prim)                                                          { return kFALSE ; }
+
+ return kTRUE ;
+}
+//-----------------------------------------------------------------------
+Bool_t AliFlowKineMaker::CheckEvent(TTree* fKTree) const
+{
+ // applies event cuts (dummy)
+ if(!fKTree) { return kFALSE ; }
+
+ return kTRUE ;
+}
+//-----------------------------------------------------------------------
+void AliFlowKineMaker::PrintCutList()
+{ 
+ // Prints the list of Cuts 
+
+ cout << " * ESD cuts list * " << endl ; 
+ if(fLabel[0]<fLabel[1])
+ { 
+  cout << "  Track Label [ " << fLabel[0] << " , " << fLabel[1] << " ] " << endl ; 
+ }
+ if(fAbsEta)
+ { 
+  cout << "  |eta| < " << fAbsEta << endl ; 
+ }
+ if(fElow<fEup)
+ { 
+  cout << "  Track Energy (P_total) [ " << fElow << " , " << fEup << " ] " << endl ; 
+ }
+ if(fPrimary)                                                   
+ { 
+  cout << "  Primary Particles " << endl ;
+ }
+ cout << " *               * " << endl ;
+}
+//-----------------------------------------------------------------------
+
+//-----------------------------------------------------------------------
+//*** USEFULL METHODS for a 3-array of double (~ TVector3) ***
+//-----------------------------------------------------------------------
+Double_t AliFlowKineMaker::Norm(Double_t nu[3])
+{ 
+ // returns the norm of a double[3] 
+
+ Double_t norm2 = nu[0]*nu[0] + nu[1]*nu[1] + nu[2]*nu[2] ;
+ return TMath::Sqrt(norm2) ; 
+}
+//-----------------------------------------------------------------------
+Double_t AliFlowKineMaker::Phi(Double_t nu[3])
+{
+ // returns the azimuthal angle of a double[3] 
+
+ if(nu[0]==0 && nu[1]==0) { return 0. ; }
+ else                    { return TMath::ATan2(nu[1],nu[0]) ; }
+}
+//-----------------------------------------------------------------------
+Double_t AliFlowKineMaker::Pt(Double_t nu[3])
+{
+ // returns the transvers momentum of a double[3] 
+
+ Double_t trans = nu[0]*nu[0] + nu[1]*nu[1] ;
+ return TMath::Sqrt(trans) ; 
+}
+//-----------------------------------------------------------------------
+Double_t AliFlowKineMaker::Eta(Double_t nu[3])
+{
+ // returns the PseudoRapidity of a double[3] 
+ // if transvers momentum = 0 --> returns +/- 1.000
+
+ Double_t m = Norm(nu) ;
+ if(nu[0]!=0 || nu[1]!=0) { return 0.5*TMath::Log((m+nu[2])/(m-nu[2])) ; }
+ else                    { return TMath::Sign((Double_t)1000.,nu[2]) ; }
+}
+//-----------------------------------------------------------------------
+
+#endif
+
+//////////////////////////////////////////////////////////////////////
+// - one way to open the alice KineTree -
+//////////////////////////////////////////////////////////////////////
+//
+//   TString fileName = "galice.root" ;
+//   rl = AliRunLoader::Open(fileName.Data(),"MyEvent","read");
+//   rl->LoadgAlice();
+//   gAlice = rl->GetAliRun();
+//   rl->LoadHeader();
+//   rl->LoadKinematics();
+//   fNumberOfEvents = rl->GetNumberOfEvents() ;
+//
+//   Int_t exitStatus = rl->GetEvent(getEv) ; if(exitStatus!=0) { return kFALSE ; }
+//
+//   TTree* pKTree = (TTree*)rl->TreeK();        // Particles' TTree (KineTree)
+//   AliStack* pStack = gAlice->Stack();         // Particles' Stack - "Label()" to get the number in the stack 
+//
+// // else if(rl)      // opens files one by one (unload and reload)
+// // {
+// //  rl->UnloadgAlice() ;
+// //  rl->UnloadHeader() ;
+// //  rl->UnloadKinematics() ;
+// //  delete rl ; rl = 0 ;
+// // }
+//
+//   fNumberOfParticles = pKTree->GetEntries() ;
+//   nPart = pStack->GetNtrack() ;
+//   cout << " Event : " << evtN << " :  particles : " << fNumberOfParticles << "  (stack: " << nPart << ") . " << endl ; }
+//
+//////////////////////////////////////////////////////////////////////
+
diff --git a/PWG2/FLOW/AliFlowKineMaker.h b/PWG2/FLOW/AliFlowKineMaker.h
new file mode 100644 (file)
index 0000000..f8fd38e
--- /dev/null
@@ -0,0 +1,131 @@
+//////////////////////////////////////////////////////////////////////
+//
+// $Id$
+//
+// Author: Emanuele Simili
+//
+//////////////////////////////////////////////////////////////////////
+//
+// Description: parser class from KineTree (pure MC simulation) to 
+//  AliFlowEvent . Nothing else to say at this point, but 3 lines of 
+//  comments (at least) must be there to make the code checker happy.
+//
+//////////////////////////////////////////////////////////////////////
+
+#ifndef ALIFLOWKINEMAKER_H
+#define ALIFLOWKINEMAKER_H
+
+class AliFlowEvent ;
+class AliFlowTrack ;
+class AliFlowV0 ;
+class TClonesArray ;
+
+class TParticle ;
+class TParticlePDG ;
+class TTree ;
+//class AliRun ;
+//class AliRunLoader ;
+//class AliStack ;
+
+class AliFlowKineMaker { 
+
+  public:
+  
+    AliFlowKineMaker();
+    virtual ~AliFlowKineMaker();
+
+  // FLOW SPECIFIC METHODS (to fill the flowEvents)
+    AliFlowEvent*  FillFlowEvent(TTree* fKTree) ;         // fills up evt quantities 
+    AliFlowTrack*  FillFlowTrack(TParticle* fParticle) ;   // fills up track quantities ; 
+    AliFlowV0*    FillFlowV0(TParticle* fParticle) ;      // fills up v0 quantities ;
+
+  // USEFULL METHODS
+    Double_t      Norm(Double_t nu[3]) ;                   // norm of a non-vector 3 array      
+    Double_t      Phi(Double_t nu[3]) ;                    // phi of a non-vector 3 array       
+    Double_t      Pt(Double_t nu[3]) ;                     // pt of a non-vector 3 array       
+    Double_t      Eta(Double_t nu[3]) ;                    // eta of a non-vector 3 array       
+
+  // Cut METHODS
+    Bool_t         CheckTrack(TParticle* fParticle) const ; // checks the particle (applies particle cuts, returns the particle charge)
+    Bool_t         CheckEvent(TTree* fKTree) const ;       // checks the KineTree (dummy)
+    void          PrintCutList() ;                         // prints the list of cuts
+
+    void          SetAbsEtaCut(Float_t aEta)               { fAbsEta = aEta ; }                        // exclude tracks with eta > aEta 
+    void          SetECut(Float_t eLow, Float_t eUp)       { fElow = eLow ; fEup = eUp ; }             // exclude tracks below and above .. GeV 
+    void          SetLabelCut(Int_t labLo, Int_t labHi)    { fLabel[0] = labLo ; fLabel[1] = labHi ; } // exclude tracks outside label interval
+    void          SetPrimaryCut(Bool_t prim = kTRUE)       { fPrimary = prim ; }                       // exclude secundaries 
+
+  // Get METHODS
+    Int_t          GetNgoodTracks() const                  { return fGoodTracks ; }      
+    Int_t          GetNgoodV0s() const                     { return fGoodV0s ; }         
+    Int_t          GetNgoodTracksEta() const               { return fGoodTracksEta ; } 
+    Int_t          GetNposiTracks() const                  { return fPosiTracks ; }    
+    Int_t          GetNnegaTracks() const                  { return fNegaTracks ; }    
+    Int_t          GetNunconstrained()  const              { return fUnconstrained ; } 
+    Int_t          GetBayesian(Int_t i = 2) const          { return fBayesianAll[i] ; }
+    Float_t        GetBayesianNorm(Int_t i = 2) const      { return (Float_t)fBayesianAll[i] / (Float_t)fSumAll ; }
+
+
+ protected:
+  // enumerators                           
+    Int_t            fEventNumber ;                        //! progressive enumeration of KineTree events
+    Int_t            fPartNumber ;                         //! progressive enumeration of TParticle
+
+    Int_t            fGoodTracks ;                         //! enumerator for good tracks
+    Int_t            fGoodV0s ;                            //! enumerator for good v0s
+    Int_t            fGoodTracksEta ;                      //! enumerator for good tracks in the good eta range (-0.9..0.9)
+    Int_t            fPosiTracks ;                         //! enumerator for positive tracks
+    Int_t            fNegaTracks ;                         //! enumerator for negative tracks
+    Int_t            fUnconstrained ;                      //! enumerator for tracks not constrainable
+    Int_t            fBayesianAll[6] ;                     //! final particles abundance -> AliFlowEvent (see Bayesian P.Id.)
+    Int_t            fSumAll ;                             //! total particles abundance (all kind)
+
+    Int_t            fCutEvts ;                            //! total enumerator for discarded events
+    Int_t            fCutParts ;                           //! total enumerator for discarded particles
+
+  // Flags
+    Bool_t           fNewAli ;                             //! enables the new features (since AliRoot 12/2006) 
+    Bool_t           fLoopParts ;                          //! flag to loop over tracks 
+
+    Int_t           fCounter ;                             //! number of processed events
+
+
+ private:
+
+  // to make the code checker happy
+   AliFlowKineMaker(const AliFlowKineMaker &flowMak) ;                   // Copy Constructor (dummy)
+   AliFlowKineMaker &operator=(const AliFlowKineMaker &flowMak) ; // Assignment Operator (dummy)
+
+  // KineTree stuff
+    TTree*                  fKTree ;                               //! KineTree
+    TParticle*              fParticle ;                            //! TParticle (momentum, decay, etc.)
+    TParticlePDG*    fParticlePDG ;                                //! TParticlePDG (type, charge, etc.)
+    Int_t            fCharge ;                             //! charge of the TParticlePDG
+    Float_t          fVertex[3] ;                                  //! primary vertex  
+//    AliStack*        fStack ;                            //! particle stack
+//    AliRunLoader*    fRunLoader ;                        //! AliRunLoader
+//    AliRun*          gAlice ;                                    //! pointer to the AliRun (gAlice)
+
+    Int_t            fRunID;                               //! last run ID
+    Int_t            fNumberOfEvents ;                     //! total number of KineTree events in file
+    Int_t            fNumberOfParticles ;                  //! total number of TParticles in the current event
+    Float_t         fMagField ;                            //! magnetic field from the ESD
+
+  // Flow
+    AliFlowEvent*    fFlowEvent ;                          //! pointer to flow event
+    AliFlowTrack*    fFlowTrack;                           //! pointer to flow track
+    AliFlowV0*       fFlowV0;                                      //! pointer to flow V0
+
+  // Tracks cuts
+    Float_t fAbsEta;                                       // exclude tracks with |eta| bigger than 
+    Float_t fElow ;                                        // exclude tracks below .. GeV (~total Momentum)
+    Float_t fEup ;                                         // exclude tracks above .. GeV (~total Momentum)
+    Int_t   fLabel[2] ;                                            // exclude tracks outside label interval
+    Bool_t  fPrimary ;                                     // exclude secundary tracks 
+    
+  ClassDef(AliFlowKineMaker,1);
+};                     
+
+#endif
+
index 59859df..8e5b526 100644 (file)
@@ -45,7 +45,9 @@ using namespace std; //required for resolving the 'cout' symbol
 
 ClassImp(AliFlowMaker) 
 //-----------------------------------------------------------------------
-AliFlowMaker::AliFlowMaker()
+AliFlowMaker::AliFlowMaker():
+  fESD(0x0), fTrack(0x0), fV0(0x0), fVertex(0x0),
+  fFlowEvent(0x0), fFlowTrack(0x0), fFlowV0(0x0)
 {
  // default constructor 
  // resets counters , sets defaults
@@ -262,6 +264,10 @@ AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
  Int_t idxr[130] ; // used for Cluster Map ( see AliESDtrack::GetTRDclusters() )    // old:    UInt
  Int_t nClus = 0 ;     
  Int_t fNFound = 0 ;                                   // *!* fNFoundable (in AliTPCtrack) ... added by M.Ianov 
+ // -
+ Double_t detPid[5] ;          
+ Float_t  detPid6[AliFlowConstants::kPid] ; 
+
  Double_t pVecAt[3] ; 
  for(Int_t gg=0;gg<3;gg++) { pVecAt[gg] = gD[gg] ; }
  Bool_t boh ; Float_t pAt = 0 ;                        // to get p at each detector
@@ -271,6 +277,12 @@ AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
  nClus = fTrack->GetTPCclusters(idXt) ;
  fNFound = fTrack->GetTPCNclsF() ;  // was 160
+ if( (fTrack->GetStatus() & AliESDtrack::kTPCpid) != 0 )
+ {
+  fTrack->GetTPCpid(detPid) ; 
+  for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
+  fFlowTrack->SetRespFunTPC(detPid6) ; 
+ }
  fFlowTrack->SetMaxPtsTPC(fNFound) ;                           
  fFlowTrack->SetFitPtsTPC(nClus) ;                             
  fFlowTrack->SetDedxTPC(fTrack->GetTPCsignal()) ;              
@@ -281,7 +293,13 @@ AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
  else       { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgITSx, fMagField, pVecAt) ; }
  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
  nClus = fTrack->GetITSclusters(idX) ;
- fNFound = 6 ;
+ fNFound = 6 ; // ? fixed
+ if( (fTrack->GetStatus() & AliESDtrack::kITSpid) != 0 )
+ {
+  fTrack->GetITSpid(detPid) ; 
+  for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
+  fFlowTrack->SetRespFunITS(detPid6) ; 
+ } 
  fFlowTrack->SetMaxPtsITS(fNFound) ;                           
  fFlowTrack->SetFitPtsITS(nClus) ;                             
  fFlowTrack->SetDedxITS(fTrack->GetITSsignal()) ;              
@@ -291,8 +309,14 @@ AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
  if(fNewAli) { boh = fTrack->GetPxPyPzAt(AliFlowConstants::fgTRDx, fMagField, pVecAt) ; } 
  else       { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgTRDx, fMagField, pVecAt) ; } 
  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
- fNFound = fTrack->GetTRDncls() ;  // was 130
  nClus = fTrack->GetTRDclusters(idxr) ;
+ fNFound = fTrack->GetTRDncls() ;  // was 130
+ if( (fTrack->GetStatus() & AliESDtrack::kTRDpid) != 0 )
+ {
+  fTrack->GetTRDpid(detPid) ; 
+  for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
+  fFlowTrack->SetRespFunTRD(detPid6) ; 
+ }
  fFlowTrack->SetMaxPtsTRD(fNFound) ;                           
  fFlowTrack->SetNhitsTRD(nClus) ;                              
  fFlowTrack->SetSigTRD(fTrack->GetTRDsignal()) ;               
@@ -302,8 +326,14 @@ AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
  if(fNewAli) { boh = fTrack->GetPxPyPzAt(AliFlowConstants::fgTOFx, fMagField, pVecAt) ; }
  else       { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgTOFx, fMagField, pVecAt) ; }
  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
- fNFound = 0 ; if(fTrack->GetTOFCalChannel() > 0) { fNFound = 1 ; }
  nClus = fTrack->GetTOFcluster() ;
+ fNFound = 0 ; if(fTrack->GetTOFCalChannel() > 0) { fNFound = 1 ; }
+ if( (fTrack->GetStatus() & AliESDtrack::kTOFpid) != 0 )
+ {
+  fTrack->GetTOFpid(detPid) ; 
+  for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
+  fFlowTrack->SetRespFunTOF(detPid6) ; 
+ }
  fFlowTrack->SetMaxPtsTOF(fNFound) ;                           
  fFlowTrack->SetNhitsTOF(nClus) ;                              
  fFlowTrack->SetTofTOF(fTrack->GetTOFsignal()) ;               
@@ -398,7 +428,7 @@ AliFlowV0* AliFlowMaker::FillFlowV0(AliESDv0* fV0)
  //fFlowV0->SetDca((Float_t)fV0->GetD()) ;    // GetDistNorm 
  //fFlowV0->SetSigma((Float_t)fV0->GetDistSigma()) ;
  //fFlowV0->SetCosPointingAngle((Float_t)fV0->GetV0CosineOfPointingAngle()) ;  
- //fFlowV0->SetDaughtersDCA(fV0->GetDcaV0Daughters()) ;
+ //fFlowV0->SetDaughtersDca(fV0->GetDcaV0Daughters()) ;
  //fFlowV0->SetChi2((Float_t)fV0->GetChi2V0()) ;     // AliRoot v4-04-Release (December 2006)  
  //fFlowV0->SetChi2((Float_t)fV0->GetChi2()) ;       // AliRoot v4-04-Release (old)
  // ...when they'll stop changing the methods I'll enable the above lines. For now:
index 2edba63..d94486a 100644 (file)
@@ -7,8 +7,8 @@
 //////////////////////////////////////////////////////////////////////
 //
 // Description: parser class from AliESD to AliFlowEvent . 
-//  Nothing else to say at thispoint, but 3 lines of comments must be 
-//  there to make the code checker happy. I hope now is ok!
+//  Nothing else to say at this point, but 3 lines of comments must be 
+//  there (at least) to make the code checker happy. I hope now is ok!
 //
 //////////////////////////////////////////////////////////////////////
 
@@ -28,6 +28,7 @@ class AliESDv0 ;
 class AliFlowMaker  {
 
   public:
+  
     AliFlowMaker();
     virtual ~AliFlowMaker();
 
@@ -95,7 +96,7 @@ class AliFlowMaker  {
    AliFlowMaker &operator=(const AliFlowMaker &flowMak) ;   // Assignment Operator (dummy)
 
   // ESDs
-    AliESD*          fESD;                                 //! "ESD" branch in fChain
+    AliESD*          fESD;                                 //! "ESD event"
     AliESDtrack*     fTrack;                               //! "ESD track" 
     AliESDv0*        fV0;                                  //! "ESD v0" 
     AliESDVertex*    fVertex;                              //! "ESD primary vertex"  
index 520900d..6dd66a1 100644 (file)
@@ -167,7 +167,7 @@ Bool_t AliFlowSelection::SelectPart(AliFlowTrack* pFlowTrack) const
  } 
  
  // PID probability
- float pidProb = pFlowTrack->MostLikelihoodProb() ;
+ float pidProb = pFlowTrack->MostLikelihoodRespFunc() ;
  if(fPidProbPart[1] > fPidProbPart[0] &&  (pidProb < fPidProbPart[0] || pidProb > fPidProbPart[1])) return kFALSE;
  
  // Constrainable
index 45adf8b..1866619 100644 (file)
@@ -175,7 +175,6 @@ class AliFlowSelection : public TObject {
   static Int_t    fgTPChits[AliFlowConstants::kSels];          // minimum number of TPC hits
   static Bool_t   fgConstrainable;                             // cut un-constrainable tracks 
 
-
   ClassDef(AliFlowSelection,1)                                 // macro for rootcint
 }; 
 
index 0f96a4e..b338c42 100644 (file)
 // Description: 
 //         an array of AliFlowTrack is the core of the AliFlowEvent. 
 // The object AliFlowTrack contains data members wich summarize the track 
-// information most useful for flow study (such as Pt, eta, phi angle, 
-// p.id hypothesis, some detector informations like fitpoints, chi2, 
-// energy loss, t.o.f., ecc.). 
+// information most useful for flow study (such as Pt, eta, phi angle), some 
+// detector signals (fitpoints, chi2 of the track fit, energy loss, time of 
+// flight, ecc.), the p.id is stored as the detector response function, for 
+// each detector (with the possibility to costomly combine them), and the 
+// combined one (see bayesian P.Id. - chap.5 of the ALICE PPR). 
 // This class is optimized for reaction plane calculation and sub-event 
 // selection, through the appropriate methods in AliFlowEvent.  
 // Two arrays of flags in the AliFlowTrack object (fSelection[][] and 
@@ -53,13 +55,13 @@ AliFlowTrack::AliFlowTrack()
  fZLastPoint = 0. ;                          
  fTrackLength = 0. ;
  fMostLikelihoodPID = 0 ;                    
- for(Int_t ii=0;ii<2;ii++)                     { fDcaSigned[ii] = 0. ; }                             
- for(Int_t ii=0;ii<AliFlowConstants::kPid;ii++) { fPidProb[ii]   = 0. ; }
- for(Int_t ii=0;ii<4;ii++)
+ for(Int_t dd=0;dd<2;dd++)                     { fDcaSigned[dd] = 0. ; }                             
+ for(Int_t ii=0;ii<AliFlowConstants::kPid;ii++) { fCombRespFun[ii]   = 0. ; }
+ for(Int_t det=0;det<4;det++)
  {
-  fFitPts[4]  = 0  ; fMaxPts[4]  = 0  ; 
-  fFitChi2[4] = 0. ; fDedx[4]    = 0. ;
-  fMom[4] = 0. ;
+  fFitPts[det]  = 0  ; fMaxPts[det]   = 0  ; 
+  fFitChi2[det] = 0. ; fDedx[det]     = 0. ; fMom[det] = 0. ; 
+  for(Int_t ii=0;ii<AliFlowConstants::kPid;ii++) { fRespFun[det][ii] = -1. ; }
  }
  ResetSelection() ;                          
 }
@@ -208,65 +210,129 @@ const Char_t* AliFlowTrack::Pid() const
  return pId.Data() ; 
 }
 //-------------------------------------------------------------
-void AliFlowTrack::PidProbs(Float_t pidN[AliFlowConstants::kPid]) const 
+void AliFlowTrack::PidProbs(Float_t *pidN) const 
 { 
- // Returns the normalized probability for the given track to be [e,mu,pi,k,p,d] 
- // The detector response is weighted by the bayesian vector of particles 
- // abundances, stored in AliFlowConstants::fgBayesian[] .
+ // Returns the normalized probability (the "bayesian weights") for the given track to be [e,mu,pi,k,p,d] .  
+ // The COMBINED detector response function is scaled by the a priori probabilities 
+ // (the normalised particle abundances is stored in AliFlowConstants::fgBayesian[]) .
 
  Double_t sum = 0 ; 
- for(Int_t n=0;n<AliFlowConstants::kPid;n++)  { sum += fPidProb[n] * AliFlowConstants::fgBayesian[n] ; }
- if(sum)
- {
-  for(Int_t n=0;n<AliFlowConstants::kPid;n++) { pidN[n] = fPidProb[n] * AliFlowConstants::fgBayesian[n] / sum ; }
+ for(Int_t n=0;n<AliFlowConstants::kPid;n++) { pidN[n] = PidProb(n) ; sum += pidN[n] ; }
+ if(sum) 
+ {  
+  for(Int_t n=0;n<AliFlowConstants::kPid;n++)  { pidN[n] /=  sum ; }
  }
  else { cout << " ERROR - Empty Bayesian Vector !!! " << endl ; }
 } 
 //-------------------------------------------------------------
-Float_t  AliFlowTrack::PidProb(Int_t nn)       const
+void AliFlowTrack::PidProbsC(Float_t *pidN) const 
+{ 
+ // Returns the normalized probability (the "bayesian weights") for the given track to be [e,mu,pi,k,p,d] .  
+ // The CUSTOM detector response function (see AliFlowTrack) is scaled by the a priori probabilities 
+ // (the normalised particle abundances is stored in AliFlowConstants::fgBayesian[]) .
+ Double_t sum = 0 ; 
+ for(Int_t n=0;n<AliFlowConstants::kPid;n++) { pidN[n] = PidProbC(n) ; sum += pidN[n] ; }
+ if(sum) 
+ {  
+  for(Int_t n=0;n<AliFlowConstants::kPid;n++)  { pidN[n] /=  sum ; }
+ }
+ else { cout << " ERROR - Empty Bayesian Vector !!! " << endl ; }
+} 
+//-------------------------------------------------------------
+Float_t  AliFlowTrack::PidProb(Int_t nPid) const
 {
- // Returns the normalized probability of the track to be [nn] (e,mu,pi,k,pi,d).
+ // Returns the bayesian weight of the track to be [nPid = e,mu,pi,k,pi,d].
+ // The detector response function in use is the combined one (from the ESD)
+
+ if(nPid > AliFlowConstants::kPid) { return 0. ; } 
 
- Float_t pidN[AliFlowConstants::kPid] ; 
- PidProbs(pidN) ;
- return pidN[nn] ;
+ return (fCombRespFun[nPid] * AliFlowConstants::fgBayesian[nPid]) ;
 }
 //-------------------------------------------------------------
-TVector AliFlowTrack::PidProbs()               const
+Float_t  AliFlowTrack::PidProbC(Int_t nPid) const
 {
- // Returns the normalized probability for the given track to be [e,mu,pi,k,p,d] 
- // as a TVector.
+ // Returns the bayesian weight of the track to be [nPid = e,mu,pi,k,pi,d].
+ // The detector response function in use is the custom one ...
 
- TVector pidNvec(AliFlowConstants::kPid) ;
- Float_t pidN[AliFlowConstants::kPid] ; 
- PidProbs(pidN) ;
- for(Int_t n=0;n<AliFlowConstants::kPid;n++)  { pidNvec[n] = pidN[n] ; }
+ if(nPid > AliFlowConstants::kPid) { return 0. ; } 
 
- return pidNvec ;
-}
-//-------------------------------------------------------------
-void AliFlowTrack::RawPidProbs(Float_t pidV[AliFlowConstants::kPid]) const 
-{ 
- // Returns the array of probabilities for the track to be [e,mu,pi,k,pi,d].
+ Float_t  customRespFun[AliFlowConstants::kPid] ;
+ GetCustomRespFun(customRespFun) ;
 
- for(Int_t ii=0;ii<AliFlowConstants::kPid;ii++) { pidV[ii] = fPidProb[ii] ; }
-} 
+ return (customRespFun[nPid] * AliFlowConstants::fgBayesian[nPid]) ;
+}
 //-------------------------------------------------------------
-Float_t AliFlowTrack::MostLikelihoodProb()         const 
+Float_t AliFlowTrack::MostLikelihoodRespFunc() const 
 { 
- // Returns the probability of the most probable P.id.
- // (Warning: THIS IS NOT WEIGHTED IN THE BAYESIAN WAY...) 
+ // Returns the detector response function for the most probable P.id. hypothesis
+ // (Warning: THIS IS NOT THE BAYESIAN WEIGHT !) 
 
  Int_t pdgCode = TMath::Abs(MostLikelihoodPID()) ;
- if(pdgCode == 11)            { return fPidProb[0] ; }
- else if(pdgCode == 13)       { return fPidProb[1] ; }
- else if(pdgCode == 211)      { return fPidProb[2] ; }
- else if(pdgCode == 321)      { return fPidProb[3] ; }
- else if(pdgCode == 2212)     { return fPidProb[4] ; }
- else if(pdgCode == 10010020) { return fPidProb[5] ; }
+ if(pdgCode == 11)            { return fCombRespFun[0] ; }
+ else if(pdgCode == 13)       { return fCombRespFun[1] ; }
+ else if(pdgCode == 211)      { return fCombRespFun[2] ; }
+ else if(pdgCode == 321)      { return fCombRespFun[3] ; }
+ else if(pdgCode == 2212)     { return fCombRespFun[4] ; }
+ else if(pdgCode == 10010020) { return fCombRespFun[5] ; }
  else { return 0. ; }
 } 
 //-------------------------------------------------------------
+void AliFlowTrack::SetRespFun(Int_t det, Float_t *r)
+{
+ // This function fills "AliFlowConstants::kPid" PID weights 
+ // (detector response functions) of the track .
+ // The method is private, cause it is called by the public 
+ // methods: SetRespFunITS, SetRespFunTPC, ...
+
+ for(Int_t i=0;i<AliFlowConstants::kPid;i++)
+ {
+  fRespFun[det][i] = r[i] ;
+ }
+}
+//-------------------------------------------------------------
+void AliFlowTrack::GetRespFun(Int_t det, Float_t *r) const
+{
+ // This function returns the response functions of the 
+ // detector [det] for the P.Id. of the track .
+ // The method is private, and it is called by the public 
+ // methods: GetRespFunITS, GetRespFunTPC, ...
+
+ for(Int_t i=0;i<AliFlowConstants::kPid;i++)
+ {
+  r[i] = fRespFun[det][i] ;
+ }
+}
+//-------------------------------------------------------------
+void AliFlowTrack::GetCombinedRespFun(Float_t *r) const 
+{
+ // This function returns the combined response functions for  
+ // [e,mu,pi,k,pi,d] as it is stored in the ESD .
+ for(Int_t i=0;i<AliFlowConstants::kPid;i++)
+ {
+  r[i] = fCombRespFun[i] ;
+ }
+}
+//------------------------------------------------------------- 
+void AliFlowTrack::GetCustomRespFun(Float_t *r) const
+{
+ // This function returns a combined response functions as setted 
+ // by the user ... at the moment just the sum of single responses ...  
+ // for the track to be [e,mu,pi,k,pi,d] .
+ Int_t sum ;
+ for(Int_t i=0;i<AliFlowConstants::kPid;i++)
+ {
+  r[i] = 0. ; sum = 0 ;
+  for(Int_t det=0;det<4;det++) 
+  { 
+   if(fRespFun[det][i]>0) { r[i] += fRespFun[det][i] ; sum += 1 ; } 
+  }
+  if(sum) { r[i] /= sum ; }
+ }
+}
+//-------------------------------------------------------------
 void AliFlowTrack::SetPid(const Char_t* pid)           
 { 
  // Sets the P.Id. hypotesis of the track from a String imput (that should be 
@@ -285,6 +351,14 @@ void AliFlowTrack::SetPid(const Char_t* pid)
  else if(strchr(pid,'-'))  { fMostLikelihoodPID = TMath::Abs(fMostLikelihoodPID) * -1 ; } 
 }
 //-------------------------------------------------------------
+void AliFlowTrack::SetPid(Int_t pdgCode)               
+{ 
+ // Sets the P.Id. hypotesis of the track from the PDG code. 
+ // Sign can be given as well. 
+ fMostLikelihoodPID = pdgCode ;
+}
+//-------------------------------------------------------------
 Bool_t AliFlowTrack::IsConstrainable()             const 
 { 
  // Returns kTRUE if the track is constrainable. 
@@ -370,8 +444,9 @@ void AliFlowTrack::ResetSelection()
 //-------------------------------------------------------------
 void AliFlowTrack::SetConstrainable()                        
 { 
- // fills the constrained parameters with the unconstrained ones, making it a constrainable track.
- //                                                   !!! TRICKY METHOD !!!
+ // fills the constrained parameters with the unconstrained ones,
+ // making the track a constrainable track.
+ //                                       !!! TRICKY METHOD !!!
 
  if(!IsConstrainable()) 
  { 
@@ -383,8 +458,9 @@ void AliFlowTrack::SetConstrainable()
 //-------------------------------------------------------------
 void AliFlowTrack::SetUnConstrainable()                              
 { 
- // deletes the constrained parameters making it  an unconstrainable track. 
- //                                                   !!! TRICKY METHOD !!!
+ // deletes the constrained parameters making the track an 
+ // unconstrainable track. 
+ //                                       !!! TRICKY METHOD !!!
 
  if(IsConstrainable()) 
  { 
index 589a5ec..7341e95 100644 (file)
@@ -65,43 +65,51 @@ public:
   Float_t       ZedDcaError()        const { return TMath::Abs(fDcaError[1]) ; }
 
  // Gets - P.id.
-  Float_t       MostLikelihoodProb()           const;
-  Float_t       PidProb(Int_t nn)              const;                   // normalized detector response (using AliFlowConstants::fgBayesian[] )
-  TVector      PidProbs()                      const; 
-  void          PidProbs(Float_t pidN[AliFlowConstants::kPid])    const; 
-  void          RawPidProbs(Float_t pidV[AliFlowConstants::kPid]) const; // raw ESD detector responses
-  Float_t       ElectronPositronProb()         const { return PidProb(0) ; } 
-  Float_t       MuonPlusMinusProb()            const { return PidProb(1) ; }
-  Float_t       PionPlusMinusProb()            const { return PidProb(2) ; }
-  Float_t       KaonPlusMinusProb()            const { return PidProb(3) ; }
-  Float_t       ProtonPbarProb()               const { return PidProb(4) ; }
-  Float_t       DeuteriumAntiDeuteriumProb()   const { return PidProb(5) ; }
-  Int_t         MostLikelihoodPID()            const { return fMostLikelihoodPID ; }
+  Float_t       MostLikelihoodRespFunc()       const; // detector Response Function for the most probable P.id. hypothesis
+  Float_t       PidProb(Int_t nPid)           const; // normalized probabilities or bayesian weights (using COMBINED response function & "a priori" AliFlowConstants::fgBayesian[])
+  void          PidProbs(Float_t *pidN)        const; // ... same for all the 6 hypothesis .
+  Float_t       PidProbC(Int_t nPid)                  const; // normalized probabilities or bayesian weights (using CUSTOM response function & "a priori" AliFlowConstants::fgBayesian[])
+  void          PidProbsC(Float_t *pidN)       const; // ... same for all the 6 hypothesis . 
+  Float_t       ElectronPositronProb()         const { return PidProb(0) ; } 
+  Float_t       MuonPlusMinusProb()            const { return PidProb(1) ; }
+  Float_t       PionPlusMinusProb()            const { return PidProb(2) ; }
+  Float_t       KaonPlusMinusProb()            const { return PidProb(3) ; }
+  Float_t       ProtonPbarProb()               const { return PidProb(4) ; }
+  Float_t       DeuteriumAntiDeuteriumProb()   const { return PidProb(5) ; }
+  Int_t         MostLikelihoodPID()            const { return fMostLikelihoodPID ; }
 
  // Gets - Detector response
-  Float_t       Chi2TPC()     const { return fFitChi2[0] ; }
-  Float_t       Chi2ITS()     const { return fFitChi2[1] ; }
-  Float_t       Chi2TRD()     const { return fFitChi2[2] ; }
-  Float_t       Chi2TOF()     const { return fFitChi2[3] ; }
-  Int_t         FitPtsTPC()   const { return fFitPts[0]  ; }
-  Int_t         FitPtsITS()   const { return fFitPts[1]  ; }
-  Int_t         NhitsTRD()    const { return fFitPts[2]  ; }
-  Int_t         NhitsTOF()    const { return fFitPts[3]  ; }
-  Int_t         MaxPtsTPC()   const { return fMaxPts[0]  ; }
-  Int_t         MaxPtsITS()   const { return fMaxPts[1]  ; }
-  Int_t         MaxPtsTRD()   const { return fMaxPts[2]  ; }
-  Int_t         MaxPtsTOF()   const { return fMaxPts[3]  ; }
-  Float_t       DedxTPC()     const { return fDedx[0]   ; }
-  Float_t       DedxITS()     const { return fDedx[1]   ; }
-  Float_t       SigTRD()      const { return fDedx[2]   ; }
-  Float_t       TofTOF()      const { return fDedx[3]   ; }
-  Float_t       PatTPC()      const { return fMom[0] ; }
-  Float_t       PatITS()      const { return fMom[1] ; }
-  Float_t       PatTRD()      const { return fMom[2] ; }
-  Float_t       PatTOF()      const { return fMom[3] ; }
+  Float_t       Chi2TPC()                     const { return fFitChi2[0] ; }
+  Float_t       Chi2ITS()                     const { return fFitChi2[1] ; }
+  Float_t       Chi2TRD()                     const { return fFitChi2[2] ; }
+  Float_t       Chi2TOF()                     const { return fFitChi2[3] ; }
+  Int_t         FitPtsTPC()                   const { return fFitPts[0]  ; }
+  Int_t         FitPtsITS()                   const { return fFitPts[1]  ; }
+  Int_t         NhitsTRD()                    const { return fFitPts[2]  ; }
+  Int_t         NhitsTOF()                    const { return fFitPts[3]  ; }
+  Int_t         MaxPtsTPC()                   const { return fMaxPts[0]  ; }
+  Int_t         MaxPtsITS()                   const { return fMaxPts[1]  ; }
+  Int_t         MaxPtsTRD()                   const { return fMaxPts[2]  ; }
+  Int_t         MaxPtsTOF()                   const { return fMaxPts[3]  ; }
+  Float_t       DedxTPC()                     const { return fDedx[0]    ; }
+  Float_t       DedxITS()                     const { return fDedx[1]    ; }
+  Float_t       SigTRD()                      const { return fDedx[2]    ; }
+  Float_t       TofTOF()                      const { return fDedx[3]    ; }
+  Float_t       PatTPC()                      const { return fMom[0] ; }
+  Float_t       PatITS()                      const { return fMom[1] ; }
+  Float_t       PatTRD()                      const { return fMom[2] ; }
+  Float_t       PatTOF()                      const { return fMom[3] ; }
+  void         GetRespFunTPC(Float_t *r)      const { GetRespFun(0,r) ; }
+  void         GetRespFunITS(Float_t *r)      const { GetRespFun(1,r) ; }
+  void         GetRespFunTRD(Float_t *r)      const { GetRespFun(2,r) ; }
+  void         GetRespFunTOF(Float_t *r)      const { GetRespFun(3,r) ; }
+
+  void          GetCombinedRespFun(Float_t *r) const ;
+  void          GetCustomRespFun(Float_t *r)   const ;
 
  // Sets - Track variables
-  void SetPid(const Char_t* pid);
+  void SetPid(const Char_t* pid) ;
+  void SetPid(Int_t pdgCode = 211) ;
   void SetCharge(Int_t charge)                           { fMostLikelihoodPID = TMath::Sign(TMath::Abs(fMostLikelihoodPID),charge) ; }
   void SetPhi(Float_t phi)                       { fPhi = phi; }
   void SetEta(Float_t eta)                       { fEta = eta; }
@@ -119,12 +127,12 @@ public:
 
  // Sets - P.id.
   void SetMostLikelihoodPID(Int_t val)           { fMostLikelihoodPID = val ; }
-  void SetElectronPositronProb(Float_t val)      { fPidProb[0] = TMath::Abs(val) ; } 
-  void SetMuonPlusMinusProb(Float_t val)         { fPidProb[1] = TMath::Abs(val) ; } 
-  void SetPionPlusMinusProb(Float_t val)         { fPidProb[2] = TMath::Abs(val) ; } 
-  void SetKaonPlusMinusProb(Float_t val)         { fPidProb[3] = TMath::Abs(val) ; } 
-  void SetProtonPbarProb(Float_t val)            { fPidProb[4] = TMath::Abs(val) ; } 
-  void SetDeuteriumAntiDeuteriumProb(Float_t val) { fPidProb[5] = TMath::Abs(val) ; }
+  void SetElectronPositronProb(Float_t val)      { fCombRespFun[0] = TMath::Abs(val) ; } 
+  void SetMuonPlusMinusProb(Float_t val)         { fCombRespFun[1] = TMath::Abs(val) ; } 
+  void SetPionPlusMinusProb(Float_t val)         { fCombRespFun[2] = TMath::Abs(val) ; } 
+  void SetKaonPlusMinusProb(Float_t val)         { fCombRespFun[3] = TMath::Abs(val) ; } 
+  void SetProtonPbarProb(Float_t val)            { fCombRespFun[4] = TMath::Abs(val) ; } 
+  void SetDeuteriumAntiDeuteriumProb(Float_t val) { fCombRespFun[5] = TMath::Abs(val) ; }
 
  // Sets - Detector response
   void SetChi2TPC(Float_t chi2)                          { fFitChi2[0] = chi2   ; }
@@ -132,21 +140,25 @@ public:
   void SetChi2TRD(Float_t chi2)                          { fFitChi2[2] = chi2   ; }
   void SetChi2TOF(Float_t chi2)                          { fFitChi2[3] = chi2   ; }
   void SetFitPtsTPC(Int_t fitPts)                { fFitPts[0]  = fitPts ; }
-  void SetFitPtsITS(Int_t nhits)                 { fFitPts[1]  = nhits ; }
+  void SetFitPtsITS(Int_t nhits)                 { fFitPts[1]  = nhits  ; }
   void SetNhitsTRD(Int_t fitPts)                 { fFitPts[2]  = fitPts ; }
-  void SetNhitsTOF(Int_t hits)                   { fFitPts[3]  = hits ; }
+  void SetNhitsTOF(Int_t hits)                   { fFitPts[3]  = hits   ; }
   void SetMaxPtsTPC(Int_t maxPts)                { fMaxPts[0]  = maxPts ; }
   void SetMaxPtsITS(Int_t maxPts)                { fMaxPts[1]  = maxPts ; }
   void SetMaxPtsTRD(Int_t maxPts)                { fMaxPts[2]  = maxPts ; }
   void SetMaxPtsTOF(Int_t maxPts)                { fMaxPts[3]  = maxPts ; }
   void SetDedxTPC(Float_t dedx)                          { fDedx[0]    = dedx   ; }
   void SetDedxITS(Float_t dedx)                          { fDedx[1]    = dedx   ; }
-  void SetSigTRD(Float_t trd)                    { fDedx[2]    = trd   ; }
-  void SetTofTOF(Float_t tof)                    { fDedx[3]    = tof   ; }
+  void SetSigTRD(Float_t trd)                    { fDedx[2]    = trd    ; }
+  void SetTofTOF(Float_t tof)                    { fDedx[3]    = tof    ; }
   void SetPatTPC(Float_t p)                      { fMom[0] = p ; }
   void SetPatITS(Float_t p)                      { fMom[1] = p ; }
   void SetPatTRD(Float_t p)                      { fMom[2] = p ; }
   void SetPatTOF(Float_t p)                      { fMom[3] = p ; }
+  void SetRespFunTPC(Float_t *r)                 { SetRespFun(0,r) ; }
+  void SetRespFunITS(Float_t *r)                 { SetRespFun(1,r) ; }
+  void SetRespFunTRD(Float_t *r)                 { SetRespFun(2,r) ; }
+  void SetRespFunTOF(Float_t *r)                 { SetRespFun(3,r) ; }
 
  // Selection methods
   void  SetSelect(Int_t harmonic, Int_t selection)                   { fSelection[harmonic][selection] = kTRUE ; }
@@ -165,34 +177,39 @@ private:
  // to make the code checker happy
   AliFlowTrack &operator=(const AliFlowTrack &flowTrack) ; // Assignment Operator (dummy)
 
+ // Used internally to set/get the detector response functions (for single detectors)
+  void     SetRespFun(Int_t det, Float_t *r) ;
+  void     GetRespFun(Int_t det, Float_t *r) const ;
+
  // Data Members
-  Float_t  fPhi;                               // azimuthal angle of the constrained track
-  Float_t  fEta;                               // pseudorapidity of the constrained track
-  Float_t  fPt;                                        // transverse momentum of the constrained track 
-  Float_t  fZFirstPoint;                       // Z position at beginning of TPC
-  Float_t  fZLastPoint;                        // Z position at end of PHOS
-  Float_t  fChi2 ;                             // constrained chi2
-  Float_t  fTrackLength;                       // lenght of the track
-  Int_t    fMostLikelihoodPID;                 // PDG code of the most probable P.Id.
-  Double_t fDcaSigned[2] ;                     // positive or negative tracks -> P changes differently including vertex ...
-  Double_t fDcaError[3] ;                      // error over the DCA
-  Float_t  fPhiGlobal;                         // azimuthal angle of the unconstrained track
-  Float_t  fEtaGlobal;                         // pseudorapidity of the unconstrained track
-  Float_t  fPtGlobal;                          // transverse momentum of the unconstrained track 
-  Int_t    fLabel ;                            // Label of the track (link: KineTree-ESD) 
+  Float_t  fPhi;                                // azimuthal angle of the constrained track
+  Float_t  fEta;                                // pseudorapidity of the constrained track
+  Float_t  fPt;                                         // transverse momentum of the constrained track 
+  Float_t  fZFirstPoint;                        // Z position at beginning of TPC
+  Float_t  fZLastPoint;                         // Z position at end of PHOS
+  Float_t  fChi2 ;                              // constrained chi2
+  Float_t  fTrackLength;                        // lenght of the track
+  Int_t    fMostLikelihoodPID;                  // PDG code of the most probable P.Id.
+  Double_t fDcaSigned[2] ;                      // positive or negative tracks -> P changes differently including vertex
+  Double_t fDcaError[3] ;                       // error over the DCA
+  Float_t  fPhiGlobal;                          // azimuthal angle of the unconstrained track
+  Float_t  fEtaGlobal;                          // pseudorapidity of the unconstrained track
+  Float_t  fPtGlobal;                           // transverse momentum of the unconstrained track 
+  Int_t    fLabel ;                             // Label of the track (link: KineTree-ESD) 
  // -
-  Float_t  fPidProb[AliFlowConstants::kPid] ;  // Array of probabilities to be (e,mu,pi,k,p,d)
-  Int_t    fFitPts[4] ;                        // Fit Points (in: TPC,ITS,TRD,TOF)
-  Int_t    fMaxPts[4] ;                        // Foundable Clusters (in: TPC,ITS,TRD,TOF)
-  Float_t  fFitChi2[4] ;                       // Chi2 (in: TPC,ITS,TRD,TOF)
-  Float_t  fDedx[4] ;                          // dE/dx from TPC and ITS , TRD signal , time of flight (TOF)
-  Float_t  fMom[4] ;                           // Track momentum at the entrance of TPC,ITS,TRD,TOF     
+  Float_t  fCombRespFun[AliFlowConstants::kPid]; // Array of probabilities to be (e,mu,pi,k,p,d)
+  Int_t    fFitPts[4] ;                         // Fit Points (in: TPC,ITS,TRD,TOF)
+  Int_t    fMaxPts[4] ;                         // Foundable Clusters (in: TPC,ITS,TRD,TOF)
+  Float_t  fFitChi2[4] ;                        // Chi2 (in: TPC,ITS,TRD,TOF)
+  Float_t  fDedx[4] ;                           // dE/dx from TPC and ITS , TRD signal , time of flight (TOF)
+  Float_t  fMom[4] ;                            // Track momentum at the entrance of TPC,ITS,TRD,TOF
+  Float_t  fRespFun[4][AliFlowConstants::kPid];  // Detector response function (single detectors)
 
  // Selection's array for R.P. & sub-event selection (not stored)  
   Bool_t   fSelection[AliFlowConstants::kHars][AliFlowConstants::kSels]; //! 
   Short_t  fSubevent[AliFlowConstants::kHars][AliFlowConstants::kSels];         //! 
 
-  ClassDef(AliFlowTrack,5) ;                   // macro for rootcint
+  ClassDef(AliFlowTrack,5) ;                    // macro for rootcint
 };
 
 #endif
index e1c5810..9572859 100644 (file)
@@ -25,7 +25,8 @@ using namespace std; //required for resolving the 'cout' symbol
 
 ClassImp(AliFlowV0) 
 //////////////////////////////////////////////////////////////////////////////
-AliFlowV0::AliFlowV0()
+AliFlowV0::AliFlowV0():
+  fDaughterP(0x0), fDaughterN(0x0)
 {
  // default constructor  
  
@@ -38,9 +39,9 @@ AliFlowV0::AliFlowV0()
  fCrossDCA = 0. ;
  fSigma = 1. ;
  fLabel = 0 ;
- for(Int_t dd=0;dd<3;dd++) { fCrossPoint[dd] = 0. ; }
- for(Int_t dd=0;dd<2;dd++) { fDaughters[dd] = 0  ; }
  fPointAngle = 0. ;
+ for(Int_t dd=0;dd<3;dd++) { fCrossPoint[dd] = 0. ; }
+ // fDaughterP = 0  ; fDaughterN = 0  ;
 }
 //////////////////////////////////////////////////////////////////////////////
 AliFlowV0::AliFlowV0(const Char_t* name) 
index 9ec46d5..4a79852 100644 (file)
@@ -44,8 +44,8 @@ public:
   Float_t       Dca()          const { return fDca ; }
   Float_t       Chi2()         const { return fChi2; }
   Float_t       Mass()         const { return fMass ; }
-  AliFlowTrack* DaughterP()    const { return fDaughters[0] ; } 
-  AliFlowTrack* DaughterN()    const { return fDaughters[1] ; } 
+  AliFlowTrack* DaughterP()    const { return fDaughterP ; } 
+  AliFlowTrack* DaughterN()    const { return fDaughterN ; } 
   Float_t       V0Lenght()     const { return TMath::Sqrt(fCrossPoint[0]*fCrossPoint[0] + fCrossPoint[1]*fCrossPoint[1] + fCrossPoint[2]*fCrossPoint[2]) ; }
   Float_t       Sigma()        const { return fSigma ; }
   Float_t       CrossPointX()  const { return fCrossPoint[0] ; }
@@ -71,7 +71,7 @@ public:
   void SetCosPointingAngle(Float_t cos) { fPointAngle = cos ; }
   void SetMostLikelihoodPID(Int_t pdgCode)               { fMostLikelihoodPID = pdgCode ; }
   void SetCrossPoint(Float_t pox,Float_t poy,Float_t poz) { fCrossPoint[0] = pox ; fCrossPoint[1] = poy ; fCrossPoint[2] = poz ; }
-  void SetDaughters(AliFlowTrack* pos, AliFlowTrack* neg) { fDaughters[0] = pos ; fDaughters[1] = neg ; }
+  void SetDaughters(AliFlowTrack* pos, AliFlowTrack* neg) { fDaughterP = pos ; fDaughterN = neg ; }
 
 
 private:
@@ -81,19 +81,20 @@ private:
   AliFlowV0 &operator=(const AliFlowV0 &flowV0) ; // Assignment Operator (dummy)
 
  // Data Members
-  Float_t  fPhi;                               // reconstructed azimuthal angle of the v0
-  Float_t  fPt;                                        // reconstructed transverse momentum of the v0 
-  Float_t  fEta;                               // reconstructed pseudorapidity of the v0
-  Float_t  fChi2;                              // chi2 of the reconstructed v0
-  Float_t  fMass;                              // reconstructed v0 mass 
-  Float_t  fDca;                               // distance of closest approach of the reconstructed v0 to the main vertex
-  Float_t  fCrossPoint[3] ;                    // crossing point coordinates of the two daughter tracks
-  Float_t  fCrossDCA ;                         // DCA between the 2 daughter tracks at the crossing point 
-  Float_t  fSigma ;                            // sigma of the DCA of the 2 daughter tracks at the crossing point 
-  Int_t    fLabel ;                            // Label of the V0 (link: KineTree-ESD) 
-  AliFlowTrack* fDaughters[2] ;                        // daughter tracks (pointers to the TracksCollection())
-  Int_t    fMostLikelihoodPID;                 // most probable P.Id. hypotesis
-  Float_t  fPointAngle;                        // cosine of the pointing angle
+  Float_t      fPhi;                           // reconstructed azimuthal angle of the v0
+  Float_t      fPt;                            // reconstructed transverse momentum of the v0 
+  Float_t      fEta;                           // reconstructed pseudorapidity of the v0
+  Float_t      fChi2;                          // chi2 of the reconstructed v0
+  Float_t      fMass;                          // reconstructed v0 mass 
+  Float_t      fDca;                           // distance of closest approach of the reconstructed v0 to the main vertex
+  Float_t      fCrossPoint[3] ;                // crossing point coordinates of the two daughter tracks
+  Float_t      fCrossDCA ;                     // DCA between the 2 daughter tracks at the crossing point 
+  Float_t      fSigma ;                        // sigma of the DCA of the 2 daughter tracks at the crossing point 
+  Int_t        fLabel ;                        // Label of the V0 (link: KineTree-ESD) 
+  Int_t        fMostLikelihoodPID;             // most probable P.Id. hypotesis
+  Float_t      fPointAngle;                    // cosine of the pointing angle
+  AliFlowTrack* fDaughterP ;                   // daughter tracks (pointers to the TracksCollection())
+  AliFlowTrack* fDaughterN ;                   // daughter tracks (pointers to the TracksCollection())
   
   ClassDef(AliFlowV0,2) ;                      // macro for rootcint
 };
index e1472c9..76d1c37 100644 (file)
@@ -53,7 +53,9 @@ using namespace std; //required for resolving the 'cout' symbol
 
 ClassImp(AliFlowWeighter) 
 //-----------------------------------------------------------------------
-AliFlowWeighter::AliFlowWeighter(const AliFlowSelection* flowSelect)
+AliFlowWeighter::AliFlowWeighter(const AliFlowSelection* flowSelect):
+  fFlowEvent(0x0), fFlowTrack(0x0), fFlowSelect(0x0), fFlowTracks(0x0),
+  fWgtFile(0x0), fPhiWgtHistList(0x0)
 {
  // default constructor (selection given or default selection)
 
index 7d3c207..9efe934 100644 (file)
@@ -230,7 +230,7 @@ void AliSelectorFoF::Begin(TTree*)
    cout << " . Writing Analysis Histograms on  : " << fFlowAnal->GetHistFileName()   << "  . " << endl ;
 
   // Analysis settings
-   fFlowAnal->SetFlowForV0(kFALSE) ;         // default kTRUE.
+   fFlowAnal->SetFlowForV0() ;         // default kTRUE.
    fFlowAnal->SetEtaSub() ;            // default kFALSE
    //fFlowAnal->SetV1Ep1Ep2() ;          // default kFALSE.
    //fFlowAnal->SetShuffle() ;           // default kFALSE. shuffles track array
@@ -240,6 +240,7 @@ void AliSelectorFoF::Begin(TTree*)
    fFlowAnal->SetUseBayWgt(kFALSE) ;    // default kFALSE. uses bayesian weights in P.id.
    //fFlowAnal->SetUsePtWgt(); // default kFALSE. uses pT as a weight for RP determination
    //fFlowAnal->SetUseEtaWgt();        // default kFALSE. uses eta as a weight for RP determination
+   //fFlowAnal->SetCustomRespFunc();   // default kFALSE. uses the combined response function from the ESD
    fFlowAnal->Init() ;
   }
  }
@@ -395,7 +396,7 @@ Bool_t AliSelectorFoF::Process(Long64_t entry)
   if(fOnFlyAnalysis)  
   {
    cout << " doing analysis :| " << endl ;
-   bool done = fFlowAnal->Analyse(fFlowEvent) ;
+   bool done = fFlowAnal->Analyze(fFlowEvent) ;
    fFlowAnal->PrintEventQuantities() ;
    if(done) { cout << "# analysis done :) " << entry << "                # ok ! #       " << endl ; } 
   }
index 2d083a8..0682cb1 100755 (executable)
@@ -24,6 +24,7 @@ rootcint -f AliFlowConstants_root.cxx -c -I$ALICE_ROOT/include -I$ROOTSYS/includ
 rootcint -f AliFlowMaker_root.cxx -c -I$ALICE_ROOT/include -I$ROOTSYS/include AliFlowMaker.h 
 rootcint -f AliFlowAnalyser_root.cxx -c -I$ALICE_ROOT/include -I$ROOTSYS/include AliFlowAnalyser.h
 rootcint -f AliFlowWeighter_root.cxx -c -I$ALICE_ROOT/include -I$ROOTSYS/include AliFlowWeighter.h
+rootcint -f AliFlowKineMaker_root.cxx -c -I$ALICE_ROOT/include -I$ROOTSYS/include AliFlowKineMaker.h 
 
 echo "                                                          code CINTed"
 
@@ -34,7 +35,7 @@ echo ""
 #echo "                                <ALIFLOW_MAK.SO>        for AliROOT compiled"
 #g++ -Wall -shared -g -I$ROOTSYS/include -o AliFlow_Ana.so AliFlowTrack*.cxx AliFlowV0*.cxx AliFlowSelection*.cxx AliFlowEvent*.cxx AliFlowAnalyser*.cxx AliFlowWeighter*.cxx AliFlowConstants*.cxx 
 #echo "                                <ALIFLOW_ANA.SO>           for ROOT compiled"
-g++ -Wall -shared -g -I$ALICE_ROOT/include -I$ROOTSYS/include -o AliFlow_All.so AliFlowTrack*.cxx AliFlowV0*.cxx AliFlowSelection*.cxx AliFlowEvent*.cxx AliFlowMaker*.cxx AliFlowAnalyser*.cxx AliFlowWeighter*.cxx AliFlowConstants*.cxx 
+g++ -Wall -shared -g -I$ALICE_ROOT/include -I$ROOTSYS/include -o AliFlow_All.so AliFlowTrack*.cxx AliFlowV0*.cxx AliFlowSelection*.cxx AliFlowEvent*.cxx AliFlowMaker*.cxx AliFlowAnalyser*.cxx AliFlowWeighter*.cxx AliFlowConstants*.cxx AliFlowKineMaker*.cxx 
 echo "                         <ALIFLOW_ALL.SO>        for AliROOT compiled"
 
 rm AliFlowTrack_root.*
@@ -43,6 +44,7 @@ rm AliFlowSelection_root.*
 rm AliFlowEvent_root.*      
 rm AliFlowConstants_root.* 
 rm AliFlowMaker_root.*     
+rm AliFlowKineMaker_root.*     
 rm AliFlowAnalyser_root.* 
 rm AliFlowWeighter_root.* 
 
diff --git a/PWG2/FLOW/macros/kOpen.C b/PWG2/FLOW/macros/kOpen.C
new file mode 100644 (file)
index 0000000..6b28a58
--- /dev/null
@@ -0,0 +1,35 @@
+//////////////////////////////////////////////////////////////////////
+// - old way -
+//////////////////////////////////////////////////////////////////////
+
+bool kOpen(int evtN=0)
+{
+ TString fileName = "galice.root" ;
+ AliRunLoader* rl = AliRunLoader::Open(fileName.Data(),"MyEvent","read");
+ rl->LoadgAlice();
+ AliRun* gAlice = rl->GetAliRun();
+ rl->LoadHeader();
+ rl->LoadKinematics();
+ Int_t fNumberOfEvents = rl->GetNumberOfEvents() ;
+ cout << " Found :  " << fNumberOfEvents << "  event(s) ... " << endl ; 
+
+ Int_t exitStatus = rl->GetEvent(evtN) ; if(exitStatus!=0) { return kFALSE ; }
+
+ TTree* pKTree = (TTree*)rl->TreeK();        // Particles TTree (KineTree)
+ AliStack* pStack = gAlice->Stack();         // Particles Stack (use "Label()" to get the number in the stack)
+
+ // else if(rl)      // opens files one by one (unload and reload)
+ // {
+ //  rl->UnloadgAlice() ;
+ //  rl->UnloadHeader() ;
+ //  rl->UnloadKinematics() ;
+ //  delete rl ; rl = 0 ;
+ // }
+
+ Int_t fNumberOfParticles = pKTree->GetEntries() ;
+ Int_t nPart = pStack->GetNtrack() ;
+ cout << " Event n. " << evtN << "  contains :  " << fNumberOfParticles << "  particles in the TTree  ( =  " << nPart << "  in the stack ) . " << endl ; 
+
+ return kTRUE ; 
+}
+
index 00a9e9e..aaa2b4d 100644 (file)
@@ -123,6 +123,7 @@ void testAnal(TString output = "flowEvtsAnal.root")
  flow->SetUseBayWgt(kFALSE) ;   // default kFALSE. uses bayesian weights in P.id.
  //flow->SetUsePtWgt();                // default kFALSE. uses pT as a weight for RP determination
  //flow->SetUseEtaWgt();       // default kFALSE. uses eta as a weight for RP determination
+ //flow->SetCustomRespFunc();  // default kFALSE. uses the combined response function from the ESD
 
  // init (make histograms, start the analysis)
  cout << endl ;
@@ -145,7 +146,7 @@ void testAnal(TString output = "flowEvtsAnal.root")
   TString evtName = flowEventsList->At(ie)->GetName() ;
   flowEventsFile->GetObject(evtName.Data(),flowEvt) ;
   // cout << "dumping event " << ie << " : " << endl ; flowEvt->Dump() ; cout << endl ; 
-  Bool_t succ = flow->Analyse(flowEvt) ;  
+  Bool_t succ = flow->Analyze(flowEvt) ;  
   flow->PrintEventQuantities() ;
   if(succ) { cout << ie << " done ... " << endl ; } 
  }
index e3187f2..01c797e 100644 (file)
@@ -1,10 +1,13 @@
 // from CreateESDChain.C - instead of  #include "CreateESDChain.C"
 TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 20, Int_t offset = 0) ;
 void LookupWrite(TChain* chain, const char* target) ;
-//
+// -
 
 void testFoF(const Char_t* dataDir=".", Int_t nRuns = 10, Int_t offset = 0)
 {
+ // include path (to find the .h files when compiling)
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include") ;
+
  // load needed libraries
  gSystem->Load("libESD");
  gSystem->Load("libPhysics.so");
diff --git a/PWG2/FLOW/macros/testKiner.C b/PWG2/FLOW/macros/testKiner.C
new file mode 100644 (file)
index 0000000..35b8133
--- /dev/null
@@ -0,0 +1,248 @@
+/////////////////////////////////////////////////////////////
+//
+// $Id$
+//
+// Author: Emanuele Simili
+//
+/////////////////////////////////////////////////////////////
+//
+// Description: AliRoot macro to make AliFlowEvents from KineTree (new way) 
+//
+/////////////////////////////////////////////////////////////
+
+#include <vector>
+#include <iostream>
+#include <fstream>
+#include "TMath.h"
+#include "TFile.h"
+#include "TObjArray"
+#include "TStopwatch.h"
+
+using namespace std; //required for resolving the 'cout' symbol
+
+TTree* kOpen(int evtN=0) ;
+//TChain* CreateKineChain(const char* aDataDir = "MCfiles.txt", Int_t aRuns = 10, Int_t offset = 0) ;
+//void LookupWrite(TChain* chain, const char* target) ;
+
+void testKiner(TString output = "flowKevts.root")
+{
+ cout << " . Here the new flow kinemaker (2007b) ... " << endl ;
+ cout << endl ;
+
+ bool kOne = kFALSE ;
+
+ TStopwatch timer;
+ timer.Start();
+
+ gSystem->Load("AliFlow_All.so");
+
+ // output file //
+
+ TFile * fFlowfile = new TFile(output.Data(),"RECREATE") ;
+ //fFlowfile->cd() ; 
+
+ // esd chain //
+
+// // TString fESDfileName = "AliESDs.root" ; TString fESDtree = "esdTree" ; 
+// TChain* pMCchain = CreateESDChain(".",10,0);
+// Int_t fNumberOfEvents = (Int_t)pMCchain->GetEntries() ;
+// cout << " tot. " << fNumberOfEvents << " events in the TChain ... " << endl ; cout << endl ;
+
+// TString fKineBranch = "TreeK" ;
+// TTree treeK = 0 ;
+// pMCchain->SetBranchAddress(fKineBranch.Data(),&treeK) ;
+
+ Int_t fNumberOfEvents = 1 ;
+ TTree* treeK = 0 ;
+
+ // flow maker //
+
+ AliFlowKineMaker * flowKiner = new  AliFlowKineMaker() ;
+ // cuts, etc.
+ flowKiner->SetAbsEtaCut(2.1) ;
+ flowKiner->SetECut(0.01,100.) ;
+ //flowKiner->SetLabelCut(..,..) ;
+ flowKiner->SetPrimaryCut(Bool_t prim = kTRUE) ;
+ flowKiner->PrintCutList() ;
+
+ // loop //
+
+ Int_t evtN = 0 ;
+ AliFlowEvent * flowEvt = 0 ;
+ for(evtN=0;evtN<fNumberOfEvents;evtN++)
+ {
+  treeK = kOpen(0) ;
+  
+  Int_t evtNN = -1 ;
+  Int_t nTrk = treeK->GetEntries() ;
+
+  cout << endl ; cout << " Event " << evtN << "  ( " << evtNN << " )  : " << nTrk << " tracks  ." << endl ;
+
+  flowEvt = flowKiner->FillFlowEvent(treeK) ;
+
+  cout << " Event filled " << flowEvt << " ... " << endl ;
+  // cout << endl ; cout << " trks : " << flowEvt->TrackCollection()->GetEntries() << endl ;
+  // flowEvt->Dump() ; cout << endl ;
+
+  TString evtID = "" ; evtID += evtN ; 
+  fFlowfile->cd() ; 
+  flowEvt->Write(evtID.Data()) ;
+  cout <<  " Event " << evtN << "  ( " << evtID.Data() << " )  -  written on disk (" << output << ") ." << endl;
+  delete flowEvt ;
+ }
+ fFlowfile->Close() ; 
+ cout <<  endl ;
+ cout << " Finished ... " << endl ;
+ cout << "  nTracks:  " << flowKiner->GetNgoodTracks() << endl ;   
+ cout << "  nV0s:  " << flowKiner->GetNgoodV0s()  << endl ;         
+ cout << "  nTracks (|eta|<0.5):  " << flowKiner->GetNgoodTracksEta() << endl ; 
+ cout << "  nTracks+:  " << flowKiner->GetNposiTracks() << endl ;           
+ cout << "  nTracks-:  " << flowKiner->GetNnegaTracks() << endl ;           
+ cout << "  nTracks unconstrained:  " << flowKiner->GetNunconstrained() << endl ;       
+ //cout << "  Bayesian :  " ; 
+ //for(int ii=0;ii<5;ii++) { cout << flowKiner->GetBayesianNorm(ii) << "   " ; } 
+ cout << " . " << endl ; 
+
+ timer.Stop() ;
+ cout << endl ;
+ timer.Print() ;
+ cout << " . here it was (kiner) ... " << endl ;  //juice!
+ cout << endl ;
+
+ // break ;
+
+}
+
+TTree* kOpen(int evtN)
+{
+ TString fileName = "./0/galice.root" ;
+ AliRunLoader* rl = AliRunLoader::Open(fileName.Data(),"MyEvent","read");
+ rl->LoadgAlice();
+ AliRun* gAlice = rl->GetAliRun();
+ rl->LoadHeader();
+ rl->LoadKinematics();
+ Int_t fNumberOfEvents = rl->GetNumberOfEvents() ;
+ cout << " Found :  " << fNumberOfEvents << "  event(s) ... " << endl ; 
+
+ Int_t exitStatus = rl->GetEvent(evtN) ; if(exitStatus!=0) { return 0 ; }
+
+ TTree* kTree = (TTree*)rl->TreeK();         // Particles TTree (KineTree)
+ AliStack* kStack = gAlice->Stack();         // Particles Stack (use "Label()" to get the number in the stack)
+
+ Int_t fNumberOfParticles = kTree->GetEntries() ;
+ Int_t nPart = kStack->GetNtrack() ;
+ cout << " Event n. " << evtN << "  contains :  " << fNumberOfParticles << "  particles in the TTree  ( =  " << nPart << "  in the stack ) . " << endl ; 
+
+ return kTree ; 
+}
+
+
+// // Helper macros for creating chains (from: CreateESDChain.C,v 1.10 jgrosseo Exp)
+// TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
+// {
+//   // creates chain of files in a given directory or file containing a list.
+//   // In case of directory the structure is expected as:
+//   // <aDataDir>/<dir0>/AliESDs.root
+//   // <aDataDir>/<dir1>/AliESDs.root
+//   // ...
+// 
+//   if (!aDataDir)
+//     return 0;
+// 
+//   Long_t id, size, flags, modtime;
+//   if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
+//   {
+//     printf("%s not found.\n", aDataDir);
+//     return 0;
+//   }
+// 
+//   TChain* chain = new TChain("esdTree");
+//   TChain* chaingAlice = 0;
+// 
+//   if (flags & 2)
+//   {
+//     TString execDir(gSystem->pwd());
+//     TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
+//     TList* dirList            = baseDir->GetListOfFiles();
+//     Int_t nDirs               = dirList->GetEntries();
+//     gSystem->cd(execDir);
+// 
+//     Int_t count = 0;
+// 
+//     for (Int_t iDir=0; iDir<nDirs; ++iDir)
+//     {
+//       TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
+//       if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
+//         continue;
+// 
+//       if (offset > 0)
+//       {
+//         --offset;
+//         continue;
+//       }
+// 
+//       if (count++ == aRuns)
+//         break;
+// 
+//       TString presentDirName(aDataDir);
+//       presentDirName += "/";
+//       presentDirName += presentDir->GetName();
+// 
+//       chain->Add(presentDirName + "/AliESDs.root/esdTree");
+//     }
+//   }
+//   else
+//   {
+//     // Open the input stream
+//     ifstream in;
+//     in.open(aDataDir);
+// 
+//     Int_t count = 0;
+// 
+//     // Read the input list of files and add them to the chain
+//     TString esdfile;
+//     while(in.good()) {
+//       in >> esdfile;
+//       if (!esdfile.Contains("root")) continue; // protection
+// 
+//       if (offset > 0)
+//       {
+//         --offset;
+//         continue;
+//       }
+// 
+//       if (count++ == aRuns)
+//         break;
+// 
+//         // add esd file
+//       chain->Add(esdfile);
+//     }
+// 
+//     in.close();
+//   }
+// 
+//   return chain;
+// }
+// 
+// void LookupWrite(TChain* chain, const char* target)
+// {
+//   // looks up the chain and writes the remaining files to the text file target
+// 
+//   chain->Lookup();
+// 
+//   TObjArray* list = chain->GetListOfFiles();
+//   TIterator* iter = list->MakeIterator();
+//   TObject* obj = 0;
+// 
+//   ofstream outfile;
+//   outfile.open(target);
+// 
+//   while ((obj = iter->Next()))
+//     outfile << obj->GetTitle() << "#AliESDs.root" << endl;
+// 
+//   outfile.close();
+// 
+//   delete iter;
+// }
index 1ff23a0..b536162 100644 (file)
 
 using namespace std; //required for resolving the 'cout' symbol
 
+TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 10, Int_t offset = 0) ;
+void LookupWrite(TChain* chain, const char* target) ;
+
 void testMaker(TString output = "flowEvts.root")
 {
- cout << " . Here the new flow maker (2007) ... " << endl ;
+ cout << " . Here the new flow maker (2007b) ... " << endl ;
  cout << endl ;
 
  bool kOne = kFALSE ;
@@ -39,32 +42,12 @@ void testMaker(TString output = "flowEvts.root")
 
  // esd chain //
 
- TString fESDfileName = "AliESDs.root" ;
- TString fESDtree = "esdTree" ;
- TString fESDbranch = "ESD" ;
-
- TChain * pESDchain = new TChain(fESDtree.Data()) ;
-
- if(kOne)
- {
-  pESDchain->Add(fESDfileName.Data()) ; // nFiles++ ;
- }
- else
- {
-  pESDchain->Add("1/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("2/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("3/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("4/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("5/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("6/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("7/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("8/AliESDs.root") ; // nFiles++ ;
-  pESDchain->Add("9/AliESDs.root") ; // nFiles++ ;
- }
-
+ // TString fESDfileName = "AliESDs.root" ; TString fESDtree = "esdTree" ; 
+ TChain* pESDchain = CreateESDChain(".",10,0);
  Int_t fNumberOfEvents = (Int_t)pESDchain->GetEntries() ;
  cout << " tot. " << fNumberOfEvents << " events in the TChain ... " << endl ; cout << endl ;
 
+ TString fESDbranch = "ESD" ;
  AliESD * pEsd = 0 ;
  pESDchain->SetBranchAddress(fESDbranch.Data(),&pEsd) ;
 
@@ -126,3 +109,112 @@ void testMaker(TString output = "flowEvts.root")
  // break ;
 
 }
+
+
+// Helper macros for creating chains (from: CreateESDChain.C,v 1.10 jgrosseo Exp)
+TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
+{
+  // creates chain of files in a given directory or file containing a list.
+  // In case of directory the structure is expected as:
+  // <aDataDir>/<dir0>/AliESDs.root
+  // <aDataDir>/<dir1>/AliESDs.root
+  // ...
+
+  if (!aDataDir)
+    return 0;
+
+  Long_t id, size, flags, modtime;
+  if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
+  {
+    printf("%s not found.\n", aDataDir);
+    return 0;
+  }
+
+  TChain* chain = new TChain("esdTree");
+  TChain* chaingAlice = 0;
+
+  if (flags & 2)
+  {
+    TString execDir(gSystem->pwd());
+    TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
+    TList* dirList            = baseDir->GetListOfFiles();
+    Int_t nDirs               = dirList->GetEntries();
+    gSystem->cd(execDir);
+
+    Int_t count = 0;
+
+    for (Int_t iDir=0; iDir<nDirs; ++iDir)
+    {
+      TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
+      if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
+        continue;
+
+      if (offset > 0)
+      {
+        --offset;
+        continue;
+      }
+
+      if (count++ == aRuns)
+        break;
+
+      TString presentDirName(aDataDir);
+      presentDirName += "/";
+      presentDirName += presentDir->GetName();
+
+      chain->Add(presentDirName + "/AliESDs.root/esdTree");
+    }
+  }
+  else
+  {
+    // Open the input stream
+    ifstream in;
+    in.open(aDataDir);
+
+    Int_t count = 0;
+
+    // Read the input list of files and add them to the chain
+    TString esdfile;
+    while(in.good()) {
+      in >> esdfile;
+      if (!esdfile.Contains("root")) continue; // protection
+
+      if (offset > 0)
+      {
+        --offset;
+        continue;
+      }
+
+      if (count++ == aRuns)
+        break;
+
+        // add esd file
+      chain->Add(esdfile);
+    }
+
+    in.close();
+  }
+
+  return chain;
+}
+
+void LookupWrite(TChain* chain, const char* target)
+{
+  // looks up the chain and writes the remaining files to the text file target
+
+  chain->Lookup();
+
+  TObjArray* list = chain->GetListOfFiles();
+  TIterator* iter = list->MakeIterator();
+  TObject* obj = 0;
+
+  ofstream outfile;
+  outfile.open(target);
+
+  while ((obj = iter->Next()))
+    outfile << obj->GetTitle() << "#AliESDs.root" << endl;
+
+  outfile.close();
+
+  delete iter;
+}