Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
RGeomData.cxx
Go to the documentation of this file.
1// Author: Sergey Linev, 14.12.2018
2
3/*************************************************************************
4 * Copyright (C) 1995-2023, Rene Brun and Fons Rademakers. *
5 * All rights reserved. *
6 * *
7 * For the licensing terms see $ROOTSYS/LICENSE. *
8 * For the list of contributors see $ROOTSYS/README/CREDITS. *
9 *************************************************************************/
10
11#include <ROOT/RGeomData.hxx>
12
15#include <ROOT/RLogger.hxx>
16#include "CsgOps.h"
17
18#include "TMath.h"
19#include "TColor.h"
20#include "TROOT.h"
21#include "TGeoNode.h"
22#include "TGeoVolume.h"
23#include "TGeoBBox.h"
24#include "TGeoSphere.h"
25#include "TGeoCone.h"
26#include "TGeoTube.h"
27#include "TGeoEltu.h"
28#include "TGeoTorus.h"
29#include "TGeoPcon.h"
30#include "TGeoPgon.h"
31#include "TGeoXtru.h"
32#include "TGeoParaboloid.h"
33#include "TGeoHype.h"
34#include "TGeoTessellated.h"
35#include "TGeoScaledShape.h"
36#include "TGeoCompositeShape.h"
37#include "TGeoManager.h"
38#include "TGeoMatrix.h"
39#include "TGeoMedium.h"
40#include "TGeoMaterial.h"
41#include "TGeoBoolNode.h"
42#include "TBuffer3D.h"
43#include "TBufferJSON.h"
44#include "TRegexp.h"
45
46#include <algorithm>
47#include <array>
48
50{
51 static ROOT::RLogChannel sLog("ROOT.Geom");
52 return sLog;
53}
54
55
56namespace ROOT {
57
58/** Iterator of hierarchical geometry structures */
59
61
63 int fParentId{-1};
64 unsigned fChild{0};
65 int fNodeId{0};
66
67 std::vector<int> fStackParents;
68 std::vector<int> fStackChilds;
69
70public:
72
73 const std::string &GetName() const { return fDesc.fDesc[fNodeId].name; }
74
75 const std::string &GetColor() const { return fDesc.fDesc[fNodeId].color; }
76
77 const std::string &GetMaterial() const { return fDesc.fDesc[fNodeId].material; }
78
79 int GetVisible() const { return fDesc.fDesc[fNodeId].vis; }
80
81 bool IsValid() const { return fNodeId >= 0; }
82
83 int GetNodeId() const { return fNodeId; }
84
85 bool HasChilds() const { return (fNodeId < 0) ? true : !fDesc.fDesc[fNodeId].chlds.empty(); }
86
87 int NumChilds() const { return (fNodeId < 0) ? 1 : fDesc.fDesc[fNodeId].chlds.size(); }
88
89 bool Enter()
90 {
91 if (fNodeId < 0) {
92 Reset();
93 fNodeId = 0;
94 return true;
95 }
96
97 if (fNodeId >= (int)fDesc.fDesc.size())
98 return false;
99
100 auto &node = fDesc.fDesc[fNodeId];
101 if (node.chlds.empty())
102 return false;
103 fStackParents.emplace_back(fParentId);
104 fStackChilds.emplace_back(fChild);
106 fChild = 0;
107 fNodeId = node.chlds[fChild];
108 return true;
109 }
110
111 bool Leave()
112 {
113 if (fStackParents.empty()) {
114 fNodeId = -1;
115 return false;
116 }
117 fParentId = fStackParents.back();
118 fChild = fStackChilds.back();
119
120 fStackParents.pop_back();
121 fStackChilds.pop_back();
122
123 if (fParentId < 0) {
124 fNodeId = 0;
125 } else {
127 }
128 return true;
129 }
130
131 bool Next()
132 {
133 // does not have parents
134 if ((fNodeId <= 0) || (fParentId < 0)) {
135 Reset();
136 return false;
137 }
138
139 auto &prnt = fDesc.fDesc[fParentId];
140 if (++fChild >= prnt.chlds.size()) {
141 fNodeId = -1; // not valid node, only Leave can be called
142 return false;
143 }
144
145 fNodeId = prnt.chlds[fChild];
146 return true;
147 }
148
149 bool Reset()
150 {
151 fParentId = -1;
152 fNodeId = -1;
153 fChild = 0;
154 fStackParents.clear();
155 fStackChilds.clear();
156
157 return true;
158 }
159
160 bool NextNode()
161 {
162 if (Enter())
163 return true;
164
165 if (Next())
166 return true;
167
168 while (Leave()) {
169 if (Next())
170 return true;
171 }
172
173 return false;
174 }
175
176 /** Navigate to specified path - path specified as string and should start with "/" */
177 bool Navigate(const std::string &path)
178 {
179 size_t pos = path.find('/');
180 if (pos != 0)
181 return false;
182
183 Reset(); // set to the top of element
184
185 while (++pos < path.length()) {
186 auto last = pos;
187
188 pos = path.find('/', last);
189
190 if (pos == std::string::npos)
191 pos = path.length();
192
193 std::string folder = path.substr(last, pos - last);
194
195 if (!Enter())
196 return false;
197
198 bool find = false;
199
200 do {
201 find = (folder.compare(GetName()) == 0);
202 } while (!find && Next());
203
204 if (!find)
205 return false;
206 }
207
208 return true;
209 }
210
211 /** Navigate to specified path */
212 bool Navigate(const std::vector<std::string> &path)
213 {
214 Reset(); // set to the top of element
215
216 for (auto &folder : path) {
217
218 if (!Enter())
219 return false;
220
221 bool find = false;
222
223 do {
224 find = (folder.compare(GetName()) == 0);
225 } while (!find && Next());
226
227 if (!find)
228 return false;
229 }
230
231 return true;
232 }
233
234 /** Navigate to specified volume - find first occurrence */
236 {
237 Reset();
238
239 while (NextNode()) {
240 if (vol == fDesc.GetVolume(GetNodeId()))
241 return true;
242 }
243
244 return false;
245 }
246
247 /// Returns array of ids to currently selected node
248 std::vector<int> CurrentIds() const
249 {
250 std::vector<int> res;
251 if (IsValid()) {
252 for (unsigned n = 1; n < fStackParents.size(); ++n)
253 res.emplace_back(fStackParents[n]);
254 if (fParentId >= 0)
255 res.emplace_back(fParentId);
256 res.emplace_back(fNodeId);
257 }
258 return res;
259 }
260};
261
262} // namespace ROOT
263
264using namespace ROOT;
265
266using namespace std::string_literals;
267
268namespace {
269
270int compare_stacks(const std::vector<int> &stack1, const std::vector<int> &stack2)
271{
272 unsigned len1 = stack1.size(), len2 = stack2.size(), len = (len1 < len2) ? len1 : len2, indx = 0;
273 while (indx < len) {
274 if (stack1[indx] < stack2[indx])
275 return -1;
276 if (stack1[indx] > stack2[indx])
277 return 1;
278 ++indx;
279 }
280
281 if (len1 < len2)
282 return -1;
283 if (len1 > len2)
284 return 1;
285
286 return 0;
287}
288} // namespace
289
290/////////////////////////////////////////////////////////////////////
291/// Issue signal, which distributed on all handlers - excluding source handler
292
293void RGeomDescription::IssueSignal(const void *handler, const std::string &kind)
294{
295 std::vector<RGeomSignalFunc_t> funcs;
296
297 {
298 TLockGuard lock(fMutex);
299 for (auto &pair : fSignals)
300 if (!handler || (pair.first != handler))
301 funcs.emplace_back(pair.second);
302 }
303
304 // invoke signal outside locked mutex to avoid any locking
305 for (auto func : funcs)
306 func(kind);
307}
308
309/////////////////////////////////////////////////////////////////////
310/// Add signal handler
311
313{
314 TLockGuard lock(fMutex);
315 fSignals.emplace_back(handler, func);
316}
317
318/////////////////////////////////////////////////////////////////////
319/// Remove signal handler
320
322{
323 TLockGuard lock(fMutex);
324
325 for (auto iter = fSignals.begin(); iter != fSignals.end(); ++iter)
326 if (handler == iter->first) {
327 fSignals.erase(iter);
328 return;
329 }
330}
331
332/////////////////////////////////////////////////////////////////////
333/// Pack matrix into vector, which can be send to client
334/// Following sizes can be used for vector:
335/// 0 - Identity matrix
336/// 3 - Translation
337/// 4 - Scale (last element always 1)
338/// 9 - Rotation
339/// 16 - Full size
340
341void RGeomDescription::PackMatrix(std::vector<float> &vect, TGeoMatrix *matr)
342{
343 vect.clear();
344
345 if (!matr || matr->IsIdentity())
346 return;
347
348 auto trans = matr->GetTranslation();
349 auto scale = matr->GetScale();
350 auto rotate = matr->GetRotationMatrix();
351
352 bool is_translate = matr->IsA() == TGeoTranslation::Class(), is_scale = matr->IsA() == TGeoScale::Class(),
353 is_rotate = matr->IsA() == TGeoRotation::Class();
354
355 if (!is_translate && !is_scale && !is_rotate) {
356 // check if trivial matrix
357
358 auto test = [](double val, double chk) { return (val == chk) || (TMath::Abs(val - chk) < 1e-20); };
359
360 bool no_scale = test(scale[0], 1) && test(scale[1], 1) && test(scale[2], 1);
361 bool no_trans = test(trans[0], 0) && test(trans[1], 0) && test(trans[2], 0);
362 bool no_rotate = test(rotate[0], 1) && test(rotate[1], 0) && test(rotate[2], 0) && test(rotate[3], 0) &&
363 test(rotate[4], 1) && test(rotate[5], 0) && test(rotate[6], 0) && test(rotate[7], 0) &&
364 test(rotate[8], 1);
365
366 if (no_scale && no_trans && no_rotate)
367 return;
368
369 if (no_scale && no_trans && !no_rotate) {
370 is_rotate = true;
371 } else if (no_scale && !no_trans && no_rotate) {
372 is_translate = true;
373 } else if (!no_scale && no_trans && no_rotate) {
374 is_scale = true;
375 }
376 }
377
378 if (is_translate) {
379 vect.resize(3);
380 vect[0] = trans[0];
381 vect[1] = trans[1];
382 vect[2] = trans[2];
383 return;
384 }
385
386 if (is_scale) {
387 vect.resize(4);
388 vect[0] = scale[0];
389 vect[1] = scale[1];
390 vect[2] = scale[2];
391 vect[3] = 1;
392 return;
393 }
394
395 if (is_rotate) {
396 vect.resize(9);
397 for (int n = 0; n < 9; ++n)
398 vect[n] = rotate[n];
399 return;
400 }
401
402 vect.resize(16);
403 vect[0] = rotate[0];
404 vect[4] = rotate[1];
405 vect[8] = rotate[2];
406 vect[12] = trans[0];
407 vect[1] = rotate[3];
408 vect[5] = rotate[4];
409 vect[9] = rotate[5];
410 vect[13] = trans[1];
411 vect[2] = rotate[6];
412 vect[6] = rotate[7];
413 vect[10] = rotate[8];
414 vect[14] = trans[2];
415 vect[3] = 0;
416 vect[7] = 0;
417 vect[11] = 0;
418 vect[15] = 1;
419}
420
421/////////////////////////////////////////////////////////////////////
422/// Collect information about geometry hierarchy into flat list
423/// like it done in JSROOT ClonedNodes.createClones
424
426{
428 if (!mgr)
429 return;
430
431 TLockGuard lock(fMutex);
432
433 // by top node visibility always enabled and harm logic
434 // later visibility can be controlled by other means
435 // mgr->GetTopNode()->GetVolume()->SetVisibility(kFALSE);
436
437 int maxnodes = mgr->GetMaxVisNodes();
438
439 SetNSegments(mgr->GetNsegments());
440 SetVisLevel(mgr->GetVisLevel());
442 SetMaxVisFaces((maxnodes > 5000 ? 5000 : (maxnodes < 1000 ? 1000 : maxnodes)) * 100);
443
444 auto topnode = mgr->GetTopNode();
445
446 BuildDescription(topnode, topnode->GetVolume());
447
448 if (!volname.empty()) {
449 auto vol = mgr->GetVolume(volname.c_str());
450 RGeomBrowserIter iter(*this);
451 if (vol && (vol != topnode->GetVolume()) && iter.Navigate(vol))
453 }
454}
455
456/////////////////////////////////////////////////////////////////////
457/// Collect information about geometry from single volume
458/// like it done in JSROOT ClonedNodes.createClones
459
461{
463 if (!vol)
464 return;
465
466 TLockGuard lock(fMutex);
467
468 fDrawVolume = vol;
469
470 fSelectedStack.clear();
471
473}
474
475/////////////////////////////////////////////////////////////////////
476/// Clear geometry description
477
479{
480 TLockGuard lock(fMutex);
481
482 fDesc.clear();
483 fNodes.clear();
484 fSortMap.clear();
486 fDrawIdCut = 0;
487 fDrawVolume = nullptr;
488 fSelectedStack.clear();
489}
490
491/////////////////////////////////////////////////////////////////////
492/// Build geometry description
493
495{
496 // vector to remember numbers
497 std::vector<int> numbers;
498 int offset = 1000000000;
499
500 // try to build flat list of all nodes
503 do {
504 if (!snode) {
505 numbers.emplace_back(offset);
506 fNodes.emplace_back(nullptr);
507 } else if (snode->GetNumber() >= offset) {
508 // artificial offset already applied, used as identifier
509 iter.Skip(); // no need to look inside
510 } else {
511 numbers.emplace_back(snode->GetNumber());
512 snode->SetNumber(offset + fNodes.size()); // use id with shift 1e9
513 fNodes.emplace_back(snode);
514 }
515 } while ((snode = iter()) != nullptr);
516
517 fDesc.reserve(fNodes.size());
518 fSortMap.reserve(fNodes.size());
519
520 // array for sorting
521 std::vector<RGeomNode *> sortarr;
522 sortarr.reserve(fNodes.size());
523
524 // create vector of desc and childs
525 int cnt = 0;
526 for (auto node : fNodes) {
527
528 fDesc.emplace_back(node ? node->GetNumber() - offset : 0);
529 TGeoVolume *vol = node ? node->GetVolume() : topvolume;
530
531 auto &desc = fDesc[cnt++];
532
533 sortarr.emplace_back(&desc);
534
535 desc.name = node ? node->GetName() : vol->GetName();
536
537 auto shape = dynamic_cast<TGeoBBox *>(vol->GetShape());
538 if (shape) {
539 desc.vol = TMath::Sqrt(shape->GetDX() * shape->GetDX() + shape->GetDY() * shape->GetDY() +
540 shape->GetDZ() * shape->GetDZ());
541 desc.nfaces = CountShapeFaces(shape);
542 }
543
544 CopyMaterialProperties(vol, desc);
545
546 auto chlds = node ? node->GetNodes() : vol->GetNodes();
547
548 PackMatrix(desc.matr, node ? node->GetMatrix() : nullptr);
549
550 if (chlds)
551 for (int n = 0; n <= chlds->GetLast(); ++n) {
552 auto chld = dynamic_cast<TGeoNode *>(chlds->At(n));
553 desc.chlds.emplace_back(chld->GetNumber() - offset);
554 }
555 }
556
557 // recover numbers
558 cnt = 0;
559 for (auto node : fNodes) {
560 auto number = numbers[cnt++];
561 if (node)
562 node->SetNumber(number);
563 }
564
565 // sort in volume descent order
566 std::sort(sortarr.begin(), sortarr.end(), [](RGeomNode *a, RGeomNode *b) { return a->vol > b->vol; });
567
568 cnt = 0;
569 for (auto &elem : sortarr) {
570 fSortMap.emplace_back(elem->id);
571 elem->sortid = cnt++; // keep place in sorted array to correctly apply cut
572 }
573
574 MarkVisible(); // set visibility flags
575
577}
578
579/////////////////////////////////////////////////////////////////////
580/// Get volume for specified nodeid
581/// If specific volume was configured, it will be returned for nodeid==0
582
584{
585 auto node = fNodes[nodeid];
586 if (node)
587 return node->GetVolume();
588 return nodeid == 0 ? fDrawVolume : nullptr;
589}
590
591/////////////////////////////////////////////////////////////////////
592/// Set visibility flag for each nodes
593
595{
596 int res = 0;
597 for (int nodeid = 0; nodeid < (int)fNodes.size(); nodeid++) {
598
599 auto node = fNodes[nodeid];
600 auto vol = GetVolume(nodeid);
601 auto &desc = fDesc[nodeid];
602 desc.vis = 0;
603 desc.nochlds = false;
604
605 if (on_screen) {
606 if (!node || node->IsOnScreen())
607 desc.vis = 99;
608 } else {
609 if (vol->IsVisible() && !vol->TestAttBit(TGeoAtt::kVisNone))
610 desc.vis = 99;
611
612 if (node && !node->IsVisDaughters())
613 desc.nochlds = true;
614
615 if ((desc.vis > 0) && (!desc.chlds.empty()) && !desc.nochlds)
616 desc.vis = 1;
617 }
618
619 if (desc.IsVisible() && desc.CanDisplay())
620 res++;
621 }
622
623 return res;
624}
625
626/////////////////////////////////////////////////////////////////////
627/// Count total number of visible childs under each node
628
630{
631 for (auto &node : fDesc)
632 node.idshift = -1;
633
634 using ScanFunc_t = std::function<int(RGeomNode &)>;
635
636 ScanFunc_t scan_func = [&, this](RGeomNode &node) {
637 if (node.idshift < 0) {
638 node.idshift = 0;
639 for (auto id : node.chlds)
640 node.idshift += scan_func(fDesc[id]);
641 }
642
643 return node.idshift + 1;
644 };
645
646 if (!fDesc.empty())
647 scan_func(fDesc[0]);
648}
649
650/////////////////////////////////////////////////////////////////////
651/// Iterate over all nodes and call function for visible
652
654{
655 if (fDesc.empty())
656 return 0;
657
658 std::vector<int> stack;
659 stack.reserve(25); // reserve enough space for most use-cases
660 int counter = 0;
661 auto viter = fVisibility.begin();
662
663 using ScanFunc_t = std::function<int(int, int, bool)>;
664
665 ScanFunc_t scan_func = [&, this](int nodeid, int lvl, bool is_inside) {
666 if (!is_inside && (fSelectedStack == stack))
667 is_inside = true;
668
669 auto &desc = fDesc[nodeid];
670 auto desc_vis = desc.vis;
671 int res = 0;
672
673 if (desc.nochlds && (lvl > 0))
674 lvl = 0;
675
676 bool can_display = desc.CanDisplay(), scan_childs = true;
677
678 if ((viter != fVisibility.end()) && (compare_stacks(viter->stack, stack) == 0)) {
679 can_display = scan_childs = viter->visible;
680 desc_vis = !viter->visible ? 0 : (!desc.chlds.empty() ? 1 : 99);
681 viter++;
682 }
683
684 // same logic as in JSROOT ClonedNodes.scanVisible
685 bool is_visible = (lvl >= 0) && (desc_vis > lvl) && can_display && is_inside;
686
687 if (is_visible || !only_visible)
688 if (func(desc, stack, is_visible, counter))
689 res++;
690
691 counter++; // count sequence id of current position in scan, will be used later for merging drawing lists
692
693 if ((!desc.chlds.empty()) && (((lvl > 0) && scan_childs) || !only_visible)) {
694 auto pos = stack.size();
695 stack.emplace_back(0);
696 for (unsigned k = 0; k < desc.chlds.size(); ++k) {
697 stack[pos] = k; // stack provides index in list of childs
698 res += scan_func(desc.chlds[k], is_inside ? lvl - 1 : lvl, is_inside);
699 }
700 stack.pop_back();
701 } else {
702 counter += desc.idshift;
703 }
704
705 return res;
706 };
707
708 if (!maxlvl && (GetVisLevel() > 0))
710 if (!maxlvl)
711 maxlvl = 4;
712 if (maxlvl > 97)
713 maxlvl = 97; // check while vis property of node is 99 normally
714
715 return scan_func(0, maxlvl, false);
716}
717
718/////////////////////////////////////////////////////////////////////
719/// Collect nodes which are used in visibles
720
722{
723 drawing.cfg = &fCfg;
724
725 drawing.numnodes = fDesc.size();
726
727 if (all_nodes) {
728 for (auto &node : fDesc)
729 drawing.nodes.emplace_back(&node);
730 return;
731 }
732
733 // TODO: for now reset all flags, later can be kept longer
734 for (auto &node : fDesc)
735 node.useflag = false;
736
737 for (auto &item : drawing.visibles) {
738 int nodeid = 0;
739 for (auto &chindx : item.stack) {
740 auto &node = fDesc[nodeid];
741 if (!node.useflag) {
742 node.useflag = true;
743 drawing.nodes.emplace_back(&node);
744 }
745 if (chindx >= (int)node.chlds.size())
746 break;
747 nodeid = node.chlds[chindx];
748 }
749
750 if (nodeid != item.nodeid)
751 printf("Nodeid mismatch %d != %d when extracting nodes for visibles\n", nodeid, item.nodeid);
752
753 auto &node = fDesc[nodeid];
754 if (!node.useflag) {
755 node.useflag = true;
756 drawing.nodes.emplace_back(&node);
757 }
758 }
759
760 // printf("SELECT NODES %d\n", (int) drawing.nodes.size());
761}
762
763/////////////////////////////////////////////////////////////////////
764/// Find description object for requested shape
765/// If not exists - will be created
766
767std::string RGeomDescription::ProcessBrowserRequest(const std::string &msg)
768{
769 TLockGuard lock(fMutex);
770
771 std::string res;
772
773 auto request = TBufferJSON::FromJSON<RBrowserRequest>(msg);
774
775 if (msg.empty()) {
776 request = std::make_unique<RBrowserRequest>();
777 request->first = 0;
778 request->number = 100;
779 }
780
781 if (!request)
782 return res;
783
784 if (request->path.empty() && (request->first == 0) && (GetNumNodes() < (IsPreferredOffline() ? 1000000 : 1000))) {
785
786 std::vector<RGeomNodeBase *> vect(fDesc.size(), nullptr);
787
788 int cnt = 0;
789 for (auto &item : fDesc)
790 vect[cnt++] = &item;
791
792 res = "DESCR:"s + TBufferJSON::ToJSON(&vect, GetJsonComp()).Data();
793
794 if (!fVisibility.empty()) {
795 res += ":__PHYSICAL_VISIBILITY__:";
797 }
798
799 res += ":__SELECTED_STACK__:";
801
802 } else {
803 std::vector<RGeoItem> temp_nodes;
804 bool toplevel = request->path.empty();
805
806 // create temporary object for the short time
808 reply.path = request->path;
809 reply.first = request->first;
810
811 RGeomBrowserIter iter(*this);
812 if (iter.Navigate(request->path)) {
813
814 reply.nchilds = iter.NumChilds();
815 // scan childs of selected nodes
816 if (iter.Enter()) {
817
818 while ((request->first > 0) && iter.Next()) {
819 request->first--;
820 }
821
822 // first element
823 auto stack = MakeStackByIds(iter.CurrentIds());
824
825 while (iter.IsValid() && (request->number > 0)) {
826 int pvis = IsPhysNodeVisible(stack);
827 temp_nodes.emplace_back(iter.GetName(), iter.NumChilds(), iter.GetNodeId(), iter.GetColor(),
828 iter.GetMaterial(), iter.GetVisible(), pvis < 0 ? iter.GetVisible() : pvis);
829 if (toplevel)
830 temp_nodes.back().SetExpanded(true);
831 if (stack == fSelectedStack)
832 temp_nodes.back().SetTop(true);
833 request->number--;
834
835 if (!stack.empty())
836 stack[stack.size() - 1]++;
837
838 if (!iter.Next())
839 break;
840 }
841 }
842 }
843
844 for (auto &n : temp_nodes)
845 reply.nodes.emplace_back(&n);
846
847 res = "BREPL:"s + TBufferJSON::ToJSON(&reply, GetJsonComp()).Data();
848 }
849
850 return res;
851}
852
853/////////////////////////////////////////////////////////////////////
854/// Find description object for requested shape
855/// If not exists - will be created
856
858{
859 for (auto &descr : fShapes)
860 if (descr.fShape == shape)
861 return descr;
862
863 fShapes.emplace_back(shape);
864 auto &elem = fShapes.back();
865 elem.id = fShapes.size() - 1;
866 return elem;
867}
868
869////////////////////////////////////////////////////////////////////////
870/// Function produces mesh for provided shape, applying matrix to the result
871
872std::unique_ptr<RootCsg::TBaseMesh> MakeGeoMesh(TGeoMatrix *matr, TGeoShape *shape)
873{
874 TGeoCompositeShape *comp = dynamic_cast<TGeoCompositeShape *>(shape);
875
876 std::unique_ptr<RootCsg::TBaseMesh> res;
877
878 if (!comp) {
879 std::unique_ptr<TBuffer3D> b3d(shape->MakeBuffer3D());
880
881 if (matr) {
882 Double_t *v = b3d->fPnts;
883 Double_t buf[3];
884 for (UInt_t i = 0; i < b3d->NbPnts(); ++i) {
885 buf[0] = v[i * 3];
886 buf[1] = v[i * 3 + 1];
887 buf[2] = v[i * 3 + 2];
888 matr->LocalToMaster(buf, &v[i * 3]);
889 }
890 }
891
892 res.reset(RootCsg::ConvertToMesh(*b3d.get()));
893 } else {
894 auto node = comp->GetBoolNode();
895
896 TGeoHMatrix mleft, mright;
897 if (matr) {
898 mleft = *matr;
899 mright = *matr;
900 }
901
902 mleft.Multiply(node->GetLeftMatrix());
903 auto left = MakeGeoMesh(&mleft, node->GetLeftShape());
904
905 mright.Multiply(node->GetRightMatrix());
906 auto right = MakeGeoMesh(&mright, node->GetRightShape());
907
908 if (node->IsA() == TGeoUnion::Class())
909 res.reset(RootCsg::BuildUnion(left.get(), right.get()));
910 if (node->IsA() == TGeoIntersection::Class())
911 res.reset(RootCsg::BuildIntersection(left.get(), right.get()));
912 if (node->IsA() == TGeoSubtraction::Class())
913 res.reset(RootCsg::BuildDifference(left.get(), right.get()));
914 }
915
916 return res;
917}
918
919/////////////////////////////////////////////////////////////////////
920/// Returns really used number of cylindrical segments
921
923{
924 int nsegm = 0;
925
926 if (GetNSegments() > 0)
927 nsegm = GetNSegments();
928 else if (gGeoManager && (gGeoManager->GetNsegments() > 0))
929 nsegm = gGeoManager->GetNsegments();
930
931 return nsegm > min ? nsegm : min;
932}
933
934/////////////////////////////////////////////////////////////////////
935/// Count number of faces for the shape
936
938{
939 if (!shape)
940 return 0;
941
942 auto countTubeFaces = [this](const std::array<Double_t, 2> &outerR, const std::array<Double_t, 2> &innerR,
943 Double_t thetaLength = 360.) -> int {
944 auto hasrmin = (innerR[0] > 0) || (innerR[1] > 0);
945
947
948 // external surface
949 int numfaces = radiusSegments * (((outerR[0] <= 0) || (outerR[1] <= 0)) ? 1 : 2);
950
951 // internal surface
952 if (hasrmin)
953 numfaces += radiusSegments * (((innerR[0] <= 0) || (innerR[1] <= 0)) ? 1 : 2);
954
955 // upper cap
956 if (outerR[0] > 0)
957 numfaces += radiusSegments * ((innerR[0] > 0) ? 2 : 1);
958 // bottom cup
959 if (outerR[1] > 0)
960 numfaces += radiusSegments * ((innerR[1] > 0) ? 2 : 1);
961
962 if (thetaLength < 360)
963 numfaces += ((outerR[0] > innerR[0]) ? 2 : 0) + ((outerR[1] > innerR[1]) ? 2 : 0);
964
965 return numfaces;
966 };
967
968 if (shape->IsA() == TGeoSphere::Class()) {
969 TGeoSphere *sphere = (TGeoSphere *)shape;
970 auto widthSegments = sphere->GetNumberOfDivisions();
971 auto heightSegments = sphere->GetNz();
972 auto phiLength = sphere->GetPhi2() - sphere->GetPhi1();
973 auto noInside = sphere->GetRmin() <= 0;
974
976 auto numtop = widthSegments * (noInside ? 1 : 2);
977 auto numbottom = widthSegments * (noInside ? 1 : 2);
978 auto numcut = (phiLength == 360.) ? 0 : heightSegments * (noInside ? 2 : 4);
979
980 return numoutside * (noInside ? 1 : 2) + numtop + numbottom + numcut;
981 } else if (shape->IsA() == TGeoCone::Class()) {
982 auto cone = (TGeoCone *)shape;
983 return countTubeFaces({cone->GetRmax2(), cone->GetRmax1()}, {cone->GetRmin2(), cone->GetRmin1()});
984 } else if (shape->IsA() == TGeoConeSeg::Class()) {
985 auto cone = (TGeoConeSeg *)shape;
986 return countTubeFaces({cone->GetRmax2(), cone->GetRmax1()}, {cone->GetRmin2(), cone->GetRmin1()},
987 cone->GetPhi2() - cone->GetPhi1());
988 } else if (shape->IsA() == TGeoTube::Class()) {
989 auto tube = (TGeoTube *)shape;
990 return countTubeFaces({tube->GetRmax(), tube->GetRmax()}, {tube->GetRmin(), tube->GetRmin()});
991 } else if (shape->IsA() == TGeoTubeSeg::Class()) {
992 auto tube = (TGeoTubeSeg *)shape;
993 return countTubeFaces({tube->GetRmax(), tube->GetRmax()}, {tube->GetRmin(), tube->GetRmin()},
994 tube->GetPhi2() - tube->GetPhi1());
995 } else if (shape->IsA() == TGeoCtub::Class()) {
996 auto tube = (TGeoCtub *)shape;
997 return countTubeFaces({tube->GetRmax(), tube->GetRmax()}, {tube->GetRmin(), tube->GetRmin()},
998 tube->GetPhi2() - tube->GetPhi1());
999 } else if (shape->IsA() == TGeoEltu::Class()) {
1000 return GetUsedNSegments(4) * 4;
1001 } else if (shape->IsA() == TGeoTorus::Class()) {
1002 auto torus = (TGeoTorus *)shape;
1004 auto tubularSegments = TMath::Max(8, TMath::Nint(torus->GetDphi() / 360. * GetUsedNSegments()));
1005 return (torus->GetRmin() > 0 ? 4 : 2) * radialSegments * (tubularSegments + (torus->GetDphi() != 360. ? 1 : 0));
1006 } else if (shape->IsA() == TGeoPcon::Class()) {
1007 auto pcon = (TGeoPcon *)shape;
1008
1009 bool hasrmin = false;
1010 int radiusSegments = TMath::Max(5, TMath::Nint(pcon->GetDphi() / 360 * GetUsedNSegments()));
1011 for (int layer = 0; layer < pcon->GetNz(); ++layer)
1012 if (pcon->GetRmin(layer) > 0.)
1013 hasrmin = true;
1014 return (hasrmin ? 4 : 2) * radiusSegments * (pcon->GetNz() - 1);
1015 } else if (shape->IsA() == TGeoPgon::Class()) {
1016 auto pgon = (TGeoPgon *)shape;
1017
1018 bool hasrmin = false;
1019 int radiusSegments = TMath::Max(5, TMath::Nint(pgon->GetDphi() / 360 * GetUsedNSegments()));
1020 for (int layer = 0; layer < pgon->GetNz(); ++layer)
1021 if (pgon->GetRmin(layer) > 0.)
1022 hasrmin = true;
1023 return (hasrmin ? 4 : 2) * radiusSegments * (pgon->GetNz() - 1);
1024 } else if (shape->IsA() == TGeoXtru::Class()) {
1025 auto xtru = (TGeoXtru *)shape;
1026 return (xtru->GetNz() - 1) * xtru->GetNvert() * 2 + xtru->GetNvert() * 3;
1027 } else if (shape->IsA() == TGeoParaboloid::Class()) {
1028 auto para = (TGeoParaboloid *)shape;
1030 int numfaces = (heightSegments + 1) * radiusSegments * 2;
1031 if (para->GetRlo() == 0.)
1032 numfaces -= radiusSegments * 2; // complete layer
1033 if (para->GetRhi() == 0.)
1034 numfaces -= radiusSegments * 2; // complete layer
1035 return numfaces;
1036 } else if (shape->IsA() == TGeoHype::Class()) {
1037 TGeoHype *hype = (TGeoHype *)shape;
1038 if ((hype->GetStIn() == 0) && (hype->GetStOut() == 0))
1039 return countTubeFaces({hype->GetRmax(), hype->GetRmax()}, {hype->GetRmin(), hype->GetRmin()});
1041 return radiusSegments * (heightSegments + 1) * ((hype->GetRmin() > 0.) ? 4 : 2);
1042 } else if (shape->IsA() == TGeoTessellated::Class()) {
1043 auto tess = (TGeoTessellated *)shape;
1044 int numfaces = 0;
1045 for (int i = 0; i < tess->GetNfacets(); ++i) {
1046 if (tess->GetFacet(i).GetNvert() == 4)
1047 numfaces += 2;
1048 else
1049 numfaces += 1;
1050 }
1051 return numfaces;
1052 } else if (shape->IsA() == TGeoScaledShape::Class()) {
1053 auto scaled = (TGeoScaledShape *)shape;
1054 return CountShapeFaces(scaled->GetShape());
1055 } else if (shape->IsA() == TGeoCompositeShape::Class()) {
1056 auto comp = (TGeoCompositeShape *)shape;
1057 if (!comp->GetBoolNode())
1058 return 0;
1059 return CountShapeFaces(comp->GetBoolNode()->GetLeftShape()) +
1060 CountShapeFaces(comp->GetBoolNode()->GetRightShape());
1061 }
1062
1063 // many of simple shapes have 12 faces
1064 return 12;
1065}
1066
1067/////////////////////////////////////////////////////////////////////
1068/// Find description object and create render information
1069
1071{
1072 auto &elem = FindShapeDescr(shape);
1073
1074 if (elem.nfaces == 0) {
1075
1076 int boundary = 3; //
1077 if (shape->IsComposite()) {
1078 // composite is most complex for client, therefore by default build on server
1079 boundary = 1;
1080 } else if (!shape->IsCylType()) {
1081 // simple box geometry is compact and can be delivered as raw
1082 boundary = 2;
1083 }
1084
1085 if (IsBuildShapes() < boundary) {
1086 elem.nfaces = 1;
1087 elem.fShapeInfo.shape = shape;
1088 } else {
1089
1090 int old_nsegm = -1;
1091 if (fCfg.nsegm > 0 && gGeoManager) {
1094 }
1095
1096 auto mesh = MakeGeoMesh(nullptr, shape);
1097
1098 if (old_nsegm > 0 && gGeoManager)
1100
1101 Int_t num_vertices = mesh->NumberOfVertices(), num_polynoms = 0;
1102
1103 for (unsigned polyIndex = 0; polyIndex < mesh->NumberOfPolys(); ++polyIndex) {
1104
1105 auto size_of_polygon = mesh->SizeOfPoly(polyIndex);
1106
1107 if (size_of_polygon >= 3)
1109 }
1110
1111 Int_t index_buffer_size = num_polynoms * 3, // triangle indexes
1112 vertex_buffer_size = num_vertices * 3; // X,Y,Z array
1113
1114 elem.nfaces = num_polynoms;
1115
1116 std::vector<float> vertices(vertex_buffer_size);
1117
1118 for (Int_t i = 0; i < num_vertices; ++i) {
1119 auto v = mesh->GetVertex(i);
1120 vertices[i * 3] = v[0];
1121 vertices[i * 3 + 1] = v[1];
1122 vertices[i * 3 + 2] = v[2];
1123 }
1124
1125 elem.fRawInfo.raw.resize(vertices.size() * sizeof(float));
1126
1127 memcpy(reinterpret_cast<char *>(elem.fRawInfo.raw.data()), vertices.data(), vertices.size() * sizeof(float));
1128
1129 auto &indexes = elem.fRawInfo.idx;
1130
1131 indexes.resize(index_buffer_size);
1132 int pos = 0;
1133
1134 for (unsigned polyIndex = 0; polyIndex < mesh->NumberOfPolys(); ++polyIndex) {
1135 auto size_of_polygon = mesh->SizeOfPoly(polyIndex);
1136
1137 // add first triangle
1138 if (size_of_polygon >= 3)
1139 for (int i = 0; i < 3; ++i)
1140 indexes[pos++] = mesh->GetVertexIndex(polyIndex, i);
1141
1142 // add following triangles
1143 if (size_of_polygon > 3)
1144 for (unsigned vertex = 3; vertex < size_of_polygon; vertex++) {
1145 indexes[pos++] = mesh->GetVertexIndex(polyIndex, 0);
1146 indexes[pos++] = mesh->GetVertexIndex(polyIndex, vertex - 1);
1147 indexes[pos++] = mesh->GetVertexIndex(polyIndex, vertex);
1148 }
1149 }
1150 }
1151 }
1152
1153 return elem;
1154}
1155
1156/////////////////////////////////////////////////////////////////////
1157/// Copy material properties
1158
1160{
1161 if (!volume)
1162 return;
1163
1164 TColor *col = nullptr;
1165
1166 if ((volume->GetFillColor() > 1) && (volume->GetLineColor() == 1))
1167 col = gROOT->GetColor(volume->GetFillColor());
1168 else if (volume->GetLineColor() >= 0)
1169 col = gROOT->GetColor(volume->GetLineColor());
1170
1171 if (volume->GetMedium() && (volume->GetMedium() != TGeoVolume::DummyMedium()) &&
1172 volume->GetMedium()->GetMaterial()) {
1173 auto material = volume->GetMedium()->GetMaterial();
1174
1175 node.material = material->GetName();
1176
1177 auto fillstyle = material->GetFillStyle();
1178 if ((fillstyle >= 3000) && (fillstyle <= 3100))
1179 node.opacity = (3100 - fillstyle) / 100.;
1180 if (!col)
1181 col = gROOT->GetColor(material->GetFillColor());
1182 } else {
1183 node.material.clear();
1184 }
1185
1186 if (col) {
1188 colbuf.Form("#%02x%02x%02x", (int)(col->GetRed() * 255), (int)(col->GetGreen() * 255),
1189 (int)(col->GetBlue() * 255));
1190 node.color = colbuf.Data();
1191 if (node.opacity == 1.)
1192 node.opacity = col->GetAlpha();
1193 } else {
1194 node.color.clear();
1195 }
1196}
1197
1198/////////////////////////////////////////////////////////////////////
1199/// Reset shape info, which used to pack binary data
1200
1202{
1203 for (auto &s : fShapes)
1204 s.reset();
1205}
1206
1207/////////////////////////////////////////////////////////////////////
1208/// Produce JSON string which can be directly used with `build`
1209/// function from JSROOT to create three.js model of configured geometry
1210///
1211/// Collect all information required to draw geometry on the client
1212/// This includes list of each visible nodes, meshes and matrixes
1213/// If @param all_nodes is true, all existing nodes will be provided,
1214/// which allows to create complete nodes hierarchy on client side
1215///
1216/// Example of usage:
1217///
1218/// void geom() {
1219/// auto f = TFile::Open("file_name.root");
1220/// auto vol = f->Get<TGeoVolume>("object_name");
1221/// ROOT::RGeomDescription desc;
1222/// desc.Build(vol);
1223/// std::ofstream fout("geom.json");
1224/// fout << desc.ProduceJson();
1225/// }
1226///
1227/// In JSROOT one loads data from JSON file and call `build` function to
1228/// produce three.js model. Also see example in tutorials/visualisation/webgui/geom/ folder
1229
1231{
1232 TLockGuard lock(fMutex);
1233
1234 std::vector<int> viscnt(fDesc.size(), 0);
1235
1236 int level = GetVisLevel();
1237
1238 // first count how many times each individual node appears
1239 int numnodes = ScanNodes(true, level, [&viscnt](RGeomNode &node, std::vector<int> &, bool, int) {
1240 viscnt[node.id]++;
1241 return true;
1242 });
1243
1244 if (GetMaxVisNodes() > 0) {
1245 while ((numnodes > GetMaxVisNodes()) && (level > 1)) {
1246 level--;
1247 viscnt.assign(viscnt.size(), 0);
1248 numnodes = ScanNodes(true, level, [&viscnt](RGeomNode &node, std::vector<int> &, bool, int) {
1249 viscnt[node.id]++;
1250 return true;
1251 });
1252 }
1253 }
1254
1255 fActualLevel = level;
1256 fDrawIdCut = 0;
1257
1258 int totalnumfaces = 0, totalnumnodes = 0;
1259
1260 // for (auto &node : fDesc)
1261 // node.SetDisplayed(false);
1262
1263 // build all shapes in volume decreasing order
1264 for (auto &sid : fSortMap) {
1265 fDrawIdCut++; //
1266 auto &desc = fDesc[sid];
1267
1268 if ((viscnt[sid] <= 0) || (desc.vol <= 0))
1269 continue;
1270
1271 auto shape = GetVolume(sid)->GetShape();
1272 if (!shape)
1273 continue;
1274
1275 // now we need to create TEveGeoPolyShape, which can provide all rendering data
1276 auto &shape_descr = MakeShapeDescr(shape);
1277
1278 // should not happen, but just in case
1279 if (shape_descr.nfaces <= 0) {
1280 R__LOG_ERROR(RGeomLog()) << "No faces for the shape " << shape->GetName() << " class " << shape->ClassName();
1281 continue;
1282 }
1283
1284 // check how many faces are created
1285 totalnumfaces += shape_descr.nfaces * viscnt[sid];
1286 if ((GetMaxVisFaces() > 0) && (totalnumfaces > GetMaxVisFaces()))
1287 break;
1288
1289 // also avoid too many nodes
1291 if ((GetMaxVisNodes() > 0) && (totalnumnodes > GetMaxVisNodes()))
1292 break;
1293
1294 // desc.SetDisplayed(true);
1295 }
1296
1297 // finally we should create data for streaming to the client
1298 // it includes list of visible nodes and rawdata
1299
1300 RGeomDrawing drawing;
1302 bool has_shape = false;
1303
1304 ScanNodes(true, level, [&, this](RGeomNode &node, std::vector<int> &stack, bool, int seqid) {
1305 if ((node.sortid < fDrawIdCut) && (viscnt[node.id] > 0)) {
1306 drawing.visibles.emplace_back(node.id, seqid, stack);
1307
1308 auto &item = drawing.visibles.back();
1309 item.color = node.color;
1310 item.opacity = node.opacity;
1311
1312 auto volume = GetVolume(node.id);
1313
1314 auto &sd = MakeShapeDescr(volume->GetShape());
1315
1316 item.ri = sd.rndr_info();
1317 if (sd.has_shape())
1318 has_shape = true;
1319 }
1320 return true;
1321 });
1322
1323 CollectNodes(drawing, all_nodes);
1324
1325 return MakeDrawingJson(drawing, has_shape);
1326}
1327
1328/////////////////////////////////////////////////////////////////////
1329/// Check if there is draw data available
1330
1332{
1333 TLockGuard lock(fMutex);
1334 return (fDrawJson.length() > 0) && (fDrawIdCut > 0);
1335}
1336
1337/////////////////////////////////////////////////////////////////////
1338/// Produces search data if necessary
1339
1341{
1342 TLockGuard lock(fMutex);
1343
1344 if (fSearch.empty() || !fSearchJson.empty())
1345 return;
1346
1347 std::string hjson;
1348
1350
1351 (void)hjson; // not used here
1352}
1353
1354/////////////////////////////////////////////////////////////////////
1355/// Collect all information required to draw geometry on the client
1356/// This includes list of each visible nodes, meshes and matrixes
1357
1359{
1360 auto json = ProduceJson();
1361
1362 TLockGuard lock(fMutex);
1363
1364 fDrawJson = "GDRAW:"s + json;
1365}
1366
1367/////////////////////////////////////////////////////////////////////
1368/// Clear drawing data.
1369/// Will be rebuild when next connection established or new message need to be send
1370
1372{
1373 TLockGuard lock(fMutex);
1374
1376}
1377
1378/////////////////////////////////////////////////////////////////////
1379/// Clear cached data, need to be clear when connection broken
1380
1382{
1383 TLockGuard lock(fMutex);
1384
1385 fShapes.clear();
1386 fSearch.clear();
1387
1389}
1390
1391/////////////////////////////////////////////////////////////////////
1392/// return true when node used in main geometry drawing and does not have childs
1393/// for such nodes one could provide optimize toggling of visibility flags
1394
1396{
1397 TLockGuard lock(fMutex);
1398
1399 if ((nodeid < 0) || (nodeid >= (int)fDesc.size()))
1400 return false;
1401
1402 auto &desc = fDesc[nodeid];
1403
1404 return (desc.sortid < fDrawIdCut) && desc.IsVisible() && desc.CanDisplay() && (desc.chlds.empty());
1405}
1406
1407/////////////////////////////////////////////////////////////////////
1408/// Search visible nodes for provided name
1409/// If number of found elements less than 100, create description and shapes for them
1410/// Returns number of match elements
1411
1412int RGeomDescription::SearchVisibles(const std::string &find, std::string &hjson, std::string &json)
1413{
1414 TLockGuard lock(fMutex);
1415
1416 hjson.clear();
1417 json.clear();
1418
1419 if (find.empty()) {
1420 hjson = "FOUND:RESET";
1421 return 0;
1422 }
1423
1424 std::vector<int> nodescnt(fDesc.size(), 0), viscnt(fDesc.size(), 0);
1425
1426 int nmatches = 0;
1427 std::string test = find;
1428 int kind = 0;
1429 if (test.compare(0, 2, "c:") == 0) {
1430 test.erase(0, 2);
1431 kind = 1;
1432 } else if (test.compare(0, 2, "m:") == 0) {
1433 test.erase(0, 2);
1434 kind = 2;
1435 }
1436
1437 TRegexp regexp(test.c_str());
1438
1439 auto match_func = [&regexp, kind](RGeomNode &node) {
1440 return (node.vol > 0) && (TString(node.GetArg(kind)).Index(regexp) >= 0);
1441 };
1442
1443 // first count how many times each individual node appears
1444 ScanNodes(false, 0,
1445 [&nodescnt, &viscnt, &match_func, &nmatches](RGeomNode &node, std::vector<int> &, bool is_vis, int) {
1446 if (match_func(node)) {
1447 nmatches++;
1448 nodescnt[node.id]++;
1449 if (is_vis)
1450 viscnt[node.id]++;
1451 };
1452 return true;
1453 });
1454
1455 // do not send too much data, limit could be made configurable later
1456 if (nmatches == 0) {
1457 hjson = "FOUND:NO";
1458 return nmatches;
1459 }
1460
1461 if ((GetMaxVisNodes() > 0) && (nmatches > 10 * GetMaxVisNodes())) {
1462 hjson = "FOUND:Too many " + std::to_string(nmatches);
1463 return nmatches;
1464 }
1465
1466 // now build all necessary shapes and check number of faces - not too many
1467
1468 int totalnumfaces = 0, totalnumnodes = 0, scnt = 0;
1469 bool send_rawdata = true;
1470
1471 // build all shapes in volume decreasing order
1472 for (auto &sid : fSortMap) {
1473 if (scnt++ < fDrawIdCut)
1474 continue; // no need to send most significant shapes
1475
1476 if (viscnt[sid] == 0)
1477 continue; // this node is not used at all
1478
1479 auto &desc = fDesc[sid];
1480 if ((viscnt[sid] <= 0) && (desc.vol <= 0))
1481 continue;
1482
1483 auto shape = GetVolume(sid)->GetShape();
1484 if (!shape)
1485 continue;
1486
1487 // create shape raw data
1488 auto &shape_descr = MakeShapeDescr(shape);
1489
1490 // should not happen, but just in case
1491 if (shape_descr.nfaces <= 0) {
1492 R__LOG_ERROR(RGeomLog()) << "No faces for the shape " << shape->GetName() << " class " << shape->ClassName();
1493 continue;
1494 }
1495
1496 // check how many faces are created
1497 totalnumfaces += shape_descr.nfaces * viscnt[sid];
1498 if ((GetMaxVisFaces() > 0) && (totalnumfaces > GetMaxVisFaces())) {
1499 send_rawdata = false;
1500 break;
1501 }
1502
1503 // also avoid too many nodes
1505 if ((GetMaxVisNodes() > 0) && (totalnumnodes > GetMaxVisNodes())) {
1506 send_rawdata = false;
1507 break;
1508 }
1509 }
1510
1511 // only for debug purposes - remove later
1512 // send_rawdata = false;
1513
1514 // finally we should create data for streaming to the client
1515 // it includes list of visible nodes and rawdata (if there is enough space)
1516
1517 std::vector<RGeomNodeBase> found_desc; ///<! hierarchy of nodes, used for search
1518 std::vector<int> found_map(fDesc.size(), -1); ///<! mapping between nodeid - > foundid
1519
1520 // these are only selected nodes to produce hierarchy
1521
1522 found_desc.emplace_back(0);
1523 found_desc[0].vis = fDesc[0].vis;
1524 found_desc[0].name = fDesc[0].name;
1525 found_desc[0].color = fDesc[0].color;
1526 found_map[0] = 0;
1527
1529
1530 RGeomDrawing drawing;
1531 bool has_shape = true;
1532
1533 ScanNodes(false, 0, [&, this](RGeomNode &node, std::vector<int> &stack, bool is_vis, int seqid) {
1534 // select only nodes which should match
1535 if (!match_func(node))
1536 return true;
1537
1538 // add entries into hierarchy of found elements
1539 int prntid = 0;
1540 for (auto &s : stack) {
1541 int chldid = fDesc[prntid].chlds[s];
1542 if (found_map[chldid] <= 0) {
1543 int newid = found_desc.size();
1544 found_desc.emplace_back(newid); // potentially original id can be used here
1545 found_map[chldid] = newid; // re-map into reduced hierarchy
1546
1547 found_desc.back().vis = fDesc[chldid].vis;
1548 found_desc.back().name = fDesc[chldid].name;
1549 found_desc.back().color = fDesc[chldid].color;
1550 found_desc.back().material = fDesc[chldid].material;
1551 }
1552
1553 auto pid = found_map[prntid];
1554 auto cid = found_map[chldid];
1555
1556 // now add entry into childs lists
1557 auto &pchlds = found_desc[pid].chlds;
1558 if (std::find(pchlds.begin(), pchlds.end(), cid) == pchlds.end())
1559 pchlds.emplace_back(cid);
1560
1561 prntid = chldid;
1562 }
1563
1564 // no need to add visibles
1565 if (!is_vis)
1566 return true;
1567
1568 drawing.visibles.emplace_back(node.id, seqid, stack);
1569
1570 // no need to transfer shape if it provided with main drawing list
1571 // also no binary will be transported when too many matches are there
1572 if (!send_rawdata || (node.sortid < fDrawIdCut)) {
1573 // do not include render data
1574 return true;
1575 }
1576
1577 auto &item = drawing.visibles.back();
1578 auto volume = GetVolume(node.id);
1579
1580 item.color = node.color;
1581 item.opacity = node.opacity;
1582
1583 auto &sd = MakeShapeDescr(volume->GetShape());
1584
1585 item.ri = sd.rndr_info();
1586 if (sd.has_shape())
1587 has_shape = true;
1588 return true;
1589 });
1590
1591 hjson = "FESCR:"s + TBufferJSON::ToJSON(&found_desc, GetJsonComp()).Data();
1592
1593 CollectNodes(drawing);
1594
1595 json = "FDRAW:"s + MakeDrawingJson(drawing, has_shape);
1596
1597 return nmatches;
1598}
1599
1600/////////////////////////////////////////////////////////////////////////////////
1601/// Returns nodeid for given stack array, returns -1 in case of failure
1602
1603int RGeomDescription::FindNodeId(const std::vector<int> &stack)
1604{
1605 TLockGuard lock(fMutex);
1606
1607 int nodeid = 0;
1608
1609 for (auto &chindx : stack) {
1610 auto &node = fDesc[nodeid];
1611 if (chindx >= (int)node.chlds.size())
1612 return -1;
1613 nodeid = node.chlds[chindx];
1614 }
1615
1616 return nodeid;
1617}
1618
1619/////////////////////////////////////////////////////////////////////////////////
1620/// Creates stack for given array of ids, first element always should be 0
1621
1622std::vector<int> RGeomDescription::MakeStackByIds(const std::vector<int> &ids)
1623{
1624 TLockGuard lock(fMutex);
1625
1626 std::vector<int> stack;
1627
1628 if (ids.empty())
1629 return stack;
1630
1631 if (ids[0] != 0) {
1632 printf("Wrong first id\n");
1633 return stack;
1634 }
1635
1636 int nodeid = 0;
1637
1638 for (unsigned k = 1; k < ids.size(); ++k) {
1639
1640 int prntid = nodeid;
1641 nodeid = ids[k];
1642
1643 if (nodeid >= (int)fDesc.size()) {
1644 printf("Wrong node id %d\n", nodeid);
1645 stack.clear();
1646 return stack;
1647 }
1648 auto &chlds = fDesc[prntid].chlds;
1649 auto pos = std::find(chlds.begin(), chlds.end(), nodeid);
1650 if (pos == chlds.end()) {
1651 printf("Wrong id %d not a child of %d - fail to find stack num %d\n", nodeid, prntid, (int)chlds.size());
1652 stack.clear();
1653 return stack;
1654 }
1655
1656 stack.emplace_back(std::distance(chlds.begin(), pos));
1657 }
1658
1659 return stack;
1660}
1661
1662/////////////////////////////////////////////////////////////////////////////////
1663/// Produce stack based on string path
1664/// Used to highlight geo volumes by browser hover event
1665
1666std::vector<int> RGeomDescription::MakeStackByPath(const std::vector<std::string> &path)
1667{
1668 TLockGuard lock(fMutex);
1669
1670 std::vector<int> res;
1671
1672 RGeomBrowserIter iter(*this);
1673
1674 if (iter.Navigate(path))
1675 res = MakeStackByIds(iter.CurrentIds());
1676
1677 return res;
1678}
1679
1680/////////////////////////////////////////////////////////////////////////////////
1681/// Produce list of node ids for given stack
1682/// If found nodes preselected - use their ids
1683
1684std::vector<int> RGeomDescription::MakeIdsByStack(const std::vector<int> &stack)
1685{
1686 TLockGuard lock(fMutex);
1687
1688 std::vector<int> ids;
1689
1690 ids.emplace_back(0);
1691 int nodeid = 0;
1692 bool failure = false;
1693
1694 for (auto s : stack) {
1695 auto &chlds = fDesc[nodeid].chlds;
1696 if (s >= (int)chlds.size()) {
1697 failure = true;
1698 break;
1699 }
1700
1701 ids.emplace_back(chlds[s]);
1702
1703 nodeid = chlds[s];
1704 }
1705
1706 if (failure) {
1707 printf("Fail to convert stack into list of nodes\n");
1708 ids.clear();
1709 }
1710
1711 return ids;
1712}
1713
1714/////////////////////////////////////////////////////////////////////////////////
1715/// Returns path string for provided stack
1716
1717std::vector<std::string> RGeomDescription::MakePathByStack(const std::vector<int> &stack)
1718{
1719 TLockGuard lock(fMutex);
1720
1721 std::vector<std::string> path;
1722
1723 auto ids = MakeIdsByStack(stack);
1724 path.reserve(ids.size());
1725for (auto &id : ids)
1726 path.emplace_back(fDesc[id].name);
1727
1728 return path;
1729}
1730
1731/////////////////////////////////////////////////////////////////////////////////
1732/// Return string with only part of nodes description which were modified
1733/// Checks also volume
1734
1736{
1737 TLockGuard lock(fMutex);
1738
1739 std::vector<RGeomNodeBase *> nodes;
1740 auto vol = GetVolume(nodeid);
1741
1742 // we take not only single node, but all there same volume is referenced
1743 // nodes.push_back(&fDesc[nodeid]);
1744
1745 int id = 0;
1746 for (auto &desc : fDesc)
1747 if (GetVolume(id++) == vol)
1748 nodes.emplace_back(&desc);
1749
1750 return "MODIF:"s + TBufferJSON::ToJSON(&nodes, GetJsonComp()).Data();
1751}
1752
1753/////////////////////////////////////////////////////////////////////////////////
1754/// Produce shape rendering data for given stack
1755/// All nodes, which are referencing same shape will be transferred
1756/// Returns true if new render information provided
1757
1758bool RGeomDescription::ProduceDrawingFor(int nodeid, std::string &json, bool check_volume)
1759{
1760 TLockGuard lock(fMutex);
1761
1762 // only this shape is interesting
1763
1764 TGeoVolume *vol = (nodeid < 0) ? nullptr : GetVolume(nodeid);
1765
1766 if (!vol || !vol->GetShape()) {
1767 json.append("NO");
1768 return false;
1769 }
1770
1771 RGeomDrawing drawing;
1772
1773 ScanNodes(true, 0, [&, this](RGeomNode &node, std::vector<int> &stack, bool, int seq_id) {
1774 // select only nodes which reference same shape
1775
1776 if (check_volume) {
1777 if (GetVolume(node.id) != vol)
1778 return true;
1779 } else {
1780 if (node.id != nodeid)
1781 return true;
1782 }
1783
1784 drawing.visibles.emplace_back(node.id, seq_id, stack);
1785
1786 auto &item = drawing.visibles.back();
1787
1788 item.color = node.color;
1789 item.opacity = node.opacity;
1790 return true;
1791 });
1792
1793 // no any visible nodes were done
1794 if (drawing.visibles.empty()) {
1795 json.append("NO");
1796 return false;
1797 }
1798
1800
1801 bool has_shape = false, has_raw = false;
1802
1803 auto &sd = MakeShapeDescr(vol->GetShape());
1804
1805 // assign shape data
1806 for (auto &item : drawing.visibles) {
1807 item.ri = sd.rndr_info();
1808 if (sd.has_shape())
1809 has_shape = true;
1810 if (sd.has_raw())
1811 has_raw = true;
1812 }
1813
1814 CollectNodes(drawing);
1815
1816 json.append(MakeDrawingJson(drawing, has_shape));
1817
1818 return has_raw || has_shape;
1819}
1820
1821/////////////////////////////////////////////////////////////////////////////////
1822/// Produce JSON for the drawing
1823/// If TGeoShape appears in the drawing, one has to keep typeinfo
1824/// But in this case one can exclude several classes which are not interesting,
1825/// but appears very often
1826
1828{
1829 int comp = GetJsonComp();
1830
1832 return TBufferJSON::ToJSON(&drawing, comp).Data();
1833
1834 comp = comp % TBufferJSON::kSkipTypeInfo; // no typeinfo skipping
1835
1838 json.SetSkipClassInfo(TClass::GetClass<RGeomDrawing>());
1839 json.SetSkipClassInfo(TClass::GetClass<RGeomNode>());
1840 json.SetSkipClassInfo(TClass::GetClass<RGeomVisible>());
1841 json.SetSkipClassInfo(TClass::GetClass<RGeomShapeRenderInfo>());
1842 json.SetSkipClassInfo(TClass::GetClass<RGeomRawRenderInfo>());
1843
1844 return json.StoreObject(&drawing, TClass::GetClass<RGeomDrawing>()).Data();
1845}
1846
1847/////////////////////////////////////////////////////////////////////////////////
1848/// Change visibility for specified element
1849/// Returns true if changes was performed
1850
1851bool RGeomDescription::ChangeNodeVisibility(const std::vector<std::string> &path, bool selected)
1852{
1853 TLockGuard lock(fMutex);
1854
1855 RGeomBrowserIter giter(*this);
1856 if (!giter.Navigate(path))
1857 return false;
1858
1859 auto nodeid = giter.GetNodeId();
1860
1861 auto &dnode = fDesc[nodeid];
1862
1863 auto vol = GetVolume(nodeid);
1864
1865 // nothing changed
1866 if (vol->IsVisible() == selected)
1867 return false;
1868
1869 dnode.vis = selected ? 99 : 0;
1870 vol->SetVisibility(selected);
1871 if (!dnode.chlds.empty()) {
1872 if (selected)
1873 dnode.vis = 1; // visibility disabled when any child
1874 vol->SetVisDaughters(selected);
1875 }
1876
1877 int id = 0;
1878 for (auto &desc : fDesc)
1879 if (GetVolume(id++) == vol)
1880 desc.vis = dnode.vis;
1881
1882 auto stack = MakeStackByIds(giter.CurrentIds());
1883
1884 // any change in logical node visibility erase individual physical node settings
1885 for (auto iter = fVisibility.begin(); iter != fVisibility.end(); iter++)
1886 if (compare_stacks(iter->stack, stack) == 0) {
1887 fVisibility.erase(iter);
1888 break;
1889 }
1890
1891 _ClearDrawData(); // after change raw data is no longer valid
1892
1893 return true;
1894}
1895
1896/////////////////////////////////////////////////////////////////////////////////
1897/// Change visibility for specified element
1898/// Returns true if changes was performed
1899
1900std::unique_ptr<RGeomNodeInfo> RGeomDescription::MakeNodeInfo(const std::vector<int> &stack)
1901{
1902 auto path = MakePathByStack(stack);
1903
1904 TLockGuard lock(fMutex);
1905
1906 std::unique_ptr<RGeomNodeInfo> res;
1907
1908 RGeomBrowserIter iter(*this);
1909
1910 if (iter.Navigate(path)) {
1911
1912 auto node = fNodes[iter.GetNodeId()];
1913
1914 auto &desc = fDesc[iter.GetNodeId()];
1915
1916 res = std::make_unique<RGeomNodeInfo>();
1917
1918 res->path = path;
1919 res->node_name = node ? node->GetName() : "node_name";
1920 res->node_type = node ? node->ClassName() : "no class";
1921
1922 auto vol = GetVolume(iter.GetNodeId());
1923
1924 TGeoShape *shape = vol ? vol->GetShape() : nullptr;
1925
1926 if (shape) {
1927 res->shape_name = shape->GetName();
1928 res->shape_type = shape->ClassName();
1929 }
1930
1931 if (shape && desc.CanDisplay()) {
1932
1933 auto &shape_descr = MakeShapeDescr(shape);
1934
1935 res->ri = shape_descr.rndr_info(); // temporary pointer, can be used preserved for short time
1936 }
1937 }
1938
1939 return res;
1940}
1941
1942/////////////////////////////////////////////////////////////////////////////////
1943/// Select top node by path
1944/// Used by the client to change active node
1945/// Returns true if selected node was changed
1946
1947bool RGeomDescription::SelectTop(const std::vector<std::string> &path)
1948{
1949 TLockGuard lock(fMutex);
1950
1951 RGeomBrowserIter iter(*this);
1952
1953 if (!iter.Navigate(path))
1954 return false;
1955
1956 auto stack = MakeStackByIds(iter.CurrentIds());
1957 if (stack == fSelectedStack)
1958 return false;
1959
1960 fSelectedStack = stack;
1961
1963
1964 return true;
1965}
1966
1967/////////////////////////////////////////////////////////////////////////////////
1968/// Set visibility of physical node by path
1969/// It overrules TGeo visibility flags - but only for specific physical node
1970
1971bool RGeomDescription::SetPhysNodeVisibility(const std::vector<std::string> &path, bool on)
1972{
1973 TLockGuard lock(fMutex);
1974
1975 RGeomBrowserIter giter(*this);
1976
1977 if (!giter.Navigate(path))
1978 return false;
1979
1980 auto stack = MakeStackByIds(giter.CurrentIds());
1981
1982 auto nodeid = giter.GetNodeId();
1983
1984 for (auto iter = fVisibility.begin(); iter != fVisibility.end(); iter++) {
1985 auto res = compare_stacks(iter->stack, stack);
1986
1987 if (res == 0) {
1988 bool changed = iter->visible != on;
1989 if (changed) {
1990 iter->visible = on;
1992
1993 // no need for custom settings if match with description
1994 if ((fDesc[nodeid].vis > 0) == on)
1995 fVisibility.erase(iter);
1996 }
1997
1998 return changed;
1999 }
2000
2001 if (res > 0) {
2002 fVisibility.emplace(iter, stack, on);
2004 return true;
2005 }
2006 }
2007
2008 fVisibility.emplace_back(stack, on);
2010 return true;
2011}
2012
2013/////////////////////////////////////////////////////////////////////////////////
2014/// Set visibility of physical node by itemname
2015/// itemname in string with path like "/TOP_1/SUB_2/NODE_3"
2016
2018{
2019 std::vector<std::string> path;
2020 std::string::size_type p1 = 0;
2021
2022 while (p1 < itemname.length()) {
2023 if (itemname[p1] == '/') {
2024 p1++;
2025 continue;
2026 }
2027 auto p = itemname.find('/', p1);
2028 if (p == std::string::npos) {
2029 path.emplace_back(itemname.substr(p1));
2030 p1 = itemname.length();
2031 } else {
2032 path.emplace_back(itemname.substr(p1, p - p1));
2033 p1 = p + 1;
2034 }
2035 }
2036
2037 return SetPhysNodeVisibility(path, on);
2038}
2039
2040/////////////////////////////////////////////////////////////////////////////////
2041/// Check if there special settings for specified physical node
2042/// returns -1 if nothing is found
2043
2044int RGeomDescription::IsPhysNodeVisible(const std::vector<int> &stack)
2045{
2046 for (auto &item : fVisibility) {
2047 unsigned sz = item.stack.size();
2048 if (stack.size() < sz)
2049 continue;
2050 bool match = true;
2051 for (unsigned n = 0; n < sz; ++n)
2052 if (stack[n] != item.stack[n]) {
2053 match = false;
2054 break;
2055 }
2056
2057 if (match)
2058 return item.visible ? 1 : 0;
2059 }
2060 return -1;
2061}
2062
2063/////////////////////////////////////////////////////////////////////////////////
2064/// Reset custom visibility of physical node by path
2065
2066bool RGeomDescription::ClearPhysNodeVisibility(const std::vector<std::string> &path)
2067{
2068 TLockGuard lock(fMutex);
2069
2070 RGeomBrowserIter giter(*this);
2071
2072 if (!giter.Navigate(path))
2073 return false;
2074
2075 auto stack = MakeStackByIds(giter.CurrentIds());
2076
2077 for (auto iter = fVisibility.begin(); iter != fVisibility.end(); iter++)
2078 if (compare_stacks(iter->stack, stack) == 0) {
2079 fVisibility.erase(iter);
2081 return true;
2082 }
2083
2084 return false;
2085}
2086
2087/////////////////////////////////////////////////////////////////////////////////
2088/// Reset all custom visibility settings
2089
2091{
2092 TLockGuard lock(fMutex);
2093
2094 if (fVisibility.empty())
2095 return false;
2096
2097 fVisibility.clear();
2099 return true;
2100}
2101
2102/////////////////////////////////////////////////////////////////////////////////
2103/// Change configuration by client
2104/// Returns true if any parameter was really changed
2105
2107{
2108 auto cfg = TBufferJSON::FromJSON<RGeomConfig>(json);
2109 if (!cfg)
2110 return false;
2111
2112 TLockGuard lock(fMutex);
2113
2114 auto json1 = TBufferJSON::ToJSON(cfg.get());
2116
2117 if (json1 == json2)
2118 return false;
2119
2120 fCfg = *cfg; // use assign
2121
2123
2124 return true;
2125}
2126
2127/////////////////////////////////////////////////////////////////////////////////
2128/// Change search query and belongs to it json string
2129/// Returns true if any parameter was really changed
2130
2131bool RGeomDescription::SetSearch(const std::string &query, const std::string &json)
2132{
2133 TLockGuard lock(fMutex);
2134
2135 bool changed = (fSearch != query) || (fSearchJson != json);
2136 fSearch = query;
2137 fSearchJson = json;
2138 return changed;
2139}
2140
2141/////////////////////////////////////////////////////////////////////////////////
2142/// Save geometry configuration as C++ macro
2143
2144void RGeomDescription::SavePrimitive(std::ostream &fs, const std::string &name)
2145{
2146 std::string prefix = " ";
2147
2148 if (fCfg.vislevel != 0)
2149 fs << prefix << name << "SetVisLevel(" << fCfg.vislevel << ");" << std::endl;
2150 if (fCfg.maxnumnodes != 0)
2151 fs << prefix << name << "SetMaxVisNodes(" << fCfg.maxnumnodes << ");" << std::endl;
2152 if (fCfg.maxnumfaces != 0)
2153 fs << prefix << name << "SetMaxVisFaces(" << fCfg.maxnumfaces << ");" << std::endl;
2154 if (fCfg.showtop)
2155 fs << prefix << name << "SetTopVisible(true);" << std::endl;
2156 if (fCfg.build_shapes != 1)
2157 fs << prefix << name << "SetBuildShapes(" << fCfg.build_shapes << ");" << std::endl;
2158 if (fCfg.nsegm != 0)
2159 fs << prefix << name << "SetNSegments(" << fCfg.nsegm << ");" << std::endl;
2160 if (!fCfg.drawopt.empty())
2161 fs << prefix << name << "SetDrawOptions(\"" << fCfg.drawopt << "\");" << std::endl;
2162 if (fJsonComp != 0)
2163 fs << prefix << name << "SetJsonComp(" << fJsonComp << ");" << std::endl;
2164
2165 // store custom visibility flags
2166 for (auto &item : fVisibility) {
2167 auto path = MakePathByStack(item.stack);
2168 fs << prefix << name << "SetPhysNodeVisibility(";
2169 for (int i = 0; i < (int)path.size(); ++i)
2170 fs << (i == 0 ? "{\"" : ", \"") << path[i] << "\"";
2171 fs << "}, " << (item.visible ? "true" : "false") << ");" << std::endl;
2172 }
2173}
nlohmann::json json
std::unique_ptr< RootCsg::TBaseMesh > MakeGeoMesh(TGeoMatrix *matr, TGeoShape *shape)
Function produces mesh for provided shape, applying matrix to the result.
#define R__LOG_ERROR(...)
Definition RLogger.hxx:357
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
#define e(i)
Definition RSha256.hxx:103
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
winID h TVirtualViewer3D TVirtualGLPainter p
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 offset
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void on
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 Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t UChar_t len
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void funcs
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize fs
char name[80]
Definition TGX11.cxx:110
R__EXTERN TGeoManager * gGeoManager
#define gROOT
Definition TROOT.h:406
Reply on browser request.
Iterator of hierarchical geometry structures.
Definition RGeomData.cxx:60
RGeomBrowserIter(RGeomDescription &desc)
Definition RGeomData.cxx:71
std::vector< int > CurrentIds() const
Returns array of ids to currently selected node.
bool Navigate(TGeoVolume *vol)
Navigate to specified volume - find first occurrence.
std::vector< int > fStackChilds
Definition RGeomData.cxx:68
bool Navigate(const std::string &path)
Navigate to specified path - path specified as string and should start with "/".
const std::string & GetMaterial() const
Definition RGeomData.cxx:77
std::vector< int > fStackParents
Definition RGeomData.cxx:67
bool HasChilds() const
Definition RGeomData.cxx:85
RGeomDescription & fDesc
Definition RGeomData.cxx:62
const std::string & GetColor() const
Definition RGeomData.cxx:75
const std::string & GetName() const
Definition RGeomData.cxx:73
bool Navigate(const std::vector< std::string > &path)
Navigate to specified path
bool showtop
show geometry top volume, off by default
int maxnumfaces
maximal number of faces
int vislevel
visible level
int maxnumnodes
maximal number of nodes
std::string drawopt
draw options for TGeoPainter
int build_shapes
when shapes build on server 0 - never, 1 - TGeoComposite, 2 - plus non-cylindrical,...
int nsegm
number of segments for cylindrical shapes
RGeomConfig fCfg
! configuration parameter editable from GUI
std::vector< std::pair< const void *, RGeomSignalFunc_t > > fSignals
! registered signals
int IsPhysNodeVisible(const std::vector< int > &stack)
Check if there special settings for specified physical node returns -1 if nothing is found.
std::vector< int > fSelectedStack
! selected branch of geometry by stack
void SetMaxVisNodes(int cnt)
Set maximal number of nodes which should be selected for drawing.
std::string ProcessBrowserRequest(const std::string &req="")
Find description object for requested shape If not exists - will be created.
std::vector< RGeomNode > fDesc
! converted description, send to client
void PackMatrix(std::vector< float > &arr, TGeoMatrix *matr)
Pack matrix into vector, which can be send to client Following sizes can be used for vector: 0 - Iden...
bool ProduceDrawingFor(int nodeid, std::string &json, bool check_volume=false)
Produce shape rendering data for given stack All nodes, which are referencing same shape will be tran...
std::unique_ptr< RGeomNodeInfo > MakeNodeInfo(const std::vector< int > &stack)
Change visibility for specified element Returns true if changes was performed.
bool HasDrawData() const
Check if there is draw data available.
std::vector< int > MakeIdsByStack(const std::vector< int > &stack)
Produce list of node ids for given stack If found nodes preselected - use their ids.
int MarkVisible(bool on_screen=false)
Set visibility flag for each nodes.
void SetVisLevel(int lvl=3)
Set maximal visible level.
void IssueSignal(const void *handler, const std::string &kind)
Issue signal, which distributed on all handlers - excluding source handler.
int GetUsedNSegments(int min=20)
Returns really used number of cylindrical segments.
bool IsPrincipalEndNode(int nodeid)
return true when node used in main geometry drawing and does not have childs for such nodes one could...
bool SetSearch(const std::string &query, const std::string &json)
Change search query and belongs to it json string Returns true if any parameter was really changed.
std::vector< RGeomNodeVisibility > fVisibility
! custom visibility flags for physical nodes
bool SelectTop(const std::vector< std::string > &path)
Select top node by path Used by the client to change active node Returns true if selected node was ch...
int GetMaxVisNodes() const
Returns maximal visible number of nodes, ignored when non-positive.
int GetVisLevel() const
Returns maximal visible level.
int GetMaxVisFaces() const
Returns maximal visible number of faces, ignored when non-positive.
void ClearCache()
Clear cached data, need to be clear when connection broken.
void ClearDescription()
Clear geometry description.
std::vector< int > MakeStackByIds(const std::vector< int > &ids)
Creates stack for given array of ids, first element always should be 0.
void SetMaxVisFaces(int cnt)
Set maximal number of faces which should be selected for drawing.
bool IsPreferredOffline() const
Is offline operations preferred.
std::vector< ShapeDescr > fShapes
! shapes with created descriptions
int fJsonComp
! default JSON compression
bool ChangeNodeVisibility(const std::vector< std::string > &path, bool on)
Change visibility for specified element Returns true if changes was performed.
std::string fSearch
! search string in hierarchy
std::string fSearchJson
! drawing json for search
void SavePrimitive(std::ostream &fs, const std::string &name)
Save geometry configuration as C++ macro.
bool ClearAllPhysVisibility()
Reset all custom visibility settings.
std::string MakeDrawingJson(RGeomDrawing &drawing, bool has_shapes=false)
Produce JSON for the drawing If TGeoShape appears in the drawing, one has to keep typeinfo But in thi...
int fActualLevel
! level can be reduced when selecting nodes
TGeoVolume * GetVolume(int nodeid)
Get volume for specified nodeid If specific volume was configured, it will be returned for nodeid==0.
int GetNumNodes() const
Number of unique nodes in the geometry.
void ProduceDrawData()
Collect all information required to draw geometry on the client This includes list of each visible no...
void SetNSegments(int n=0)
Set number of segments for cylindrical shapes, if 0 - default value will be used.
bool SetPhysNodeVisibility(const std::vector< std::string > &path, bool on=true)
Set visibility of physical node by path It overrules TGeo visibility flags - but only for specific ph...
bool ClearPhysNodeVisibility(const std::vector< std::string > &path)
Reset custom visibility of physical node by path.
void BuildDescription(TGeoNode *topnode, TGeoVolume *topvolume)
Build geometry description.
int SearchVisibles(const std::string &find, std::string &hjson, std::string &json)
Search visible nodes for provided name If number of found elements less than 100, create description ...
std::vector< int > fSortMap
! nodes in order large -> smaller volume
std::string ProduceJson(bool all_nodes=false)
Produce JSON string which can be directly used with build function from JSROOT to create three....
void ClearDrawData()
Clear drawing data.
std::vector< std::string > MakePathByStack(const std::vector< int > &stack)
Returns path string for provided stack.
TVirtualMutex * fMutex
! external mutex used to protect all data
void AddSignalHandler(const void *handler, RGeomSignalFunc_t func)
Add signal handler.
bool ChangeConfiguration(const std::string &json)
Change configuration by client Returns true if any parameter was really changed.
void CopyMaterialProperties(TGeoVolume *vol, RGeomNode &node)
Copy material properties.
std::vector< int > MakeStackByPath(const std::vector< std::string > &path)
Produce stack based on string path Used to highlight geo volumes by browser hover event.
int ScanNodes(bool only_visible, int maxlvl, RGeomScanFunc_t func)
Iterate over all nodes and call function for visible.
int IsBuildShapes() const
Returns true if binary 3D model build already by C++ server (default)
void ProduceSearchData()
Produces search data if necessary.
void CollectNodes(RGeomDrawing &drawing, bool all_nodes=false)
Collect nodes which are used in visibles.
int fDrawIdCut
! sortid used for selection of most-significant nodes
TGeoVolume * fDrawVolume
! select volume independent from TGeoManager
int CountShapeFaces(TGeoShape *shape)
Count number of faces for the shape.
ShapeDescr & MakeShapeDescr(TGeoShape *shape)
Find description object and create render information.
void _ClearDrawData()
clear drawing data without locking mutex
ShapeDescr & FindShapeDescr(TGeoShape *shape)
Find description object for requested shape If not exists - will be created.
std::string fDrawJson
! JSON with main nodes drawn by client
void RemoveSignalHandler(const void *handler)
Remove signal handler.
int FindNodeId(const std::vector< int > &stack)
Returns nodeid for given stack array, returns -1 in case of failure.
void ProduceIdShifts()
Count total number of visible childs under each node.
std::string ProduceModifyReply(int nodeid)
Return string with only part of nodes description which were modified Checks also volume.
void ResetRndrInfos()
Reset shape info, which used to pack binary data.
int GetNSegments() const
Return of segments for cylindrical shapes, if 0 - default value will be used.
std::vector< TGeoNode * > fNodes
! flat list of all nodes
int GetJsonComp() const
Returns JSON compression level for data transfer.
void Build(TGeoManager *mgr, const std::string &volname="")
Collect information about geometry hierarchy into flat list like it done in JSROOT ClonedNodes....
Object with full description for drawing geometry It includes list of visible items and list of nodes...
int numnodes
total number of nodes in description
std::vector< RGeomVisible > visibles
all visible items
RGeomConfig * cfg
current configurations
std::vector< RGeomNode * > nodes
all used nodes to display visible items and not known for client
std::string material
name of the material
Definition RGeomData.hxx:50
int sortid
! place in sorted array, to check cuts, or id of original node when used search structures
Definition RGeomData.hxx:51
std::string color
rgb code in hex format
Definition RGeomData.hxx:49
int id
node id, index in array
Definition RGeomData.hxx:43
Full node description including matrices and other attributes.
Definition RGeomData.hxx:68
float opacity
! opacity of the color
Definition RGeomData.hxx:75
A log configuration for a channel, e.g.
Definition RLogger.hxx:98
const_iterator begin() const
const_iterator end() const
virtual Color_t GetFillColor() const
Return the fill area color.
Definition TAttFill.h:31
virtual Color_t GetLineColor() const
Return the line color.
Definition TAttLine.h:35
Class for serializing object to and from JavaScript Object Notation (JSON) format.
Definition TBufferJSON.h:30
static TString ToJSON(const T *obj, Int_t compact=0, const char *member_name=nullptr)
Definition TBufferJSON.h:75
@ kSkipTypeInfo
do not store typenames in JSON
Definition TBufferJSON.h:48
void SetCompact(int level)
Set level of space/newline/array compression Lower digit of compact parameter define formatting rules...
The color creation and management class.
Definition TColor.h:22
Float_t GetRed() const
Definition TColor.h:61
static Int_t GetColor(const char *hexcolor)
Static method returning color number for color specified by hex color string of form: "#rrggbb",...
Definition TColor.cxx:1924
Float_t GetAlpha() const
Definition TColor.h:67
Float_t GetBlue() const
Definition TColor.h:63
Float_t GetGreen() const
Definition TColor.h:62
@ kVisNone
Definition TGeoAtt.h:25
Box class.
Definition TGeoBBox.h:17
Composite shapes are Boolean combinations of two or more shape components.
static TClass * Class()
A cone segment is a cone having a range in phi.
Definition TGeoCone.h:99
static TClass * Class()
The cones are defined by 5 parameters:
Definition TGeoCone.h:17
static TClass * Class()
The cut tubes constructor has the form:
Definition TGeoTube.h:173
static TClass * Class()
static TClass * Class()
Matrix class used for computing global transformations Should NOT be used for node definition.
Definition TGeoMatrix.h:458
void Multiply(const TGeoMatrix *right)
multiply to the right with an other transformation if right is identity matrix, just return
A hyperboloid is represented as a solid limited by two planes perpendicular to the Z axis (top and bo...
Definition TGeoHype.h:17
static TClass * Class()
static TClass * Class()
A geometry iterator.
Definition TGeoNode.h:248
void Skip()
Stop iterating the current branch.
The manager class for any TGeo geometry.
Definition TGeoManager.h:44
void SetNsegments(Int_t nseg)
Set number of segments for approximating circles in drawing.
Int_t GetNsegments() const
Get number of segments approximating circles.
Geometrical transformation package.
Definition TGeoMatrix.h:38
virtual const Double_t * GetTranslation() const =0
TClass * IsA() const override
Definition TGeoMatrix.h:104
virtual void LocalToMaster(const Double_t *local, Double_t *master) const
convert a point by multiplying its column vector (x, y, z, 1) to matrix inverse
virtual const Double_t * GetScale() const =0
Bool_t IsIdentity() const
Definition TGeoMatrix.h:63
virtual const Double_t * GetRotationMatrix() const =0
TGeoMaterial * GetMaterial() const
Definition TGeoMedium.h:49
A node represent a volume positioned inside another.They store links to both volumes and to the TGeoM...
Definition TGeoNode.h:39
A paraboloid is defined by the revolution surface generated by a parabola and is bounded by two plane...
static TClass * Class()
A polycone is represented by a sequence of tubes/cones, glued together at defined Z planes.
Definition TGeoPcon.h:17
static TClass * Class()
Polygons are defined in the same way as polycones, the difference being just that the segments betwee...
Definition TGeoPgon.h:20
static TClass * Class()
static TClass * Class()
static TClass * Class()
A shape scaled by a TGeoScale transformation.
static TClass * Class()
Base abstract class for all shapes.
Definition TGeoShape.h:25
virtual Bool_t IsComposite() const
Definition TGeoShape.h:138
virtual Bool_t IsCylType() const =0
const char * GetName() const override
Get the shape name.
TClass * IsA() const override
Definition TGeoShape.h:179
virtual TBuffer3D * MakeBuffer3D() const
Definition TGeoShape.h:154
TGeoSphere are not just balls having internal and external radii, but sectors of a sphere having defi...
Definition TGeoSphere.h:17
static TClass * Class()
static TClass * Class()
Tessellated solid class.
static TClass * Class()
The torus is defined by its axial radius, its inner and outer radius.
Definition TGeoTorus.h:17
static TClass * Class()
static TClass * Class()
A tube segment is a tube having a range in phi.
Definition TGeoTube.h:94
static TClass * Class()
Cylindrical tube class.
Definition TGeoTube.h:17
static TClass * Class()
static TClass * Class()
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition TGeoVolume.h:43
TGeoMedium * GetMedium() const
Definition TGeoVolume.h:175
TObjArray * GetNodes()
Definition TGeoVolume.h:169
static TGeoMedium * DummyMedium()
TGeoShape * GetShape() const
Definition TGeoVolume.h:190
A TGeoXtru shape is represented by the extrusion of an arbitrary polygon with fixed outline between s...
Definition TGeoXtru.h:22
static TClass * Class()
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:225
Regular expression class.
Definition TRegexp.h:31
Basic string class.
Definition TString.h:139
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:651
const Int_t n
Definition legend1.C:16
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
RLogChannel & RGeomLog()
Log channel for Geomviewer diagnostics.
Definition RGeomData.cxx:49
std::function< bool(RGeomNode &, std::vector< int > &, bool, int)> RGeomScanFunc_t
std::function< void(const std::string &)> RGeomSignalFunc_t
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
Definition TMath.h:697
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:250
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:666
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:123