Added IsInFiducualAcceptance() (Alessandro)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliRDHFCutsDStartoKpipi.cxx
CommitLineData
45bc9115 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 DStar->Kpipi
19//
20// Author: A.Grelli, alessandro.grelli@uu.nl
21/////////////////////////////////////////////////////////////
22
23#include <TDatabasePDG.h>
24#include <Riostream.h>
25
26#include "AliAODRecoDecayHF2Prong.h"
27#include "AliAODRecoCascadeHF.h"
28#include "AliRDHFCutsD0toKpi.h"
29#include "AliRDHFCutsDStartoKpipi.h"
30#include "AliAODTrack.h"
31#include "AliESDtrack.h"
da6fefc3 32#include "AliAODVertex.h"
33#include "AliESDVertex.h"
45bc9115 34
ee1128c6 35
36
45bc9115 37ClassImp(AliRDHFCutsDStartoKpipi)
38
39//--------------------------------------------------------------------------
da6fefc3 40AliRDHFCutsDStartoKpipi::AliRDHFCutsDStartoKpipi(const char* name) :
41AliRDHFCuts(name),
42fTrackCutsSoftPi(0)
45bc9115 43{
44 //
45 // Default Constructor
46 //
47
48 Int_t nvars=14;
49 SetNVars(nvars);
50 TString varNames[14]={
51 "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 "inv. mass half width of D* [GeV]",
61 "half width of (M_Kpipi-M_D0) [GeV]",
62 "PtMin of pi_s [GeV/c]",
63 "PtMax of pi_s [GeV/c]",
64 "theta, angle between the pi_s and decay plane of the D0 [rad]"};
65 Bool_t isUpperCut[14]={
66 kTRUE,
67 kTRUE,
68 kTRUE,
69 kFALSE,
70 kFALSE,
71 kTRUE,
72 kTRUE,
73 kTRUE,
74 kFALSE,
75 kTRUE,
76 kTRUE,
77 kTRUE,
78 kTRUE,
79 kFALSE};
80 SetVarNames(nvars,varNames,isUpperCut);
81 Bool_t forOpt[14]={
82 kFALSE,
83 kTRUE,
84 kTRUE,
85 kFALSE,
86 kFALSE,
87 kFALSE,
88 kFALSE,
89 kTRUE,
90 kTRUE,
91 kFALSE,
92 kTRUE,
93 kFALSE,
94 kFALSE,
95 kFALSE};
96 SetVarsForOpt(5,forOpt);
97 Float_t limits[2]={0,999999999.};
98 SetPtBins(2,limits);
99}
100//--------------------------------------------------------------------------
101AliRDHFCutsDStartoKpipi::AliRDHFCutsDStartoKpipi(const AliRDHFCutsDStartoKpipi &source) :
da6fefc3 102 AliRDHFCuts(source),
103 fTrackCutsSoftPi(0)
45bc9115 104{
105 //
106 // Copy constructor
107 //
ee1128c6 108
da6fefc3 109 if(source.GetTrackCutsSoftPi()) AddTrackCutsSoftPi(source.GetTrackCutsSoftPi());
ee1128c6 110
45bc9115 111}
112//--------------------------------------------------------------------------
113AliRDHFCutsDStartoKpipi &AliRDHFCutsDStartoKpipi::operator=(const AliRDHFCutsDStartoKpipi &source)
114{
115 //
116 // assignment operator
117 //
118 if(&source == this) return *this;
119
120 AliRDHFCuts::operator=(source);
121
da6fefc3 122 if(source.GetTrackCutsSoftPi()) AddTrackCutsSoftPi(source.GetTrackCutsSoftPi());
45bc9115 123 return *this;
124}
125
126
127//---------------------------------------------------------------------------
128void AliRDHFCutsDStartoKpipi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
129 //
130 // Fills in vars the values of the variables
131 //
132 if(nvars!=fnVarsForOpt) {
133 printf("AliRDHFCutsDStartoKpipi::GetCutsVarsForOpt: wrong number of variables\n");
134 return;
135 }
136
137
138 AliAODRecoCascadeHF* dstarD0pi = (AliAODRecoCascadeHF*)d;
139
140 AliAODTrack *softPi = (AliAODTrack*)dstarD0pi->GetBachelor();
141
142 AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();
143
144 Int_t iter=-1;
145 if(fVarsForOpt[0]){
146 iter++;
147 if(TMath::Abs(pdgdaughters[0])==211) {
148 vars[iter]=dd->InvMassD0();
149 } else {
150 vars[iter]=dd->InvMassD0bar();
151 }
152 }
153 if(fVarsForOpt[1]){
154 iter++;
155 vars[iter]=dd->GetDCA();
156 }
157 if(fVarsForOpt[2]){
158 iter++;
159 if(TMath::Abs(pdgdaughters[0])==211) {
160 vars[iter] = dd->CosThetaStarD0();
161 } else {
162 vars[iter] = dd->CosThetaStarD0bar();
163 }
164 }
165 if(fVarsForOpt[3]){
166 iter++;
167 if(TMath::Abs(pdgdaughters[0])==321) {
168 vars[iter]=dd->PtProng(0);
169 }
170 else{
171 vars[iter]=dd->PtProng(1);
172 }
173 }
174 if(fVarsForOpt[4]){
175 iter++;
176 if(TMath::Abs(pdgdaughters[0])==211) {
177 vars[iter]=dd->PtProng(0);
178 }
179 else{
180 vars[iter]=dd->PtProng(1);
181 }
182 }
183 if(fVarsForOpt[5]){
184 iter++;
185 if(TMath::Abs(pdgdaughters[0])==321) {
186 vars[iter]=dd->Getd0Prong(0);
187 }
188 else{
189 vars[iter]=dd->Getd0Prong(1);
190 }
191 }
192 if(fVarsForOpt[6]){
193 iter++;
194 if(TMath::Abs(pdgdaughters[0])==211) {
195 vars[iter]=dd->Getd0Prong(0);
196 }
197 else{
198 vars[iter]=dd->Getd0Prong(1);
199 }
200 }
201 if(fVarsForOpt[7]){
202 iter++;
203 vars[iter]= dd->Prodd0d0();
204 }
205 if(fVarsForOpt[8]){
206 iter++;
207 vars[iter]=dd->CosPointingAngle();
208 }
209 if(fVarsForOpt[9]){
210 iter++;
211 vars[iter]=dstarD0pi->InvMassDstarKpipi();
212 }
213 if(fVarsForOpt[10]){
214 iter++;
215 vars[iter]=dstarD0pi->DeltaInvMass();
216 }
217 if(fVarsForOpt[11]){
218 iter++;
219 vars[iter] = softPi->Pt();
220 }
221 if(fVarsForOpt[12]){
222 iter++;
223 vars[iter] = softPi->Pt();
224 }
225 if(fVarsForOpt[13]){
226 iter++;
227 vars[iter] =dstarD0pi->AngleD0dkpPisoft();
228 }
229
230 return;
231}
232//---------------------------------------------------------------------------
233Int_t AliRDHFCutsDStartoKpipi::IsSelected(TObject* obj,Int_t selectionLevel) {
234 //
235 // Apply selection for D*.
236 //
237 if(!fCutsRD){
238 cout<<"Cut matrice not inizialized. Exit..."<<endl;
239 return 0;
240 }
241
242 AliAODRecoCascadeHF* d = (AliAODRecoCascadeHF*)obj;
243 if(!d){
244 cout<<"AliAODRecoCascadeHF null"<<endl;
245 return 0;
246 }
247
248 AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)d->Get2Prong();
249 if(!dd){
250 cout<<"AliAODRecoDecayHF2Prong null"<<endl;
251 return 0;
252 }
253
254 AliAODTrack *b = (AliAODTrack*)d->GetBachelor();
255
256 // selection on daughter tracks
257 if(selectionLevel==AliRDHFCuts::kAll ||
258 selectionLevel==AliRDHFCuts::kTracks) {
da6fefc3 259 if(!AreDaughtersSelected(dd)) return 0;
260 if(fTrackCutsSoftPi) {
261 AliAODVertex *vAOD = d->GetPrimaryVtx();
262 Double_t pos[3],cov[6];
263 vAOD->GetXYZ(pos);
264 vAOD->GetCovarianceMatrix(cov);
265 const AliESDVertex vESD(pos,cov,100.,100);
266 if(!IsDaughterSelected(b,&vESD,fTrackCutsSoftPi)) return 0;
267 }
45bc9115 268 }
269
e4f1bb0a 270 Int_t returnvalue=1;
45bc9115 271
272 // selection on candidate
273 if(selectionLevel==AliRDHFCuts::kAll ||
274 selectionLevel==AliRDHFCuts::kCandidate) {
275
276 Double_t pt=d->Pt();
277 Int_t ptbin=PtBin(pt);
278
279 // select D0 that passes D* cuts
280 returnvalue = IsD0FromDStarSelected(pt,dd,selectionLevel);
ee1128c6 281 //if(retunvalue==0) return 0;
282
283 if((b->Charge()==+1 && returnvalue==2) || (b->Charge()==-1 && returnvalue==1)) return 0;
45bc9115 284
285 // DStarMass and D0mass
286 Double_t mDSPDG = TDatabasePDG::Instance()->GetParticle(413)->Mass();
287 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
288 // delta mass PDG
289 Double_t deltaPDG = mDSPDG-mD0PDG;
45bc9115 290
45bc9115 291 // Half width DStar mass
ee1128c6 292 if(TMath::Abs(mDSPDG - (d->InvMassDstarKpipi()))>fCutsRD[GetGlobalIndex(9,ptbin)]) return 0;
45bc9115 293 // Half width Delta mass
294
ee1128c6 295 if(TMath::Abs(deltaPDG-(d->DeltaInvMass())) > fCutsRD[GetGlobalIndex(10,ptbin)]) return 0;
45bc9115 296
297 // cut on soft pion pt
298 if(b->Pt() < fCutsRD[GetGlobalIndex(11,ptbin)] || b->Pt() > fCutsRD[GetGlobalIndex(12,ptbin)]) return 0;
299 // cut on the angle between D0 decay plane and soft pion
300 if(d->AngleD0dkpPisoft() > fCutsRD[GetGlobalIndex(13,ptbin)]) return 0;
301
302 }
303
304 return returnvalue;
305}
306//_________________________________________________________________________________________________
da6fefc3 307Int_t AliRDHFCutsDStartoKpipi::IsD0FromDStarSelected(Double_t pt, TObject* obj,Int_t selectionLevel) const {
45bc9115 308 //
309 // Apply selection for D0 from D*. The selection in on D0 prongs
310 //
311
312 if(!fCutsRD){
313 cout<<"Cut matrice not inizialized. Exit..."<<endl;
314 return 0;
315 }
316
317 AliAODRecoDecayHF2Prong* dd = (AliAODRecoDecayHF2Prong*)obj;
318
319 if(!dd){
320 cout<<"AliAODRecoDecayHF2Prong null"<<endl;
321 return 0;
322 }
323
da6fefc3 324 // selection on daughter tracks is done in IsSelected()
45bc9115 325
326 Int_t returnvalue=1;
327
328 // selection on candidate
329 if(selectionLevel==AliRDHFCuts::kAll ||
330 selectionLevel==AliRDHFCuts::kCandidate) {
331
332 // D0 mass
333 Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
334 // delta mass PDG
335
336 Int_t ptbin=PtBin(pt);
337
338 Double_t mD0,mD0bar,ctsD0,ctsD0bar;
339
340 Int_t okD0 =0;
341 Int_t okD0bar =0;
342 okD0=1; okD0bar=1;
343
344 if(dd->PtProng(1) < fCutsRD[GetGlobalIndex(3,ptbin)] || dd->PtProng(0) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
345 if(dd->PtProng(0) < fCutsRD[GetGlobalIndex(3,ptbin)] || dd->PtProng(1) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;
346
347 if(!okD0 && !okD0bar) return 0;
348
349 if(TMath::Abs(dd->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] ||
350 TMath::Abs(dd->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
351 if(TMath::Abs(dd->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||
352 TMath::Abs(dd->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;
353 if(!okD0 && !okD0bar) return 0;
354
355 if(dd->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) return 0;
356
357 dd->InvMassD0(mD0,mD0bar);
358 if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
359 if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0bar = 0;
360 if(!okD0 && !okD0bar) return 0;
361
362 dd->CosThetaStarD0(ctsD0,ctsD0bar);
363 if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0;
364 if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;
365 if(!okD0 && !okD0bar) return 0;
366
367 if(dd->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)]) return 0;
368
369 if(dd->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)]) return 0;
370
371 if (okD0) returnvalue=1; //cuts passed as D0
372 if (okD0bar) returnvalue=2; //cuts passed as D0bar
373 if (okD0 && okD0bar) returnvalue=3; //both
374 }
375
ee1128c6 376 //if((Charge()==+1 && !okD0) || (Charge()==-1 && !okD0bar)) return kFALSE;
45bc9115 377 return returnvalue;
378}
96462a8c 379//----------------------------------------------------------------------------------
380Bool_t AliRDHFCutsDStartoKpipi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
381{
382 //
383 // D* fiducial acceptance region
384 //
45bc9115 385
96462a8c 386 if(pt > 5.) {
387 // applying cut for pt > 5 GeV
388 AliDebug(4,Form("pt of D* = %f (> 5), cutting at |y| < 0.8\n",pt));
389 if (TMath::Abs(y) > 0.8){
390 return kFALSE;
391 }
392 } else {
393 // appliying smooth cut for pt < 5 GeV
394 Double_t maxFiducialY = -0.151/15*pt*pt+1.9/15*pt+0.4;
395 Double_t minFiducialY = 0.151/15*pt*pt-1.9/15*pt-0.4;
396 AliDebug(4,Form("pt of D* = %f (< 5), cutting according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY));
397 if (y < minFiducialY || y > maxFiducialY){
398 return kFALSE;
399 }
400 }
45bc9115 401
96462a8c 402 return kTRUE;
403}