AliPhysics  a0db429 (a0db429)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros
AliRDHFCutsLctopK0sfromAODtracks.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2010, 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 
19 //
20 // Class for cuts on AOD reconstructed Lc->p+K0s
21 //
22 // Modified by Y.S Watanabe - wyosuke@cns.s.u-tokyo.ac.jp
23 //
25 
26 #include <Riostream.h>
27 
28 #include <TDatabasePDG.h>
29 #include <TMath.h>
30 #include <TLorentzVector.h>
31 #include <TVector3.h>
32 
33 #include "AliAnalysisManager.h"
34 #include "AliInputEventHandler.h"
35 #include "AliPIDResponse.h"
37 #include "AliAODRecoCascadeHF.h"
38 #include "AliAODTrack.h"
39 #include "AliESDtrack.h"
40 #include "AliESDVertex.h"
41 #include "AliAODVertex.h"
42 #include "AliAODv0.h"
43 #include "AliAODcascade.h"
44 #include "AliESDv0.h"
45 
46 using std::cout;
47 using std::endl;
48 
52 
53 
54 //--------------------------------------------------------------------------
56 AliRDHFCuts(name),
57  fPIDStrategy(kNSigmaCuts),
58  fCombinedPIDThreshold(0.),
59  fUseOnTheFlyV0(kFALSE),
60  fProdTrackTPCNclsPIDMin(0),
61  fProdTrackTPCNclsRatioMin(0.0),
62  fProdUseAODFilterBit(kTRUE),
63  fProdV0MassTolK0s(0.01),
64  fProdV0MassRejLambda(0.00),
65  fProdV0MassRejPhoton(0.00),
66  fProdV0PtMin(0.5),
67  fProdV0CosPointingAngleToPrimVtxMin(0.99),
68  fProdV0DcaDaughtersMax(1.5),
69  fProdV0DaughterEtaRange(0.8),
70  fProdV0DaughterPtMin(0.0),
71  fProdV0DaughterTPCClusterMin(70),
72  fProdV0EtaMin(-9999.),
73  fProdV0EtaMax(9999.),
74  fProdV0RapMin(-9999.),
75  fProdV0RapMax(9999.),
76  fProdRoughMassTol(0.25),
77  fProdRoughPtMin(0.0),
78  fNWeightingProtonBinLimits(0),
79  fWeightingProtonBins(0),
80  fNWeightingK0sBinLimits(0),
81  fWeightingK0sBins(0),
82  fNWeightingBins(0),
83  fWeight_p0(0),
84  fWeight_p1(0),
85  fWeight_p2(0),
86  fWeight_p3(0),
87  fTagV0MassTol(0)
88 {
89  //
90  // Default Constructor
91  //
92 
93  const Int_t nvars=9;
94  SetNVars(nvars);
95  TString varNames[nvars]={"Lc inv. mass [GeV/c2]", // 0
96  "Opening angle [rad]",//1
97  "n#sigma_{TPC} max",//2
98  "n#sigma_{TOF} min",//3
99  "Decay Length min [cm]",//4
100  "cos #theta min",//5
101  "cos #theta max",//6
102  "Proton d0 max",//7
103  "V0 d0 max"//8
104  };
105 
106  Bool_t isUpperCut[nvars]={kTRUE, // 0
107  kFALSE, //1: opening angle
108  kTRUE, //2: nsigma_tpc max
109  kFALSE, //3: nsigma tof min
110  kFALSE, //4: decay length min
111  kFALSE, //5: cos the min
112  kTRUE, //6: cos the max
113  kTRUE, //7:
114  kTRUE //8:
115  };
116  SetVarNames(nvars,varNames,isUpperCut);
117  Bool_t forOpt[nvars]={kFALSE, // 0
118  kFALSE, //1
119  kTRUE, //2
120  kTRUE, //3
121  kTRUE, //4
122  kTRUE, //5
123  kTRUE, //6
124  kTRUE, //7
125  kTRUE //8
126  };
127  SetVarsForOpt(nvars,forOpt);
128 
129  Float_t limits[2]={0,999999999.};
130  SetPtBins(2,limits);
131 }
132 //--------------------------------------------------------------------------
134  AliRDHFCuts(source),
135  fPIDStrategy(source.fPIDStrategy),
136  fCombinedPIDThreshold(source.fCombinedPIDThreshold),
137  fUseOnTheFlyV0(source.fUseOnTheFlyV0),
138  fProdTrackTPCNclsPIDMin(source.fProdTrackTPCNclsPIDMin),
139  fProdTrackTPCNclsRatioMin(source.fProdTrackTPCNclsRatioMin),
140  fProdUseAODFilterBit(source.fProdUseAODFilterBit),
141  fProdV0MassTolK0s(source.fProdV0MassTolK0s),
142  fProdV0MassRejLambda(source.fProdV0MassRejLambda),
143  fProdV0MassRejPhoton(source.fProdV0MassRejPhoton),
144  fProdV0PtMin(source.fProdV0PtMin),
145  fProdV0CosPointingAngleToPrimVtxMin(source.fProdV0CosPointingAngleToPrimVtxMin),
146  fProdV0DcaDaughtersMax(source.fProdV0DcaDaughtersMax),
147  fProdV0DaughterEtaRange(source.fProdV0DaughterEtaRange),
148  fProdV0DaughterPtMin(source.fProdV0DaughterPtMin),
149  fProdV0DaughterTPCClusterMin(source.fProdV0DaughterTPCClusterMin),
150  fProdV0EtaMin(source.fProdV0EtaMin),
151  fProdV0EtaMax(source.fProdV0EtaMax),
152  fProdV0RapMin(source.fProdV0RapMin),
153  fProdV0RapMax(source.fProdV0RapMax),
154  fProdRoughMassTol(source.fProdRoughMassTol),
155  fProdRoughPtMin(source.fProdRoughPtMin),
156  fNWeightingProtonBinLimits(source.fNWeightingProtonBinLimits),
157  fWeightingProtonBins(NULL),
158  fNWeightingK0sBinLimits(source.fNWeightingK0sBinLimits),
159  fWeightingK0sBins(NULL),
160  fNWeightingBins(source.fNWeightingBins),
161  fWeight_p0(NULL),
162  fWeight_p1(NULL),
163  fWeight_p2(NULL),
164  fWeight_p3(NULL),
165  fTagV0MassTol(source.fTagV0MassTol)
166 {
167  //
168  // Copy constructor
169  //
172  fWeight_p0 = new Double_t[fNWeightingBins];
173  fWeight_p1 = new Double_t[fNWeightingBins];
174  fWeight_p2 = new Double_t[fNWeightingBins];
175  fWeight_p3 = new Double_t[fNWeightingBins];
176  for(Int_t i=0;i<fNWeightingProtonBinLimits;i++){
178  }
179  for(Int_t i=0;i<fNWeightingK0sBinLimits;i++){
180  fWeightingK0sBins[i] = source.fWeightingK0sBins[i];
181  }
182  for(Int_t i=0;i<fNWeightingBins;i++){
183  fWeight_p0[i] = source.fWeight_p0[i];
184  fWeight_p1[i] = source.fWeight_p1[i];
185  fWeight_p2[i] = source.fWeight_p2[i];
186  fWeight_p3[i] = source.fWeight_p3[i];
187  }
188 }
189 //--------------------------------------------------------------------------
191 {
192  //
193  // assignment operator
194  //
195 
196  if (this != &source) {
197  AliRDHFCuts::operator=(source);
198  }
199 
200  fPIDStrategy = source.fPIDStrategy;
209  fProdV0PtMin = source.fProdV0PtMin;
215  fProdV0EtaMin = source.fProdV0EtaMin;
216  fProdV0EtaMax = source.fProdV0EtaMax;
217  fProdV0RapMin = source.fProdV0RapMin;
218  fProdV0RapMax = source.fProdV0RapMax;
224  fTagV0MassTol = source.fTagV0MassTol;
225 
228  fWeight_p0 = new Double_t[fNWeightingBins];
229  fWeight_p1 = new Double_t[fNWeightingBins];
230  fWeight_p2 = new Double_t[fNWeightingBins];
231  fWeight_p3 = new Double_t[fNWeightingBins];
232  for(Int_t i=0;i<fNWeightingProtonBinLimits;i++){
234  }
235  for(Int_t i=0;i<fNWeightingK0sBinLimits;i++){
236  fWeightingK0sBins[i] = source.fWeightingK0sBins[i];
237  }
238  for(Int_t i=0;i<fNWeightingBins;i++){
239  fWeight_p0[i] = source.fWeight_p0[i];
240  fWeight_p1[i] = source.fWeight_p1[i];
241  fWeight_p2[i] = source.fWeight_p2[i];
242  fWeight_p3[i] = source.fWeight_p3[i];
243  }
244 
245 
246  return *this;
247 }
248 
249 //---------------------------------------------------------------------------
251  //
252  // Default Destructor
253  //
254 
255  delete [] fWeightingProtonBins;
256  delete [] fWeightingK0sBins;
257  delete [] fWeight_p0;
258  delete [] fWeight_p1;
259  delete [] fWeight_p2;
260  delete [] fWeight_p3;
261 
262 }
263 
264 //---------------------------------------------------------------------------
265 void AliRDHFCutsLctopK0sfromAODtracks::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
266  //
267  // Fills in vars the values of the variables
268  //
269 
270  if (pdgdaughters[0]==-9999) return; // dummy
271 
273  if(!dd){
274  AliError("No AliAODRecoCascadeHF object found\n");
275  return;
276  }
277 
278  if (nvars!=fnVarsForOpt) {
279  AliError("Cut object seems to have the wrong number of variables\n");
280  return;
281  }
282 
283  //Double_t ptD=d->Pt();
284  //Int_t ptbin=PtBin(ptD);
285  Int_t iter=-1;
286 
287  if(fVarsForOpt[0]){
288  iter++;
289  Double_t mlcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
290  vars[iter]= TMath::Abs(dd->InvMassLctoK0sP()-mlcPDG) ;
291  }
292  if(fVarsForOpt[1]){
293  iter++;
294  vars[iter]= dd->Pt();
295  }
296 
297  return;
298 }
299 
300 //---------------------------------------------------------------------------
301 Int_t AliRDHFCutsLctopK0sfromAODtracks::IsSelected(TObject* obj, Int_t selectionLevel) {
302  //
303  // Apply selection
304  //
305 
306  if (!fCutsRD) {
307  AliFatal("Cut matrix not inizialized. Exit...");
308  return 0;
309  }
310 
312  if(!d){
313  AliDebug(2,"AliAODRecoCascadeHF null");
314  return 0;
315  }
316 
317  Double_t ptD=d->Pt();
318  if(ptD<fMinPtCand) return 0;
319  if(ptD>fMaxPtCand) return 0;
320 
321  Double_t pt=d->Pt();
322  Int_t ptbin=PtBin(pt);
323  if (ptbin==-1) {
324  return 0;
325  }
326 
327  if (selectionLevel==AliRDHFCuts::kAll ||
328  selectionLevel==AliRDHFCuts::kTracks) {
329  //Performed in production stage
330  }
331 
332  Int_t returnvalueCuts=1;
333  // selection on candidate
334  if (selectionLevel==AliRDHFCuts::kAll ||
335  selectionLevel==AliRDHFCuts::kCandidate) {
336 
337  Bool_t okcand=kTRUE;
338 
339  Double_t mlcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
340  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
341  Double_t mk0PDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
342  Double_t v0px = d->PxProng(1);
343  Double_t v0py = d->PyProng(1);
344  Double_t v0pz = d->PzProng(1);
345  Double_t epx = d->PxProng(0);
346  Double_t epy = d->PyProng(0);
347  Double_t epz = d->PzProng(0);
348  Double_t cosoa = (v0px*epx+v0py*epy+v0pz*epz)/sqrt(v0px*v0px+v0py*v0py+v0pz*v0pz)/sqrt(epx*epx+epy*epy+epz*epz);
349 
350  TLorentzVector vpr, vk0s,vlc;
351  vpr.SetXYZM(epx,epy,epz,mprPDG);
352  vk0s.SetXYZM(v0px,v0py,v0pz,mk0PDG);
353  vlc = vpr + vk0s;
354  TVector3 vboost = vlc.BoostVector();
355  vpr.Boost(-vboost);
356  Double_t bachcosthe = cos(vpr.Angle(vlc.Vect()));
357 
358  if(TMath::Abs(d->InvMassLctoK0sP()-mlcPDG) > fCutsRD[GetGlobalIndex(0,ptbin)])
359  {
360  okcand = kFALSE;
361  }
362  if(cosoa < fCutsRD[GetGlobalIndex(1,ptbin)])
363  {
364  okcand = kFALSE;
365  }
366  if(d->DecayLengthXY() < fCutsRD[GetGlobalIndex(4,ptbin)])
367  {
368  okcand = kFALSE;
369  }
370  if(bachcosthe < fCutsRD[GetGlobalIndex(5,ptbin)])
371  {
372  okcand = kFALSE;
373  }
374  if(bachcosthe > fCutsRD[GetGlobalIndex(6,ptbin)])
375  {
376  okcand = kFALSE;
377  }
378  if(fabs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(7,ptbin)])
379  {
380  okcand = kFALSE;
381  }
382  if(fabs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(8,ptbin)])
383  {
384  okcand = kFALSE;
385  }
386 
387  if(!okcand) return 0;
388  returnvalueCuts = 1;
389  }
390 
391  Int_t returnvaluePID=1;
392  if(selectionLevel==AliRDHFCuts::kAll ||
393  selectionLevel==AliRDHFCuts::kCandidate||
394  selectionLevel==AliRDHFCuts::kPID) {
395 
396  AliAODTrack *part = (AliAODTrack*)d->GetSecondaryVtx()->GetDaughter(0);
397 
398  Double_t nSigmaTPCpr = fPidHF->GetPidResponse()->NumberOfSigmasTPC(part,AliPID::kProton);
399  if(nSigmaTPCpr>fCutsRD[GetGlobalIndex(2,ptbin)]){
400  returnvaluePID = -1;
401  }
402 
403  Double_t nSigmaTOFpr = fPidHF->GetPidResponse()->NumberOfSigmasTOF(part,AliPID::kProton);
404  if(nSigmaTOFpr<fCutsRD[GetGlobalIndex(3,ptbin)]){
405  returnvaluePID = -1;
406  }
407 
408  }
409 
410  Int_t returnvalue = 0;
411  if(returnvalueCuts==1 && returnvaluePID==1) returnvalue=1;
412 
413  return returnvalue;
414 }
415 
416 //---------------------------------------------------------------------------
417 Int_t AliRDHFCutsLctopK0sfromAODtracks::IsSelected(TLorentzVector* vtrk, TLorentzVector *vv0, Double_t *cutvars, Int_t selectionLevel) {
418  //
419  // Apply selection on mixed event tracks
420  //
421 
422  if (!fCutsRD) {
423  AliFatal("Cut matrix not inizialized. Exit...");
424  return 0;
425  }
426 
427  Double_t ptD=cutvars[1];
428  if(ptD<fMinPtCand) return 0;
429  if(ptD>fMaxPtCand) return 0;
430 
431  Double_t pt=cutvars[1];
432  Int_t ptbin=PtBin(pt);
433  if (ptbin==-1) {
434  return 0;
435  }
436 
437  if (selectionLevel==AliRDHFCuts::kAll ||
438  selectionLevel==AliRDHFCuts::kTracks) {
439  //Performed in production stage
440  }
441 
442  Int_t returnvalueCuts=1;
443  // selection on candidate
444  if (selectionLevel==AliRDHFCuts::kAll ||
445  selectionLevel==AliRDHFCuts::kCandidate) {
446 
447  Bool_t okcand=kTRUE;
448 
449  Double_t mlcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
450  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
451  Double_t mk0PDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
452  Double_t v0px = vv0->Px();
453  Double_t v0py = vv0->Py();
454  Double_t v0pz = vv0->Pz();
455  Double_t epx = vtrk->Px();
456  Double_t epy = vtrk->Py();
457  Double_t epz = vtrk->Pz();
458  Double_t cosoa = (v0px*epx+v0py*epy+v0pz*epz)/sqrt(v0px*v0px+v0py*v0py+v0pz*v0pz)/sqrt(epx*epx+epy*epy+epz*epz);
459 
460  TLorentzVector vpr, vk0s,vlc;
461  vpr.SetXYZM(epx,epy,epz,mprPDG);
462  vk0s.SetXYZM(v0px,v0py,v0pz,mk0PDG);
463  vlc = vpr + vk0s;
464  TVector3 vboost = vlc.BoostVector();
465  vpr.Boost(-vboost);
466  Double_t bachcosthe = cos(vpr.Angle(vlc.Vect()));
467 
468  if(TMath::Abs(cutvars[0]-mlcPDG) > fCutsRD[GetGlobalIndex(0,ptbin)])
469  {
470  okcand = kFALSE;
471  }
472  if(cosoa < fCutsRD[GetGlobalIndex(1,ptbin)])
473  {
474  okcand = kFALSE;
475  }
476  if(cutvars[4] < fCutsRD[GetGlobalIndex(4,ptbin)])
477  {
478  okcand = kFALSE;
479  }
480  if(bachcosthe < fCutsRD[GetGlobalIndex(5,ptbin)])
481  {
482  okcand = kFALSE;
483  }
484  if(bachcosthe > fCutsRD[GetGlobalIndex(6,ptbin)])
485  {
486  okcand = kFALSE;
487  }
488  if(cutvars[5] > fCutsRD[GetGlobalIndex(7,ptbin)])
489  {
490  okcand = kFALSE;
491  }
492  if(cutvars[6] > fCutsRD[GetGlobalIndex(8,ptbin)])
493  {
494  okcand = kFALSE;
495  }
496 
497  if(!okcand) return 0;
498  returnvalueCuts = 1;
499  }
500 
501  Int_t returnvaluePID=1;
502  if(selectionLevel==AliRDHFCuts::kAll ||
503  selectionLevel==AliRDHFCuts::kCandidate||
504  selectionLevel==AliRDHFCuts::kPID) {
505 
506  Double_t nSigmaTPCpr = cutvars[2];
507  if(fabs(nSigmaTPCpr)>fCutsRD[GetGlobalIndex(2,ptbin)]){
508  returnvaluePID = -1;
509  }
510 
511  Double_t nSigmaTOFpr = cutvars[3];
512  if(fabs(nSigmaTOFpr)<fCutsRD[GetGlobalIndex(3,ptbin)]){
513  returnvaluePID = -1;
514  }
515 
516  }
517 
518  Int_t returnvalue = 0;
519  if(returnvalueCuts==1 && returnvaluePID==1) returnvalue=1;
520 
521  return returnvalue;
522 }
523 
524 //---------------------------------------------------------------------------
526 {
527  //
528  // IsSelectedPID ( not used)
529  //
530 
531  if(!fUsePID || !obj) return 1;
532 
534  AliAODTrack *part = dd->GetBachelor();
535 
536  Int_t returnvalue=1;
537 
538  if(fPidHF->GetPidResponse()==0x0){
539  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
540  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
541  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
542  fPidHF->SetPidResponse(pidResp);
543  }
544 
545  Int_t isProton=fPidHF->MakeRawPid(part,4);
546  if(isProton<1) returnvalue = 0;
547 
548  return returnvalue;
549 }
550 
551 //---------------------------------------------------------------------------
553 {
554  //
555  // IsSelectedProtonID
556  //
557 
558  if(fPidHF->GetPidResponse()==0x0){
559  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
560  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
561  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
562  fPidHF->SetPidResponse(pidResp);
563  }
564 
565  //Int_t isProton=fPidHF->MakeRawPid(part,4);
566  //if(isProton<1) return kFALSE;
567  Double_t nsigmatpc_proton = fPidHF->GetSigma(0);
568  Double_t nsigmatof_proton = fPidHF->GetSigma(3);
569  Double_t nSigmaTPCpr = fPidHF->GetPidResponse()->NumberOfSigmasTPC(part,AliPID::kProton);
570  Double_t nSigmaTOFpr = fPidHF->GetPidResponse()->NumberOfSigmasTOF(part,AliPID::kProton);
571 
572  if(fabs(nSigmaTPCpr)<nsigmatpc_proton && fabs(nSigmaTOFpr)<nsigmatof_proton){
573  return kTRUE;
574  }
575 
576  return kFALSE;
577 }
578 //---------------------------------------------------------------------------
580 {
581  //
582  // IsSelectedKaonID
583  //
584 
585  if(fPidHF->GetPidResponse()==0x0){
586  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
587  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
588  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
589  fPidHF->SetPidResponse(pidResp);
590  }
591 
592  //Int_t isKaon=fPidHF->MakeRawPid(part,3);
593  //if(isKaon<1) return kFALSE;
594  Double_t nsigmatpc_kaon = fPidHF->GetSigma(0);
595  Double_t nsigmatof_kaon = fPidHF->GetSigma(3);
596  Double_t nSigmaTPCka = fPidHF->GetPidResponse()->NumberOfSigmasTPC(part,AliPID::kKaon);
597  Double_t nSigmaTOFka = fPidHF->GetPidResponse()->NumberOfSigmasTOF(part,AliPID::kKaon);
598 
599  if(fabs(nSigmaTPCka)<nsigmatpc_kaon && fabs(nSigmaTOFka)<nsigmatof_kaon){
600  return kTRUE;
601  }
602 
603  return kFALSE;
604 }
605 
606 //---------------------------------------------------------------------------
608  //
609  // IsSelectedCombinedPID
610  //
611 
612  if(!fUsePID || !obj) {return 1;}
613 
615  AliAODTrack *part = dd->GetBachelor();
616  if(!part) return 0;
617 
618  Int_t returnvalue=1;
619  Double_t probProton = GetProtonProbabilityTPCTOF(part);
620  if(probProton<fCombinedPIDThreshold) returnvalue = 0;
621 
622  return returnvalue;
623 }
624 
625 //________________________________________________________________________
627 {
628  //
629  // Get Proton probability
630  //
631  if(!fPidHF->GetUseCombined()) return -9999.;
632  fPidHF->GetPidCombined()->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);
633  Double_t prob1[AliPID::kSPECIES];
634  UInt_t detUsed1 = fPidHF->GetPidCombined()->ComputeProbabilities(trk, fPidHF->GetPidResponse(), prob1);
635  if (detUsed1 != (UInt_t)fPidHF->GetPidCombined()->GetDetectorMask() )
636  {
637  fPidHF->GetPidCombined()->SetDetectorMask(AliPIDResponse::kDetTPC);
638  detUsed1 = fPidHF->GetPidCombined()->ComputeProbabilities(trk, fPidHF->GetPidResponse(), prob1);
639  }
640  return prob1[AliPID::kProton];
641 }
642 
643 //________________________________________________________________________
644 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SingleTrkCuts(AliAODTrack *trk, AliAODVertex *primVert)
645 {
646  //
647  // Single Track Cut to be applied before object creation
648  //
649 
650  if(trk->GetStatus()&AliESDtrack::kITSpureSA) return kFALSE;
651  if(!(trk->GetStatus()&AliESDtrack::kITSin)) return kFALSE;
652  if(fProdUseAODFilterBit && !trk->TestFilterMask(BIT(4))) return kFALSE;
653  Double_t pos[3]; primVert->GetXYZ(pos);
654  Double_t cov[6]; primVert->GetCovarianceMatrix(cov);
655  const AliESDVertex vESD(pos,cov,100.,100);
656  if(fTrackCuts&&!IsDaughterSelected(trk,&vESD,fTrackCuts)) return kFALSE;
657 
658  if(trk->GetTPCsignalN()<fProdTrackTPCNclsPIDMin) return kFALSE;
659  if(trk->GetTPCNclsF()>0){
660  Float_t tpcratio = (Float_t)trk->GetTPCncls()/(Float_t)trk->GetTPCNclsF();
661  if(tpcratio<fProdTrackTPCNclsRatioMin) return kFALSE;
662  }
663 
664  if(fUsePID)
665  {
666  if(fPidHF->GetPidResponse()==0x0){
667  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
668  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
669  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
670  fPidHF->SetPidResponse(pidResp);
671  }
672 
673  switch(fPIDStrategy){
674  case kNSigmaCuts:
675  return IsSelectedProtonID(trk);
676  break;
677  }
678  }
679 
680  return kTRUE;
681 }
682 //________________________________________________________________________
683 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SingleKaonCuts(AliAODTrack *trk, AliAODVertex *primVert)
684 {
685  //
686  // Single Track Cut to be applied before object creation
687  //
688 
689  if(trk->GetStatus()&AliESDtrack::kITSpureSA) return kFALSE;
690  if(!(trk->GetStatus()&AliESDtrack::kITSin)) return kFALSE;
691  if(fProdUseAODFilterBit && !trk->TestFilterMask(BIT(4))) return kFALSE;
692  Double_t pos[3]; primVert->GetXYZ(pos);
693  Double_t cov[6]; primVert->GetCovarianceMatrix(cov);
694  const AliESDVertex vESD(pos,cov,100.,100);
695  if(fTrackCuts&&!IsDaughterSelected(trk,&vESD,fTrackCuts)) return kFALSE;
696 
697  if(trk->GetTPCsignalN()<fProdTrackTPCNclsPIDMin) return kFALSE;
698  if(trk->GetTPCNclsF()>0){
699  Float_t tpcratio = (Float_t)trk->GetTPCncls()/(Float_t)trk->GetTPCNclsF();
700  if(tpcratio<fProdTrackTPCNclsRatioMin) return kFALSE;
701  }
702 
703  if(fUsePID)
704  {
705  if(fPidHF->GetPidResponse()==0x0){
706  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
707  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
708  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
709  fPidHF->SetPidResponse(pidResp);
710  }
711 
712  switch(fPIDStrategy){
713  case kNSigmaCuts:
714  return IsSelectedKaonID(trk);
715  break;
716  }
717  }
718 
719  return kTRUE;
720 }
721 //________________________________________________________________________
723 {
725  if(!lcobj){
726  AliError("No AliAODRecoCascadeHF object found\n");
727  return -9999.;
728  }
729  Double_t dvertx = lcobj->GetSecondaryVtx()->GetX()-lcobj->GetOwnPrimaryVtx()->GetX();
730  Double_t dverty = lcobj->GetSecondaryVtx()->GetY()-lcobj->GetOwnPrimaryVtx()->GetY();
731  Double_t px = lcobj->Px();
732  Double_t py = lcobj->Py();
733  if(TMath::Sqrt(dvertx*dvertx+dverty*dverty)==0) return -9999.;
734  if(TMath::Sqrt(px*px+py*py)==0) return -9999.;
735  Double_t cospaxy = (dvertx*px+dverty*py)/TMath::Sqrt(dvertx*dvertx+dverty*dverty)/TMath::Sqrt(px*px+py*py);
736  return (Float_t)cospaxy;
737 }
738 //________________________________________________________________________
739 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SingleV0Cuts(AliAODv0 *v0, AliAODVertex *primVert)
740 {
741  //
742  // Single V0 Cut to be applied before object creation
743  //
744 
745  Bool_t onFlyV0 = v0->GetOnFlyStatus(); // on-the-flight V0s
746  if ( onFlyV0 && !fUseOnTheFlyV0 ) return kFALSE;
747 
748  AliAODTrack *cptrack = (AliAODTrack*)(v0->GetDaughter(0));
749  AliAODTrack *cntrack = (AliAODTrack*)(v0->GetDaughter(1));
750  if(!cptrack || !cntrack) return kFALSE;
751  if ( cptrack->Charge() == cntrack->Charge() ) return kFALSE;
752  if(!(cptrack->GetStatus() & AliESDtrack::kTPCrefit) ||
753  !(cntrack->GetStatus() & AliESDtrack::kTPCrefit)) return kFALSE;
754  AliAODVertex *maybeKinkPos = (AliAODVertex*)cptrack->GetProdVertex();
755  AliAODVertex *maybeKinkNeg = (AliAODVertex*)cntrack->GetProdVertex();
756  if (maybeKinkPos->GetType()==AliAODVertex::kKink || maybeKinkNeg->GetType()==AliAODVertex::kKink)
757  return kFALSE;
758 
759  if ( ( ( cptrack->GetTPCClusterInfo(2,1) ) < (Float_t)fProdV0DaughterTPCClusterMin ) ||
760  ( ( cntrack->GetTPCClusterInfo(2,1) ) < (Float_t)fProdV0DaughterTPCClusterMin) ) return kFALSE;
761 
762  Double_t massK0S = v0->MassK0Short();
763  Double_t mk0sPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
764  if(fabs(massK0S-mk0sPDG)>fProdV0MassTolK0s) return kFALSE;
765 
766  Double_t massLambda = v0->MassLambda();
767  Double_t massAntiLambda = v0->MassAntiLambda();
768  Double_t mlamPDG = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
769  if((fabs(massAntiLambda-mlamPDG)<fProdV0MassRejLambda) || (fabs(massLambda-mlamPDG)<fProdV0MassRejLambda)) return kFALSE;
770 
771  Double_t pxe1 = v0->MomPosX();
772  Double_t pye1 = v0->MomPosY();
773  Double_t pze1 = v0->MomPosZ();
774  Double_t Ee1 = sqrt(pxe1*pxe1+pye1*pye1+pze1*pze1+0.000511*0.000511);
775  Double_t pxe2 = v0->MomNegX();
776  Double_t pye2 = v0->MomNegY();
777  Double_t pze2 = v0->MomNegZ();
778  Double_t Ee2 = sqrt(pxe2*pxe2+pye2*pye2+pze2*pze2+0.000511*0.000511);
779  Double_t mphoton = sqrt(pow(Ee1+Ee2,2)-pow(pxe1+pxe2,2)-pow(pye1+pye2,2)-pow(pze1+pze2,2));
780  if(mphoton<fProdV0MassRejPhoton) return kFALSE;
781 
782 
783  if(TMath::Abs(v0->DcaV0Daughters())>fProdV0DcaDaughtersMax) return kFALSE;
784  Double_t posVtx[3] = {0.,0.,0.};
785  primVert->GetXYZ(posVtx);
786  Double_t cospav0 = v0->CosPointingAngle(posVtx);
787  if(cospav0<fProdV0CosPointingAngleToPrimVtxMin) return kFALSE;
788  if(v0->Pt()<fProdV0PtMin) return kFALSE;
789  if(fabs(cptrack->Eta())>fProdV0DaughterEtaRange) return kFALSE;
790  if(fabs(cntrack->Eta())>fProdV0DaughterEtaRange) return kFALSE;
791  if(cptrack->Pt()<fProdV0DaughterPtMin) return kFALSE;
792  if(cntrack->Pt()<fProdV0DaughterPtMin) return kFALSE;
793 
794  Double_t RapK0s = v0->RapK0Short();
795  if(RapK0s<fProdV0RapMin || RapK0s>fProdV0RapMax) return kFALSE;
796 
797  Double_t EtaK0s = v0->PseudoRapV0();
798  if(EtaK0s<fProdV0EtaMin || EtaK0s>fProdV0EtaMax) return kFALSE;
799 
800  return kTRUE;
801 }
802 
803 //________________________________________________________________________
804 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SelectWithRoughCuts(AliAODv0 *v0, AliAODTrack *part)
805 {
806  //
807  // Mass and pT Cut to be applied before object creation
808  //
809 
810  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
811  Double_t mLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
812 
813  Double_t pxpr_init = part->Px();
814  Double_t pypr_init = part->Py();
815  Double_t pzpr_init = part->Pz();
816  Double_t Epr_init = sqrt(pxpr_init*pxpr_init+pypr_init*pypr_init+pzpr_init*pzpr_init+mprPDG*mprPDG);
817 
818  Double_t pxv0_init = v0->Px();
819  Double_t pyv0_init = v0->Py();
820  Double_t pzv0_init = v0->Pz();
821  Double_t massv0_init = v0->MassK0Short();
822  Double_t Ev0_init = sqrt(pxv0_init*pxv0_init+pyv0_init*pyv0_init+pzv0_init*pzv0_init+massv0_init*massv0_init);
823 
824  Double_t pxlc_init = pxpr_init+pxv0_init;
825  Double_t pylc_init = pypr_init+pyv0_init;
826  Double_t pzlc_init = pzpr_init+pzv0_init;
827  Double_t Elc_init = Epr_init+Ev0_init;
828  Double_t lcmass_init = sqrt(Elc_init*Elc_init-pxlc_init*pxlc_init-pylc_init*pylc_init-pzlc_init*pzlc_init);
829 
830  if(lcmass_init<mLcPDG-fProdRoughMassTol || lcmass_init>mLcPDG+fProdRoughMassTol) return kFALSE;
831  if(sqrt(pxlc_init*pxlc_init+pylc_init*pylc_init)<fProdRoughPtMin) return kFALSE;
832 
833  return kTRUE;
834 }
835 //________________________________________________________________________
836 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SelectWithRoughCutsWS(AliAODTrack *vka, AliAODTrack *part)
837 {
838  //
839  // Mass and pT Cut to be applied before object creation
840  //
841 
842  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
843  Double_t mLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
844 
845  Double_t pxpr_init = part->Px();
846  Double_t pypr_init = part->Py();
847  Double_t pzpr_init = part->Pz();
848  Double_t Epr_init = sqrt(pxpr_init*pxpr_init+pypr_init*pypr_init+pzpr_init*pzpr_init+mprPDG*mprPDG);
849 
850  Double_t pxv0_init = vka->Px();
851  Double_t pyv0_init = vka->Py();
852  Double_t pzv0_init = vka->Pz();
853  Double_t massv0_init = 0.493677;
854  Double_t Ev0_init = sqrt(pxv0_init*pxv0_init+pyv0_init*pyv0_init+pzv0_init*pzv0_init+massv0_init*massv0_init);
855 
856  Double_t pxlc_init = pxpr_init+pxv0_init;
857  Double_t pylc_init = pypr_init+pyv0_init;
858  Double_t pzlc_init = pzpr_init+pzv0_init;
859  Double_t Elc_init = Epr_init+Ev0_init;
860  Double_t lcmass_init = sqrt(Elc_init*Elc_init-pxlc_init*pxlc_init-pylc_init*pylc_init-pzlc_init*pzlc_init);
861 
862  if(lcmass_init<mLcPDG-fProdRoughMassTol || lcmass_init>mLcPDG+fProdRoughMassTol) return kFALSE;
863  if(sqrt(pxlc_init*pxlc_init+pylc_init*pylc_init)<fProdRoughPtMin) return kFALSE;
864 
865  return kTRUE;
866 }
867 //________________________________________________________________________
868 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SelectWithRoughCuts(TLorentzVector *v0, TLorentzVector *part)
869 {
870  //
871  // Mass and pT Cut to be applied before object creation
872  //
873 
874  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
875  Double_t mLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
876 
877  Double_t pxpr_init = part->Px();
878  Double_t pypr_init = part->Py();
879  Double_t pzpr_init = part->Pz();
880  Double_t Epr_init = sqrt(pxpr_init*pxpr_init+pypr_init*pypr_init+pzpr_init*pzpr_init+mprPDG*mprPDG);
881 
882  Double_t pxv0_init = v0->Px();
883  Double_t pyv0_init = v0->Py();
884  Double_t pzv0_init = v0->Pz();
885  Double_t massv0_init = v0->M();
886  Double_t Ev0_init = sqrt(pxv0_init*pxv0_init+pyv0_init*pyv0_init+pzv0_init*pzv0_init+massv0_init*massv0_init);
887 
888  Double_t pxlc_init = pxpr_init+pxv0_init;
889  Double_t pylc_init = pypr_init+pyv0_init;
890  Double_t pzlc_init = pzpr_init+pzv0_init;
891  Double_t Elc_init = Epr_init+Ev0_init;
892  Double_t lcmass_init = sqrt(Elc_init*Elc_init-pxlc_init*pxlc_init-pylc_init*pylc_init-pzlc_init*pzlc_init);
893 
894  if(lcmass_init<mLcPDG-fProdRoughMassTol || lcmass_init>mLcPDG+fProdRoughMassTol) return kFALSE;
895  if(sqrt(pxlc_init*pxlc_init+pylc_init*pylc_init)<fProdRoughPtMin) return kFALSE;
896 
897  return kTRUE;
898 }
899 //________________________________________________________________________
900 Bool_t AliRDHFCutsLctopK0sfromAODtracks::SelectWithRoughCutsWS(TLorentzVector *vka, TLorentzVector *part)
901 {
902  //
903  // Mass and pT Cut to be applied before object creation
904  //
905 
906  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
907  Double_t mLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
908 
909  Double_t pxpr_init = part->Px();
910  Double_t pypr_init = part->Py();
911  Double_t pzpr_init = part->Pz();
912  Double_t Epr_init = sqrt(pxpr_init*pxpr_init+pypr_init*pypr_init+pzpr_init*pzpr_init+mprPDG*mprPDG);
913 
914  Double_t pxv0_init = vka->Px();
915  Double_t pyv0_init = vka->Py();
916  Double_t pzv0_init = vka->Pz();
917  Double_t massv0_init = 0.493677;
918  Double_t Ev0_init = sqrt(pxv0_init*pxv0_init+pyv0_init*pyv0_init+pzv0_init*pzv0_init+massv0_init*massv0_init);
919 
920  Double_t pxlc_init = pxpr_init+pxv0_init;
921  Double_t pylc_init = pypr_init+pyv0_init;
922  Double_t pzlc_init = pzpr_init+pzv0_init;
923  Double_t Elc_init = Epr_init+Ev0_init;
924  Double_t lcmass_init = sqrt(Elc_init*Elc_init-pxlc_init*pxlc_init-pylc_init*pylc_init-pzlc_init*pzlc_init);
925 
926  if(lcmass_init<mLcPDG-fProdRoughMassTol || lcmass_init>mLcPDG+fProdRoughMassTol) return kFALSE;
927  if(sqrt(pxlc_init*pxlc_init+pylc_init*pylc_init)<fProdRoughPtMin) return kFALSE;
928 
929  return kTRUE;
930 }
931 //________________________________________________________________________
932 void AliRDHFCutsLctopK0sfromAODtracks::SetMixingWeights(Int_t nbinpr, Double_t *binspr, Int_t nbink0s, Double_t *binsk0s, Double_t *p0val, Double_t *p1val, Double_t *p2val, Double_t *p3val)
933 {
934  //
935  // Set weighting factor for mixing
936  //
937  fNWeightingProtonBinLimits = nbinpr+1;
938  fNWeightingK0sBinLimits = nbink0s+1;
939  fNWeightingBins = nbinpr*nbink0s;
940 
943  fWeight_p0 = new Double_t[fNWeightingBins];
944  fWeight_p1 = new Double_t[fNWeightingBins];
945  fWeight_p2 = new Double_t[fNWeightingBins];
946  fWeight_p3 = new Double_t[fNWeightingBins];
947 
948  for(Int_t i=0;i<fNWeightingProtonBinLimits;i++){
949  fWeightingProtonBins[i] = binspr[i];
950  }
951  for(Int_t i=0;i<fNWeightingK0sBinLimits;i++){
952  fWeightingK0sBins[i] = binsk0s[i];
953  }
954  for(Int_t i=0;i<fNWeightingBins;i++){
955  fWeight_p0[i] = p0val[i];
956  fWeight_p1[i] = p1val[i];
957  fWeight_p2[i] = p2val[i];
958  fWeight_p3[i] = p3val[i];
959  }
960  return;
961 }
962 //________________________________________________________________________
963 Double_t AliRDHFCutsLctopK0sfromAODtracks::GetMixingWeight(Double_t dphi, Double_t deta, Double_t pt_pr, Double_t pt_k0s)
964 {
965  //
966  // Get weighting factor for mixing
967  //
968  if(fNWeightingBins==0) return 1.;
969  if(dphi>M_PI/2.) return 1.;//does not support away side
970 
971  Int_t ibin_pr = -9999;
972  for(Int_t i=0;i<fNWeightingProtonBinLimits-1;i++){
973  if(fWeightingProtonBins[i]<pt_pr && fWeightingProtonBins[i+1]>pt_pr){
974  ibin_pr = i;
975  break;
976  }
977  }
978  Int_t ibin_k0s = -9999;
979  for(Int_t i=0;i<fNWeightingK0sBinLimits-1;i++){
980  if(fWeightingK0sBins[i]<pt_k0s && fWeightingK0sBins[i+1]>pt_k0s){
981  ibin_k0s = i;
982  break;
983  }
984  }
985  if(ibin_pr<0 || ibin_pr > fNWeightingProtonBinLimits-1) return 1.;
986  if(ibin_k0s<0 || ibin_k0s > fNWeightingK0sBinLimits-1) return 1.;
987 
988  Int_t ibin_comb = ibin_pr*(fNWeightingK0sBinLimits-1)+ibin_k0s;
989 
990  Double_t p0 = fWeight_p0[ibin_comb];
991  Double_t p1 = fWeight_p1[ibin_comb];
992  Double_t p2 = fWeight_p2[ibin_comb];
993  Double_t p3 = fWeight_p3[ibin_comb];
994  Double_t weight = p0 + p1 *TMath::Gaus(dphi,0.,p2)*TMath::Gaus(deta,0.,p3);
995 
996  return weight;
997 }
998 //________________________________________________________________________
999 Bool_t AliRDHFCutsLctopK0sfromAODtracks::TagV0(AliAODTrack *ptrk, AliAODEvent *evt, Int_t ntrk, Double_t &minmass)
1000 {
1001  //
1002  // Tag conversion electron tracks
1003  //
1004  if(fTagV0MassTol<0.) return kFALSE;
1005 
1006  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
1007  Double_t mpiPDG = TDatabasePDG::Instance()->GetParticle(211)->Mass();
1008  Double_t mlamPDG = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
1009 
1010  Bool_t isv0 = kFALSE;
1011  Bool_t islam = kFALSE;
1012  minmass = 9999.;
1013 
1014  Int_t trkid = ptrk->GetID();
1015  Double_t px1 = ptrk->Px();
1016  Double_t py1 = ptrk->Py();
1017  Double_t pz1 = ptrk->Pz();
1018  Double_t Epr1 = sqrt(px1*px1+py1*py1+pz1*pz1+mprPDG*mprPDG);
1019 
1020  for(Int_t it=0;it<ntrk;it++){
1021  AliAODTrack *trk2 = (AliAODTrack*) evt->GetTrack(it);
1022  if(!trk2) continue;
1023  Int_t trkid2 = trk2->GetID();
1024  if(trkid==trkid2) continue;
1025  if(ptrk->Charge()*trk2->Charge()>0) continue;
1026  if(!ptrk->TestFilterMask(BIT(4))) continue;
1027 
1028  Double_t px2 = trk2->Px();
1029  Double_t py2 = trk2->Py();
1030  Double_t pz2 = trk2->Pz();
1031  Double_t E2 = sqrt(px2*px2+py2*py2+pz2*pz2+mpiPDG*mpiPDG);
1032 
1033  Double_t mass_lam = sqrt(pow(Epr1+E2,2)-pow(px1+px2,2)-pow(py1+py2,2)-pow(pz1+pz2,2));
1034  Double_t dlam = mass_lam-mlamPDG;
1035  if(fabs(dlam)<fabs(minmass)){
1036  minmass = dlam;
1037  islam = kTRUE;
1038  }
1039  }
1040 
1041  if(fabs(minmass)<fTagV0MassTol) isv0 = kTRUE;
1042 
1043  if(islam) minmass += mlamPDG;
1044 
1045  return isv0;
1046 }
1047 
1048 //________________________________________________________________________
1049 Bool_t AliRDHFCutsLctopK0sfromAODtracks::TagV0SameSign(AliAODTrack *ptrk, AliAODEvent *evt, Int_t ntrk, Double_t &minmass)
1050 {
1051  //
1052  // Tag conversion electron tracks
1053  //
1054  if(fTagV0MassTol<0.) return kFALSE;
1055 
1056  Double_t mprPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
1057  Double_t mpiPDG = TDatabasePDG::Instance()->GetParticle(211)->Mass();
1058  Double_t mlamPDG = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
1059 
1060  Bool_t isv0 = kFALSE;
1061  Bool_t islam = kFALSE;
1062  minmass = 9999.;
1063 
1064  Int_t trkid = ptrk->GetID();
1065  Double_t px1 = ptrk->Px();
1066  Double_t py1 = ptrk->Py();
1067  Double_t pz1 = ptrk->Pz();
1068  Double_t Epr1 = sqrt(px1*px1+py1*py1+pz1*pz1+mprPDG*mprPDG);
1069 
1070  for(Int_t it=0;it<ntrk;it++){
1071  AliAODTrack *trk2 = (AliAODTrack*) evt->GetTrack(it);
1072  if(!trk2) continue;
1073  Int_t trkid2 = trk2->GetID();
1074  if(trkid==trkid2) continue;
1075  if(ptrk->Charge()*trk2->Charge()<0) continue;
1076  if(!ptrk->TestFilterMask(BIT(4))) continue;
1077 
1078  Double_t px2 = trk2->Px();
1079  Double_t py2 = trk2->Py();
1080  Double_t pz2 = trk2->Pz();
1081  Double_t E2 = sqrt(px2*px2+py2*py2+pz2*pz2+mpiPDG*mpiPDG);
1082 
1083  Double_t mass_lam = sqrt(pow(Epr1+E2,2)-pow(px1+px2,2)-pow(py1+py2,2)-pow(pz1+pz2,2));
1084  Double_t dlam = mass_lam-mlamPDG;
1085  if(fabs(dlam)<fabs(minmass)){
1086  minmass = dlam;
1087  islam = kTRUE;
1088  }
1089  }
1090 
1091  if(fabs(minmass)<fTagV0MassTol) isv0 = kTRUE;
1092 
1093  if(islam) minmass += mlamPDG;
1094 
1095  return isv0;
1096 }
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Double_t * fWeight_p0
Number of bins for mixing weight should be proton x k0s.
Double_t fProdV0MassRejLambda
K0s mass selection used before object creation.
Double_t fProdV0EtaMin
V0 daughter Minimum TPC cluster pT used before object creation.
Double_t GetMixingWeight(Double_t dphi, Double_t deta, Double_t pt_pr, Double_t pt_k0s)
Double_t GetSigma(Int_t idet) const
Definition: AliAODPidHF.h:134
Double_t fProdV0DaughterTPCClusterMin
V0 Daughter pT min used before object creation.
Double_t fProdTrackTPCNclsRatioMin
Min. Number of TPC PID cluster.
Bool_t GetUseCombined()
Definition: AliAODPidHF.h:165
Int_t fProdTrackTPCNclsPIDMin
Flag to check if we use on-the-fly v0.
void SetNVars(Int_t nVars)
Definition: AliRDHFCuts.h:362
AliRDHFCutsLctopK0sfromAODtracks & operator=(const AliRDHFCutsLctopK0sfromAODtracks &source)
Double_t fProdRoughPtMin
Mass cut for Lc used before object creation.
Double_t fProdV0RapMax
Minimum rapidity of cascade.
AliRDHFCuts & operator=(const AliRDHFCuts &source)
Double_t * fWeightingK0sBins
Number of bins for k0s.
Bool_t fUsePID
Definition: AliRDHFCuts.h:389
Double_t fProdV0DaughterEtaRange
Max DCA between V0 daughters used before object creation.
Double_t fProdV0MassTolK0s
Flag for AOD filter Bit used before object creation.
Bool_t SingleV0Cuts(AliAODv0 *v0, AliAODVertex *vert)
Double_t InvMassLctoK0sP() const
Double_t fProdV0CosPointingAngleToPrimVtxMin
Minimum K0s pT used before object creation.
Double_t fProdRoughMassTol
Maximum rapidity of cascade.
Double_t fProdV0PtMin
photon mass rejection used before object creation
Int_t fnVarsForOpt
Definition: AliRDHFCuts.h:384
AliPIDCombined * GetPidCombined() const
Definition: AliAODPidHF.h:161
Double_t fMaxPtCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:411
Double_t CalculateLcCosPAXY(AliAODRecoDecayHF *obj)
AliESDtrackCuts * fTrackCuts
quality cuts on the daughter tracks
Definition: AliRDHFCuts.h:377
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Bool_t fUseOnTheFlyV0
Threshold used in IsSelectedCombinedPID.
Bool_t IsDaughterSelected(AliAODTrack *track, const AliESDVertex *primary, AliESDtrackCuts *cuts) const
AliAODTrack * GetBachelor() const
AliAODVertex * GetOwnPrimaryVtx() const
Bool_t SingleTrkCuts(AliAODTrack *trk, AliAODVertex *vtx)
Double_t fProdV0DcaDaughtersMax
V0 pointing angle used before object creation.
Double_t fProdV0RapMin
Maximum eta of cascade.
virtual Int_t IsSelectedPID(AliAODRecoDecayHF *obj)
Int_t MakeRawPid(AliAODTrack *track, Int_t specie)
AliPIDResponse * GetPidResponse() const
Definition: AliAODPidHF.h:160
Bool_t SelectWithRoughCutsWS(AliAODTrack *vka, AliAODTrack *trk1)
Float_t * fCutsRD
fnVars*fnPtBins
Definition: AliRDHFCuts.h:387
Bool_t fProdUseAODFilterBit
Min. Number of TPC PID cluster.
Double_t fProdV0DaughterPtMin
V0Daughter eta range used before object creation.
void SetVarsForOpt(Int_t nVars, Bool_t *forOpt)
Double_t DecayLengthXY() const
void SetVarNames(Int_t nVars, TString *varNames, Bool_t *isUpperCut)
Bool_t * fVarsForOpt
number of cut vars to be optimized for candidates
Definition: AliRDHFCuts.h:385
Double_t fProdV0EtaMax
Minimum eta of cascade.
Double_t * fWeightingProtonBins
Number of bins for proton.
Bool_t SingleKaonCuts(AliAODTrack *trk, AliAODVertex *vtx)
virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d, Float_t *vars, Int_t nvars, Int_t *pdgdaughters)
AliRDHFCutsLctopK0sfromAODtracks(const char *name="CutsLctopK0s")
Int_t fNWeightingProtonBinLimits
pT cut for Lc used before object creation
Bool_t TagV0(AliAODTrack *etrk, AliAODEvent *evt, Int_t ntrk, Double_t &minmass)
void SetPtBins(Int_t nPtBinLimits, Float_t *ptBinLimits)
void SetMixingWeights(Int_t nbinpr, Double_t *bins_pr, Int_t nbink0s, Double_t *bins_k0s, Double_t *p0val, Double_t *p1val, Double_t *p2val, Double_t *p3val)
Bool_t SelectWithRoughCuts(AliAODv0 *v0, AliAODTrack *trk1)
Bool_t TagV0SameSign(AliAODTrack *etrk, AliAODEvent *evt, Int_t ntrk, Double_t &minmass)
AliAODPidHF * fPidHF
enable AOD049 centrality cleanup
Definition: AliRDHFCuts.h:391
Int_t PtBin(Double_t pt) const
Int_t GetGlobalIndex(Int_t iVar, Int_t iPtBin) const
Double_t fMinPtCand
outcome of PID selection
Definition: AliRDHFCuts.h:410
void SetPidResponse(AliPIDResponse *pidResp)
Definition: AliAODPidHF.h:111
Double_t fProdV0MassRejPhoton
lambda mass rejection used before object creation