Logo ROOT   master
Reference Guide
event_demo.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_eve7
3 /// This example display geometry, tracks and hits in web browser
4 ///
5 /// \macro_code
6 ///
7 
8 #include <vector>
9 #include <string>
10 #include <iostream>
11 
12 #include "TClass.h"
13 #include "TRandom.h"
14 #include "TGeoTube.h"
15 #include "TGeoSphere.h"
16 #include "TParticle.h"
17 #include "TApplication.h"
18 #include "TMatrixDSym.h"
19 #include "TVector.h"
20 #include "TMatrixDEigen.h"
21 
22 #include <ROOT/REveGeoShape.hxx>
23 #include <ROOT/REveScene.hxx>
24 #include <ROOT/REveViewer.hxx>
25 #include <ROOT/REveElement.hxx>
26 #include <ROOT/REveManager.hxx>
27 #include <ROOT/REveUtil.hxx>
28 #include <ROOT/REveGeoShape.hxx>
31 #include <ROOT/REvePointSet.hxx>
32 #include <ROOT/REveJetCone.hxx>
33 #include <ROOT/REveTrans.hxx>
34 
35 #include <ROOT/REveTrack.hxx>
37 #include <ROOT/REveEllipsoid.hxx>
38 
39 namespace REX = ROOT::Experimental;
40 
41 // globals
49 
50 const Double_t kR_min = 240;
51 const Double_t kR_max = 250;
52 const Double_t kZ_d = 300;
53 
54 
55 REX::REvePointSet *getPointSet(int npoints = 2, float s=2, int color=28)
56 {
57  TRandom &r = *gRandom;
58 
59  auto ps = new REX::REvePointSet("fu", "", npoints);
60 
61  for (Int_t i=0; i<npoints; ++i)
62  ps->SetNextPoint(r.Uniform(-s,s), r.Uniform(-s,s), r.Uniform(-s,s));
63 
64  ps->SetMarkerColor(color);
65  ps->SetMarkerSize(3+r.Uniform(1, 7));
66  ps->SetMarkerStyle(4);
67  return ps;
68 }
69 
70 void addPoints()
71 {
73 
74  auto pntHolder = new REX::REveElement("Hits");
75 
76  auto ps1 = getPointSet(20, 100);
77  ps1->SetName("Points_1");
78  ps1->SetTitle("Points_1 title"); // used as tooltip
79 
80  pntHolder->AddElement(ps1);
81 
82  auto ps2 = getPointSet(10, 200, 4);
83  ps2->SetName("Points_2");
84  ps2->SetTitle("Points_2 title"); // used as tooltip
85  pntHolder->AddElement(ps2);
86 
87  event->AddElement(pntHolder);
88 }
89 
90 void addTracks()
91 {
92  TRandom &r = *gRandom;
93 
95  auto prop = new REX::REveTrackPropagator();
96  prop->SetMagFieldObj(new REX::REveMagFieldDuo(350, -3.5, 2.0));
97  prop->SetMaxR(300);
98  prop->SetMaxZ(600);
99  prop->SetMaxOrbs(6);
100 
101  auto trackHolder = new REX::REveElement("Tracks");
102 
103  double v = 0.2;
104  double m = 5;
105 
106  int N_Tracks = 10 + r.Integer(20);
107  for (int i = 0; i < N_Tracks; i++)
108  {
109  TParticle* p = new TParticle();
110 
111  int pdg = 11* (r.Integer(2) -1);
112  p->SetPdgCode(pdg);
113 
114  p->SetProductionVertex(r.Uniform(-v,v), r.Uniform(-v,v), r.Uniform(-v,v), 1);
115  p->SetMomentum(r.Uniform(-m,m), r.Uniform(-m,m), r.Uniform(-m,m)*r.Uniform(1, 3), 1);
116  auto track = new REX::REveTrack(p, 1, prop);
117  track->MakeTrack();
118  if (i % 4 == 3) track->SetLineStyle(2); // enabled dashed style for some tracks
119  track->SetMainColor(kBlue);
120  track->SetName(Form("RandomTrack_%d", i));
121  track->SetTitle(Form("RandomTrack_%d title", i)); // used as tooltip
122  trackHolder->AddElement(track);
123  }
124 
125  event->AddElement(trackHolder);
126 }
127 
128 void addJets()
129 {
130  TRandom &r = *gRandom;
131 
133  auto jetHolder = new REX::REveElement("Jets");
134 
135  int N_Jets = 5 + r.Integer(5);
136  for (int i = 0; i < N_Jets; i++)
137  {
138  auto jet = new REX::REveJetCone(Form("Jet_%d", i));
139  jet->SetTitle(Form("Jet_%d title", i)); // used as tooltip
140  jet->SetCylinder(2*kR_max, 2*kZ_d);
141  jet->AddEllipticCone(r.Uniform(-3.5, 3.5), r.Uniform(0, TMath::TwoPi()),
142  r.Uniform(0.02, 0.2), r.Uniform(0.02, 0.3));
143  jet->SetFillColor(kPink - 8);
144  jet->SetLineColor(kViolet - 7);
145 
146  jetHolder->AddElement(jet);
147  }
148  event->AddElement(jetHolder);
149 }
150 
151 void addVertex()
152 {
153  float pos[3] = {1.46589e-06,-1.30522e-05,-1.98267e-05};
154 
155  // symnetric matrix
156 
157  double a[16] = {1.46589e-01,-1.30522e-02,-1.98267e-02, 0,
158  -1.30522e-02, 4.22955e-02,-5.86628e-03, 0,
159  -1.98267e-02,-5.86628e-03, 2.12836e-01, 0,
160  0, 0, 0, 1};
161 
162  REX::REveTrans t;
163  t.SetFrom(a);
164  TMatrixDSym xxx(3);
165  for(int i = 0; i < 3; i++)
166  for(int j = 0; j < 3; j++)
167  {
168  xxx(i,j) = t(i+1,j+1);
169  }
170 
171  TMatrixDEigen eig(xxx);
172  TVectorD xxxEig ( eig.GetEigenValues() );
173  xxxEig = xxxEig.Sqrt();
174 
175  TMatrixD vecEig = eig.GetEigenVectors();
176  REX::REveVector v[3]; int ei = 0;
177  for (int i = 0; i < 3; ++i)
178  {
179  v[i].Set(vecEig(0,i), vecEig(1,i), vecEig(2,i));
180  v[i] *= xxxEig(i);
181  }
182  REX::REveEllipsoid* ell = new REX::REveEllipsoid("VertexError");
183  ell->InitMainTrans();
184  ell->SetMainColor(kGreen + 10);
185  ell->SetLineWidth(2);
186  ell->SetBaseVectors(v[0], v[1], v[2]);
187  ell->Outline();
189  event->AddElement(ell);
190  return;
191  //center
192  auto ps = new REX::REvePointSet();
193  ps->SetMainColor(kGreen + 10);
194  ps->SetNextPoint(pos[0], pos[1], pos[2]);
195  ps->SetMarkerStyle(4);
196  ps->SetMarkerSize(4);
197  event->AddElement(ps);
198 }
199 
200 
202 {
203  addPoints();
204  addTracks();
205  addJets();
206  addVertex();
207 }
208 
210 {
211  auto b1 = new REX::REveGeoShape("Barrel 1");
212  b1->SetShape(new TGeoTube(kR_min, kR_max, kZ_d));
213  b1->SetMainColor(kCyan);
215 
216  // Debug of surface fill in RPhi (index buffer screwed).
217  // b1->SetNSegments(3);
218  b1->SetNSegments(40);
219 }
220 
221 
223 {
224  // project RhoPhi
225  rPhiGeomScene = eveMng->SpawnNewScene("RPhi Geometry","RPhi");
226  rPhiEventScene = eveMng->SpawnNewScene("RPhi Event Data","RPhi");
227 
229 
230  rphiView = eveMng->SpawnNewViewer("RPhi View", "");
233 
234  // ----------------------------------------------------------------
235 
236  rhoZGeomScene = eveMng->SpawnNewScene("RhoZ Geometry", "RhoZ");
237  rhoZEventScene = eveMng->SpawnNewScene("RhoZ Event Data","RhoZ");
238 
240 
241  rhoZView = eveMng->SpawnNewViewer("RhoZ View", "");
244 }
245 
246 void projectScenes(bool geomp, bool eventp)
247 {
248  if (geomp)
249  {
250  for (auto &ie : eveMng->GetGlobalScene()->RefChildren())
251  {
254  }
255  }
256  if (eventp)
257  {
258  for (auto &ie : eveMng->GetEventScene()->RefChildren())
259  {
262  }
263  }
264 
265  // auto t0 = eveMng->GetEventScene()->FindChild("Tracks")->FirstChild();
266  // printf("t0=%p, %s %s\n", t0, t0->GetElementName(), t0->IsA()->GetName());
267  // dynamic_cast<REX::REveTrack*>(t0)->Print("all");
268 
269  // auto t1 = rPhiEventScene->FindChild("Tracks [P]")->FirstChild();
270  // printf("t1=%p, %s %s\n", t1, t1->GetElementName(), t1->IsA()->GetName());
271  // dynamic_cast<REX::REveTrack*>(t1)->Print("all");
272 }
273 
274 //==============================================================================
275 
276 class EventManager : public REX::REveElement
277 {
278 public:
279  EventManager() = default;
280 
281  virtual ~EventManager() {}
282 
283  virtual void NextEvent()
284  {
286  auto scene = eveMng->GetEventScene();
287  scene->DestroyElements();
288  makeEventScene();
289  for (auto &ie : scene->RefChildren())
290  {
291  if (mngRhoPhi)
293  if (mngRhoZ)
295  }
296  eveMng->EnableRedraw();
297  eveMng->DoRedraw3D();
298  }
299 
300  virtual void QuitRoot()
301  {
302  printf("Quit ROOT\n");
304  }
305 
306 };
307 
309 {
310  // disable browser cache - all scripts and html files will be loaded every time, useful for development
311  // gEnv->SetValue("WebGui.HttpMaxAge", 0);
312 
313  gRandom->SetSeed(0); // make random seed
314 
316 
317  auto eventMng = new EventManager();
318  eventMng->SetName("EventManager");
319  eveMng->GetWorld()->AddElement(eventMng);
320 
321  eveMng->GetWorld()->AddCommand("QuitRoot", "sap-icon://log", eventMng, "QuitRoot()");
322 
323  eveMng->GetWorld()->AddCommand("NextEvent", "sap-icon://step", eventMng, "NextEvent()");
324 
326  makeEventScene();
327 
328  if (1) {
330  projectScenes(true, true);
331  }
332 
333  eveMng->Show();
334 }
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
REX::REveScene * rPhiGeomScene
Definition: event_demo.C:45
Cylindrical tube class.
Definition: TGeoTube.h:17
REveTrackPropagator Calculates path of a particle taking into account special path-marks and imposed ...
auto * m
Definition: textangle.C:8
REveViewer Reve representation of TGLViewer.
Definition: REveViewer.hxx:27
void SetProductionVertex(Double_t vx, Double_t vy, Double_t vz, Double_t t)
Definition: TParticle.h:168
virtual void InitMainTrans(Bool_t can_edit=kTRUE)
Initialize the main transformation to identity matrix.
TMatrixDEigen.
Definition: TMatrixDEigen.h:26
constexpr Double_t TwoPi()
Definition: TMath.h:45
Description of the dynamic properties of a particle.
Definition: TParticle.h:26
REX::REvePointSet * getPointSet(int npoints=2, float s=2, int color=28)
Definition: event_demo.C:55
TVectorT.
Definition: TMatrixTBase.h:79
void projectScenes(bool geomp, bool eventp)
Definition: event_demo.C:246
void addJets()
Definition: event_demo.C:128
void createProjectionStuff()
Definition: event_demo.C:222
static constexpr double ps
int Int_t
Definition: RtypesCore.h:41
void event_demo()
Definition: event_demo.C:308
void SetFrom(Double_t *carr)
Definition: REveTrans.cxx:982
void AddCommand(const std::string &name, const std::string &icon, const REveElement *element, const std::string &action)
Definition: REveScene.hxx:121
Definition: Rtypes.h:64
R__EXTERN TApplication * gApplication
Definition: TApplication.h:166
virtual void Outline()
Draw archade around base vectors.
TMatrixT.
Definition: TMatrixDfwd.h:22
Definition: Rtypes.h:65
This is the base class for the ROOT Random number generators.
Definition: TRandom.h:27
virtual void SetSeed(ULong_t seed=0)
Set the random generator seed.
Definition: TRandom.cxx:597
REveScene * GetEventScene() const
REX::REveScene * rPhiEventScene
Definition: event_demo.C:45
static constexpr double s
virtual void AddElement(REveElement *el)
Add el to the list of children.
TVectorT< Element > & Sqrt()
Take square root of all elements.
Definition: TVectorT.cxx:498
const Double_t kR_max
Definition: event_demo.C:51
REveScene * SpawnNewScene(const char *name, const char *title="")
Create a new scene.
const TMatrixD GetEigenValues() const
Computes the block diagonal eigenvalue matrix.
REX::REveScene * rhoZGeomScene
Definition: event_demo.C:46
REveMagFieldDuo Interface to magnetic field with two different values depending on radius...
void SetBaseVectors(REveVector &v0, REveVector &v1, REveVector &v3)
Three defining base vectors of ellipse.
const Double_t kR_min
Definition: event_demo.C:50
ROOT::R::TRInterface & r
Definition: Object.C:4
virtual void SetMainColor(Color_t color)
Set main color of the element.
void Show(const RWebDisplayArgs &args="")
Show eve manager in specified browser.
auto * a
Definition: textangle.C:12
static REveManager * Create()
If global REveManager* REX::gEve is not set initialize it.
REveProjectionManager Manager class for steering of projections and managing projected objects...
char * Form(const char *fmt,...)
virtual void AddScene(REveScene *scene)
Add &#39;scene&#39; to the list of scenes.
Definition: REveViewer.cxx:59
R__EXTERN TRandom * gRandom
Definition: TRandom.h:62
const TMatrixD & GetEigenVectors() const
Definition: TMatrixDEigen.h:55
void addTracks()
Definition: event_demo.C:90
double Double_t
Definition: RtypesCore.h:55
REveScene * GetGlobalScene() const
void addPoints()
Definition: event_demo.C:70
void SetMomentum(Double_t px, Double_t py, Double_t pz, Double_t e)
Definition: TParticle.h:166
Definition: Rtypes.h:64
virtual void Terminate(Int_t status=0)
REveViewer * SpawnNewViewer(const char *name, const char *title="")
Create a new GL viewer.
virtual REveElement * ImportElements(REveElement *el, REveElement *ext_list=nullptr)
Recursively import elements and apply projection to the newly imported objects.
void makeGeometryScene()
Definition: event_demo.C:209
REX::REveManager * eveMng
Definition: event_demo.C:42
REX::REveViewer * rhoZView
Definition: event_demo.C:48
void SetPdgCode(Int_t pdg)
Change the PDG code for this particle.
Definition: TParticle.cxx:353
REveTrack Track with given vertex, momentum and optional referece-points (path-marks) along its path...
Definition: REveTrack.hxx:39
Definition: Rtypes.h:65
REX::REveProjectionManager * mngRhoPhi
Definition: event_demo.C:43
virtual void DestroyElements()
Destroy all children of this element.
REX::REveViewer * rphiView
Definition: event_demo.C:47
void makeEventScene()
Definition: event_demo.C:201
const Double_t kZ_d
Definition: event_demo.C:52
void DoRedraw3D()
Perform 3D redraw of scenes and viewers whose contents has changed.
void addVertex()
Definition: event_demo.C:151
Definition: Rtypes.h:64
REX::REveScene * rhoZEventScene
Definition: event_demo.C:46
REX::REveProjectionManager * mngRhoZ
Definition: event_demo.C:44
REveScene * GetWorld() const