Record changes.
[u/mrichter/AliRoot.git] / EVE / Reve / Track.cxx
CommitLineData
5a5a1232 1// $Header$
2
3#include "Track.h"
4#include "MCHelixLine.hi"
7cbaaebd 5#include "PointSet.h"
5a5a1232 6
7#include <TPolyLine3D.h>
8#include <TPolyMarker3D.h>
9#include <TColor.h>
10
11// Updates
12#include <Reve/RGTopFrame.h>
13#include <TCanvas.h>
14
15#include <vector>
16
17using namespace Reve;
18
19//______________________________________________________________________
20// Track
21//
22
23ClassImp(Reve::Track)
24
265ecb21 25Track::Track() :
eadaa89b 26 Line(),
5a5a1232 27
265ecb21 28 fV(),
29 fP(),
30 fBeta(0),
31 fCharge(0),
7cbaaebd 32 fLabel(-1),
33 fIndex(-1),
265ecb21 34 fPathMarks(),
35
eadaa89b 36 fRnrStyle(0)
265ecb21 37{}
38
3d75306d 39Track::Track(TParticle* t, Int_t label, TrackRnrStyle* rs):
eadaa89b 40 Line(),
3d75306d 41
42 fV(t->Vx(), t->Vy(), t->Vz()),
43 fP(t->Px(), t->Py(), t->Pz()),
44 fBeta(t->P()/t->Energy()),
45 fCharge(0),
46 fLabel(label),
7cbaaebd 47 fIndex(-1),
3d75306d 48 fPathMarks(),
49
eadaa89b 50 fRnrStyle(rs)
3d75306d 51{
52 fLineColor = fRnrStyle->GetColor();
53 fMainColorPtr = &fLineColor;
54
55 TParticlePDG* pdgp = t->GetPDG();
56 if (pdgp)
57 fCharge = (Int_t) TMath::Nint(pdgp->Charge()/3);
eadaa89b 58
59 SetName(t->GetName());
3d75306d 60}
61
265ecb21 62Track::Track(Reve::MCTrack* t, TrackRnrStyle* rs):
eadaa89b 63 Line(),
5a5a1232 64
265ecb21 65 fV(t->Vx(), t->Vy(), t->Vz()),
66 fP(t->Px(), t->Py(), t->Pz()),
67 fBeta(t->P()/t->Energy()),
68 fCharge(0),
69 fLabel(t->label),
fa446a68 70 fIndex(t->index),
265ecb21 71 fPathMarks(),
72
eadaa89b 73 fRnrStyle(rs)
265ecb21 74{
5a5a1232 75 fLineColor = fRnrStyle->GetColor();
76 fMainColorPtr = &fLineColor;
77
5a5a1232 78 TParticlePDG* pdgp = t->GetPDG();
79 if(pdgp == 0) {
80 t->ResetPdgCode(); pdgp = t->GetPDG();
81 }
48dc973d 82 fCharge = (Int_t) TMath::Nint(pdgp->Charge()/3);
eadaa89b 83
84 SetName(t->GetName());
5a5a1232 85}
86
265ecb21 87Track::Track(Reve::RecTrack* t, TrackRnrStyle* rs) :
eadaa89b 88 Line(),
265ecb21 89
90 fV(t->V),
91 fP(t->P),
92 fBeta(t->beta),
93 fCharge(t->sign),
94 fLabel(t->label),
fa446a68 95 fIndex(t->index),
265ecb21 96 fPathMarks(),
97
eadaa89b 98 fRnrStyle(rs)
5a5a1232 99{
5a5a1232 100 fLineColor = fRnrStyle->GetColor();
101 fMainColorPtr = &fLineColor;
eadaa89b 102
103 SetName(t->GetName());
5a5a1232 104}
105
106Track::~Track()
265ecb21 107{
108 for (vpPathMark_i i=fPathMarks.begin(); i!=fPathMarks.end(); ++i)
109 delete *i;
110}
5a5a1232 111
eadaa89b 112/*
5a5a1232 113void Track::Reset(Int_t n_points)
114{
115 delete [] TPolyLine3D::fP; TPolyLine3D::fP = 0;
116 fN = n_points;
117 if(fN) TPolyLine3D::fP = new Float_t [3*fN];
118 memset(TPolyLine3D::fP, 0, 3*fN*sizeof(Float_t));
119 fLastPoint = -1;
120}
eadaa89b 121*/
5a5a1232 122
7cbaaebd 123 /**************************************************************************/
5a5a1232 124
125void Track::MakeTrack()
126{
127
128 TrackRnrStyle& RS((fRnrStyle != 0) ? *fRnrStyle : TrackRnrStyle::fgDefStyle);
129
130 Float_t px = fP.x, py = fP.y, pz = fP.z;
131
132 MCVertex mc_v0;
133 mc_v0.x = fV.x;
134 mc_v0.y = fV.y;
135 mc_v0.z = fV.z;
136 mc_v0.t = 0;
137
138 std::vector<MCVertex> track_points;
48dc973d 139 Bool_t decay = kFALSE;
5a5a1232 140
141 if ((TMath::Abs(fV.z) > RS.fMaxZ) || (fV.x*fV.x + fV.y*fV.y > RS.fMaxR*RS.fMaxR))
142 goto make_polyline;
143
01024c63 144 if (fCharge != 0 && TMath::Abs(RS.fMagField) > 1e-5) {
145
146 // Charged particle in magnetic field
5a5a1232 147
48dc973d 148 Float_t a = RS.fgkB2C * RS.fMagField * fCharge;
5a5a1232 149
150 MCHelix helix(fRnrStyle, &mc_v0, TMath::C()*fBeta, &track_points, a); //m->cm
151 helix.Init(TMath::Sqrt(px*px+py*py), pz);
152
7cbaaebd 153 if(!fPathMarks.empty())
154 {
155 for(std::vector<Reve::PathMark*>::iterator i=fPathMarks.begin(); i!=fPathMarks.end(); ++i)
156 {
5a5a1232 157 Reve::PathMark* pm = *i;
158
7cbaaebd 159 if (RS.fFitReferences && pm->type == Reve::PathMark::Reference)
160 {
5a5a1232 161 if(TMath::Abs(pm->V.z) > RS.fMaxZ
162 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
163 goto helix_bounds;
164
7cbaaebd 165 // printf("%s fit reference \n", fName.Data());
166 helix.LoopToVertex(px, py, pz, pm->V.x, pm->V.y, pm->V.z);
167 px = pm->P.x;
168 py = pm->P.y;
169 pz = pm->P.z;
5a5a1232 170 }
7cbaaebd 171 else if(RS.fFitDaughters && pm->type == Reve::PathMark::Daughter)
172 {
5a5a1232 173 if(TMath::Abs(pm->V.z) > RS.fMaxZ
174 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
175 goto helix_bounds;
7cbaaebd 176
177 // printf("%s fit daughter \n", fName.Data());
178 helix.LoopToVertex(px, py, pz, pm->V.x, pm->V.y, pm->V.z);
179 px -= pm->P.x;
180 py -= pm->P.y;
181 pz -= pm->P.z;
182 }
183 else if(RS.fFitDecay && pm->type == Reve::PathMark::Decay)
184 {
185 if(TMath::Abs(pm->V.z) > RS.fMaxZ
186 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
187 goto helix_bounds;
188 helix.LoopToVertex(px, py, pz, pm->V.x, pm->V.y, pm->V.z);
5a5a1232 189 decay = true;
190 break;
191 }
192 }
193 }
194 helix_bounds:
7cbaaebd 195 // go to bounds
48dc973d 196 if(!decay || RS.fFitDecay == kFALSE){
5a5a1232 197 helix.LoopToBounds(px,py,pz);
198 // printf("%s loop to bounds \n",fName.Data() );
199 }
200
01024c63 201 } else {
202
203 // Neutral particle or no field
5a5a1232 204
205 MCLine line(fRnrStyle, &mc_v0, TMath::C()*fBeta, &track_points);
206
207 if(!fPathMarks.empty()){
208 for(std::vector<Reve::PathMark*>::iterator i=fPathMarks.begin(); i!=fPathMarks.end(); ++i) {
209 Reve::PathMark* pm = *i;
210
211 if(RS.fFitDaughters && pm->type == Reve::PathMark::Daughter){
212 if(TMath::Abs(pm->V.z) > RS.fMaxZ
213 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
214 goto line_bounds;
215 line.GotoVertex(pm->V.x, pm->V.y, pm->V.z);
216 fP.x -= pm->P.x;
217 fP.y -= pm->P.y;
218 fP.z -= pm->P.z;
219 }
220
221 if(RS.fFitDecay && pm->type == Reve::PathMark::Decay){
222 if(TMath::Abs(pm->V.z) > RS.fMaxZ
223 || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
224 goto line_bounds;
225 line.GotoVertex(pm->V.x, pm->V.y, pm->V.z);
226 decay = true;
227 break;
228 }
229 }
230 }
231
232 line_bounds:
48dc973d 233 if(!decay || RS.fFitDecay == kFALSE)
5a5a1232 234 line.GotoBounds(px,py,pz);
235
236 }
237make_polyline:
238 Reset(track_points.size());
239 for(std::vector<MCVertex>::iterator i=track_points.begin(); i!=track_points.end(); ++i)
240 SetNextPoint(i->x, i->y, i->z);
241}
242
243/**************************************************************************/
244
48dc973d 245void Track::ImportHits()
246{
247 Reve::LoadMacro("hits_from_label.C");
248 gROOT->ProcessLine(Form("hits_from_label(%d);", fLabel));
249}
250
251void Track::ImportClusters()
252{
253 Reve::LoadMacro("clusters_from_label.C");
254 gROOT->ProcessLine(Form("clusters_from_label(%d);", fLabel));
255}
256
fa446a68 257void Track::ImportClustersFromIndex()
258{
7cbaaebd 259 static const Exc_t eH("Track::ImportClustersFromIndex ");
260
261 if (fIndex < 0)
262 throw(eH + "index not set.");
263
fa446a68 264 Reve::LoadMacro("clusters_from_index.C");
265 gROOT->ProcessLine(Form("clusters_from_index(%d);", fIndex));
266}
267
a525aa60 268void Track::ImportDaughters()
269{
270 static const Exc_t eH("Track::ImportDaughters ");
271
272 if (fLabel < 0)
273 throw(eH + "label not set.");
274
be0e42e3 275 Reve::LoadMacro("kine_daughter_tracks.C");
276 gROOT->ProcessLine(Form("kine_daughter_tracks(%d);", fLabel));
a525aa60 277}
278
279/**************************************************************************/
280
b1f08706 281void Track::PrintKineStack()
282{
283 Reve::LoadMacro("print_kine_from_label.C");
284 gROOT->ProcessLine(Form("print_kine_from_label(%d);", fLabel));
285}
286
a525aa60 287
288void Track::PrintPathMarks()
289{
290 static const Exc_t eH("Track::PrintPathMarks ");
291
292 if (fLabel < 0)
293 throw(eH + "label not set.");
294
295 printf("Number of path marks %d label %d\n",
296 fPathMarks.size(), fLabel);
297
298 PathMark* pm;
299 for(vpPathMark_i i=fPathMarks.begin(); i!=fPathMarks.end(); i++)
300 {
301 pm = *i;
302 printf("Reve::PathMark: %-9s p: %8f %8f %8f Vertex: %8e %8e %8e %g \n",
303 pm->type_name(),
304 pm->P.x, pm->P.y, pm->P.z,
305 pm->V.x, pm->V.y, pm->V.z,
306 pm->time);
307 }
308}
309
9c39ede9 310/**************************************************************************/
311
312void Track::CtrlClicked(Reve::Track* track)
313{
314 Emit("CtrlClicked(Reve::Track*)", (Long_t)track);
315}
316
48dc973d 317
318/**************************************************************************/
319/**************************************************************************/
320
5a5a1232 321//______________________________________________________________________
322// TrackRnrStyle
323//
324
325ClassImp(Reve::TrackRnrStyle)
326
48dc973d 327Float_t TrackRnrStyle::fgDefMagField = 5;
328const Float_t TrackRnrStyle::fgkB2C = 0.299792458e-3;
5a5a1232 329TrackRnrStyle TrackRnrStyle::fgDefStyle;
330
265ecb21 331TrackRnrStyle::TrackRnrStyle() :
332 TObject(),
5a5a1232 333
265ecb21 334 fColor(1),
2074deef 335 fWidth(1),
265ecb21 336 fMagField(fgDefMagField),
5a5a1232 337
265ecb21 338 fMaxR (350),
339 fMaxZ (450),
5a5a1232 340
265ecb21 341 fMaxOrbs (0.5),
342 fMinAng (45),
343 fDelta (0.1),
5a5a1232 344
22df2a83 345 fMinPt (0.1),
346 fMaxPt (10),
347
7cbaaebd 348 fFitDaughters (kTRUE),
349 fFitReferences (kTRUE),
350 fFitDecay (kTRUE),
5a5a1232 351
7cbaaebd 352 fRnrDaughters (kTRUE),
353 fRnrReferences (kTRUE),
354 fRnrDecay (kTRUE)
355{}
5a5a1232 356/**************************************************************************/
357/**************************************************************************/
358
359//______________________________________________________________________
360// TrackList
361//
362
363ClassImp(Reve::TrackList)
364
365void TrackList::Init()
366{
7cbaaebd 367 fMarkerStyle = 5;
5a5a1232 368 fMarkerColor = 5;
369 // fMarker->SetMarkerSize(0.05);
370
a8600b56 371 if (fRnrStyle== 0) fRnrStyle = new TrackRnrStyle;
372 SetMainColorPtr(&fRnrStyle->fColor);
5a5a1232 373}
374
a8600b56 375TrackList::TrackList(Int_t n_tracks, TrackRnrStyle* rs) :
7d42b6c2 376 RenderElementListBase(),
265ecb21 377 TPolyMarker3D(n_tracks),
378
379 fTitle(),
380
7cbaaebd 381 fRnrStyle (rs),
382 fRnrTracks (kTRUE),
383 fEditPathMarks (kFALSE)
5a5a1232 384{
385 Init();
386}
387
a8600b56 388TrackList::TrackList(const Text_t* name, Int_t n_tracks, TrackRnrStyle* rs) :
7d42b6c2 389 RenderElementListBase(),
265ecb21 390 TPolyMarker3D(n_tracks),
a8600b56 391
265ecb21 392 fTitle(),
393
a8600b56 394 fRnrStyle (rs),
265ecb21 395 fRnrTracks (kTRUE)
5a5a1232 396{
397 Init();
398 SetName(name);
399}
400
401void TrackList::Reset(Int_t n_tracks)
402{
403 delete [] fP; fP = 0;
404 fN = n_tracks;
405 if(fN) fP = new Float_t [3*fN];
406 memset(fP, 0, 3*fN*sizeof(Float_t));
407 fLastPoint = -1;
408}
409
410/**************************************************************************/
411
412void TrackList::Paint(Option_t* option)
413{
414 if(fRnrElement) {
415 if(fRnrMarkers) {
416 TPolyMarker3D::Paint(option);
417 }
418 if(fRnrTracks) {
9c39ede9 419 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
5a5a1232 420 if((*i)->GetRnrElement())
421 (*i)->GetObject()->Paint(option);
422 }
423 }
424 }
425}
426
427/**************************************************************************/
428
429void TrackList::AddElement(RenderElement* el)
430{
431 static const Exc_t eH("TrackList::AddElement ");
432 if (dynamic_cast<Track*>(el) == 0)
433 throw(eH + "new element not a Track.");
434 RenderElementListBase::AddElement(el);
435}
436
437/**************************************************************************/
438
5a5a1232 439void TrackList::MakeTracks()
440{
9c39ede9 441 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
5a5a1232 442 ((Track*)(*i))->MakeTrack();
443 }
444 gReve->Redraw3D();
445}
446
447
448void TrackList::MakeMarkers()
449{
7d42b6c2 450 Reset(fChildren.size());
9c39ede9 451 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
5a5a1232 452 Track& t = *((Track*)(*i));
453 if(t.GetN() > 0)
454 SetNextPoint(t.fV.x, t.fV.y, t.fV.z);
455 }
456 gReve->Redraw3D();
457}
458
459/**************************************************************************/
460/*************************************************************************/
461
2074deef 462void TrackList::SetWidth(Width_t w)
463{
464 Width_t oldw = fRnrStyle->fWidth;
465 fRnrStyle->fWidth = w;
9c39ede9 466 for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
2074deef 467 Track& t = *((Track*)(*i));
468 if (t.GetLineWidth() == oldw)
469 t.SetLineWidth(w);
470 }
471}
472
5a5a1232 473void TrackList::SetMaxR(Float_t x)
474{
a8600b56 475 fRnrStyle->fMaxR = x;
5a5a1232 476 MakeTracks();
477 MakeMarkers();
478}
479
480void TrackList::SetMaxZ(Float_t x)
481{
a8600b56 482 fRnrStyle->fMaxZ = x;
5a5a1232 483 MakeTracks();
484 MakeMarkers();
485}
486
487void TrackList::SetMaxOrbs(Float_t x)
488{
a8600b56 489 fRnrStyle->fMaxOrbs = x;
5a5a1232 490 MakeTracks();
491}
492
493void TrackList::SetMinAng(Float_t x)
494{
a8600b56 495 fRnrStyle->fMinAng = x;
5a5a1232 496 MakeTracks();
497}
498
499void TrackList::SetDelta(Float_t x)
500{
a8600b56 501 fRnrStyle->fDelta = x;
5a5a1232 502 MakeTracks();
503}
504
505void TrackList::SetFitDaughters(Bool_t x)
506{
a8600b56 507 fRnrStyle->fFitDaughters = x;
5a5a1232 508 MakeTracks();
509}
510
7cbaaebd 511void TrackList::SetFitReferences(Bool_t x)
512{
513 fRnrStyle->fFitReferences = x;
514 MakeTracks();
515}
516
5a5a1232 517void TrackList::SetFitDecay(Bool_t x)
518{
a8600b56 519 fRnrStyle->fFitDecay = x;
5a5a1232 520 MakeTracks();
521}
522
7cbaaebd 523void TrackList::SetRnrDecay(Bool_t rnr)
524{
525 fRnrStyle->fRnrDecay = rnr;
526 MakeTracks();
527}
528
529void TrackList::SetRnrDaughters(Bool_t rnr)
530{
531 fRnrStyle->fRnrDaughters = rnr;
532 MakeTracks();
533}
534
535void TrackList::SetRnrReferences(Bool_t rnr)
536{
537 fRnrStyle->fRnrReferences = rnr;
538 MakeTracks();
539}
540
541void TrackList::SetRnrMarkers(Bool_t rnr)
542{
543 fRnrMarkers = rnr;
544 gReve->Redraw3D();
545}
546
547void TrackList::SetRnrTracks(Bool_t rnr)
548{
549
550 fRnrTracks = rnr;
551 gReve->Redraw3D();
552}
553
5a5a1232 554/**************************************************************************/
555/**************************************************************************/
556
557void TrackList::SelectByPt(Float_t min_pt, Float_t max_pt)
558{
22df2a83 559 fRnrStyle->fMinPt = min_pt;
560 fRnrStyle->fMaxPt = max_pt;
561
5a5a1232 562 Float_t minptsq = min_pt*min_pt;
563 Float_t maxptsq = max_pt*max_pt;
564 Float_t ptsq;
565
9c39ede9 566 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
5a5a1232 567 ptsq = ((Track*)(*i))->fP.Perp2();
568 (*i)->SetRnrElement(ptsq >= minptsq && ptsq <= maxptsq);
569 }
570}
571
572/**************************************************************************/
573
b99aed53 574void TrackList::ImportHits()
5a5a1232 575{
9c39ede9 576 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
b99aed53 577 ((Track*)(*i))->ImportHits();
578 }
5a5a1232 579}
580
b99aed53 581void TrackList::ImportClusters()
5a5a1232 582{
9c39ede9 583 for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
b99aed53 584 ((Track*)(*i))->ImportClusters();
585 }
5a5a1232 586}
9c39ede9 587
588/**************************************************************************/
589/**************************************************************************/
590/**************************************************************************/
591
592#include "RGEditor.h"
593
594//______________________________________________________________________
595// TrackCounter
596//
597
598ClassImp(TrackCounter)
599
b1f08706 600TrackCounter* TrackCounter::fgInstance = 0;
601
9c39ede9 602TrackCounter::TrackCounter(const Text_t* name, const Text_t* title) :
603 RenderElement(),
604 TNamed(name, title),
605
606 fBadLineStyle (6),
607 fClickAction (CA_ToggleTrack),
608 fAllTracks (0),
609 fGoodTracks (0),
610 fTrackLists ()
611{
b1f08706 612 if (fgInstance == 0) fgInstance = this;
9c39ede9 613 TQObject::Connect("Reve::Track", "CtrlClicked(Reve::Track*)",
614 "Reve::TrackCounter", this, "DoTrackAction(Reve::Track*)");
615}
616
617TrackCounter::~TrackCounter()
618{
619 TQObject::Disconnect("Reve::Track", "DoTrackAction(Reve::Track*)");
b1f08706 620 if (fgInstance == this) fgInstance = 0;
9c39ede9 621}
622
623/**************************************************************************/
624
625void TrackCounter::Reset()
626{
627 printf("TrackCounter::Reset()\n");
628 fAllTracks = 0;
629 fGoodTracks = 0;
630 TIter next(&fTrackLists);
631 TrackList* tlist;
632 while ((tlist = dynamic_cast<TrackList*>(next())))
633 tlist->RemoveParent(this);
634 fTrackLists.Clear();
635}
636
637void TrackCounter::RegisterTracks(TrackList* tlist, Bool_t goodTracks)
638{
639 // printf("TrackCounter::RegisterTracks '%s', %s\n",
640 // tlist->GetObject()->GetName(), goodTracks ? "good" : "bad");
641
642 tlist->AddParent(this);
643 fTrackLists.Add(tlist);
644
645 List_i i = tlist->BeginChildren();
646 while (i != tlist->EndChildren())
647 {
648 Track* t = dynamic_cast<Track*>(*i);
649 if (t != 0)
650 {
651 if (goodTracks)
652 {
653 ++fGoodTracks;
654 } else {
655 t->SetLineStyle(fBadLineStyle);
656 }
657 ++fAllTracks;
658 }
659 ++i;
660 }
661}
662
663void TrackCounter::DoTrackAction(Track* track)
664{
665 // !!!! No check done if ok.
666 // !!!! Should also override RemoveElementLocal
7af62191 667 // !!!! But then ... should also sotre local information if track is ok.
9c39ede9 668
669 switch (fClickAction)
670 {
671
672 case CA_PrintTrackInfo:
673 {
674 printf("Track '%s'\n", track->GetObject()->GetName());
675 Vector &v = track->fV, &p = track->fP;
676 printf(" Vx=%f, Vy=%f, Vz=%f; Pt=%f, Pz=%f, phi=%f)\n",
677 v.x, v.y, v.z, p.Perp(), p.z, TMath::RadToDeg()*p.Phi());
678 printf(" <other information should be printed ... full AliESDtrack>\n");
679 break;
680 }
681
682 case CA_ToggleTrack:
683 {
684 if (track->GetLineStyle() == 1)
685 {
686 track->SetLineStyle(fBadLineStyle);
687 --fGoodTracks;
688 } else {
689 track->SetLineStyle(1);
690 ++fGoodTracks;
691 }
692 gReve->Redraw3D();
693
694 printf("TrackCounter::CountTrack All=%d, Good=%d, Bad=%d\n",
695 fAllTracks, fGoodTracks, fAllTracks-fGoodTracks);
696
697 if (gReve->GetEditor()->GetModel() == GetObject())
698 gReve->EditRenderElement(this);
699
700 break;
701 }
702
703 } // end switch fClickAction
704}
7af62191 705
706/**************************************************************************/
707
708void TrackCounter::OutputEventTracks(FILE* out)
709{
710 if (out == 0)
711 {
712 out = stdout;
713 fprintf(out, "TrackCounter::FinalizeEvent()\n");
714 }
715
716 fprintf(out, "Event = %d Ntracks = %d\n", fEventId, fGoodTracks);
717
718 TIter tlists(&fTrackLists);
719 TrackList* tlist;
720 Int_t cnt = 0;
721 while ((tlist = (TrackList*) tlists()) != 0)
722 {
723 List_i i = tlist->BeginChildren();
724 while (i != tlist->EndChildren())
725 {
726 Track* t = dynamic_cast<Track*>(*i);
727 if (t != 0 && t->GetLineStyle() == 1)
728 {
729 ++cnt;
730 fprintf(out, " %2d: chg=%+2d pt=%8.5f eta=%+8.5f\n",
731 cnt, t->fCharge, t->fP.Perp(), t->fP.Eta());
732 }
733 ++i;
734 }
7af62191 735 }
736}