Added new cuts for (TPC pid & TOF match) and (TPC pid & TOF veto) for K* pA analysis
[u/mrichter/AliRoot.git] / PWGLF / RESONANCES / AliRsnCutV0.cxx
CommitLineData
8082e538 1//
2// Class AliRsnCutV0
3//
4// General implementation of a single cut strategy, which can be:
5// - a value contained in a given interval [--> IsBetween() ]
6// - a value equal to a given reference [--> MatchesValue()]
7//
8// In all cases, the reference value(s) is (are) given as data members
9// and each kind of cut requires a given value type (Int, UInt, Double),
10// but the cut check procedure is then automatized and chosen thanks to
11// an enumeration of the implemented cut types.
12// At the end, the user (or any other point which uses this object) has
13// to use the method IsSelected() to check if this cut has been passed.
14//
76d89202 15// authors: Massimo Venaruzzo (massimo.venaruzzo@ts.infn.it)
16// modified: Enrico Fragiacomo (enrico.fragiacomo@ts.infn.it)
8082e538 17//
18
19#include <Riostream.h>
20#include <TFormula.h>
21#include <TBits.h>
22
23#include "AliLog.h"
24#include "AliESDtrackCuts.h"
25
26#include "AliRsnEvent.h"
27#include "AliRsnDaughter.h"
28#include "AliRsnCutV0.h"
29
30ClassImp(AliRsnCutV0)
31
32//_________________________________________________________________________________________________
3da8cef7 33AliRsnCutV0::AliRsnCutV0(const char *name, Int_t hypothesis, AliPID::EParticleType pid, AliPID::EParticleType pid2) :
8082e538 34 AliRsnCut(name, AliRsnTarget::kDaughter),
35 fHypothesis(0),
36 fMass(0.0),
3da8cef7 37 fTolerance(0.01),
8082e538 38 fMaxDCAVertex(0.3),
39 fMinCosPointAngle(0.95),
3da8cef7 40 fMaxDaughtersDCA(0.5),
9d4bb6d8 41 fMinTPCcluster(70),
3da8cef7 42 fMaxRapidity(0.8),
92932daf 43 fDCARptFormula(""),
3da8cef7 44 fPID(pid),
45 fPID2(pid2),
92932daf 46 fPIDCutProton(0),
47 fPIDCutPion(0),
76d89202 48 fESDtrackCuts(0x0),
49 fCutQuality(Form("%sDaughtersQuality", name)),
50 fAODTestFilterBit(5)
8082e538 51{
52//
53// Default constructor.
54// Initializes all cuts in such a way that all of them are disabled.
55//
56
57 SetHypothesis(hypothesis);
58}
59
60//_________________________________________________________________________________________________
61AliRsnCutV0::AliRsnCutV0(const AliRsnCutV0 &copy) :
62 AliRsnCut(copy),
63 fHypothesis(copy.fHypothesis),
64 fMass(copy.fMass),
65 fTolerance(copy.fTolerance),
66 fMaxDCAVertex(copy.fMaxDCAVertex),
67 fMinCosPointAngle(copy.fMinCosPointAngle),
68 fMaxDaughtersDCA(copy.fMaxDaughtersDCA),
9d4bb6d8 69 fMinTPCcluster(copy.fMinTPCcluster),
e6952ec7 70 fMaxRapidity(copy.fMaxRapidity),
92932daf 71 fDCARptFormula(copy.fDCARptFormula),
3da8cef7 72 fPID(copy.fPID),
73 fPID2(copy.fPID2),
92932daf 74 fPIDCutProton(copy.fPIDCutProton),
75 fPIDCutPion(copy.fPIDCutPion),
76d89202 76 fESDtrackCuts(copy.fESDtrackCuts),
77 fCutQuality(copy.fCutQuality),
78 fAODTestFilterBit(copy.fAODTestFilterBit)
8082e538 79{
80//
81// Copy constructor.
76d89202 82// Just copy all data member values.:IsSelected: Object is not a V0 (RESONANCES/AliRsnCutV0.cxx:149)
83
8082e538 84//
76d89202 85 fCutQuality.SetPtRange(0.15, 1E+20);
86 fCutQuality.SetEtaRange(-0.8, 0.8);
92932daf 87 fCutQuality.SetDCARPtFormula(fDCARptFormula);
76d89202 88 fCutQuality.SetDCAZmax(2.0);
89 fCutQuality.SetSPDminNClusters(1);
90 fCutQuality.SetITSminNClusters(0);
91 fCutQuality.SetITSmaxChi2(1E+20);
9d4bb6d8 92 fCutQuality.SetTPCminNClusters(fMinTPCcluster);
76d89202 93 fCutQuality.SetTPCmaxChi2(4.0);
94 fCutQuality.SetRejectKinkDaughters();
95 fCutQuality.SetAODTestFilterBit(5);
96
8082e538 97}
98
99//_________________________________________________________________________________________________
e6952ec7 100AliRsnCutV0 &AliRsnCutV0::operator=(const AliRsnCutV0 &copy)
8082e538 101{
102//
103// Assignment operator.
104// Just copy all data member values.
105//
4cc7af7d 106 if (this == &copy)
107 return *this;
61f275d1 108 fHypothesis = copy.fHypothesis;
109 fMass = copy.fMass;
110 fTolerance = copy.fTolerance;
111 fMaxDCAVertex = copy.fMaxDCAVertex;
112 fMinCosPointAngle = copy.fMinCosPointAngle;
113 fMaxDaughtersDCA = copy.fMaxDaughtersDCA;
9d4bb6d8 114 fMinTPCcluster = copy.fMinTPCcluster;
3da8cef7 115 fMaxRapidity = copy.fMaxRapidity;
92932daf 116 fDCARptFormula = copy.fDCARptFormula;
117 fCutQuality = copy.fCutQuality;
3da8cef7 118 fPID = copy.fPID;
119 fPID2 = copy.fPID2;
92932daf 120 fPIDCutProton = copy.fPIDCutProton;
121 fPIDCutPion = copy.fPIDCutPion;
61f275d1 122 fESDtrackCuts = copy.fESDtrackCuts;
76d89202 123 fCutQuality = copy.fCutQuality;
124 fAODTestFilterBit = copy.fAODTestFilterBit;
61f275d1 125
8082e538 126 return (*this);
e6952ec7 127}
8082e538 128
129//_________________________________________________________________________________________________
130Bool_t AliRsnCutV0::IsSelected(TObject *object)
131{
76d89202 132//:IsSelected: Object is not a V0 (RESONANCES/AliRsnCutV0.cxx:149)
133
8082e538 134// Cut checker.
135// Checks the type of object being evaluated
136// and then calls the appropriate sub-function (for ESD or AOD)
137//
138
139 // coherence check
140 if (!TargetOK(object)) return kFALSE;
e6952ec7 141
8082e538 142 // check cast
143 AliESDv0 *v0esd = fDaughter->Ref2ESDv0();
144 AliAODv0 *v0aod = fDaughter->Ref2AODv0();
145 //cout << fDaughter->GetRef()->ClassName() << ' ' << v0esd << ' ' << v0aod << endl;
e6952ec7 146
76d89202 147 // operate depending on cast:IsSelected: Object is not a V0 (RESONANCES/AliRsnCutV0.cxx:149)
61f275d1 148
8082e538 149 if (v0esd) {
150 return CheckESD(v0esd);
151 } else if (v0aod) {
152 return CheckAOD(v0aod);
153 } else {
154 AliDebugClass(1, "Object is not a V0");
155 return kFALSE;
156 }
157}
158
159//_________________________________________________________________________________________________
160Bool_t AliRsnCutV0::CheckESD(AliESDv0 *v0)
161{
162//
163// Check an ESD V0.
164// This is done using the default track checker for ESD.
165// It is declared static, not to recreate it every time.
166//
167
168 AliDebugClass(1, "Check ESD");
169 if (v0->GetOnFlyStatus()) {
170 AliDebugClass(1, "Rejecting V0 in 'on fly' status");
171 return kFALSE; // if kTRUE, then this V0 is recontructed
172 }
e6952ec7 173
8082e538 174 // retrieve pointer to owner event
175 AliESDEvent *lESDEvent = fEvent->GetRefESD();
176 Double_t xPrimaryVertex = lESDEvent->GetPrimaryVertex()->GetX();
177 Double_t yPrimaryVertex = lESDEvent->GetPrimaryVertex()->GetY();
178 Double_t zPrimaryVertex = lESDEvent->GetPrimaryVertex()->GetZ();
179 AliDebugClass(2, Form("Primary vertex: %f %f %f", xPrimaryVertex, yPrimaryVertex, zPrimaryVertex));
e6952ec7 180
8082e538 181 // retrieve the V0 daughters
182 UInt_t lIdxPos = (UInt_t) TMath::Abs(v0->GetPindex());
183 UInt_t lIdxNeg = (UInt_t) TMath::Abs(v0->GetNindex());
184 AliESDtrack *pTrack = lESDEvent->GetTrack(lIdxPos);
185 AliESDtrack *nTrack = lESDEvent->GetTrack(lIdxNeg);
e6952ec7 186
8082e538 187 // check quality cuts
188 if (fESDtrackCuts) {
189 AliDebugClass(2, "Checking quality cuts");
190 if (!fESDtrackCuts->IsSelected(pTrack)) {
191 AliDebugClass(2, "Positive daughter failed quality cuts");
192 return kFALSE;
193 }
194 if (!fESDtrackCuts->IsSelected(nTrack)) {
195 AliDebugClass(2, "Negative daughter failed quality cuts");
196 return kFALSE;
197 }
198 }
e6952ec7 199
8082e538 200 // filter like-sign V0
e6952ec7 201 if ( TMath::Abs( ((pTrack->GetSign()) - (nTrack->GetSign())) ) < 0.1) {
202 AliDebugClass(2, "Failed like-sign V0 check");
203 return kFALSE;
204 }
205
61f275d1 206
8082e538 207 // check compatibility with expected species hypothesis
208 v0->ChangeMassHypothesis(fHypothesis);
209 if ((TMath::Abs(v0->GetEffMass() - fMass)) > fTolerance) {
210 AliDebugClass(2, "V0 is not in the expected inv mass range");
211 return kFALSE;
212 }
e6952ec7 213
8082e538 214 // topological checks
215 if (TMath::Abs(v0->GetD(xPrimaryVertex, yPrimaryVertex, zPrimaryVertex)) > fMaxDCAVertex) {
216 AliDebugClass(2, "Failed check on DCA to primary vertes");
217 return kFALSE;
218 }
219 if (TMath::Abs(v0->GetV0CosineOfPointingAngle()) < fMinCosPointAngle) {
220 AliDebugClass(2, "Failed check on cosine of pointing angle");
221 return kFALSE;
222 }
223 if (TMath::Abs(v0->GetDcaV0Daughters()) > fMaxDaughtersDCA) {
224 AliDebugClass(2, "Failed check on DCA between daughters");
225 return kFALSE;
226 }
3da8cef7 227 if (TMath::Abs(v0->Y(fHypothesis)) > fMaxRapidity) {
e6952ec7 228 AliDebugClass(2, "Failed check on V0 rapidity");
229 return kFALSE;
76d89202 230 }
e6952ec7 231
232
76d89202 233 // check PID on proton or antiproton from V0
e6952ec7 234
235 // check initialization of PID object
236 AliPIDResponse *pid = fEvent->GetPIDResponse();
237 if (!pid) {
238 AliFatal("NULL PID response");
239 return kFALSE;
240 }
241
242 // check if TOF is matched
243 // and computes all values used in the PID cut
244 //Bool_t isTOFpos = MatchTOF(ptrack);
245 //Bool_t isTOFneg = MatchTOF(ntrack);
92932daf 246 //Double_t pospTPC = pTrack->GetTPCmomentum();
247 //Double_t negpTPC = nTrack->GetTPCmomentum();
e6952ec7 248 //Double_t posp = pTrack->P();
249 //Double_t negp = nTrack->P();
250 Double_t posnsTPC = TMath::Abs(pid->NumberOfSigmasTPC(pTrack, fPID));
251 Double_t posnsTPC2 = TMath::Abs(pid->NumberOfSigmasTPC(pTrack, fPID2));
252 //Double_t posnsTOF = TMath::Abs(pid->NumberOfSigmasTOF(ptrack, fPID));
253 Double_t negnsTPC = TMath::Abs(pid->NumberOfSigmasTPC(nTrack, fPID));
254 Double_t negnsTPC2 = TMath::Abs(pid->NumberOfSigmasTPC(nTrack, fPID2));
255 //Double_t negnsTOF = TMath::Abs(pid->NumberOfSigmasTOF(ntrack, fPID));
256 Double_t maxTPC = 1E20;
257 Double_t maxTPC2 = 1E20;
258 //Double_t maxTOF = 1E20;
259
260 // applies the cut differently depending on the PID and the momentum
261
262 if(fHypothesis==kLambda0) {
263 //if (isTOFpos) {
264 // TPC: 5sigma cut for all
265 //if (posnsTPC > 5.0) return kFALSE;
266 // TOF: 3sigma
267 // maxTOF = 3.0;
268 //return (posnsTOF <= maxTOF);
269 //} else {
270 // TPC:
e6952ec7 271
92932daf 272 maxTPC = fPIDCutProton;
273 maxTPC2 = fPIDCutPion;
e6952ec7 274
275 if (! ((posnsTPC <= maxTPC) && (negnsTPC2 <= maxTPC2)) ) {
276 AliDebugClass(2, "Failed check on V0 PID");
277 return kFALSE;
278 }
279 }
280
281 //}
282
283 if(fHypothesis==kLambda0Bar) {
284 //if (isTOFneg) {
285 // TPC: 5sigma cut for all
286 //if (negnsTPC > 5.0) return kFALSE;
287 // TOF: 3sigma
288 // maxTOF = 3.0;
289 //return (negnsTOF <= maxTOF);
290 //} else {
291 // TPC:
e6952ec7 292
92932daf 293
294 maxTPC = fPIDCutProton;
295 maxTPC2 = fPIDCutPion;
e6952ec7 296
297 if(! ((negnsTPC <= maxTPC) && (posnsTPC2 <= maxTPC2)) ) {
298 AliDebugClass(2, "Failed check on V0 PID");
299 return kFALSE;
300 }
301 }
302 //}
303
304
76d89202 305 // if we reach this point, all checks were successful
306 AliDebugClass(2, "Good V0 (hallelujah)");
e6952ec7 307 return kTRUE;
76d89202 308}
309
310//_________________________________________________________________________________________________
311Bool_t AliRsnCutV0::CheckAOD(AliAODv0 *v0)
312{
313//
314// Check an AOD V0.
315// This is done doing directly all checks, since there is not
316// an equivalend checker for AOD tracks
317//
318
319 AliDebugClass(2, "Check AOD");
320 if (v0->GetOnFlyStatus()) {
321 AliDebugClass(2, "Rejecting V0 in 'on fly' status");
322 return kFALSE; // if kTRUE, then this V0 is recontructed
323 }
324
e6952ec7 325
326
76d89202 327 // retrieve pointer to owner event
328 AliAODEvent *lAODEvent = fEvent->GetRefAOD();
329 Double_t xPrimaryVertex = lAODEvent->GetPrimaryVertex()->GetX();
330 Double_t yPrimaryVertex = lAODEvent->GetPrimaryVertex()->GetY();
331 Double_t zPrimaryVertex = lAODEvent->GetPrimaryVertex()->GetZ();
332 AliDebugClass(2, Form("Primary vertex: %f %f %f", xPrimaryVertex, yPrimaryVertex, zPrimaryVertex));
333
334 // retrieve the V0 daughters
335 AliAODTrack *pTrack = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
336 AliAODTrack *nTrack = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));
337
338
339 // check quality cuts
340 UInt_t filtermapP = 9999;
341 UInt_t filtermapN = 9999;
342 filtermapP = pTrack->GetFilterMap();
343 filtermapN = nTrack->GetFilterMap();
e6952ec7 344
345
76d89202 346 if ( !pTrack->TestFilterBit(fAODTestFilterBit) ) {
e6952ec7 347 AliDebugClass(2, Form("Positive daughter failed quality cuts filtermapP=%d",filtermapP));
348 return kFALSE;
76d89202 349 }
350 if ( !nTrack->TestFilterBit(fAODTestFilterBit) ) {
e6952ec7 351 AliDebugClass(2, Form("Negative daughter failed quality cuts filtermapN=%d",filtermapN));
352 return kFALSE;
76d89202 353 }
e6952ec7 354
355
76d89202 356
1c7bcf66 357 /*AliDebugClass(1, Form("fESDtrackCuts=%p",fESDtrackCuts));
76d89202 358 if (fESDtrackCuts) { // use fESDtrackCuts to retrieve cuts values
359
e6952ec7 360
361 AliDebugClass(2, "Checking quality cuts");
362
363
364
365 const Bool_t CutAcceptKinkDaughters = fESDtrackCuts->GetAcceptKinkDaughters(); // 0 = kFalse
366 const Float_t CutMaxChi2PerClusterTPC = fESDtrackCuts->GetMaxChi2PerClusterTPC();
367 const Int_t CutMinNClusterTPC = fESDtrackCuts->GetMinNClusterTPC();
368 const Bool_t CutRequireTPCRefit = fESDtrackCuts->GetRequireTPCRefit();
369 //AliDebugClass(2, Form("accept kink=%d maxchi2=%f minnclSTPC=%d requireTPCrefit=%d", CutAcceptKinkDaughters,
370 // CutMaxChi2PerClusterTPC, CutMinNClusterTPC, CutRequireTPCRefit));
371 //AliDebugClass(2, Form("pTrack->TPCNcls=%d", pTrack->GetTPCNcls() ));
372 //AliDebugClass(2, Form("nTrack->TPCNcls=%d", nTrack->GetTPCNcls() ));
373
374
76d89202 375 if(pTrack->GetTPCNcls() < CutMinNClusterTPC) {AliDebugClass(2, "Positive daughter not MinNclsTPC"); return kFALSE;}
376 if(nTrack->GetTPCNcls() < CutMinNClusterTPC) {AliDebugClass(2, "Negative daughter not MinNclsTPC"); return kFALSE;}
377
1c7bcf66 378 }*/
76d89202 379
380 // filter like-sign V0
e6952ec7 381 if ( TMath::Abs( ((pTrack->Charge()) - (nTrack->Charge())) ) < 0.1) {
382 AliDebugClass(2, "Failed like-sign V0 check");
383 return kFALSE;
384 }
76d89202 385
386 // check compatibility with expected species hypothesis
e6952ec7 387 Double_t mass = 0.0;
388 if(fHypothesis==kLambda0) {
389 mass = v0->MassLambda();
390 }
391 else if (fHypothesis==kLambda0Bar) {
392 mass = v0->MassAntiLambda();
393 }
394 if ((TMath::Abs(mass - fMass)) > fTolerance) {
395 AliDebugClass(2, Form("V0 is not in the expected inv mass range Mass: %d %f %f", fHypothesis, fMass, mass));
396 return kFALSE;
397 }
76d89202 398 AliDebugClass(2, Form("Mass: %d %f %f", fHypothesis, fMass, mass));
e6952ec7 399
400
401 // topological checks
402 if (TMath::Abs(v0->DcaV0ToPrimVertex()) > fMaxDCAVertex) {
403 AliDebugClass(2, Form("Failed check on DCA to primary vertes dca=%f maxdca=%f",TMath::Abs(v0->DcaV0ToPrimVertex()),fMaxDCAVertex));
404 return kFALSE;
405 }
406 if (TMath::Abs( TMath::Cos(v0->OpenAngleV0()) ) < fMinCosPointAngle) {
76d89202 407 AliDebugClass(2, "Failed check on cosine of pointing angle");
3da8cef7 408 return kFALSE;
409 }
76d89202 410 if (TMath::Abs(v0->DcaV0Daughters()) > fMaxDaughtersDCA) {
411 AliDebugClass(2, "Failed check on DCA between daughters");
412 return kFALSE;
413 }
414 if (TMath::Abs(v0->RapLambda()) > fMaxRapidity) {
e6952ec7 415 AliDebugClass(2, "Failed check on V0 rapidity");
416 return kFALSE;
76d89202 417 }
3da8cef7 418
419
76d89202 420
e6952ec7 421
3da8cef7 422
423 // check initialization of PID object
424 AliPIDResponse *pid = fEvent->GetPIDResponse();
425 if (!pid) {
e6952ec7 426 AliFatal("NULL PID response");
427 return kFALSE;
3da8cef7 428 }
429
430 // check if TOF is matched
431 // and computes all values used in the PID cut
432 //Bool_t isTOFpos = MatchTOF(ptrack);
433 //Bool_t isTOFneg = MatchTOF(ntrack);
92932daf 434 //Double_t pospTPC = pTrack->GetTPCmomentum();
435 //Double_t negpTPC = nTrack->GetTPCmomentum();
3da8cef7 436 //Double_t posp = pTrack->P();
437 //Double_t negp = nTrack->P();
438 Double_t posnsTPC = TMath::Abs(pid->NumberOfSigmasTPC(pTrack, fPID));
439 Double_t posnsTPC2 = TMath::Abs(pid->NumberOfSigmasTPC(pTrack, fPID2));
440 //Double_t posnsTOF = TMath::Abs(pid->NumberOfSigmasTOF(ptrack, fPID));
441 Double_t negnsTPC = TMath::Abs(pid->NumberOfSigmasTPC(nTrack, fPID));
442 Double_t negnsTPC2 = TMath::Abs(pid->NumberOfSigmasTPC(nTrack, fPID2));
443 //Double_t negnsTOF = TMath::Abs(pid->NumberOfSigmasTOF(ntrack, fPID));
444 Double_t maxTPC = 1E20;
445 Double_t maxTPC2 = 1E20;
446 //Double_t maxTOF = 1E20;
e6952ec7 447
3da8cef7 448 // applies the cut differently depending on the PID and the momentum
e6952ec7 449
450 if(fHypothesis==kLambda0) {
451
452
453 //if (isTOFpos) {
454 // TPC: 5sigma cut for all
455 //if (posnsTPC > 5.0) return kFALSE;
456 // TOF: 3sigma
457 // maxTOF = 3.0;
458 //return (posnsTOF <= maxTOF);
459 //} else {
460 // TPC:
92932daf 461
462
463 maxTPC = fPIDCutProton;
464 maxTPC2 = fPIDCutPion;
e6952ec7 465
466 if (! ((posnsTPC <= maxTPC) && (negnsTPC2 <= maxTPC2)) ) {
467 AliDebugClass(2, "Failed check on V0 PID");
468 return kFALSE;
469 }
3da8cef7 470 }
e6952ec7 471
472 if(fHypothesis==kLambda0Bar) {
473
474 //if (isTOFneg) {
475 // TPC: 5sigma cut for all
476 //if (negnsTPC > 5.0) return kFALSE;
477 // TOF: 3sigma
478 // maxTOF = 3.0;
479 //return (negnsTOF <= maxTOF);
480 //} else {
481 // TPC:
e6952ec7 482
92932daf 483
484 maxTPC = fPIDCutProton;
485 maxTPC2 = fPIDCutPion;
e6952ec7 486
487 if(! ((negnsTPC <= maxTPC) && (posnsTPC2 <= maxTPC2)) ) {
488 AliDebugClass(2, "Failed check on V0 PID");
489 return kFALSE;
490 }
3da8cef7 491 }
e6952ec7 492
493
61f275d1 494
8082e538 495 // if we reach this point, all checks were successful
76d89202 496 AliDebugClass(1, "Good AOD V0 (hallelujah)");
497 AliDebugClass(1, Form("Mass: %d %f %f %d %d", fHypothesis, fMass, mass, filtermapP, filtermapN));
e6952ec7 498 return kTRUE;
61f275d1 499
8082e538 500}
501
502//_________________________________________________________________________________________________
503void AliRsnCutV0::Print(const Option_t *) const
504{
505//
506// Print information on this cut
507//
508}