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