silvermy@ornl.gov - basic interface to calibration info for testbeam (CERN 2007)...
[u/mrichter/AliRoot.git] / EVE / EveBase / AliEveKineTools.cxx
CommitLineData
d810d0de 1// $Id$
2// Main authors: Matevz Tadel & Alja Mrak-Tadel: 2006, 2007
c63fcc1d 3
d810d0de 4/**************************************************************************
5 * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. *
6 * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for *
51346b82 7 * full copyright notice. *
d810d0de 8 **************************************************************************/
9
10#include "AliEveKineTools.h"
c63fcc1d 11
12#include <TObject.h>
13#include <TTree.h>
14#include <TBranchElement.h>
7704f096 15#include <TClonesArray.h>
16
c63fcc1d 17#include <AliStack.h>
18#include <AliTrackReference.h>
19
84aff7a4 20#include <TEveTrack.h>
21#include <TEveElement.h>
c63fcc1d 22
7704f096 23#include <map>
c63fcc1d 24
57ffa5fb 25//______________________________________________________________________________
d810d0de 26// AliEveKineTools
c63fcc1d 27//
40790e5b 28// Tools for import of kinematics. Preliminary version.
29//
d810d0de 30
d810d0de 31ClassImp(AliEveKineTools)
c63fcc1d 32
40790e5b 33namespace {
34
35 typedef std::map<Int_t, TEveTrack*> TrackMap_t;
36
37 void MapTracks(TrackMap_t& map, TEveElement* cont, Bool_t recurse)
38 {
39 TEveElement::List_i i = cont->BeginChildren();
40 while (i != cont->EndChildren()) {
41 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
42 map[track->GetLabel()] = track;
43 if (recurse)
44 MapTracks(map, track, recurse);
45 ++i;
46 }
47 }
48}
49
50/**************************************************************************/
76461dad 51
d810d0de 52void AliEveKineTools::SetDaughterPathMarks(TEveElement* cont, AliStack* stack, Bool_t recurse)
c63fcc1d 53{
76461dad 54 // Import daughters birth points.
55
84aff7a4 56 TEveElement::List_i iter = cont->BeginChildren();
c63fcc1d 57 while(iter != cont->EndChildren())
58 {
51346b82 59 TEveTrack* track = dynamic_cast<TEveTrack*>(*iter);
c63fcc1d 60 TParticle* p = stack->Particle(track->GetLabel());
a15e6d7d 61 if (p->GetNDaughters())
62 {
c63fcc1d 63 Int_t d0 = p->GetDaughter(0), d1 = p->GetDaughter(1);
a15e6d7d 64 for(int d = d0; d > 0 && d <= d1; ++d)
51346b82 65 {
c63fcc1d 66 TParticle* dp = stack->Particle(d);
a15e6d7d 67 track->AddPathMark(TEvePathMark(TEvePathMark::kDaughter,
68 TEveVector(dp->Vx(), dp->Vy(), dp->Vz()),
69 TEveVector(dp->Px(), dp->Py(), dp->Pz()),
70 dp->T()));
c63fcc1d 71 }
27aa96b0 72 if (recurse)
73 SetDaughterPathMarks(track, stack, recurse);
c63fcc1d 74 }
75 ++iter;
76 }
77}
78
40790e5b 79/**************************************************************************/
974767e6 80
d810d0de 81void AliEveKineTools::SetTrackReferences(TEveElement* cont, TTree* treeTR, Bool_t recurse)
974767e6 82{
40790e5b 83 // Set decay and track reference path-marks.
974767e6 84
a15e6d7d 85 static const TEveException kEH("AliEveKineTools::ImportPathMarks");
974767e6 86
40790e5b 87 TrackMap_t map;
88 MapTracks(map, cont, recurse);
51346b82 89
c63fcc1d 90 Int_t nPrimaries = (Int_t) treeTR->GetEntries();
91 TIter next(treeTR->GetListOfBranches());
92 TBranchElement* el;
93 Bool_t isRef = kTRUE;
c63fcc1d 94
95 while ((el = (TBranchElement*) next()))
96 {
97 if (strcmp("AliRun",el->GetName()) == 0)
98 isRef = kFALSE;
99
1d74e2b0 100 TClonesArray* arr = 0;
101 el->SetAddress(&arr);
51346b82 102 for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++)
c63fcc1d 103 {
1d74e2b0 104 el->GetEntry(iPrimPart);
105
106 Int_t last_label = -1;
40790e5b 107 TrackMap_t::iterator iter = map.end();
1d74e2b0 108 Int_t Nent = arr->GetEntriesFast();
51346b82 109 for (Int_t iTrackRef = 0; iTrackRef < Nent; iTrackRef++)
c63fcc1d 110 {
1d74e2b0 111 AliTrackReference* atr = (AliTrackReference*)arr->UncheckedAt(iTrackRef);
112
113 Int_t label = atr->GetTrack();
114 if (label < 0)
a15e6d7d 115 throw(kEH + Form("negative label for entry %d in branch %s.",
1d74e2b0 116 iTrackRef, el->GetName()));
51346b82 117
a15e6d7d 118 if (label != last_label)
119 {
120 iter = map.find(label);
1d74e2b0 121 last_label = label;
122 }
123
a15e6d7d 124 if (iter != map.end())
125 {
126 iter->second->AddPathMark
127 (TEvePathMark(isRef ? TEvePathMark::kReference : TEvePathMark::kDecay,
128 TEveVector(atr->X(), atr->Y(), atr->Z()),
129 TEveVector(atr->Px(), atr->Py(), atr->Pz()),
130 atr->GetTime()));
c63fcc1d 131 }
51346b82 132 } // loop track refs
40790e5b 133 } // loop primaries in clones arrays
1d74e2b0 134 delete arr;
c63fcc1d 135 } // end loop through top branches
27aa96b0 136}
137
40790e5b 138/**************************************************************************/
139
140void AliEveKineTools::SortPathMarks(TEveElement* el, Bool_t recurse)
27aa96b0 141{
142 // Sort path-marks for all tracks by time.
143
a15e6d7d 144 // !!!!! MT ... this is one level deep only!
145
40790e5b 146 TEveTrack* track = dynamic_cast<TEveTrack*>(el);
a15e6d7d 147 if (track) track->SortPathMarksByTime();
c63fcc1d 148
40790e5b 149 TEveElement::List_i i = el->BeginChildren();
a15e6d7d 150 while (i != el->EndChildren() && recurse)
151 {
40790e5b 152 track = dynamic_cast<TEveTrack*>(el);
153 if (track) track->SortPathMarksByTime();
a15e6d7d 154 ++i;
c63fcc1d 155 }
156}