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