Using of the normalized energy deposition in default AliESDtrack.
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowCommon / AliFlowEventSimple.cxx
CommitLineData
f1d945a1 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
f1d945a1 16#include "Riostream.h"
17#include "TObjArray.h"
26c4cbb9 18#include "TFile.h"
03a02aca 19#include "TList.h"
f1d945a1 20#include "TMath.h"
26c4cbb9 21#include "TH1F.h"
22#include "TH1D.h"
23#include "TProfile.h"
f1d945a1 24#include "AliFlowVector.h"
25#include "AliFlowTrackSimple.h"
26#include "AliFlowEventSimple.h"
27
26c4cbb9 28/**************************************
29 * AliFlowEventSimple: A simple event *
30 * for flow analysis *
31 * *
32 * authors: Naomi van der Kolk *
33 * (kolk@nikhef.nl) *
34 * Ante Bilandzic *
35 * (anteb@nikhef.nl) *
36 * ***********************************/
37
f1d945a1 38ClassImp(AliFlowEventSimple)
39
40//-----------------------------------------------------------------------
41
e35ddff0 42 AliFlowEventSimple::AliFlowEventSimple(Int_t aLenght):
d29ba078 43 fTrackCollection(NULL),
f1d945a1 44 fNumberOfTracks(0),
03a02aca 45 fEventNSelTracksIntFlow(0)
f1d945a1 46{
47 //constructor
e35ddff0 48 fTrackCollection = new TObjArray(aLenght) ;
f1d945a1 49}
50
51//-----------------------------------------------------------------------
52
e35ddff0 53AliFlowEventSimple::AliFlowEventSimple(const AliFlowEventSimple& anEvent):
bc6b015e 54 TObject(),
e35ddff0 55 fTrackCollection(anEvent.fTrackCollection),
56 fNumberOfTracks(anEvent.fNumberOfTracks),
03a02aca 57 fEventNSelTracksIntFlow(anEvent.fEventNSelTracksIntFlow)
f1d945a1 58{
59 //copy constructor
f1d945a1 60}
61
62//-----------------------------------------------------------------------
63
e35ddff0 64AliFlowEventSimple& AliFlowEventSimple::operator=(const AliFlowEventSimple& anEvent)
f1d945a1 65{
26c4cbb9 66 *fTrackCollection = *anEvent.fTrackCollection ;
e35ddff0 67 fNumberOfTracks = anEvent.fNumberOfTracks;
68 fEventNSelTracksIntFlow = anEvent.fEventNSelTracksIntFlow;
26c4cbb9 69
f1d945a1 70 return *this;
f1d945a1 71}
72
f1d945a1 73//-----------------------------------------------------------------------
74
75AliFlowEventSimple::~AliFlowEventSimple()
76{
77 //destructor
78 fTrackCollection->Delete() ; delete fTrackCollection ;
79}
80
81//-----------------------------------------------------------------------
82
83AliFlowTrackSimple* AliFlowEventSimple::GetTrack(Int_t i)
84{
85 //get track i from collection
e35ddff0 86 AliFlowTrackSimple* pTrack = (AliFlowTrackSimple*)TrackCollection()->At(i) ;
87 return pTrack;
f1d945a1 88}
89
90//-----------------------------------------------------------------------
03a02aca 91AliFlowVector AliFlowEventSimple::GetQ(Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
f1d945a1 92{
ae733b3b 93 // calculate Q-vector in harmonic n without weights (default harmonic n=2)
e35ddff0 94 Double_t dQX = 0.;
95 Double_t dQY = 0.;
96 AliFlowVector vQ;
97 vQ.Set(0.,0.);
9825d4a9 98
99 Int_t iOrder = n;
ae733b3b 100 Double_t iUsedTracks = 0;
26c4cbb9 101 Double_t dPhi=0.;
ae733b3b 102 Double_t dPt=0.;
103 Double_t dEta=0.;
26c4cbb9 104
105 AliFlowTrackSimple* pTrack = NULL;
106
107 Int_t nBinsPhi=0;
108 Double_t dBinWidthPt=0.;
ae733b3b 109 Double_t dPtMin=0.;
26c4cbb9 110 Double_t dBinWidthEta=0.;
ae733b3b 111 Double_t dEtaMin=0.;
26c4cbb9 112
ae733b3b 113 Double_t wPhi=1.; // weight Phi
114 Double_t wPt=1.; // weight Pt
115 Double_t wEta=1.; // weight Eta
26c4cbb9 116
03a02aca 117 TH1F *phiWeights = NULL;
ae733b3b 118 TH1D *ptWeights = NULL;
03a02aca 119 TH1D *etaWeights = NULL;
120
ae733b3b 121 Double_t dSumOfWeightsToPower2 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 2)
122 Double_t dSumOfWeightsToPower3 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 3)
123 Double_t dSumOfWeightsToPower4 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 4)
124 Double_t dSumOfWeightsToPower5 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 5)
125 Double_t dSumOfWeightsToPower6 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 6)
126 Double_t dSumOfWeightsToPower7 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 7)
127 Double_t dSumOfWeightsToPower8 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 8)
128
03a02aca 129 if(weightsList)
26c4cbb9 130 {
03a02aca 131 if(usePhiWeights)
26c4cbb9 132 {
03a02aca 133 phiWeights = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights"));
134 if(phiWeights) nBinsPhi = phiWeights->GetNbinsX();
135 }
136 if(usePtWeights)
26c4cbb9 137 {
03a02aca 138 ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
139 if(ptWeights)
140 {
141 dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
ae733b3b 142 dPtMin = (ptWeights->GetXaxis())->GetXmin();
03a02aca 143 }
144 }
145 if(useEtaWeights)
26c4cbb9 146 {
03a02aca 147 etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
148 if(etaWeights)
149 {
150 dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
ae733b3b 151 dEtaMin = (etaWeights->GetXaxis())->GetXmin();
03a02aca 152 }
153 }
154 } // end of if(weightsList)
26c4cbb9 155
03a02aca 156 // loop over tracks
26c4cbb9 157 for(Int_t i=0;i<fNumberOfTracks;i++)
158 {
159 pTrack = (AliFlowTrackSimple*)TrackCollection()->At(i);
160 if(pTrack)
161 {
162 if(pTrack->UseForIntegratedFlow())
f1d945a1 163 {
26c4cbb9 164 dPhi = pTrack->Phi();
165 dPt = pTrack->Pt();
166 dEta = pTrack->Eta();
ae733b3b 167
168 // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
169 if(phiWeights && nBinsPhi)
26c4cbb9 170 {
ae733b3b 171 wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
26c4cbb9 172 }
ae733b3b 173 // determine v'(pt) weight:
174 if(ptWeights && dBinWidthPt)
26c4cbb9 175 {
ae733b3b 176 wPt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
26c4cbb9 177 }
ae733b3b 178 // determine v'(eta) weight:
179 if(etaWeights && dBinWidthEta)
26c4cbb9 180 {
ae733b3b 181 wEta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
03a02aca 182 }
ae733b3b 183
184 // building up the weighted Q-vector:
26c4cbb9 185 dQX += wPhi*wPt*wEta*TMath::Cos(iOrder*dPhi);
186 dQY += wPhi*wPt*wEta*TMath::Sin(iOrder*dPhi);
ae733b3b 187
188 // weighted multiplicity:
189 iUsedTracks+=wPhi*wPt*wEta;
190
191 // weights raised to various powers are summed up:
192 dSumOfWeightsToPower2+=pow(wPhi*wPt*wEta, 2);
193 dSumOfWeightsToPower3+=pow(wPhi*wPt*wEta, 3);
194 dSumOfWeightsToPower4+=pow(wPhi*wPt*wEta, 4);
195 dSumOfWeightsToPower5+=pow(wPhi*wPt*wEta, 5);
196 dSumOfWeightsToPower6+=pow(wPhi*wPt*wEta, 6);
197 dSumOfWeightsToPower7+=pow(wPhi*wPt*wEta, 7);
198 dSumOfWeightsToPower8+=pow(wPhi*wPt*wEta, 8);
03a02aca 199
ae733b3b 200 } // end of if (pTrack->UseForIntegratedFlow())
201 } // end of if (pTrack)
26c4cbb9 202 else {cerr << "no particle!!!"<<endl;}
ae733b3b 203 } // loop over particles
26c4cbb9 204
e35ddff0 205 vQ.Set(dQX,dQY);
206 vQ.SetMult(iUsedTracks);
ae733b3b 207 vQ.SetSumOfWeightsToPower2(dSumOfWeightsToPower2);
208 vQ.SetSumOfWeightsToPower3(dSumOfWeightsToPower3);
209 vQ.SetSumOfWeightsToPower4(dSumOfWeightsToPower4);
210 vQ.SetSumOfWeightsToPower5(dSumOfWeightsToPower5);
211 vQ.SetSumOfWeightsToPower6(dSumOfWeightsToPower6);
212 vQ.SetSumOfWeightsToPower7(dSumOfWeightsToPower7);
213 vQ.SetSumOfWeightsToPower8(dSumOfWeightsToPower8);
26c4cbb9 214
e35ddff0 215 return vQ;
f1d945a1 216
5fef318d 217}
218
219
26c4cbb9 220
221
222
223
224
225
226
227