Add full support for SPD tracklets.
[u/mrichter/AliRoot.git] / EVE / EveBase / AliEveTrackCounter.cxx
1 // @(#)root/eve:$Id$
2 // Author: Matevz Tadel 2007
3
4 /**************************************************************************
5  * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. *
6  * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for          *
7  * full copyright notice.                                                 *
8  **************************************************************************/
9
10 #include "AliEveTrackCounter.h"
11
12 #include "TEveManager.h"
13 #include "AliEveTrack.h"
14 #include "AliEveTracklet.h"
15
16 #include <TEveGedEditor.h>
17
18 //==============================================================================
19 // AliEveTrackCounter
20 //==============================================================================
21
22 //______________________________________________________________________________
23 //
24 // Provides event-based method for tagging of good / bad (or primary /
25 // secondary) tracks. A report can be written into a text file.
26 //
27 // AliEveTrack status is toggled by using secondary-selection / ctrl-click
28 // functionality of the GL viewer.
29 //
30 // Some of the functionality is implemented in AliEveTrackCounterEditor
31 // class.
32
33 ClassImp(AliEveTrackCounter)
34
35 //______________________________________________________________________________
36 AliEveTrackCounter* AliEveTrackCounter::fgInstance = 0;
37
38 //______________________________________________________________________________
39 AliEveTrackCounter::AliEveTrackCounter(const Text_t* name, const Text_t* title) :
40   TEveElement(),
41   TNamed(name, title),
42
43   fBadLineStyle (6),
44   fClickAction  (kCA_ToggleTrack),
45   fEventId      (-1),
46   fAllTracks    (0), fGoodTracks   (0),
47   fAllTracklets (0), fGoodTracklets(0),
48   fTrackLists   (),  fTrackletLists()
49 {
50   // Constructor.
51   // Connects to global signal "AliEveTrack", "SecSelected(AliEveTrack*)".
52
53   if (fgInstance == 0) fgInstance = this;
54
55   TQObject::Connect("AliEveTrack", "SecSelectedTrack(AliEveTrack*)",
56                     "AliEveTrackCounter", this, "DoTrackAction(AliEveTrack*)");
57   TQObject::Connect("AliEveTracklet", "SecSelectedTracklet(AliEveTracklet*)",
58                     "AliEveTrackCounter", this, "DoTrackletAction(AliEveTracklet*)");
59 }
60
61 //______________________________________________________________________________
62 AliEveTrackCounter::~AliEveTrackCounter()
63 {
64   // Destructor.
65   // Disconnect from the global track signals.
66
67   TQObject::Disconnect("AliEveTrack", "DoTrackAction(AliEveTrack*)");
68   TQObject::Disconnect("AliEveTracklet", "DoTrackletAction(AliEveTracklet*)");
69   if (fgInstance == this) fgInstance = 0;
70 }
71
72 /******************************************************************************/
73
74 //______________________________________________________________________________
75 void AliEveTrackCounter::Reset()
76 {
77   // Reset internal track-counters and track-list.
78
79   fAllTracks     = fGoodTracks    = 0;
80   fAllTracklets  = fGoodTracklets = 0;
81   {
82     TIter next(&fTrackLists);
83     TEveTrackList* tlist;
84     while ((tlist = dynamic_cast<TEveTrackList*>(next())))
85       tlist->DecDenyDestroy();
86     fTrackLists.Clear("nodelete");
87   }
88   {
89     TIter next(&fTrackletLists);
90     TEveTrackList* tlist;
91     while ((tlist = dynamic_cast<TEveTrackList*>(next())))
92       tlist->DecDenyDestroy();
93     fTrackletLists.Clear("nodelete");
94   }
95 }
96
97 //______________________________________________________________________________
98 void AliEveTrackCounter::RegisterTracks(TEveTrackList* tlist, Bool_t goodTracks)
99 {
100   // Register tracks from tlist and tlist itself.
101   // If goodTracks is true, they are considered as primary/good
102   // tracks.
103
104   tlist->IncDenyDestroy();
105   fTrackLists.Add(tlist);
106
107   List_i i = tlist->BeginChildren();
108   while (i != tlist->EndChildren())
109   {
110     AliEveTrack* t = dynamic_cast<AliEveTrack*>(*i);
111     if (t != 0)
112     {
113       if (goodTracks)
114       {
115         ++fGoodTracks;
116       } else {
117         t->SetLineStyle(fBadLineStyle);
118       }
119       ++fAllTracks;
120     }
121     ++i;
122   }
123 }
124
125 //______________________________________________________________________________
126 void AliEveTrackCounter::RegisterTracklets(TEveTrackList* tlist, Bool_t goodTracks)
127 {
128   // Register tracklets from tlist and tlist itself.
129   // If goodTracks is true, they are considered as primary/good
130   // tracks.
131
132   tlist->IncDenyDestroy();
133   fTrackletLists.Add(tlist);
134
135   List_i i = tlist->BeginChildren();
136   while (i != tlist->EndChildren())
137   {
138     AliEveTracklet* t = dynamic_cast<AliEveTracklet*>(*i);
139     if (t != 0)
140     {
141       if (goodTracks)
142       {
143         ++fGoodTracklets;
144       } else {
145         t->SetLineStyle(fBadLineStyle);
146       }
147       ++fAllTracklets;
148     }
149     ++i;
150   }
151 }
152
153 //______________________________________________________________________________
154 void AliEveTrackCounter::DoTrackAction(AliEveTrack* track)
155 {
156    // Slot called when track is secondary selected.
157    //
158    // No check is done if track actually belongs to one of the
159    // registered track-lists.
160    //
161    // Probably it would be safer to copy good/bad tracks into special
162    // sub-containers.
163    // In this case one should also override RemoveElementLocal.
164
165    static const TEveException eh("AliEveTrackCounter::DoTrackAction ");
166
167    switch (fClickAction)
168    {
169
170       case kCA_PrintTrackInfo:
171       {
172          printf("AliEveTrack '%s'\n", track->GetObject(eh)->GetName());
173          const TEveVector &v = track->GetVertex();
174          const TEveVector &p = track->GetMomentum();;
175          printf("  Vx=%f, Vy=%f, Vz=%f; Pt=%f, Pz=%f, phi=%f)\n",
176                 v.fX, v.fY, v.fZ, p.Perp(), p.fZ, TMath::RadToDeg()*p.Phi());
177          printf("  <other information should be printed ... full AliESDtrack>\n");
178          break;
179       }
180
181       case kCA_ToggleTrack:
182       {
183          if (track->GetLineStyle() == 1)
184          {
185             track->SetLineStyle(fBadLineStyle);
186             --fGoodTracks;
187          } else {
188             track->SetLineStyle(1);
189             ++fGoodTracks;
190          }
191          track->ElementChanged();
192          gEve->Redraw3D();
193
194          printf("AliEveTrackCounter::DoTrackAction All=%d, Good=%d, Bad=%d\n",
195                 fAllTracks, fGoodTracks, fAllTracks-fGoodTracks);
196
197          if (gEve->GetEditor()->GetModel() == GetObject(eh))
198             gEve->EditElement(this);
199
200          break;
201       }
202
203    } // end switch fClickAction
204 }
205
206 //______________________________________________________________________________
207 void AliEveTrackCounter::DoTrackletAction(AliEveTracklet* track)
208 {
209    // Slot called when tracklet is secondary selected.
210    //
211    // No check is done if track actually belongs to one of the
212    // registered track-lists.
213    //
214    // Probably it would be safer to copy good/bad tracks into special
215    // sub-containers.
216    // In this case one should also override RemoveElementLocal.
217
218    static const TEveException eh("AliEveTrackCounter::DoTrackletAction ");
219
220    switch (fClickAction)
221    {
222
223       case kCA_PrintTrackInfo:
224       {
225          printf("AliEveTracklet '%s'\n", track->GetObject(eh)->GetName());
226          const TEveVector &v = track->GetVertex();
227          const TEveVector &p = track->GetMomentum();;
228          printf("  Vx=%f, Vy=%f, Vz=%f; Pt=%f, Pz=%f, phi=%f)\n",
229                 v.fX, v.fY, v.fZ, p.Perp(), p.fZ, TMath::RadToDeg()*p.Phi());
230          printf("  <other information should be printed ... full AliESDtrack>\n");
231          break;
232       }
233
234       case kCA_ToggleTrack:
235       {
236          if (track->GetLineStyle() == 1)
237          {
238             track->SetLineStyle(fBadLineStyle);
239             --fGoodTracklets;
240          } else {
241             track->SetLineStyle(1);
242             ++fGoodTracklets;
243          }
244          track->ElementChanged();
245          gEve->Redraw3D();
246
247          printf("AliEveTrackCounter::DoTrackletAction All=%d, Good=%d, Bad=%d\n",
248                 fAllTracklets, fGoodTracklets, fAllTracklets-fGoodTracklets);
249
250          if (gEve->GetEditor()->GetModel() == GetObject(eh))
251             gEve->EditElement(this);
252
253          break;
254       }
255
256    } // end switch fClickAction
257 }
258
259 /******************************************************************************/
260
261 //______________________________________________________________________________
262 void AliEveTrackCounter::OutputEventTracks(FILE* out)
263 {
264   // Print good-track summary into a plain-text file by iteration
265   // through all registered track-lists.
266   // State of each track is determined by its line-style, it is
267   // considered a good track if it's line style is solid.
268
269   if (out == 0)
270   {
271     out = stdout;
272     fprintf(out, "AliEveTrackCounter::OutputEventTracks()\n");
273   }
274
275   fprintf(out, "Event=%d\n", fEventId);
276   fprintf(out, "GoodTracks=%d  AllTracks=%d\n", fGoodTracks, fAllTracks);
277
278   {
279     TIter tlists(&fTrackLists);
280     TEveTrackList* tlist;
281     Int_t cnt = 0;
282     while ((tlist = (TEveTrackList*) tlists()) != 0)
283     {
284       List_i i = tlist->BeginChildren();
285       while (i != tlist->EndChildren())
286       {
287         AliEveTrack* t = dynamic_cast<AliEveTrack*>(*i);
288         if (t != 0 && t->GetLineStyle() == 1)
289         {
290           ++cnt;
291           fprintf(out, " %2d: chg=%+2d  pt=%8.5f  eta=%+8.5f  phi=%+8.5f\n",
292                   cnt, t->GetCharge(), t->GetMomentum().Perp(),
293                   t->GetMomentum().Eta(), t->GetMomentum().Phi());
294         }
295         ++i;
296       }
297     }
298   }
299
300   fprintf(out, "GoodTracklets=%d  AllTracklets=%d\n", fGoodTracklets, fAllTracklets);
301   {
302     TIter tlists(&fTrackletLists);
303     TEveTrackList* tlist;
304     Int_t cnt = 0;
305     while ((tlist = (TEveTrackList*) tlists()) != 0)
306     {
307       List_i i = tlist->BeginChildren();
308       while (i != tlist->EndChildren())
309       {
310         AliEveTracklet* t = dynamic_cast<AliEveTracklet*>(*i);
311         if (t != 0 && t->GetLineStyle() == 1)
312         {
313           ++cnt;
314           fprintf(out, " %2d: theta=%+8.5f  eta=%+8.5f  phi=%+8.5f\n",
315                   cnt, t->GetMomentum().Theta(), t->GetMomentum().Eta(), t->GetMomentum().Phi());
316         }
317         ++i;
318       }
319     }
320   }
321 }