]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliESDv0KineCuts.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Added commet to explain AliAODTrack::fType, AliAODTrack::XAtDCA(), YAtDCA(), ZAtDCA...
[u/mrichter/AliRoot.git] / ANALYSIS / AliESDv0KineCuts.cxx
CommitLineData
a65a7e70 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/*
17 * author: M.Kalisky@gsi.de
18 * 08/Dec/2010
19 *
20 * Description: This class allows with purely kinematical cuts
21 * to select clean samples of electrons, pions and protons from the
22 * V0 online finder ESD V0 candidates for PID and dectector resonse
23 * studies.
24 */
25
26#include <TVector3.h>
27#include <TDatabasePDG.h>
28
29#include "AliESDv0.h"
30#include "AliESDtrack.h"
31#include "AliESDEvent.h"
32#include "AliVEvent.h"
33#include "AliLog.h"
34#include "AliKFParticle.h"
35#include "AliVTrack.h"
36#include "AliKFVertex.h"
37
38#include "AliESDv0KineCuts.h"
39
40ClassImp(AliESDv0KineCuts)
41
42//____________________________________________________________________
43AliESDv0KineCuts::AliESDv0KineCuts() :
44 fEvent(0x0)
45 , fPrimaryVertex(0x0)
46 , fType(0)
47 , fMode(0)
48 , fTPCNcls(1)
49 , fTPCrefit(kTRUE)
50 , fTPCchi2perCls(4.0)
51 , fTPCclsRatio(0.6)
52 , fNoKinks(kTRUE)
53 , fGcutChi2NDF(10)
54 , fGcutInvMass(0.05)
55 , fK0cutChi2NDF(10)
56 , fLcutChi2NDF(10)
57 , fUseExternalVertex(kFALSE)
58 , fDeleteVertex(kFALSE)
59{
60 //
61 // Default constructor
62 //
63
64 // default single track cuts
65 fTPCNcls = 1; // minimal number of the TPC clusters
66 fTPCrefit = kTRUE; // TPC refit
67 fTPCchi2perCls = 4.0; // chi2 per TPC cluster
68 fTPCclsRatio = 0.6; // minimal foun/findable TPC cluster ratio
69 fNoKinks = kTRUE; // kinks - no [kTRUE] or do not care [kFalse]
70
71
72 // default gamma cuts values
73 fGcutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle gamma
74 fGcutCosPoint[0] = 0; // cos of the pointing angle [min, max]
75 fGcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
76 fGcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
77 fGcutDCA[1] = 0.25; // DCA between the daughter tracks [min, max]
78 fGcutVertexR[0] = 3.; // radius of the conversion point [min, max]
79 fGcutVertexR[1] = 90.; // radius of the conversion point [min, max]
80 fGcutPsiPair[0] = 0.; // value of the psi pair cut [min, max]
81 fGcutPsiPair[1] = 0.05; // value of the psi pair cut [min, max]
82 fGcutInvMass = 0.05; // upper value on the gamma invariant mass
83 // default K0 cuts
84 fK0cutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle K0
85 fK0cutCosPoint[0] = 0.; // cos of the pointing angle [min, max]
86 fK0cutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
87 fK0cutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
88 fK0cutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]
89 fK0cutVertexR[0] = 2.0; // radius of the decay point [min, max]
90 fK0cutVertexR[1] = 30.0; // radius of the decay point [min, max]
91 fK0cutInvMass[0] = 0.486; // invariant mass window
92 fK0cutInvMass[1] = 0.508; // invariant mass window
93 // Lambda & anti-Lambda cut values
94 fLcutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle K0
95 fLcutCosPoint[0] = 0.; // cos of the pointing angle [min, max]
96 fLcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
97 fLcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
98 fLcutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]
99 fLcutVertexR[0] = 2.0; // radius of the decay point [min, max]
100 fLcutVertexR[1] = 40.0; // radius of the decay point [min, max]
101 fLcutInvMass[0] = 1.11; // invariant mass window
102 fLcutInvMass[1] = 1.12; // invariant mass window
103
104}
105//____________________________________________________________________
106AliESDv0KineCuts::~AliESDv0KineCuts(){
107 //
108 // Destructor
109 //
110
111
112}
113//____________________________________________________________________
114AliESDv0KineCuts::AliESDv0KineCuts(const AliESDv0KineCuts &ref):
115 TObject(ref)
116 , fEvent(0x0)
117 , fPrimaryVertex(0x0)
118 , fType(0)
119 , fMode(0)
120 , fTPCNcls(1)
121 , fTPCrefit(kTRUE)
122 , fTPCchi2perCls(4.0)
123 , fTPCclsRatio(0.6)
124 , fNoKinks(kTRUE)
125 , fGcutChi2NDF(10)
126 , fGcutInvMass(0.05)
127 , fK0cutChi2NDF(10)
128 , fLcutChi2NDF(10)
129 , fUseExternalVertex(kFALSE)
130 , fDeleteVertex(kFALSE)
131{
132 //
133 // Copy operator
134 //
135
136 ref.Copy(*this);
137}
138//____________________________________________________________________
139AliESDv0KineCuts &AliESDv0KineCuts::operator=(const AliESDv0KineCuts &ref){
140 //
141 // assignment operator
142 //
143 if(this != &ref)
144 ref.Copy(*this);
145 return *this;
146}
147//____________________________________________________________________
148void AliESDv0KineCuts::Copy(TObject &ref) const {
149 //
150 // Performs the copying of the object
151 //
152
153 TObject::Copy(ref);
154
155 AliESDv0KineCuts &target = dynamic_cast<AliESDv0KineCuts &>(ref);
156
157 // default single track cuts
158 target.fTPCNcls = fTPCNcls;
159 target.fTPCrefit = fTPCrefit;
160 target.fTPCchi2perCls = fTPCchi2perCls;
161 target.fTPCclsRatio = fTPCclsRatio;
162 target.fNoKinks = fNoKinks;
163 target.fUseExternalVertex = fUseExternalVertex; //added december 2nd 2011
164 target.fDeleteVertex = fDeleteVertex; //added december 2nd 2011
165
166 // default gamma cuts values
167 target.fGcutChi2NDF = fGcutChi2NDF;
168 memcpy(target.fGcutCosPoint, fGcutCosPoint, sizeof(Float_t) * 2);
169 memcpy(target.fGcutDCA, fGcutDCA, sizeof(Float_t) * 2);
170 memcpy(target.fGcutVertexR, fGcutVertexR, sizeof(Float_t) * 2);
171 memcpy(target.fGcutPsiPair, fGcutPsiPair, sizeof(Float_t) * 2);
172 target.fGcutInvMass = fGcutInvMass;
173 // default K0 cuts
174 target.fK0cutChi2NDF = fK0cutChi2NDF;
175 memcpy(target.fK0cutCosPoint, fK0cutCosPoint, sizeof(Float_t) * 2);
176 memcpy(target.fK0cutDCA, fK0cutDCA, sizeof(Float_t) * 2);
177 memcpy(target.fK0cutVertexR, fK0cutVertexR, sizeof(Float_t) * 2);
178 memcpy(target.fK0cutInvMass, fK0cutInvMass, sizeof(Float_t) * 2);
179 // Lambda & anti-Lambda cut values
180 target.fLcutChi2NDF = fLcutChi2NDF;
181 memcpy(target.fLcutCosPoint, fLcutCosPoint, sizeof(Float_t) * 2);
182 memcpy(target.fLcutDCA, fLcutDCA, sizeof(Float_t) * 2);
183 memcpy(target.fLcutVertexR, fLcutVertexR, sizeof(Float_t) * 2);
184 memcpy(target.fLcutInvMass, fLcutInvMass, sizeof(Float_t) * 2);
185
186}
187//____________________________________________________________________
188Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const
189{
190 //
191 // main user function
192 //
193
194 if(!v0) return kFALSE;
195 if(!fEvent){
196 AliErrorClass("No valid Event pointer available, provide it first");
197 return kFALSE;
198 }
199
200 if(!V0CutsCommon(v0)) return kFALSE;
201
202 const Int_t id = PreselectV0(v0);
203
204 if(!SingleTrackCuts(v0)) return kFALSE;
205
206 switch(id){
207 case kUndef:
208 return kFALSE;
209 case kGamma:
210 return CaseGamma(v0, pdgV0, pdgP, pdgN);
211 case kK0:
212 return CaseK0(v0, pdgV0, pdgP, pdgN);
213 case kLambda:
214 return CaseLambda(v0, pdgV0, pdgP, pdgN, 0);
215 case kALambda:
216 return CaseLambda(v0, pdgV0, pdgP, pdgN, 1);
217 default:
218 return kFALSE;
219 }
220
221 return kFALSE;
222}
223//____________________________________________________________________
224Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgP, Int_t &pdgN) const
225{
226 //
227 // main user function, simplified if the V0 identity is not necessary
228 //
229
230 if(!v0) return kFALSE;
231 if(!fEvent){
232 AliErrorClass("No valid Event pointer available, provide it first");
233 return kFALSE;
234 }
235
236 Int_t idV0 = -1;
237 return ProcessV0(v0, idV0, pdgP, pdgN);
238
239}
240//____________________________________________________________________
241Int_t AliESDv0KineCuts::PreselectV0(AliESDv0* const v0) const
242{
243 //
244 // Make a preselection (exclusive) of the V0 cadidates based on
245 // Armenteros plot
246 // the armenteros cut values are currently fixed and user is not able to set them via
247 // set funcions. The reason is that these cuts are optimized and furneter changes should
248 // not be necessary. To prove otherwise please study in detail before changing the values
249 //
250
251 Float_t ap[2] = {-1., -1.};
252 Armenteros(v0, ap);
253 // for clarity
254 const Float_t alpha = ap[0];
255 const Float_t qt = ap[1];
256
257 // selection cuts
258 // - the reagions for different candidates must not overlap
259
260 // Gamma cuts
261 const Double_t cutAlphaG = 0.35;
262 const Double_t cutQTG = 0.05;
263 const Double_t cutAlphaG2[2] = {0.6, 0.8};
264 const Double_t cutQTG2 = 0.04;
265
266 // K0 cuts
267 const Float_t cutQTK0[2] = {0.1075, 0.215};
268 const Float_t cutAPK0[2] = {0.199, 0.8}; // parameters for curved QT cut
269
270 // Lambda & A-Lambda cuts
271 const Float_t cutQTL = 0.03;
272 const Float_t cutAlphaL[2] = {0.35, 0.7};
273 const Float_t cutAlphaAL[2] = {-0.7, -0.35};
274 const Float_t cutAPL[3] = {0.107, -0.69, 0.5}; // parameters fir curved QT cut
275
276
277 if(kPurity == fMode){
278 // Check for Gamma candidates
279 if(qt < cutQTG){
280 if( (TMath::Abs(alpha) < cutAlphaG) ) return kGamma;
281 }
282 // additional region - should help high pT gammas
283 if(qt < cutQTG2){
284 if( (TMath::Abs(alpha) > cutAlphaG2[0]) && (TMath::Abs(alpha) < cutAlphaG2[1]) ) return kGamma;
285 }
286 }
287 if(kEffGamma == fMode){
288 if(qt < cutQTG) return kGamma;
289 }
290
291
292 // Check for K0 candidates
293 Float_t q = cutAPK0[0] * TMath::Sqrt(TMath::Abs(1 - alpha*alpha/(cutAPK0[1]*cutAPK0[1])));
294 if( (qt > cutQTK0[0]) && (qt < cutQTK0[1]) && (qt > q) ){
295 return kK0;
296 }
297
298 // Check for Lambda candidates
299 q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha + cutAPL[1]) * (alpha + cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));
300 if( (alpha > cutAlphaL[0]) && (alpha < cutAlphaL[1]) && (qt > cutQTL) && (qt < q) ){
301 return kLambda;
302 }
303
304 // Check for A-Lambda candidates
305 q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha - cutAPL[1]) * (alpha - cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));
306 if( (alpha > cutAlphaAL[0]) && (alpha < cutAlphaAL[1]) && (qt > cutQTL) && (qt < q) ){
307 return kALambda;
308 }
309
310 return kUndef;
311}
312//____________________________________________________________________
313Bool_t AliESDv0KineCuts::SingleTrackCuts(AliESDv0 * const v0) const
314{
315 //
316 // apply single track cuts
317 // correct sign not relevat here
318 //
319
320 if(!v0) return kFALSE;
321
322 Int_t pIndex = 0, nIndex = 0;
323 pIndex = v0->GetPindex();
324 nIndex = v0->GetNindex();
325 AliESDtrack* d[2];
326 d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
327 d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
328
329 for(Int_t i=0; i<2; ++i){
330 if(!d[i]) return kFALSE;
331
332 // status word
333 ULong_t status = d[i]->GetStatus();
334
335 // No. of TPC clusters leave to the users
336 if(d[i]->GetTPCNcls() < 1) return kFALSE;
337
338 // TPC refit
339 if(!(status & AliESDtrack::kTPCrefit)) return kFALSE;
340
341 // Chi2 per TPC cluster
342 Int_t nTPCclusters = d[i]->GetTPCNcls();
343 Float_t chi2perTPCcluster = d[i]->GetTPCchi2()/Float_t(nTPCclusters);
344 if(chi2perTPCcluster > 4) return kFALSE;
345
346 // TPC cluster ratio
347 Float_t cRatioTPC = d[i]->GetTPCNclsF() > 0. ? static_cast<Float_t>(d[i]->GetTPCNcls())/static_cast<Float_t> (d[i]->GetTPCNclsF()) : 1.;
348 if(cRatioTPC < 0.6) return kFALSE;
349
350 // kinks
351 if(d[i]->GetKinkIndex(0) != 0) return kFALSE;
352
353 }
354
355 return kTRUE;
356}
357//____________________________________________________________________
358Bool_t AliESDv0KineCuts::CaseGamma(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const
359{
360 //
361 // process the gamma conversion candidate
362 //
363
364 if(!v0) return kFALSE;
365
366 AliVTrack* daughter[2];
367 Int_t pIndex = 0, nIndex = 0;
368
369 Bool_t sign = CheckSigns(v0);
370 if(sign){
371 pIndex = v0->GetPindex();
372 nIndex = v0->GetNindex();
373 }
374 else{
375 pIndex = v0->GetNindex();
376 nIndex = v0->GetPindex();
377 }
378 daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
379 daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
380 if(!daughter[0] || !daughter[1]) return kFALSE;
381
382 AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kElectron), TMath::Abs(kElectron));
383 if(!kfMother) return kFALSE;
384
385 AliESDtrack* d[2];
386 d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
387 d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
388
389 Float_t iMass = v0->GetEffMass(0, 0);
390
391 // cos pointing angle
392 Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
393 cosPoint = TMath::ACos(cosPoint);
394
395 // DCA between daughters
396 Double_t dca = v0->GetDcaV0Daughters();
397
398 // Production vertex
399 Double_t x, y, z;
400 v0->GetXYZ(x,y,z);
401 Double_t r = TMath::Sqrt(x*x + y*y);
402
403 Double_t xy[2];
404 Double_t r2 = -1.;
405 if ( GetConvPosXY(d[0], d[1], xy) ){
406 r2 = TMath::Sqrt(xy[0]*xy[0] + xy[1]*xy[1]);
407 }
408
409 // psi pair
410 Double_t psiPair = PsiPair(v0);
411
412 // V0 chi2/ndf
413 Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
414
415 if(kfMother) delete kfMother;
416
417 // apply the cuts
418
419 if(iMass > fGcutInvMass) return kFALSE;
420
421 if(chi2ndf > fGcutChi2NDF) return kFALSE;
422
423 if(cosPoint < fGcutCosPoint[0] || cosPoint > fGcutCosPoint[1]) return kFALSE;
424
425 if(dca < fGcutDCA[0] || dca > fGcutDCA[1]) return kFALSE;
426
427 if(r < fGcutVertexR[0] || r > fGcutVertexR[1]) return kFALSE;
428
429 if(psiPair < fGcutPsiPair[0] || psiPair > fGcutPsiPair[1]) return kFALSE;
430
431 // all cuts passed
432
433 pdgV0 = 22;
434 if(sign){
435 pdgP = -11;
436 pdgN = 11;
437 }
438 else{
439 pdgP = 11;
440 pdgN = -11;
441 }
442
443 return kTRUE;
444}
445//____________________________________________________________________
446Bool_t AliESDv0KineCuts::CaseK0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const {
447 //
448 // process the K0 candidate
449 //
450
451 if(!v0) return kFALSE;
452
453 AliVTrack* daughter[2];
454 Int_t pIndex = 0, nIndex = 0;
455 Bool_t sign = CheckSigns(v0);
456 if(sign){
457 pIndex = v0->GetPindex();
458 nIndex = v0->GetNindex();
459 }
460 else{
461 pIndex = v0->GetNindex();
462 nIndex = v0->GetPindex();
463 }
464
465 daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
466 daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
467 if(!daughter[0] || !daughter[1]) return kFALSE;
468
469 AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kPiPlus));
470 if(!kfMother) return kFALSE;
471
472 AliESDtrack* d[2];
473 d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
474 d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
475
476 Float_t iMass = v0->GetEffMass(2, 2);
477
478 // cos pointing angle
479 Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
480 cosPoint = TMath::ACos(cosPoint);
481
482 // DCA between daughters
483 Double_t dca = v0->GetDcaV0Daughters();
484
485 // Production vertex
486 Double_t x, y, z;
487 v0->GetXYZ(x,y,z);
488
489 Double_t r = TMath::Sqrt(x*x + y*y);
490
491 // V0 chi2/ndf
492 Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
493
494 if(kfMother) delete kfMother;
495
496 //
497 // apply the cuts
498 //
499 if(iMass < fK0cutInvMass[0] || iMass > fK0cutInvMass[1]) return kFALSE;
500
501 if(chi2ndf > fK0cutChi2NDF) return kFALSE;
502
503 if(cosPoint < fK0cutCosPoint[0] || cosPoint > fK0cutCosPoint[1]) return kFALSE;
504
505 if(dca < fK0cutDCA[0] || dca > fK0cutDCA[1]) return kFALSE;
506
507 if(r < fK0cutVertexR[0] || r > fK0cutVertexR[1]) return kFALSE;
508
509 // all cuts passed
510 pdgV0 = 310;
511 if(sign){
512 pdgP = 211;
513 pdgN = -211;
514 }
515 else{
516 pdgP = -211;
517 pdgN = 211;
518 }
519
520 return kTRUE;
521}
522//____________________________________________________________________
523Bool_t AliESDv0KineCuts::CaseLambda(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN, Int_t id) const {
524 //
525 // process teh Lambda and Anti-Lambda candidate
526 //
527
528 if(!v0) return kFALSE;
529
530 const Double_t cL0mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(); // PDG lambda mass
531
532 AliVTrack* daughter[2];
533 Int_t pIndex = 0, nIndex = 0;
534 Float_t mMass[2] = {-1., -1.};
535 Bool_t sign = CheckSigns(v0);
536 if(sign){
537 pIndex = v0->GetPindex();
538 nIndex = v0->GetNindex();
539 mMass[0] = v0->GetEffMass(4, 2);
540 mMass[1] = v0->GetEffMass(2, 4);
541 }
542 else{
543 pIndex = v0->GetNindex();
544 nIndex = v0->GetPindex();
545 mMass[0] = v0->GetEffMass(2, 4);
546 mMass[1] = v0->GetEffMass(4, 2);
547 }
548
549 daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
550 daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
551 if(!daughter[0] || !daughter[1]) return kFALSE;
552
553 AliKFParticle *kfMother[2] = {0x0, 0x0};
554 // Lambda
555 kfMother[0] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kProton), TMath::Abs(kPiPlus));
556 if(!kfMother[0]) return kFALSE;
557
558 // Anti-Lambda
559 kfMother[1] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kProton));
560 if(!kfMother[1]) return kFALSE;
561
2942f542 562 Float_t dMass[2] = {static_cast<Float_t>(TMath::Abs(mMass[0] - cL0mass)), static_cast<Float_t>(TMath::Abs(mMass[1] - cL0mass))};
a65a7e70 563
564 AliESDtrack* d[2];
565 d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
566 d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
567 if(!d[0] || !d[1]) return kFALSE;
568
2942f542 569 Float_t p[2] = {static_cast<Float_t>(d[0]->GetP()), static_cast<Float_t>(d[1]->GetP())};
a65a7e70 570
571 // check the 3 lambda - antilambda variables
572 Int_t check[2] = {-1, -1}; // 0 : lambda, 1 : antilambda
573 // 1) momentum of the daughter particles - proton is expected to have higher momentum than pion
574 check[0] = (p[0] > p[1]) ? 0 : 1;
575 // 2) mass of the mother particle
576 check[1] = (dMass[0] < dMass[1]) ? 0 : 1;
577
578 // require positive correlation of (1) and (2)
579 if(check[0] != check[1]){
580 if(kfMother[0]) delete kfMother[0];
581 if(kfMother[1]) delete kfMother[1];
582 return kFALSE;
583 }
584
585 // now that the check[0] == check[1]
586 const Int_t type = check[0];
587
588 // require that the input armenteros preselection agree:
589 if(type != id) return kFALSE;
590
591 Float_t iMass =0.;
592 if(sign){
593 iMass = (type == 0) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);
594 }
595 else{
596 iMass = (type == 0) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
597 }
598
599 // cos pointing angle
600 Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
601 cosPoint = TMath::ACos(cosPoint);
602
603 // DCA between daughters
604 Double_t dca = v0->GetDcaV0Daughters();
605
606 // Production vertex
607 Double_t x, y, z;
608 v0->GetXYZ(x,y,z);
609 Double_t r = TMath::Sqrt(x*x + y*y);
610
611 // proton - pion indices
612 Int_t ix[2] = {0, 1};
613 if(1 == type){
614 ix[0] = 1;
615 ix[1] = 0;
616 }
617
618 // V0 chi2/ndf
619 Double_t chi2ndf = kfMother[type]->GetChi2()/kfMother[type]->GetNDF();
620
621 if(kfMother[0]) delete kfMother[0];
622 if(kfMother[1]) delete kfMother[1];
623
624 //
625 // apply the cuts
626 //
627
628 if(iMass < fLcutInvMass[0] || iMass > fLcutInvMass[1]) return kFALSE;
629
630 if(chi2ndf > fLcutChi2NDF) return kFALSE;
631
632 if(cosPoint < fLcutCosPoint[0] || cosPoint > fLcutCosPoint[1]) return kFALSE;
633
634 if(dca < fLcutDCA[0] || dca > fLcutDCA[1]) return kFALSE;
635
636 if(r < fLcutVertexR[0] || r > fLcutVertexR[1]) return kFALSE;
637
638 // all cuts passed
639
640 if(0 == type){
641 pdgV0 = 3122;
642 if(sign){
643 pdgP = 2212;
644 pdgN = -211;
645 }
646 else{
647 pdgP = -211;
648 pdgN = 2212;
649 }
650 }
651 else{
652 pdgV0 = -3122;
653 if(sign){
654 pdgP = 211;
655 pdgN = -2212;
656 }
657 else{
658 pdgP = -2212;
659 pdgN = 211;
660 }
661 }
662
663 return kTRUE;
664}
665//____________________________________________________________________
666Bool_t AliESDv0KineCuts::V0CutsCommon(const AliESDv0 * const v0) const
667{
668 //
669 // V0 cuts common to all V0s
670 //
671
672 AliESDtrack* dN, *dP;
673
674 dP = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetPindex()));
675 dN = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetNindex()));
676
677 if(!dN || !dP) return kFALSE;
678
679 Int_t qP = dP->Charge();
680 Int_t qN = dN->Charge();
681
682 if((qP*qN) != -1) return kFALSE;
683
684 return kTRUE;
685}
686//____________________________________________________________________
687void AliESDv0KineCuts::Armenteros(AliESDv0* const v0, Float_t val[2]) const
688{
689 //
690 // computes the Armenteros variables for given V0
691 // fills the histogram
692 // returns the values via "val"
693 //
694
695 Double_t mn[3] = {0,0,0};
696 Double_t mp[3] = {0,0,0};
697 Double_t mm[3] = {0,0,0};
698
699 if(CheckSigns(v0)){
700 v0->GetNPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
701 v0->GetPPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
702 }
703 else{
704 v0->GetPPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
705 v0->GetNPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
706 }
707 v0->GetPxPyPz(mm[0],mm[1],mm[2]); //reconstructed cartesian momentum components of mother
708
709 TVector3 vecN(mn[0],mn[1],mn[2]);
710 TVector3 vecP(mp[0],mp[1],mp[2]);
711 TVector3 vecM(mm[0],mm[1],mm[2]);
712
713 Double_t thetaP = acos((vecP * vecM)/(vecP.Mag() * vecM.Mag()));
714 Double_t thetaN = acos((vecN * vecM)/(vecN.Mag() * vecM.Mag()));
715
716 Double_t alfa = ((vecP.Mag())*cos(thetaP)-(vecN.Mag())*cos(thetaN))/
717 ((vecP.Mag())*cos(thetaP)+(vecN.Mag())*cos(thetaN)) ;
718 Double_t qt = vecP.Mag()*sin(thetaP);
719
720 val[0] = alfa;
721 val[1] = qt;
722}
723//____________________________________________________________________
724Bool_t AliESDv0KineCuts::CheckSigns(AliESDv0* const v0) const
725{
726 //
727 // check wheter the sign was correctly applied to
728 // V0 daughter tracks
729 //
730
731 Bool_t correct = kFALSE;
732
733 Int_t pIndex = 0, nIndex = 0;
734 pIndex = v0->GetPindex();
735 nIndex = v0->GetNindex();
736
737 AliESDtrack* d[2];
738 d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
739 d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
740
741 Int_t sign[2];
742 sign[0] = (int)d[0]->GetSign();
743 sign[1] = (int)d[1]->GetSign();
744
745 if(-1 == sign[0] && 1 == sign[1]){
746 correct = kFALSE;
747 }
748 else{
749 correct = kTRUE;
750 }
751
752 return correct;
753}
754//________________________________________________________________
755Double_t AliESDv0KineCuts::PsiPair(AliESDv0* const v0) const
756{
757 //
758 // Angle between daughter momentum plane and plane
759 //
760
761 if(!fEvent) return -1.;
762
763 Float_t magField = fEvent->GetMagneticField();
764
765 Int_t pIndex = -1;
766 Int_t nIndex = -1;
767 if(CheckSigns(v0)){
768 pIndex = v0->GetPindex();
769 nIndex = v0->GetNindex();
770 }
771 else{
772 pIndex = v0->GetNindex();
773 nIndex = v0->GetPindex();
774 }
775
776
777 AliESDtrack* daughter[2];
778
779 daughter[0] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(pIndex));
780 daughter[1] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(nIndex));
781
782 Double_t x, y, z;
783 v0->GetXYZ(x,y,z);//Reconstructed coordinates of V0; to be replaced by Markus Rammler's method in case of conversions!
784
785 Double_t mn[3] = {0,0,0};
786 Double_t mp[3] = {0,0,0};
787
788
789 v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
790 v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
791
792
793 Double_t deltat = 1.;
794 deltat = TMath::ATan(mp[2]/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1])+1.e-13)) - TMath::ATan(mn[2]/(TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1])+1.e-13));//difference of angles of the two daughter tracks with z-axis
795
796 Double_t radiussum = TMath::Sqrt(x*x + y*y) + 50;//radius to which tracks shall be propagated
797
798 Double_t momPosProp[3];
799 Double_t momNegProp[3];
800
801 AliExternalTrackParam pt(*daughter[0]), nt(*daughter[1]);
802
803 Double_t psiPair = 4.;
804
805 if(nt.PropagateTo(radiussum,magField) == 0)//propagate tracks to the outside
806 psiPair = -5.;
807 if(pt.PropagateTo(radiussum,magField) == 0)
808 psiPair = -5.;
809 pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
810 nt.GetPxPyPz(momNegProp);
811
812 Double_t pEle =
813 TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
814 Double_t pPos =
815 TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
816
817 Double_t scalarproduct =
818 momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
819
820 Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
821
822 psiPair = TMath::Abs(TMath::ASin(deltat/chipair));
823
824 return psiPair;
825}
826//___________________________________________________________________
827Bool_t AliESDv0KineCuts::GetConvPosXY(AliESDtrack * const ptrack, AliESDtrack * const ntrack, Double_t convpos[2]) const
828{
829 //
830 // recalculate the gamma conversion XY postition
831 //
832
833 const Double_t b = fEvent->GetMagneticField();
834
835 Double_t helixcenterpos[2];
836 GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);
837
838 Double_t helixcenterneg[2];
839 GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);
840
841 Double_t poshelix[6];
842 ptrack->GetHelixParameters(poshelix,b);
843 Double_t posradius = TMath::Abs(1./poshelix[4]);
844
845 Double_t neghelix[6];
846 ntrack->GetHelixParameters(neghelix,b);
847 Double_t negradius = TMath::Abs(1./neghelix[4]);
848
849 Double_t xpos = helixcenterpos[0];
850 Double_t ypos = helixcenterpos[1];
851 Double_t xneg = helixcenterneg[0];
852 Double_t yneg = helixcenterneg[1];
853
854 convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
855 convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);
856
857 return 1;
858}
859//___________________________________________________________________
860Bool_t AliESDv0KineCuts::GetHelixCenter(AliESDtrack * const track, Double_t b,Int_t charge, Double_t center[2]) const
861{
862 //
863 // computes the center of the track helix
864 //
865
866 Double_t pi = TMath::Pi();
867
868 Double_t helix[6];
869 track->GetHelixParameters(helix,b);
870
871 Double_t xpos = helix[5];
872 Double_t ypos = helix[0];
873 Double_t radius = TMath::Abs(1./helix[4]);
874 Double_t phi = helix[2];
875
876 if(phi < 0){
877 phi = phi + 2*pi;
878 }
879
880 phi -= pi/2.;
881 Double_t xpoint = radius * TMath::Cos(phi);
882 Double_t ypoint = radius * TMath::Sin(phi);
883
884 if(b<0){
885 if(charge > 0){
886 xpoint = - xpoint;
887 ypoint = - ypoint;
888 }
889 /* avoid self assignment
890 if(charge < 0){
891 xpoint = xpoint;
892 ypoint = ypoint;
893 }
894 */
895 }
896 if(b>0){
897 /* avoid self assignment
898 if(charge > 0){
899 xpoint = xpoint;
900 ypoint = ypoint;
901 }
902 */
903 if(charge < 0){
904 xpoint = - xpoint;
905 ypoint = - ypoint;
906 }
907 }
908 center[0] = xpos + xpoint;
909 center[1] = ypos + ypoint;
910
911 return 1;
912}
913//___________________________________________________________________
914AliKFParticle *AliESDv0KineCuts::CreateMotherParticle(const AliVTrack* const pdaughter, const AliVTrack* const ndaughter, Int_t pspec, Int_t nspec) const
915{
916 //
917 // Creates a mother particle
918 //
919 AliKFParticle pkfdaughter(*pdaughter, pspec);
920 AliKFParticle nkfdaughter(*ndaughter, nspec);
921
922
923 // Create the mother particle
924 AliKFParticle *m = new AliKFParticle(pkfdaughter, nkfdaughter);
925 m->SetField(fEvent->GetMagneticField());
926 if(TMath::Abs(kElectron) == pspec && TMath::Abs(kElectron) == nspec) m->SetMassConstraint(0, 0.001);
927 else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(), 0.);
928 else if(TMath::Abs(kProton) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
929 else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kProton) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
930 else{
931 AliErrorClass("Wrong daughter ID - mass constraint can not be set");
932 }
933
934 AliKFVertex improvedVertex = *fPrimaryVertex;
935 improvedVertex += *m;
936 m->SetProductionVertex(improvedVertex);
937
938 // update 15/06/2010
939 // mother particle will not be added to primary vertex but only to its copy
940 // as this confilcts with calling
941 // m->SetPrimaryVertex() function and
942 // subsequently removing the mother particle afterwards
943 // Source: Sergey Gorbunov
944
945 return m;
946}
947//____________________________________________________________________
948void AliESDv0KineCuts::SetEvent(AliESDEvent* const event){
949 //
950 // direct setter of ESD event
951 //
952 fEvent = event;
953 if(!fEvent){
954 AliErrorClass("Invalid input event pointer");
955 return;
956 }
957if (fUseExternalVertex) return;
958else{
959 if(fPrimaryVertex && fDeleteVertex){
960 delete fPrimaryVertex;
961 fPrimaryVertex=0x0;
962 }
963 fPrimaryVertex = new AliKFVertex(*(fEvent->GetPrimaryVertex()));
964 fDeleteVertex=kTRUE;
965 }
966
967
968
969}
970//____________________________________________________________________
971void AliESDv0KineCuts::SetEvent(AliVEvent* const event){
972 //
973 // direct setter of ESD event
974 //
975
976 fEvent = dynamic_cast<AliESDEvent*>(event);
977 if(!fEvent){
978 AliErrorClass("Invalid input event pointer");
979 return;
980 }
981
982 if (fUseExternalVertex) return;
983 else{
984 if(fPrimaryVertex && fDeleteVertex){
985 delete fPrimaryVertex;
986 fPrimaryVertex=0x0;
987 }
988 fPrimaryVertex = new AliKFVertex(*(fEvent->GetPrimaryVertex()));
989 fDeleteVertex=kTRUE;
990 }
991}
992
993
994//________________________________________________________________
995void AliESDv0KineCuts::UseExternalVertex(Bool_t use_external){
996 //
997 // Reenable primary Vertex from ESD event
998 //
999 if (use_external) fUseExternalVertex =kTRUE;
1000 else fUseExternalVertex =kFALSE;
1001}
1002
1003
1004
1005
1006//________________________________________________________________
1007void AliESDv0KineCuts::SetPrimaryVertex(AliKFVertex* const v){
1008 //
1009 // set the primary vertex of the event
1010 //
1011 if(fPrimaryVertex && fDeleteVertex){
1012 delete fPrimaryVertex;
1013 fPrimaryVertex =0x0;
1014 fDeleteVertex = kFALSE;
1015 }
1016 fUseExternalVertex=kTRUE;
1017 fPrimaryVertex = v; // set primary Vertex
1018 if(!fPrimaryVertex){
1019 AliErrorClass("Failed to initialize the primary vertex");
1020 return;
1021 }
1022}
1023//___________________________________________________________________
1024void AliESDv0KineCuts::SetMode(Int_t mode, Int_t type){
1025 //
1026 // this function allows the user to select (prior running the 'ProcessV0' function)
1027 // to select different approaches to V0 selection - the 'mode'
1028 // - and -
1029 // different systems (pp, PbPb) - 'type'
1030 //
1031 // To see the cut values for different modes please refer to the
1032 // function SetCuts()
1033 //
1034 // Important notice: based on the parameters particular sets of cuts will
1035 // be activated for teh V0 selection. If some additional changes to single
1036 // cuts are needed please us the SetXXXcut function (see the header file)
1037 //
1038
1039 switch(mode){
1040 case kPurity:
1041 fMode = kPurity; // used to obtain highest purity possible - the efficiency may be low
1042 break;
1043 case kEffGamma:
1044 fMode = kEffGamma; // used to obtain highes efficiency possible - the purity may be worse
1045 break;
1046 default:
1047 AliError("V0 selection mode not recognozed, setting 'kPurity'");
1048 fMode = kPurity;
1049 }
1050
1051 switch(type){
1052 case kPP:
1053 fType = kPP; // cuts optimized for low multiplicity
1054 break;
1055 case kPbPb:
1056 fType = kPbPb; // cuts optimized for high multiplicity
1057 break;
1058 }
1059
1060 // setup the cut values for selected mode & type
1061 SetCuts();
1062
1063}
1064//___________________________________________________________________
1065void AliESDv0KineCuts::SetMode(Int_t mode, const char* type){
1066 //
1067 // overloaded function - please see above
1068 //
1069
1070 Int_t t = -1;
1071
1072 if(!strcmp("pp", type)) t = kPP;
1073 else if(!(strcmp("PbPb", type))) t = kPbPb;
1074 else{
1075 AliError("data type not recognized, setting 'pp'");
1076 t = kPP;
1077 }
1078
1079 SetMode(mode, t);
1080
1081}
1082//___________________________________________________________________
1083void AliESDv0KineCuts::SetCuts(){
1084 //
1085 // this funciton sets the default cut values based on the selected
1086 // fMode and fType.
1087 // please note that only the cuts that have different values than the default
1088 // cuts are updated here
1089 //
1090
1091 // last update: 14/02/2011
1092 // as a very preliminary - the only change to default cuts is to apply
1093 // less restricting gamma conversion selection in PreselectV0() function
1094
1095
1096
1097}