AliPhysics  v5-07-15-01 (b3d7633)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliRDHFCutsD0toKpi.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 D0->Kpi
21 //
22 // Author: A.Dainese, andrea.dainese@pd.infn.it
24 
25 #include <TDatabasePDG.h>
26 #include <Riostream.h>
27 
28 #include "AliRDHFCutsD0toKpi.h"
30 #include "AliAODTrack.h"
31 #include "AliESDtrack.h"
32 #include "AliAODPid.h"
33 #include "AliTPCPIDResponse.h"
34 #include "AliAODVertex.h"
35 #include "AliKFVertex.h"
36 #include "AliKFParticle.h"
37 
38 using std::cout;
39 using std::endl;
40 
41 
45 
46 
47 //--------------------------------------------------------------------------
49  AliRDHFCuts(name),
50  fUseSpecialCuts(kFALSE),
51  fLowPt(kTRUE),
52  fDefaultPID(kFALSE),
53  fUseKF(kFALSE),
54  fPtLowPID(2.),
55  fPtMaxSpecialCuts(9999.),
56  fmaxPtrackForPID(9999),
57  fCombPID(kFALSE),
58  fnSpecies(AliPID::kSPECIES),
59  fWeightsPositive(0),
60  fWeightsNegative(0),
61  fProbThreshold(0.5),
62  fBayesianStrategy(0),
63  fBayesianCondition(0)
64 {
65  //
66  // Default Constructor
67  //
68  Int_t nvars=11;
69  SetNVars(nvars);
70  TString varNames[11]={"inv. mass [GeV]",
71  "dca [cm]",
72  "cosThetaStar",
73  "pTK [GeV/c]",
74  "pTPi [GeV/c]",
75  "d0K [cm]",
76  "d0Pi [cm]",
77  "d0d0 [cm^2]",
78  "cosThetaPoint",
79  "|cosThetaPointXY|",
80  "NormDecayLenghtXY"};
81  Bool_t isUpperCut[11]={kTRUE,
82  kTRUE,
83  kTRUE,
84  kFALSE,
85  kFALSE,
86  kTRUE,
87  kTRUE,
88  kTRUE,
89  kFALSE,
90  kFALSE,
91  kFALSE};
92  SetVarNames(nvars,varNames,isUpperCut);
93  Bool_t forOpt[11]={kFALSE,
94  kTRUE,
95  kTRUE,
96  kFALSE,
97  kFALSE,
98  kFALSE,
99  kFALSE,
100  kTRUE,
101  kTRUE,
102  kFALSE,
103  kFALSE};
104  SetVarsForOpt(4,forOpt);
105  Float_t limits[2]={0,999999999.};
106  SetPtBins(2,limits);
107 
108  fWeightsNegative = new Double_t[AliPID::kSPECIES];
109  fWeightsPositive = new Double_t[AliPID::kSPECIES];
110 
111  for (Int_t i = 0; i<AliPID::kSPECIES; i++) {
112  fWeightsPositive[i] = 0.;
113  fWeightsNegative[i] = 0.;
114  }
115 
116 }
117 //--------------------------------------------------------------------------
119  AliRDHFCuts(source),
120  fUseSpecialCuts(source.fUseSpecialCuts),
121  fLowPt(source.fLowPt),
122  fDefaultPID(source.fDefaultPID),
123  fUseKF(source.fUseKF),
124  fPtLowPID(source.fPtLowPID),
125  fPtMaxSpecialCuts(source.fPtMaxSpecialCuts),
126  fmaxPtrackForPID(source.fmaxPtrackForPID),
127  fCombPID(source.fCombPID),
128  fnSpecies(source.fnSpecies),
129  fWeightsPositive(source.fWeightsPositive),
130  fWeightsNegative(source.fWeightsNegative),
131  fProbThreshold(source.fProbThreshold),
132  fBayesianStrategy(source.fBayesianStrategy),
133  fBayesianCondition(source.fBayesianCondition)
134 {
135  //
136  // Copy constructor
137  //
138 
139 }
140 //--------------------------------------------------------------------------
142 {
143  //
144  // assignment operator
145  //
146  if(&source == this) return *this;
147 
148  AliRDHFCuts::operator=(source);
150  fLowPt=source.fLowPt;
151  fDefaultPID=source.fDefaultPID;
152  fUseKF=source.fUseKF;
153  fPtLowPID=source.fPtLowPID;
156  fCombPID = source.fCombPID;
157  fnSpecies = source.fnSpecies;
163  return *this;
164 }
165 
166 
167 
168 //---------------------------------------------------------------------------
170 {
171 //
172 // Destructor
173 //
174  if (fWeightsNegative) {
175  delete[] fWeightsNegative;
176  fWeightsNegative = 0;
177  }
178  if (fWeightsPositive) {
179  delete[] fWeightsNegative;
180  fWeightsNegative = 0;
181  }
182 }
183 
184 //---------------------------------------------------------------------------
185 void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent *aod) {
186  //
187  // Fills in vars the values of the variables
188  //
189 
190  if(nvars!=fnVarsForOpt) {
191  printf("AliRDHFCutsD0toKpi::GetCutsVarsForOpt: wrong number of variables\n");
192  return;
193  }
194 
196 
197  //recalculate vertex w/o daughters
198  Bool_t cleanvtx=kFALSE;
199  AliAODVertex *origownvtx=0x0;
201  if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());
202  cleanvtx=kTRUE;
203  if(!RecalcOwnPrimaryVtx(dd,aod)) {
204  CleanOwnPrimaryVtx(dd,aod,origownvtx);
205  cleanvtx=kFALSE;
206  }
207  }
208 
209  Int_t iter=-1;
210  if(fVarsForOpt[0]){
211  iter++;
212  if(TMath::Abs(pdgdaughters[0])==211) {
213  vars[iter]=dd->InvMassD0();
214  } else {
215  vars[iter]=dd->InvMassD0bar();
216  }
217  }
218  if(fVarsForOpt[1]){
219  iter++;
220  vars[iter]=dd->GetDCA();
221  }
222  if(fVarsForOpt[2]){
223  iter++;
224  if(TMath::Abs(pdgdaughters[0])==211) {
225  vars[iter] = dd->CosThetaStarD0();
226  } else {
227  vars[iter] = dd->CosThetaStarD0bar();
228  }
229  }
230  if(fVarsForOpt[3]){
231  iter++;
232  if(TMath::Abs(pdgdaughters[0])==321) {
233  vars[iter]=dd->PtProng(0);
234  }
235  else{
236  vars[iter]=dd->PtProng(1);
237  }
238  }
239  if(fVarsForOpt[4]){
240  iter++;
241  if(TMath::Abs(pdgdaughters[0])==211) {
242  vars[iter]=dd->PtProng(0);
243  }
244  else{
245  vars[iter]=dd->PtProng(1);
246  }
247  }
248  if(fVarsForOpt[5]){
249  iter++;
250  if(TMath::Abs(pdgdaughters[0])==321) {
251  vars[iter]=dd->Getd0Prong(0);
252  }
253  else{
254  vars[iter]=dd->Getd0Prong(1);
255  }
256  }
257  if(fVarsForOpt[6]){
258  iter++;
259  if(TMath::Abs(pdgdaughters[0])==211) {
260  vars[iter]=dd->Getd0Prong(0);
261  }
262  else{
263  vars[iter]=dd->Getd0Prong(1);
264  }
265  }
266  if(fVarsForOpt[7]){
267  iter++;
268  vars[iter]= dd->Prodd0d0();
269  }
270  if(fVarsForOpt[8]){
271  iter++;
272  vars[iter]=dd->CosPointingAngle();
273  }
274 
275  if(fVarsForOpt[9]){
276  iter++;
277  vars[iter]=TMath::Abs(dd->CosPointingAngleXY());
278  }
279 
280  if(fVarsForOpt[10]){
281  iter++;
282  vars[iter]=dd->NormalizedDecayLengthXY();
283  }
284 
285  if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
286 
287  return;
288 }
289 //---------------------------------------------------------------------------
290 Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* aod) {
291  //
292  // Apply selection
293  //
294 
295 
296  fIsSelectedCuts=0;
297  fIsSelectedPID=0;
298 
299  if(!fCutsRD){
300  cout<<"Cut matrice not inizialized. Exit..."<<endl;
301  return 0;
302  }
303  //PrintAll();
305  if(!d){
306  cout<<"AliAODRecoDecayHF2Prong null"<<endl;
307  return 0;
308  }
309 
310  if(fKeepSignalMC) if(IsSignalMC(d,aod,421)) return 3;
311 
312  Double_t ptD=d->Pt();
313  if(ptD<fMinPtCand) return 0;
314  if(ptD>fMaxPtCand) return 0;
315 
317 
318  // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
319  Int_t returnvaluePID=3;
320  Int_t returnvalueCuts=3;
321 
322  // selection on candidate
323  if(selectionLevel==AliRDHFCuts::kAll ||
324  selectionLevel==AliRDHFCuts::kCandidate) {
325 
326  if(!fUseKF) {
327 
328  //recalculate vertex w/o daughters
329  AliAODVertex *origownvtx=0x0;
331  if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
332  if(!RecalcOwnPrimaryVtx(d,aod)) {
333  CleanOwnPrimaryVtx(d,aod,origownvtx);
334  return 0;
335  }
336  }
337 
338  if(fUseMCVertex) {
339  if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
340  if(!SetMCPrimaryVtx(d,aod)) {
341  CleanOwnPrimaryVtx(d,aod,origownvtx);
342  return 0;
343  }
344  }
345 
346  Double_t pt=d->Pt();
347 
348  Int_t okD0=0,okD0bar=0;
349 
350  Int_t ptbin=PtBin(pt);
351  if (ptbin==-1) {
352  CleanOwnPrimaryVtx(d,aod,origownvtx);
353  return 0;
354  }
355 
356  Double_t mD0,mD0bar,ctsD0,ctsD0bar;
357  okD0=1; okD0bar=1;
358 
359  Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
360 
361  d->InvMassD0(mD0,mD0bar);
362  if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
363  if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;
364  if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
365 
366  if(d->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
367 
368 
369  if(d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(3,ptbin)]*fCutsRD[GetGlobalIndex(3,ptbin)] || d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(4,ptbin)]*fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
370  if(d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(3,ptbin)]*fCutsRD[GetGlobalIndex(3,ptbin)] || d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(4,ptbin)]*fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;
371  if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
372 
373 
374  if(TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
375  TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
376  if(TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||
377  TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;
378  if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
379 
380  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
381 
382 
383  d->CosThetaStarD0(ctsD0,ctsD0bar);
384  if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
385  if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;
386  if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
387 
388  if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
389 
390 
391  if(TMath::Abs(d->CosPointingAngleXY()) < fCutsRD[GetGlobalIndex(9,ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
392 
393  Double_t normalDecayLengXY=d->NormalizedDecayLengthXY();
394  if (normalDecayLengXY < fCutsRD[GetGlobalIndex(10, ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
395 
396  if (returnvalueCuts!=0) {
397  if (okD0) returnvalueCuts=1; //cuts passed as D0
398  if (okD0bar) returnvalueCuts=2; //cuts passed as D0bar
399  if (okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar
400  }
401 
402  // call special cuts
403  Int_t special=1;
405  if(!special) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
406 
407  // unset recalculated primary vertex when not needed any more
408  CleanOwnPrimaryVtx(d,aod,origownvtx);
409 
410  } else {
411  // go to selection with Kalman vertexing, if requested
412  returnvalueCuts = IsSelectedKF(d,aod);
413  }
414 
415  fIsSelectedCuts=returnvalueCuts;
416  if(!returnvalueCuts) return 0;
417  }
418 
419 
420  // selection on PID
421  if(selectionLevel==AliRDHFCuts::kAll ||
422  selectionLevel==AliRDHFCuts::kCandidate ||
423  selectionLevel==AliRDHFCuts::kPID) {
424  if (!fCombPID) {
425  returnvaluePID = IsSelectedPID(d);
426  fIsSelectedPID=returnvaluePID;
427  if(!returnvaluePID) return 0;
428  } else {
429  returnvaluePID = IsSelectedCombPID(d);
430  if(!returnvaluePID) return 0;
431  }
432  }
433 
434 
435 
436  Int_t returnvalueComb=CombineSelectionLevels(3,returnvalueCuts,returnvaluePID);
437 
438  if(!returnvalueComb) return 0;
439 
440  // selection on daughter tracks
441  if(selectionLevel==AliRDHFCuts::kAll ||
442  selectionLevel==AliRDHFCuts::kTracks) {
443  if(!AreDaughtersSelected(d)) return 0;
444  }
445 
446  // cout<<"Pid = "<<returnvaluePID<<endl;
447  return returnvalueComb;
448 }
449 
450 //------------------------------------------------------------------------------------------
452  AliAODEvent* aod) const {
453  //
454  // Apply selection using KF-vertexing
455  //
456 
457  AliAODTrack *track0 = (AliAODTrack*)d->GetDaughter(0);
458  AliAODTrack *track1 = (AliAODTrack*)d->GetDaughter(1);
459 
460  if(!track0 || !track1) {
461  cout<<"one or two D0 daughters missing!"<<endl;
462  return 0;
463  }
464 
465  // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
466  Int_t returnvalueCuts=3;
467 
468  // candidate selection with AliKF
469  AliKFParticle::SetField(aod->GetMagneticField()); // set the magnetic field
470 
471  Int_t okD0=0,okD0bar=0;
472  okD0=1; okD0bar=1;
473 
474  // convert tracks into AliKFParticles
475 
476  AliKFParticle negPiKF(*track1,-211); // neg pion kandidate
477  AliKFParticle negKKF(*track1,-321); // neg kaon kandidate
478  AliKFParticle posPiKF(*track0,211); // pos pion kandidate
479  AliKFParticle posKKF(*track0,321); // pos kaon kandidate
480 
481  // build D0 candidates
482 
483  AliKFParticle d0c(negKKF,posPiKF); // D0 candidate
484  AliKFParticle ad0c(posKKF,negPiKF); // D0bar candidate
485 
486  // create kf primary vertices
487 
488  AliAODVertex *vtx1 = aod->GetPrimaryVertex();
489  AliKFVertex primVtx(*vtx1);
490  AliKFVertex aprimVtx(*vtx1);
491 
492  if(primVtx.GetNContributors()<=0) okD0 = 0;
493  if(aprimVtx.GetNContributors()<=0) okD0bar = 0;
494  if(!okD0 && !okD0bar) returnvalueCuts=0;
495 
496  // calculate mass
497 
498  Double_t d0mass = d0c.GetMass();
499  Double_t ad0mass = ad0c.GetMass();
500 
501  // calculate P of D0 and D0bar
502  Double_t d0P = d0c.GetP();
503  Double_t d0Px = d0c.GetPx();
504  Double_t d0Py = d0c.GetPy();
505  Double_t d0Pz = d0c.GetPz();
506  Double_t ad0P = ad0c.GetP();
507  Double_t ad0Px = ad0c.GetPx();
508  Double_t ad0Py = ad0c.GetPy();
509  Double_t ad0Pz = ad0c.GetPz();
510 
511  //calculate Pt of D0 and D0bar
512 
513  Double_t pt=d0c.GetPt();
514  Double_t apt=ad0c.GetPt();
515 
516  // remove D0 daughters from primary vertices (if used in vertex fit) and add D0-candidates
517 
518  if(track0->GetUsedForPrimVtxFit()) {
519  primVtx -= posPiKF;
520  aprimVtx -= posKKF;
521  }
522 
523  if(track1->GetUsedForPrimVtxFit()) {
524  primVtx -= negKKF;
525  aprimVtx -= negPiKF;
526  }
527 
528  primVtx += d0c;
529  aprimVtx += ad0c;
530 
531  if(primVtx.GetNContributors()<=0) okD0 = 0;
532  if(aprimVtx.GetNContributors()<=0) okD0bar = 0;
533  if(!okD0 && !okD0bar) returnvalueCuts=0;
534 
535  //calculate cut variables
536 
537  // calculate impact params of daughters w.r.t recalculated vertices
538 
539  Double_t impactPi = posPiKF.GetDistanceFromVertexXY(primVtx);
540  Double_t aimpactPi = negPiKF.GetDistanceFromVertexXY(aprimVtx);
541  Double_t impactKa = negKKF.GetDistanceFromVertexXY(primVtx);
542  Double_t aimpactKa = posKKF.GetDistanceFromVertexXY(aprimVtx);
543 
544  // calculate Product of Impact Params
545 
546  Double_t prodParam = impactPi*impactKa;
547  Double_t aprodParam = aimpactPi*aimpactKa;
548 
549  // calculate cosine of pointing angles
550 
551  TVector3 mom(d0c.GetPx(),d0c.GetPy(),d0c.GetPz());
552  TVector3 fline(d0c.GetX()-primVtx.GetX(),
553  d0c.GetY()-primVtx.GetY(),
554  d0c.GetZ()-primVtx.GetZ());
555  Double_t pta = mom.Angle(fline);
556  Double_t cosP = TMath::Cos(pta); // cosine of pta for D0 candidate
557 
558  TVector3 amom(ad0c.GetPx(),ad0c.GetPy(),ad0c.GetPz());
559  TVector3 afline(ad0c.GetX()-aprimVtx.GetX(),
560  ad0c.GetY()-aprimVtx.GetY(),
561  ad0c.GetZ()-aprimVtx.GetZ());
562  Double_t apta = amom.Angle(afline);
563  Double_t acosP = TMath::Cos(apta); // cosine of pta for D0bar candidate
564 
565  // calculate P of Pions at Decay Position of D0 and D0bar candidates
566  negKKF.TransportToParticle(d0c);
567  posPiKF.TransportToParticle(d0c);
568  posKKF.TransportToParticle(ad0c);
569  negPiKF.TransportToParticle(ad0c);
570 
571  Double_t pxPi = posPiKF.GetPx();
572  Double_t pyPi = posPiKF.GetPy();
573  Double_t pzPi = posPiKF.GetPz();
574  Double_t ptPi = posPiKF.GetPt();
575 
576  Double_t apxPi = negPiKF.GetPx();
577  Double_t apyPi = negPiKF.GetPy();
578  Double_t apzPi = negPiKF.GetPz();
579  Double_t aptPi = negPiKF.GetPt();
580 
581  // calculate Pt of Kaons at Decay Position of D0 and D0bar candidates
582 
583  Double_t ptK = negKKF.GetPt();
584  Double_t aptK = posKKF.GetPt();
585 
586  //calculate cos(thetastar)
587  Double_t massvtx = TDatabasePDG::Instance()->GetParticle(421)->Mass();
588  Double_t massp[2];
589  massp[0] = TDatabasePDG::Instance()->GetParticle(321)->Mass();
590  massp[1] = TDatabasePDG::Instance()->GetParticle(211)->Mass();
591  Double_t pStar = TMath::Sqrt(TMath::Power(massvtx*massvtx-massp[0]*massp[0]-massp[1]*massp[1],2.)
592  -4.*massp[0]*massp[0]*massp[1]*massp[1])/(2.*massvtx);
593 
594  // cos(thetastar) for D0 and Pion
595 
596  Double_t d0E = TMath::Sqrt(massvtx*massvtx + d0P*d0P);
597  Double_t beta = d0P/d0E;
598  Double_t gamma = d0E/massvtx;
599  TVector3 momPi(pxPi,pyPi,pzPi);
600  TVector3 momTot(d0Px,d0Py,d0Pz);
601  Double_t q1 = momPi.Dot(momTot)/momTot.Mag();
602  Double_t cts = (q1/gamma-beta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar;
603 
604  // cos(thetastar) for D0bar and Pion
605 
606  Double_t ad0E = TMath::Sqrt(massvtx*massvtx + ad0P*ad0P);
607  Double_t abeta = ad0P/ad0E;
608  Double_t agamma = ad0E/massvtx;
609  TVector3 amomPi(apxPi,apyPi,apzPi);
610  TVector3 amomTot(ad0Px,ad0Py,ad0Pz);
611  Double_t aq1 = amomPi.Dot(amomTot)/amomTot.Mag();
612  Double_t acts = (aq1/agamma-abeta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar;
613 
614  // calculate reduced Chi2 for the full D0 fit
615  d0c.SetProductionVertex(primVtx);
616  ad0c.SetProductionVertex(aprimVtx);
617  negKKF.SetProductionVertex(d0c);
618  posPiKF.SetProductionVertex(d0c);
619  posKKF.SetProductionVertex(ad0c);
620  negPiKF.SetProductionVertex(ad0c);
621  d0c.TransportToProductionVertex();
622  ad0c.TransportToProductionVertex();
623 
624  // calculate the decay length
625  Double_t decayLengthD0 = d0c.GetDecayLength();
626  Double_t adecayLengthD0 = ad0c.GetDecayLength();
627 
628  Double_t chi2D0 = 50.;
629  if(d0c.GetNDF() > 0 && d0c.GetChi2() >= 0) {
630  chi2D0 = d0c.GetChi2()/d0c.GetNDF();
631  }
632 
633  Double_t achi2D0 = 50.;
634  if(ad0c.GetNDF() > 0 && ad0c.GetChi2() >= 0) {
635  achi2D0 = ad0c.GetChi2()/ad0c.GetNDF();
636  }
637 
638  // Get the Pt-bins
639  Int_t ptbin=PtBin(pt);
640  Int_t aptbin=PtBin(apt);
641 
642  if(ptbin < 0) okD0 = 0;
643  if(aptbin < 0) okD0bar = 0;
644  if(!okD0 && !okD0bar) returnvalueCuts=0;
645 
646  if(ptK < fCutsRD[GetGlobalIndex(3,ptbin)] || ptPi < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
647  if(aptK < fCutsRD[GetGlobalIndex(3,aptbin)] || aptPi < fCutsRD[GetGlobalIndex(4,aptbin)]) okD0bar = 0;
648  if(!okD0 && !okD0bar) returnvalueCuts=0;
649 
650 
651  if(TMath::Abs(impactKa) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
652  TMath::Abs(impactPi) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
653 
654  if(TMath::Abs(aimpactKa) > fCutsRD[GetGlobalIndex(5,aptbin)] ||
655  TMath::Abs(aimpactPi) > fCutsRD[GetGlobalIndex(6,aptbin)]) okD0bar = 0;
656 
657  if(!okD0 && !okD0bar) returnvalueCuts=0;
658 
659  // for the moment via the standard method due to bug in AliKF
660  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) okD0 = 0;
661  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,aptbin)]) okD0bar = 0;
662  if(!okD0 && !okD0bar) returnvalueCuts=0;
663 
664 
665  if(TMath::Abs(d0mass-massvtx) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
666  if(TMath::Abs(ad0mass-massvtx) > fCutsRD[GetGlobalIndex(0,aptbin)]) okD0bar = 0;
667  if(!okD0 && !okD0bar) returnvalueCuts=0;
668 
669 
670  if(TMath::Abs(cts) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
671  if(TMath::Abs(acts) > fCutsRD[GetGlobalIndex(2,aptbin)]) okD0bar = 0;
672  if(!okD0 && !okD0bar) returnvalueCuts=0;
673 
674  if(prodParam > fCutsRD[GetGlobalIndex(7,ptbin)]) okD0 = 0;
675  if(aprodParam > fCutsRD[GetGlobalIndex(7,aptbin)]) okD0bar = 0;
676  if(!okD0 && !okD0bar) returnvalueCuts=0;
677 
678  if(cosP < fCutsRD[GetGlobalIndex(8,ptbin)]) okD0 = 0;
679  if(acosP < fCutsRD[GetGlobalIndex(8,aptbin)]) okD0bar = 0;
680  if(!okD0 && !okD0bar) returnvalueCuts=0;
681 
682  if(chi2D0 > fCutsRD[GetGlobalIndex(10,ptbin)]) okD0 = 0;
683  if(achi2D0 > fCutsRD[GetGlobalIndex(10,aptbin)]) okD0bar = 0;
684  if(!okD0 && !okD0bar) returnvalueCuts=0;
685 
686  if(decayLengthD0 < fCutsRD[GetGlobalIndex(9,ptbin)]) okD0 = 0;
687  if(adecayLengthD0 < fCutsRD[GetGlobalIndex(9,aptbin)]) okD0bar = 0;
688  if(!okD0 && !okD0bar) returnvalueCuts=0;
689 
690  if(returnvalueCuts!=0) {
691  if(okD0) returnvalueCuts=1; //cuts passed as D0
692  if(okD0bar) returnvalueCuts=2; //cuts passed as D0bar
693  if(okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar
694  }
695 
696  return returnvalueCuts;
697 }
698 
699 //---------------------------------------------------------------------------
700 
701 Bool_t AliRDHFCutsD0toKpi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
702 {
703  //
704  // Checking if D0 is in fiducial acceptance region
705  //
706 
707  if(fMaxRapidityCand>-998.){
708  if(TMath::Abs(y) > fMaxRapidityCand) return kFALSE;
709  else return kTRUE;
710  }
711 
712  if(pt > 5.) {
713  // applying cut for pt > 5 GeV
714  AliDebug(2,Form("pt of D0 = %f (> 5), cutting at |y| < 0.8\n",pt));
715  if (TMath::Abs(y) > 0.8){
716  return kFALSE;
717  }
718  } else {
719  // appliying smooth cut for pt < 5 GeV
720  Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
721  Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
722  AliDebug(2,Form("pt of D0 = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
723  if (y < minFiducialY || y > maxFiducialY){
724  return kFALSE;
725  }
726  }
727 
728  return kTRUE;
729 }
730 //---------------------------------------------------------------------------
732 {
733  // ############################################################
734  //
735  // Apply PID selection
736  //
737  //
738  // ############################################################
739 
740  if(!fUsePID) return 3;
741  if(fDefaultPID) return IsSelectedPIDdefault(d);
742  if (fCombPID) return IsSelectedCombPID(d); //to use Bayesian method if applicable
743  fWhyRejection=0;
744  Int_t isD0D0barPID[2]={1,2};
745  Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value;
746  // same for prong 2
747  // values convention -1 = discarded
748  // 0 = not identified (but compatible) || No PID (->hasPID flag)
749  // 1 = identified
750  // PID search: pion (TPC) or not K (TOF), Kaon hypothesis for both
751  // Initial hypothesis: unknwon (but compatible)
752  combinedPID[0][0]=0; // prima figlia, Kaon
753  combinedPID[0][1]=0; // prima figlia, pione
754  combinedPID[1][0]=0; // seconda figlia, Kaon
755  combinedPID[1][1]=0; // seconda figlia, pion
756 
757  Bool_t checkPIDInfo[2]={kTRUE,kTRUE};
758  Double_t sigma_tmp[3]={fPidHF->GetSigma(0),fPidHF->GetSigma(1),fPidHF->GetSigma(2)};
759  Bool_t isTOFused=fPidHF->GetTOF(),isCompat=fPidHF->GetCompat();
760  for(Int_t daught=0;daught<2;daught++){
761  //Loop con prongs
762  AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);
763  if(fPidHF->IsTOFPiKexcluded(aodtrack,5.)) return 0;
764 
765  if(!(fPidHF->CheckStatus(aodtrack,"TPC")) && !(fPidHF->CheckStatus(aodtrack,"TOF"))) {
766  checkPIDInfo[daught]=kFALSE;
767  continue;
768  }
769  Double_t pProng=aodtrack->P();
770 
771  // identify kaon
772  if(pProng<fmaxPtrackForPID){
773  combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
774  }
775  // identify pion
776  if(pProng<fmaxPtrackForPID){
777  if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {
778  combinedPID[daught][1]=0;
779  }else{
780  fPidHF->SetTOF(kFALSE);
781  combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
782  if(isTOFused)fPidHF->SetTOF(kTRUE);
783  if(isCompat)fPidHF->SetCompat(kTRUE);
784  }
785  }
786 
787  if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){ // if not a K- and not a pi- both D0 and D0bar excluded
788  isD0D0barPID[0]=0;
789  isD0D0barPID[1]=0;
790  }
791  else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){
792  if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded
793  else isD0D0barPID[0]=0;// if K+ D0 excluded
794  }
795  /* else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){
796  isD0D0barPID[0]=0;
797  isD0D0barPID[1]=0;
798  }
799  */
800  else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){
801  if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded
802  else isD0D0barPID[0]=0;// not a D0 if K+ or if pi+ excluded
803  }
804  else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){
805  if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded
806  else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded
807  }
808 
809  if(fLowPt && d->Pt()<fPtLowPID){
810  Double_t sigmaTPC[3]={3.,2.,0.};
811  fPidHF->SetSigmaForTPC(sigmaTPC);
812  // identify kaon
813  combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
814 
815  if(pProng<0.6){
816  fPidHF->SetCompat(kFALSE);
817  combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
818  fPidHF->SetCompat(kTRUE);
819  }
820 
821  if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {
822  combinedPID[daught][1]=0;
823  }else{
824  fPidHF->SetTOF(kFALSE);
825  Double_t sigmaTPCpi[3]={3.,3.,0.};
826  fPidHF->SetSigmaForTPC(sigmaTPCpi);
827  combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
828  fPidHF->SetTOF(kTRUE);
829  if(pProng<0.8){
830  Bool_t isTPCpion=fPidHF->IsPionRaw(aodtrack,"TPC");
831  if(isTPCpion){
832  combinedPID[daught][1]=1;
833  }else{
834  combinedPID[daught][1]=-1;
835  }
836  }
837  }
838 
839  }
840  fPidHF->SetSigmaForTPC(sigma_tmp);
841  }// END OF LOOP ON DAUGHTERS
842 
843  if(!checkPIDInfo[0] && !checkPIDInfo[1]) {
844  if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
845  return 0;
846  }
847 
848 
849  // FURTHER PID REQUEST (both daughter info is needed)
850  if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){
851  fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found
852  if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
853  return 0;
854  }
855 
856  if(fLowPt && d->Pt()<fPtLowPID){
857  if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
858  fWhyRejection=32;// reject cases where the Kaon is not identified
859  fPidHF->SetSigmaForTPC(sigma_tmp);
860  return 0;
861  }
862  }
863  if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
864 
865  // cout<<"Why? "<<fWhyRejection<<endl;
866  return isD0D0barPID[0]+isD0D0barPID[1];
867 }
868 //---------------------------------------------------------------------------
870 {
871  // ############################################################
872  //
873  // Apply PID selection
874  //
875  //
876  // temporary selection: PID AS USED FOR D0 by Andrea Rossi (up to 28/06/2010)
877  //
878  // d must be a AliAODRecoDecayHF2Prong object
879  // returns 0 if both D0 and D0bar are rejectecd
880  // 1 if D0 is accepted while D0bar is rejected
881  // 2 if D0bar is accepted while D0 is rejected
882  // 3 if both are accepted
883  // fWhyRejection variable (not returned for the moment, print it if needed)
884  // keeps some information on why a candidate has been
885  // rejected according to the following (unfriendly?) scheme
886  // if more rejection cases are considered interesting, just add numbers
887  //
888  // TO BE CONSIDERED WITH A GRAIN OF SALT (the order in which cut are applied is relevant)
889  // from 20 to 30: "detector" selection (PID acceptance)
890  // 26: TPC refit
891  // 27: ITS refit
892  // 28: no (TOF||TPC) pid information (no kTOFpid,kTOFout,kTIME,kTPCpid,...)
893  //
894  // from 30 to 40: PID selection
895  // 31: no Kaon compatible tracks found between daughters
896  // 32: no Kaon identified tracks found (strong sel. at low momenta)
897  // 33: both mass hypotheses are rejected
898  //
899  // ############################################################
900 
901  if(!fUsePID) return 3;
902  fWhyRejection=0;
903  Int_t isD0D0barPID[2]={1,2};
904  Double_t nsigmaTPCpi=-1., nsigmaTPCK=-1.; //used for TPC pid
905  Double_t tofSig,times[5];// used fot TOF pid
906  Int_t hasPID[2]={2,2};// flag to count how many detectors give PID info for the daughters
907  Int_t isKaonPionTOF[2][2],isKaonPionTPC[2][2];
908  Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value;
909  // same for prong 2
910  // values convention -1 = discarded
911  // 0 = not identified (but compatible) || No PID (->hasPID flag)
912  // 1 = identified
913  // PID search: pion (TPC) or not K (TOF), Kaon hypothesis for both
914  // Initial hypothesis: unknwon (but compatible)
915  isKaonPionTOF[0][0]=0;
916  isKaonPionTOF[0][1]=0;
917  isKaonPionTOF[1][0]=0;
918  isKaonPionTOF[1][1]=0;
919 
920  isKaonPionTPC[0][0]=0;
921  isKaonPionTPC[0][1]=0;
922  isKaonPionTPC[1][0]=0;
923  isKaonPionTPC[1][1]=0;
924 
925  combinedPID[0][0]=0;
926  combinedPID[0][1]=0;
927  combinedPID[1][0]=0;
928  combinedPID[1][1]=0;
929 
930 
931 
932  for(Int_t daught=0;daught<2;daught++){
933  //Loop con prongs
934 
935  // ########### Step 0- CHECKING minimal PID "ACCEPTANCE" ####################
936 
937  AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);
938 
939  if(!(aodtrack->GetStatus()&AliESDtrack::kTPCrefit)){
940  fWhyRejection=26;
941  return 0;
942  }
943  if(!(aodtrack->GetStatus()&AliESDtrack::kITSrefit)){
944  fWhyRejection=27;
945  return 0;
946  }
947 
948  AliAODPid *pid=aodtrack->GetDetPid();
949  if(!pid) {
950  //delete esdtrack;
951  hasPID[daught]--;
952  continue;
953  }
954 
955  // ########### Step 1- Check of TPC and TOF response ####################
956 
957  Double_t ptrack=aodtrack->P();
958  //#################### TPC PID #######################
959  if (!(aodtrack->GetStatus()&AliESDtrack::kTPCpid )){
960  // NO TPC PID INFO FOR THIS TRACK
961  hasPID[daught]--;
962  }
963  else {
964  static AliTPCPIDResponse theTPCpid;
965  AliAODPid *pidObj = aodtrack->GetDetPid();
966  Double_t ptProng=pidObj->GetTPCmomentum();
967  nsigmaTPCpi = theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kPion);
968  nsigmaTPCK = theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kKaon);
969  //if(ptrack<0.6){
970  if(ptProng<0.6){
971  if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;
972  else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
973  if(TMath::Abs(nsigmaTPCpi)<2.)isKaonPionTPC[daught][1]=1;
974  else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
975  }
976  //else if(ptrack<.8){
977  else if(ptProng<.8){
978  if(TMath::Abs(nsigmaTPCK)<1.)isKaonPionTPC[daught][0]=1;
979  else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
980  if(TMath::Abs(nsigmaTPCpi)<1.)isKaonPionTPC[daught][1]=1;
981  else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
982  }
983  else {
984  // if(nsigmaTPCK>-2.&&nsigmaTPCK<1.)isKaonPionTPC[daught][0]=1;
985  if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
986  //if(nsigmaTPCpi>-1.&&nsigmaTPCpi<2.)isKaonPionTPC[daught][1]=1;
987  if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
988  }
989  }
990 
991 
992  // ##### TOF PID: do not ask nothing for pion/protons ############
993  if(!((aodtrack->GetStatus()&AliESDtrack::kTOFpid)&&(aodtrack->GetStatus()&AliESDtrack::kTOFout)&&(aodtrack->GetStatus()&AliESDtrack::kTIME))){
994  // NO TOF PID INFO FOR THIS TRACK
995  hasPID[daught]--;
996  }
997  else{
998  tofSig=pid->GetTOFsignal();
999  pid->GetIntegratedTimes(times);
1000  if((tofSig-times[3])>5.*160.)return 0;// PROTON REJECTION
1001  if(TMath::Abs(tofSig-times[3])>3.*160.){
1002  isKaonPionTOF[daught][0]=-1;
1003  }
1004  else {
1005  if(ptrack<1.5){
1006  isKaonPionTOF[daught][0]=1;
1007  }
1008  }
1009  }
1010 
1011  //######### Step 2: COMBINE TOF and TPC PID ###############
1012  // we apply the following convention: if TPC and TOF disagree (discarded Vs identified) -> unknown
1013  combinedPID[daught][0]=isKaonPionTOF[daught][0]+isKaonPionTPC[daught][0];
1014  combinedPID[daught][1]=isKaonPionTOF[daught][1]+isKaonPionTPC[daught][1];
1015 
1016 
1017  //######### Step 3: USE PID INFO
1018 
1019  if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){// if not a K- and not a pi- both D0 and D0bar excluded
1020  isD0D0barPID[0]=0;
1021  isD0D0barPID[1]=0;
1022  }
1023  else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){// if in conflict (both pi- and K-), if k for both TPC and TOF -> is K
1024  if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded
1025  else isD0D0barPID[0]=0;// if K+ D0 excluded
1026  }
1027  else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){// if in conflict (both pi- and K-) and k- only for TPC or TOF -> reject
1028  isD0D0barPID[0]=0;
1029  isD0D0barPID[1]=0;
1030  }
1031  else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){
1032  if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded
1033  else isD0D0barPID[0]=0;// not a D0 if K+ or if pi+ excluded
1034  }
1035  else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){
1036  if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded
1037  else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded
1038  }
1039 
1040  // ########## ALSO DIFFERENT TPC PID REQUEST FOR LOW pt D0: request of K identification ###############################
1041  // ########## more tolerant criteria for single particle ID-> more selective criteria for D0 ##############################
1042  // ############### NOT OPTIMIZED YET ###################################
1043  if(d->Pt()<2.){
1044  isKaonPionTPC[daught][0]=0;
1045  isKaonPionTPC[daught][1]=0;
1046  AliAODPid *pidObj = aodtrack->GetDetPid();
1047  Double_t ptProng=pidObj->GetTPCmomentum();
1048  //if(ptrack<0.6){
1049  if(ptProng<0.6){
1050  if(TMath::Abs(nsigmaTPCK)<3.)isKaonPionTPC[daught][0]=1;
1051  else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
1052  if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;
1053  else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
1054  }
1055  //else if(ptrack<.8){
1056  else if(ptProng<.8){
1057  if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;
1058  else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
1059  if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;
1060  else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
1061  }
1062  else {
1063  if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
1064  if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
1065  }
1066  }
1067 
1068  }// END OF LOOP ON DAUGHTERS
1069 
1070  // FURTHER PID REQUEST (both daughter info is needed)
1071  if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){
1072  fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found
1073  return 0;
1074  }
1075  else if(hasPID[0]==0&&hasPID[1]==0){
1076  fWhyRejection=28;// reject cases in which no PID info is available
1077  return 0;
1078  }
1079  if(d->Pt()<2.){
1080  // request of K identification at low D0 pt
1081  combinedPID[0][0]=0;
1082  combinedPID[0][1]=0;
1083  combinedPID[1][0]=0;
1084  combinedPID[1][1]=0;
1085 
1086  combinedPID[0][0]=isKaonPionTOF[0][0]+isKaonPionTPC[0][0];
1087  combinedPID[0][1]=isKaonPionTOF[0][1]+isKaonPionTPC[0][1];
1088  combinedPID[1][0]=isKaonPionTOF[1][0]+isKaonPionTPC[1][0];
1089  combinedPID[1][1]=isKaonPionTOF[1][1]+isKaonPionTPC[1][1];
1090 
1091  if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
1092  fWhyRejection=32;// reject cases where the Kaon is not identified
1093  return 0;
1094  }
1095  }
1096 
1097  // cout<<"Why? "<<fWhyRejection<<endl;
1098  return isD0D0barPID[0]+isD0D0barPID[1];
1099 }
1100 
1101 
1102 
1103 //---------------------------------------------------------------------------
1104 Int_t AliRDHFCutsD0toKpi::CombineSelectionLevels(Int_t selectionvalTrack,
1105  Int_t selectionvalCand,
1106  Int_t selectionvalPID) const
1107 {
1108  //
1109  // This method combines the tracks, PID and cuts selection results
1110  //
1111  if(selectionvalTrack==0) return 0;
1112 
1113  Int_t returnvalue;
1114 
1115  switch(selectionvalPID) {
1116  case 0:
1117  returnvalue=0;
1118  break;
1119  case 1:
1120  returnvalue=((selectionvalCand==1 || selectionvalCand==3) ? 1 : 0);
1121  break;
1122  case 2:
1123  returnvalue=((selectionvalCand==2 || selectionvalCand==3) ? 2 : 0);
1124  break;
1125  case 3:
1126  returnvalue=selectionvalCand;
1127  break;
1128  default:
1129  returnvalue=0;
1130  break;
1131  }
1132 
1133  return returnvalue;
1134 }
1135 //----------------------------------------------------------------------------
1137 {
1138  //
1139  // Note: this method is temporary
1140  // Additional cuts on decay lenght and lower cut for d0 norm are applied using vertex without candidate's daughters
1141  //
1142 
1143  //apply cuts
1144 
1145  Float_t normDecLengthCut=1.,decLengthCut=TMath::Min(d->P()*0.0066+0.01,0.06/*cm*/), normd0Cut=0.5;
1146  // "decay length" expo law with tau' = beta*gamma*ctau= p/m*ctau =p*0.0123/1.864~p*0.0066
1147  // decay lenght > ctau' implies to retain (1-1/e) (for signal without considering detector resolution),
1148 
1149  Int_t returnvalue=3; //cut passed
1150  for(Int_t i=0;i<2/*prongs*/;i++){
1151  if(TMath::Abs(d->Normalizedd0Prong(i))<normd0Cut) return 0; //normd0Cut not passed
1152  }
1153  if(d->DecayLength2()<decLengthCut*decLengthCut) return 0; //decLengthCut not passed
1154  if(d->NormalizedDecayLength2()<normDecLengthCut*normDecLengthCut) return 0; //decLengthCut not passed
1155 
1156  return returnvalue;
1157 }
1158 
1159 //----------------------------------------------
1161 {
1162  //switch on candidate selection via AliKFparticle
1163  if(!useKF) return;
1164  if(useKF){
1165  fUseKF=useKF;
1166  Int_t nvarsKF=11;
1167  SetNVars(nvarsKF);
1168  TString varNamesKF[11]={"inv. mass [GeV]",
1169  "dca [cm]",
1170  "cosThetaStar",
1171  "pTK [GeV/c]",
1172  "pTPi [GeV/c]",
1173  "d0K [cm]",
1174  "d0Pi [cm]",
1175  "d0d0 [cm^2]",
1176  "cosThetaPoint"
1177  "DecayLength[cm]",
1178  "RedChi2"};
1179  Bool_t isUpperCutKF[11]={kTRUE,
1180  kTRUE,
1181  kTRUE,
1182  kFALSE,
1183  kFALSE,
1184  kTRUE,
1185  kTRUE,
1186  kTRUE,
1187  kFALSE,
1188  kFALSE,
1189  kTRUE};
1190  SetVarNames(nvarsKF,varNamesKF,isUpperCutKF);
1191  SetGlobalIndex();
1192  Bool_t forOpt[11]={kFALSE,
1193  kTRUE,
1194  kTRUE,
1195  kFALSE,
1196  kFALSE,
1197  kFALSE,
1198  kFALSE,
1199  kTRUE,
1200  kTRUE,
1201  kFALSE,
1202  kFALSE};
1203  SetVarsForOpt(4,forOpt);
1204  }
1205  return;
1206 }
1207 
1208 
1209 //---------------------------------------------------------------------------
1211  //
1212  //STANDARD CUTS USED FOR 2010 pp analysis
1213  //dca cut will be enlarged soon to 400 micron
1214  //
1215 
1216  SetName("D0toKpiCutsStandard");
1217  SetTitle("Standard Cuts for D0 analysis");
1218 
1219  // PILE UP REJECTION
1221 
1222  // EVENT CUTS
1223  SetMinVtxContr(1);
1224  // MAX Z-VERTEX CUT
1225  SetMaxVtxZ(10.);
1226 
1227  // TRACKS ON SINGLE TRACKS
1228  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
1229  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1230  esdTrackCuts->SetRequireTPCRefit(kTRUE);
1231  esdTrackCuts->SetRequireITSRefit(kTRUE);
1232  // esdTrackCuts->SetMinNClustersITS(4);
1233  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
1234  esdTrackCuts->SetMinDCAToVertexXY(0.);
1235  esdTrackCuts->SetEtaRange(-0.8,0.8);
1236  esdTrackCuts->SetPtRange(0.3,1.e10);
1237 
1238  AddTrackCuts(esdTrackCuts);
1239  delete esdTrackCuts;
1240  esdTrackCuts=NULL;
1241 
1242  const Int_t nptbins =14;
1243  const Double_t ptmax = 9999.;
1244  const Int_t nvars=11;
1245  Float_t ptbins[nptbins+1];
1246  ptbins[0]=0.;
1247  ptbins[1]=0.5;
1248  ptbins[2]=1.;
1249  ptbins[3]=2.;
1250  ptbins[4]=3.;
1251  ptbins[5]=4.;
1252  ptbins[6]=5.;
1253  ptbins[7]=6.;
1254  ptbins[8]=7.;
1255  ptbins[9]=8.;
1256  ptbins[10]=12.;
1257  ptbins[11]=16.;
1258  ptbins[12]=20.;
1259  ptbins[13]=24.;
1260  ptbins[14]=ptmax;
1261 
1262  SetGlobalIndex(nvars,nptbins);
1263  SetPtBins(nptbins+1,ptbins);
1264 
1265  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* pt<0.5*/
1266  {0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* 0.5<pt<1*/
1267  {0.400,300.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.80,0.,0.},/* 1<pt<2 */
1268  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.,0.},/* 2<pt<3 */
1269  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 3<pt<4 */
1270  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 4<pt<5 */
1271  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 5<pt<6 */
1272  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 6<pt<7 */
1273  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-7000.*1E-8,0.85,0.,0.},/* 7<pt<8 */
1274  {0.400,300.*1E-4,0.9,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.,0.},/* 8<pt<12 */
1275  {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,10000.*1E-8,0.85,0.,0.},/* 12<pt<16 */
1276  {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 16<pt<20 */
1277  {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 20<pt<24 */
1278  {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.}};/* pt>24 */
1279 
1280 
1281  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1282  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1283  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1284 
1285  for (Int_t ibin=0;ibin<nptbins;ibin++){
1286  for (Int_t ivar = 0; ivar<nvars; ivar++){
1287  cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1288  }
1289  }
1290 
1291  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1292  SetUseSpecialCuts(kTRUE);
1294 
1295  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1296  delete [] cutsMatrixTransposeStand;
1297  cutsMatrixTransposeStand=NULL;
1298 
1299  // PID SETTINGS
1300  AliAODPidHF* pidObj=new AliAODPidHF();
1301  //pidObj->SetName("pid4D0");
1302  Int_t mode=1;
1303  const Int_t nlims=2;
1304  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1305  Bool_t compat=kTRUE; //effective only for this mode
1306  Bool_t asym=kTRUE;
1307  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1308  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1309  pidObj->SetMatch(mode);
1310  pidObj->SetPLimit(plims,nlims);
1311  pidObj->SetSigma(sigmas);
1312  pidObj->SetCompat(compat);
1313  pidObj->SetTPC(kTRUE);
1314  pidObj->SetTOF(kTRUE);
1315  pidObj->SetPCompatTOF(1.5);
1316  pidObj->SetSigmaForTPCCompat(3.);
1317  pidObj->SetSigmaForTOFCompat(3.);
1318  pidObj->SetOldPid(kTRUE);
1319 
1320  SetPidHF(pidObj);
1321  SetUsePID(kTRUE);
1322  SetUseDefaultPID(kFALSE);
1323 
1324  SetLowPt(kFALSE);
1325 
1326  PrintAll();
1327 
1328  delete pidObj;
1329  pidObj=NULL;
1330 
1331  return;
1332 
1333 }
1334 
1335 //___________________________________________________________________________
1337  //
1338  // Cuts for 2010 pp 7 TeV data analysis in Ntracklets bins (vs multiplicity)
1339  //
1340  SetName("D0toKpiCuts");
1341  SetTitle("Cuts for D0 analysis in 2010-data pp 7 TeV vs multiplicity");
1342 
1343  //
1344  // Track cuts
1345  //
1346  AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
1347  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1348  //default
1349  esdTrackCuts->SetRequireTPCRefit(kTRUE);
1350  esdTrackCuts->SetRequireITSRefit(kTRUE);
1351  esdTrackCuts->SetEtaRange(-0.8,0.8);
1352  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
1353  AliESDtrackCuts::kAny);
1354  // default is kBoth, otherwise kAny
1355  esdTrackCuts->SetMinDCAToVertexXY(0.);
1356  esdTrackCuts->SetPtRange(0.3,1.e10);
1357 
1358  AddTrackCuts(esdTrackCuts);
1359  delete esdTrackCuts;
1360  esdTrackCuts=NULL;
1361 
1362 
1363  //
1364  // Cut values per pt bin
1365  //
1366  const Int_t nvars=11;
1367  const Int_t nptbins=14;
1368  Float_t* ptbins;
1369  ptbins=new Float_t[nptbins+1];
1370  ptbins[0]=0.;
1371  ptbins[1]=0.5;
1372  ptbins[2]=1.;
1373  ptbins[3]=2.;
1374  ptbins[4]=3.;
1375  ptbins[5]=4.;
1376  ptbins[6]=5.;
1377  ptbins[7]=6.;
1378  ptbins[8]=7.;
1379  ptbins[9]=8.;
1380  ptbins[10]=12.;
1381  ptbins[11]=16.;
1382  ptbins[12]=20.;
1383  ptbins[13]=24.;
1384  ptbins[14]=9999.;
1385 
1386  SetPtBins(nptbins+1,ptbins);
1387 
1388  //setting cut values
1389  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,350.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,3.},/* pt<0.5*/
1390  {0.400,350.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,3.},/* 0.5<pt<1*/
1391  {0.400,300.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-30000.*1E-8,0.80,0.,4.},/* 1<pt<2 */
1392  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-20000.*1E-8,0.85,0.,4.},/* 2<pt<3 */
1393  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.,4.},/* 3<pt<4 */
1394  {0.400,300.*1E-4,0.75,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.875,0.,0.},/* 4<pt<5 */
1395  {0.400,300.*1E-4,0.75,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.875,0.,0.},/* 5<pt<6 */
1396  {0.400,400.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 6<pt<7 */
1397  {0.400,400.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 7<pt<8 */
1398  {0.400,0.06,0.85,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.00001,0.85,0.,0.},/* 8<pt<12 */
1399  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 12<pt<16 */
1400  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 16<pt<20 */
1401  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 20<pt<24 */
1402  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.}};/* pt>24 */
1403 
1404 
1405  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1406  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1407  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1408 
1409  for (Int_t ibin=0;ibin<nptbins;ibin++){
1410  for (Int_t ivar = 0; ivar<nvars; ivar++){
1411  cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1412  }
1413  }
1414 
1415  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1416  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1417  delete [] cutsMatrixTransposeStand;
1418 
1419  SetUseSpecialCuts(kTRUE);
1421 
1422  //Do recalculate the vertex
1424 
1425  //
1426  // Pid settings
1427  //
1428  Bool_t pidflag=kTRUE;
1429  SetUsePID(pidflag);
1430  if(pidflag) cout<<"PID is used"<<endl;
1431  else cout<<"PID is not used"<<endl;
1432  //
1433  AliAODPidHF* pidObj=new AliAODPidHF();
1434  Int_t mode=1;
1435  const Int_t nlims=2;
1436  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1437  Bool_t compat=kTRUE; //effective only for this mode
1438  Bool_t asym=kTRUE;
1439  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1440  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1441  pidObj->SetMatch(mode);
1442  pidObj->SetPLimit(plims,nlims);
1443  pidObj->SetSigma(sigmas);
1444  pidObj->SetCompat(compat);
1445  pidObj->SetTPC(kTRUE);
1446  pidObj->SetTOF(kTRUE);
1447  pidObj->SetPCompatTOF(1.5);
1448  pidObj->SetSigmaForTPCCompat(3.);
1449  pidObj->SetSigmaForTOFCompat(3.);
1450  pidObj->SetOldPid(kTRUE);
1451  // pidObj->SetOldPid(kFALSE);
1452  SetPidHF(pidObj);
1453 
1454  SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
1455 
1456  SetLowPt(kFALSE);
1457 
1458  //activate pileup rejection (for pp)
1460 
1461 
1462  PrintAll();
1463 
1464  delete pidObj;
1465  pidObj=NULL;
1466  delete [] ptbins;
1467  ptbins=NULL;
1468 
1469  return;
1470 }
1471 
1472 //---------------------------------------------------------------------------
1474  //
1475  // STANDARD CUTS USED FOR 2011 pp analysis at 2.76TeV
1476  //
1477 
1478  SetName("D0toKpiCutsStandard");
1479  SetTitle("Standard Cuts for D0 analysis in pp2011 at 2.76TeV run");
1480 
1481  //
1482  // Track cuts
1483  //
1484  AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
1485  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1486  //default
1487  esdTrackCuts->SetRequireTPCRefit(kTRUE);
1488  esdTrackCuts->SetRequireITSRefit(kTRUE);
1489  esdTrackCuts->SetEtaRange(-0.8,0.8);
1490  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
1491  AliESDtrackCuts::kAny);
1492  // default is kBoth, otherwise kAny
1493  esdTrackCuts->SetMinDCAToVertexXY(0.);
1494  esdTrackCuts->SetPtRange(0.3,1.e10);
1495 
1496  esdTrackCuts->SetMaxDCAToVertexXY(1.);
1497  esdTrackCuts->SetMaxDCAToVertexZ(1.);
1498  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0075*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");
1499 
1500  AddTrackCuts(esdTrackCuts);
1501  delete esdTrackCuts;
1502  esdTrackCuts=NULL;
1503 
1504 
1505  const Int_t nvars=11;
1506  const Int_t nptbins=13;
1507  Float_t ptbins[nptbins+1];
1508  ptbins[0]=0.;
1509  ptbins[1]=0.5;
1510  ptbins[2]=1.;
1511  ptbins[3]=2.;
1512  ptbins[4]=3.;
1513  ptbins[5]=4.;
1514  ptbins[6]=5.;
1515  ptbins[7]=6.;
1516  ptbins[8]=8.;
1517  ptbins[9]=12.;
1518  ptbins[10]=16.;
1519  ptbins[11]=20.;
1520  ptbins[12]=24.;
1521  ptbins[13]=9999.;
1522 
1523  SetPtBins(nptbins+1,ptbins);
1524 
1525 
1526  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,0.04,0.75,0.3,0.3,1000.*1E-4,1000.*1E-4,0.,0.85,0.,0.},/* pt<0.5*/
1527  {0.400,0.04,0.75,0.3,0.3,1000.*1E-4,1000.*1E-4,0.,0.85,0.,0.},/* 0.5<pt<1*/
1528  {0.400,0.03,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.8,0.,0.},/* 1<pt<2 */
1529  {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.0003,0.9,0.,0.},/* 2<pt<3 */
1530  {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.0002,0.9,0.,0.},/* 3<pt<4 */
1531  {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.00015,0.9,0.,0.},/* 4<pt<5 */
1532  {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.0001,0.9,0.,0.},/* 5<pt<6 */
1533  {0.400,0.09,0.85,0.7,0.7,1000.*1E-4,1000.*1E-4,0.,0.85,0.,0.},/* 6<pt<8 */
1534  {0.400,0.06,0.85,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.00001,0.85,0.,0.},/* 8<pt<12 */
1535  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 12<pt<16 */
1536  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 16<pt<20 */
1537  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 20<pt<24 */
1538  {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.}};/* pt>24 */
1539 
1540  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1541  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1542  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1543  for (Int_t ibin=0;ibin<nptbins;ibin++){
1544  for (Int_t ivar = 0; ivar<nvars; ivar++){
1545  cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1546  }
1547  }
1548  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1549  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1550  delete [] cutsMatrixTransposeStand;
1551 
1552 
1553  //pid settings
1554  AliAODPidHF* pidObj=new AliAODPidHF();
1555  Int_t mode=1;
1556  const Int_t nlims=2;
1557  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1558  Bool_t compat=kTRUE; //effective only for this mode
1559  Bool_t asym=kTRUE;
1560  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1561  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1562  pidObj->SetMatch(mode);
1563  pidObj->SetPLimit(plims,nlims);
1564  pidObj->SetSigma(sigmas);
1565  pidObj->SetCompat(compat);
1566  pidObj->SetTPC(kTRUE);
1567  pidObj->SetTOF(kTRUE);
1568  pidObj->SetOldPid(kTRUE);
1569  SetPidHF(pidObj);
1570 
1571  SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
1572 
1573  SetUsePID(kTRUE);
1574 
1575  //activate pileup rejection (for pp)
1577 
1578  //Do recalculate the vertex
1580 
1581  // Cut on the zvtx
1582  SetMaxVtxZ(10.);
1583 
1584  // Use the kFirst selection for tracks with candidate D with pt<2
1585  SetSelectCandTrackSPDFirst(kTRUE,2.);
1586 
1587  // Use special cuts for D candidates with pt<2
1588  SetUseSpecialCuts(kTRUE);
1590 
1591  PrintAll();
1592 
1593  delete pidObj;
1594  pidObj=NULL;
1595 
1596  return;
1597 }
1598 
1599 
1600 //---------------------------------------------------------------------------
1602  //
1603  // Final CUTS USED FOR 2010 PbPb analysis of central events
1604  // REMEMBER TO EVENTUALLY SWITCH ON
1605  // SetUseAOD049(kTRUE);
1606  //
1607 
1608  SetName("D0toKpiCutsStandard");
1609  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run");
1610 
1611 
1612  // CENTRALITY SELECTION
1613  SetMinCentrality(0.);
1614  SetMaxCentrality(80.);
1616 
1617 
1618 
1619  // EVENT CUTS
1620  SetMinVtxContr(1);
1621  // MAX Z-VERTEX CUT
1622  SetMaxVtxZ(10.);
1623  // PILE UP
1625  // PILE UP REJECTION
1626  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
1627 
1628  // TRACKS ON SINGLE TRACKS
1629  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
1630  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1631  esdTrackCuts->SetRequireTPCRefit(kTRUE);
1632  esdTrackCuts->SetRequireITSRefit(kTRUE);
1633  // esdTrackCuts->SetMinNClustersITS(4);
1634  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
1635  esdTrackCuts->SetMinDCAToVertexXY(0.);
1636  esdTrackCuts->SetEtaRange(-0.8,0.8);
1637  esdTrackCuts->SetPtRange(0.7,1.e10);
1638 
1639  esdTrackCuts->SetMaxDCAToVertexXY(1.);
1640  esdTrackCuts->SetMaxDCAToVertexZ(1.);
1641  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0075*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");
1642 
1643 
1644  AddTrackCuts(esdTrackCuts);
1645  delete esdTrackCuts;
1646  esdTrackCuts=NULL;
1647  // SPD k FIRST for Pt<3 GeV/c
1648  SetSelectCandTrackSPDFirst(kTRUE, 3);
1649 
1650  // CANDIDATE CUTS
1651  const Int_t nptbins =13;
1652  const Double_t ptmax = 9999.;
1653  const Int_t nvars=11;
1654  Float_t ptbins[nptbins+1];
1655  ptbins[0]=0.;
1656  ptbins[1]=0.5;
1657  ptbins[2]=1.;
1658  ptbins[3]=2.;
1659  ptbins[4]=3.;
1660  ptbins[5]=4.;
1661  ptbins[6]=5.;
1662  ptbins[7]=6.;
1663  ptbins[8]=8.;
1664  ptbins[9]=12.;
1665  ptbins[10]=16.;
1666  ptbins[11]=20.;
1667  ptbins[12]=24.;
1668  ptbins[13]=ptmax;
1669 
1670  SetGlobalIndex(nvars,nptbins);
1671  SetPtBins(nptbins+1,ptbins);
1672  SetMinPtCandidate(2.);
1673 
1674  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.85,0.99,2.},/* pt<0.5*/
1675  {0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-35000.*1E-8,0.9,0.99,2.},/* 0.5<pt<1*/
1676  {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-43000.*1E-8,0.85,0.995,8.},/* 1<pt<2 */
1677  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-45000.*1E-8,0.95,0.998,7.},/* 2<pt<3 */
1678  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-36000.*1E-8,0.95,0.998,5.},/* 3<pt<4 */
1679  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-27000.*1E-8,0.95,0.998,5.},/* 4<pt<5 */
1680  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-21000.*1E-8,0.92,0.998,5.},/* 5<pt<6 */
1681  {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-14000.*1E-8,0.88,0.998,5.},/* 6<pt<8 */
1682  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.998,5.},/* 8<pt<12 */
1683  {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.83,0.998,8.},/* 12<pt<16 */
1684  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.82,0.995,6.},/* 16<pt<20 */
1685  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.81,0.995,6.},/* 20<pt<24 */
1686  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.8,0.99,2.}};/* pt>24 */
1687 
1688 
1689  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1690  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1691  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1692 
1693  for (Int_t ibin=0;ibin<nptbins;ibin++){
1694  for (Int_t ivar = 0; ivar<nvars; ivar++){
1695  cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1696  }
1697  }
1698 
1699  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1700  SetUseSpecialCuts(kTRUE);
1701  SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb
1702  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1703  delete [] cutsMatrixTransposeStand;
1704  cutsMatrixTransposeStand=NULL;
1705 
1706  // PID SETTINGS
1707  AliAODPidHF* pidObj=new AliAODPidHF();
1708  //pidObj->SetName("pid4D0");
1709  Int_t mode=1;
1710  const Int_t nlims=2;
1711  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1712  Bool_t compat=kTRUE; //effective only for this mode
1713  Bool_t asym=kTRUE;
1714  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1715  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1716  pidObj->SetMatch(mode);
1717  pidObj->SetPLimit(plims,nlims);
1718  pidObj->SetSigma(sigmas);
1719  pidObj->SetCompat(compat);
1720  pidObj->SetTPC(kTRUE);
1721  pidObj->SetTOF(kTRUE);
1722  pidObj->SetPCompatTOF(2.);
1723  pidObj->SetSigmaForTPCCompat(3.);
1724  pidObj->SetSigmaForTOFCompat(3.);
1725  pidObj->SetOldPid(kTRUE);
1726 
1727 
1728  SetPidHF(pidObj);
1729  SetUsePID(kTRUE);
1730  SetUseDefaultPID(kFALSE);
1731 
1732 
1733  PrintAll();
1734 
1735 
1736  delete pidObj;
1737  pidObj=NULL;
1738 
1739  return;
1740 
1741 }
1742 
1743 //---------------------------------------------------------------------------
1745  // CUTS USED FOR D0 RAA for the CENTRALITY RANGE 20-80%
1746 
1747 
1748  SetName("D0toKpiCutsStandard");
1749  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run, for peripheral events");
1750 
1751 
1752  // CENTRALITY SELECTION
1753  SetMinCentrality(40.);
1754  SetMaxCentrality(80.);
1756 
1757  // EVENT CUTS
1758  SetMinVtxContr(1);
1759  // MAX Z-VERTEX CUT
1760  SetMaxVtxZ(10.);
1761  // PILE UP
1763  // PILE UP REJECTION
1764  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
1765 
1766  // TRACKS ON SINGLE TRACKS
1767  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
1768  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1769  esdTrackCuts->SetRequireTPCRefit(kTRUE);
1770  esdTrackCuts->SetRequireITSRefit(kTRUE);
1771  // esdTrackCuts->SetMinNClustersITS(4);
1772  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
1773  esdTrackCuts->SetMinDCAToVertexXY(0.);
1774  esdTrackCuts->SetEtaRange(-0.8,0.8);
1775  esdTrackCuts->SetPtRange(0.5,1.e10);
1776 
1777  esdTrackCuts->SetMaxDCAToVertexXY(1.);
1778  esdTrackCuts->SetMaxDCAToVertexZ(1.);
1779  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0025*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");
1780 
1781 
1782  AddTrackCuts(esdTrackCuts);
1783  delete esdTrackCuts;
1784  esdTrackCuts=NULL;
1785  // SPD k FIRST for Pt<3 GeV/c
1786  SetSelectCandTrackSPDFirst(kTRUE, 3);
1787 
1788  // CANDIDATE CUTS
1789  const Int_t nptbins =13;
1790  const Double_t ptmax = 9999.;
1791  const Int_t nvars=11;
1792  Float_t ptbins[nptbins+1];
1793  ptbins[0]=0.;
1794  ptbins[1]=0.5;
1795  ptbins[2]=1.;
1796  ptbins[3]=2.;
1797  ptbins[4]=3.;
1798  ptbins[5]=4.;
1799  ptbins[6]=5.;
1800  ptbins[7]=6.;
1801  ptbins[8]=8.;
1802  ptbins[9]=12.;
1803  ptbins[10]=16.;
1804  ptbins[11]=20.;
1805  ptbins[12]=24.;
1806  ptbins[13]=ptmax;
1807 
1808  SetGlobalIndex(nvars,nptbins);
1809  SetPtBins(nptbins+1,ptbins);
1810  SetMinPtCandidate(0.);
1811 
1812  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-20000.*1E-8,0.7,0.993,2.},/* pt<0.5*/
1813  {0.400,400.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.85,0.993,2.},/* 0.5<pt<1*/
1814  {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.85,0.995,6.},/* 1<pt<2 */
1815  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.95,0.991,5.},/* 2<pt<3 */
1816  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-36000.*1E-8,0.95,0.993,5.},/* 3<pt<4 */
1817  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-27000.*1E-8,0.95,0.995,5.},/* 4<pt<5 */
1818  {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-21000.*1E-8,0.92,0.998,5.},/* 5<pt<6 */
1819  {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-14000.*1E-8,0.88,0.998,5.},/* 6<pt<8 */
1820  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.995,5.},/* 8<pt<12 */
1821  {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.83,0.995,4.},/* 12<pt<16 */
1822  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.82,0.995,4.},/* 16<pt<20 */
1823  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.81,0.995,4.},/* 20<pt<24 */
1824  {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.8,0.995,4.}};/* pt>24 */
1825 
1826 
1827  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1828  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1829  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1830 
1831  for (Int_t ibin=0;ibin<nptbins;ibin++){
1832  for (Int_t ivar = 0; ivar<nvars; ivar++){
1833  cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1834  }
1835  }
1836 
1837  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1838  SetUseSpecialCuts(kTRUE);
1839  SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb
1840  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1841  delete [] cutsMatrixTransposeStand;
1842  cutsMatrixTransposeStand=NULL;
1843 
1844  // PID SETTINGS
1845  AliAODPidHF* pidObj=new AliAODPidHF();
1846  //pidObj->SetName("pid4D0");
1847  Int_t mode=1;
1848  const Int_t nlims=2;
1849  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1850  Bool_t compat=kTRUE; //effective only for this mode
1851  Bool_t asym=kTRUE;
1852  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1853  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1854  pidObj->SetMatch(mode);
1855  pidObj->SetPLimit(plims,nlims);
1856  pidObj->SetSigma(sigmas);
1857  pidObj->SetCompat(compat);
1858  pidObj->SetTPC(kTRUE);
1859  pidObj->SetTOF(kTRUE);
1860  pidObj->SetPCompatTOF(2.);
1861  pidObj->SetSigmaForTPCCompat(3.);
1862  pidObj->SetSigmaForTOFCompat(3.);
1863  pidObj->SetOldPid(kTRUE);
1864 
1865  SetPidHF(pidObj);
1866  SetUsePID(kTRUE);
1867  SetUseDefaultPID(kFALSE);
1868 
1869  SetLowPt(kTRUE,2.);
1870 
1871  PrintAll();
1872 
1873 
1874  delete pidObj;
1875  pidObj=NULL;
1876 
1877  return;
1878 
1879 }
1880 
1881 
1882 //---------------------------------------------------------------------------
1884 
1885  // Default 2010 PbPb cut object
1887  AliAODPidHF *pidobj=GetPidHF();
1888 
1889  pidobj->SetOldPid(kFALSE);
1890 
1891  //
1892  // Enable all 2011 PbPb run triggers
1893  //
1894  SetTriggerClass("");
1898 
1899 }
1900 
1901 //---------------------------------------------------------------------------
1903 {
1904  // ############################################################
1905  //
1906  // Apply Bayesian PID selection
1907  // Implementation of Bayesian PID: Jeremy Wilkinson
1908  //
1909  // ############################################################
1910 
1911  if (!fUsePID || !d) return 3;
1912 
1913 
1915  //WeightNoFilter: Accept all particles (no PID cut) but fill mass histos with weights in task
1917  return 3;
1918  }
1919 
1920 
1921 
1922  Int_t isD0 = 0;
1923  Int_t isD0bar = 0;
1924  Int_t returnvalue = 0;
1925 
1926  Int_t isPosKaon = 0, isNegKaon = 0, isPosPion = 0, isNegPion = 0;
1927 
1928  //Bayesian methods used here check for ID of kaon, and whether it is positive or negative.
1929 
1930 
1931  Bool_t checkPIDInfo[2] = {kTRUE, kTRUE};
1932  AliAODTrack *aodtrack[2] = {(AliAODTrack*)d->GetDaughter(0), (AliAODTrack*)d->GetDaughter(1)};
1933 
1934  if ((aodtrack[0]->Charge()*aodtrack[1]->Charge()) != -1) {
1935  return 0; //reject if charges not opposite
1936  }
1937  Double_t momentumpositive=0., momentumnegative=0.; //Used in "prob > prior" method
1938  for (Int_t daught = 0; daught < 2; daught++) {
1939  //Loop over prongs
1940 
1941  if (aodtrack[daught]->Charge() == -1) {
1942  momentumnegative = aodtrack[daught]->P();
1943  }
1944  if (aodtrack[daught]->Charge() == +1) {
1945  momentumpositive = aodtrack[daught]->P();
1946  }
1947  if (!(fPidHF->CheckStatus(aodtrack[daught], "TPC")) && !(fPidHF->CheckStatus(aodtrack[daught], "TOF"))) {
1948  checkPIDInfo[daught] = kFALSE;
1949  continue;
1950  }
1951 
1952  }
1953 
1954  //Loop over daughters ends here
1955 
1956  if (!checkPIDInfo[0] && !checkPIDInfo[1]) {
1957  return 0; //Reject if both daughters lack both TPC and TOF info
1958  }
1959 
1960 
1961  CalculateBayesianWeights(d); //Calculates all Bayesian probabilities for both positive and negative tracks
1962  // Double_t prob0[AliPID::kSPECIES];
1963  // fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), prob0);
1965  switch (fBayesianCondition) {
1967  case kMaxProb:
1968  if (TMath::MaxElement(AliPID::kSPECIES, fWeightsPositive) == fWeightsPositive[AliPID::kKaon]) { //If highest probability lies with kaon
1969  isPosKaon = 1; //flag [daught] as a kaon
1970  }
1971 
1972  if (TMath::MaxElement(AliPID::kSPECIES, fWeightsPositive) == fWeightsPositive[AliPID::kPion]) { //If highest probability lies with pion
1973  isPosPion = 1; //flag [daught] as a pion
1974  }
1975 
1976  if (TMath::MaxElement(AliPID::kSPECIES, fWeightsNegative) == fWeightsNegative[AliPID::kKaon]) { //If highest probability lies with kaon
1977  isNegKaon = 1; //flag [daught] as a kaon
1978  }
1979 
1980  if (TMath::MaxElement(AliPID::kSPECIES, fWeightsNegative) == fWeightsNegative[AliPID::kPion]) { //If highest probability lies with kaon
1981  isNegPion = 1; //flag [daught] as a pion
1982  }
1983 
1984 
1985  break;
1987  case kAbovePrior:
1988 
1989  if (fWeightsNegative[AliPID::kKaon] > (fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->
1990  GetBinContent(fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->FindBin(momentumnegative)))) { //Retrieves relevant prior, gets value at momentum
1991  isNegKaon = 1;
1992  }
1993  if (fWeightsPositive[AliPID::kKaon] > (fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->
1994  GetBinContent(fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->FindBin(momentumpositive)))) { //Retrieves relevant prior, gets value at momentum
1995  isPosKaon = 1;
1996  }
1997 
1998  break;
1999 
2001  case kThreshold:
2002  if (fWeightsNegative[AliPID::kKaon] > fProbThreshold) {
2003  isNegKaon = 1;
2004  }
2005  if (fWeightsNegative[AliPID::kPion] > fProbThreshold) {
2006  isNegPion = 1;
2007  }
2008 
2009  if (fWeightsPositive[AliPID::kKaon] > fProbThreshold) {
2010  isPosKaon = 1;
2011  }
2012 
2013  if (fWeightsPositive[AliPID::kPion] > fProbThreshold) {
2014  isPosPion = 1;
2015  }
2016 
2017  break;
2018  }
2019 
2020 
2021  //Momentum-based selection (also applied to filtered weighted method)
2022 
2024  if (isNegKaon && isPosKaon) { // If both are kaons, reject
2025  isD0 = 1;
2026  isD0bar = 1;
2027  } else if (isNegKaon && isPosPion) { //If negative kaon present, D0
2028  isD0 = 1;
2029  } else if (isPosKaon && isNegPion) { //If positive kaon present, D0bar
2030  isD0bar = 1;
2031  } else { //If neither ID'd as kaon, subject to extra tests
2032  isD0 = 1;
2033  isD0bar = 1;
2034  if (aodtrack[0]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[0], "TOF"))) { //If it's low momentum and not ID'd as a kaon, we assume it must be a pion
2035  if (aodtrack[0]->Charge() == -1) {
2036  isD0 = 0; //Check charge and reject based on pion hypothesis
2037  }
2038  if (aodtrack[0]->Charge() == 1) {
2039  isD0bar = 0;
2040  }
2041  }
2042  if (aodtrack[1]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[1], "TOF"))) {
2043  if (aodtrack[1]->Charge() == -1) {
2044  isD0 = 0;
2045  }
2046  if (aodtrack[1]->Charge() == 1) {
2047  isD0bar = 0;
2048  }
2049  }
2050  }
2051 
2052 
2053 
2054  if (isD0 && isD0bar) {
2055  returnvalue = 3;
2056  }
2057  if (isD0 && !isD0bar) {
2058  returnvalue = 1;
2059  }
2060  if (!isD0 && isD0bar) {
2061  returnvalue = 2;
2062  }
2063  if (!isD0 && !isD0bar) {
2064  returnvalue = 0;
2065  }
2066  }
2067 
2068  //Simple Bayesian
2070 
2071  if (isPosKaon && isNegKaon) { //If both are ID'd as kaons, accept as possible
2072  returnvalue = 3;
2073  } else if (isNegKaon && isPosPion) { //If negative kaon, D0
2074  returnvalue = 1;
2075  } else if (isPosKaon && isNegPion) { //If positive kaon, D0-bar
2076  returnvalue = 2;
2077  } else if (isPosPion && isNegPion) {
2078  returnvalue = 0; //If neither kaon, reject
2079  } else {returnvalue = 0;} //default
2080 
2081  }
2082 
2083  return returnvalue;
2084 
2085 
2086 
2087 }
2088 
2089 
2090 
2091 //---------------------------------------------------------------------------
2093 
2094 {
2095  //Function to compute weights for Bayesian method
2096 
2097  AliAODTrack *aodtrack[2] = {(AliAODTrack*)d->GetDaughter(0), (AliAODTrack*)d->GetDaughter(1)};
2098  if ((aodtrack[0]->Charge() * aodtrack[1]->Charge()) != -1) {
2099  return; //Reject if charges do not oppose
2100  }
2101  for (Int_t daught = 0; daught < 2; daught++) {
2102  //Loop over prongs
2103 
2104  if (!(fPidHF->CheckStatus(aodtrack[daught], "TPC")) && !(fPidHF->CheckStatus(aodtrack[daught], "TOF"))) {
2105  continue;
2106  }
2107 
2108 
2109  // identify kaon, define weights
2110  if (aodtrack[daught]->Charge() == +1) {
2111  fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), fWeightsPositive);
2112  }
2113 
2114  if (aodtrack[daught]->Charge() == -1) {
2115  fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), fWeightsNegative);
2116  }
2117  }
2118 }
2119 
2120 
2121 
Double_t NormalizedDecayLengthXY() const
Int_t fIsSelectedCuts
fix the daughter track references
Definition: AliRDHFCuts.h:407
Double_t fPtMaxSpecialCuts
transverse momentum below which the strong PID is applied
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
void SetAsym(Bool_t asym)
Definition: AliAODPidHF.h:88
Int_t fnSpecies
switch for Bayesian
Int_t fWhyRejection
PID for heavy flavours manager.
Definition: AliRDHFCuts.h:391
void SetUseDefaultPID(Bool_t defPID)
Bool_t IsSignalMC(AliAODRecoDecay *d, AliAODEvent *aod, Int_t pdg) const
Bool_t IsTOFPiKexcluded(AliAODTrack *track, Double_t nsigmaK)
general method to perform PID using raw signals
Bool_t fUseKF
flag to switch on/off the default pid
Bool_t fUseMCVertex
flag to switch on the removal of duaghters from the primary vertex computation
Definition: AliRDHFCuts.h:394
void SetSigmaForTPC(Double_t *sigma)
Definition: AliAODPidHF.h:43
Double_t GetSigma(Int_t idet) const
Definition: AliAODPidHF.h:134
void SetMaxVtxZ(Float_t z=1e6)
Definition: AliRDHFCuts.h:60
Bool_t SetMCPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod) const
Bool_t fCombPID
max momentum for applying PID
Bool_t fRemoveDaughtersFromPrimary
Definition: AliRDHFCuts.h:393
void SetUseCentrality(Int_t flag=1)
virtual Int_t IsSelectedCombPID(AliAODRecoDecayHF *d)
Double_t NormalizedDecayLength2() const
Double_t fmaxPtrackForPID
transverse momentum below which the special cuts are applied
void SetPCompatTOF(Double_t pTOF)
Definition: AliAODPidHF.h:106
virtual void SetStandardCutsPP2011_276TeV()
void SetNVars(Int_t nVars)
Definition: AliRDHFCuts.h:361
virtual void SetStandardCutsPP2010()
Double_t CosPointingAngleXY() const
Double_t fMaxRapidityCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:411
void EnableSemiCentralTrigger()
Definition: AliRDHFCuts.h:93
AliRDHFCutsD0toKpi(const char *name="CutsD0toKpi")
void SetGlobalIndex()
Definition: AliRDHFCuts.h:196
AliRDHFCuts & operator=(const AliRDHFCuts &source)
Bool_t fLowPt
flag to switch on/off special cuts
virtual Bool_t IsInFiducialAcceptance(Double_t pt, Double_t y) const
Bool_t fUsePID
Definition: AliRDHFCuts.h:388
void SetSelectCandTrackSPDFirst(Bool_t flag, Double_t ptmax)
Flag and pt-maximum to check if the candidate daughters fulfill the kFirst criteria.
Definition: AliRDHFCuts.h:349
void SetPidHF(AliAODPidHF *pidObj)
see enum below
Definition: AliRDHFCuts.h:210
AliAODPidHF * GetPidHF() const
Definition: AliRDHFCuts.h:231
void SetTOF(Bool_t tof)
Definition: AliAODPidHF.h:95
Bool_t HasBadDaughters() const
Bool_t fDefaultPID
flag to switch on/off different pid for low pt D0
void SetMinVtxContr(Int_t contr=1)
Definition: AliRDHFCuts.h:58
Bool_t fUseTrackSelectionWithFilterBits
flag to reject kink daughters
Definition: AliRDHFCuts.h:421
void SetCuts(Int_t nVars, Int_t nPtBins, Float_t **cutsRD)
Int_t fnVarsForOpt
Definition: AliRDHFCuts.h:383
void SetUseKF(Bool_t useKF)
void ResetMaskAndEnableMBTrigger()
Definition: AliRDHFCuts.h:73
void SetMaximumPtSpecialCuts(Double_t pt)
void SetMinCentrality(Float_t minCentrality=0.)
Definition: AliRDHFCuts.h:51
AliPIDCombined * GetPidCombined() const
Definition: AliAODPidHF.h:161
Double_t CosThetaStarD0bar() const
angle of K
const Double_t ptmax
Double_t fMaxPtCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:410
Bool_t CheckStatus(AliAODTrack *track, TString detectors) const
Double_t sigmas[nPtBins]
AliAODVertex * GetOwnPrimaryVtx() const
Bool_t fKeepSignalMC
max rapidity of candidate (if !=-999 overrides IsInFiducialAcceptance)
Definition: AliRDHFCuts.h:412
Bool_t AreDaughtersSelected(AliAODRecoDecayHF *rd) const
virtual void CalculateBayesianWeights(AliAODRecoDecayHF *d)
Int_t fIsSelectedPID
outcome of cuts selection
Definition: AliRDHFCuts.h:408
Bool_t GetTOF() const
Definition: AliAODPidHF.h:142
Int_t MakeRawPid(AliAODTrack *track, Int_t specie)
Int_t mode
Definition: anaM.C:40
AliPIDResponse * GetPidResponse() const
Definition: AliAODPidHF.h:160
void SetSigma(Double_t *sigma)
Definition: AliAODPidHF.h:39
void SetSigmaForTOFCompat(Double_t sigma)
Definition: AliAODPidHF.h:45
Float_t * fCutsRD
fnVars*fnPtBins
Definition: AliRDHFCuts.h:386
Int_t IsSelectedKF(AliAODRecoDecayHF2Prong *d, AliAODEvent *aod) const
Double_t DecayLength2() const
kinematics & topology
void SetSigmaForTPCCompat(Double_t sigma)
Definition: AliAODPidHF.h:44
void SetMaxCentrality(Float_t maxCentrality=100.)
Definition: AliRDHFCuts.h:52
AliRDHFCutsD0toKpi & operator=(const AliRDHFCutsD0toKpi &source)
Int_t IsSelectedSpecialCuts(AliAODRecoDecayHF *d) const
Int_t CombineSelectionLevels(Int_t selectionvalTrack, Int_t selectionvalCand, Int_t selectionvalPID) const
void SetVarsForOpt(Int_t nVars, Bool_t *forOpt)
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:384
void EnableCentralTrigger()
Definition: AliRDHFCuts.h:81
Bool_t IsPionRaw(AliAODTrack *track, TString detector) const
void SetUsePID(Bool_t flag=kTRUE)
Definition: AliRDHFCuts.h:204
Double_t Normalizedd0Prong(Int_t ip) const
virtual void PrintAll() const
virtual Int_t IsSelectedPID(AliAODRecoDecayHF *rd)
void CleanOwnPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod, AliAODVertex *origownvtx) const
void SetRemoveDaughtersFromPrim(Bool_t removeDaughtersPrim)
Definition: AliRDHFCuts.h:214
Bool_t GetCompat() const
Definition: AliAODPidHF.h:147
void SetUseSpecialCuts(Bool_t useSpecialCuts)
void SetOldPid(Bool_t oldPid)
Definition: AliAODPidHF.h:108
void SetPtBins(Int_t nPtBinLimits, Float_t *ptBinLimits)
virtual void SetStandardCutsPP2010vsMult()
void SetLowPt(Bool_t lowpt, Double_t ptlow=2.)
Double_t * fWeightsPositive
number of species (used only for array declaration)
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Int_t fBayesianCondition
Switch for which Bayesian PID strategy to use.
void SetMatch(Int_t match)
Definition: AliAODPidHF.h:98
void AddTrackCuts(const AliESDtrackCuts *cuts)
Definition: AliRDHFCuts.h:202
void SetPLimit(Double_t *plim, Int_t npLim)
void SetTPC(Bool_t tpc)
Definition: AliAODPidHF.h:94
void SetMinPtCandidate(Double_t ptCand=-1.)
Definition: AliRDHFCuts.h:215
Double_t CosPointingAngle() const
void SetCompat(Bool_t comp)
Definition: AliAODPidHF.h:100
void SetTriggerClass(TString trclass0, TString trclass1="")
Definition: AliRDHFCuts.h:192
AliAODPidHF * fPidHF
enable AOD049 centrality cleanup
Definition: AliRDHFCuts.h:390
Bool_t RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod) const
virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d, Float_t *vars, Int_t nvars, Int_t *pdgdaughters)
Int_t PtBin(Double_t pt) const
void SetOptPileup(Int_t opt=0)
Definition: AliRDHFCuts.h:218
virtual void SetStandardCutsPbPb2010()
Int_t GetGlobalIndex(Int_t iVar, Int_t iPtBin) const
Double_t fPtLowPID
flag to switch on/off D0 selection via KF
Int_t nptbins
Double_t fMinPtCand
outcome of PID selection
Definition: AliRDHFCuts.h:409
Int_t IsSelectedPIDdefault(AliAODRecoDecayHF *rd)
virtual void SetStandardCutsPbPb2011()