]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG/FLOW/Tasks/AliFlowTrackCuts.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge remote-tracking branch 'origin/master' into TPCdev
[u/mrichter/AliRoot.git] / PWG / FLOW / Tasks / AliFlowTrackCuts.cxx
CommitLineData
daf66719 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. *
86a97121 14 **********************************************************i****************/
daf66719 15
16/* $Id$ */
17
18// AliFlowTrackCuts:
3c67c769 19// Data selection for flow framework
daf66719 20//
21// origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
86a97121 22// mods: Redmer A. Bertens (rbertens@cern.ch)
5559ce24 23//
24// This class gurantees consistency of cut methods, trackparameter
25// selection (global tracks, TPC only, etc..) and parameter mixing
26// in the flow framework. Transparently handles different input types:
27// ESD, MC, AOD.
28// This class works in 2 steps: first the requested track parameters are
29// constructed (to be set by SetParamType() ), then cuts are applied.
30// the constructed track can be requested AFTER checking the cuts by
31// calling GetTrack(), in this case the cut object stays in control,
32// caller does not have to delete the track.
33// Additionally caller can request an AliFlowTrack object to be constructed
34// according the parameter mixing scenario requested by SetParamMix().
35// AliFlowTrack is made using MakeFlowTrack() method, its an 'object factory'
36// so caller needs to take care of the freshly created object.
daf66719 37
38#include <limits.h>
39#include <float.h>
32b846cd 40#include <TMatrix.h>
875ee12a 41#include "TMCProcess.h"
9a0783cc 42#include "TParticle.h"
1a80f9f6 43#include "TH2F.h"
9a0783cc 44#include "AliStack.h"
499fe731 45#include "TBrowser.h"
daf66719 46#include "AliMCEvent.h"
9a0783cc 47#include "AliESDEvent.h"
8eca5d19 48#include "AliAODEvent.h"
daf66719 49#include "AliVParticle.h"
86a97121 50#include "AliVVZERO.h"
daf66719 51#include "AliMCParticle.h"
dd294bb6
MK		 52#include "AliESDkink.h"
53#include "AliESDv0.h"
daf66719 54#include "AliESDtrack.h"
2e5052c5 55#include "AliESDMuonTrack.h" // XZhang 20120604
12b2b8bc 56#include "AliMultiplicity.h"
daf66719 57#include "AliAODTrack.h"
2e5052c5 58#include "AliAODTracklets.h" // XZhang 20120615
8eca5d19 59#include "AliFlowTrackSimple.h"
daf66719 60#include "AliFlowTrack.h"
61#include "AliFlowTrackCuts.h"
62#include "AliLog.h"
32b846cd 63#include "AliESDpid.h"
d148af7e 64#include "AliESDPmdTrack.h"
92b507bb 65#include "AliESDUtils.h" //TODO
cc0afcfc 66#include "AliFlowBayesianPID.h"
dd294bb6
MK		 67#include "AliFlowCandidateTrack.h"
68#include "AliKFParticle.h"
69#include "AliESDVZERO.h"
70#include "AliFlowCommonConstants.h"
daf66719 71
72ClassImp(AliFlowTrackCuts)
73
74//-----------------------------------------------------------------------
75AliFlowTrackCuts::AliFlowTrackCuts():
441ea1cf 76 AliFlowTrackSimpleCuts(),
77 fAliESDtrackCuts(NULL),
2e5052c5 78 fMuonTrackCuts(NULL), // XZhang 20120604
441ea1cf 79 fQA(NULL),
80 fCutMC(kFALSE),
f21bdf48 81 fCutMChasTrackReferences(kFALSE),
441ea1cf 82 fCutMCprocessType(kFALSE),
83 fMCprocessType(kPNoProcess),
84 fCutMCPID(kFALSE),
85 fMCPID(0),
2ed97fda 86 fCutMCfirstMotherPID(kFALSE),
87 fMCfirstMotherPID(0),
a14b8f3c 88 fIgnoreSignInMCPID(kFALSE),
441ea1cf 89 fCutMCisPrimary(kFALSE),
90 fRequireTransportBitForPrimaries(kTRUE),
91 fMCisPrimary(kFALSE),
92 fRequireCharge(kFALSE),
93 fFakesAreOK(kTRUE),
94 fCutSPDtrackletDeltaPhi(kFALSE),
95 fSPDtrackletDeltaPhiMax(FLT_MAX),
96 fSPDtrackletDeltaPhiMin(-FLT_MAX),
97 fIgnoreTPCzRange(kFALSE),
98 fIgnoreTPCzRangeMax(FLT_MAX),
99 fIgnoreTPCzRangeMin(-FLT_MAX),
100 fCutChi2PerClusterTPC(kFALSE),
101 fMaxChi2PerClusterTPC(FLT_MAX),
102 fMinChi2PerClusterTPC(-FLT_MAX),
103 fCutNClustersTPC(kFALSE),
104 fNClustersTPCMax(INT_MAX),
105 fNClustersTPCMin(INT_MIN),
42655d16 106 fCutNClustersITS(kFALSE),
107 fNClustersITSMax(INT_MAX),
108 fNClustersITSMin(INT_MIN),
549b5f40
RAB		 109 fUseAODFilterBit(kTRUE),
110 fAODFilterBit(1),
a63303bd 111 fCutDCAToVertexXY(kFALSE),
112 fCutDCAToVertexZ(kFALSE),
2b1eaa10 113 fCutMinimalTPCdedx(kFALSE),
114 fMinimalTPCdedx(0.),
92b507bb 115 fLinearizeVZEROresponse(kFALSE),
1a80f9f6 116 fCutPmdDet(kFALSE),
117 fPmdDet(0),
118 fCutPmdAdc(kFALSE),
119 fPmdAdc(0.),
120 fCutPmdNcell(kFALSE),
121 fPmdNcell(0.),
dd294bb6
MK		 122 fMinKinkAngle(TMath::DegToRad()*2.),
123 fMinKinkRadius(130.),
124 fMaxKinkRadius(200.),
125 fMinKinkQt(.05),
126 fMaxKinkQt(.5),
127 fMinKinkInvMassKmu(0.),
128 fMaxKinkInvMassKmu(0.6),
129 fForceTPCstandalone(kFALSE),
130 fRequireKinkDaughters(kFALSE),
441ea1cf 131 fParamType(kGlobal),
132 fParamMix(kPure),
dd294bb6
MK		 133 fKink(NULL),
134 fV0(NULL),
441ea1cf 135 fTrack(NULL),
dd294bb6
MK		 136 fTrackMass(0.),
137 fTrackPt(0.),
441ea1cf 138 fTrackPhi(0.),
139 fTrackEta(0.),
2281e7c5 140 fTrackWeight(1.),
441ea1cf 141 fTrackLabel(INT_MIN),
142 fMCevent(NULL),
143 fMCparticle(NULL),
144 fEvent(NULL),
145 fTPCtrack(),
aab6527a 146 fESDpid(),
9126fd0c 147 fBayesianResponse(NULL),
2b1eaa10 148 fPIDsource(kTOFpid),
441ea1cf 149 fTPCpidCuts(NULL),
150 fTOFpidCuts(NULL),
a14b8f3c 151 fParticleID(AliPID::kUnknown),
499fe731 152 fParticleProbability(.9),
875ee12a 153 fAllowTOFmismatchFlag(kFALSE),
0b14d0f4 154 fRequireStrictTOFTPCagreement(kFALSE),
9126fd0c 155 fCutRejectElectronsWithTPCpid(kFALSE),
8b649d8c 156 fProbBayes(0.0),
86a97121 157 fCurrCentr(0.0),
c27333c2 158 fVZEROgainEqualization(NULL),
86a97121 159 fApplyRecentering(kFALSE),
c27333c2 160 fVZEROgainEqualizationPerRing(kFALSE)
441ea1cf 161{
162 //io constructor
422c61c7 163 SetPriors(); //init arrays
b147dd2f 164
165 // New PID procedure (Bayesian Combined PID)
166 // allocating here is necessary because we don't
167 // stream this member
168 // TODO: fix streaming problems AliFlowBayesianPID
169 fBayesianResponse = new AliFlowBayesianPID();
170 fBayesianResponse->SetNewTrackParam();
86a97121 171 for(Int_t i(0); i < 4; i++) {
c27333c2
MK		 172 fVZEROApol[i] = 0;
173 fVZEROCpol[i] = 0;
86a97121 174 }
175 for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = kTRUE;
b147dd2f 176
441ea1cf 177}
178
179//-----------------------------------------------------------------------
180AliFlowTrackCuts::AliFlowTrackCuts(const char* name):
cc0afcfc 181 AliFlowTrackSimpleCuts(name),
1a80f9f6 182 fAliESDtrackCuts(NULL),
2e5052c5 183 fMuonTrackCuts(NULL), // XZhang 20120604
a1c43d26 184 fQA(NULL),
1ff4bda1 185 fCutMC(kFALSE),
f21bdf48 186 fCutMChasTrackReferences(kFALSE),
daf66719 187 fCutMCprocessType(kFALSE),
188 fMCprocessType(kPNoProcess),
189 fCutMCPID(kFALSE),
190 fMCPID(0),
2ed97fda 191 fCutMCfirstMotherPID(kFALSE),
192 fMCfirstMotherPID(0),
a14b8f3c 193 fIgnoreSignInMCPID(kFALSE),
daf66719 194 fCutMCisPrimary(kFALSE),
441ea1cf 195 fRequireTransportBitForPrimaries(kTRUE),
daf66719 196 fMCisPrimary(kFALSE),
957517fa 197 fRequireCharge(kFALSE),
127a5825 198 fFakesAreOK(kTRUE),
9a0783cc 199 fCutSPDtrackletDeltaPhi(kFALSE),
200 fSPDtrackletDeltaPhiMax(FLT_MAX),
201 fSPDtrackletDeltaPhiMin(-FLT_MAX),
1ff4bda1 202 fIgnoreTPCzRange(kFALSE),
203 fIgnoreTPCzRangeMax(FLT_MAX),
204 fIgnoreTPCzRangeMin(-FLT_MAX),
2948ac5a 205 fCutChi2PerClusterTPC(kFALSE),
206 fMaxChi2PerClusterTPC(FLT_MAX),
207 fMinChi2PerClusterTPC(-FLT_MAX),
32b846cd 208 fCutNClustersTPC(kFALSE),
209 fNClustersTPCMax(INT_MAX),
210 fNClustersTPCMin(INT_MIN),
42655d16 211 fCutNClustersITS(kFALSE),
212 fNClustersITSMax(INT_MAX),
5ba02b32 213 fNClustersITSMin(INT_MIN),
549b5f40
RAB		 214 fUseAODFilterBit(kTRUE),
215 fAODFilterBit(1),
a63303bd 216 fCutDCAToVertexXY(kFALSE),
217 fCutDCAToVertexZ(kFALSE),
2b1eaa10 218 fCutMinimalTPCdedx(kFALSE),
219 fMinimalTPCdedx(0.),
92b507bb 220 fLinearizeVZEROresponse(kFALSE),
1a80f9f6 221 fCutPmdDet(kFALSE),
222 fPmdDet(0),
223 fCutPmdAdc(kFALSE),
224 fPmdAdc(0.),
225 fCutPmdNcell(kFALSE),
226 fPmdNcell(0.),
dd294bb6
MK		 227 fMinKinkAngle(TMath::DegToRad()*2.),
228 fMinKinkRadius(130.),
229 fMaxKinkRadius(200.),
230 fMinKinkQt(0.05),
231 fMaxKinkQt(0.5),
232 fMinKinkInvMassKmu(0.0),
233 fMaxKinkInvMassKmu(0.6),
234 fForceTPCstandalone(kFALSE),
235 fRequireKinkDaughters(kFALSE),
12b2b8bc 236 fParamType(kGlobal),
daf66719 237 fParamMix(kPure),
dd294bb6
MK		 238 fKink(NULL),
239 fV0(NULL),
daf66719 240 fTrack(NULL),
dd294bb6
MK		 241 fTrackMass(0.),
242 fTrackPt(0.),
12b2b8bc 243 fTrackPhi(0.),
244 fTrackEta(0.),
2281e7c5 245 fTrackWeight(1.),
127a5825 246 fTrackLabel(INT_MIN),
957517fa 247 fMCevent(NULL),
9a0783cc 248 fMCparticle(NULL),
1ff4bda1 249 fEvent(NULL),
32b846cd 250 fTPCtrack(),
aab6527a 251 fESDpid(),
9126fd0c 252 fBayesianResponse(NULL),
2b1eaa10 253 fPIDsource(kTOFpid),
32b846cd 254 fTPCpidCuts(NULL),
255 fTOFpidCuts(NULL),
a14b8f3c 256 fParticleID(AliPID::kUnknown),
499fe731 257 fParticleProbability(.9),
875ee12a 258 fAllowTOFmismatchFlag(kFALSE),
0b14d0f4 259 fRequireStrictTOFTPCagreement(kFALSE),
9126fd0c 260 fCutRejectElectronsWithTPCpid(kFALSE),
8b649d8c 261 fProbBayes(0.0),
86a97121 262 fCurrCentr(0.0),
c27333c2 263 fVZEROgainEqualization(NULL),
86a97121 264 fApplyRecentering(kFALSE),
c27333c2 265 fVZEROgainEqualizationPerRing(kFALSE)
daf66719 266{
267 //constructor
441ea1cf 268 SetTitle("AliFlowTrackCuts");
aab6527a 269 fESDpid.GetTPCResponse().SetBetheBlochParameters( 0.0283086,
270 2.63394e+01,
271 5.04114e-11,
272 2.12543e+00,
273 4.88663e+00 );
422c61c7 274 SetPriors(); //init arrays
1a80f9f6 275
9126fd0c 276 // New PID procedure (Bayesian Combined PID)
62f2ea04 277 fBayesianResponse = new AliFlowBayesianPID();
278 fBayesianResponse->SetNewTrackParam();
86a97121 279 for(Int_t i(0); i < 4; i++) {
c27333c2
MK		 280 fVZEROApol[i] = 0;
281 fVZEROCpol[i] = 0;
86a97121 282 }
283 for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = kTRUE;
daf66719 284}
285
286//-----------------------------------------------------------------------
ee242db3 287AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
288 AliFlowTrackSimpleCuts(that),
441ea1cf 289 fAliESDtrackCuts(NULL),
2e5052c5 290 fMuonTrackCuts(NULL), // XZhang 20120604
a1c43d26 291 fQA(NULL),
1ff4bda1 292 fCutMC(that.fCutMC),
f21bdf48 293 fCutMChasTrackReferences(that.fCutMChasTrackReferences),
ee242db3 294 fCutMCprocessType(that.fCutMCprocessType),
295 fMCprocessType(that.fMCprocessType),
296 fCutMCPID(that.fCutMCPID),
297 fMCPID(that.fMCPID),
2ed97fda 298 fCutMCfirstMotherPID(that.fCutMCfirstMotherPID),
299 fMCfirstMotherPID(that.fMCfirstMotherPID),
a14b8f3c 300 fIgnoreSignInMCPID(that.fIgnoreSignInMCPID),
ee242db3 301 fCutMCisPrimary(that.fCutMCisPrimary),
441ea1cf 302 fRequireTransportBitForPrimaries(that.fRequireTransportBitForPrimaries),
ee242db3 303 fMCisPrimary(that.fMCisPrimary),
304 fRequireCharge(that.fRequireCharge),
305 fFakesAreOK(that.fFakesAreOK),
306 fCutSPDtrackletDeltaPhi(that.fCutSPDtrackletDeltaPhi),
307 fSPDtrackletDeltaPhiMax(that.fSPDtrackletDeltaPhiMax),
308 fSPDtrackletDeltaPhiMin(that.fSPDtrackletDeltaPhiMin),
1ff4bda1 309 fIgnoreTPCzRange(that.fIgnoreTPCzRange),
310 fIgnoreTPCzRangeMax(that.fIgnoreTPCzRangeMax),
311 fIgnoreTPCzRangeMin(that.fIgnoreTPCzRangeMin),
2948ac5a 312 fCutChi2PerClusterTPC(that.fCutChi2PerClusterTPC),
313 fMaxChi2PerClusterTPC(that.fMaxChi2PerClusterTPC),
314 fMinChi2PerClusterTPC(that.fMinChi2PerClusterTPC),
32b846cd 315 fCutNClustersTPC(that.fCutNClustersTPC),
316 fNClustersTPCMax(that.fNClustersTPCMax),
317 fNClustersTPCMin(that.fNClustersTPCMin),
42655d16 318 fCutNClustersITS(that.fCutNClustersITS),
319 fNClustersITSMax(that.fNClustersITSMax),
320 fNClustersITSMin(that.fNClustersITSMin),
5ba02b32 321 fUseAODFilterBit(that.fUseAODFilterBit),
322 fAODFilterBit(that.fAODFilterBit),
a63303bd 323 fCutDCAToVertexXY(that.fCutDCAToVertexXY),
324 fCutDCAToVertexZ(that.fCutDCAToVertexZ),
2b1eaa10 325 fCutMinimalTPCdedx(that.fCutMinimalTPCdedx),
326 fMinimalTPCdedx(that.fMinimalTPCdedx),
92b507bb 327 fLinearizeVZEROresponse(that.fLinearizeVZEROresponse),
1a80f9f6 328 fCutPmdDet(that.fCutPmdDet),
329 fPmdDet(that.fPmdDet),
330 fCutPmdAdc(that.fCutPmdAdc),
331 fPmdAdc(that.fPmdAdc),
332 fCutPmdNcell(that.fCutPmdNcell),
333 fPmdNcell(that.fPmdNcell),
dd294bb6
MK		 334 fMinKinkAngle(that.fMinKinkAngle),
335 fMinKinkRadius(that.fMinKinkRadius),
336 fMaxKinkRadius(that.fMaxKinkRadius),
337 fMinKinkQt(that.fMinKinkQt),
338 fMaxKinkQt(that.fMaxKinkQt),
339 fMinKinkInvMassKmu(that.fMinKinkInvMassKmu),
340 fMaxKinkInvMassKmu(that.fMaxKinkInvMassKmu),
341 fForceTPCstandalone(that.fForceTPCstandalone),
342 fRequireKinkDaughters(that.fRequireKinkDaughters),
ee242db3 343 fParamType(that.fParamType),
344 fParamMix(that.fParamMix),
dd294bb6
MK		 345 fKink(NULL),
346 fV0(NULL),
daf66719 347 fTrack(NULL),
dd294bb6
MK		 348 fTrackMass(0.),
349 fTrackPt(0.),
ee242db3 350 fTrackPhi(0.),
351 fTrackEta(0.),
2281e7c5 352 fTrackWeight(1.),
127a5825 353 fTrackLabel(INT_MIN),
957517fa 354 fMCevent(NULL),
9a0783cc 355 fMCparticle(NULL),
1ff4bda1 356 fEvent(NULL),
32b846cd 357 fTPCtrack(),
358 fESDpid(that.fESDpid),
9126fd0c 359 fBayesianResponse(NULL),
32b846cd 360 fPIDsource(that.fPIDsource),
361 fTPCpidCuts(NULL),
362 fTOFpidCuts(NULL),
2b1eaa10 363 fParticleID(that.fParticleID),
499fe731 364 fParticleProbability(that.fParticleProbability),
875ee12a 365 fAllowTOFmismatchFlag(that.fAllowTOFmismatchFlag),
0b14d0f4 366 fRequireStrictTOFTPCagreement(that.fRequireStrictTOFTPCagreement),
9126fd0c 367 fCutRejectElectronsWithTPCpid(that.fCutRejectElectronsWithTPCpid),
8b649d8c 368 fProbBayes(0.0),
86a97121 369 fCurrCentr(0.0),
c27333c2 370 fVZEROgainEqualization(NULL),
86a97121 371 fApplyRecentering(that.fApplyRecentering),
c27333c2 372 fVZEROgainEqualizationPerRing(that.fVZEROgainEqualizationPerRing)
daf66719 373{
374 //copy constructor
86a97121 375 printf(" \n\n claling copy ctor \n\n" );
32b846cd 376 if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
377 if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
441ea1cf 378 if (that.fAliESDtrackCuts) fAliESDtrackCuts = new AliESDtrackCuts(*(that.fAliESDtrackCuts));
2e5052c5 379 if (that.fMuonTrackCuts) fMuonTrackCuts = new AliMuonTrackCuts(*(that.fMuonTrackCuts)); // XZhang 20120604
422c61c7 380 SetPriors(); //init arrays
1a80f9f6 381 if (that.fQA) DefineHistograms();
9126fd0c 382
383 // New PID procedure (Bayesian Combined PID)
62f2ea04 384 fBayesianResponse = new AliFlowBayesianPID();
385 fBayesianResponse->SetNewTrackParam();
86a97121 386
c27333c2
MK		 387 // VZERO gain calibration
388 // no reason to init fVZEROgainEqualizationPerRing, will be initialized on node if necessary
86a97121 389 // pointer is set to NULL in initialization list of this constructor
c27333c2 390// if (that.fVZEROgainEqualization) fVZEROgainEqualization = new TH1(*(that.fVZEROgainEqualization));
86a97121 391 for(Int_t i(0); i < 4; i++) { // no use to copy these guys since they're only initialized on worked node
c27333c2
MK		 392 fVZEROApol[i] = 0.;
393 fVZEROCpol[i] = 0.;
86a97121 394 }
395 for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = that.fUseVZERORing[i];
396
daf66719 397}
398
399//-----------------------------------------------------------------------
ee242db3 400AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that)
daf66719 401{
402 //assignment
cdc20344 403 if (this==&that) return *this;
404
ee242db3 405 AliFlowTrackSimpleCuts::operator=(that);
1a80f9f6 406 //the following may seem excessive but if AliESDtrackCuts properly does copy and clone
407 //this approach is better memory-fragmentation-wise in some cases
408 if (that.fAliESDtrackCuts && fAliESDtrackCuts) *fAliESDtrackCuts=*(that.fAliESDtrackCuts);
409 if (that.fAliESDtrackCuts && !fAliESDtrackCuts) fAliESDtrackCuts=new AliESDtrackCuts(*(that.fAliESDtrackCuts));
410 if (!that.fAliESDtrackCuts) delete fAliESDtrackCuts; fAliESDtrackCuts=NULL;
2e5052c5 411
412 if ( that.fMuonTrackCuts && fMuonTrackCuts) *fMuonTrackCuts = *(that.fMuonTrackCuts); // XZhang 20120604
52990716 413 if ( that.fMuonTrackCuts && !fMuonTrackCuts) fMuonTrackCuts = new AliMuonTrackCuts(*(that.fMuonTrackCuts)); // XZhang 20120604
2e5052c5 414 if (!that.fMuonTrackCuts) delete fMuonTrackCuts; fMuonTrackCuts = NULL; // XZhang 20120604
c27333c2 415 if (!that.fVZEROgainEqualization) delete fVZEROgainEqualization; fVZEROgainEqualization = NULL;
1a80f9f6 416 //these guys we don't need to copy, just reinit
417 if (that.fQA) {fQA->Delete(); delete fQA; fQA=NULL; DefineHistograms();}
1ff4bda1 418 fCutMC=that.fCutMC;
f21bdf48 419 fCutMChasTrackReferences=that.fCutMChasTrackReferences;
ee242db3 420 fCutMCprocessType=that.fCutMCprocessType;
421 fMCprocessType=that.fMCprocessType;
422 fCutMCPID=that.fCutMCPID;
423 fMCPID=that.fMCPID;
2ed97fda 424 fCutMCfirstMotherPID=that.fCutMCfirstMotherPID;
425 fMCfirstMotherPID=that.fMCfirstMotherPID;
a14b8f3c 426 fIgnoreSignInMCPID=that.fIgnoreSignInMCPID,
ee242db3 427 fCutMCisPrimary=that.fCutMCisPrimary;
441ea1cf 428 fRequireTransportBitForPrimaries=that.fRequireTransportBitForPrimaries;
ee242db3 429 fMCisPrimary=that.fMCisPrimary;
430 fRequireCharge=that.fRequireCharge;
431 fFakesAreOK=that.fFakesAreOK;
432 fCutSPDtrackletDeltaPhi=that.fCutSPDtrackletDeltaPhi;
433 fSPDtrackletDeltaPhiMax=that.fSPDtrackletDeltaPhiMax;
434 fSPDtrackletDeltaPhiMin=that.fSPDtrackletDeltaPhiMin;
1ff4bda1 435 fIgnoreTPCzRange=that.fIgnoreTPCzRange;
436 fIgnoreTPCzRangeMax=that.fIgnoreTPCzRangeMax;
437 fIgnoreTPCzRangeMin=that.fIgnoreTPCzRangeMin;
2948ac5a 438 fCutChi2PerClusterTPC=that.fCutChi2PerClusterTPC;
439 fMaxChi2PerClusterTPC=that.fMaxChi2PerClusterTPC;
440 fMinChi2PerClusterTPC=that.fMinChi2PerClusterTPC;
32b846cd 441 fCutNClustersTPC=that.fCutNClustersTPC;
442 fNClustersTPCMax=that.fNClustersTPCMax;
443 fNClustersTPCMin=that.fNClustersTPCMin;
42655d16 444 fCutNClustersITS=that.fCutNClustersITS;
445 fNClustersITSMax=that.fNClustersITSMax;
446 fNClustersITSMin=that.fNClustersITSMin;
a63303bd 447 fUseAODFilterBit=that.fUseAODFilterBit;
448 fAODFilterBit=that.fAODFilterBit;
449 fCutDCAToVertexXY=that.fCutDCAToVertexXY;
450 fCutDCAToVertexZ=that.fCutDCAToVertexZ;
2b1eaa10 451 fCutMinimalTPCdedx=that.fCutMinimalTPCdedx;
452 fMinimalTPCdedx=that.fMinimalTPCdedx;
92b507bb 453 fLinearizeVZEROresponse=that.fLinearizeVZEROresponse;
1a80f9f6 454 fCutPmdDet=that.fCutPmdDet;
455 fPmdDet=that.fPmdDet;
456 fCutPmdAdc=that.fCutPmdAdc;
457 fPmdAdc=that.fPmdAdc;
458 fCutPmdNcell=that.fCutPmdNcell;
459 fPmdNcell=that.fPmdNcell;
dd294bb6
MK		 460 fMinKinkAngle=that.fMinKinkAngle;
461 fMinKinkRadius=that.fMinKinkRadius;
462 fMaxKinkRadius=that.fMaxKinkRadius;
463 fMinKinkQt=that.fMinKinkQt;
464 fMaxKinkQt=that.fMaxKinkQt;
465 fMinKinkInvMassKmu=that.fMinKinkInvMassKmu;
466 fMaxKinkInvMassKmu=that.fMaxKinkInvMassKmu;
467 fForceTPCstandalone=that.fForceTPCstandalone;
468 fRequireKinkDaughters=that.fRequireKinkDaughters;
1a80f9f6 469
ee242db3 470 fParamType=that.fParamType;
471 fParamMix=that.fParamMix;
daf66719 472
dd294bb6
MK		 473 fKink=NULL;
474 fV0=NULL;
daf66719 475 fTrack=NULL;
dd294bb6
MK		 476 fTrackMass=0.;
477 fTrackPt=0.;
1a80f9f6 478 fTrackEta=0.;
ee242db3 479 fTrackPhi=0.;
2281e7c5 480 fTrackWeight=1.;
127a5825 481 fTrackLabel=INT_MIN;
957517fa 482 fMCevent=NULL;
daf66719 483 fMCparticle=NULL;
9a0783cc 484 fEvent=NULL;
daf66719 485
32b846cd 486 fESDpid = that.fESDpid;
487 fPIDsource = that.fPIDsource;
488
cdc20344 489 delete fTPCpidCuts;
490 delete fTOFpidCuts;
32b846cd 491 if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
492 if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
32b846cd 493
2b1eaa10 494 fParticleID=that.fParticleID;
495 fParticleProbability=that.fParticleProbability;
875ee12a 496 fAllowTOFmismatchFlag=that.fAllowTOFmismatchFlag;
497 fRequireStrictTOFTPCagreement=that.fRequireStrictTOFTPCagreement;
0b14d0f4 498 fCutRejectElectronsWithTPCpid=that.fCutRejectElectronsWithTPCpid;
8b649d8c 499 fProbBayes = that.fProbBayes;
9126fd0c 500 fCurrCentr = that.fCurrCentr;
86a97121 501
502 fApplyRecentering = that.fApplyRecentering;
c27333c2 503 fVZEROgainEqualizationPerRing = that.fVZEROgainEqualizationPerRing;
93ec5f2f 504#if ROOT_VERSION_CODE < ROOT_VERSION(5,99,0)
c27333c2 505 if (that.fVZEROgainEqualization) fVZEROgainEqualization = new TH1(*(that.fVZEROgainEqualization));
93ec5f2f 506#else
507 //PH Lets try Clone, however the result might be wrong
c27333c2 508 if (that.fVZEROgainEqualization) fVZEROgainEqualization = (TH1*)that.fVZEROgainEqualization->Clone();
93ec5f2f 509#endif
86a97121 510 for(Int_t i(0); i < 4; i++) { // no use to copy these guys since they're only initialized on worked node
c27333c2
MK		 511 fVZEROApol[i] = that.fVZEROApol[i];
512 fVZEROCpol[i] = that.fVZEROCpol[i];
86a97121 513 }
514 for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = that.fUseVZERORing[i];
9126fd0c 515
daf66719 516 return *this;
517}
518
519//-----------------------------------------------------------------------
520AliFlowTrackCuts::~AliFlowTrackCuts()
521{
522 //dtor
daf66719 523 delete fAliESDtrackCuts;
32b846cd 524 delete fTPCpidCuts;
525 delete fTOFpidCuts;
2e5052c5 526 if (fMuonTrackCuts) delete fMuonTrackCuts; // XZhang 20120604
499fe731 527 if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
c27333c2 528 if (fVZEROgainEqualization) delete fVZEROgainEqualization;
daf66719 529}
530
aab6527a 531//-----------------------------------------------------------------------
532void AliFlowTrackCuts::SetEvent(AliVEvent* event, AliMCEvent* mcEvent)
533{
534 //set the event
535 Clear();
536 fEvent=event;
537 fMCevent=mcEvent;
538
539 //do the magic for ESD
540 AliESDEvent* myESD = dynamic_cast<AliESDEvent*>(event);
0d16f553 541 AliAODEvent* myAOD = dynamic_cast<AliAODEvent*>(event);
aab6527a 542 if (fCutPID && myESD)
543 {
544 //TODO: maybe call it only for the TOF options?
545 // Added by F. Noferini for TOF PID
9126fd0c 546 // old procedure now implemented inside fBayesianResponse
547 // fESDpid.MakePID(myESD,kFALSE);
548 // new procedure
549 fBayesianResponse->SetDetResponse(myESD, fCurrCentr,AliESDpid::kTOF_T0,kTRUE); // centrality = PbPb centrality class (0-100%) or -1 for pp collisions
aab6527a 550 fESDpid.SetTOFResponse(myESD,AliESDpid::kTOF_T0);
aab6527a 551 // End F. Noferini added part
552 }
0d16f553 553 if (fCutPID && myAOD){
554 fBayesianResponse->SetDetResponse(myAOD, fCurrCentr,AliESDpid::kTOF_T0); // centrality = PbPb centrality class (0-100%) or -1 for pp collisions
555 if(myAOD->GetTOFHeader()){
556 fESDpid.SetTOFResponse(myAOD,AliESDpid::kTOF_T0);
557 }
558 else{ // corrected on the fly track by track if tof header is not present (old AODs)
559 }
560 // End F. Noferini added part
561 }
562
563 if(fPIDsource==kTOFbayesian) fBayesianResponse->SetDetAND(1);
564 else if(fPIDsource==kTPCbayesian) fBayesianResponse->ResetDetOR(1);
aab6527a 565}
566
567//-----------------------------------------------------------------------
568void AliFlowTrackCuts::SetCutMC( Bool_t b )
569{
570 //will we be cutting on MC information?
571 fCutMC=b;
572
573 //if we cut on MC info then also the Bethe Bloch should be the one tuned for MC
574 if (fCutMC)
575 {
576 fESDpid.GetTPCResponse().SetBetheBlochParameters( 2.15898e+00/50.,
577 1.75295e+01,
578 3.40030e-09,
579 1.96178e+00,
580 3.91720e+00);
581 }
582}
583
daf66719 584//-----------------------------------------------------------------------
12b2b8bc 585Bool_t AliFlowTrackCuts::IsSelected(TObject* obj, Int_t id)
daf66719 586{
cac040c5 587 //check cuts
daf66719 588 AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj);
2e5052c5 589//if (vparticle) return PassesCuts(vparticle); // XZhang 20120604
590 if (vparticle) { // XZhang 20120604
591 if (fParamType==kMUON) return PassesMuonCuts(vparticle); // XZhang 20120604
592 return PassesCuts(vparticle); // XZhang 20120604
593 } // XZhang 20120604
4aae2a93 594
daf66719 595 AliFlowTrackSimple* flowtrack = dynamic_cast<AliFlowTrackSimple*>(obj);
596 if (flowtrack) return PassesCuts(flowtrack);
12b2b8bc 597 AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj);
598 if (tracklets) return PassesCuts(tracklets,id);
2e5052c5 599 AliAODTracklets* trkletAOD = dynamic_cast<AliAODTracklets*>(obj); // XZhang 20120615
600 if (trkletAOD) return PassesCuts(trkletAOD,id); // XZhang 20120615
d148af7e 601 AliESDPmdTrack* pmdtrack = dynamic_cast<AliESDPmdTrack*>(obj);
602 if (pmdtrack) return PassesPMDcuts(pmdtrack);
647a09d3 603 AliVEvent* vvzero = dynamic_cast<AliVEvent*>(obj); // should be removed; left for protection only
c27333c2 604 if (vvzero) return PassesVZEROcuts(id);
dd294bb6
MK		 605 AliESDkink* kink = dynamic_cast<AliESDkink*>(obj);
606 if (kink) return PassesCuts(kink);
607 //AliESDv0* v0 = dynamic_cast<AliESDv0*>(obj);
608 //if (v0) return PassesCuts(v0);
609
daf66719 610 return kFALSE; //default when passed wrong type of object
611}
612
1ff4bda1 613//-----------------------------------------------------------------------
614Bool_t AliFlowTrackCuts::IsSelectedMCtruth(TObject* obj, Int_t id)
615{
616 //check cuts
617 AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj);
618 if (vparticle)
619 {
2281e7c5 620 return PassesMCcuts(fMCevent,(fFakesAreOK)?TMath::Abs(vparticle->GetLabel()):vparticle->GetLabel());
1ff4bda1 621 }
622 AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj);
623 if (tracklets)
624 {
625 Int_t label0 = tracklets->GetLabel(id,0);
626 Int_t label1 = tracklets->GetLabel(id,1);
2281e7c5 627 Int_t label = (fFakesAreOK)?TMath::Abs(label0):label0;
628 if (label0!=label1) label = -666;
1ff4bda1 629 return PassesMCcuts(fMCevent,label);
630 }
631 return kFALSE; //default when passed wrong type of object
632}
633
daf66719 634//-----------------------------------------------------------------------
1a80f9f6 635Bool_t AliFlowTrackCuts::PassesCuts(const AliFlowTrackSimple* track)
daf66719 636{
637 //check cuts on a flowtracksimple
5559ce24 638
639 //clean up from last iteration
dd294bb6 640 ClearTrack();
daf66719 641 return AliFlowTrackSimpleCuts::PassesCuts(track);
642}
643
12b2b8bc 644//-----------------------------------------------------------------------
1a80f9f6 645Bool_t AliFlowTrackCuts::PassesCuts(const AliMultiplicity* tracklet, Int_t id)
12b2b8bc 646{
647 //check cuts on a tracklets
648
b35580e6 649 if (id<0) return kFALSE;
650
9a0783cc 651 //clean up from last iteration, and init label
dd294bb6 652 ClearTrack();
12b2b8bc 653
654 fTrackPhi = tracklet->GetPhi(id);
655 fTrackEta = tracklet->GetEta(id);
656 fTrackWeight = 1.0;
657 if (fCutEta) {if ( fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) return kFALSE;}
658 if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) return kFALSE;}
659
660 //check MC info if available
9a0783cc 661 //if the 2 clusters have different label track cannot be good
662 //and should therefore not pass the mc cuts
663 Int_t label0 = tracklet->GetLabel(id,0);
664 Int_t label1 = tracklet->GetLabel(id,1);
d0471ea0 665 //if possible get label and mcparticle
9a0783cc 666 fTrackLabel = (label0==label1)?tracklet->GetLabel(id,1):-1;
7afa829c 667 if (!fFakesAreOK && fTrackLabel<0) return kFALSE;
d0471ea0 668 if (fTrackLabel>=0 && fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
669 //check MC cuts
1ff4bda1 670 if (fCutMC && !PassesMCcuts()) return kFALSE;
12b2b8bc 671 return kTRUE;
672}
673
2e5052c5 674//-----------------------------------------------------------------------
675Bool_t AliFlowTrackCuts::PassesCuts(const AliAODTracklets* tracklet, Int_t id)
676{
677 // XZhang 20120615
678 //check cuts on a tracklets
679
680 if (id<0) return kFALSE;
681
682 //clean up from last iteration, and init label
dd294bb6 683 ClearTrack();
2e5052c5 684
685 fTrackPhi = tracklet->GetPhi(id);
686//fTrackEta = tracklet->GetEta(id);
687 fTrackEta = -1.*TMath::Log(TMath::Tan(tracklet->GetTheta(id)/2.));
688 fTrackWeight = 1.0;
689 if (fCutEta) {if ( fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) return kFALSE;}
690 if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) return kFALSE;}
691
692 //check MC info if available
693 //if the 2 clusters have different label track cannot be good
694 //and should therefore not pass the mc cuts
695 Int_t label0 = tracklet->GetLabel(id,0);
696 Int_t label1 = tracklet->GetLabel(id,1);
697 //if possible get label and mcparticle
698 fTrackLabel = (label0==label1)?tracklet->GetLabel(id,1):-1;
699 if (!fFakesAreOK && fTrackLabel<0) return kFALSE;
700 if (fTrackLabel>=0 && fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
701 //check MC cuts
702 if (fCutMC && !PassesMCcuts()) return kFALSE;
703 return kTRUE;
704}
705
12b2b8bc 706//-----------------------------------------------------------------------
1ff4bda1 707Bool_t AliFlowTrackCuts::PassesMCcuts(AliMCEvent* mcEvent, Int_t label)
12b2b8bc 708{
709 //check the MC info
1ff4bda1 710 if (!mcEvent) return kFALSE;
711 if (label<0) return kFALSE;//otherwise AliCMevent prints a warning before returning NULL
712 AliMCParticle* mcparticle = static_cast<AliMCParticle*>(mcEvent->GetTrack(label));
713 if (!mcparticle) {AliError("no MC track"); return kFALSE;}
12b2b8bc 714
715 if (fCutMCisPrimary)
716 {
441ea1cf 717 if (IsPhysicalPrimary(mcEvent,label,fRequireTransportBitForPrimaries) != fMCisPrimary) return kFALSE;
12b2b8bc 718 }
719 if (fCutMCPID)
720 {
1ff4bda1 721 Int_t pdgCode = mcparticle->PdgCode();
a14b8f3c 722 if (fIgnoreSignInMCPID)
4cbcbead 723 {
724 if (TMath::Abs(fMCPID) != TMath::Abs(pdgCode)) return kFALSE;
725 }
726 else
727 {
728 if (fMCPID != pdgCode) return kFALSE;
729 }
12b2b8bc 730 }
2ed97fda 731 if (fCutMCfirstMotherPID)
732 {
733
734 TParticle* tparticle=mcparticle->Particle();
735 Int_t firstMotherLabel = 0;
736 if (tparticle) { firstMotherLabel = tparticle->GetFirstMother(); }
737 AliVParticle* firstMotherParticle = mcEvent->GetTrack(firstMotherLabel);
738 Int_t pdgcodeFirstMother = 0;
739 if (firstMotherParticle) { pdgcodeFirstMother = firstMotherParticle->PdgCode(); }
740 if (pdgcodeFirstMother != fMCfirstMotherPID) return kFALSE;
741 }
12b2b8bc 742 if ( fCutMCprocessType )
743 {
1ff4bda1 744 TParticle* particle = mcparticle->Particle();
12b2b8bc 745 Int_t processID = particle->GetUniqueID();
746 if (processID != fMCprocessType ) return kFALSE;
747 }
f21bdf48 748 if (fCutMChasTrackReferences)
749 {
750 if (mcparticle->GetNumberOfTrackReferences()<1) return kFALSE;
751 }
12b2b8bc 752 return kTRUE;
753}
1a80f9f6 754
1ff4bda1 755//-----------------------------------------------------------------------
756Bool_t AliFlowTrackCuts::PassesMCcuts()
757{
1a80f9f6 758 //check MC info
1ff4bda1 759 if (!fMCevent) return kFALSE;
760 if (fTrackLabel<0) return kFALSE;//otherwise AliCMevent prints a warning before returning NULL
761 fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
762 return PassesMCcuts(fMCevent,fTrackLabel);
763}
12b2b8bc 764
dd294bb6
MK		 765//-----------------------------------------------------------------------
766Bool_t AliFlowTrackCuts::PassesCuts(const AliESDv0* v0)
767{
768 //check if the v0 passes cuts
769
770 //clean up from last iteration
771 ClearTrack();
772
773 fV0 = const_cast<AliESDv0*>(v0);
774
775 Bool_t pass = kTRUE;
776
777 //first, extract/prepare the v0
778 if (!v0->GetOnFlyStatus()) return kFALSE; //skip offline V0 finder
779 const AliExternalTrackParam *negHelix=v0->GetParamN();
780 const AliExternalTrackParam *posHelix=v0->GetParamP();
781 AliVParticle *v0tracks[2];
782 v0tracks[0] = fEvent->GetTrack(v0->GetNindex());
783 v0tracks[1] = fEvent->GetTrack(v0->GetPindex());
784 if( v0tracks[1]->Charge() < 0)
785 {
786 v0tracks[1] = fEvent->GetTrack(v0->GetNindex());
787 v0tracks[0] = fEvent->GetTrack(v0->GetPindex());
788 negHelix=v0->GetParamP();
789 posHelix=v0->GetParamN();
790 }
791
792 int KalmanPidPairs[4][2] =
793 {
794 {-11,11}, // e-,e+ (gamma)
795 {-211,2212}, // pi- p (lambda)
796 {-2212,211}, // anti-p pi+ (anti-lambda)
797 {-211,211} // pi- pi+ (Kshort)
798 // {-321,321} // K- K+ (phi)
799 };
800
801 Int_t id = 3;
802 //refit using a mass hypothesis
803 AliKFParticle v0trackKFneg(*(negHelix),KalmanPidPairs[id][0]);
804 AliKFParticle v0trackKFpos(*(posHelix),KalmanPidPairs[id][1]);
805 AliKFParticle v0particleRefit;
806 v0particleRefit += v0trackKFneg;
807 v0particleRefit += v0trackKFpos;
808 Double_t invMassErr= -999;
809 v0particleRefit.GetMass(fTrackMass,invMassErr);
810 //Double_t openAngle = v0trackKFneg.GetAngle(v0trackKFpos);
811 fTrackEta = v0particleRefit.GetEta();
812 fTrackPt = v0particleRefit.GetPt();
813 fTrackPhi = TMath::Pi()+v0particleRefit.GetPhi();
814 ////find out which mass bin and put the number in fPOItype
815 //Int_t massBins = AliFlowCommonConstants::GetMaster()->GetNbinsMass();
816 //Double_t massMin = AliFlowCommonConstants::GetMaster()->GetMassMin();
817 //Double_t massMax = AliFlowCommonConstants::GetMaster()->GetMassMax();
818 //fPOItype = 1 + int(massBins*(fTrackMass-massMin)/(massMax-massMin));
819
820 /////////////////////////////////////////////////////////////////////////////
821 //apply cuts
822 if ( v0tracks[0]->Charge() == v0tracks[1]->Charge() ) pass=kFALSE;
823 if ( v0tracks[0]->Pt()<0.15 || v0tracks[1]->Pt()<0.15 ) pass=kFALSE;
824
825 return pass;
826}
827
828//-----------------------------------------------------------------------
829Bool_t AliFlowTrackCuts::PassesCuts(const AliESDkink* kink)
830{
831 //check if the kink passes cuts
832
833 //clean up from last iteration
834 ClearTrack();
835 fKink=const_cast<AliESDkink*>(kink);
836
837 Bool_t pass = kTRUE;
838
839 Float_t kinkAngle = kink->GetAngle(2);
840 if (kinkAngle<fMinKinkAngle) pass = kFALSE;
841 Double_t kinkRadius = kink->GetR();
842 if (kinkRadius<fMinKinkRadius || kinkRadius>fMaxKinkRadius) pass = kFALSE;
843
844 //Qt
845 const TVector3 motherMfromKink(kink->GetMotherP());
846 const TVector3 daughterMfromKink(kink->GetDaughterP());
847 Float_t qt=kink->GetQt();
848 if ( qt < fMinKinkQt || qt > fMaxKinkQt) pass = kFALSE;
849
850 //invariant mass
851 Float_t energyDaughterMu = TMath::Sqrt( daughterMfromKink.Mag()*daughterMfromKink.Mag()+
852 0.105658*0.105658 );
853 Float_t p1XM = motherMfromKink.Px();
854 Float_t p1YM = motherMfromKink.Py();
855 Float_t p1ZM = motherMfromKink.Pz();
856 Float_t p2XM = daughterMfromKink.Px();
857 Float_t p2YM = daughterMfromKink.Py();
858 Float_t p2ZM = daughterMfromKink.Pz();
859 Float_t p3Daughter = TMath::Sqrt( ((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+
860 ((p1ZM-p2ZM)*(p1ZM-p2ZM)) );
861 Double_t invariantMassKmu = TMath::Sqrt( (energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-
862 motherMfromKink.Mag()*motherMfromKink.Mag() );
863 if (invariantMassKmu > fMaxKinkInvMassKmu) pass=kFALSE;
864 if (invariantMassKmu < fMinKinkInvMassKmu) pass=kFALSE;
865 fTrackMass=invariantMassKmu;
866
867 if (fQA)
868 {
869 QAbefore(13)->Fill(qt);
870 if (pass) QAafter(13)->Fill(qt);
871 QAbefore(14)->Fill(invariantMassKmu);
872 if (pass) QAafter(14)->Fill(invariantMassKmu);
873 const Double_t* kinkPosition = kink->GetPosition();
874 QAbefore(15)->Fill(kinkPosition[0],kinkPosition[1]);
875 if (pass) QAafter(15)->Fill(kinkPosition[0],kinkPosition[1]);
876 QAbefore(16)->Fill(motherMfromKink.Mag(),kinkAngle*TMath::RadToDeg());
877 if (pass) QAafter(16)->Fill(motherMfromKink.Mag(),kinkAngle*TMath::RadToDeg());
878 }
879
880 //mother track cuts
881 Int_t indexKinkMother = kink->GetIndex(0);
882 AliESDtrack* motherTrack = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(indexKinkMother));
883 if (!motherTrack) return kFALSE;
884 if (!PassesCuts(motherTrack)) pass = kFALSE;
885
886 return pass;
887}
888
daf66719 889//-----------------------------------------------------------------------
890Bool_t AliFlowTrackCuts::PassesCuts(AliVParticle* vparticle)
891{
892 //check cuts for an ESD vparticle
893
127a5825 894 ////////////////////////////////////////////////////////////////
895 // start by preparing the track parameters to cut on //////////
896 ////////////////////////////////////////////////////////////////
5559ce24 897 //clean up from last iteration
dd294bb6
MK		 898 ClearTrack();
899 Bool_t pass=kTRUE;
5559ce24 900
957517fa 901 //get the label and the mc particle
127a5825 902 fTrackLabel = (fFakesAreOK)?TMath::Abs(vparticle->GetLabel()):vparticle->GetLabel();
903 if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
daf66719 904 else fMCparticle=NULL;
905
957517fa 906 Bool_t isMCparticle = kFALSE; //some things are different for MC particles, check!
daf66719 907 AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*>(vparticle);
42655d16 908 AliAODTrack* aodTrack = NULL;
daf66719 909 if (esdTrack)
7d27a354 910 {
42655d16 911 //for an ESD track we do some magic sometimes like constructing TPC only parameters
912 //or doing some hybrid, handle that here
daf66719 913 HandleESDtrack(esdTrack);
7d27a354 914 }
daf66719 915 else
957517fa 916 {
daf66719 917 HandleVParticle(vparticle);
957517fa 918 //now check if produced particle is MC
919 isMCparticle = (dynamic_cast<AliMCParticle*>(fTrack))!=NULL;
42655d16 920 aodTrack = dynamic_cast<AliAODTrack*>(vparticle); //keep the additional dynamic cast out of the way for ESDs
957517fa 921 }
127a5825 922 ////////////////////////////////////////////////////////////////
923 ////////////////////////////////////////////////////////////////
924
1ff4bda1 925 if (!fTrack) return kFALSE;
42655d16 926 //because it may be different from global, not needed for aodTrack because we dont do anything funky there
927 if (esdTrack) esdTrack = static_cast<AliESDtrack*>(fTrack);
2948ac5a 928
9a0783cc 929 //check the common cuts for the current particle fTrack (MC,AOD,ESD)
32b846cd 930 Double_t pt = fTrack->Pt();
f21bdf48 931 Double_t p = fTrack->P();
7afa829c 932 if (!fFakesAreOK) {if (fTrackLabel<0) pass=kFALSE;}
32b846cd 933 if (fCutPt) {if (pt < fPtMin || pt >= fPtMax ) pass=kFALSE;}
924b02b0 934 if (fCutEta) {if (fTrack->Eta() < fEtaMin || fTrack->Eta() >= fEtaMax ) pass=kFALSE;}
935 if (fCutPhi) {if (fTrack->Phi() < fPhiMin || fTrack->Phi() >= fPhiMax ) pass=kFALSE;}
936 if (fRequireCharge) {if (fTrack->Charge() == 0) pass=kFALSE;}
937 if (fCutCharge && !isMCparticle) {if (fTrack->Charge() != fCharge) pass=kFALSE;}
957517fa 938 if (fCutCharge && isMCparticle)
939 {
940 //in case of an MC particle the charge is stored in units of 1/3|e|
941 Int_t charge = TMath::Nint(fTrack->Charge()/3.0); //mc particles have charge in units of 1/3e
32b846cd 942 if (charge!=fCharge) pass=kFALSE;
957517fa 943 }
924b02b0 944 //if(fCutPID) {if (fTrack->PID() != fPID) pass=kFALSE;}
daf66719 945
957517fa 946 //when additionally MC info is required
1ff4bda1 947 if (fCutMC && !PassesMCcuts()) pass=kFALSE;
daf66719 948
42655d16 949 //the case of ESD or AOD
bb6ca457 950 if (esdTrack) { if (!PassesESDcuts(esdTrack)) { pass=kFALSE; } }
60875c3c 951 if (aodTrack) { if (!PassesAODcuts(aodTrack,pass)) { pass=kFALSE; } }
42655d16 952
f21bdf48 953 if (fQA)
954 {
c3eb0b6c 955 if (fMCparticle)
f21bdf48 956 {
875ee12a 957 TParticle* tparticle=fMCparticle->Particle();
958 Int_t processID = tparticle->GetUniqueID();
8eca5d19 959 Int_t firstMotherLabel = tparticle->GetFirstMother();
3c67c769 960 Bool_t primary = IsPhysicalPrimary();
875ee12a 961 //TLorentzVector v;
962 //mcparticle->Particle()->ProductionVertex(v);
963 //Double_t prodvtxX = v.X();
964 //Double_t prodvtxY = v.Y();
965
c3eb0b6c 966 Int_t pdgcode = fMCparticle->PdgCode();
483a17da 967 AliVParticle* firstMotherParticle = fMCevent->GetTrack(firstMotherLabel);
968 Int_t pdgcodeFirstMother = 0;
969 if (firstMotherParticle) {pdgcodeFirstMother = firstMotherParticle->PdgCode();}
8eca5d19 970
483a17da 971 Float_t pdg = 0;
c3eb0b6c 972 switch (TMath::Abs(pdgcode))
973 {
974 case 11:
975 pdg = AliPID::kElectron + 0.5; break;
976 case 13:
977 pdg = AliPID::kMuon + 0.5; break;
978 case 211:
979 pdg = AliPID::kPion + 0.5; break;
980 case 321:
981 pdg = AliPID::kKaon + 0.5; break;
982 case 2212:
983 pdg = AliPID::kProton + 0.5; break;
984 default:
985 pdg = AliPID::kUnknown + 0.5; break;
986 }
987 pdg = TMath::Sign(pdg,static_cast<Float_t>(pdgcode));
8eca5d19 988
2281e7c5 989 Float_t geantCode = 0;
990 switch (pdgcodeFirstMother)
8eca5d19 991 {
2281e7c5 992 case 22: //#gamma
993 geantCode=1;
0284e4ef 994 break;
2281e7c5 995 case -11: //e^{+}
996 geantCode=2;
0284e4ef 997 break;
2281e7c5 998 case 11: //e^{-}
999 geantCode=3;
0284e4ef 1000 break;
2281e7c5 1001 case 12: case 14: case 16: //#nu
1002 geantCode=4;
1003 break;
1004 case -13:
1005 geantCode=5; //#mu^{+}
1006 break;
1007 case 13:
1008 geantCode=6; //#mu^{-}
1009 break;
1010 case 111: //#pi^{0}
1011 geantCode=7;
1012 break;
1013 case 211: //#pi^{+}
1014 geantCode=8;
1015 break;
1016 case -211: //#pi^{-}
1017 geantCode=9;
1018 break;
1019 case 130: //K^{0}_{L}
1020 geantCode=10;
1021 break;
1022 case 321: //K^{+}
1023 geantCode=11;
1024 break;
1025 case -321: //K^{-}
1026 geantCode=12;
1027 break;
1028 case 2112:
1029 geantCode = 13; //n
1030 break;
1031 case 2212:
1032 geantCode = 14; //p
1033 break;
1034 case -2212:
1035 geantCode=15; //#bar{p}
1036 break;
1037 case 310:
1038 geantCode=16; //K^{0}_{S}
1039 break;
1040 case 221:
1041 geantCode=17; //#eta
1042 break;
1043 case 3122:
1044 geantCode=18; //#Lambda
1045 break;
1046 case 3222:
1047 geantCode=19; //#Sigma^{+}
1048 break;
1049 case 3212:
1050 geantCode=20; //#Sigma^{0}
1051 break;
1052 case 3112:
1053 geantCode=21; //#Sigma^{-}
1054 break;
1055 case 3322:
1056 geantCode=22; //#Theta^{0}
1057 break;
1058 case 3312:
1059 geantCode=23; //#Theta^{-}
1060 break;
1061 case 3332:
1062 geantCode=24; //#Omega^{-}
1063 break;
1064 case -2112:
1065 geantCode=25; //#bar{n}
1066 break;
1067 case -3122:
1068 geantCode=26; //#bar{#Lambda}
1069 break;
1070 case -3112:
1071 geantCode=27; //#bar{#Sigma}^{-}
1072 break;
1073 case -3212:
1074 geantCode=28; //#bar{#Sigma}^{0}
1075 break;
1076 case -3222:
1077 geantCode=29; //#bar{#Sigma}^{+}
1078 break;
1079 case -3322:
1080 geantCode=30; //#bar{#Theta}^{0}
1081 break;
1082 case -3312:
1083 geantCode=31; //#bar{#Theta}^{+}
1084 break;
1085 case -3332:
1086 geantCode=32; //#bar{#Omega}^{+}
1087 break;
1088 case -15:
1089 geantCode=33; //#tau^{+}
1090 break;
1091 case 15:
1092 geantCode=34; //#tau^{-}
1093 break;
1094 case 411:
1095 geantCode=35; //D^{+}
1096 break;
1097 case -411:
1098 geantCode=36; //D^{-}
1099 break;
1100 case 421:
1101 geantCode=37; //D^{0}
1102 break;
1103 case -421:
1104 geantCode=38; //#bar{D}^{0}
1105 break;
1106 case 431:
1107 geantCode=39; //D_{s}^{+}
1108 break;
1109 case -431:
1110 geantCode=40; //#bar{D_{s}}^{-}
1111 break;
1112 case 4122:
1113 geantCode=41; //#Lambda_{c}^{+}
1114 break;
1115 case 24:
1116 geantCode=42; //W^{+}
1117 break;
1118 case -24:
1119 geantCode=43; //W^{-}
1120 break;
1121 case 23:
1122 geantCode=44; //Z^{0}
0284e4ef 1123 break;
8eca5d19 1124 default:
2281e7c5 1125 geantCode = 100;
0284e4ef 1126 break;
8eca5d19 1127 }
8eca5d19 1128
c3eb0b6c 1129 QAbefore(2)->Fill(p,pdg);
3c67c769 1130 QAbefore(3)->Fill(p,primary?0.5:-0.5);
875ee12a 1131 QAbefore(4)->Fill(p,static_cast<Float_t>(processID));
2281e7c5 1132 QAbefore(7)->Fill(p,geantCode+0.5);
c3eb0b6c 1133 if (pass) QAafter(2)->Fill(p,pdg);
3c67c769 1134 if (pass) QAafter(3)->Fill(p,primary?0.5:-0.5);
875ee12a 1135 if (pass) QAafter(4)->Fill(p,static_cast<Float_t>(processID));
2281e7c5 1136 if (pass) QAafter(7)->Fill(p,geantCode);
f21bdf48 1137 }
f21bdf48 1138 }
1139
42655d16 1140 //true by default, if we didn't set any cuts
1141 return pass;
1142}
1143
1144//_______________________________________________________________________
60875c3c 1145Bool_t AliFlowTrackCuts::PassesAODcuts(const AliAODTrack* track, Bool_t passedFid)
42655d16 1146{
1a80f9f6 1147 //check cuts for AOD
60875c3c 1148 Bool_t pass = passedFid;
42655d16 1149
1150 if (fCutNClustersTPC)
1ff4bda1 1151 {
42655d16 1152 Int_t ntpccls = track->GetTPCNcls();
1153 if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
1154 }
1155
1156 if (fCutNClustersITS)
1157 {
1158 Int_t nitscls = track->GetITSNcls();
1159 if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE;
1160 }
5ba02b32 1161
1162 if (fCutChi2PerClusterTPC)
1163 {
1164 Double_t chi2tpc = track->Chi2perNDF();
1165 if (chi2tpc < fMinChi2PerClusterTPC || chi2tpc > fMaxChi2PerClusterTPC) pass=kFALSE;
1166 }
1167
1168 if (GetRequireTPCRefit() && !(track->GetStatus() & AliESDtrack::kTPCrefit) ) pass=kFALSE;
1169 if (GetRequireITSRefit() && !(track->GetStatus() & AliESDtrack::kITSrefit) ) pass=kFALSE;
1170
1171 if (fUseAODFilterBit && !track->TestFilterBit(fAODFilterBit)) pass=kFALSE;
1172
a63303bd 1173 if (fCutDCAToVertexXY && track->DCA()>GetMaxDCAToVertexXY()) pass=kFALSE;
42655d16 1174
8eca5d19 1175 if (fCutDCAToVertexZ && track->ZAtDCA()>GetMaxDCAToVertexZ()) pass=kFALSE;
1176
60875c3c 1177 Double_t dedx = track->GetTPCsignal();
1178 if (dedx < fMinimalTPCdedx) pass=kFALSE;
80efa607 1179 Double_t time[9];
60875c3c 1180 track->GetIntegratedTimes(time);
1181 if (fQA) {
1182 Double_t momTPC = track->GetTPCmomentum();
1183 QAbefore( 1)->Fill(momTPC,dedx);
1184 QAbefore( 5)->Fill(track->Pt(),track->DCA());
1185 QAbefore( 6)->Fill(track->Pt(),track->ZAtDCA());
1186 if (pass) QAafter( 1)->Fill(momTPC,dedx);
1187 if (pass) QAafter( 5)->Fill(track->Pt(),track->DCA());
1188 if (pass) QAafter( 6)->Fill(track->Pt(),track->ZAtDCA());
1189 QAbefore( 8)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kElectron]));
1190 if (pass) QAafter( 8)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kElectron]));
1191 QAbefore( 9)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kMuon]));
1192 if (pass) QAafter( 9)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kMuon]));
1193 QAbefore(10)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kPion]));
1194 if (pass) QAafter( 10)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kPion]));
1195 QAbefore(11)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kKaon]));
1196 if (pass) QAafter( 11)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kKaon]));
1197 QAbefore(12)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kProton]));
1198 if (pass) QAafter( 12)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kProton]));
1199 }
8eca5d19 1200
0d16f553 1201 if (fCutPID && (fParticleID!=AliPID::kUnknown)) //if kUnknown don't cut on PID
1202 {
1203 if (!PassesAODpidCut(track)) pass=kFALSE;
1204 }
1205
42655d16 1206 return pass;
1207}
1208
1209//_______________________________________________________________________
1210Bool_t AliFlowTrackCuts::PassesESDcuts(AliESDtrack* track)
1211{
1a80f9f6 1212 //check cuts on ESD tracks
42655d16 1213 Bool_t pass=kTRUE;
cd62a2a8 1214 Float_t dcaxy=0.0;
1215 Float_t dcaz=0.0;
1216 track->GetImpactParameters(dcaxy,dcaz);
499fe731 1217 const AliExternalTrackParam* pout = track->GetOuterParam();
1218 const AliExternalTrackParam* pin = track->GetInnerParam();
42655d16 1219 if (fIgnoreTPCzRange)
1220 {
42655d16 1221 if (pin&&pout)
441ea1cf 1222 {
42655d16 1223 Double_t zin = pin->GetZ();
1224 Double_t zout = pout->GetZ();
1225 if (zin*zout<0) pass=kFALSE; //reject if cross the membrane
1226 if (zin < fIgnoreTPCzRangeMin || zin > fIgnoreTPCzRangeMax) pass=kFALSE;
1227 if (zout < fIgnoreTPCzRangeMin || zout > fIgnoreTPCzRangeMax) pass=kFALSE;
441ea1cf 1228 }
42655d16 1229 }
32b846cd 1230
1a80f9f6 1231 Int_t ntpccls = ( fParamType==kTPCstandalone )?
42655d16 1232 track->GetTPCNclsIter1():track->GetTPCNcls();
1233 if (fCutChi2PerClusterTPC)
1234 {
1a80f9f6 1235 Float_t tpcchi2 = (fParamType==kTPCstandalone)?
42655d16 1236 track->GetTPCchi2Iter1():track->GetTPCchi2();
1237 tpcchi2 = (ntpccls>0)?tpcchi2/ntpccls:-FLT_MAX;
1238 if (tpcchi2<fMinChi2PerClusterTPC || tpcchi2 >=fMaxChi2PerClusterTPC)
1239 pass=kFALSE;
1240 }
1241
2b1eaa10 1242 if (fCutMinimalTPCdedx)
1243 {
1244 if (track->GetTPCsignal() < fMinimalTPCdedx) pass=kFALSE;
1245 }
1246
42655d16 1247 if (fCutNClustersTPC)
1248 {
1249 if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
1250 }
1251
1252 Int_t nitscls = track->GetNcls(0);
1253 if (fCutNClustersITS)
1254 {
1255 if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE;
1256 }
32b846cd 1257
1a80f9f6 1258 //some stuff is still handled by AliESDtrackCuts class - delegate
1259 if (fAliESDtrackCuts)
1260 {
1261 if (!fAliESDtrackCuts->IsSelected(track)) pass=kFALSE;
1262 }
1263
cd62a2a8 1264 //PID part with pid QA
499fe731 1265 Double_t beta = GetBeta(track);
1266 Double_t dedx = Getdedx(track);
1267 if (fQA)
1268 {
1269 if (pass) QAbefore(0)->Fill(track->GetP(),beta);
1270 if (pass) QAbefore(1)->Fill(pin->GetP(),dedx);
1271 }
a14b8f3c 1272 if (fCutPID && (fParticleID!=AliPID::kUnknown)) //if kUnknown don't cut on PID
42655d16 1273 {
3c67c769 1274 if (!PassesESDpidCut(track)) pass=kFALSE;
1275 }
0b14d0f4 1276 if (fCutRejectElectronsWithTPCpid)
1277 {
1278 //reject electrons using standard TPC pid
3c67c769 1279 //TODO this should be rethought....
0b14d0f4 1280 Double_t pidTPC[AliPID::kSPECIES];
1281 track->GetTPCpid(pidTPC);
1282 if (pidTPC[AliPID::kElectron]<fParticleProbability) pass=kFALSE;
1283 }
499fe731 1284 if (fQA)
1285 {
1286 if (pass) QAafter(0)->Fill(track->GetP(),beta);
1287 if (pass) QAafter(1)->Fill(pin->GetP(),dedx);
1288 }
cd62a2a8 1289 //end pid part with qa
3c67c769 1290
cd62a2a8 1291 if (fQA)
1292 {
1293 Double_t pt = track->Pt();
1294 QAbefore(5)->Fill(pt,dcaxy);
1295 QAbefore(6)->Fill(pt,dcaz);
1296 if (pass) QAafter(5)->Fill(pt,dcaxy);
1297 if (pass) QAafter(6)->Fill(pt,dcaz);
dd294bb6
MK		 1298 QAbefore(17)->Fill(Float_t(track->GetKinkIndex(0)));
1299 if (pass) QAafter(17)->Fill(Float_t(track->GetKinkIndex(0)));
cd62a2a8 1300 }
32b846cd 1301
42655d16 1302 return pass;
daf66719 1303}
1304
1305//-----------------------------------------------------------------------
1306void AliFlowTrackCuts::HandleVParticle(AliVParticle* track)
1307{
1308 //handle the general case
daf66719 1309 switch (fParamType)
1310 {
daf66719 1311 default:
daf66719 1312 fTrack = track;
32b846cd 1313 break;
daf66719 1314 }
1315}
1316
1317//-----------------------------------------------------------------------
1318void AliFlowTrackCuts::HandleESDtrack(AliESDtrack* track)
1319{
1320 //handle esd track
daf66719 1321 switch (fParamType)
1322 {
12b2b8bc 1323 case kGlobal:
daf66719 1324 fTrack = track;
daf66719 1325 break;
1a80f9f6 1326 case kTPCstandalone:
2b1eaa10 1327 if (!track->FillTPCOnlyTrack(fTPCtrack))
1ff4bda1 1328 {
1329 fTrack=NULL;
1330 fMCparticle=NULL;
2281e7c5 1331 fTrackLabel=-998;
1ff4bda1 1332 return;
1333 }
1334 fTrack = &fTPCtrack;
957517fa 1335 //recalculate the label and mc particle, they may differ as TPClabel != global label
127a5825 1336 fTrackLabel = (fFakesAreOK)?TMath::Abs(fTrack->GetLabel()):fTrack->GetLabel();
1337 if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
957517fa 1338 else fMCparticle=NULL;
daf66719 1339 break;
daf66719 1340 default:
dd294bb6
MK		 1341 if (fForceTPCstandalone)
1342 {
1343 if (!track->FillTPCOnlyTrack(fTPCtrack))
1344 {
1345 fTrack=NULL;
1346 fMCparticle=NULL;
1347 fTrackLabel=-998;
1348 return;
1349 }
1350 fTrack = &fTPCtrack;
1351 //recalculate the label and mc particle, they may differ as TPClabel != global label
1352 fTrackLabel = (fFakesAreOK)?TMath::Abs(fTrack->GetLabel()):fTrack->GetLabel();
1353 if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
1354 else fMCparticle=NULL;
1355 }
1356 else
1357 fTrack = track;
32b846cd 1358 break;
daf66719 1359 }
1360}
1361
a0241c3a 1362//-----------------------------------------------------------------------
1363Int_t AliFlowTrackCuts::Count(AliVEvent* event)
1364{
1365 //calculate the number of track in given event.
1366 //if argument is NULL(default) take the event attached
1367 //by SetEvent()
1368 Int_t multiplicity = 0;
1369 if (!event)
1370 {
1371 for (Int_t i=0; i<GetNumberOfInputObjects(); i++)
1372 {
1373 if (IsSelected(GetInputObject(i))) multiplicity++;
1374 }
1375 }
1376 else
1377 {
1378 for (Int_t i=0; i<event->GetNumberOfTracks(); i++)
1379 {
1380 if (IsSelected(event->GetTrack(i))) multiplicity++;
1381 }
1382 }
1383 return multiplicity;
1384}
1385
1a80f9f6 1386//-----------------------------------------------------------------------
1387AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts()
54dd223b 1388{
1389 //returns the lhc10h vzero track cuts, this function
1390 //is left here for backward compatibility
1391 //if a run is recognized as 11h, the calibration method will
1392 //switch to 11h calbiration, which means that the cut
1393 //object is updated but not replaced.
1394 //calibratin is only available for PbPb runs
1395 return GetStandardVZEROOnlyTrackCuts2010();
1396}
1397//-----------------------------------------------------------------------
1398AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts2010()
1a80f9f6 1399{
c27333c2 1400 //get standard VZERO cuts
54dd223b 1401 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard vzero flow cuts 2010");
c27333c2 1402 cuts->SetParamType(kVZERO);
1a80f9f6 1403 cuts->SetEtaRange( -10, +10 );
1404 cuts->SetPhiMin( 0 );
1405 cuts->SetPhiMax( TMath::TwoPi() );
86a97121 1406 // options for the reweighting
c27333c2 1407 cuts->SetVZEROgainEqualizationPerRing(kFALSE);
54dd223b 1408 cuts->SetApplyRecentering(kTRUE);
86a97121 1409 // to exclude a ring , do e.g.
54dd223b 1410 // cuts->SetUseVZERORing(7, kFALSE);
1411 return cuts;
1412}
1413//-----------------------------------------------------------------------
1414AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts2011()
1415{
c27333c2 1416 //get standard VZERO cuts for 2011 data
54dd223b 1417 //in this case, the vzero segments will be weighted by
1418 //VZEROEqMultiplicity,
1419 //if recentering is enableded, the sub-q vectors
1420 //will be taken from the event header, so make sure to run
1421 //the VZERO event plane selection task before this task !
1422 //recentering replaces the already evaluated q-vectors, so
1423 //when chosen, additional settings (e.g. excluding rings)
1424 //have no effect. recentering is true by default
1425 //
1426 //NOTE user is responsible for running the vzero event plane
1427 //selection task in advance, e.g. add to your launcher macro
1428 //
1429 // gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskVZEROEPSelection.C");
1430 // AddTaskVZEROEPSelection();
1431 //
1432 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard vzero flow cuts 2011");
c27333c2 1433 cuts->SetParamType(kVZERO);
54dd223b 1434 cuts->SetEtaRange( -10, +10 );
1435 cuts->SetPhiMin( 0 );
1436 cuts->SetPhiMax( TMath::TwoPi() );
1437 cuts->SetApplyRecentering(kTRUE);
1a80f9f6 1438 return cuts;
1439}
a0241c3a 1440//-----------------------------------------------------------------------
1441AliFlowTrackCuts* AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()
1442{
1443 //get standard cuts
1a80f9f6 1444 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard Global tracks");
a0241c3a 1445 cuts->SetParamType(kGlobal);
1446 cuts->SetPtRange(0.2,5.);
1447 cuts->SetEtaRange(-0.8,0.8);
1448 cuts->SetMinNClustersTPC(70);
1449 cuts->SetMinChi2PerClusterTPC(0.1);
1450 cuts->SetMaxChi2PerClusterTPC(4.0);
1451 cuts->SetMinNClustersITS(2);
1452 cuts->SetRequireITSRefit(kTRUE);
1453 cuts->SetRequireTPCRefit(kTRUE);
1454 cuts->SetMaxDCAToVertexXY(0.3);
1455 cuts->SetMaxDCAToVertexZ(0.3);
1456 cuts->SetAcceptKinkDaughters(kFALSE);
2b1eaa10 1457 cuts->SetMinimalTPCdedx(10.);
a0241c3a 1458
1459 return cuts;
1460}
1461
1462//-----------------------------------------------------------------------
1a80f9f6 1463AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010()
a0241c3a 1464{
1465 //get standard cuts
1a80f9f6 1466 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPC standalone 2010");
1467 cuts->SetParamType(kTPCstandalone);
a0241c3a 1468 cuts->SetPtRange(0.2,5.);
1469 cuts->SetEtaRange(-0.8,0.8);
1470 cuts->SetMinNClustersTPC(70);
1471 cuts->SetMinChi2PerClusterTPC(0.2);
1472 cuts->SetMaxChi2PerClusterTPC(4.0);
1473 cuts->SetMaxDCAToVertexXY(3.0);
1474 cuts->SetMaxDCAToVertexZ(3.0);
1475 cuts->SetDCAToVertex2D(kTRUE);
1476 cuts->SetAcceptKinkDaughters(kFALSE);
2b1eaa10 1477 cuts->SetMinimalTPCdedx(10.);
a0241c3a 1478
1479 return cuts;
1480}
1481
daf66719 1482//-----------------------------------------------------------------------
1a80f9f6 1483AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts()
daf66719 1484{
1485 //get standard cuts
1a80f9f6 1486 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPC standalone");
1487 cuts->SetParamType(kTPCstandalone);
a0241c3a 1488 cuts->SetPtRange(0.2,5.);
1489 cuts->SetEtaRange(-0.8,0.8);
1490 cuts->SetMinNClustersTPC(70);
1491 cuts->SetMinChi2PerClusterTPC(0.2);
1492 cuts->SetMaxChi2PerClusterTPC(4.0);
1493 cuts->SetMaxDCAToVertexXY(3.0);
1494 cuts->SetMaxDCAToVertexZ(3.0);
1495 cuts->SetDCAToVertex2D(kTRUE);
1496 cuts->SetAcceptKinkDaughters(kFALSE);
2b1eaa10 1497 cuts->SetMinimalTPCdedx(10.);
a0241c3a 1498
daf66719 1499 return cuts;
1500}
1501
1502//-----------------------------------------------------------------------
1503AliFlowTrackCuts* AliFlowTrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
1504{
1505 //get standard cuts
441ea1cf 1506 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard global track cuts 2009");
daf66719 1507 delete cuts->fAliESDtrackCuts;
1508 cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selPrimaries);
12b2b8bc 1509 cuts->SetParamType(kGlobal);
daf66719 1510 return cuts;
1511}
1512
2e5052c5 1513//-----------------------------------------------------------------------------
84714904 1514AliFlowTrackCuts* AliFlowTrackCuts::GetStandardMuonTrackCuts(Bool_t isMC, Int_t passN)
2e5052c5 1515{
1516// XZhang 20120604
1517 AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard muon track cuts");
1518 cuts->SetParamType(kMUON);
1519 cuts->SetStandardMuonTrackCuts();
1520 cuts->SetIsMuonMC(isMC);
84714904 1521 cuts->SetMuonPassNumber(passN);
2e5052c5 1522 return cuts;
1523}
1524
dd294bb6
MK		 1525//-----------------------------------------------------------------------
1526//AliFlowTrack* AliFlowTrackCuts::FillFlowTrackV0(TObjArray* trackCollection, Int_t trackIndex) const
1527//{
1528// //fill flow track from a reconstructed V0 (topological)
1529// AliFlowTrack* flowtrack = static_cast<AliFlowTrack*>(trackCollection->At(trackIndex));
1530// if (flowtrack)
1531// {
1532// flowtrack->Clear();
1533// }
1534// else
1535// {
1536// flowtrack = new AliFlowCandidateTrack();
1537// trackCollection->AddAtAndExpand(flowtrack,trackIndex);
1538// }
1539//
1540// TParticle *tmpTParticle=NULL;
1541// AliMCParticle* tmpAliMCParticle=NULL;
1542// AliExternalTrackParam* externalParams=NULL;
1543// AliESDtrack* esdtrack=NULL;
1544// switch(fParamMix)
1545// {
1546// case kPure:
1547// flowtrack->Set(fTrack);
1548// break;
1549// case kTrackWithMCkine:
1550// flowtrack->Set(fMCparticle);
1551// break;
1552// case kTrackWithMCPID:
1553// flowtrack->Set(fTrack);
1554// //flowtrack->setPID(...) from mc, when implemented
1555// break;
1556// case kTrackWithMCpt:
1557// if (!fMCparticle) return NULL;
1558// flowtrack->Set(fTrack);
1559// flowtrack->SetPt(fMCparticle->Pt());
1560// break;
1561// case kTrackWithPtFromFirstMother:
1562// if (!fMCparticle) return NULL;
1563// flowtrack->Set(fTrack);
1564// tmpTParticle = fMCparticle->Particle();
1565// tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1566// flowtrack->SetPt(tmpAliMCParticle->Pt());
1567// break;
1568// case kTrackWithTPCInnerParams:
1569// esdtrack = dynamic_cast<AliESDtrack*>(fTrack);
1570// if (!esdtrack) return NULL;
1571// externalParams = const_cast<AliExternalTrackParam*>(esdtrack->GetTPCInnerParam());
1572// if (!externalParams) return NULL;
1573// flowtrack->Set(externalParams);
1574// break;
1575// default:
1576// flowtrack->Set(fTrack);
1577// break;
1578// }
1579// if (fParamType==kMC)
1580// {
1581// flowtrack->SetSource(AliFlowTrack::kFromMC);
1582// flowtrack->SetID(fTrackLabel);
1583// }
1584// else if (dynamic_cast<AliAODTrack*>(fTrack))
1585// {
1586// if (fParamType==kMUON) // XZhang 20120604
1587// flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
1588// else // XZhang 20120604
1589// flowtrack->SetSource(AliFlowTrack::kFromAOD); // XZhang 20120604
1590// flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1591// }
1592// else if (dynamic_cast<AliMCParticle*>(fTrack))
1593// {
1594// flowtrack->SetSource(AliFlowTrack::kFromMC);
1595// flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1596// }
1597//
1598// if (fV0)
1599// {
1600// //add daughter indices
1601// }
1602//
1603// return flowtrack;
1604//}
1605
1606//-----------------------------------------------------------------------
1607AliFlowTrack* AliFlowTrackCuts::FillFlowTrackKink(TObjArray* trackCollection, Int_t trackIndex) const
1608{
1609 //fill flow track from AliVParticle (ESD,AOD,MC)
1610 AliFlowTrack* flowtrack = static_cast<AliFlowTrack*>(trackCollection->At(trackIndex));
1611 if (flowtrack)
1612 {
1613 flowtrack->Clear();
1614 }
1615 else
1616 {
1617 flowtrack = new AliFlowCandidateTrack();
1618 trackCollection->AddAtAndExpand(flowtrack,trackIndex);
1619 }
1620
1621 TParticle *tmpTParticle=NULL;
1622 AliMCParticle* tmpAliMCParticle=NULL;
1623 AliExternalTrackParam* externalParams=NULL;
1624 AliESDtrack* esdtrack=NULL;
1625 switch(fParamMix)
1626 {
1627 case kPure:
1628 flowtrack->Set(fTrack);
1629 break;
1630 case kTrackWithMCkine:
1631 flowtrack->Set(fMCparticle);
1632 break;
1633 case kTrackWithMCPID:
1634 flowtrack->Set(fTrack);
1635 //flowtrack->setPID(...) from mc, when implemented
1636 break;
1637 case kTrackWithMCpt:
1638 if (!fMCparticle) return NULL;
1639 flowtrack->Set(fTrack);
1640 flowtrack->SetPt(fMCparticle->Pt());
1641 break;
1642 case kTrackWithPtFromFirstMother:
1643 if (!fMCparticle) return NULL;
1644 flowtrack->Set(fTrack);
1645 tmpTParticle = fMCparticle->Particle();
1646 tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1647 flowtrack->SetPt(tmpAliMCParticle->Pt());
1648 break;
1649 case kTrackWithTPCInnerParams:
1650 esdtrack = dynamic_cast<AliESDtrack*>(fTrack);
1651 if (!esdtrack) return NULL;
1652 externalParams = const_cast<AliExternalTrackParam*>(esdtrack->GetTPCInnerParam());
1653 if (!externalParams) return NULL;
1654 flowtrack->Set(externalParams);
1655 break;
1656 default:
1657 flowtrack->Set(fTrack);
1658 break;
1659 }
1660 if (fParamType==kMC)
1661 {
1662 flowtrack->SetSource(AliFlowTrack::kFromMC);
1663 flowtrack->SetID(fTrackLabel);
1664 }
1665 else if (dynamic_cast<AliESDtrack*>(fTrack))
1666 {
1667 flowtrack->SetSource(AliFlowTrack::kFromESD);
1668 flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1669 }
1670 else if (dynamic_cast<AliESDMuonTrack*>(fTrack)) // XZhang 20120604
1671 { // XZhang 20120604
1672 flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
1673 flowtrack->SetID((Int_t)static_cast<AliESDMuonTrack*>(fTrack)->GetUniqueID()); // XZhang 20120604
1674 } // XZhang 20120604
1675 else if (dynamic_cast<AliAODTrack*>(fTrack))
1676 {
1677 if (fParamType==kMUON) // XZhang 20120604
1678 flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
1679 else // XZhang 20120604
1680 flowtrack->SetSource(AliFlowTrack::kFromAOD); // XZhang 20120604
1681 flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1682 }
1683 else if (dynamic_cast<AliMCParticle*>(fTrack))
1684 {
1685 flowtrack->SetSource(AliFlowTrack::kFromMC);
1686 flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1687 }
1688
1689 if (fKink)
1690 {
1691 Int_t indexMother = fKink->GetIndex(0);
1692 Int_t indexDaughter = fKink->GetIndex(1);
1693 flowtrack->SetID(indexMother);
1694 flowtrack->AddDaughter(indexDaughter);
1695 flowtrack->SetMass(1.);
1696 flowtrack->SetSource(AliFlowTrack::kFromKink);
1697 }
1698
1699 flowtrack->SetMass(fTrackMass);
1700
1701 return flowtrack;
1702}
1703
1704//-----------------------------------------------------------------------
1705AliFlowTrack* AliFlowTrackCuts::FillFlowTrackGeneric(TObjArray* trackCollection, Int_t trackIndex) const
1706{
1707 //fill a flow track from tracklet,vzero,pmd,...
1708 AliFlowTrack* flowtrack = static_cast<AliFlowTrack*>(trackCollection->At(trackIndex));
1709 if (flowtrack)
1710 {
1711 flowtrack->Clear();
1712 }
1713 else
1714 {
1715 flowtrack = new AliFlowTrack();
1716 trackCollection->AddAtAndExpand(flowtrack,trackIndex);
1717 }
1718
1719 if (FillFlowTrackGeneric(flowtrack)) return flowtrack;
1720 else
1721 {
1722 trackCollection->RemoveAt(trackIndex);
1723 delete flowtrack;
1724 return NULL;
1725 }
1726}
1727
7d27a354 1728//-----------------------------------------------------------------------
1729Bool_t AliFlowTrackCuts::FillFlowTrackGeneric(AliFlowTrack* flowtrack) const
1730{
1731 //fill a flow track from tracklet,vzero,pmd,...
1732 TParticle *tmpTParticle=NULL;
1733 AliMCParticle* tmpAliMCParticle=NULL;
1734 switch (fParamMix)
1735 {
1736 case kPure:
1737 flowtrack->SetPhi(fTrackPhi);
1738 flowtrack->SetEta(fTrackEta);
2281e7c5 1739 flowtrack->SetWeight(fTrackWeight);
dd294bb6
MK		 1740 flowtrack->SetPt(fTrackPt);
1741 flowtrack->SetMass(fTrackMass);
7d27a354 1742 break;
1743 case kTrackWithMCkine:
1744 if (!fMCparticle) return kFALSE;
1745 flowtrack->SetPhi( fMCparticle->Phi() );
1746 flowtrack->SetEta( fMCparticle->Eta() );
1747 flowtrack->SetPt( fMCparticle->Pt() );
dd294bb6 1748 flowtrack->SetMass(fTrackMass);
7d27a354 1749 break;
1750 case kTrackWithMCpt:
1751 if (!fMCparticle) return kFALSE;
1752 flowtrack->SetPhi(fTrackPhi);
1753 flowtrack->SetEta(fTrackEta);
2281e7c5 1754 flowtrack->SetWeight(fTrackWeight);
7d27a354 1755 flowtrack->SetPt(fMCparticle->Pt());
dd294bb6 1756 flowtrack->SetMass(fTrackMass);
7d27a354 1757 break;
1758 case kTrackWithPtFromFirstMother:
1759 if (!fMCparticle) return kFALSE;
1760 flowtrack->SetPhi(fTrackPhi);
1761 flowtrack->SetEta(fTrackEta);
2281e7c5 1762 flowtrack->SetWeight(fTrackWeight);
7d27a354 1763 tmpTParticle = fMCparticle->Particle();
1764 tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1765 flowtrack->SetPt(tmpAliMCParticle->Pt());
dd294bb6 1766 flowtrack->SetMass(fTrackMass);
7d27a354 1767 break;
1768 default:
1769 flowtrack->SetPhi(fTrackPhi);
1770 flowtrack->SetEta(fTrackEta);
2281e7c5 1771 flowtrack->SetWeight(fTrackWeight);
dd294bb6
MK		 1772 flowtrack->SetPt(fTrackPt);
1773 flowtrack->SetMass(fTrackMass);
7d27a354 1774 break;
1775 }
1776 flowtrack->SetSource(AliFlowTrack::kFromTracklet);
1777 return kTRUE;
1778}
1779
dd294bb6
MK		 1780//-----------------------------------------------------------------------
1781AliFlowTrack* AliFlowTrackCuts::FillFlowTrackVParticle(TObjArray* trackCollection, Int_t trackIndex) const
1782{
1783 //fill flow track from AliVParticle (ESD,AOD,MC)
1784
1785 AliFlowTrack* flowtrack = static_cast<AliFlowTrack*>(trackCollection->At(trackIndex));
1786 if (flowtrack)
1787 {
1788 flowtrack->Clear();
1789 }
1790 else
1791 {
1792 flowtrack = new AliFlowTrack();
1793 trackCollection->AddAtAndExpand(flowtrack,trackIndex);
1794 }
1795
1796 if (FillFlowTrackVParticle(flowtrack)) return flowtrack;
1797 else
1798 {
1799 trackCollection->RemoveAt(trackIndex);
1800 delete flowtrack;
1801 return NULL;
1802 }
1803}
1804
7d27a354 1805//-----------------------------------------------------------------------
1806Bool_t AliFlowTrackCuts::FillFlowTrackVParticle(AliFlowTrack* flowtrack) const
1807{
1808 //fill flow track from AliVParticle (ESD,AOD,MC)
1809 if (!fTrack) return kFALSE;
1810 TParticle *tmpTParticle=NULL;
1811 AliMCParticle* tmpAliMCParticle=NULL;
1812 AliExternalTrackParam* externalParams=NULL;
1813 AliESDtrack* esdtrack=NULL;
1814 switch(fParamMix)
1815 {
1816 case kPure:
1817 flowtrack->Set(fTrack);
1818 break;
1819 case kTrackWithMCkine:
1820 flowtrack->Set(fMCparticle);
1821 break;
1822 case kTrackWithMCPID:
1823 flowtrack->Set(fTrack);
1824 //flowtrack->setPID(...) from mc, when implemented
1825 break;
1826 case kTrackWithMCpt:
1827 if (!fMCparticle) return kFALSE;
1828 flowtrack->Set(fTrack);
1829 flowtrack->SetPt(fMCparticle->Pt());
1830 break;
1831 case kTrackWithPtFromFirstMother:
1832 if (!fMCparticle) return kFALSE;
1833 flowtrack->Set(fTrack);
1834 tmpTParticle = fMCparticle->Particle();
1835 tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1836 flowtrack->SetPt(tmpAliMCParticle->Pt());
1837 break;
1838 case kTrackWithTPCInnerParams:
1839 esdtrack = dynamic_cast<AliESDtrack*>(fTrack);
1840 if (!esdtrack) return kFALSE;
1841 externalParams = const_cast<AliExternalTrackParam*>(esdtrack->GetTPCInnerParam());
1842 if (!externalParams) return kFALSE;
1843 flowtrack->Set(externalParams);
1844 break;
1845 default:
1846 flowtrack->Set(fTrack);
1847 break;
1848 }
499fe731 1849 if (fParamType==kMC)
1850 {
1851 flowtrack->SetSource(AliFlowTrack::kFromMC);
2281e7c5 1852 flowtrack->SetID(fTrackLabel);
499fe731 1853 }
1854 else if (dynamic_cast<AliESDtrack*>(fTrack))
1855 {
1856 flowtrack->SetSource(AliFlowTrack::kFromESD);
1857 flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1858 }
2e5052c5 1859 else if (dynamic_cast<AliESDMuonTrack*>(fTrack)) // XZhang 20120604
1860 { // XZhang 20120604
1861 flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
1862 flowtrack->SetID((Int_t)static_cast<AliESDMuonTrack*>(fTrack)->GetUniqueID()); // XZhang 20120604
1863 } // XZhang 20120604
1864 else if (dynamic_cast<AliAODTrack*>(fTrack))
1865 {
1866 if (fParamType==kMUON) // XZhang 20120604
1867 flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
1868 else // XZhang 20120604
1869 flowtrack->SetSource(AliFlowTrack::kFromAOD); // XZhang 20120604
499fe731 1870 flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1871 }
1872 else if (dynamic_cast<AliMCParticle*>(fTrack))
1873 {
1874 flowtrack->SetSource(AliFlowTrack::kFromMC);
1875 flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1876 }
dd294bb6 1877 flowtrack->SetMass(fTrackMass);
7d27a354 1878 return kTRUE;
1879}
1880
dd294bb6
MK		 1881//-----------------------------------------------------------------------
1882AliFlowTrack* AliFlowTrackCuts::FillFlowTrack(TObjArray* trackCollection, Int_t trackIndex) const
1883{
1884 //fill a flow track constructed from whatever we applied cuts on
1885 //return true on success
1886 switch (fParamType)
1887 {
1888 case kSPDtracklet:
1889 return FillFlowTrackGeneric(trackCollection, trackIndex);
1890 case kPMD:
1891 return FillFlowTrackGeneric(trackCollection, trackIndex);
1892 case kVZERO:
1893 return FillFlowTrackGeneric(trackCollection, trackIndex);
1894 case kKink:
1895 return FillFlowTrackKink(trackCollection, trackIndex);
1896 //case kV0:
1897 // return FillFlowTrackV0(trackCollection, trackIndex);
1898 default:
1899 return FillFlowTrackVParticle(trackCollection, trackIndex);
1900 }
1901}
1902
7d27a354 1903//-----------------------------------------------------------------------
1904Bool_t AliFlowTrackCuts::FillFlowTrack(AliFlowTrack* track) const
1905{
1906 //fill a flow track constructed from whatever we applied cuts on
1907 //return true on success
1908 switch (fParamType)
1909 {
1910 case kSPDtracklet:
1911 return FillFlowTrackGeneric(track);
1912 case kPMD:
1913 return FillFlowTrackGeneric(track);
c27333c2 1914 case kVZERO:
7d27a354 1915 return FillFlowTrackGeneric(track);
1916 default:
1917 return FillFlowTrackVParticle(track);
1918 }
1919}
1920
cac040c5
MK		 1921////-----------------------------------------------------------------------
1922//AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackSPDtracklet() const
1923//{
1924// //make a flow track from tracklet
1925// AliFlowTrack* flowtrack=NULL;
1926// TParticle *tmpTParticle=NULL;
1927// AliMCParticle* tmpAliMCParticle=NULL;
1928// switch (fParamMix)
1929// {
1930// case kPure:
1931// flowtrack = new AliFlowTrack();
1932// flowtrack->SetPhi(fTrackPhi);
1933// flowtrack->SetEta(fTrackEta);
1934// flowtrack->SetWeight(fTrackWeight);
1935// break;
1936// case kTrackWithMCkine:
1937// if (!fMCparticle) return NULL;
1938// flowtrack = new AliFlowTrack();
1939// flowtrack->SetPhi( fMCparticle->Phi() );
1940// flowtrack->SetEta( fMCparticle->Eta() );
1941// flowtrack->SetPt( fMCparticle->Pt() );
1942// break;
1943// case kTrackWithMCpt:
1944// if (!fMCparticle) return NULL;
1945// flowtrack = new AliFlowTrack();
1946// flowtrack->SetPhi(fTrackPhi);
1947// flowtrack->SetEta(fTrackEta);
1948// flowtrack->SetWeight(fTrackWeight);
1949// flowtrack->SetPt(fMCparticle->Pt());
1950// break;
1951// case kTrackWithPtFromFirstMother:
1952// if (!fMCparticle) return NULL;
1953// flowtrack = new AliFlowTrack();
1954// flowtrack->SetPhi(fTrackPhi);
1955// flowtrack->SetEta(fTrackEta);
1956// flowtrack->SetWeight(fTrackWeight);
1957// tmpTParticle = fMCparticle->Particle();
1958// tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1959// flowtrack->SetPt(tmpAliMCParticle->Pt());
1960// break;
1961// default:
1962// flowtrack = new AliFlowTrack();
1963// flowtrack->SetPhi(fTrackPhi);
1964// flowtrack->SetEta(fTrackEta);
1965// flowtrack->SetWeight(fTrackWeight);
1966// break;
1967// }
1968// flowtrack->SetSource(AliFlowTrack::kFromTracklet);
1969// return flowtrack;
1970//}
1971//
1972////-----------------------------------------------------------------------
1973//AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackVParticle() const
1974//{
1975// //make flow track from AliVParticle (ESD,AOD,MC)
1976// if (!fTrack) return NULL;
1977// AliFlowTrack* flowtrack=NULL;
1978// TParticle *tmpTParticle=NULL;
1979// AliMCParticle* tmpAliMCParticle=NULL;
1980// AliExternalTrackParam* externalParams=NULL;
1981// AliESDtrack* esdtrack=NULL;
1982// switch(fParamMix)
1983// {
1984// case kPure:
1985// flowtrack = new AliFlowTrack(fTrack);
1986// break;
1987// case kTrackWithMCkine:
1988// flowtrack = new AliFlowTrack(fMCparticle);
1989// break;
1990// case kTrackWithMCPID:
1991// flowtrack = new AliFlowTrack(fTrack);
1992// //flowtrack->setPID(...) from mc, when implemented
1993// break;
1994// case kTrackWithMCpt:
1995// if (!fMCparticle) return NULL;
1996// flowtrack = new AliFlowTrack(fTrack);
1997// flowtrack->SetPt(fMCparticle->Pt());
1998// break;
1999// case kTrackWithPtFromFirstMother:
2000// if (!fMCparticle) return NULL;
2001// flowtrack = new AliFlowTrack(fTrack);
2002// tmpTParticle = fMCparticle->Particle();
2003// tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
2004// flowtrack->SetPt(tmpAliMCParticle->Pt());
2005// break;
2006// case kTrackWithTPCInnerParams:
2007// esdtrack = dynamic_cast<AliESDtrack*>(fTrack);
2008// if (!esdtrack) return NULL;
2009// externalParams = const_cast<AliExternalTrackParam*>(esdtrack->GetTPCInnerParam());
2010// if (!externalParams) return NULL;
2011// flowtrack = new AliFlowTrack(externalParams);
2012// break;
2013// default:
2014// flowtrack = new AliFlowTrack(fTrack);
2015// break;
2016// }
2017// if (fParamType==kMC)
2018// {
2019// flowtrack->SetSource(AliFlowTrack::kFromMC);
2020// flowtrack->SetID(fTrackLabel);
2021// }
2022// else if (dynamic_cast<AliESDtrack*>(fTrack))
2023// {
2024// flowtrack->SetSource(AliFlowTrack::kFromESD);
2025// flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
2026// }
2027// else if (dynamic_cast<AliAODTrack*>(fTrack))
2028// {
2029// flowtrack->SetSource(AliFlowTrack::kFromAOD);
2030// flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
2031// }
2032// else if (dynamic_cast<AliMCParticle*>(fTrack))
2033// {
2034// flowtrack->SetSource(AliFlowTrack::kFromMC);
2035// flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
2036// }
2037// return flowtrack;
2038//}
2039//
2040////-----------------------------------------------------------------------
2041//AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackPMDtrack() const
2042//{
2043// //make a flow track from PMD track
2044// AliFlowTrack* flowtrack=NULL;
2045// TParticle *tmpTParticle=NULL;
2046// AliMCParticle* tmpAliMCParticle=NULL;
2047// switch (fParamMix)
2048// {
2049// case kPure:
2050// flowtrack = new AliFlowTrack();
2051// flowtrack->SetPhi(fTrackPhi);
2052// flowtrack->SetEta(fTrackEta);
2053// flowtrack->SetWeight(fTrackWeight);
2054// break;
2055// case kTrackWithMCkine:
2056// if (!fMCparticle) return NULL;
2057// flowtrack = new AliFlowTrack();
2058// flowtrack->SetPhi( fMCparticle->Phi() );
2059// flowtrack->SetEta( fMCparticle->Eta() );
2060// flowtrack->SetWeight(fTrackWeight);
2061// flowtrack->SetPt( fMCparticle->Pt() );
2062// break;
2063// case kTrackWithMCpt:
2064// if (!fMCparticle) return NULL;
2065// flowtrack = new AliFlowTrack();
2066// flowtrack->SetPhi(fTrackPhi);
2067// flowtrack->SetEta(fTrackEta);
2068// flowtrack->SetWeight(fTrackWeight);
2069// flowtrack->SetPt(fMCparticle->Pt());
2070// break;
2071// case kTrackWithPtFromFirstMother:
2072// if (!fMCparticle) return NULL;
2073// flowtrack = new AliFlowTrack();
2074// flowtrack->SetPhi(fTrackPhi);
2075// flowtrack->SetEta(fTrackEta);
2076// flowtrack->SetWeight(fTrackWeight);
2077// tmpTParticle = fMCparticle->Particle();
2078// tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
2079// flowtrack->SetPt(tmpAliMCParticle->Pt());
2080// break;
2081// default:
2082// flowtrack = new AliFlowTrack();
2083// flowtrack->SetPhi(fTrackPhi);
2084// flowtrack->SetEta(fTrackEta);
2085// flowtrack->SetWeight(fTrackWeight);
2086// break;
2087// }
2088//
2089// flowtrack->SetSource(AliFlowTrack::kFromPMD);
2090// return flowtrack;
2091//}
2092//
2093////-----------------------------------------------------------------------
2094//AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackVZERO() const
2095//{
2096// //make a flow track from VZERO
2097// AliFlowTrack* flowtrack=NULL;
2098// TParticle *tmpTParticle=NULL;
2099// AliMCParticle* tmpAliMCParticle=NULL;
2100// switch (fParamMix)
2101// {
2102// case kPure:
2103// flowtrack = new AliFlowTrack();
2104// flowtrack->SetPhi(fTrackPhi);
2105// flowtrack->SetEta(fTrackEta);
2106// flowtrack->SetWeight(fTrackWeight);
2107// break;
2108// case kTrackWithMCkine:
2109// if (!fMCparticle) return NULL;
2110// flowtrack = new AliFlowTrack();
2111// flowtrack->SetPhi( fMCparticle->Phi() );
2112// flowtrack->SetEta( fMCparticle->Eta() );
2113// flowtrack->SetWeight(fTrackWeight);
2114// flowtrack->SetPt( fMCparticle->Pt() );
2115// break;
2116// case kTrackWithMCpt:
2117// if (!fMCparticle) return NULL;
2118// flowtrack = new AliFlowTrack();
2119// flowtrack->SetPhi(fTrackPhi);
2120// flowtrack->SetEta(fTrackEta);
2121// flowtrack->SetWeight(fTrackWeight);
2122// flowtrack->SetPt(fMCparticle->Pt());
2123// break;
2124// case kTrackWithPtFromFirstMother:
2125// if (!fMCparticle) return NULL;
2126// flowtrack = new AliFlowTrack();
2127// flowtrack->SetPhi(fTrackPhi);
2128// flowtrack->SetEta(fTrackEta);
2129// flowtrack->SetWeight(fTrackWeight);
2130// tmpTParticle = fMCparticle->Particle();
2131// tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
2132// flowtrack->SetPt(tmpAliMCParticle->Pt());
2133// break;
2134// default:
2135// flowtrack = new AliFlowTrack();
2136// flowtrack->SetPhi(fTrackPhi);
2137// flowtrack->SetEta(fTrackEta);
2138// flowtrack->SetWeight(fTrackWeight);
2139// break;
2140// }
2141//
2142// flowtrack->SetSource(AliFlowTrack::kFromVZERO);
2143// return flowtrack;
2144//}
2145//
2146////-----------------------------------------------------------------------
2147//AliFlowTrack* AliFlowTrackCuts::MakeFlowTrack() const
2148//{
2149// //get a flow track constructed from whatever we applied cuts on
2150// //caller is resposible for deletion
2151// //if construction fails return NULL
2152// //TODO: for tracklets, PMD and VZERO we probably need just one method,
2153// //something like MakeFlowTrackGeneric(), wait with this until
2154// //requirements quirks are known.
2155// switch (fParamType)
2156// {
2157// case kSPDtracklet:
2158// return MakeFlowTrackSPDtracklet();
2159// case kPMD:
2160// return MakeFlowTrackPMDtrack();
2161// case kVZERO:
2162// return MakeFlowTrackVZERO();
2163// default:
2164// return MakeFlowTrackVParticle();
2165// }
2166//}
dd294bb6 2167
127a5825 2168//-----------------------------------------------------------------------
2169Bool_t AliFlowTrackCuts::IsPhysicalPrimary() const
2170{
2171 //check if current particle is a physical primary
9a0783cc 2172 if (!fMCevent) return kFALSE;
2173 if (fTrackLabel<0) return kFALSE;
441ea1cf 2174 return IsPhysicalPrimary(fMCevent, fTrackLabel, fRequireTransportBitForPrimaries);
9a0783cc 2175}
2176
2177//-----------------------------------------------------------------------
441ea1cf 2178Bool_t AliFlowTrackCuts::IsPhysicalPrimary(AliMCEvent* mcEvent, Int_t label, Bool_t requiretransported)
9a0783cc 2179{
2180 //check if current particle is a physical primary
2181 Bool_t physprim=mcEvent->IsPhysicalPrimary(label);
9a0783cc 2182 AliMCParticle* track = static_cast<AliMCParticle*>(mcEvent->GetTrack(label));
2183 if (!track) return kFALSE;
2184 TParticle* particle = track->Particle();
2185 Bool_t transported = particle->TestBit(kTransportBit);
441ea1cf 2186 //printf("label: %i prim: %s, transp: %s, pass: %s\n",label, (physprim)?"YES":"NO ",(transported)?"YES":"NO ",
2187 //(physprim && (transported || !requiretransported))?"YES":"NO" );
2188 return (physprim && (transported || !requiretransported));
127a5825 2189}
12b2b8bc 2190
924b02b0 2191//-----------------------------------------------------------------------
2192void AliFlowTrackCuts::DefineHistograms()
2193{
2a745a5f 2194 //define qa histograms
1a80f9f6 2195 if (fQA) return;
2196
2281e7c5 2197 const Int_t kNbinsP=200;
499fe731 2198 Double_t binsP[kNbinsP+1];
2199 binsP[0]=0.0;
6a043318 2200 for(int i=1; i<kNbinsP+1; i++)
1a80f9f6 2201 {
2281e7c5 2202 //if(binsP[i-1]+0.05<1.01)
2203 // binsP[i]=binsP[i-1]+0.05;
2204 //else
2205 binsP[i]=binsP[i-1]+0.05;
1a80f9f6 2206 }
2207
cd62a2a8 2208 const Int_t nBinsDCA=1000;
2209 Double_t binsDCA[nBinsDCA+1];
b9cf8f8e 2210 for(int i=0; i<nBinsDCA+1; i++) {binsDCA[i]=0.01*i-5.;}
cd62a2a8 2211 //for(int i=1; i<41; i++) {binsDCA[i+100]=0.1*i+1.0;}
2212
499fe731 2213 Bool_t adddirstatus = TH1::AddDirectoryStatus();
2214 TH1::AddDirectory(kFALSE);
a14b8f3c 2215 fQA=new TList(); fQA->SetOwner();
499fe731 2216 fQA->SetName(Form("%s QA",GetName()));
a14b8f3c 2217 TList* before = new TList(); before->SetOwner();
499fe731 2218 before->SetName("before");
a14b8f3c 2219 TList* after = new TList(); after->SetOwner();
499fe731 2220 after->SetName("after");
2221 fQA->Add(before);
2222 fQA->Add(after);
f21bdf48 2223 before->Add(new TH2F("TOFbeta",";p [GeV/c];#beta",kNbinsP,binsP,1000,0.4,1.1)); //0
2224 after->Add(new TH2F("TOFbeta",";p [GeV/c];#beta",kNbinsP,binsP,1000,0.4,1.1)); //0
2225 before->Add(new TH2F("TPCdedx",";p [GeV/c];dEdx",kNbinsP,binsP,500,0,500)); //1
2226 after->Add(new TH2F("TPCdedx",";p [GeV/c];dEdx",kNbinsP,binsP,500,0,500)); //1
2227 before->Add(new TH2F("MC pid",";p[GeV/c];species",kNbinsP,binsP,10,-5, 5)); //2
2228 after->Add(new TH2F("MC pid",";p[GeV/c];species",kNbinsP,binsP,10,-5, 5)); //2
483a17da 2229 //primary
2230 TH2F* hb = new TH2F("MC primary",";p[GeV/c];primary",kNbinsP,binsP,2,-1,1);
2231 TH2F* ha = new TH2F("MC primary",";p[GeV/c];primary",kNbinsP,binsP,2,-1,1);
2232 TAxis* axis = NULL;
2233 axis = hb->GetYaxis();
2234 axis->SetBinLabel(1,"secondary");
2235 axis->SetBinLabel(2,"primary");
2236 axis = ha->GetYaxis();
2237 axis->SetBinLabel(1,"secondary");
2238 axis->SetBinLabel(2,"primary");
2239 before->Add(hb); //3
2240 after->Add(ha); //3
875ee12a 2241 //production process
483a17da 2242 hb = new TH2F("MC production process",";p[GeV/c];",kNbinsP,binsP,kMaxMCProcess,
875ee12a 2243 -0.5, kMaxMCProcess-0.5);
483a17da 2244 ha = new TH2F("MC production process",";p[GeV/c];",kNbinsP,binsP,kMaxMCProcess,
875ee12a 2245 -0.5, kMaxMCProcess-0.5);
483a17da 2246 axis = hb->GetYaxis();
875ee12a 2247 for (Int_t i=0; i<kMaxMCProcess; i++)
2248 {
2249 axis->SetBinLabel(i+1,TMCProcessName[i]);
2250 }
644f1368 2251 axis = ha->GetYaxis();
875ee12a 2252 for (Int_t i=0; i<kMaxMCProcess; i++)
2253 {
2254 axis->SetBinLabel(i+1,TMCProcessName[i]);
2255 }
2256 before->Add(hb); //4
2257 after->Add(ha); //4
cd62a2a8 2258 //DCA
2259 before->Add(new TH2F("DCAxy",";p_{t}[GeV/c];DCAxy[cm]", 100, 0., 10., nBinsDCA, binsDCA));//5
2260 after->Add(new TH2F("DCAxy",";p_{t}[GeV/c];DCAxy[cm]", 100, 0., 10., nBinsDCA, binsDCA));//5
2261 before->Add(new TH2F("DCAz",";p_{t}[GeV/c];DCAz[cm]", 100, 0., 10., nBinsDCA, binsDCA));//6
2262 after->Add(new TH2F("DCAz",";p_{t}[GeV/c];DCAz[cm]", 100, 0., 10., nBinsDCA, binsDCA));//6
8eca5d19 2263 //first mother
2281e7c5 2264 hb = new TH2F("MC first mother",";p[GeV/c];",kNbinsP,binsP,44,1.,45.);
2265 ha = new TH2F("MC first mother",";p[GeV/c];",kNbinsP,binsP,44,1.,45.);
2266 hb->GetYaxis()->SetBinLabel(1, "#gamma");
2267 ha->GetYaxis()->SetBinLabel(1, "#gamma");
2268 hb->GetYaxis()->SetBinLabel(2, "e^{+}");
2269 ha->GetYaxis()->SetBinLabel(2, "e^{+}");
2270 hb->GetYaxis()->SetBinLabel(3, "e^{-}");
2271 ha->GetYaxis()->SetBinLabel(3, "e^{-}");
2272 hb->GetYaxis()->SetBinLabel(4, "#nu");
2273 ha->GetYaxis()->SetBinLabel(4, "#nu");
2274 hb->GetYaxis()->SetBinLabel(5, "#mu^{+}");
2275 ha->GetYaxis()->SetBinLabel(5, "#mu^{+}");
2276 hb->GetYaxis()->SetBinLabel(6, "#mu^{-}");
2277 ha->GetYaxis()->SetBinLabel(6, "#mu^{-}");
2278 hb->GetYaxis()->SetBinLabel(7, "#pi^{0}");
2279 ha->GetYaxis()->SetBinLabel(7, "#pi^{0}");
2280 hb->GetYaxis()->SetBinLabel(8, "#pi^{+}");
2281 ha->GetYaxis()->SetBinLabel(8, "#pi^{+}");
2282 hb->GetYaxis()->SetBinLabel(9, "#pi^{-}");
2283 ha->GetYaxis()->SetBinLabel(9, "#pi^{-}");
2284 hb->GetYaxis()->SetBinLabel(10, "K^{0}_{L}");
2285 ha->GetYaxis()->SetBinLabel(10, "K^{0}_{L}");
2286 hb->GetYaxis()->SetBinLabel(11, "K^{+}");
2287 ha->GetYaxis()->SetBinLabel(11, "K^{+}");
2288 hb->GetYaxis()->SetBinLabel(12, "K^{-}");
2289 ha->GetYaxis()->SetBinLabel(12, "K^{-}");
2290 hb->GetYaxis()->SetBinLabel( 13, "n");
2291 ha->GetYaxis()->SetBinLabel( 13, "n");
2292 hb->GetYaxis()->SetBinLabel( 14, "p");
2293 ha->GetYaxis()->SetBinLabel( 14, "p");
2294 hb->GetYaxis()->SetBinLabel(15, "#bar{p}");
2295 ha->GetYaxis()->SetBinLabel(15, "#bar{p}");
2296 hb->GetYaxis()->SetBinLabel(16, "K^{0}_{S}");
2297 ha->GetYaxis()->SetBinLabel(16, "K^{0}_{S}");
2298 hb->GetYaxis()->SetBinLabel(17, "#eta");
2299 ha->GetYaxis()->SetBinLabel(17, "#eta");
2300 hb->GetYaxis()->SetBinLabel(18, "#Lambda");
2301 ha->GetYaxis()->SetBinLabel(18, "#Lambda");
2302 hb->GetYaxis()->SetBinLabel(19, "#Sigma^{+}");
2303 ha->GetYaxis()->SetBinLabel(19, "#Sigma^{+}");
2304 hb->GetYaxis()->SetBinLabel(20, "#Sigma^{0}");
2305 ha->GetYaxis()->SetBinLabel(20, "#Sigma^{0}");
2306 hb->GetYaxis()->SetBinLabel(21, "#Sigma^{-}");
2307 ha->GetYaxis()->SetBinLabel(21, "#Sigma^{-}");
2308 hb->GetYaxis()->SetBinLabel(22, "#Theta^{0}");
2309 ha->GetYaxis()->SetBinLabel(22, "#Theta^{0}");
2310 hb->GetYaxis()->SetBinLabel(23, "#Theta^{-}");
2311 ha->GetYaxis()->SetBinLabel(23, "#Theta^{-}");
2312 hb->GetYaxis()->SetBinLabel(24, "#Omega^{-}");
2313 ha->GetYaxis()->SetBinLabel(24, "#Omega^{-}");
2314 hb->GetYaxis()->SetBinLabel(25, "#bar{n}");
2315 ha->GetYaxis()->SetBinLabel(25, "#bar{n}");
2316 hb->GetYaxis()->SetBinLabel(26, "#bar{#Lambda}");
2317 ha->GetYaxis()->SetBinLabel(26, "#bar{#Lambda}");
2318 hb->GetYaxis()->SetBinLabel(27, "#bar{#Sigma}^{-}");
2319 ha->GetYaxis()->SetBinLabel(27, "#bar{#Sigma}^{-}");
2320 hb->GetYaxis()->SetBinLabel(28, "#bar{#Sigma}^{0}");
2321 ha->GetYaxis()->SetBinLabel(28, "#bar{#Sigma}^{0}");
2322 hb->GetYaxis()->SetBinLabel(29, "#bar{#Sigma}^{+}");
2323 ha->GetYaxis()->SetBinLabel(29, "#bar{#Sigma}^{+}");
2324 hb->GetYaxis()->SetBinLabel(30, "#bar{#Theta}^{0}");
2325 ha->GetYaxis()->SetBinLabel(30, "#bar{#Theta}^{0}");
2326 hb->GetYaxis()->SetBinLabel(31, "#bar{#Theta}^{+}");
2327 ha->GetYaxis()->SetBinLabel(31, "#bar{#Theta}^{+}");
2328 hb->GetYaxis()->SetBinLabel(32, "#bar{#Omega}^{+}");
2329 ha->GetYaxis()->SetBinLabel(32, "#bar{#Omega}^{+}");
2330 hb->GetYaxis()->SetBinLabel(33, "#tau^{+}");
2331 ha->GetYaxis()->SetBinLabel(33, "#tau^{+}");
2332 hb->GetYaxis()->SetBinLabel(34, "#tau^{-}");
2333 ha->GetYaxis()->SetBinLabel(34, "#tau^{-}");
2334 hb->GetYaxis()->SetBinLabel(35, "D^{+}");
2335 ha->GetYaxis()->SetBinLabel(35, "D^{+}");
2336 hb->GetYaxis()->SetBinLabel(36, "D^{-}");
2337 ha->GetYaxis()->SetBinLabel(36, "D^{-}");
2338 hb->GetYaxis()->SetBinLabel(37, "D^{0}");
2339 ha->GetYaxis()->SetBinLabel(37, "D^{0}");
2340 hb->GetYaxis()->SetBinLabel(38, "#bar{D}^{0}");
2341 ha->GetYaxis()->SetBinLabel(38, "#bar{D}^{0}");
2342 hb->GetYaxis()->SetBinLabel(39, "D_{s}^{+}");
2343 ha->GetYaxis()->SetBinLabel(39, "D_{s}^{+}");
2344 hb->GetYaxis()->SetBinLabel(40, "#bar{D_{s}}^{-}");
2345 ha->GetYaxis()->SetBinLabel(40, "#bar{D_{s}}^{-}");
2346 hb->GetYaxis()->SetBinLabel(41, "#Lambda_{c}^{+}");
2347 ha->GetYaxis()->SetBinLabel(41, "#Lambda_{c}^{+}");
2348 hb->GetYaxis()->SetBinLabel(42, "W^{+}");
2349 ha->GetYaxis()->SetBinLabel(42, "W^{+}");
2350 hb->GetYaxis()->SetBinLabel(43, "W^{-}");
2351 ha->GetYaxis()->SetBinLabel(43, "W^{-}");
2352 hb->GetYaxis()->SetBinLabel(44, "Z^{0}");
2353 ha->GetYaxis()->SetBinLabel(44, "Z^{0}");
8eca5d19 2354 before->Add(hb);//7
2355 after->Add(ha);//7
2281e7c5 2356
60875c3c 2357 before->Add(new TH2F("TOFkElectron",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//8
2358 after->Add( new TH2F("TOFkElectron",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//8
2359
2360 before->Add(new TH2F("TOFkMuon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//9
2361 after->Add( new TH2F("TOFkMuon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//9
2362
2363 before->Add(new TH2F("TOFkPion",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//10
2364 after->Add( new TH2F("TOFkPion",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//10
2365
2366 before->Add(new TH2F("TOFkKaon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//11
2367 after->Add( new TH2F("TOFkKaon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//11
2368
2369 before->Add(new TH2F("TOFkProton",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//12
2370 after->Add( new TH2F("TOFkProton",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//12
2371
dd294bb6
MK		 2372 //kink stuff
2373 before->Add(new TH1F("KinkQt",";q_{t}[GeV/c];counts", 200, 0., 0.3));//13
2374 after->Add(new TH1F("KinkQt",";q_{t}[GeV/c];counts", 200, 0., 0.3));//13
2375
2376 before->Add(new TH1F("KinkMinv",";m_{inv}(#mu#nu)[GeV/c^{2}];counts;Kink M_{inv}", 200, 0., 0.7));//14
2377 after->Add(new TH1F("KinkMinv",";m_{inv}(#mu#nu)[GeV/c^{2}];counts;Kink M_{inv}", 200, 0., 0.7));//14
2378
2379 before->Add(new TH2F("KinkVertex",";x[cm];y[cm];Kink vertex position",250,-250.,250., 250,-250.,250.));//15
2380 after->Add(new TH2F("KinkVertex",";x[cm];y[cm];Kink vertex position",250,-250.,250., 250,-250.,250.));//15
2381
2382 before->Add(new TH2F("KinkAngleMp",";p_{mother}[GeV/c];Kink decay angle [deg];", 100,0.,6., 100,0.,80.));//16
2383 after->Add(new TH2F("KinkAngleMp",";p_{mother}[GeV/c];Kink decay angle [deg];", 100,0.,6., 100,0.,80.));//16
2384
2385 before->Add(new TH1F("KinkIndex",";Kink index;counts", 2, 0., 2.));//17
2386 after->Add(new TH1F("KinkIndex",";Kink index;counts", 2, 0., 2.));//17
499fe731 2387 TH1::AddDirectory(adddirstatus);
924b02b0 2388}
9a0783cc 2389
2390//-----------------------------------------------------------------------
2391Int_t AliFlowTrackCuts::GetNumberOfInputObjects() const
2392{
2393 //get the number of tracks in the input event according source
2394 //selection (ESD tracks, tracklets, MC particles etc.)
2395 AliESDEvent* esd=NULL;
2e5052c5 2396 AliAODEvent* aod=NULL; // XZhang 20120615
9a0783cc 2397 switch (fParamType)
2398 {
1a80f9f6 2399 case kSPDtracklet:
2e5052c5 2400 if (!fEvent) return 0; // XZhang 20120615
9a0783cc 2401 esd = dynamic_cast<AliESDEvent*>(fEvent);
2e5052c5 2402 aod = dynamic_cast<AliAODEvent*>(fEvent); // XZhang 20120615
2403// if (!esd) return 0; // XZhang 20120615
2404// return esd->GetMultiplicity()->GetNumberOfTracklets(); // XZhang 20120615
2405 if (esd) return esd->GetMultiplicity()->GetNumberOfTracklets(); // XZhang 20120615
2406 if (aod) return aod->GetTracklets()->GetNumberOfTracklets(); // XZhang 20120615
9a0783cc 2407 case kMC:
2408 if (!fMCevent) return 0;
2409 return fMCevent->GetNumberOfTracks();
d148af7e 2410 case kPMD:
2411 esd = dynamic_cast<AliESDEvent*>(fEvent);
2412 if (!esd) return 0;
2413 return esd->GetNumberOfPmdTracks();
c27333c2
MK		 2414 case kVZERO:
2415 return fgkNumberOfVZEROtracks;
2e5052c5 2416 case kMUON: // XZhang 20120604
2417 if (!fEvent) return 0; // XZhang 20120604
2418 esd = dynamic_cast<AliESDEvent*>(fEvent); // XZhang 20120604
2419 if (esd) return esd->GetNumberOfMuonTracks(); // XZhang 20120604
2420 return fEvent->GetNumberOfTracks(); // if AOD // XZhang 20120604
dd294bb6
MK		 2421 case kKink:
2422 esd = dynamic_cast<AliESDEvent*>(fEvent);
2423 if (!esd) return 0;
2424 return esd->GetNumberOfKinks();
2425 case kV0:
2426 esd = dynamic_cast<AliESDEvent*>(fEvent);
2427 if (!esd) return 0;
2428 return esd->GetNumberOfV0s();
9a0783cc 2429 default:
2430 if (!fEvent) return 0;
2431 return fEvent->GetNumberOfTracks();
2432 }
2433 return 0;
2434}
2435
2436//-----------------------------------------------------------------------
2437TObject* AliFlowTrackCuts::GetInputObject(Int_t i)
2438{
2439 //get the input object according the data source selection:
2440 //(esd tracks, traclets, mc particles,etc...)
2441 AliESDEvent* esd=NULL;
2e5052c5 2442 AliAODEvent* aod=NULL; // XZhang 20120615
9a0783cc 2443 switch (fParamType)
2444 {
1a80f9f6 2445 case kSPDtracklet:
2e5052c5 2446 if (!fEvent) return NULL; // XZhang 20120615
9a0783cc 2447 esd = dynamic_cast<AliESDEvent*>(fEvent);
2e5052c5 2448 aod = dynamic_cast<AliAODEvent*>(fEvent); // XZhang 20120615
2449// if (!esd) return NULL; // XZhang 20120615
2450// return const_cast<AliMultiplicity*>(esd->GetMultiplicity()); // XZhang 20120615
2451 if (esd) return const_cast<AliMultiplicity*>(esd->GetMultiplicity()); // XZhang 20120615
2452 if (aod) return const_cast<AliAODTracklets*>(aod->GetTracklets()); // XZhang 20120615
9a0783cc 2453 case kMC:
2454 if (!fMCevent) return NULL;
2455 return fMCevent->GetTrack(i);
d148af7e 2456 case kPMD:
2457 esd = dynamic_cast<AliESDEvent*>(fEvent);
2458 if (!esd) return NULL;
2459 return esd->GetPmdTrack(i);
c27333c2 2460 case kVZERO:
dd294bb6
MK		 2461 esd = dynamic_cast<AliESDEvent*>(fEvent);
2462 if (!esd) //contributed by G.Ortona
2463 {
2464 aod = dynamic_cast<AliAODEvent*>(fEvent);
2465 if(!aod)return NULL;
2466 return aod->GetVZEROData();
2467 }
2468 return esd->GetVZEROData();
2e5052c5 2469 case kMUON: // XZhang 20120604
2470 if (!fEvent) return NULL; // XZhang 20120604
2471 esd = dynamic_cast<AliESDEvent*>(fEvent); // XZhang 20120604
2472 if (esd) return esd->GetMuonTrack(i); // XZhang 20120604
2473 return fEvent->GetTrack(i); // if AOD // XZhang 20120604
dd294bb6
MK		 2474 case kKink:
2475 esd = dynamic_cast<AliESDEvent*>(fEvent);
2476 if (!esd) return NULL;
2477 return esd->GetKink(i);
2478 case kV0:
2479 esd = dynamic_cast<AliESDEvent*>(fEvent);
2480 if (!esd) return NULL;
2481 return esd->GetV0(i);
9a0783cc 2482 default:
2483 if (!fEvent) return NULL;
2484 return fEvent->GetTrack(i);
2485 }
2486}
a1c43d26 2487
2488//-----------------------------------------------------------------------
2489void AliFlowTrackCuts::Clear(Option_t*)
2490{
2491 //clean up
a1c43d26 2492 fMCevent=NULL;
dd294bb6
MK		 2493 fEvent=NULL;
2494 ClearTrack();
2495}
2496
2497//-----------------------------------------------------------------------
2498void AliFlowTrackCuts::ClearTrack(Option_t*)
2499{
2500 //clean up last track
2501 fKink=NULL;
2502 fV0=NULL;
2503 fTrack=NULL;
a1c43d26 2504 fMCparticle=NULL;
2281e7c5 2505 fTrackLabel=-997;
2506 fTrackWeight=1.0;
a1c43d26 2507 fTrackEta=0.0;
2508 fTrackPhi=0.0;
dd294bb6
MK		 2509 fTrackPt=0.0;
2510 fPOItype=1;
2511 fTrackMass=0.;
a1c43d26 2512}
0d16f553 2513//-----------------------------------------------------------------------
2514Bool_t AliFlowTrackCuts::PassesAODpidCut(const AliAODTrack* track )
2515{
2516 if(!track->GetAODEvent()->GetTOFHeader()){
2517 AliAODPid *pidObj = track->GetDetPid();
2518 if (!pidObj) fESDpid.GetTOFResponse().SetTimeResolution(84.);
2519 else{
2520 Double_t sigmaTOFPidInAOD[10];
2521 pidObj->GetTOFpidResolution(sigmaTOFPidInAOD);
2522 if(sigmaTOFPidInAOD[0] > 84.){
2523 fESDpid.GetTOFResponse().SetTimeResolution(sigmaTOFPidInAOD[0]); // use the electron TOF PID sigma as time resolution (including the T0 used)
2524 }
2525 }
2526 }
2527
2528 //check if passes the selected pid cut for ESDs
2529 Bool_t pass = kTRUE;
2530 switch (fPIDsource)
2531 {
2532 case kTOFbeta:
2533 if (!PassesTOFbetaCut(track)) pass=kFALSE;
2534 break;
2535 case kTOFbayesian:
2536 if (!PassesTOFbayesianCut(track)) pass=kFALSE;
2537 break;
2538 case kTPCbayesian:
2539 if (!PassesTPCbayesianCut(track)) pass=kFALSE;
2540 break;
2541 default:
2542 return kTRUE;
2543 break;
2544 }
2545 return pass;
32b846cd 2546
0d16f553 2547}
3c67c769 2548//-----------------------------------------------------------------------
2549Bool_t AliFlowTrackCuts::PassesESDpidCut(const AliESDtrack* track )
2550{
2551 //check if passes the selected pid cut for ESDs
2552 Bool_t pass = kTRUE;
2553 switch (fPIDsource)
2554 {
2555 case kTPCpid:
2556 if (!PassesTPCpidCut(track)) pass=kFALSE;
2557 break;
2558 case kTPCdedx:
2559 if (!PassesTPCdedxCut(track)) pass=kFALSE;
2560 break;
2561 case kTOFpid:
2562 if (!PassesTOFpidCut(track)) pass=kFALSE;
2563 break;
2564 case kTOFbeta:
2565 if (!PassesTOFbetaCut(track)) pass=kFALSE;
2566 break;
2567 case kTOFbetaSimple:
2568 if (!PassesTOFbetaSimpleCut(track)) pass=kFALSE;
2569 break;
2570 case kTPCbayesian:
2571 if (!PassesTPCbayesianCut(track)) pass=kFALSE;
2572 break;
2573 // part added by F. Noferini
2574 case kTOFbayesian:
2575 if (!PassesTOFbayesianCut(track)) pass=kFALSE;
2576 break;
2577 // end part added by F. Noferini
2578
2579 //part added by Natasha
2580 case kTPCNuclei:
2581 if (!PassesNucleiSelection(track)) pass=kFALSE;
2582 break;
2583 //end part added by Natasha
2584
2585 default:
2586 printf("AliFlowTrackCuts::PassesCuts() this should never be called!\n");
2587 pass=kFALSE;
2588 break;
2589 }
2590 return pass;
2591}
2592
32b846cd 2593//-----------------------------------------------------------------------
1a80f9f6 2594Bool_t AliFlowTrackCuts::PassesTOFbetaSimpleCut(const AliESDtrack* track )
32b846cd 2595{
2596 //check if passes PID cut using timing in TOF
8b649d8c 2597 Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
b092f828 2598 (track->GetTOFsignal() > 12000) &&
2599 (track->GetTOFsignal() < 100000) &&
499fe731 2600 (track->GetIntegratedLength() > 365);
2601
875ee12a 2602 if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
2603
2604 Bool_t statusMatchingHard = TPCTOFagree(track);
2605 if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
2606 return kFALSE;
aab6527a 2607
2608 if (!goodtrack) return kFALSE;
2609
a0241c3a 2610 const Float_t c = 2.99792457999999984e-02;
b092f828 2611 Float_t p = track->GetP();
1a80f9f6 2612 Float_t l = track->GetIntegratedLength();
aab6527a 2613 Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p);
b092f828 2614 Float_t timeTOF = track->GetTOFsignal()- trackT0;
1a80f9f6 2615 Float_t beta = l/timeTOF/c;
80efa607 2616 Double_t integratedTimes[9] = {-1.0,-1.0,-1.0,-1.0,-1.0, -1.0, -1.0, -1.0, -1.0};
1a80f9f6 2617 track->GetIntegratedTimes(integratedTimes);
2618 Float_t betaHypothesis[5] = {0.0,0.0,0.0,0.0,0.0};
2619 Float_t s[5] = {0.0,0.0,0.0,0.0,0.0};
2620 for (Int_t i=0;i<5;i++)
2621 {
2622 betaHypothesis[i] = l/integratedTimes[i]/c;
2623 s[i] = beta-betaHypothesis[i];
2624 }
2625
2626 switch (fParticleID)
2627 {
2628 case AliPID::kPion:
2629 return ( (s[2]<0.015) && (s[2]>-0.015) &&
7288f660 2630 (s[3]>0.025) &&
2631 (s[4]>0.03) );
1a80f9f6 2632 case AliPID::kKaon:
2633 return ( (s[3]<0.015) && (s[3]>-0.015) &&
7288f660 2634 (s[2]<-0.03) &&
2635 (s[4]>0.03) );
1a80f9f6 2636 case AliPID::kProton:
499fe731 2637 return ( (s[4]<0.015) && (s[4]>-0.015) &&
7288f660 2638 (s[3]<-0.025) &&
2639 (s[2]<-0.025) );
1a80f9f6 2640 default:
2641 return kFALSE;
2642 }
2643 return kFALSE;
2644}
2645
499fe731 2646//-----------------------------------------------------------------------
0d16f553 2647Float_t AliFlowTrackCuts::GetBeta(const AliVTrack* track)
499fe731 2648{
2649 //get beta
80efa607 2650 Double_t integratedTimes[9] = {-1.0,-1.0,-1.0,-1.0,-1.0, -1.0, -1.0, -1.0, -1.0};
0d16f553 2651 track->GetIntegratedTimes(integratedTimes);
2652
499fe731 2653 const Float_t c = 2.99792457999999984e-02;
0d16f553 2654 Float_t p = track->P();
2655 Float_t l = integratedTimes[0]*c;
499fe731 2656 Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p);
2657 Float_t timeTOF = track->GetTOFsignal()- trackT0;
2658 return l/timeTOF/c;
2659}
0d16f553 2660//-----------------------------------------------------------------------
2661Bool_t AliFlowTrackCuts::PassesTOFbetaCut(const AliAODTrack* track )
2662{
2663 //check if track passes pid selection with an asymmetric TOF beta cut
2664 if (!fTOFpidCuts)
2665 {
2666 //printf("no TOFpidCuts\n");
2667 return kFALSE;
2668 }
2669
2670 //check if passes PID cut using timing in TOF
2671 Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
2672 (track->GetTOFsignal() > 12000) &&
2673 (track->GetTOFsignal() < 100000);
2674
2675 if (!goodtrack) return kFALSE;
2676
2677 const Float_t c = 2.99792457999999984e-02;
80efa607 2678 Double_t integratedTimes[9] = {-1.0,-1.0,-1.0,-1.0,-1.0, -1.0, -1.0, -1.0, -1.0};
0d16f553 2679 track->GetIntegratedTimes(integratedTimes);
2680 Float_t l = integratedTimes[0]*c;
499fe731 2681
0d16f553 2682 goodtrack = goodtrack && (l > 365);
2683
2684 if (!goodtrack) return kFALSE;
2685
2686 if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
2687
2688 Bool_t statusMatchingHard = TPCTOFagree(track);
2689 if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
2690 return kFALSE;
2691
2692
2693 Float_t beta = GetBeta(track);
2694
2695 //construct the pid index because it's not AliPID::EParticleType
2696 Int_t pid = 0;
2697 switch (fParticleID)
2698 {
2699 case AliPID::kPion:
2700 pid=2;
2701 break;
2702 case AliPID::kKaon:
2703 pid=3;
2704 break;
2705 case AliPID::kProton:
2706 pid=4;
2707 break;
2708 default:
2709 return kFALSE;
2710 }
2711
2712 //signal to cut on
2713 Float_t p = track->P();
2714 Float_t betahypothesis = l/integratedTimes[pid]/c;
2715 Float_t betadiff = beta-betahypothesis;
2716
2717 Float_t* arr = fTOFpidCuts->GetMatrixArray();
2718 Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(p));
2719 if (col<0) return kFALSE;
2720 Float_t min = (*fTOFpidCuts)(1,col);
2721 Float_t max = (*fTOFpidCuts)(2,col);
2722
2723 Bool_t pass = (betadiff>min && betadiff<max);
2724
2725 return pass;
2726}
1a80f9f6 2727//-----------------------------------------------------------------------
2728Bool_t AliFlowTrackCuts::PassesTOFbetaCut(const AliESDtrack* track )
2729{
3c67c769 2730 //check if track passes pid selection with an asymmetric TOF beta cut
9b692328 2731 if (!fTOFpidCuts)
2732 {
2733 //printf("no TOFpidCuts\n");
2734 return kFALSE;
2735 }
2736
1a80f9f6 2737 //check if passes PID cut using timing in TOF
a14b8f3c 2738 Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
1a80f9f6 2739 (track->GetTOFsignal() > 12000) &&
2740 (track->GetTOFsignal() < 100000) &&
499fe731 2741 (track->GetIntegratedLength() > 365);
2742
875ee12a 2743 if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
2744
0d16f553 2745 if (!goodtrack) return kFALSE;
2746
875ee12a 2747 Bool_t statusMatchingHard = TPCTOFagree(track);
2748 if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
2749 return kFALSE;
1a80f9f6 2750
499fe731 2751 Float_t beta = GetBeta(track);
80efa607 2752 Double_t integratedTimes[9] = {-1.0,-1.0,-1.0,-1.0,-1.0, -1.0, -1.0, -1.0, -1.0};
a0241c3a 2753 track->GetIntegratedTimes(integratedTimes);
2754
2755 //construct the pid index because it's not AliPID::EParticleType
32b846cd 2756 Int_t pid = 0;
2b1eaa10 2757 switch (fParticleID)
32b846cd 2758 {
2759 case AliPID::kPion:
2760 pid=2;
2761 break;
2762 case AliPID::kKaon:
2763 pid=3;
2764 break;
2765 case AliPID::kProton:
2766 pid=4;
2767 break;
2768 default:
2769 return kFALSE;
2770 }
a0241c3a 2771
2772 //signal to cut on
499fe731 2773 const Float_t c = 2.99792457999999984e-02;
2774 Float_t l = track->GetIntegratedLength();
2775 Float_t p = track->GetP();
1a80f9f6 2776 Float_t betahypothesis = l/integratedTimes[pid]/c;
2777 Float_t betadiff = beta-betahypothesis;
32b846cd 2778
2779 Float_t* arr = fTOFpidCuts->GetMatrixArray();
a0241c3a 2780 Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(p));
32b846cd 2781 if (col<0) return kFALSE;
2782 Float_t min = (*fTOFpidCuts)(1,col);
2783 Float_t max = (*fTOFpidCuts)(2,col);
2784
1a80f9f6 2785 Bool_t pass = (betadiff>min && betadiff<max);
2786
1a80f9f6 2787 return pass;
32b846cd 2788}
2789
2b1eaa10 2790//-----------------------------------------------------------------------
7d27a354 2791Bool_t AliFlowTrackCuts::PassesTOFpidCut(const AliESDtrack* track) const
2b1eaa10 2792{
2793 //check if passes PID cut using default TOF pid
7d27a354 2794 Double_t pidTOF[AliPID::kSPECIES];
2795 track->GetTOFpid(pidTOF);
2796 if (pidTOF[fParticleID]>=fParticleProbability) return kTRUE;
2b1eaa10 2797 return kFALSE;
2798}
2799
32b846cd 2800//-----------------------------------------------------------------------
1a80f9f6 2801Bool_t AliFlowTrackCuts::PassesTPCpidCut(const AliESDtrack* track) const
2b1eaa10 2802{
2803 //check if passes PID cut using default TPC pid
2804 Double_t pidTPC[AliPID::kSPECIES];
2805 track->GetTPCpid(pidTPC);
2806 Double_t probablity = 0.;
2807 switch (fParticleID)
2808 {
2809 case AliPID::kPion:
2810 probablity = pidTPC[AliPID::kPion] + pidTPC[AliPID::kMuon];
2811 break;
2812 default:
2813 probablity = pidTPC[fParticleID];
2814 }
2815 if (probablity >= fParticleProbability) return kTRUE;
2816 return kFALSE;
2817}
2818
499fe731 2819//-----------------------------------------------------------------------
3c67c769 2820Float_t AliFlowTrackCuts::Getdedx(const AliESDtrack* track) const
499fe731 2821{
2822 //get TPC dedx
2823 return track->GetTPCsignal();
2824}
2825
2b1eaa10 2826//-----------------------------------------------------------------------
1a80f9f6 2827Bool_t AliFlowTrackCuts::PassesTPCdedxCut(const AliESDtrack* track)
32b846cd 2828{
2829 //check if passes PID cut using dedx signal in the TPC
32b846cd 2830 if (!fTPCpidCuts)
2831 {
9b692328 2832 //printf("no TPCpidCuts\n");
32b846cd 2833 return kFALSE;
2834 }
2835
3c75bbd5 2836 const AliExternalTrackParam* tpcparam = track->GetInnerParam(); //tpc only params at the inner wall
2837 if (!tpcparam) return kFALSE;
d79d61fb 2838 Double_t p = tpcparam->GetP();
2839 Float_t sigExp = fESDpid.GetTPCResponse().GetExpectedSignal(p, fParticleID);
32b846cd 2840 Float_t sigTPC = track->GetTPCsignal();
2841 Float_t s = (sigTPC-sigExp)/sigExp;
32b846cd 2842
2843 Float_t* arr = fTPCpidCuts->GetMatrixArray();
1a80f9f6 2844 Int_t arrSize = fTPCpidCuts->GetNcols();
2845 Int_t col = TMath::BinarySearch( arrSize, arr, static_cast<Float_t>(p));
32b846cd 2846 if (col<0) return kFALSE;
2847 Float_t min = (*fTPCpidCuts)(1,col);
2848 Float_t max = (*fTPCpidCuts)(2,col);
2849
2850 //printf("------------TPC pid cut %s\n",(s>min && s<max)?"PASS":"FAIL");
2851 return (s>min && s<max);
2852}
2853
2854//-----------------------------------------------------------------------
2855void AliFlowTrackCuts::InitPIDcuts()
2856{
2857 //init matrices with PID cuts
2858 TMatrixF* t = NULL;
2859 if (!fTPCpidCuts)
2860 {
2b1eaa10 2861 if (fParticleID==AliPID::kPion)
32b846cd 2862 {
d79d61fb 2863 t = new TMatrixF(3,15);
2864 (*t)(0,0) = 0.20; (*t)(1,0) = -0.4; (*t)(2,0) = 0.0;
2865 (*t)(0,1) = 0.25; (*t)(1,1) = -0.4; (*t)(2,1) = 0.1;
2866 (*t)(0,2) = 0.30; (*t)(1,2) = -0.4; (*t)(2,2) = 0.2;
2867 (*t)(0,3) = 0.35; (*t)(1,3) = -0.4; (*t)(2,3) = 0.2;
32b846cd 2868 (*t)(0,4) = 0.40; (*t)(1,4) = -0.4; (*t)(2,4) = 0.3;
2869 (*t)(0,5) = 0.45; (*t)(1,5) = -0.4; (*t)(2,5) = 0.3;
2870 (*t)(0,6) = 0.50; (*t)(1,6) = -0.4; (*t)(2,6) = 0.25;
2871 (*t)(0,7) = 0.55; (*t)(1,7) = -0.4; (*t)(2,7) = 0.15;
2872 (*t)(0,8) = 0.60; (*t)(1,8) = -0.4; (*t)(2,8) = 0.1;
d79d61fb 2873 (*t)(0,9) = 0.65; (*t)(1,9) = -0.4; (*t)(2,9) = 0.05;
2874 (*t)(0,10) = 0.70; (*t)(1,10) = -0.4; (*t)(2,10) = 0;
2875 (*t)(0,11) = 0.75; (*t)(1,11) = -0.4; (*t)(2,11) = 0;
2876 (*t)(0,12) = 0.80; (*t)(1,12) = -0.4; (*t)(2,12) = -0.05;
2877 (*t)(0,13) = 0.85; (*t)(1,13) = -0.4; (*t)(2,13) = -0.1;
2878 (*t)(0,14) = 0.90; (*t)(1,14) = 0; (*t)(2,14) = 0;
32b846cd 2879 }
2880 else
2b1eaa10 2881 if (fParticleID==AliPID::kKaon)
3ca4688a 2882 {
d79d61fb 2883 t = new TMatrixF(3,12);
2884 (*t)(0,0) = 0.20; (*t)(1,0) = -0.2; (*t)(2,0) = 0.2;
2885 (*t)(0,1) = 0.25; (*t)(1,1) = -0.2; (*t)(2,1) = 0.2;
2886 (*t)(0,2) = 0.30; (*t)(1,2) = -0.2; (*t)(2,2) = 0.2;
2887 (*t)(0,3) = 0.35; (*t)(1,3) = -0.2; (*t)(2,3) = 0.2;
2888 (*t)(0,4) = 0.40; (*t)(1,4) = -0.1; (*t)(2,4) = 0.2;
2889 (*t)(0,5) = 0.45; (*t)(1,5) = -0.1; (*t)(2,5) = 0.2;
2890 (*t)(0,6) = 0.50; (*t)(1,6) =-0.05; (*t)(2,6) = 0.2;
2891 (*t)(0,7) = 0.55; (*t)(1,7) = -0.1; (*t)(2,7) = 0.1;
2892 (*t)(0,8) = 0.60; (*t)(1,8) =-0.05; (*t)(2,8) = 0.1;
2893 (*t)(0,9) = 0.65; (*t)(1,9) = 0; (*t)(2,9) = 0.15;
2894 (*t)(0,10) = 0.70; (*t)(1,10) = 0.05; (*t)(2,10) = 0.2;
1a80f9f6 2895 (*t)(0,11) = 0.75; (*t)(1,11) = 0; (*t)(2,11) = 0;
3ca4688a 2896 }
2897 else
2b1eaa10 2898 if (fParticleID==AliPID::kProton)
32b846cd 2899 {
d79d61fb 2900 t = new TMatrixF(3,9);
2901 (*t)(0,0) = 0.20; (*t)(1,0) = -0.1; (*t)(2,0) = 0.1;
2902 (*t)(0,1) = 0.25; (*t)(1,1) = -0.2; (*t)(2,1) = 0.2;
2903 (*t)(0,2) = 0.80; (*t)(1,2) = -0.1; (*t)(2,2) = 0.2;
2904 (*t)(0,3) = 0.85; (*t)(1,3) =-0.05; (*t)(2,3) = 0.2;
2905 (*t)(0,4) = 0.90; (*t)(1,4) =-0.05; (*t)(2,4) = 0.25;
2906 (*t)(0,5) = 0.95; (*t)(1,5) =-0.05; (*t)(2,5) = 0.25;
2907 (*t)(0,6) = 1.00; (*t)(1,6) = -0.1; (*t)(2,6) = 0.25;
2908 (*t)(0,7) = 1.10; (*t)(1,7) =-0.05; (*t)(2,7) = 0.3;
2909 (*t)(0,8) = 1.20; (*t)(1,8) = 0; (*t)(2,8) = 0;
32b846cd 2910 }
d79d61fb 2911 delete fTPCpidCuts;
32b846cd 2912 fTPCpidCuts=t;
2913 }
7d9ab4fb 2914 t = NULL;
32b846cd 2915 if (!fTOFpidCuts)
2916 {
2b1eaa10 2917 if (fParticleID==AliPID::kPion)
32b846cd 2918 {
a0241c3a 2919 //TOF pions, 0.9 purity
2920 t = new TMatrixF(3,61);
d79d61fb 2921 (*t)(0,0) = 0.000; (*t)(1,0) = 0.000; (*t)(2,0) = 0.000;
2922 (*t)(0,1) = 0.050; (*t)(1,1) = 0.000; (*t)(2,1) = 0.000;
2923 (*t)(0,2) = 0.100; (*t)(1,2) = 0.000; (*t)(2,2) = 0.000;
2924 (*t)(0,3) = 0.150; (*t)(1,3) = 0.000; (*t)(2,3) = 0.000;
1a80f9f6 2925 (*t)(0,4) = 0.200; (*t)(1,4) = -0.030; (*t)(2,4) = 0.030;
d79d61fb 2926 (*t)(0,5) = 0.250; (*t)(1,5) = -0.036; (*t)(2,5) = 0.032;
2927 (*t)(0,6) = 0.300; (*t)(1,6) = -0.038; (*t)(2,6) = 0.032;
2928 (*t)(0,7) = 0.350; (*t)(1,7) = -0.034; (*t)(2,7) = 0.032;
2929 (*t)(0,8) = 0.400; (*t)(1,8) = -0.032; (*t)(2,8) = 0.020;
2930 (*t)(0,9) = 0.450; (*t)(1,9) = -0.030; (*t)(2,9) = 0.020;
2931 (*t)(0,10) = 0.500; (*t)(1,10) = -0.030; (*t)(2,10) = 0.020;
2932 (*t)(0,11) = 0.550; (*t)(1,11) = -0.030; (*t)(2,11) = 0.020;
2933 (*t)(0,12) = 0.600; (*t)(1,12) = -0.030; (*t)(2,12) = 0.020;
2934 (*t)(0,13) = 0.650; (*t)(1,13) = -0.030; (*t)(2,13) = 0.020;
2935 (*t)(0,14) = 0.700; (*t)(1,14) = -0.030; (*t)(2,14) = 0.020;
2936 (*t)(0,15) = 0.750; (*t)(1,15) = -0.030; (*t)(2,15) = 0.020;
2937 (*t)(0,16) = 0.800; (*t)(1,16) = -0.030; (*t)(2,16) = 0.020;
2938 (*t)(0,17) = 0.850; (*t)(1,17) = -0.030; (*t)(2,17) = 0.020;
2939 (*t)(0,18) = 0.900; (*t)(1,18) = -0.030; (*t)(2,18) = 0.020;
2940 (*t)(0,19) = 0.950; (*t)(1,19) = -0.028; (*t)(2,19) = 0.028;
2941 (*t)(0,20) = 1.000; (*t)(1,20) = -0.028; (*t)(2,20) = 0.028;
2942 (*t)(0,21) = 1.100; (*t)(1,21) = -0.028; (*t)(2,21) = 0.028;
2943 (*t)(0,22) = 1.200; (*t)(1,22) = -0.026; (*t)(2,22) = 0.028;
2944 (*t)(0,23) = 1.300; (*t)(1,23) = -0.024; (*t)(2,23) = 0.028;
2945 (*t)(0,24) = 1.400; (*t)(1,24) = -0.020; (*t)(2,24) = 0.028;
2946 (*t)(0,25) = 1.500; (*t)(1,25) = -0.018; (*t)(2,25) = 0.028;
2947 (*t)(0,26) = 1.600; (*t)(1,26) = -0.016; (*t)(2,26) = 0.028;
2948 (*t)(0,27) = 1.700; (*t)(1,27) = -0.014; (*t)(2,27) = 0.028;
2949 (*t)(0,28) = 1.800; (*t)(1,28) = -0.012; (*t)(2,28) = 0.026;
2950 (*t)(0,29) = 1.900; (*t)(1,29) = -0.010; (*t)(2,29) = 0.026;
2951 (*t)(0,30) = 2.000; (*t)(1,30) = -0.008; (*t)(2,30) = 0.026;
2952 (*t)(0,31) = 2.100; (*t)(1,31) = -0.008; (*t)(2,31) = 0.024;
2953 (*t)(0,32) = 2.200; (*t)(1,32) = -0.006; (*t)(2,32) = 0.024;
2954 (*t)(0,33) = 2.300; (*t)(1,33) = -0.004; (*t)(2,33) = 0.024;
2955 (*t)(0,34) = 2.400; (*t)(1,34) = -0.004; (*t)(2,34) = 0.024;
2956 (*t)(0,35) = 2.500; (*t)(1,35) = -0.002; (*t)(2,35) = 0.024;
2957 (*t)(0,36) = 2.600; (*t)(1,36) = -0.002; (*t)(2,36) = 0.024;
2958 (*t)(0,37) = 2.700; (*t)(1,37) = 0.000; (*t)(2,37) = 0.024;
2959 (*t)(0,38) = 2.800; (*t)(1,38) = 0.000; (*t)(2,38) = 0.026;
2960 (*t)(0,39) = 2.900; (*t)(1,39) = 0.000; (*t)(2,39) = 0.024;
2961 (*t)(0,40) = 3.000; (*t)(1,40) = 0.002; (*t)(2,40) = 0.026;
2962 (*t)(0,41) = 3.100; (*t)(1,41) = 0.002; (*t)(2,41) = 0.026;
2963 (*t)(0,42) = 3.200; (*t)(1,42) = 0.002; (*t)(2,42) = 0.026;
2964 (*t)(0,43) = 3.300; (*t)(1,43) = 0.002; (*t)(2,43) = 0.026;
2965 (*t)(0,44) = 3.400; (*t)(1,44) = 0.002; (*t)(2,44) = 0.026;
2966 (*t)(0,45) = 3.500; (*t)(1,45) = 0.002; (*t)(2,45) = 0.026;
2967 (*t)(0,46) = 3.600; (*t)(1,46) = 0.002; (*t)(2,46) = 0.026;
2968 (*t)(0,47) = 3.700; (*t)(1,47) = 0.002; (*t)(2,47) = 0.026;
2969 (*t)(0,48) = 3.800; (*t)(1,48) = 0.002; (*t)(2,48) = 0.026;
2970 (*t)(0,49) = 3.900; (*t)(1,49) = 0.004; (*t)(2,49) = 0.024;
2971 (*t)(0,50) = 4.000; (*t)(1,50) = 0.004; (*t)(2,50) = 0.026;
2972 (*t)(0,51) = 4.100; (*t)(1,51) = 0.004; (*t)(2,51) = 0.026;
2973 (*t)(0,52) = 4.200; (*t)(1,52) = 0.004; (*t)(2,52) = 0.024;
2974 (*t)(0,53) = 4.300; (*t)(1,53) = 0.006; (*t)(2,53) = 0.024;
2975 (*t)(0,54) = 4.400; (*t)(1,54) = 0.000; (*t)(2,54) = 0.000;
2976 (*t)(0,55) = 4.500; (*t)(1,55) = 0.000; (*t)(2,55) = 0.000;
2977 (*t)(0,56) = 4.600; (*t)(1,56) = 0.000; (*t)(2,56) = 0.000;
2978 (*t)(0,57) = 4.700; (*t)(1,57) = 0.000; (*t)(2,57) = 0.000;
2979 (*t)(0,58) = 4.800; (*t)(1,58) = 0.000; (*t)(2,58) = 0.000;
2980 (*t)(0,59) = 4.900; (*t)(1,59) = 0.000; (*t)(2,59) = 0.000;
2981 (*t)(0,60) = 5.900; (*t)(1,60) = 0.000; (*t)(2,60) = 0.000;
32b846cd 2982 }
2983 else
2b1eaa10 2984 if (fParticleID==AliPID::kProton)
32b846cd 2985 {
a0241c3a 2986 //TOF protons, 0.9 purity
2987 t = new TMatrixF(3,61);
d79d61fb 2988 (*t)(0,0) = 0.000; (*t)(1,0) = 0.000; (*t)(2,0) = 0.000;
2989 (*t)(0,1) = 0.050; (*t)(1,1) = 0.000; (*t)(2,1) = 0.000;
2990 (*t)(0,2) = 0.100; (*t)(1,2) = 0.000; (*t)(2,2) = 0.000;
2991 (*t)(0,3) = 0.150; (*t)(1,3) = 0.000; (*t)(2,3) = 0.000;
1a80f9f6 2992 (*t)(0,4) = 0.200; (*t)(1,4) = -0.07; (*t)(2,4) = 0.07;
2993 (*t)(0,5) = 0.200; (*t)(1,5) = -0.07; (*t)(2,5) = 0.07;
2994 (*t)(0,6) = 0.200; (*t)(1,6) = -0.07; (*t)(2,6) = 0.07;
2995 (*t)(0,7) = 0.200; (*t)(1,7) = -0.07; (*t)(2,7) = 0.07;
2996 (*t)(0,8) = 0.200; (*t)(1,8) = -0.07; (*t)(2,8) = 0.07;
2997 (*t)(0,9) = 0.200; (*t)(1,9) = -0.07; (*t)(2,9) = 0.07;
2998 (*t)(0,10) = 0.200; (*t)(1,10) = -0.07; (*t)(2,10) = 0.07;
2999 (*t)(0,11) = 0.200; (*t)(1,11) = -0.07; (*t)(2,11) = 0.07;
3000 (*t)(0,12) = 0.200; (*t)(1,12) = -0.07; (*t)(2,12) = 0.07;
3001 (*t)(0,13) = 0.200; (*t)(1,13) = -0.07; (*t)(2,13) = 0.07;
3002 (*t)(0,14) = 0.200; (*t)(1,14) = -0.07; (*t)(2,14) = 0.07;
3003 (*t)(0,15) = 0.200; (*t)(1,15) = -0.07; (*t)(2,15) = 0.07;
3004 (*t)(0,16) = 0.200; (*t)(1,16) = -0.07; (*t)(2,16) = 0.07;
d79d61fb 3005 (*t)(0,17) = 0.850; (*t)(1,17) = -0.070; (*t)(2,17) = 0.070;
3006 (*t)(0,18) = 0.900; (*t)(1,18) = -0.072; (*t)(2,18) = 0.072;
3007 (*t)(0,19) = 0.950; (*t)(1,19) = -0.072; (*t)(2,19) = 0.072;
3008 (*t)(0,20) = 1.000; (*t)(1,20) = -0.074; (*t)(2,20) = 0.074;
3009 (*t)(0,21) = 1.100; (*t)(1,21) = -0.032; (*t)(2,21) = 0.032;
3010 (*t)(0,22) = 1.200; (*t)(1,22) = -0.026; (*t)(2,22) = 0.026;
3011 (*t)(0,23) = 1.300; (*t)(1,23) = -0.026; (*t)(2,23) = 0.026;
3012 (*t)(0,24) = 1.400; (*t)(1,24) = -0.024; (*t)(2,24) = 0.024;
3013 (*t)(0,25) = 1.500; (*t)(1,25) = -0.024; (*t)(2,25) = 0.024;
3014 (*t)(0,26) = 1.600; (*t)(1,26) = -0.026; (*t)(2,26) = 0.026;
3015 (*t)(0,27) = 1.700; (*t)(1,27) = -0.026; (*t)(2,27) = 0.026;
3016 (*t)(0,28) = 1.800; (*t)(1,28) = -0.026; (*t)(2,28) = 0.026;
3017 (*t)(0,29) = 1.900; (*t)(1,29) = -0.026; (*t)(2,29) = 0.026;
3018 (*t)(0,30) = 2.000; (*t)(1,30) = -0.026; (*t)(2,30) = 0.026;
3019 (*t)(0,31) = 2.100; (*t)(1,31) = -0.026; (*t)(2,31) = 0.026;
3020 (*t)(0,32) = 2.200; (*t)(1,32) = -0.026; (*t)(2,32) = 0.024;
3021 (*t)(0,33) = 2.300; (*t)(1,33) = -0.028; (*t)(2,33) = 0.022;
3022 (*t)(0,34) = 2.400; (*t)(1,34) = -0.028; (*t)(2,34) = 0.020;
3023 (*t)(0,35) = 2.500; (*t)(1,35) = -0.028; (*t)(2,35) = 0.018;
3024 (*t)(0,36) = 2.600; (*t)(1,36) = -0.028; (*t)(2,36) = 0.016;
3025 (*t)(0,37) = 2.700; (*t)(1,37) = -0.028; (*t)(2,37) = 0.016;
3026 (*t)(0,38) = 2.800; (*t)(1,38) = -0.030; (*t)(2,38) = 0.014;
3027 (*t)(0,39) = 2.900; (*t)(1,39) = -0.030; (*t)(2,39) = 0.012;
3028 (*t)(0,40) = 3.000; (*t)(1,40) = -0.030; (*t)(2,40) = 0.012;
3029 (*t)(0,41) = 3.100; (*t)(1,41) = -0.030; (*t)(2,41) = 0.010;
3030 (*t)(0,42) = 3.200; (*t)(1,42) = -0.030; (*t)(2,42) = 0.010;
3031 (*t)(0,43) = 3.300; (*t)(1,43) = -0.030; (*t)(2,43) = 0.010;
3032 (*t)(0,44) = 3.400; (*t)(1,44) = -0.030; (*t)(2,44) = 0.008;
3033 (*t)(0,45) = 3.500; (*t)(1,45) = -0.030; (*t)(2,45) = 0.008;
3034 (*t)(0,46) = 3.600; (*t)(1,46) = -0.030; (*t)(2,46) = 0.008;
3035 (*t)(0,47) = 3.700; (*t)(1,47) = -0.030; (*t)(2,47) = 0.006;
3036 (*t)(0,48) = 3.800; (*t)(1,48) = -0.030; (*t)(2,48) = 0.006;
3037 (*t)(0,49) = 3.900; (*t)(1,49) = -0.030; (*t)(2,49) = 0.006;
3038 (*t)(0,50) = 4.000; (*t)(1,50) = -0.028; (*t)(2,50) = 0.004;
3039 (*t)(0,51) = 4.100; (*t)(1,51) = -0.030; (*t)(2,51) = 0.004;
3040 (*t)(0,52) = 4.200; (*t)(1,52) = -0.030; (*t)(2,52) = 0.004;
3041 (*t)(0,53) = 4.300; (*t)(1,53) = -0.028; (*t)(2,53) = 0.002;
3042 (*t)(0,54) = 4.400; (*t)(1,54) = -0.030; (*t)(2,54) = 0.002;
3043 (*t)(0,55) = 4.500; (*t)(1,55) = -0.028; (*t)(2,55) = 0.002;
3044 (*t)(0,56) = 4.600; (*t)(1,56) = -0.028; (*t)(2,56) = 0.002;
3045 (*t)(0,57) = 4.700; (*t)(1,57) = -0.028; (*t)(2,57) = 0.000;
3046 (*t)(0,58) = 4.800; (*t)(1,58) = -0.028; (*t)(2,58) = 0.002;
3047 (*t)(0,59) = 4.900; (*t)(1,59) = 0.000; (*t)(2,59) = 0.000;
3048 (*t)(0,60) = 5.900; (*t)(1,60) = 0.000; (*t)(2,60) = 0.000;
32b846cd 3049 }
3050 else
2b1eaa10 3051 if (fParticleID==AliPID::kKaon)
32b846cd 3052 {
a0241c3a 3053 //TOF kaons, 0.9 purity
3054 t = new TMatrixF(3,61);
d79d61fb 3055 (*t)(0,0) = 0.000; (*t)(1,0) = 0.000; (*t)(2,0) = 0.000;
3056 (*t)(0,1) = 0.050; (*t)(1,1) = 0.000; (*t)(2,1) = 0.000;
3057 (*t)(0,2) = 0.100; (*t)(1,2) = 0.000; (*t)(2,2) = 0.000;
3058 (*t)(0,3) = 0.150; (*t)(1,3) = 0.000; (*t)(2,3) = 0.000;
1a80f9f6 3059 (*t)(0,4) = 0.200; (*t)(1,4) = -0.05; (*t)(2,4) = 0.05;
3060 (*t)(0,5) = 0.200; (*t)(1,5) = -0.05; (*t)(2,5) = 0.05;
3061 (*t)(0,6) = 0.200; (*t)(1,6) = -0.05; (*t)(2,6) = 0.05;
3062 (*t)(0,7) = 0.200; (*t)(1,7) = -0.05; (*t)(2,7) = 0.05;
3063 (*t)(0,8) = 0.200; (*t)(1,8) = -0.05; (*t)(2,8) = 0.05;
3064 (*t)(0,9) = 0.200; (*t)(1,9) = -0.05; (*t)(2,9) = 0.05;
3065 (*t)(0,10) = 0.200; (*t)(1,10) = -0.05; (*t)(2,10) = 0.05;
d79d61fb 3066 (*t)(0,11) = 0.550; (*t)(1,11) = -0.026; (*t)(2,11) = 0.026;
3067 (*t)(0,12) = 0.600; (*t)(1,12) = -0.026; (*t)(2,12) = 0.026;
3068 (*t)(0,13) = 0.650; (*t)(1,13) = -0.026; (*t)(2,13) = 0.026;
3069 (*t)(0,14) = 0.700; (*t)(1,14) = -0.026; (*t)(2,14) = 0.026;
3070 (*t)(0,15) = 0.750; (*t)(1,15) = -0.026; (*t)(2,15) = 0.026;
3071 (*t)(0,16) = 0.800; (*t)(1,16) = -0.026; (*t)(2,16) = 0.026;
3072 (*t)(0,17) = 0.850; (*t)(1,17) = -0.024; (*t)(2,17) = 0.024;
3073 (*t)(0,18) = 0.900; (*t)(1,18) = -0.024; (*t)(2,18) = 0.024;
3074 (*t)(0,19) = 0.950; (*t)(1,19) = -0.024; (*t)(2,19) = 0.024;
3075 (*t)(0,20) = 1.000; (*t)(1,20) = -0.024; (*t)(2,20) = 0.024;
3076 (*t)(0,21) = 1.100; (*t)(1,21) = -0.024; (*t)(2,21) = 0.024;
3077 (*t)(0,22) = 1.200; (*t)(1,22) = -0.024; (*t)(2,22) = 0.022;
3078 (*t)(0,23) = 1.300; (*t)(1,23) = -0.024; (*t)(2,23) = 0.020;
3079 (*t)(0,24) = 1.400; (*t)(1,24) = -0.026; (*t)(2,24) = 0.016;
3080 (*t)(0,25) = 1.500; (*t)(1,25) = -0.028; (*t)(2,25) = 0.014;
3081 (*t)(0,26) = 1.600; (*t)(1,26) = -0.028; (*t)(2,26) = 0.012;
3082 (*t)(0,27) = 1.700; (*t)(1,27) = -0.028; (*t)(2,27) = 0.010;
3083 (*t)(0,28) = 1.800; (*t)(1,28) = -0.028; (*t)(2,28) = 0.010;
3084 (*t)(0,29) = 1.900; (*t)(1,29) = -0.028; (*t)(2,29) = 0.008;
3085 (*t)(0,30) = 2.000; (*t)(1,30) = -0.028; (*t)(2,30) = 0.006;
3086 (*t)(0,31) = 2.100; (*t)(1,31) = -0.026; (*t)(2,31) = 0.006;
3087 (*t)(0,32) = 2.200; (*t)(1,32) = -0.024; (*t)(2,32) = 0.004;
3088 (*t)(0,33) = 2.300; (*t)(1,33) = -0.020; (*t)(2,33) = 0.002;
3089 (*t)(0,34) = 2.400; (*t)(1,34) = -0.020; (*t)(2,34) = 0.002;
3090 (*t)(0,35) = 2.500; (*t)(1,35) = -0.018; (*t)(2,35) = 0.000;
3091 (*t)(0,36) = 2.600; (*t)(1,36) = -0.016; (*t)(2,36) = 0.000;
3092 (*t)(0,37) = 2.700; (*t)(1,37) = -0.014; (*t)(2,37) = -0.002;
3093 (*t)(0,38) = 2.800; (*t)(1,38) = -0.014; (*t)(2,38) = -0.004;
3094 (*t)(0,39) = 2.900; (*t)(1,39) = -0.012; (*t)(2,39) = -0.004;
3095 (*t)(0,40) = 3.000; (*t)(1,40) = -0.010; (*t)(2,40) = -0.006;
3096 (*t)(0,41) = 3.100; (*t)(1,41) = 0.000; (*t)(2,41) = 0.000;
3097 (*t)(0,42) = 3.200; (*t)(1,42) = 0.000; (*t)(2,42) = 0.000;
3098 (*t)(0,43) = 3.300; (*t)(1,43) = 0.000; (*t)(2,43) = 0.000;
3099 (*t)(0,44) = 3.400; (*t)(1,44) = 0.000; (*t)(2,44) = 0.000;
3100 (*t)(0,45) = 3.500; (*t)(1,45) = 0.000; (*t)(2,45) = 0.000;
3101 (*t)(0,46) = 3.600; (*t)(1,46) = 0.000; (*t)(2,46) = 0.000;
3102 (*t)(0,47) = 3.700; (*t)(1,47) = 0.000; (*t)(2,47) = 0.000;
3103 (*t)(0,48) = 3.800; (*t)(1,48) = 0.000; (*t)(2,48) = 0.000;
3104 (*t)(0,49) = 3.900; (*t)(1,49) = 0.000; (*t)(2,49) = 0.000;
3105 (*t)(0,50) = 4.000; (*t)(1,50) = 0.000; (*t)(2,50) = 0.000;
3106 (*t)(0,51) = 4.100; (*t)(1,51) = 0.000; (*t)(2,51) = 0.000;
3107 (*t)(0,52) = 4.200; (*t)(1,52) = 0.000; (*t)(2,52) = 0.000;
3108 (*t)(0,53) = 4.300; (*t)(1,53) = 0.000; (*t)(2,53) = 0.000;
3109 (*t)(0,54) = 4.400; (*t)(1,54) = 0.000; (*t)(2,54) = 0.000;
3110 (*t)(0,55) = 4.500; (*t)(1,55) = 0.000; (*t)(2,55) = 0.000;
3111 (*t)(0,56) = 4.600; (*t)(1,56) = 0.000; (*t)(2,56) = 0.000;
3112 (*t)(0,57) = 4.700; (*t)(1,57) = 0.000; (*t)(2,57) = 0.000;
3113 (*t)(0,58) = 4.800; (*t)(1,58) = 0.000; (*t)(2,58) = 0.000;
3114 (*t)(0,59) = 4.900; (*t)(1,59) = 0.000; (*t)(2,59) = 0.000;
3115 (*t)(0,60) = 5.900; (*t)(1,60) = 0.000; (*t)(2,60) = 0.000;
32b846cd 3116 }
d79d61fb 3117 delete fTOFpidCuts;
32b846cd 3118 fTOFpidCuts=t;
3119 }
3120}
3abf7ecc 3121
0d16f553 3122//-----------------------------------------------------------------------
3123//-----------------------------------------------------------------------
3124Bool_t AliFlowTrackCuts::PassesTPCbayesianCut(const AliAODTrack* track)
3125{
3126 fBayesianResponse->ComputeProb(track,track->GetAODEvent()); // fCurrCentr is needed for mismatch fraction
3127 Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mism$
3128
3129 Int_t kTPC = fBayesianResponse->GetCurrentMask(0); // is TPC on
3130
3131 if(! kTPC) return kFALSE;
3132
3133 fProbBayes = 0.0;
3134
3135switch (fParticleID)
3136 {
3137 case AliPID::kPion:
3138 fProbBayes = probabilities[2];
3139 break;
3140 case AliPID::kKaon:
3141 fProbBayes = probabilities[3];
3142 break;
3143 case AliPID::kProton:
3144 fProbBayes = probabilities[4];
3145 break;
3146 case AliPID::kElectron:
3147 fProbBayes = probabilities[0];
3148 break;
3149 case AliPID::kMuon:
3150 fProbBayes = probabilities[1];
3151 break;
3152 case AliPID::kDeuteron:
3153 fProbBayes = probabilities[5];
3154 break;
3155 case AliPID::kTriton:
3156 fProbBayes = probabilities[6];
3157 break;
3158 case AliPID::kHe3:
3159 fProbBayes = probabilities[7];
3160 break;
3161 default:
3162 return kFALSE;
3163 }
3164
3165 if(fProbBayes > fParticleProbability){
3166 if(!fCutCharge)
3167 return kTRUE;
3168 else if (fCutCharge && fCharge * track->Charge() > 0)
3169 return kTRUE;
3170 }
3171 return kFALSE;
3172}
cd62a2a8 3173//-----------------------------------------------------------------------
3174Bool_t AliFlowTrackCuts::PassesTPCbayesianCut(const AliESDtrack* track)
3175{
3176 //cut on TPC bayesian pid
cd62a2a8 3177 //if (! fAllowTOFmismatchFlag) {if (track->GetStatus() & AliESDtrack::kTOFmismatch) return kFALSE;}
3178
3179 //Bool_t statusMatchingHard = TPCTOFagree(track);
3180 //if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
3181 // return kFALSE;
8b649d8c 3182 fBayesianResponse->ComputeProb(track,fCurrCentr); // fCurrCentr is needed for mismatch fraction
9126fd0c 3183 Int_t kTPC = fBayesianResponse->GetCurrentMask(0); // is TPC on
0d16f553 3184 //Int_t kTOF = fBayesianResponse->GetCurrentMask(1); // is TOF on
cd62a2a8 3185
9126fd0c 3186 if(! kTPC) return kFALSE;
3187
0d16f553 3188 // Bool_t statusMatchingHard = 1;
3189 // Float_t mismProb = 0;
3190 // if(kTOF){
3191 // statusMatchingHard = TPCTOFagree(track);
3192 // mismProb = fBayesianResponse->GetTOFMismProb();
3193 // }
3194 // if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
3195 // return kFALSE;
9126fd0c 3196
9126fd0c 3197 Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mismatch parameterization
3198
8b649d8c 3199 fProbBayes = 0.0;
cd62a2a8 3200
cd62a2a8 3201 switch (fParticleID)
3202 {
3203 case AliPID::kPion:
8b649d8c 3204 fProbBayes = probabilities[2];
cd62a2a8 3205 break;
3206 case AliPID::kKaon:
8b649d8c 3207 fProbBayes = probabilities[3];
cd62a2a8 3208 break;
3209 case AliPID::kProton:
8b649d8c 3210 fProbBayes = probabilities[4];
cd62a2a8 3211 break;
3212 case AliPID::kElectron:
8b649d8c 3213 fProbBayes = probabilities[0];
3214 break;
3215 case AliPID::kMuon:
3216 fProbBayes = probabilities[1];
3217 break;
3218 case AliPID::kDeuteron:
3219 fProbBayes = probabilities[5];
3220 break;
3221 case AliPID::kTriton:
3222 fProbBayes = probabilities[6];
3223 break;
3224 case AliPID::kHe3:
3225 fProbBayes = probabilities[7];
cd62a2a8 3226 break;
3227 default:
3228 return kFALSE;
3229 }
3230
0d16f553 3231 if(fProbBayes > fParticleProbability)
9126fd0c 3232 {
3233 if(!fCutCharge)
3234 return kTRUE;
3235 else if (fCutCharge && fCharge * track->GetSign() > 0)
3236 return kTRUE;
3237 }
cd62a2a8 3238 return kFALSE;
3239}
0d16f553 3240//-----------------------------------------------------------------------
3241Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(const AliAODTrack* track)
3242{
3243//check is track passes bayesian combined TOF+TPC pid cut
3244 Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFout) &&
3245 (track->GetStatus() & AliESDtrack::kTIME) &&
3246 (track->GetTOFsignal() > 12000) &&
3247 (track->GetTOFsignal() < 100000);
cd62a2a8 3248
0d16f553 3249 if (! goodtrack)
3250 return kFALSE;
3251
3252 Bool_t statusMatchingHard = TPCTOFagree(track);
3253//ciao
3254 if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
3255 return kFALSE;
3256
3257 fBayesianResponse->ComputeProb(track,track->GetAODEvent()); // fCurrCentr is needed for mismatch fraction
3258 Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mism$
3259
3260 Float_t mismProb = fBayesianResponse->GetTOFMismProb(); // mismatch Bayesian probabilities
3261
3262 fProbBayes = 0.0;
3263
3264 switch (fParticleID)
3265 {
3266 case AliPID::kPion:
3267 fProbBayes = probabilities[2];
3268 break;
3269 case AliPID::kKaon:
3270 fProbBayes = probabilities[3];
3271 break;
3272 case AliPID::kProton:
3273 fProbBayes = probabilities[4];
3274 break;
3275 case AliPID::kElectron:
3276 fProbBayes = probabilities[0];
3277 break;
3278 case AliPID::kMuon:
3279 fProbBayes = probabilities[1];
3280 break;
3281 case AliPID::kDeuteron:
3282 fProbBayes = probabilities[5];
3283 break;
3284 case AliPID::kTriton:
3285 fProbBayes = probabilities[6];
3286 break;
3287 case AliPID::kHe3:
3288 fProbBayes = probabilities[7];
3289 break;
3290 default:
3291 return kFALSE;
3292 }
3293
3294 if(fProbBayes > fParticleProbability && mismProb < 0.5){
3295 if(!fCutCharge)
3296 return kTRUE;
3297 else if (fCutCharge && fCharge * track->Charge() > 0)
3298 return kTRUE;
3299 }
3300 return kFALSE;
3301
3302}
3abf7ecc 3303//-----------------------------------------------------------------------
3304// part added by F. Noferini (some methods)
cd62a2a8 3305Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(const AliESDtrack* track)
875ee12a 3306{
3c67c769 3307 //check is track passes bayesian combined TOF+TPC pid cut
8b649d8c 3308 Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFout) &&
3309 (track->GetStatus() & AliESDtrack::kTIME) &&
875ee12a 3310 (track->GetTOFsignal() > 12000) &&
3311 (track->GetTOFsignal() < 100000) &&
3312 (track->GetIntegratedLength() > 365);
3abf7ecc 3313
3abf7ecc 3314 if (! goodtrack)
3315 return kFALSE;
3316
6a043318 3317 if (! fAllowTOFmismatchFlag) {if (track->GetStatus() & AliESDtrack::kTOFmismatch) return kFALSE;}
3318
875ee12a 3319 Bool_t statusMatchingHard = TPCTOFagree(track);
3320 if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
3321 return kFALSE;
3322
9126fd0c 3323 fBayesianResponse->ComputeProb(track,fCurrCentr); // fCurrCentr is needed for mismatch fraction
3324 Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mismatch parameterization
3325
9126fd0c 3326 Float_t mismProb = fBayesianResponse->GetTOFMismProb(); // mismatch Bayesian probabilities
3abf7ecc 3327
8b649d8c 3328 fProbBayes = 0.0;
3329
2b1eaa10 3330 switch (fParticleID)
3abf7ecc 3331 {
3332 case AliPID::kPion:
8b649d8c 3333 fProbBayes = probabilities[2];
3abf7ecc 3334 break;
3335 case AliPID::kKaon:
8b649d8c 3336 fProbBayes = probabilities[3];
3abf7ecc 3337 break;
3338 case AliPID::kProton:
8b649d8c 3339 fProbBayes = probabilities[4];
3abf7ecc 3340 break;
3341 case AliPID::kElectron:
8b649d8c 3342 fProbBayes = probabilities[0];
3abf7ecc 3343 break;
8b649d8c 3344 case AliPID::kMuon:
3345 fProbBayes = probabilities[1];
3346 break;
3347 case AliPID::kDeuteron:
3348 fProbBayes = probabilities[5];
3349 break;
3350 case AliPID::kTriton:
3351 fProbBayes = probabilities[6];
3352 break;
3353 case AliPID::kHe3:
3354 fProbBayes = probabilities[7];
3355 break;
3356 default:
3abf7ecc 3357 return kFALSE;
3358 }
3359
3360 // printf("pt = %f -- all prob = [%4.2f,%4.2f,%4.2f,%4.2f,%4.2f] -- prob = %f\n",track->Pt(),fProbBayes[0],fProbBayes[1],fProbBayes[2],fProbBayes[3],fProbBayes[4],prob);
8b649d8c 3361 if(fProbBayes > fParticleProbability && mismProb < 0.5){
3abf7ecc 3362 if(!fCutCharge)
3363 return kTRUE;
3364 else if (fCutCharge && fCharge * track->GetSign() > 0)
3365 return kTRUE;
3366 }
3367 return kFALSE;
3368}
875ee12a 3369
8eca5d19 3370
3371//-----------------------------------------------------------------------
3372 // part added by Natasha
3373Bool_t AliFlowTrackCuts::PassesNucleiSelection(const AliESDtrack* track)
3374{
3c67c769 3375 //pid selection for heavy nuclei
8eca5d19 3376 Bool_t select=kFALSE;
3377
3378 //if (!track) continue;
6887199d 3379
8eca5d19 3380 if (!track->GetInnerParam())
3381 return kFALSE; //break;
6887199d 3382
8eca5d19 3383 const AliExternalTrackParam* tpcTrack = track->GetInnerParam();
6887199d 3384
8eca5d19 3385 Double_t ptotTPC = tpcTrack->GetP();
3386 Double_t sigTPC = track->GetTPCsignal();
3387 Double_t dEdxBBA = 0.;
3388 Double_t dSigma = 0.;
3389
3390 switch (fParticleID)
6887199d 3391 {
8eca5d19 3392 case AliPID::kDeuteron:
3393 //pid=10;
3394 dEdxBBA = AliExternalTrackParam::BetheBlochAleph(ptotTPC/1.8756,
6887199d 3395 4.60e+00,
3396 8.9684e+00,
3397 1.640e-05,
3398 2.35e+00,
3399 2.35e+00);
8eca5d19 3400 dSigma = (sigTPC - dEdxBBA)/dEdxBBA;
6887199d 3401
8eca5d19 3402 if( ptotTPC<=1.1 && (dSigma < (0.5 - (0.1818*ptotTPC)) ) && (dSigma > ( (0.218*ptotTPC - 0.4) ) ) )
6887199d 3403 {select=kTRUE;}
8eca5d19 3404 break;
6887199d 3405
8eca5d19 3406 case AliPID::kTriton:
3407 //pid=11;
3408 select=kFALSE;
3409 break;
6887199d 3410
8eca5d19 3411 case AliPID::kHe3:
3412 //pid=12;
6887199d 3413 // ----- Pass 2 -------
3414 dEdxBBA = 4.0 * AliExternalTrackParam::BetheBlochAleph( (2.*ptotTPC)/2.8084,
3415 1.74962,
3416 27.4992,
3417 4.00313e-15,
3418 2.42485,
3419 8.31768);
3420 dSigma = (sigTPC - dEdxBBA)/dEdxBBA;
3421 if(ptotTPC<=5.0 && (dSigma >= (-0.03968*ptotTPC - 0.1)) && (dSigma <= (0.31 - 0.0217*ptotTPC)))
3422 {select=kTRUE;}
8eca5d19 3423 break;
6887199d 3424
8eca5d19 3425 case AliPID::kAlpha:
3426 //pid=13;
3427 select=kFALSE;
3428 break;
6887199d 3429
8eca5d19 3430 default:
3431 return kFALSE;
6887199d 3432 }
3433
8eca5d19 3434 return select;
3435}
3436// end part added by Natasha
3abf7ecc 3437//-----------------------------------------------------------------------
875ee12a 3438void AliFlowTrackCuts::SetPriors(Float_t centrCur){
cc0afcfc 3439 //set priors for the bayesian pid selection
9126fd0c 3440 fCurrCentr = centrCur;
3441
875ee12a 3442 fBinLimitPID[0] = 0.300000;
3443 fBinLimitPID[1] = 0.400000;
3444 fBinLimitPID[2] = 0.500000;
3445 fBinLimitPID[3] = 0.600000;
3446 fBinLimitPID[4] = 0.700000;
3447 fBinLimitPID[5] = 0.800000;
3448 fBinLimitPID[6] = 0.900000;
3449 fBinLimitPID[7] = 1.000000;
3450 fBinLimitPID[8] = 1.200000;
3451 fBinLimitPID[9] = 1.400000;
3452 fBinLimitPID[10] = 1.600000;
3453 fBinLimitPID[11] = 1.800000;
3454 fBinLimitPID[12] = 2.000000;
3455 fBinLimitPID[13] = 2.200000;
3456 fBinLimitPID[14] = 2.400000;
3457 fBinLimitPID[15] = 2.600000;
3458 fBinLimitPID[16] = 2.800000;
3459 fBinLimitPID[17] = 3.000000;
3460
3461 // 0-10%
3462 if(centrCur < 10){
3463 fC[0][0] = 0.005;
3464 fC[0][1] = 0.005;
3465 fC[0][2] = 1.0000;
3466 fC[0][3] = 0.0010;
3467 fC[0][4] = 0.0010;
3468
3469 fC[1][0] = 0.005;
3470 fC[1][1] = 0.005;
3471 fC[1][2] = 1.0000;
3472 fC[1][3] = 0.0168;
3473 fC[1][4] = 0.0122;
3474
3475 fC[2][0] = 0.005;
3476 fC[2][1] = 0.005;
3477 fC[2][2] = 1.0000;
3478 fC[2][3] = 0.0272;
3479 fC[2][4] = 0.0070;
3480
3481 fC[3][0] = 0.005;
3482 fC[3][1] = 0.005;
3483 fC[3][2] = 1.0000;
3484 fC[3][3] = 0.0562;
3485 fC[3][4] = 0.0258;
3486
3487 fC[4][0] = 0.005;
3488 fC[4][1] = 0.005;
3489 fC[4][2] = 1.0000;
3490 fC[4][3] = 0.0861;
3491 fC[4][4] = 0.0496;
3492
3493 fC[5][0] = 0.005;
3494 fC[5][1] = 0.005;
3495 fC[5][2] = 1.0000;
3496 fC[5][3] = 0.1168;
3497 fC[5][4] = 0.0740;
3498
3499 fC[6][0] = 0.005;
3500 fC[6][1] = 0.005;
3501 fC[6][2] = 1.0000;
3502 fC[6][3] = 0.1476;
3503 fC[6][4] = 0.0998;
3504
3505 fC[7][0] = 0.005;
3506 fC[7][1] = 0.005;
3507 fC[7][2] = 1.0000;
3508 fC[7][3] = 0.1810;
3509 fC[7][4] = 0.1296;
3510
3511 fC[8][0] = 0.005;
3512 fC[8][1] = 0.005;
3513 fC[8][2] = 1.0000;
3514 fC[8][3] = 0.2240;
3515 fC[8][4] = 0.1827;
3516
3517 fC[9][0] = 0.005;
3518 fC[9][1] = 0.005;
3519 fC[9][2] = 1.0000;
3520 fC[9][3] = 0.2812;
3521 fC[9][4] = 0.2699;
3522
3523 fC[10][0] = 0.005;
3524 fC[10][1] = 0.005;
3525 fC[10][2] = 1.0000;
3526 fC[10][3] = 0.3328;
3527 fC[10][4] = 0.3714;
3528
3529 fC[11][0] = 0.005;
3530 fC[11][1] = 0.005;
3531 fC[11][2] = 1.0000;
3532 fC[11][3] = 0.3780;
3533 fC[11][4] = 0.4810;
3534
3535 fC[12][0] = 0.005;
3536 fC[12][1] = 0.005;
3537 fC[12][2] = 1.0000;
3538 fC[12][3] = 0.4125;
3539 fC[12][4] = 0.5771;
3540
3541 fC[13][0] = 0.005;
3542 fC[13][1] = 0.005;
3543 fC[13][2] = 1.0000;
3544 fC[13][3] = 0.4486;
3545 fC[13][4] = 0.6799;
3546
3547 fC[14][0] = 0.005;
3548 fC[14][1] = 0.005;
3549 fC[14][2] = 1.0000;
3550 fC[14][3] = 0.4840;
3551 fC[14][4] = 0.7668;
3552
3553 fC[15][0] = 0.005;
3554 fC[15][1] = 0.005;
3555 fC[15][2] = 1.0000;
3556 fC[15][3] = 0.4971;
3557 fC[15][4] = 0.8288;
3558
3559 fC[16][0] = 0.005;
3560 fC[16][1] = 0.005;
3561 fC[16][2] = 1.0000;
3562 fC[16][3] = 0.4956;
3563 fC[16][4] = 0.8653;
3564
3565 fC[17][0] = 0.005;
3566 fC[17][1] = 0.005;
3567 fC[17][2] = 1.0000;
3568 fC[17][3] = 0.5173;
3569 fC[17][4] = 0.9059;
3570 }
3571 // 10-20%
3572 else if(centrCur < 20){
3573 fC[0][0] = 0.005;
3574 fC[0][1] = 0.005;
3575 fC[0][2] = 1.0000;
3576 fC[0][3] = 0.0010;
3577 fC[0][4] = 0.0010;
3578
3579 fC[1][0] = 0.005;
3580 fC[1][1] = 0.005;
3581 fC[1][2] = 1.0000;
3582 fC[1][3] = 0.0132;
3583 fC[1][4] = 0.0088;
3584
3585 fC[2][0] = 0.005;
3586 fC[2][1] = 0.005;
3587 fC[2][2] = 1.0000;
3588 fC[2][3] = 0.0283;
3589 fC[2][4] = 0.0068;
3590
3591 fC[3][0] = 0.005;
3592 fC[3][1] = 0.005;
3593 fC[3][2] = 1.0000;
3594 fC[3][3] = 0.0577;
3595 fC[3][4] = 0.0279;
3596
3597 fC[4][0] = 0.005;
3598 fC[4][1] = 0.005;
3599 fC[4][2] = 1.0000;
3600 fC[4][3] = 0.0884;
3601 fC[4][4] = 0.0534;
3602
3603 fC[5][0] = 0.005;
3604 fC[5][1] = 0.005;
3605 fC[5][2] = 1.0000;
3606 fC[5][3] = 0.1179;
3607 fC[5][4] = 0.0794;
3608
3609 fC[6][0] = 0.005;
3610 fC[6][1] = 0.005;
3611 fC[6][2] = 1.0000;
3612 fC[6][3] = 0.1480;
3613 fC[6][4] = 0.1058;
3614
3615 fC[7][0] = 0.005;
3616 fC[7][1] = 0.005;
3617 fC[7][2] = 1.0000;
3618 fC[7][3] = 0.1807;
3619 fC[7][4] = 0.1366;
3620
3621 fC[8][0] = 0.005;
3622 fC[8][1] = 0.005;
3623 fC[8][2] = 1.0000;
3624 fC[8][3] = 0.2219;
3625 fC[8][4] = 0.1891;
3626
3627 fC[9][0] = 0.005;
3628 fC[9][1] = 0.005;
3629 fC[9][2] = 1.0000;
3630 fC[9][3] = 0.2804;
3631 fC[9][4] = 0.2730;
3632
3633 fC[10][0] = 0.005;
3634 fC[10][1] = 0.005;
3635 fC[10][2] = 1.0000;
3636 fC[10][3] = 0.3283;
3637 fC[10][4] = 0.3660;
3638
3639 fC[11][0] = 0.005;
3640 fC[11][1] = 0.005;
3641 fC[11][2] = 1.0000;
3642 fC[11][3] = 0.3710;
3643 fC[11][4] = 0.4647;
3644
3645 fC[12][0] = 0.005;
3646 fC[12][1] = 0.005;
3647 fC[12][2] = 1.0000;
3648 fC[12][3] = 0.4093;
3649 fC[12][4] = 0.5566;
3650
3651 fC[13][0] = 0.005;
3652 fC[13][1] = 0.005;
3653 fC[13][2] = 1.0000;
3654 fC[13][3] = 0.4302;
3655 fC[13][4] = 0.6410;
3656
3657 fC[14][0] = 0.005;
3658 fC[14][1] = 0.005;
3659 fC[14][2] = 1.0000;
3660 fC[14][3] = 0.4649;
3661 fC[14][4] = 0.7055;
3662
3663 fC[15][0] = 0.005;
3664 fC[15][1] = 0.005;
3665 fC[15][2] = 1.0000;
3666 fC[15][3] = 0.4523;
3667 fC[15][4] = 0.7440;
3668
3669 fC[16][0] = 0.005;
3670 fC[16][1] = 0.005;
3671 fC[16][2] = 1.0000;
3672 fC[16][3] = 0.4591;
3673 fC[16][4] = 0.7799;
3674
3675 fC[17][0] = 0.005;
3676 fC[17][1] = 0.005;
3677 fC[17][2] = 1.0000;
3678 fC[17][3] = 0.4804;
3679 fC[17][4] = 0.8218;
3680 }
3681 // 20-30%
3682 else if(centrCur < 30){
3683 fC[0][0] = 0.005;
3684 fC[0][1] = 0.005;
3685 fC[0][2] = 1.0000;
3686 fC[0][3] = 0.0010;
3687 fC[0][4] = 0.0010;
3688
3689 fC[1][0] = 0.005;
3690 fC[1][1] = 0.005;
3691 fC[1][2] = 1.0000;
3692 fC[1][3] = 0.0102;
3693 fC[1][4] = 0.0064;
3694
3695 fC[2][0] = 0.005;
3696 fC[2][1] = 0.005;
3697 fC[2][2] = 1.0000;
3698 fC[2][3] = 0.0292;
3699 fC[2][4] = 0.0066;
3700
3701 fC[3][0] = 0.005;
3702 fC[3][1] = 0.005;
3703 fC[3][2] = 1.0000;
3704 fC[3][3] = 0.0597;
3705 fC[3][4] = 0.0296;
3706
3707 fC[4][0] = 0.005;
3708 fC[4][1] = 0.005;
3709 fC[4][2] = 1.0000;
3710 fC[4][3] = 0.0900;
3711 fC[4][4] = 0.0589;
3712
3713 fC[5][0] = 0.005;
3714 fC[5][1] = 0.005;
3715 fC[5][2] = 1.0000;
3716 fC[5][3] = 0.1199;
3717 fC[5][4] = 0.0859;
3718
3719 fC[6][0] = 0.005;
3720 fC[6][1] = 0.005;
3721 fC[6][2] = 1.0000;
3722 fC[6][3] = 0.1505;
3723 fC[6][4] = 0.1141;
3724
3725 fC[7][0] = 0.005;
3726 fC[7][1] = 0.005;
3727 fC[7][2] = 1.0000;
3728 fC[7][3] = 0.1805;
3729 fC[7][4] = 0.1454;
3730
3731 fC[8][0] = 0.005;
3732 fC[8][1] = 0.005;
3733 fC[8][2] = 1.0000;
3734 fC[8][3] = 0.2221;
3735 fC[8][4] = 0.2004;
3736
3737 fC[9][0] = 0.005;
3738 fC[9][1] = 0.005;
3739 fC[9][2] = 1.0000;
3740 fC[9][3] = 0.2796;
3741 fC[9][4] = 0.2838;
3742
3743 fC[10][0] = 0.005;
3744 fC[10][1] = 0.005;
3745 fC[10][2] = 1.0000;
3746 fC[10][3] = 0.3271;
3747 fC[10][4] = 0.3682;
3748
3749 fC[11][0] = 0.005;
3750 fC[11][1] = 0.005;
3751 fC[11][2] = 1.0000;
3752 fC[11][3] = 0.3648;
3753 fC[11][4] = 0.4509;
3754
3755 fC[12][0] = 0.005;
3756 fC[12][1] = 0.005;
3757 fC[12][2] = 1.0000;
3758 fC[12][3] = 0.3988;
3759 fC[12][4] = 0.5339;
3760
3761 fC[13][0] = 0.005;
3762 fC[13][1] = 0.005;
3763 fC[13][2] = 1.0000;
3764 fC[13][3] = 0.4315;
3765 fC[13][4] = 0.5995;
3766
3767 fC[14][0] = 0.005;
3768 fC[14][1] = 0.005;
3769 fC[14][2] = 1.0000;
3770 fC[14][3] = 0.4548;
3771 fC[14][4] = 0.6612;
3772
3773 fC[15][0] = 0.005;
3774 fC[15][1] = 0.005;
3775 fC[15][2] = 1.0000;
3776 fC[15][3] = 0.4744;
3777 fC[15][4] = 0.7060;
3778
3779 fC[16][0] = 0.005;
3780 fC[16][1] = 0.005;
3781 fC[16][2] = 1.0000;
3782 fC[16][3] = 0.4899;
3783 fC[16][4] = 0.7388;
3784
3785 fC[17][0] = 0.005;
3786 fC[17][1] = 0.005;
3787 fC[17][2] = 1.0000;
3788 fC[17][3] = 0.4411;
3789 fC[17][4] = 0.7293;
3790 }
3791 // 30-40%
3792 else if(centrCur < 40){
3793 fC[0][0] = 0.005;
3794 fC[0][1] = 0.005;
3795 fC[0][2] = 1.0000;
3796 fC[0][3] = 0.0010;
3797 fC[0][4] = 0.0010;
3798
3799 fC[1][0] = 0.005;
3800 fC[1][1] = 0.005;
3801 fC[1][2] = 1.0000;
3802 fC[1][3] = 0.0102;
3803 fC[1][4] = 0.0048;
3804
3805 fC[2][0] = 0.005;
3806 fC[2][1] = 0.005;
3807 fC[2][2] = 1.0000;
3808 fC[2][3] = 0.0306;
3809 fC[2][4] = 0.0079;
3810
3811 fC[3][0] = 0.005;
3812 fC[3][1] = 0.005;
3813 fC[3][2] = 1.0000;
3814 fC[3][3] = 0.0617;
3815 fC[3][4] = 0.0338;
3816
3817 fC[4][0] = 0.005;
3818 fC[4][1] = 0.005;
3819 fC[4][2] = 1.0000;
3820 fC[4][3] = 0.0920;
3821 fC[4][4] = 0.0652;
3822
3823 fC[5][0] = 0.005;
3824 fC[5][1] = 0.005;
3825 fC[5][2] = 1.0000;
3826 fC[5][3] = 0.1211;
3827 fC[5][4] = 0.0955;
3828
3829 fC[6][0] = 0.005;
3830 fC[6][1] = 0.005;
3831 fC[6][2] = 1.0000;
3832 fC[6][3] = 0.1496;
3833 fC[6][4] = 0.1242;
3834
3835 fC[7][0] = 0.005;
3836 fC[7][1] = 0.005;
3837 fC[7][2] = 1.0000;
3838 fC[7][3] = 0.1807;
3839 fC[7][4] = 0.1576;
3840
3841 fC[8][0] = 0.005;
3842 fC[8][1] = 0.005;
3843 fC[8][2] = 1.0000;
3844 fC[8][3] = 0.2195;
3845 fC[8][4] = 0.2097;
3846
3847 fC[9][0] = 0.005;
3848 fC[9][1] = 0.005;
3849 fC[9][2] = 1.0000;
3850 fC[9][3] = 0.2732;
3851 fC[9][4] = 0.2884;
3852
3853 fC[10][0] = 0.005;
3854 fC[10][1] = 0.005;
3855 fC[10][2] = 1.0000;
3856 fC[10][3] = 0.3204;
3857 fC[10][4] = 0.3679;
3858
3859 fC[11][0] = 0.005;
3860 fC[11][1] = 0.005;
3861 fC[11][2] = 1.0000;
3862 fC[11][3] = 0.3564;
3863 fC[11][4] = 0.4449;
3864
3865 fC[12][0] = 0.005;
3866 fC[12][1] = 0.005;
3867 fC[12][2] = 1.0000;
3868 fC[12][3] = 0.3791;
3869 fC[12][4] = 0.5052;
3870
3871 fC[13][0] = 0.005;
3872 fC[13][1] = 0.005;
3873 fC[13][2] = 1.0000;
3874 fC[13][3] = 0.4062;
3875 fC[13][4] = 0.5647;
3876
3877 fC[14][0] = 0.005;
3878 fC[14][1] = 0.005;
3879 fC[14][2] = 1.0000;
3880 fC[14][3] = 0.4234;
3881 fC[14][4] = 0.6203;
3882
3883 fC[15][0] = 0.005;
3884 fC[15][1] = 0.005;
3885 fC[15][2] = 1.0000;
3886 fC[15][3] = 0.4441;
3887 fC[15][4] = 0.6381;
3888
3889 fC[16][0] = 0.005;
3890 fC[16][1] = 0.005;
3891 fC[16][2] = 1.0000;
3892 fC[16][3] = 0.4629;
3893 fC[16][4] = 0.6496;
3894
3895 fC[17][0] = 0.005;
3896 fC[17][1] = 0.005;
3897 fC[17][2] = 1.0000;
3898 fC[17][3] = 0.4293;
3899 fC[17][4] = 0.6491;
3900 }
3901 // 40-50%
3902 else if(centrCur < 50){
3903 fC[0][0] = 0.005;
3904 fC[0][1] = 0.005;
3905 fC[0][2] = 1.0000;
3906 fC[0][3] = 0.0010;
3907 fC[0][4] = 0.0010;
3908
3909 fC[1][0] = 0.005;
3910 fC[1][1] = 0.005;
3911 fC[1][2] = 1.0000;
3912 fC[1][3] = 0.0093;
3913 fC[1][4] = 0.0057;
3914
3915 fC[2][0] = 0.005;
3916 fC[2][1] = 0.005;
3917 fC[2][2] = 1.0000;
3918 fC[2][3] = 0.0319;
3919 fC[2][4] = 0.0075;
3920
3921 fC[3][0] = 0.005;
3922 fC[3][1] = 0.005;
3923 fC[3][2] = 1.0000;
3924 fC[3][3] = 0.0639;
3925 fC[3][4] = 0.0371;
3926
3927 fC[4][0] = 0.005;
3928 fC[4][1] = 0.005;
3929 fC[4][2] = 1.0000;
3930 fC[4][3] = 0.0939;
3931 fC[4][4] = 0.0725;
3932
3933 fC[5][0] = 0.005;
3934 fC[5][1] = 0.005;
3935 fC[5][2] = 1.0000;
3936 fC[5][3] = 0.1224;
3937 fC[5][4] = 0.1045;
3938
3939 fC[6][0] = 0.005;
3940 fC[6][1] = 0.005;
3941 fC[6][2] = 1.0000;
3942 fC[6][3] = 0.1520;
3943 fC[6][4] = 0.1387;
3944
3945 fC[7][0] = 0.005;
3946 fC[7][1] = 0.005;
3947 fC[7][2] = 1.0000;
3948 fC[7][3] = 0.1783;
3949 fC[7][4] = 0.1711;
3950
3951 fC[8][0] = 0.005;
3952 fC[8][1] = 0.005;
3953 fC[8][2] = 1.0000;
3954 fC[8][3] = 0.2202;
3955 fC[8][4] = 0.2269;
3956
3957 fC[9][0] = 0.005;
3958 fC[9][1] = 0.005;
3959 fC[9][2] = 1.0000;
3960 fC[9][3] = 0.2672;
3961 fC[9][4] = 0.2955;
3962
3963 fC[10][0] = 0.005;
3964 fC[10][1] = 0.005;
3965 fC[10][2] = 1.0000;
3966 fC[10][3] = 0.3191;
3967 fC[10][4] = 0.3676;
3968
3969 fC[11][0] = 0.005;
3970 fC[11][1] = 0.005;
3971 fC[11][2] = 1.0000;
3972 fC[11][3] = 0.3434;
3973 fC[11][4] = 0.4321;
3974
3975 fC[12][0] = 0.005;
3976 fC[12][1] = 0.005;
3977 fC[12][2] = 1.0000;
3978 fC[12][3] = 0.3692;
3979 fC[12][4] = 0.4879;
3980
3981 fC[13][0] = 0.005;
3982 fC[13][1] = 0.005;
3983 fC[13][2] = 1.0000;
3984 fC[13][3] = 0.3993;
3985 fC[13][4] = 0.5377;
3986
3987 fC[14][0] = 0.005;
3988 fC[14][1] = 0.005;
3989 fC[14][2] = 1.0000;
3990 fC[14][3] = 0.3818;
3991 fC[14][4] = 0.5547;
3992
3993 fC[15][0] = 0.005;
3994 fC[15][1] = 0.005;
3995 fC[15][2] = 1.0000;
3996 fC[15][3] = 0.4003;
3997 fC[15][4] = 0.5484;
3998
3999 fC[16][0] = 0.005;
4000 fC[16][1] = 0.005;
4001 fC[16][2] = 1.0000;
4002 fC[16][3] = 0.4281;
4003 fC[16][4] = 0.5383;
4004
4005 fC[17][0] = 0.005;
4006 fC[17][1] = 0.005;
4007 fC[17][2] = 1.0000;
4008 fC[17][3] = 0.3960;
4009 fC[17][4] = 0.5374;
4010 }
4011 // 50-60%
4012 else if(centrCur < 60){
4013 fC[0][0] = 0.005;
4014 fC[0][1] = 0.005;
4015 fC[0][2] = 1.0000;
4016 fC[0][3] = 0.0010;
4017 fC[0][4] = 0.0010;
4018
4019 fC[1][0] = 0.005;
4020 fC[1][1] = 0.005;
4021 fC[1][2] = 1.0000;
4022 fC[1][3] = 0.0076;
4023 fC[1][4] = 0.0032;
4024
4025 fC[2][0] = 0.005;
4026 fC[2][1] = 0.005;
4027 fC[2][2] = 1.0000;
4028 fC[2][3] = 0.0329;
4029 fC[2][4] = 0.0085;
4030
4031 fC[3][0] = 0.005;
4032 fC[3][1] = 0.005;
4033 fC[3][2] = 1.0000;
4034 fC[3][3] = 0.0653;
4035 fC[3][4] = 0.0423;
4036
4037 fC[4][0] = 0.005;
4038 fC[4][1] = 0.005;
4039 fC[4][2] = 1.0000;
4040 fC[4][3] = 0.0923;
4041 fC[4][4] = 0.0813;
4042
4043 fC[5][0] = 0.005;
4044 fC[5][1] = 0.005;
4045 fC[5][2] = 1.0000;
4046 fC[5][3] = 0.1219;
4047 fC[5][4] = 0.1161;
4048
4049 fC[6][0] = 0.005;
4050 fC[6][1] = 0.005;
4051 fC[6][2] = 1.0000;
4052 fC[6][3] = 0.1519;
4053 fC[6][4] = 0.1520;
4054
4055 fC[7][0] = 0.005;
4056 fC[7][1] = 0.005;
4057 fC[7][2] = 1.0000;
4058 fC[7][3] = 0.1763;
4059 fC[7][4] = 0.1858;
4060
4061 fC[8][0] = 0.005;
4062 fC[8][1] = 0.005;
4063 fC[8][2] = 1.0000;
4064 fC[8][3] = 0.2178;
4065 fC[8][4] = 0.2385;
4066
4067 fC[9][0] = 0.005;
4068 fC[9][1] = 0.005;
4069 fC[9][2] = 1.0000;
4070 fC[9][3] = 0.2618;
4071 fC[9][4] = 0.3070;
4072
4073 fC[10][0] = 0.005;
4074 fC[10][1] = 0.005;
4075 fC[10][2] = 1.0000;
4076 fC[10][3] = 0.3067;
4077 fC[10][4] = 0.3625;
4078
4079 fC[11][0] = 0.005;
4080 fC[11][1] = 0.005;
4081 fC[11][2] = 1.0000;
4082 fC[11][3] = 0.3336;
4083 fC[11][4] = 0.4188;
4084
4085 fC[12][0] = 0.005;
4086 fC[12][1] = 0.005;
4087 fC[12][2] = 1.0000;
4088 fC[12][3] = 0.3706;
4089 fC[12][4] = 0.4511;
4090
4091 fC[13][0] = 0.005;
4092 fC[13][1] = 0.005;
4093 fC[13][2] = 1.0000;
4094 fC[13][3] = 0.3765;
4095 fC[13][4] = 0.4729;
4096
4097 fC[14][0] = 0.005;
4098 fC[14][1] = 0.005;
4099 fC[14][2] = 1.0000;
4100 fC[14][3] = 0.3942;
4101 fC[14][4] = 0.4855;
4102
4103 fC[15][0] = 0.005;
4104 fC[15][1] = 0.005;
4105 fC[15][2] = 1.0000;
4106 fC[15][3] = 0.4051;
4107 fC[15][4] = 0.4762;
4108
4109 fC[16][0] = 0.005;
4110 fC[16][1] = 0.005;
4111 fC[16][2] = 1.0000;
4112 fC[16][3] = 0.3843;
4113 fC[16][4] = 0.4763;
4114
4115 fC[17][0] = 0.005;
4116 fC[17][1] = 0.005;
4117 fC[17][2] = 1.0000;
4118 fC[17][3] = 0.4237;
4119 fC[17][4] = 0.4773;
4120 }
4121 // 60-70%
4122 else if(centrCur < 70){
4123 fC[0][0] = 0.005;
4124 fC[0][1] = 0.005;
4125 fC[0][2] = 1.0000;
4126 fC[0][3] = 0.0010;
4127 fC[0][4] = 0.0010;
4128
4129 fC[1][0] = 0.005;
4130 fC[1][1] = 0.005;
4131 fC[1][2] = 1.0000;
4132 fC[1][3] = 0.0071;
4133 fC[1][4] = 0.0012;
4134
4135 fC[2][0] = 0.005;
4136 fC[2][1] = 0.005;
4137 fC[2][2] = 1.0000;
4138 fC[2][3] = 0.0336;
4139 fC[2][4] = 0.0097;
4140
4141 fC[3][0] = 0.005;
4142 fC[3][1] = 0.005;
4143 fC[3][2] = 1.0000;
4144 fC[3][3] = 0.0662;
4145 fC[3][4] = 0.0460;
4146
4147 fC[4][0] = 0.005;
4148 fC[4][1] = 0.005;
4149 fC[4][2] = 1.0000;
4150 fC[4][3] = 0.0954;
4151 fC[4][4] = 0.0902;
4152
4153 fC[5][0] = 0.005;
4154 fC[5][1] = 0.005;
4155 fC[5][2] = 1.0000;
4156 fC[5][3] = 0.1181;
4157 fC[5][4] = 0.1306;
4158
4159 fC[6][0] = 0.005;
4160 fC[6][1] = 0.005;
4161 fC[6][2] = 1.0000;
4162 fC[6][3] = 0.1481;
4163 fC[6][4] = 0.1662;
4164
4165 fC[7][0] = 0.005;
4166 fC[7][1] = 0.005;
4167 fC[7][2] = 1.0000;
4168 fC[7][3] = 0.1765;
4169 fC[7][4] = 0.1963;
4170
4171 fC[8][0] = 0.005;
4172 fC[8][1] = 0.005;
4173 fC[8][2] = 1.0000;
4174 fC[8][3] = 0.2155;
4175 fC[8][4] = 0.2433;
4176
4177 fC[9][0] = 0.005;
4178 fC[9][1] = 0.005;
4179 fC[9][2] = 1.0000;
4180 fC[9][3] = 0.2580;
4181 fC[9][4] = 0.3022;
4182
4183 fC[10][0] = 0.005;
4184 fC[10][1] = 0.005;
4185 fC[10][2] = 1.0000;
4186 fC[10][3] = 0.2872;
4187 fC[10][4] = 0.3481;
4188
4189 fC[11][0] = 0.005;
4190 fC[11][1] = 0.005;
4191 fC[11][2] = 1.0000;
4192 fC[11][3] = 0.3170;
4193 fC[11][4] = 0.3847;
4194
4195 fC[12][0] = 0.005;
4196 fC[12][1] = 0.005;
4197 fC[12][2] = 1.0000;
4198 fC[12][3] = 0.3454;
4199 fC[12][4] = 0.4258;
4200
4201 fC[13][0] = 0.005;
4202 fC[13][1] = 0.005;
4203 fC[13][2] = 1.0000;
4204 fC[13][3] = 0.3580;
4205 fC[13][4] = 0.4299;
4206
4207 fC[14][0] = 0.005;
4208 fC[14][1] = 0.005;
4209 fC[14][2] = 1.0000;
4210 fC[14][3] = 0.3903;
4211 fC[14][4] = 0.4326;
4212
4213 fC[15][0] = 0.005;
4214 fC[15][1] = 0.005;
4215 fC[15][2] = 1.0000;
4216 fC[15][3] = 0.3690;
4217 fC[15][4] = 0.4491;
4218
4219 fC[16][0] = 0.005;
4220 fC[16][1] = 0.005;
4221 fC[16][2] = 1.0000;
4222 fC[16][3] = 0.4716;
4223 fC[16][4] = 0.4298;
4224
4225 fC[17][0] = 0.005;
4226 fC[17][1] = 0.005;
4227 fC[17][2] = 1.0000;
4228 fC[17][3] = 0.3875;
4229 fC[17][4] = 0.4083;
4230 }
4231 // 70-80%
4232 else if(centrCur < 80){
4233 fC[0][0] = 0.005;
4234 fC[0][1] = 0.005;
4235 fC[0][2] = 1.0000;
4236 fC[0][3] = 0.0010;
4237 fC[0][4] = 0.0010;
4238
4239 fC[1][0] = 0.005;
4240 fC[1][1] = 0.005;
4241 fC[1][2] = 1.0000;
4242 fC[1][3] = 0.0075;
4243 fC[1][4] = 0.0007;
4244
4245 fC[2][0] = 0.005;
4246 fC[2][1] = 0.005;
4247 fC[2][2] = 1.0000;
4248 fC[2][3] = 0.0313;
4249 fC[2][4] = 0.0124;
4250
4251 fC[3][0] = 0.005;
4252 fC[3][1] = 0.005;
4253 fC[3][2] = 1.0000;
4254 fC[3][3] = 0.0640;
4255 fC[3][4] = 0.0539;
4256
4257 fC[4][0] = 0.005;
4258 fC[4][1] = 0.005;
4259 fC[4][2] = 1.0000;
4260 fC[4][3] = 0.0923;
4261 fC[4][4] = 0.0992;
4262
4263 fC[5][0] = 0.005;
4264 fC[5][1] = 0.005;
4265 fC[5][2] = 1.0000;
4266 fC[5][3] = 0.1202;
4267 fC[5][4] = 0.1417;
4268
4269 fC[6][0] = 0.005;
4270 fC[6][1] = 0.005;
4271 fC[6][2] = 1.0000;
4272 fC[6][3] = 0.1413;
4273 fC[6][4] = 0.1729;
4274
4275 fC[7][0] = 0.005;
4276 fC[7][1] = 0.005;
4277 fC[7][2] = 1.0000;
4278 fC[7][3] = 0.1705;
4279 fC[7][4] = 0.1999;
4280
4281 fC[8][0] = 0.005;
4282 fC[8][1] = 0.005;
4283 fC[8][2] = 1.0000;
4284 fC[8][3] = 0.2103;
4285 fC[8][4] = 0.2472;
4286
4287 fC[9][0] = 0.005;
4288 fC[9][1] = 0.005;
4289 fC[9][2] = 1.0000;
4290 fC[9][3] = 0.2373;
4291 fC[9][4] = 0.2916;
4292
4293 fC[10][0] = 0.005;
4294 fC[10][1] = 0.005;
4295 fC[10][2] = 1.0000;
4296 fC[10][3] = 0.2824;
4297 fC[10][4] = 0.3323;
4298
4299 fC[11][0] = 0.005;
4300 fC[11][1] = 0.005;
4301 fC[11][2] = 1.0000;
4302 fC[11][3] = 0.3046;
4303 fC[11][4] = 0.3576;
4304
4305 fC[12][0] = 0.005;
4306 fC[12][1] = 0.005;
4307 fC[12][2] = 1.0000;
4308 fC[12][3] = 0.3585;
4309 fC[12][4] = 0.4003;
4310
4311 fC[13][0] = 0.005;
4312 fC[13][1] = 0.005;
4313 fC[13][2] = 1.0000;
4314 fC[13][3] = 0.3461;
4315 fC[13][4] = 0.3982;
4316
4317 fC[14][0] = 0.005;
4318 fC[14][1] = 0.005;
4319 fC[14][2] = 1.0000;
4320 fC[14][3] = 0.3362;
4321 fC[14][4] = 0.3776;
4322
4323 fC[15][0] = 0.005;
4324 fC[15][1] = 0.005;
4325 fC[15][2] = 1.0000;
4326 fC[15][3] = 0.3071;
4327 fC[15][4] = 0.3500;
4328
4329 fC[16][0] = 0.005;
4330 fC[16][1] = 0.005;
4331 fC[16][2] = 1.0000;
4332 fC[16][3] = 0.2914;
4333 fC[16][4] = 0.3937;
4334
4335 fC[17][0] = 0.005;
4336 fC[17][1] = 0.005;
4337 fC[17][2] = 1.0000;
4338 fC[17][3] = 0.3727;
4339 fC[17][4] = 0.3877;
4340 }
4341 // 80-100%
4342 else{
4343 fC[0][0] = 0.005;
4344 fC[0][1] = 0.005;
4345 fC[0][2] = 1.0000;
4346 fC[0][3] = 0.0010;
4347 fC[0][4] = 0.0010;
4348
4349 fC[1][0] = 0.005;
4350 fC[1][1] = 0.005;
4351 fC[1][2] = 1.0000;
4352 fC[1][3] = 0.0060;
4353 fC[1][4] = 0.0035;
4354
4355 fC[2][0] = 0.005;
4356 fC[2][1] = 0.005;
4357 fC[2][2] = 1.0000;
4358 fC[2][3] = 0.0323;
4359 fC[2][4] = 0.0113;
4360
4361 fC[3][0] = 0.005;
4362 fC[3][1] = 0.005;
4363 fC[3][2] = 1.0000;
4364 fC[3][3] = 0.0609;
4365 fC[3][4] = 0.0653;
4366
4367 fC[4][0] = 0.005;
4368 fC[4][1] = 0.005;
4369 fC[4][2] = 1.0000;
4370 fC[4][3] = 0.0922;
4371 fC[4][4] = 0.1076;
4372
4373 fC[5][0] = 0.005;
4374 fC[5][1] = 0.005;
4375 fC[5][2] = 1.0000;
4376 fC[5][3] = 0.1096;
4377 fC[5][4] = 0.1328;
4378
4379 fC[6][0] = 0.005;
4380 fC[6][1] = 0.005;
4381 fC[6][2] = 1.0000;
4382 fC[6][3] = 0.1495;
4383 fC[6][4] = 0.1779;
4384
4385 fC[7][0] = 0.005;
4386 fC[7][1] = 0.005;
4387 fC[7][2] = 1.0000;
4388 fC[7][3] = 0.1519;
4389 fC[7][4] = 0.1989;
4390
4391 fC[8][0] = 0.005;
4392 fC[8][1] = 0.005;
4393 fC[8][2] = 1.0000;
4394 fC[8][3] = 0.1817;
4395 fC[8][4] = 0.2472;
4396
4397 fC[9][0] = 0.005;
4398 fC[9][1] = 0.005;
4399 fC[9][2] = 1.0000;
4400 fC[9][3] = 0.2429;
4401 fC[9][4] = 0.2684;
4402
4403 fC[10][0] = 0.005;
4404 fC[10][1] = 0.005;
4405 fC[10][2] = 1.0000;
4406 fC[10][3] = 0.2760;
4407 fC[10][4] = 0.3098;
4408
4409 fC[11][0] = 0.005;
4410 fC[11][1] = 0.005;
4411 fC[11][2] = 1.0000;
4412 fC[11][3] = 0.2673;
4413 fC[11][4] = 0.3198;
4414
4415 fC[12][0] = 0.005;
4416 fC[12][1] = 0.005;
4417 fC[12][2] = 1.0000;
4418 fC[12][3] = 0.3165;
4419 fC[12][4] = 0.3564;
4420
4421 fC[13][0] = 0.005;
4422 fC[13][1] = 0.005;
4423 fC[13][2] = 1.0000;
4424 fC[13][3] = 0.3526;
4425 fC[13][4] = 0.3011;
4426
4427 fC[14][0] = 0.005;
4428 fC[14][1] = 0.005;
4429 fC[14][2] = 1.0000;
4430 fC[14][3] = 0.3788;
4431 fC[14][4] = 0.3011;
4432
4433 fC[15][0] = 0.005;
4434 fC[15][1] = 0.005;
4435 fC[15][2] = 1.0000;
4436 fC[15][3] = 0.3788;
4437 fC[15][4] = 0.3011;
4438
4439 fC[16][0] = 0.005;
4440 fC[16][1] = 0.005;
4441 fC[16][2] = 1.0000;
4442 fC[16][3] = 0.3788;
4443 fC[16][4] = 0.3011;
4444
4445 fC[17][0] = 0.005;
4446 fC[17][1] = 0.005;
4447 fC[17][2] = 1.0000;
4448 fC[17][3] = 0.3788;
4449 fC[17][4] = 0.3011;
4450 }
4451
4452 for(Int_t i=18;i<fgkPIDptBin;i++){
4453 fBinLimitPID[i] = 3.0 + 0.2 * (i-18);
4454 fC[i][0] = fC[17][0];
4455 fC[i][1] = fC[17][1];
4456 fC[i][2] = fC[17][2];
4457 fC[i][3] = fC[17][3];
4458 fC[i][4] = fC[17][4];
4459 }
4460}
4461
4462//---------------------------------------------------------------//
0d16f553 4463Bool_t AliFlowTrackCuts::TPCTOFagree(const AliVTrack *track)
875ee12a 4464{
3c67c769 4465 //check pid agreement between TPC and TOF
875ee12a 4466 Bool_t status = kFALSE;
0d16f553 4467
4468 const Float_t c = 2.99792457999999984e-02;
4469
875ee12a 4470 Float_t mass[5] = {5.10998909999999971e-04,1.05658000000000002e-01,1.39570000000000000e-01,4.93676999999999977e-01,9.38271999999999995e-01};
4471
4472
80efa607 4473 Double_t exptimes[9];
875ee12a 4474 track->GetIntegratedTimes(exptimes);
4475
4476 Float_t dedx = track->GetTPCsignal();
4477
4478 Float_t p = track->P();
4479 Float_t time = track->GetTOFsignal()- fESDpid.GetTOFResponse().GetStartTime(p);
0d16f553 4480 Float_t tl = exptimes[0]*c; // to work both for ESD and AOD
875ee12a 4481
4482 Float_t betagammares = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[4], mass[4]);
4483
4484 Float_t betagamma1 = tl/(time-5 *betagammares) * 33.3564095198152043;
4485
4486// printf("betagamma1 = %f\n",betagamma1);
4487
4488 if(betagamma1 < 0.1) betagamma1 = 0.1;
4489
4490 if(betagamma1 < 0.99999) betagamma1 /= TMath::Sqrt(1-betagamma1*betagamma1);
4491 else betagamma1 = 100;
4492
4493 Float_t betagamma2 = tl/(time+5 *betagammares) * 33.3564095198152043;
4494// printf("betagamma2 = %f\n",betagamma2);
4495
4496 if(betagamma2 < 0.1) betagamma2 = 0.1;
4497
4498 if(betagamma2 < 0.99999) betagamma2 /= TMath::Sqrt(1-betagamma2*betagamma2);
4499 else betagamma2 = 100;
4500
4501
0d16f553 4502 Float_t momtpc=track->GetTPCmomentum();
875ee12a 4503
4504 for(Int_t i=0;i < 5;i++){
4505 Float_t resolutionTOF = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[i], mass[i]);
4506 if(TMath::Abs(exptimes[i] - time) < 5 * resolutionTOF){
4507 Float_t dedxExp = 0;
4508 if(i==0) dedxExp = fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kElectron);
4509 else if(i==1) dedxExp = fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kMuon);
4510 else if(i==2) dedxExp = fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kPion);
4511 else if(i==3) dedxExp = fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kKaon);
4512 else if(i==4) dedxExp = fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kProton);
4513
4514 Float_t resolutionTPC = 2;
4515 if(i==0) resolutionTPC = fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kElectron);
4516 else if(i==1) resolutionTPC = fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kMuon);
4517 else if(i==2) resolutionTPC = fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kPion);
4518 else if(i==3) resolutionTPC = fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kKaon);
4519 else if(i==4) resolutionTPC = fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kProton);
4520
4521 if(TMath::Abs(dedx - dedxExp) < 3 * resolutionTPC){
4522 status = kTRUE;
4523 }
4524 }
4525 }
4526
4527 Float_t bb1 = fESDpid.GetTPCResponse().Bethe(betagamma1);
4528 Float_t bb2 = fESDpid.GetTPCResponse().Bethe(betagamma2);
4529 Float_t bbM = fESDpid.GetTPCResponse().Bethe((betagamma1+betagamma2)*0.5);
4530
4531
4532 // status = kFALSE;
4533 // for nuclei
4534 Float_t resolutionTOFpr = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[4], mass[4]);
4535 Float_t resolutionTPCpr = fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kProton);
4536 if(TMath::Abs(dedx-bb1) < resolutionTPCpr*3 && exptimes[4] < time-7*resolutionTOFpr){
4537 status = kTRUE;
4538 }
4539 else if(TMath::Abs(dedx-bb2) < resolutionTPCpr*3 && exptimes[4] < time-7*resolutionTOFpr){
4540 status = kTRUE;
4541 }
4542 else if(TMath::Abs(dedx-bbM) < resolutionTPCpr*3 && exptimes[4] < time-7*resolutionTOFpr){
4543 status = kTRUE;
4544 }
4545
4546 return status;
3abf7ecc 4547}
4548// end part added by F. Noferini
4549//-----------------------------------------------------------------------
2b1eaa10 4550
2b1eaa10 4551//-----------------------------------------------------------------------
4552const char* AliFlowTrackCuts::PIDsourceName(PIDsource s)
4553{
4554 //get the name of the particle id source
4555 switch (s)
4556 {
4557 case kTPCdedx:
4558 return "TPCdedx";
b9cf8f8e 4559 case kTPCbayesian:
4560 return "TPCbayesian";
2b1eaa10 4561 case kTOFbeta:
4562 return "TOFbeta";
4563 case kTPCpid:
4564 return "TPCpid";
4565 case kTOFpid:
4566 return "TOFpid";
4567 case kTOFbayesian:
4568 return "TOFbayesianPID";
1a80f9f6 4569 case kTOFbetaSimple:
4570 return "TOFbetaSimple";
8eca5d19 4571 case kTPCNuclei:
4572 return "TPCnuclei";
2b1eaa10 4573 default:
4574 return "NOPID";
4575 }
4576}
4577
4578//-----------------------------------------------------------------------
4579const char* AliFlowTrackCuts::GetParamTypeName(trackParameterType type)
4580{
4581 //return the name of the selected parameter type
4582 switch (type)
4583 {
4584 case kMC:
4585 return "MC";
4586 case kGlobal:
1a80f9f6 4587 return "Global";
4588 case kTPCstandalone:
4589 return "TPCstandalone";
4590 case kSPDtracklet:
4591 return "SPDtracklets";
d148af7e 4592 case kPMD:
4593 return "PMD";
c27333c2
MK		 4594 case kVZERO:
4595 return "VZERO";
2e5052c5 4596 case kMUON: // XZhang 20120604
4597 return "MUON"; // XZhang 20120604
dd294bb6
MK		 4598 case kKink:
4599 return "Kink";
4600 case kV0:
4601 return "V0";
2b1eaa10 4602 default:
4603 return "unknown";
4604 }
4605}
4606
d148af7e 4607//-----------------------------------------------------------------------
3c67c769 4608Bool_t AliFlowTrackCuts::PassesPMDcuts(const AliESDPmdTrack* track )
d148af7e 4609{
4610 //check PMD specific cuts
4611 //clean up from last iteration, and init label
1a80f9f6 4612 Int_t det = track->GetDetector();
4613 //Int_t smn = track->GetSmn();
4614 Float_t clsX = track->GetClusterX();
4615 Float_t clsY = track->GetClusterY();
4616 Float_t clsZ = track->GetClusterZ();
4617 Float_t ncell = track->GetClusterCells();
4618 Float_t adc = track->GetClusterADC();
d148af7e 4619
4620 fTrack = NULL;
4621 fMCparticle=NULL;
2281e7c5 4622 fTrackLabel=-996;
d148af7e 4623
1a80f9f6 4624 fTrackEta = GetPmdEta(clsX,clsY,clsZ);
4625 fTrackPhi = GetPmdPhi(clsX,clsY);
d148af7e 4626 fTrackWeight = 1.0;
4627
4628 Bool_t pass=kTRUE;
4629 if (fCutEta) {if ( fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) pass = kFALSE;}
4630 if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) pass = kFALSE;}
1a80f9f6 4631 if (fCutPmdDet) {if(det != fPmdDet) pass = kFALSE;}
4632 if (fCutPmdAdc) {if(adc < fPmdAdc) pass = kFALSE;}
4633 if (fCutPmdNcell) {if(ncell < fPmdNcell) pass = kFALSE;}
d148af7e 4634
1a80f9f6 4635 return pass;
d148af7e 4636}
22289738 4637
4638//-----------------------------------------------------------------------
c27333c2 4639Bool_t AliFlowTrackCuts::PassesVZEROcuts(Int_t id)
22289738 4640{
c27333c2 4641 //check VZERO cuts
b35580e6 4642 if (id<0) return kFALSE;
4643
22289738 4644 //clean up from last iter
dd294bb6 4645 ClearTrack();
22289738 4646
1a80f9f6 4647 fTrackPhi = TMath::PiOver4()*(0.5+id%8);
647a09d3 4648
86a97121 4649 // 10102013 weighting vzero tiles - rbertens@cern.ch
c27333c2 4650 if(!fVZEROgainEqualization) {
86a97121 4651 // if for some reason the equalization is not initialized (e.g. 2011 data)
c27333c2 4652 // the fVZEROxpol[] weights are used to enable or disable vzero rings
54dd223b 4653 if(id<32) { // v0c side
86a97121 4654 fTrackEta = -3.45+0.5*(id/8);
c27333c2
MK		 4655 if(id < 8) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROCpol[0];
4656 else if (id < 16 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROCpol[1];
4657 else if (id < 24 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROCpol[2];
4658 else if (id < 32 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROCpol[3];
86a97121 4659 } else { // v0a side
4660 fTrackEta = +4.8-0.6*((id/8)-4);
c27333c2
MK		 4661 if( id < 40) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROApol[0];
4662 else if ( id < 48 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROApol[1];
4663 else if ( id < 56 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROApol[2];
4664 else if ( id < 64 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fVZEROApol[3];
86a97121 4665 }
4666 } else { // the equalization is initialized
4667 // note that disabled rings have already been excluded on calibration level in
c27333c2 4668 // AliFlowEvent (so for a disabled ring, fVZEROxpol is zero
86a97121 4669 if(id<32) { // v0c side
4670 fTrackEta = -3.45+0.5*(id/8);
c27333c2
MK		 4671 if(id < 8) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROCpol[0]/fVZEROgainEqualization->GetBinContent(1+id);
4672 else if (id < 16 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROCpol[1]/fVZEROgainEqualization->GetBinContent(1+id);
4673 else if (id < 24 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROCpol[2]/fVZEROgainEqualization->GetBinContent(1+id);
4674 else if (id < 32 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROCpol[3]/fVZEROgainEqualization->GetBinContent(1+id);
86a97121 4675 } else { // v0a side
4676 fTrackEta = +4.8-0.6*((id/8)-4);
c27333c2
MK		 4677 if( id < 40) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROApol[0]/fVZEROgainEqualization->GetBinContent(1+id);
4678 else if ( id < 48 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROApol[1]/fVZEROgainEqualization->GetBinContent(1+id);
4679 else if ( id < 56 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROApol[2]/fVZEROgainEqualization->GetBinContent(1+id);
4680 else if ( id < 64 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fVZEROApol[3]/fVZEROgainEqualization->GetBinContent(1+id);
86a97121 4681 }
4682 // printf ( " tile %i and weight %.2f \n", id, fTrackWeight);
4683 }
647a09d3 4684
1d78aabb 4685 if (fLinearizeVZEROresponse && id < 64)
92b507bb 4686 {
4687 //this is only needed in pass1 of LHC10h
c27333c2 4688 Float_t multVZERO[fgkNumberOfVZEROtracks];
92b507bb 4689 Float_t dummy=0.0;
cac040c5 4690 AliESDUtils::GetCorrV0((AliESDEvent*)fEvent,dummy,multVZERO);
c27333c2 4691 fTrackWeight = multVZERO[id];
92b507bb 4692 }
22289738 4693
4694 Bool_t pass=kTRUE;
4695 if (fCutEta) {if ( fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) pass = kFALSE;}
4696 if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) pass = kFALSE;}
4697
1a80f9f6 4698 return pass;
22289738 4699}
4700
2e5052c5 4701//-----------------------------------------------------------------------------
4702Bool_t AliFlowTrackCuts::PassesMuonCuts(AliVParticle* vparticle)
4703{
4704// XZhang 20120604
dd294bb6 4705 ClearTrack();
2e5052c5 4706 fTrackLabel = (fFakesAreOK)?TMath::Abs(vparticle->GetLabel()):vparticle->GetLabel();
4707 if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
4708 else fMCparticle=NULL;
4709
4710 AliESDMuonTrack *esdTrack = dynamic_cast<AliESDMuonTrack*>(vparticle);
4711 AliAODTrack *aodTrack = dynamic_cast<AliAODTrack*>(vparticle);
4712 if ((!esdTrack) && (!aodTrack)) return kFALSE;
4713 if (!fMuonTrackCuts->IsSelected(vparticle)) return kFALSE;
4714 HandleVParticle(vparticle); if (!fTrack) return kFALSE;
4715 return kTRUE;
4716}
4717
1a80f9f6 4718//----------------------------------------------------------------------------//
4719Double_t AliFlowTrackCuts::GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos)
4720{
3c67c769 4721 //get PMD "track" eta coordinate
1a80f9f6 4722 Float_t rpxpy, theta, eta;
4723 rpxpy = TMath::Sqrt(xPos*xPos + yPos*yPos);
4724 theta = TMath::ATan2(rpxpy,zPos);
4725 eta = -TMath::Log(TMath::Tan(0.5*theta));
4726 return eta;
4727}
4728
4729//--------------------------------------------------------------------------//
4730Double_t AliFlowTrackCuts::GetPmdPhi(Float_t xPos, Float_t yPos)
4731{
3c67c769 4732 //Get PMD "track" phi coordinate
1a80f9f6 4733 Float_t pybypx, phi = 0., phi1;
4734 if(xPos==0)
4735 {
4736 if(yPos>0) phi = 90.;
4737 if(yPos<0) phi = 270.;
4738 }
4739 if(xPos != 0)
4740 {
4741 pybypx = yPos/xPos;
4742 if(pybypx < 0) pybypx = - pybypx;
4743 phi1 = TMath::ATan(pybypx)*180./3.14159;
4744
4745 if(xPos > 0 && yPos > 0) phi = phi1; // 1st Quadrant
4746 if(xPos < 0 && yPos > 0) phi = 180 - phi1; // 2nd Quadrant
4747 if(xPos < 0 && yPos < 0) phi = 180 + phi1; // 3rd Quadrant
4748 if(xPos > 0 && yPos < 0) phi = 360 - phi1; // 4th Quadrant
4749
4750 }
4751 phi = phi*3.14159/180.;
4752 return phi;
4753}
499fe731 4754
4755//---------------------------------------------------------------//
4756void AliFlowTrackCuts::Browse(TBrowser* b)
4757{
4758 //some browsing capabilities
4759 if (fQA) b->Add(fQA);
4760}
4761
1a80f9f6 4762//---------------------------------------------------------------//
499fe731 4763Long64_t AliFlowTrackCuts::Merge(TCollection* list)
4764{
4765 //merge
8eca5d19 4766 if (!fQA || !list) return 0;
4767 if (list->IsEmpty()) return 0;
4768 AliFlowTrackCuts* obj=NULL;
4769 TList tmplist;
499fe731 4770 TIter next(list);
4771 while ( (obj = dynamic_cast<AliFlowTrackCuts*>(next())) )
4772 {
4773 if (obj==this) continue;
8eca5d19 4774 tmplist.Add(obj->GetQA());
499fe731 4775 }
8eca5d19 4776 return fQA->Merge(&tmplist);
499fe731 4777}
f21bdf48 4778
9126fd0c 4779
2281e7c5 4780