227{
230
232 EAxis targetHeaderAxis;
233 G4double targetHeaderMin, targetHeaderMax, targetHeaderNodeWidth;
234 G4int targetHeaderNoSlices;
236
237 G4double minSafety= previousMinSafety;
239 unsigned int checkedNum= 0;
240
241 fVoxelDepth++;
242
243
244 targetHeaderAxis = targetVoxelHeader->
GetAxis();
245 targetHeaderNoSlices = targetVoxelHeader->
GetNoSlices();
248
249 targetHeaderNodeWidth = (targetHeaderMax-targetHeaderMin)
250 / targetHeaderNoSlices;
251
252 G4double localCrd= localPoint(targetHeaderAxis);
253
254 const G4int candNodeNo=
G4int( (localCrd-targetHeaderMin)
255 / targetHeaderNodeWidth );
256
257
258#ifdef G4DEBUG_VOXELISATION
259 if( candNodeNo < 0 || candNodeNo > targetHeaderNoSlices-1 )
260 {
262 ed << " Potential ERROR."
263 <<
" Point is outside range of Voxel in current coordinate" <<
G4endl;
264 ed << " Node number of point " << localPoint
265 <<
"is outside the range. " <<
G4endl;
266 ed << " Voxel node Num= " << candNodeNo << " versus minimum= " << 0
267 <<
" and maximum= " << targetHeaderNoSlices-1 <<
G4endl;
268 ed << " Axis = " << targetHeaderAxis
269 <<
" No of slices = " << targetHeaderNoSlices <<
G4endl;
270 ed << " Local coord = " << localCrd
271 << " Voxel Min = " << targetHeaderMin
272 <<
" Max = " << targetHeaderMax <<
G4endl;
274 ed <<
" Current volume (physical) = " << currentPhysical.
GetName()
279 G4Exception(
"G4VoxelSafety::SafetyForVoxelHeader()",
"GeomNav1003",
281 "Point is outside range of Voxel in current coordinate");
282 }
283#endif
284
285 const G4int pointNodeNo =
286 std::max( 0, std::min( candNodeNo, targetHeaderNoSlices-1 ) );
287
288#ifdef G4VERBOSE
289 if( fVerbose > 2 )
290 {
292 G4cout <<
"**** G4VoxelSafety::SafetyForVoxelHeader " <<
G4endl;
293 G4cout <<
" Called at Depth = " << fVoxelDepth;
294 G4cout <<
" Calculated pointNodeNo= " << pointNodeNo
295 << " from position= " << localPoint(targetHeaderAxis)
296 << " min= " << targetHeaderMin
297 << " max= " << targetHeaderMax
298 << " width= " << targetHeaderNodeWidth
299 << " no-slices= " << targetHeaderNoSlices
300 <<
" axis= " << targetHeaderAxis <<
G4endl;
301 }
302 else if (fVerbose == 1)
303 {
304 G4cout <<
" VoxelSafety: Depth = " << fVoxelDepth
305 << " Number of Slices = " << targetHeaderNoSlices
306 <<
" Header (address) = " << targetVoxelHeader <<
G4endl;
307 }
308#endif
309
310
311
312 fVoxelAxisStack[fVoxelDepth] = targetHeaderAxis;
313 fVoxelNoSlicesStack[fVoxelDepth] = targetHeaderNoSlices;
314 fVoxelSliceWidthStack[fVoxelDepth] = targetHeaderNodeWidth;
315
316 fVoxelHeaderStack[fVoxelDepth] = pHeader;
317
318 G4int trialUp= -1, trialDown= -1;
320
321
322
323 G4int nextUp= pointNodeNo+1;
324 G4int nextDown= pointNodeNo-1;
325
326 G4int nextNodeNo= pointNodeNo;
328 distAxis= 0.0;
329
330 G4bool nextIsInside=
false;
331
332 G4double distMaxInterest= std::min( previousMinSafety, maxLength);
333
334
335
336
337 targetNodeNo= pointNodeNo;
338 do
339 {
341 fVoxelNodeNoStack[fVoxelDepth] = targetNodeNo;
342
343 checkedNum++;
344
345 sampleProxy = targetVoxelHeader->
GetSlice(targetNodeNo);
346
347#ifdef G4DEBUG_NAVIGATION
348 if( fVerbose > 2 )
349 {
350 G4cout <<
" -Checking node " << targetNodeNo
351 <<
" is proxy with address " << sampleProxy <<
G4endl;
352 }
353#endif
354
355 if ( sampleProxy == 0 )
356 {
358 ed << " Problem for node number= " << targetNodeNo
359 << " Number of slides = " << targetHeaderNoSlices
361 G4Exception(
"G4VoxelSafety::SafetyForVoxelHeader()",
"GeomNav0003",
363 "Problem sampleProxy is Zero. Failure in loop.");
364 }
365 else if ( sampleProxy->
IsNode() )
366 {
367 targetVoxelNode = sampleProxy->
GetNode();
368
369
370
372#ifdef G4DEBUG_NAVIGATION
373 if( fVerbose > 2 )
374 {
375 G4cout <<
" -- It is a Node ";
376 G4cout <<
" its safety= " << nodeSafety
377 << " our level Saf = " << ourSafety
378 <<
" MinSaf= " << minSafety <<
G4endl;
379 }
380#endif
381 ourSafety= std::min( ourSafety, nodeSafety );
382
385 }
386 else
387 {
389
391 distCombined_Sq = distUpperDepth_Sq + distAxis*distAxis;
392
393#ifdef G4DEBUG_NAVIGATION
394 if( fVerbose > 2 )
395 {
396 G4double distCombined= std::sqrt( distCombined_Sq );
397 G4double distUpperDepth= std::sqrt ( distUpperDepth_Sq );
399 G4cout <<
" Recurse to deal with next level, fVoxelDepth= "
400 << fVoxelDepth+1 <<
G4endl;
401 G4cout <<
" Distances: Upper= " << distUpperDepth
402 << " Axis= " << distAxis
403 <<
" Combined= " << distCombined <<
G4endl;
404 }
405#endif
406
407
408
410 maxLength, currentPhysical,
411 distCombined_Sq, minSafety);
412 ourSafety= std::min( ourSafety, headerSafety );
413
414#ifdef G4DEBUG_NAVIGATION
415 if( fVerbose > 2 )
416 {
417 G4cout <<
" >> Header safety = " << headerSafety
418 <<
" our level Saf = " << ourSafety <<
G4endl;
419 }
420#endif
423 }
424 minSafety= std::min( minSafety, ourSafety );
425
426
427
428
429
430 if( targetNodeNo >= pointNodeNo )
431 {
432 nextUp = trialUp;
433 distUp = targetHeaderMax-localCrd ;
434 if( distUp < 0.0 )
435 {
437 }
438#ifdef G4DEBUG_NAVIGATION
439 if( fVerbose > 2 )
440 {
442 }
443#endif
444 }
445
446 if( targetNodeNo <= pointNodeNo )
447 {
448 nextDown = trialDown;
449 distDown = localCrd-targetHeaderMin;
450 if( distDown < 0.0 )
451 {
453 }
454#ifdef G4DEBUG_NAVIGATION
455 if( fVerbose > 2 )
456 {
457 G4cout <<
" > Updated nextDown= " << nextDown <<
G4endl;
458 }
459#endif
460 }
461
462#ifdef G4DEBUG_NAVIGATION
463 if( fVerbose > 2 )
464 {
465 G4cout <<
" Node= " << pointNodeNo
466 << " Up: next= " << nextUp << " d# "
467 << nextUp - pointNodeNo
468 << " trialUp= " << trialUp << " d# "
469 << trialUp - pointNodeNo
470 << " Down: next= " << nextDown << " d# "
471 << targetNodeNo - nextDown
472 << " trialDwn= " << trialDown << " d# "
473 << targetNodeNo - trialDown
474 << " condition (next is Inside)= " << nextIsInside
476 }
477#endif
478
480 UpIsClosest= distUp < distDown;
481
482 if( UpIsClosest || (nextDown < 0) )
483 {
484 nextNodeNo=nextUp;
485 distAxis = distUp;
486 ++nextUp;
487#ifdef G4VERBOSE
488 if( fVerbose > 2 )
489 {
490 G4cout <<
" > Chose Up. Depth= " << fVoxelDepth
491 << " Nodes: next= " << nextNodeNo
492 << " new nextUp= " << nextUp
493 <<
" Dist = " << distAxis <<
G4endl;
494 }
495#endif
496 }
497 else
498 {
499 nextNodeNo=nextDown;
500 distAxis = distDown;
501 --nextDown;
502#ifdef G4VERBOSE
503 if( fVerbose > 2 )
504 {
505 G4cout <<
" > Chose Down. Depth= " << fVoxelDepth
506 << " Nodes: next= " << nextNodeNo
507 << " new nextDown= " << nextDown
508 <<
" Dist = " << distAxis <<
G4endl;
509 }
510#endif
511 }
512
513 nextIsInside = (nextNodeNo >= 0) && (nextNodeNo < targetHeaderNoSlices);
514 if( nextIsInside )
515 {
516 targetNodeNo= nextNodeNo;
517
518#ifdef G4DEBUG_NAVIGATION
519 G4bool bContinue= (distAxis<minSafety);
520 if( !bContinue )
521 {
522 if( fVerbose > 2 )
523 {
524 G4cout <<
" Can skip remaining at depth " << targetHeaderAxis
525 << " >> distAxis= " << distAxis
526 <<
" minSafety= " << minSafety <<
G4endl;
527 }
528 }
529#endif
530 }
531 else
532 {
533#ifdef G4DEBUG_NAVIGATION
534 if( fVerbose > 2)
535 {
537 G4cout <<
" VoxSaf> No more candidates: nodeDown= " << nextDown
538 << " nodeUp= " << nextUp
539 <<
" lastSlice= " << targetHeaderNoSlices <<
G4endl;
540 }
541#endif
542 }
543
544
545
546 distMaxInterest = std::min( minSafety, maxLength );
547
548 } while ( nextIsInside && ( distAxis*distAxis + distUpperDepth_Sq
549 < distMaxInterest*distMaxInterest ) );
550
551#ifdef G4VERBOSE
552 if( fVerbose > 0 )
553 {
554 G4cout <<
" Ended for targetNodeNo -- checked " << checkedNum <<
" out of "
555 << targetHeaderNoSlices <<
" slices." <<
G4endl;
556 G4cout <<
" ===== Returning from SafetyForVoxelHeader "
557 <<
" Depth= " << fVoxelDepth <<
G4endl
559 }
560#endif
561
562
563 fVoxelDepth--;
564
565 return ourSafety;
566}
G4int GetMaxEquivalentSliceNo() const
G4int GetMinEquivalentSliceNo() const
G4SmartVoxelNode * GetNode() const
G4SmartVoxelHeader * GetHeader() const
virtual G4GeometryType GetEntityType() const =0
G4double SafetyForVoxelNode(const G4SmartVoxelNode *curVoxelNode, const G4ThreeVector &localPoint)
std::ostringstream G4ExceptionDescription