]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/AliRDHFCutsD0toKpipipi.cxx
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliRDHFCutsD0toKpipipi.cxx
CommitLineData
ba27eecd 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
af5d687e 16/* $Id$ */
17
ba27eecd 18/////////////////////////////////////////////////////////////
19//
20// Class for cuts on AOD reconstructed D0->Kpipipi
21//
601736df 22// Author: r.romita@gsi.de, andrea.dainese@pd.infn.it,
ba27eecd 23// fabio.colamaria@ba.infn.it
24/////////////////////////////////////////////////////////////
25
26#include <TDatabasePDG.h>
27#include <Riostream.h>
28
29#include "AliRDHFCutsD0toKpipipi.h"
30#include "AliAODRecoDecayHF4Prong.h"
31#include "AliAODTrack.h"
32#include "AliESDtrack.h"
601736df 33#include "AliAODPidHF.h"
ba27eecd 34
c64cb1f6 35using std::cout;
36using std::endl;
37
ba27eecd 38ClassImp(AliRDHFCutsD0toKpipipi)
39
40//--------------------------------------------------------------------------
41AliRDHFCutsD0toKpipipi::AliRDHFCutsD0toKpipipi(const char* name) :
42AliRDHFCuts(name)
43{
44 //
45 // Default Constructor
46 //
47 Int_t nvars=9;
48 SetNVars(nvars);
49 TString varNames[9]={"inv. mass [GeV]",
50 "dca [cm]",
51 "Dist 2-trk Vtx to PrimVtx [cm]",
52 "Dist 3-trk Vtx to PrimVtx [cm]",
53 "Dist 4-trk Vtx to PrimVtx [cm]",
54 "cosThetaPoint",
55 "pt [GeV/c]",
56 "rho mass [GeV]",
57 "PID cut"};
58 Bool_t isUpperCut[9]={kTRUE,
59 kTRUE,
60 kFALSE,
61 kFALSE,
62 kFALSE,
63 kFALSE,
64 kFALSE,
65 kTRUE,
66 kFALSE};
67 SetVarNames(nvars,varNames,isUpperCut);
68 Bool_t forOpt[9]={kFALSE,
69 kTRUE,
70 kTRUE,
71 kTRUE,
72 kTRUE,
73 kTRUE,
74 kFALSE,
75 kFALSE,
76 kFALSE};
77 SetVarsForOpt(5,forOpt);
78 Float_t limits[2]={0,999999999.};
79 SetPtBins(2,limits);
80}
81//--------------------------------------------------------------------------
82AliRDHFCutsD0toKpipipi::AliRDHFCutsD0toKpipipi(const AliRDHFCutsD0toKpipipi &source) :
83 AliRDHFCuts(source)
84{
85 //
86 // Copy constructor
87 //
88
89}
90//--------------------------------------------------------------------------
91AliRDHFCutsD0toKpipipi &AliRDHFCutsD0toKpipipi::operator=(const AliRDHFCutsD0toKpipipi &source)
92{
93 //
94 // assignment operator
95 //
96 if(&source == this) return *this;
97
98 AliRDHFCuts::operator=(source);
99
100 return *this;
101}
102
103
104//---------------------------------------------------------------------------
105void AliRDHFCutsD0toKpipipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent* aod) {
106 //
107 // Fills in vars the values of the variables
108 //
109
110 if(nvars!=fnVarsForOpt) {
111 printf("AliRDHFCutsD0toKpipipi::GetCutsVarsForOpt: wrong number of variables\n");
112 return;
113 }
114
115 AliAODRecoDecayHF4Prong *dd = (AliAODRecoDecayHF4Prong*)d;
116
a6003e0a 117 //recalculate vertex w/o daughters
118 Bool_t cleanvtx=kFALSE;
119 AliAODVertex *origownvtx=0x0;
120 if(fRemoveDaughtersFromPrimary) {
121 if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());
122 cleanvtx=kTRUE;
123 if(!RecalcOwnPrimaryVtx(dd,aod)) {
124 CleanOwnPrimaryVtx(dd,aod,origownvtx);
125 cleanvtx=kFALSE;
126 }
127 }
128
ba27eecd 129 Int_t iter=-1;
130
131 if(fVarsForOpt[0]) {
132 iter++;
133 Double_t mD0[2],mD0bar[2];
134 if(TMath::Abs(pdgdaughters[1])==321 || TMath::Abs(pdgdaughters[3])==321) {
135 dd->InvMassD0(mD0);
136 if(TMath::Abs(pdgdaughters[1])==321) {
137 vars[iter]=mD0[0];
138 }else{
139 vars[iter]=mD0[1];
140 }
141 } else {
142 dd->InvMassD0bar(mD0bar);
143 if(TMath::Abs(pdgdaughters[0])==321) {
144 vars[iter]=mD0bar[0];
145 }else{
146 vars[iter]=mD0bar[1];
147 }
148 }
149 }
150
151 if(fVarsForOpt[1]){
152 iter++;
153 vars[iter]=dd->GetDCA();
154 }
155
156 if(fVarsForOpt[2]){
157 iter++;
158 vars[iter]=dd->GetDist12toPrim();
159 }
160 if(fVarsForOpt[3]){
161 iter++;
162 vars[iter]=dd->GetDist3toPrim();
163 }
164 if(fVarsForOpt[4]){
165 iter++;
166 vars[iter]=dd->GetDist4toPrim();
167 }
168 if(fVarsForOpt[5]){
169 iter++;
170 vars[iter]=dd->CosPointingAngle();
171 }
172 if(fVarsForOpt[6]){
173 iter++;
174 vars[iter]=dd->Pt();
175 }
176 if(fVarsForOpt[7]){
177 iter++;
178 vars[iter]=999999999.;
179 printf("ERROR: optmization for rho mass cut not implemented\n");
180 }
181 if(fVarsForOpt[8]){
182 iter++;
183 vars[iter]=999999999.;
184 printf("ERROR: optmization for PID cut not implemented\n");
185 }
186
a6003e0a 187 if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
ba27eecd 188 return;
189}
190//---------------------------------------------------------------------------
191Int_t AliRDHFCutsD0toKpipipi::IsSelected(TObject* obj,Int_t selectionLevel) {
192 //
193 // Apply selection
194 //
195
196 if(!fCutsRD){
197 cout<<"Cut matrix not inizialized. Exit..."<<endl;
198 return 0;
199 }
200 //PrintAll();
201 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;
202
203 if(!d){
204 cout<<"AliAODRecoDecayHF4Prong null"<<endl;
205 return 0;
206 }
207
208 Double_t ptD=d->Pt();
209 if(ptD<fMinPtCand) return 0;
601736df 210 if(ptD>fMaxPtCand) return 0;
ba27eecd 211
3e075b37 212 if(fUseTrackSelectionWithFilterBits && d->HasBadDaughters()) return 0;
ba27eecd 213
214 // selection on daughter tracks
215 if(selectionLevel==AliRDHFCuts::kAll ||
216 selectionLevel==AliRDHFCuts::kTracks) {
217 if(!AreDaughtersSelected(d)) return 0;
218 }
219
220
221 Int_t returnvalue=1;
222
223 // selection on candidate
224 if(selectionLevel==AliRDHFCuts::kAll ||
225 selectionLevel==AliRDHFCuts::kCandidate) {
226
227 Int_t ptbin=PtBin(d->Pt());
601736df 228
ba27eecd 229 Int_t okD0=1,okD0bar=1;
230 Double_t mD0[2],mD0bar[2];
231 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
232
233 d->InvMassD0(mD0);
234 if(TMath::Abs(mD0[0]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)] &&
235 TMath::Abs(mD0[1]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
236 d->InvMassD0bar(mD0bar);
237 if(TMath::Abs(mD0bar[0]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)] &&
238 TMath::Abs(mD0bar[1]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;
239 if(!okD0 && !okD0bar) return 0;
240
601736df 241 if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]
242 || d->GetDCA(3) > fCutsRD[GetGlobalIndex(1,ptbin)]
243 || d->GetDCA(2) > fCutsRD[GetGlobalIndex(1,ptbin)]
ba27eecd 244 || d->GetDCA(5) > fCutsRD[GetGlobalIndex(1,ptbin)]) return 0;
245 if(d->GetDist12toPrim() < fCutsRD[GetGlobalIndex(2,ptbin)]) return 0;
246 if(d->GetDist3toPrim() < fCutsRD[GetGlobalIndex(3,ptbin)]) return 0;
247 if(d->GetDist4toPrim() < fCutsRD[GetGlobalIndex(4,ptbin)]) return 0;
248 if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(5,ptbin)]) return 0;
249 if(d->Pt() < fCutsRD[GetGlobalIndex(6,ptbin)]) return 0;
250 if(!d->CutRhoMass(mD0,mD0bar,fCutsRD[GetGlobalIndex(0,ptbin)],fCutsRD[GetGlobalIndex(7,ptbin)])) return 0;
251
252 if (okD0) returnvalue=1; //cuts passed as D0
253 if (okD0bar) returnvalue=2; //cuts passed as D0bar
254 if (okD0 && okD0bar) returnvalue=3; //cuts passed as D0 and D0bar
255 }
256
257 // selection on PID (from AliAODPidHF)
258 if(selectionLevel==AliRDHFCuts::kAll ||
259 selectionLevel==AliRDHFCuts::kPID) {
601736df 260
261 Int_t selD01 = D01Selected(d,AliRDHFCuts::kCandidate);
262 Int_t selD02 = D02Selected(d,AliRDHFCuts::kCandidate);
263 Int_t selD0bar1 = D0bar1Selected(d,AliRDHFCuts::kCandidate);
264 Int_t selD0bar2 = D0bar2Selected(d,AliRDHFCuts::kCandidate);
265
266 Int_t d01PID = 0, d02PID = 0, d0bar1PID = 0, d0bar2PID = 0;
267 returnvalue = IsSelectedFromPID(d, &d01PID, &d02PID, &d0bar1PID, &d0bar2PID); //This returnvalue is dummy! Now it's modified as it must be!
268
269returnvalue = 0;
270
ba27eecd 271 if((selD01 == 1 && d01PID == 1)||(selD02 == 1 && d02PID == 1)||(selD0bar1 == 1 && d0bar1PID == 1)||(selD0bar2 == 1 && d0bar2PID == 1)) returnvalue = 1;
601736df 272 }
273
ba27eecd 274 return returnvalue;
601736df 275}
ba27eecd 276
277//---------------------------------------------------------------------------
278Int_t AliRDHFCutsD0toKpipipi::IsSelectedFromPID(AliAODRecoDecayHF4Prong *d, Int_t *hyp1, Int_t *hyp2, Int_t *hyp3, Int_t *hyp4) {
279 //
601736df 280 // Apply selection (using AliAODPidHF methods)
ba27eecd 281 // Mass hypothesis true if each particle is at least compatible with specie of hypothesis
282 //
601736df 283
284 Int_t output=0;
285
286 Int_t matchK[4], matchPi[4];
ba27eecd 287 Double_t ptlimit[2] = {0.6,0.8};
601736df 288 AliAODTrack* trk[4];
289 trk[0] = (AliAODTrack*)d->GetDaughter(0);
290 trk[1] = (AliAODTrack*)d->GetDaughter(1);
291 trk[2] = (AliAODTrack*)d->GetDaughter(2);
292 trk[3] = (AliAODTrack*)d->GetDaughter(3);
293// if(!trk[0] || !trk[1] || !trk[2] || !trk[3]) { //REMOVED (needs also AliAODEvent to be passed, here and in IsSelected method)
294// trk[0]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(0)]);
295// trk[1]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(1)]);
296// trk[2]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(2)]);
297// trk[3]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(3)]);}
298
70e398d3 299 AliAODPidHF* pidObj = new AliAODPidHF();
300 pidObj->SetAsym(kTRUE);
301 pidObj->SetPLimit(ptlimit);
302 pidObj->SetSigma(0,2.); //TPC sigma, in three pT ranges
303 pidObj->SetSigma(1,1.);
304 pidObj->SetSigma(2,0.);
305 pidObj->SetSigma(3,2.); //TOF sigma, whole pT range
306 pidObj->SetTPC(kTRUE);
307 pidObj->SetTOF(kTRUE);
308 pidObj->SetMatch(1);
309 pidObj->SetCompat(kTRUE);
601736df 310 for(Int_t ii=0; ii<4; ii++) {
70e398d3 311 pidObj->SetSigma(0,2.);
312 matchK[ii] = pidObj->MatchTPCTOF(trk[ii],3); //arguments: track, mode, particle#, compatibility allowed
313 pidObj->SetSigma(0,2.);
314 matchPi[ii] = pidObj->MatchTPCTOF(trk[ii],2);
601736df 315 }
316
317 //Rho invariant mass under various hypotheses (to be matched with PID infos in order to accet the candidate)
318 Int_t d01rho03 = 0, d01rho23 = 0, d02rho01 = 0, d02rho12 = 0, d0bar1rho12 = 0, d0bar1rho23 = 0, d0bar2rho01 = 0, d0bar2rho03 = 0;
319 if(TMath::Abs(0.775 - d->InvMassRho(0,3))<0.1) {d01rho03 = 1; d0bar2rho03 = 1;}
320 if(TMath::Abs(0.775 - d->InvMassRho(2,3))<0.1) {d01rho23 = 1; d0bar1rho23 = 1;}
321 if(TMath::Abs(0.775 - d->InvMassRho(0,1))<0.1) {d02rho01 = 1; d0bar2rho01 = 1;}
322 if(TMath::Abs(0.775 - d->InvMassRho(1,2))<0.1) {d02rho12 = 1; d0bar1rho12 = 1;}
323 Int_t d01rho = 0, d02rho = 0, d0bar1rho = 0, d0bar2rho = 0;
324 if(d01rho03==1||d01rho23==1) d01rho = 1;
325 if(d02rho01==1||d02rho12==1) d02rho = 1;
326 if(d0bar1rho12==1||d0bar1rho23==1) d0bar1rho = 1;
327 if(d0bar2rho01==1||d0bar2rho03==1) d0bar2rho = 1;
328
ba27eecd 329 //This way because there could be multiple hypotheses accepted
601736df 330 if(d01rho==1 && (matchK[1]>=0 && matchPi[0]>=0 && matchPi[2]>=0 && matchPi[3]>=0)) {*hyp1 = 1; output = 1;} //d01 hyp
331 if(d02rho==1 && (matchK[3]>=0 && matchPi[0]>=0 && matchPi[1]>=0 && matchPi[2]>=0)) {*hyp2 = 1; output = 1;} //d02 hyp
332 if(d0bar1rho==1 && (matchK[0]>=0 && matchPi[1]>=0 && matchPi[2]>=0 && matchPi[3]>=0)) {*hyp3 = 1; output = 1;} //d0bar1 hyp
333 if(d0bar2rho==1 && (matchK[2]>=0 && matchPi[0]>=0 && matchPi[1]>=0 && matchPi[3]>=0)) {*hyp4 = 1; output = 1;} //d0bar2 hyp
334
ba27eecd 335 return output;
601736df 336}
ba27eecd 337//---------------------------------------------------------------------------
338Int_t AliRDHFCutsD0toKpipipi::D01Selected(TObject* obj,Int_t selectionLevel) {
339 //
340 // Apply selection
341 //
342
343 if(!fCutsRD){
344 cout<<"Cut matrix not inizialized. Exit..."<<endl;
345 return 0;
346 }
347 //PrintAll();
348 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;
349
350 if(!d){
351 cout<<"AliAODRecoDecayHF4Prong null"<<endl;
352 return 0;
353 }
354
355 Int_t returnvalue=0;
356
357 // selection on candidate
358 if(selectionLevel==AliRDHFCuts::kAll ||
359 selectionLevel==AliRDHFCuts::kCandidate) {
360
361 Int_t ptbin=PtBin(d->Pt());
362
363 Double_t mD0[2];
364 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
365
366 d->InvMassD0(mD0);
367 if(TMath::Abs(mD0[0]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;
368 }
369
370 return returnvalue;
371}
372//---------------------------------------------------------------------------
373Int_t AliRDHFCutsD0toKpipipi::D02Selected(TObject* obj,Int_t selectionLevel) {
374 //
375 // Apply selection
376 //
377
378 if(!fCutsRD){
379 cout<<"Cut matrix not inizialized. Exit..."<<endl;
380 return 0;
381 }
382 //PrintAll();
383 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;
384
385 if(!d){
386 cout<<"AliAODRecoDecayHF4Prong null"<<endl;
387 return 0;
388 }
389
390 Int_t returnvalue=0;
391
392 // selection on candidate
393 if(selectionLevel==AliRDHFCuts::kAll ||
394 selectionLevel==AliRDHFCuts::kCandidate) {
395
396 Int_t ptbin=PtBin(d->Pt());
397
398 Double_t mD0[2];
399 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
400
401 d->InvMassD0(mD0);
601736df 402 if(TMath::Abs(mD0[1]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;
ba27eecd 403 }
404
405 return returnvalue;
406}//---------------------------------------------------------------------------
407Int_t AliRDHFCutsD0toKpipipi::D0bar1Selected(TObject* obj,Int_t selectionLevel) {
408 //
409 // Apply selection
410 //
411
412 if(!fCutsRD){
413 cout<<"Cut matrix not inizialized. Exit..."<<endl;
414 return 0;
415 }
416 //PrintAll();
417 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;
418
419 if(!d){
420 cout<<"AliAODRecoDecayHF4Prong null"<<endl;
421 return 0;
422 }
423
424 Int_t returnvalue=0;
425
426 // selection on candidate
427 if(selectionLevel==AliRDHFCuts::kAll ||
428 selectionLevel==AliRDHFCuts::kCandidate) {
429
430 Int_t ptbin=PtBin(d->Pt());
431
432 Double_t mD0bar[2];
433 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
434
435 d->InvMassD0bar(mD0bar);
436 if(TMath::Abs(mD0bar[0]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;
437 }
438
439 return returnvalue;
601736df 440}
ba27eecd 441//---------------------------------------------------------------------------
442Int_t AliRDHFCutsD0toKpipipi::D0bar2Selected(TObject* obj,Int_t selectionLevel) {
443 //
444 // Apply selection
445 //
446
447 if(!fCutsRD){
448 cout<<"Cut matrix not inizialized. Exit..."<<endl;
449 return 0;
450 }
451 //PrintAll();
452 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;
453
454 if(!d){
455 cout<<"AliAODRecoDecayHF4Prong null"<<endl;
456 return 0;
457 }
458
459 Int_t returnvalue=0;
460
461 // selection on candidate
462 if(selectionLevel==AliRDHFCuts::kAll ||
463 selectionLevel==AliRDHFCuts::kCandidate) {
464
465 Int_t ptbin=PtBin(d->Pt());
466
467 Double_t mD0bar[2];
468 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
469
470 d->InvMassD0bar(mD0bar);
471 if(TMath::Abs(mD0bar[1]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;
472 }
473
474 return returnvalue;
601736df 475}
476//---------------------------------------------------------------------------
477Bool_t AliRDHFCutsD0toKpipipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
478{
479 //
480 // Checking if D0 is in fiducial acceptance region
481 //
482
af5d687e 483 if(fMaxRapidityCand>-998.){
484 if(TMath::Abs(y) > fMaxRapidityCand) return kFALSE;
485 else return kTRUE;
486 }
487
601736df 488 if(pt > 5.) {
489 // applying cut for pt > 5 GeV
490 AliDebug(4,Form("pt of D0 = %f (> 5), cutting at |y| < 0.8\n",pt));
491 if (TMath::Abs(y) > 0.8){
492 return kFALSE;
493 }
494 } else {
495 // appliying smooth cut for pt < 5 GeV
496 Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
497 Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
498 AliDebug(4,Form("pt of D0 = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
499 if (y < minFiducialY || y > maxFiducialY){
500 return kFALSE;
501 }
502 }
503
504 return kTRUE;
505}
506