Logo ROOT  

Reference Guide
 
Loading...
Searching...
No Matches
TGeoNavigator.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Mihaela Gheata 30/05/07
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/** \class TGeoNavigator
13\ingroup Geometry_classes
14
15 Class providing navigation API for TGeo geometries. Several instances are
16allowed for a single geometry.
17A default navigator is provided for any geometry but one may add several
18others for parallel navigation:
19
20~~~ {.cpp}
21TGeoNavigator *navig = new TGeoNavigator(gGeoManager);
22Int_t inav = gGeoManager->AddNavigator(navig);
23gGeoManager->SetCurrentNavigator(inav);
24~~~
25
26.... and then switch back to the default navigator:
27
28~~~ {.cpp}
29gGeoManager->SetCurrentNavigator(0);
30~~~
31
32*/
33
34#include "TGeoNavigator.h"
35
36#include "TGeoManager.h"
37#include "TGeoMatrix.h"
38#include "TGeoNode.h"
39#include "TGeoVolume.h"
40#include "TGeoPatternFinder.h"
41#include "TGeoVoxelFinder.h"
42#include "TMath.h"
43#include "TGeoParallelWorld.h"
44#include "TGeoPhysicalNode.h"
45
46static Double_t gTolerance = TGeoShape::Tolerance();
47const char *kGeoOutsidePath = " ";
48const Int_t kN3 = 3 * sizeof(Double_t);
49
50Bool_t TGeoNavigator::fgUsePWSafetyCaching = kFALSE;
51
52
53////////////////////////////////////////////////////////////////////////////////
54/// Constructor
55
56TGeoNavigator::TGeoNavigator()
57 : fStep(0.),
58 fSafety(0.),
59 fLastSafety(0.),
60 fLastPWSafety(-1.),
61 fThreadId(0),
62 fLevel(0),
63 fNmany(0),
64 fNextDaughterIndex(0),
65 fOverlapSize(0),
66 fOverlapMark(0),
67 fOverlapClusters(nullptr),
68 fSearchOverlaps(kFALSE),
69 fCurrentOverlapping(kFALSE),
70 fStartSafe(kFALSE),
71 fIsEntering(kFALSE),
72 fIsExiting(kFALSE),
73 fIsStepEntering(kFALSE),
74 fIsStepExiting(kFALSE),
75 fIsOutside(kFALSE),
76 fIsOnBoundary(kFALSE),
77 fIsSameLocation(kFALSE),
78 fIsNullStep(kFALSE),
79 fGeometry(nullptr),
80 fCache(nullptr),
81 fCurrentVolume(nullptr),
82 fCurrentNode(nullptr),
83 fTopNode(nullptr),
84 fLastNode(nullptr),
85 fNextNode(nullptr),
86 fForcedNode(nullptr),
87 fBackupState(nullptr),
88 fCurrentMatrix(nullptr),
89 fGlobalMatrix(nullptr),
90 fDivMatrix(nullptr),
91 fPath()
92
93{
94 // dummy constructor
95 for (Int_t i = 0; i < 3; i++) {
96 fNormal[i] = 0.;
97 fCldir[i] = 0.;
98 fCldirChecked[i] = 0.;
99 fPoint[i] = 0.;
100 fDirection[i] = 0.;
101 fLastPoint[i] = 0.;
102 }
103}
104
105////////////////////////////////////////////////////////////////////////////////
106/// Constructor
107
108TGeoNavigator::TGeoNavigator(TGeoManager *geom)
109 : fStep(0.),
110 fSafety(0.),
111 fLastSafety(0.),
112 fLastPWSafety(-1.),
113 fThreadId(0),
114 fLevel(0),
115 fNmany(0),
116 fNextDaughterIndex(-2),
117 fOverlapSize(1000),
118 fOverlapMark(0),
119 fOverlapClusters(nullptr),
120 fSearchOverlaps(kFALSE),
121 fCurrentOverlapping(kFALSE),
122 fStartSafe(kTRUE),
123 fIsEntering(kFALSE),
124 fIsExiting(kFALSE),
125 fIsStepEntering(kFALSE),
126 fIsStepExiting(kFALSE),
127 fIsOutside(kFALSE),
128 fIsOnBoundary(kFALSE),
129 fIsSameLocation(kTRUE),
130 fIsNullStep(kFALSE),
131 fGeometry(geom),
132 fCache(nullptr),
133 fCurrentVolume(nullptr),
134 fCurrentNode(nullptr),
135 fTopNode(nullptr),
136 fLastNode(nullptr),
137 fNextNode(nullptr),
138 fForcedNode(nullptr),
139 fBackupState(nullptr),
140 fCurrentMatrix(nullptr),
141 fGlobalMatrix(nullptr),
142 fDivMatrix(nullptr),
143 fPath()
144
145{
146 // Default constructor.
147 fThreadId = TGeoManager::ThreadId();
148 // printf("Navigator: threadId=%d\n", fThreadId);
149 for (Int_t i = 0; i < 3; i++) {
150 fNormal[i] = 0.;
151 fCldir[i] = 0.;
152 fCldirChecked[i] = 0;
153 fPoint[i] = 0.;
154 fDirection[i] = 0.;
155 fLastPoint[i] = 0.;
156 }
157 fCurrentMatrix = new TGeoHMatrix();
158 fCurrentMatrix->RegisterYourself();
159 fDivMatrix = new TGeoHMatrix();
160 fDivMatrix->RegisterYourself();
161 fOverlapClusters = new Int_t[fOverlapSize];
162 ResetAll();
163}
164
165////////////////////////////////////////////////////////////////////////////////
166/// Destructor.
167
168TGeoNavigator::~TGeoNavigator()
169{
170 if (fCache)
171 delete fCache;
172 if (fBackupState)
173 delete fBackupState;
174 if (fOverlapClusters)
175 delete[] fOverlapClusters;
176}
177
178////////////////////////////////////////////////////////////////////////////////
179/// Builds the cache for physical nodes and global matrices.
180
181void TGeoNavigator::BuildCache(Bool_t /*dummy*/, Bool_t nodeid)
182{
183 static Bool_t first = kTRUE;
184 Int_t verbose = TGeoManager::GetVerboseLevel();
185 Int_t nlevel = fGeometry->GetMaxLevel();
186 if (nlevel <= 0)
187 nlevel = 100;
188 if (!fCache) {
189 if (nlevel == 100) {
190 if (first && verbose > 0)
191 Info("BuildCache", "--- Maximum geometry depth set to 100");
192 } else {
193 if (first && verbose > 0)
194 Info("BuildCache", "--- Maximum geometry depth is %i", nlevel);
195 }
196 // build cache
197 fCache = new TGeoNodeCache(fGeometry->GetTopNode(), nodeid, nlevel + 1);
198 fGlobalMatrix = fCache->GetCurrentMatrix();
199 fBackupState = new TGeoCacheState(nlevel + 1);
200 }
201 first = kFALSE;
202}
203
204////////////////////////////////////////////////////////////////////////////////
205/// Browse the tree of nodes starting from top node according to pathname.
206/// Changes the path accordingly. The path is changed to point to the top node
207/// in case of failure.
208
209Bool_t TGeoNavigator::cd(const char *path)
210{
211 CdTop();
212 if (!path[0])
213 return kTRUE;
214 TString spath = path;
215 TGeoVolume *vol;
216 Int_t length = spath.Length();
217 Int_t ind1 = spath.Index("/");
218 if (ind1 == length - 1)
219 ind1 = -1;
220 Int_t ind2 = 0;
221 Bool_t end = kFALSE;
222 Bool_t first = kTRUE;
223 TString name;
224 TGeoNode *node;
225 while (!end) {
226 ind2 = spath.Index("/", ind1 + 1);
227 if (ind2 < 0 || ind2 == length - 1) {
228 if (ind2 < 0)
229 ind2 = length;
230 end = kTRUE;
231 }
232 name = spath(ind1 + 1, ind2 - ind1 - 1);
233 vol = fCurrentNode->GetVolume();
234 if (first) {
235 first = kFALSE;
236 if (name.BeginsWith(vol->GetName())) {
237 ind1 = ind2;
238 continue;
239 }
240 }
241 node = vol->GetNode(name.Data());
242 if (!node) {
243 Error("cd", "Path %s not valid", path);
244 return kFALSE;
245 }
246 CdDown(vol->GetIndex(node));
247 ind1 = ind2;
248 }
249 return kTRUE;
250}
251
252////////////////////////////////////////////////////////////////////////////////
253/// Check if a geometry path is valid without changing the state of the navigator.
254
255Bool_t TGeoNavigator::CheckPath(const char *path) const
256{
257 if (!path[0])
258 return kTRUE;
259 TGeoNode *crtnode = fGeometry->GetTopNode();
260 TString spath = path;
261 TGeoVolume *vol;
262 Int_t length = spath.Length();
263 Int_t ind1 = spath.Index("/");
264 if (ind1 == length - 1)
265 ind1 = -1;
266 Int_t ind2 = 0;
267 Bool_t end = kFALSE;
268 Bool_t first = kTRUE;
269 TString name;
270 TGeoNode *node;
271 while (!end) {
272 ind2 = spath.Index("/", ind1 + 1);
273 if (ind2 < 0 || ind2 == length - 1) {
274 if (ind2 < 0)
275 ind2 = length;
276 end = kTRUE;
277 }
278 name = spath(ind1 + 1, ind2 - ind1 - 1);
279 vol = crtnode->GetVolume();
280 if (first) {
281 first = kFALSE;
282 if (name.BeginsWith(vol->GetName())) {
283 ind1 = ind2;
284 continue;
285 }
286 }
287 node = vol->GetNode(name.Data());
288 if (!node)
289 return kFALSE;
290 crtnode = node;
291 ind1 = ind2;
292 }
293 return kTRUE;
294}
295
296////////////////////////////////////////////////////////////////////////////////
297/// Change current path to point to the node having this id.
298/// Node id has to be in range : 0 to fNNodes-1 (no check for performance reasons)
299
300void TGeoNavigator::CdNode(Int_t nodeid)
301{
302 if (fCache) {
303 fCache->CdNode(nodeid);
304 fGlobalMatrix = fCache->GetCurrentMatrix();
305 }
306}
307
308////////////////////////////////////////////////////////////////////////////////
309/// Make a daughter of current node current. Can be called only with a valid
310/// daughter index (no check). Updates cache accordingly.
311
312void TGeoNavigator::CdDown(Int_t index)
313{
314 TGeoNode *node = fCurrentNode->GetDaughter(index);
315 Bool_t is_offset = node->IsOffset();
316 if (is_offset)
317 node->cd();
318 else
319 fCurrentOverlapping = node->IsOverlapping();
320 fCache->CdDown(index);
321 fCurrentNode = node;
322 fGlobalMatrix = fCache->GetCurrentMatrix();
323 if (fCurrentOverlapping)
324 fNmany++;
325 fLevel++;
326}
327
328////////////////////////////////////////////////////////////////////////////////
329/// Make a daughter of current node current. Can be called only with a valid
330/// daughter node (no check). Updates cache accordingly.
331
332void TGeoNavigator::CdDown(TGeoNode *node)
333{
334 Bool_t is_offset = node->IsOffset();
335 if (is_offset)
336 node->cd();
337 else
338 fCurrentOverlapping = node->IsOverlapping();
339 fCache->CdDown(node);
340 fCurrentNode = node;
341 fGlobalMatrix = fCache->GetCurrentMatrix();
342 if (fCurrentOverlapping)
343 fNmany++;
344 fLevel++;
345}
346
347////////////////////////////////////////////////////////////////////////////////
348/// Go one level up in geometry. Updates cache accordingly.
349/// Determine the overlapping state of current node.
350
351void TGeoNavigator::CdUp()
352{
353 if (!fLevel || !fCache)
354 return;
355 fLevel--;
356 if (!fLevel) {
357 CdTop();
358 return;
359 }
360 fCache->CdUp();
361 if (fCurrentOverlapping) {
362 fLastNode = fCurrentNode;
363 fNmany--;
364 }
365 fCurrentNode = fCache->GetNode();
366 fGlobalMatrix = fCache->GetCurrentMatrix();
367 if (!fCurrentNode->IsOffset()) {
368 fCurrentOverlapping = fCurrentNode->IsOverlapping();
369 } else {
370 Int_t up = 1;
371 Bool_t offset = kTRUE;
372 TGeoNode *mother = nullptr;
373 while (offset) {
374 mother = GetMother(up++);
375 offset = mother->IsOffset();
376 }
377 fCurrentOverlapping = mother->IsOverlapping();
378 }
379}
380
381////////////////////////////////////////////////////////////////////////////////
382/// Make top level node the current node. Updates the cache accordingly.
383/// Determine the overlapping state of current node.
384
400
401////////////////////////////////////////////////////////////////////////////////
402/// Do a cd to the node found next by FindNextBoundary
403
404void TGeoNavigator::CdNext()
405{
406 if (fNextDaughterIndex == -2 || !fCache)
407 return;
408 if (fNextDaughterIndex == -3) {
409 // Next node is a many - restore it
410 DoRestoreState();
411 fNextDaughterIndex = -2;
412 return;
413 }
414 if (fNextDaughterIndex == -1) {
415 CdUp();
416 while (fCurrentNode->GetVolume()->IsAssembly())
417 CdUp();
418 fNextDaughterIndex--;
419 return;
420 }
421 if (fCurrentNode && fNextDaughterIndex < fCurrentNode->GetNdaughters()) {
422 CdDown(fNextDaughterIndex);
423 Int_t nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
424 while (nextindex >= 0) {
425 CdDown(nextindex);
426 nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
427 }
428 }
429 fNextDaughterIndex = -2;
430}
431
432////////////////////////////////////////////////////////////////////////////////
433/// Fill volume names of current branch into an array.
434
435void TGeoNavigator::GetBranchNames(Int_t *names) const
436{
437 fCache->GetBranchNames(names);
438}
439
440////////////////////////////////////////////////////////////////////////////////
441/// Fill node copy numbers of current branch into an array.
442
447
448////////////////////////////////////////////////////////////////////////////////
449/// Fill node copy numbers of current branch into an array.
450
455
456////////////////////////////////////////////////////////////////////////////////
457/// Cross a division cell. Distance to exit contained in fStep, current node
458/// points to the cell node.
459
460TGeoNode *TGeoNavigator::CrossDivisionCell()
461{
462 TGeoPatternFinder *finder = fCurrentNode->GetFinder();
463 if (!finder) {
464 Fatal("CrossDivisionCell", "Volume has no pattern finder");
465 return nullptr;
466 }
467 // Mark current node and go up to the level of the divided volume
468 TGeoNode *skip = fCurrentNode;
469 CdUp();
470 Double_t point[3], newpoint[3], dir[3];
471 fGlobalMatrix->MasterToLocal(fPoint, newpoint);
472 fGlobalMatrix->MasterToLocalVect(fDirection, dir);
473 // Does step cross a boundary along division axis ?
474 Bool_t onbound = finder->IsOnBoundary(newpoint);
475 if (onbound) {
476 // Work along division axis
477 // Get the starting point
478 point[0] = newpoint[0] - dir[0] * fStep * (1. - gTolerance);
479 point[1] = newpoint[1] - dir[1] * fStep * (1. - gTolerance);
480 point[2] = newpoint[2] - dir[2] * fStep * (1. - gTolerance);
481 // Find which is the next crossed cell.
482 finder->FindNode(point, dir);
483 Int_t inext = finder->GetNext();
484 if (inext < 0) {
485 // step fully exits the division along the division axis
486 // Do step exits in a mother cell ?
487 if (fCurrentNode->IsOffset()) {
488 Double_t dist = fCurrentNode->GetVolume()->GetShape()->DistFromInside(point, dir, 3);
489 // Do step exit also from mother cell ?
490 if (dist < fStep + 2. * gTolerance) {
491 // Step exits mother on its own division axis
492 return CrossDivisionCell();
493 }
494 // We end up here
495 return fCurrentNode;
496 }
497 // Exiting in a non-divided volume
498 while (fCurrentNode->GetVolume()->IsAssembly()) {
499 // Move always to mother for assemblies
500 skip = fCurrentNode;
501 if (!fLevel)
502 break;
503 CdUp();
504 }
505 return CrossBoundaryAndLocate(kFALSE, skip);
506 }
507 // step enters a new cell
508 CdDown(inext + finder->GetDivIndex());
509 skip = fCurrentNode;
510 return CrossBoundaryAndLocate(kTRUE, skip);
511 }
512 // step exits on an axis other than the division axis -> get next slice
513 if (fCurrentNode->IsOffset())
514 return CrossDivisionCell();
515 return CrossBoundaryAndLocate(kFALSE, skip);
516}
517
518////////////////////////////////////////////////////////////////////////////////
519/// Cross next boundary and locate within current node
520/// The current point must be on the boundary of fCurrentNode.
521
522TGeoNode *TGeoNavigator::CrossBoundaryAndLocate(Bool_t downwards, TGeoNode *skipnode)
523{
524 // Extrapolate current point with estimated error.
525 Double_t *tr = fGlobalMatrix->GetTranslation();
526 Double_t trmax = 1. + TMath::Abs(tr[0]) + TMath::Abs(tr[1]) + TMath::Abs(tr[2]);
527 Double_t extra = 100. * (trmax + fStep) * gTolerance;
528 const Int_t idebug = TGeoManager::GetVerboseLevel();
529 TGeoNode *crtstate[10];
530 Int_t level = fLevel + 1;
531 Bool_t samepath = kFALSE;
532 for (Int_t i = 0; i < 10; ++i)
533 crtstate[i] = nullptr;
534
535 if (!downwards) {
536 for (Int_t i = 0; i < fLevel; ++i) {
537 crtstate[i] = GetMother(i);
538 if (i == 9)
539 break;
540 }
541 }
542 fPoint[0] += extra * fDirection[0];
543 fPoint[1] += extra * fDirection[1];
544 fPoint[2] += extra * fDirection[2];
545 TGeoNode *current = SearchNode(downwards, skipnode);
546 fForcedNode = nullptr;
547 fPoint[0] -= extra * fDirection[0];
548 fPoint[1] -= extra * fDirection[1];
549 fPoint[2] -= extra * fDirection[2];
550 if (!current)
551 return nullptr;
552 if (downwards) {
553 Int_t nextindex = current->GetVolume()->GetNextNodeIndex();
554 while (nextindex >= 0) {
555 CdDown(nextindex);
556 current = fCurrentNode;
557 nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
558 }
559 if (idebug > 4) {
560 printf("CrossBoundaryAndLocate: entered %s\n", GetPath());
561 }
562 return current;
563 }
564
565 if (skipnode) {
566 if (current == skipnode) {
567 samepath = kTRUE;
568 if (!downwards) {
569 level = TMath::Min(level, 10);
570 for (Int_t i = 1; i < level; i++) {
571 if (crtstate[i - 1] != GetMother(i)) {
572 samepath = kFALSE;
573 break;
574 }
575 }
576 }
577 }
578 }
579
580 if (samepath || current->GetVolume()->IsAssembly()) {
581 if (!fLevel) {
582 fIsOutside = kTRUE;
583 if (idebug > 4) {
584 printf("CrossBoundaryAndLocate: Exited geometry\n");
585 }
586 return fGeometry->GetCurrentNode();
587 }
588 CdUp();
589 while (fLevel && fCurrentNode->GetVolume()->IsAssembly())
590 CdUp();
591 if (!fLevel && fCurrentNode->GetVolume()->IsAssembly()) {
592 fIsOutside = kTRUE;
593 if (idebug > 4) {
594 printf("CrossBoundaryAndLocate: Exited geometry\n");
595 }
596 if (idebug > 4) {
597 printf("CrossBoundaryAndLocate: entered %s\n", GetPath());
598 }
599 return fCurrentNode;
600 }
601 return fCurrentNode;
602 }
603 if (idebug > 4) {
604 printf("CrossBoundaryAndLocate: entered %s\n", GetPath());
605 }
606 return current;
607}
608
609////////////////////////////////////////////////////////////////////////////////
610/// Find distance to next boundary and store it in fStep. Returns node to which this
611/// boundary belongs. If PATH is specified, compute only distance to the node to which
612/// PATH points. If STEPMAX is specified, compute distance only in case fSafety is smaller
613/// than this value. STEPMAX represent the step to be made imposed by other reasons than
614/// geometry (usually physics processes). Therefore in this case this method provides the
615/// answer to the question : "Is STEPMAX a safe step ?" returning a NULL node and filling
616/// fStep with a big number.
617/// In case frombdr=kTRUE, the isotropic safety is set to zero.
618///
619/// Note : safety distance for the current point is computed ONLY in case STEPMAX is
620/// specified, otherwise users have to call explicitly TGeoManager::Safety() if
621/// they want this computed for the current point.
622
623TGeoNode *TGeoNavigator::FindNextBoundary(Double_t stepmax, const char *path, Bool_t frombdr)
624{
625 // convert current point and direction to local reference
626 Int_t iact = 3;
627 Int_t idebug = TGeoManager::GetVerboseLevel();
628 fNextDaughterIndex = -2;
629 fStep = TGeoShape::Big();
630 fIsStepEntering = kFALSE;
631 fIsStepExiting = kFALSE;
632 fForcedNode = nullptr;
633 Bool_t computeGlobal = kFALSE;
634 fIsOnBoundary = frombdr;
635 fSafety = 0.;
636 TGeoNode *top_node = fGeometry->GetTopNode();
637 TGeoVolume *top_volume = top_node->GetVolume();
638 // If inside an assembly, go logically up in the hierarchy
639 while (fCurrentNode->GetVolume()->IsAssembly() && fLevel)
640 CdUp();
641 if (stepmax < 1E29) {
642 if (stepmax <= 0) {
643 stepmax = -stepmax;
644 computeGlobal = kTRUE;
645 }
646 // if (fLastSafety>0 && IsSamePoint(fPoint[0], fPoint[1], fPoint[2])) fSafety = fLastSafety;
647 fSafety = Safety();
648 // Try to get out easy if proposed step within safe region
649 if (!frombdr && (fSafety > 0) && IsSafeStep(stepmax + gTolerance, fSafety)) {
650 fStep = stepmax;
651 fNextNode = fCurrentNode;
652 return fCurrentNode;
653 }
654 fSafety = TMath::Abs(fSafety);
655 memcpy(fLastPoint, fPoint, kN3);
656 fLastSafety = fSafety;
657 if (fSafety < gTolerance)
658 fIsOnBoundary = kTRUE;
659 else
660 fIsOnBoundary = kFALSE;
661 fStep = stepmax;
662 if (stepmax + gTolerance < fSafety) {
663 fNextNode = fCurrentNode;
664 return fCurrentNode;
665 }
666 }
667 if (computeGlobal)
668 fCurrentMatrix->CopyFrom(fGlobalMatrix);
669 Double_t snext = TGeoShape::Big();
670 Double_t safe;
671 Double_t point[3];
672 Double_t dir[3];
673 if (idebug > 4) {
674 printf("TGeoManager::FindNextBoundary: point=(%19.16f, %19.16f, %19.16f)\n", fPoint[0], fPoint[1], fPoint[2]);
675 printf(" dir= (%19.16f, %19.16f, %19.16f)\n", fDirection[0], fDirection[1],
676 fDirection[2]);
677 printf(" pstep=%9.6g path=%s\n", stepmax, GetPath());
678 }
679 if (path[0]) {
680 PushPath();
681 if (!cd(path)) {
682 PopPath();
683 return nullptr;
684 }
685 if (computeGlobal)
686 fCurrentMatrix->CopyFrom(fGlobalMatrix);
687 fNextNode = fCurrentNode;
688 TGeoVolume *tvol = fCurrentNode->GetVolume();
689 fGlobalMatrix->MasterToLocal(fPoint, &point[0]);
690 fGlobalMatrix->MasterToLocalVect(fDirection, &dir[0]);
691 if (idebug > 4) {
692 printf("=== To path: %s\n", path);
693 printf("=== local to %s: (%19.16f, %19.16f, %19.16f, %19.16f, %19.16f, %19.16f)\n", tvol->GetName(), point[0],
694 point[1], point[2], dir[0], dir[1], dir[2]);
695 }
696 if (tvol->Contains(point)) {
697 if (idebug > 4)
698 printf("=== volume %s contains point\n", tvol->GetName());
699 fStep = tvol->GetShape()->DistFromInside(&point[0], &dir[0], iact, fStep, &safe);
700 } else {
701 fStep = tvol->GetShape()->DistFromOutside(&point[0], &dir[0], iact, fStep, &safe);
702 if (idebug > 4) {
703 printf("=== volume %s does not contain point\n", tvol->GetName());
704 printf("=== distance to path: %g\n", fStep);
705 tvol->InspectShape();
706 if (fStep < 1.E20) {
707 Double_t newpt[3];
708 newpt[0] = point[0] + fStep * dir[0];
709 newpt[1] = point[1] + fStep * dir[1];
710 newpt[2] = point[2] + fStep * dir[2];
711 printf("=== Propagated point: (%19.16f, %19.16f, %19.16f)", newpt[0], newpt[1], newpt[2]);
712 }
713 while (fLevel) {
714 CdUp();
715 tvol = fCurrentNode->GetVolume();
716 fGlobalMatrix->MasterToLocal(fPoint, &point[0]);
717 fGlobalMatrix->MasterToLocalVect(fDirection, &dir[0]);
718 printf("=== local to %s: (%19.16f, %19.16f, %19.16f, %19.16f, %19.16f, %19.16f)\n", tvol->GetName(),
719 point[0], point[1], point[2], dir[0], dir[1], dir[2]);
720 if (tvol->Contains(point)) {
721 printf("=== volume %s contains point\n", tvol->GetName());
722 } else {
723 printf("=== volume %s does not contain point\n", tvol->GetName());
724 snext = tvol->GetShape()->DistFromOutside(&point[0], &dir[0], iact, 1.E30, &safe);
725 }
726 }
727 }
728 }
729 PopPath();
730 return fNextNode;
731 }
732 // compute distance to exit point from current node and the distance to its
733 // closest boundary
734 // if point is outside, just check the top node
735 if (fIsOutside) {
736 snext = top_volume->GetShape()->DistFromOutside(fPoint, fDirection, iact, fStep, &safe);
737 fNextNode = top_node;
738 if (snext < fStep - gTolerance) {
739 fIsStepEntering = kTRUE;
740 fStep = snext;
741 Int_t indnext = fNextNode->GetVolume()->GetNextNodeIndex();
742 fNextDaughterIndex = indnext;
743 while (indnext >= 0) {
744 fNextNode = fNextNode->GetDaughter(indnext);
745 if (computeGlobal)
746 fCurrentMatrix->Multiply(fNextNode->GetMatrix());
747 indnext = fNextNode->GetVolume()->GetNextNodeIndex();
748 }
749 return fNextNode;
750 }
751 return nullptr;
752 }
753 fGlobalMatrix->MasterToLocal(fPoint, &point[0]);
754 fGlobalMatrix->MasterToLocalVect(fDirection, &dir[0]);
755 TGeoVolume *vol = fCurrentNode->GetVolume();
756 if (idebug > 4) {
757 printf(" -> from local=(%19.16f, %19.16f, %19.16f)\n", point[0], point[1], point[2]);
758 printf(" ldir =(%19.16f, %19.16f, %19.16f)\n", dir[0], dir[1], dir[2]);
759 }
760 // find distance to exiting current node
761 snext = vol->GetShape()->DistFromInside(&point[0], &dir[0], iact, fStep, &safe);
762 if (idebug > 4) {
763 printf(" exiting %s shape %s at snext=%g\n", vol->GetName(), vol->GetShape()->ClassName(), snext);
764 }
765 if (snext < fStep - gTolerance) {
766 fNextNode = fCurrentNode;
767 fNextDaughterIndex = -1;
768 fIsStepExiting = kTRUE;
769 fStep = snext;
770 fIsStepEntering = kFALSE;
771 if (fStep < 1E-6)
772 return fCurrentNode;
773 }
774 fNextNode = (fStep < 1E20) ? fCurrentNode : nullptr;
775 // Find next daughter boundary for the current volume
776 Int_t idaughter = -1;
777 FindNextDaughterBoundary(point, dir, idaughter, computeGlobal);
778 if (idaughter >= 0)
779 fNextDaughterIndex = idaughter;
780 TGeoNode *current = nullptr;
781 TGeoNode *dnode = nullptr;
782 TGeoVolume *mother = nullptr;
783 // if we are in an overlapping node, check also the mother(s)
784 if (fNmany) {
785 Double_t mothpt[3];
786 Double_t vecpt[3];
787 Double_t dpt[3], dvec[3];
788 Int_t novlps;
789 Int_t idovlp = -1;
790 Int_t safelevel = GetSafeLevel();
791 PushPath(safelevel + 1);
792 while (fCurrentOverlapping) {
793 Int_t *ovlps = fCurrentNode->GetOverlaps(novlps);
794 CdUp();
795 mother = fCurrentNode->GetVolume();
796 fGlobalMatrix->MasterToLocal(fPoint, &mothpt[0]);
797 fGlobalMatrix->MasterToLocalVect(fDirection, &vecpt[0]);
798 // check distance to out
799 snext = TGeoShape::Big();
800 if (!mother->IsAssembly())
801 snext = mother->GetShape()->DistFromInside(&mothpt[0], &vecpt[0], iact, fStep, &safe);
802 if (snext < fStep - gTolerance) {
803 fIsStepExiting = kTRUE;
804 fIsStepEntering = kFALSE;
805 fStep = snext;
806 if (computeGlobal)
807 fCurrentMatrix->CopyFrom(fGlobalMatrix);
808 fNextNode = fCurrentNode;
809 fNextDaughterIndex = -3;
810 DoBackupState();
811 }
812 // check overlapping nodes
813 for (Int_t i = 0; i < novlps; i++) {
814 current = mother->GetNode(ovlps[i]);
815 if (!current->IsOverlapping()) {
816 current->cd();
817 current->MasterToLocal(&mothpt[0], &dpt[0]);
818 current->MasterToLocalVect(&vecpt[0], &dvec[0]);
819 // Current point may be inside the other node - geometry error that we ignore
820 snext = current->GetVolume()->GetShape()->DistFromOutside(&dpt[0], &dvec[0], iact, fStep, &safe);
821 if (snext < fStep - gTolerance) {
822 if (computeGlobal) {
823 fCurrentMatrix->CopyFrom(fGlobalMatrix);
824 fCurrentMatrix->Multiply(current->GetMatrix());
825 }
826 fIsStepExiting = kTRUE;
827 fIsStepEntering = kFALSE;
828 fStep = snext;
829 fNextNode = current;
830 fNextDaughterIndex = -3;
831 CdDown(ovlps[i]);
832 DoBackupState();
833 CdUp();
834 }
835 } else {
836 // another many - check if point is in or out
837 current->cd();
838 current->MasterToLocal(&mothpt[0], &dpt[0]);
839 current->MasterToLocalVect(&vecpt[0], &dvec[0]);
840 if (current->GetVolume()->Contains(dpt)) {
841 if (current->GetVolume()->GetNdaughters()) {
842 CdDown(ovlps[i]);
843 fIsStepEntering = kFALSE;
844 fIsStepExiting = kTRUE;
845 dnode = FindNextDaughterBoundary(dpt, dvec, idovlp, computeGlobal);
846 if (dnode) {
847 if (computeGlobal) {
848 fCurrentMatrix->CopyFrom(fGlobalMatrix);
849 fCurrentMatrix->Multiply(dnode->GetMatrix());
850 }
851 fNextNode = dnode;
852 fNextDaughterIndex = -3;
853 CdDown(idovlp);
854 Int_t indnext = fCurrentNode->GetVolume()->GetNextNodeIndex();
855 Int_t iup = 0;
856 while (indnext >= 0) {
857 CdDown(indnext);
858 iup++;
859 indnext = fCurrentNode->GetVolume()->GetNextNodeIndex();
860 }
861 DoBackupState();
862 while (iup > 0) {
863 CdUp();
864 iup--;
865 }
866 CdUp();
867 }
868 CdUp();
869 }
870 } else {
871 snext = current->GetVolume()->GetShape()->DistFromOutside(&dpt[0], &dvec[0], iact, fStep, &safe);
872 if (snext < fStep - gTolerance) {
873 if (computeGlobal) {
874 fCurrentMatrix->CopyFrom(fGlobalMatrix);
875 fCurrentMatrix->Multiply(current->GetMatrix());
876 }
877 fIsStepExiting = kTRUE;
878 fIsStepEntering = kFALSE;
879 fStep = snext;
880 fNextNode = current;
881 fNextDaughterIndex = -3;
882 CdDown(ovlps[i]);
883 DoBackupState();
884 CdUp();
885 }
886 }
887 }
888 }
889 }
890 // Now we are in a non-overlapping node
891 if (fNmany) {
892 // We have overlaps up in the branch, check distance to exit
893 Int_t up = 1;
894 Int_t imother;
895 Int_t nmany = fNmany;
896 Bool_t ovlp = kFALSE;
897 Bool_t nextovlp = kFALSE;
898 Bool_t offset = kFALSE;
899 TGeoNode *currentnode = fCurrentNode;
900 TGeoNode *mothernode, *mup;
901 TGeoHMatrix *matrix;
902 while (nmany) {
903 mothernode = GetMother(up);
904 if (!mothernode) {
905 Fatal("FindNextBoundary", "Cannot find mother node");
906 return nullptr;
907 }
908 mup = mothernode;
909 imother = up + 1;
910 offset = kFALSE;
911 while (mup->IsOffset()) {
912 mup = GetMother(imother++);
913 offset = kTRUE;
914 }
915 nextovlp = mup->IsOverlapping();
916 if (offset) {
917 mothernode = mup;
918 if (nextovlp)
919 nmany -= imother - up;
920 up = imother - 1;
921 } else {
922 if (ovlp)
923 nmany--;
924 }
925 if (ovlp || nextovlp) {
926 matrix = GetMotherMatrix(up);
927 if (!matrix) {
928 Fatal("FindNextBoundary", "Cannot find mother matrix");
929 return nullptr;
930 }
931 matrix->MasterToLocal(fPoint, dpt);
932 matrix->MasterToLocalVect(fDirection, dvec);
933 // If the node below this mother is a MANY node, its overlap
934 // candidates have to be checked even when the current node is
935 // one of its ONLY descendants.
936 if (ovlp && !offset) {
937 Int_t *ovlps = currentnode->GetOverlaps(novlps);
938 for (Int_t i = 0; i < novlps; i++) {
939 current = mothernode->GetVolume()->GetNode(ovlps[i]);
940 if (current->IsOverlapping())
941 continue;
942 current->cd();
943 current->MasterToLocal(dpt, mothpt);
944 current->MasterToLocalVect(dvec, vecpt);
945 snext = current->GetVolume()->GetShape()->DistFromOutside(mothpt, vecpt, iact, fStep, &safe);
946 if (snext < fStep - gTolerance) {
947 if (computeGlobal) {
948 fCurrentMatrix->CopyFrom(matrix);
949 fCurrentMatrix->Multiply(current->GetMatrix());
950 }
951 fIsStepExiting = kTRUE;
952 fIsStepEntering = kFALSE;
953 fStep = snext;
954 fNextNode = current;
955 fNextDaughterIndex = -3;
956 PushPath();
957 Int_t iup = up;
958 while (iup--)
959 CdUp();
960 CdDown(ovlps[i]);
961 DoBackupState();
962 PopPath();
963 }
964 }
965 }
966 snext = TGeoShape::Big();
967 if (!mothernode->GetVolume()->IsAssembly())
968 snext = mothernode->GetVolume()->GetShape()->DistFromInside(dpt, dvec, iact, fStep);
969 if (snext < fStep - gTolerance) {
970 fIsStepEntering = kFALSE;
971 fIsStepExiting = kTRUE;
972 fStep = snext;
973 fNextNode = mothernode;
974 fNextDaughterIndex = -3;
975 if (computeGlobal)
976 fCurrentMatrix->CopyFrom(matrix);
977 while (up--)
978 CdUp();
979 DoBackupState();
980 up = 1;
981 currentnode = fCurrentNode;
982 ovlp = currentnode->IsOverlapping();
983 continue;
984 }
985 }
986 currentnode = mothernode;
987 ovlp = nextovlp;
988 up++;
989 }
990 }
991 PopPath();
992 }
993 // Compute now the distance in case we have a parallel world
994 Double_t parstep = TGeoShape::Big();
995 if (fGeometry->IsParallelWorldNav()) {
996 // printf("path: %s next node %s at %g\n", GetPath(), fNextNode->GetName(), fStep);
997 TGeoPhysicalNode *pnode = fGeometry->GetParallelWorld()->FindNextBoundary(fPoint, fDirection, parstep, fStep);
998 if (pnode) {
999 // A boundary is hit at less than fPStep
1000 fStep = parstep;
1001 fNextNode = pnode->GetNode();
1002 fNextDaughterIndex = -2; // No way to store it for CdNext
1003 fIsStepEntering = kTRUE;
1004 fIsStepExiting = kFALSE;
1005 Int_t nextindex = fNextNode->GetVolume()->GetNextNodeIndex();
1006 while (nextindex >= 0) {
1007 fNextNode = fNextNode->GetDaughter(nextindex);
1008 nextindex = fNextNode->GetVolume()->GetNextNodeIndex();
1009 }
1010 }
1011 }
1012 return fNextNode;
1013}
1014
1015////////////////////////////////////////////////////////////////////////////////
1016/// Computes as fStep the distance to next daughter of the current volume.
1017/// The point and direction must be converted in the coordinate system of the current volume.
1018/// The proposed step limit is fStep.
1019
1020TGeoNode *TGeoNavigator::FindNextDaughterBoundary(Double_t *point, Double_t *dir, Int_t &idaughter, Bool_t compmatrix)
1021{
1022 Double_t snext = TGeoShape::Big();
1023 Int_t idebug = TGeoManager::GetVerboseLevel();
1024 idaughter = -1; // nothing crossed
1025 TGeoNode *nodefound = nullptr;
1026 // Get number of daughters. If no daughters we are done.
1027
1028 TGeoVolume *vol = fCurrentNode->GetVolume();
1029 Int_t nd = vol->GetNdaughters();
1030 if (!nd)
1031 return nullptr; // No daughter
1032 if (fGeometry->IsActivityEnabled() && !vol->IsActiveDaughters())
1033 return nullptr;
1034 Double_t lpoint[3], ldir[3];
1035 TGeoNode *current = nullptr;
1036 Int_t i = 0;
1037 // if current volume is divided, we are in the non-divided region. We
1038 // check first if we are inside a cell in which case compute distance to next cell
1039 TGeoPatternFinder *finder = vol->GetFinder();
1040 if (finder) {
1041 Int_t ifirst = finder->GetDivIndex();
1042 Int_t ilast = ifirst + finder->GetNdiv() - 1;
1043 current = finder->FindNode(point);
1044 if (current) {
1045 // Point inside a cell: find distance to next cell
1046 Int_t index = current->GetIndex();
1047 if ((index - 1) >= ifirst)
1048 ifirst = index - 1;
1049 else
1050 ifirst = -1;
1051 if ((index + 1) <= ilast)
1052 ilast = index + 1;
1053 else
1054 ilast = -1;
1055 }
1056 if (ifirst >= 0) {
1057 current = vol->GetNode(ifirst);
1058 current->cd();
1059 current->MasterToLocal(&point[0], lpoint);
1060 current->MasterToLocalVect(&dir[0], ldir);
1061 snext = current->GetVolume()->GetShape()->DistFromOutside(lpoint, ldir, 3, fStep);
1062 if (snext < fStep - gTolerance) {
1063 if (compmatrix) {
1064 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1065 fCurrentMatrix->Multiply(current->GetMatrix());
1066 }
1067 fIsStepExiting = kFALSE;
1068 fIsStepEntering = kTRUE;
1069 fStep = snext;
1070 fNextNode = current;
1071 nodefound = current;
1072 idaughter = ifirst;
1073 }
1074 }
1075 if (ilast == ifirst)
1076 return nodefound;
1077 if (ilast >= 0) {
1078 current = vol->GetNode(ilast);
1079 current->cd();
1080 current->MasterToLocal(&point[0], lpoint);
1081 current->MasterToLocalVect(&dir[0], ldir);
1082 snext = current->GetVolume()->GetShape()->DistFromOutside(lpoint, ldir, 3, fStep);
1083 if (snext < fStep - gTolerance) {
1084 if (compmatrix) {
1085 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1086 fCurrentMatrix->Multiply(current->GetMatrix());
1087 }
1088 fIsStepExiting = kFALSE;
1089 fIsStepEntering = kTRUE;
1090 fStep = snext;
1091 fNextNode = current;
1092 nodefound = current;
1093 idaughter = ilast;
1094 }
1095 }
1096 return nodefound;
1097 }
1098 // if only few daughters, check all and exit
1099 TGeoVoxelFinder *voxels = vol->GetVoxels();
1100 Int_t indnext;
1101 if (idebug > 4)
1102 printf(" Checking distance to %d daughters...\n", nd);
1103 if (nd < 5 || !voxels) {
1104 for (i = 0; i < nd; i++) {
1105 current = vol->GetNode(i);
1106 if (fGeometry->IsActivityEnabled() && !current->GetVolume()->IsActive())
1107 continue;
1108 current->cd();
1109 if (voxels && voxels->IsSafeVoxel(point, i, fStep))
1110 continue;
1111 current->MasterToLocal(point, lpoint);
1112 current->MasterToLocalVect(dir, ldir);
1113 if (current->IsOverlapping() && current->GetVolume()->Contains(lpoint) &&
1114 current->GetVolume()->GetShape()->Safety(lpoint, kTRUE) > gTolerance)
1115 continue;
1116 snext = current->GetVolume()->GetShape()->DistFromOutside(lpoint, ldir, 3, fStep);
1117 if (snext < fStep - gTolerance) {
1118 if (idebug > 4) {
1119 printf(" -> from local=(%19.16f, %19.16f, %19.16f)\n", lpoint[0], lpoint[1], lpoint[2]);
1120 printf(" ldir =(%19.16f, %19.16f, %19.16f)\n", ldir[0], ldir[1], ldir[2]);
1121 printf(" -> to: %s shape %s snext=%g\n", current->GetName(),
1122 current->GetVolume()->GetShape()->ClassName(), snext);
1123 }
1124 indnext = current->GetVolume()->GetNextNodeIndex();
1125 if (compmatrix) {
1126 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1127 fCurrentMatrix->Multiply(current->GetMatrix());
1128 }
1129 fIsStepExiting = kFALSE;
1130 fIsStepEntering = kTRUE;
1131 fStep = snext;
1132 fNextNode = current;
1133 nodefound = fNextNode;
1134 idaughter = i;
1135 while (indnext >= 0) {
1136 current = current->GetDaughter(indnext);
1137 if (compmatrix)
1138 fCurrentMatrix->Multiply(current->GetMatrix());
1139 fNextNode = current;
1140 nodefound = current;
1141 indnext = current->GetVolume()->GetNextNodeIndex();
1142 }
1143 }
1144 }
1145 if (vol->IsAssembly())
1146 ((TGeoVolumeAssembly *)vol)->SetNextNodeIndex(idaughter);
1147 return nodefound;
1148 }
1149 // if current volume is voxelized, first get current voxel
1150 Int_t ncheck = 0;
1151 Int_t sumchecked = 0;
1152 Int_t *vlist = nullptr;
1153 TGeoStateInfo &info = *fCache->GetInfo();
1154 voxels->SortCrossedVoxels(point, dir, info);
1155 while ((sumchecked < nd) && (vlist = voxels->GetNextVoxel(point, dir, ncheck, info))) {
1156 for (i = 0; i < ncheck; i++) {
1157 current = vol->GetNode(vlist[i]);
1158 if (fGeometry->IsActivityEnabled() && !current->GetVolume()->IsActive())
1159 continue;
1160 current->cd();
1161 current->MasterToLocal(point, lpoint);
1162 current->MasterToLocalVect(dir, ldir);
1163 if (current->IsOverlapping() && current->GetVolume()->Contains(lpoint) &&
1164 current->GetVolume()->GetShape()->Safety(lpoint, kTRUE) > gTolerance)
1165 continue;
1166 snext = current->GetVolume()->GetShape()->DistFromOutside(lpoint, ldir, 3, fStep);
1167 sumchecked++;
1168 // printf("checked %d from %d : snext=%g\n", sumchecked, nd, snext);
1169 if (snext < fStep - gTolerance) {
1170 if (idebug > 4) {
1171 printf(" -> from local=(%19.16f, %19.16f, %19.16f)\n", lpoint[0], lpoint[1], lpoint[2]);
1172 printf(" ldir =(%19.16f, %19.16f, %19.16f)\n", ldir[0], ldir[1], ldir[2]);
1173 printf(" -> to: %s shape %s snext=%g\n", current->GetName(),
1174 current->GetVolume()->GetShape()->ClassName(), snext);
1175 }
1176 indnext = current->GetVolume()->GetNextNodeIndex();
1177 if (compmatrix) {
1178 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1179 fCurrentMatrix->Multiply(current->GetMatrix());
1180 }
1181 fIsStepExiting = kFALSE;
1182 fIsStepEntering = kTRUE;
1183 fStep = snext;
1184 fNextNode = current;
1185 nodefound = fNextNode;
1186 idaughter = vlist[i];
1187 while (indnext >= 0) {
1188 current = current->GetDaughter(indnext);
1189 if (compmatrix)
1190 fCurrentMatrix->Multiply(current->GetMatrix());
1191 fNextNode = current;
1192 nodefound = current;
1193 indnext = current->GetVolume()->GetNextNodeIndex();
1194 }
1195 }
1196 }
1197 }
1198 fCache->ReleaseInfo();
1199 if (vol->IsAssembly())
1200 ((TGeoVolumeAssembly *)vol)->SetNextNodeIndex(idaughter);
1201 return nodefound;
1202}
1203
1204////////////////////////////////////////////////////////////////////////////////
1205/// Compute distance to next boundary within STEPMAX. If no boundary is found,
1206/// propagate current point along current direction with fStep=STEPMAX. Otherwise
1207/// propagate with fStep=SNEXT (distance to boundary) and locate/return the next
1208/// node.
1209
1210TGeoNode *TGeoNavigator::FindNextBoundaryAndStep(Double_t stepmax, Bool_t compsafe)
1211{
1212 Int_t iact = 3;
1213 Int_t idebug = TGeoManager::GetVerboseLevel();
1214 Int_t nextindex;
1215 Bool_t is_assembly;
1216 fForcedNode = nullptr;
1217 fIsStepExiting = kFALSE;
1218 TGeoNode *skip;
1219 fIsStepEntering = kFALSE;
1220 fStep = stepmax;
1221 Double_t snext = TGeoShape::Big();
1222 // If inside an assembly, go logically up in the hierarchy
1223 while (fCurrentNode->GetVolume()->IsAssembly() && fLevel)
1224 CdUp();
1225 if (compsafe) {
1226 // Try to get out easy if proposed step within safe region
1227 fIsOnBoundary = kFALSE;
1228 if (IsSafeStep(stepmax + gTolerance, fSafety)) {
1229 fPoint[0] += stepmax * fDirection[0];
1230 fPoint[1] += stepmax * fDirection[1];
1231 fPoint[2] += stepmax * fDirection[2];
1232 return fCurrentNode;
1233 }
1234 Safety();
1235 fLastSafety = fSafety;
1236 memcpy(fLastPoint, fPoint, kN3);
1237 // If proposed step less than safety, nothing to check
1238 if (fSafety > stepmax + gTolerance) {
1239 fPoint[0] += stepmax * fDirection[0];
1240 fPoint[1] += stepmax * fDirection[1];
1241 fPoint[2] += stepmax * fDirection[2];
1242 return fCurrentNode;
1243 }
1244 }
1245 Double_t extra = (fIsOnBoundary) ? gTolerance : 0.0;
1246 fIsOnBoundary = kFALSE;
1247 fPoint[0] += extra * fDirection[0];
1248 fPoint[1] += extra * fDirection[1];
1249 fPoint[2] += extra * fDirection[2];
1250 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1251 if (idebug > 4) {
1252 printf("TGeoManager::FindNextBAndStep: point=(%19.16f, %19.16f, %19.16f)\n", fPoint[0], fPoint[1], fPoint[2]);
1253 printf(" dir= (%19.16f, %19.16f, %19.16f)\n", fDirection[0], fDirection[1],
1254 fDirection[2]);
1255 printf(" pstep=%9.6g path=%s\n", stepmax, GetPath());
1256 }
1257
1258 if (fIsOutside) {
1259 snext = fGeometry->GetTopVolume()->GetShape()->DistFromOutside(fPoint, fDirection, iact, fStep);
1260 if (snext < fStep - gTolerance) {
1261 if (snext <= 0) {
1262 snext = 0.0;
1263 fStep = snext;
1264 fPoint[0] -= extra * fDirection[0];
1265 fPoint[1] -= extra * fDirection[1];
1266 fPoint[2] -= extra * fDirection[2];
1267 } else {
1268 fStep = snext + extra;
1269 }
1270 fIsStepEntering = kTRUE;
1271 fNextNode = fGeometry->GetTopNode();
1272 nextindex = fNextNode->GetVolume()->GetNextNodeIndex();
1273 while (nextindex >= 0) {
1274 CdDown(nextindex);
1275 fNextNode = fCurrentNode;
1276 nextindex = fNextNode->GetVolume()->GetNextNodeIndex();
1277 if (nextindex < 0)
1278 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1279 }
1280 // Update global point
1281 fPoint[0] += snext * fDirection[0];
1282 fPoint[1] += snext * fDirection[1];
1283 fPoint[2] += snext * fDirection[2];
1284 fIsOnBoundary = kTRUE;
1285 fIsOutside = kFALSE;
1286 fForcedNode = fCurrentNode;
1287 return CrossBoundaryAndLocate(kTRUE, fCurrentNode);
1288 }
1289 if (snext < TGeoShape::Big()) {
1290 // New point still outside, but the top node is reachable
1291 fNextNode = fGeometry->GetTopNode();
1292 fPoint[0] += (fStep - extra) * fDirection[0];
1293 fPoint[1] += (fStep - extra) * fDirection[1];
1294 fPoint[2] += (fStep - extra) * fDirection[2];
1295 return fNextNode;
1296 }
1297 // top node not reachable from current point/direction
1298 fNextNode = nullptr;
1299 fIsOnBoundary = kFALSE;
1300 return nullptr;
1301 }
1302 Double_t point[3], dir[3];
1303 Int_t icrossed = -2;
1304 fGlobalMatrix->MasterToLocal(fPoint, &point[0]);
1305 fGlobalMatrix->MasterToLocalVect(fDirection, &dir[0]);
1306 TGeoVolume *vol = fCurrentNode->GetVolume();
1307 // find distance to exiting current node
1308 if (idebug > 4) {
1309 printf(" -> from local=(%19.16f, %19.16f, %19.16f)\n", point[0], point[1], point[2]);
1310 printf(" ldir =(%19.16f, %19.16f, %19.16f)\n", dir[0], dir[1], dir[2]);
1311 }
1312 // find distance to exiting current node
1313 snext = vol->GetShape()->DistFromInside(point, dir, iact, fStep);
1314 if (idebug > 4) {
1315 printf(" exiting %s shape %s at snext=%g\n", vol->GetName(), vol->GetShape()->ClassName(), snext);
1316 }
1317 fNextNode = fCurrentNode;
1318 if (snext <= gTolerance) {
1319 // Current point on the boundary while track exiting
1320 snext = gTolerance;
1321 fStep = snext;
1322 fIsOnBoundary = kTRUE;
1323 fIsStepEntering = kFALSE;
1324 fIsStepExiting = kTRUE;
1325 skip = fCurrentNode;
1326 fPoint[0] += fStep * fDirection[0];
1327 fPoint[1] += fStep * fDirection[1];
1328 fPoint[2] += fStep * fDirection[2];
1329 is_assembly = fCurrentNode->GetVolume()->IsAssembly();
1330 if (!fLevel && !is_assembly) {
1331 fIsOutside = kTRUE;
1332 return nullptr;
1333 }
1334 if (fCurrentNode->IsOffset())
1335 return CrossDivisionCell();
1336 if (fLevel)
1337 CdUp();
1338 else
1339 skip = nullptr;
1340 return CrossBoundaryAndLocate(kFALSE, skip);
1341 }
1342
1343 if (snext < fStep - gTolerance) {
1344 // Currently the minimum step chosen is the exiting one
1345 icrossed = -1;
1346 fStep = snext;
1347 fIsStepEntering = kFALSE;
1348 fIsStepExiting = kTRUE;
1349 }
1350 // Find next daughter boundary for the current volume
1351 Int_t idaughter = -1;
1352 TGeoNode *crossed = FindNextDaughterBoundary(point, dir, idaughter, kTRUE);
1353 if (crossed) {
1354 fIsStepExiting = kFALSE;
1355 icrossed = idaughter;
1356 fIsStepEntering = kTRUE;
1357 }
1358 TGeoNode *current = nullptr;
1359 TGeoNode *dnode = nullptr;
1360 TGeoVolume *mother = nullptr;
1361 // if we are in an overlapping node, check also the mother(s)
1362 if (fNmany) {
1363 Double_t mothpt[3];
1364 Double_t vecpt[3];
1365 Double_t dpt[3], dvec[3];
1366 Int_t novlps;
1367 Int_t safelevel = GetSafeLevel();
1368 PushPath(safelevel + 1);
1369 while (fCurrentOverlapping) {
1370 Int_t *ovlps = fCurrentNode->GetOverlaps(novlps);
1371 CdUp();
1372 mother = fCurrentNode->GetVolume();
1373 fGlobalMatrix->MasterToLocal(fPoint, &mothpt[0]);
1374 fGlobalMatrix->MasterToLocalVect(fDirection, &vecpt[0]);
1375 // check distance to out
1376 snext = TGeoShape::Big();
1377 if (!mother->IsAssembly())
1378 snext = mother->GetShape()->DistFromInside(mothpt, vecpt, iact, fStep);
1379 if (snext < fStep - gTolerance) {
1380 // exiting mother first (extrusion)
1381 icrossed = -1;
1382 PopDummy();
1383 PushPath(safelevel + 1);
1384 fIsStepEntering = kFALSE;
1385 fIsStepExiting = kTRUE;
1386 fStep = snext;
1387 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1388 fNextNode = fCurrentNode;
1389 }
1390 // check overlapping nodes
1391 for (Int_t i = 0; i < novlps; i++) {
1392 current = mother->GetNode(ovlps[i]);
1393 if (!current->IsOverlapping()) {
1394 current->cd();
1395 current->MasterToLocal(&mothpt[0], &dpt[0]);
1396 current->MasterToLocalVect(&vecpt[0], &dvec[0]);
1397 snext = current->GetVolume()->GetShape()->DistFromOutside(dpt, dvec, iact, fStep);
1398 if (snext < fStep - gTolerance) {
1399 PopDummy();
1400 PushPath(safelevel + 1);
1401 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1402 fCurrentMatrix->Multiply(current->GetMatrix());
1403 fIsStepEntering = kFALSE;
1404 fIsStepExiting = kTRUE;
1405 icrossed = ovlps[i];
1406 fStep = snext;
1407 fNextNode = current;
1408 }
1409 } else {
1410 // another many - check if point is in or out
1411 current->cd();
1412 current->MasterToLocal(&mothpt[0], &dpt[0]);
1413 current->MasterToLocalVect(&vecpt[0], &dvec[0]);
1414 if (current->GetVolume()->Contains(dpt)) {
1415 if (current->GetVolume()->GetNdaughters()) {
1416 CdDown(ovlps[i]);
1417 dnode = FindNextDaughterBoundary(dpt, dvec, idaughter, kFALSE);
1418 if (dnode) {
1419 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1420 fCurrentMatrix->Multiply(dnode->GetMatrix());
1421 icrossed = idaughter;
1422 PopDummy();
1423 PushPath(safelevel + 1);
1424 fIsStepEntering = kFALSE;
1425 fIsStepExiting = kTRUE;
1426 fNextNode = dnode;
1427 }
1428 CdUp();
1429 }
1430 } else {
1431 snext = current->GetVolume()->GetShape()->DistFromOutside(dpt, dvec, iact, fStep);
1432 if (snext < fStep - gTolerance) {
1433 fCurrentMatrix->CopyFrom(fGlobalMatrix);
1434 fCurrentMatrix->Multiply(current->GetMatrix());
1435 fIsStepEntering = kFALSE;
1436 fIsStepExiting = kTRUE;
1437 fStep = snext;
1438 fNextNode = current;
1439 icrossed = ovlps[i];
1440 PopDummy();
1441 PushPath(safelevel + 1);
1442 }
1443 }
1444 }
1445 }
1446 }
1447 // Now we are in a non-overlapping node
1448 if (fNmany) {
1449 // We have overlaps up in the branch, check distance to exit
1450 Int_t up = 1;
1451 Int_t imother;
1452 Int_t nmany = fNmany;
1453 Bool_t ovlp = kFALSE;
1454 Bool_t nextovlp = kFALSE;
1455 Bool_t offset = kFALSE;
1456 TGeoNode *currentnode = fCurrentNode;
1457 TGeoNode *mothernode, *mup;
1458 TGeoHMatrix *matrix;
1459 while (nmany) {
1460 mothernode = GetMother(up);
1461 mup = mothernode;
1462 imother = up + 1;
1463 offset = kFALSE;
1464 while (mup->IsOffset()) {
1465 mup = GetMother(imother++);
1466 offset = kTRUE;
1467 }
1468 nextovlp = mup->IsOverlapping();
1469 if (offset) {
1470 mothernode = mup;
1471 if (nextovlp)
1472 nmany -= imother - up;
1473 up = imother - 1;
1474 } else {
1475 if (ovlp)
1476 nmany--;
1477 }
1478 if (ovlp || nextovlp) {
1479 matrix = GetMotherMatrix(up);
1480 matrix->MasterToLocal(fPoint, dpt);
1481 matrix->MasterToLocalVect(fDirection, dvec);
1482 /*
1483 * The live leaf can be ONLY while currentnode is a MANY ancestor.
1484 * Check its ONLY overlap candidates in mothernode. If one wins,
1485 * commit the climb to mothernode and restart because the path,
1486 * MANY count, and coordinate frame have changed.
1487 */
1488 Bool_t restart = kFALSE;
1489 if (ovlp && !offset) {
1490 Int_t *ovlps = currentnode->GetOverlaps(novlps);
1491 Int_t icandidate = -1;
1492 TGeoNode *candidate = nullptr;
1493 for (Int_t i = 0; i < novlps; i++) {
1494 current = mothernode->GetVolume()->GetNode(ovlps[i]);
1495 if (current->IsOverlapping())
1496 continue;
1497 current->cd();
1498 current->MasterToLocal(dpt, mothpt);
1499 current->MasterToLocalVect(dvec, vecpt);
1500 snext = current->GetVolume()->GetShape()->DistFromOutside(mothpt, vecpt, iact, fStep);
1501 if (snext < fStep - gTolerance) {
1502 fCurrentMatrix->CopyFrom(matrix);
1503 fCurrentMatrix->Multiply(current->GetMatrix());
1504 fIsStepEntering = kFALSE;
1505 fIsStepExiting = kTRUE;
1506 fStep = snext;
1507 fNextNode = current;
1508 icandidate = ovlps[i];
1509 candidate = current;
1510 }
1511 }
1512 if (icandidate >= 0) {
1513 icrossed = icandidate;
1514 current = candidate;
1515 Int_t iup = up;
1516 while (iup--)
1517 CdUp();
1518 PopDummy();
1519 PushPath(fLevel + 1);
1520 nmany = fNmany;
1521 up = 1;
1522 currentnode = fCurrentNode;
1523 ovlp = currentnode->IsOverlapping();
1524 restart = kTRUE;
1525 }
1526 if (restart)
1527 continue;
1528 }
1529 snext = TGeoShape::Big();
1530 if (!mothernode->GetVolume()->IsAssembly())
1531 snext = mothernode->GetVolume()->GetShape()->DistFromInside(dpt, dvec, iact, fStep);
1532 fIsStepEntering = kFALSE;
1533 fIsStepExiting = kTRUE;
1534 if (snext < fStep - gTolerance) {
1535 fNextNode = mothernode;
1536 fCurrentMatrix->CopyFrom(matrix);
1537 fStep = snext;
1538 while (up--)
1539 CdUp();
1540 PopDummy();
1541 PushPath();
1542 icrossed = -1;
1543 up = 1;
1544 currentnode = fCurrentNode;
1545 ovlp = currentnode->IsOverlapping();
1546 continue;
1547 }
1548 }
1549 currentnode = mothernode;
1550 ovlp = nextovlp;
1551 up++;
1552 }
1553 }
1554 PopPath();
1555 }
1556 // Compute now the distance in case we have a parallel world
1557 Double_t parstep = TGeoShape::Big();
1558 TGeoPhysicalNode *pnode = nullptr;
1559 if (fGeometry->IsParallelWorldNav()) {
1560 pnode = fGeometry->GetParallelWorld()->FindNextBoundary(fPoint, fDirection, parstep, fStep);
1561 if (pnode) {
1562 // A boundary is hit at less than fPStep
1563 fStep = parstep;
1564 fPoint[0] += fStep * fDirection[0];
1565 fPoint[1] += fStep * fDirection[1];
1566 fPoint[2] += fStep * fDirection[2];
1567 fNextNode = pnode->GetNode();
1568 // icrossed = -4; //
1569 fIsStepEntering = kTRUE;
1570 fIsStepExiting = kFALSE;
1571 cd(pnode->GetName());
1572 nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
1573 while (nextindex >= 0) {
1574 current = fCurrentNode;
1575 CdDown(nextindex);
1576 nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
1577 }
1578 return fCurrentNode;
1579 }
1580 }
1581 fPoint[0] += fStep * fDirection[0];
1582 fPoint[1] += fStep * fDirection[1];
1583 fPoint[2] += fStep * fDirection[2];
1584 fStep += extra;
1585 if (icrossed == -2) {
1586 // Nothing crossed within stepmax -> propagate and return same location
1587 fIsOnBoundary = kFALSE;
1588 return fCurrentNode;
1589 }
1590 fIsOnBoundary = kTRUE;
1591 if (icrossed == -1) {
1592 // Exiting current node.
1593 skip = fCurrentNode;
1594 is_assembly = fCurrentNode->GetVolume()->IsAssembly();
1595 if (!fLevel && !is_assembly) {
1596 fIsOutside = kTRUE;
1597 return nullptr;
1598 }
1599 if (fCurrentNode->IsOffset())
1600 return CrossDivisionCell();
1601 if (fLevel)
1602 CdUp();
1603 else
1604 skip = nullptr;
1605 return CrossBoundaryAndLocate(kFALSE, skip);
1606 }
1607
1608 CdDown(icrossed);
1609 nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
1610 while (nextindex >= 0) {
1611 current = fCurrentNode;
1612 CdDown(nextindex);
1613 nextindex = fCurrentNode->GetVolume()->GetNextNodeIndex();
1614 }
1615 fForcedNode = fCurrentNode;
1616 return CrossBoundaryAndLocate(kTRUE, current);
1617}
1618
1619////////////////////////////////////////////////////////////////////////////////
1620/// Returns deepest node containing current point.
1621
1622TGeoNode *TGeoNavigator::FindNode(Bool_t safe_start)
1623{
1624 fSafety = 0;
1625 fSearchOverlaps = kFALSE;
1626 fIsOutside = kFALSE;
1627 fIsEntering = fIsExiting = kFALSE;
1628 fIsOnBoundary = kFALSE;
1629 fStartSafe = safe_start;
1630 fIsSameLocation = kTRUE;
1631 TGeoNode *last = fCurrentNode;
1632 TGeoNode *found = SearchNode();
1633 if (found != last) {
1634 fIsSameLocation = kFALSE;
1635 } else {
1636 if (last->IsOverlapping())
1637 fIsSameLocation = kTRUE;
1638 }
1639 return found;
1640}
1641
1642////////////////////////////////////////////////////////////////////////////////
1643/// Returns deepest node containing current point.
1644
1645TGeoNode *TGeoNavigator::FindNode(Double_t x, Double_t y, Double_t z)
1646{
1647 fPoint[0] = x;
1648 fPoint[1] = y;
1649 fPoint[2] = z;
1650 fSafety = 0;
1651 fSearchOverlaps = kFALSE;
1652 fIsOutside = kFALSE;
1653 fIsEntering = fIsExiting = kFALSE;
1654 fIsOnBoundary = kFALSE;
1655 fStartSafe = kTRUE;
1656 fIsSameLocation = kTRUE;
1657 TGeoNode *last = fCurrentNode;
1658 TGeoNode *found = SearchNode();
1659 if (found != last) {
1660 fIsSameLocation = kFALSE;
1661 } else {
1662 if (last->IsOverlapping())
1663 fIsSameLocation = kTRUE;
1664 }
1665 return found;
1666}
1667
1668////////////////////////////////////////////////////////////////////////////////
1669/// Computes fast normal to next crossed boundary, assuming that the current point
1670/// is close enough to the boundary. Works only after calling FindNextBoundary.
1671
1685
1686////////////////////////////////////////////////////////////////////////////////
1687/// Computes normal vector to the next surface that will be or was already
1688/// crossed when propagating on a straight line from a given point/direction.
1689/// Returns the normal vector cosines in the MASTER coordinate system. The dot
1690/// product of the normal and the current direction is positive defined.
1691
1692Double_t *TGeoNavigator::FindNormal(Bool_t /*forward*/)
1693{
1694 return FindNormalFast();
1695}
1696
1697////////////////////////////////////////////////////////////////////////////////
1698/// Initialize current point and current direction vector (normalized)
1699/// in MARS. Return corresponding node.
1700
1701TGeoNode *TGeoNavigator::InitTrack(const Double_t *point, const Double_t *dir)
1702{
1703 SetCurrentPoint(point);
1704 SetCurrentDirection(dir);
1705 return FindNode();
1706}
1707
1708////////////////////////////////////////////////////////////////////////////////
1709/// Initialize current point and current direction vector (normalized)
1710/// in MARS. Return corresponding node.
1711
1718
1719////////////////////////////////////////////////////////////////////////////////
1720/// Reset current state flags.
1721
1730
1731//////////////////////////////////////////////////////////////////////////////////
1732/// Wrapper for getting the safety from the parallel world. Takes care of
1733/// caching mechanics + talking to the parallel world.
1734
1735double TGeoNavigator::GetPWSafety(Double_t cpoint[3], Double_t saf_max)
1736{
1737 if (!IsPWSafetyCaching()) {
1738 return fGeometry->GetParallelWorld()->Safety(cpoint, saf_max);
1739 }
1740 const auto cached = GetPWSafetyEstimateFromCache(cpoint);
1741 if (cached > 0) {
1742 // if cache is valid, just use it
1743 return cached;
1744 }
1745 // otherwise we need to evaluate it and update the cache
1746 // we evaluate this with saf_max = infinity to get the best
1747 // possible safety value
1748 auto pw = fGeometry->GetParallelWorld();
1749 const auto newsafety = pw->Safety(cpoint /*saf_max*/);
1750
1751 // we need to be a bit careful: A returned safety value of TGeoShape::Big()
1752 // is not the actual safety and should not be cached
1753 if (newsafety < TGeoShape::Big()) {
1754 fLastPWSafety = newsafety;
1755 fLastPWSaftyPnt[0] = cpoint[0];
1756 fLastPWSaftyPnt[1] = cpoint[1];
1757 fLastPWSaftyPnt[2] = cpoint[2];
1758 } else {
1759 fLastPWSafety = -1;
1760 }
1761 return newsafety;
1762}
1763
1764////////////////////////////////////////////////////////////////////////////////
1765/// Compute safe distance from the current point. This represent the distance
1766/// from POINT to the closest boundary.
1767
1768Double_t TGeoNavigator::Safety(Bool_t inside)
1769{
1770 if (fIsOnBoundary) {
1771 fSafety = 0;
1772 return fSafety;
1773 }
1774 Double_t point[3];
1775 Double_t safpar = TGeoShape::Big(); // safety from parallel world
1776 if (!inside)
1777 fSafety = TGeoShape::Big();
1778
1779 // Check if parallel navigation is enabled
1780 const bool have_PW = fGeometry->IsParallelWorldNav();
1781
1782 if (fIsOutside) {
1783 fSafety = fGeometry->GetTopVolume()->GetShape()->Safety(fPoint, kFALSE);
1784 if (fSafety < gTolerance) {
1785 fSafety = 0;
1786 fIsOnBoundary = kTRUE;
1787 return fSafety;
1788 }
1789
1790 // cross-check against the parallel world safety, using fSafety as limit
1791 if (have_PW) {
1792 safpar = GetPWSafety(fPoint, fSafety);
1793 }
1794 return TMath::Min(fSafety, safpar);
1795 }
1796 //---> convert point to local reference frame of current node
1797 fGlobalMatrix->MasterToLocal(fPoint, point);
1798
1799 //---> compute safety to current node
1800 TGeoVolume *vol = fCurrentNode->GetVolume();
1801 if (!inside) {
1802 fSafety = vol->GetShape()->Safety(point, kTRUE);
1803 //---> if we were just entering, return this safety
1804 if (fSafety < gTolerance) {
1805 fSafety = 0;
1806 fIsOnBoundary = kTRUE;
1807 return fSafety;
1808 }
1809 }
1810
1811 //---> Check against the parallel geometry safety
1812 // cross-check against the parallel world safety, using fSafety as limit
1813 if (have_PW) {
1814 safpar = GetPWSafety(fPoint, fSafety);
1815 }
1816 if (safpar < fSafety)
1817 fSafety = safpar;
1818
1819 //---> if we were just exiting, return this safety
1820 TObjArray *nodes = vol->GetNodes();
1821 Int_t nd = fCurrentNode->GetNdaughters();
1822 if (!nd && !fCurrentOverlapping)
1823 return fSafety;
1824 TGeoNode *node;
1825 Double_t safe;
1826 Int_t id;
1827
1828 // if current volume is divided, we are in the non-divided region. We
1829 // check only the first and the last cell
1830 TGeoPatternFinder *finder = vol->GetFinder();
1831 if (finder) {
1832 Int_t ifirst = finder->GetDivIndex();
1833 node = (TGeoNode *)nodes->UncheckedAt(ifirst);
1834 node->cd();
1835 safe = node->Safety(point, kFALSE);
1836 if (safe < gTolerance) {
1837 fSafety = 0;
1838 fIsOnBoundary = kTRUE;
1839 return fSafety;
1840 }
1841 if (safe < fSafety)
1842 fSafety = safe;
1843 Int_t ilast = ifirst + finder->GetNdiv() - 1;
1844 if (ilast == ifirst)
1845 return fSafety;
1846 node = (TGeoNode *)nodes->UncheckedAt(ilast);
1847 node->cd();
1848 safe = node->Safety(point, kFALSE);
1849 if (safe < gTolerance) {
1850 fSafety = 0;
1851 fIsOnBoundary = kTRUE;
1852 return fSafety;
1853 }
1854 if (safe < fSafety)
1855 fSafety = safe;
1856 if (fCurrentOverlapping && !inside)
1857 SafetyOverlaps();
1858 return fSafety;
1859 }
1860
1861 //---> If no voxels just loop daughters
1862 TGeoVoxelFinder *voxels = vol->GetVoxels();
1863 if (!voxels) {
1864 for (id = 0; id < nd; id++) {
1865 node = (TGeoNode *)nodes->UncheckedAt(id);
1866 safe = node->Safety(point, kFALSE);
1867 if (safe < gTolerance) {
1868 fSafety = 0;
1869 fIsOnBoundary = kTRUE;
1870 return fSafety;
1871 }
1872 if (safe < fSafety)
1873 fSafety = safe;
1874 }
1875 if (fNmany && !inside)
1876 SafetyOverlaps();
1877 return fSafety;
1878 } else {
1879 if (voxels->NeedRebuild()) {
1880 voxels->Voxelize();
1881 vol->FindOverlaps();
1882 }
1883 }
1884
1885 //---> check fast unsafe voxels
1886 Double_t *boxes = voxels->GetBoxes();
1887 for (id = 0; id < nd; id++) {
1888 Int_t ist = 6 * id;
1889 Double_t dxyz = 0.;
1890 Double_t dxyz0 = TMath::Abs(point[0] - boxes[ist + 3]) - boxes[ist];
1891 if (dxyz0 > fSafety)
1892 continue;
1893 Double_t dxyz1 = TMath::Abs(point[1] - boxes[ist + 4]) - boxes[ist + 1];
1894 if (dxyz1 > fSafety)
1895 continue;
1896 Double_t dxyz2 = TMath::Abs(point[2] - boxes[ist + 5]) - boxes[ist + 2];
1897 if (dxyz2 > fSafety)
1898 continue;
1899 if (dxyz0 > 0)
1900 dxyz += dxyz0 * dxyz0;
1901 if (dxyz1 > 0)
1902 dxyz += dxyz1 * dxyz1;
1903 if (dxyz2 > 0)
1904 dxyz += dxyz2 * dxyz2;
1905 if (dxyz >= fSafety * fSafety)
1906 continue;
1907 node = (TGeoNode *)nodes->UncheckedAt(id);
1908 safe = node->Safety(point, kFALSE);
1909 if (safe < gTolerance) {
1910 fSafety = 0;
1911 fIsOnBoundary = kTRUE;
1912 return fSafety;
1913 }
1914 if (safe < fSafety)
1915 fSafety = safe;
1916 }
1917 if (fNmany && !inside)
1918 SafetyOverlaps();
1919 return fSafety;
1920}
1921
1922////////////////////////////////////////////////////////////////////////////////
1923/// Compute safe distance from the current point within an overlapping node
1924
1925void TGeoNavigator::SafetyOverlaps()
1926{
1927 Double_t point[3], local[3];
1928 Double_t safe;
1929 Bool_t contains;
1930 TGeoNode *nodeovlp;
1931 TGeoVolume *vol;
1932 Int_t novlp, io;
1933 Int_t *ovlp;
1934 Int_t safelevel = GetSafeLevel();
1935 PushPath(safelevel + 1);
1936 while (fCurrentOverlapping) {
1937 ovlp = fCurrentNode->GetOverlaps(novlp);
1938 CdUp();
1939 vol = fCurrentNode->GetVolume();
1940 fGeometry->MasterToLocal(fPoint, point);
1941 contains = fCurrentNode->GetVolume()->Contains(point);
1942 safe = fCurrentNode->GetVolume()->GetShape()->Safety(point, contains);
1943 if (safe < fSafety && safe >= 0)
1944 fSafety = safe;
1945 if (!novlp || !contains)
1946 continue;
1947 // we are now in the container, check safety to all candidates
1948 for (io = 0; io < novlp; io++) {
1949 nodeovlp = vol->GetNode(ovlp[io]);
1950 nodeovlp->GetMatrix()->MasterToLocal(point, local);
1951 contains = nodeovlp->GetVolume()->Contains(local);
1952 if (contains) {
1953 CdDown(ovlp[io]);
1954 safe = Safety(kTRUE);
1955 CdUp();
1956 } else {
1957 safe = nodeovlp->GetVolume()->GetShape()->Safety(local, kFALSE);
1958 }
1959 if (safe < fSafety && safe >= 0)
1960 fSafety = safe;
1961 }
1962 }
1963 if (fNmany) {
1964 // We have overlaps up in the branch, check distance to exit
1965 Int_t up = 1;
1966 Int_t imother;
1967 Int_t nmany = fNmany;
1968 Bool_t crtovlp = kFALSE;
1969 Bool_t nextovlp = kFALSE;
1970 TGeoNode *mother, *mup;
1971 TGeoHMatrix *matrix;
1972 while (nmany) {
1973 mother = GetMother(up);
1974 mup = mother;
1975 imother = up + 1;
1976 while (mup->IsOffset())
1977 mup = GetMother(imother++);
1978 nextovlp = mup->IsOverlapping();
1979 if (crtovlp)
1980 nmany--;
1981 if (crtovlp || nextovlp) {
1982 matrix = GetMotherMatrix(up);
1983 matrix->MasterToLocal(fPoint, local);
1984 safe = mother->GetVolume()->GetShape()->Safety(local, kTRUE);
1985 if (safe < fSafety)
1986 fSafety = safe;
1987 crtovlp = nextovlp;
1988 }
1989 up++;
1990 }
1991 }
1992 PopPath();
1993 if (fSafety < gTolerance) {
1994 fSafety = 0.;
1995 fIsOnBoundary = kTRUE;
1996 }
1997}
1998
1999////////////////////////////////////////////////////////////////////////////////
2000/// Returns the deepest node containing fPoint, which must be set a priori.
2001/// Check if parallel world navigation is enabled
2002
2003TGeoNode *TGeoNavigator::SearchNode(Bool_t downwards, const TGeoNode *skipnode)
2004{
2005 if (fGeometry->IsParallelWorldNav()) {
2006 TGeoPhysicalNode *pnode = fGeometry->GetParallelWorld()->FindNode(fPoint);
2007 if (pnode) {
2008 // A node from the parallel world contains the point -> stop the search
2009 // and synchronize with navigation state
2010 pnode->cd();
2011 Int_t crtindex = fCurrentNode->GetVolume()->GetCurrentNodeIndex();
2012 while (crtindex >= 0) {
2013 // Make sure we did not end up in an assembly.
2014 CdDown(crtindex);
2015 crtindex = fCurrentNode->GetVolume()->GetCurrentNodeIndex();
2016 }
2017 return fCurrentNode;
2018 }
2019 }
2020 Double_t point[3];
2021 fNextDaughterIndex = -2;
2022 TGeoVolume *vol = nullptr;
2023 Int_t idebug = TGeoManager::GetVerboseLevel();
2024 Bool_t inside_current = (fCurrentNode == skipnode) ? kTRUE : kFALSE;
2025 if (!downwards) {
2026 // we are looking upwards until inside current node or exit
2027 if (fGeometry->IsActivityEnabled() && !fCurrentNode->GetVolume()->IsActive()) {
2028 // We are inside an inactive volume-> go upwards
2029 CdUp();
2030 fIsSameLocation = kFALSE;
2031 return SearchNode(kFALSE, skipnode);
2032 }
2033 // Check if the current point is still inside the current volume
2034 vol = fCurrentNode->GetVolume();
2035 if (vol->IsAssembly())
2036 inside_current = kTRUE;
2037 // If the current node is not to be skipped
2038 if (!inside_current) {
2039 fGlobalMatrix->MasterToLocal(fPoint, point);
2040 inside_current = vol->Contains(point);
2041 }
2042 // Point might be inside an overlapping node
2043 if (fNmany) {
2044 inside_current = GotoSafeLevel();
2045 }
2046 if (!inside_current) {
2047 // If not, go upwards
2048 fIsSameLocation = kFALSE;
2049 TGeoNode *skip = fCurrentNode; // skip current node at next search
2050 // check if we can go up
2051 if (!fLevel) {
2052 fIsOutside = kTRUE;
2053 return nullptr;
2054 }
2055 CdUp();
2056 return SearchNode(kFALSE, skip);
2057 }
2058 }
2059 vol = fCurrentNode->GetVolume();
2060 fGlobalMatrix->MasterToLocal(fPoint, point);
2061 if (!inside_current && downwards) {
2062 // we are looking downwards
2063 if (fCurrentNode == fForcedNode)
2064 inside_current = kTRUE;
2065 else
2066 inside_current = vol->Contains(point);
2067 if (!inside_current) {
2068 fIsSameLocation = kFALSE;
2069 return nullptr;
2070 } else {
2071 if (fIsOutside) {
2072 fIsOutside = kFALSE;
2073 fIsSameLocation = kFALSE;
2074 }
2075 if (idebug > 4) {
2076 printf("Search node local=(%19.16f, %19.16f, %19.16f) -> %s\n", point[0], point[1], point[2],
2077 fCurrentNode->GetName());
2078 }
2079 }
2080 }
2081 // point inside current (safe) node -> search downwards
2082 TGeoNode *node;
2083 Int_t ncheck = 0;
2084 // if inside an non-overlapping node, reset overlap searches
2085 if (!fCurrentOverlapping) {
2086 fSearchOverlaps = kFALSE;
2087 }
2088
2089 Int_t crtindex = vol->GetCurrentNodeIndex();
2090 while (crtindex >= 0 && downwards) {
2091 // Make sure we did not end up in an assembly.
2092 CdDown(crtindex);
2093 vol = fCurrentNode->GetVolume();
2094 crtindex = vol->GetCurrentNodeIndex();
2095 if (crtindex < 0)
2096 fGlobalMatrix->MasterToLocal(fPoint, point);
2097 }
2098
2099 Int_t nd = vol->GetNdaughters();
2100 // in case there are no daughters
2101 if (!nd)
2102 return fCurrentNode;
2103 if (fGeometry->IsActivityEnabled() && !vol->IsActiveDaughters())
2104 return fCurrentNode;
2105
2106 TGeoPatternFinder *finder = vol->GetFinder();
2107 // point is inside the current node
2108 // first check if inside a division
2109 if (finder) {
2110 node = finder->FindNode(point);
2111 if (!node && fForcedNode) {
2112 // Point *HAS* to be inside a cell
2113 Double_t dir[3];
2114 fGlobalMatrix->MasterToLocalVect(fDirection, dir);
2115 finder->FindNode(point, dir);
2116 node = finder->CdNext();
2117 if (!node)
2118 return fCurrentNode; // inside divided volume but not in a cell
2119 }
2120 if (node && node != skipnode) {
2121 // go inside the division cell and search downwards
2122 fIsSameLocation = kFALSE;
2123 CdDown(node->GetIndex());
2124 fForcedNode = nullptr;
2125 return SearchNode(kTRUE, node);
2126 }
2127 // point is not inside the division, but might be in other nodes
2128 // at the same level (NOT SUPPORTED YET)
2129 while (fCurrentNode && fCurrentNode->IsOffset())
2130 CdUp();
2131 return fCurrentNode;
2132 }
2133 // second, look if current volume is voxelized
2134 TGeoVoxelFinder *voxels = vol->GetVoxels();
2135 Int_t *check_list = nullptr;
2136 Int_t id;
2137 if (voxels) {
2138 // get the list of nodes passing thorough the current voxel
2139 check_list = voxels->GetCheckList(&point[0], ncheck, *fCache->GetInfo());
2140 // if none in voxel, see if this is the last one
2141 if (!check_list) {
2142 if (!fCurrentNode->GetVolume()->IsAssembly()) {
2143 fCache->ReleaseInfo();
2144 return fCurrentNode;
2145 }
2146 // Point in assembly - go up
2147 node = fCurrentNode;
2148 if (!fLevel) {
2149 fIsOutside = kTRUE;
2150 fCache->ReleaseInfo();
2151 return nullptr;
2152 }
2153 CdUp();
2154 fCache->ReleaseInfo();
2155 return SearchNode(kFALSE, node);
2156 }
2157 // loop all nodes in voxel
2158 for (id = 0; id < ncheck; id++) {
2159 node = vol->GetNode(check_list[id]);
2160 if (node == skipnode)
2161 continue;
2162 if (fGeometry->IsActivityEnabled() && !node->GetVolume()->IsActive())
2163 continue;
2164 if ((id < (ncheck - 1)) && node->IsOverlapping()) {
2165 // make the cluster of overlaps
2166 if (ncheck + fOverlapMark > fOverlapSize) {
2167 fOverlapSize = 2 * (ncheck + fOverlapMark);
2168 delete[] fOverlapClusters;
2169 fOverlapClusters = new Int_t[fOverlapSize];
2170 }
2171 Int_t *cluster = fOverlapClusters + fOverlapMark;
2172 Int_t nc = GetTouchedCluster(id, &point[0], check_list, ncheck, cluster);
2173 if (nc > 1) {
2174 fOverlapMark += nc;
2175 node = FindInCluster(cluster, nc);
2176 fOverlapMark -= nc;
2177 fCache->ReleaseInfo();
2178 return node;
2179 }
2180 }
2181 CdDown(check_list[id]);
2182 fForcedNode = nullptr;
2183 node = SearchNode(kTRUE);
2184 if (node) {
2185 fIsSameLocation = kFALSE;
2186 fCache->ReleaseInfo();
2187 return node;
2188 }
2189 CdUp();
2190 }
2191 if (!fCurrentNode->GetVolume()->IsAssembly()) {
2192 fCache->ReleaseInfo();
2193 return fCurrentNode;
2194 }
2195 node = fCurrentNode;
2196 if (!fLevel) {
2197 fIsOutside = kTRUE;
2198 fCache->ReleaseInfo();
2199 return nullptr;
2200 }
2201 CdUp();
2202 fCache->ReleaseInfo();
2203 return SearchNode(kFALSE, node);
2204 }
2205 // if there are no voxels just loop all daughters
2206 for (id = 0; id < nd; id++) {
2207 node = fCurrentNode->GetDaughter(id);
2208 if (node == skipnode)
2209 continue;
2210 if (fGeometry->IsActivityEnabled() && !node->GetVolume()->IsActive())
2211 continue;
2212 CdDown(id);
2213 fForcedNode = nullptr;
2214 node = SearchNode(kTRUE);
2215 if (node) {
2216 fIsSameLocation = kFALSE;
2217 return node;
2218 }
2219 CdUp();
2220 }
2221 // point is not inside one of the daughters, so it is in the current vol
2222 if (fCurrentNode->GetVolume()->IsAssembly()) {
2223 node = fCurrentNode;
2224 if (!fLevel) {
2225 fIsOutside = kTRUE;
2226 return nullptr;
2227 }
2228 CdUp();
2229 return SearchNode(kFALSE, node);
2230 }
2231 return fCurrentNode;
2232}
2233
2234////////////////////////////////////////////////////////////////////////////////
2235/// Find a node inside a cluster of overlapping nodes. Current node must
2236/// be on top of all the nodes in cluster. Always nc>1.
2237
2238TGeoNode *TGeoNavigator::FindInCluster(Int_t *cluster, Int_t nc)
2239{
2240 TGeoNode *clnode = nullptr;
2241 TGeoNode *priority = fLastNode;
2242 // save current node
2243 TGeoNode *current = fCurrentNode;
2244 TGeoNode *found = nullptr;
2245 // save path
2246 Int_t ipop = PushPath();
2247 // mark this search
2248 fSearchOverlaps = kTRUE;
2249 Int_t deepest = fLevel;
2250 Int_t deepest_virtual = fLevel - GetVirtualLevel();
2251 Int_t found_virtual = 0;
2252 Bool_t replace = kFALSE;
2253 Bool_t added = kFALSE;
2254 Int_t i;
2255 for (i = 0; i < nc; i++) {
2256 clnode = current->GetDaughter(cluster[i]);
2257 CdDown(cluster[i]);
2258 Bool_t max_priority = (clnode == fNextNode) ? kTRUE : kFALSE;
2259 found = SearchNode(kTRUE, clnode);
2260 if (!fSearchOverlaps || max_priority) {
2261 // an only was found during the search -> exiting
2262 // The node given by FindNextBoundary returned -> exiting
2263 PopDummy(ipop);
2264 return found;
2265 }
2266 found_virtual = fLevel - GetVirtualLevel();
2267 if (added) {
2268 // we have put something in stack -> check it
2269 if (found_virtual > deepest_virtual) {
2270 replace = kTRUE;
2271 } else {
2272 if (found_virtual == deepest_virtual) {
2273 if (fLevel > deepest) {
2274 replace = kTRUE;
2275 } else {
2276 if ((fLevel == deepest) && (clnode == priority))
2277 replace = kTRUE;
2278 else
2279 replace = kFALSE;
2280 }
2281 } else
2282 replace = kFALSE;
2283 }
2284 // if this was the last checked node
2285 if (i == (nc - 1)) {
2286 if (replace) {
2287 PopDummy(ipop);
2288 return found;
2289 } else {
2290 fCurrentOverlapping = PopPath();
2291 PopDummy(ipop);
2292 return fCurrentNode;
2293 }
2294 }
2295 // we still have to go on
2296 if (replace) {
2297 // reset stack
2298 PopDummy();
2299 PushPath();
2300 deepest = fLevel;
2301 deepest_virtual = found_virtual;
2302 }
2303 // restore top of cluster
2304 fCurrentOverlapping = PopPath(ipop);
2305 } else {
2306 // the stack was clean, push new one
2307 PushPath();
2308 added = kTRUE;
2309 deepest = fLevel;
2310 deepest_virtual = found_virtual;
2311 // restore original path
2312 fCurrentOverlapping = PopPath(ipop);
2313 }
2314 }
2315 PopDummy(ipop);
2316 return fCurrentNode;
2317}
2318
2319////////////////////////////////////////////////////////////////////////////////
2320/// Make the cluster of overlapping nodes in a voxel, containing point in reference
2321/// of the mother. Returns number of nodes containing the point. Nodes should not be
2322/// offsets.
2323
2324Int_t TGeoNavigator::GetTouchedCluster(Int_t start, Double_t *point, Int_t *check_list, Int_t ncheck, Int_t *result)
2325{
2326 // we are in the mother reference system
2327 TGeoNode *current = fCurrentNode->GetDaughter(check_list[start]);
2328 Int_t novlps = 0;
2329 Int_t *ovlps = current->GetOverlaps(novlps);
2330 if (!ovlps)
2331 return 0;
2332 Double_t local[3];
2333 // intersect check list with overlap list
2334 Int_t ntotal = 0;
2335 current->MasterToLocal(point, &local[0]);
2336 if (current->GetVolume()->Contains(&local[0])) {
2337 result[ntotal++] = check_list[start];
2338 }
2339
2340 Int_t jst = 0, i, j;
2341 while ((jst < novlps) && (ovlps[jst] <= check_list[start]))
2342 jst++;
2343 if (jst == novlps)
2344 return 0;
2345 for (i = start; i < ncheck; i++) {
2346 for (j = jst; j < novlps; j++) {
2347 if (check_list[i] == ovlps[j]) {
2348 // overlapping node in voxel -> check if touched
2349 current = fCurrentNode->GetDaughter(check_list[i]);
2350 if (fGeometry->IsActivityEnabled() && !current->GetVolume()->IsActive())
2351 continue;
2352 current->MasterToLocal(point, &local[0]);
2353 if (current->GetVolume()->Contains(&local[0])) {
2354 result[ntotal++] = check_list[i];
2355 }
2356 }
2357 }
2358 }
2359 return ntotal;
2360}
2361
2362////////////////////////////////////////////////////////////////////////////////
2363/// Make a rectiliniar step of length fStep from current point (fPoint) on current
2364/// direction (fDirection). If the step is imposed by geometry, is_geom flag
2365/// must be true (default). The cross flag specifies if the boundary should be
2366/// crossed in case of a geometry step (default true). Returns new node after step.
2367/// Set also on boundary condition.
2368
2369TGeoNode *TGeoNavigator::Step(Bool_t is_geom, Bool_t cross)
2370{
2371 Double_t epsil = 0;
2372 if (fStep < 1E-6) {
2373 fIsNullStep = kTRUE;
2374 if (fStep < 0)
2375 fStep = 0.;
2376 } else {
2377 fIsNullStep = kFALSE;
2378 }
2379 if (is_geom)
2380 epsil = (cross) ? 1E-6 : -1E-6;
2381 TGeoNode *old = fCurrentNode;
2382 Int_t idold = GetNodeId();
2383 if (fIsOutside)
2384 old = nullptr;
2385 fStep += epsil;
2386 for (Int_t i = 0; i < 3; i++)
2387 fPoint[i] += fStep * fDirection[i];
2388 TGeoNode *current = FindNode();
2389 if (is_geom) {
2390 fIsEntering = (current == old) ? kFALSE : kTRUE;
2391 if (!fIsEntering) {
2392 Int_t id = GetNodeId();
2393 fIsEntering = (id == idold) ? kFALSE : kTRUE;
2394 }
2395 fIsExiting = !fIsEntering;
2396 if (fIsEntering && fIsNullStep)
2397 fIsNullStep = kFALSE;
2398 fIsOnBoundary = kTRUE;
2399 } else {
2400 fIsEntering = fIsExiting = kFALSE;
2401 fIsOnBoundary = kFALSE;
2402 }
2403 return current;
2404}
2405
2406////////////////////////////////////////////////////////////////////////////////
2407/// Find level of virtuality of current overlapping node (number of levels
2408/// up having the same tracking media.
2409
2410Int_t TGeoNavigator::GetVirtualLevel()
2411{
2412 // return if the current node is ONLY
2413 if (!fCurrentOverlapping)
2414 return 0;
2415 Int_t new_media = 0;
2416 TGeoMedium *medium = fCurrentNode->GetMedium();
2417 Int_t virtual_level = 1;
2418 TGeoNode *mother = nullptr;
2419
2420 while ((mother = GetMother(virtual_level))) {
2421 if (!mother->IsOverlapping() && !mother->IsOffset()) {
2422 if (!new_media)
2423 new_media = (mother->GetMedium() == medium) ? 0 : virtual_level;
2424 break;
2425 }
2426 if (!new_media)
2427 new_media = (mother->GetMedium() == medium) ? 0 : virtual_level;
2428 virtual_level++;
2429 }
2430 return (new_media == 0) ? virtual_level : (new_media - 1);
2431}
2432
2433////////////////////////////////////////////////////////////////////////////////
2434/// Go upwards the tree until a non-overlapping node
2435
2436Bool_t TGeoNavigator::GotoSafeLevel()
2437{
2438 while (fCurrentOverlapping && fLevel)
2439 CdUp();
2440 Double_t point[3];
2441 fGlobalMatrix->MasterToLocal(fPoint, point);
2442 if (!fCurrentNode->GetVolume()->Contains(point))
2443 return kFALSE;
2444 if (fNmany) {
2445 // We still have overlaps on the branch
2446 Int_t up = 1;
2447 Int_t imother;
2448 Int_t nmany = fNmany;
2449 Bool_t ovlp = kFALSE;
2450 Bool_t nextovlp = kFALSE;
2451 TGeoNode *mother, *mup;
2452 TGeoHMatrix *matrix;
2453 while (nmany) {
2454 mother = GetMother(up);
2455 if (!mother)
2456 return kTRUE;
2457 mup = mother;
2458 imother = up + 1;
2459 while (mup->IsOffset())
2460 mup = GetMother(imother++);
2461 nextovlp = mup->IsOverlapping();
2462 if (ovlp)
2463 nmany--;
2464 if (ovlp || nextovlp) {
2465 // check if the point is in the next node up
2466 matrix = GetMotherMatrix(up);
2467 matrix->MasterToLocal(fPoint, point);
2468 if (!mother->GetVolume()->Contains(point)) {
2469 up++;
2470 while (up--)
2471 CdUp();
2472 return GotoSafeLevel();
2473 }
2474 }
2475 ovlp = nextovlp;
2476 up++;
2477 }
2478 }
2479 return kTRUE;
2480}
2481
2482////////////////////////////////////////////////////////////////////////////////
2483/// Go upwards the tree until a non-overlapping node
2484
2485Int_t TGeoNavigator::GetSafeLevel() const
2486{
2487 Bool_t overlapping = fCurrentOverlapping;
2488 if (!overlapping)
2489 return fLevel;
2490 Int_t level = fLevel;
2491 TGeoNode *node;
2492 while (overlapping && level) {
2493 level--;
2494 node = GetMother(fLevel - level);
2495 if (!node->IsOffset())
2496 overlapping = node->IsOverlapping();
2497 }
2498 return level;
2499}
2500
2501////////////////////////////////////////////////////////////////////////////////
2502/// Inspects path and all flags for the current state.
2503
2504void TGeoNavigator::InspectState() const
2505{
2506 Info("InspectState", "Current path is: %s", GetPath());
2507 Int_t level;
2508 TGeoNode *node;
2509 Bool_t is_offset, is_overlapping;
2510 for (level = 0; level < fLevel + 1; level++) {
2511 node = GetMother(fLevel - level);
2512 if (!node)
2513 continue;
2514 is_offset = node->IsOffset();
2515 is_overlapping = node->IsOverlapping();
2516 Info("InspectState", "level %i: %s div=%i many=%i", level, node->GetName(), is_offset, is_overlapping);
2517 }
2518 Info("InspectState", "on_bound=%i entering=%i", fIsOnBoundary, fIsEntering);
2519}
2520
2521////////////////////////////////////////////////////////////////////////////////
2522/// Checks if point (x,y,z) is still in the current node.
2523/// check if this is an overlapping node
2524
2525Bool_t TGeoNavigator::IsSameLocation(Double_t x, Double_t y, Double_t z, Bool_t change)
2526{
2527 Double_t oldpt[3];
2528 if (fLastSafety > 0) {
2529 Double_t dx = (x - fLastPoint[0]);
2530 Double_t dy = (y - fLastPoint[1]);
2531 Double_t dz = (z - fLastPoint[2]);
2532 Double_t dsq = dx * dx + dy * dy + dz * dz;
2533 if (dsq < fLastSafety * fLastSafety) {
2534 if (change) {
2535 fPoint[0] = x;
2536 fPoint[1] = y;
2537 fPoint[2] = z;
2538 memcpy(fLastPoint, fPoint, 3 * sizeof(Double_t));
2539 fLastSafety -= TMath::Sqrt(dsq);
2540 }
2541 return kTRUE;
2542 }
2543 if (change)
2544 fLastSafety = 0;
2545 }
2546 if (fCurrentOverlapping) {
2547 // TGeoNode *current = fCurrentNode;
2548 Int_t cid = GetCurrentNodeId();
2549 if (!change)
2550 PushPoint();
2551 memcpy(oldpt, fPoint, kN3);
2552 SetCurrentPoint(x, y, z);
2553 SearchNode();
2554 memcpy(fPoint, oldpt, kN3);
2555 Bool_t same = (cid == GetCurrentNodeId()) ? kTRUE : kFALSE;
2556 if (!change)
2557 PopPoint();
2558 return same;
2559 }
2560
2561 Double_t point[3];
2562 point[0] = x;
2563 point[1] = y;
2564 point[2] = z;
2565 if (change)
2566 memcpy(fPoint, point, kN3);
2567 TGeoVolume *vol = fCurrentNode->GetVolume();
2568 if (fIsOutside) {
2569 if (vol->GetShape()->Contains(point)) {
2570 if (!change)
2571 return kFALSE;
2572 FindNode(x, y, z);
2573 return kFALSE;
2574 }
2575 return kTRUE;
2576 }
2577 Double_t local[3];
2578 // convert to local frame
2579 fGlobalMatrix->MasterToLocal(point, local);
2580 // check if still in current volume.
2581 if (!vol->GetShape()->Contains(local)) {
2582 if (!change)
2583 return kFALSE;
2584 CdUp();
2585 FindNode(x, y, z);
2586 return kFALSE;
2587 }
2588
2589 // Check if the point is in a parallel world volume
2590 if (fGeometry->IsParallelWorldNav()) {
2591 TGeoPhysicalNode *pnode = fGeometry->GetParallelWorld()->FindNode(fPoint);
2592 if (pnode) {
2593 if (!change)
2594 return kFALSE;
2595 pnode->cd();
2596 Int_t crtindex = fCurrentNode->GetVolume()->GetCurrentNodeIndex();
2597 while (crtindex >= 0) {
2598 // Make sure we did not end up in an assembly.
2599 CdDown(crtindex);
2600 crtindex = fCurrentNode->GetVolume()->GetCurrentNodeIndex();
2601 }
2602 return kFALSE;
2603 }
2604 }
2605 // check if there are daughters
2606 Int_t nd = vol->GetNdaughters();
2607 if (!nd)
2608 return kTRUE;
2609
2610 TGeoNode *node;
2611 TGeoPatternFinder *finder = vol->GetFinder();
2612 if (finder) {
2613 node = finder->FindNode(local);
2614 if (node) {
2615 if (!change)
2616 return kFALSE;
2617 CdDown(node->GetIndex());
2618 SearchNode(kTRUE, node);
2619 return kFALSE;
2620 }
2621 return kTRUE;
2622 }
2623 // if we are not allowed to do changes, save the current path
2624 TGeoVoxelFinder *voxels = vol->GetVoxels();
2625 Int_t *check_list = nullptr;
2626 Int_t ncheck = 0;
2627 Double_t local1[3];
2628 if (voxels) {
2629 check_list = voxels->GetCheckList(local, ncheck, *fCache->GetInfo());
2630 if (!check_list) {
2631 fCache->ReleaseInfo();
2632 return kTRUE;
2633 }
2634 if (!change)
2635 PushPath();
2636 for (Int_t id = 0; id < ncheck; id++) {
2637 // node = vol->GetNode(check_list[id]);
2638 CdDown(check_list[id]);
2639 fGlobalMatrix->MasterToLocal(point, local1);
2640 if (fCurrentNode->GetVolume()->GetShape()->Contains(local1)) {
2641 if (!change) {
2642 PopPath();
2643 fCache->ReleaseInfo();
2644 return kFALSE;
2645 }
2646 SearchNode(kTRUE);
2647 fCache->ReleaseInfo();
2648 return kFALSE;
2649 }
2650 CdUp();
2651 }
2652 if (!change)
2653 PopPath();
2654 fCache->ReleaseInfo();
2655 return kTRUE;
2656 }
2657 Int_t id = 0;
2658 if (!change)
2659 PushPath();
2660 while (fCurrentNode && fCurrentNode->GetDaughter(id++)) {
2661 CdDown(id - 1);
2662 fGlobalMatrix->MasterToLocal(point, local1);
2663 if (fCurrentNode->GetVolume()->GetShape()->Contains(local1)) {
2664 if (!change) {
2665 PopPath();
2666 return kFALSE;
2667 }
2668 SearchNode(kTRUE);
2669 return kFALSE;
2670 }
2671 CdUp();
2672 if (id == nd) {
2673 if (!change)
2674 PopPath();
2675 return kTRUE;
2676 }
2677 }
2678 if (!change)
2679 PopPath();
2680 return kTRUE;
2681}
2682
2683////////////////////////////////////////////////////////////////////////////////
2684/// In case a previous safety value was computed, check if the safety region is
2685/// still safe for the current point and proposed step. Return value changed only
2686/// if proposed distance is safe.
2687
2688Bool_t TGeoNavigator::IsSafeStep(Double_t proposed, Double_t &newsafety) const
2689{
2690 // Last safety not computed.
2691 if (fLastSafety < gTolerance)
2692 return kFALSE;
2693 // Proposed step too small
2694 if (proposed < gTolerance) {
2695 newsafety = fLastSafety - proposed;
2696 return kTRUE;
2697 }
2698 // Normal step
2699 Double_t dist = (fPoint[0] - fLastPoint[0]) * (fPoint[0] - fLastPoint[0]) +
2700 (fPoint[1] - fLastPoint[1]) * (fPoint[1] - fLastPoint[1]) +
2701 (fPoint[2] - fLastPoint[2]) * (fPoint[2] - fLastPoint[2]);
2702 dist = TMath::Sqrt(dist);
2703 Double_t safe = fLastSafety - dist;
2704 if (safe < proposed)
2705 return kFALSE;
2706 newsafety = safe;
2707 return kTRUE;
2708}
2709
2710////////////////////////////////////////////////////////////////////////////////
2711/// Check if a new point with given coordinates is the same as the last located one.
2712
2713Bool_t TGeoNavigator::IsSamePoint(Double_t x, Double_t y, Double_t z) const
2714{
2715 if (TMath::Abs(x - fLastPoint[0]) < 1.E-20) {
2716 if (TMath::Abs(y - fLastPoint[1]) < 1.E-20) {
2717 if (TMath::Abs(z - fLastPoint[2]) < 1.E-20)
2718 return kTRUE;
2719 }
2720 }
2721 return kFALSE;
2722}
2723
2724////////////////////////////////////////////////////////////////////////////////
2725/// Backup the current state without affecting the cache stack.
2726
2732
2733////////////////////////////////////////////////////////////////////////////////
2734/// Restore a backed-up state without affecting the cache stack.
2735
2745
2746////////////////////////////////////////////////////////////////////////////////
2747/// Return stored current matrix (global matrix of the next touched node).
2748
2757
2758////////////////////////////////////////////////////////////////////////////////
2759/// Get path to the current node in the form /node0/node1/...
2760
2761const char *TGeoNavigator::GetPath() const
2762{
2763 if (fIsOutside)
2764 return kGeoOutsidePath;
2765 return fCache->GetPath();
2766}
2767
2768////////////////////////////////////////////////////////////////////////////////
2769/// Convert coordinates from master volume frame to top.
2770
2775
2776////////////////////////////////////////////////////////////////////////////////
2777/// Convert coordinates from top volume frame to master.
2778
2783
2784////////////////////////////////////////////////////////////////////////////////
2785/// Reset the navigator.
2786
2787void TGeoNavigator::ResetAll()
2788{
2789 GetHMatrix();
2790 *fCurrentMatrix = gGeoIdentity;
2791 fCurrentNode = fGeometry->GetTopNode();
2792 ResetState();
2793 fStep = 0.;
2794 fSafety = 0.;
2795 fLastSafety = 0.;
2796 fLevel = 0;
2797 fNmany = 0;
2798 fNextDaughterIndex = -2;
2799 fCurrentOverlapping = kFALSE;
2800 fStartSafe = kFALSE;
2801 fIsSameLocation = kFALSE;
2802 fIsNullStep = kFALSE;
2803 fCurrentVolume = fGeometry->GetTopVolume();
2804 fCurrentNode = fGeometry->GetTopNode();
2805 fLastNode = nullptr;
2806 fNextNode = nullptr;
2807 fPath = "";
2808 if (fCache) {
2809 Bool_t dummy = fCache->IsDummy();
2810 Bool_t nodeid = fCache->HasIdArray();
2811 delete fCache;
2812 fCache = nullptr;
2813 delete fBackupState;
2814 fBackupState = nullptr;
2815 BuildCache(dummy, nodeid);
2816 }
2817}
2818
2819
2820////////////////////////////////////////////////////////////////////////////////
2821/// Add a new navigator to the array.
2822
Bool_t
bool Bool_t
Boolean (0=false, 1=true) (bool)
Definition RtypesCore.h:77
Int_t
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
kFALSE
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
Double_t
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
kTRUE
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
TRangeDynCast
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Definition TCollection.h:360
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 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
Definition TGWin32VirtualXProxy.cxx:245
result
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 result
Definition TGWin32VirtualXProxy.cxx:174
index
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
Definition TGWin32VirtualXProxy.cxx:168
length
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 length
Definition TGWin32VirtualXProxy.cxx:245
id
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize id
Definition TGWin32VirtualXProxy.cxx:94
name
char name[80]
Definition TGX11.cxx:145
kN3
const Int_t kN3
Definition TGeoMatrix.cxx:221
gGeoIdentity
R__EXTERN TGeoIdentity * gGeoIdentity
Definition TGeoMatrix.h:538
kN3
const Int_t kN3
Definition TGeoNavigator.cxx:48
kGeoOutsidePath
const char * kGeoOutsidePath
Definition TGeoNavigator.cxx:47
gTolerance
static Double_t gTolerance
Definition TGeoNavigator.cxx:46
TCollection::SetOwner
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
Definition TCollection.cxx:777
TGeoCacheState
Class storing the state of the cache at a given moment.
Definition TGeoCache.h:28
TGeoCacheState::SetState
void SetState(Int_t level, Int_t startlevel, Int_t nmany, Bool_t ovlp, Double_t *point=nullptr)
Fill current modeller state.
Definition TGeoCache.cxx:591
TGeoHMatrix
Matrix class used for computing global transformations Should NOT be used for node definition.
Definition TGeoMatrix.h:459
TGeoHMatrix::CopyFrom
void CopyFrom(const TGeoMatrix *other)
Fast copy method.
Definition TGeoMatrix.cxx:2520
TGeoHMatrix::GetTranslation
const Double_t * GetTranslation() const override
Definition TGeoMatrix.h:528
TGeoHMatrix::Multiply
void Multiply(const TGeoMatrix *right)
multiply to the right with an other transformation if right is identity matrix, just return
Definition TGeoMatrix.cxx:2625
TGeoManager
The manager class for any TGeo geometry.
Definition TGeoManager.h:46
TGeoManager::GetCurrentNode
TGeoNode * GetCurrentNode() const
Definition TGeoManager.h:511
TGeoManager::GetParallelWorld
TGeoParallelWorld * GetParallelWorld() const
Definition TGeoManager.h:613
TGeoManager::GetMaxLevel
Int_t GetMaxLevel() const
Definition TGeoManager.h:520
TGeoManager::IsParallelWorldNav
Bool_t IsParallelWorldNav() const
Definition TGeoManager.h:615
TGeoManager::IsActivityEnabled
Bool_t IsActivityEnabled() const
Definition TGeoManager.h:421
TGeoManager::GetTopNode
TGeoNode * GetTopNode() const
Definition TGeoManager.h:525
TGeoManager::GetVerboseLevel
static Int_t GetVerboseLevel()
Set verbosity level (static function).
Definition TGeoManager.cxx:4192
TGeoManager::MasterToLocal
void MasterToLocal(const Double_t *master, Double_t *local) const
Definition TGeoManager.h:550
TGeoManager::GetTopVolume
TGeoVolume * GetTopVolume() const
Definition TGeoManager.h:524
TGeoManager::ThreadId
static Int_t ThreadId()
Translates the current thread id to an ordinal number.
Definition TGeoManager.cxx:1058
TGeoMatrix::LocalToMasterVect
virtual void LocalToMasterVect(const Double_t *local, Double_t *master) const
convert a vector by multiplying its column vector (x, y, z, 1) to matrix inverse
Definition TGeoMatrix.cxx:388
TGeoMatrix::MasterToLocal
virtual void MasterToLocal(const Double_t *master, Double_t *local) const
convert a point by multiplying its column vector (x, y, z, 1) to matrix
Definition TGeoMatrix.cxx:427
TGeoMatrix::MasterToLocalVect
virtual void MasterToLocalVect(const Double_t *master, Double_t *local) const
convert a point by multiplying its column vector (x, y, z, 1) to matrix
Definition TGeoMatrix.cxx:452
TGeoMatrix::RegisterYourself
virtual void RegisterYourself()
Register the matrix in the current manager, which will become the owner.
Definition TGeoMatrix.cxx:540
TGeoMatrix::LocalToMaster
virtual void LocalToMaster(const Double_t *local, Double_t *master) const
convert a point by multiplying its column vector (x, y, z, 1) to matrix inverse
Definition TGeoMatrix.cxx:366
TGeoMedium
Media are used to store properties related to tracking and which are useful only when using geometry ...
Definition TGeoMedium.h:23
TGeoNavigatorArray::AddNavigator
TGeoNavigator * AddNavigator()
Add a new navigator to the array.
Definition TGeoNavigator.cxx:2823
TGeoNavigatorArray::fGeoManager
TGeoManager * fGeoManager
Definition TGeoNavigator.h:299
TGeoNavigatorArray::SetCurrentNavigator
TGeoNavigator * SetCurrentNavigator(Int_t inav)
Definition TGeoNavigator.h:311
TGeoNavigator
Class providing navigation API for TGeo geometries.
Definition TGeoNavigator.h:34
TGeoNavigator::CdUp
void CdUp()
Go one level up in geometry.
Definition TGeoNavigator.cxx:351
TGeoNavigator::DoBackupState
void DoBackupState()
Backup the current state without affecting the cache stack.
Definition TGeoNavigator.cxx:2727
TGeoNavigator::GetMother
TGeoNode * GetMother(Int_t up=1) const
Definition TGeoNavigator.h:157
TGeoNavigator::fLastPWSaftyPnt
Double_t fLastPWSaftyPnt[3]
! last point for which parallel world safety was "evaluated"
Definition TGeoNavigator.h:54
TGeoNavigator::DoRestoreState
void DoRestoreState()
Restore a backed-up state without affecting the cache stack.
Definition TGeoNavigator.cxx:2736
TGeoNavigator::fPoint
Double_t fPoint[3]
! current point
Definition TGeoNavigator.h:51
TGeoNavigator::fSearchOverlaps
Bool_t fSearchOverlaps
! flag set when an overlapping cluster is searched
Definition TGeoNavigator.h:63
TGeoNavigator::fIsExiting
Bool_t fIsExiting
! flag that current track is about to leave current node
Definition TGeoNavigator.h:67
TGeoNavigator::fPath
TString fPath
! path to current node
Definition TGeoNavigator.h:86
TGeoNavigator::fDivMatrix
TGeoHMatrix * fDivMatrix
! current local matrix of the selected division cell
Definition TGeoNavigator.h:85
TGeoNavigator::CrossBoundaryAndLocate
TGeoNode * CrossBoundaryAndLocate(Bool_t downwards, TGeoNode *skipnode)
Cross next boundary and locate within current node The current point must be on the boundary of fCurr...
Definition TGeoNavigator.cxx:522
TGeoNavigator::fLastPWSafety
Double_t fLastPWSafety
! last safety returned from parallel world (negative if invalid)
Definition TGeoNavigator.h:55
TGeoNavigator::GetHMatrix
TGeoHMatrix * GetHMatrix()
Return stored current matrix (global matrix of the next touched node).
Definition TGeoNavigator.cxx:2749
TGeoNavigator::fCache
TGeoNodeCache * fCache
! cache of states
Definition TGeoNavigator.h:75
TGeoNavigator::fStartSafe
Bool_t fStartSafe
! flag a safe start for point classification
Definition TGeoNavigator.h:65
TGeoNavigator::CdNext
void CdNext()
Do a cd to the node found next by FindNextBoundary.
Definition TGeoNavigator.cxx:404
TGeoNavigator::Safety
Double_t Safety(Bool_t inside=kFALSE)
Compute safe distance from the current point.
Definition TGeoNavigator.cxx:1768
TGeoNavigator::GotoSafeLevel
Bool_t GotoSafeLevel()
Go upwards the tree until a non-overlapping node.
Definition TGeoNavigator.cxx:2436
TGeoNavigator::fNormal
Double_t fNormal[3]
! cosine of incident angle on current checked surface
Definition TGeoNavigator.h:48
TGeoNavigator::cd
Bool_t cd(const char *path="")
Browse the tree of nodes starting from top node according to pathname.
Definition TGeoNavigator.cxx:209
TGeoNavigator::fLastPoint
Double_t fLastPoint[3]
! last point for which safety was computed
Definition TGeoNavigator.h:53
TGeoNavigator::IsSameLocation
Bool_t IsSameLocation() const
Definition TGeoNavigator.h:135
TGeoNavigator::fCldir
Double_t fCldir[3]
! unit vector to current closest shape
Definition TGeoNavigator.h:49
TGeoNavigator::fIsStepEntering
Bool_t fIsStepEntering
! flag that next geometric step will enter new volume
Definition TGeoNavigator.h:68
TGeoNavigator::GetNodeId
Int_t GetNodeId() const
Definition TGeoNavigator.h:154
TGeoNavigator::GetVirtualLevel
Int_t GetVirtualLevel()
Find level of virtuality of current overlapping node (number of levels up having the same tracking me...
Definition TGeoNavigator.cxx:2410
TGeoNavigator::PopPoint
Bool_t PopPoint()
Definition TGeoNavigator.h:266
TGeoNavigator::fOverlapSize
Int_t fOverlapSize
! current size of fOverlapClusters
Definition TGeoNavigator.h:60
TGeoNavigator::InitTrack
TGeoNode * InitTrack(const Double_t *point, const Double_t *dir)
Initialize current point and current direction vector (normalized) in MARS.
Definition TGeoNavigator.cxx:1701
TGeoNavigator::InspectState
void InspectState() const
Inspects path and all flags for the current state.
Definition TGeoNavigator.cxx:2504
TGeoNavigator::PushPoint
Int_t PushPoint(Int_t startlevel=0)
Definition TGeoNavigator.h:265
TGeoNavigator::Step
TGeoNode * Step(Bool_t is_geom=kTRUE, Bool_t cross=kTRUE)
Make a rectiliniar step of length fStep from current point (fPoint) on current direction (fDirection)...
Definition TGeoNavigator.cxx:2369
TGeoNavigator::FindInCluster
TGeoNode * FindInCluster(Int_t *cluster, Int_t nc)
Find a node inside a cluster of overlapping nodes.
Definition TGeoNavigator.cxx:2238
TGeoNavigator::fStep
Double_t fStep
! step to be done from current point and direction
Definition TGeoNavigator.h:45
TGeoNavigator::SafetyOverlaps
void SafetyOverlaps()
Compute safe distance from the current point within an overlapping node.
Definition TGeoNavigator.cxx:1925
TGeoNavigator::CrossDivisionCell
TGeoNode * CrossDivisionCell()
Cross a division cell.
Definition TGeoNavigator.cxx:460
TGeoNavigator::ResetState
void ResetState()
Reset current state flags.
Definition TGeoNavigator.cxx:1722
TGeoNavigator::IsPWSafetyCaching
static Bool_t IsPWSafetyCaching()
Definition TGeoNavigator.h:234
TGeoNavigator::FindNextDaughterBoundary
TGeoNode * FindNextDaughterBoundary(Double_t *point, Double_t *dir, Int_t &idaughter, Bool_t compmatrix=kFALSE)
Computes as fStep the distance to next daughter of the current volume.
Definition TGeoNavigator.cxx:1020
TGeoNavigator::fIsSameLocation
Bool_t fIsSameLocation
! flag that a new point is in the same node as previous
Definition TGeoNavigator.h:72
TGeoNavigator::GetBranchNumbers
void GetBranchNumbers(Int_t *copyNumbers, Int_t *volumeNumbers) const
Fill node copy numbers of current branch into an array.
Definition TGeoNavigator.cxx:443
TGeoNavigator::CheckPath
Bool_t CheckPath(const char *path) const
Check if a geometry path is valid without changing the state of the navigator.
Definition TGeoNavigator.cxx:255
TGeoNavigator::GetMotherMatrix
TGeoHMatrix * GetMotherMatrix(Int_t up=1) const
Definition TGeoNavigator.h:158
TGeoNavigator::fCurrentVolume
TGeoVolume * fCurrentVolume
! current volume
Definition TGeoNavigator.h:76
TGeoNavigator::fLastNode
TGeoNode * fLastNode
! last searched node
Definition TGeoNavigator.h:79
TGeoNavigator::fThreadId
Int_t fThreadId
! thread id for this navigator
Definition TGeoNavigator.h:56
TGeoNavigator::GetPWSafetyEstimateFromCache
Double_t GetPWSafetyEstimateFromCache(Double_t cpoint[3]) const
Definition TGeoNavigator.h:210
TGeoNavigator::fDirection
Double_t fDirection[3]
! current direction
Definition TGeoNavigator.h:52
TGeoNavigator::PopDummy
void PopDummy(Int_t ipop=9999)
Definition TGeoNavigator.h:282
TGeoNavigator::GetTouchedCluster
Int_t GetTouchedCluster(Int_t start, Double_t *point, Int_t *check_list, Int_t ncheck, Int_t *result)
Make the cluster of overlapping nodes in a voxel, containing point in reference of the mother.
Definition TGeoNavigator.cxx:2324
TGeoNavigator::FindNextBoundary
TGeoNode * FindNextBoundary(Double_t stepmax=TGeoShape::Big(), const char *path="", Bool_t frombdr=kFALSE)
Find distance to next boundary and store it in fStep.
Definition TGeoNavigator.cxx:623
TGeoNavigator::FindNode
TGeoNode * FindNode(Bool_t safe_start=kTRUE)
Returns deepest node containing current point.
Definition TGeoNavigator.cxx:1622
TGeoNavigator::fOverlapMark
Int_t fOverlapMark
! current recursive position in fOverlapClusters
Definition TGeoNavigator.h:61
TGeoNavigator::FindNextBoundaryAndStep
TGeoNode * FindNextBoundaryAndStep(Double_t stepmax=TGeoShape::Big(), Bool_t compsafe=kFALSE)
Compute distance to next boundary within STEPMAX.
Definition TGeoNavigator.cxx:1210
TGeoNavigator::CdTop
void CdTop()
Make top level node the current node.
Definition TGeoNavigator.cxx:385
TGeoNavigator::fNmany
Int_t fNmany
! number of overlapping nodes on current branch
Definition TGeoNavigator.h:58
TGeoNavigator::fGeometry
TGeoManager * fGeometry
! current geometry
Definition TGeoNavigator.h:74
TGeoNavigator::fGlobalMatrix
TGeoHMatrix * fGlobalMatrix
! current pointer to cached global matrix
Definition TGeoNavigator.h:84
TGeoNavigator::MasterToTop
void MasterToTop(const Double_t *master, Double_t *top) const
Convert coordinates from master volume frame to top.
Definition TGeoNavigator.cxx:2771
TGeoNavigator::GetCurrentNodeId
Int_t GetCurrentNodeId() const
Definition TGeoNavigator.h:162
TGeoNavigator::FindNormalFast
Double_t * FindNormalFast()
Computes fast normal to next crossed boundary, assuming that the current point is close enough to the...
Definition TGeoNavigator.cxx:1672
TGeoNavigator::PushPath
Int_t PushPath(Int_t startlevel=0)
Definition TGeoNavigator.h:248
TGeoNavigator::fIsStepExiting
Bool_t fIsStepExiting
! flag that next geometric step will exit current volume
Definition TGeoNavigator.h:69
TGeoNavigator::fIsOnBoundary
Bool_t fIsOnBoundary
! flag that current point is on some boundary
Definition TGeoNavigator.h:71
TGeoNavigator::fgUsePWSafetyCaching
static Bool_t fgUsePWSafetyCaching
! global mode is caching enabled for parallel world safety calls
Definition TGeoNavigator.h:88
TGeoNavigator::GetBranchOnlys
void GetBranchOnlys(Int_t *isonly) const
Fill node copy numbers of current branch into an array.
Definition TGeoNavigator.cxx:451
TGeoNavigator::TopToMaster
void TopToMaster(const Double_t *top, Double_t *master) const
Convert coordinates from top volume frame to master.
Definition TGeoNavigator.cxx:2779
TGeoNavigator::SetCurrentPoint
void SetCurrentPoint(const Double_t *point)
Definition TGeoNavigator.h:174
TGeoNavigator::fCurrentMatrix
TGeoHMatrix * fCurrentMatrix
! current stored global matrix
Definition TGeoNavigator.h:83
TGeoNavigator::fOverlapClusters
Int_t * fOverlapClusters
! internal array for overlaps
Definition TGeoNavigator.h:62
TGeoNavigator::IsSafeStep
Bool_t IsSafeStep(Double_t proposed, Double_t &newsafety) const
In case a previous safety value was computed, check if the safety region is still safe for the curren...
Definition TGeoNavigator.cxx:2688
TGeoNavigator::SearchNode
TGeoNode * SearchNode(Bool_t downwards=kFALSE, const TGeoNode *skipnode=nullptr)
Returns the deepest node containing fPoint, which must be set a priori.
Definition TGeoNavigator.cxx:2003
TGeoNavigator::fLastSafety
Double_t fLastSafety
! last computed safety radius
Definition TGeoNavigator.h:47
TGeoNavigator::~TGeoNavigator
~TGeoNavigator() override
Destructor.
Definition TGeoNavigator.cxx:168
TGeoNavigator::TGeoNavigator
TGeoNavigator()
Constructor.
Definition TGeoNavigator.cxx:56
TGeoNavigator::SetCurrentDirection
void SetCurrentDirection(const Double_t *dir)
Definition TGeoNavigator.h:187
TGeoNavigator::BuildCache
void BuildCache(Bool_t dummy=kFALSE, Bool_t nodeid=kFALSE)
Builds the cache for physical nodes and global matrices.
Definition TGeoNavigator.cxx:181
TGeoNavigator::fNextDaughterIndex
Int_t fNextDaughterIndex
! next daughter index after FindNextBoundary
Definition TGeoNavigator.h:59
TGeoNavigator::fIsNullStep
Bool_t fIsNullStep
! flag that last geometric step was null
Definition TGeoNavigator.h:73
TGeoNavigator::CdNode
void CdNode(Int_t nodeid)
Change current path to point to the node having this id.
Definition TGeoNavigator.cxx:300
TGeoNavigator::fNextNode
TGeoNode * fNextNode
! next node that will be crossed
Definition TGeoNavigator.h:80
TGeoNavigator::fCldirChecked
Double_t fCldirChecked[3]
! unit vector to current checked shape
Definition TGeoNavigator.h:50
TGeoNavigator::fLevel
Int_t fLevel
! current geometry level;
Definition TGeoNavigator.h:57
TGeoNavigator::GetPWSafety
Double_t GetPWSafety(Double_t cpoint[3], Double_t saf_max)
Wrapper for getting the safety from the parallel world.
Definition TGeoNavigator.cxx:1735
TGeoNavigator::ResetAll
void ResetAll()
Reset the navigator.
Definition TGeoNavigator.cxx:2787
TGeoNavigator::fBackupState
TGeoCacheState * fBackupState
! backup state
Definition TGeoNavigator.h:82
TGeoNavigator::IsSamePoint
Bool_t IsSamePoint(Double_t x, Double_t y, Double_t z) const
Check if a new point with given coordinates is the same as the last located one.
Definition TGeoNavigator.cxx:2713
TGeoNavigator::fCurrentOverlapping
Bool_t fCurrentOverlapping
! flags the type of the current node
Definition TGeoNavigator.h:64
TGeoNavigator::fIsOutside
Bool_t fIsOutside
! flag that current point is outside geometry
Definition TGeoNavigator.h:70
TGeoNavigator::CdDown
void CdDown(Int_t index)
Make a daughter of current node current.
Definition TGeoNavigator.cxx:312
TGeoNavigator::fIsEntering
Bool_t fIsEntering
! flag if current step just got into a new node
Definition TGeoNavigator.h:66
TGeoNavigator::fForcedNode
TGeoNode * fForcedNode
! current point is supposed to be inside this node
Definition TGeoNavigator.h:81
TGeoNavigator::GetPath
const char * GetPath() const
Get path to the current node in the form /node0/node1/...
Definition TGeoNavigator.cxx:2761
TGeoNavigator::GetSafeLevel
Int_t GetSafeLevel() const
Go upwards the tree until a non-overlapping node.
Definition TGeoNavigator.cxx:2485
TGeoNavigator::fCurrentNode
TGeoNode * fCurrentNode
! current node
Definition TGeoNavigator.h:77
TGeoNavigator::fSafety
Double_t fSafety
! safety radius from current point
Definition TGeoNavigator.h:46
TGeoNavigator::FindNormal
Double_t * FindNormal(Bool_t forward=kTRUE)
Computes normal vector to the next surface that will be or was already crossed when propagating on a ...
Definition TGeoNavigator.cxx:1692
TGeoNavigator::GetBranchNames
void GetBranchNames(Int_t *names) const
Fill volume names of current branch into an array.
Definition TGeoNavigator.cxx:435
TGeoNodeCache
Special pool of reusable nodes.
Definition TGeoCache.h:56
TGeoNodeCache::IsDummy
Bool_t IsDummy() const
Definition TGeoCache.h:126
TGeoNodeCache::GetNode
TGeoNode * GetNode() const
Definition TGeoCache.h:116
TGeoNodeCache::CdNode
void CdNode(Int_t nodeid)
Change current path to point to the node having this id.
Definition TGeoCache.cxx:175
TGeoNodeCache::CdTop
void CdTop()
Definition TGeoCache.h:92
TGeoNodeCache::GetBranchOnlys
void GetBranchOnlys(Int_t *isonly) const
Fill copy numbers of current branch nodes.
Definition TGeoCache.cxx:328
TGeoNodeCache::GetPath
const char * GetPath()
Returns the current geometry path.
Definition TGeoCache.cxx:366
TGeoNodeCache::GetBranchNumbers
void GetBranchNumbers(Int_t *copyNumbers, Int_t *volumeNumbers) const
Fill copy numbers of current branch nodes.
Definition TGeoCache.cxx:317
TGeoNodeCache::GetInfo
TGeoStateInfo * GetInfo()
Get next state info pointer.
Definition TGeoCache.cxx:341
TGeoNodeCache::HasIdArray
Bool_t HasIdArray() const
Definition TGeoCache.h:125
TGeoNodeCache::CdDown
Bool_t CdDown(Int_t index)
Make daughter INDEX of current node the active state. Compute global matrix.
Definition TGeoCache.cxx:218
TGeoNodeCache::RestoreState
Bool_t RestoreState(Int_t &nmany, TGeoCacheState *state, Double_t *point=nullptr)
Pop next state/point from a backed-up state.
Definition TGeoCache.cxx:417
TGeoNodeCache::CdUp
void CdUp()
Make mother of current node the active state.
Definition TGeoCache.cxx:270
TGeoNodeCache::GetBranchNames
void GetBranchNames(Int_t *names) const
Fill names with current branch volume names (4 char - used by GEANT3 interface).
Definition TGeoCache.cxx:305
TGeoNodeCache::GetCurrentMatrix
TGeoHMatrix * GetCurrentMatrix() const
Definition TGeoCache.h:109
TGeoNodeCache::ReleaseInfo
void ReleaseInfo()
Release last used state info pointer.
Definition TGeoCache.cxx:358
TGeoNodeCache::GetLevel
Int_t GetLevel() const
Definition TGeoCache.h:121
TGeoNode
A node represent a volume positioned inside another.They store links to both volumes and to the TGeoM...
Definition TGeoNode.h:39
TGeoNode::GetMedium
TGeoMedium * GetMedium() const
Definition TGeoNode.h:90
TGeoNode::IsOverlapping
Bool_t IsOverlapping() const
Definition TGeoNode.h:108
TGeoNode::GetVolume
TGeoVolume * GetVolume() const
Definition TGeoNode.h:100
TGeoNode::IsOffset
Bool_t IsOffset() const
Definition TGeoNode.h:106
TGeoNode::GetNdaughters
Int_t GetNdaughters() const
Definition TGeoNode.h:92
TGeoNode::GetDaughter
TGeoNode * GetDaughter(Int_t ind) const
Definition TGeoNode.h:84
TGeoNode::GetMatrix
virtual TGeoMatrix * GetMatrix() const =0
TGeoNode::cd
virtual void cd() const
Definition TGeoNode.h:71
TGeoNode::GetOverlaps
Int_t * GetOverlaps(Int_t &novlp) const
Definition TGeoNode.h:95
TGeoNode::MasterToLocal
virtual void MasterToLocal(const Double_t *master, Double_t *local) const
Convert the point coordinates from mother reference to local reference system.
Definition TGeoNode.cxx:692
TGeoNode::GetIndex
virtual Int_t GetIndex() const
Definition TGeoNode.h:88
TGeoNode::GetFinder
virtual TGeoPatternFinder * GetFinder() const
Definition TGeoNode.h:89
TGeoNode::MasterToLocalVect
virtual void MasterToLocalVect(const Double_t *master, Double_t *local) const
Convert a vector from mother reference to local reference system.
Definition TGeoNode.cxx:700
TGeoNode::Safety
Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
computes the closest distance from given point to this shape
Definition TGeoNode.cxx:808
TGeoParallelWorld::Safety
Double_t Safety(Double_t point[3], Double_t safmax=1.E30)
Definition TGeoParallelWorld.h:118
TGeoParallelWorld::FindNextBoundary
TGeoPhysicalNode * FindNextBoundary(Double_t point[3], Double_t dir[3], Double_t &step, Double_t stepmax=1.E30)
Definition TGeoParallelWorld.h:130
TGeoParallelWorld::FindNode
TGeoPhysicalNode * FindNode(Double_t point[3])
Definition TGeoParallelWorld.h:106
TGeoPatternFinder
base finder class for patterns. A pattern is specifying a division type
Definition TGeoPatternFinder.h:25
TGeoPhysicalNode
Physical nodes are the actual 'touchable' objects in the geometry, representing a path of positioned ...
Definition TGeoPhysicalNode.h:29
TGeoShape::Big
static Double_t Big()
Definition TGeoShape.h:95
TGeoShape::DistFromInside
virtual Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const =0
TGeoShape::Safety
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const =0
TGeoShape::ComputeNormal
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) const =0
TGeoShape::DistFromOutside
virtual Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const =0
TGeoShape::Contains
virtual Bool_t Contains(const Double_t *point) const =0
TGeoShape::Tolerance
static Double_t Tolerance()
Definition TGeoShape.h:98
TGeoVolumeAssembly
Volume assemblies.
Definition TGeoVolume.h:317
TGeoVolume
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition TGeoVolume.h:43
TGeoVolume::GetNextNodeIndex
virtual Int_t GetNextNodeIndex() const
Definition TGeoVolume.h:169
TGeoVolume::IsActiveDaughters
Bool_t IsActiveDaughters() const
Definition TGeoVolume.h:147
TGeoVolume::Contains
Bool_t Contains(const Double_t *point) const
Definition TGeoVolume.h:105
TGeoVolume::GetNdaughters
Int_t GetNdaughters() const
Definition TGeoVolume.h:363
TGeoVolume::GetNodes
TObjArray * GetNodes()
Definition TGeoVolume.h:170
TGeoVolume::FindOverlaps
void FindOverlaps() const
loop all nodes marked as overlaps and find overlapping brothers
Definition TGeoVolume.cxx:2343
TGeoVolume::GetNode
TGeoNode * GetNode(const char *name) const
get the pointer to a daughter node
Definition TGeoVolume.cxx:2309
TGeoVolume::GetIndex
Int_t GetIndex(const TGeoNode *node) const
get index number for a given daughter
Definition TGeoVolume.cxx:1859
TGeoVolume::GetFinder
TGeoPatternFinder * GetFinder() const
Definition TGeoVolume.h:178
TGeoVolume::GetVoxels
TGeoVoxelFinder * GetVoxels() const
Getter for optimization structure.
Definition TGeoVolume.cxx:1918
TGeoVolume::GetShape
TGeoShape * GetShape() const
Definition TGeoVolume.h:191
TGeoVolume::GetCurrentNodeIndex
virtual Int_t GetCurrentNodeIndex() const
Definition TGeoVolume.h:168
TGeoVolume::IsAssembly
virtual Bool_t IsAssembly() const
Returns true if the volume is an assembly or a scaled assembly.
Definition TGeoVolume.cxx:1938
TGeoVolume::IsActive
Bool_t IsActive() const
Definition TGeoVolume.h:146
TGeoVoxelFinder
Finder class handling voxels.
Definition TGeoVoxelFinder.h:20
TNamed::GetName
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
TObjArray
An array of TObjects.
Definition TObjArray.h:31
TObjArray::GetEntriesFast
Int_t GetEntriesFast() const
Definition TObjArray.h:58
TObjArray::UncheckedAt
TObject * UncheckedAt(Int_t i) const
Definition TObjArray.h:90
TObjArray::Add
void Add(TObject *obj) override
Definition TObjArray.h:68
TObject::ClassName
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:224
TObject::Error
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1095
TObject::Fatal
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition TObject.cxx:1123
TObject::Info
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:1069
TString
Basic string class.
Definition TString.h:138
y
Double_t y[n]
Definition legend1.C:17
x
Double_t x[n]
Definition legend1.C:17
TMath::Sqrt
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:675
TMath::Min
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:197
TMath::Abs
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:122
TGeoStateInfo
Statefull info for the current geometry level.
Definition TGeoStateInfo.h:21