]> 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
Updates D* cuts for peripheral PbPb events
[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
a6003e0a 100void AliRDHFCutsD0toKpipipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent* aod) {\r
4755453e 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
a6003e0a 112 //recalculate vertex w/o daughters
113 Bool_t cleanvtx=kFALSE;
114 AliAODVertex *origownvtx=0x0;
115 if(fRemoveDaughtersFromPrimary) {
116 if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());
117 cleanvtx=kTRUE;
118 if(!RecalcOwnPrimaryVtx(dd,aod)) {
119 CleanOwnPrimaryVtx(dd,aod,origownvtx);
120 cleanvtx=kFALSE;
121 }
122 }
123
4755453e 124 Int_t iter=-1;\r
125\r
126 if(fVarsForOpt[0]) {\r
127 iter++;\r
128 Double_t mD0[2],mD0bar[2];\r
129 if(TMath::Abs(pdgdaughters[1])==321 || TMath::Abs(pdgdaughters[3])==321) {\r
130 dd->InvMassD0(mD0);\r
131 if(TMath::Abs(pdgdaughters[1])==321) {\r
132 vars[iter]=mD0[0];\r
133 }else{\r
134 vars[iter]=mD0[1];\r
135 }\r
136 } else {\r
137 dd->InvMassD0bar(mD0bar);\r
138 if(TMath::Abs(pdgdaughters[0])==321) {\r
e11ae259 139 vars[iter]=mD0bar[0];\r
4755453e 140 }else{\r
e11ae259 141 vars[iter]=mD0bar[1];\r
4755453e 142 }\r
143 }\r
144 }\r
145\r
146 if(fVarsForOpt[1]){\r
147 iter++;\r
148 vars[iter]=dd->GetDCA();\r
149 }\r
150\r
151 if(fVarsForOpt[2]){\r
152 iter++;\r
153 vars[iter]=dd->GetDist12toPrim();\r
154 }\r
155 if(fVarsForOpt[3]){\r
156 iter++;\r
157 vars[iter]=dd->GetDist3toPrim();\r
158 }\r
159 if(fVarsForOpt[4]){\r
160 iter++;\r
161 vars[iter]=dd->GetDist4toPrim();\r
162 }\r
163 if(fVarsForOpt[5]){\r
164 iter++;\r
165 vars[iter]=dd->CosPointingAngle();\r
166 }\r
167 if(fVarsForOpt[6]){\r
168 iter++;\r
169 vars[iter]=dd->Pt();\r
170 }\r
171 if(fVarsForOpt[7]){\r
172 iter++;\r
173 vars[iter]=999999999.;\r
174 printf("ERROR: optmization for rho mass cut not implemented\n");\r
175 }\r
176 if(fVarsForOpt[8]){\r
177 iter++;\r
178 vars[iter]=999999999.;\r
179 printf("ERROR: optmization for PID cut not implemented\n");\r
180 }\r
181 \r
a6003e0a 182 if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
4755453e 183 return;\r
184}\r
185//---------------------------------------------------------------------------\r
186Int_t AliRDHFCutsD0toKpipipi::IsSelected(TObject* obj,Int_t selectionLevel) {\r
187 //\r
188 // Apply selection\r
189 //\r
190\r
191 if(!fCutsRD){\r
192 cout<<"Cut matrix not inizialized. Exit..."<<endl;\r
193 return 0;\r
194 }\r
195 //PrintAll();\r
196 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;\r
197\r
198 if(!d){\r
199 cout<<"AliAODRecoDecayHF4Prong null"<<endl;\r
200 return 0;\r
201 }\r
202\r
6af2d88a 203 Double_t ptD=d->Pt();\r
204 if(ptD<fMinPtCand) return 0;\r
601736df 205 if(ptD>fMaxPtCand) return 0;
4755453e 206\r
207 // selection on daughter tracks \r
208 if(selectionLevel==AliRDHFCuts::kAll || \r
209 selectionLevel==AliRDHFCuts::kTracks) {\r
210 if(!AreDaughtersSelected(d)) return 0;\r
211 }\r
212\r
213\r
214 Int_t returnvalue=1;\r
215\r
216 // selection on candidate\r
217 if(selectionLevel==AliRDHFCuts::kAll || \r
218 selectionLevel==AliRDHFCuts::kCandidate) {\r
219\r
220 Int_t ptbin=PtBin(d->Pt());\r
601736df 221
222 Int_t okD0=1,okD0bar=1; \r
4755453e 223 Double_t mD0[2],mD0bar[2];\r
224 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
225\r
226 d->InvMassD0(mD0);\r
227 if(TMath::Abs(mD0[0]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)] &&\r
228 TMath::Abs(mD0[1]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;\r
229 d->InvMassD0bar(mD0bar);\r
230 if(TMath::Abs(mD0bar[0]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)] &&\r
231 TMath::Abs(mD0bar[1]-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;\r
232 if(!okD0 && !okD0bar) return 0;\r
233 \r
601736df 234 if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]
235 || d->GetDCA(3) > fCutsRD[GetGlobalIndex(1,ptbin)]
236 || d->GetDCA(2) > fCutsRD[GetGlobalIndex(1,ptbin)]
237 || d->GetDCA(5) > fCutsRD[GetGlobalIndex(1,ptbin)]) return 0;\r
4755453e 238 if(d->GetDist12toPrim() < fCutsRD[GetGlobalIndex(2,ptbin)]) return 0;\r
239 if(d->GetDist3toPrim() < fCutsRD[GetGlobalIndex(3,ptbin)]) return 0;\r
240 if(d->GetDist4toPrim() < fCutsRD[GetGlobalIndex(4,ptbin)]) return 0;\r
241 if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(5,ptbin)]) return 0;\r
242 if(d->Pt() < fCutsRD[GetGlobalIndex(6,ptbin)]) return 0;\r
243 if(!d->CutRhoMass(mD0,mD0bar,fCutsRD[GetGlobalIndex(0,ptbin)],fCutsRD[GetGlobalIndex(7,ptbin)])) return 0;\r
601736df 244
4755453e 245 if (okD0) returnvalue=1; //cuts passed as D0\r
246 if (okD0bar) returnvalue=2; //cuts passed as D0bar\r
247 if (okD0 && okD0bar) returnvalue=3; //cuts passed as D0 and D0bar\r
248 }\r
601736df 249
250 // selection on PID (from AliAODPidHF)\r
251 if(selectionLevel==AliRDHFCuts::kAll || \r
252 selectionLevel==AliRDHFCuts::kPID) {\r
253
254 Int_t selD01 = D01Selected(d,AliRDHFCuts::kCandidate);
255 Int_t selD02 = D02Selected(d,AliRDHFCuts::kCandidate);
256 Int_t selD0bar1 = D0bar1Selected(d,AliRDHFCuts::kCandidate);
257 Int_t selD0bar2 = D0bar2Selected(d,AliRDHFCuts::kCandidate);
258
259 Int_t d01PID = 0, d02PID = 0, d0bar1PID = 0, d0bar2PID = 0;
260 returnvalue = IsSelectedFromPID(d, &d01PID, &d02PID, &d0bar1PID, &d0bar2PID); //This returnvalue is dummy! Now it's modified as it must be!
261
262returnvalue = 0;
263
264 if((selD01 == 1 && d01PID == 1)||(selD02 == 1 && d02PID == 1)||(selD0bar1 == 1 && d0bar1PID == 1)||(selD0bar2 == 1 && d0bar2PID == 1)) returnvalue = 1;\r
265 }
266
4755453e 267 return returnvalue;\r
601736df 268}
269\r
270//---------------------------------------------------------------------------\r
271Int_t AliRDHFCutsD0toKpipipi::IsSelectedFromPID(AliAODRecoDecayHF4Prong *d, Int_t *hyp1, Int_t *hyp2, Int_t *hyp3, Int_t *hyp4) {\r
272 //\r
273 // Apply selection (using AliAODPidHF methods)
274 // Mass hypothesis true if each particle is at least compatible with specie of hypothesis\r
275 // \r
276
277 Int_t output=0;
278
279 Int_t matchK[4], matchPi[4];
280 Double_t ptlimit[2] = {0.6,0.8};\r
281 AliAODTrack* trk[4];
282 trk[0] = (AliAODTrack*)d->GetDaughter(0);
283 trk[1] = (AliAODTrack*)d->GetDaughter(1);
284 trk[2] = (AliAODTrack*)d->GetDaughter(2);
285 trk[3] = (AliAODTrack*)d->GetDaughter(3);
286// if(!trk[0] || !trk[1] || !trk[2] || !trk[3]) { //REMOVED (needs also AliAODEvent to be passed, here and in IsSelected method)
287// trk[0]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(0)]);
288// trk[1]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(1)]);
289// trk[2]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(2)]);
290// trk[3]=aodIn->GetTrack(trkIDtoEntry[d->GetProngID(3)]);}
291
292 AliAODPidHF* PidObj = new AliAODPidHF();
293 PidObj->SetAsym(kTRUE);
294 PidObj->SetPLimit(ptlimit);
295 PidObj->SetSigma(0,2.); //TPC sigma, in three pT ranges
296 PidObj->SetSigma(1,1.);
297 PidObj->SetSigma(2,0.);
298 PidObj->SetSigma(3,2.); //TOF sigma, whole pT range
299 PidObj->SetTPC(kTRUE);
300 PidObj->SetTOF(kTRUE);
301
302 for(Int_t ii=0; ii<4; ii++) {
303 PidObj->SetSigma(0,2.);
304 matchK[ii] = PidObj->MatchTPCTOF(trk[ii],1,3,kTRUE); //arguments: track, mode, particle#, compatibility allowed
305 PidObj->SetSigma(0,2.);
306 matchPi[ii] = PidObj->MatchTPCTOF(trk[ii],1,2,kTRUE);
307 }
308
309 //Rho invariant mass under various hypotheses (to be matched with PID infos in order to accet the candidate)
310 Int_t d01rho03 = 0, d01rho23 = 0, d02rho01 = 0, d02rho12 = 0, d0bar1rho12 = 0, d0bar1rho23 = 0, d0bar2rho01 = 0, d0bar2rho03 = 0;
311 if(TMath::Abs(0.775 - d->InvMassRho(0,3))<0.1) {d01rho03 = 1; d0bar2rho03 = 1;}
312 if(TMath::Abs(0.775 - d->InvMassRho(2,3))<0.1) {d01rho23 = 1; d0bar1rho23 = 1;}
313 if(TMath::Abs(0.775 - d->InvMassRho(0,1))<0.1) {d02rho01 = 1; d0bar2rho01 = 1;}
314 if(TMath::Abs(0.775 - d->InvMassRho(1,2))<0.1) {d02rho12 = 1; d0bar1rho12 = 1;}
315 Int_t d01rho = 0, d02rho = 0, d0bar1rho = 0, d0bar2rho = 0;
316 if(d01rho03==1||d01rho23==1) d01rho = 1;
317 if(d02rho01==1||d02rho12==1) d02rho = 1;
318 if(d0bar1rho12==1||d0bar1rho23==1) d0bar1rho = 1;
319 if(d0bar2rho01==1||d0bar2rho03==1) d0bar2rho = 1;
320
321 //This way because there could be multiple hypotheses accepted\r
322 if(d01rho==1 && (matchK[1]>=0 && matchPi[0]>=0 && matchPi[2]>=0 && matchPi[3]>=0)) {*hyp1 = 1; output = 1;} //d01 hyp
323 if(d02rho==1 && (matchK[3]>=0 && matchPi[0]>=0 && matchPi[1]>=0 && matchPi[2]>=0)) {*hyp2 = 1; output = 1;} //d02 hyp
324 if(d0bar1rho==1 && (matchK[0]>=0 && matchPi[1]>=0 && matchPi[2]>=0 && matchPi[3]>=0)) {*hyp3 = 1; output = 1;} //d0bar1 hyp
325 if(d0bar2rho==1 && (matchK[2]>=0 && matchPi[0]>=0 && matchPi[1]>=0 && matchPi[3]>=0)) {*hyp4 = 1; output = 1;} //d0bar2 hyp
326
327 return output;\r
4755453e 328}\r
329//---------------------------------------------------------------------------\r
601736df 330Int_t AliRDHFCutsD0toKpipipi::D01Selected(TObject* obj,Int_t selectionLevel) {\r
331 //\r
332 // Apply selection\r
333 //\r
334\r
335 if(!fCutsRD){\r
336 cout<<"Cut matrix not inizialized. Exit..."<<endl;\r
337 return 0;\r
338 }\r
339 //PrintAll();\r
340 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;\r
341\r
342 if(!d){\r
343 cout<<"AliAODRecoDecayHF4Prong null"<<endl;\r
344 return 0;\r
345 }\r
346\r
347 Int_t returnvalue=0;\r
348\r
349 // selection on candidate\r
350 if(selectionLevel==AliRDHFCuts::kAll || \r
351 selectionLevel==AliRDHFCuts::kCandidate) {\r
352\r
353 Int_t ptbin=PtBin(d->Pt());\r
354 \r
355 Double_t mD0[2];\r
356 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
357\r
358 d->InvMassD0(mD0);\r
359 if(TMath::Abs(mD0[0]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;\r
360 }\r
361\r
362 return returnvalue;\r
363}
364//---------------------------------------------------------------------------\r
365Int_t AliRDHFCutsD0toKpipipi::D02Selected(TObject* obj,Int_t selectionLevel) {\r
366 //\r
367 // Apply selection\r
368 //\r
369\r
370 if(!fCutsRD){\r
371 cout<<"Cut matrix not inizialized. Exit..."<<endl;\r
372 return 0;\r
373 }\r
374 //PrintAll();\r
375 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;\r
376\r
377 if(!d){\r
378 cout<<"AliAODRecoDecayHF4Prong null"<<endl;\r
379 return 0;\r
380 }\r
381\r
382 Int_t returnvalue=0;\r
383\r
384 // selection on candidate\r
385 if(selectionLevel==AliRDHFCuts::kAll || \r
386 selectionLevel==AliRDHFCuts::kCandidate) {\r
387\r
388 Int_t ptbin=PtBin(d->Pt());\r
389 \r
390 Double_t mD0[2];\r
391 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
392\r
393 d->InvMassD0(mD0);\r
394 if(TMath::Abs(mD0[1]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;
395 }\r
396\r
397 return returnvalue;\r
398}\r//---------------------------------------------------------------------------\r
399Int_t AliRDHFCutsD0toKpipipi::D0bar1Selected(TObject* obj,Int_t selectionLevel) {\r
400 //\r
401 // Apply selection\r
402 //\r
403\r
404 if(!fCutsRD){\r
405 cout<<"Cut matrix not inizialized. Exit..."<<endl;\r
406 return 0;\r
407 }\r
408 //PrintAll();\r
409 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;\r
410\r
411 if(!d){\r
412 cout<<"AliAODRecoDecayHF4Prong null"<<endl;\r
413 return 0;\r
414 }\r
415\r
416 Int_t returnvalue=0;\r
417\r
418 // selection on candidate\r
419 if(selectionLevel==AliRDHFCuts::kAll || \r
420 selectionLevel==AliRDHFCuts::kCandidate) {\r
421\r
422 Int_t ptbin=PtBin(d->Pt());\r
423 \r
424 Double_t mD0bar[2];\r
425 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
426
427 d->InvMassD0bar(mD0bar);\r
428 if(TMath::Abs(mD0bar[0]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;\r
429 }\r
430\r
431 return returnvalue;\r
432}
433//---------------------------------------------------------------------------\r
434Int_t AliRDHFCutsD0toKpipipi::D0bar2Selected(TObject* obj,Int_t selectionLevel) {\r
435 //\r
436 // Apply selection\r
437 //\r
438\r
439 if(!fCutsRD){\r
440 cout<<"Cut matrix not inizialized. Exit..."<<endl;\r
441 return 0;\r
442 }\r
443 //PrintAll();\r
444 AliAODRecoDecayHF4Prong* d=(AliAODRecoDecayHF4Prong*)obj;\r
445\r
446 if(!d){\r
447 cout<<"AliAODRecoDecayHF4Prong null"<<endl;\r
448 return 0;\r
449 }\r
450\r
451 Int_t returnvalue=0;\r
452\r
453 // selection on candidate\r
454 if(selectionLevel==AliRDHFCuts::kAll || \r
455 selectionLevel==AliRDHFCuts::kCandidate) {\r
456\r
457 Int_t ptbin=PtBin(d->Pt());\r
458 \r
459 Double_t mD0bar[2];\r
460 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
461
462 d->InvMassD0bar(mD0bar);\r
463 if(TMath::Abs(mD0bar[1]-mD0PDG) < fCutsRD[GetGlobalIndex(0,ptbin)]) returnvalue = 1;\r
464 }\r
465\r
466 return returnvalue;\r
467}
468//---------------------------------------------------------------------------
469Bool_t AliRDHFCutsD0toKpipipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
470{
471 //
472 // Checking if D0 is in fiducial acceptance region
473 //
474
475 if(pt > 5.) {
476 // applying cut for pt > 5 GeV
477 AliDebug(4,Form("pt of D0 = %f (> 5), cutting at |y| < 0.8\n",pt));
478 if (TMath::Abs(y) > 0.8){
479 return kFALSE;
480 }
481 } else {
482 // appliying smooth cut for pt < 5 GeV
483 Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5;
484 Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;
485 AliDebug(4,Form("pt of D0 = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
486 if (y < minFiducialY || y > maxFiducialY){
487 return kFALSE;
488 }
489 }
490
491 return kTRUE;
492}
493