small mods for separate task
[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
46bec39c 42AliFlowEventSimple::AliFlowEventSimple():
43 fTrackCollection(NULL),
44 fNumberOfTracks(0),
45 fEventNSelTracksIntFlow(0)
46{
47 cout << "AliFlowEventSimple: Default constructor to be used only by root for io" << endl;
48}
49
50//-----------------------------------------------------------------------
51
52AliFlowEventSimple::AliFlowEventSimple(Int_t aLenght):
53 fTrackCollection(NULL),
54 fNumberOfTracks(0),
55 fEventNSelTracksIntFlow(0)
f1d945a1 56{
57 //constructor
e35ddff0 58 fTrackCollection = new TObjArray(aLenght) ;
f1d945a1 59}
60
61//-----------------------------------------------------------------------
62
e35ddff0 63AliFlowEventSimple::AliFlowEventSimple(const AliFlowEventSimple& anEvent):
bc6b015e 64 TObject(),
e35ddff0 65 fTrackCollection(anEvent.fTrackCollection),
66 fNumberOfTracks(anEvent.fNumberOfTracks),
03a02aca 67 fEventNSelTracksIntFlow(anEvent.fEventNSelTracksIntFlow)
f1d945a1 68{
69 //copy constructor
f1d945a1 70}
71
72//-----------------------------------------------------------------------
73
e35ddff0 74AliFlowEventSimple& AliFlowEventSimple::operator=(const AliFlowEventSimple& anEvent)
f1d945a1 75{
26c4cbb9 76 *fTrackCollection = *anEvent.fTrackCollection ;
e35ddff0 77 fNumberOfTracks = anEvent.fNumberOfTracks;
78 fEventNSelTracksIntFlow = anEvent.fEventNSelTracksIntFlow;
26c4cbb9 79
f1d945a1 80 return *this;
f1d945a1 81}
82
f1d945a1 83//-----------------------------------------------------------------------
84
85AliFlowEventSimple::~AliFlowEventSimple()
86{
87 //destructor
46bec39c 88 if (fTrackCollection) {
89 fTrackCollection->Delete() ; delete fTrackCollection ;
90 }
91 else {
92 cout << "AliFlowEventSimple: Warning trying to delete track collections NULL pointer" << endl;
93 }
f1d945a1 94}
95
96//-----------------------------------------------------------------------
97
98AliFlowTrackSimple* AliFlowEventSimple::GetTrack(Int_t i)
99{
100 //get track i from collection
e35ddff0 101 AliFlowTrackSimple* pTrack = (AliFlowTrackSimple*)TrackCollection()->At(i) ;
102 return pTrack;
f1d945a1 103}
104
105//-----------------------------------------------------------------------
03a02aca 106AliFlowVector AliFlowEventSimple::GetQ(Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
f1d945a1 107{
ae733b3b 108 // calculate Q-vector in harmonic n without weights (default harmonic n=2)
e35ddff0 109 Double_t dQX = 0.;
110 Double_t dQY = 0.;
111 AliFlowVector vQ;
112 vQ.Set(0.,0.);
9825d4a9 113
114 Int_t iOrder = n;
ae733b3b 115 Double_t iUsedTracks = 0;
26c4cbb9 116 Double_t dPhi=0.;
ae733b3b 117 Double_t dPt=0.;
118 Double_t dEta=0.;
26c4cbb9 119
120 AliFlowTrackSimple* pTrack = NULL;
121
122 Int_t nBinsPhi=0;
123 Double_t dBinWidthPt=0.;
ae733b3b 124 Double_t dPtMin=0.;
26c4cbb9 125 Double_t dBinWidthEta=0.;
ae733b3b 126 Double_t dEtaMin=0.;
26c4cbb9 127
ae733b3b 128 Double_t wPhi=1.; // weight Phi
129 Double_t wPt=1.; // weight Pt
130 Double_t wEta=1.; // weight Eta
26c4cbb9 131
03a02aca 132 TH1F *phiWeights = NULL;
ae733b3b 133 TH1D *ptWeights = NULL;
03a02aca 134 TH1D *etaWeights = NULL;
135
ae733b3b 136 Double_t dSumOfWeightsToPower2 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 2)
137 Double_t dSumOfWeightsToPower3 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 3)
138 Double_t dSumOfWeightsToPower4 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 4)
139 Double_t dSumOfWeightsToPower5 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 5)
140 Double_t dSumOfWeightsToPower6 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 6)
141 Double_t dSumOfWeightsToPower7 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 7)
142 Double_t dSumOfWeightsToPower8 = 0.; // sum_{i=1}^{n} pow((wPhi*wPt*wEta)_i, 8)
143
03a02aca 144 if(weightsList)
26c4cbb9 145 {
03a02aca 146 if(usePhiWeights)
26c4cbb9 147 {
03a02aca 148 phiWeights = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights"));
149 if(phiWeights) nBinsPhi = phiWeights->GetNbinsX();
150 }
151 if(usePtWeights)
26c4cbb9 152 {
03a02aca 153 ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
154 if(ptWeights)
155 {
156 dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
ae733b3b 157 dPtMin = (ptWeights->GetXaxis())->GetXmin();
03a02aca 158 }
159 }
160 if(useEtaWeights)
26c4cbb9 161 {
03a02aca 162 etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
163 if(etaWeights)
164 {
165 dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
ae733b3b 166 dEtaMin = (etaWeights->GetXaxis())->GetXmin();
03a02aca 167 }
168 }
169 } // end of if(weightsList)
26c4cbb9 170
03a02aca 171 // loop over tracks
26c4cbb9 172 for(Int_t i=0;i<fNumberOfTracks;i++)
173 {
174 pTrack = (AliFlowTrackSimple*)TrackCollection()->At(i);
175 if(pTrack)
176 {
177 if(pTrack->UseForIntegratedFlow())
f1d945a1 178 {
26c4cbb9 179 dPhi = pTrack->Phi();
180 dPt = pTrack->Pt();
181 dEta = pTrack->Eta();
ae733b3b 182
183 // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
184 if(phiWeights && nBinsPhi)
26c4cbb9 185 {
ae733b3b 186 wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
26c4cbb9 187 }
ae733b3b 188 // determine v'(pt) weight:
189 if(ptWeights && dBinWidthPt)
26c4cbb9 190 {
ae733b3b 191 wPt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
26c4cbb9 192 }
ae733b3b 193 // determine v'(eta) weight:
194 if(etaWeights && dBinWidthEta)
26c4cbb9 195 {
ae733b3b 196 wEta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
03a02aca 197 }
ae733b3b 198
199 // building up the weighted Q-vector:
26c4cbb9 200 dQX += wPhi*wPt*wEta*TMath::Cos(iOrder*dPhi);
201 dQY += wPhi*wPt*wEta*TMath::Sin(iOrder*dPhi);
ae733b3b 202
203 // weighted multiplicity:
204 iUsedTracks+=wPhi*wPt*wEta;
205
206 // weights raised to various powers are summed up:
207 dSumOfWeightsToPower2+=pow(wPhi*wPt*wEta, 2);
208 dSumOfWeightsToPower3+=pow(wPhi*wPt*wEta, 3);
209 dSumOfWeightsToPower4+=pow(wPhi*wPt*wEta, 4);
210 dSumOfWeightsToPower5+=pow(wPhi*wPt*wEta, 5);
211 dSumOfWeightsToPower6+=pow(wPhi*wPt*wEta, 6);
212 dSumOfWeightsToPower7+=pow(wPhi*wPt*wEta, 7);
213 dSumOfWeightsToPower8+=pow(wPhi*wPt*wEta, 8);
03a02aca 214
ae733b3b 215 } // end of if (pTrack->UseForIntegratedFlow())
216 } // end of if (pTrack)
26c4cbb9 217 else {cerr << "no particle!!!"<<endl;}
ae733b3b 218 } // loop over particles
26c4cbb9 219
e35ddff0 220 vQ.Set(dQX,dQY);
221 vQ.SetMult(iUsedTracks);
ae733b3b 222 vQ.SetSumOfWeightsToPower2(dSumOfWeightsToPower2);
223 vQ.SetSumOfWeightsToPower3(dSumOfWeightsToPower3);
224 vQ.SetSumOfWeightsToPower4(dSumOfWeightsToPower4);
225 vQ.SetSumOfWeightsToPower5(dSumOfWeightsToPower5);
226 vQ.SetSumOfWeightsToPower6(dSumOfWeightsToPower6);
227 vQ.SetSumOfWeightsToPower7(dSumOfWeightsToPower7);
228 vQ.SetSumOfWeightsToPower8(dSumOfWeightsToPower8);
26c4cbb9 229
e35ddff0 230 return vQ;
f1d945a1 231
5fef318d 232}
233
234
26c4cbb9 235
236
237
238
239
240
241
242