]>
Commit | Line | Data |
---|---|---|
4292b3b6 | 1 | #include <TParticle.h> |
2 | #include <TLorentzVector.h> | |
3 | #include <TClonesArray.h> | |
4 | ||
5 | #include "AliMCEventHandler.h" | |
6 | #include "AliMCEvent.h" | |
7 | #include "AliStack.h" | |
8 | #include "AliAODMCParticle.h" | |
9 | #include "AliESDEvent.h" | |
10 | #include "AliESDMuonTrack.h" | |
11 | #include "AliESDMuonCluster.h" | |
12 | #include "AliAODTrack.h" | |
13 | #include "AliAODMuonTrack.h" | |
14 | #include "AliMCMuonTrack.h" | |
15 | ||
16 | ClassImp(AliMCMuonTrack) | |
17 | ||
18 | const Double_t AliMCMuonTrack::fgkSigmaCut = 10.; | |
19 | ||
20 | //----------------------------------------------------------------------------- | |
21 | AliMCMuonTrack::AliMCMuonTrack() : | |
22 | AliAODMuonTrack(), | |
23 | fIsFull(kFALSE), | |
24 | fPGen(), | |
25 | fTrackIndex(-1), | |
26 | fTrackPDGCode(0), | |
27 | fSource(-1), | |
28 | fNParents(0), | |
29 | fOscillation(kFALSE), | |
30 | fWeight(0.) | |
31 | { | |
32 | // | |
33 | // default constructor | |
34 | // | |
35 | for (Int_t i=0; i<fgkNParentsMax; i++) { | |
36 | fParentIndex[i] = -1; | |
37 | fParentPDGCode[i] = 0; | |
38 | } | |
39 | for (Int_t i=0; i<4; i++) { | |
40 | fQuarkIndex[i] = -1; | |
41 | fQuarkPDGCode[i] = 0; | |
42 | } | |
43 | } | |
44 | ||
45 | //----------------------------------------------------------------------------- | |
46 | AliMCMuonTrack::AliMCMuonTrack(AliAODTrack *trkAOD, TClonesArray *mcClArr, Bool_t full) : | |
47 | AliAODMuonTrack(trkAOD), | |
48 | fIsFull(full), | |
49 | fPGen(), | |
50 | fTrackIndex(-1), | |
51 | fTrackPDGCode(0), | |
52 | fSource(-1), | |
53 | fNParents(0), | |
54 | fOscillation(kFALSE), | |
55 | fWeight(0.) | |
56 | { | |
57 | // | |
58 | // default constructor | |
59 | // | |
60 | AliAODMCParticle *pMC = this->FindTrackRef(trkAOD, mcClArr); | |
61 | if (pMC) this->SetMCInfo(pMC, mcClArr); | |
62 | } | |
63 | ||
64 | //----------------------------------------------------------------------------- | |
65 | AliMCMuonTrack::AliMCMuonTrack(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH, Bool_t full) : | |
66 | AliAODMuonTrack(trkESD), | |
67 | fIsFull(full), | |
68 | fPGen(), | |
69 | fTrackIndex(-1), | |
70 | fTrackPDGCode(0), | |
71 | fSource(-1), | |
72 | fNParents(0), | |
73 | fOscillation(kFALSE), | |
74 | fWeight(0.) | |
75 | { | |
76 | // | |
77 | // default constructor | |
78 | // | |
79 | TParticle *pMC = this->FindTrackRef(trkESD, esd, mcH); // do nothing when running on official train | |
80 | if (pMC) this->SetMCInfo(pMC, mcH); | |
81 | // MC infomation is not implement when running on the official train for they depend on the MUON module. | |
82 | // In the case of private running, one can get the information from ESD by uncommenting the definitions | |
83 | // of FindTrackRef() and CovESDtoMuonTrack() | |
84 | // and including the following head files in both AliMCMuonTrack.h and this file, | |
85 | // #include "AliMUONRecoCheck.h" | |
86 | // #include "AliMUONVTrackStore.h" | |
87 | // #include "AliMUONVClusterStore.h" | |
88 | // #include "AliMUONTrackParam.h" | |
89 | // #include "AliMUONTrack.h" | |
90 | // #include "AliMUONVCluster.h" | |
91 | // #include "AliMUONESDInterface.h" | |
92 | } | |
93 | ||
94 | ||
95 | ||
96 | //----------------------------------------------------------------------------- | |
97 | AliMCMuonTrack::~AliMCMuonTrack() | |
98 | { | |
99 | // | |
100 | // destructor | |
101 | // | |
102 | } | |
103 | ||
104 | //----------------------------------------------------------------------------- | |
105 | AliAODMCParticle* AliMCMuonTrack::FindTrackRef(AliAODTrack *trkAOD, TClonesArray *mcClArr) | |
106 | { | |
107 | AliAODMCParticle *pMC = 0; | |
108 | fTrackIndex = trkAOD->GetLabel(); | |
109 | if (fTrackIndex>=0) | |
110 | pMC = (AliAODMCParticle*)mcClArr->At(fTrackIndex); | |
111 | return pMC; | |
112 | } | |
113 | ||
114 | //----------------------------------------------------------------------------- | |
115 | TParticle* AliMCMuonTrack::FindTrackRef(AliESDMuonTrack *trkESD, AliESDEvent *esd, AliMCEventHandler *mcH) | |
116 | { | |
117 | TParticle *pMCRef = 0; | |
118 | ||
119 | /*AliMUONTrack *trkMuon = CovESDtoMuonTrack(*trkESD); | |
120 | AliMUONRecoCheck rc(esd,mcH); | |
121 | AliMUONVTrackStore* trkRefArr = rc.TrackRefs(-1); | |
122 | TIter next(trkRefArr->CreateIterator()); | |
123 | AliMUONTrack* trkRef = 0; | |
124 | while ((trkRef=static_cast<AliMUONTrack*>(next()))) { | |
125 | Bool_t trkCompArr[10]; | |
126 | Int_t nMatchClusters = trkMuon->CompatibleTrack(trkRef, fgkSigmaCut, trkCompArr); | |
127 | Double_t matchClusterFrac = ((Double_t)nMatchClusters) / ((Double_t)trkMuon->GetNClusters()); | |
128 | if ((trkCompArr[0] || trkCompArr[1] || trkCompArr[2] || trkCompArr[3]) && | |
129 | (trkCompArr[6] || trkCompArr[7] || trkCompArr[8] || trkCompArr[9]) && matchClusterFrac>0.5) { | |
130 | fTrackIndex = trkRef->GetUniqueID(); | |
131 | trkRefArr->Remove(*trkRef); | |
132 | break; | |
133 | } | |
134 | } | |
135 | if (fTrackIndex>=0) pMCRef = mcH->MCEvent()->Stack()->Particle(fTrackIndex);*/ | |
136 | ||
137 | return pMCRef; | |
138 | } | |
139 | ||
140 | //----------------------------------------------------------------------------- | |
141 | void AliMCMuonTrack::SetMCInfo(AliAODMCParticle *pMC, TClonesArray *mcClArr) | |
142 | { | |
143 | fPGen.SetPxPyPzE(pMC->Px(), pMC->Py(), pMC->Pz(), pMC->E()); | |
144 | fTrackPDGCode = pMC->GetPdgCode(); | |
145 | if (TMath::Abs(fTrackPDGCode)!=13) { | |
146 | fSource = 4; | |
147 | return; | |
148 | } | |
149 | ||
150 | Int_t lineM = pMC->GetMother(); | |
151 | if (lineM<0) { | |
152 | fSource = 2; | |
153 | return; | |
154 | } | |
155 | ||
156 | Bool_t isPrimary = | |
157 | ((AliAODMCParticle*)mcClArr->At(lineM))->IsPrimary(); | |
158 | if (!isPrimary) { | |
159 | fSource = 3; | |
160 | return; | |
161 | } | |
162 | ||
163 | this->FillHistoryParents(lineM, mcClArr); | |
164 | fSource = this->SelectHFMuon(); | |
165 | return; | |
166 | } | |
167 | ||
168 | //----------------------------------------------------------------------------- | |
169 | void AliMCMuonTrack::SetMCInfo(TParticle *pMC, AliMCEventHandler *mcH) | |
170 | { | |
171 | fPGen.SetPxPyPzE(pMC->Px(), pMC->Py(), pMC->Pz(), pMC->Energy()); | |
172 | fTrackPDGCode = pMC->GetPdgCode(); | |
173 | if (TMath::Abs(fTrackPDGCode)!=13) { | |
174 | fSource = 4; | |
175 | return; | |
176 | } | |
177 | ||
178 | Int_t lineM = pMC->GetFirstMother(); | |
179 | if (lineM<0) { | |
180 | fSource = 2; | |
181 | return; | |
182 | } | |
183 | ||
184 | AliStack *stack = mcH->MCEvent()->Stack(); | |
185 | if (lineM>=stack->GetNprimary()) { | |
186 | fSource = 3; | |
187 | return; | |
188 | } | |
189 | ||
190 | this->FillHistoryParents(lineM, stack); | |
191 | fSource = this->SelectHFMuon(); | |
192 | return; | |
193 | } | |
194 | ||
195 | //----------------------------------------------------------------------------- | |
196 | void AliMCMuonTrack::FillHistoryParents(Int_t lineM, TClonesArray *mcClArr) | |
197 | { | |
198 | Int_t countP=-1, pdg=0; | |
199 | Int_t parents[10], parLine[10]; | |
200 | AliAODMCParticle *mother = 0; | |
201 | while(lineM>=0){ | |
202 | mother = (AliAODMCParticle*)mcClArr->At(lineM); | |
203 | pdg = mother->GetPdgCode(); | |
204 | if(pdg==92 || pdg==21 || TMath::Abs(pdg)<10 || IsDiquark(pdg)) break; | |
205 | parents[++countP] = pdg; | |
206 | parLine[countP] = lineM; | |
207 | lineM = mother->GetMother(); | |
208 | } | |
209 | for(Int_t i=0; i<=countP; i++){ | |
210 | fParentIndex[i] = parLine[countP-i]; | |
211 | fParentPDGCode[i] = parents[countP-i]; | |
212 | } | |
213 | fNParents = countP + 1; | |
214 | ||
215 | if (fIsFull && lineM>=0) this->FillHistoryQuarks(lineM, mcClArr); | |
216 | return; | |
217 | } | |
218 | ||
219 | //----------------------------------------------------------------------------- | |
220 | void AliMCMuonTrack::FillHistoryParents(Int_t lineM, AliStack *stack) | |
221 | { | |
222 | Int_t countP=-1, pdg=0; | |
223 | Int_t parents[10], parLine[10]; | |
224 | TParticle *mother = 0; | |
225 | while(lineM>=0){ | |
226 | mother = stack->Particle(lineM); | |
227 | pdg = mother->GetPdgCode(); | |
228 | if(pdg==92 || pdg==21 || TMath::Abs(pdg)<10 || IsDiquark(pdg)) break; | |
229 | parents[++countP] = pdg; | |
230 | parLine[countP] = lineM; | |
231 | lineM = mother->GetFirstMother(); | |
232 | } | |
233 | for(Int_t i=0; i<=countP; i++){ | |
234 | fParentIndex[i] = parLine[countP-i]; | |
235 | fParentPDGCode[i] = parents[countP-i]; | |
236 | } | |
237 | fNParents = countP + 1; | |
238 | ||
239 | if (fIsFull && lineM>=0) this->FillHistoryQuarks(lineM, stack); | |
240 | return; | |
241 | } | |
242 | ||
243 | //----------------------------------------------------------------------------- | |
244 | void AliMCMuonTrack::FillHistoryQuarks(Int_t lineM, TClonesArray *mcClArr) | |
245 | { | |
246 | // method in $ALICE_ROOT/MUON/AliMUONTrackLight.cxx | |
247 | if (lineM<0) return; | |
248 | Int_t countP=-1, pdg=0; | |
249 | AliAODMCParticle *mother = 0; | |
250 | while(lineM>=0){ | |
251 | mother = (AliAODMCParticle*)mcClArr->At(lineM); | |
252 | pdg = mother->GetPdgCode(); | |
253 | fQuarkIndex[++countP] = lineM; | |
254 | fQuarkPDGCode[countP] = pdg; | |
255 | lineM = mother->GetMother(); | |
256 | } | |
257 | ||
258 | // for PYTHIA checking | |
259 | countP = 1; | |
260 | for(Int_t par=0; par<4; par++) { | |
261 | if(TMath::Abs(this->GetQuarkPDGCode(par))<6) { countP=par; break; } | |
262 | } | |
263 | if(this->GetQuarkIndex(countP)>-1 && (this->GetParentFlavour(0)==4 || this->GetParentFlavour(0)==5)) { | |
264 | if(this->GetParentFlavour(0)!=TMath::Abs(this->GetQuarkPDGCode(countP))) { | |
265 | AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n", | |
266 | this->GetParentFlavour(0), TMath::Abs(this->GetQuarkPDGCode(countP)))); | |
267 | ||
268 | pdg = this->GetQuarkPDGCode(countP); | |
269 | Int_t line = this->GetQuarkIndex(countP); | |
270 | this->ResetQuarkInfo(); | |
271 | while(TMath::Abs(pdg)!=this->GetParentFlavour(0)) { | |
272 | pdg = ((AliAODMCParticle*)mcClArr->At(++line))->GetPdgCode(); | |
273 | } | |
274 | while(line>=0){ | |
275 | mother = (AliAODMCParticle*)mcClArr->At(line); | |
276 | pdg = mother->GetPdgCode(); | |
277 | fQuarkIndex[countP] = line; | |
278 | fQuarkPDGCode[countP++] = pdg; | |
279 | line = mother->GetMother(); | |
280 | } | |
281 | } | |
282 | } | |
283 | return; | |
284 | } | |
285 | ||
286 | //----------------------------------------------------------------------------- | |
287 | void AliMCMuonTrack::FillHistoryQuarks(Int_t lineM, AliStack *stack) | |
288 | { | |
289 | // method in $ALICE_ROOT/MUON/AliMUONTrackLight.cxx | |
290 | if (lineM<0) return; | |
291 | Int_t countP=-1, pdg=0; | |
292 | TParticle *mother = 0; | |
293 | while(lineM>=0){ | |
294 | mother = stack->Particle(lineM); | |
295 | pdg = mother->GetPdgCode(); | |
296 | fQuarkIndex[++countP] = lineM; | |
297 | fQuarkPDGCode[countP] = pdg; | |
298 | lineM = mother->GetFirstMother(); | |
299 | } | |
300 | ||
301 | // for PYTHIA checking | |
302 | countP = 1; | |
303 | for(Int_t par=0; par<4; par++) { | |
304 | if(TMath::Abs(this->GetQuarkPDGCode(par))<6) { countP=par; break; } | |
305 | } | |
306 | if(this->GetQuarkIndex(countP)>-1 && (this->GetParentFlavour(0)==4 || this->GetParentFlavour(0)==5)) { | |
307 | if(this->GetParentFlavour(0)!=TMath::Abs(this->GetQuarkPDGCode(countP))) { | |
308 | AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n", | |
309 | this->GetParentFlavour(0), TMath::Abs(this->GetQuarkPDGCode(countP)))); | |
310 | ||
311 | pdg = this->GetQuarkPDGCode(countP); | |
312 | Int_t line = this->GetQuarkIndex(countP); | |
313 | this->ResetQuarkInfo(); | |
314 | while(TMath::Abs(pdg)!=this->GetParentFlavour(0)) { | |
315 | pdg = stack->Particle(++line)->GetPdgCode(); | |
316 | } | |
317 | while(line>=0){ | |
318 | mother = stack->Particle(line); | |
319 | pdg = mother->GetPdgCode(); | |
320 | fQuarkIndex[countP] = line; | |
321 | fQuarkPDGCode[countP++] = pdg; | |
322 | line = mother->GetFirstMother(); | |
323 | } | |
324 | } | |
325 | } | |
326 | return; | |
327 | } | |
328 | ||
329 | //----------------------------------------------------------------------------- | |
330 | /*AliMUONTrack* AliMCMuonTrack::CovESDtoMuonTrack(AliESDMuonTrack &trkESD) | |
331 | { | |
332 | // method in $ALICE_ROOT/PWG3/muondep/AliAnalysisTaskESDMCLabelAddition.cxx | |
333 | AliMUONTrack *trkMuon = new AliMUONTrack(); | |
334 | if(!trkESD.ClustersStored()) return trkMuon; | |
335 | ||
336 | AliMUONTrackParam param; | |
337 | AliMUONESDInterface::GetParamAtFirstCluster(trkESD, param); | |
338 | AliMUONESDInterface::GetParamCov(trkESD, param); | |
339 | AliMUONVClusterStore *cStore = AliMUONESDInterface::NewClusterStore(); | |
340 | AliMUONVCluster *cluster = cStore->CreateCluster(0,0,0); | |
341 | ||
342 | AliESDMuonCluster *esdCluster = (AliESDMuonCluster*)trkESD.GetClusters().First(); | |
343 | while (esdCluster) { | |
344 | AliMUONESDInterface::ESDToMUON(*esdCluster, *cluster); | |
345 | param.SetZ(cluster->GetZ()); | |
346 | trkMuon->AddTrackParamAtCluster(param, *cluster, kTRUE); | |
347 | esdCluster = (AliESDMuonCluster*)trkESD.GetClusters().After(esdCluster); | |
348 | } | |
349 | ||
350 | delete cluster; cluster = 0; | |
351 | delete cStore; cStore = 0; | |
352 | return trkMuon; | |
353 | }*/ | |
354 | ||
355 | //----------------------------------------------------------------------------- | |
356 | Int_t AliMCMuonTrack::SelectHFMuon() | |
357 | { | |
358 | Int_t flv = GetParentFlavour(0); | |
359 | if (flv!=4 && flv!=5) return 2; | |
360 | ||
361 | Bool_t isRes = kFALSE; | |
362 | Int_t i=0, nparents=this->GetNParents(); | |
363 | while (i<nparents && !isRes) { | |
364 | isRes = IsMotherAResonance(i++); | |
365 | } | |
366 | ||
367 | if (isRes) return 2; | |
368 | if (flv==5) return 0; | |
369 | else return 1; | |
370 | } | |
371 | ||
372 | //----------------------------------------------------------------------------- | |
373 | Bool_t AliMCMuonTrack::IsDiquark(Int_t pdg) | |
374 | { | |
375 | // copy from $ALICE_ROOT/MUON/AliMUONTrackLight.cxx | |
376 | pdg = TMath::Abs(pdg); | |
377 | if(pdg>1000 && (pdg%100)<10) return kTRUE; | |
378 | else return kFALSE; | |
379 | } | |
380 | ||
381 | //----------------------------------------------------------------------------- | |
382 | void AliMCMuonTrack::ResetQuarkInfo() | |
383 | { | |
384 | for(Int_t pos=1; pos<4; pos++) { | |
385 | fQuarkIndex[pos] = -1; | |
386 | fQuarkPDGCode[pos] = 0; | |
387 | } | |
388 | return; | |
389 | } | |
390 | ||
391 | //----------------------------------------------------------------------------- | |
392 | Int_t AliMCMuonTrack::GetParentFlavour(Int_t i) const | |
393 | { | |
394 | Int_t pdg = GetParentPDGCode(i); | |
395 | pdg = TMath::Abs(pdg/100); | |
396 | if(pdg>9) pdg /= 10; | |
397 | return pdg; | |
398 | } | |
399 | ||
400 | //----------------------------------------------------------------------------- | |
401 | Bool_t AliMCMuonTrack::IsMotherAResonance(Int_t i) const | |
402 | { | |
403 | // copy from $ALICE_ROOT/MUON/AliMUONTrackLight.cxx | |
404 | Int_t pdg = GetParentPDGCode(i); | |
405 | Int_t id=pdg%100000; | |
406 | return (!((id-id%10)%110)); | |
407 | } |