]>
Commit | Line | Data |
---|---|---|
259c3296 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, 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 | * QA class of primary vertex study for Heavy Flavor electrons * | |
18 | * * | |
19 | * Authors: * | |
20 | * MinJung Kweon <minjung@physi.uni-heidelberg.de> * | |
21 | * * | |
22 | **************************************************************************/ | |
23 | ||
24 | ||
25 | #include <iostream> | |
26 | #include <TH2F.h> | |
27 | #include <TH3F.h> | |
28 | #include "TLorentzVector.h" | |
29 | #include <TNtuple.h> | |
30 | #include <TParticle.h> | |
31 | ||
32 | #include "AliESDEvent.h" | |
33 | #include <AliESDtrack.h> | |
34 | #include <AliStack.h> | |
35 | ||
36 | #include <AliLog.h> | |
37 | #include <AliKFParticle.h> | |
38 | #include <AliKFVertex.h> | |
39 | ||
40 | #include "AliHFEpriVtx.h" | |
41 | ||
42 | ||
43 | ClassImp(AliHFEpriVtx) | |
44 | ||
45 | //_______________________________________________________________________________________________ | |
46 | AliHFEpriVtx::AliHFEpriVtx(): | |
47 | fESD1(0x0) | |
48 | ,fStack(0x0) | |
49 | ,fNtrackswoPid(0) | |
50 | ,fHNtrackswoPid(0x0) | |
51 | ,fNESDprimVtxContributor(0x0) | |
52 | ,fNESDprimVtxIndices(0x0) | |
53 | ,fDiffDCAvsPt(0x0) | |
54 | ,fDiffDCAvsNt(0x0) | |
55 | { | |
56 | // | |
57 | // Default constructor | |
58 | // | |
59 | ||
60 | Init(); | |
61 | ||
62 | } | |
63 | ||
64 | //_______________________________________________________________________________________________ | |
65 | AliHFEpriVtx::AliHFEpriVtx(const AliHFEpriVtx &p): | |
66 | TObject(p) | |
67 | ,fESD1(0x0) | |
68 | ,fStack(0x0) | |
69 | ,fNtrackswoPid(p.fNtrackswoPid) | |
70 | ,fHNtrackswoPid(0x0) | |
71 | ,fNESDprimVtxContributor(0x0) | |
72 | ,fNESDprimVtxIndices(0x0) | |
73 | ,fDiffDCAvsPt(0x0) | |
74 | ,fDiffDCAvsNt(0x0) | |
75 | { | |
76 | // | |
77 | // Copy constructor | |
78 | // | |
79 | } | |
80 | ||
81 | //_______________________________________________________________________________________________ | |
82 | AliHFEpriVtx& | |
83 | AliHFEpriVtx::operator=(const AliHFEpriVtx &) | |
84 | { | |
85 | // | |
86 | // Assignment operator | |
87 | // | |
88 | ||
89 | AliInfo("Not yet implemented."); | |
90 | return *this; | |
91 | } | |
92 | ||
93 | //_______________________________________________________________________________________________ | |
94 | AliHFEpriVtx::~AliHFEpriVtx() | |
95 | { | |
96 | // | |
97 | // Destructor | |
98 | // | |
99 | ||
100 | AliInfo("Analysis Done."); | |
101 | } | |
102 | ||
103 | //__________________________________________ | |
104 | void AliHFEpriVtx::Init() | |
105 | { | |
106 | // | |
107 | // initialize counters | |
108 | // | |
109 | ||
110 | fNtrackswoPid = 0; | |
111 | for (int i=0; i<10; i++){ | |
112 | fPrimVtx[i].fNtrackCount = 0 ; | |
113 | fPrimVtx[i].fNprimVtxContributorCount = 0 ; | |
114 | } | |
115 | } | |
116 | ||
117 | //_______________________________________________________________________________________________ | |
118 | void AliHFEpriVtx::CreateHistograms(TString hnopt) | |
119 | { | |
120 | // | |
121 | // create histograms | |
122 | // | |
123 | ||
124 | fkSourceLabel[kAll]="all"; | |
125 | fkSourceLabel[kDirectCharm]="directCharm"; | |
126 | fkSourceLabel[kDirectBeauty]="directBeauty"; | |
127 | fkSourceLabel[kBeautyCharm]="beauty2charm"; | |
128 | fkSourceLabel[kGamma]="gamma"; | |
129 | fkSourceLabel[kPi0]="pi0"; | |
130 | fkSourceLabel[kElse]="others"; | |
131 | fkSourceLabel[kBeautyGamma]="beauty22gamma"; | |
132 | fkSourceLabel[kBeautyPi0]="beauty22pi0"; | |
133 | fkSourceLabel[kBeautyElse]="beauty22others"; | |
134 | ||
135 | ||
136 | TString hname; | |
137 | for (Int_t isource = 0; isource < 10; isource++ ){ | |
138 | ||
139 | hname=hnopt+"ntracks_"+fkSourceLabel[isource]; | |
140 | fPrimVtx[isource].fNtracks = new TH1F(hname,hname,50,0,50); | |
141 | hname=hnopt+"nPrimVtxContributor_"+fkSourceLabel[isource]; | |
142 | fPrimVtx[isource].fNprimVtxContributor = new TH1F(hname,hname,100,0,100); | |
143 | hname=hnopt+"PtElec_"+fkSourceLabel[isource]; | |
144 | fPrimVtx[isource].fPtElec = new TH1F(hname,hname,150,0,30); | |
145 | hname=hnopt+"PtElecContributor_"+fkSourceLabel[isource]; | |
146 | fPrimVtx[isource].fPtElecContributor = new TH1F(hname,hname,150,0,30); | |
147 | ||
148 | } | |
149 | ||
150 | hname=hnopt+"ntrackswopid"; | |
151 | fHNtrackswoPid = new TH1F(hname,hname,50,0,50); | |
152 | hname=hnopt+"nESDprimVtxContributor"; | |
153 | fNESDprimVtxContributor = new TH1I(hname,hname,100,0,100); | |
154 | hname=hnopt+"nESDprimVtxIndices"; | |
155 | fNESDprimVtxIndices= new TH1I(hname,hname,100,0,100); | |
156 | hname=hnopt+"diffDCAvsPt"; | |
157 | fDiffDCAvsPt = new TH2F(hname,hname,150,0,30,500,0,1); | |
158 | hname=hnopt+"diffDCAvsNt"; | |
159 | fDiffDCAvsNt = new TH2F(hname,hname,100,0,100,500,0,1); | |
160 | ||
161 | } | |
162 | ||
163 | //_______________________________________________________________________________________________ | |
164 | void AliHFEpriVtx::CountNtracks(Int_t sourcePart, Int_t recpid, Double_t recprob) | |
165 | { | |
166 | // | |
167 | // count number of tracks passed certain cuts | |
168 | // | |
169 | ||
170 | fNtrackswoPid++; | |
171 | ||
172 | if (!recpid && recprob>0.5) | |
173 | fPrimVtx[kAll].fNtrackCount++; | |
174 | if(sourcePart<0) return; | |
175 | fPrimVtx[sourcePart].fNtrackCount++; | |
176 | ||
177 | } | |
178 | ||
179 | //_______________________________________________________________________________________________ | |
180 | void AliHFEpriVtx::FillNtracks() | |
181 | { | |
182 | // | |
183 | // count number of tracks passed certain cuts | |
184 | // | |
185 | ||
186 | fHNtrackswoPid->Fill(fNtrackswoPid); | |
187 | for (int i=0; i<10; i++){ | |
188 | fPrimVtx[i].fNtracks->Fill(fPrimVtx[i].fNtrackCount); | |
189 | } | |
190 | ||
191 | } | |
192 | ||
193 | //_______________________________________________________________________________________________ | |
194 | Int_t AliHFEpriVtx::GetMCPID(AliESDtrack *track) | |
195 | { | |
196 | // | |
197 | // get MC pid | |
198 | // | |
199 | ||
200 | Int_t label = TMath::Abs(track->GetLabel()); | |
201 | TParticle* mcpart = fStack->Particle(label); | |
202 | if ( !mcpart ) return 0; | |
203 | Int_t pdgCode = mcpart->GetPdgCode(); | |
204 | ||
205 | return pdgCode; | |
206 | } | |
207 | ||
208 | //_______________________________________________________________________________________________ | |
209 | void AliHFEpriVtx::GetNPriVxtContributor() | |
210 | { | |
211 | // | |
212 | // count number of primary vertex contributor | |
213 | // | |
214 | ||
215 | const AliESDVertex *primvtx = fESD1->GetPrimaryVertex(); | |
216 | fNESDprimVtxContributor->Fill(primvtx->GetNContributors()); | |
217 | fNESDprimVtxIndices->Fill(primvtx->GetNIndices()); | |
218 | } | |
219 | ||
220 | //_______________________________________________________________________________________________ | |
221 | void AliHFEpriVtx::CountPriVxtElecContributor(AliESDtrack *ESDelectron, Int_t sourcePart, Int_t recpid, Double_t recprob) | |
222 | { | |
223 | // | |
224 | // count number of electrons contributing to the primary vertex | |
225 | // | |
226 | ||
227 | ||
228 | if (recpid || recprob<0.5) return; | |
229 | ||
230 | // get track id of our selected electron | |
231 | Int_t elecTrkID = ESDelectron->GetID(); | |
232 | ||
233 | Int_t label = TMath::Abs(ESDelectron->GetLabel()); | |
234 | TParticle* mcpart = fStack->Particle(label); | |
235 | ||
236 | AliKFParticle::SetField(fESD1->GetMagneticField()); | |
237 | AliKFParticle kfElectron(*ESDelectron,11); | |
238 | ||
239 | // prepare kfprimary vertex | |
240 | AliKFVertex kfESDprimary; | |
241 | ||
242 | // Reconstructed Primary Vertex (with ESD tracks) | |
243 | const AliESDVertex *primvtx = fESD1->GetPrimaryVertex(); | |
244 | // Int_t nt = primvtx->GetNContributors(); | |
245 | Int_t n=primvtx->GetNIndices(); | |
246 | ||
247 | if (n>0 && primvtx->GetStatus()){ | |
248 | ||
249 | kfESDprimary = AliKFVertex(*primvtx); | |
250 | Double_t dcaBFrmElec = kfElectron.GetDistanceFromVertexXY(kfESDprimary); | |
251 | ||
252 | UShort_t *priIndex = primvtx->GetIndices(); | |
253 | ||
254 | fPrimVtx[kAll].fPtElec->Fill(mcpart->Pt()); | |
255 | if(sourcePart>=0) fPrimVtx[sourcePart].fPtElec->Fill(mcpart->Pt()); | |
256 | ||
257 | for (Int_t i=0;i<n;i++){ | |
258 | ||
259 | Int_t idx = Int_t(priIndex[i]); | |
260 | if (idx == elecTrkID){ | |
261 | fPrimVtx[kAll].fNprimVtxContributorCount++; | |
262 | fPrimVtx[kAll].fPtElecContributor->Fill(mcpart->Pt()); | |
263 | if(sourcePart<0) continue; | |
264 | fPrimVtx[sourcePart].fNprimVtxContributorCount++; | |
265 | fPrimVtx[sourcePart].fPtElecContributor->Fill(mcpart->Pt()); | |
266 | ||
267 | kfESDprimary -= kfElectron; | |
268 | Double_t dcaAFrmElec = kfElectron.GetDistanceFromVertexXY(kfESDprimary); | |
269 | fDiffDCAvsPt->Fill(mcpart->Pt(),dcaBFrmElec-dcaAFrmElec); | |
270 | fDiffDCAvsNt->Fill(n,dcaBFrmElec-dcaAFrmElec); | |
271 | } | |
272 | } | |
273 | } | |
274 | ||
275 | } | |
276 | ||
277 | //_______________________________________________________________________________________________ | |
278 | void AliHFEpriVtx::FillNprimVtxContributor() | |
279 | { | |
280 | // | |
281 | // Fill histogram with number of electrons contributing to the primary vertex | |
282 | // | |
283 | ||
284 | for (int i=0; i<10; i++){ | |
285 | fPrimVtx[i].fNprimVtxContributor->Fill(fPrimVtx[i].fNprimVtxContributorCount); | |
286 | } | |
287 | ||
288 | } |