]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliRDHFCutsD0toKpi.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
possibility to cut on the pt of candidate (Rossella)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliRDHFCutsD0toKpi.cxx
CommitLineData
650b3ced 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/////////////////////////////////////////////////////////////
17//
18// Class for cuts on AOD reconstructed D0->Kpi
19//
20// Author: A.Dainese, andrea.dainese@pd.infn.it
21/////////////////////////////////////////////////////////////
22
23#include <TDatabasePDG.h>
24#include <Riostream.h>
25
26#include "AliRDHFCutsD0toKpi.h"
27#include "AliAODRecoDecayHF2Prong.h"
28#include "AliAODTrack.h"
29#include "AliESDtrack.h"
937e1290 30#include "AliAODPid.h"
31#include "AliTPCPIDResponse.h"
32#include "AliAODVertex.h"
5f25117b 33#include "AliKFVertex.h"
34#include "AliKFParticle.h"
650b3ced 35
36ClassImp(AliRDHFCutsD0toKpi)
37
38//--------------------------------------------------------------------------
a9b75906 39AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const char* name) :
937e1290 40AliRDHFCuts(name),
08ea5c32 41fUseSpecialCuts(kFALSE),
42fLowPt(kTRUE),
5f25117b 43fDefaultPID(kTRUE),
44fUseKF(kFALSE)
650b3ced 45{
46 //
47 // Default Constructor
48 //
49 Int_t nvars=9;
50 SetNVars(nvars);
51 TString varNames[9]={"inv. mass [GeV]",
52 "dca [cm]",
53 "cosThetaStar",
54 "pTK [GeV/c]",
55 "pTPi [GeV/c]",
56 "d0K [cm]",
57 "d0Pi [cm]",
58 "d0d0 [cm^2]",
59 "cosThetaPoint"};
60 Bool_t isUpperCut[9]={kTRUE,
61 kTRUE,
62 kTRUE,
63 kFALSE,
64 kFALSE,
65 kTRUE,
66 kTRUE,
67 kTRUE,
68 kFALSE};
69 SetVarNames(nvars,varNames,isUpperCut);
70 Bool_t forOpt[9]={kFALSE,
71 kTRUE,
72 kTRUE,
73 kFALSE,
74 kFALSE,
75 kFALSE,
76 kFALSE,
77 kTRUE,
78 kTRUE};
79 SetVarsForOpt(4,forOpt);
80 Float_t limits[2]={0,999999999.};
81 SetPtBins(2,limits);
937e1290 82
650b3ced 83}
84//--------------------------------------------------------------------------
85AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const AliRDHFCutsD0toKpi &source) :
08ea5c32 86 AliRDHFCuts(source),
87 fUseSpecialCuts(source.fUseSpecialCuts),
88 fLowPt(source.fLowPt),
5f25117b 89 fDefaultPID(source.fDefaultPID),
90 fUseKF(source.fUseKF)
650b3ced 91{
92 //
93 // Copy constructor
94 //
95
96}
97//--------------------------------------------------------------------------
98AliRDHFCutsD0toKpi &AliRDHFCutsD0toKpi::operator=(const AliRDHFCutsD0toKpi &source)
99{
100 //
101 // assignment operator
102 //
103 if(&source == this) return *this;
104
08ea5c32 105 AliRDHFCuts::operator=(source);
937e1290 106 fUseSpecialCuts=source.fUseSpecialCuts;
08ea5c32 107 fLowPt=source.fLowPt;
108 fDefaultPID=source.fDefaultPID;
650b3ced 109
110 return *this;
111}
112
113
114//---------------------------------------------------------------------------
115void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
116 //
117 // Fills in vars the values of the variables
118 //
119
650b3ced 120 if(nvars!=fnVarsForOpt) {
121 printf("AliRDHFCutsD0toKpi::GetCutsVarsForOpt: wrong number of variables\n");
122 return;
123 }
124
125 AliAODRecoDecayHF2Prong *dd = (AliAODRecoDecayHF2Prong*)d;
650b3ced 126
650b3ced 127 Int_t iter=-1;
128 if(fVarsForOpt[0]){
129 iter++;
130 if(TMath::Abs(pdgdaughters[0])==211) {
131 vars[iter]=dd->InvMassD0();
132 } else {
133 vars[iter]=dd->InvMassD0bar();
134 }
135 }
136 if(fVarsForOpt[1]){
137 iter++;
138 vars[iter]=dd->GetDCA();
139 }
140 if(fVarsForOpt[2]){
141 iter++;
142 if(TMath::Abs(pdgdaughters[0])==211) {
143 vars[iter] = dd->CosThetaStarD0();
144 } else {
145 vars[iter] = dd->CosThetaStarD0bar();
146 }
147 }
148 if(fVarsForOpt[3]){
149 iter++;
150 if(TMath::Abs(pdgdaughters[0])==321) {
151 vars[iter]=dd->PtProng(0);
152 }
153 else{
154 vars[iter]=dd->PtProng(1);
155 }
156 }
157 if(fVarsForOpt[4]){
158 iter++;
159 if(TMath::Abs(pdgdaughters[0])==211) {
160 vars[iter]=dd->PtProng(0);
161 }
162 else{
163 vars[iter]=dd->PtProng(1);
164 }
165 }
166 if(fVarsForOpt[5]){
167 iter++;
168 if(TMath::Abs(pdgdaughters[0])==321) {
169 vars[iter]=dd->Getd0Prong(0);
170 }
171 else{
172 vars[iter]=dd->Getd0Prong(1);
173 }
174 }
175 if(fVarsForOpt[6]){
176 iter++;
177 if(TMath::Abs(pdgdaughters[0])==211) {
178 vars[iter]=dd->Getd0Prong(0);
179 }
180 else{
181 vars[iter]=dd->Getd0Prong(1);
182 }
183 }
184 if(fVarsForOpt[7]){
185 iter++;
186 vars[iter]= dd->Prodd0d0();
187 }
188 if(fVarsForOpt[8]){
189 iter++;
190 vars[iter]=dd->CosPointingAngle();
191 }
192
193 return;
194}
195//---------------------------------------------------------------------------
937e1290 196Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* aod) {
650b3ced 197 //
198 // Apply selection
199 //
47aa3d55 200
201
d7688946 202 fIsSelectedCuts=0;
203 fIsSelectedPID=0;
204
650b3ced 205 if(!fCutsRD){
206 cout<<"Cut matrice not inizialized. Exit..."<<endl;
3cbc09cd 207 return 0;
650b3ced 208 }
209 //PrintAll();
3cbc09cd 210 AliAODRecoDecayHF2Prong* d=(AliAODRecoDecayHF2Prong*)obj;
650b3ced 211
650b3ced 212 if(!d){
213 cout<<"AliAODRecoDecayHF2Prong null"<<endl;
214 return 0;
215 }
216
47aa3d55 217 Double_t ptD=d->Pt();
218 if(ptD<fMinPtCand) return 0;
219 if(ptD>fMaxPtCand) return 0;
937e1290 220
c96634a2 221 // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
222 Int_t returnvaluePID=3;
223 Int_t returnvalueCuts=3;
c96634a2 224
650b3ced 225 // selection on candidate
226 if(selectionLevel==AliRDHFCuts::kAll ||
227 selectionLevel==AliRDHFCuts::kCandidate) {
937e1290 228
d7688946 229 if(!fUseKF) {
5f25117b 230
231 //recalculate vertex w/o daughters
232 AliAODVertex *origownvtx=0x0;
233 AliAODVertex *recvtx=0x0;
5f25117b 234 if(fRemoveDaughtersFromPrimary) {
d7688946 235 if(!RecalcOwnPrimaryVtx(d,aod,origownvtx,recvtx)) return 0;
5f25117b 236 }
237
238
239 Double_t pt=d->Pt();
240
241 Int_t okD0=0,okD0bar=0;
242
243 Int_t ptbin=PtBin(pt);
244 if (ptbin==-1) {
d7688946 245 CleanOwnPrimaryVtx(d,origownvtx);
937e1290 246 return 0;
247 }
d7688946 248
5f25117b 249 Double_t mD0,mD0bar,ctsD0,ctsD0bar;
250 okD0=1; okD0bar=1;
251
252 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
d7688946 253
254 d->InvMassD0(mD0,mD0bar);
255 if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
256 if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;
257 if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
258
259 if(d->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)]) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
260
5f25117b 261
be9cad0c 262 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;
263 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;
d7688946 264 if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
5f25117b 265
266
267 if(TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
268 TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
269 if(TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||
270 TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;
d7688946 271 if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
5f25117b 272
d7688946 273 if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
5f25117b 274
650b3ced 275
5f25117b 276 d->CosThetaStarD0(ctsD0,ctsD0bar);
277 if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
278 if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;
d7688946 279 if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
5f25117b 280
be9cad0c 281 if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)]) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
5f25117b 282
283 if (returnvalueCuts!=0) {
284 if (okD0) returnvalueCuts=1; //cuts passed as D0
285 if (okD0bar) returnvalueCuts=2; //cuts passed as D0bar
286 if (okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar
287 }
288
289 // call special cuts
290 Int_t special=1;
291 if(fUseSpecialCuts) special=IsSelectedSpecialCuts(d);
d7688946 292 if(!special) {CleanOwnPrimaryVtx(d,origownvtx); return 0;}
5f25117b 293
294 // unset recalculated primary vertex when not needed any more
d7688946 295 CleanOwnPrimaryVtx(d,origownvtx);
650b3ced 296
d7688946 297 } else {
298 // go to selection with Kalman vertexing, if requested
299 returnvalueCuts = IsSelectedKF(d,aod);
300 }
dcc2ade0 301
d7688946 302 fIsSelectedCuts=returnvalueCuts;
dcc2ade0 303 if(!returnvalueCuts) return 0;
5f25117b 304 }
305
650b3ced 306
dcc2ade0 307 // selection on PID
308 if(selectionLevel==AliRDHFCuts::kAll ||
309 selectionLevel==AliRDHFCuts::kCandidate ||
310 selectionLevel==AliRDHFCuts::kPID) {
311 returnvaluePID = IsSelectedPID(d);
d7688946 312 fIsSelectedPID=returnvaluePID;
dcc2ade0 313 if(!returnvaluePID) return 0;
314 }
315
316 Int_t returnvalueComb=CombineSelectionLevels(3,returnvalueCuts,returnvaluePID);
317
318 if(!returnvalueComb) return 0;
319
320 // selection on daughter tracks
321 if(selectionLevel==AliRDHFCuts::kAll ||
322 selectionLevel==AliRDHFCuts::kTracks) {
323 if(!AreDaughtersSelected(d)) return 0;
324 }
325
5f25117b 326 // cout<<"Pid = "<<returnvaluePID<<endl;
e2b51a19 327 return returnvalueComb;
5f25117b 328}
329
330//------------------------------------------------------------------------------------------
331Int_t AliRDHFCutsD0toKpi::IsSelectedKF(AliAODRecoDecayHF2Prong *d,
332 AliAODEvent* aod) const {
333 //
334 // Apply selection using KF-vertexing
335 //
336
337 AliAODTrack *track0 = (AliAODTrack*)d->GetDaughter(0);
338 AliAODTrack *track1 = (AliAODTrack*)d->GetDaughter(1);
339
340 if(!track0 || !track1) {
341 cout<<"one or two D0 daughters missing!"<<endl;
342 return 0;
343 }
344
345 // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
346 Int_t returnvalueCuts=3;
347
348 // candidate selection with AliKF
349 AliKFParticle::SetField(aod->GetMagneticField()); // set the magnetic field
937e1290 350
5f25117b 351 Int_t okD0=0,okD0bar=0;
352 okD0=1; okD0bar=1;
353
354 // convert tracks into AliKFParticles
355
356 AliKFParticle negPiKF(*track1,-211); // neg pion kandidate
357 AliKFParticle negKKF(*track1,-321); // neg kaon kandidate
358 AliKFParticle posPiKF(*track0,211); // pos pion kandidate
359 AliKFParticle posKKF(*track0,321); // pos kaon kandidate
360
361 // build D0 candidates
362
363 AliKFParticle d0c(negKKF,posPiKF); // D0 candidate
364 AliKFParticle ad0c(posKKF,negPiKF); // D0bar candidate
365
366 // create kf primary vertices
367
368 AliAODVertex *vtx1 = aod->GetPrimaryVertex();
369 AliKFVertex primVtx(*vtx1);
370 AliKFVertex aprimVtx(*vtx1);
371
372 if(primVtx.GetNContributors()<=0) okD0 = 0;
373 if(aprimVtx.GetNContributors()<=0) okD0bar = 0;
374 if(!okD0 && !okD0bar) returnvalueCuts=0;
375
376 // calculate mass
377
378 Double_t d0mass = d0c.GetMass();
379 Double_t ad0mass = ad0c.GetMass();
380
381 // calculate P of D0 and D0bar
382 Double_t d0P = d0c.GetP();
383 Double_t d0Px = d0c.GetPx();
384 Double_t d0Py = d0c.GetPy();
385 Double_t d0Pz = d0c.GetPz();
386 Double_t ad0P = ad0c.GetP();
387 Double_t ad0Px = ad0c.GetPx();
388 Double_t ad0Py = ad0c.GetPy();
389 Double_t ad0Pz = ad0c.GetPz();
390
391 //calculate Pt of D0 and D0bar
392
393 Double_t pt=d0c.GetPt();
394 Double_t apt=ad0c.GetPt();
395
396 // remove D0 daughters from primary vertices (if used in vertex fit) and add D0-candidates
397
398 if(track0->GetUsedForPrimVtxFit()) {
399 primVtx -= posPiKF;
400 aprimVtx -= posKKF;
401 }
402
403 if(track1->GetUsedForPrimVtxFit()) {
404 primVtx -= negKKF;
405 aprimVtx -= negPiKF;
406 }
407
408 primVtx += d0c;
409 aprimVtx += ad0c;
410
411 if(primVtx.GetNContributors()<=0) okD0 = 0;
412 if(aprimVtx.GetNContributors()<=0) okD0bar = 0;
413 if(!okD0 && !okD0bar) returnvalueCuts=0;
414
415 //calculate cut variables
416
417 // calculate impact params of daughters w.r.t recalculated vertices
418
419 Double_t impactPi = posPiKF.GetDistanceFromVertexXY(primVtx);
420 Double_t aimpactPi = negPiKF.GetDistanceFromVertexXY(aprimVtx);
421 Double_t impactKa = negKKF.GetDistanceFromVertexXY(primVtx);
422 Double_t aimpactKa = posKKF.GetDistanceFromVertexXY(aprimVtx);
423
424 // calculate Product of Impact Params
425
426 Double_t prodParam = impactPi*impactKa;
427 Double_t aprodParam = aimpactPi*aimpactKa;
428
429 // calculate cosine of pointing angles
430
431 TVector3 mom(d0c.GetPx(),d0c.GetPy(),d0c.GetPz());
432 TVector3 fline(d0c.GetX()-primVtx.GetX(),
433 d0c.GetY()-primVtx.GetY(),
434 d0c.GetZ()-primVtx.GetZ());
435 Double_t pta = mom.Angle(fline);
436 Double_t cosP = TMath::Cos(pta); // cosine of pta for D0 candidate
437
438 TVector3 amom(ad0c.GetPx(),ad0c.GetPy(),ad0c.GetPz());
439 TVector3 afline(ad0c.GetX()-aprimVtx.GetX(),
440 ad0c.GetY()-aprimVtx.GetY(),
441 ad0c.GetZ()-aprimVtx.GetZ());
442 Double_t apta = amom.Angle(afline);
443 Double_t acosP = TMath::Cos(apta); // cosine of pta for D0bar candidate
444
445 // calculate P of Pions at Decay Position of D0 and D0bar candidates
446 negKKF.TransportToParticle(d0c);
447 posPiKF.TransportToParticle(d0c);
448 posKKF.TransportToParticle(ad0c);
449 negPiKF.TransportToParticle(ad0c);
450
451 Double_t pxPi = posPiKF.GetPx();
452 Double_t pyPi = posPiKF.GetPy();
453 Double_t pzPi = posPiKF.GetPz();
454 Double_t ptPi = posPiKF.GetPt();
455
456 Double_t apxPi = negPiKF.GetPx();
457 Double_t apyPi = negPiKF.GetPy();
458 Double_t apzPi = negPiKF.GetPz();
459 Double_t aptPi = negPiKF.GetPt();
460
461 // calculate Pt of Kaons at Decay Position of D0 and D0bar candidates
462
463 Double_t ptK = negKKF.GetPt();
464 Double_t aptK = posKKF.GetPt();
465
466 //calculate cos(thetastar)
467 Double_t massvtx = TDatabasePDG::Instance()->GetParticle(421)->Mass();
468 Double_t massp[2];
469 massp[0] = TDatabasePDG::Instance()->GetParticle(321)->Mass();
470 massp[1] = TDatabasePDG::Instance()->GetParticle(211)->Mass();
471 Double_t pStar = TMath::Sqrt(TMath::Power(massvtx*massvtx-massp[0]*massp[0]-massp[1]*massp[1],2.)
472 -4.*massp[0]*massp[0]*massp[1]*massp[1])/(2.*massvtx);
473
474 // cos(thetastar) for D0 and Pion
475
476 Double_t d0E = TMath::Sqrt(massvtx*massvtx + d0P*d0P);
477 Double_t beta = d0P/d0E;
478 Double_t gamma = d0E/massvtx;
479 TVector3 momPi(pxPi,pyPi,pzPi);
480 TVector3 momTot(d0Px,d0Py,d0Pz);
481 Double_t q1 = momPi.Dot(momTot)/momTot.Mag();
482 Double_t cts = (q1/gamma-beta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar;
483
484 // cos(thetastar) for D0bar and Pion
485
486 Double_t ad0E = TMath::Sqrt(massvtx*massvtx + ad0P*ad0P);
487 Double_t abeta = ad0P/ad0E;
488 Double_t agamma = ad0E/massvtx;
489 TVector3 amomPi(apxPi,apyPi,apzPi);
490 TVector3 amomTot(ad0Px,ad0Py,ad0Pz);
491 Double_t aq1 = amomPi.Dot(amomTot)/amomTot.Mag();
492 Double_t acts = (aq1/agamma-abeta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar;
493
494 // calculate reduced Chi2 for the full D0 fit
495 d0c.SetProductionVertex(primVtx);
496 ad0c.SetProductionVertex(aprimVtx);
497 negKKF.SetProductionVertex(d0c);
498 posPiKF.SetProductionVertex(d0c);
499 posKKF.SetProductionVertex(ad0c);
500 negPiKF.SetProductionVertex(ad0c);
501 d0c.TransportToProductionVertex();
502 ad0c.TransportToProductionVertex();
503
504 // calculate the decay length
505 Double_t decayLengthD0 = d0c.GetDecayLength();
506 Double_t adecayLengthD0 = ad0c.GetDecayLength();
507
508 Double_t chi2D0 = 50.;
509 if(d0c.GetNDF() > 0 && d0c.GetChi2() >= 0) {
510 chi2D0 = d0c.GetChi2()/d0c.GetNDF();
511 }
512
513 Double_t achi2D0 = 50.;
514 if(ad0c.GetNDF() > 0 && ad0c.GetChi2() >= 0) {
515 achi2D0 = ad0c.GetChi2()/ad0c.GetNDF();
516 }
517
518 // Get the Pt-bins
519 Int_t ptbin=PtBin(pt);
520 Int_t aptbin=PtBin(apt);
521
522 if(ptbin < 0) okD0 = 0;
523 if(aptbin < 0) okD0bar = 0;
524 if(!okD0 && !okD0bar) returnvalueCuts=0;
525
526 if(ptK < fCutsRD[GetGlobalIndex(3,ptbin)] || ptPi < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
527 if(aptK < fCutsRD[GetGlobalIndex(3,aptbin)] || aptPi < fCutsRD[GetGlobalIndex(4,aptbin)]) okD0bar = 0;
528 if(!okD0 && !okD0bar) returnvalueCuts=0;
529
530
531 if(TMath::Abs(impactKa) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
532 TMath::Abs(impactPi) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
533
534 if(TMath::Abs(aimpactKa) > fCutsRD[GetGlobalIndex(5,aptbin)] ||
535 TMath::Abs(aimpactPi) > fCutsRD[GetGlobalIndex(6,aptbin)]) okD0bar = 0;
536
537 if(!okD0 && !okD0bar) returnvalueCuts=0;
538
539 // for the moment via the standard method due to bug in AliKF
540 if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) okD0 = 0;
541 if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,aptbin)]) okD0bar = 0;
542 if(!okD0 && !okD0bar) returnvalueCuts=0;
650b3ced 543
650b3ced 544
5f25117b 545 if(TMath::Abs(d0mass-massvtx) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
546 if(TMath::Abs(ad0mass-massvtx) > fCutsRD[GetGlobalIndex(0,aptbin)]) okD0bar = 0;
547 if(!okD0 && !okD0bar) returnvalueCuts=0;
548
549
550 if(TMath::Abs(cts) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
551 if(TMath::Abs(acts) > fCutsRD[GetGlobalIndex(2,aptbin)]) okD0bar = 0;
552 if(!okD0 && !okD0bar) returnvalueCuts=0;
553
554 if(prodParam > fCutsRD[GetGlobalIndex(7,ptbin)]) okD0 = 0;
555 if(aprodParam > fCutsRD[GetGlobalIndex(7,aptbin)]) okD0bar = 0;
556 if(!okD0 && !okD0bar) returnvalueCuts=0;
650b3ced 557
5f25117b 558 if(cosP < fCutsRD[GetGlobalIndex(8,ptbin)]) okD0 = 0;
559 if(acosP < fCutsRD[GetGlobalIndex(8,aptbin)]) okD0bar = 0;
560 if(!okD0 && !okD0bar) returnvalueCuts=0;
561
562 if(chi2D0 > fCutsRD[GetGlobalIndex(10,ptbin)]) okD0 = 0;
563 if(achi2D0 > fCutsRD[GetGlobalIndex(10,aptbin)]) okD0bar = 0;
564 if(!okD0 && !okD0bar) returnvalueCuts=0;
565
566 if(decayLengthD0 < fCutsRD[GetGlobalIndex(9,ptbin)]) okD0 = 0;
567 if(adecayLengthD0 < fCutsRD[GetGlobalIndex(9,aptbin)]) okD0bar = 0;
568 if(!okD0 && !okD0bar) returnvalueCuts=0;
650b3ced 569
5f25117b 570 if(returnvalueCuts!=0) {
571 if(okD0) returnvalueCuts=1; //cuts passed as D0
572 if(okD0bar) returnvalueCuts=2; //cuts passed as D0bar
573 if(okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar
650b3ced 574 }
575
5f25117b 576 return returnvalueCuts;
650b3ced 577}
5f25117b 578
650b3ced 579//---------------------------------------------------------------------------
5f25117b 580
aa8d25fe 581Bool_t AliRDHFCutsD0toKpi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
582{
583 //
584 // Checking if D0 is in fiducial acceptance region
585 //
586
c96634a2 587 if(pt > 5.) {
aa8d25fe 588 // applying cut for pt > 5 GeV
937e1290 589 AliDebug(2,Form("pt of D0 = %f (> 5), cutting at |y| < 0.8\n",pt));
aa8d25fe 590 if (TMath::Abs(y) > 0.8){
591 return kFALSE;
592 }
593 } else {
594 // appliying smooth cut for pt < 5 GeV
595 Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
596 Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
937e1290 597 AliDebug(2,Form("pt of D0 = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
aa8d25fe 598 if (y < minFiducialY || y > maxFiducialY){
599 return kFALSE;
600 }
601 }
c96634a2 602
aa8d25fe 603 return kTRUE;
aa8d25fe 604}
c96634a2 605//---------------------------------------------------------------------------
937e1290 606Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d)
08ea5c32 607{
608 // ############################################################
609 //
610 // Apply PID selection
611 //
612 //
613 // ############################################################
614
615 if(!fUsePID) return 3;
616 if(fDefaultPID) return IsSelectedPIDdefault(d);
617 fWhyRejection=0;
618 Int_t isD0D0barPID[2]={1,2};
619 Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value;
620 // same for prong 2
621 // values convention -1 = discarded
622 // 0 = not identified (but compatible) || No PID (->hasPID flag)
623 // 1 = identified
624 // PID search: pion (TPC) or not K (TOF), Kaon hypothesis for both
625 // Initial hypothesis: unknwon (but compatible)
626 combinedPID[0][0]=0; // prima figlia, Kaon
627 combinedPID[0][1]=0; // prima figlia, pione
628 combinedPID[1][0]=0; // seconda figlia, Kaon
629 combinedPID[1][1]=0; // seconda figlia, pion
630
631 Bool_t checkPIDInfo[2]={kTRUE,kTRUE};
632 Double_t sigma_tmp[3]={fPidHF->GetSigma(0),fPidHF->GetSigma(1),fPidHF->GetSigma(2)};
633 for(Int_t daught=0;daught<2;daught++){
634 //Loop con prongs
635 AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);
636
637 if(!(fPidHF->CheckStatus(aodtrack,"TPC")) && !(fPidHF->CheckStatus(aodtrack,"TOF"))) {
638 checkPIDInfo[daught]=kFALSE;
639 continue;
640 }
641
642 // identify kaon
643 combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
644
645 // identify pion
646
647 if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {
648 combinedPID[daught][1]=0;
649 }else{
650 fPidHF->SetTOF(kFALSE);
651 combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
652 fPidHF->SetTOF(kTRUE);
653 fPidHF->SetCompat(kTRUE);
654 }
655
656
657 if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){ // if not a K- and not a pi- both D0 and D0bar excluded
658 isD0D0barPID[0]=0;
659 isD0D0barPID[1]=0;
660 }
661 else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){
662 if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded
663 else isD0D0barPID[0]=0;// if K+ D0 excluded
664 }
665 else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){
666 isD0D0barPID[0]=0;
667 isD0D0barPID[1]=0;
668 }
669 else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){
670 if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded
671 else isD0D0barPID[0]=0;// not a D0 if K+ or if pi+ excluded
672 }
673 else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){
674 if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded
675 else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded
676 }
677
678 if(fLowPt && d->Pt()<2.){
679 Double_t sigmaTPC[3]={3.,2.,0.};
680 fPidHF->SetSigmaForTPC(sigmaTPC);
681 // identify kaon
682 combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
683
684 Double_t ptProng=aodtrack->P();
685
686 if(ptProng<0.6){
687 fPidHF->SetCompat(kFALSE);
688 combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
689 fPidHF->SetCompat(kTRUE);
690 }
691
692 if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {
693 combinedPID[daught][1]=0;
694 }else{
695 fPidHF->SetTOF(kFALSE);
696 Double_t sigmaTPCpi[3]={3.,3.,0.};
697 fPidHF->SetSigmaForTPC(sigmaTPCpi);
698 combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
699 fPidHF->SetTOF(kTRUE);
700 if(ptProng<0.8){
701 Bool_t isTPCpion=fPidHF->IsPionRaw(aodtrack,"TPC");
702 if(isTPCpion){
703 combinedPID[daught][1]=1;
704 }else{
705 combinedPID[daught][1]=-1;
706 }
707 }
708 }
709
710 }
711 fPidHF->SetSigmaForTPC(sigma_tmp);
712 }// END OF LOOP ON DAUGHTERS
713
714 if(!checkPIDInfo[0] && !checkPIDInfo[1]) {
715 if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
716 return 0;
717 }
718
719
720 // FURTHER PID REQUEST (both daughter info is needed)
721 if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){
722 fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found
723 if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
724 return 0;
725 }
726
727 if(d->Pt()<2.){
728 if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
729 if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
730 fWhyRejection=32;// reject cases where the Kaon is not identified
731 return 0;
732 }
733 }
734 if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
735
736 // cout<<"Why? "<<fWhyRejection<<endl;
737 return isD0D0barPID[0]+isD0D0barPID[1];
738}
739//---------------------------------------------------------------------------
740Int_t AliRDHFCutsD0toKpi::IsSelectedPIDdefault(AliAODRecoDecayHF* d)
937e1290 741{
742 // ############################################################
c96634a2 743 //
744 // Apply PID selection
c96634a2 745 //
937e1290 746 //
747 // temporary selection: PID AS USED FOR D0 by Andrea Rossi (up to 28/06/2010)
748 //
749 // d must be a AliAODRecoDecayHF2Prong object
750 // returns 0 if both D0 and D0bar are rejectecd
751 // 1 if D0 is accepted while D0bar is rejected
752 // 2 if D0bar is accepted while D0 is rejected
753 // 3 if both are accepted
754 // fWhyRejection variable (not returned for the moment, print it if needed)
755 // keeps some information on why a candidate has been
756 // rejected according to the following (unfriendly?) scheme
757 // if more rejection cases are considered interesting, just add numbers
758 //
759 // TO BE CONSIDERED WITH A GRAIN OF SALT (the order in which cut are applied is relevant)
760 // from 20 to 30: "detector" selection (PID acceptance)
761 // 26: TPC refit
762 // 27: ITS refit
763 // 28: no (TOF||TPC) pid information (no kTOFpid,kTOFout,kTIME,kTPCpid,...)
764 //
765 // from 30 to 40: PID selection
766 // 31: no Kaon compatible tracks found between daughters
767 // 32: no Kaon identified tracks found (strong sel. at low momenta)
768 // 33: both mass hypotheses are rejected
769 //
770 // ############################################################
771
772 if(!fUsePID) return 3;
773 fWhyRejection=0;
774 Int_t isD0D0barPID[2]={1,2};
775 Double_t nsigmaTPCpi=-1., nsigmaTPCK=-1.; //used for TPC pid
776 Double_t tofSig,times[5];// used fot TOF pid
777 Int_t hasPID[2]={2,2};// flag to count how many detectors give PID info for the daughters
778 Int_t isKaonPionTOF[2][2],isKaonPionTPC[2][2];
779 Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value;
780 // same for prong 2
781 // values convention -1 = discarded
782 // 0 = not identified (but compatible) || No PID (->hasPID flag)
783 // 1 = identified
784 // PID search: pion (TPC) or not K (TOF), Kaon hypothesis for both
785 // Initial hypothesis: unknwon (but compatible)
786 isKaonPionTOF[0][0]=0;
787 isKaonPionTOF[0][1]=0;
788 isKaonPionTOF[1][0]=0;
789 isKaonPionTOF[1][1]=0;
790
791 isKaonPionTPC[0][0]=0;
792 isKaonPionTPC[0][1]=0;
793 isKaonPionTPC[1][0]=0;
794 isKaonPionTPC[1][1]=0;
795
796 combinedPID[0][0]=0;
797 combinedPID[0][1]=0;
798 combinedPID[1][0]=0;
799 combinedPID[1][1]=0;
800
801
802
803 for(Int_t daught=0;daught<2;daught++){
804 //Loop con prongs
805
806 // ########### Step 0- CHECKING minimal PID "ACCEPTANCE" ####################
807
808 AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);
809
810 if(!(aodtrack->GetStatus()&AliESDtrack::kTPCrefit)){
811 fWhyRejection=26;
812 return 0;
813 }
814 if(!(aodtrack->GetStatus()&AliESDtrack::kITSrefit)){
815 fWhyRejection=27;
816 return 0;
817 }
818
819 AliAODPid *pid=aodtrack->GetDetPid();
820 if(!pid) {
821 //delete esdtrack;
822 hasPID[daught]--;
823 continue;
824 }
825
826 // ########### Step 1- Check of TPC and TOF response ####################
827
828 Double_t ptrack=aodtrack->P();
829 //#################### TPC PID #######################
830 if (!(aodtrack->GetStatus()&AliESDtrack::kTPCpid )){
831 // NO TPC PID INFO FOR THIS TRACK
832 hasPID[daught]--;
833 }
834 else {
835 static AliTPCPIDResponse theTPCpid;
08ea5c32 836 AliAODPid *pidObj = aodtrack->GetDetPid();
837 Double_t ptProng=pidObj->GetTPCmomentum();
838 nsigmaTPCpi = theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kPion);
839 nsigmaTPCK = theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kKaon);
840 //if(ptrack<0.6){
841 if(ptProng<0.6){
937e1290 842 if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;
843 else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
844 if(TMath::Abs(nsigmaTPCpi)<2.)isKaonPionTPC[daught][1]=1;
845 else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
846 }
08ea5c32 847 //else if(ptrack<.8){
848 else if(ptProng<.8){
937e1290 849 if(TMath::Abs(nsigmaTPCK)<1.)isKaonPionTPC[daught][0]=1;
850 else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
851 if(TMath::Abs(nsigmaTPCpi)<1.)isKaonPionTPC[daught][1]=1;
852 else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
853 }
854 else {
855 // if(nsigmaTPCK>-2.&&nsigmaTPCK<1.)isKaonPionTPC[daught][0]=1;
856 if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
857 //if(nsigmaTPCpi>-1.&&nsigmaTPCpi<2.)isKaonPionTPC[daught][1]=1;
858 if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
859 }
860 }
861
862
863 // ##### TOF PID: do not ask nothing for pion/protons ############
864 if(!((aodtrack->GetStatus()&AliESDtrack::kTOFpid)&&(aodtrack->GetStatus()&AliESDtrack::kTOFout)&&(aodtrack->GetStatus()&AliESDtrack::kTIME))){
865 // NO TOF PID INFO FOR THIS TRACK
866 hasPID[daught]--;
867 }
868 else{
869 tofSig=pid->GetTOFsignal();
870 pid->GetIntegratedTimes(times);
b14fe7f1 871 if((tofSig-times[3])>5.*160.)return 0;// PROTON REJECTION
937e1290 872 if(TMath::Abs(tofSig-times[3])>3.*160.){
873 isKaonPionTOF[daught][0]=-1;
874 }
875 else {
876 if(ptrack<1.5){
877 isKaonPionTOF[daught][0]=1;
878 }
879 }
880 }
881
882 //######### Step 2: COMBINE TOF and TPC PID ###############
883 // we apply the following convention: if TPC and TOF disagree (discarded Vs identified) -> unknown
884 combinedPID[daught][0]=isKaonPionTOF[daught][0]+isKaonPionTPC[daught][0];
885 combinedPID[daught][1]=isKaonPionTOF[daught][1]+isKaonPionTPC[daught][1];
886
887
888 //######### Step 3: USE PID INFO
889
890 if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){// if not a K- and not a pi- both D0 and D0bar excluded
891 isD0D0barPID[0]=0;
892 isD0D0barPID[1]=0;
893 }
894 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
895 if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded
896 else isD0D0barPID[0]=0;// if K+ D0 excluded
897 }
898 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
899 isD0D0barPID[0]=0;
900 isD0D0barPID[1]=0;
901 }
902 else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){
903 if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded
904 else isD0D0barPID[0]=0;// not a D0 if K+ or if pi+ excluded
905 }
906 else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){
907 if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded
908 else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded
909 }
910
911 // ########## ALSO DIFFERENT TPC PID REQUEST FOR LOW pt D0: request of K identification ###############################
912 // ########## more tolerant criteria for single particle ID-> more selective criteria for D0 ##############################
913 // ############### NOT OPTIMIZED YET ###################################
914 if(d->Pt()<2.){
915 isKaonPionTPC[daught][0]=0;
916 isKaonPionTPC[daught][1]=0;
08ea5c32 917 AliAODPid *pidObj = aodtrack->GetDetPid();
918 Double_t ptProng=pidObj->GetTPCmomentum();
919 //if(ptrack<0.6){
920 if(ptProng<0.6){
937e1290 921 if(TMath::Abs(nsigmaTPCK)<3.)isKaonPionTPC[daught][0]=1;
922 else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
923 if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;
924 else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
925 }
08ea5c32 926 //else if(ptrack<.8){
927 else if(ptProng<.8){
937e1290 928 if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;
929 else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
930 if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;
931 else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
932 }
933 else {
934 if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
935 if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
936 }
937 }
938
939 }// END OF LOOP ON DAUGHTERS
940
941 // FURTHER PID REQUEST (both daughter info is needed)
942 if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){
943 fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found
944 return 0;
945 }
946 else if(hasPID[0]==0&&hasPID[1]==0){
947 fWhyRejection=28;// reject cases in which no PID info is available
948 return 0;
949 }
950 if(d->Pt()<2.){
951 // request of K identification at low D0 pt
952 combinedPID[0][0]=0;
953 combinedPID[0][1]=0;
954 combinedPID[1][0]=0;
955 combinedPID[1][1]=0;
956
957 combinedPID[0][0]=isKaonPionTOF[0][0]+isKaonPionTPC[0][0];
958 combinedPID[0][1]=isKaonPionTOF[0][1]+isKaonPionTPC[0][1];
aaf00f06 959 combinedPID[1][0]=isKaonPionTOF[1][0]+isKaonPionTPC[1][0];
960 combinedPID[1][1]=isKaonPionTOF[1][1]+isKaonPionTPC[1][1];
937e1290 961
962 if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
963 fWhyRejection=32;// reject cases where the Kaon is not identified
964 return 0;
965 }
966 }
967
968 // cout<<"Why? "<<fWhyRejection<<endl;
969 return isD0D0barPID[0]+isD0D0barPID[1];
970}
971
08ea5c32 972
973
937e1290 974//---------------------------------------------------------------------------
975Int_t AliRDHFCutsD0toKpi::CombineSelectionLevels(Int_t selectionvalTrack,
976 Int_t selectionvalCand,
977 Int_t selectionvalPID) const
978{
979 //
980 // This method combines the tracks, PID and cuts selection results
981 //
982 if(selectionvalTrack==0) return 0;
983
984 Int_t returnvalue;
c96634a2 985
937e1290 986 switch(selectionvalPID) {
987 case 0:
988 returnvalue=0;
989 break;
990 case 1:
991 returnvalue=((selectionvalCand==1 || selectionvalCand==3) ? 1 : 0);
992 break;
993 case 2:
994 returnvalue=((selectionvalCand==2 || selectionvalCand==3) ? 2 : 0);
995 break;
996 case 3:
997 returnvalue=selectionvalCand;
998 break;
999 default:
1000 returnvalue=0;
1001 break;
1002 }
1003
1004 return returnvalue;
1005}
1006//----------------------------------------------------------------------------
1007Int_t AliRDHFCutsD0toKpi::IsSelectedSpecialCuts(AliAODRecoDecayHF *d) const
1008{
1009 //
1010 // Note: this method is temporary
1011 // Additional cuts on decay lenght and lower cut for d0 norm are applied using vertex without candidate's daughters
1012 //
c96634a2 1013
937e1290 1014 //apply cuts
aa8d25fe 1015
937e1290 1016 Float_t normDecLengthCut=1.,decLengthCut=TMath::Min(d->P()*0.0066+0.01,0.06/*cm*/), normd0Cut=0.5;
1017 // "decay length" expo law with tau' = beta*gamma*ctau= p/m*ctau =p*0.0123/1.864~p*0.0066
1018 // decay lenght > ctau' implies to retain (1-1/e) (for signal without considering detector resolution),
c96634a2 1019
937e1290 1020 Int_t returnvalue=3; //cut passed
1021 for(Int_t i=0;i<2/*prongs*/;i++){
1022 if(TMath::Abs(d->Normalizedd0Prong(i))<normd0Cut) return 0; //normd0Cut not passed
1023 }
be9cad0c 1024 if(d->DecayLength2()<decLengthCut*decLengthCut) return 0; //decLengthCut not passed
1025 if(d->NormalizedDecayLength2()<normDecLengthCut*normDecLengthCut) return 0; //decLengthCut not passed
937e1290 1026
c96634a2 1027
1028 return returnvalue;
1029}
08ea5c32 1030
5f25117b 1031//----------------------------------------------
5f25117b 1032void AliRDHFCutsD0toKpi::SetUseKF(Bool_t useKF)
1033{
1034 //switch on candidate selection via AliKFparticle
1035 if(!useKF) return;
1036 if(useKF){
1037 fUseKF=useKF;
1038 Int_t nvarsKF=11;
1039 SetNVars(nvarsKF);
1040 TString varNamesKF[11]={"inv. mass [GeV]",
1041 "dca [cm]",
1042 "cosThetaStar",
1043 "pTK [GeV/c]",
1044 "pTPi [GeV/c]",
1045 "d0K [cm]",
1046 "d0Pi [cm]",
1047 "d0d0 [cm^2]",
1048 "cosThetaPoint"
1049 "DecayLength[cm]",
1050 "RedChi2"};
1051 Bool_t isUpperCutKF[11]={kTRUE,
1052 kTRUE,
1053 kTRUE,
1054 kFALSE,
1055 kFALSE,
1056 kTRUE,
1057 kTRUE,
1058 kTRUE,
1059 kFALSE,
1060 kFALSE,
1061 kTRUE};
1062 SetVarNames(nvarsKF,varNamesKF,isUpperCutKF);
1063 SetGlobalIndex();
1064 Bool_t forOpt[11]={kFALSE,
1065 kTRUE,
1066 kTRUE,
1067 kFALSE,
1068 kFALSE,
1069 kFALSE,
1070 kFALSE,
1071 kTRUE,
1072 kTRUE,
1073 kFALSE,
1074 kFALSE};
1075 SetVarsForOpt(4,forOpt);
1076 }
1077 return;
1078}
b14fe7f1 1079
1080
73f4ccee 1081void AliRDHFCutsD0toKpi::SetStandardCutsPP2010() {
1082 //
1083 //STANDARD CUTS USED FOR 2010 pp analysis
1084 //dca cut will be enlarged soon to 400 micron
1085 //
b14fe7f1 1086
1087 SetName("D0toKpiCutsStandard");
1088 SetTitle("Standard Cuts for D0 analysis");
1089
1090 // PILE UP REJECTION
1091 SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
1092
1093 // EVENT CUTS
1094 SetMinVtxContr(1);
1095
1096
1097 // TRACKS ON SINGLE TRACKS
1098 AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
1099 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1100 esdTrackCuts->SetRequireTPCRefit(kTRUE);
1101 esdTrackCuts->SetRequireITSRefit(kTRUE);
1102 // esdTrackCuts->SetMinNClustersITS(4);
1103 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
1104 esdTrackCuts->SetMinDCAToVertexXY(0.);
1105 esdTrackCuts->SetEtaRange(-0.8,0.8);
1106 esdTrackCuts->SetPtRange(0.3,1.e10);
1107
1108 AddTrackCuts(esdTrackCuts);
1109
1110 const Int_t nptbins =13;
1111 const Double_t ptmax = 9999.;
1112 const Int_t nvars=9;
1113 Float_t ptbins[nptbins+1];
1114 ptbins[0]=0.;
1115 ptbins[1]=0.5;
1116 ptbins[2]=1.;
1117 ptbins[3]=2.;
1118 ptbins[4]=3.;
1119 ptbins[5]=4.;
1120 ptbins[6]=5.;
1121 ptbins[7]=6.;
1122 ptbins[8]=8.;
1123 ptbins[9]=12.;
1124 ptbins[10]=16.;
1125 ptbins[11]=20.;
1126 ptbins[12]=24.;
1127 ptbins[13]=ptmax;
1128
1129 SetGlobalIndex(nvars,nptbins);
1130 SetPtBins(nptbins+1,ptbins);
1131
1132 Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,300.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-35000.*1E-8,0.73},/* pt<0.5*/
1133 {0.400,300.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-35000.*1E-8,0.73},/* 0.5<pt<1*/
1134 {0.400,200.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.75},/* 1<pt<2 */
1135 {0.400,200.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.8},/* 2<pt<3 */
1136 {0.400,200.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85},/* 3<pt<4 */
1137 {0.400,200.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85},/* 4<pt<5 */
1138 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85},/* 5<pt<6 */
1139 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85},/* 6<pt<8 */
1140 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85},/* 8<pt<12 */
1141 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.*1E-8,0.85},/* 12<pt<16 */
1142 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.85},/* 16<pt<20 */
1143 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.85},/* 20<pt<24 */
1144 {0.400,150.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.85}};/* pt>24 */
1145
1146
1147 //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1148 Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1149 for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1150
1151 for (Int_t ibin=0;ibin<nptbins;ibin++){
1152 for (Int_t ivar = 0; ivar<nvars; ivar++){
1153 cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1154 }
1155 }
1156
1157 SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1158 SetUseSpecialCuts(kTRUE);
1159 SetRemoveDaughtersFromPrim(kTRUE);
1160
e11ae259 1161 for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1162 delete [] cutsMatrixTransposeStand;
1163 cutsMatrixTransposeStand=NULL;
1164
b14fe7f1 1165 // PID SETTINGS
1166 AliAODPidHF* pidObj=new AliAODPidHF();
1167 //pidObj->SetName("pid4D0");
1168 Int_t mode=1;
1169 const Int_t nlims=2;
1170 Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1171 Bool_t compat=kTRUE; //effective only for this mode
1172 Bool_t asym=kTRUE;
1173 Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1174 pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1175 pidObj->SetMatch(mode);
1176 pidObj->SetPLimit(plims,nlims);
1177 pidObj->SetSigma(sigmas);
1178 pidObj->SetCompat(compat);
1179 pidObj->SetTPC(kTRUE);
1180 pidObj->SetTOF(kTRUE);
1181
1182 SetPidHF(pidObj);
1183 SetUsePID(kTRUE);
1184 SetUseDefaultPID(kFALSE);
1185
1186
1187 PrintAll();
1188
e11ae259 1189 delete pidObj;
1190 pidObj=NULL;
1191
b14fe7f1 1192 return;
1193
1194}
1195
1196
73f4ccee 1197void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010() {
1198 //
1199 //PRELIMINARY CUTS USED FOR 2010 PbPb analysis
1200 //... EVOLVING SOON
1201 //
b14fe7f1 1202
1203 SetName("D0toKpiCutsStandard");
1204 SetTitle("Standard Cuts for D0 analysis in PbPb2010 run");
1205
1206 // PILE UP REJECTION
1207 //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
1208
1209 // EVENT CUTS
1210 SetMinVtxContr(1);
1211
1212
1213 // TRACKS ON SINGLE TRACKS
1214 AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
1215 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
1216 esdTrackCuts->SetRequireTPCRefit(kTRUE);
1217 esdTrackCuts->SetRequireITSRefit(kTRUE);
1218 // esdTrackCuts->SetMinNClustersITS(4);
1219 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
1220 esdTrackCuts->SetMinDCAToVertexXY(0.);
1221 esdTrackCuts->SetEtaRange(-0.8,0.8);
1222 esdTrackCuts->SetPtRange(0.3,1.e10);
1223
1224 AddTrackCuts(esdTrackCuts);
1225
1226 const Int_t nptbins =13;
1227 const Double_t ptmax = 9999.;
1228 const Int_t nvars=9;
1229 Float_t ptbins[nptbins+1];
1230 ptbins[0]=0.;
1231 ptbins[1]=0.5;
1232 ptbins[2]=1.;
1233 ptbins[3]=2.;
1234 ptbins[4]=3.;
1235 ptbins[5]=4.;
1236 ptbins[6]=5.;
1237 ptbins[7]=6.;
1238 ptbins[8]=8.;
1239 ptbins[9]=12.;
1240 ptbins[10]=16.;
1241 ptbins[11]=20.;
1242 ptbins[12]=24.;
1243 ptbins[13]=ptmax;
1244
1245 SetGlobalIndex(nvars,nptbins);
1246 SetPtBins(nptbins+1,ptbins);
1247
1248 Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,300.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.8},/* pt<0.5*/
1249 {0.400,300.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.8},/* 0.5<pt<1*/
1250 {0.400,250.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-32000.*1E-8,0.8},/* 1<pt<2 */
1251 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-26000.*1E-8,0.94},/* 2<pt<3 */
1252 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-1500.*1E-8,0.88},/* 3<pt<4 */
1253 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-1500.*1E-8,0.88},/* 4<pt<5 */
1254 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.90},/* 5<pt<6 */
1255 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.90},/* 6<pt<8 */
1256 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.90},/* 8<pt<12 */
1257 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.90},/* 12<pt<16 */
1258 {0.400,350.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85},/* 16<pt<20 */
1259 {0.400,350.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.*1E-8,0.85},/* 20<pt<24 */
1260 {0.400,350.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.*1E-8,0.82}};/* pt>24 */
1261
1262
1263 //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
1264 Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
1265 for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
1266
1267 for (Int_t ibin=0;ibin<nptbins;ibin++){
1268 for (Int_t ivar = 0; ivar<nvars; ivar++){
1269 cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
1270 }
1271 }
1272
1273 SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
1274 SetUseSpecialCuts(kTRUE);
1275 SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb
e11ae259 1276 for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
1277 delete [] cutsMatrixTransposeStand;
1278 cutsMatrixTransposeStand=NULL;
b14fe7f1 1279
1280 // PID SETTINGS
1281 AliAODPidHF* pidObj=new AliAODPidHF();
1282 //pidObj->SetName("pid4D0");
1283 Int_t mode=1;
1284 const Int_t nlims=2;
1285 Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
1286 Bool_t compat=kTRUE; //effective only for this mode
1287 Bool_t asym=kTRUE;
1288 Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
1289 pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
1290 pidObj->SetMatch(mode);
1291 pidObj->SetPLimit(plims,nlims);
1292 pidObj->SetSigma(sigmas);
1293 pidObj->SetCompat(compat);
1294 pidObj->SetTPC(kTRUE);
1295 pidObj->SetTOF(kTRUE);
1296
1297 SetPidHF(pidObj);
1298 SetUsePID(kTRUE);
1299 SetUseDefaultPID(kTRUE);// TEMPORARY: PROTON EXCLUSION SET ONLY IN DEFAULT PID
1300
1301
1302 PrintAll();
1303
e11ae259 1304
1305 delete pidObj;
1306 pidObj=NULL;
1307
b14fe7f1 1308 return;
1309
1310}