]>
Commit | Line | Data |
---|---|---|
c63fcc1d | 1 | // $Header$ |
2 | ||
3 | #include "KineTools.h" | |
4 | ||
5 | #include <TObject.h> | |
6 | #include <TTree.h> | |
7 | #include <TBranchElement.h> | |
7704f096 | 8 | #include <TClonesArray.h> |
9 | ||
c63fcc1d | 10 | #include <AliStack.h> |
11 | #include <AliTrackReference.h> | |
12 | ||
13 | #include "Reve/Track.h" | |
14 | #include "Reve/RenderElement.h" | |
15 | ||
16 | #include <algorithm> | |
7704f096 | 17 | #include <map> |
c63fcc1d | 18 | |
19 | //______________________________________________________________________ | |
20 | // KineTools | |
21 | // | |
22 | ||
23 | using namespace Reve; | |
24 | using namespace Alieve; | |
25 | using namespace std; | |
26 | ||
27 | ClassImp(KineTools) | |
28 | ||
29 | KineTools::KineTools() | |
76461dad | 30 | {} |
c63fcc1d | 31 | |
32 | /**************************************************************************/ | |
76461dad | 33 | |
974767e6 | 34 | void KineTools::SetDaughterPathMarks(RenderElement* cont, AliStack* stack, Bool_t recurse) |
c63fcc1d | 35 | { |
76461dad | 36 | // Import daughters birth points. |
37 | ||
c63fcc1d | 38 | RenderElement::List_i iter = cont->BeginChildren(); |
39 | ||
40 | while(iter != cont->EndChildren()) | |
41 | { | |
42 | Track* track = dynamic_cast<Track*>(*iter); | |
43 | TParticle* p = stack->Particle(track->GetLabel()); | |
44 | if(p->GetNDaughters()) { | |
45 | Int_t d0 = p->GetDaughter(0), d1 = p->GetDaughter(1); | |
46 | for(int d=d0; d>0 && d<=d1;++d) | |
47 | { | |
48 | TParticle* dp = stack->Particle(d); | |
1d74e2b0 | 49 | Reve::PathMark* pm = new PathMark(PathMark::Daughter); |
c63fcc1d | 50 | pm->V.Set(dp->Vx(),dp->Vy(), dp->Vz()); |
51 | pm->P.Set(dp->Px(),dp->Py(), dp->Pz()); | |
52 | pm->time = dp->T(); | |
53 | track->AddPathMark(pm); | |
974767e6 | 54 | if (recurse) |
55 | SetDaughterPathMarks(track, stack, recurse); | |
c63fcc1d | 56 | } |
57 | } | |
58 | ++iter; | |
59 | } | |
60 | } | |
61 | ||
62 | /**************************************************************************/ | |
63 | ||
64 | namespace { | |
974767e6 | 65 | struct cmp_pathmark |
66 | { | |
c63fcc1d | 67 | bool operator()(PathMark* const & a, PathMark* const & b) |
68 | { return a->time < b->time; } | |
69 | }; | |
c63fcc1d | 70 | |
974767e6 | 71 | void slurp_tracks(map<Int_t, Track*>& tracks, RenderElement* cont, Bool_t recurse) |
c63fcc1d | 72 | { |
974767e6 | 73 | RenderElement::List_i citer = cont->BeginChildren(); |
1d74e2b0 | 74 | while(citer != cont->EndChildren()) |
75 | { | |
76 | Track* track = dynamic_cast<Track*>(*citer); | |
77 | tracks[track->GetLabel()] = track; | |
974767e6 | 78 | if (recurse) |
79 | slurp_tracks(tracks, track, recurse); | |
1d74e2b0 | 80 | ++citer; |
c63fcc1d | 81 | } |
974767e6 | 82 | } |
83 | ||
84 | } | |
85 | ||
86 | void KineTools::SetTrackReferences(RenderElement* cont, TTree* treeTR, Bool_t recurse) | |
87 | { | |
88 | // set decay and reference points | |
89 | ||
90 | static const Exc_t eH("KineTools::ImportPathMarks"); | |
91 | ||
92 | // Fill map | |
93 | map<Int_t, Track*> tracks; | |
94 | slurp_tracks(tracks, cont, recurse); | |
1d74e2b0 | 95 | |
c63fcc1d | 96 | Int_t nPrimaries = (Int_t) treeTR->GetEntries(); |
97 | TIter next(treeTR->GetListOfBranches()); | |
98 | TBranchElement* el; | |
99 | Bool_t isRef = kTRUE; | |
c63fcc1d | 100 | |
101 | while ((el = (TBranchElement*) next())) | |
102 | { | |
103 | if (strcmp("AliRun",el->GetName()) == 0) | |
104 | isRef = kFALSE; | |
105 | ||
1d74e2b0 | 106 | TClonesArray* arr = 0; |
107 | el->SetAddress(&arr); | |
c63fcc1d | 108 | for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++) |
109 | { | |
1d74e2b0 | 110 | el->GetEntry(iPrimPart); |
111 | ||
112 | Int_t last_label = -1; | |
113 | map<Int_t, Track*>::iterator iter = tracks.end(); | |
114 | Int_t Nent = arr->GetEntriesFast(); | |
115 | for (Int_t iTrackRef = 0; iTrackRef < Nent; iTrackRef++) | |
c63fcc1d | 116 | { |
1d74e2b0 | 117 | AliTrackReference* atr = (AliTrackReference*)arr->UncheckedAt(iTrackRef); |
118 | ||
119 | Int_t label = atr->GetTrack(); | |
120 | if (label < 0) | |
121 | throw(eH + Form("negative label for entry %d in branch %s.", | |
122 | iTrackRef, el->GetName())); | |
123 | ||
124 | if(label != last_label) { | |
125 | iter = tracks.find(label); | |
126 | last_label = label; | |
127 | } | |
128 | ||
129 | if (iter != tracks.end()) { | |
130 | PathMark* pm = new PathMark(isRef ? PathMark::Reference : PathMark::Decay); | |
131 | pm->V.Set(atr->X(),atr->Y(), atr->Z()); | |
132 | pm->P.Set(atr->Px(),atr->Py(), atr->Pz()); | |
133 | pm->time = atr->GetTime(); | |
134 | Track* track = iter->second; | |
135 | track->AddPathMark(pm); | |
c63fcc1d | 136 | } |
137 | } // loop track refs | |
c63fcc1d | 138 | } // loop primaries, clones arrays |
1d74e2b0 | 139 | delete arr; |
c63fcc1d | 140 | } // end loop through top branches |
c63fcc1d | 141 | |
1d74e2b0 | 142 | // sort |
974767e6 | 143 | for(map<Int_t, Track*>::iterator j=tracks.begin(); j!=tracks.end(); ++j) |
144 | { | |
1d74e2b0 | 145 | (j->second)->SortPathMarksByTime(); |
c63fcc1d | 146 | } |
147 | } |