]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliAODRecoCascadeHF.cxx
Addition of Lc->V0+bachelor analysis (Zaida)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAODRecoCascadeHF.cxx
CommitLineData
05d80dd6 1/**************************************************************************
2 * Copyright(c) 1998-2008, 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 AOD reconstructed heavy-flavour cascades
19//
c3fe1eed 20// Author: X-M. Zhang, zhangxm@iopp.ccnu.edu.cn
05d80dd6 21/////////////////////////////////////////////////////////////
22
23#include <TVector3.h>
24#include <TDatabasePDG.h>
b3999999 25#include <TClonesArray.h>
26#include "AliAODMCParticle.h"
05d80dd6 27#include "AliAODRecoDecay.h"
28#include "AliAODVertex.h"
29#include "AliAODRecoDecayHF2Prong.h"
30#include "AliAODRecoCascadeHF.h"
31
32ClassImp(AliAODRecoCascadeHF)
33//-----------------------------------------------------------------------------
34
35AliAODRecoCascadeHF::AliAODRecoCascadeHF() :
dcfa35b3 36 AliAODRecoDecayHF2Prong()
05d80dd6 37{
38 //
39 // Default Constructor
40 //
41}
42//-----------------------------------------------------------------------------
43AliAODRecoCascadeHF::AliAODRecoCascadeHF(AliAODVertex *vtx2, Short_t charge,
44 Double_t *px, Double_t *py, Double_t *pz,
45 Double_t *d0, Double_t *d0err, Double_t dca) :
dcfa35b3 46 AliAODRecoDecayHF2Prong(vtx2, px, py, pz, d0, d0err, dca)
05d80dd6 47{
48 //
49 // Constructor with AliAODVertex for decay vertex
50 //
51 SetCharge(charge);
52}
53//-----------------------------------------------------------------------------
54AliAODRecoCascadeHF::AliAODRecoCascadeHF(AliAODVertex *vtx2, Short_t charge,
55 Double_t *d0, Double_t *d0err, Double_t dca) :
dcfa35b3 56 AliAODRecoDecayHF2Prong(vtx2, d0, d0err, dca)
05d80dd6 57{
58 //
59 // Constructor with decay vertex and without prongs momenta
60 //
61 SetCharge(charge);
62}
63//-----------------------------------------------------------------------------
64AliAODRecoCascadeHF::AliAODRecoCascadeHF(const AliAODRecoCascadeHF &source) :
dcfa35b3 65 AliAODRecoDecayHF2Prong(source)
05d80dd6 66{
67 //
68 // Copy constructor
69 //
70}
71//-----------------------------------------------------------------------------
72AliAODRecoCascadeHF &AliAODRecoCascadeHF::operator=(const AliAODRecoCascadeHF &source)
73{
74 //
75 // assignment operator
76 //
77 if(&source == this) return *this;
78
79 AliAODRecoDecayHF2Prong::operator=(source);
80
05d80dd6 81 return *this;
82}
83//-----------------------------------------------------------------------------
84AliAODRecoCascadeHF::~AliAODRecoCascadeHF()
85{
86 //
87 // Default Destructor
88 //
89}
90//-----------------------------------------------------------------------------
91Double_t AliAODRecoCascadeHF::InvMassDstarKpipi() const
92{
93 //
94 // 3 prong invariant mass of the D0 daughters and the soft pion
95 //
96
97 Double_t px[3],py[3],pz[3];
98 UInt_t pdg[3]={321,211,211};
99 pdg[0] = (Charge()>0 ? 211 : 321); // positive daughter of D0
100 px[0] = Get2Prong()->PxProng(0);
101 py[0] = Get2Prong()->PyProng(0);
102 pz[0] = Get2Prong()->PzProng(0);
103 pdg[1] = (Charge()>0 ? 321 : 211); // negative daughter of D0
104 px[1] = Get2Prong()->PxProng(1);
105 py[1] = Get2Prong()->PyProng(1);
106 pz[1] = Get2Prong()->PzProng(1);
107 pdg[2] = 211; // soft pion
108 px[2] = PxProng(0);
109 py[2] = PyProng(0);
110 pz[2] = PzProng(0);
111 Short_t dummycharge=0;
112 Double_t dummyd0[3]={0,0,0};
113 AliAODRecoDecay *rd = new AliAODRecoDecay(0x0,3,dummycharge,px,py,pz,dummyd0);
114
115 Double_t minv = rd->InvMass(3,pdg);
116
117 delete rd; rd=NULL;
118
119 return minv;
120}
b3999999 121//----------------------------------------------------------------------------
122Int_t AliAODRecoCascadeHF::MatchToMC(Int_t pdgabs,Int_t pdgabs2prong,
383937f0 123 Int_t *pdgDg,Int_t *pdgDg2prong,
b3999999 124 TClonesArray *mcArray) const
125{
126 //
127 // Check if this candidate is matched to a MC signal
128 // If no, return -1
129 // If yes, return label (>=0) of the AliAODMCParticle
130 //
131
8862e6bf 132 Int_t ndg=GetNDaughters();
133 if(!ndg) {
b3999999 134 AliError("No daughters available");
135 return -1;
136 }
137
138 AliAODRecoDecayHF2Prong *the2Prong = Get2Prong();
383937f0 139 Int_t lab2Prong = the2Prong->MatchToMC(pdgabs2prong,mcArray,2,pdgDg2prong);
b3999999 140
141 if(lab2Prong<0) return -1;
142
8862e6bf 143 Int_t dgLabels[10];
b3999999 144
145 // loop on daughters and write labels
8862e6bf 146 for(Int_t i=0; i<ndg; i++) {
b3999999 147 AliVTrack *trk = (AliVTrack*)GetDaughter(i);
148 Int_t lab = trk->GetLabel();
149 if(lab==-1) { // this daughter is the 2prong
150 lab=lab2Prong;
151 } else if(lab<-1) {
152 printf("daughter with negative label\n");
153 continue;
154 }
155 dgLabels[i] = lab;
156 }
157
383937f0 158 return AliAODRecoDecay::MatchToMC(pdgabs,mcArray,dgLabels,2,2,pdgDg);
b3999999 159}
05d80dd6 160//-----------------------------------------------------------------------------
161Bool_t AliAODRecoCascadeHF::SelectDstar(const Double_t *cutsDstar,
162 const Double_t *cutsD0,
163 Bool_t testD0) const
164{
165 //
166 // cutsDstar[0] = inv. mass half width of D* [GeV]
167 // cutsDstar[1] = half width of (M_Kpipi-M_D0) [GeV]
168 // cutsDstar[2] = PtMin of pi_s [GeV/c]
169 // cutsDstar[3] = PtMax of pi_s [GeV/c]
170 // cutsDstar[4] = theta, angle between the pi_s and decay plane of the D0 [rad]
171 //
172 // cutsD0[0] = inv. mass half width [GeV]
173 // cutsD0[1] = dca [cm]
174 // cutsD0[2] = cosThetaStar
175 // cutsD0[3] = pTK [GeV/c]
176 // cutsD0[4] = pTPi [GeV/c]
177 // cutsD0[5] = d0K [cm] upper limit!
178 // cutsD0[6] = d0Pi [cm] upper limit!
179 // cutsD0[7] = d0d0 [cm^2]
180 // cutsD0[8] = cosThetaPoint
181
182
183 // check that the D0 passes the cuts
184 // (if we have a D*+, it has to pass as D0,
185 // if we have a D*-, it has to pass as D0bar)
186
187 if(testD0) {
188 Int_t okD0=0,okD0bar=0;
189 Get2Prong()->SelectD0(cutsD0,okD0,okD0bar);
190 if((Charge()==+1 && !okD0) || (Charge()==-1 && !okD0bar)) return kFALSE;
191 }
192
193 if( (PtProng(0)<cutsDstar[2]) || (PtProng(0)>cutsDstar[3]) ) return kFALSE;
194
195 Double_t mDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();
196 Double_t invmDstar = InvMassDstarKpipi();
197 if(TMath::Abs(mDstar-invmDstar)>cutsDstar[0]) return kFALSE;
198
199 Double_t mD0 = TDatabasePDG::Instance()->GetParticle(421)->Mass();
200 if(TMath::Abs((mDstar-mD0)-DeltaInvMass())>cutsDstar[1]) return kFALSE;
201
202 TVector3 p3Trk0(Get2Prong()->PxProng(0),Get2Prong()->PyProng(0),Get2Prong()->PzProng(0)); // from D0
203 TVector3 p3Trk1(Get2Prong()->PxProng(1),Get2Prong()->PyProng(1),Get2Prong()->PzProng(1)); // from D0
204 TVector3 p3Trk2(PxProng(0),PyProng(0),PzProng(0)); // pi_s
205
206 TVector3 perp = p3Trk0.Cross(p3Trk1);
207 Double_t theta = p3Trk2.Angle(perp);
208 if(theta>(TMath::Pi()-theta)) theta = TMath::Pi() - theta;
209 theta = TMath::Pi()/2. - theta;
210
211 if(theta>cutsDstar[4]) return kFALSE;
212
213 Double_t alpha = p3Trk0.Angle(p3Trk2);
214 Double_t belta = p3Trk1.Angle(p3Trk2);
215
216 Double_t cosphi01 = TMath::Cos(alpha) / TMath::Cos(theta);
217 Double_t cosphi02 = TMath::Cos(belta) / TMath::Cos(theta);
218
219 Double_t phi01 = TMath::ACos(cosphi01);
220 Double_t phi02 = TMath::ACos(cosphi02);
221 Double_t phi00 = p3Trk0.Angle(p3Trk1);
222
223 if((phi01>phi00) || (phi02>phi00)) return kFALSE;
224
225 return kTRUE;
226}
227//-----------------------------------------------------------------------------
a07ad8e0 228
229//-----------------------------------------------------------------------------
230Bool_t AliAODRecoCascadeHF::SelectLctoV0(const Double_t *cutsLctoV0,
231 Bool_t okLck0sp, Bool_t okLcLpi) const
232{
233 // cuts on Lambdac candidates to V0+bachelor
234 // (to be passed to AliAODRecoDecayHF3Prong::SelectLctoV0())
235 // 0 = inv. mass half width in K0s hypothesis [GeV]
236 // 1 = inv. mass half width in Lambda hypothesis [GeV]
237 // 2 = inv. mass V0 in K0s hypothesis half width [GeV]
238 // 3 = inv. mass V0 in Lambda hypothesis half width [GeV]
239 // 4 = pT min Bachelor track [GeV/c]
240 // 5 = pT min V0-Positive track [GeV/c]
241 // 6 = pT min V0-Negative track [GeV/c]
242 // 7 = dca cut on the V0 (cm)
243 // 8 = dca cut on the cascade (cm)
244
245// if ( !Getv0() || !Getv0PositiveTrack() || !Getv0NegativeTrack() )
246// { AliInfo(Form("Not adapted for ESDv0s, return true...")); return false; }
247
248 Double_t mLck0sp,mLcLpi;
249 okLck0sp=1; okLcLpi=1;
250
251 Double_t mLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
252 Double_t mk0sPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
253 Double_t mLPDG = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
254
255 // k0s + p
256 double mk0s = Getv0()->MassK0Short();
257 mLck0sp = InvMassLctoK0sP();
258
259 // lambda + pi
260 double mlambda = Getv0()->MassLambda();
261 double malambda = Getv0()->MassAntiLambda();
262 mLcLpi = InvMassLctoLambdaPi();
263
264 // cut on Lc mass
265 // with k0s p hypothesis
266 if(TMath::Abs(mLck0sp-mLcPDG)>cutsLctoV0[0]) okLck0sp = 0;
267 // with Lambda pi hypothesis
268 if(TMath::Abs(mLcLpi-mLcPDG)>cutsLctoV0[1]) okLcLpi = 0;
269
270 // cuts on the v0 mass
271 if(TMath::Abs(mk0s-mk0sPDG)>cutsLctoV0[2]) okLck0sp = 0;
272 if( TMath::Abs(mlambda-mLPDG)>cutsLctoV0[3] &&
273 TMath::Abs(malambda-mLPDG)>cutsLctoV0[3] ) okLcLpi = 0;
274
275 if(!okLck0sp && !okLcLpi) return 0;
276
277 // cuts on the minimum pt of the tracks
278 if(TMath::Abs(GetBachelor()->Pt()) < cutsLctoV0[4]) return 0;
279 if(TMath::Abs(Getv0PositiveTrack()->Pt()) < cutsLctoV0[5]) return 0;
280 if(TMath::Abs(Getv0NegativeTrack()->Pt()) < cutsLctoV0[6]) return 0;
281
282 // cut on the v0 dca
283 if(TMath::Abs(Getv0()->DcaV0Daughters()) > cutsLctoV0[7]) return 0;
284
285 // cut on the cascade dca
286 if( TMath::Abs(GetDCA(0))>cutsLctoV0[8] ||
287 TMath::Abs(Getv0()->DcaPosToPrimVertex())>cutsLctoV0[8] ||
288 TMath::Abs(Getv0()->DcaNegToPrimVertex())>cutsLctoV0[8] ) return 0;
289
290 return true;
291
292}
293//-----------------------------------------------------------------------------