| | 1 | /************************************************************************** |
| | 2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * |
| | 3 | * * |
| | 4 | * Author: The ALICE Off-line Project. * |
| | 5 | * Contributors are mentioned in the code where appropriate. * |
| | 6 | * * |
| | 7 | * Permission to use, copy, modify and distribute this software and its * |
| | 8 | * documentation strictly for non-commercial purposes is hereby granted * |
| | 9 | * without fee, provided that the above copyright notice appears in all * |
| | 10 | * copies and that both the copyright notice and this permission notice * |
| | 11 | * appear in the supporting documentation. The authors make no claims * |
| | 12 | * about the suitability of this software for any purpose. It is * |
| | 13 | * provided "as is" without express or implied warranty. * |
| | 14 | **************************************************************************/ |
| | 15 | |
| | 16 | /* $Id$ */ |
| | 17 | |
| | 18 | // AliFlowTrackCuts: |
| | 19 | // ESD track cuts for flow framework |
| | 20 | // |
| | 21 | // origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch) |
| | 22 | // |
| | 23 | // This class gurantees consistency of cut methods, trackparameter |
| | 24 | // selection (global tracks, TPC only, etc..) and parameter mixing |
| | 25 | // in the flow framework. Transparently handles different input types: |
| | 26 | // ESD, MC, AOD. |
| | 27 | // This class works in 2 steps: first the requested track parameters are |
| | 28 | // constructed (to be set by SetParamType() ), then cuts are applied. |
| | 29 | // the constructed track can be requested AFTER checking the cuts by |
| | 30 | // calling GetTrack(), in this case the cut object stays in control, |
| | 31 | // caller does not have to delete the track. |
| | 32 | // Additionally caller can request an AliFlowTrack object to be constructed |
| | 33 | // according the parameter mixing scenario requested by SetParamMix(). |
| | 34 | // AliFlowTrack is made using MakeFlowTrack() method, its an 'object factory' |
| | 35 | // so caller needs to take care of the freshly created object. |
| | 36 | |
| | 37 | #include <limits.h> |
| | 38 | #include <float.h> |
| | 39 | #include <TMatrix.h> |
| | 40 | #include "TParticle.h" |
| | 41 | #include "TObjArray.h" |
| | 42 | #include "AliStack.h" |
| | 43 | #include "AliMCEvent.h" |
| | 44 | #include "AliESDEvent.h" |
| | 45 | #include "AliVParticle.h" |
| | 46 | #include "AliMCParticle.h" |
| | 47 | #include "AliESDtrack.h" |
| | 48 | #include "AliMultiplicity.h" |
| | 49 | #include "AliAODTrack.h" |
| | 50 | #include "AliFlowTrack.h" |
| | 51 | #include "AliFlowTrackCuts.h" |
| | 52 | #include "AliLog.h" |
| | 53 | #include "AliESDpid.h" |
| | 54 | |
| | 55 | ClassImp(AliFlowTrackCuts) |
| | 56 | |
| | 57 | //----------------------------------------------------------------------- |
| | 58 | AliFlowTrackCuts::AliFlowTrackCuts(): |
| | 59 | AliFlowTrackSimpleCuts(), |
| | 60 | fAliESDtrackCuts(NULL), |
| | 61 | fQA(NULL), |
| | 62 | fCutMC(kFALSE), |
| | 63 | fCutMCprocessType(kFALSE), |
| | 64 | fMCprocessType(kPNoProcess), |
| | 65 | fCutMCPID(kFALSE), |
| | 66 | fMCPID(0), |
| | 67 | fIgnoreSignInPID(kFALSE), |
| | 68 | fCutMCisPrimary(kFALSE), |
| | 69 | fRequireTransportBitForPrimaries(kTRUE), |
| | 70 | fMCisPrimary(kFALSE), |
| | 71 | fRequireCharge(kFALSE), |
| | 72 | fFakesAreOK(kTRUE), |
| | 73 | fCutSPDtrackletDeltaPhi(kFALSE), |
| | 74 | fSPDtrackletDeltaPhiMax(FLT_MAX), |
| | 75 | fSPDtrackletDeltaPhiMin(-FLT_MAX), |
| | 76 | fIgnoreTPCzRange(kFALSE), |
| | 77 | fIgnoreTPCzRangeMax(FLT_MAX), |
| | 78 | fIgnoreTPCzRangeMin(-FLT_MAX), |
| | 79 | fCutChi2PerClusterTPC(kFALSE), |
| | 80 | fMaxChi2PerClusterTPC(FLT_MAX), |
| | 81 | fMinChi2PerClusterTPC(-FLT_MAX), |
| | 82 | fCutNClustersTPC(kFALSE), |
| | 83 | fNClustersTPCMax(INT_MAX), |
| | 84 | fNClustersTPCMin(INT_MIN), |
| | 85 | fCutNClustersITS(kFALSE), |
| | 86 | fNClustersITSMax(INT_MAX), |
| | 87 | fNClustersITSMin(INT_MIN), |
| | 88 | fUseAODFilterBit(kFALSE), |
| | 89 | fAODFilterBit(0), |
| | 90 | fCutDCAToVertexXY(kFALSE), |
| | 91 | fCutDCAToVertexZ(kFALSE), |
| | 92 | fParamType(kGlobal), |
| | 93 | fParamMix(kPure), |
| | 94 | fTrack(NULL), |
| | 95 | fTrackPhi(0.), |
| | 96 | fTrackEta(0.), |
| | 97 | fTrackWeight(0.), |
| | 98 | fTrackLabel(INT_MIN), |
| | 99 | fMCevent(NULL), |
| | 100 | fMCparticle(NULL), |
| | 101 | fEvent(NULL), |
| | 102 | fTPCtrack(), |
| | 103 | fESDpid(), |
| | 104 | fPIDsource(kTPCTOFpid), |
| | 105 | fTPCpidCuts(NULL), |
| | 106 | fTOFpidCuts(NULL), |
| | 107 | fTPCTOFpidCrossOverPt(0.4), |
| | 108 | fAliPID(AliPID::kPion) |
| | 109 | { |
| | 110 | //io constructor |
| | 111 | SetPriors(); //init arrays |
| | 112 | } |
| | 113 | |
| | 114 | //----------------------------------------------------------------------- |
| | 115 | AliFlowTrackCuts::AliFlowTrackCuts(const char* name): |
| | 116 | AliFlowTrackSimpleCuts(), |
| | 117 | fAliESDtrackCuts(new AliESDtrackCuts()), |
| | 118 | fQA(NULL), |
| | 119 | fCutMC(kFALSE), |
| | 120 | fCutMCprocessType(kFALSE), |
| | 121 | fMCprocessType(kPNoProcess), |
| | 122 | fCutMCPID(kFALSE), |
| | 123 | fMCPID(0), |
| | 124 | fIgnoreSignInPID(kFALSE), |
| | 125 | fCutMCisPrimary(kFALSE), |
| | 126 | fRequireTransportBitForPrimaries(kTRUE), |
| | 127 | fMCisPrimary(kFALSE), |
| | 128 | fRequireCharge(kFALSE), |
| | 129 | fFakesAreOK(kTRUE), |
| | 130 | fCutSPDtrackletDeltaPhi(kFALSE), |
| | 131 | fSPDtrackletDeltaPhiMax(FLT_MAX), |
| | 132 | fSPDtrackletDeltaPhiMin(-FLT_MAX), |
| | 133 | fIgnoreTPCzRange(kFALSE), |
| | 134 | fIgnoreTPCzRangeMax(FLT_MAX), |
| | 135 | fIgnoreTPCzRangeMin(-FLT_MAX), |
| | 136 | fCutChi2PerClusterTPC(kFALSE), |
| | 137 | fMaxChi2PerClusterTPC(FLT_MAX), |
| | 138 | fMinChi2PerClusterTPC(-FLT_MAX), |
| | 139 | fCutNClustersTPC(kFALSE), |
| | 140 | fNClustersTPCMax(INT_MAX), |
| | 141 | fNClustersTPCMin(INT_MIN), |
| | 142 | fCutNClustersITS(kFALSE), |
| | 143 | fNClustersITSMax(INT_MAX), |
| | 144 | fNClustersITSMin(INT_MIN), |
| | 145 | fUseAODFilterBit(kFALSE), |
| | 146 | fAODFilterBit(0), |
| | 147 | fCutDCAToVertexXY(kFALSE), |
| | 148 | fCutDCAToVertexZ(kFALSE), |
| | 149 | fParamType(kGlobal), |
| | 150 | fParamMix(kPure), |
| | 151 | fTrack(NULL), |
| | 152 | fTrackPhi(0.), |
| | 153 | fTrackEta(0.), |
| | 154 | fTrackWeight(0.), |
| | 155 | fTrackLabel(INT_MIN), |
| | 156 | fMCevent(NULL), |
| | 157 | fMCparticle(NULL), |
| | 158 | fEvent(NULL), |
| | 159 | fTPCtrack(), |
| | 160 | fESDpid(), |
| | 161 | fPIDsource(kTPCTOFpid), |
| | 162 | fTPCpidCuts(NULL), |
| | 163 | fTOFpidCuts(NULL), |
| | 164 | fTPCTOFpidCrossOverPt(0.4), |
| | 165 | fAliPID(AliPID::kPion) |
| | 166 | { |
| | 167 | //constructor |
| | 168 | SetName(name); |
| | 169 | SetTitle("AliFlowTrackCuts"); |
| | 170 | fESDpid.GetTPCResponse().SetBetheBlochParameters( 0.0283086, |
| | 171 | 2.63394e+01, |
| | 172 | 5.04114e-11, |
| | 173 | 2.12543e+00, |
| | 174 | 4.88663e+00 ); |
| | 175 | SetPriors(); //init arrays |
| | 176 | } |
| | 177 | |
| | 178 | //----------------------------------------------------------------------- |
| | 179 | AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that): |
| | 180 | AliFlowTrackSimpleCuts(that), |
| | 181 | fAliESDtrackCuts(NULL), |
| | 182 | fQA(NULL), |
| | 183 | fCutMC(that.fCutMC), |
| | 184 | fCutMCprocessType(that.fCutMCprocessType), |
| | 185 | fMCprocessType(that.fMCprocessType), |
| | 186 | fCutMCPID(that.fCutMCPID), |
| | 187 | fMCPID(that.fMCPID), |
| | 188 | fIgnoreSignInPID(that.fIgnoreSignInPID), |
| | 189 | fCutMCisPrimary(that.fCutMCisPrimary), |
| | 190 | fRequireTransportBitForPrimaries(that.fRequireTransportBitForPrimaries), |
| | 191 | fMCisPrimary(that.fMCisPrimary), |
| | 192 | fRequireCharge(that.fRequireCharge), |
| | 193 | fFakesAreOK(that.fFakesAreOK), |
| | 194 | fCutSPDtrackletDeltaPhi(that.fCutSPDtrackletDeltaPhi), |
| | 195 | fSPDtrackletDeltaPhiMax(that.fSPDtrackletDeltaPhiMax), |
| | 196 | fSPDtrackletDeltaPhiMin(that.fSPDtrackletDeltaPhiMin), |
| | 197 | fIgnoreTPCzRange(that.fIgnoreTPCzRange), |
| | 198 | fIgnoreTPCzRangeMax(that.fIgnoreTPCzRangeMax), |
| | 199 | fIgnoreTPCzRangeMin(that.fIgnoreTPCzRangeMin), |
| | 200 | fCutChi2PerClusterTPC(that.fCutChi2PerClusterTPC), |
| | 201 | fMaxChi2PerClusterTPC(that.fMaxChi2PerClusterTPC), |
| | 202 | fMinChi2PerClusterTPC(that.fMinChi2PerClusterTPC), |
| | 203 | fCutNClustersTPC(that.fCutNClustersTPC), |
| | 204 | fNClustersTPCMax(that.fNClustersTPCMax), |
| | 205 | fNClustersTPCMin(that.fNClustersTPCMin), |
| | 206 | fCutNClustersITS(that.fCutNClustersITS), |
| | 207 | fNClustersITSMax(that.fNClustersITSMax), |
| | 208 | fNClustersITSMin(that.fNClustersITSMin), |
| | 209 | fUseAODFilterBit(that.fUseAODFilterBit), |
| | 210 | fAODFilterBit(that.fAODFilterBit), |
| | 211 | fCutDCAToVertexXY(that.fCutDCAToVertexXY), |
| | 212 | fCutDCAToVertexZ(that.fCutDCAToVertexZ), |
| | 213 | fParamType(that.fParamType), |
| | 214 | fParamMix(that.fParamMix), |
| | 215 | fTrack(NULL), |
| | 216 | fTrackPhi(0.), |
| | 217 | fTrackEta(0.), |
| | 218 | fTrackWeight(0.), |
| | 219 | fTrackLabel(INT_MIN), |
| | 220 | fMCevent(NULL), |
| | 221 | fMCparticle(NULL), |
| | 222 | fEvent(NULL), |
| | 223 | fTPCtrack(), |
| | 224 | fESDpid(that.fESDpid), |
| | 225 | fPIDsource(that.fPIDsource), |
| | 226 | fTPCpidCuts(NULL), |
| | 227 | fTOFpidCuts(NULL), |
| | 228 | fTPCTOFpidCrossOverPt(that.fTPCTOFpidCrossOverPt), |
| | 229 | fAliPID(that.fAliPID) |
| | 230 | { |
| | 231 | //copy constructor |
| | 232 | if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts)); |
| | 233 | if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts)); |
| | 234 | if (that.fAliESDtrackCuts) fAliESDtrackCuts = new AliESDtrackCuts(*(that.fAliESDtrackCuts)); |
| | 235 | SetPriors(); //init arrays |
| | 236 | } |
| | 237 | |
| | 238 | //----------------------------------------------------------------------- |
| | 239 | AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that) |
| | 240 | { |
| | 241 | //assignment |
| | 242 | AliFlowTrackSimpleCuts::operator=(that); |
| | 243 | if (that.fAliESDtrackCuts) *fAliESDtrackCuts=*(that.fAliESDtrackCuts); |
| | 244 | fQA=NULL; |
| | 245 | fCutMC=that.fCutMC; |
| | 246 | fCutMCprocessType=that.fCutMCprocessType; |
| | 247 | fMCprocessType=that.fMCprocessType; |
| | 248 | fCutMCPID=that.fCutMCPID; |
| | 249 | fMCPID=that.fMCPID; |
| | 250 | fIgnoreSignInPID=that.fIgnoreSignInPID, |
| | 251 | fCutMCisPrimary=that.fCutMCisPrimary; |
| | 252 | fRequireTransportBitForPrimaries=that.fRequireTransportBitForPrimaries; |
| | 253 | fMCisPrimary=that.fMCisPrimary; |
| | 254 | fRequireCharge=that.fRequireCharge; |
| | 255 | fFakesAreOK=that.fFakesAreOK; |
| | 256 | fCutSPDtrackletDeltaPhi=that.fCutSPDtrackletDeltaPhi; |
| | 257 | fSPDtrackletDeltaPhiMax=that.fSPDtrackletDeltaPhiMax; |
| | 258 | fSPDtrackletDeltaPhiMin=that.fSPDtrackletDeltaPhiMin; |
| | 259 | fIgnoreTPCzRange=that.fIgnoreTPCzRange; |
| | 260 | fIgnoreTPCzRangeMax=that.fIgnoreTPCzRangeMax; |
| | 261 | fIgnoreTPCzRangeMin=that.fIgnoreTPCzRangeMin; |
| | 262 | fCutChi2PerClusterTPC=that.fCutChi2PerClusterTPC; |
| | 263 | fMaxChi2PerClusterTPC=that.fMaxChi2PerClusterTPC; |
| | 264 | fMinChi2PerClusterTPC=that.fMinChi2PerClusterTPC; |
| | 265 | fCutNClustersTPC=that.fCutNClustersTPC; |
| | 266 | fNClustersTPCMax=that.fNClustersTPCMax; |
| | 267 | fNClustersTPCMin=that.fNClustersTPCMin; |
| | 268 | fCutNClustersITS=that.fCutNClustersITS; |
| | 269 | fNClustersITSMax=that.fNClustersITSMax; |
| | 270 | fNClustersITSMin=that.fNClustersITSMin; |
| | 271 | fUseAODFilterBit=that.fUseAODFilterBit; |
| | 272 | fAODFilterBit=that.fAODFilterBit; |
| | 273 | fCutDCAToVertexXY=that.fCutDCAToVertexXY; |
| | 274 | fCutDCAToVertexZ=that.fCutDCAToVertexZ; |
| | 275 | fParamType=that.fParamType; |
| | 276 | fParamMix=that.fParamMix; |
| | 277 | |
| | 278 | fTrack=NULL; |
| | 279 | fTrackPhi=0.; |
| | 280 | fTrackPhi=0.; |
| | 281 | fTrackWeight=0.; |
| | 282 | fTrackLabel=INT_MIN; |
| | 283 | fMCevent=NULL; |
| | 284 | fMCparticle=NULL; |
| | 285 | fEvent=NULL; |
| | 286 | |
| | 287 | fESDpid = that.fESDpid; |
| | 288 | fPIDsource = that.fPIDsource; |
| | 289 | |
| | 290 | if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts)); |
| | 291 | if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts)); |
| | 292 | fTPCTOFpidCrossOverPt=that.fTPCTOFpidCrossOverPt; |
| | 293 | |
| | 294 | fAliPID=that.fAliPID; |
| | 295 | |
| | 296 | return *this; |
| | 297 | } |
| | 298 | |
| | 299 | //----------------------------------------------------------------------- |
| | 300 | AliFlowTrackCuts::~AliFlowTrackCuts() |
| | 301 | { |
| | 302 | //dtor |
| | 303 | delete fAliESDtrackCuts; |
| | 304 | delete fTPCpidCuts; |
| | 305 | delete fTOFpidCuts; |
| | 306 | } |
| | 307 | |
| | 308 | //----------------------------------------------------------------------- |
| | 309 | void AliFlowTrackCuts::SetEvent(AliVEvent* event, AliMCEvent* mcEvent) |
| | 310 | { |
| | 311 | //set the event |
| | 312 | Clear(); |
| | 313 | fEvent=event; |
| | 314 | fMCevent=mcEvent; |
| | 315 | |
| | 316 | //do the magic for ESD |
| | 317 | AliESDEvent* myESD = dynamic_cast<AliESDEvent*>(event); |
| | 318 | if (fCutPID && myESD) |
| | 319 | { |
| | 320 | //TODO: maybe call it only for the TOF options? |
| | 321 | // Added by F. Noferini for TOF PID |
| | 322 | fESDpid.SetTOFResponse(myESD,AliESDpid::kTOF_T0); |
| | 323 | fESDpid.MakePID(myESD,kFALSE); |
| | 324 | // End F. Noferini added part |
| | 325 | } |
| | 326 | |
| | 327 | //TODO: AOD |
| | 328 | } |
| | 329 | |
| | 330 | //----------------------------------------------------------------------- |
| | 331 | void AliFlowTrackCuts::SetCutMC( Bool_t b ) |
| | 332 | { |
| | 333 | //will we be cutting on MC information? |
| | 334 | fCutMC=b; |
| | 335 | |
| | 336 | //if we cut on MC info then also the Bethe Bloch should be the one tuned for MC |
| | 337 | if (fCutMC) |
| | 338 | { |
| | 339 | fESDpid.GetTPCResponse().SetBetheBlochParameters( 2.15898e+00/50., |
| | 340 | 1.75295e+01, |
| | 341 | 3.40030e-09, |
| | 342 | 1.96178e+00, |
| | 343 | 3.91720e+00); |
| | 344 | } |
| | 345 | } |
| | 346 | |
| | 347 | //----------------------------------------------------------------------- |
| | 348 | Bool_t AliFlowTrackCuts::IsSelected(TObject* obj, Int_t id) |
| | 349 | { |
| | 350 | //check cuts |
| | 351 | AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj); |
| | 352 | if (vparticle) return PassesCuts(vparticle); |
| | 353 | AliFlowTrackSimple* flowtrack = dynamic_cast<AliFlowTrackSimple*>(obj); |
| | 354 | if (flowtrack) return PassesCuts(flowtrack); |
| | 355 | AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj); |
| | 356 | if (tracklets) return PassesCuts(tracklets,id); |
| | 357 | return kFALSE; //default when passed wrong type of object |
| | 358 | } |
| | 359 | |
| | 360 | //----------------------------------------------------------------------- |
| | 361 | Bool_t AliFlowTrackCuts::IsSelectedMCtruth(TObject* obj, Int_t id) |
| | 362 | { |
| | 363 | //check cuts |
| | 364 | AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj); |
| | 365 | if (vparticle) |
| | 366 | { |
| | 367 | return PassesMCcuts(fMCevent,vparticle->GetLabel()); |
| | 368 | } |
| | 369 | AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj); |
| | 370 | if (tracklets) |
| | 371 | { |
| | 372 | Int_t label0 = tracklets->GetLabel(id,0); |
| | 373 | Int_t label1 = tracklets->GetLabel(id,1); |
| | 374 | Int_t label = (label0==label1)?tracklets->GetLabel(id,1):-666; |
| | 375 | return PassesMCcuts(fMCevent,label); |
| | 376 | } |
| | 377 | return kFALSE; //default when passed wrong type of object |
| | 378 | } |
| | 379 | |
| | 380 | //----------------------------------------------------------------------- |
| | 381 | Bool_t AliFlowTrackCuts::PassesCuts(AliFlowTrackSimple* track) |
| | 382 | { |
| | 383 | //check cuts on a flowtracksimple |
| | 384 | |
| | 385 | //clean up from last iteration |
| | 386 | fTrack = NULL; |
| | 387 | return AliFlowTrackSimpleCuts::PassesCuts(track); |
| | 388 | } |
| | 389 | |
| | 390 | //----------------------------------------------------------------------- |
| | 391 | Bool_t AliFlowTrackCuts::PassesCuts(AliMultiplicity* tracklet, Int_t id) |
| | 392 | { |
| | 393 | //check cuts on a tracklets |
| | 394 | |
| | 395 | //clean up from last iteration, and init label |
| | 396 | fTrack = NULL; |
| | 397 | fMCparticle=NULL; |
| | 398 | fTrackLabel=-1; |
| | 399 | |
| | 400 | fTrackPhi = tracklet->GetPhi(id); |
| | 401 | fTrackEta = tracklet->GetEta(id); |
| | 402 | fTrackWeight = 1.0; |
| | 403 | if (fCutEta) {if ( fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) return kFALSE;} |
| | 404 | if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) return kFALSE;} |
| | 405 | |
| | 406 | //check MC info if available |
| | 407 | //if the 2 clusters have different label track cannot be good |
| | 408 | //and should therefore not pass the mc cuts |
| | 409 | Int_t label0 = tracklet->GetLabel(id,0); |
| | 410 | Int_t label1 = tracklet->GetLabel(id,1); |
| | 411 | //if possible get label and mcparticle |
| | 412 | fTrackLabel = (label0==label1)?tracklet->GetLabel(id,1):-1; |
| | 413 | if (!fFakesAreOK && fTrackLabel<0) return kFALSE; |
| | 414 | if (fTrackLabel>=0 && fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel)); |
| | 415 | //check MC cuts |
| | 416 | if (fCutMC && !PassesMCcuts()) return kFALSE; |
| | 417 | return kTRUE; |
| | 418 | } |
| | 419 | |
| | 420 | //----------------------------------------------------------------------- |
| | 421 | Bool_t AliFlowTrackCuts::PassesMCcuts(AliMCEvent* mcEvent, Int_t label) |
| | 422 | { |
| | 423 | //check the MC info |
| | 424 | if (!mcEvent) return kFALSE; |
| | 425 | if (label<0) return kFALSE;//otherwise AliCMevent prints a warning before returning NULL |
| | 426 | AliMCParticle* mcparticle = static_cast<AliMCParticle*>(mcEvent->GetTrack(label)); |
| | 427 | if (!mcparticle) {AliError("no MC track"); return kFALSE;} |
| | 428 | |
| | 429 | if (fCutMCisPrimary) |
| | 430 | { |
| | 431 | if (IsPhysicalPrimary(mcEvent,label,fRequireTransportBitForPrimaries) != fMCisPrimary) return kFALSE; |
| | 432 | } |
| | 433 | if (fCutMCPID) |
| | 434 | { |
| | 435 | Int_t pdgCode = mcparticle->PdgCode(); |
| | 436 | if (fIgnoreSignInPID) |
| | 437 | { |
| | 438 | if (TMath::Abs(fMCPID) != TMath::Abs(pdgCode)) return kFALSE; |
| | 439 | } |
| | 440 | else |
| | 441 | { |
| | 442 | if (fMCPID != pdgCode) return kFALSE; |
| | 443 | } |
| | 444 | } |
| | 445 | if ( fCutMCprocessType ) |
| | 446 | { |
| | 447 | TParticle* particle = mcparticle->Particle(); |
| | 448 | Int_t processID = particle->GetUniqueID(); |
| | 449 | if (processID != fMCprocessType ) return kFALSE; |
| | 450 | } |
| | 451 | return kTRUE; |
| | 452 | } |
| | 453 | //----------------------------------------------------------------------- |
| | 454 | Bool_t AliFlowTrackCuts::PassesMCcuts() |
| | 455 | { |
| | 456 | if (!fMCevent) return kFALSE; |
| | 457 | if (fTrackLabel<0) return kFALSE;//otherwise AliCMevent prints a warning before returning NULL |
| | 458 | fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel)); |
| | 459 | return PassesMCcuts(fMCevent,fTrackLabel); |
| | 460 | } |
| | 461 | |
| | 462 | //----------------------------------------------------------------------- |
| | 463 | Bool_t AliFlowTrackCuts::PassesCuts(AliVParticle* vparticle) |
| | 464 | { |
| | 465 | //check cuts for an ESD vparticle |
| | 466 | |
| | 467 | //////////////////////////////////////////////////////////////// |
| | 468 | // start by preparing the track parameters to cut on ////////// |
| | 469 | //////////////////////////////////////////////////////////////// |
| | 470 | //clean up from last iteration |
| | 471 | fTrack=NULL; |
| | 472 | |
| | 473 | //get the label and the mc particle |
| | 474 | fTrackLabel = (fFakesAreOK)?TMath::Abs(vparticle->GetLabel()):vparticle->GetLabel(); |
| | 475 | if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel)); |
| | 476 | else fMCparticle=NULL; |
| | 477 | |
| | 478 | Bool_t isMCparticle = kFALSE; //some things are different for MC particles, check! |
| | 479 | AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*>(vparticle); |
| | 480 | AliAODTrack* aodTrack = NULL; |
| | 481 | if (esdTrack) |
| | 482 | //for an ESD track we do some magic sometimes like constructing TPC only parameters |
| | 483 | //or doing some hybrid, handle that here |
| | 484 | HandleESDtrack(esdTrack); |
| | 485 | else |
| | 486 | { |
| | 487 | HandleVParticle(vparticle); |
| | 488 | //now check if produced particle is MC |
| | 489 | isMCparticle = (dynamic_cast<AliMCParticle*>(fTrack))!=NULL; |
| | 490 | aodTrack = dynamic_cast<AliAODTrack*>(vparticle); //keep the additional dynamic cast out of the way for ESDs |
| | 491 | } |
| | 492 | //////////////////////////////////////////////////////////////// |
| | 493 | //////////////////////////////////////////////////////////////// |
| | 494 | |
| | 495 | if (!fTrack) return kFALSE; |
| | 496 | //because it may be different from global, not needed for aodTrack because we dont do anything funky there |
| | 497 | if (esdTrack) esdTrack = static_cast<AliESDtrack*>(fTrack); |
| | 498 | |
| | 499 | Bool_t pass=kTRUE; |
| | 500 | //check the common cuts for the current particle fTrack (MC,AOD,ESD) |
| | 501 | Double_t pt = fTrack->Pt(); |
| | 502 | if (!fFakesAreOK) {if (fTrackLabel<0) pass=kFALSE;} |
| | 503 | if (fCutPt) {if (pt < fPtMin || pt >= fPtMax ) pass=kFALSE;} |
| | 504 | if (fCutEta) {if (fTrack->Eta() < fEtaMin || fTrack->Eta() >= fEtaMax ) pass=kFALSE;} |
| | 505 | if (fCutPhi) {if (fTrack->Phi() < fPhiMin || fTrack->Phi() >= fPhiMax ) pass=kFALSE;} |
| | 506 | if (fRequireCharge) {if (fTrack->Charge() == 0) pass=kFALSE;} |
| | 507 | if (fCutCharge && !isMCparticle) {if (fTrack->Charge() != fCharge) pass=kFALSE;} |
| | 508 | if (fCutCharge && isMCparticle) |
| | 509 | { |
| | 510 | //in case of an MC particle the charge is stored in units of 1/3|e| |
| | 511 | Int_t charge = TMath::Nint(fTrack->Charge()/3.0); //mc particles have charge in units of 1/3e |
| | 512 | if (charge!=fCharge) pass=kFALSE; |
| | 513 | } |
| | 514 | //if(fCutPID) {if (fTrack->PID() != fPID) pass=kFALSE;} |
| | 515 | |
| | 516 | //when additionally MC info is required |
| | 517 | if (fCutMC && !PassesMCcuts()) pass=kFALSE; |
| | 518 | |
| | 519 | //the case of ESD or AOD |
| | 520 | if (esdTrack) PassesESDcuts(esdTrack); |
| | 521 | if (aodTrack) PassesAODcuts(aodTrack); |
| | 522 | |
| | 523 | //true by default, if we didn't set any cuts |
| | 524 | return pass; |
| | 525 | } |
| | 526 | |
| | 527 | //_______________________________________________________________________ |
| | 528 | Bool_t AliFlowTrackCuts::PassesAODcuts(AliAODTrack* track) |
| | 529 | { |
| | 530 | Bool_t pass = kTRUE; |
| | 531 | |
| | 532 | if (fCutNClustersTPC) |
| | 533 | { |
| | 534 | Int_t ntpccls = track->GetTPCNcls(); |
| | 535 | if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE; |
| | 536 | } |
| | 537 | |
| | 538 | if (fCutNClustersITS) |
| | 539 | { |
| | 540 | Int_t nitscls = track->GetITSNcls(); |
| | 541 | if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE; |
| | 542 | } |
| | 543 | |
| | 544 | if (fCutChi2PerClusterTPC) |
| | 545 | { |
| | 546 | Double_t chi2tpc = track->Chi2perNDF(); |
| | 547 | if (chi2tpc < fMinChi2PerClusterTPC || chi2tpc > fMaxChi2PerClusterTPC) pass=kFALSE; |
| | 548 | } |
| | 549 | |
| | 550 | if (GetRequireTPCRefit() && !(track->GetStatus() & AliESDtrack::kTPCrefit) ) pass=kFALSE; |
| | 551 | if (GetRequireITSRefit() && !(track->GetStatus() & AliESDtrack::kITSrefit) ) pass=kFALSE; |
| | 552 | |
| | 553 | if (fUseAODFilterBit && !track->TestFilterBit(fAODFilterBit)) pass=kFALSE; |
| | 554 | |
| | 555 | if (fCutDCAToVertexXY && track->DCA()>GetMaxDCAToVertexXY()) pass=kFALSE; |
| | 556 | |
| | 557 | |
| | 558 | return pass; |
| | 559 | } |
| | 560 | |
| | 561 | //_______________________________________________________________________ |
| | 562 | Bool_t AliFlowTrackCuts::PassesESDcuts(AliESDtrack* track) |
| | 563 | { |
| | 564 | Bool_t pass=kTRUE; |
| | 565 | if (fIgnoreTPCzRange) |
| | 566 | { |
| | 567 | const AliExternalTrackParam* pin = track->GetOuterParam(); |
| | 568 | const AliExternalTrackParam* pout = track->GetInnerParam(); |
| | 569 | if (pin&&pout) |
| | 570 | { |
| | 571 | Double_t zin = pin->GetZ(); |
| | 572 | Double_t zout = pout->GetZ(); |
| | 573 | if (zin*zout<0) pass=kFALSE; //reject if cross the membrane |
| | 574 | if (zin < fIgnoreTPCzRangeMin || zin > fIgnoreTPCzRangeMax) pass=kFALSE; |
| | 575 | if (zout < fIgnoreTPCzRangeMin || zout > fIgnoreTPCzRangeMax) pass=kFALSE; |
| | 576 | } |
| | 577 | } |
| | 578 | |
| | 579 | Int_t ntpccls = ( fParamType==kESD_TPConly )? |
| | 580 | track->GetTPCNclsIter1():track->GetTPCNcls(); |
| | 581 | if (fCutChi2PerClusterTPC) |
| | 582 | { |
| | 583 | Float_t tpcchi2 = (fParamType==kESD_TPConly)? |
| | 584 | track->GetTPCchi2Iter1():track->GetTPCchi2(); |
| | 585 | tpcchi2 = (ntpccls>0)?tpcchi2/ntpccls:-FLT_MAX; |
| | 586 | if (tpcchi2<fMinChi2PerClusterTPC || tpcchi2 >=fMaxChi2PerClusterTPC) |
| | 587 | pass=kFALSE; |
| | 588 | } |
| | 589 | |
| | 590 | if (fCutNClustersTPC) |
| | 591 | { |
| | 592 | if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE; |
| | 593 | } |
| | 594 | |
| | 595 | Int_t nitscls = track->GetNcls(0); |
| | 596 | if (fCutNClustersITS) |
| | 597 | { |
| | 598 | if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE; |
| | 599 | } |
| | 600 | |
| | 601 | Double_t pt = fTrack->Pt(); |
| | 602 | if (fCutPID) |
| | 603 | { |
| | 604 | switch (fPIDsource) |
| | 605 | { |
| | 606 | case kTPCpid: |
| | 607 | if (!PassesTPCpidCut(track)) pass=kFALSE; |
| | 608 | break; |
| | 609 | case kTOFpid: |
| | 610 | if (!PassesTOFpidCut(track)) pass=kFALSE; |
| | 611 | break; |
| | 612 | case kTPCTOFpid: |
| | 613 | if (pt< fTPCTOFpidCrossOverPt) |
| | 614 | { |
| | 615 | if (!PassesTPCpidCut(track)) pass=kFALSE; |
| | 616 | } |
| | 617 | else //if (pt>=fTPCTOFpidCrossOverPt) |
| | 618 | { |
| | 619 | if (!PassesTOFpidCut(track)) pass=kFALSE; |
| | 620 | } |
| | 621 | break; |
| | 622 | // part added by F. Noferini |
| | 623 | case kTOFbayesian: |
| | 624 | if (!PassesTOFbayesianCut(track)) pass=kFALSE; |
| | 625 | break; |
| | 626 | // end part added by F. Noferini |
| | 627 | default: |
| | 628 | printf("AliFlowTrackCuts::PassesCuts() this should never be called!\n"); |
| | 629 | pass=kFALSE; |
| | 630 | break; |
| | 631 | } |
| | 632 | } |
| | 633 | |
| | 634 | //some stuff is still handled by AliESDtrackCuts class - delegate |
| | 635 | if (fAliESDtrackCuts) |
| | 636 | { |
| | 637 | if (!fAliESDtrackCuts->IsSelected(track)) pass=kFALSE; |
| | 638 | } |
| | 639 | |
| | 640 | return pass; |
| | 641 | } |
| | 642 | |
| | 643 | //----------------------------------------------------------------------- |
| | 644 | void AliFlowTrackCuts::HandleVParticle(AliVParticle* track) |
| | 645 | { |
| | 646 | //handle the general case |
| | 647 | switch (fParamType) |
| | 648 | { |
| | 649 | default: |
| | 650 | fTrack = track; |
| | 651 | break; |
| | 652 | } |
| | 653 | } |
| | 654 | |
| | 655 | //----------------------------------------------------------------------- |
| | 656 | void AliFlowTrackCuts::HandleESDtrack(AliESDtrack* track) |
| | 657 | { |
| | 658 | //handle esd track |
| | 659 | AliExternalTrackParam* ip=NULL; |
| | 660 | switch (fParamType) |
| | 661 | { |
| | 662 | case kGlobal: |
| | 663 | fTrack = track; |
| | 664 | break; |
| | 665 | case kESD_TPConly: |
| | 666 | if (!track->FillTPCOnlyTrack(fTPCtrack)) |
| | 667 | { |
| | 668 | fTrack=NULL; |
| | 669 | fMCparticle=NULL; |
| | 670 | fTrackLabel=-1; |
| | 671 | return; |
| | 672 | } |
| | 673 | ip = const_cast<AliExternalTrackParam*>(fTPCtrack.GetInnerParam()); |
| | 674 | if (!ip) { ip = new AliExternalTrackParam(*(track->GetInnerParam())); } |
| | 675 | else { *ip = *(track->GetInnerParam()); } |
| | 676 | fTrack = &fTPCtrack; |
| | 677 | //recalculate the label and mc particle, they may differ as TPClabel != global label |
| | 678 | fTrackLabel = (fFakesAreOK)?TMath::Abs(fTrack->GetLabel()):fTrack->GetLabel(); |
| | 679 | if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel)); |
| | 680 | else fMCparticle=NULL; |
| | 681 | break; |
| | 682 | default: |
| | 683 | fTrack = track; |
| | 684 | break; |
| | 685 | } |
| | 686 | } |
| | 687 | |
| | 688 | //----------------------------------------------------------------------- |
| | 689 | AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCOnlyTrackCuts() |
| | 690 | { |
| | 691 | //get standard cuts |
| | 692 | AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPConly cuts"); |
| | 693 | delete cuts->fAliESDtrackCuts; |
| | 694 | cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts(); |
| | 695 | cuts->SetParamType(kESD_TPConly); |
| | 696 | return cuts; |
| | 697 | } |
| | 698 | |
| | 699 | //----------------------------------------------------------------------- |
| | 700 | AliFlowTrackCuts* AliFlowTrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries) |
| | 701 | { |
| | 702 | //get standard cuts |
| | 703 | AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard global track cuts 2009"); |
| | 704 | delete cuts->fAliESDtrackCuts; |
| | 705 | cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selPrimaries); |
| | 706 | cuts->SetParamType(kGlobal); |
| | 707 | return cuts; |
| | 708 | } |
| | 709 | |
| | 710 | //----------------------------------------------------------------------- |
| | 711 | AliFlowTrack* AliFlowTrackCuts::MakeFlowTrack() const |
| | 712 | { |
| | 713 | //get a flow track constructed from whatever we applied cuts on |
| | 714 | //caller is resposible for deletion |
| | 715 | //if construction fails return NULL |
| | 716 | AliFlowTrack* flowtrack=NULL; |
| | 717 | TParticle *tmpTParticle=NULL; |
| | 718 | AliMCParticle* tmpAliMCParticle=NULL; |
| | 719 | if (fParamType==kESD_SPDtracklet) |
| | 720 | { |
| | 721 | switch (fParamMix) |
| | 722 | { |
| | 723 | case kPure: |
| | 724 | flowtrack = new AliFlowTrack(); |
| | 725 | flowtrack->SetPhi(fTrackPhi); |
| | 726 | flowtrack->SetEta(fTrackEta); |
| | 727 | break; |
| | 728 | case kTrackWithMCkine: |
| | 729 | if (!fMCparticle) return NULL; |
| | 730 | flowtrack = new AliFlowTrack(); |
| | 731 | flowtrack->SetPhi( fMCparticle->Phi() ); |
| | 732 | flowtrack->SetEta( fMCparticle->Eta() ); |
| | 733 | flowtrack->SetPt( fMCparticle->Pt() ); |
| | 734 | break; |
| | 735 | case kTrackWithMCpt: |
| | 736 | if (!fMCparticle) return NULL; |
| | 737 | flowtrack = new AliFlowTrack(); |
| | 738 | flowtrack->SetPhi(fTrackPhi); |
| | 739 | flowtrack->SetEta(fTrackEta); |
| | 740 | flowtrack->SetPt(fMCparticle->Pt()); |
| | 741 | break; |
| | 742 | case kTrackWithPtFromFirstMother: |
| | 743 | if (!fMCparticle) return NULL; |
| | 744 | flowtrack = new AliFlowTrack(); |
| | 745 | flowtrack->SetPhi(fTrackPhi); |
| | 746 | flowtrack->SetEta(fTrackEta); |
| | 747 | tmpTParticle = fMCparticle->Particle(); |
| | 748 | tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother())); |
| | 749 | flowtrack->SetPt(tmpAliMCParticle->Pt()); |
| | 750 | break; |
| | 751 | default: |
| | 752 | flowtrack = new AliFlowTrack(); |
| | 753 | flowtrack->SetPhi(fTrackPhi); |
| | 754 | flowtrack->SetEta(fTrackEta); |
| | 755 | break; |
| | 756 | } |
| | 757 | flowtrack->SetSource(AliFlowTrack::kFromTracklet); |
| | 758 | } |
| | 759 | else |
| | 760 | { |
| | 761 | if (!fTrack) return NULL; |
| | 762 | switch(fParamMix) |
| | 763 | { |
| | 764 | case kPure: |
| | 765 | flowtrack = new AliFlowTrack(fTrack); |
| | 766 | break; |
| | 767 | case kTrackWithMCkine: |
| | 768 | flowtrack = new AliFlowTrack(fMCparticle); |
| | 769 | break; |
| | 770 | case kTrackWithMCPID: |
| | 771 | flowtrack = new AliFlowTrack(fTrack); |
| | 772 | //flowtrack->setPID(...) from mc, when implemented |
| | 773 | break; |
| | 774 | case kTrackWithMCpt: |
| | 775 | if (!fMCparticle) return NULL; |
| | 776 | flowtrack = new AliFlowTrack(fTrack); |
| | 777 | flowtrack->SetPt(fMCparticle->Pt()); |
| | 778 | break; |
| | 779 | case kTrackWithPtFromFirstMother: |
| | 780 | if (!fMCparticle) return NULL; |
| | 781 | flowtrack = new AliFlowTrack(fTrack); |
| | 782 | tmpTParticle = fMCparticle->Particle(); |
| | 783 | tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother())); |
| | 784 | flowtrack->SetPt(tmpAliMCParticle->Pt()); |
| | 785 | break; |
| | 786 | default: |
| | 787 | flowtrack = new AliFlowTrack(fTrack); |
| | 788 | break; |
| | 789 | } |
| | 790 | if (fParamType==kMC) flowtrack->SetSource(AliFlowTrack::kFromMC); |
| | 791 | else if (dynamic_cast<AliESDtrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromESD); |
| | 792 | else if (dynamic_cast<AliAODTrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromAOD); |
| | 793 | else if (dynamic_cast<AliMCParticle*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromMC); |
| | 794 | } |
| | 795 | |
| | 796 | return flowtrack; |
| | 797 | } |
| | 798 | |
| | 799 | //----------------------------------------------------------------------- |
| | 800 | Bool_t AliFlowTrackCuts::IsPhysicalPrimary() const |
| | 801 | { |
| | 802 | //check if current particle is a physical primary |
| | 803 | if (!fMCevent) return kFALSE; |
| | 804 | if (fTrackLabel<0) return kFALSE; |
| | 805 | return IsPhysicalPrimary(fMCevent, fTrackLabel, fRequireTransportBitForPrimaries); |
| | 806 | } |
| | 807 | |
| | 808 | //----------------------------------------------------------------------- |
| | 809 | Bool_t AliFlowTrackCuts::IsPhysicalPrimary(AliMCEvent* mcEvent, Int_t label, Bool_t requiretransported) |
| | 810 | { |
| | 811 | //check if current particle is a physical primary |
| | 812 | Bool_t physprim=mcEvent->IsPhysicalPrimary(label); |
| | 813 | AliMCParticle* track = static_cast<AliMCParticle*>(mcEvent->GetTrack(label)); |
| | 814 | if (!track) return kFALSE; |
| | 815 | TParticle* particle = track->Particle(); |
| | 816 | Bool_t transported = particle->TestBit(kTransportBit); |
| | 817 | //printf("label: %i prim: %s, transp: %s, pass: %s\n",label, (physprim)?"YES":"NO ",(transported)?"YES":"NO ", |
| | 818 | //(physprim && (transported || !requiretransported))?"YES":"NO" ); |
| | 819 | return (physprim && (transported || !requiretransported)); |
| | 820 | } |
| | 821 | |
| | 822 | //----------------------------------------------------------------------- |
| | 823 | const char* AliFlowTrackCuts::GetParamTypeName(trackParameterType type) |
| | 824 | { |
| | 825 | //return the name of the selected parameter type |
| | 826 | switch (type) |
| | 827 | { |
| | 828 | case kMC: |
| | 829 | return "MC"; |
| | 830 | case kGlobal: |
| | 831 | return "ESD global"; |
| | 832 | case kESD_TPConly: |
| | 833 | return "TPC only"; |
| | 834 | case kESD_SPDtracklet: |
| | 835 | return "SPD tracklet"; |
| | 836 | default: |
| | 837 | return "unknown"; |
| | 838 | } |
| | 839 | } |
| | 840 | |
| | 841 | //----------------------------------------------------------------------- |
| | 842 | void AliFlowTrackCuts::DefineHistograms() |
| | 843 | { |
| | 844 | //define qa histograms |
| | 845 | } |
| | 846 | |
| | 847 | //----------------------------------------------------------------------- |
| | 848 | Int_t AliFlowTrackCuts::GetNumberOfInputObjects() const |
| | 849 | { |
| | 850 | //get the number of tracks in the input event according source |
| | 851 | //selection (ESD tracks, tracklets, MC particles etc.) |
| | 852 | AliESDEvent* esd=NULL; |
| | 853 | switch (fParamType) |
| | 854 | { |
| | 855 | case kESD_SPDtracklet: |
| | 856 | esd = dynamic_cast<AliESDEvent*>(fEvent); |
| | 857 | if (!esd) return 0; |
| | 858 | return esd->GetMultiplicity()->GetNumberOfTracklets(); |
| | 859 | case kMC: |
| | 860 | if (!fMCevent) return 0; |
| | 861 | return fMCevent->GetNumberOfTracks(); |
| | 862 | default: |
| | 863 | if (!fEvent) return 0; |
| | 864 | return fEvent->GetNumberOfTracks(); |
| | 865 | } |
| | 866 | return 0; |
| | 867 | } |
| | 868 | |
| | 869 | //----------------------------------------------------------------------- |
| | 870 | TObject* AliFlowTrackCuts::GetInputObject(Int_t i) |
| | 871 | { |
| | 872 | //get the input object according the data source selection: |
| | 873 | //(esd tracks, traclets, mc particles,etc...) |
| | 874 | AliESDEvent* esd=NULL; |
| | 875 | switch (fParamType) |
| | 876 | { |
| | 877 | case kESD_SPDtracklet: |
| | 878 | esd = dynamic_cast<AliESDEvent*>(fEvent); |
| | 879 | if (!esd) return NULL; |
| | 880 | return const_cast<AliMultiplicity*>(esd->GetMultiplicity()); |
| | 881 | case kMC: |
| | 882 | if (!fMCevent) return NULL; |
| | 883 | return fMCevent->GetTrack(i); |
| | 884 | default: |
| | 885 | if (!fEvent) return NULL; |
| | 886 | return fEvent->GetTrack(i); |
| | 887 | } |
| | 888 | } |
| | 889 | |
| | 890 | //----------------------------------------------------------------------- |
| | 891 | void AliFlowTrackCuts::Clear(Option_t*) |
| | 892 | { |
| | 893 | //clean up |
| | 894 | fTrack=NULL; |
| | 895 | fMCevent=NULL; |
| | 896 | fMCparticle=NULL; |
| | 897 | fTrackLabel=0; |
| | 898 | fTrackWeight=0.0; |
| | 899 | fTrackEta=0.0; |
| | 900 | fTrackPhi=0.0; |
| | 901 | } |
| | 902 | |
| | 903 | //----------------------------------------------------------------------- |
| | 904 | Bool_t AliFlowTrackCuts::PassesTOFpidCut(AliESDtrack* track ) |
| | 905 | { |
| | 906 | //check if passes PID cut using timing in TOF |
| | 907 | Bool_t goodtrack = (track) && |
| | 908 | (track->GetStatus() & AliESDtrack::kTOFpid) && |
| | 909 | (track->GetTOFsignal() > 12000) && |
| | 910 | (track->GetTOFsignal() < 100000) && |
| | 911 | (track->GetIntegratedLength() > 365) && |
| | 912 | !(track->GetStatus() & AliESDtrack::kTOFmismatch); |
| | 913 | |
| | 914 | if (!goodtrack) return kFALSE; |
| | 915 | |
| | 916 | Float_t pt = track->Pt(); |
| | 917 | Float_t p = track->GetP(); |
| | 918 | Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p); |
| | 919 | Float_t timeTOF = track->GetTOFsignal()- trackT0; |
| | 920 | //2=pion 3=kaon 4=protons |
| | 921 | Double_t inttimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0}; |
| | 922 | track->GetIntegratedTimes(inttimes); |
| | 923 | //construct the pid index because it's screwed up in TOF |
| | 924 | Int_t pid = 0; |
| | 925 | switch (fAliPID) |
| | 926 | { |
| | 927 | case AliPID::kPion: |
| | 928 | pid=2; |
| | 929 | break; |
| | 930 | case AliPID::kKaon: |
| | 931 | pid=3; |
| | 932 | break; |
| | 933 | case AliPID::kProton: |
| | 934 | pid=4; |
| | 935 | break; |
| | 936 | default: |
| | 937 | return kFALSE; |
| | 938 | } |
| | 939 | Float_t s = timeTOF-inttimes[pid]; |
| | 940 | |
| | 941 | Float_t* arr = fTOFpidCuts->GetMatrixArray(); |
| | 942 | Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(pt)); |
| | 943 | if (col<0) return kFALSE; |
| | 944 | Float_t min = (*fTOFpidCuts)(1,col); |
| | 945 | Float_t max = (*fTOFpidCuts)(2,col); |
| | 946 | |
| | 947 | //printf("--------------TOF pid cut %s\n",(s>min && s<max)?"PASS":"FAIL"); |
| | 948 | return (s>min && s<max); |
| | 949 | } |
| | 950 | |
| | 951 | //----------------------------------------------------------------------- |
| | 952 | Bool_t AliFlowTrackCuts::PassesTPCpidCut(AliESDtrack* track) |
| | 953 | { |
| | 954 | //check if passes PID cut using dedx signal in the TPC |
| | 955 | if (!fTPCpidCuts) |
| | 956 | { |
| | 957 | printf("no TPCpidCuts\n"); |
| | 958 | return kFALSE; |
| | 959 | } |
| | 960 | |
| | 961 | const AliExternalTrackParam* tpcparam = track->GetInnerParam(); //tpc only params at the inner wall |
| | 962 | if (!tpcparam) return kFALSE; |
| | 963 | Float_t sigExp = fESDpid.GetTPCResponse().GetExpectedSignal(tpcparam->GetP(), fAliPID); |
| | 964 | Float_t sigTPC = track->GetTPCsignal(); |
| | 965 | Float_t s = (sigTPC-sigExp)/sigExp; |
| | 966 | Double_t pt = track->Pt(); |
| | 967 | |
| | 968 | Float_t* arr = fTPCpidCuts->GetMatrixArray(); |
| | 969 | Int_t col = TMath::BinarySearch(fTPCpidCuts->GetNcols(),arr,static_cast<Float_t>(pt)); |
| | 970 | if (col<0) return kFALSE; |
| | 971 | Float_t min = (*fTPCpidCuts)(1,col); |
| | 972 | Float_t max = (*fTPCpidCuts)(2,col); |
| | 973 | |
| | 974 | //printf("------------TPC pid cut %s\n",(s>min && s<max)?"PASS":"FAIL"); |
| | 975 | return (s>min && s<max); |
| | 976 | } |
| | 977 | |
| | 978 | //----------------------------------------------------------------------- |
| | 979 | void AliFlowTrackCuts::InitPIDcuts() |
| | 980 | { |
| | 981 | //init matrices with PID cuts |
| | 982 | TMatrixF* t = NULL; |
| | 983 | if (!fTPCpidCuts) |
| | 984 | { |
| | 985 | if (fAliPID==AliPID::kPion) |
| | 986 | { |
| | 987 | t = new TMatrixF(3,10); |
| | 988 | (*t)(0,0) = 0.20; (*t)(1,0) = -0.4; (*t)(2,0) = 0.2; |
| | 989 | (*t)(0,1) = 0.25; (*t)(1,1) = -0.4; (*t)(2,1) = 0.2; |
| | 990 | (*t)(0,2) = 0.30; (*t)(1,2) = -0.4; (*t)(2,2) = 0.25; |
| | 991 | (*t)(0,3) = 0.35; (*t)(1,3) = -0.4; (*t)(2,3) = 0.25; |
| | 992 | (*t)(0,4) = 0.40; (*t)(1,4) = -0.4; (*t)(2,4) = 0.3; |
| | 993 | (*t)(0,5) = 0.45; (*t)(1,5) = -0.4; (*t)(2,5) = 0.3; |
| | 994 | (*t)(0,6) = 0.50; (*t)(1,6) = -0.4; (*t)(2,6) = 0.25; |
| | 995 | (*t)(0,7) = 0.55; (*t)(1,7) = -0.4; (*t)(2,7) = 0.15; |
| | 996 | (*t)(0,8) = 0.60; (*t)(1,8) = -0.4; (*t)(2,8) = 0.1; |
| | 997 | (*t)(0,9) = 0.65; (*t)(1,9) = 0; (*t)(2,9) = 0; |
| | 998 | } |
| | 999 | else |
| | 1000 | if (fAliPID==AliPID::kKaon) |
| | 1001 | { |
| | 1002 | t = new TMatrixF(3,7); |
| | 1003 | (*t)(0,0) = 0.20; (*t)(1,0) = -0.2; (*t)(2,0) = 0.4; |
| | 1004 | (*t)(0,1) = 0.25; (*t)(1,1) =-0.15; (*t)(2,1) = 0.4; |
| | 1005 | (*t)(0,2) = 0.30; (*t)(1,2) = -0.1; (*t)(2,2) = 0.4; |
| | 1006 | (*t)(0,3) = 0.35; (*t)(1,3) = -0.1; (*t)(2,3) = 0.4; |
| | 1007 | (*t)(0,4) = 0.40; (*t)(1,4) = -0.1; (*t)(2,4) = 0.6; |
| | 1008 | (*t)(0,5) = 0.45; (*t)(1,5) = -0.1; (*t)(2,5) = 0.6; |
| | 1009 | (*t)(0,6) = 0.50; (*t)(1,6) = 0; (*t)(2,6) = 0; |
| | 1010 | } |
| | 1011 | else |
| | 1012 | if (fAliPID==AliPID::kProton) |
| | 1013 | { |
| | 1014 | t = new TMatrixF(3,16); |
| | 1015 | (*t)(0,0) = 0.20; (*t)(1,0) = 0; (*t)(2,0) = 0; |
| | 1016 | (*t)(0,1) = 0.25; (*t)(1,1) = -0.2; (*t)(2,1) = 0.3; |
| | 1017 | (*t)(0,2) = 0.30; (*t)(1,2) = -0.2; (*t)(2,2) = 0.6; |
| | 1018 | (*t)(0,3) = 0.35; (*t)(1,3) = -0.2; (*t)(2,3) = 0.6; |
| | 1019 | (*t)(0,4) = 0.40; (*t)(1,4) = -0.2; (*t)(2,4) = 0.6; |
| | 1020 | (*t)(0,5) = 0.45; (*t)(1,5) = -0.15; (*t)(2,5) = 0.6; |
| | 1021 | (*t)(0,6) = 0.50; (*t)(1,6) = -0.1; (*t)(2,6) = 0.6; |
| | 1022 | (*t)(0,7) = 0.55; (*t)(1,7) = -0.05; (*t)(2,7) = 0.6; |
| | 1023 | (*t)(0,8) = 0.60; (*t)(1,8) = -0.05; (*t)(2,8) = 0.45; |
| | 1024 | (*t)(0,9) = 0.65; (*t)(1,9) = -0.05; (*t)(2,9) = 0.45; |
| | 1025 | (*t)(0,10) = 0.70; (*t)(1,10) = -0.05; (*t)(2,10) = 0.45; |
| | 1026 | (*t)(0,11) = 0.75; (*t)(1,11) = -0.05; (*t)(2,11) = 0.45; |
| | 1027 | (*t)(0,12) = 0.80; (*t)(1,12) = 0; (*t)(2,12) = 0.45; |
| | 1028 | (*t)(0,13) = 0.85; (*t)(1,13) = 0; (*t)(2,13) = 0.45; |
| | 1029 | (*t)(0,14) = 0.90; (*t)(1,14) = 0; (*t)(2,14) = 0.45; |
| | 1030 | (*t)(0,15) = 0.95; (*t)(1,15) = 0; (*t)(2,15) = 0; |
| | 1031 | } |
| | 1032 | fTPCpidCuts=t; |
| | 1033 | } |
| | 1034 | t = NULL; |
| | 1035 | if (!fTOFpidCuts) |
| | 1036 | { |
| | 1037 | if (fAliPID==AliPID::kPion) |
| | 1038 | { |
| | 1039 | t = new TMatrixF(3,31); |
| | 1040 | (*t)(0,0) = 0.3; (*t)(1,0) = -700; (*t)(2,0) = 700; |
| | 1041 | (*t)(0,1) = 0.35; (*t)(1,1) = -800; (*t)(2,1) = 800; |
| | 1042 | (*t)(0,2) = 0.40; (*t)(1,2) = -600; (*t)(2,2) = 800; |
| | 1043 | (*t)(0,3) = 0.45; (*t)(1,3) = -500; (*t)(2,3) = 700; |
| | 1044 | (*t)(0,4) = 0.50; (*t)(1,4) = -400; (*t)(2,4) = 700; |
| | 1045 | (*t)(0,5) = 0.55; (*t)(1,5) = -400; (*t)(2,5) = 700; |
| | 1046 | (*t)(0,6) = 0.60; (*t)(1,6) = -400; (*t)(2,6) = 700; |
| | 1047 | (*t)(0,7) = 0.65; (*t)(1,7) = -400; (*t)(2,7) = 700; |
| | 1048 | (*t)(0,8) = 0.70; (*t)(1,8) = -400; (*t)(2,8) = 700; |
| | 1049 | (*t)(0,9) = 0.75; (*t)(1,9) = -400; (*t)(2,9) = 700; |
| | 1050 | (*t)(0,10) = 0.80; (*t)(1,10) = -400; (*t)(2,10) = 600; |
| | 1051 | (*t)(0,11) = 0.85; (*t)(1,11) = -400; (*t)(2,11) = 600; |
| | 1052 | (*t)(0,12) = 0.90; (*t)(1,12) = -400; (*t)(2,12) = 600; |
| | 1053 | (*t)(0,13) = 0.95; (*t)(1,13) = -400; (*t)(2,13) = 600; |
| | 1054 | (*t)(0,14) = 1.00; (*t)(1,14) = -400; (*t)(2,14) = 550; |
| | 1055 | (*t)(0,15) = 1.10; (*t)(1,15) = -400; (*t)(2,15) = 450; |
| | 1056 | (*t)(0,16) = 1.20; (*t)(1,16) = -400; (*t)(2,16) = 400; |
| | 1057 | (*t)(0,17) = 1.30; (*t)(1,17) = -400; (*t)(2,17) = 300; |
| | 1058 | (*t)(0,18) = 1.40; (*t)(1,18) = -400; (*t)(2,18) = 300; |
| | 1059 | (*t)(0,19) = 1.50; (*t)(1,19) = -400; (*t)(2,19) = 250; |
| | 1060 | (*t)(0,20) = 1.60; (*t)(1,20) = -400; (*t)(2,20) = 200; |
| | 1061 | (*t)(0,21) = 1.70; (*t)(1,21) = -400; (*t)(2,21) = 150; |
| | 1062 | (*t)(0,22) = 1.80; (*t)(1,22) = -400; (*t)(2,22) = 100; |
| | 1063 | (*t)(0,23) = 1.90; (*t)(1,23) = -400; (*t)(2,23) = 70; |
| | 1064 | (*t)(0,24) = 2.00; (*t)(1,24) = -400; (*t)(2,24) = 50; |
| | 1065 | (*t)(0,25) = 2.10; (*t)(1,25) = -400; (*t)(2,25) = 0; |
| | 1066 | (*t)(0,26) = 2.20; (*t)(1,26) = -400; (*t)(2,26) = 0; |
| | 1067 | (*t)(0,27) = 2.30; (*t)(1,26) = -400; (*t)(2,26) = 0; |
| | 1068 | (*t)(0,28) = 2.40; (*t)(1,26) = -400; (*t)(2,26) = -50; |
| | 1069 | (*t)(0,29) = 2.50; (*t)(1,26) = -400; (*t)(2,26) = -50; |
| | 1070 | (*t)(0,30) = 2.60; (*t)(1,26) = 0; (*t)(2,26) = 0; |
| | 1071 | } |
| | 1072 | else |
| | 1073 | if (fAliPID==AliPID::kProton) |
| | 1074 | { |
| | 1075 | t = new TMatrixF(3,39); |
| | 1076 | (*t)(0,0) = 0.3; (*t)(1,0) = 0; (*t)(2,0) = 0; |
| | 1077 | (*t)(0,1) = 0.35; (*t)(1,1) = 0; (*t)(2,1) = 0; |
| | 1078 | (*t)(0,2) = 0.40; (*t)(1,2) = 0; (*t)(2,2) = 0; |
| | 1079 | (*t)(0,3) = 0.45; (*t)(1,3) = 0; (*t)(2,3) = 0; |
| | 1080 | (*t)(0,4) = 0.50; (*t)(1,4) = 0; (*t)(2,4) = 0; |
| | 1081 | (*t)(0,5) = 0.55; (*t)(1,5) = -900; (*t)(2,5) = 600; |
| | 1082 | (*t)(0,6) = 0.60; (*t)(1,6) = -800; (*t)(2,6) = 600; |
| | 1083 | (*t)(0,7) = 0.65; (*t)(1,7) = -800; (*t)(2,7) = 600; |
| | 1084 | (*t)(0,8) = 0.70; (*t)(1,8) = -800; (*t)(2,8) = 600; |
| | 1085 | (*t)(0,9) = 0.75; (*t)(1,9) = -700; (*t)(2,9) = 500; |
| | 1086 | (*t)(0,10) = 0.80; (*t)(1,10) = -700; (*t)(2,10) = 500; |
| | 1087 | (*t)(0,11) = 0.85; (*t)(1,11) = -700; (*t)(2,11) = 500; |
| | 1088 | (*t)(0,12) = 0.90; (*t)(1,12) = -600; (*t)(2,12) = 500; |
| | 1089 | (*t)(0,13) = 0.95; (*t)(1,13) = -600; (*t)(2,13) = 500; |
| | 1090 | (*t)(0,14) = 1.00; (*t)(1,14) = -600; (*t)(2,14) = 500; |
| | 1091 | (*t)(0,15) = 1.10; (*t)(1,15) = -600; (*t)(2,15) = 500; |
| | 1092 | (*t)(0,16) = 1.20; (*t)(1,16) = -500; (*t)(2,16) = 500; |
| | 1093 | (*t)(0,17) = 1.30; (*t)(1,17) = -500; (*t)(2,17) = 500; |
| | 1094 | (*t)(0,18) = 1.40; (*t)(1,18) = -500; (*t)(2,18) = 500; |
| | 1095 | (*t)(0,19) = 1.50; (*t)(1,19) = -500; (*t)(2,19) = 500; |
| | 1096 | (*t)(0,20) = 1.60; (*t)(1,20) = -400; (*t)(2,20) = 500; |
| | 1097 | (*t)(0,21) = 1.70; (*t)(1,21) = -400; (*t)(2,21) = 500; |
| | 1098 | (*t)(0,22) = 1.80; (*t)(1,22) = -400; (*t)(2,22) = 500; |
| | 1099 | (*t)(0,23) = 1.90; (*t)(1,23) = -400; (*t)(2,23) = 500; |
| | 1100 | (*t)(0,24) = 2.00; (*t)(1,24) = -400; (*t)(2,24) = 500; |
| | 1101 | (*t)(0,25) = 2.10; (*t)(1,25) = -350; (*t)(2,25) = 500; |
| | 1102 | (*t)(0,26) = 2.20; (*t)(1,26) = -350; (*t)(2,26) = 500; |
| | 1103 | (*t)(0,27) = 2.30; (*t)(1,27) = -300; (*t)(2,27) = 500; |
| | 1104 | (*t)(0,28) = 2.40; (*t)(1,28) = -300; (*t)(2,28) = 500; |
| | 1105 | (*t)(0,29) = 2.50; (*t)(1,29) = -300; (*t)(2,29) = 500; |
| | 1106 | (*t)(0,30) = 2.60; (*t)(1,30) = -250; (*t)(2,30) = 500; |
| | 1107 | (*t)(0,31) = 2.70; (*t)(1,31) = -200; (*t)(2,31) = 500; |
| | 1108 | (*t)(0,32) = 2.80; (*t)(1,32) = -150; (*t)(2,32) = 500; |
| | 1109 | (*t)(0,33) = 2.90; (*t)(1,33) = -150; (*t)(2,33) = 500; |
| | 1110 | (*t)(0,34) = 3.00; (*t)(1,34) = -100; (*t)(2,34) = 400; |
| | 1111 | (*t)(0,35) = 3.10; (*t)(1,35) = -100; (*t)(2,35) = 400; |
| | 1112 | (*t)(0,36) = 3.50; (*t)(1,36) = -100; (*t)(2,36) = 400; |
| | 1113 | (*t)(0,37) = 3.60; (*t)(1,37) = 0; (*t)(2,37) = 0; |
| | 1114 | (*t)(0,38) = 3.70; (*t)(1,38) = 0; (*t)(2,38) = 0; |
| | 1115 | } |
| | 1116 | else |
| | 1117 | if (fAliPID==AliPID::kKaon) |
| | 1118 | { |
| | 1119 | t = new TMatrixF(3,26); |
| | 1120 | (*t)(0,0) = 0.3; (*t)(1,0) = 0; (*t)(2,0) = 0; |
| | 1121 | (*t)(0,1) = 0.35; (*t)(1,1) = 0; (*t)(2,1) = 0; |
| | 1122 | (*t)(0,2) = 0.40; (*t)(1,2) = -800; (*t)(2,2) = 600; |
| | 1123 | (*t)(0,3) = 0.45; (*t)(1,3) = -800; (*t)(2,3) = 600; |
| | 1124 | (*t)(0,4) = 0.50; (*t)(1,4) = -800; (*t)(2,4) = 600; |
| | 1125 | (*t)(0,5) = 0.55; (*t)(1,5) = -800; (*t)(2,5) = 600; |
| | 1126 | (*t)(0,6) = 0.60; (*t)(1,6) = -800; (*t)(2,6) = 600; |
| | 1127 | (*t)(0,7) = 0.65; (*t)(1,7) = -700; (*t)(2,7) = 600; |
| | 1128 | (*t)(0,8) = 0.70; (*t)(1,8) = -600; (*t)(2,8) = 600; |
| | 1129 | (*t)(0,9) = 0.75; (*t)(1,9) = -600; (*t)(2,9) = 500; |
| | 1130 | (*t)(0,10) = 0.80; (*t)(1,10) = -500; (*t)(2,10) = 500; |
| | 1131 | (*t)(0,11) = 0.85; (*t)(1,11) = -500; (*t)(2,11) = 500; |
| | 1132 | (*t)(0,12) = 0.90; (*t)(1,12) = -400; (*t)(2,12) = 500; |
| | 1133 | (*t)(0,13) = 0.95; (*t)(1,13) = -400; (*t)(2,13) = 500; |
| | 1134 | (*t)(0,14) = 1.00; (*t)(1,14) = -400; (*t)(2,14) = 500; |
| | 1135 | (*t)(0,15) = 1.10; (*t)(1,15) = -350; (*t)(2,15) = 450; |
| | 1136 | (*t)(0,16) = 1.20; (*t)(1,16) = -300; (*t)(2,16) = 400; |
| | 1137 | (*t)(0,17) = 1.30; (*t)(1,17) = -300; (*t)(2,17) = 400; |
| | 1138 | (*t)(0,18) = 1.40; (*t)(1,18) = -250; (*t)(2,18) = 400; |
| | 1139 | (*t)(0,19) = 1.50; (*t)(1,19) = -200; (*t)(2,19) = 400; |
| | 1140 | (*t)(0,20) = 1.60; (*t)(1,20) = -150; (*t)(2,20) = 400; |
| | 1141 | (*t)(0,21) = 1.70; (*t)(1,21) = -100; (*t)(2,21) = 400; |
| | 1142 | (*t)(0,22) = 1.80; (*t)(1,22) = -50; (*t)(2,22) = 400; |
| | 1143 | (*t)(0,23) = 1.90; (*t)(1,22) = 0; (*t)(2,22) = 400; |
| | 1144 | (*t)(0,24) = 2.00; (*t)(1,22) = 50; (*t)(2,22) = 400; |
| | 1145 | (*t)(0,25) = 2.10; (*t)(1,22) = 0; (*t)(2,22) = 0; |
| | 1146 | } |
| | 1147 | fTOFpidCuts=t; |
| | 1148 | } |
| | 1149 | } |
| | 1150 | |
| | 1151 | //----------------------------------------------------------------------- |
| | 1152 | // part added by F. Noferini (some methods) |
| | 1153 | Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(AliESDtrack* track){ |
| | 1154 | |
| | 1155 | Bool_t goodtrack = track && (track->GetStatus() & AliESDtrack::kTOFpid) && (track->GetTOFsignal() > 12000) && (track->GetTOFsignal() < 100000) && (track->GetIntegratedLength() > 365) && !(track->GetStatus() & AliESDtrack::kTOFmismatch); |
| | 1156 | |
| | 1157 | if (! goodtrack) |
| | 1158 | return kFALSE; |
| | 1159 | |
| | 1160 | Int_t pdg = GetESDPdg(track,"bayesianALL"); |
| | 1161 | // printf("pdg set to %i\n",pdg); |
| | 1162 | |
| | 1163 | Int_t pid = 0; |
| | 1164 | Float_t prob = 0; |
| | 1165 | switch (fAliPID) |
| | 1166 | { |
| | 1167 | case AliPID::kPion: |
| | 1168 | pid=211; |
| | 1169 | prob = fProbBayes[2]; |
| | 1170 | break; |
| | 1171 | case AliPID::kKaon: |
| | 1172 | pid=321; |
| | 1173 | prob = fProbBayes[3]; |
| | 1174 | break; |
| | 1175 | case AliPID::kProton: |
| | 1176 | pid=2212; |
| | 1177 | prob = fProbBayes[4]; |
| | 1178 | break; |
| | 1179 | case AliPID::kElectron: |
| | 1180 | pid=-11; |
| | 1181 | prob = fProbBayes[0]; |
| | 1182 | break; |
| | 1183 | default: |
| | 1184 | return kFALSE; |
| | 1185 | } |
| | 1186 | |
| | 1187 | // 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); |
| | 1188 | if(TMath::Abs(pdg) == TMath::Abs(pid) && prob > 0.8){ |
| | 1189 | if(!fCutCharge) |
| | 1190 | return kTRUE; |
| | 1191 | else if (fCutCharge && fCharge * track->GetSign() > 0) |
| | 1192 | return kTRUE; |
| | 1193 | } |
| | 1194 | return kFALSE; |
| | 1195 | } |
| | 1196 | //----------------------------------------------------------------------- |
| | 1197 | Int_t AliFlowTrackCuts::GetESDPdg(AliESDtrack *track,Option_t *option,Int_t ipart,Float_t cPi,Float_t cKa,Float_t cPr){ |
| | 1198 | Int_t pdg = 0; |
| | 1199 | Int_t pdgvalues[5] = {-11,-13,211,321,2212}; |
| | 1200 | Float_t mass[5] = {5.10998909999999971e-04,1.05658000000000002e-01,1.39570000000000000e-01,4.93676999999999977e-01,9.38271999999999995e-01}; |
| | 1201 | |
| | 1202 | if(strstr(option,"bayesianTOF")){ // Bayesian TOF PID |
| | 1203 | Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05}; |
| | 1204 | Double_t rcc=0.; |
| | 1205 | |
| | 1206 | Float_t pt = track->Pt(); |
| | 1207 | |
| | 1208 | Int_t iptesd = 0; |
| | 1209 | while(pt > fBinLimitPID[iptesd] && iptesd < fnPIDptBin-1) iptesd++; |
| | 1210 | |
| | 1211 | if(cPi < 0){ |
| | 1212 | c[0] = fC[iptesd][0]; |
| | 1213 | c[1] = fC[iptesd][1]; |
| | 1214 | c[2] = fC[iptesd][2]; |
| | 1215 | c[3] = fC[iptesd][3]; |
| | 1216 | c[4] = fC[iptesd][4]; |
| | 1217 | } |
| | 1218 | else{ |
| | 1219 | c[0] = 0.0; |
| | 1220 | c[1] = 0.0; |
| | 1221 | c[2] = cPi; |
| | 1222 | c[3] = cKa; |
| | 1223 | c[4] = cPr; |
| | 1224 | } |
| | 1225 | |
| | 1226 | Double_t r1[10]; track->GetTOFpid(r1); |
| | 1227 | |
| | 1228 | Int_t i; |
| | 1229 | for (i=0; i<5; i++) rcc+=(c[i]*r1[i]); |
| | 1230 | |
| | 1231 | Double_t w[10]; |
| | 1232 | for (i=0; i<5; i++){ |
| | 1233 | w[i]=c[i]*r1[i]/rcc; |
| | 1234 | fProbBayes[i] = w[i]; |
| | 1235 | } |
| | 1236 | if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion |
| | 1237 | pdg = 211*Int_t(track->GetSign()); |
| | 1238 | } |
| | 1239 | else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton |
| | 1240 | pdg = 2212*Int_t(track->GetSign()); |
| | 1241 | } |
| | 1242 | else if (w[3]>=w[1] && w[3]>=w[0]){//kaon |
| | 1243 | pdg = 321*Int_t(track->GetSign()); |
| | 1244 | } |
| | 1245 | else if (w[0]>=w[1]) { //electrons |
| | 1246 | pdg = -11*Int_t(track->GetSign()); |
| | 1247 | } |
| | 1248 | else{ // muon |
| | 1249 | pdg = -13*Int_t(track->GetSign()); |
| | 1250 | } |
| | 1251 | } |
| | 1252 | |
| | 1253 | else if(strstr(option,"bayesianTPC")){ // Bayesian TPC PID |
| | 1254 | Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05}; |
| | 1255 | Double_t rcc=0.; |
| | 1256 | |
| | 1257 | Float_t pt = track->Pt(); |
| | 1258 | |
| | 1259 | Int_t iptesd = 0; |
| | 1260 | while(pt > fBinLimitPID[iptesd] && iptesd < fnPIDptBin-1) iptesd++; |
| | 1261 | |
| | 1262 | if(cPi < 0){ |
| | 1263 | c[0] = fC[iptesd][0]; |
| | 1264 | c[1] = fC[iptesd][1]; |
| | 1265 | c[2] = fC[iptesd][2]; |
| | 1266 | c[3] = fC[iptesd][3]; |
| | 1267 | c[4] = fC[iptesd][4]; |
| | 1268 | } |
| | 1269 | else{ |
| | 1270 | c[0] = 0.0; |
| | 1271 | c[1] = 0.0; |
| | 1272 | c[2] = cPi; |
| | 1273 | c[3] = cKa; |
| | 1274 | c[4] = cPr; |
| | 1275 | } |
| | 1276 | |
| | 1277 | Double_t r1[10]; track->GetTPCpid(r1); |
| | 1278 | |
| | 1279 | Int_t i; |
| | 1280 | for (i=0; i<5; i++) rcc+=(c[i]*r1[i]); |
| | 1281 | |
| | 1282 | Double_t w[10]; |
| | 1283 | for (i=0; i<5; i++){ |
| | 1284 | w[i]=c[i]*r1[i]/rcc; |
| | 1285 | fProbBayes[i] = w[i]; |
| | 1286 | } |
| | 1287 | if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion |
| | 1288 | pdg = 211*Int_t(track->GetSign()); |
| | 1289 | } |
| | 1290 | else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton |
| | 1291 | pdg = 2212*Int_t(track->GetSign()); |
| | 1292 | } |
| | 1293 | else if (w[3]>=w[1] && w[3]>=w[0]){//kaon |
| | 1294 | pdg = 321*Int_t(track->GetSign()); |
| | 1295 | } |
| | 1296 | else if (w[0]>=w[1]) { //electrons |
| | 1297 | pdg = -11*Int_t(track->GetSign()); |
| | 1298 | } |
| | 1299 | else{ // muon |
| | 1300 | pdg = -13*Int_t(track->GetSign()); |
| | 1301 | } |
| | 1302 | } |
| | 1303 | |
| | 1304 | else if(strstr(option,"bayesianALL")){ |
| | 1305 | Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05}; |
| | 1306 | Double_t rcc=0.; |
| | 1307 | |
| | 1308 | Float_t pt = track->Pt(); |
| | 1309 | |
| | 1310 | Int_t iptesd = 0; |
| | 1311 | while(pt > fBinLimitPID[iptesd] && iptesd < fnPIDptBin-1) iptesd++; |
| | 1312 | |
| | 1313 | if(cPi < 0){ |
| | 1314 | c[0] = fC[iptesd][0]; |
| | 1315 | c[1] = fC[iptesd][1]; |
| | 1316 | c[2] = fC[iptesd][2]; |
| | 1317 | c[3] = fC[iptesd][3]; |
| | 1318 | c[4] = fC[iptesd][4]; |
| | 1319 | } |
| | 1320 | else{ |
| | 1321 | c[0] = 0.0; |
| | 1322 | c[1] = 0.0; |
| | 1323 | c[2] = cPi; |
| | 1324 | c[3] = cKa; |
| | 1325 | c[4] = cPr; |
| | 1326 | } |
| | 1327 | |
| | 1328 | Double_t r1[10]; track->GetTOFpid(r1); |
| | 1329 | Double_t r2[10]; track->GetTPCpid(r2); |
| | 1330 | |
| | 1331 | Int_t i; |
| | 1332 | for (i=0; i<5; i++) rcc+=(c[i]*r1[i]*r2[i]); |
| | 1333 | |
| | 1334 | |
| | 1335 | Double_t w[10]; |
| | 1336 | for (i=0; i<5; i++){ |
| | 1337 | w[i]=c[i]*r1[i]*r2[i]/rcc; |
| | 1338 | fProbBayes[i] = w[i]; |
| | 1339 | } |
| | 1340 | |
| | 1341 | if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion |
| | 1342 | pdg = 211*Int_t(track->GetSign()); |
| | 1343 | } |
| | 1344 | else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton |
| | 1345 | pdg = 2212*Int_t(track->GetSign()); |
| | 1346 | } |
| | 1347 | else if (w[3]>=w[1] && w[3]>=w[0]){//kaon |
| | 1348 | pdg = 321*Int_t(track->GetSign()); |
| | 1349 | } |
| | 1350 | else if (w[0]>=w[1]) { //electrons |
| | 1351 | pdg = -11*Int_t(track->GetSign()); |
| | 1352 | } |
| | 1353 | else{ // muon |
| | 1354 | pdg = -13*Int_t(track->GetSign()); |
| | 1355 | } |
| | 1356 | } |
| | 1357 | |
| | 1358 | else if(strstr(option,"sigmacutTOF")){ |
| | 1359 | printf("PID not implemented yet: %s\nNO PID!!!!\n",option); |
| | 1360 | Float_t p = track->P(); |
| | 1361 | |
| | 1362 | // Take expected times |
| | 1363 | Double_t exptimes[5]; |
| | 1364 | track->GetIntegratedTimes(exptimes); |
| | 1365 | |
| | 1366 | // Take resolution for TOF response |
| | 1367 | // like fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[ipart], mass[ipart]); |
| | 1368 | Float_t resolution = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[ipart], mass[ipart]); |
| | 1369 | |
| | 1370 | if(TMath::Abs(exptimes[ipart] - track->GetTOFsignal()) < 3 * resolution){ |
| | 1371 | pdg = pdgvalues[ipart] * Int_t(track->GetSign()); |
| | 1372 | } |
| | 1373 | } |
| | 1374 | |
| | 1375 | else{ |
| | 1376 | printf("Invalid PID option: %s\nNO PID!!!!\n",option); |
| | 1377 | } |
| | 1378 | |
| | 1379 | return pdg; |
| | 1380 | } |
| | 1381 | //----------------------------------------------------------------------- |
| | 1382 | void AliFlowTrackCuts::SetPriors(){ |
| | 1383 | // set abbundancies |
| | 1384 | fBinLimitPID[0] = 0.30; |
| | 1385 | fC[0][0] = 0.015; |
| | 1386 | fC[0][1] = 0.015; |
| | 1387 | fC[0][2] = 1; |
| | 1388 | fC[0][3] = 0.0025; |
| | 1389 | fC[0][4] = 0.000015; |
| | 1390 | fBinLimitPID[1] = 0.35; |
| | 1391 | fC[1][0] = 0.015; |
| | 1392 | fC[1][1] = 0.015; |
| | 1393 | fC[1][2] = 1; |
| | 1394 | fC[1][3] = 0.01; |
| | 1395 | fC[1][4] = 0.001; |
| | 1396 | fBinLimitPID[2] = 0.40; |
| | 1397 | fC[2][0] = 0.015; |
| | 1398 | fC[2][1] = 0.015; |
| | 1399 | fC[2][2] = 1; |
| | 1400 | fC[2][3] = 0.026; |
| | 1401 | fC[2][4] = 0.004; |
| | 1402 | fBinLimitPID[3] = 0.45; |
| | 1403 | fC[3][0] = 0.015; |
| | 1404 | fC[3][1] = 0.015; |
| | 1405 | fC[3][2] = 1; |
| | 1406 | fC[3][3] = 0.026; |
| | 1407 | fC[3][4] = 0.004; |
| | 1408 | fBinLimitPID[4] = 0.50; |
| | 1409 | fC[4][0] = 0.015; |
| | 1410 | fC[4][1] = 0.015; |
| | 1411 | fC[4][2] = 1.000000; |
| | 1412 | fC[4][3] = 0.05; |
| | 1413 | fC[4][4] = 0.01; |
| | 1414 | fBinLimitPID[5] = 0.60; |
| | 1415 | fC[5][0] = 0.012; |
| | 1416 | fC[5][1] = 0.012; |
| | 1417 | fC[5][2] = 1; |
| | 1418 | fC[5][3] = 0.085; |
| | 1419 | fC[5][4] = 0.022; |
| | 1420 | fBinLimitPID[6] = 0.70; |
| | 1421 | fC[6][0] = 0.01; |
| | 1422 | fC[6][1] = 0.01; |
| | 1423 | fC[6][2] = 1; |
| | 1424 | fC[6][3] = 0.12; |
| | 1425 | fC[6][4] = 0.036; |
| | 1426 | fBinLimitPID[7] = 0.80; |
| | 1427 | fC[7][0] = 0.0095; |
| | 1428 | fC[7][1] = 0.0095; |
| | 1429 | fC[7][2] = 1; |
| | 1430 | fC[7][3] = 0.15; |
| | 1431 | fC[7][4] = 0.05; |
| | 1432 | fBinLimitPID[8] = 0.90; |
| | 1433 | fC[8][0] = 0.0085; |
| | 1434 | fC[8][1] = 0.0085; |
| | 1435 | fC[8][2] = 1; |
| | 1436 | fC[8][3] = 0.18; |
| | 1437 | fC[8][4] = 0.074; |
| | 1438 | fBinLimitPID[9] = 1; |
| | 1439 | fC[9][0] = 0.008; |
| | 1440 | fC[9][1] = 0.008; |
| | 1441 | fC[9][2] = 1; |
| | 1442 | fC[9][3] = 0.22; |
| | 1443 | fC[9][4] = 0.1; |
| | 1444 | fBinLimitPID[10] = 1.20; |
| | 1445 | fC[10][0] = 0.007; |
| | 1446 | fC[10][1] = 0.007; |
| | 1447 | fC[10][2] = 1; |
| | 1448 | fC[10][3] = 0.28; |
| | 1449 | fC[10][4] = 0.16; |
| | 1450 | fBinLimitPID[11] = 1.40; |
| | 1451 | fC[11][0] = 0.0066; |
| | 1452 | fC[11][1] = 0.0066; |
| | 1453 | fC[11][2] = 1; |
| | 1454 | fC[11][3] = 0.35; |
| | 1455 | fC[11][4] = 0.23; |
| | 1456 | fBinLimitPID[12] = 1.60; |
| | 1457 | fC[12][0] = 0.0075; |
| | 1458 | fC[12][1] = 0.0075; |
| | 1459 | fC[12][2] = 1; |
| | 1460 | fC[12][3] = 0.40; |
| | 1461 | fC[12][4] = 0.31; |
| | 1462 | fBinLimitPID[13] = 1.80; |
| | 1463 | fC[13][0] = 0.0062; |
| | 1464 | fC[13][1] = 0.0062; |
| | 1465 | fC[13][2] = 1; |
| | 1466 | fC[13][3] = 0.45; |
| | 1467 | fC[13][4] = 0.39; |
| | 1468 | fBinLimitPID[14] = 2.00; |
| | 1469 | fC[14][0] = 0.005; |
| | 1470 | fC[14][1] = 0.005; |
| | 1471 | fC[14][2] = 1; |
| | 1472 | fC[14][3] = 0.46; |
| | 1473 | fC[14][4] = 0.47; |
| | 1474 | fBinLimitPID[15] = 2.20; |
| | 1475 | fC[15][0] = 0.0042; |
| | 1476 | fC[15][1] = 0.0042; |
| | 1477 | fC[15][2] = 1; |
| | 1478 | fC[15][3] = 0.5; |
| | 1479 | fC[15][4] = 0.55; |
| | 1480 | fBinLimitPID[16] = 2.40; |
| | 1481 | fC[16][0] = 0.007; |
| | 1482 | fC[16][1] = 0.007; |
| | 1483 | fC[16][2] = 1; |
| | 1484 | fC[16][3] = 0.5; |
| | 1485 | fC[16][4] = 0.6; |
| | 1486 | |
| | 1487 | for(Int_t i=17;i<fnPIDptBin;i++){ |
| | 1488 | fBinLimitPID[i] = 2.0 + 0.2 * (i-14); |
| | 1489 | fC[i][0] = fC[13][0]; |
| | 1490 | fC[i][1] = fC[13][1]; |
| | 1491 | fC[i][2] = fC[13][2]; |
| | 1492 | fC[i][3] = fC[13][3]; |
| | 1493 | fC[i][4] = fC[13][4]; |
| | 1494 | } |
| | 1495 | } |
| | 1496 | // end part added by F. Noferini |
| | 1497 | //----------------------------------------------------------------------- |
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental