Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoVolume.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Andrei Gheata 30/05/02
3// Divide(), CheckOverlaps() implemented by Mihaela Gheata
4
5/*************************************************************************
6 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13/** \class TGeoVolume
14\ingroup Shapes_classes
15
16TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes
17
18 Volumes are the basic objects used in building the geometrical hierarchy.
19They represent unpositioned objects but store all information about the
20placement of the other volumes they may contain. Therefore a volume can
21be replicated several times in the geometry. In order to create a volume, one
22has to put together a shape and a medium which are already defined. Volumes
23have to be named by users at creation time. Every different name may represent a
24an unique volume object, but may also represent more general a family (class)
25of volume objects having the same shape type and medium, but possibly
26different shape parameters. It is the user's task to provide different names
27for different volume families in order to avoid ambiguities at tracking time.
28A generic family rather than a single volume is created only in two cases :
29when a generic shape is provided to the volume constructor or when a division
30operation is applied. Each volume in the geometry stores an unique
31ID corresponding to its family. In order to ease-up their creation, the manager
32class is providing an API that allows making a shape and a volume in a single step.
33
34 Volumes are objects that can be visualized, therefore having visibility,
35colour, line and fill attributes that can be defined or modified any time after
36the volume creation. It is advisable however to define these properties just
37after the first creation of a volume namespace, since in case of volume families
38any new member created by the modeler inherits these properties.
39
40 In order to provide navigation features, volumes have to be able to find
41the proper container of any point defined in the local reference frame. This
42can be the volume itself, one of its positioned daughter volumes or none if
43the point is actually outside. On the other hand, volumes have to provide also
44other navigation methods such as finding the distances to its shape boundaries
45or which daughter will be crossed first. The implementation of these features
46is done at shape level, but the local mother-daughters management is handled
47by volumes that builds additional optimisation structures upon geometry closure.
48In order to have navigation features properly working one has to follow the
49general rules for building a valid geometry (see TGeoManager class).
50
51 Now let's make a simple volume representing a copper wire. We suppose that
52a medium is already created (see TGeoMedium class on how to create media).
53We will create a TUBE shape for our wire, having Rmin=0cm, Rmax=0.01cm
54and a half-length dZ=1cm :
55
56~~~ {.cpp}
57 TGeoTube *tube = new TGeoTube("wire_tube", 0, 0.01, 1);
58~~~
59
60One may omit the name for the shape if no retrieving by name is further needed
61during geometry building. The same shape can be shared by different volumes
62having different names and materials. Now let's make the volume for our wire.
63The prototype for volumes constructor looks like :
64
65 TGeoVolume::TGeoVolume(const char *name, TGeoShape *shape, TGeoMedium *med)
66
67Since TGeoTube derives from the base shape class, we can provide it to the volume
68constructor :
69
70~~~ {.cpp}
71 TGeoVolume *wire_co = new TGeoVolume("WIRE_CO", tube, ptrCOPPER);
72~~~
73
74Do not bother to delete neither the media, shapes or volumes that you have
75created since all will be automatically cleaned on exit by the manager class.
76If we would have taken a look inside TGeoManager::MakeTube() method, we would
77have been able to create our wire with a single line :
78
79~~~ {.cpp}
80 TGeoVolume *wire_co = gGeoManager->MakeTube("WIRE_CO", ptrCOPPER, 0, 0.01, 1);
81~~~
82
83The same applies for all primitive shapes, for which there can be found
84corresponding MakeSHAPE() methods. Their usage is much more convenient unless
85a shape has to be shared between more volumes. Let's make now an aluminium wire
86having the same shape, supposing that we have created the copper wire with the
87line above :
88
89~~~ {.cpp}
90 TGeoVolume *wire_al = new TGeoVolume("WIRE_AL", wire_co->GetShape(), ptrAL);
91~~~
92
93Now that we have learned how to create elementary volumes, let's see how we
94can create a geometrical hierarchy.
95
96
97### Positioning volumes
98
99 When creating a volume one does not specify if this will contain or not other
100volumes. Adding daughters to a volume implies creating those and adding them
101one by one to the list of daughters. Since the volume has to know the position
102of all its daughters, we will have to supply at the same time a geometrical
103transformation with respect to its local reference frame for each of them.
104The objects referencing a volume and a transformation are called NODES and
105their creation is fully handled by the modeler. They represent the link
106elements in the hierarchy of volumes. Nodes are unique and distinct geometrical
107objects ONLY from their container point of view. Since volumes can be replicated
108in the geometry, the same node may be found on different branches.
109
110\image html geom_t_example.png width=600px
111
112 An important observation is that volume objects are owned by the TGeoManager
113class. This stores a list of all volumes in the geometry, that is cleaned
114upon destruction.
115
116 Let's consider positioning now our wire in the middle of a gas chamber. We
117need first to define the gas chamber :
118
119~~~ {.cpp}
120 TGeoVolume *chamber = gGeoManager->MakeTube("CHAMBER", ptrGAS, 0, 1, 1);
121~~~
122
123Now we can put the wire inside :
124
125~~~ {.cpp}
126 chamber->AddNode(wire_co, 1);
127~~~
128
129If we inspect now the chamber volume in a browser, we will notice that it has
130one daughter. Of course the gas has some container also, but let's keep it like
131that for the sake of simplicity. The full prototype of AddNode() is :
132
133~~~ {.cpp}
134 TGeoVolume::AddNode(TGeoVolume *daughter, Int_t usernumber,
135 TGeoMatrix *matrix=gGeoIdentity)
136~~~
137
138Since we did not supplied the third argument, the wire will be positioned with
139an identity transformation inside the chamber. One will notice that the inner
140radii of the wire and chamber are both zero - therefore, aren't the two volumes
141overlapping ? The answer is no, the modeler is even relaying on the fact that
142any daughter is fully contained by its mother. On the other hand, neither of
143the nodes positioned inside a volume should overlap with each other. We will
144see that there are allowed some exceptions to those rules.
145
146### Overlapping volumes
147
148 Positioning volumes that does not overlap their neighbours nor extrude
149their container is sometimes quite strong constraint. Some parts of the geometry
150might overlap naturally, e.g. two crossing tubes. The modeller supports such
151cases only if the overlapping nodes are declared by the user. In order to do
152that, one should use TGeoVolume::AddNodeOverlap() instead of TGeoVolume::AddNode().
153 When 2 or more positioned volumes are overlapping, not all of them have to
154be declared so, but at least one. A point inside an overlapping region equally
155belongs to all overlapping nodes, but the way these are defined can enforce
156the modeler to give priorities.
157 The general rule is that the deepest node in the hierarchy containing a point
158have the highest priority. For the same geometry level, non-overlapping is
159prioritised over overlapping. In order to illustrate this, we will consider
160few examples. We will designate non-overlapping nodes as ONLY and the others
161MANY as in GEANT3, where this concept was introduced:
162 1. The part of a MANY node B extruding its container A will never be "seen"
163during navigation, as if B was in fact the result of the intersection of A and B.
164 2. If we have two nodes A (ONLY) and B (MANY) inside the same container, all
165points in the overlapping region of A and B will be designated as belonging to A.
166 3. If A an B in the above case were both MANY, points in the overlapping
167part will be designated to the one defined first. Both nodes must have the
168same medium.
169 4. The slices of a divided MANY will be as well MANY.
170
171One needs to know that navigation inside geometry parts MANY nodes is much
172slower. Any overlapping part can be defined based on composite shapes - this
173is always recommended.
174
175### Replicating volumes
176
177 What can we do if our chamber contains two identical wires instead of one ?
178What if then we would need 1000 chambers in our detector ? Should we create
1792000 wires and 1000 chamber volumes ? No, we will just need to replicate the
180ones that we have already created.
181
182~~~ {.cpp}
183 chamber->AddNode(wire_co, 1, new TGeoTranslation(-0.2,0,0));
184 chamber->AddNode(wire_co, 2, new TGeoTranslation(0.2,0,0));
185~~~
186
187 The 2 nodes that we have created inside chamber will both point to a wire_co
188object, but will be completely distinct : WIRE_CO_1 and WIRE_CO_2. We will
189want now to place symmetrically 1000 chambers on a pad, following a pattern
190of 20 rows and 50 columns. One way to do this will be to replicate our chamber
191by positioning it 1000 times in different positions of the pad. Unfortunately,
192this is far from being the optimal way of doing what we want.
193Imagine that we would like to find out which of the 1000 chambers is containing
194a (x,y,z) point defined in the pad reference. You will never have to do that,
195since the modeller will take care of it for you, but let's guess what it has
196to do. The most simple algorithm will just loop over all daughters, convert
197the point from mother to local reference and check if the current chamber
198contains the point or not. This might be efficient for pads with few chambers,
199but definitely not for 1000. Fortunately the modeler is smarter than that and
200create for each volume some optimization structures called voxels (see Voxelization)
201to minimize the penalty having too many daughters, but if you have 100 pads like
202this in your geometry you will anyway loose a lot in your tracking performance.
203
204 The way out when volumes can be arranged according to simple patterns is the
205usage of divisions. We will describe them in detail later on. Let's think now
206at a different situation : instead of 1000 chambers of the same type, we may
207have several types of chambers. Let's say all chambers are cylindrical and have
208a wire inside, but their dimensions are different. However, we would like all
209to be represented by a single volume family, since they have the same properties.
210*/
211
212/** \class TGeoVolumeMulti
213\ingroup Geometry_classes
214
215Volume families
216
217A volume family is represented by the class TGeoVolumeMulti. It represents
218a class of volumes having the same shape type and each member will be
219identified by the same name and volume ID. Any operation applied to a
220TGeoVolume equally affects all volumes in that family. The creation of a
221family is generally not a user task, but can be forced in particular cases:
222
223~~~ {.cpp}
224 TGeoManager::Volume(const char *vname, const char *shape, Int_t nmed);
225~~~
226
227where VNAME is the family name, NMED is the medium number and SHAPE is the
228shape type that can be:
229
230~~~ {.cpp}
231 box - for TGeoBBox
232 trd1 - for TGeoTrd1
233 trd2 - for TGeoTrd2
234 trap - for TGeoTrap
235 gtra - for TGeoGtra
236 para - for TGeoPara
237 tube, tubs - for TGeoTube, TGeoTubeSeg
238 cone, cons - for TGeoCone, TgeoCons
239 eltu - for TGeoEltu
240 ctub - for TGeoCtub
241 pcon - for TGeoPcon
242 pgon - for TGeoPgon
243~~~
244
245Volumes are then added to a given family upon adding the generic name as node
246inside other volume:
247
248~~~ {.cpp}
249 TGeoVolume *box_family = gGeoManager->Volume("BOXES", "box", nmed);
250 ...
251 gGeoManager->Node("BOXES", Int_t copy_no, "mother_name",
252 Double_t x, Double_t y, Double_t z, Int_t rot_index,
253 Bool_t is_only, Double_t *upar, Int_t npar);
254~~~
255
256here:
257
258~~~ {.cpp}
259 BOXES - name of the family of boxes
260 copy_no - user node number for the created node
261 mother_name - name of the volume to which we want to add the node
262 x,y,z - translation components
263 rot_index - indx of a rotation matrix in the list of matrices
264 upar - array of actual shape parameters
265 npar - number of parameters
266~~~
267
268The parameters order and number are the same as in the corresponding shape
269constructors.
270
271 Another particular case where volume families are used is when we want
272that a volume positioned inside a container to match one ore more container
273limits. Suppose we want to position the same box inside 2 different volumes
274and we want the Z size to match the one of each container:
275
276~~~ {.cpp}
277 TGeoVolume *container1 = gGeoManager->MakeBox("C1", imed, 10,10,30);
278 TGeoVolume *container2 = gGeoManager->MakeBox("C2", imed, 10,10,20);
279 TGeoVolume *pvol = gGeoManager->MakeBox("PVOL", jmed, 3,3,-1);
280 container1->AddNode(pvol, 1);
281 container2->AddNode(pvol, 1);
282~~~
283
284 Note that the third parameter of PVOL is negative, which does not make sense
285as half-length on Z. This is interpreted as: when positioned, create a box
286replacing all invalid parameters with the corresponding dimensions of the
287container. This is also internally handled by the TGeoVolumeMulti class, which
288does not need to be instantiated by users.
289
290### Dividing volumes
291
292 Volumes can be divided according a pattern. The most simple division can
293be done along one axis, that can be: X, Y, Z, Phi, Rxy or Rxyz. Let's take
294the most simple case: we would like to divide a box in N equal slices along X
295coordinate, representing a new volume family. Supposing we already have created
296the initial box, this can be done like:
297
298~~~ {.cpp}
299 TGeoVolume *slicex = box->Divide("SLICEX", 1, N);
300~~~
301
302where SLICE is the name of the new family representing all slices and 1 is the
303slicing axis. The meaning of the axis index is the following: for all volumes
304having shapes like box, trd1, trd2, trap, gtra or para - 1,2,3 means X,Y,Z; for
305tube, tubs, cone, cons - 1 means Rxy, 2 means phi and 3 means Z; for pcon and
306pgon - 2 means phi and 3 means Z; for spheres 1 means R and 2 means phi.
307 In fact, the division operation has the same effect as positioning volumes
308in a given order inside the divided container - the advantage being that the
309navigation in such a structure is much faster. When a volume is divided, a
310volume family corresponding to the slices is created. In case all slices can
311be represented by a single shape, only one volume is added to the family and
312positioned N times inside the divided volume, otherwise, each slice will be
313represented by a distinct volume in the family.
314 Divisions can be also performed in a given range of one axis. For that, one
315have to specify also the starting coordinate value and the step:
316
317~~~ {.cpp}
318 TGeoVolume *slicex = box->Divide("SLICEX", 1, N, start, step);
319~~~
320
321A check is always done on the resulting division range : if not fitting into
322the container limits, an error message is posted. If we will browse the divided
323volume we will notice that it will contain N nodes starting with index 1 upto
324N. The first one has the lower X limit at START position, while the last one
325will have the upper X limit at START+N*STEP. The resulting slices cannot
326be positioned inside an other volume (they are by default positioned inside the
327divided one) but can be further divided and may contain other volumes:
328
329~~~ {.cpp}
330 TGeoVolume *slicey = slicex->Divide("SLICEY", 2, N1);
331 slicey->AddNode(other_vol, index, some_matrix);
332~~~
333
334 When doing that, we have to remember that SLICEY represents a family, therefore
335all members of the family will be divided on Y and the other volume will be
336added as node inside all.
337 In the example above all the resulting slices had the same shape as the
338divided volume (box). This is not always the case. For instance, dividing a
339volume with TUBE shape on PHI axis will create equal slices having TUBESEG
340shape. Other divisions can also create slices having shapes with different
341dimensions, e.g. the division of a TRD1 volume on Z.
342 When positioning volumes inside slices, one can do it using the generic
343volume family (e.g. slicey). This should be done as if the coordinate system
344of the generic slice was the same as the one of the divided volume. The generic
345slice in case of PHI division is centered with respect to X axis. If the
346family contains slices of different sizes, any volume positioned inside should
347fit into the smallest one.
348 Examples for specific divisions according to shape types can be found inside
349shape classes.
350
351~~~ {.cpp}
352 TGeoVolume::Divide(N, Xmin, Xmax, "X");
353~~~
354
355 The GEANT3 option MANY is supported by TGeoVolumeOverlap class. An overlapping
356volume is in fact a virtual container that does not represent a physical object.
357It contains a list of nodes that are not its daughters but that must be checked
358always before the container itself. This list must be defined by users and it
359is checked and resolved in a priority order. Note that the feature is non-standard
360to geometrical modelers and it was introduced just to support conversions of
361GEANT3 geometries, therefore its extensive usage should be avoided.
362*/
363
364/** \class TGeoVolumeAssembly
365\ingroup Geometry_classes
366
367Volume assemblies
368
369Assemblies a volumes that have neither a shape or a material/medium. Assemblies
370behave exactly like normal volumes grouping several daughters together, but
371the daughters can never extrude the assembly since this has no shape. However,
372a bounding box and a voxelization structure are built for assemblies as for
373normal volumes, so that navigation is still optimized. Assemblies are useful
374for grouping hierarchically volumes which are otherwise defined in a flat
375manner, but also to avoid clashes between container shapes.
376To define an assembly one should just input a name, then start adding other
377volumes (or volume assemblies) as content.
378*/
379
380#include <fstream>
381#include <iomanip>
382
383#include "TString.h"
384#include "TBuffer.h"
385#include "TBrowser.h"
386#include "TStyle.h"
387#include "TH2F.h"
388#include "TROOT.h"
389#include "TEnv.h"
390#include "TMap.h"
391#include "TFile.h"
392#include "TKey.h"
393
394#include "TGeoManager.h"
395#include "TGeoNode.h"
396#include "TGeoMatrix.h"
397#include "TVirtualGeoPainter.h"
398#include "TVirtualGeoChecker.h"
399#include "TGeoVolume.h"
400#include "TGeoShapeAssembly.h"
401#include "TGeoScaledShape.h"
402#include "TGeoCompositeShape.h"
403#include "TGeoVoxelFinder.h"
404#include "TGeoExtension.h"
405
407
409
410////////////////////////////////////////////////////////////////////////////////
411/// Create a dummy medium
412
414{
415 if (fgDummyMedium)
416 return;
418 fgDummyMedium->SetName("dummy");
420 dummyMaterial->SetName("dummy");
421 fgDummyMedium->SetMaterial(dummyMaterial);
422}
423
424////////////////////////////////////////////////////////////////////////////////
425
427{
428 if (fFinder)
430 if (fShape)
432}
433
434////////////////////////////////////////////////////////////////////////////////
435
443
444////////////////////////////////////////////////////////////////////////////////
445
450
451////////////////////////////////////////////////////////////////////////////////
452/// dummy constructor
453
455{
456 fNodes = nullptr;
457 fShape = nullptr;
458 fMedium = nullptr;
459 fFinder = nullptr;
460 fVoxels = nullptr;
462 fField = nullptr;
463 fOption = "";
464 fNumber = 0;
465 fNtotal = 0;
466 fRefCount = 0;
467 fUserExtension = nullptr;
468 fFWExtension = nullptr;
469 fTransparency = -1;
471}
472
473////////////////////////////////////////////////////////////////////////////////
474/// default constructor
475
476TGeoVolume::TGeoVolume(const char *name, const TGeoShape *shape, const TGeoMedium *med) : TNamed(name, "")
477{
478 fName = fName.Strip();
479 fNodes = nullptr;
480 fShape = (TGeoShape *)shape;
481 if (fShape) {
483 Warning("Ctor", "volume %s has invalid shape", name);
484 }
485 if (!fShape->IsValid()) {
486 Fatal("ctor", "Shape of volume %s invalid. Aborting!", fName.Data());
487 }
488 }
490 if (fMedium && fMedium->GetMaterial())
492 fFinder = nullptr;
493 fVoxels = nullptr;
495 fField = nullptr;
496 fOption = "";
497 fNumber = 0;
498 fNtotal = 0;
499 fRefCount = 0;
500 fUserExtension = nullptr;
501 fFWExtension = nullptr;
502 fTransparency = -1;
503 if (fGeoManager)
506}
507
508////////////////////////////////////////////////////////////////////////////////
509/// Destructor
510
512{
513 if (fNodes) {
515 fNodes->Delete();
516 }
517 delete fNodes;
518 }
520 delete fFinder;
521 if (fVoxels)
522 delete fVoxels;
523 if (fUserExtension) {
525 fUserExtension = nullptr;
526 }
527 if (fFWExtension) {
529 fFWExtension = nullptr;
530 }
531}
532
533////////////////////////////////////////////////////////////////////////////////
534/// How to browse a volume
535
537{
538 if (!b)
539 return;
540
541 // if (!GetNdaughters()) b->Add(this, GetName(), IsVisible());
543 TString title;
544 for (Int_t i = 0; i < GetNdaughters(); i++) {
545 daughter = GetNode(i)->GetVolume();
546 if (daughter->GetTitle()[0]) {
547 if (daughter->IsAssembly())
548 title.TString::Format("Assembly with %d daughter(s)", daughter->GetNdaughters());
549 else if (daughter->GetFinder()) {
550 TString s1 = daughter->GetFinder()->ClassName();
551 s1.ReplaceAll("TGeoPattern", "");
552 title.TString::Format("Volume having %s shape divided in %d %s slices", daughter->GetShape()->ClassName(),
553 daughter->GetNdaughters(), s1.Data());
554
555 } else
556 title.TString::Format("Volume with %s shape having %d daughter(s)", daughter->GetShape()->ClassName(),
557 daughter->GetNdaughters());
558 daughter->SetTitle(title.Data());
559 }
560 b->Add(daughter, daughter->GetName(), daughter->IsVisible());
561 // if (IsVisDaughters())
562 // b->AddCheckBox(daughter, daughter->IsVisible());
563 // else
564 // b->AddCheckBox(daughter, kFALSE);
565 }
566}
567
568////////////////////////////////////////////////////////////////////////////////
569/// Computes the capacity of this [cm^3] as the capacity of its shape.
570/// In case of assemblies, the capacity is computed as the sum of daughter's capacities.
571
573{
574 if (!IsAssembly())
575 return fShape->Capacity();
576 Double_t capacity = 0.0;
577 Int_t nd = GetNdaughters();
578 Int_t i;
579 for (i = 0; i < nd; i++)
580 capacity += GetNode(i)->GetVolume()->Capacity();
581 return capacity;
582}
583
584////////////////////////////////////////////////////////////////////////////////
585/// Shoot nrays with random directions from starting point (startx, starty, startz)
586/// in the reference frame of this volume. Track each ray until exiting geometry, then
587/// shoot backwards from exiting point and compare boundary crossing points.
588
600
601////////////////////////////////////////////////////////////////////////////////
602/// Overlap checking tool. Check for illegal overlaps within a limit OVLP.
603/// Use option="s[number]" to force overlap checking by sampling volume with
604/// [number] points.
605///
606/// Ex:
607/// ~~~ {.cpp}
608/// myVol->CheckOverlaps(0.01, "s10000000"); // shoot 10000000 points
609/// myVol->CheckOverlaps(0.01, "s"); // shoot the default value of 1e6 points
610/// ~~~
611
613{
614 if (!GetNdaughters() || fFinder)
615 return;
617 TString opt(option);
618 opt.ToLower();
619 if (opt.Contains("s"))
620 sampling = kTRUE;
622 if (!sampling)
626 // Info("CheckOverlaps", "=== Checking overlaps for volume %s ===\n", GetName());
627 }
628 checker->CheckOverlaps(this, ovlp, option);
629 // if (sampling) return;
633 Int_t novlps = overlaps->GetEntriesFast();
634 TNamed *obj;
636 for (Int_t i = 0; i < novlps; i++) {
637 obj = (TNamed *)overlaps->At(i);
638 if (novlps < 1000)
639 name = TString::Format("ov%03d", i);
640 else
641 name = TString::Format("ov%06d", i);
642 obj->SetName(name);
643 }
644 if (novlps)
645 Info("CheckOverlaps", "Number of illegal overlaps/extrusions for volume %s: %d\n", GetName(), novlps);
646 }
647}
648
649////////////////////////////////////////////////////////////////////////////////
650/// Tests for checking the shape navigation algorithms. See TGeoShape::CheckShape()
651
656
657////////////////////////////////////////////////////////////////////////////////
658/// Clean data of the volume.
659
661{
662 ClearNodes();
663 ClearShape();
664}
665
666////////////////////////////////////////////////////////////////////////////////
667/// Clear the shape of this volume from the list held by the current manager.
668
673
674////////////////////////////////////////////////////////////////////////////////
675/// check for negative parameters in shapes.
676
678{
679 if (fShape->IsRunTimeShape()) {
680 Error("CheckShapes", "volume %s has run-time shape", GetName());
681 InspectShape();
682 return;
683 }
684 if (!fNodes)
685 return;
687 TGeoNode *node = nullptr;
689 const TGeoShape *shape = nullptr;
691 for (Int_t i = 0; i < nd; i++) {
692 node = (TGeoNode *)fNodes->At(i);
693 // check if node has name
694 if (!node->GetName()[0])
695 printf("Daughter %i of volume %s - NO NAME!!!\n", i, GetName());
696 old_vol = node->GetVolume();
697 shape = old_vol->GetShape();
698 if (shape->IsRunTimeShape()) {
699 // printf(" Node %s/%s has shape with negative parameters. \n",
700 // GetName(), node->GetName());
701 // old_vol->InspectShape();
702 // make a copy of the node
703 new_node = node->MakeCopyNode();
704 if (!new_node) {
705 Fatal("CheckShapes", "Cannot make copy node for %s", node->GetName());
706 return;
707 }
709 if (!new_shape) {
710 Error("CheckShapes", "cannot resolve runtime shape for volume %s/%s\n", GetName(), old_vol->GetName());
711 continue;
712 }
713 TGeoVolume *new_volume = old_vol->MakeCopyVolume(new_shape);
714 // printf(" new volume %s shape params :\n", new_volume->GetName());
715 // new_volume->InspectShape();
716 new_node->SetVolume(new_volume);
717 // decouple the old node and put the new one instead
718 fNodes->AddAt(new_node, i);
719 // new_volume->CheckShapes();
720 }
721 }
722}
723
724////////////////////////////////////////////////////////////////////////////////
725/// Count total number of subnodes starting from this volume, nlevels down
726/// - option = 0 (default) - count only once per volume
727/// - option = 1 - count every time
728/// - option = 2 - count volumes on visible branches
729/// - option = 3 - return maximum level counted already with option = 0
730
732{
733 static Int_t maxlevel = 0;
734 static Int_t nlev = 0;
735
737 option = 0;
738 Int_t visopt = 0;
739 Int_t nd = GetNdaughters();
740 Bool_t last = (!nlevels || !nd) ? kTRUE : kFALSE;
741 switch (option) {
742 case 0:
743 if (fNtotal)
744 return fNtotal;
745 case 1: fNtotal = 1; break;
746 case 2:
748 if (!IsVisDaughters())
749 last = kTRUE;
750 switch (visopt) {
751 case TVirtualGeoPainter::kGeoVisDefault: fNtotal = (IsVisible()) ? 1 : 0; break;
752 case TVirtualGeoPainter::kGeoVisLeaves: fNtotal = (IsVisible() && last) ? 1 : 0;
753 }
754 if (!IsVisibleDaughters())
755 return fNtotal;
756 break;
757 case 3: return maxlevel;
758 }
759 if (last)
760 return fNtotal;
761 if (gGeoManager->GetTopVolume() == this) {
762 maxlevel = 0;
763 nlev = 0;
764 }
765 if (nlev > maxlevel)
766 maxlevel = nlev;
767 TGeoNode *node;
768 TGeoVolume *vol;
769 nlev++;
770 for (Int_t i = 0; i < nd; i++) {
771 node = GetNode(i);
772 vol = node->GetVolume();
773 fNtotal += vol->CountNodes(nlevels - 1, option);
774 }
775 nlev--;
776 return fNtotal;
777}
778
779////////////////////////////////////////////////////////////////////////////////
780/// Return TRUE if volume and all daughters are invisible.
781
783{
784 if (IsVisible())
785 return kFALSE;
786 Int_t nd = GetNdaughters();
787 for (Int_t i = 0; i < nd; i++)
788 if (GetNode(i)->GetVolume()->IsVisible())
789 return kFALSE;
790 return kTRUE;
791}
792
793////////////////////////////////////////////////////////////////////////////////
794/// Make volume and each of it daughters (in)visible.
795
797{
799 Int_t nd = GetNdaughters();
800 TObjArray *list = new TObjArray(nd + 1);
801 list->Add(this);
802 TGeoVolume *vol;
803 for (Int_t i = 0; i < nd; i++) {
804 vol = GetNode(i)->GetVolume();
805 vol->SetAttVisibility(!flag);
806 list->Add(vol);
807 }
808 TIter next(gROOT->GetListOfBrowsers());
809 TBrowser *browser = nullptr;
810 while ((browser = (TBrowser *)next())) {
811 for (Int_t i = 0; i < nd + 1; i++) {
812 vol = (TGeoVolume *)list->At(i);
813 browser->CheckObjectItem(vol, !flag);
814 }
815 browser->Refresh();
816 }
817 delete list;
819}
820
821////////////////////////////////////////////////////////////////////////////////
822/// Return TRUE if volume contains nodes
823
825{
826 return kTRUE;
827}
828
829////////////////////////////////////////////////////////////////////////////////
830/// check if the visibility and attributes are the default ones
831
833{
834 if (!IsVisible())
835 return kFALSE;
836 if (GetLineColor() != gStyle->GetLineColor())
837 return kFALSE;
838 if (GetLineStyle() != gStyle->GetLineStyle())
839 return kFALSE;
840 if (GetLineWidth() != gStyle->GetLineWidth())
841 return kFALSE;
842 return kTRUE;
843}
844
845////////////////////////////////////////////////////////////////////////////////
846/// True if this is the top volume of the geometry
847
849{
850 if (fGeoManager->GetTopVolume() == this)
851 return kTRUE;
852 return kFALSE;
853}
854
855////////////////////////////////////////////////////////////////////////////////
856/// Check if the painter is currently ray-tracing the content of this volume.
857
862
863////////////////////////////////////////////////////////////////////////////////
864/// Inspect the material for this volume.
865
867{
868 GetMaterial()->Print();
869}
870
871////////////////////////////////////////////////////////////////////////////////
872/// Import a volume from a file.
873
874TGeoVolume *TGeoVolume::Import(const char *filename, const char *name, Option_t * /*option*/)
875{
876 if (!gGeoManager)
877 gGeoManager = new TGeoManager("geometry", "");
878 if (!filename)
879 return nullptr;
880 TGeoVolume *volume = nullptr;
881 if (strstr(filename, ".gdml")) {
882 // import from a gdml file
883 } else {
884 // import from a root file
887 if (!f || f->IsZombie()) {
888 printf("Error: TGeoVolume::Import : Cannot open file %s\n", filename);
889 return nullptr;
890 }
891 if (name && name[0]) {
892 volume = (TGeoVolume *)f->Get(name);
893 } else {
894 TIter next(f->GetListOfKeys());
895 TKey *key;
896 while ((key = (TKey *)next())) {
897 if (strcmp(key->GetClassName(), "TGeoVolume") != 0)
898 continue;
899 volume = (TGeoVolume *)key->ReadObj();
900 break;
901 }
902 }
903 delete f;
904 }
905 if (!volume)
906 return nullptr;
907 volume->RegisterYourself();
908 return volume;
909}
910
911////////////////////////////////////////////////////////////////////////////////
912/// Export this volume to a file.
913///
914/// - Case 1: root file or root/xml file
915/// if filename end with ".root". The key will be named name
916/// if filename end with ".xml" a root/xml file is produced.
917///
918/// - Case 2: C++ script
919/// if filename end with ".C"
920///
921/// - Case 3: gdml file
922/// if filename end with ".gdml"
923///
924/// NOTE that to use this option, the PYTHONPATH must be defined like
925/// export PYTHONPATH=$ROOTSYS/lib:$ROOTSYS/gdml
926///
927
929{
931 if (sfile.Contains(".C")) {
932 // Save volume as a C++ script
933 Info("Export", "Exporting volume %s as C++ code", GetName());
934 SaveAs(filename, "");
935 return 1;
936 }
937 if (sfile.Contains(".gdml")) {
938 // Save geometry as a gdml file
939 Info("Export", "Exporting %s as gdml code - not implemented yet", GetName());
940 return 0;
941 }
942 if (sfile.Contains(".root") || sfile.Contains(".xml")) {
943 // Save volume in a root file
944 Info("Export", "Exporting %s as root file.", GetName());
945 TString opt(option);
946 if (!opt.Length())
947 opt = "recreate";
948 TFile *f = TFile::Open(filename, opt.Data());
949 if (!f || f->IsZombie()) {
950 Error("Export", "Cannot open file");
951 return 0;
952 }
954 if (keyname.IsNull())
955 keyname = GetName();
957 delete f;
958 return nbytes;
959 }
960 return 0;
961}
962
963////////////////////////////////////////////////////////////////////////////////
964/// Actualize matrix of node indexed `<inode>`
965
967{
968 if (fFinder)
970}
971
972////////////////////////////////////////////////////////////////////////////////
973/// Add a TGeoNode to the list of nodes. This is the usual method for adding
974/// daughters inside the container volume.
975
977{
979 if (matrix == nullptr)
981 else
982 matrix->RegisterYourself();
983 if (!vol) {
984 Error("AddNode", "Volume is NULL");
985 return nullptr;
986 }
987 if (!vol->IsValid()) {
988 Error("AddNode", "Won't add node with invalid shape");
989 printf("### invalid volume was : %s\n", vol->GetName());
990 return nullptr;
991 }
992 if (!fNodes)
993 fNodes = new TObjArray();
994
995 if (fFinder) {
996 // volume already divided.
997 Error("AddNode", "Cannot add node %s_%i into divided volume %s", vol->GetName(), copy_no, GetName());
998 return nullptr;
999 }
1000
1001 TGeoNodeMatrix *node = nullptr;
1002 node = new TGeoNodeMatrix(vol, matrix);
1003 node->SetMotherVolume(this);
1004 fNodes->Add(node);
1005 TString name = TString::Format("%s_%d", vol->GetName(), copy_no);
1006 // if (fNodes->FindObject(name))
1007 // Warning("AddNode", "Volume %s : added node %s with same name", GetName(), name.Data());
1008 node->SetName(name);
1009 node->SetNumber(copy_no);
1010 fRefCount++;
1011 vol->Grab();
1012 return node;
1013}
1014
1015////////////////////////////////////////////////////////////////////////////////
1016/// Add a division node to the list of nodes. The method is called by
1017/// TGeoVolume::Divide() for creating the division nodes.
1018
1020{
1021 if (!vol) {
1022 Error("AddNodeOffset", "invalid volume");
1023 return;
1024 }
1025 if (!vol->IsValid()) {
1026 Error("AddNode", "Won't add node with invalid shape");
1027 printf("### invalid volume was : %s\n", vol->GetName());
1028 return;
1029 }
1030 if (!fNodes)
1031 fNodes = new TObjArray();
1032 TGeoNode *node = new TGeoNodeOffset(vol, copy_no, offset);
1033 node->SetMotherVolume(this);
1034 fNodes->Add(node);
1035 TString name = TString::Format("%s_%d", vol->GetName(), copy_no + 1);
1036 node->SetName(name);
1037 node->SetNumber(copy_no + 1);
1038 vol->Grab();
1039}
1040
1041////////////////////////////////////////////////////////////////////////////////
1042/// Add a TGeoNode to the list of nodes. This is the usual method for adding
1043/// daughters inside the container volume.
1044
1046{
1047 if (!vol) {
1048 Error("AddNodeOverlap", "Volume is NULL");
1049 return;
1050 }
1051 if (!vol->IsValid()) {
1052 Error("AddNodeOverlap", "Won't add node with invalid shape");
1053 printf("### invalid volume was : %s\n", vol->GetName());
1054 return;
1055 }
1056 if (vol->IsAssembly()) {
1057 Warning("AddNodeOverlap",
1058 "Declaring assembly %s as possibly overlapping inside %s not allowed. Using AddNode instead !",
1059 vol->GetName(), GetName());
1060 AddNode(vol, copy_no, mat, option);
1061 return;
1062 }
1064 if (matrix == nullptr)
1066 else
1067 matrix->RegisterYourself();
1068 if (!fNodes)
1069 fNodes = new TObjArray();
1070
1071 if (fFinder) {
1072 // volume already divided.
1073 Error("AddNodeOverlap", "Cannot add node %s_%i into divided volume %s", vol->GetName(), copy_no, GetName());
1074 return;
1075 }
1076
1077 TGeoNodeMatrix *node = new TGeoNodeMatrix(vol, matrix);
1078 node->SetMotherVolume(this);
1079 fNodes->Add(node);
1080 TString name = TString::Format("%s_%d", vol->GetName(), copy_no);
1081 if (fNodes->FindObject(name))
1082 Warning("AddNode", "Volume %s : added node %s with same name", GetName(), name.Data());
1083 node->SetName(name);
1084 node->SetNumber(copy_no);
1085 node->SetOverlapping();
1086 if (vol->GetMedium() == fMedium)
1087 node->SetVirtual();
1088 vol->Grab();
1089}
1090
1091////////////////////////////////////////////////////////////////////////////////
1092/// Division a la G3. The volume will be divided along IAXIS (see shape classes), in NDIV
1093/// slices, from START with given STEP. The division volumes will have medium number NUMED.
1094/// If NUMED=0 they will get the medium number of the divided volume (this). If NDIV<=0,
1095/// all range of IAXIS will be divided and the resulting number of divisions will be centered on
1096/// IAXIS. If STEP<=0, the real STEP will be computed as the full range of IAXIS divided by NDIV.
1097/// Options (case insensitive):
1098/// - N - divide all range in NDIV cells (same effect as STEP<=0) (GSDVN in G3)
1099/// - NX - divide range starting with START in NDIV cells (GSDVN2 in G3)
1100/// - S - divide all range with given STEP. NDIV is computed and divisions will be centered
1101/// in full range (same effect as NDIV<=0) (GSDVS, GSDVT in G3)
1102/// - SX - same as DVS, but from START position. (GSDVS2, GSDVT2 in G3)
1103
1106{
1107 if (fFinder) {
1108 // volume already divided.
1109 Fatal("Divide", "volume %s already divided", GetName());
1110 return nullptr;
1111 }
1112 TString opt(option);
1113 opt.ToLower();
1114 TString stype = fShape->ClassName();
1115 if (!fNodes)
1116 fNodes = new TObjArray();
1117 Double_t xlo, xhi, range;
1118 range = fShape->GetAxisRange(iaxis, xlo, xhi);
1119 // for phi divisions correct the range
1120 if (!strcmp(fShape->GetAxisName(iaxis), "PHI")) {
1121 if ((start - xlo) < -1E-3)
1122 start += 360.;
1124 xlo = start;
1125 xhi = start + range;
1126 }
1127 }
1128 if (range <= 0) {
1129 InspectShape();
1130 Fatal("Divide", "cannot divide volume %s (%s) on %s axis", GetName(), stype.Data(), fShape->GetAxisName(iaxis));
1131 return nullptr;
1132 }
1133 if (ndiv <= 0 || opt.Contains("s")) {
1134 if (step <= 0) {
1135 Fatal("Divide", "invalid division type for volume %s : ndiv=%i, step=%g", GetName(), ndiv, step);
1136 return nullptr;
1137 }
1138 if (opt.Contains("x")) {
1139 if ((xlo - start) > 1E-3 || (xhi - start) < -1E-3) {
1140 Fatal("Divide", "invalid START=%g for division on axis %s of volume %s. Range is (%g, %g)", start,
1141 fShape->GetAxisName(iaxis), GetName(), xlo, xhi);
1142 return nullptr;
1143 }
1144 xlo = start;
1145 range = xhi - xlo;
1146 }
1147 ndiv = Int_t((range + 0.1 * step) / step);
1148 Double_t ddx = range - ndiv * step;
1149 // always center the division in this case
1150 if (ddx > 1E-3)
1151 Warning("Divide", "division of volume %s on %s axis (ndiv=%d) will be centered in the full range", GetName(),
1152 fShape->GetAxisName(iaxis), ndiv);
1153 start = xlo + 0.5 * ddx;
1154 }
1155 if (step <= 0 || opt.Contains("n")) {
1156 if (opt.Contains("x")) {
1157 if ((xlo - start) > 1E-3 || (xhi - start) < -1E-3) {
1158 Fatal("Divide", "invalid START=%g for division on axis %s of volume %s. Range is (%g, %g)", start,
1159 fShape->GetAxisName(iaxis), GetName(), xlo, xhi);
1160 return nullptr;
1161 }
1162 xlo = start;
1163 range = xhi - xlo;
1164 }
1165 step = range / ndiv;
1166 start = xlo;
1167 }
1168
1169 Double_t end = start + ndiv * step;
1170 if (((start - xlo) < -1E-3) || ((end - xhi) > 1E-3)) {
1171 Fatal("Divide", "division of volume %s on axis %s exceed range (%g, %g)", GetName(), fShape->GetAxisName(iaxis),
1172 xlo, xhi);
1173 return nullptr;
1174 }
1175 TGeoVolume *voldiv = fShape->Divide(this, divname, iaxis, ndiv, start, step);
1176 if (numed) {
1178 if (!medium) {
1179 Fatal("Divide", "invalid medium number %d for division volume %s", numed, divname);
1180 return voldiv;
1181 }
1182 voldiv->SetMedium(medium);
1183 if (medium->GetMaterial())
1184 medium->GetMaterial()->SetUsed();
1185 }
1186 return voldiv;
1187}
1188
1189////////////////////////////////////////////////////////////////////////////////
1190/// compute the closest distance of approach from point px,py to this volume
1191
1193{
1194 if (gGeoManager != fGeoManager)
1197 Int_t dist = 9999;
1198 if (!painter)
1199 return dist;
1200 dist = painter->DistanceToPrimitiveVol(this, px, py);
1201 return dist;
1202}
1203
1204////////////////////////////////////////////////////////////////////////////////
1205/// draw top volume according to option
1206
1208{
1209 if (gGeoManager != fGeoManager)
1214 if (!IsVisContainers())
1215 SetVisLeaves();
1216 if (option && option[0] > 0) {
1217 painter->DrawVolume(this, option);
1218 } else {
1219 painter->DrawVolume(this, gEnv->GetValue("Viewer3D.DefaultDrawOption", ""));
1220 }
1221}
1222
1223////////////////////////////////////////////////////////////////////////////////
1224/// draw only this volume
1225
1227{
1228 if (IsAssembly()) {
1229 Info("DrawOnly", "Volume assemblies do not support this option.");
1230 return;
1231 }
1232 if (gGeoManager != fGeoManager)
1234 SetVisOnly();
1237 if (option && option[0] > 0) {
1238 painter->DrawVolume(this, option);
1239 } else {
1240 painter->DrawVolume(this, gEnv->GetValue("Viewer3D.DefaultDrawOption", ""));
1241 }
1242}
1243
1244////////////////////////////////////////////////////////////////////////////////
1245/// Perform an extensive sampling to find which type of voxelization is
1246/// most efficient.
1247
1249{
1250 printf("Optimizing volume %s ...\n", GetName());
1252 return checker->TestVoxels(this);
1253}
1254
1255////////////////////////////////////////////////////////////////////////////////
1256/// Print volume info
1257
1259{
1260 printf("== Volume: %s type %s positioned %d times\n", GetName(), ClassName(), fRefCount);
1261 InspectShape();
1263}
1264
1265////////////////////////////////////////////////////////////////////////////////
1266/// paint volume
1267
1269{
1271 painter->SetTopVolume(this);
1272 // painter->Paint(option);
1273 if (option && option[0] > 0) {
1274 painter->Paint(option);
1275 } else {
1276 painter->Paint(gEnv->GetValue("Viewer3D.DefaultDrawOption", ""));
1277 }
1278}
1279
1280////////////////////////////////////////////////////////////////////////////////
1281/// Print the voxels for this volume.
1282
1284{
1285 if (fVoxels)
1286 fVoxels->Print();
1287}
1288
1289////////////////////////////////////////////////////////////////////////////////
1290/// Recreate the content of the other volume without pointer copying. Voxels are
1291/// ignored and supposed to be created in a later step via Voxelize.
1292
1294{
1295 Int_t nd = other->GetNdaughters();
1296 if (!nd)
1297 return;
1298 TGeoPatternFinder *finder = other->GetFinder();
1299 if (finder) {
1300 Int_t iaxis = finder->GetDivAxis();
1301 Int_t ndiv = finder->GetNdiv();
1302 Double_t start = finder->GetStart();
1303 Double_t step = finder->GetStep();
1304 Int_t numed = other->GetNode(0)->GetVolume()->GetMedium()->GetId();
1305 TGeoVolume *voldiv = Divide(other->GetNode(0)->GetVolume()->GetName(), iaxis, ndiv, start, step, numed);
1306 voldiv->ReplayCreation(other->GetNode(0)->GetVolume());
1307 return;
1308 }
1309 for (Int_t i = 0; i < nd; i++) {
1310 TGeoNode *node = other->GetNode(i);
1311 if (node->IsOverlapping())
1312 AddNodeOverlap(node->GetVolume(), node->GetNumber(), node->GetMatrix());
1313 else
1314 AddNode(node->GetVolume(), node->GetNumber(), node->GetMatrix());
1315 }
1316}
1317
1318////////////////////////////////////////////////////////////////////////////////
1319/// print nodes
1320
1322{
1323 Int_t nd = GetNdaughters();
1324 for (Int_t i = 0; i < nd; i++) {
1325 printf("%s\n", GetNode(i)->GetName());
1326 cd(i);
1327 GetNode(i)->GetMatrix()->Print();
1328 }
1329}
1330////////////////////////////////////////////////////////////////////////////////
1331/// Generate a lego plot fot the top volume, according to option.
1332
1335{
1338 if (old_vol != this)
1340 else
1341 old_vol = nullptr;
1342 TH2F *hist = checker->LegoPlot(ntheta, themin, themax, nphi, phimin, phimax, rmin, rmax, option);
1343 hist->Draw("lego1sph");
1344 return hist;
1345}
1346
1347////////////////////////////////////////////////////////////////////////////////
1348/// Register the volume and all materials/media/matrices/shapes to the manager.
1349
1351{
1353 return;
1354 // Register volume
1355 fGeoManager->AddVolume(this);
1356 // Register shape
1358 if (fShape->IsComposite()) {
1360 comp->RegisterYourself();
1361 } else {
1363 }
1364 }
1365 // Register medium/material
1370 }
1371 // Register matrices for nodes.
1373 TGeoNode *node;
1374 Int_t nd = GetNdaughters();
1375 Int_t i;
1376 for (i = 0; i < nd; i++) {
1377 node = GetNode(i);
1378 matrix = node->GetMatrix();
1379 if (!matrix->IsRegistered())
1380 matrix->RegisterYourself();
1383 }
1384 }
1385 // Call RegisterYourself recursively
1386 for (i = 0; i < nd; i++)
1388}
1389
1390////////////////////////////////////////////////////////////////////////////////
1391/// Draw random points in the bounding box of this volume.
1392
1394{
1395 if (gGeoManager != fGeoManager)
1398 if (old_vol != this)
1400 else
1401 old_vol = nullptr;
1403 if (old_vol)
1405}
1406
1407////////////////////////////////////////////////////////////////////////////////
1408/// Random raytracing method.
1409
1424
1425////////////////////////////////////////////////////////////////////////////////
1426/// Draw this volume with current settings and perform raytracing in the pad.
1427
1429{
1431 if (gGeoManager != fGeoManager)
1434 Bool_t drawn = (painter->GetDrawnVolume() == this) ? kTRUE : kFALSE;
1435 if (!drawn) {
1436 painter->DrawVolume(this, "");
1438 painter->ModifiedPad();
1439 return;
1440 }
1442 painter->ModifiedPad();
1443}
1444
1445////////////////////////////////////////////////////////////////////////////////
1446/// Save geometry having this as top volume as a C++ macro.
1447
1449{
1450 if (!filename)
1451 return;
1452 std::ofstream out;
1453 out.open(filename, std::ios::out);
1454 if (out.bad()) {
1455 Error("SavePrimitive", "Bad file name: %s", filename);
1456 return;
1457 }
1458 if (fGeoManager->GetTopVolume() != this)
1460
1462 Int_t ind = fname.Index(".");
1463 if (ind > 0)
1464 fname.Remove(ind);
1465 out << "void " << fname << "() {" << std::endl;
1466 out << " gSystem->Load(\"libGeom\");" << std::endl;
1468 out << std::setprecision(prec);
1469 ((TGeoVolume *)this)->SavePrimitive(out, option);
1470 out << "}" << std::endl;
1471}
1472
1473////////////////////////////////////////////////////////////////////////////////
1474/// Connect user-defined extension to the volume. The volume "grabs" a copy, so
1475/// the original object can be released by the producer. Release the previously
1476/// connected extension if any.
1477///
1478/// NOTE: This interface is intended for user extensions and is guaranteed not
1479/// to be used by TGeo
1480
1482{
1484 fUserExtension = nullptr;
1485 if (ext)
1486 fUserExtension = ext->Grab();
1487 if (tmp)
1488 tmp->Release();
1489}
1490
1491////////////////////////////////////////////////////////////////////////////////
1492/// Connect framework defined extension to the volume. The volume "grabs" a copy,
1493/// so the original object can be released by the producer. Release the previously
1494/// connected extension if any.
1495///
1496/// NOTE: This interface is intended for the use by TGeo and the users should
1497/// NOT connect extensions using this method
1498
1500{
1502 fFWExtension = nullptr;
1503 if (ext)
1504 fFWExtension = ext->Grab();
1505 if (tmp)
1506 tmp->Release();
1507}
1508
1509////////////////////////////////////////////////////////////////////////////////
1510/// Get a copy of the user extension pointer. The user must call Release() on
1511/// the copy pointer once this pointer is not needed anymore (equivalent to
1512/// delete() after calling new())
1513
1515{
1516 if (fUserExtension)
1517 return fUserExtension->Grab();
1518 return nullptr;
1519}
1520
1521////////////////////////////////////////////////////////////////////////////////
1522/// Get a copy of the framework extension pointer. The user must call Release() on
1523/// the copy pointer once this pointer is not needed anymore (equivalent to
1524/// delete() after calling new())
1525
1527{
1528 if (fFWExtension)
1529 return fFWExtension->Grab();
1530 return nullptr;
1531}
1532
1533////////////////////////////////////////////////////////////////////////////////
1534/// Save a primitive as a C++ statement(s) on output stream "out".
1535
1536void TGeoVolume::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
1537{
1538 Int_t i, icopy;
1539 Int_t nd = GetNdaughters();
1541 TGeoNode *dnode;
1543
1544 // check if we need to save shape/volume
1546 if (fGeoManager->GetGeomPainter()->GetTopVolume() == this)
1547 mustDraw = kTRUE;
1548 if (!option[0]) {
1550 out << " new TGeoManager(\"" << fGeoManager->GetName() << "\", \"" << fGeoManager->GetTitle() << "\");"
1551 << std::endl
1552 << std::endl;
1553 // if (mustDraw) out << " Bool_t mustDraw = kTRUE;" << std::endl;
1554 // else out << " Bool_t mustDraw = kFALSE;" << std::endl;
1555 out << " Double_t dx, dy, dz;" << std::endl;
1556 out << " Double_t dx1, dx2, dy1, dy2;" << std::endl;
1557 out << " Double_t vert[20], par[20];" << std::endl;
1558 out << " Double_t theta, phi, h1, bl1, tl1, alpha1, h2, bl2, tl2, alpha2;" << std::endl;
1559 out << " Double_t twist;" << std::endl;
1560 out << " Double_t origin[3];" << std::endl;
1561 out << " Double_t rmin, rmax, rmin1, rmax1, rmin2, rmax2;" << std::endl;
1562 out << " Double_t r, rlo, rhi;" << std::endl;
1563 out << " Double_t a, b;" << std::endl;
1564 out << " Double_t point[3], norm[3];" << std::endl;
1565 out << " Double_t rin, stin, rout, stout;" << std::endl;
1566 out << " Double_t thx, phx, thy, phy, thz, phz;" << std::endl;
1567 out << " Double_t alpha, theta1, theta2, phi1, phi2, dphi;" << std::endl;
1568 out << " Double_t tr[3], rot[9];" << std::endl;
1569 out << " Double_t z, density, radl, absl, w;" << std::endl;
1570 out << " Double_t lx, ly, lz, tx, ty, tz;" << std::endl;
1571 out << " Double_t xvert[50], yvert[50];" << std::endl;
1572 out << " Double_t zsect, x0, y0, scale0;" << std::endl;
1573 out << " Int_t nel, numed, nz, nedges, nvert;" << std::endl;
1574 out << " TGeoBoolNode *pBoolNode = nullptr;" << std::endl << std::endl;
1575 // first save materials/media
1576 out << " // MATERIALS, MIXTURES AND TRACKING MEDIA" << std::endl;
1577 SavePrimitive(out, "m");
1578 // then, save matrices
1579 out << std::endl << " // TRANSFORMATION MATRICES" << std::endl;
1580 SavePrimitive(out, "x");
1581 // save this volume and shape
1582 SavePrimitive(out, "s");
1583 out << std::endl << " // SET TOP VOLUME OF GEOMETRY" << std::endl;
1584 out << " gGeoManager->SetTopVolume(" << GetPointerName() << ");" << std::endl;
1585 // save daughters
1586 out << std::endl << " // SHAPES, VOLUMES AND GEOMETRICAL HIERARCHY" << std::endl;
1587 SavePrimitive(out, "d");
1588 out << std::endl << " // CLOSE GEOMETRY" << std::endl;
1589 out << " gGeoManager->CloseGeometry();" << std::endl;
1590 if (mustDraw) {
1591 if (!IsRaytracing())
1592 out << " gGeoManager->GetTopVolume()->Draw();" << std::endl;
1593 else
1594 out << " gGeoManager->GetTopVolume()->Raytrace();" << std::endl;
1595 }
1596 return;
1597 }
1598 // check if we need to save shape/volume
1599 if (!strcmp(option, "s")) {
1600 // create the shape for this volume
1602 return;
1603 if (!IsAssembly()) {
1605 out << " // Volume: " << GetName() << std::endl;
1606 out << " TGeoVolume *" << GetPointerName() << " = new TGeoVolume(\"" << GetName() << "\","
1607 << fShape->GetPointerName();
1608 if (fMedium)
1609 out << ", " << fMedium->GetPointerName();
1610 out << ");" << std::endl;
1611 } else {
1612 out << " // Assembly: " << GetName() << std::endl;
1613 out << " " << GetPointerName() << " = new TGeoVolumeAssembly(\"" << GetName() << "\""
1614 << ");" << std::endl;
1615 }
1616 SaveLineAttributes(out, GetPointerName(), 1, 1, 1);
1617 if (!IsVisible() && !IsAssembly())
1618 out << " " << GetPointerName() << "->SetVisibility(kFALSE);" << std::endl;
1619 if (!IsVisibleDaughters())
1620 out << " " << GetPointerName() << "->VisibleDaughters(kFALSE);" << std::endl;
1621 if (IsVisContainers())
1622 out << " " << GetPointerName() << "->SetVisContainers(kTRUE);" << std::endl;
1623 if (IsVisLeaves())
1624 out << " " << GetPointerName() << "->SetVisLeaves(kTRUE);" << std::endl;
1626 }
1627 // check if we need to save the media
1628 if (!strcmp(option, "m")) {
1629 if (fMedium)
1631 for (i = 0; i < nd; i++) {
1632 dvol = GetNode(i)->GetVolume();
1633 dvol->SavePrimitive(out, option);
1634 }
1635 return;
1636 }
1637 // check if we need to save the matrices
1638 if (!strcmp(option, "x")) {
1639 if (fFinder) {
1640 dvol = GetNode(0)->GetVolume();
1641 dvol->SavePrimitive(out, option);
1642 return;
1643 }
1644 for (i = 0; i < nd; i++) {
1645 dnode = GetNode(i);
1646 matrix = dnode->GetMatrix();
1647 if (!matrix->IsIdentity())
1648 matrix->SavePrimitive(out, option);
1649 dnode->GetVolume()->SavePrimitive(out, option);
1650 }
1651 return;
1652 }
1653 // check if we need to save volume daughters
1654 if (!strcmp(option, "d")) {
1655 if (!nd)
1656 return;
1658 return;
1660 if (fFinder) {
1661 // volume divided: generate volume->Divide()
1662 dnode = GetNode(0);
1663 dvol = dnode->GetVolume();
1664 out << " TGeoVolume *" << dvol->GetPointerName() << " = ";
1665 out << GetPointerName() << "->Divide(\"" << dvol->GetName() << "\", ";
1667 if (fMedium != dvol->GetMedium())
1668 out << ", " << dvol->GetMedium()->GetId();
1669 out << ");" << std::endl;
1670 dvol->SavePrimitive(out, "d");
1671 return;
1672 }
1673 for (i = 0; i < nd; i++) {
1674 dnode = GetNode(i);
1675 dvol = dnode->GetVolume();
1676 dvol->SavePrimitive(out, "s");
1677 matrix = dnode->GetMatrix();
1678 icopy = dnode->GetNumber();
1679 // generate AddNode()
1680 out << " " << GetPointerName() << "->AddNode";
1681 if (dnode->IsOverlapping())
1682 out << "Overlap";
1683 out << "(" << dvol->GetPointerName() << ", " << icopy;
1684 if (!matrix->IsIdentity())
1685 out << ", " << matrix->GetPointerName();
1686 out << ");" << std::endl;
1687 }
1688 // Recursive loop to daughters
1689 for (i = 0; i < nd; i++) {
1690 dnode = GetNode(i);
1691 dvol = dnode->GetVolume();
1692 dvol->SavePrimitive(out, "d");
1693 }
1694 }
1695}
1696
1697////////////////////////////////////////////////////////////////////////////////
1698/// Reset SavePrimitive bits.
1699
1707
1708////////////////////////////////////////////////////////////////////////////////
1709/// Execute mouse actions on this volume.
1710
1712{
1714 if (!painter)
1715 return;
1716 painter->ExecuteVolumeEvent(this, event, px, py);
1717}
1718
1719////////////////////////////////////////////////////////////////////////////////
1720/// search a daughter inside the list of nodes
1721
1723{
1724 return ((TGeoNode *)fNodes->FindObject(name));
1725}
1726
1727////////////////////////////////////////////////////////////////////////////////
1728/// Get the index of a daughter within check_list by providing the node pointer.
1729
1731{
1732 TGeoNode *current = nullptr;
1733 for (Int_t i = 0; i < ncheck; i++) {
1734 current = (TGeoNode *)fNodes->At(check_list[i]);
1735 if (current == node)
1736 return check_list[i];
1737 }
1738 return -1;
1739}
1740
1741////////////////////////////////////////////////////////////////////////////////
1742/// get index number for a given daughter
1743
1745{
1746 TGeoNode *current = nullptr;
1747 Int_t nd = GetNdaughters();
1748 if (!nd)
1749 return -1;
1750 for (Int_t i = 0; i < nd; i++) {
1751 current = (TGeoNode *)fNodes->At(i);
1752 if (current == node)
1753 return i;
1754 }
1755 return -1;
1756}
1757
1758////////////////////////////////////////////////////////////////////////////////
1759/// Get volume info for the browser.
1760
1762{
1763 TGeoVolume *vol = (TGeoVolume *)this;
1765 if (!painter)
1766 return nullptr;
1767 return (char *)painter->GetVolumeInfo(vol, px, py);
1768}
1769
1770////////////////////////////////////////////////////////////////////////////////
1771/// Returns true if cylindrical voxelization is optimal.
1772
1774{
1775 Int_t nd = GetNdaughters();
1776 if (!nd)
1777 return kFALSE;
1778 Int_t id;
1779 Int_t ncyl = 0;
1780 TGeoNode *node;
1781 for (id = 0; id < nd; id++) {
1782 node = (TGeoNode *)fNodes->At(id);
1783 ncyl += node->GetOptimalVoxels();
1784 }
1785 if (ncyl > (nd / 2))
1786 return kTRUE;
1787 return kFALSE;
1788}
1789
1790////////////////////////////////////////////////////////////////////////////////
1791/// Provide a pointer name containing uid.
1792
1794{
1795 static TString name;
1796 name.Form("p%s_%zx", GetName(), (size_t)this);
1797 return name.Data();
1798}
1799
1800////////////////////////////////////////////////////////////////////////////////
1801/// Getter for optimization structure.
1802
1804{
1805 if (fVoxels && !fVoxels->IsInvalid())
1806 return fVoxels;
1807 return nullptr;
1808}
1809
1810////////////////////////////////////////////////////////////////////////////////
1811/// Move perspective view focus to this volume
1812
1814{
1816 if (painter)
1817 painter->GrabFocus();
1818}
1819
1820////////////////////////////////////////////////////////////////////////////////
1821/// Returns true if the volume is an assembly or a scaled assembly.
1822
1824{
1825 return fShape->IsAssembly();
1826}
1827
1828////////////////////////////////////////////////////////////////////////////////
1829/// Clone this volume.
1830/// build a volume with same name, shape and medium
1831
1833{
1834 TGeoVolume *vol = new TGeoVolume(GetName(), fShape, fMedium);
1835 Int_t i;
1836 // copy volume attributes
1837 vol->SetTitle(GetTitle());
1838 vol->SetLineColor(GetLineColor());
1839 vol->SetLineStyle(GetLineStyle());
1840 vol->SetLineWidth(GetLineWidth());
1841 vol->SetFillColor(GetFillColor());
1842 vol->SetFillStyle(GetFillStyle());
1843 // copy other attributes
1844 Int_t nbits = 8 * sizeof(UInt_t);
1845 for (i = 0; i < nbits; i++)
1846 vol->SetAttBit(1 << i, TGeoAtt::TestAttBit(1 << i));
1847 for (i = 14; i < 24; i++)
1848 vol->SetBit(1 << i, TestBit(1 << i));
1849
1850 // copy field
1851 vol->SetField(fField);
1852 // Set bits
1853 for (i = 0; i < nbits; i++)
1854 vol->SetBit(1 << i, TObject::TestBit(1 << i));
1855 vol->SetBit(kVolumeClone);
1856 // copy nodes
1857 // CloneNodesAndConnect(vol);
1858 vol->MakeCopyNodes(this);
1859 // if volume is divided, copy finder
1860 vol->SetFinder(fFinder);
1861 // copy voxels
1862 TGeoVoxelFinder *voxels = nullptr;
1863 if (fVoxels) {
1864 voxels = new TGeoVoxelFinder(vol);
1865 vol->SetVoxelFinder(voxels);
1866 }
1867 // copy option, uid
1868 vol->SetOption(fOption);
1869 vol->SetNumber(fNumber);
1870 vol->SetNtotal(fNtotal);
1871 // copy extensions
1875 return vol;
1876}
1877
1878////////////////////////////////////////////////////////////////////////////////
1879/// Clone the array of nodes.
1880
1882{
1883 if (!fNodes)
1884 return;
1885 TGeoNode *node;
1886 Int_t nd = fNodes->GetEntriesFast();
1887 if (!nd)
1888 return;
1889 // create new list of nodes
1890 TObjArray *list = new TObjArray(nd);
1891 // attach it to new volume
1892 newmother->SetNodes(list);
1893 // ((TObject*)newmother)->SetBit(kVolumeImportNodes);
1894 for (Int_t i = 0; i < nd; i++) {
1895 // create copies of nodes and add them to list
1896 node = GetNode(i)->MakeCopyNode();
1897 if (!node) {
1898 Fatal("CloneNodesAndConnect", "cannot make copy node");
1899 return;
1900 }
1902 list->Add(node);
1903 }
1904}
1905
1906////////////////////////////////////////////////////////////////////////////////
1907/// make a new list of nodes and copy all nodes of other volume inside
1908
1910{
1911 Int_t nd = other->GetNdaughters();
1912 if (!nd)
1913 return;
1914 if (fNodes) {
1916 fNodes->Delete();
1917 delete fNodes;
1918 }
1919 fNodes = new TObjArray();
1920 for (Int_t i = 0; i < nd; i++)
1921 fNodes->Add(other->GetNode(i));
1923}
1924
1925////////////////////////////////////////////////////////////////////////////////
1926/// make a copy of this volume
1927/// build a volume with same name, shape and medium
1928
1930{
1932 // copy volume attributes
1933 vol->SetVisibility(IsVisible());
1934 vol->SetLineColor(GetLineColor());
1935 vol->SetLineStyle(GetLineStyle());
1936 vol->SetLineWidth(GetLineWidth());
1937 vol->SetFillColor(GetFillColor());
1938 vol->SetFillStyle(GetFillStyle());
1939 // copy field
1940 vol->SetField(fField);
1941 // if divided, copy division object
1942 if (fFinder) {
1943 // Error("MakeCopyVolume", "volume %s divided", GetName());
1944 vol->SetFinder(fFinder);
1945 }
1946 // Copy extensions
1950 // ((TObject*)vol)->SetBit(kVolumeImportNodes);
1951 ((TObject *)vol)->SetBit(kVolumeClone);
1953 return vol;
1954}
1955
1956////////////////////////////////////////////////////////////////////////////////
1957/// Make a copy of this volume which is reflected with respect to XY plane.
1958
1960{
1961 static TMap map(100);
1962 if (!fGeoManager->IsClosed()) {
1963 Error("MakeReflectedVolume", "Geometry must be closed.");
1964 return nullptr;
1965 }
1966 TGeoVolume *vol = (TGeoVolume *)map.GetValue(this);
1967 if (vol) {
1968 if (newname && newname[0])
1969 vol->SetName(newname);
1970 return vol;
1971 }
1972 // printf("Making reflection for volume: %s\n", GetName());
1973 vol = CloneVolume();
1974 if (!vol) {
1975 Fatal("MakeReflectedVolume", "Cannot clone volume %s\n", GetName());
1976 return nullptr;
1977 }
1978 map.Add((TObject *)this, vol);
1979 if (newname && newname[0])
1980 vol->SetName(newname);
1981 delete vol->GetNodes();
1982 vol->SetNodes(nullptr);
1985 // The volume is now properly cloned, but with the same shape.
1986 // Reflect the shape (if any) and connect it.
1987 if (fShape) {
1991 }
1992 // Reflect the daughters.
1993 Int_t nd = vol->GetNdaughters();
1994 if (!nd)
1995 return vol;
1996 TGeoNodeMatrix *node;
1999 if (!vol->GetFinder()) {
2000 for (Int_t i = 0; i < nd; i++) {
2001 node = (TGeoNodeMatrix *)vol->GetNode(i);
2002 local = node->GetMatrix();
2003 // printf("%s before\n", node->GetName());
2004 // local->Print();
2005 Bool_t reflected = local->IsReflection();
2007 local_cloned->RegisterYourself();
2008 node->SetMatrix(local_cloned);
2009 if (!reflected) {
2010 // We need to reflect only the translation and propagate to daughters.
2011 // H' = Sz * H * Sz
2012 local_cloned->ReflectZ(kTRUE);
2013 local_cloned->ReflectZ(kFALSE);
2014 // printf("%s after\n", node->GetName());
2015 // node->GetMatrix()->Print();
2017 node->SetVolume(new_vol);
2018 continue;
2019 }
2020 // The next daughter is already reflected, so reflect on Z everything and stop
2021 local_cloned->ReflectZ(kTRUE); // rot + tr
2022 // printf("%s already reflected... After:\n", node->GetName());
2023 // node->GetMatrix()->Print();
2024 }
2025 if (vol->GetVoxels())
2026 vol->GetVoxels()->Voxelize();
2027 return vol;
2028 }
2029 // Volume is divided, so we have to reflect the division.
2030 // printf(" ... divided %s\n", fFinder->ClassName());
2032 if (!new_finder) {
2033 Fatal("MakeReflectedVolume", "Could not copy finder for volume %s", GetName());
2034 return nullptr;
2035 }
2036 new_finder->SetVolume(vol);
2037 vol->SetFinder(new_finder);
2039 new_vol = nullptr;
2040 for (Int_t i = 0; i < nd; i++) {
2041 nodeoff = (TGeoNodeOffset *)vol->GetNode(i);
2042 nodeoff->SetFinder(new_finder);
2043 new_vol = nodeoff->GetVolume()->MakeReflectedVolume();
2044 nodeoff->SetVolume(new_vol);
2045 }
2046 return vol;
2047}
2048
2049////////////////////////////////////////////////////////////////////////////////
2050/// Set this volume as the TOP one (the whole geometry starts from here)
2051
2056
2057////////////////////////////////////////////////////////////////////////////////
2058/// Set the current tracking point.
2059
2064
2065////////////////////////////////////////////////////////////////////////////////
2066/// set the shape associated with this volume
2067
2069{
2070 if (!shape) {
2071 Error("SetShape", "No shape");
2072 return;
2073 }
2074 fShape = (TGeoShape *)shape;
2075}
2076
2077////////////////////////////////////////////////////////////////////////////////
2078/// sort nodes by decreasing volume of the bounding box. ONLY nodes comes first,
2079/// then overlapping nodes and finally division nodes.
2080
2082{
2083 if (!Valid()) {
2084 Error("SortNodes", "Bounding box not valid");
2085 return;
2086 }
2087 Int_t nd = GetNdaughters();
2088 // printf("volume : %s, nd=%i\n", GetName(), nd);
2089 if (!nd)
2090 return;
2091 if (fFinder)
2092 return;
2093 // printf("Nodes for %s\n", GetName());
2094 Int_t id = 0;
2095 TGeoNode *node = nullptr;
2096 TObjArray *nodes = new TObjArray(nd);
2097 Int_t inode = 0;
2098 // first put ONLY's
2099 for (id = 0; id < nd; id++) {
2100 node = GetNode(id);
2101 if (node->InheritsFrom(TGeoNodeOffset::Class()) || node->IsOverlapping())
2102 continue;
2103 nodes->Add(node);
2104 // printf("inode %i ONLY\n", inode);
2105 inode++;
2106 }
2107 // second put overlapping nodes
2108 for (id = 0; id < nd; id++) {
2109 node = GetNode(id);
2110 if (node->InheritsFrom(TGeoNodeOffset::Class()) || (!node->IsOverlapping()))
2111 continue;
2112 nodes->Add(node);
2113 // printf("inode %i MANY\n", inode);
2114 inode++;
2115 }
2116 // third put the divided nodes
2117 if (fFinder) {
2119 for (id = 0; id < nd; id++) {
2120 node = GetNode(id);
2121 if (!node->InheritsFrom(TGeoNodeOffset::Class()))
2122 continue;
2123 nodes->Add(node);
2124 // printf("inode %i DIV\n", inode);
2125 inode++;
2126 }
2127 }
2128 if (inode != nd)
2129 printf(" volume %s : number of nodes does not match!!!\n", GetName());
2130 delete fNodes;
2131 fNodes = nodes;
2132}
2133
2134////////////////////////////////////////////////////////////////////////////////
2135/// Stream an object of class TGeoVolume.
2136
2138{
2139 if (R__b.IsReading()) {
2140 R__b.ReadClassBuffer(TGeoVolume::Class(), this);
2141 if (fVoxels && fVoxels->IsInvalid())
2142 Voxelize("");
2143 } else {
2144 if (!fVoxels) {
2145 R__b.WriteClassBuffer(TGeoVolume::Class(), this);
2146 } else {
2149 fVoxels = nullptr;
2150 R__b.WriteClassBuffer(TGeoVolume::Class(), this);
2151 fVoxels = voxels;
2152 } else {
2153 R__b.WriteClassBuffer(TGeoVolume::Class(), this);
2154 }
2155 }
2156 }
2157}
2158
2159////////////////////////////////////////////////////////////////////////////////
2160/// Set the current options (none implemented)
2161
2163{
2164 fOption = option;
2165}
2166
2167////////////////////////////////////////////////////////////////////////////////
2168/// Set the line color.
2169
2174
2175////////////////////////////////////////////////////////////////////////////////
2176/// Set the line style.
2177
2182
2183////////////////////////////////////////////////////////////////////////////////
2184/// Set the line width.
2185
2190
2191////////////////////////////////////////////////////////////////////////////////
2192/// get the pointer to a daughter node
2193
2195{
2196 if (!fNodes)
2197 return nullptr;
2198 TGeoNode *node = (TGeoNode *)fNodes->FindObject(name);
2199 return node;
2200}
2201
2202////////////////////////////////////////////////////////////////////////////////
2203/// get the total size in bytes for this volume
2204
2206{
2207 Int_t count = 28 + 2 + 6 + 4 + 0; // TNamed+TGeoAtt+TAttLine+TAttFill+TAtt3D
2208 count += fName.Capacity() + fTitle.Capacity(); // name+title
2209 count += 7 * sizeof(char *); // fShape + fMedium + fFinder + fField + fNodes + 2 extensions
2210 count += fOption.Capacity(); // fOption
2211 if (fShape)
2212 count += fShape->GetByteCount();
2213 if (fFinder)
2214 count += fFinder->GetByteCount();
2215 if (fNodes) {
2216 count += 32 + 4 * fNodes->GetEntries(); // TObjArray
2217 TIter next(fNodes);
2218 TGeoNode *node;
2219 while ((node = (TGeoNode *)next()))
2220 count += node->GetByteCount();
2221 }
2222 return count;
2223}
2224
2225////////////////////////////////////////////////////////////////////////////////
2226/// loop all nodes marked as overlaps and find overlapping brothers
2227
2229{
2230 if (!Valid()) {
2231 Error("FindOverlaps", "Bounding box not valid");
2232 return;
2233 }
2234 if (!fVoxels)
2235 return;
2236 Int_t nd = GetNdaughters();
2237 if (!nd)
2238 return;
2239 TGeoNode *node = nullptr;
2240 Int_t inode = 0;
2241 for (inode = 0; inode < nd; inode++) {
2242 node = GetNode(inode);
2243 if (!node->IsOverlapping())
2244 continue;
2246 }
2247}
2248
2249////////////////////////////////////////////////////////////////////////////////
2250/// Remove an existing daughter.
2251
2253{
2254 if (!fNodes || !fNodes->GetEntriesFast())
2255 return;
2256 if (!fNodes->Remove(node))
2257 return;
2258 fNodes->Compress();
2259 if (fVoxels)
2261 if (IsAssembly())
2263}
2264
2265////////////////////////////////////////////////////////////////////////////////
2266/// Replace an existing daughter with a new volume having the same name but
2267/// possibly a new shape, position or medium. Not allowed for positioned assemblies.
2268/// For division cells, the new shape/matrix are ignored.
2269
2271{
2273 if (ind < 0)
2274 return nullptr;
2275 TGeoVolume *oldvol = nodeorig->GetVolume();
2276 if (oldvol->IsAssembly()) {
2277 Error("ReplaceNode", "Cannot replace node %s since it is an assembly", nodeorig->GetName());
2278 return nullptr;
2279 }
2280 TGeoShape *shape = oldvol->GetShape();
2281 if (newshape && !nodeorig->IsOffset())
2282 shape = newshape;
2283 TGeoMedium *med = oldvol->GetMedium();
2284 if (newmed)
2285 med = newmed;
2286 // Make a new volume
2287 TGeoVolume *vol = new TGeoVolume(oldvol->GetName(), shape, med);
2288 // copy volume attributes
2289 vol->SetVisibility(oldvol->IsVisible());
2290 vol->SetLineColor(oldvol->GetLineColor());
2291 vol->SetLineStyle(oldvol->GetLineStyle());
2292 vol->SetLineWidth(oldvol->GetLineWidth());
2293 vol->SetFillColor(oldvol->GetFillColor());
2294 vol->SetFillStyle(oldvol->GetFillStyle());
2295 // copy field
2296 vol->SetField(oldvol->GetField());
2297 // Make a copy of the node
2298 TGeoNode *newnode = nodeorig->MakeCopyNode();
2299 if (!newnode) {
2300 Fatal("ReplaceNode", "Cannot make copy node for %s", nodeorig->GetName());
2301 return nullptr;
2302 }
2303 // Change the volume for the new node
2304 newnode->SetVolume(vol);
2305 // Replace the matrix
2306 if (newpos && !nodeorig->IsOffset()) {
2308 nodemat->SetMatrix(newpos);
2309 }
2310 // Replace nodeorig with new one
2313 if (fVoxels)
2315 if (IsAssembly())
2317 return newnode;
2318}
2319
2320////////////////////////////////////////////////////////////////////////////////
2321/// Select this volume as matching an arbitrary criteria. The volume is added to
2322/// a static list and the flag TGeoVolume::kVolumeSelected is set. All flags need
2323/// to be reset at the end by calling the method with CLEAR=true. This will also clear
2324/// the list.
2325
2327{
2328 static TObjArray array(256);
2329 static Int_t len = 0;
2330 Int_t i;
2331 TObject *vol;
2332 if (clear) {
2333 for (i = 0; i < len; i++) {
2334 vol = array.At(i);
2336 }
2337 array.Clear();
2338 len = 0;
2339 return;
2340 }
2342 array.AddAtAndExpand(this, len++);
2343}
2344
2345////////////////////////////////////////////////////////////////////////////////
2346/// set visibility of this volume
2347
2349{
2351 if (fGeoManager->IsClosed())
2354 TSeqCollection *brlist = gROOT->GetListOfBrowsers();
2355 TIter next(brlist);
2356 TBrowser *browser = nullptr;
2357 while ((browser = (TBrowser *)next())) {
2358 browser->CheckObjectItem(this, vis);
2359 browser->Refresh();
2360 }
2361}
2362
2363////////////////////////////////////////////////////////////////////////////////
2364/// Set visibility for containers.
2365
2376
2377////////////////////////////////////////////////////////////////////////////////
2378/// Set visibility for leaves.
2379
2390
2391////////////////////////////////////////////////////////////////////////////////
2392/// Set visibility for leaves.
2393
2406
2407////////////////////////////////////////////////////////////////////////////////
2408/// Check if the shape of this volume is valid.
2409
2411{
2412 return fShape->IsValidBox();
2413}
2414
2415////////////////////////////////////////////////////////////////////////////////
2416/// Find a daughter node having VOL as volume and fill TGeoManager::fHMatrix
2417/// with its global matrix.
2418
2420{
2421 if (vol == this)
2422 return kTRUE;
2423 Int_t nd = GetNdaughters();
2424 if (!nd)
2425 return kFALSE;
2427 if (!global)
2428 return kFALSE;
2429 TGeoNode *dnode;
2432 Int_t i;
2433 for (i = 0; i < nd; i++) {
2434 dnode = GetNode(i);
2435 dvol = dnode->GetVolume();
2436 if (dvol == vol) {
2437 local = dnode->GetMatrix();
2438 global->MultiplyLeft(local);
2439 return kTRUE;
2440 }
2441 }
2442 for (i = 0; i < nd; i++) {
2443 dnode = GetNode(i);
2444 dvol = dnode->GetVolume();
2445 if (dvol->FindMatrixOfDaughterVolume(vol))
2446 return kTRUE;
2447 }
2448 return kFALSE;
2449}
2450
2451////////////////////////////////////////////////////////////////////////////////
2452/// set visibility for daughters
2453
2461
2462////////////////////////////////////////////////////////////////////////////////
2463/// build the voxels for this volume
2464
2466{
2467 if (!Valid()) {
2468 Error("Voxelize", "Bounding box not valid");
2469 return;
2470 }
2471 // do not voxelize divided volumes
2472 if (fFinder)
2473 return;
2474 // or final leaves
2475 Int_t nd = GetNdaughters();
2476 if (!nd)
2477 return;
2478 // If this is an assembly, re-compute bounding box
2479 if (IsAssembly())
2481 // delete old voxelization if any
2482 if (fVoxels) {
2484 delete fVoxels;
2485 fVoxels = nullptr;
2486 }
2487 // Create the voxels structure
2488 fVoxels = new TGeoVoxelFinder(this);
2490 if (fVoxels) {
2491 if (fVoxels->IsInvalid()) {
2492 delete fVoxels;
2493 fVoxels = nullptr;
2494 }
2495 }
2496}
2497
2498////////////////////////////////////////////////////////////////////////////////
2499/// Estimate the weight of a volume (in kg) with SIGMA(M)/M better than PRECISION.
2500/// Option can contain : v - verbose, a - analytical (default)
2501
2503{
2505 if (top != this)
2507 else
2508 top = nullptr;
2509 Double_t weight = fGeoManager->Weight(precision, option);
2510 if (top)
2512 return weight;
2513}
2514
2515////////////////////////////////////////////////////////////////////////////////
2516/// Analytical computation of the weight.
2517
2519{
2520 Double_t capacity = Capacity();
2521 Double_t weight = 0.0;
2522 Int_t i;
2523 Int_t nd = GetNdaughters();
2525 for (i = 0; i < nd; i++) {
2526 daughter = GetNode(i)->GetVolume();
2527 weight += daughter->WeightA();
2528 capacity -= daughter->Capacity();
2529 }
2530 Double_t density = 0.0;
2531 if (!IsAssembly()) {
2532 if (fMedium)
2534 if (density < 0.01)
2535 density = 0.0; // do not weight gases
2536 }
2537 weight += 0.001 * capacity * density; //[kg]
2538 return weight;
2539}
2540
2542
2543////////////////////////////////////////////////////////////////////////////////
2544/// dummy constructor
2545
2547{
2548 fVolumes = nullptr;
2549 fDivision = nullptr;
2550 fNumed = 0;
2551 fNdiv = 0;
2552 fAxis = 0;
2553 fStart = 0;
2554 fStep = 0;
2555 fAttSet = kFALSE;
2557}
2558
2559////////////////////////////////////////////////////////////////////////////////
2560/// default constructor
2561
2563{
2564 fVolumes = new TObjArray();
2565 fDivision = nullptr;
2566 fNumed = 0;
2567 fNdiv = 0;
2568 fAxis = 0;
2569 fStart = 0;
2570 fStep = 0;
2571 fAttSet = kFALSE;
2573 SetName(name);
2574 SetMedium(med);
2575 fGeoManager->AddVolume(this);
2576 // printf("--- volume multi %s created\n", name);
2577}
2578
2579////////////////////////////////////////////////////////////////////////////////
2580/// Destructor
2581
2583{
2584 if (fVolumes)
2585 delete fVolumes;
2586}
2587
2588////////////////////////////////////////////////////////////////////////////////
2589/// Add a volume with valid shape to the list of volumes. Copy all existing nodes
2590/// to this volume
2591
2593{
2594 Int_t idx = fVolumes->GetEntriesFast();
2595 fVolumes->AddAtAndExpand(vol, idx);
2596 vol->SetUniqueID(idx + 1);
2599 if (fDivision) {
2601 if (!div) {
2602 Fatal("AddVolume", "Cannot divide volume %s", vol->GetName());
2603 return;
2604 }
2605 for (Int_t i = 0; i < div->GetNvolumes(); i++) {
2606 cell = div->GetVolume(i);
2608 }
2609 }
2610 if (fNodes) {
2611 Int_t nd = fNodes->GetEntriesFast();
2612 for (Int_t id = 0; id < nd; id++) {
2613 TGeoNode *node = (TGeoNode *)fNodes->At(id);
2614 Bool_t many = node->IsOverlapping();
2615 if (many)
2616 vol->AddNodeOverlap(node->GetVolume(), node->GetNumber(), node->GetMatrix());
2617 else
2618 vol->AddNode(node->GetVolume(), node->GetNumber(), node->GetMatrix());
2619 }
2620 }
2621 // vol->MakeCopyNodes(this);
2622}
2623
2624////////////////////////////////////////////////////////////////////////////////
2625/// Add a new node to the list of nodes. This is the usual method for adding
2626/// daughters inside the container volume.
2627
2629{
2632 TGeoVolume *volume = nullptr;
2633 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2634 volume = GetVolume(ivo);
2635 volume->SetLineColor(GetLineColor());
2636 volume->SetLineStyle(GetLineStyle());
2637 volume->SetLineWidth(GetLineWidth());
2638 volume->SetVisibility(IsVisible());
2639 volume->AddNode(vol, copy_no, mat, option);
2640 }
2641 // printf("--- vmulti %s : node %s added to %i components\n", GetName(), vol->GetName(), nvolumes);
2642 return n;
2643}
2644
2645////////////////////////////////////////////////////////////////////////////////
2646/// Add a new node to the list of nodes, This node is possibly overlapping with other
2647/// daughters of the volume or extruding the volume.
2648
2650{
2653 TGeoVolume *volume = nullptr;
2654 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2655 volume = GetVolume(ivo);
2656 volume->SetLineColor(GetLineColor());
2657 volume->SetLineStyle(GetLineStyle());
2658 volume->SetLineWidth(GetLineWidth());
2659 volume->SetVisibility(IsVisible());
2660 volume->AddNodeOverlap(vol, copy_no, mat, option);
2661 }
2662 // printf("--- vmulti %s : node ovlp %s added to %i components\n", GetName(), vol->GetName(), nvolumes);
2663}
2664
2665////////////////////////////////////////////////////////////////////////////////
2666/// Returns the last shape.
2667
2669{
2671 if (!vol)
2672 return nullptr;
2673 return vol->GetShape();
2674}
2675
2676////////////////////////////////////////////////////////////////////////////////
2677/// division of multiple volumes
2678
2680 Int_t numed, const char *option)
2681{
2682 if (fDivision) {
2683 Error("Divide", "volume %s already divided", GetName());
2684 return nullptr;
2685 }
2688 if (numed) {
2690 if (!medium) {
2691 Error("Divide", "Invalid medium number %d for division volume %s", numed, divname);
2692 medium = fMedium;
2693 }
2694 }
2695 if (!nvolumes) {
2696 // this is a virtual volume
2698 fNumed = medium->GetId();
2699 fOption = option;
2700 fAxis = iaxis;
2701 fNdiv = ndiv;
2702 fStart = start;
2703 fStep = step;
2704 // nothing else to do at this stage
2705 return fDivision;
2706 }
2707 TGeoVolume *vol = nullptr;
2709 if (medium)
2710 fNumed = medium->GetId();
2711 fOption = option;
2712 fAxis = iaxis;
2713 fNdiv = ndiv;
2714 fStart = start;
2715 fStep = step;
2716 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2717 vol = GetVolume(ivo);
2718 vol->SetLineColor(GetLineColor());
2719 vol->SetLineStyle(GetLineStyle());
2720 vol->SetLineWidth(GetLineWidth());
2721 vol->SetVisibility(IsVisible());
2722 fDivision->AddVolume(vol->Divide(divname, iaxis, ndiv, start, step, numed, option));
2723 }
2724 // printf("--- volume multi %s (%i volumes) divided\n", GetName(), nvolumes);
2725 if (numed)
2727 return fDivision;
2728}
2729
2730////////////////////////////////////////////////////////////////////////////////
2731/// Make a copy of this volume
2732/// build a volume with same name, shape and medium
2733
2735{
2737 Int_t i = 0;
2738 // copy volume attributes
2739 vol->SetVisibility(IsVisible());
2740 vol->SetLineColor(GetLineColor());
2741 vol->SetLineStyle(GetLineStyle());
2742 vol->SetLineWidth(GetLineWidth());
2743 vol->SetFillColor(GetFillColor());
2744 vol->SetFillStyle(GetFillStyle());
2745 // copy field
2746 vol->SetField(fField);
2747 // Copy extensions
2750 // if divided, copy division object
2751 // if (fFinder) {
2752 // Error("MakeCopyVolume", "volume %s divided", GetName());
2753 // vol->SetFinder(fFinder);
2754 // }
2755 if (fDivision) {
2759 if (!div) {
2760 Fatal("MakeCopyVolume", "Cannot divide volume %s", vol->GetName());
2761 return nullptr;
2762 }
2763 for (i = 0; i < div->GetNvolumes(); i++) {
2764 cell = div->GetVolume(i);
2766 }
2767 }
2768
2769 if (!fNodes)
2770 return vol;
2771 TGeoNode *node;
2772 Int_t nd = fNodes->GetEntriesFast();
2773 if (!nd)
2774 return vol;
2775 // create new list of nodes
2776 TObjArray *list = new TObjArray();
2777 // attach it to new volume
2778 vol->SetNodes(list);
2779 ((TObject *)vol)->SetBit(kVolumeImportNodes);
2780 for (i = 0; i < nd; i++) {
2781 // create copies of nodes and add them to list
2782 node = GetNode(i)->MakeCopyNode();
2783 if (!node) {
2784 Fatal("MakeCopyNode", "cannot make copy node for daughter %d of %s", i, GetName());
2785 return nullptr;
2786 }
2787 node->SetMotherVolume(vol);
2788 list->Add(node);
2789 }
2790 return vol;
2791}
2792
2793////////////////////////////////////////////////////////////////////////////////
2794/// Set the line color for all components.
2795
2797{
2800 TGeoVolume *vol = nullptr;
2801 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2802 vol = GetVolume(ivo);
2803 vol->SetLineColor(lcolor);
2804 }
2805}
2806
2807////////////////////////////////////////////////////////////////////////////////
2808/// Set the line style for all components.
2809
2811{
2814 TGeoVolume *vol = nullptr;
2815 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2816 vol = GetVolume(ivo);
2817 vol->SetLineStyle(lstyle);
2818 }
2819}
2820
2821////////////////////////////////////////////////////////////////////////////////
2822/// Set the line width for all components.
2823
2825{
2828 TGeoVolume *vol = nullptr;
2829 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2830 vol = GetVolume(ivo);
2831 vol->SetLineWidth(lwidth);
2832 }
2833}
2834
2835////////////////////////////////////////////////////////////////////////////////
2836/// Set medium for a multiple volume.
2837
2839{
2842 TGeoVolume *vol = nullptr;
2843 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2844 vol = GetVolume(ivo);
2845 vol->SetMedium(med);
2846 }
2847}
2848
2849////////////////////////////////////////////////////////////////////////////////
2850/// Set visibility for all components.
2851
2853{
2856 TGeoVolume *vol = nullptr;
2857 for (Int_t ivo = 0; ivo < nvolumes; ivo++) {
2858 vol = GetVolume(ivo);
2859 vol->SetVisibility(vis);
2860 }
2861}
2862
2864
2865////////////////////////////////////////////////////////////////////////////////
2866/// Constructor.
2867
2869
2870////////////////////////////////////////////////////////////////////////////////
2871/// Destructor.
2872
2874
2875////////////////////////////////////////////////////////////////////////////////
2876
2882
2883////////////////////////////////////////////////////////////////////////////////
2884
2886{
2887 std::lock_guard<std::mutex> guard(fMutex);
2889 std::vector<ThreadData_t *>::iterator i = fThreadData.begin();
2890 while (i != fThreadData.end()) {
2891 delete *i;
2892 ++i;
2893 }
2894 fThreadData.clear();
2895 fThreadSize = 0;
2896}
2897
2898////////////////////////////////////////////////////////////////////////////////
2899
2901{
2902 std::lock_guard<std::mutex> guard(fMutex);
2903 // Create assembly thread data here
2904 fThreadData.resize(nthreads);
2906 for (Int_t tid = 0; tid < nthreads; tid++) {
2907 if (fThreadData[tid] == nullptr) {
2909 }
2910 }
2912}
2913
2914////////////////////////////////////////////////////////////////////////////////
2915
2920
2921////////////////////////////////////////////////////////////////////////////////
2922
2927
2928////////////////////////////////////////////////////////////////////////////////
2929
2934
2935////////////////////////////////////////////////////////////////////////////////
2936
2941
2942////////////////////////////////////////////////////////////////////////////////
2943/// Default constructor
2944
2950
2951////////////////////////////////////////////////////////////////////////////////
2952/// Constructor. Just the name has to be provided. Assemblies does not have their own
2953/// shape or medium.
2954
2956{
2957 fName = name;
2958 fName = fName.Strip();
2959 fShape = new TGeoShapeAssembly(this);
2960 if (fGeoManager)
2961 fNumber = fGeoManager->AddVolume(this);
2962 fThreadSize = 0;
2964}
2965
2966////////////////////////////////////////////////////////////////////////////////
2967/// Destructor. The assembly is owner of its "shape".
2968
2970{
2972 if (fShape)
2973 delete fShape;
2974}
2975
2976////////////////////////////////////////////////////////////////////////////////
2977/// Add a component to the assembly.
2978
2980{
2982 // ((TGeoShapeAssembly*)fShape)->RecomputeBoxLast();
2983 ((TGeoShapeAssembly *)fShape)->NeedsBBoxRecompute();
2984 return node;
2985}
2986
2987////////////////////////////////////////////////////////////////////////////////
2988/// Add an overlapping node - not allowed for assemblies.
2989
2991{
2992 Warning("AddNodeOverlap",
2993 "Declaring assembly %s as possibly overlapping inside %s not allowed. Using AddNode instead !",
2994 vol->GetName(), GetName());
2995 AddNode(vol, copy_no, mat, option);
2996}
2997
2998////////////////////////////////////////////////////////////////////////////////
2999/// Clone this volume.
3000/// build a volume with same name, shape and medium
3001
3003{
3005 Int_t i;
3006 // copy other attributes
3007 Int_t nbits = 8 * sizeof(UInt_t);
3008 for (i = 0; i < nbits; i++)
3009 vol->SetAttBit(1 << i, TGeoAtt::TestAttBit(1 << i));
3010 for (i = 14; i < 24; i++)
3011 vol->SetBit(1 << i, TestBit(1 << i));
3012
3013 // copy field
3014 vol->SetField(fField);
3015 // Set bits
3016 for (i = 0; i < nbits; i++)
3017 vol->SetBit(1 << i, TObject::TestBit(1 << i));
3018 vol->SetBit(kVolumeClone);
3019 // make copy nodes
3020 vol->MakeCopyNodes(this);
3021 // CloneNodesAndConnect(vol);
3022 ((TGeoShapeAssembly *)vol->GetShape())->NeedsBBoxRecompute();
3023 // copy voxels
3024 TGeoVoxelFinder *voxels = nullptr;
3025 if (fVoxels) {
3026 voxels = new TGeoVoxelFinder(vol);
3027 vol->SetVoxelFinder(voxels);
3028 }
3029 // copy option, uid
3030 vol->SetOption(fOption);
3031 vol->SetNumber(fNumber);
3032 vol->SetNtotal(fNtotal);
3033 vol->SetTitle(GetTitle());
3034 // copy extensions
3037 return vol;
3038}
3039
3040////////////////////////////////////////////////////////////////////////////////
3041/// Division makes no sense for assemblies.
3042
3044{
3045 Error("Divide", "Assemblies cannot be divided");
3046 return nullptr;
3047}
3048
3049////////////////////////////////////////////////////////////////////////////////
3050/// Assign to the assembly a collection of identical volumes positioned according
3051/// a predefined pattern. The option can be spaced out or touching depending on the empty
3052/// space between volumes.
3053
3055{
3056 if (fNodes) {
3057 Error("Divide", "Cannot divide assembly %s since it has nodes", GetName());
3058 return nullptr;
3059 }
3060 if (fFinder) {
3061 Error("Divide", "Assembly %s already divided", GetName());
3062 return nullptr;
3063 }
3064 Int_t ncells = pattern->GetNdiv();
3065 if (!ncells || pattern->GetStep() <= 0) {
3066 Error("Divide", "Pattern finder for dividing assembly %s not initialized. Use SetRange() method.", GetName());
3067 return nullptr;
3068 }
3069 fFinder = pattern;
3070 TString opt(option);
3071 opt.ToLower();
3072 if (opt.Contains("spacedout"))
3074 else
3076 // Position volumes
3077 for (Int_t i = 0; i < ncells; i++) {
3078 fFinder->cd(i);
3079 TGeoNodeOffset *node = new TGeoNodeOffset(cell, i, 0.);
3080 node->SetFinder(fFinder);
3081 fNodes->Add(node);
3082 }
3083 return cell;
3084}
3085
3086////////////////////////////////////////////////////////////////////////////////
3087/// Make a clone of volume VOL but which is an assembly.
3088
3090{
3091 if (volorig->IsAssembly() || volorig->IsVolumeMulti())
3092 return nullptr;
3093 Int_t nd = volorig->GetNdaughters();
3094 if (!nd)
3095 return nullptr;
3096 TGeoVolumeAssembly *vol = new TGeoVolumeAssembly(volorig->GetName());
3097 Int_t i;
3098 // copy other attributes
3099 Int_t nbits = 8 * sizeof(UInt_t);
3100 for (i = 0; i < nbits; i++)
3101 vol->SetAttBit(1 << i, volorig->TestAttBit(1 << i));
3102 for (i = 14; i < 24; i++)
3103 vol->SetBit(1 << i, volorig->TestBit(1 << i));
3104
3105 // copy field
3106 vol->SetField(volorig->GetField());
3107 // Set bits
3108 for (i = 0; i < nbits; i++)
3109 vol->SetBit(1 << i, volorig->TestBit(1 << i));
3110 vol->SetBit(kVolumeClone);
3111 // make copy nodes
3112 vol->MakeCopyNodes(volorig);
3113 // volorig->CloneNodesAndConnect(vol);
3114 vol->GetShape()->ComputeBBox();
3115 // copy voxels
3116 TGeoVoxelFinder *voxels = nullptr;
3117 if (volorig->GetVoxels()) {
3118 voxels = new TGeoVoxelFinder(vol);
3119 vol->SetVoxelFinder(voxels);
3120 }
3121 // copy option, uid
3122 vol->SetOption(volorig->GetOption());
3123 vol->SetNumber(volorig->GetNumber());
3124 vol->SetNtotal(volorig->GetNtotal());
3125 return vol;
3126}
#define b(i)
Definition RSha256.hxx:100
#define f(i)
Definition RSha256.hxx:104
#define s1(x)
Definition RSha256.hxx:91
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
unsigned int UInt_t
Unsigned integer 4 bytes (unsigned int)
Definition RtypesCore.h:60
short Width_t
Line width (short)
Definition RtypesCore.h:98
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
#define ClassImp(name)
Definition Rtypes.h:376
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
R__EXTERN TEnv * gEnv
Definition TEnv.h:170
Option_t Option_t option
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 filename
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 char Point_t Rectangle_t WindowAttributes_t index
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize id
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
char name[80]
Definition TGX11.cxx:110
R__EXTERN TGeoManager * gGeoManager
R__EXTERN TGeoIdentity * gGeoIdentity
Definition TGeoMatrix.h:537
#define gROOT
Definition TROOT.h:411
R__EXTERN TStyle * gStyle
Definition TStyle.h:442
virtual Color_t GetFillColor() const
Return the fill area color.
Definition TAttFill.h:31
virtual Style_t GetFillStyle() const
Return the fill area style.
Definition TAttFill.h:32
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition TAttFill.h:38
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition TAttFill.h:40
virtual Color_t GetLineColor() const
Return the line color.
Definition TAttLine.h:35
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition TAttLine.h:44
virtual Width_t GetLineWidth() const
Return the line width.
Definition TAttLine.h:37
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition TAttLine.h:45
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:42
virtual Style_t GetLineStyle() const
Return the line style.
Definition TAttLine.h:36
virtual void SaveLineAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1, Int_t widdef=1)
Save line attributes as C++ statement(s) on output stream out.
Definition TAttLine.cxx:275
Using a TBrowser one can browse all ROOT objects.
Definition TBrowser.h:37
Buffer base class used for serializing objects.
Definition TBuffer.h:43
TDirectory::TContext keeps track and restore the current directory.
Definition TDirectory.h:89
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition TEnv.cxx:491
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition TFile.h:131
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:3760
Bool_t IsVisRaytrace() const
Definition TGeoAtt.h:82
virtual void SetVisOnly(Bool_t flag=kTRUE)
Set branch type visibility.
Definition TGeoAtt.cxx:94
Bool_t TestAttBit(UInt_t f) const
Definition TGeoAtt.h:64
virtual void SetVisLeaves(Bool_t flag=kTRUE)
Set branch type visibility.
Definition TGeoAtt.cxx:84
@ kSaveNodesAtt
Definition TGeoAtt.h:50
@ kSavePrimitiveAtt
Definition TGeoAtt.h:50
void SetVisDaughters(Bool_t vis=kTRUE)
Set visibility for the daughters.
Definition TGeoAtt.cxx:116
void ResetAttBit(UInt_t f)
Definition TGeoAtt.h:63
void SetVisRaytrace(Bool_t flag=kTRUE)
Definition TGeoAtt.h:66
Bool_t IsVisDaughters() const
Definition TGeoAtt.h:84
virtual void SetVisibility(Bool_t vis=kTRUE)
Set visibility for this object.
Definition TGeoAtt.cxx:104
void SetAttBit(UInt_t f)
Definition TGeoAtt.h:61
void SetVisTouched(Bool_t vis=kTRUE)
Mark visualization attributes as "modified".
Definition TGeoAtt.cxx:138
virtual void SetVisContainers(Bool_t flag=kTRUE)
Set branch type visibility.
Definition TGeoAtt.cxx:76
Class describing rotation + translation.
Definition TGeoMatrix.h:317
Composite shapes are Boolean combinations of two or more shape components.
ABC for user objects attached to TGeoVolume or TGeoNode.
virtual TGeoExtension * Grab()=0
virtual void Release() const =0
Matrix class used for computing global transformations Should NOT be used for node definition.
Definition TGeoMatrix.h:458
The manager class for any TGeo geometry.
Definition TGeoManager.h:45
TObjArray * GetListOfOverlaps()
TList * GetListOfMedia() const
void SetUserPaintVolume(TGeoVolume *vol)
TVirtualGeoChecker * GetGeomChecker()
Make a default checker if none present. Returns pointer to it.
TObjArray * GetListOfVolumes() const
void ClearOverlaps()
Clear the list of overlaps.
TObjArray * GetListOfMatrices() const
Bool_t IsClosed() const
TVirtualGeoPainter * GetGeomPainter()
Make a default painter if none present. Returns pointer to it.
void SetVisOption(Int_t option=0)
set drawing mode :
Int_t AddMaterial(const TGeoMaterial *material)
Add a material to the list. Returns index of the material in list.
TGeoHMatrix * GetHMatrix()
Return stored current matrix (global matrix of the next touched node).
Int_t AddVolume(TGeoVolume *volume)
Add a volume to the list. Returns index of the volume in list.
Bool_t IsStreamingVoxels() const
void SetCurrentPoint(Double_t *point)
Int_t GetVisOption() const
Returns current depth to which geometry is drawn.
void SetTopVolume(TGeoVolume *vol)
Set the top volume and corresponding node as starting point of the geometry.
void ClearShape(const TGeoShape *shape)
Remove a shape from the list of shapes.
TGeoMedium * GetMedium(const char *medium) const
Search for a named tracking medium. All trailing blanks stripped.
Double_t Weight(Double_t precision=0.01, Option_t *option="va")
Estimate weight of volume VOL with a precision SIGMA(W)/W better than PRECISION.
void SetNsegments(Int_t nseg)
Set number of segments for approximating circles in drawing.
static UInt_t GetExportPrecision()
TVirtualGeoPainter * GetPainter() const
Bool_t IsCheckingOverlaps() const
void RandomPoints(const TGeoVolume *vol, Int_t npoints=10000, Option_t *option="")
Draw random points in the bounding box of a volume.
TList * GetListOfMaterials() const
void RandomRays(Int_t nrays=1000, Double_t startx=0, Double_t starty=0, Double_t startz=0, const char *target_vol=nullptr, Bool_t check_norm=kFALSE)
Randomly shoot nrays and plot intersections with surfaces for current top node.
void SetAllIndex()
Assigns uid's for all materials,media and matrices.
TObjArray * GetListOfShapes() const
TGeoVolume * GetTopVolume() const
static Int_t ThreadId()
Translates the current thread id to an ordinal number.
Int_t AddShape(const TGeoShape *shape)
Add a shape to the list. Returns index of the shape in list.
void SortOverlaps()
Sort overlaps by decreasing overlap distance. Extrusions comes first.
Base class describing materials.
void SetUsed(Bool_t flag=kTRUE)
void Print(const Option_t *option="") const override
print characteristics of this material
virtual Double_t GetDensity() const
Geometrical transformation package.
Definition TGeoMatrix.h:38
Media are used to store properties related to tracking and which are useful only when using geometry ...
Definition TGeoMedium.h:23
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
TGeoMaterial * GetMaterial() const
Definition TGeoMedium.h:49
const char * GetPointerName() const
Provide a pointer name containing uid.
A node containing local transformation.
Definition TGeoNode.h:154
void SetMatrix(const TGeoMatrix *matrix)
Matrix setter.
Definition TGeoNode.cxx:830
TGeoMatrix * GetMatrix() const override
Definition TGeoNode.h:171
Node containing an offset.
Definition TGeoNode.h:184
void SetFinder(TGeoPatternFinder *finder)
Definition TGeoNode.h:210
static TClass * Class()
A node represent a volume positioned inside another.They store links to both volumes and to the TGeoM...
Definition TGeoNode.h:39
Bool_t IsOverlapping() const
Definition TGeoNode.h:107
TGeoVolume * GetVolume() const
Definition TGeoNode.h:99
void SetVolume(TGeoVolume *volume)
Definition TGeoNode.h:117
virtual Int_t GetByteCount() const
Definition TGeoNode.h:82
void SetOverlapping(Bool_t flag=kTRUE)
Definition TGeoNode.h:120
virtual Int_t GetOptimalVoxels() const
Definition TGeoNode.h:101
virtual TGeoMatrix * GetMatrix() const =0
void SetMotherVolume(TGeoVolume *mother)
Definition TGeoNode.h:125
virtual TGeoNode * MakeCopyNode() const
Definition TGeoNode.h:113
void SetVirtual()
Definition TGeoNode.h:121
Int_t GetNumber() const
Definition TGeoNode.h:93
void SetNumber(Int_t number)
Definition TGeoNode.h:118
base finder class for patterns. A pattern is specifying a division type
virtual void cd(Int_t)
void SetSpacedOut(Bool_t flag)
void SetDivIndex(Int_t index)
virtual TGeoPatternFinder * MakeCopy(Bool_t reflect=kFALSE)=0
Int_t GetNdiv() const
Double_t GetStep() const
void ClearThreadData() const
virtual Int_t GetByteCount() const
void CreateThreadData(Int_t nthreads)
Create thread data for n threads max.
Class describing scale transformations.
Definition TGeoMatrix.h:253
static TGeoShape * MakeScaledShape(const char *name, TGeoShape *shape, TGeoScale *scale)
Create a scaled shape starting from a non-scaled one.
The shape encapsulating an assembly (union) of volumes.
Base abstract class for all shapes.
Definition TGeoShape.h:25
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const =0
virtual void CreateThreadData(Int_t)
Definition TGeoShape.h:74
Bool_t IsValid() const
Definition TGeoShape.h:151
virtual const char * GetAxisName(Int_t iaxis) const =0
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)=0
virtual Bool_t IsComposite() const
Definition TGeoShape.h:138
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
Bool_t IsRunTimeShape() const
Definition TGeoShape.h:150
virtual void ClearThreadData() const
Definition TGeoShape.h:73
const char * GetPointerName() const
Provide a pointer name containing uid.
void CheckShape(Int_t testNo, Int_t nsamples=10000, Option_t *option="")
Test for shape navigation methods.
virtual Bool_t IsValidBox() const =0
virtual Int_t GetByteCount() const =0
const char * GetName() const override
Get the shape name.
virtual void ComputeBBox()=0
virtual Double_t Capacity() const =0
@ kGeoSavePrimitive
Definition TGeoShape.h:64
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const =0
virtual Bool_t IsAssembly() const
Definition TGeoShape.h:137
Bool_t TestShapeBit(UInt_t f) const
Definition TGeoShape.h:175
Volume assemblies.
Definition TGeoVolume.h:316
static TGeoVolumeAssembly * MakeAssemblyFromVolume(TGeoVolume *vol)
Make a clone of volume VOL but which is an assembly.
~TGeoVolumeAssembly() override
Destructor. The assembly is owner of its "shape".
Int_t GetNextNodeIndex() const override
void CreateThreadData(Int_t nthreads) override
TGeoNode * AddNode(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat=nullptr, Option_t *option="") override
Add a component to the assembly.
void ClearThreadData() const override
TGeoVolume * CloneVolume() const override
Clone this volume.
void AddNodeOverlap(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat, Option_t *option) override
Add an overlapping node - not allowed for assemblies.
std::vector< ThreadData_t * > fThreadData
Definition TGeoVolume.h:331
std::mutex fMutex
Thread vector size.
Definition TGeoVolume.h:333
Int_t fThreadSize
Thread specific data vector.
Definition TGeoVolume.h:332
void SetNextNodeIndex(Int_t index)
Int_t GetCurrentNodeIndex() const override
void SetCurrentNodeIndex(Int_t index)
TGeoVolumeAssembly()
Default constructor.
TGeoVolume * Divide(const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step, Int_t numed=0, Option_t *option="") override
Division makes no sense for assemblies.
ThreadData_t & GetThreadData() const
Volume families.
Definition TGeoVolume.h:266
TGeoVolume * Divide(const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step, Int_t numed=0, Option_t *option="") override
division of multiple volumes
void SetLineColor(Color_t lcolor) override
Set the line color for all components.
Double_t fStart
Definition TGeoVolume.h:273
void AddVolume(TGeoVolume *vol)
Add a volume with valid shape to the list of volumes.
TGeoVolume * MakeCopyVolume(TGeoShape *newshape) override
Make a copy of this volume build a volume with same name, shape and medium.
void SetMedium(TGeoMedium *medium) override
Set medium for a multiple volume.
TGeoNode * AddNode(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat, Option_t *option="") override
Add a new node to the list of nodes.
TGeoVolume * GetVolume(Int_t id) const
Definition TGeoVolume.h:286
void SetVisibility(Bool_t vis=kTRUE) override
Set visibility for all components.
~TGeoVolumeMulti() override
Destructor.
TGeoVolumeMulti * fDivision
Definition TGeoVolume.h:269
TGeoVolumeMulti()
dummy constructor
void SetLineWidth(Width_t lwidth) override
Set the line width for all components.
TGeoShape * GetLastShape() const
Returns the last shape.
void AddNodeOverlap(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat, Option_t *option="") override
Add a new node to the list of nodes, This node is possibly overlapping with other daughters of the vo...
void SetLineStyle(Style_t lstyle) override
Set the line style for all components.
TObjArray * fVolumes
Definition TGeoVolume.h:268
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition TGeoVolume.h:43
Double_t WeightA() const
Analytical computation of the weight.
void AddNodeOffset(TGeoVolume *vol, Int_t copy_no, Double_t offset=0, Option_t *option="")
Add a division node to the list of nodes.
void SetVisContainers(Bool_t flag=kTRUE) override
Set visibility for containers.
virtual void cd(Int_t inode) const
Actualize matrix of node indexed <inode>
virtual void ClearThreadData() const
void SetVisibility(Bool_t vis=kTRUE) override
set visibility of this volume
Bool_t IsVisContainers() const
Definition TGeoVolume.h:157
void SetVoxelFinder(TGeoVoxelFinder *finder)
Definition TGeoVolume.h:243
void RemoveNode(TGeoNode *node)
Remove an existing daughter.
Int_t GetNodeIndex(const TGeoNode *node, Int_t *check_list, Int_t ncheck) const
Get the index of a daughter within check_list by providing the node pointer.
Bool_t Valid() const
Check if the shape of this volume is valid.
Bool_t IsAllInvisible() const
Return TRUE if volume and all daughters are invisible.
Int_t fNtotal
Definition TGeoVolume.h:56
void MakeCopyNodes(const TGeoVolume *other)
make a new list of nodes and copy all nodes of other volume inside
void SetUserExtension(TGeoExtension *ext)
Connect user-defined extension to the volume.
TGeoExtension * GrabFWExtension() const
Get a copy of the framework extension pointer.
void SetNumber(Int_t number)
Definition TGeoVolume.h:245
void ClearNodes()
Definition TGeoVolume.h:95
void SetLineWidth(Width_t lwidth) override
Set the line width.
void Raytrace(Bool_t flag=kTRUE)
Draw this volume with current settings and perform raytracing in the pad.
void RandomRays(Int_t nrays=10000, Double_t startx=0, Double_t starty=0, Double_t startz=0, const char *target_vol=nullptr, Bool_t check_norm=kFALSE)
Random raytracing method.
TGeoVolume()
dummy constructor
TGeoMedium * GetMedium() const
Definition TGeoVolume.h:175
char * GetObjectInfo(Int_t px, Int_t py) const override
Get volume info for the browser.
void Print(Option_t *option="") const override
Print volume info.
void CloneNodesAndConnect(TGeoVolume *newmother) const
Clone the array of nodes.
void SortNodes()
sort nodes by decreasing volume of the bounding box.
Bool_t FindMatrixOfDaughterVolume(TGeoVolume *vol) const
Find a daughter node having VOL as volume and fill TGeoManager::fHMatrix with its global matrix.
void Voxelize(Option_t *option)
build the voxels for this volume
Double_t Capacity() const
Computes the capacity of this [cm^3] as the capacity of its shape.
virtual TGeoVolume * MakeCopyVolume(TGeoShape *newshape)
make a copy of this volume build a volume with same name, shape and medium
void ReplayCreation(const TGeoVolume *other)
Recreate the content of the other volume without pointer copying.
Double_t Weight(Double_t precision=0.01, Option_t *option="va")
Estimate the weight of a volume (in kg) with SIGMA(M)/M better than PRECISION.
Int_t fNumber
option - if any
Definition TGeoVolume.h:55
virtual void CreateThreadData(Int_t nthreads)
virtual Int_t GetByteCount() const
get the total size in bytes for this volume
virtual TGeoNode * AddNode(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat=nullptr, Option_t *option="")
Add a TGeoNode to the list of nodes.
Bool_t OptimizeVoxels()
Perform an extensive sampling to find which type of voxelization is most efficient.
void Browse(TBrowser *b) override
How to browse a volume.
void SetCurrentPoint(Double_t x, Double_t y, Double_t z)
Set the current tracking point.
void Paint(Option_t *option="") override
paint volume
void SetVisOnly(Bool_t flag=kTRUE) override
Set visibility for leaves.
TGeoManager * fGeoManager
Definition TGeoVolume.h:51
TH2F * LegoPlot(Int_t ntheta=20, Double_t themin=0., Double_t themax=180., Int_t nphi=60, Double_t phimin=0., Double_t phimax=360., Double_t rmin=0., Double_t rmax=9999999, Option_t *option="")
Generate a lego plot fot the top volume, according to option.
void Draw(Option_t *option="") override
draw top volume according to option
TGeoVoxelFinder * fVoxels
Definition TGeoVolume.h:50
TGeoMaterial * GetMaterial() const
Definition TGeoVolume.h:174
TGeoExtension * GrabUserExtension() const
Get a copy of the user extension pointer.
@ kVolumeSelected
Definition TGeoVolume.h:73
@ kVolumeImportNodes
Definition TGeoVolume.h:76
Int_t CountNodes(Int_t nlevels=1000, Int_t option=0)
Count total number of subnodes starting from this volume, nlevels down.
void GrabFocus()
Move perspective view focus to this volume.
void UnmarkSaved()
Reset SavePrimitive bits.
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Execute mouse actions on this volume.
virtual TGeoVolume * CloneVolume() const
Clone this volume.
void SetFinder(TGeoPatternFinder *finder)
Definition TGeoVolume.h:244
Int_t GetNdaughters() const
Definition TGeoVolume.h:362
Bool_t IsValid() const
Definition TGeoVolume.h:154
void Grab()
Definition TGeoVolume.h:136
void CheckGeometry(Int_t nrays=1, Double_t startx=0, Double_t starty=0, Double_t startz=0) const
Shoot nrays with random directions from starting point (startx, starty, startz) in the reference fram...
void SelectVolume(Bool_t clear=kFALSE)
Select this volume as matching an arbitrary criteria.
const char * GetPointerName() const
Provide a pointer name containing uid.
static TClass * Class()
TObjArray * GetNodes()
Definition TGeoVolume.h:169
void ClearShape()
Clear the shape of this volume from the list held by the current manager.
void SetFWExtension(TGeoExtension *ext)
Connect framework defined extension to the volume.
void VisibleDaughters(Bool_t vis=kTRUE)
set visibility for daughters
void FindOverlaps() const
loop all nodes marked as overlaps and find overlapping brothers
virtual void AddNodeOverlap(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat=nullptr, Option_t *option="")
Add a TGeoNode to the list of nodes.
TGeoNode * GetNode(const char *name) const
get the pointer to a daughter node
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
compute the closest distance of approach from point px,py to this volume
void RandomPoints(Int_t npoints=1000000, Option_t *option="")
Draw random points in the bounding box of this volume.
void CheckShapes()
check for negative parameters in shapes.
void SetNtotal(Int_t ntotal)
Definition TGeoVolume.h:246
Bool_t GetOptimalVoxels() const
Returns true if cylindrical voxelization is optimal.
TGeoNode * ReplaceNode(TGeoNode *nodeorig, TGeoShape *newshape=nullptr, TGeoMatrix *newpos=nullptr, TGeoMedium *newmed=nullptr)
Replace an existing daughter with a new volume having the same name but possibly a new shape,...
void InvisibleAll(Bool_t flag=kTRUE)
Make volume and each of it daughters (in)visible.
Bool_t IsVisibleDaughters() const
Definition TGeoVolume.h:156
TString fOption
just a hook for now
Definition TGeoVolume.h:54
Int_t GetIndex(const TGeoNode *node) const
get index number for a given daughter
void SetNodes(TObjArray *nodes)
Definition TGeoVolume.h:223
TGeoPatternFinder * GetFinder() const
Definition TGeoVolume.h:177
void PrintVoxels() const
Print the voxels for this volume.
TGeoExtension * fUserExtension
Definition TGeoVolume.h:59
virtual void SetMedium(TGeoMedium *medium)
Definition TGeoVolume.h:242
TGeoVoxelFinder * GetVoxels() const
Getter for optimization structure.
void SetAttVisibility(Bool_t vis)
Definition TGeoVolume.h:233
~TGeoVolume() override
Destructor.
void SetShape(const TGeoShape *shape)
set the shape associated with this volume
static TGeoMedium * DummyMedium()
TObject * fField
pointer to TGeoManager owning this volume
Definition TGeoVolume.h:53
void SetLineColor(Color_t lcolor) override
Set the line color.
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
void CleanAll()
Clean data of the volume.
Bool_t IsTopVolume() const
True if this is the top volume of the geometry.
TGeoMedium * fMedium
Definition TGeoVolume.h:47
TGeoShape * GetShape() const
Definition TGeoVolume.h:190
void InspectMaterial() const
Inspect the material for this volume.
void PrintNodes() const
print nodes
static TGeoMedium * fgDummyMedium
Definition TGeoVolume.h:48
void RegisterYourself(Option_t *option="")
Register the volume and all materials/media/matrices/shapes to the manager.
virtual TGeoVolume * Divide(const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step, Int_t numed=0, Option_t *option="")
Division a la G3.
Bool_t IsRaytracing() const
Check if the painter is currently ray-tracing the content of this volume.
void SaveAs(const char *filename="", Option_t *option="") const override
Save geometry having this as top volume as a C++ macro.
TGeoShape * fShape
Definition TGeoVolume.h:46
void SetField(TObject *field)
Definition TGeoVolume.h:231
static TGeoVolume * Import(const char *filename, const char *name="", Option_t *option="")
Import a volume from a file.
Bool_t IsStyleDefault() const
check if the visibility and attributes are the default ones
Char_t fTransparency
Definition TGeoVolume.h:58
static void CreateDummyMedium()
Create a dummy medium.
TGeoExtension * fFWExtension
Transient user-defined extension to volumes.
Definition TGeoVolume.h:60
void SetAsTopVolume()
Set this volume as the TOP one (the whole geometry starts from here)
Bool_t IsVisLeaves() const
Definition TGeoVolume.h:158
TGeoPatternFinder * fFinder
dummy medium
Definition TGeoVolume.h:49
Int_t Export(const char *filename, const char *name="", Option_t *option="")
Export this volume to a file.
virtual void DrawOnly(Option_t *option="")
draw only this volume
void SetLineStyle(Style_t lstyle) override
Set the line style.
void SetOption(const char *option)
Set the current options (none implemented)
virtual Bool_t IsAssembly() const
Returns true if the volume is an assembly or a scaled assembly.
TGeoVolume * MakeReflectedVolume(const char *newname="") const
Make a copy of this volume which is reflected with respect to XY plane.
TObjArray * fNodes
Definition TGeoVolume.h:45
virtual Bool_t IsVisible() const
Definition TGeoVolume.h:155
void SetVisLeaves(Bool_t flag=kTRUE) override
Set visibility for leaves.
void InspectShape() const
Definition TGeoVolume.h:195
Bool_t IsFolder() const override
Return TRUE if volume contains nodes.
TGeoNode * FindNode(const char *name) const
search a daughter inside the list of nodes
void CheckOverlaps(Double_t ovlp=0.1, Option_t *option="") const
Overlap checking tool.
void SetOverlappingCandidate(Bool_t flag)
Definition TGeoVolume.h:228
Bool_t IsOverlappingCandidate() const
Definition TGeoVolume.h:148
Int_t fRefCount
Definition TGeoVolume.h:57
void Streamer(TBuffer &) override
Stream an object of class TGeoVolume.
void CheckShape(Int_t testNo, Int_t nsamples=10000, Option_t *option="")
Tests for checking the shape navigation algorithms. See TGeoShape::CheckShape()
Finder class handling voxels.
Bool_t IsInvalid() const
void Print(Option_t *option="") const override
Print the voxels.
void SetNeedRebuild(Bool_t flag=kTRUE)
virtual void Voxelize(Option_t *option="")
Voxelize attached volume according to option If the volume is an assembly, make sure the bbox is comp...
virtual void FindOverlaps(Int_t inode) const
create the list of nodes for which the bboxes overlap with inode's bbox
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3038
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:307
Book space in a file, create I/O buffers, to fill them, (un)compress them.
Definition TKey.h:28
virtual const char * GetClassName() const
Definition TKey.h:75
virtual TObject * ReadObj()
To read a TObject* from the file.
Definition TKey.cxx:761
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition TList.cxx:576
void Add(TObject *obj) override
Definition TList.h:81
TMap implements an associative array of (key,value) pairs using a THashTable for efficient retrieval ...
Definition TMap.h:40
void Add(TObject *obj) override
This function may not be used (but we need to provide it since it is a pure virtual in TCollection).
Definition TMap.cxx:54
TObject * GetValue(const char *keyname) const
Returns a pointer to the value associated with keyname as name of the key.
Definition TMap.cxx:236
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:174
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:50
TString fTitle
Definition TNamed.h:33
TString fName
Definition TNamed.h:32
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition TNamed.cxx:150
An array of TObjects.
Definition TObjArray.h:31
Int_t GetEntriesFast() const
Definition TObjArray.h:58
void AddAt(TObject *obj, Int_t idx) override
Add object at position ids.
void Clear(Option_t *option="") override
Remove all objects from the array.
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
virtual void Compress()
Remove empty slots from array.
Int_t GetEntries() const override
Return the number of objects in array (i.e.
void Delete(Option_t *option="") override
Remove all objects from the array AND delete all heap based objects.
TObject * At(Int_t idx) const override
Definition TObjArray.h:164
TObject * Remove(TObject *obj) override
Remove object from array.
TObject * RemoveAt(Int_t idx) override
Remove object at index idx.
TObject * FindObject(const char *name) const override
Find an object in this collection using its name.
void Add(TObject *obj) override
Definition TObjArray.h:68
Mother of all ROOT objects.
Definition TObject.h:41
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:202
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:226
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition TObject.cxx:835
virtual Int_t Write(const char *name=nullptr, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
Definition TObject.cxx:964
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:864
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition TObject.cxx:543
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition TObject.cxx:1099
virtual void SetUniqueID(UInt_t uid)
Set the unique object id.
Definition TObject.cxx:875
void ResetBit(UInt_t f)
Definition TObject.h:201
Sequenceable collection abstract base class.
Basic string class.
Definition TString.h:139
Ssiz_t Length() const
Definition TString.h:417
void ToLower()
Change string to lower-case.
Definition TString.cxx:1182
TSubString Strip(EStripType s=kTrailing, char c=' ') const
Return a substring of self stripped at beginning and/or end.
Definition TString.cxx:1163
const char * Data() const
Definition TString.h:376
Ssiz_t Capacity() const
Definition TString.h:364
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition TString.cxx:2378
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:632
Abstract class for geometry painters.
virtual TGeoVolume * GetTopVolume() const =0
std::ostream & Info()
Definition hadd.cxx:177
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
ThreadData_t()
index of next node to be entered