// 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.
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)
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 ;
//-----------------------------------------------------------------------
// ###
//-----------------------------------------------------------------------
-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
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
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
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
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
// - 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.
//-------------------------------------------------------------
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) ;
}
}
//-------------------------------------------------------------
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
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)
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; }
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
};
--- /dev/null
+//////////////////////////////////////////////////////////////////////
+//
+// $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 ; }
+//
+//////////////////////////////////////////////////////////////////////
+
--- /dev/null
+//////////////////////////////////////////////////////////////////////
+//
+// $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
+
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
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
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()) ;
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()) ;
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()) ;
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()) ;
//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:
//////////////////////////////////////////////////////////////////////
//
// 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!
//
//////////////////////////////////////////////////////////////////////
class AliFlowMaker {
public:
+
AliFlowMaker();
virtual ~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"
}
// PID probability
- float pidProb = pFlowTrack->MostLikelihoodProb() ;
+ float pidProb = pFlowTrack->MostLikelihoodRespFunc() ;
if(fPidProbPart[1] > fPidProbPart[0] && (pidProb < fPidProbPart[0] || pidProb > fPidProbPart[1])) return kFALSE;
// Constrainable
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
};
// 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
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() ;
}
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
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.
//-------------------------------------------------------------
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())
{
//-------------------------------------------------------------
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())
{
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; }
// 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 ; }
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 ; }
// 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
ClassImp(AliFlowV0)
//////////////////////////////////////////////////////////////////////////////
-AliFlowV0::AliFlowV0()
+AliFlowV0::AliFlowV0():
+ fDaughterP(0x0), fDaughterN(0x0)
{
// default constructor
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)
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] ; }
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:
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
};
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)
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
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() ;
}
}
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 ; }
}
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"
#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.*
rm AliFlowEvent_root.*
rm AliFlowConstants_root.*
rm AliFlowMaker_root.*
+rm AliFlowKineMaker_root.*
rm AliFlowAnalyser_root.*
rm AliFlowWeighter_root.*
--- /dev/null
+//////////////////////////////////////////////////////////////////////
+// - 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 ;
+}
+
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 ;
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 ; }
}
// 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");
--- /dev/null
+/////////////////////////////////////////////////////////////
+//
+// $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;
+// }
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 ;
// 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) ;
// 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;
+}