| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 //
27 //
28 //
29 // John Allison 31st December 1997.
30 // Model for physical volumes.
31
32 #include "G4PhysicalVolumeModel.hh"
33
34 #include "G4VGraphicsScene.hh"
35 #include "G4VPhysicalVolume.hh"
36 #include "G4PhysicalVolumeStore.hh"
37 #include "G4VPVParameterisation.hh"
38 #include "G4LogicalVolume.hh"
39 #include "G4VSolid.hh"
40 #include "G4SubtractionSolid.hh"
41 #include "G4IntersectionSolid.hh"
42 #include "G4Material.hh"
43 #include "G4VisAttributes.hh"
44 #include "G4BoundingExtentScene.hh"
45 #include "G4TransportationManager.hh"
46 #include "G4Polyhedron.hh"
47 #include "HepPolyhedronProcessor.h"
48 #include "G4AttDefStore.hh"
49 #include "G4AttDef.hh"
50 #include "G4AttValue.hh"
51 #include "G4UnitsTable.hh"
52 #include "G4Vector3D.hh"
53 #include "G4Mesh.hh"
54
55 #include <sstream>
56 #include <iomanip>
57
58 #define G4warn G4cout
59
60 G4PhysicalVolumeModel::G4PhysicalVolumeModel
61 (G4VPhysicalVolume* pVPV
62 , G4int requestedDepth
63 , const G4Transform3D& modelTransform
64 , const G4ModelingParameters* pMP
65 , G4bool useFullExtent
66 , const std::vector<G4PhysicalVolumeNodeID>&
67 : G4VModel (pMP)
68 , fpTopPV (pVPV)
69 , fTopPVCopyNo (pVPV? pVPV->GetCopyNo():
70 , fRequestedDepth (requestedDepth)
71 , fUseFullExtent (useFullExtent)
72 , fTransform (modelTransform)
73 , fCurrentDepth (0)
74 , fpCurrentPV (fpTopPV)
75 , fCurrentPVCopyNo (fpTopPV? fpTopPV->GetCop
76 , fpCurrentLV (fpTopPV? fpTopPV->GetLog
77 , fpCurrentMaterial (fpCurrentLV? fpCurrentLV
78 , fCurrentTransform (modelTransform)
79 , fBaseFullPVPath (baseFullPVPath)
80 , fFullPVPath (fBaseFullPVPath)
81 , fAbort (false)
82 , fCurtailDescent (false)
83 , fpClippingSolid (0)
84 , fClippingMode (subtraction)
85 , fNClippers (0)
86 , fTotalTouchables (0)
87 {
88 fType = "G4PhysicalVolumeModel";
89
90 if (!fpTopPV) {
91
92 // In some circumstances creating an "empt
93 // allowed, so I have supressed the G4Exce
94 // be a problem we might have to re-instat
95 // be used except by visualisation experts
96 // where it is used simply to get a list o
97 // G4Exception
98 // ("G4PhysicalVolumeModel::G4PhysicalV
99 // "modeling0010", FatalException, "Nu
100
101 fTopPVName = "NULL";
102 fGlobalTag = "Empty";
103 fGlobalDescription = "G4PhysicalVolumeMode
104
105 } else {
106
107 fTopPVName = fpTopPV -> GetName ();
108 std::ostringstream oss;
109 oss << fpTopPV->GetName() << ':' << fpTopP
110 << " BasePath:" << fBaseFullPVPath;
111 fGlobalTag = oss.str();
112 fGlobalDescription = "G4PhysicalVolumeMode
113 CalculateExtent ();
114 }
115 }
116
117 G4PhysicalVolumeModel::~G4PhysicalVolumeModel
118 {
119 delete fpClippingSolid;
120 }
121
122 G4ModelingParameters::PVNameCopyNoPath G4Physi
123 (const std::vector<G4PhysicalVolumeNodeID>& pa
124 {
125 G4ModelingParameters::PVNameCopyNoPath PVNam
126 for (const auto& node: path) {
127 PVNameCopyNoPath.push_back
128 (G4ModelingParameters::PVNameCopyNo
129 (node.GetPhysicalVolume()->GetName(),node
130 }
131 return PVNameCopyNoPath;
132 }
133
134 G4String G4PhysicalVolumeModel::GetPVNamePathS
135 (const std::vector<G4PhysicalVolumeNodeID>& pa
136 // Converts to path string, e.g., " World 0 En
137 // Note leading space character.
138 {
139 std::ostringstream oss;
140 oss << path;
141 return oss.str();
142 }
143
144 void G4PhysicalVolumeModel::CalculateExtent ()
145 {
146 // To handle paramaterisations, set copy num
147 // to get extent right
148 G4VPVParameterisation* pP = fpTopPV -> GetPa
149 if (pP) {
150 fpTopPV -> SetCopyNo (fTopPVCopyNo);
151 G4VSolid* solid = pP -> ComputeSolid (fTop
152 solid -> ComputeDimensions (pP, fTopPVCopy
153 }
154 if (fUseFullExtent) {
155 fExtent = fpTopPV -> GetLogicalVolume () -
156 } else {
157 // Calculate extent of *drawn* volumes, i.
158 // invisible volumes, by traversing the wh
159 // this physical volume.
160 G4BoundingExtentScene beScene(this);
161 const G4int tempRequestedDepth = fRequeste
162 const G4Transform3D tempTransform = fTrans
163 const G4ModelingParameters* tempMP = fpMP;
164 fRequestedDepth = -1; // Always search to
165 fTransform = G4Transform3D(); // Extent i
166 G4ModelingParameters mParams
167 (0, // No default vis attributes ne
168 G4ModelingParameters::wf, // wireframe
169 true, // Global culling.
170 true, // Cull invisible volumes.
171 false, // Density culling.
172 0., // Density (not relevant if den
173 true, // Cull daughters of opaque mot
174 24); // No of sides (not relevant fo
175 mParams.SetSpecialMeshRendering(true); //
176 fpMP = &mParams;
177 DescribeYourselfTo (beScene);
178 fExtent = beScene.GetBoundingExtent();
179 fpMP = tempMP;
180 fTransform = tempTransform;
181 fRequestedDepth = tempRequestedDepth;
182 }
183 G4double radius = fExtent.GetExtentRadius();
184 if (radius < 0.) { // Nothing in the scene
185 fExtent = fpTopPV -> GetLogicalVolume () -
186 }
187 fExtent.Transform(fTransform);
188 }
189
190 void G4PhysicalVolumeModel::DescribeYourselfTo
191 (G4VGraphicsScene& sceneHandler)
192 {
193 if (!fpTopPV) {
194 G4Exception
195 ("G4PhysicalVolumeModel::DescribeYourselfT
196 "modeling0012", FatalException, "No model
197 return; // Should never reach here, but k
198 }
199
200 if (!fpMP) {
201 G4Exception
202 ("G4PhysicalVolumeModel::DescribeYourselfT
203 "modeling0013", FatalException, "No model
204 return; // Should never reach here, but k
205 }
206
207 fNClippers = 0;
208 G4DisplacedSolid* pSectionSolid = fpMP->GetS
209 G4DisplacedSolid* pCutawaySolid = fpMP->GetC
210 if (fpClippingSolid) fNClippers++;
211 if (pSectionSolid) fNClippers++;
212 if (pCutawaySolid) fNClippers++;
213 if (fNClippers > 1) {
214 G4ExceptionDescription ed;
215 ed << "More than one solid cutter/clipper:
216 if (fpClippingSolid) ed << "\nclipper in f
217 if (pSectionSolid) ed << "\nsectioner in
218 if (pCutawaySolid) ed << "\ncutaway in f
219 G4Exception("G4PhysicalVolumeModel::Descri
220 }
221
222 G4Transform3D startingTransformation = fTran
223
224 fNTouchables.clear(); // Keeps count of tou
225
226 VisitGeometryAndGetVisReps
227 (fpTopPV,
228 fRequestedDepth,
229 startingTransformation,
230 sceneHandler);
231
232 fTotalTouchables = 0;
233 for (const auto& entry : fNTouchables) {
234 fTotalTouchables += entry.second;
235 }
236
237 // Reset or clear data...
238 fCurrentDepth = 0;
239 fpCurrentPV = fpTopPV;
240 fCurrentPVCopyNo = fpTopPV->GetCopyNo();
241 fpCurrentLV = fpTopPV->GetLogicalVolum
242 fpCurrentMaterial = fpCurrentLV? fpCurrentLV
243 fFullPVPath = fBaseFullPVPath;
244 fDrawnPVPath.clear();
245 fAbort = false;
246 fCurtailDescent = false;
247 }
248
249 G4String G4PhysicalVolumeModel::GetCurrentTag
250 {
251 if (fpCurrentPV) {
252 std::ostringstream o;
253 o << fpCurrentPV -> GetCopyNo ();
254 return fpCurrentPV -> GetName () + ":" + o
255 }
256 else {
257 return "WARNING: NO CURRENT VOLUME - globa
258 }
259 }
260
261 G4String G4PhysicalVolumeModel::GetCurrentDesc
262 {
263 return "G4PhysicalVolumeModel " + GetCurrent
264 }
265
266 void G4PhysicalVolumeModel::VisitGeometryAndGe
267 (G4VPhysicalVolume* pVPV,
268 G4int requestedDepth,
269 const G4Transform3D& theAT,
270 G4VGraphicsScene& sceneHandler)
271 {
272 // Visits geometry structure to a given dept
273 // at given physical volume with given sta
274 // describes volumes to the scene handler.
275 // requestedDepth < 0 (default) implies full
276 // theAT is the Accumulated Transformation.
277
278 // Find corresponding logical volume and (la
279 // local variables to preserve re-entrancy.
280 G4LogicalVolume* pLV = pVPV -> GetLogicalVo
281 G4VSolid* pSol = nullptr;
282 G4Material* pMaterial = nullptr;
283
284 if (!(pVPV -> IsReplicated ())) {
285 // Non-replicated physical volume.
286 pSol = pLV -> GetSolid ();
287 pMaterial = pLV -> GetMaterial ();
288 DescribeAndDescend (pVPV, requestedDepth,
289 theAT, sceneHandler);
290 }
291 else {
292 // Replicated or parametrised physical vol
293 EAxis axis;
294 G4int nReplicas;
295 G4double width;
296 G4double offset;
297 G4bool consuming;
298 pVPV -> GetReplicationData (axis, nReplica
299 G4int nBegin = 0;
300 G4int nEnd = nReplicas;
301 if (fCurrentDepth == 0) { // i.e., top vol
302 nBegin = fTopPVCopyNo; // Describe only
303 nEnd = nBegin + 1; // specified by
304 }
305 G4VPVParameterisation* pP = pVPV -> GetPar
306 if (pP) { // Parametrised volume.
307 for (int n = nBegin; n < nEnd; n++) {
308 pSol = pP -> ComputeSolid (n, pVPV);
309 pP -> ComputeTransformation (n, pVPV);
310 pSol -> ComputeDimensions (pP, n, pVPV
311 pVPV -> SetCopyNo (n);
312 fCurrentPVCopyNo = n;
313 // Create a touchable of current paren
314 // fFullPVPath has not been updated ye
315 // corresponds to the parent.
316 G4PhysicalVolumeModelTouchable parentT
317 pMaterial = pP -> ComputeMaterial (n,
318 DescribeAndDescend (pVPV, requestedDep
319 theAT, sceneHandle
320 }
321 }
322 else { // Plain replicated volume. From
323 // The replica's positions are claculate
324 // Replication may occur along:
325 //
326 // o Cartesian axes (kXAxis,kYAxis,kZAxi
327 //
328 // The replications, of specified widt
329 // form (-width*(nReplicas-1)*0.5+n*wi
330 // for the case of kXAxis, and are unr
331 //
332 // o Radial axis (cylindrical polar) (kR
333 //
334 // The replications are cons/tubs sect
335 // and are unrotated.
336 // They have radii of width*n+offset t
337 // where n=0..nRepl
338 //
339 // o Phi axis (cylindrical polar) (kPhi)
340 // The replications are `phi sections'
341 // They have phi of offset+n*width to
342 // n=0..nReplicas-1
343 //
344 pSol = pLV -> GetSolid ();
345 pMaterial = pLV -> GetMaterial ();
346 G4ThreeVector originalTranslation = pVPV
347 G4RotationMatrix* pOriginalRotation = pV
348 G4double originalRMin = 0., originalRMax
349 if (axis == kRho && pSol->GetEntityType(
350 originalRMin = ((G4Tubs*)pSol)->GetInn
351 originalRMax = ((G4Tubs*)pSol)->GetOut
352 }
353 G4bool visualisable = true;
354 for (int n = nBegin; n < nEnd; n++) {
355 G4ThreeVector translation; // Identit
356 G4RotationMatrix rotation; // Identit
357 G4RotationMatrix* pRotation = 0;
358 switch (axis) {
359 default:
360 case kXAxis:
361 translation = G4ThreeVector (-widt
362 break;
363 case kYAxis:
364 translation = G4ThreeVector (0,-wi
365 break;
366 case kZAxis:
367 translation = G4ThreeVector (0,0,-
368 break;
369 case kRho:
370 if (pSol->GetEntityType() == "G4Tu
371 ((G4Tubs*)pSol)->SetInnerRadius(
372 ((G4Tubs*)pSol)->SetOuterRadius(
373 } else {
374 if (fpMP->IsWarning())
375 G4warn <<
376 "G4PhysicalVolumeModel::VisitG
377 "\n built-in replicated volum
378 << pSol->GetEntityType() <<
379 "-type\n solids (your solid \
380 << pSol->GetName() <<
381 "\") are not visualisable."
382 << G4endl;
383 visualisable = false;
384 }
385 break;
386 case kPhi:
387 rotation.rotateZ (-(offset+(n+0.5)
388 // Minus Sign because for the phys
389 // coordinate system rotation.
390 pRotation = &rotation;
391 break;
392 }
393 pVPV -> SetTranslation (translation);
394 pVPV -> SetRotation (pRotation);
395 pVPV -> SetCopyNo (n);
396 fCurrentPVCopyNo = n;
397 if (visualisable) {
398 DescribeAndDescend (pVPV, requestedD
399 theAT, sceneHand
400 }
401 }
402 // Restore originals...
403 pVPV -> SetTranslation (originalTranslat
404 pVPV -> SetRotation (pOriginalRotatio
405 if (axis == kRho && pSol->GetEntityType(
406 ((G4Tubs*)pSol)->SetInnerRadius(origin
407 ((G4Tubs*)pSol)->SetOuterRadius(origin
408 }
409 }
410 }
411 }
412
413 void G4PhysicalVolumeModel::DescribeAndDescend
414 (G4VPhysicalVolume* pVPV,
415 G4int requestedDepth,
416 G4LogicalVolume* pLV,
417 G4VSolid* pSol,
418 G4Material* pMaterial,
419 const G4Transform3D& theAT,
420 G4VGraphicsScene& sceneHandler)
421 {
422 // Maintain useful data members...
423 fpCurrentPV = pVPV;
424 fCurrentPVCopyNo = pVPV->GetCopyNo();
425 fpCurrentLV = pLV;
426 fpCurrentMaterial = pMaterial;
427
428 // Create a nodeID for use below - note the
429 G4int copyNo = fpCurrentPV->GetCopyNo();
430 auto nodeID = G4PhysicalVolumeNodeID
431 (fpCurrentPV,copyNo,fCurrentDepth,fCurrentTr
432
433 // Update full path of physical volumes...
434 fFullPVPath.push_back(nodeID);
435
436 const G4RotationMatrix objectRotation = pVPV
437 const G4ThreeVector& translation = pVPV
438 G4Transform3D theLT (G4Transform3D (objectRo
439
440 // Compute the accumulated transformation...
441 // Note that top volume's transformation rel
442 // coordinate system is specified in theAT =
443 // = fTransform (see DescribeYourselfTo), so
444 // volume's own transformation, which is onl
445 // mother, i.e., not relative to the world c
446 // not be accumulated.
447 G4Transform3D theNewAT (theAT);
448 if (fCurrentDepth != 0) theNewAT = theAT * t
449 fCurrentTransform = theNewAT;
450
451 const G4VisAttributes* pVisAttribs = pLV->Ge
452 // If the volume does not have any vis attr
453 G4VisAttributes* tempVisAtts = nullptr;
454 if (!pVisAttribs) {
455 if (fpMP->GetDefaultVisAttributes()) {
456 tempVisAtts = new G4VisAttributes(*fpMP-
457 } else {
458 tempVisAtts = new G4VisAttributes;
459 }
460 // The user may request /vis/viewer/set/co
461 if (fpMP->GetCBDAlgorithmNumber() == 1) {
462 // Algorithm 1: 3 parameters: Simple rai
463 if (fpMP->GetCBDParameters().size() != 3
464 G4Exception("G4PhysicalVolumeModelTouc
465 "modeling0014",
466 FatalErrorInArgument,
467 "Algorithm-parameter misma
468 } else {
469 const G4double d = pMaterial? pMateria
470 const G4double d0 = fpMP->GetCBDParame
471 const G4double d1 = fpMP->GetCBDParame
472 const G4double d2 = fpMP->GetCBDParame
473 if (d < d0) { // Density < d0 is invis
474 tempVisAtts->SetVisibility(false);
475 } else { // Intermediate densities are
476 G4double red, green, blue;
477 if (d < d1) {
478 red = (d1-d)/(d1-d0); green = (d-d
479 } else if (d < d2) {
480 red = 0.; green = (d2-d)/(d2-d1);
481 } else { // Density >= d2 is blue.
482 red = 0.; green = 0.; blue = 1.;
483 }
484 tempVisAtts->SetColour(G4Colour(red,
485 }
486 }
487 } else if (fpMP->GetCBDAlgorithmNumber() =
488 // Algorithm 2
489 // ...etc.
490 }
491 pVisAttribs = tempVisAtts;
492 }
493 // From here, can assume pVisAttribs is a va
494 // because PreAddSolid needs a vis attribute
495
496 // Check if vis attributes are to be modifie
497 const auto& vams = fpMP->GetVisAttributesMod
498 if (vams.size()) {
499 // OK, we have some VAMs (Vis Attributes M
500 for (const auto& vam: vams) {
501 const auto& vamPath = vam.GetPVNameCopyN
502 if (vamPath.size() == fFullPVPath.size()
503 // OK, we have a size match.
504 // Check the volume name/copy number p
505 auto iVAMNameCopyNo = vamPath.begin();
506 auto iPVNodeId = fFullPVPath.begin();
507 for (; iVAMNameCopyNo != vamPath.end()
508 if (!(
509 iVAMNameCopyNo->GetName() ==
510 iPVNodeId->GetPhysicalVolume()
511 iVAMNameCopyNo->GetCopyNo() ==
512 iPVNodeId->GetPhysicalVolume()
513 )) {
514 // This path element does NOT matc
515 break;
516 }
517 }
518 if (iVAMNameCopyNo == vamPath.end()) {
519 // OK, the paths match (the above lo
520 // Create a vis atts object for the
521 // It is static so that we may retur
522 static G4VisAttributes modifiedVisAt
523 // Initialise it with the current vi
524 modifiedVisAtts = *pVisAttribs;
525 pVisAttribs = &modifiedVisAtts;
526 const G4VisAttributes& transVisAtts
527 switch (vam.GetVisAttributesSignifie
528 case G4ModelingParameters::VASVisi
529 modifiedVisAtts.SetVisibility(tr
530 break;
531 case G4ModelingParameters::VASDaug
532 modifiedVisAtts.SetDaughtersInvi
533 (transVisAtts.IsDaughtersInvisib
534 break;
535 case G4ModelingParameters::VASColo
536 modifiedVisAtts.SetColour(transV
537 break;
538 case G4ModelingParameters::VASLine
539 modifiedVisAtts.SetLineStyle(tra
540 break;
541 case G4ModelingParameters::VASLine
542 modifiedVisAtts.SetLineWidth(tra
543 break;
544 case G4ModelingParameters::VASForc
545 if (transVisAtts.IsForceDrawingS
546 if (transVisAtts.GetForcedDraw
547 G4VisAttributes::wireframe
548 modifiedVisAtts.SetForceWire
549 }
550 }
551 break;
552 case G4ModelingParameters::VASForc
553 if (transVisAtts.IsForceDrawingS
554 if (transVisAtts.GetForcedDraw
555 G4VisAttributes::solid) {
556 modifiedVisAtts.SetForceSoli
557 }
558 }
559 break;
560 case G4ModelingParameters::VASForc
561 if (transVisAtts.IsForceDrawingS
562 if (transVisAtts.GetForcedDraw
563 G4VisAttributes::cloud) {
564 modifiedVisAtts.SetForceClou
565 }
566 }
567 break;
568 case G4ModelingParameters::VASForc
569 modifiedVisAtts.SetForceNumberOf
570 (transVisAtts.GetForcedNumberOfC
571 break;
572 case G4ModelingParameters::VASForc
573 if (transVisAtts.IsForceAuxEdgeV
574 modifiedVisAtts.SetForceAuxEdg
575 (transVisAtts.IsForcedAuxEdgeV
576 }
577 break;
578 case G4ModelingParameters::VASForc
579 modifiedVisAtts.SetForceLineSegm
580 (transVisAtts.GetForcedLineSegme
581 break;
582 }
583 }
584 }
585 }
586 }
587
588 // Check for special mesh rendering
589 if (fpMP->IsSpecialMeshRendering()) {
590 G4bool potentialG4Mesh = false;
591 if (fpMP->GetSpecialMeshVolumes().empty())
592 // No volumes specified - all are potent
593 potentialG4Mesh = true;
594 } else {
595 // Name and (optionally) copy number of
596 for (const auto& pvNameCopyNo: fpMP->Get
597 if (pVPV->GetName() == pvNameCopyNo.Ge
598 // We have a name match
599 if (pvNameCopyNo.GetCopyNo() < 0) {
600 // Any copy number is OK
601 potentialG4Mesh = true;
602 } else {
603 if (pVPV->GetCopyNo() == pvNameCop
604 // We have a name and copy numbe
605 potentialG4Mesh = true;
606 }
607 }
608 }
609 }
610 }
611 if (potentialG4Mesh) {
612 // Create - or at least attempt to creat
613 // out of this pVPV the type will be "in
614 G4Mesh mesh(pVPV,theNewAT);
615 if (mesh.GetMeshType() != G4Mesh::invali
616 // Create "artificial" nodeID to repre
617 G4int artCopyNo = 0;
618 auto artPV = mesh.GetParameterisedVolu
619 auto artDepth = fCurrentDepth + 1;
620 auto artNodeID = G4PhysicalVolumeNodeI
621 fFullPVPath.push_back(artNodeID);
622 fDrawnPVPath.push_back(artNodeID);
623 sceneHandler.AddCompound(mesh);
624 fFullPVPath.pop_back();
625 fDrawnPVPath.pop_back();
626 delete tempVisAtts; // Needs cleaning
627 return; // Mesh found and processed -
628 } // else continue processing
629 }
630 }
631
632 // Make decision to draw...
633 G4bool thisToBeDrawn = true;
634
635 // There are various reasons why this volume
636 // might not be drawn...
637 G4bool culling = fpMP->IsCulling();
638 G4bool cullingInvisible = fpMP->IsCullingInv
639 G4bool markedVisible
640 = pVisAttribs->IsVisible() && pVisAttribs->G
641 G4bool cullingLowDensity = fpMP->IsDensityCu
642 G4double density = pMaterial? pMaterial->Get
643 G4double densityCut = fpMP -> GetVisibleDens
644
645 // 1) Global culling is on....
646 if (culling) {
647 // 2) Culling of invisible volumes is on..
648 if (cullingInvisible) {
649 // 3) ...and the volume is marked not vi
650 if (!markedVisible) thisToBeDrawn = fals
651 }
652 // 4) Or culling of low density volumes is
653 if (cullingLowDensity) {
654 // 5) ...and density is less than cut va
655 if (density < densityCut) thisToBeDrawn
656 }
657 }
658 // 6) The user has asked for all further tra
659 if (fAbort) thisToBeDrawn = false;
660
661 // Set "drawn" flag (it was true by default)
662 nodeID.SetDrawn(thisToBeDrawn);
663
664 if (thisToBeDrawn) {
665
666 // Update path of drawn physical volumes..
667 fDrawnPVPath.push_back(nodeID);
668
669 if (fpMP->IsExplode() && fDrawnPVPath.size
670 // For top-level drawn volumes, explode
671 G4Transform3D centering = G4Translate3D(
672 G4Transform3D centred = centering.invers
673 G4Scale3D oldScale;
674 G4Rotate3D oldRotation;
675 G4Translate3D oldTranslation;
676 centred.getDecomposition(oldScale, oldRo
677 G4double explodeFactor = fpMP->GetExplod
678 G4Translate3D newTranslation =
679 G4Translate3D(explodeFactor * oldTranslation
680 explodeFactor * oldTranslation.dy(),
681 explodeFactor * oldTranslation.dz())
682 theNewAT = centering * newTranslation *
683 }
684
685 auto fullDepth = fCurrentDepth + (G4int)fB
686 fNTouchables[fullDepth]++; // Increment f
687
688 DescribeSolid (theNewAT, pSol, pVisAttribs
689
690 }
691
692 // Make decision to draw daughters, if any.
693 // reasons why daughters might not be drawn.
694
695 // First, reasons that do not depend on cull
696 G4int nDaughters = (G4int)pLV->GetNoDaughter
697 G4bool daughtersToBeDrawn = true;
698 // 1) There are no daughters...
699 if (!nDaughters) daughtersToBeDrawn = false;
700 // 2) We are at the limit if requested depth
701 else if (requestedDepth == 0) daughtersToBeD
702 // 3) The user has asked for all further tra
703 else if (fAbort) daughtersToBeDrawn = false;
704 // 4) The user has asked that the descent be
705 else if (fCurtailDescent) daughtersToBeDrawn
706
707 // Now, reasons that depend on culling polic
708 else {
709 G4bool daughtersInvisible = pVisAttribs->I
710 // Culling of covered daughters request.
711 // G4VSceneHandler::CreateModelingParamete
712 // parameters...
713 G4bool cullingCovered = fpMP->IsCullingCov
714 G4bool surfaceDrawing =
715 fpMP->GetDrawingStyle() == G4ModelingPar
716 fpMP->GetDrawingStyle() == G4ModelingPar
717 if (pVisAttribs->IsForceDrawingStyle()) {
718 switch (pVisAttribs->GetForcedDrawingSty
719 default:
720 case G4VisAttributes::wireframe: surface
721 case G4VisAttributes::solid: surfaceDraw
722 }
723 }
724 G4bool opaque = pVisAttribs->GetColour().G
725 // 5) Global culling is on....
726 if (culling) {
727 // 6) ..and culling of invisible volumes
728 if (cullingInvisible) {
729 // 7) ...and the mother requests daughters i
730 if (daughtersInvisible) daughtersToBeDrawn =
731 }
732 // 8) Or culling of covered daughters is
733 if (cullingCovered) {
734 // 9) ...and surface drawing is operating...
735 if (surfaceDrawing) {
736 // 10) ...but only if mother is visible...
737 if (thisToBeDrawn) {
738 // 11) ...and opaque...
739 if (opaque) daughtersToBeDrawn = false
740 }
741 }
742 }
743 }
744 }
745
746 if (daughtersToBeDrawn) {
747 for (G4int iDaughter = 0; iDaughter < nDau
748 // Store daughter pVPV in local variable
749 G4VPhysicalVolume* pDaughterVPV = pLV ->
750 // Descend the geometry structure recurs
751 fCurrentDepth++;
752 VisitGeometryAndGetVisReps
753 (pDaughterVPV, requestedDepth - 1, theNewAT,
754 fCurrentDepth--;
755 }
756 }
757
758 // Clean up
759 delete tempVisAtts;
760
761 // Reset for normal descending of next volum
762 fCurtailDescent = false;
763
764 // Pop item from paths physical volumes...
765 fFullPVPath.pop_back();
766 if (thisToBeDrawn) {
767 fDrawnPVPath.pop_back();
768 }
769 }
770
771 namespace
772 {
773 G4bool SubtractionBoundingLimits(const G4VSo
774 {
775 // Algorithm from G4SubtractionSolid::Boun
776 // Since it is unclear how the shape of th
777 // after subtraction, just return its orig
778 G4ThreeVector pMin, pMax;
779 const auto& pSolA = target;
780 pSolA->BoundingLimits(pMin,pMax);
781 // Check correctness of the bounding box
782 if (pMin.x() >= pMax.x() || pMin.y() >= pM
783 // Bad bounding box (min >= max)
784 // This signifies a subtraction of non-i
785 return false;
786 }
787 return true;
788 }
789
790 G4bool IntersectionBoundingLimits(const G4VS
791 {
792 // Algorithm from G4IntersectionSolid::Bou
793 G4ThreeVector pMin, pMax;
794 G4ThreeVector minA,maxA, minB,maxB;
795 const auto& pSolA = target;
796 const auto& pSolB = intersector;
797 pSolA->BoundingLimits(minA,maxA);
798 pSolB->BoundingLimits(minB,maxB);
799 pMin.set(std::max(minA.x(),minB.x()),
800 std::max(minA.y(),minB.y()),
801 std::max(minA.z(),minB.z()));
802 pMax.set(std::min(maxA.x(),maxB.x()),
803 std::min(maxA.y(),maxB.y()),
804 std::min(maxA.z(),maxB.z()));
805 if (pMin.x() >= pMax.x() || pMin.y() >= pM
806 // Bad bounding box (min >= max)
807 // This signifies a subtraction of non-i
808 return false;
809 }
810 return true;
811 }
812 }
813
814 void G4PhysicalVolumeModel::DescribeSolid
815 (const G4Transform3D& theAT,
816 G4VSolid* pSol,
817 const G4VisAttributes* pVisAttribs,
818 G4VGraphicsScene& sceneHandler)
819 {
820 G4DisplacedSolid* pSectionSolid = fpMP->GetS
821 G4DisplacedSolid* pCutawaySolid = fpMP->GetC
822
823 if (fNClippers <= 0 || fNClippers > 1) {
824
825 // Normal case - no clipping, etc. - or, i
826 sceneHandler.PreAddSolid (theAT, *pVisAttr
827 pSol -> DescribeYourselfTo (sceneHandler);
828 sceneHandler.PostAddSolid ();
829
830 } else { // fNClippers == 1
831
832 G4VSolid* pResultantSolid = nullpt
833 G4DisplacedSolid* pDisplacedSolid = nullpt
834
835 if (fpClippingSolid) {
836 pDisplacedSolid = new G4DisplacedSolid("
837 switch (fClippingMode) {
838 case subtraction:
839 if (SubtractionBoundingLimits(pSol))
840 pResultantSolid = new G4Subtractio
841 ("subtracted_clipped_solid", pSol,
842 }
843 break;
844 case intersection:
845 if (IntersectionBoundingLimits(pSol,
846 pResultantSolid = new G4Intersecti
847 ("intersected_clipped_solid", pSol
848 }
849 break;
850 }
851
852 } else if (pSectionSolid) {
853 pDisplacedSolid = new G4DisplacedSolid("
854 if (IntersectionBoundingLimits(pSol, pDi
855 pResultantSolid = new G4IntersectionSo
856 }
857
858 } else if (pCutawaySolid) {
859 pDisplacedSolid = new G4DisplacedSolid("
860 switch (fpMP->GetCutawayMode()) {
861 case G4ModelingParameters::cutawayUnio
862 if (SubtractionBoundingLimits(pSol))
863 pResultantSolid = new G4Subtractio
864 }
865 break;
866 case G4ModelingParameters::cutawayInte
867 if (IntersectionBoundingLimits(pSol,
868 pResultantSolid = new G4Intersecti
869 }
870 break;
871 }
872 }
873
874 if (pResultantSolid) {
875 sceneHandler.PreAddSolid (theAT, *pVisAt
876 pResultantSolid -> DescribeYourselfTo (s
877 sceneHandler.PostAddSolid ();
878 }
879
880 delete pResultantSolid;
881 delete pDisplacedSolid;
882 }
883 }
884
885 G4bool G4PhysicalVolumeModel::Validate (G4bool
886 {
887 // Not easy to see how to validate this sort o
888 // a check that a volume of the same name (fTo
889 // the geometry tree but under some circumstan
890 // time. Instead, let us simply check that the
891 // physical volume store.
892 const auto& pvStore = G4PhysicalVolumeStore:
893 auto iterator = find(pvStore->begin(),pvStor
894 if (iterator == pvStore->end()) {
895 if (warn) {
896 G4ExceptionDescription ed;
897 ed << "Attempt to validate a volume that
898 G4Exception("G4PhysicalVolumeModel::Vali
899 }
900 return false;
901 } else {
902 return true;
903 }
904 }
905
906 const std::map<G4String,G4AttDef>* G4PhysicalV
907 {
908 G4bool isNew;
909 std::map<G4String,G4AttDef>* store
910 = G4AttDefStore::GetInstance("G4Physical
911 if (isNew) {
912 (*store)["PVPath"] =
913 G4AttDef("PVPath","Physical Volume Path"
914 (*store)["BasePVPath"] =
915 G4AttDef("BasePVPath","Base Physical Vol
916 (*store)["LVol"] =
917 G4AttDef("LVol","Logical Volume","Physic
918 (*store)["Solid"] =
919 G4AttDef("Solid","Solid Name","Physics",
920 (*store)["EType"] =
921 G4AttDef("EType","Entity Type","Physics"
922 (*store)["DmpSol"] =
923 G4AttDef("DmpSol","Dump of Solid propert
924 (*store)["LocalTrans"] =
925 G4AttDef("LocalTrans","Local transformat
926 (*store)["LocalExtent"] =
927 G4AttDef("LocalExtent","Local extent of
928 (*store)["GlobalTrans"] =
929 G4AttDef("GlobalTrans","Global transform
930 (*store)["GlobalExtent"] =
931 G4AttDef("GlobalExtent","Global extent o
932 (*store)["Material"] =
933 G4AttDef("Material","Material Name","Phy
934 (*store)["Density"] =
935 G4AttDef("Density","Material Density","P
936 (*store)["State"] =
937 G4AttDef("State","Material State (enum u
938 (*store)["Radlen"] =
939 G4AttDef("Radlen","Material Radiation Le
940 (*store)["Region"] =
941 G4AttDef("Region","Cuts Region","Physics
942 (*store)["RootRegion"] =
943 G4AttDef("RootRegion","Root Region (0/1
944 }
945 return store;
946 }
947
948 static std::ostream& operator<< (std::ostream&
949 {
950 using namespace std;
951
952 G4Scale3D sc;
953 G4Rotate3D r;
954 G4Translate3D tl;
955 t.getDecomposition(sc, r, tl);
956
957 const int w = 10;
958
959 // Transformation itself
960 o << setw(w) << t.xx() << setw(w) << t.xy()
961 o << setw(w) << t.yx() << setw(w) << t.yy()
962 o << setw(w) << t.zx() << setw(w) << t.zy()
963
964 // Translation
965 o << "= translation:" << endl;
966 o << setw(w) << tl.dx() << setw(w) << tl.dy(
967
968 // Rotation
969 o << "* rotation:" << endl;
970 o << setw(w) << r.xx() << setw(w) << r.xy()
971 o << setw(w) << r.yx() << setw(w) << r.yy()
972 o << setw(w) << r.zx() << setw(w) << r.zy()
973
974 // Scale
975 o << "* scale:" << endl;
976 o << setw(w) << sc.xx() << setw(w) << sc.yy(
977
978 // Transformed axes
979 o << "Transformed axes:" << endl;
980 o << "x': " << r * G4Vector3D(1., 0., 0.) <<
981 o << "y': " << r * G4Vector3D(0., 1., 0.) <<
982 o << "z': " << r * G4Vector3D(0., 0., 1.) <<
983
984 return o;
985 }
986
987 std::vector<G4AttValue>* G4PhysicalVolumeModel
988 {
989 std::vector<G4AttValue>* values = new std::v
990
991 if (!fpCurrentLV) {
992 G4Exception
993 ("G4PhysicalVolumeModel::CreateCurrent
994 "modeling0004",
995 JustWarning,
996 "Current logical volume not defined."
997 return values;
998 }
999
1000 std::ostringstream oss; oss << fFullPVPath;
1001 values->push_back(G4AttValue("PVPath", oss.
1002
1003 oss.str(""); oss << fBaseFullPVPath;
1004 values->push_back(G4AttValue("BasePVPath",
1005
1006 values->push_back(G4AttValue("LVol", fpCurr
1007 G4VSolid* pSol = fpCurrentLV->GetSolid();
1008
1009 values->push_back(G4AttValue("Solid", pSol-
1010
1011 values->push_back(G4AttValue("EType", pSol-
1012
1013 oss.str(""); oss << '\n' << *pSol;
1014 values->push_back(G4AttValue("DmpSol", oss.
1015
1016 const G4RotationMatrix localRotation = fpCu
1017 const G4ThreeVector& localTranslation = fpC
1018 oss.str(""); oss << '\n' << G4Transform3D(l
1019 values->push_back(G4AttValue("LocalTrans",
1020
1021 oss.str(""); oss << '\n' << pSol->GetExtent
1022 values->push_back(G4AttValue("LocalExtent",
1023
1024 oss.str(""); oss << '\n' << fCurrentTransfo
1025 values->push_back(G4AttValue("GlobalTrans",
1026
1027 oss.str(""); oss << '\n' << (pSol->GetExten
1028 values->push_back(G4AttValue("GlobalExtent"
1029
1030 G4String matName = fpCurrentMaterial? fpCur
1031 values->push_back(G4AttValue("Material", ma
1032
1033 G4double matDensity = fpCurrentMaterial? fp
1034 values->push_back(G4AttValue("Density", G4B
1035
1036 G4State matState = fpCurrentMaterial? fpCur
1037 oss.str(""); oss << matState;
1038 values->push_back(G4AttValue("State", oss.s
1039
1040 G4double matRadlen = fpCurrentMaterial? fpC
1041 values->push_back(G4AttValue("Radlen", G4Be
1042
1043 G4Region* region = fpCurrentLV->GetRegion()
1044 G4String regionName = region? region->GetNa
1045 values->push_back(G4AttValue("Region", regi
1046
1047 oss.str(""); oss << fpCurrentLV->IsRootRegi
1048 values->push_back(G4AttValue("RootRegion",
1049
1050 return values;
1051 }
1052
1053 G4bool G4PhysicalVolumeModel::G4PhysicalVolum
1054 (const G4PhysicalVolumeModel::G4PhysicalVol
1055 {
1056 if (fpPV < right.fpPV) return true;
1057 if (fpPV == right.fpPV) {
1058 if (fCopyNo < right.fCopyNo) return true;
1059 if (fCopyNo == right.fCopyNo)
1060 return fNonCulledDepth < right.fNonCull
1061 }
1062 return false;
1063 }
1064
1065 G4bool G4PhysicalVolumeModel::G4PhysicalVolum
1066 (const G4PhysicalVolumeModel::G4PhysicalVol
1067 {
1068 if (fpPV != right.fpPV ||
1069 fCopyNo != right.fCopyNo ||
1070 fNonCulledDepth != right.fNonCulledDept
1071 fTransform != right.fTransform ||
1072 fDrawn != right.fDrawn) return
1073 return false;
1074 }
1075
1076 std::ostream& operator<<
1077 (std::ostream& os, const G4PhysicalVolumeMo
1078 {
1079 G4VPhysicalVolume* pPV = node.GetPhysicalVo
1080 if (pPV) {
1081 os << pPV->GetName()
1082 << ' ' << node.GetCopyNo()
1083 // << '[' << node.GetNonCulledDepth() <
1084 // << ':' << node.GetTransform()
1085 ;
1086 // os << " (";
1087 // if (!node.GetDrawn()) os << "not ";
1088 // os << "drawn)";
1089 } else {
1090 os << " (Null PV node)";
1091 }
1092 return os;
1093 }
1094
1095 std::ostream& operator<<
1096 (std::ostream& os, const std::vector<G4Physic
1097 {
1098 if (path.empty()) {
1099 os << " TOP";
1100 } else {
1101 for (const auto& nodeID: path) {
1102 os << ' ' << nodeID;
1103 }
1104 }
1105 return os;
1106 }
1107
1108 G4PhysicalVolumeModel::G4PhysicalVolumeModelT
1109 (const std::vector<G4PhysicalVolumeNodeID>& f
1110 fFullPVPath(fullPVPath) {}
1111
1112 const G4ThreeVector& G4PhysicalVolumeModel::G
1113 {
1114 size_t i = fFullPVPath.size() - depth - 1;
1115 if (i >= fFullPVPath.size()) {
1116 G4Exception("G4PhysicalVolumeModelTouchab
1117 "modeling0005",
1118 FatalErrorInArgument,
1119 "Index out of range. Asking for non-exist
1120 }
1121 static G4ThreeVector tempTranslation;
1122 tempTranslation = fFullPVPath[i].GetTransfo
1123 return tempTranslation;
1124 }
1125
1126 const G4RotationMatrix* G4PhysicalVolumeModel
1127 {
1128 size_t i = fFullPVPath.size() - depth - 1;
1129 if (i >= fFullPVPath.size()) {
1130 G4Exception("G4PhysicalVolumeModelTouchab
1131 "modeling0006",
1132 FatalErrorInArgument,
1133 "Index out of range. Asking for non-exist
1134 }
1135 static G4RotationMatrix tempRotation;
1136 tempRotation = fFullPVPath[i].GetTransform(
1137 return &tempRotation;
1138 }
1139
1140 G4VPhysicalVolume* G4PhysicalVolumeModel::G4P
1141 {
1142 size_t i = fFullPVPath.size() - depth - 1;
1143 if (i >= fFullPVPath.size()) {
1144 G4Exception("G4PhysicalVolumeModelTouchab
1145 "modeling0007",
1146 FatalErrorInArgument,
1147 "Index out of range. Asking for non-exist
1148 }
1149 return fFullPVPath[i].GetPhysicalVolume();
1150 }
1151
1152 G4VSolid* G4PhysicalVolumeModel::G4PhysicalVo
1153 {
1154 size_t i = fFullPVPath.size() - depth - 1;
1155 if (i >= fFullPVPath.size()) {
1156 G4Exception("G4PhysicalVolumeModelTouchab
1157 "modeling0008",
1158 FatalErrorInArgument,
1159 "Index out of range. Asking for non-exist
1160 }
1161 return fFullPVPath[i].GetPhysicalVolume()->
1162 }
1163
1164 G4int G4PhysicalVolumeModel::G4PhysicalVolume
1165 {
1166 size_t i = fFullPVPath.size() - depth - 1;
1167 if (i >= fFullPVPath.size()) {
1168 G4Exception("G4PhysicalVolumeModelTouchab
1169 "modeling0009",
1170 FatalErrorInArgument,
1171 "Index out of range. Asking for non-exist
1172 }
1173 return fFullPVPath[i].GetCopyNo();
1174 }
1175