1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 // AliFlowTrackSimple:
19 // A simple track class to the the AliFlowEventSimple for flow analysis
22 // author: N. van der Kolk (kolk@nikhef.nl)
23 // mods: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
26 #include "TParticle.h"
27 #include "TParticlePDG.h"
28 #include "AliFlowTrackSimple.h"
32 ClassImp(AliFlowTrackSimple)
34 //-----------------------------------------------------------------------
35 AliFlowTrackSimple::AliFlowTrackSimple():
50 //-----------------------------------------------------------------------
51 AliFlowTrackSimple::AliFlowTrackSimple(Double_t phi, Double_t eta, Double_t pt, Double_t weight, Int_t charge, Double_t mass):
66 //-----------------------------------------------------------------------
67 AliFlowTrackSimple::AliFlowTrackSimple( TParticle* p ):
80 TParticlePDG* ppdg = p->GetPDG();
81 fCharge = TMath::Nint(ppdg->Charge()/3.0);
85 //-----------------------------------------------------------------------
86 void AliFlowTrackSimple::Set(TParticle* p)
88 //set from a TParticle
93 TParticlePDG* ppdg = p->GetPDG();
94 fCharge = TMath::Nint(ppdg->Charge()/3.0);
98 //-----------------------------------------------------------------------
99 AliFlowTrackSimple::AliFlowTrackSimple(const AliFlowTrackSimple& aTrack):
104 fTrackWeight(aTrack.fTrackWeight),
105 fCharge(aTrack.fCharge),
107 fFlowBits(aTrack.fFlowBits),
108 fSubEventBits(aTrack.fSubEventBits),
114 //-----------------------------------------------------------------------
115 AliFlowTrackSimple* AliFlowTrackSimple::Clone(const char* /*option*/) const
117 //clone "constructor"
118 return new AliFlowTrackSimple(*this);
121 //-----------------------------------------------------------------------
122 AliFlowTrackSimple& AliFlowTrackSimple::operator=(const AliFlowTrackSimple& aTrack)
125 if (&aTrack==this) return *this; //handle self assignmnet
129 fTrackWeight = aTrack.fTrackWeight;
130 fCharge = aTrack.fCharge;
131 fMass = aTrack.fMass;
132 fFlowBits = aTrack.fFlowBits;
133 fSubEventBits = aTrack.fSubEventBits;
139 //-----------------------------------------------------------------------
140 AliFlowTrackSimple::~AliFlowTrackSimple()
146 //-----------------------------------------------------------------------
147 void AliFlowTrackSimple::ResolutionPt(Double_t res)
149 //smear the pt by a gaussian with sigma=res
150 fPt += gRandom->Gaus(0.,res);
153 //-----------------------------------------------------------------------
154 void AliFlowTrackSimple::AddV1( Double_t v1,
155 Double_t reactionPlaneAngle,
156 Double_t precisionPhi,
157 Int_t maxNumberOfIterations )
159 //afterburner, adds v1, uses Newton-Raphson iteration
165 for (Int_t i=0; i<maxNumberOfIterations; i++)
167 phiprev=fPhi; //store last value for comparison
168 f = fPhi-phi0+2.0*v1*TMath::Sin(fPhi-reactionPlaneAngle);
169 fp = 1.0+2.0*v1*TMath::Cos(fPhi-reactionPlaneAngle); //first derivative
171 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
175 //-----------------------------------------------------------------------
176 void AliFlowTrackSimple::AddV2( Double_t v2,
177 Double_t reactionPlaneAngle,
178 Double_t precisionPhi,
179 Int_t maxNumberOfIterations )
181 //afterburner, adds v2, uses Newton-Raphson iteration
187 for (Int_t i=0; i<maxNumberOfIterations; i++)
189 phiprev=fPhi; //store last value for comparison
190 f = fPhi-phi0+v2*TMath::Sin(2.*(fPhi-reactionPlaneAngle));
191 fp = 1.0+2.0*v2*TMath::Cos(2.*(fPhi-reactionPlaneAngle)); //first derivative
193 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
197 //-----------------------------------------------------------------------
198 void AliFlowTrackSimple::AddV3( Double_t v3,
199 Double_t reactionPlaneAngle,
200 Double_t precisionPhi,
201 Int_t maxNumberOfIterations )
203 //afterburner, adds v3, uses Newton-Raphson iteration
209 for (Int_t i=0; i<maxNumberOfIterations; i++)
211 phiprev=fPhi; //store last value for comparison
212 f = fPhi-phi0+2./3.*v3*TMath::Sin(3.*(fPhi-reactionPlaneAngle));
213 fp = 1.0+2.0*v3*TMath::Cos(3.*(fPhi-reactionPlaneAngle)); //first derivative
215 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
219 //-----------------------------------------------------------------------
220 void AliFlowTrackSimple::AddV4( Double_t v4,
221 Double_t reactionPlaneAngle,
222 Double_t precisionPhi,
223 Int_t maxNumberOfIterations )
225 //afterburner, adds v4, uses Newton-Raphson iteration
231 for (Int_t i=0; i<maxNumberOfIterations; i++)
233 phiprev=fPhi; //store last value for comparison
234 f = fPhi-phi0+0.5*v4*TMath::Sin(4.*(fPhi-reactionPlaneAngle));
235 fp = 1.0+2.0*v4*TMath::Cos(4.*(fPhi-reactionPlaneAngle)); //first derivative
237 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
241 //-----------------------------------------------------------------------
242 void AliFlowTrackSimple::AddV5( Double_t v5,
243 Double_t reactionPlaneAngle,
244 Double_t precisionPhi,
245 Int_t maxNumberOfIterations )
247 //afterburner, adds v4, uses Newton-Raphson iteration
253 for (Int_t i=0; i<maxNumberOfIterations; i++)
255 phiprev=fPhi; //store last value for comparison
256 f = fPhi-phi0+0.4*v5*TMath::Sin(5.*(fPhi-reactionPlaneAngle));
257 fp = 1.0+2.0*v5*TMath::Cos(5.*(fPhi-reactionPlaneAngle)); //first derivative
259 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
263 //______________________________________________________________________________
264 void AliFlowTrackSimple::AddFlow( Double_t v1,
274 Double_t precisionPhi,
275 Int_t maxNumberOfIterations )
277 //afterburner, adds v1,v2,v4 uses Newton-Raphson iteration
283 for (Int_t i=0; i<maxNumberOfIterations; i++)
285 phiprev=fPhi; //store last value for comparison
287 +2.0* v1*TMath::Sin( fPhi-rp1)
288 + v2*TMath::Sin(2.*(fPhi-rp2))
289 +2./3.*v3*TMath::Sin(3.*(fPhi-rp3))
290 +0.5* v4*TMath::Sin(4.*(fPhi-rp4))
291 +0.4* v5*TMath::Sin(5.*(fPhi-rp5))
295 +v1*TMath::Cos( fPhi-rp1)
296 +v2*TMath::Cos(2.*(fPhi-rp2))
297 +v3*TMath::Cos(3.*(fPhi-rp3))
298 +v4*TMath::Cos(4.*(fPhi-rp4))
299 +v5*TMath::Cos(5.*(fPhi-rp5))
300 ); //first derivative
302 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
306 //______________________________________________________________________________
307 void AliFlowTrackSimple::AddFlow( Double_t v1,
313 Double_t precisionPhi,
314 Int_t maxNumberOfIterations )
316 //afterburner, adds v1,v2,v4 uses Newton-Raphson iteration
317 AddFlow(v1,v2,v3,v4,v5,rp,rp,rp,rp,rp,precisionPhi,maxNumberOfIterations);
320 //______________________________________________________________________________
321 void AliFlowTrackSimple::AddFlow( Double_t v1,
325 Double_t reactionPlaneAngle,
326 Double_t precisionPhi,
327 Int_t maxNumberOfIterations )
329 //afterburner, adds v1,v2,v4 uses Newton-Raphson iteration
335 for (Int_t i=0; i<maxNumberOfIterations; i++)
337 phiprev=fPhi; //store last value for comparison
339 +2.0* v1*TMath::Sin( fPhi-reactionPlaneAngle)
340 + v2*TMath::Sin(2.*(fPhi-reactionPlaneAngle))
341 +2./3.*v3*TMath::Sin(3.*(fPhi-reactionPlaneAngle))
342 +0.5* v4*TMath::Sin(4.*(fPhi-reactionPlaneAngle))
346 +v1*TMath::Cos( fPhi-reactionPlaneAngle)
347 +v2*TMath::Cos(2.*(fPhi-reactionPlaneAngle))
348 +v3*TMath::Cos(3.*(fPhi-reactionPlaneAngle))
349 +v4*TMath::Cos(4.*(fPhi-reactionPlaneAngle))
350 ); //first derivative
352 if (TMath::AreEqualAbs(phiprev,fPhi,precisionPhi)) break;
356 //______________________________________________________________________________
357 void AliFlowTrackSimple::Print( Option_t* /*option*/ ) const
360 printf("Phi: %.3f, Eta: %+.3f, Pt: %.3f, weight: %.3f",fPhi,fEta,fPt,fTrackWeight);
361 if (InRPSelection()) printf(", RP");
362 if (InPOISelection()) printf(", POI");
363 for (Int_t i=0; i<2; i++)
365 if (InSubevent(i)) printf(", subevent %i",i);
370 //______________________________________________________________________________
371 void AliFlowTrackSimple::Clear(Option_t*)
380 fFlowBits.ResetAllBits();
381 fSubEventBits.ResetAllBits();