]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/AliRDHFCutsDstoKKpi.cxx
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliRDHFCutsDstoKKpi.cxx
CommitLineData
42828f21 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
18/////////////////////////////////////////////////////////////
19//
20// Class for cuts on AOD reconstructed Ds->KKpi
21//
22// Author: A.Dainese, andrea.dainese@pd.infn.it
23/////////////////////////////////////////////////////////////
24
25#include <TDatabasePDG.h>
26#include <Riostream.h>
27
28#include "AliRDHFCutsDstoKKpi.h"
29#include "AliAODRecoDecayHF3Prong.h"
30#include "AliAODTrack.h"
31#include "AliESDtrack.h"
32
c64cb1f6 33using std::cout;
34using std::endl;
35
42828f21 36ClassImp(AliRDHFCutsDstoKKpi)
37
38//--------------------------------------------------------------------------
39AliRDHFCutsDstoKKpi::AliRDHFCutsDstoKKpi(const char* name) :
40AliRDHFCuts(name),
addffa76 41 fCutOnResonances(kTRUE),
42 fPidOption(0),
43 fMaxPtStrongPid(0.),
44 fMaxPStrongPidK(0.),
45 fMaxPStrongPidpi(0.),
46 fDistToMaxProb(0.01),
47 fBayesThreshold(0.05),
48 fWeightKKpi(1.),
49 fWeightpiKK(1.)
42828f21 50{
51 //
52 // Default Constructor
53 //
ac8cab2d 54 Int_t nvars=20;
42828f21 55 SetNVars(nvars);
ac8cab2d 56 TString varNames[20]={"inv. mass [GeV]",
42828f21 57 "pTK [GeV/c]",
58 "pTPi [GeV/c]",
59 "d0K [cm]",
60 "d0Pi [cm]",
61 "dist12 [cm]",
62 "sigmavert [cm]",
63 "decLen [cm]",
64 "ptMax [GeV/c]",
65 "cosThetaPoint",
66 "Sum d0^2 (cm^2)",
67 "dca [cm]",
68 "inv. mass (Mphi-MKK) [GeV]",
69 "inv. mass (MKo*-MKpi) [GeV]",
70 "Abs(CosineKpiPhiRFrame)^3",
ac8cab2d 71 "CosPiDsLabFrame",
72 "decLenXY [cm]"
73 "NormdecLen",
74 "NormdecLenXY [cm]",
75 "cosThetaPointXY"};
42828f21 76
ac8cab2d 77 Bool_t isUpperCut[20]={kTRUE,
42828f21 78 kFALSE,
79 kFALSE,
80 kFALSE,
81 kFALSE,
82 kFALSE,
83 kTRUE,
84 kFALSE,
85 kFALSE,
86 kFALSE,
87 kFALSE,
88 kTRUE,
89 kTRUE,
90 kTRUE,
91 kFALSE,
ac8cab2d 92 kTRUE,
93 kFALSE,
94 kFALSE,
95 kFALSE,
96 kFALSE};
97 SetVarNames(20,varNames,isUpperCut);
98 Bool_t forOpt[20]={kFALSE,
42828f21 99 kFALSE,
100 kFALSE,
101 kFALSE,
102 kFALSE,
103 kFALSE,
104 kTRUE,
105 kTRUE,
106 kTRUE,
107 kTRUE,
108 kTRUE,
109 kFALSE,
110 kTRUE,
111 kTRUE,
112 kFALSE,
ac8cab2d 113 kFALSE,
114 kTRUE,
115 kTRUE,
116 kTRUE,
117 kTRUE};
118
119 SetVarsForOpt(11,forOpt);
42828f21 120 Float_t limits[2]={0,999999999.};
121 SetPtBins(2,limits);
122 if(fPidHF)delete fPidHF;
123 fPidHF=new AliAODPidHF();
124 Double_t plim[2]={0.6,0.8};
125 Double_t nsigma[5]={2.,1.,2.,3.,0.};
126
127 fPidHF->SetPLimit(plim);
128 fPidHF->SetAsym(kTRUE);
129 fPidHF->SetSigma(nsigma);
130 fPidHF->SetMatch(1);
131 fPidHF->SetTPC(1);
132 fPidHF->SetTOF(1);
133 fPidHF->SetITS(0);
134 fPidHF->SetTRD(0);
135 fPidHF->SetCompat(kTRUE);
136
137}
138//--------------------------------------------------------------------------
139AliRDHFCutsDstoKKpi::AliRDHFCutsDstoKKpi(const AliRDHFCutsDstoKKpi &source) :
140 AliRDHFCuts(source),
addffa76 141 fCutOnResonances(source.fCutOnResonances),
ac8cab2d 142 fPidOption(source.fPidOption),
143 fMaxPtStrongPid(source.fMaxPtStrongPid),
144 fMaxPStrongPidK(source.fMaxPStrongPidK),
addffa76 145 fMaxPStrongPidpi(source.fMaxPStrongPidpi),
146 fDistToMaxProb(source.fDistToMaxProb),
147 fBayesThreshold(source.fBayesThreshold),
148 fWeightKKpi(source.fWeightKKpi),
149 fWeightpiKK(source.fWeightpiKK)
42828f21 150{
151 //
152 // Copy constructor
153 //
154
155}
156//--------------------------------------------------------------------------
157AliRDHFCutsDstoKKpi &AliRDHFCutsDstoKKpi::operator=(const AliRDHFCutsDstoKKpi &source)
158{
159 //
160 // assignment operator
161 //
162 if(&source == this) return *this;
163
164 AliRDHFCuts::operator=(source);
addffa76 165
166 fCutOnResonances=source.fCutOnResonances;
42828f21 167 fPidOption=source.fPidOption;
ac8cab2d 168 fMaxPtStrongPid=source.fMaxPtStrongPid;
169 fMaxPStrongPidK=source.fMaxPStrongPidK;
170 fMaxPStrongPidpi=source.fMaxPStrongPidpi;
addffa76 171 fDistToMaxProb=source.fDistToMaxProb;
172 fBayesThreshold=source.fBayesThreshold;
173 fWeightKKpi=source.fWeightKKpi;
174 fWeightpiKK=source.fWeightpiKK;
42828f21 175
176 return *this;
177}
178
179
180//---------------------------------------------------------------------------
181void AliRDHFCutsDstoKKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent *aod) {
182 //
183 // Fills in vars the values of the variables
184 //
185
186 if(nvars!=fnVarsForOpt) {
187 printf("AliRDHFCutsDstoKKpi::GetCutsVarsForOpt: wrong number of variables\n");
188 return;
189 }
190
191 AliAODRecoDecayHF3Prong *dd = (AliAODRecoDecayHF3Prong*)d;
192
193 //recalculate vertex w/o daughters
194 Bool_t cleanvtx=kFALSE;
195 AliAODVertex *origownvtx=0x0;
196 if(fRemoveDaughtersFromPrimary) {
197 if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());
198 cleanvtx=kTRUE;
199 if(!RecalcOwnPrimaryVtx(dd,aod)) {
200 CleanOwnPrimaryVtx(dd,aod,origownvtx);
201 cleanvtx=kFALSE;
202 }
203 }
204
205 Int_t iter=-1;
206 if(fVarsForOpt[0]){
207 iter++;
208 if(TMath::Abs(pdgdaughters[0])==321){
209 vars[iter]=dd->InvMassDsKKpi();
210 }else{
211 vars[iter]=dd->InvMassDspiKK();
212 }
213 }
214 if(fVarsForOpt[1]){
215 iter++;
216 Float_t minPtDau=99999.;
217 for(Int_t iprong=0;iprong<3;iprong++){
218 if(TMath::Abs(pdgdaughters[iprong])==321 &&
219 dd->PtProng(iprong)<minPtDau) minPtDau=dd->PtProng(iprong);
220 }
221 vars[iter]=minPtDau;
222 }
223 if(fVarsForOpt[2]){
224 iter++;
225 for(Int_t iprong=0;iprong<3;iprong++){
226 if(TMath::Abs(pdgdaughters[iprong])==211) {
227 vars[iter]=dd->PtProng(iprong);
228 }
229 }
230 }
231 if(fVarsForOpt[3]){
232 iter++;
233 Float_t minImpParDau=99999.;
234 for(Int_t iprong=0;iprong<3;iprong++){
235 if(TMath::Abs(pdgdaughters[iprong])==321 &&
236 dd->Getd0Prong(iprong)<minImpParDau) minImpParDau=dd->Getd0Prong(iprong);
237 }
238 vars[iter]=minImpParDau;
239 }
240 if(fVarsForOpt[4]){
241 iter++;
242 for(Int_t iprong=0;iprong<3;iprong++){
243 if(TMath::Abs(pdgdaughters[iprong])==211) {
244 vars[iter]=dd->Getd0Prong(iprong);
245 }
246 }
247 }
248 if(fVarsForOpt[5]){
249 iter++;
250 Float_t minDistPair=TMath::Min(dd->GetDist12toPrim(),dd->GetDist23toPrim());
251 vars[iter]=minDistPair;
252 }
253 if(fVarsForOpt[6]){
254 iter++;
255 vars[iter]=dd->GetSigmaVert(aod);
256 }
257 if(fVarsForOpt[7]){
258 iter++;
259 vars[iter] = dd->DecayLength();
260 }
261 if(fVarsForOpt[8]){
262 iter++;
263 Float_t ptmax=0;
264 for(Int_t i=0;i<3;i++){
265 if(dd->PtProng(i)>ptmax)ptmax=dd->PtProng(i);
266 }
267 vars[iter]=ptmax;
268 }
269 if(fVarsForOpt[9]){
270 iter++;
271 vars[iter]=dd->CosPointingAngle();
272 }
273 if(fVarsForOpt[10]){
274 iter++;
275 vars[iter]=dd->Getd0Prong(0)*dd->Getd0Prong(0)+dd->Getd0Prong(1)*dd->Getd0Prong(1)+dd->Getd0Prong(2)*dd->Getd0Prong(2);
276 }
277 if(fVarsForOpt[11]){
278 iter++;
279 Float_t maxDCA=0.;
280 for(Int_t i=0;i<3;i++){
281 if(d->GetDCA(i)>maxDCA) maxDCA=d->GetDCA(i);
282 }
283 vars[iter]=maxDCA;
284 }
285 if(fVarsForOpt[12]){
286 iter++;
287 Double_t mPDGPhi = TDatabasePDG::Instance()->GetParticle(333)->Mass();
288 if(TMath::Abs(pdgdaughters[0])==321){
289
290 Double_t phimass01=d->InvMass2Prongs(0,1,321,321);
291 vars[iter]=TMath::Abs(phimass01-mPDGPhi);
292 // vars[iter]=dd->InvMass2Prongs(0,1,321,321);
293 }else{
294 Double_t phimass12=d->InvMass2Prongs(1,2,321,321);
295 vars[iter]=TMath::Abs(phimass12-mPDGPhi);
296 // vars[iter]=dd->InvMass2Prongs(1,2,321,321);
297 }
298 }
299 if(fVarsForOpt[13]){
300 iter++;
301 Double_t mPDGK0star = TDatabasePDG::Instance()->GetParticle(313)->Mass();
302 if(TMath::Abs(pdgdaughters[0])==321){
303
304 Double_t mass12kpi=d->InvMass2Prongs(1,2,321,211);
305 vars[iter]=TMath::Abs(mass12kpi-mPDGK0star);
306 // vars[iter]=dd->InvMass2Prongs(1,2,321,211);
307 }else{
308 Double_t mass01pik=d->InvMass2Prongs(0,1,211,321);
309 vars[iter]=TMath::Abs(mass01pik-mPDGK0star);
310 // vars[iter]=dd->InvMass2Prongs(0,1,211,321);
311 }
312 }
313 if(fVarsForOpt[14]){
314 iter++;
315 if(TMath::Abs(pdgdaughters[0])==321){
316 vars[iter]=dd->CosPiKPhiRFrameKKpi();
317 }else{
318 vars[iter]=dd->CosPiKPhiRFramepiKK();
319 }
320 }
321 if(fVarsForOpt[15]){
322 iter++;
323 if(TMath::Abs(pdgdaughters[0])==321){
324 vars[iter]=dd->CosPiDsLabFrameKKpi();
325 }else{
326 vars[iter]=dd->CosPiDsLabFramepiKK();
327 }
328 }
ac8cab2d 329
330 if(fVarsForOpt[16]){
331 iter++;
332 vars[iter]=dd->DecayLengthXY();
333 }
334
335 if(fVarsForOpt[17]){
336 iter++;
337 vars[iter]=dd->NormalizedDecayLength();
338 }
339
340 if(fVarsForOpt[18]){
341 iter++;
342 vars[iter]=dd->NormalizedDecayLengthXY();
343 }
344
345 if(fVarsForOpt[19]){
346 iter++;
347 vars[iter]=dd->CosPointingAngleXY();
348 }
42828f21 349
350 if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
351 return;
352}
353//---------------------------------------------------------------------------
354Bool_t AliRDHFCutsDstoKKpi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
355{
356 //
357 // Checking if Ds is in fiducial acceptance region
358 //
359
af5d687e 360 if(fMaxRapidityCand>-998.){
361 if(TMath::Abs(y) > fMaxRapidityCand) return kFALSE;
362 else return kTRUE;
363 }
364
42828f21 365 if(pt > 5.) {
366 // applying cut for pt > 5 GeV
367 AliDebug(2,Form("pt of Ds = %f (> 5), cutting at |y| < 0.8",pt));
368 if (TMath::Abs(y) > 0.8) return kFALSE;
369
370 } else {
371 // appliying smooth cut for pt < 5 GeV
372 Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
373 Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
374 AliDebug(2,Form("pt of Ds = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
375 if (y < minFiducialY || y > maxFiducialY) return kFALSE;
376 }
377
378 return kTRUE;
379}
380
addffa76 381//---------------------------------------------------------------------------
382Int_t AliRDHFCutsDstoKKpi::IsSelectedPIDBayes(AliAODRecoDecayHF *rd) {
383 Int_t retCode=3;
384 Bool_t okKKpi=kTRUE;
385 Bool_t okpiKK=kTRUE;
386 if(!fUsePID || !rd) return retCode;
387 if(!fPidHF){
388 AliWarning("AliAODPidHF not created!");
389 return retCode;
390 }
391 if(fPidOption!=kBayesianMaxProb && fPidOption!=kBayesianThreshold && fPidOption!=kBayesianWeights){
392 AliWarning("Wrong call to Bayesian PID");
393 return retCode;
394 }
395
396 AliPIDCombined* copid=fPidHF->GetPidCombined();
397 copid->SetDetectorMask(AliPIDResponse::kDetTPC | AliPIDResponse::kDetTOF);
398 AliPIDResponse* pidres=fPidHF->GetPidResponse();
399 Double_t bayesProb[AliPID::kSPECIES];
400 Int_t nKaons=0;
401 Int_t nNotKaons=0;
402 Int_t sign= rd->GetCharge();
403 fWeightKKpi=1.;
404 fWeightpiKK=1.;
405 for(Int_t iDaught=0; iDaught<3; iDaught++){
406 AliAODTrack *track=(AliAODTrack*)rd->GetDaughter(iDaught);
407
408 Int_t isPion=0;
409 Int_t isKaon=0;
410 Int_t isProton=0;
411 UInt_t usedDet=copid->ComputeProbabilities(track,pidres,bayesProb);
412 if(usedDet!=0){
413 if(fPidOption==kBayesianMaxProb){
414 Double_t maxProb=TMath::MaxElement(AliPID::kSPECIES,bayesProb);
415 if(TMath::Abs(maxProb-bayesProb[AliPID::kPion])<fDistToMaxProb) isPion=1;
416 else isPion=-1;
417 if(TMath::Abs(maxProb-bayesProb[AliPID::kKaon])<fDistToMaxProb) isKaon=1;
418 else isKaon=-1;
419 if(TMath::Abs(maxProb-bayesProb[AliPID::kProton])<fDistToMaxProb) isProton=1;
420 else isProton=-1;
421 }
422 if(fPidOption==kBayesianThreshold){
423 if(bayesProb[AliPID::kPion]>fBayesThreshold) isPion=1;
424 else isPion=-1;
425 if(bayesProb[AliPID::kKaon]>fBayesThreshold) isKaon=1;
426 else isKaon=-1;
427 if(bayesProb[AliPID::kProton]>fBayesThreshold) isProton=1;
428 else isProton=-1;
429 }
430 }
431 if(fPidOption==kBayesianWeights){ // store the probabilities in the case kBayesianWeights
432 if(iDaught==0){
433 fWeightKKpi*=bayesProb[AliPID::kKaon];
434 fWeightpiKK*=bayesProb[AliPID::kPion];
435 }else if(iDaught==1){
436 fWeightKKpi*=bayesProb[AliPID::kKaon];
437 fWeightpiKK*=bayesProb[AliPID::kKaon];
438 }else if(iDaught==2){
439 fWeightKKpi*=bayesProb[AliPID::kPion];
440 fWeightpiKK*=bayesProb[AliPID::kKaon];
441 }
442 }else{ // selection for the other cases
443 if(isProton>0 && isKaon<0 && isPion<0) return 0;
444 if(sign!=track->Charge()){// must be kaon
445 if(isKaon<0) return 0;
446 }
447 if(isKaon>0 && isPion<0) nKaons++;
448 if(isKaon<0) nNotKaons++;
449 if(iDaught==0){
450 if(isKaon<0) okKKpi=kFALSE;
451 if(isPion<0) okpiKK=kFALSE;
452 }else if(iDaught==2){
453 if(isKaon<0) okpiKK=kFALSE;
454 if(isPion<0) okKKpi=kFALSE;
455 }
456 }
457 }
458 if(fPidOption==kBayesianWeights) return retCode;
459
460 if(nKaons>2)return 0;
461 if(nNotKaons>1) return 0;
462
463 if(!okKKpi) retCode-=1;
464 if(!okpiKK) retCode-=2;
465
466 return retCode;
467}
468
42828f21 469//---------------------------------------------------------------------------
470Int_t AliRDHFCutsDstoKKpi::IsSelectedPID(AliAODRecoDecayHF *rd) {
471 // PID selection
472 // return values: 0->NOT OK, 1->OK as KKpi, 2->OK as piKK, 3->OK as both
473 Int_t retCode=3;
474 Bool_t okKKpi=kTRUE;
475 Bool_t okpiKK=kTRUE;
476 if(!fUsePID || !rd) return retCode;
477 if(!fPidHF){
478 AliWarning("AliAODPidHF not created!");
479 return retCode;
480 }
addffa76 481 if(fPidOption==kBayesianMaxProb || fPidOption==kBayesianThreshold || fPidOption==kBayesianWeights){
482 // call method for Bayesian probability
483 return IsSelectedPIDBayes(rd);
484 }
42828f21 485
486 Double_t origCompatTOF=fPidHF->GetPCompatTOF();
487 Double_t origThreshTPC=fPidHF->GetPtThresholdTPC();
488 if(fPidOption==kStrong){
489 fPidHF->SetPCompatTOF(999999.);
490 fPidHF->SetPtThresholdTPC(999999.);
491 }
492
42828f21 493 Int_t nKaons=0;
494 Int_t nNotKaons=0;
495 Int_t sign= rd->GetCharge();
496 for(Int_t iDaught=0; iDaught<3; iDaught++){
497 AliAODTrack *track=(AliAODTrack*)rd->GetDaughter(iDaught);
498
499 Int_t isPion=fPidHF->MakeRawPid(track,AliPID::kPion);
500 Int_t isKaon=fPidHF->MakeRawPid(track,AliPID::kKaon);
501 Int_t isProton=fPidHF->MakeRawPid(track,AliPID::kProton);
502
503 if(isProton>0 && isKaon<0 && isPion<0){
504 fPidHF->SetPCompatTOF(origCompatTOF);
505 fPidHF->SetPtThresholdTPC(origThreshTPC);
506 return 0;
507 }
508 if(sign!=track->Charge()){// must be kaon
509 if(isKaon<0){
510 fPidHF->SetPCompatTOF(origCompatTOF);
511 fPidHF->SetPtThresholdTPC(origThreshTPC);
512 return 0;
513 }
ac8cab2d 514 if(fPidOption==kStrong && rd->Pt()<fMaxPtStrongPid && isKaon<=0){
42828f21 515 fPidHF->SetPCompatTOF(origCompatTOF);
516 fPidHF->SetPtThresholdTPC(origThreshTPC);
517 return 0;
518 }
ac8cab2d 519 if(fPidOption==kStrongPDep && rd->Pt()<fMaxPtStrongPid){
520 if(isKaon<=0 && track->P()<fMaxPStrongPidK) return 0;
521 }
42828f21 522 }
ac8cab2d 523
42828f21 524 if(isKaon>0 && isPion<0) nKaons++;
525 if(isKaon<0) nNotKaons++;
526 if(iDaught==0){
527 if(isKaon<0) okKKpi=kFALSE;
528 if(isPion<0) okpiKK=kFALSE;
ac8cab2d 529 if(fPidOption==kStrong && rd->Pt()<fMaxPtStrongPid){
42828f21 530 if(isKaon<=0) okKKpi=kFALSE;
531 if(isPion<=0) okpiKK=kFALSE;
532 }
ac8cab2d 533 if(fPidOption==kStrongPDep && rd->Pt()<fMaxPtStrongPid){
534 if(isKaon<=0 && track->P()<fMaxPStrongPidK) okKKpi=kFALSE;
535 if(isPion<=0 && track->P()<fMaxPStrongPidpi) okpiKK=kFALSE;
536 }
42828f21 537 }
538 else if(iDaught==2){
539 if(isKaon<0) okpiKK=kFALSE;
540 if(isPion<0) okKKpi=kFALSE;
ac8cab2d 541 if(fPidOption==kStrong && rd->Pt()<fMaxPtStrongPid){
542 if(isKaon<=0) okpiKK=kFALSE;
543 if(isPion<=0) okKKpi=kFALSE;
544 }
545 if(fPidOption==kStrongPDep && rd->Pt()<fMaxPtStrongPid){
546 if(isKaon<=0 && track->P()<fMaxPStrongPidK) okpiKK=kFALSE;
547 if(isPion<=0 && track->P()<fMaxPStrongPidpi) okKKpi=kFALSE;
42828f21 548 }
549 }
550 }
551
552 fPidHF->SetPCompatTOF(origCompatTOF);
553 fPidHF->SetPtThresholdTPC(origThreshTPC);
554
555 if(nKaons>2)return 0;
556 if(nNotKaons>1) return 0;
557
558 if(!okKKpi) retCode-=1;
559 if(!okpiKK) retCode-=2;
560
561 return retCode;
562}
563
564//---------------------------------------------------------------------------
565Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel, AliAODEvent* aod) {
566 //
567 // Apply selection
568 //
569
570 if(!fCutsRD){
571 cout<<"Cut matrix not inizialized. Exit..."<<endl;
572 return 0;
573 }
574 //PrintAll();
575 AliAODRecoDecayHF3Prong* d=(AliAODRecoDecayHF3Prong*)obj;
576
577 if(!d){
578 cout<<"AliAODRecoDecayHF3Prong null"<<endl;
579 return 0;
580 }
581
582 if(fKeepSignalMC) if(IsSignalMC(d,aod,431)) return 3;
583
584 Double_t ptD=d->Pt();
585 if(ptD<fMinPtCand) return 0;
586 if(ptD>fMaxPtCand) return 0;
587
3e075b37 588 if(fUseTrackSelectionWithFilterBits && d->HasBadDaughters()) return 0;
42828f21 589
590
591 // selection on daughter tracks
592 if(selectionLevel==AliRDHFCuts::kAll ||
593 selectionLevel==AliRDHFCuts::kTracks) {
594 if(!AreDaughtersSelected(d)) return 0;
595 }
596
597
598
599
600 // selection on candidate
601 if(selectionLevel==AliRDHFCuts::kAll ||
602 selectionLevel==AliRDHFCuts::kCandidate) {
603 //recalculate vertex w/o daughters
604 AliAODVertex *origownvtx=0x0;
605 if(fRemoveDaughtersFromPrimary) {
606 if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
607 if(!RecalcOwnPrimaryVtx(d,aod)) {
608 CleanOwnPrimaryVtx(d,aod,origownvtx);
609 return 0;
610 }
611 }
612
613 Int_t okDsKKpi=1;
614 Int_t okDspiKK=1;
615 Int_t okMassPhiKKpi=0;
616 Int_t okMassPhipiKK=0;
617 Int_t okMassK0starKKpi=0;
618 Int_t okMassK0starpiKK=0;
619 Int_t okDsPhiKKpi=0;
620 Int_t okDsPhipiKK=0;
621 Int_t okDsK0starKKpi=0;
622 Int_t okDsK0starpiKK=0;
623
624 Double_t pt=d->Pt();
625 Int_t ptbin=PtBin(pt);
626 if (ptbin==-1) {
627 CleanOwnPrimaryVtx(d,aod,origownvtx);
628 return 0;
629 }
630
631 Double_t mDsPDG = TDatabasePDG::Instance()->GetParticle(431)->Mass();
632 Double_t mDsKKpi=d->InvMassDsKKpi();
633 Double_t mDspiKK=d->InvMassDspiKK();
634 if(TMath::Abs(mDsKKpi-mDsPDG)>fCutsRD[GetGlobalIndex(0,ptbin)]) okDsKKpi = 0;
635 if(TMath::Abs(mDspiKK-mDsPDG)>fCutsRD[GetGlobalIndex(0,ptbin)]) okDspiKK = 0;
636 if(!okDsKKpi && !okDspiKK){
637 CleanOwnPrimaryVtx(d,aod,origownvtx);
638 return 0;
639 }
640
641
642
643 // cuts on resonant decays (via Phi or K0*)
addffa76 644 if(fCutOnResonances){
645 Double_t mPhiPDG = TDatabasePDG::Instance()->GetParticle(333)->Mass();
646 Double_t mK0starPDG = TDatabasePDG::Instance()->GetParticle(313)->Mass();
647 if(okDsKKpi){
648 Double_t mass01phi=d->InvMass2Prongs(0,1,321,321);
649 Double_t mass12K0s=d->InvMass2Prongs(1,2,321,211);
650 if(TMath::Abs(mass01phi-mPhiPDG)<fCutsRD[GetGlobalIndex(12,ptbin)]) okMassPhiKKpi=1;
651 if(TMath::Abs(mass12K0s-mK0starPDG)<fCutsRD[GetGlobalIndex(13,ptbin)]) okMassK0starKKpi = 1;
652 if(!okMassPhiKKpi && !okMassK0starKKpi) okDsKKpi=0;
653 if(okMassPhiKKpi) okDsPhiKKpi=1;
654 if(okMassK0starKKpi) okDsK0starKKpi=1;
655 }
656 if(okDspiKK){
657 Double_t mass01K0s=d->InvMass2Prongs(0,1,211,321);
658 Double_t mass12phi=d->InvMass2Prongs(1,2,321,321);
659 if(TMath::Abs(mass01K0s-mK0starPDG)<fCutsRD[GetGlobalIndex(13,ptbin)]) okMassK0starpiKK = 1;
660 if(TMath::Abs(mass12phi-mPhiPDG)<fCutsRD[GetGlobalIndex(12,ptbin)]) okMassPhipiKK=1;
661 if(!okMassPhipiKK && !okMassK0starpiKK) okDspiKK=0;
662 if(okMassPhipiKK) okDsPhipiKK=1;
663 if(okMassK0starpiKK) okDsK0starpiKK=1;
664 }
665 if(!okDsKKpi && !okDspiKK){
666 CleanOwnPrimaryVtx(d,aod,origownvtx);
667 return 0;
668 }
42828f21 669 }
670
671 // Cuts on track pairs
672 for(Int_t i=0;i<3;i++){
673 if(d->GetDCA(i)>fCutsRD[GetGlobalIndex(11,ptbin)]){
674 CleanOwnPrimaryVtx(d,aod,origownvtx);
675 return 0;
676 }
677 }
678 if(d->GetDist12toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)] ||
679 d->GetDist23toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)]){
680 CleanOwnPrimaryVtx(d,aod,origownvtx);
681 return 0;
682 }
683
684
685
686 //single track
687 if(TMath::Abs(d->Pt2Prong(1)) < fCutsRD[GetGlobalIndex(1,ptbin)]*fCutsRD[GetGlobalIndex(1,ptbin)] ||
688 TMath::Abs(d->Getd0Prong(1))<fCutsRD[GetGlobalIndex(3,ptbin)]){
689 CleanOwnPrimaryVtx(d,aod,origownvtx);
690 return 0;
691 }
692
693 if(okDsKKpi){
694 if(TMath::Abs(d->Pt2Prong(0)) < fCutsRD[GetGlobalIndex(1,ptbin)]*fCutsRD[GetGlobalIndex(1,ptbin)] ||
695 TMath::Abs(d->Getd0Prong(0))<fCutsRD[GetGlobalIndex(3,ptbin)]) okDsKKpi=0;
696 if(TMath::Abs(d->Pt2Prong(2)) < fCutsRD[GetGlobalIndex(2,ptbin)]*fCutsRD[GetGlobalIndex(2,ptbin)] ||
697 TMath::Abs(d->Getd0Prong(2))<fCutsRD[GetGlobalIndex(4,ptbin)]) okDsKKpi=0;
698 }
699 if(okDspiKK){
700 if(TMath::Abs(d->Pt2Prong(0)) < fCutsRD[GetGlobalIndex(2,ptbin)]*fCutsRD[GetGlobalIndex(2,ptbin)] ||
701 TMath::Abs(d->Getd0Prong(0))<fCutsRD[GetGlobalIndex(4,ptbin)]) okDspiKK=0;
702 if(TMath::Abs(d->Pt2Prong(2)) < fCutsRD[GetGlobalIndex(1,ptbin)]*fCutsRD[GetGlobalIndex(1,ptbin)] ||
703 TMath::Abs(d->Getd0Prong(2))<fCutsRD[GetGlobalIndex(3,ptbin)]) okDspiKK=0;
704 }
705 if(!okDsKKpi && !okDspiKK){
706 CleanOwnPrimaryVtx(d,aod,origownvtx);
707 return 0;
708 }
709
710 // Cuts on candidate triplet
711
712
713 if(d->CosPointingAngle()< fCutsRD[GetGlobalIndex(9,ptbin)]){
714 CleanOwnPrimaryVtx(d,aod,origownvtx);
715 return 0;
716 }
717
718 if(d->Pt2Prong(0)<fCutsRD[GetGlobalIndex(8,ptbin)]*fCutsRD[GetGlobalIndex(8,ptbin)] &&
719 d->Pt2Prong(1)<fCutsRD[GetGlobalIndex(8,ptbin)]*fCutsRD[GetGlobalIndex(8,ptbin)] &&
720 d->Pt2Prong(2)<fCutsRD[GetGlobalIndex(8,ptbin)]*fCutsRD[GetGlobalIndex(8,ptbin)]) {
721 CleanOwnPrimaryVtx(d,aod,origownvtx);
722 return 0;
723 }
724
725 if(d->DecayLength2()<fCutsRD[GetGlobalIndex(7,ptbin)]*fCutsRD[GetGlobalIndex(7,ptbin)]){
726 CleanOwnPrimaryVtx(d,aod,origownvtx);
727 return 0;
728 }
729
730
731 Double_t sum2=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
732 if(sum2<fCutsRD[GetGlobalIndex(10,ptbin)]){
733 CleanOwnPrimaryVtx(d,aod,origownvtx);
734 return 0;
735 }
736
737
738 //sec vert
739 Double_t sigmavert=d->GetSigmaVert(aod);
740 if(sigmavert>fCutsRD[GetGlobalIndex(6,ptbin)]){
741 CleanOwnPrimaryVtx(d,aod,origownvtx);
742 return 0;
743 }
ac8cab2d 744
745 // decay length XY
746 if(d->DecayLengthXY()<fCutsRD[GetGlobalIndex(16,ptbin)]){
747 CleanOwnPrimaryVtx(d,aod,origownvtx);
748 return 0;
749 }
750
751 //norm decay length
752 if(d->NormalizedDecayLength()<fCutsRD[GetGlobalIndex(17,ptbin)]){
753 CleanOwnPrimaryVtx(d,aod,origownvtx);
754 return 0;
755 }
756
757 //norm decay length XY
758 if(d->NormalizedDecayLengthXY()<fCutsRD[GetGlobalIndex(18,ptbin)]){
759 CleanOwnPrimaryVtx(d,aod,origownvtx);
760 return 0;
761 }
762
763 //cos pointing XY
764 if(d->CosPointingAngleXY()<fCutsRD[GetGlobalIndex(19,ptbin)]){
765 CleanOwnPrimaryVtx(d,aod,origownvtx);
766 return 0;
767 }
768
42828f21 769
770 if(okDsKKpi){
771 Double_t cosPiKPhiRFKKpi=d->CosPiKPhiRFrameKKpi();
772 Double_t kincutPiKPhiKKpi=TMath::Abs(cosPiKPhiRFKKpi*cosPiKPhiRFKKpi*cosPiKPhiRFKKpi);
773 if(kincutPiKPhiKKpi<fCutsRD[GetGlobalIndex(14,ptbin)]) okDsKKpi=0;
774 }
775 if(okDspiKK){
776 Double_t cosPiKPhiRFpiKK=d->CosPiKPhiRFramepiKK();
777 Double_t kincutPiKPhipiKK=TMath::Abs(cosPiKPhiRFpiKK*cosPiKPhiRFpiKK*cosPiKPhiRFpiKK);
778 if(kincutPiKPhipiKK<fCutsRD[GetGlobalIndex(14,ptbin)]) okDspiKK=0;
779 }
780 if(!okDsKKpi && !okDspiKK){
781 CleanOwnPrimaryVtx(d,aod,origownvtx);
782 return 0;
783 }
784
785
786
787 if(okDsKKpi){
788 Double_t cosPiDsLabFrameKKpi=d->CosPiDsLabFrameKKpi();
789 if(cosPiDsLabFrameKKpi>fCutsRD[GetGlobalIndex(15,ptbin)]) okDsKKpi=0;
790 }
791 if(okDspiKK){
792 Double_t cosPiDsLabFramepiKK=d->CosPiDsLabFramepiKK();
793 if(cosPiDsLabFramepiKK>fCutsRD[GetGlobalIndex(15,ptbin)]) okDspiKK=0;
794 }
795 if(!okDsKKpi && !okDspiKK){
796 CleanOwnPrimaryVtx(d,aod,origownvtx);
797 return 0;
798 }
799
800 // unset recalculated primary vertex when not needed any more
801 CleanOwnPrimaryVtx(d,aod,origownvtx);
802
803
804
805 if(!okDsKKpi){
806 okDsPhiKKpi=0;
807 okDsK0starKKpi=0;
808 }
809 if(!okDspiKK){
810 okDsPhipiKK=0;
811 okDsK0starpiKK=0;
812 }
813
814 // PID selection
815 Int_t returnvaluePID=3;
816 if(selectionLevel==AliRDHFCuts::kAll ||
817 selectionLevel==AliRDHFCuts::kCandidate ||
818 selectionLevel==AliRDHFCuts::kPID) {
819 returnvaluePID = IsSelectedPID(d);
820 fIsSelectedPID=returnvaluePID;
821 }
822 if(returnvaluePID==0)return 0;
823
824 Bool_t okPidDsKKpi=returnvaluePID&1;
825 Bool_t okPidDspiKK=returnvaluePID&2;
826 if(!okPidDsKKpi){
827 okDsPhiKKpi=0;
828 okDsK0starKKpi=0;
829 }
830 if(!okPidDspiKK){
831 okDsPhipiKK=0;
832 okDsK0starpiKK=0;
833 }
834
835 if((okPidDsKKpi && okDsKKpi)||(okPidDspiKK && okDspiKK)){
836 Int_t returnvalue=0;
837 if(okDsKKpi) returnvalue+=1;
838 if(okDspiKK) returnvalue+=2;
839 if(okDsPhiKKpi) returnvalue+=4;
840 if(okDsPhipiKK) returnvalue+=8;
841 if(okDsK0starKKpi) returnvalue+=16;
842 if(okDsK0starpiKK) returnvalue+=32;
843 return returnvalue;
844 }else{
845 return 0;
846 }
847 }
848 return 15;
849
850}
851
852//--------------------------------------------------------------------------
853
854UInt_t AliRDHFCutsDstoKKpi::GetPIDTrackTPCTOFBitMap(AliAODTrack *track) const{
855
856 UInt_t bitmap=0;
857
858 Double_t sigmaTPCPionHyp=-999.;
859 Double_t sigmaTPCKaonHyp=-999.;
860 Double_t sigmaTPCProtonHyp=-999.;
861 Double_t sigmaTOFPionHyp=-999.;
862 Double_t sigmaTOFKaonHyp=-999.;
863 Double_t sigmaTOFProtonHyp=-999.;
864
865 Int_t oksigmaTPCPionHyp=fPidHF->GetnSigmaTPC(track,2,sigmaTPCPionHyp);
866 Int_t oksigmaTPCKaonHyp=fPidHF->GetnSigmaTPC(track,3,sigmaTPCKaonHyp);
867 Int_t oksigmaTPCProtonHyp=fPidHF->GetnSigmaTPC(track,4,sigmaTPCProtonHyp);
868 Int_t oksigmaTOFPionHyp=fPidHF->GetnSigmaTOF(track,2,sigmaTOFPionHyp);
869 Int_t oksigmaTOFKaonHyp=fPidHF->GetnSigmaTOF(track,3,sigmaTOFKaonHyp);
870 Int_t oksigmaTOFProtonHyp=fPidHF->GetnSigmaTOF(track,4,sigmaTOFProtonHyp);
871
4940d5bf 872 sigmaTPCPionHyp=TMath::Abs(sigmaTPCPionHyp);
873 sigmaTPCKaonHyp=TMath::Abs(sigmaTPCKaonHyp);
874 sigmaTPCProtonHyp=TMath::Abs(sigmaTPCProtonHyp);
875 sigmaTOFPionHyp=TMath::Abs(sigmaTOFPionHyp);
876 sigmaTOFKaonHyp=TMath::Abs(sigmaTOFKaonHyp);
2d719fc0 877 sigmaTOFProtonHyp=TMath::Abs(sigmaTOFProtonHyp);
4940d5bf 878
879 if (oksigmaTPCPionHyp && sigmaTPCPionHyp>0.){
42828f21 880 if (sigmaTPCPionHyp<1.) bitmap+=1<<kTPCPionLess1;
881 else{
882 if (sigmaTPCPionHyp<2.) bitmap+=1<<kTPCPionMore1Less2;
883 else {
884 if (sigmaTPCPionHyp<3.) bitmap+=1<<kTPCPionMore2Less3;
885 else bitmap+=1<<kTPCPionMore3;
886 }
887 }
888 }
889
890 if (oksigmaTPCKaonHyp && sigmaTPCKaonHyp>0.){
891 if (sigmaTPCKaonHyp<1.) bitmap+=1<<kTPCKaonLess1;
892 else{
893 if (sigmaTPCKaonHyp<2.) bitmap+=1<<kTPCKaonMore1Less2;
894 else {
895 if (sigmaTPCKaonHyp<3.) bitmap+=1<<kTPCKaonMore2Less3;
896 else bitmap+=1<<kTPCKaonMore3;
897 }
898 }
899 }
900
901 if (oksigmaTPCProtonHyp && sigmaTPCProtonHyp>0.){
902 if (sigmaTPCProtonHyp<1.) bitmap+=1<<kTPCProtonLess1;
903 else{
904 if (sigmaTPCProtonHyp<2.) bitmap+=1<<kTPCProtonMore1Less2;
905 else {
906 if (sigmaTPCProtonHyp<3.) bitmap+=1<<kTPCProtonMore2Less3;
907 else bitmap+=1<<kTPCProtonMore3;
908 }
909 }
910 }
911
912 if (oksigmaTOFPionHyp && sigmaTOFPionHyp>0.){
913 if (sigmaTOFPionHyp<1.) bitmap+=1<<kTOFPionLess1;
914 else{
915 if (sigmaTOFPionHyp<2.) bitmap+=1<<kTOFPionMore1Less2;
916 else {
917 if (sigmaTOFPionHyp<3.) bitmap+=1<<kTOFPionMore2Less3;
918 else bitmap+=1<<kTOFPionMore3;
919 }
920 }
921 }
922
923 if (oksigmaTOFKaonHyp && sigmaTOFKaonHyp>0.){
924 if (sigmaTOFKaonHyp<1.) bitmap+=1<<kTOFKaonLess1;
925 else{
926 if (sigmaTOFKaonHyp<2.) bitmap+=1<<kTOFKaonMore1Less2;
927 else {
928 if (sigmaTOFKaonHyp<3.) bitmap+=1<<kTOFKaonMore2Less3;
929 else bitmap+=1<<kTOFKaonMore3;
930 }
931 }
932 }
933
934 if (oksigmaTOFProtonHyp && sigmaTOFProtonHyp>0.){
935 if (sigmaTOFProtonHyp<1.) bitmap+=1<<kTOFProtonLess1;
936 else{
937 if (sigmaTOFProtonHyp<2.) bitmap+=1<<kTOFProtonMore1Less2;
938 else {
939 if (sigmaTOFProtonHyp<3.) bitmap+=1<<kTOFProtonMore2Less3;
940 else bitmap+=1<<kTOFProtonMore3;
941 }
942 }
943 }
944
945
946
947 return bitmap;
948
949}
950
951
952//---------------------------------------------------------------------------
953
954void AliRDHFCutsDstoKKpi::SetStandardCutsPP2010() {
955 //
956 //STANDARD CUTS USED FOR 2010 pp analysis
957 //
958
959 SetName("DstoKKpiCutsStandard");
960 SetTitle("Standard Cuts for D+s analysis");
961
962 // PILE UP REJECTION
963 SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
964
965 // EVENT CUTS
966 SetMinVtxContr(1);
967
968 AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
969 esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
970 //default
971 esdTrackCuts->SetRequireTPCRefit(kTRUE);
972 esdTrackCuts->SetRequireITSRefit(kTRUE);
973 //esdTrackCuts->SetMinNClustersITS(4); // default is 5
974 esdTrackCuts->SetMinNClustersTPC(70);
975 esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
976 AliESDtrackCuts::kAny);
977 // default is kBoth, otherwise kAny
978 esdTrackCuts->SetMinDCAToVertexXY(0.);
979 esdTrackCuts->SetPtRange(0.3,1.e10);
980
981 AddTrackCuts(esdTrackCuts);
57d15652 982 delete esdTrackCuts;
983 esdTrackCuts=NULL;
42828f21 984
985
986 const Int_t nptbins=4;
8195786d 987 Float_t ptbins[nptbins+1];
42828f21 988 ptbins[0]=2.;
989 ptbins[1]=4.;
990 ptbins[2]=6.;
991 ptbins[3]=8.;
992 ptbins[4]=12.;
993
ac8cab2d 994 const Int_t nvars=20;
42828f21 995
996 Float_t** anacutsval;
997 anacutsval=new Float_t*[nvars];
998
999 for(Int_t ic=0;ic<nvars;ic++){anacutsval[ic]=new Float_t[nptbins];}
1000 for(Int_t ipt=0;ipt<nptbins;ipt++){
1001
1002 anacutsval[0][ipt]=0.35;
1003 anacutsval[1][ipt]=0.3;
1004 anacutsval[2][ipt]=0.3;
1005 anacutsval[3][ipt]=0.;
1006 anacutsval[4][ipt]=0.;
1007 anacutsval[5][ipt]=0.005;
1008 anacutsval[8][ipt]=0.;
1009 anacutsval[10][ipt]=0.;
1010 anacutsval[11][ipt]=1000.0;
1011 anacutsval[13][ipt]=0.1;
ac8cab2d 1012 anacutsval[16][ipt]=0.;
1013 anacutsval[17][ipt]=0.;
1014 anacutsval[18][ipt]=0.;
1015 anacutsval[19][ipt]=-1.;
1016
1017
42828f21 1018 }
1019
1020 //sigmavertex
1021
1022 anacutsval[6][0]=0.020;
1023 anacutsval[6][1]=0.030;
1024 anacutsval[6][2]=0.030;
1025 anacutsval[6][3]=0.060;
1026
1027 //Decay length
1028
1029 anacutsval[7][0]=0.035;
1030 anacutsval[7][1]=0.035;
1031 anacutsval[7][2]=0.040;
1032 anacutsval[7][3]=0.040;
1033
1034 //cosThetaPoint
1035
1036 anacutsval[9][0]=0.94;
1037 anacutsval[9][1]=0.94;
1038 anacutsval[9][2]=0.94;
1039 anacutsval[9][3]=0.94;
1040
1041 //phi DeltaMass
1042
1043 anacutsval[12][0]=0.0080;
1044 anacutsval[12][1]=0.0050;
1045 anacutsval[12][2]=0.0045;
1046 anacutsval[12][3]=0.0090;
1047
1048 //Kin1
1049
1050 anacutsval[14][0]=0.10;
1051 anacutsval[14][1]=0.05;
1052 anacutsval[14][2]=0.0;
1053 anacutsval[14][3]=0.05;
1054
1055 //Kin2
1056
1057 anacutsval[15][0]=0.95;
1058 anacutsval[15][1]=0.95;
1059 anacutsval[15][2]=1.;
1060 anacutsval[15][3]=0.95;
1061
2a6fcbb7 1062 fPidHF->SetOldPid(kTRUE);
42828f21 1063 SetUsePID(kTRUE);
1064 SetPidOption(1);
ac8cab2d 1065 SetMaxPtStrongPid(9999.);
42828f21 1066 SetGlobalIndex(nvars,nptbins);
1067 SetPtBins(nptbins+1,ptbins);
1068 SetCuts(nvars,nptbins,anacutsval);
1069 SetRemoveDaughtersFromPrim(kTRUE);
1070
1071 PrintAll();
1072
1073 for(Int_t iic=0;iic<nvars;iic++){delete [] anacutsval[iic];}
1074 delete [] anacutsval;
1075 anacutsval=NULL;
1076
1077 return;
1078}
1079