]>
Commit | Line | Data |
---|---|---|
96b85d73 | 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 | /* $Id: AliAnalysisTaskESDfilter.cxx 24535 2008-03-16 22:43:30Z fca $ */ | |
17 | ||
18 | #include <TChain.h> | |
7b4cee8d | 19 | #include <TTree.h> |
ec6465fe | 20 | #include <TList.h> |
e6500713 | 21 | #include <TArrayI.h> |
d7bdc804 | 22 | #include <TRandom.h> |
da97a08a | 23 | #include <TParticle.h> |
96b85d73 | 24 | |
25 | #include "AliAnalysisTaskESDfilter.h" | |
26 | #include "AliAnalysisManager.h" | |
27 | #include "AliESDEvent.h" | |
da97a08a | 28 | #include "AliStack.h" |
96b85d73 | 29 | #include "AliAODEvent.h" |
da97a08a | 30 | #include "AliMCEvent.h" |
31 | #include "AliMCEventHandler.h" | |
96b85d73 | 32 | #include "AliESDInputHandler.h" |
33 | #include "AliAODHandler.h" | |
da97a08a | 34 | #include "AliAODMCParticle.h" |
96b85d73 | 35 | #include "AliAnalysisFilter.h" |
96b85d73 | 36 | #include "AliESDMuonTrack.h" |
37 | #include "AliESDVertex.h" | |
38 | #include "AliESDv0.h" | |
39 | #include "AliESDkink.h" | |
40 | #include "AliESDcascade.h" | |
41 | #include "AliESDPmdTrack.h" | |
42 | #include "AliESDCaloCluster.h" | |
43 | #include "AliESDCaloCells.h" | |
44 | #include "AliMultiplicity.h" | |
45 | #include "AliLog.h" | |
46 | ||
47 | ClassImp(AliAnalysisTaskESDfilter) | |
48 | ||
49 | //////////////////////////////////////////////////////////////////////// | |
50 | ||
51 | AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(): | |
52 | AliAnalysisTaskSE(), | |
53 | fTrackFilter(0x0), | |
54 | fKinkFilter(0x0), | |
d7bdc804 | 55 | fV0Filter(0x0), |
9a244178 | 56 | fCascadeFilter(0x0), |
d7bdc804 | 57 | fHighPthreshold(0), |
980e4563 | 58 | fPtshape(0x0) |
96b85d73 | 59 | { |
60 | // Default constructor | |
61 | } | |
62 | ||
980e4563 | 63 | AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name): |
96b85d73 | 64 | AliAnalysisTaskSE(name), |
65 | fTrackFilter(0x0), | |
66 | fKinkFilter(0x0), | |
d7bdc804 | 67 | fV0Filter(0x0), |
9a244178 | 68 | fCascadeFilter(0x0), |
d7bdc804 | 69 | fHighPthreshold(0), |
980e4563 | 70 | fPtshape(0x0) |
96b85d73 | 71 | { |
72 | // Constructor | |
73 | } | |
74 | ||
75 | void AliAnalysisTaskESDfilter::UserCreateOutputObjects() | |
76 | { | |
77 | // Create the output container | |
78 | OutputTree()->GetUserInfo()->Add(fTrackFilter); | |
79 | } | |
80 | ||
81 | void AliAnalysisTaskESDfilter::Init() | |
82 | { | |
83 | // Initialization | |
84 | if (fDebug > 1) AliInfo("Init() \n"); | |
85 | // Call configuration file | |
86 | } | |
87 | ||
88 | ||
89 | void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/) | |
90 | { | |
91 | // Execute analysis for current event | |
92 | // | |
93 | ||
94 | Long64_t ientry = Entry(); | |
2827d046 | 95 | if (fDebug > 0) printf("Filter: Analysing event # %5d\n", (Int_t) ientry); |
d7bdc804 | 96 | if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track"); |
97 | if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold"); | |
da97a08a | 98 | |
96b85d73 | 99 | ConvertESDtoAOD(); |
100 | } | |
101 | ||
102 | void AliAnalysisTaskESDfilter::ConvertESDtoAOD() { | |
103 | // ESD Filter analysis task executed for each event | |
da97a08a | 104 | |
96b85d73 | 105 | AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent()); |
9a244178 | 106 | AliESD* old = esd->GetAliESDOld(); |
da97a08a | 107 | |
108 | // Fetch Stack for debuggging if available | |
109 | AliStack *pStack = 0; | |
110 | AliMCEventHandler *mcH = 0; | |
111 | if(MCEvent()){ | |
112 | pStack = MCEvent()->Stack(); | |
113 | mcH = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); | |
114 | } | |
96b85d73 | 115 | // set arrays and pointers |
116 | Float_t posF[3]; | |
117 | Double_t pos[3]; | |
118 | Double_t p[3]; | |
9a244178 | 119 | Double_t momPos[3]; |
120 | Double_t momNeg[3]; | |
121 | Double_t momBach[3]; | |
122 | Double_t momPosAtV0vtx[3]; | |
123 | Double_t momNegAtV0vtx[3]; | |
124 | Double_t momBachAtCascadeVtx[3]; | |
96b85d73 | 125 | Double_t covVtx[6]; |
126 | Double_t covTr[21]; | |
127 | Double_t pid[10]; | |
128 | ||
129 | for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.; | |
130 | for (Int_t i = 0; i < 21; i++) covTr [i] = 0.; | |
131 | ||
132 | ||
da97a08a | 133 | // loop over events and fill them |
134 | ||
135 | // Multiplicity information needed by the header (to be revised!) | |
96b85d73 | 136 | Int_t nTracks = esd->GetNumberOfTracks(); |
8d69965b | 137 | // if (fDebug > 0) printf("-------------------Bo: Number of ESD tracks %d \n",nTracks); |
138 | ||
96b85d73 | 139 | Int_t nPosTracks = 0; |
3e492e99 | 140 | // for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack) |
141 | // if (esd->GetTrack(iTrack)->GetSign()> 0) nPosTracks++; | |
96b85d73 | 142 | |
143 | // Update the header | |
144 | ||
145 | AliAODHeader* header = AODEvent()->GetHeader(); | |
9a244178 | 146 | |
96b85d73 | 147 | header->SetRunNumber(esd->GetRunNumber()); |
148 | if (old) { | |
149 | header->SetBunchCrossNumber(0); | |
150 | header->SetOrbitNumber(0); | |
151 | header->SetPeriodNumber(0); | |
152 | header->SetEventType(0); | |
153 | header->SetMuonMagFieldScale(-999.); // FIXME | |
154 | header->SetCentrality(-999.); // FIXME | |
155 | } else { | |
156 | header->SetBunchCrossNumber(esd->GetBunchCrossNumber()); | |
157 | header->SetOrbitNumber(esd->GetOrbitNumber()); | |
158 | header->SetPeriodNumber(esd->GetPeriodNumber()); | |
159 | header->SetEventType(esd->GetEventType()); | |
160 | header->SetMuonMagFieldScale(-999.); // FIXME | |
161 | header->SetCentrality(-999.); // FIXME | |
162 | } | |
163 | ||
164 | header->SetTriggerMask(esd->GetTriggerMask()); | |
165 | header->SetTriggerCluster(esd->GetTriggerCluster()); | |
166 | header->SetMagneticField(esd->GetMagneticField()); | |
167 | header->SetZDCN1Energy(esd->GetZDCN1Energy()); | |
168 | header->SetZDCP1Energy(esd->GetZDCP1Energy()); | |
169 | header->SetZDCN2Energy(esd->GetZDCN2Energy()); | |
170 | header->SetZDCP2Energy(esd->GetZDCP2Energy()); | |
171 | header->SetZDCEMEnergy(esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1)); | |
613fc341 | 172 | Float_t diamxy[2]={esd->GetDiamondX(),esd->GetDiamondY()}; |
173 | Float_t diamcov[3]; esd->GetDiamondCovXY(diamcov); | |
174 | header->SetDiamond(diamxy,diamcov); | |
96b85d73 | 175 | // |
9a244178 | 176 | // |
96b85d73 | 177 | Int_t nV0s = esd->GetNumberOfV0s(); |
178 | Int_t nCascades = esd->GetNumberOfCascades(); | |
179 | Int_t nKinks = esd->GetNumberOfKinks(); | |
9a244178 | 180 | Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/; |
96b85d73 | 181 | Int_t nJets = 0; |
182 | Int_t nCaloClus = esd->GetNumberOfCaloClusters(); | |
183 | Int_t nFmdClus = 0; | |
184 | Int_t nPmdClus = esd->GetNumberOfPmdTracks(); | |
185 | ||
2827d046 | 186 | if (fDebug > 0) |
187 | printf(" NV0=%d NCASCADES=%d NKINKS=%d\n", nV0s, nCascades, nKinks); | |
9a244178 | 188 | |
189 | AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus); | |
96b85d73 | 190 | |
7dd90a3c | 191 | |
9a244178 | 192 | AliAODTrack *aodTrack = 0x0; |
193 | AliAODPid *detpid = 0x0; | |
194 | Double_t timezero = 0; //TO BE FIXED | |
195 | AliAODVertex *vV0FromCascade = 0x0; | |
196 | AliAODv0 *aodV0 = 0x0; | |
197 | AliAODcascade *aodCascade = 0x0; | |
198 | ||
d3fc6709 | 199 | // RefArray to store the mapping between esd track number and newly created AOD-Track |
9a244178 | 200 | TRefArray *aodTrackRefs = NULL; |
201 | if (nTracks > 0) aodTrackRefs = new TRefArray(nTracks); | |
202 | ||
203 | // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0 | |
204 | TRefArray *aodV0VtxRefs = NULL; | |
205 | if (nV0s > 0) aodV0VtxRefs = new TRefArray(nV0s); | |
206 | ||
207 | // RefArray to store the mapping between esd V0 number and newly created AOD-V0 | |
208 | TRefArray *aodV0Refs = NULL; | |
209 | if (nV0s > 0) aodV0Refs = new TRefArray(nV0s); | |
210 | ||
211 | ||
212 | ||
213 | ||
8d69965b | 214 | |
96b85d73 | 215 | // Array to take into account the tracks already added to the AOD |
216 | Bool_t * usedTrack = NULL; | |
217 | if (nTracks>0) { | |
218 | usedTrack = new Bool_t[nTracks]; | |
219 | for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE; | |
220 | } | |
9a244178 | 221 | // Array to take into account the V0s already added to the AOD (V0 within cascades) |
96b85d73 | 222 | Bool_t * usedV0 = NULL; |
223 | if (nV0s>0) { | |
224 | usedV0 = new Bool_t[nV0s]; | |
225 | for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE; | |
226 | } | |
227 | // Array to take into account the kinks already added to the AOD | |
228 | Bool_t * usedKink = NULL; | |
229 | if (nKinks>0) { | |
230 | usedKink = new Bool_t[nKinks]; | |
231 | for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE; | |
232 | } | |
233 | ||
234 | // Access to the AOD container of vertices | |
235 | TClonesArray &vertices = *(AODEvent()->GetVertices()); | |
236 | Int_t jVertices=0; | |
237 | ||
238 | // Access to the AOD container of tracks | |
239 | TClonesArray &tracks = *(AODEvent()->GetTracks()); | |
240 | Int_t jTracks=0; | |
241 | ||
9a244178 | 242 | // Access to the AOD container of Cascades |
243 | TClonesArray &cascades = *(AODEvent()->GetCascades()); | |
244 | Int_t jCascades=0; | |
245 | ||
96b85d73 | 246 | // Access to the AOD container of V0s |
9a244178 | 247 | TClonesArray &v0s = *(AODEvent()->GetV0s()); |
96b85d73 | 248 | Int_t jV0s=0; |
249 | ||
250 | // Add primary vertex. The primary tracks will be defined | |
251 | // after the loops on the composite objects (V0, cascades, kinks) | |
252 | const AliESDVertex *vtx = esd->GetPrimaryVertex(); | |
253 | ||
254 | vtx->GetXYZ(pos); // position | |
255 | vtx->GetCovMatrix(covVtx); //covariance matrix | |
256 | ||
257 | AliAODVertex * primary = new(vertices[jVertices++]) | |
c06d243f | 258 | AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary); |
f7361742 | 259 | primary->SetName(vtx->GetName()); |
260 | primary->SetTitle(vtx->GetTitle()); | |
261 | ||
2827d046 | 262 | if (fDebug > 0) primary->Print(); |
96b85d73 | 263 | |
264 | // Create vertices starting from the most complex objects | |
265 | Double_t chi2 = 0.; | |
266 | ||
9a244178 | 267 | // Cascades (Modified by A.Maire - February 2009) |
96b85d73 | 268 | for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) { |
96b85d73 | 269 | |
9a244178 | 270 | if (fDebug > 1) |
271 | printf("\n ******** Cascade number : %d/%d *********\n", nCascade, nCascades); | |
272 | ||
273 | // 0- Preparation | |
274 | // | |
275 | AliESDcascade *esdCascade = esd->GetCascade(nCascade); | |
276 | Int_t idxPosFromV0Dghter = esdCascade->GetPindex(); | |
277 | Int_t idxNegFromV0Dghter = esdCascade->GetNindex(); | |
278 | Int_t idxBachFromCascade = esdCascade->GetBindex(); | |
96b85d73 | 279 | |
9a244178 | 280 | AliESDtrack *esdCascadePos = esd->GetTrack( idxPosFromV0Dghter); |
281 | AliESDtrack *esdCascadeNeg = esd->GetTrack( idxNegFromV0Dghter); | |
282 | AliESDtrack *esdCascadeBach = esd->GetTrack( idxBachFromCascade); | |
283 | ||
284 | // Identification of the V0 within the esdCascade (via both daughter track indices) | |
285 | AliESDv0 * currentV0 = 0x0; | |
286 | Int_t idxV0FromCascade = -1; | |
96b85d73 | 287 | |
9a244178 | 288 | for (Int_t iV0=0; iV0<nV0s; ++iV0) { |
96b85d73 | 289 | |
9a244178 | 290 | currentV0 = esd->GetV0(iV0); |
291 | Int_t posCurrentV0 = currentV0->GetPindex(); | |
292 | Int_t negCurrentV0 = currentV0->GetNindex(); | |
96b85d73 | 293 | |
9a244178 | 294 | if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) { |
295 | idxV0FromCascade = iV0; | |
96b85d73 | 296 | break; |
9a244178 | 297 | } |
96b85d73 | 298 | } |
9a244178 | 299 | |
300 | if(idxV0FromCascade < 0){ | |
301 | printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n"); | |
302 | continue; | |
303 | }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds" | |
304 | ||
305 | if (fDebug > 1) | |
306 | printf("Cascade %d - V0fromCascade ind : %d/%d\n", nCascade, idxV0FromCascade, nV0s); | |
307 | ||
308 | AliESDv0 *esdV0FromCascade = esd->GetV0(idxV0FromCascade); | |
96b85d73 | 309 | |
9a244178 | 310 | |
311 | // 1 - Cascade selection | |
96b85d73 | 312 | |
9a244178 | 313 | // AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd->GetPrimaryVertex())); |
314 | // TList cascadeObjects; | |
315 | // cascadeObjects.AddAt(esdV0FromCascade, 0); | |
316 | // cascadeObjects.AddAt(esdCascadePos, 1); | |
317 | // cascadeObjects.AddAt(esdCascadeNeg, 2); | |
318 | // cascadeObjects.AddAt(esdCascade, 3); | |
319 | // cascadeObjects.AddAt(esdCascadeBach, 4); | |
320 | // cascadeObjects.AddAt(esdPrimVtx, 5); | |
321 | // | |
322 | // UInt_t selectCascade = 0; | |
323 | // if (fCascadeFilter) { | |
324 | // // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects); | |
325 | // // FIXME AliESDCascadeCuts to be implemented ... | |
326 | // | |
327 | // // Here we may encounter a moot point at the V0 level | |
328 | // // between the cascade selections and the V0 ones : | |
329 | // // the V0 selected along with the cascade (secondary V0) may | |
330 | // // usually be removed from the dedicated V0 selections (prim V0) ... | |
331 | // // -> To be discussed ! | |
332 | // | |
333 | // // this is a little awkward but otherwise the | |
334 | // // list wants to access the pointer (delete it) | |
335 | // // again when going out of scope | |
336 | // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct | |
337 | // esdPrimVtx = 0; | |
338 | // if (!selectCascade) | |
339 | // continue; | |
340 | // } | |
341 | // else{ | |
342 | // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct | |
343 | // esdPrimVtx = 0; | |
344 | // } | |
96b85d73 | 345 | |
9a244178 | 346 | // 2 - Add the cascade vertex |
347 | ||
348 | esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]); | |
349 | esdCascade->GetPosCovXi(covVtx); | |
350 | chi2 = esdCascade->GetChi2Xi(); | |
351 | ||
352 | AliAODVertex *vCascade = new(vertices[jVertices++]) AliAODVertex( pos, | |
353 | covVtx, | |
354 | chi2, // FIXME = Chi2/NDF will be needed | |
355 | primary, | |
356 | nCascade, // id | |
357 | AliAODVertex::kCascade); | |
358 | primary->AddDaughter(vCascade); | |
96b85d73 | 359 | |
9a244178 | 360 | if (fDebug > 2) { |
361 | printf("---- Cascade / Cascade Vertex (AOD) : \n"); | |
362 | vCascade->Print(); | |
96b85d73 | 363 | } |
9a244178 | 364 | |
365 | ||
366 | // 3 - Add the bachelor track from the cascade | |
96b85d73 | 367 | |
9a244178 | 368 | if (!usedTrack[idxBachFromCascade]) { |
369 | ||
370 | esdCascadeBach->GetPxPyPz(momBach); | |
371 | esdCascadeBach->GetXYZ(pos); | |
372 | esdCascadeBach->GetCovarianceXYZPxPyPz(covTr); | |
373 | esdCascadeBach->GetESDpid(pid); | |
374 | ||
375 | usedTrack[idxBachFromCascade] = kTRUE; | |
fae68d32 | 376 | UInt_t selectInfo = 0; |
9a244178 | 377 | if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach); |
378 | if(mcH)mcH->SelectParticle(esdCascadeBach->GetLabel()); | |
379 | aodTrack = new(tracks[jTracks++]) AliAODTrack(esdCascadeBach->GetID(), | |
380 | esdCascadeBach->GetLabel(), | |
381 | momBach, | |
382 | kTRUE, | |
383 | pos, | |
384 | kFALSE, // Why kFALSE for "isDCA" ? FIXME | |
385 | covTr, | |
386 | (Short_t)esdCascadeBach->GetSign(), | |
387 | esdCascadeBach->GetITSClusterMap(), | |
388 | pid, | |
389 | vCascade, | |
390 | kTRUE, // usedForVtxFit = kFALSE ? FIXME | |
391 | vtx->UsesTrack(esdCascadeBach->GetID()), | |
392 | AliAODTrack::kSecondary, | |
393 | selectInfo); | |
394 | aodTrackRefs->AddAt(aodTrack,idxBachFromCascade); | |
395 | ||
396 | if (esdCascadeBach->GetSign() > 0) nPosTracks++; | |
96b85d73 | 397 | aodTrack->ConvertAliPIDtoAODPID(); |
9a244178 | 398 | aodTrack->SetFlags(esdCascadeBach->GetStatus()); |
a401ba0b | 399 | SetAODPID(esdCascadeBach,aodTrack,detpid,timezero,esd->GetMagneticField()); |
96b85d73 | 400 | } |
401 | else { | |
9a244178 | 402 | aodTrack = dynamic_cast<AliAODTrack*>( aodTrackRefs->At(idxBachFromCascade) ); |
403 | } | |
404 | ||
405 | vCascade->AddDaughter(aodTrack); | |
406 | ||
407 | if (fDebug > 4) { | |
408 | printf("---- Cascade / bach dghter : \n"); | |
409 | aodTrack->Print(); | |
96b85d73 | 410 | } |
411 | ||
9a244178 | 412 | |
413 | // 4 - Add the V0 from the cascade. | |
414 | // = V0vtx + both pos and neg daughter tracks + the aodV0 itself | |
415 | // | |
416 | ||
417 | if ( !usedV0[idxV0FromCascade] ) { | |
418 | // 4.A - if VO structure hasn't been created yet | |
419 | ||
420 | // 4.A.1 - Create the V0 vertex of the cascade | |
96b85d73 | 421 | |
9a244178 | 422 | esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]); |
423 | esdV0FromCascade->GetPosCov(covVtx); | |
424 | chi2 = esdV0FromCascade->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 ? | |
425 | ||
426 | vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos, | |
427 | covVtx, | |
428 | chi2, | |
429 | vCascade, | |
430 | idxV0FromCascade, //id of ESDv0 | |
431 | AliAODVertex::kV0); | |
432 | // Note: | |
433 | // one V0 can be used by several cascades. | |
434 | // So, one AOD V0 vtx can have several parent vtx. | |
435 | // This is not directly allowed by AliAODvertex. | |
436 | // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash | |
437 | // but to a problem of consistency within AODEvent. | |
438 | // -> See below paragraph 4.B, for the proposed treatment of such a case. | |
439 | ||
440 | // Add the vV0FromCascade to the aodVOVtxRefs | |
441 | aodV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade); | |
442 | ||
443 | ||
444 | // 4.A.2 - Add the positive tracks from the V0 | |
445 | ||
446 | esdCascadePos->GetPxPyPz(momPos); | |
447 | esdCascadePos->GetXYZ(pos); | |
448 | esdCascadePos->GetCovarianceXYZPxPyPz(covTr); | |
449 | esdCascadePos->GetESDpid(pid); | |
450 | ||
451 | ||
452 | if (!usedTrack[idxPosFromV0Dghter]) { | |
453 | usedTrack[idxPosFromV0Dghter] = kTRUE; | |
96b85d73 | 454 | |
9a244178 | 455 | UInt_t selectInfo = 0; |
456 | if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos); | |
457 | if(mcH) mcH->SelectParticle(esdCascadePos->GetLabel()); | |
458 | aodTrack = new(tracks[jTracks++]) AliAODTrack( esdCascadePos->GetID(), | |
459 | esdCascadePos->GetLabel(), | |
460 | momPos, | |
461 | kTRUE, | |
462 | pos, | |
463 | kFALSE, // Why kFALSE for "isDCA" ? FIXME | |
464 | covTr, | |
465 | (Short_t)esdCascadePos->GetSign(), | |
466 | esdCascadePos->GetITSClusterMap(), | |
467 | pid, | |
468 | vV0FromCascade, | |
469 | kTRUE, // usedForVtxFit = kFALSE ? FIXME | |
470 | vtx->UsesTrack(esdCascadePos->GetID()), | |
471 | AliAODTrack::kSecondary, | |
472 | selectInfo); | |
473 | aodTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter); | |
474 | ||
475 | if (esdCascadePos->GetSign() > 0) nPosTracks++; | |
476 | aodTrack->ConvertAliPIDtoAODPID(); | |
477 | aodTrack->SetFlags(esdCascadePos->GetStatus()); | |
a401ba0b | 478 | SetAODPID(esdCascadePos,aodTrack,detpid,timezero,esd->GetMagneticField()); |
9a244178 | 479 | } |
480 | else { | |
481 | aodTrack = dynamic_cast<AliAODTrack*>(aodTrackRefs->At(idxPosFromV0Dghter)); | |
482 | } | |
483 | vV0FromCascade->AddDaughter(aodTrack); | |
484 | ||
485 | ||
486 | // 4.A.3 - Add the negative tracks from the V0 | |
487 | ||
488 | esdCascadeNeg->GetPxPyPz(momNeg); | |
489 | esdCascadeNeg->GetXYZ(pos); | |
490 | esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr); | |
491 | esdCascadeNeg->GetESDpid(pid); | |
492 | ||
493 | ||
494 | if (!usedTrack[idxNegFromV0Dghter]) { | |
495 | usedTrack[idxNegFromV0Dghter] = kTRUE; | |
496 | ||
497 | UInt_t selectInfo = 0; | |
498 | if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg); | |
499 | if(mcH)mcH->SelectParticle(esdCascadeNeg->GetLabel()); | |
500 | aodTrack = new(tracks[jTracks++]) AliAODTrack( esdCascadeNeg->GetID(), | |
501 | esdCascadeNeg->GetLabel(), | |
502 | momNeg, | |
503 | kTRUE, | |
504 | pos, | |
505 | kFALSE, // Why kFALSE for "isDCA" ? FIXME | |
506 | covTr, | |
507 | (Short_t)esdCascadeNeg->GetSign(), | |
508 | esdCascadeNeg->GetITSClusterMap(), | |
509 | pid, | |
510 | vV0FromCascade, | |
511 | kTRUE, // usedForVtxFit = kFALSE ? FIXME | |
512 | vtx->UsesTrack(esdCascadeNeg->GetID()), | |
513 | AliAODTrack::kSecondary, | |
514 | selectInfo); | |
515 | ||
516 | aodTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter); | |
517 | ||
518 | if (esdCascadeNeg->GetSign() > 0) nPosTracks++; | |
519 | aodTrack->ConvertAliPIDtoAODPID(); | |
520 | aodTrack->SetFlags(esdCascadeNeg->GetStatus()); | |
a401ba0b | 521 | SetAODPID(esdCascadeNeg,aodTrack,detpid,timezero,esd->GetMagneticField()); |
9a244178 | 522 | } |
523 | else { | |
524 | aodTrack = dynamic_cast<AliAODTrack*>(aodTrackRefs->At(idxNegFromV0Dghter)); | |
525 | } | |
526 | ||
527 | vV0FromCascade->AddDaughter(aodTrack); | |
528 | ||
529 | ||
530 | // 4.A.4 - Add the V0 from cascade to the V0 array | |
96b85d73 | 531 | |
9a244178 | 532 | Double_t dcaV0Daughters = esdV0FromCascade->GetDcaV0Daughters(); |
533 | Double_t dcaV0ToPrimVertex = esdV0FromCascade->GetD( esd->GetPrimaryVertex()->GetX(), | |
534 | esd->GetPrimaryVertex()->GetY(), | |
535 | esd->GetPrimaryVertex()->GetZ() ); | |
536 | esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] ); | |
537 | esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] ); | |
96b85d73 | 538 | |
9a244178 | 539 | Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg |
540 | dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD( esd->GetPrimaryVertex()->GetX(), | |
541 | esd->GetPrimaryVertex()->GetY(), | |
542 | esd->GetMagneticField()) ); | |
543 | dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD( esd->GetPrimaryVertex()->GetX(), | |
544 | esd->GetPrimaryVertex()->GetY(), | |
545 | esd->GetMagneticField()) ); | |
546 | ||
547 | aodV0 = new(v0s[jV0s++]) AliAODv0( vV0FromCascade, | |
548 | dcaV0Daughters, | |
549 | dcaV0ToPrimVertex, | |
550 | momPosAtV0vtx, | |
551 | momNegAtV0vtx, | |
552 | dcaDaughterToPrimVertex); | |
553 | // set the aod v0 on-the-fly status | |
554 | aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus()); | |
555 | ||
556 | // Add the aodV0 to the aodVORefs | |
557 | aodV0Refs->AddAt(aodV0,idxV0FromCascade); | |
558 | ||
559 | usedV0[idxV0FromCascade] = kTRUE; | |
560 | ||
561 | } else { | |
562 | // 4.B - if V0 structure already used | |
563 | ||
564 | // Note : | |
565 | // one V0 can be used by several cascades (frequent in PbPb evts) : | |
566 | // same V0 which used but attached to different bachelor tracks | |
567 | // -> aodVORefs and aodV0VtxRefs are needed. | |
568 | // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array. | |
96b85d73 | 569 | |
9a244178 | 570 | vV0FromCascade = dynamic_cast<AliAODVertex*>( aodV0VtxRefs->At(idxV0FromCascade) ); |
571 | aodV0 = dynamic_cast<AliAODv0*> ( aodV0Refs ->At(idxV0FromCascade) ); | |
572 | ||
573 | // - Treatment of the parent for such a "re-used" V0 : | |
574 | // Insert the cascade that reuses the V0 vertex in the lineage chain | |
575 | // Before : vV0 -> vCascade1 -> vPrimary | |
576 | // - Hyp : cascade2 uses the same V0 as cascade1 | |
577 | // After : vV0 -> vCascade2 -> vCascade1 -> vPrimary | |
578 | ||
579 | AliAODVertex *vCascadePreviousParent = dynamic_cast<AliAODVertex*> (vV0FromCascade->GetParent()); | |
580 | vV0FromCascade->SetParent(vCascade); | |
581 | vCascade ->SetParent(vCascadePreviousParent); | |
582 | ||
583 | if(fDebug > 2) | |
584 | printf("---- Cascade / Lineage insertion\n" | |
585 | "Parent of V0 vtx = Cascade vtx %p\n" | |
586 | "Parent of the cascade vtx = Cascade vtx %p\n" | |
587 | "Parent of the parent cascade vtx = Cascade vtx %p\n", | |
588 | static_cast<void*> (vV0FromCascade->GetParent()), | |
589 | static_cast<void*> (vCascade->GetParent()), | |
590 | static_cast<void*> (vCascadePreviousParent->GetParent()) ); | |
591 | ||
592 | }// end if V0 structure already used | |
593 | ||
594 | if (fDebug > 2) { | |
595 | printf("---- Cascade / V0 vertex: \n"); | |
596 | vV0FromCascade->Print(); | |
96b85d73 | 597 | } |
9a244178 | 598 | |
599 | if (fDebug > 4) { | |
600 | printf("---- Cascade / pos dghter : \n"); | |
601 | aodTrack->Print(); | |
602 | printf("---- Cascade / neg dghter : \n"); | |
603 | aodTrack->Print(); | |
604 | printf("---- Cascade / aodV0 : \n"); | |
605 | aodV0->Print(); | |
96b85d73 | 606 | } |
9a244178 | 607 | |
608 | // In any case (used V0 or not), add the V0 vertex to the cascade one. | |
609 | vCascade->AddDaughter(vV0FromCascade); | |
96b85d73 | 610 | |
9a244178 | 611 | |
612 | // 5 - Add the primary track of the cascade (if any) | |
613 | ||
614 | ||
615 | // 6 - Add the cascade to the AOD array of cascades | |
616 | ||
617 | Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd->GetPrimaryVertex()->GetX(), | |
618 | esd->GetPrimaryVertex()->GetY(), | |
619 | esd->GetMagneticField()) ); | |
620 | ||
28740a50 | 621 | esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]); |
9a244178 | 622 | |
623 | aodCascade = new(cascades[jCascades++]) AliAODcascade( vCascade, | |
624 | esdCascade->Charge(), | |
625 | esdCascade->GetDcaXiDaughters(), | |
626 | -999., | |
627 | // DCAXiToPrimVtx -> needs to be calculated ----| | |
628 | // doesn't exist at ESD level; | |
629 | // See AODcascade::DcaXiToPrimVertex(Double, Double, Double) | |
630 | dcaBachToPrimVertexXY, | |
631 | momBachAtCascadeVtx, | |
632 | *aodV0); | |
96b85d73 | 633 | |
9a244178 | 634 | if (fDebug > 3) { |
635 | printf("---- Cascade / AOD cascade : \n\n"); | |
636 | aodCascade->PrintXi(primary->GetX(), primary->GetY(), primary->GetZ()); | |
637 | } | |
638 | ||
96b85d73 | 639 | } // end of the loop on cascades |
ec6465fe | 640 | |
9a244178 | 641 | cascades.Expand(jCascades); |
642 | ||
643 | ||
ec6465fe | 644 | // |
96b85d73 | 645 | // V0s |
ec6465fe | 646 | // |
96b85d73 | 647 | |
648 | for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) { | |
649 | ||
7dd90a3c | 650 | if (usedV0[nV0]) continue; // skip if already added to the AOD |
96b85d73 | 651 | |
652 | AliESDv0 *v0 = esd->GetV0(nV0); | |
ec6465fe | 653 | Int_t posFromV0 = v0->GetPindex(); |
654 | Int_t negFromV0 = v0->GetNindex(); | |
9a244178 | 655 | |
ec6465fe | 656 | // V0 selection |
657 | // | |
9a244178 | 658 | AliESDVertex *esdVtx = new AliESDVertex(*(esd->GetPrimaryVertex())); |
659 | AliESDtrack *esdV0Pos = esd->GetTrack(posFromV0); | |
660 | AliESDtrack *esdV0Neg = esd->GetTrack(negFromV0); | |
c12f1b16 | 661 | TList v0objects; |
ec6465fe | 662 | v0objects.AddAt(v0, 0); |
663 | v0objects.AddAt(esdV0Pos, 1); | |
664 | v0objects.AddAt(esdV0Neg, 2); | |
665 | v0objects.AddAt(esdVtx, 3); | |
666 | UInt_t selectV0 = 0; | |
667 | if (fV0Filter) { | |
396d2b33 | 668 | selectV0 = fV0Filter->IsSelected(&v0objects); |
c12f1b16 | 669 | // this is a little awkward but otherwise the |
70186db5 | 670 | // list wants to access the pointer (delete it) |
671 | // again when going out of scope | |
672 | delete v0objects.RemoveAt(3); // esdVtx created via copy construct | |
673 | esdVtx = 0; | |
396d2b33 | 674 | if (!selectV0) |
675 | continue; | |
676 | } | |
677 | else{ | |
70186db5 | 678 | delete v0objects.RemoveAt(3); // esdVtx created via copy construct |
679 | esdVtx = 0; | |
ec6465fe | 680 | } |
ec6465fe | 681 | |
96b85d73 | 682 | v0->GetXYZ(pos[0], pos[1], pos[2]); |
683 | ||
684 | if (!old) { | |
685 | chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 | |
686 | v0->GetPosCov(covVtx); | |
687 | } else { | |
688 | chi2 = -999.; | |
689 | for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.; | |
690 | } | |
691 | ||
692 | ||
693 | AliAODVertex * vV0 = | |
694 | new(vertices[jVertices++]) AliAODVertex(pos, | |
695 | covVtx, | |
696 | chi2, | |
697 | primary, | |
698 | nV0, | |
699 | AliAODVertex::kV0); | |
700 | primary->AddDaughter(vV0); | |
701 | ||
7dd90a3c | 702 | |
96b85d73 | 703 | // Add the positive tracks from the V0 |
704 | ||
9a244178 | 705 | esdV0Pos->GetPxPyPz(momPos); |
ec6465fe | 706 | esdV0Pos->GetXYZ(pos); |
707 | esdV0Pos->GetCovarianceXYZPxPyPz(covTr); | |
708 | esdV0Pos->GetESDpid(pid); | |
9a244178 | 709 | |
ec6465fe | 710 | if (!usedTrack[posFromV0]) { |
96b85d73 | 711 | usedTrack[posFromV0] = kTRUE; |
fae68d32 | 712 | UInt_t selectInfo = 0; |
ec6465fe | 713 | if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos); |
da97a08a | 714 | if(mcH)mcH->SelectParticle(esdV0Pos->GetLabel()); |
ec6465fe | 715 | aodTrack = new(tracks[jTracks++]) AliAODTrack(esdV0Pos->GetID(), |
716 | esdV0Pos->GetLabel(), | |
9a244178 | 717 | momPos, |
8d69965b | 718 | kTRUE, |
719 | pos, | |
720 | kFALSE, | |
721 | covTr, | |
ec6465fe | 722 | (Short_t)esdV0Pos->GetSign(), |
723 | esdV0Pos->GetITSClusterMap(), | |
8d69965b | 724 | pid, |
725 | vV0, | |
726 | kTRUE, // check if this is right | |
f7361742 | 727 | vtx->UsesTrack(esdV0Pos->GetID()), |
3de26c63 | 728 | AliAODTrack::kSecondary, |
729 | selectInfo); | |
9a244178 | 730 | aodTrackRefs->AddAt(aodTrack,posFromV0); |
8d69965b | 731 | // if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt()); |
ec6465fe | 732 | if (esdV0Pos->GetSign() > 0) nPosTracks++; |
96b85d73 | 733 | aodTrack->ConvertAliPIDtoAODPID(); |
ec6465fe | 734 | aodTrack->SetFlags(esdV0Pos->GetStatus()); |
a401ba0b | 735 | SetAODPID(esdV0Pos,aodTrack,detpid,timezero,esd->GetMagneticField()); |
ec6465fe | 736 | } |
737 | else { | |
9a244178 | 738 | aodTrack = dynamic_cast<AliAODTrack*>(aodTrackRefs->At(posFromV0)); |
8d69965b | 739 | // if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt()); |
96b85d73 | 740 | } |
ec6465fe | 741 | vV0->AddDaughter(aodTrack); |
742 | ||
96b85d73 | 743 | // Add the negative tracks from the V0 |
744 | ||
9a244178 | 745 | esdV0Neg->GetPxPyPz(momNeg); |
ec6465fe | 746 | esdV0Neg->GetXYZ(pos); |
747 | esdV0Neg->GetCovarianceXYZPxPyPz(covTr); | |
748 | esdV0Neg->GetESDpid(pid); | |
ec6465fe | 749 | |
750 | if (!usedTrack[negFromV0]) { | |
96b85d73 | 751 | usedTrack[negFromV0] = kTRUE; |
fae68d32 | 752 | UInt_t selectInfo = 0; |
ec6465fe | 753 | if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg); |
da97a08a | 754 | if(mcH)mcH->SelectParticle(esdV0Neg->GetLabel()); |
ec6465fe | 755 | aodTrack = new(tracks[jTracks++]) AliAODTrack(esdV0Neg->GetID(), |
756 | esdV0Neg->GetLabel(), | |
9a244178 | 757 | momNeg, |
8d69965b | 758 | kTRUE, |
759 | pos, | |
760 | kFALSE, | |
761 | covTr, | |
ec6465fe | 762 | (Short_t)esdV0Neg->GetSign(), |
763 | esdV0Neg->GetITSClusterMap(), | |
8d69965b | 764 | pid, |
765 | vV0, | |
766 | kTRUE, // check if this is right | |
f7361742 | 767 | vtx->UsesTrack(esdV0Neg->GetID()), |
3de26c63 | 768 | AliAODTrack::kSecondary, |
769 | selectInfo); | |
ec6465fe | 770 | |
9a244178 | 771 | aodTrackRefs->AddAt(aodTrack,negFromV0); |
8d69965b | 772 | // if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt()); |
ec6465fe | 773 | if (esdV0Neg->GetSign() > 0) nPosTracks++; |
96b85d73 | 774 | aodTrack->ConvertAliPIDtoAODPID(); |
ec6465fe | 775 | aodTrack->SetFlags(esdV0Neg->GetStatus()); |
a401ba0b | 776 | SetAODPID(esdV0Neg,aodTrack,detpid,timezero,esd->GetMagneticField()); |
ec6465fe | 777 | } |
778 | else { | |
9a244178 | 779 | aodTrack = dynamic_cast<AliAODTrack*>(aodTrackRefs->At(negFromV0)); |
8d69965b | 780 | // if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt()); |
96b85d73 | 781 | } |
ec6465fe | 782 | vV0->AddDaughter(aodTrack); |
9a244178 | 783 | |
784 | ||
785 | // Add the V0 the V0 array as well | |
786 | ||
787 | Double_t dcaV0Daughters = v0->GetDcaV0Daughters(); | |
788 | Double_t dcaV0ToPrimVertex = v0->GetD(esd->GetPrimaryVertex()->GetX(), | |
789 | esd->GetPrimaryVertex()->GetY(), | |
790 | esd->GetPrimaryVertex()->GetZ()); | |
791 | v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]); | |
792 | v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]); | |
793 | ||
794 | Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg | |
795 | dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD( esd->GetPrimaryVertex()->GetX(), | |
796 | esd->GetPrimaryVertex()->GetY(), | |
797 | esd->GetMagneticField()) ); | |
798 | dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD( esd->GetPrimaryVertex()->GetX(), | |
799 | esd->GetPrimaryVertex()->GetY(), | |
800 | esd->GetMagneticField()) ); | |
801 | ||
802 | aodV0 = new(v0s[jV0s++]) AliAODv0(vV0, | |
803 | dcaV0Daughters, | |
804 | dcaV0ToPrimVertex, | |
805 | momPosAtV0vtx, | |
806 | momNegAtV0vtx, | |
807 | dcaDaughterToPrimVertex); | |
808 | ||
7dd90a3c | 809 | // set the aod v0 on-the-fly status |
810 | aodV0->SetOnFlyStatus(v0->GetOnFlyStatus()); | |
9a244178 | 811 | }//End of loop on V0s |
812 | ||
813 | v0s.Expand(jV0s); | |
814 | ||
815 | if (fDebug > 0) printf(" NAODCascades=%d / NAODV0s=%d\n", jCascades, jV0s); | |
816 | // end of V0 parts | |
817 | ||
818 | ||
96b85d73 | 819 | // Kinks: it is a big mess the access to the information in the kinks |
820 | // The loop is on the tracks in order to find the mother and daugther of each kink | |
821 | ||
822 | ||
823 | for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) { | |
824 | ||
825 | AliESDtrack * esdTrack = esd->GetTrack(iTrack); | |
826 | ||
827 | Int_t ikink = esdTrack->GetKinkIndex(0); | |
828 | ||
829 | if (ikink && nKinks) { | |
830 | // Negative kink index: mother, positive: daughter | |
831 | ||
832 | // Search for the second track of the kink | |
833 | ||
834 | for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) { | |
835 | ||
836 | AliESDtrack * esdTrack1 = esd->GetTrack(jTrack); | |
837 | ||
838 | Int_t jkink = esdTrack1->GetKinkIndex(0); | |
839 | ||
840 | if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) { | |
841 | ||
842 | // The two tracks are from the same kink | |
843 | ||
844 | if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks | |
845 | ||
846 | Int_t imother = -1; | |
847 | Int_t idaughter = -1; | |
848 | ||
849 | if (ikink<0 && jkink>0) { | |
850 | ||
851 | imother = iTrack; | |
852 | idaughter = jTrack; | |
853 | } | |
854 | else if (ikink>0 && jkink<0) { | |
855 | ||
856 | imother = jTrack; | |
857 | idaughter = iTrack; | |
858 | } | |
859 | else { | |
860 | // cerr << "Error: Wrong combination of kink indexes: " | |
861 | // << ikink << " " << jkink << endl; | |
862 | continue; | |
863 | } | |
864 | ||
5dddb02b | 865 | // Add the mother track if it passed primary track selection cuts |
96b85d73 | 866 | |
867 | AliAODTrack * mother = NULL; | |
868 | ||
fae68d32 | 869 | UInt_t selectInfo = 0; |
870 | if (fTrackFilter) { | |
871 | selectInfo = fTrackFilter->IsSelected(esd->GetTrack(imother)); | |
872 | if (!selectInfo) continue; | |
873 | } | |
5dddb02b | 874 | |
96b85d73 | 875 | if (!usedTrack[imother]) { |
876 | ||
877 | usedTrack[imother] = kTRUE; | |
878 | ||
03a5cc9f | 879 | AliESDtrack *esdTrackM = esd->GetTrack(imother); |
880 | esdTrackM->GetPxPyPz(p); | |
881 | esdTrackM->GetXYZ(pos); | |
882 | esdTrackM->GetCovarianceXYZPxPyPz(covTr); | |
883 | esdTrackM->GetESDpid(pid); | |
da97a08a | 884 | if(mcH)mcH->SelectParticle(esdTrackM->GetLabel()); |
96b85d73 | 885 | mother = |
03a5cc9f | 886 | new(tracks[jTracks++]) AliAODTrack(esdTrackM->GetID(), |
887 | esdTrackM->GetLabel(), | |
96b85d73 | 888 | p, |
889 | kTRUE, | |
890 | pos, | |
891 | kFALSE, | |
892 | covTr, | |
03a5cc9f | 893 | (Short_t)esdTrackM->GetSign(), |
894 | esdTrackM->GetITSClusterMap(), | |
96b85d73 | 895 | pid, |
896 | primary, | |
897 | kTRUE, // check if this is right | |
f7361742 | 898 | vtx->UsesTrack(esdTrack->GetID()), |
3de26c63 | 899 | AliAODTrack::kPrimary, |
900 | selectInfo); | |
9a244178 | 901 | aodTrackRefs->AddAt(mother, imother); |
d3fc6709 | 902 | |
03a5cc9f | 903 | if (esdTrackM->GetSign() > 0) nPosTracks++; |
904 | mother->SetFlags(esdTrackM->GetStatus()); | |
635ef56c | 905 | mother->ConvertAliPIDtoAODPID(); |
96b85d73 | 906 | primary->AddDaughter(mother); |
907 | mother->ConvertAliPIDtoAODPID(); | |
a401ba0b | 908 | SetAODPID(esdTrackM,mother,detpid,timezero,esd->GetMagneticField()); |
96b85d73 | 909 | } |
910 | else { | |
911 | // cerr << "Error: event " << esd->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1 | |
912 | // << " track " << imother << " has already been used!" << endl; | |
913 | } | |
914 | ||
915 | // Add the kink vertex | |
916 | AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1); | |
917 | ||
918 | AliAODVertex * vkink = | |
919 | new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(), | |
920 | NULL, | |
921 | 0., | |
922 | mother, | |
923 | esdTrack->GetID(), // This is the track ID of the mother's track! | |
924 | AliAODVertex::kKink); | |
925 | // Add the daughter track | |
926 | ||
927 | AliAODTrack * daughter = NULL; | |
928 | ||
929 | if (!usedTrack[idaughter]) { | |
930 | ||
931 | usedTrack[idaughter] = kTRUE; | |
932 | ||
03a5cc9f | 933 | AliESDtrack *esdTrackD = esd->GetTrack(idaughter); |
934 | esdTrackD->GetPxPyPz(p); | |
935 | esdTrackD->GetXYZ(pos); | |
936 | esdTrackD->GetCovarianceXYZPxPyPz(covTr); | |
937 | esdTrackD->GetESDpid(pid); | |
c12f1b16 | 938 | selectInfo = 0; |
fae68d32 | 939 | if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD); |
da97a08a | 940 | if(mcH)mcH->SelectParticle(esdTrackD->GetLabel()); |
96b85d73 | 941 | daughter = |
03a5cc9f | 942 | new(tracks[jTracks++]) AliAODTrack(esdTrackD->GetID(), |
943 | esdTrackD->GetLabel(), | |
96b85d73 | 944 | p, |
945 | kTRUE, | |
946 | pos, | |
947 | kFALSE, | |
948 | covTr, | |
03a5cc9f | 949 | (Short_t)esdTrackD->GetSign(), |
950 | esdTrackD->GetITSClusterMap(), | |
96b85d73 | 951 | pid, |
952 | vkink, | |
953 | kTRUE, // check if this is right | |
f7361742 | 954 | vtx->UsesTrack(esdTrack->GetID()), |
3de26c63 | 955 | AliAODTrack::kSecondary, |
956 | selectInfo); | |
d3fc6709 | 957 | |
9a244178 | 958 | aodTrackRefs->AddAt(daughter, idaughter); |
d3fc6709 | 959 | |
03a5cc9f | 960 | if (esdTrackD->GetSign() > 0) nPosTracks++; |
961 | daughter->SetFlags(esdTrackD->GetStatus()); | |
635ef56c | 962 | daughter->ConvertAliPIDtoAODPID(); |
96b85d73 | 963 | vkink->AddDaughter(daughter); |
964 | daughter->ConvertAliPIDtoAODPID(); | |
a401ba0b | 965 | SetAODPID(esdTrackD,daughter,detpid,timezero,esd->GetMagneticField()); |
96b85d73 | 966 | } |
967 | else { | |
968 | // cerr << "Error: event " << esd->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1 | |
969 | // << " track " << idaughter << " has already been used!" << endl; | |
970 | } | |
971 | } | |
972 | } | |
973 | } | |
974 | } | |
975 | ||
976 | ||
977 | // Tracks (primary and orphan) | |
978 | ||
3e492e99 | 979 | if (fDebug > 0) printf("NUMBER OF ESD TRACKS %5d\n", nTracks); |
96b85d73 | 980 | |
981 | for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) { | |
982 | ||
983 | ||
984 | if (usedTrack[nTrack]) continue; | |
985 | ||
986 | AliESDtrack *esdTrack = esd->GetTrack(nTrack); | |
987 | UInt_t selectInfo = 0; | |
988 | // | |
989 | // Track selection | |
990 | if (fTrackFilter) { | |
991 | selectInfo = fTrackFilter->IsSelected(esdTrack); | |
f7361742 | 992 | if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue; |
96b85d73 | 993 | } |
994 | ||
995 | // | |
996 | esdTrack->GetPxPyPz(p); | |
997 | esdTrack->GetXYZ(pos); | |
998 | esdTrack->GetCovarianceXYZPxPyPz(covTr); | |
999 | esdTrack->GetESDpid(pid); | |
da97a08a | 1000 | if(mcH)mcH->SelectParticle(esdTrack->GetLabel()); |
fae68d32 | 1001 | primary->AddDaughter(aodTrack = |
1002 | new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(), | |
1003 | esdTrack->GetLabel(), | |
1004 | p, | |
1005 | kTRUE, | |
1006 | pos, | |
1007 | kFALSE, | |
1008 | covTr, | |
1009 | (Short_t)esdTrack->GetSign(), | |
1010 | esdTrack->GetITSClusterMap(), | |
1011 | pid, | |
1012 | primary, | |
1013 | kTRUE, // check if this is right | |
f7361742 | 1014 | vtx->UsesTrack(esdTrack->GetID()), |
fae68d32 | 1015 | AliAODTrack::kPrimary, |
1016 | selectInfo) | |
da97a08a | 1017 | ); |
9a244178 | 1018 | aodTrackRefs->AddAt(aodTrack, nTrack); |
d3fc6709 | 1019 | |
fae68d32 | 1020 | if (esdTrack->GetSign() > 0) nPosTracks++; |
1021 | aodTrack->SetFlags(esdTrack->GetStatus()); | |
1022 | aodTrack->ConvertAliPIDtoAODPID(); | |
a401ba0b | 1023 | SetAODPID(esdTrack,aodTrack,detpid,timezero,esd->GetMagneticField()); |
96b85d73 | 1024 | } // end of loop on tracks |
1025 | ||
3e492e99 | 1026 | // Update number of AOD tracks in header at the end of track loop (M.G.) |
1027 | header->SetRefMultiplicity(jTracks); | |
1028 | header->SetRefMultiplicityPos(nPosTracks); | |
1029 | header->SetRefMultiplicityNeg(jTracks - nPosTracks); | |
1030 | if (fDebug > 0) | |
1031 | printf(" NAODTRACKS=%d NPOS=%d NNEG=%d\n", jTracks, nPosTracks, jTracks - nPosTracks); | |
1032 | // Do not shrink the array of tracks - other filters may add to it (M.G) | |
9a244178 | 1033 | // tracks.Expand(jTracks); // remove 'empty slots' due to unwritten tracks |
96b85d73 | 1034 | |
1035 | // Access to the AOD container of PMD clusters | |
1036 | TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters()); | |
1037 | Int_t jPmdClusters=0; | |
1038 | ||
1039 | for (Int_t iPmd = 0; iPmd < nPmdClus; ++iPmd) { | |
1040 | // file pmd clusters, to be revised! | |
1041 | AliESDPmdTrack *pmdTrack = esd->GetPmdTrack(iPmd); | |
1042 | Int_t nLabel = 0; | |
1043 | Int_t *label = 0x0; | |
03a5cc9f | 1044 | Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()}; |
73d219bc | 1045 | Double_t pidPmd[13] = { 0.}; // to be revised! |
96b85d73 | 1046 | // type not set! |
1047 | // assoc cluster not set | |
03a5cc9f | 1048 | new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd); |
96b85d73 | 1049 | } |
1050 | ||
1051 | // Access to the AOD container of clusters | |
1052 | TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters()); | |
1053 | Int_t jClusters=0; | |
1054 | ||
1055 | for (Int_t iClust=0; iClust<nCaloClus; ++iClust) { | |
1056 | ||
1057 | AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust); | |
1058 | ||
e6500713 | 1059 | Int_t id = cluster->GetID(); |
1060 | Int_t nLabel = cluster->GetNLabels(); | |
1061 | TArrayI* labels = cluster->GetLabels(); | |
1062 | Int_t *label = 0; | |
da97a08a | 1063 | if (labels){ |
1064 | label = (cluster->GetLabels())->GetArray(); | |
1065 | for(int i = 0;i < labels->GetSize();++i){ | |
1066 | if(mcH)mcH->SelectParticle(label[i]); | |
1067 | } | |
1068 | } | |
e6500713 | 1069 | |
96b85d73 | 1070 | Float_t energy = cluster->E(); |
1071 | cluster->GetPosition(posF); | |
e6500713 | 1072 | Char_t ttype = AliAODCluster::kUndef; |
96b85d73 | 1073 | |
1074 | if (cluster->GetClusterType() == AliESDCaloCluster::kPHOSCluster) { | |
1075 | ttype=AliAODCluster::kPHOSNeutral; | |
1076 | } | |
1077 | else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALClusterv1) { | |
1078 | ttype = AliAODCluster::kEMCALClusterv1; | |
1079 | } | |
1080 | ||
1081 | ||
1082 | AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id, | |
1083 | nLabel, | |
1084 | label, | |
1085 | energy, | |
636555d2 | 1086 | posF, |
96b85d73 | 1087 | NULL, |
1088 | ttype); | |
1089 | ||
e6500713 | 1090 | caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(), |
1091 | cluster->GetClusterDisp(), | |
78902954 | 1092 | cluster->GetM20(), cluster->GetM02(), |
1093 | cluster->GetEmcCpvDistance(), | |
1094 | cluster->GetNExMax(),cluster->GetTOF()) ; | |
e6500713 | 1095 | |
3cbe9a4d | 1096 | caloCluster->SetPIDFromESD(cluster->GetPid()); |
e6500713 | 1097 | caloCluster->SetNCells(cluster->GetNCells()); |
1098 | caloCluster->SetCellsAbsId(cluster->GetCellsAbsId()); | |
1099 | caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction()); | |
78902954 | 1100 | |
1101 | TArrayI* matchedT = cluster->GetTracksMatched(); | |
8097629b | 1102 | if (nTracks>0 && matchedT && cluster->GetTrackMatched() >= 0) { |
78902954 | 1103 | for (Int_t im = 0; im < matchedT->GetSize(); im++) { |
d3fc6709 | 1104 | Int_t iESDtrack = matchedT->At(im);; |
9a244178 | 1105 | if (aodTrackRefs->At(iESDtrack) != 0) { |
1106 | caloCluster->AddTrackMatched((AliAODTrack*)aodTrackRefs->At(iESDtrack)); | |
d3fc6709 | 1107 | } |
78902954 | 1108 | } |
1109 | } | |
d3fc6709 | 1110 | |
96b85d73 | 1111 | } |
1112 | caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters | |
1113 | // end of loop on calo clusters | |
1114 | ||
1115 | // fill EMCAL cell info | |
1116 | if (esd->GetEMCALCells()) { // protection against missing ESD information | |
1117 | AliESDCaloCells &esdEMcells = *(esd->GetEMCALCells()); | |
1118 | Int_t nEMcell = esdEMcells.GetNumberOfCells() ; | |
1119 | ||
1120 | AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells()); | |
1121 | aodEMcells.CreateContainer(nEMcell); | |
1122 | aodEMcells.SetType(AliAODCaloCells::kEMCAL); | |
1123 | for (Int_t iCell = 0; iCell < nEMcell; iCell++) { | |
1124 | aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell)); | |
1125 | } | |
1126 | aodEMcells.Sort(); | |
1127 | } | |
1128 | ||
1129 | // fill PHOS cell info | |
1130 | if (esd->GetPHOSCells()) { // protection against missing ESD information | |
1131 | AliESDCaloCells &esdPHcells = *(esd->GetPHOSCells()); | |
1132 | Int_t nPHcell = esdPHcells.GetNumberOfCells() ; | |
1133 | ||
1134 | AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells()); | |
1135 | aodPHcells.CreateContainer(nPHcell); | |
1136 | aodPHcells.SetType(AliAODCaloCells::kPHOS); | |
1137 | for (Int_t iCell = 0; iCell < nPHcell; iCell++) { | |
1138 | aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell)); | |
1139 | } | |
1140 | aodPHcells.Sort(); | |
1141 | } | |
1142 | ||
1143 | // tracklets | |
1144 | AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets()); | |
1145 | const AliMultiplicity *mult = esd->GetMultiplicity(); | |
1146 | if (mult) { | |
1147 | if (mult->GetNumberOfTracklets()>0) { | |
1148 | SPDTracklets.CreateContainer(mult->GetNumberOfTracklets()); | |
1149 | ||
1150 | for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) { | |
da97a08a | 1151 | if(mcH){ |
1152 | mcH->SelectParticle(mult->GetLabel(n, 0)); | |
1153 | mcH->SelectParticle(mult->GetLabel(n, 1)); | |
1154 | } | |
1155 | SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1)); | |
96b85d73 | 1156 | } |
1157 | } | |
1158 | } else { | |
1159 | //Printf("ERROR: AliMultiplicity could not be retrieved from ESD"); | |
1160 | } | |
1161 | ||
96b85d73 | 1162 | delete [] usedKink; |
9a244178 | 1163 | delete [] usedV0; |
1164 | delete [] usedTrack; | |
1165 | delete aodV0Refs; | |
1166 | delete aodV0VtxRefs; | |
1167 | delete aodTrackRefs; | |
96b85d73 | 1168 | |
1169 | return; | |
1170 | } | |
1171 | ||
7b4cee8d | 1172 | |
a401ba0b | 1173 | void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid, Double_t timezero, Double_t bfield) |
d7bdc804 | 1174 | { |
1175 | // | |
1176 | // Setter for the raw PID detector signals | |
1177 | // | |
1178 | ||
1179 | if(esdtrack->Pt()>fHighPthreshold) { | |
1180 | detpid = new AliAODPid(); | |
a401ba0b | 1181 | SetDetectorRawSignals(detpid,esdtrack,timezero, bfield); |
d7bdc804 | 1182 | aodtrack->SetDetPID(detpid); |
1183 | } else { | |
1184 | if(fPtshape){ | |
1185 | if(esdtrack->Pt()> fPtshape->GetXmin()){ | |
1186 | Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold); | |
1187 | if(gRandom->Rndm(0)<1./y){ | |
1188 | detpid = new AliAODPid(); | |
a401ba0b | 1189 | SetDetectorRawSignals(detpid,esdtrack,timezero, bfield); |
d7bdc804 | 1190 | aodtrack->SetDetPID(detpid); |
1191 | }//end rndm | |
1192 | }//end if p < pmin | |
1193 | }//end if p function | |
1194 | }// end else | |
1195 | } | |
1196 | ||
a401ba0b | 1197 | void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track, Double_t timezero, Double_t bfield) |
373fc041 | 1198 | { |
1199 | // | |
1200 | //assignment of the detector signals (AliXXXesdPID inspired) | |
1201 | // | |
1202 | if(!track){ | |
1203 | AliInfo("no ESD track found. .....exiting"); | |
1204 | return; | |
1205 | } | |
1206 | ||
1207 | aodpid->SetITSsignal(track->GetITSsignal()); | |
1208 | aodpid->SetTPCsignal(track->GetTPCsignal()); | |
1209 | //n TRD planes = 6 | |
1210 | ||
1211 | Int_t nslices = track->GetNumberOfTRDslices()*6; | |
1212 | Double_t *trdslices = new Double_t[nslices]; | |
1213 | for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) { | |
1214 | for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl); | |
1215 | } | |
d7bdc804 | 1216 | |
373fc041 | 1217 | |
1218 | aodpid->SetTRDsignal(track->GetNumberOfTRDslices()*6,trdslices); | |
1219 | Double_t times[AliAODPid::kSPECIES]; track->GetIntegratedTimes(times); | |
1220 | aodpid->SetIntegratedTimes(times); | |
1221 | ||
1222 | aodpid->SetTOFsignal(track->GetTOFsignal()-timezero); // to be fixed | |
1223 | aodpid->SetHMPIDsignal(track->GetHMPIDsignal()); | |
1224 | ||
a401ba0b | 1225 | //Extrapolate track to EMCAL surface for AOD-level track-cluster matching |
1226 | Double_t emcpos[3] = {0.,0.,0.}; | |
1227 | Double_t emcmom[3] = {0.,0.,0.}; | |
1228 | aodpid->SetEMCALPosition(emcpos); | |
1229 | aodpid->SetEMCALMomentum(emcmom); | |
1230 | ||
1231 | AliExternalTrackParam *outerparam = (AliExternalTrackParam*)track->GetOuterParam(); | |
1232 | if(!outerparam) return; | |
1233 | ||
1234 | //To be replaced by call to AliEMCALGeoUtils when the class becomes available | |
1235 | Double_t radius = 441.0; //[cm] EMCAL radius +13cm | |
1236 | ||
1237 | Bool_t okpos = outerparam->GetXYZAt(radius,bfield,emcpos); | |
1238 | Bool_t okmom = outerparam->GetPxPyPzAt(radius,bfield,emcmom); | |
1239 | if(!(okpos && okmom)) return; | |
1240 | ||
1241 | aodpid->SetEMCALPosition(emcpos); | |
1242 | aodpid->SetEMCALMomentum(emcmom); | |
1243 | ||
373fc041 | 1244 | } |
7b4cee8d | 1245 | |
96b85d73 | 1246 | void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/) |
1247 | { | |
1248 | // Terminate analysis | |
1249 | // | |
1250 | if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n"); | |
1251 | } | |
1252 | ||
da97a08a | 1253 | void AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){ |
1254 | if(!pStack)return; | |
1255 | label = TMath::Abs(label); | |
1256 | TParticle *part = pStack->Particle(label); | |
1257 | Printf("########################"); | |
1258 | Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P()); | |
1259 | part->Print(); | |
1260 | TParticle* mother = part; | |
1261 | Int_t imo = part->GetFirstMother(); | |
1262 | Int_t nprim = pStack->GetNprimary(); | |
1263 | // while((imo >= nprim) && (mother->GetUniqueID() == 4)) { | |
1264 | while((imo >= nprim)) { | |
1265 | mother = pStack->Particle(imo); | |
1266 | Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P()); | |
1267 | mother->Print(); | |
1268 | imo = mother->GetFirstMother(); | |
1269 | } | |
1270 | Printf("########################"); | |
1271 | } |