Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TEveTrack.cxx
Go to the documentation of this file.
1// @(#)root/eve:$Id$
2// Authors: Matevz Tadel & Alja Mrak-Tadel: 2006, 2007
3
4/*************************************************************************
5 * Copyright (C) 1995-2007, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "TEveTrack.h"
13#include "TEveTrackPropagator.h"
14#include "TEvePointSet.h"
15
16#include "TParticle.h"
17#include "TPolyMarker3D.h"
18#include "TParticlePDG.h"
19
20// Updates
21#include "TEveManager.h"
22#include "TEveBrowser.h"
23#include "TEveTrackProjected.h"
24#include "TEveVSD.h"
25
26#include <iostream>
27#include <vector>
28#include <algorithm>
29#include <functional>
30
31/** \class TEveTrack
32\ingroup TEve
33Visual representation of a track.
34
35If member fDpDs is set, the momentum is reduced on all path-marks that do
36not fix the momentum according to the distance travelled from the previous
37pathmark.
38*/
39
40
41////////////////////////////////////////////////////////////////////////////////
42/// Default constructor.
43
45 TEveLine(),
46
47 fV(),
48 fP(),
49 fPEnd(),
50 fBeta(0),
51 fDpDs(0),
52 fPdg(0),
53 fCharge(0),
54 fLabel(kMinInt),
55 fIndex(kMinInt),
56 fStatus(0),
57 fLockPoints(kFALSE),
58 fPathMarks(),
59 fLastPMIdx(0),
60 fPropagator(nullptr)
61{
62}
63
64////////////////////////////////////////////////////////////////////////////////
65/// Constructor from TParticle.
66
68 TEveLine(),
69
70 fV(t->Vx(), t->Vy(), t->Vz()),
71 fP(t->Px(), t->Py(), t->Pz()),
72 fPEnd(),
73 fBeta(t->P()/t->Energy()),
74 fDpDs(0),
75 fPdg(0),
76 fCharge(0),
77 fLabel(label),
78 fIndex(kMinInt),
79 fStatus(t->GetStatusCode()),
80 fLockPoints(kFALSE),
81 fPathMarks(),
82 fLastPMIdx(0),
83 fPropagator(nullptr)
84{
87
88 TParticlePDG* pdgp = t->GetPDG();
89 if (pdgp) {
90 fPdg = pdgp->PdgCode();
91 fCharge = (Int_t) TMath::Nint(pdgp->Charge()/3);
92 }
93
94 SetName(t->GetName());
95}
96
97////////////////////////////////////////////////////////////////////////////////
98
100 TEveLine(),
101
102 fV(t->Vx(), t->Vy(), t->Vz()),
103 fP(t->Px(), t->Py(), t->Pz()),
104 fPEnd(),
105 fBeta(t->P()/t->Energy()),
106 fDpDs(0),
107 fPdg(0),
108 fCharge(0),
109 fLabel(t->fLabel),
110 fIndex(t->fIndex),
111 fStatus(t->GetStatusCode()),
112 fLockPoints(kFALSE),
113 fPathMarks(),
114 fLastPMIdx(0),
115 fPropagator(nullptr)
116{
117 // Constructor from TEveUtil Monte Carlo track.
118
121
122 TParticlePDG* pdgp = t->GetPDG();
123 if (pdgp) {
124 fCharge = (Int_t) TMath::Nint(pdgp->Charge()/3);
125 }
126
127 SetName(t->GetName());
128}
129
130////////////////////////////////////////////////////////////////////////////////
131/// Constructor from TEveRecTrack<double> reconstructed track.
132
134 TEveLine(),
135
136 fV(t->fV),
137 fP(t->fP),
138 fPEnd(),
139 fBeta(t->fBeta),
140 fDpDs(0),
141 fPdg(0),
142 fCharge(t->fSign),
143 fLabel(t->fLabel),
144 fIndex(t->fIndex),
145 fStatus(t->fStatus),
146 fLockPoints(kFALSE),
147 fPathMarks(),
148 fLastPMIdx(0),
149 fPropagator(nullptr)
150{
153
154 SetName(t->GetName());
155}
156
157////////////////////////////////////////////////////////////////////////////////
158/// Constructor from TEveRecTrack<float> reconstructed track.
159/// It is recommended to use constructor with TEveRecTrack<double> since
160/// TEveTrackPropagator operates with double type.
161
163 TEveLine(),
164
165 fV(t->fV),
166 fP(t->fP),
167 fPEnd(),
168 fBeta(t->fBeta),
169 fDpDs(0),
170 fPdg(0),
171 fCharge(t->fSign),
172 fLabel(t->fLabel),
173 fIndex(t->fIndex),
174 fStatus(t->fStatus),
175 fLockPoints(kFALSE),
176 fPathMarks(),
177 fLastPMIdx(0),
178 fPropagator(nullptr)
179{
182
183 SetName(t->GetName());
184}
185
186////////////////////////////////////////////////////////////////////////////////
187/// Copy constructor. Track parameters are copied but the
188/// extrapolation is not performed so you should still call
189/// MakeTrack() to do that.
190/// If points of 't' are locked, they are cloned.
191
193 TEveLine(),
194 fV(t.fV),
195 fP(t.fP),
196 fPEnd(),
197 fBeta(t.fBeta),
198 fDpDs(t.fDpDs),
199 fPdg(t.fPdg),
200 fCharge(t.fCharge),
201 fLabel(t.fLabel),
202 fIndex(t.fIndex),
203 fStatus(t.fStatus),
204 fLockPoints(t.fLockPoints),
205 fPathMarks(),
206 fLastPMIdx(t.fLastPMIdx),
207 fPropagator(nullptr)
208{
209 if (fLockPoints)
210 ClonePoints(t);
211
212 SetPathMarks(t);
214
215 CopyVizParams(&t);
216}
217
218////////////////////////////////////////////////////////////////////////////////
219/// Destructor.
220
222{
223 SetPropagator(nullptr);
224}
225
226
227////////////////////////////////////////////////////////////////////////////////
228/// Returns list-tree icon for TEveTrack.
229
234
235////////////////////////////////////////////////////////////////////////////////
236/// Compute the bounding box of the track.
237
239{
240 if (Size() > 0 || ! fPathMarks.empty())
241 {
242 BBoxInit();
243 Int_t n = Size();
245 for (Int_t i = 0; i < n; ++i, p += 3)
246 {
248 }
249 for (vPathMark_ci i = fPathMarks.begin(); i != fPathMarks.end(); ++i)
250 {
251 BBoxCheckPoint(i->fV.fX, i->fV.fY,i->fV.fZ);
252 }
253 }
254 else
255 {
256 BBoxZero();
257 }
258}
259
260////////////////////////////////////////////////////////////////////////////////
261/// Set standard track title based on most data-member values.
262
264{
265 TString idx(fIndex == kMinInt ? "<undef>" : Form("%d", fIndex));
266 TString lbl(fLabel == kMinInt ? "<undef>" : Form("%d", fLabel));
267 SetTitle(Form("Index=%s, Label=%s\nChg=%d, Pdg=%d\n"
268 "pT=%.3f, pZ=%.3f\nV=(%.3f, %.3f, %.3f)",
269 idx.Data(), lbl.Data(), fCharge, fPdg,
270 fP.Perp(), fP.fZ, fV.fX, fV.fY, fV.fZ));
271}
272
273////////////////////////////////////////////////////////////////////////////////
274/// Copy track parameters from t. Track-propagator is set, too.
275/// PathMarks are cleared - you can copy them via SetPathMarks(t).
276/// If track 't' is locked, you should probably clone its points
277/// over - use TEvePointSet::ClonePoints(t);
278
280{
281 fV = t.fV;
282 fP = t.fP;
283 fBeta = t.fBeta;
284 fPdg = t.fPdg;
285 fCharge = t.fCharge;
286 fLabel = t.fLabel;
287 fIndex = t.fIndex;
288
289 fPathMarks.clear();
291}
292
293////////////////////////////////////////////////////////////////////////////////
294/// Copy path-marks from t.
295
297{
298 std::copy(t.RefPathMarks().begin(), t.RefPathMarks().end(),
299 std::back_insert_iterator<vPathMark_t>(fPathMarks));
300}
301
302////////////////////////////////////////////////////////////////////////////////
303/// Set track's render style.
304/// Reference counts of old and new propagator are updated.
305
313
314////////////////////////////////////////////////////////////////////////////////
315/// Set line and marker attributes from TEveTrackList.
316
318{
319 SetRnrLine(tl->GetRnrLine());
320 SetLineColor(tl->GetLineColor());
321 SetLineStyle(tl->GetLineStyle());
322 SetLineWidth(tl->GetLineWidth());
323
324 SetRnrPoints(tl->GetRnrPoints());
325 SetMarkerColor(tl->GetMarkerColor());
326 SetMarkerStyle(tl->GetMarkerStyle());
327 SetMarkerSize(tl->GetMarkerSize());
328}
329
330////////////////////////////////////////////////////////////////////////////////
331/// Calculate track representation based on track data and current
332/// settings of the propagator.
333/// If recurse is true, descend into children.
334
336{
337 if (!fLockPoints)
338 {
339 Reset(0);
340 fLastPMIdx = 0;
341
343
344 const Double_t maxRsq = rTP.GetMaxR() * rTP.GetMaxR();
345 const Double_t maxZ = rTP.GetMaxZ();
346
348 {
351 rTP.InitTrack(fV, fCharge);
352 for (vPathMark_i pm = fPathMarks.begin(); pm != fPathMarks.end(); ++pm, ++fLastPMIdx)
353 {
354 Int_t start_point = rTP.GetCurrentPoint();
355
356 if (rTP.GetFitReferences() && pm->fType == TEvePathMarkD::kReference)
357 {
359 break;
360 if (rTP.GoToVertex(pm->fV, currP))
361 {
362 currP.fX = pm->fP.fX; currP.fY = pm->fP.fY; currP.fZ = pm->fP.fZ;
363 }
364 else
365 {
366 break;
367 }
368 }
369 else if (rTP.GetFitDaughters() && pm->fType == TEvePathMarkD::kDaughter)
370 {
372 break;
373 if (rTP.GoToVertex(pm->fV, currP))
374 {
375 currP.fX -= pm->fP.fX; currP.fY -= pm->fP.fY; currP.fZ -= pm->fP.fZ;
376 if (fDpDs != 0)
377 {
378 Double_t dp = fDpDs * rTP.GetTrackLength(start_point);
379 Double_t p = currP.Mag();
380 if (p > dp) currP *= 1.0 - dp / p;
381 }
382 }
383 else
384 {
385 break;
386 }
387 }
388 else if (rTP.GetFitDecay() && pm->fType == TEvePathMarkD::kDecay)
389 {
391 break;
392 rTP.GoToVertex(pm->fV, currP);
393 decay = kTRUE;
394 ++fLastPMIdx;
395 break;
396 }
397 else if (rTP.GetFitCluster2Ds() && pm->fType == TEvePathMarkD::kCluster2D)
398 {
400 if (rTP.IntersectPlane(currP, pm->fV, pm->fP, itsect))
401 {
402 TEveVectorD delta = itsect - pm->fV;
403 TEveVectorD vtopass = pm->fV + pm->fE*(pm->fE.Dot(delta));
405 break;
406 if ( ! rTP.GoToVertex(vtopass, currP))
407 break;
408
409 if (fDpDs != 0)
410 {
411 Double_t dp = fDpDs * rTP.GetTrackLength(start_point);
412 Double_t p = currP.Mag();
413 if (p > dp) currP *= 1.0 - dp / p;
414 }
415 }
416 else
417 {
418 Warning("TEveTrack::MakeTrack", "Failed to intersect plane for Cluster2D. Ignoring path-mark.");
419 }
420 }
421 else if (rTP.GetFitLineSegments() && pm->fType == TEvePathMarkD::kLineSegment)
422 {
424 break;
425
426 if (rTP.GoToLineSegment(pm->fV, pm->fE, currP))
427 {
428 if (fDpDs != 0)
429 {
430 Double_t dp = fDpDs * rTP.GetTrackLength(start_point);
431 Double_t p = currP.Mag();
432 if (p > dp) currP *= 1.0 - dp / p;
433 }
434 }
435 else
436 {
437 break;
438 }
439 }
440 else
441 {
443 break;
444 }
445 } // loop path-marks
446
447 if (!decay)
448 {
449 // printf("%s loop to bounds \n",fName.Data() );
450 rTP.GoToBounds(currP);
451 }
452 fPEnd = currP;
453 // make_polyline:
454 rTP.FillPointSet(this);
455 rTP.ResetTrack();
456 }
457 }
458
459 if (recurse)
460 {
461 for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i)
462 {
463 TEveTrack* t = dynamic_cast<TEveTrack*>(*i);
464 if (t) t->MakeTrack(recurse);
465 }
466 }
467}
468
469////////////////////////////////////////////////////////////////////////////////
470/// Copy visualization parameters from element el.
471
473{
474 // No local parameters.
475
476 // const TEveTrack* t = dynamic_cast<const TEveTrack*>(el);
477 // if (t)
478 // {}
479
481}
482
483////////////////////////////////////////////////////////////////////////////////
484/// Write visualization parameters.
485
486void TEveTrack::WriteVizParams(std::ostream& out, const TString& var)
487{
488 TEveLine::WriteVizParams(out, var);
489
490 // TString t = " " + var + "->";
491}
492
493////////////////////////////////////////////////////////////////////////////////
494/// Virtual from TEveProjectable, return TEveTrackProjected class.
495
500
501namespace
502{
503 struct Cmp_pathmark_t
504 {
505 bool operator()(TEvePathMarkD const & a, TEvePathMarkD const & b)
506 { return a.fTime < b.fTime; }
507 };
508}
509
510////////////////////////////////////////////////////////////////////////////////
511/// Sort registered pat-marks by time.
512
514{
515 std::sort(fPathMarks.begin(), fPathMarks.end(), Cmp_pathmark_t());
516}
517
518////////////////////////////////////////////////////////////////////////////////
519/// Print registered path-marks.
520
522{
523 static const TEveException eh("TEveTrack::PrintPathMarks ");
524
525 printf("TEveTrack '%s', number of path marks %d, label %d\n",
526 GetName(), (Int_t)fPathMarks.size(), fLabel);
527
528 for (vPathMark_i pm = fPathMarks.begin(); pm != fPathMarks.end(); ++pm)
529 {
530 printf(" %-9s p: %8f %8f %8f Vertex: %8e %8e %8e %g Extra:%8f %8f %8f\n",
531 pm->TypeName(),
532 pm->fP.fX, pm->fP.fY, pm->fP.fZ,
533 pm->fV.fX, pm->fV.fY, pm->fV.fZ,
534 pm->fE.fX, pm->fE.fY, pm->fE.fZ,
535 pm->fTime);
536 }
537}
538
539//------------------------------------------------------------------------------
540
541////////////////////////////////////////////////////////////////////////////////
542/// Emits "SecSelected(TEveTrack*)" signal.
543/// Called from TEveTrackGL on secondary-selection.
544
546{
547 Emit("SecSelected(TEveTrack*)", (Longptr_t)track);
548}
549
550/** \class TEveTrackList
551\ingroup TEve
552A list of tracks supporting change of common attributes and
553selection based on track parameters.
554*/
555
556
557////////////////////////////////////////////////////////////////////////////////
558/// Constructor. If track-propagator argument is 0, a new default
559/// one is created.
560
563 TAttMarker(1, 20, 1),
564 TAttLine(1,1,1),
565
566 fPropagator(nullptr),
567 fRecurse(kTRUE),
568 fRnrLine(kTRUE),
569 fRnrPoints(kFALSE),
570
571 fMinPt (0), fMaxPt (0), fLimPt (0),
572 fMinP (0), fMaxP (0), fLimP (0)
573{
574
575 fChildClass = TEveTrack::Class(); // override member from base TEveElementList
576
578
579 if (prop == nullptr) prop = new TEveTrackPropagator;
581}
582
583////////////////////////////////////////////////////////////////////////////////
584/// Constructor. If track-propagator argument is 0, a new default
585/// one is created.
586
589 TAttMarker(1, 20, 1),
590 TAttLine(1,1,1),
591
592 fPropagator(nullptr),
593 fRecurse(kTRUE),
594 fRnrLine(kTRUE),
595 fRnrPoints(kFALSE),
596
597 fMinPt (0), fMaxPt (0), fLimPt (0),
598 fMinP (0), fMaxP (0), fLimP (0)
599{
600 fChildClass = TEveTrack::Class(); // override member from base TEveElementList
601
603
604 if (prop == nullptr) prop = new TEveTrackPropagator;
606}
607
608////////////////////////////////////////////////////////////////////////////////
609/// Destructor.
610
615
616////////////////////////////////////////////////////////////////////////////////
617/// Set default propagator for tracks.
618/// This is not enforced onto the tracks themselves but this is the
619/// propagator that is shown in the TEveTrackListEditor.
620
628
629////////////////////////////////////////////////////////////////////////////////
630/// Regenerate the visual representations of tracks.
631/// The momentum limits are rescanned during the same traversal.
632
634{
635 fLimPt = fLimP = 0;
636
637 for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i)
638 {
639 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
640 if (track)
641 {
642 track->MakeTrack(recurse);
643
644 fLimPt = TMath::Max(fLimPt, track->fP.Perp());
645 fLimP = TMath::Max(fLimP, track->fP.Mag());
646 }
647 if (recurse)
649 }
650
653
655}
656
657////////////////////////////////////////////////////////////////////////////////
658/// Loop over children and find highest pT and p of contained TEveTracks.
659/// These are stored in members fLimPt and fLimP.
660
662{
663 fLimPt = fLimP = 0;
664
665 if (HasChildren())
666 {
667 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
668 {
669 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
670 if (track)
671 {
672 fLimPt = TMath::Max(fLimPt, track->fP.Perp());
673 fLimP = TMath::Max(fLimP, track->fP.Mag());
674 }
675 if (recurse)
677 }
678
681 }
682
684}
685
686////////////////////////////////////////////////////////////////////////////////
687/// Loop over track elements of argument el and find highest pT and p.
688/// These are stored in members fLimPt and fLimP.
689
691{
692 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
693 {
694 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
695 if (track)
696 {
697 fLimPt = TMath::Max(fLimPt, track->fP.Perp());
698 fLimP = TMath::Max(fLimP, track->fP.Mag());
699 }
700 if (recurse)
702 }
703}
704
705////////////////////////////////////////////////////////////////////////////////
706/// Round the momentum limit up to a nice value.
707
709{
710 using namespace TMath;
711
712 if (x < 1e-3) return 1e-3;
713
714 Double_t fac = Power(10, 1 - Floor(Log10(x)));
715 return Ceil(fac*x) / fac;
716}
717
718////////////////////////////////////////////////////////////////////////////////
719/// Set Min/Max cuts so that they are within detected limits.
720
722{
723 using namespace TMath;
724
725 fMinPt = Min(fMinPt, fLimPt);
726 fMaxPt = fMaxPt == 0 ? fLimPt : Min(fMaxPt, fLimPt);
727 fMinP = Min(fMinP, fLimP);
728 fMaxP = fMaxP == 0 ? fLimP : Min(fMaxP, fLimP);
729}
730
731////////////////////////////////////////////////////////////////////////////////
732/// Set rendering of track as line for the list and the elements.
733
735{
736 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
737 {
738 TEveTrack* track = (TEveTrack*)(*i);
739 if (track->GetRnrLine() == fRnrLine)
740 track->SetRnrLine(rnr);
741 if (fRecurse)
742 SetRnrLine(rnr, *i);
743 }
744 fRnrLine = rnr;
745}
746
747////////////////////////////////////////////////////////////////////////////////
748/// Set rendering of track as line for children of el.
749
751{
753 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
754 {
755 track = dynamic_cast<TEveTrack*>(*i);
756 if (track && (track->GetRnrLine() == fRnrLine))
757 track->SetRnrLine(rnr);
758 if (fRecurse)
759 SetRnrLine(rnr, *i);
760 }
761}
762
763////////////////////////////////////////////////////////////////////////////////
764/// Set rendering of track as points for the list and the elements.
765
767{
768 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
769 {
770 TEveTrack* track = (TEveTrack*)(*i);
771 if (track->GetRnrPoints() == fRnrPoints)
772 track->SetRnrPoints(rnr);
773 if (fRecurse)
774 SetRnrPoints(rnr, *i);
775 }
776 fRnrPoints = rnr;
777}
778
779////////////////////////////////////////////////////////////////////////////////
780/// Set rendering of track as points for children of el.
781
783{
785 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
786 {
787 track = dynamic_cast<TEveTrack*>(*i);
788 if (track)
789 if (track->GetRnrPoints() == fRnrPoints)
790 track->SetRnrPoints(rnr);
791 if (fRecurse)
792 SetRnrPoints(rnr, *i);
793 }
794}
795
796////////////////////////////////////////////////////////////////////////////////
797/// Set main (line) color for the list and the elements.
798
800{
801 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
802 {
803 TEveTrack* track = (TEveTrack*)(*i);
804 if (track->GetLineColor() == fLineColor)
805 track->SetLineColor(col);
806 if (fRecurse)
807 SetLineColor(col, *i);
808 }
810}
811
812////////////////////////////////////////////////////////////////////////////////
813/// Set line color for children of el.
814
816{
818 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
819 {
820 track = dynamic_cast<TEveTrack*>(*i);
821 if (track && track->GetLineColor() == fLineColor)
822 track->SetLineColor(col);
823 if (fRecurse)
824 SetLineColor(col, *i);
825 }
826}
827
828////////////////////////////////////////////////////////////////////////////////
829/// Set line width for the list and the elements.
830
832{
833 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
834 {
835 TEveTrack* track = (TEveTrack*)(*i);
836 if (track->GetLineWidth() == fLineWidth)
837 track->SetLineWidth(width);
838 if (fRecurse)
839 SetLineWidth(width, *i);
840 }
842}
843
844////////////////////////////////////////////////////////////////////////////////
845/// Set line width for children of el.
846
848{
850 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
851 {
852 track = dynamic_cast<TEveTrack*>(*i);
853 if (track && track->GetLineWidth() == fLineWidth)
854 track->SetLineWidth(width);
855 if (fRecurse)
856 SetLineWidth(width, *i);
857 }
858}
859
860////////////////////////////////////////////////////////////////////////////////
861/// Set line style for the list and the elements.
862
864{
865 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
866 {
867 TEveTrack* track = (TEveTrack*)(*i);
868 if (track->GetLineStyle() == fLineStyle)
869 track->SetLineStyle(style);
870 if (fRecurse)
871 SetLineStyle(style, *i);
872 }
874}
875
876////////////////////////////////////////////////////////////////////////////////
877/// Set line style for children of el.
878
880{
882 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
883 {
884 track = dynamic_cast<TEveTrack*>(*i);
885 if (track && track->GetLineStyle() == fLineStyle)
886 track->SetLineStyle(style);
887 if (fRecurse)
888 SetLineStyle(style, *i);
889 }
890}
891
892////////////////////////////////////////////////////////////////////////////////
893/// Set marker style for the list and the elements.
894
896{
897 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
898 {
899 TEveTrack* track = (TEveTrack*)(*i);
900 if (track->GetMarkerStyle() == fMarkerStyle)
901 track->SetMarkerStyle(style);
902 if (fRecurse)
904 }
906}
907
908////////////////////////////////////////////////////////////////////////////////
909/// Set marker style for children of el.
910
912{
914 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
915 {
916 track = dynamic_cast<TEveTrack*>(*i);
917 if (track && track->GetMarkerStyle() == fMarkerStyle)
918 track->SetMarkerStyle(style);
919 if (fRecurse)
921 }
922}
923
924////////////////////////////////////////////////////////////////////////////////
925/// Set marker color for the list and the elements.
926
928{
929 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
930 {
931 TEveTrack* track = (TEveTrack*)(*i);
932 if (track->GetMarkerColor() == fMarkerColor)
933 track->SetMarkerColor(col);
934 if (fRecurse)
935 SetMarkerColor(col, *i);
936 }
937 fMarkerColor = col;
938}
939
940////////////////////////////////////////////////////////////////////////////////
941/// Set marker color for children of el.
942
944{
946 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
947 {
948 track = dynamic_cast<TEveTrack*>(*i);
949 if (track && track->GetMarkerColor() == fMarkerColor)
950 track->SetMarkerColor(col);
951 if (fRecurse)
952 SetMarkerColor(col, *i);
953 }
954}
955
956////////////////////////////////////////////////////////////////////////////////
957/// Set marker size for the list and the elements.
958
960{
961 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
962 {
963 TEveTrack* track = (TEveTrack*)(*i);
964 if (track->GetMarkerSize() == fMarkerSize)
965 track->SetMarkerSize(size);
966 if (fRecurse)
967 SetMarkerSize(size, *i);
968 }
970}
971
972////////////////////////////////////////////////////////////////////////////////
973/// Set marker size for children of el.
974
976{
978 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
979 {
980 track = dynamic_cast<TEveTrack*>(*i);
981 if (track && track->GetMarkerSize() == fMarkerSize)
982 track->SetMarkerSize(size);
983 if (fRecurse)
984 SetMarkerSize(size, *i);
985 }
986}
987
988////////////////////////////////////////////////////////////////////////////////
989/// Select visibility of tracks by transverse momentum.
990/// If data-member fRecurse is set, the selection is applied
991/// recursively to all children.
992
994{
995 fMinPt = min_pt;
996 fMaxPt = max_pt;
997
1000
1001 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
1002 {
1003 const Double_t ptsq = ((TEveTrack*)(*i))->fP.Perp2();
1004 Bool_t on = ptsq >= minptsq && ptsq <= maxptsq;
1005 (*i)->SetRnrState(on);
1006 if (on && fRecurse)
1007 SelectByPt(min_pt, max_pt, *i);
1008 }
1009}
1010
1011////////////////////////////////////////////////////////////////////////////////
1012/// Select visibility of el's children tracks by transverse momentum.
1013
1015{
1016 const Double_t minptsq = min_pt*min_pt;
1017 const Double_t maxptsq = max_pt*max_pt;
1018
1019 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
1020 {
1021 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
1022 if (track)
1023 {
1024 const Double_t ptsq = track->fP.Perp2();
1025 Bool_t on = ptsq >= minptsq && ptsq <= maxptsq;
1026 track->SetRnrState(on);
1027 if (on && fRecurse)
1028 SelectByPt(min_pt, max_pt, *i);
1029 }
1030 }
1031}
1032
1033////////////////////////////////////////////////////////////////////////////////
1034/// Select visibility of tracks by momentum.
1035/// If data-member fRecurse is set, the selection is applied
1036/// recursively to all children.
1037
1039{
1040 fMinP = min_p;
1041 fMaxP = max_p;
1042
1043 const Double_t minpsq = min_p*min_p;
1044 const Double_t maxpsq = max_p*max_p;
1045
1046 for (List_i i=BeginChildren(); i!=EndChildren(); ++i)
1047 {
1048 const Double_t psq = ((TEveTrack*)(*i))->fP.Mag2();
1049 Bool_t on = psq >= minpsq && psq <= maxpsq;
1050 (*i)->SetRnrState(psq >= minpsq && psq <= maxpsq);
1051 if (on && fRecurse)
1052 SelectByP(min_p, max_p, *i);
1053 }
1054}
1055
1056////////////////////////////////////////////////////////////////////////////////
1057/// Select visibility of el's children tracks by momentum.
1058
1060{
1061 const Double_t minpsq = min_p*min_p;
1062 const Double_t maxpsq = max_p*max_p;
1063
1064 for (List_i i=el->BeginChildren(); i!=el->EndChildren(); ++i)
1065 {
1066 TEveTrack* track = dynamic_cast<TEveTrack*>(*i);
1067 if (track)
1068 {
1069 const Double_t psq = ((TEveTrack*)(*i))->fP.Mag2();
1070 Bool_t on = psq >= minpsq && psq <= maxpsq;
1071 track->SetRnrState(on);
1072 if (on && fRecurse)
1073 SelectByP(min_p, max_p, *i);
1074 }
1075 }
1076}
1077
1078////////////////////////////////////////////////////////////////////////////////
1079/// Find track by label, select it and display it in the editor.
1080
1082{
1083 for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i)
1084 {
1085 if (((TEveTrack*)(*i))->GetLabel() == label)
1086 {
1088 TGListTreeItem *mlti = lt->GetSelected();
1089 if (mlti->GetUserData() != this)
1091 TGListTreeItem *tlti = (*i)->FindListTreeItem(lt, mlti);
1092 lt->HighlightItem(tlti);
1093 lt->SetSelected(tlti);
1094 gEve->EditElement(*i);
1095 return (TEveTrack*) *i;
1096 }
1097 }
1098 return nullptr;
1099}
1100
1101////////////////////////////////////////////////////////////////////////////////
1102/// Find track by index, select it and display it in the editor.
1103
1105{
1106 for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i)
1107 {
1108 if (((TEveTrack*)(*i))->GetIndex() == index)
1109 {
1111 TGListTreeItem *mlti = lt->GetSelected();
1112 if (mlti->GetUserData() != this)
1114 TGListTreeItem *tlti = (*i)->FindListTreeItem(lt, mlti);
1115 lt->HighlightItem(tlti);
1116 lt->SetSelected(tlti);
1117 gEve->EditElement(*i);
1118 return (TEveTrack*) *i;
1119 }
1120 }
1121 return nullptr;
1122}
1123
1124////////////////////////////////////////////////////////////////////////////////
1125/// Copy visualization parameters from element el.
1126
1128{
1129 const TEveTrackList* m = dynamic_cast<const TEveTrackList*>(el);
1130 if (m)
1131 {
1132 TAttMarker::operator=(*m);
1133 TAttLine::operator=(*m);
1134 fRecurse = m->fRecurse;
1135 fRnrLine = m->fRnrLine;
1136 fRnrPoints = m->fRnrPoints;
1137 fMinPt = m->fMinPt;
1138 fMaxPt = m->fMaxPt;
1139 fLimPt = m->fLimPt;
1140 fMinP = m->fMinP;
1141 fMaxP = m->fMaxP;
1142 fLimP = m->fLimP;
1143 }
1144
1146}
1147
1148////////////////////////////////////////////////////////////////////////////////
1149/// Write visualization parameters.
1150
1151void TEveTrackList::WriteVizParams(std::ostream& out, const TString& var)
1152{
1154
1155 TString t = " " + var + "->";
1157 TAttLine ::SaveLineAttributes (out, var);
1158 out << t << "SetRecurse(" << ToString(fRecurse) << ");\n";
1159 out << t << "SetRnrLine(" << ToString(fRnrLine) << ");\n";
1160 out << t << "SetRnrPoints(" << ToString(fRnrPoints) << ");\n";
1161 // These setters are not available -- need proper AND/OR mode.
1162 // out << t << "SetMinPt(" << fMinPt << ");\n";
1163 // out << t << "SetMaxPt(" << fMaxPt << ");\n";
1164 // out << t << "SetLimPt(" << fLimPt << ");\n";
1165 // out << t << "SetMinP(" << fMinP << ");\n";
1166 // out << t << "SetMaxP(" << fMaxP << ");\n";
1167 // out << t << "SetLimP(" << fLimP << ");\n";
1168}
1169
1170////////////////////////////////////////////////////////////////////////////////
1171/// Virtual from TEveProjectable, returns TEveTrackListProjected class.
1172
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
#define e(i)
Definition RSha256.hxx:103
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
short Style_t
Style number (short)
Definition RtypesCore.h:96
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
short Color_t
Color number (short)
Definition RtypesCore.h:99
float Size_t
Attribute size (float)
Definition RtypesCore.h:103
long Longptr_t
Integer large enough to hold a pointer (platform-dependent)
Definition RtypesCore.h:89
constexpr Int_t kMinInt
Definition RtypesCore.h:120
short Width_t
Line width (short)
Definition RtypesCore.h:98
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
void Warning(const char *location, const char *msgfmt,...)
Use this function in warning situations.
Definition TError.cxx:252
R__EXTERN TEveManager * gEve
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t SetLineWidth
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h prop
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void on
Option_t Option_t SetLineColor
Option_t Option_t SetMarkerStyle
Option_t Option_t width
Option_t Option_t style
char name[80]
Definition TGX11.cxx:110
TRObject operator()(const T1 &t1) const
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2495
void BBoxCheckPoint(Float_t x, Float_t y, Float_t z)
Definition TAttBBox.h:69
void BBoxZero(Float_t epsilon=0, Float_t x=0, Float_t y=0, Float_t z=0)
Create cube of volume (2*epsilon)^3 at (x,y,z).
Definition TAttBBox.cxx:41
void BBoxInit(Float_t infinity=1e6)
Dynamic Float_t[6] X(min,max), Y(min,max), Z(min,max)
Definition TAttBBox.cxx:28
Line Attributes class.
Definition TAttLine.h:20
Width_t fLineWidth
Line width.
Definition TAttLine.h:25
Style_t fLineStyle
Line style.
Definition TAttLine.h:24
Color_t fLineColor
Line color.
Definition TAttLine.h:23
Marker Attributes class.
Definition TAttMarker.h:20
virtual void SaveMarkerAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1, Int_t sizdef=1)
Save line attributes as C++ statement(s) on output stream out.
Color_t fMarkerColor
Marker color.
Definition TAttMarker.h:23
Size_t fMarkerSize
Marker size.
Definition TAttMarker.h:25
Style_t fMarkerStyle
Marker style.
Definition TAttMarker.h:24
TClass instances represent classes, structs and namespaces in the ROOT type system.
Definition TClass.h:84
A list of TEveElements.
TClass * fChildClass
Base class for TEveUtil visualization elements, providing hierarchy management, rendering control and...
Definition TEveElement.h:36
List_t fChildren
Definition TEveElement.h:81
List_i EndChildren()
Bool_t HasChildren() const
virtual TGListTreeItem * FindListTreeItem(TGListTree *ltree)
Find any list-tree-item of this element in list-tree 'ltree'.
virtual void SetMainColor(Color_t color)
Set main color of the element.
Color_t * fMainColorPtr
Definition TEveElement.h:99
virtual void CopyVizParams(const TEveElement *el)
Copy visualization parameters from element el.
static const char * ToString(Bool_t b)
Convert Bool_t to string - kTRUE or kFALSE.
static const TGPicture * fgListTreeIcons[9]
Definition TEveElement.h:54
virtual void WriteVizParams(std::ostream &out, const TString &var)
Write-out visual parameters for this object.
List_i BeginChildren()
List_t::iterator List_i
Definition TEveElement.h:72
Exception class thrown by TEve classes and macros.
Definition TEveUtil.h:102
TGListTree * GetListTree() const
An arbitrary polyline with fixed line and marker attributes.
Definition TEveLine.h:26
void SetLineStyle(Style_t lstyle) override
Set line-style of the line.
Definition TEveLine.cxx:85
void SetRnrLine(Bool_t r)
Set rendering of line. Propagate to projected lines.
Definition TEveLine.cxx:124
void WriteVizParams(std::ostream &out, const TString &var) override
Write visualization parameters.
Definition TEveLine.cxx:276
void CopyVizParams(const TEveElement *el) override
Copy visualization parameters from element el.
Definition TEveLine.cxx:259
void SetMarkerColor(Color_t col) override
Set marker color. Propagate to projected lines.
Definition TEveLine.cxx:65
void SetRnrPoints(Bool_t r)
Set rendering of points. Propagate to projected lines.
Definition TEveLine.cxx:143
TEveGListTreeEditorFrame * GetLTEFrame() const
void EditElement(TEveElement *element)
Show element in default editor.
virtual void ClonePoints(const TEvePointSet &e)
Clone points and all point-related information from point-set 'e'.
virtual void SetTitle(const char *t)
void Reset(Int_t n_points=0, Int_t n_int_ids=0)
Drop all data and set-up the data structures to recive new data.
void SetMarkerSize(Size_t msize=1) override
Set marker size, propagate to projecteds.
Base-class for non-linear projections.
virtual void IncRefCount(TEveElement *re)
Increase reference count and add re to the list of back-references.
Definition TEveUtil.cxx:571
virtual void DecRefCount(TEveElement *re)
Decrease reference count and remove re from the list of back-references.
Definition TEveUtil.cxx:580
static TClass * Class()
A list of tracks supporting change of common attributes and selection based on track parameters.
Definition TEveTrack.h:140
Double_t fMinPt
Definition TEveTrack.h:155
TClass * ProjectedClass(const TEveProjection *p) const override
Virtual from TEveProjectable, returns TEveTrackListProjected class.
void SetPropagator(TEveTrackPropagator *prop)
Set default propagator for tracks.
void SelectByP(Double_t min_p, Double_t max_p)
Select visibility of tracks by momentum.
void SetRnrPoints(Bool_t r)
Set rendering of track as points for the list and the elements.
void SetMarkerStyle(Style_t s) override
Set marker style for the list and the elements.
Double_t fMaxPt
Definition TEveTrack.h:156
void SanitizeMinMaxCuts()
Set Min/Max cuts so that they are within detected limits.
Bool_t fRecurse
Definition TEveTrack.h:150
void SetLineStyle(Style_t s) override
Set line style for the list and the elements.
Double_t fMaxP
Definition TEveTrack.h:159
Bool_t fRnrPoints
Definition TEveTrack.h:153
void SetMarkerColor(Color_t c) override
Set marker color for the list and the elements.
~TEveTrackList() override
Destructor.
TEveTrackPropagator * fPropagator
Definition TEveTrack.h:148
void SetLineWidth(Width_t w) override
Set line width for the list and the elements.
TEveTrack * FindTrackByIndex(Int_t index)
Find track by index, select it and display it in the editor.
void SetMainColor(Color_t c) override
Set main (line) color for the list and the elements.
TEveTrackList(const TEveTrackList &)
void WriteVizParams(std::ostream &out, const TString &var) override
Write visualization parameters.
void SetLineColor(Color_t c) override
Set the line color.
Definition TEveTrack.h:183
void SelectByPt(Double_t min_pt, Double_t max_pt)
Select visibility of tracks by transverse momentum.
TEveTrack * FindTrackByLabel(Int_t label)
Find track by label, select it and display it in the editor.
void SetMarkerSize(Size_t s) override
Set marker size for the list and the elements.
void MakeTracks(Bool_t recurse=kTRUE)
Regenerate the visual representations of tracks.
void CopyVizParams(const TEveElement *el) override
Copy visualization parameters from element el.
Double_t fLimP
Definition TEveTrack.h:160
Double_t RoundMomentumLimit(Double_t x)
Round the momentum limit up to a nice value.
Double_t fLimPt
Definition TEveTrack.h:157
Bool_t fRnrLine
Definition TEveTrack.h:152
Double_t fMinP
Definition TEveTrack.h:158
void SetRnrLine(Bool_t rnr)
Set rendering of track as line for the list and the elements.
void FindMomentumLimits(TEveElement *el, Bool_t recurse=kTRUE)
Loop over track elements of argument el and find highest pT and p.
static TClass * Class()
Holding structure for a number of track rendering parameters.
static Bool_t IsOutsideBounds(const TEveVectorD &point, Double_t maxRsqr, Double_t maxZ)
static TEveTrackPropagator fgDefault
Visual representation of a track.
Definition TEveTrack.h:33
vPathMark_t::iterator vPathMark_i
Definition TEveTrack.h:43
void ComputeBBox() override
Compute the bounding box of the track.
virtual void SecSelected(TEveTrack *)
Emits "SecSelected(TEveTrack*)" signal.
Int_t fCharge
Definition TEveTrack.h:56
void SetPropagator(TEveTrackPropagator *prop)
Set track's render style.
Int_t fPdg
Definition TEveTrack.h:55
vPathMark_t & RefPathMarks()
Definition TEveTrack.h:111
virtual void MakeTrack(Bool_t recurse=kTRUE)
Calculate track representation based on track data and current settings of the propagator.
TEveTrackPropagator * fPropagator
Last path-mark index tried in track-propagation.
Definition TEveTrack.h:64
TEveVectorD fP
Definition TEveTrack.h:51
vPathMark_t::const_iterator vPathMark_ci
Definition TEveTrack.h:44
TClass * ProjectedClass(const TEveProjection *p) const override
Virtual from TEveProjectable, return TEveTrackProjected class.
void SetAttLineAttMarker(TEveTrackList *tl)
Set line and marker attributes from TEveTrackList.
const TGPicture * GetListTreeIcon(Bool_t open=kFALSE) override
Returns list-tree icon for TEveTrack.
Int_t fIndex
Definition TEveTrack.h:58
virtual void SetPathMarks(const TEveTrack &t)
Copy path-marks from t.
static TClass * Class()
void SortPathMarksByTime()
Sort registered pat-marks by time.
void PrintPathMarks()
Print registered path-marks.
~TEveTrack() override
Destructor.
TEveTrack()
Default constructor.
Definition TEveTrack.cxx:44
Double_t fDpDs
Definition TEveTrack.h:54
void CopyVizParams(const TEveElement *el) override
Copy visualization parameters from element el.
Int_t fLabel
Definition TEveTrack.h:57
Bool_t fLockPoints
Definition TEveTrack.h:60
virtual void SetTrackParams(const TEveTrack &t)
Copy track parameters from t.
Double_t fBeta
Definition TEveTrack.h:53
TEveVectorD fV
Definition TEveTrack.h:50
virtual void SetStdTitle()
Set standard track title based on most data-member values.
vPathMark_t fPathMarks
Definition TEveTrack.h:61
void WriteVizParams(std::ostream &out, const TString &var) override
Write visualization parameters.
Int_t fLastPMIdx
Definition TEveTrack.h:62
TEveVectorD fPEnd
Definition TEveTrack.h:52
TT Perp() const
Definition TEveVector.h:102
A list tree is a widget that can contain a number of items arranged in a tree structure.
Definition TGListTree.h:197
The TGPicture class implements pictures and icons used in the different GUI elements and widgets.
Definition TGPicture.h:25
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:457
Description of the static properties of a particle.
Description of the dynamic properties of a particle.
Definition TParticle.h:26
TParticlePDG * GetPDG(Int_t mode=0) const
Returns a pointer to the TParticlePDG object using the pdgcode.
const char * GetName() const override
Return particle name.
const char * GetName() const override
Returns name of object.
virtual Int_t Size() const
void Emit(const char *signal, const T &arg)
Activate signal with single parameter.
Definition TQObject.h:164
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
Double_t x[n]
Definition legend1.C:17
gr SetName("gr")
const Int_t n
Definition legend1.C:16
TMath.
Definition TMathBase.h:35
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
Definition TMath.h:704
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
TMarker m
Definition textangle.C:8