Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4PhysicalVolumeModel Class Reference

#include <G4PhysicalVolumeModel.hh>

+ Inheritance diagram for G4PhysicalVolumeModel:

Classes

class  G4PhysicalVolumeModelTouchable
 
class  G4PhysicalVolumeNodeID
 
struct  TouchableProperties
 

Public Types

enum  { UNLIMITED = -1 }
 
enum  ClippingMode { subtraction , intersection }
 

Public Member Functions

 G4PhysicalVolumeModel (G4VPhysicalVolume *=0, G4int requestedDepth=UNLIMITED, const G4Transform3D &modelTransformation=G4Transform3D(), const G4ModelingParameters *=0, G4bool useFullExtent=false, const std::vector< G4PhysicalVolumeNodeID > &baseFullPVPath=std::vector< G4PhysicalVolumeNodeID >())
 
virtual ~G4PhysicalVolumeModel ()
 
void DescribeYourselfTo (G4VGraphicsScene &)
 
G4String GetCurrentDescription () const
 
G4String GetCurrentTag () const
 
G4VPhysicalVolumeGetTopPhysicalVolume () const
 
G4int GetRequestedDepth () const
 
const G4VSolidGetClippingSolid () const
 
G4int GetCurrentDepth () const
 
const G4Transform3DGetTransformation () const
 
G4VPhysicalVolumeGetCurrentPV () const
 
G4int GetCurrentPVCopyNo () const
 
G4LogicalVolumeGetCurrentLV () const
 
G4MaterialGetCurrentMaterial () const
 
const G4Transform3DGetCurrentTransform () const
 
const std::vector< G4PhysicalVolumeNodeID > & GetBaseFullPVPath () const
 
const std::vector< G4PhysicalVolumeNodeID > & GetFullPVPath () const
 
const std::vector< G4PhysicalVolumeNodeID > & GetDrawnPVPath () const
 
const std::map< G4String, G4AttDef > * GetAttDefs () const
 
std::vector< G4AttValue > * CreateCurrentAttValues () const
 
void SetRequestedDepth (G4int requestedDepth)
 
void SetClippingSolid (G4VSolid *pClippingSolid)
 
void SetClippingMode (ClippingMode mode)
 
G4bool Validate (G4bool warn)
 
void Abort () const
 
void CurtailDescent () const
 
void CalculateExtent ()
 
- Public Member Functions inherited from G4VModel
 G4VModel (const G4ModelingParameters *=0)
 
virtual ~G4VModel ()
 
virtual void DescribeYourselfTo (G4VGraphicsScene &)=0
 
const G4ModelingParametersGetModelingParameters () const
 
const G4StringGetType () const
 
virtual G4String GetCurrentDescription () const
 
virtual G4String GetCurrentTag () const
 
const G4VisExtentGetExtent () const
 
const G4StringGetGlobalDescription () const
 
const G4StringGetGlobalTag () const
 
void SetModelingParameters (const G4ModelingParameters *)
 
void SetExtent (const G4VisExtent &)
 
void SetType (const G4String &)
 
void SetGlobalDescription (const G4String &)
 
void SetGlobalTag (const G4String &)
 
virtual G4bool Validate (G4bool warn=true)
 

Static Public Member Functions

static G4ModelingParameters::PVNameCopyNoPath GetPVNameCopyNoPath (const std::vector< G4PhysicalVolumeNodeID > &)
 

Protected Member Functions

void VisitGeometryAndGetVisReps (G4VPhysicalVolume *, G4int requestedDepth, const G4Transform3D &, G4VGraphicsScene &)
 
void DescribeAndDescend (G4VPhysicalVolume *, G4int requestedDepth, G4LogicalVolume *, G4VSolid *, G4Material *, const G4Transform3D &, G4VGraphicsScene &)
 
virtual void DescribeSolid (const G4Transform3D &theAT, G4VSolid *pSol, const G4VisAttributes *pVisAttribs, G4VGraphicsScene &sceneHandler)
 

Protected Attributes

G4VPhysicalVolumefpTopPV
 
G4String fTopPVName
 
G4int fTopPVCopyNo
 
G4int fRequestedDepth
 
G4bool fUseFullExtent
 
G4Transform3D fTransform
 
G4int fCurrentDepth
 
G4VPhysicalVolumefpCurrentPV
 
G4int fCurrentPVCopyNo
 
G4LogicalVolumefpCurrentLV
 
G4MaterialfpCurrentMaterial
 
G4Transform3D fCurrentTransform
 
std::vector< G4PhysicalVolumeNodeIDfBaseFullPVPath
 
std::vector< G4PhysicalVolumeNodeIDfFullPVPath
 
std::vector< G4PhysicalVolumeNodeIDfDrawnPVPath
 
G4bool fAbort
 
G4bool fCurtailDescent
 
G4VSolidfpClippingSolid
 
ClippingMode fClippingMode
 
- Protected Attributes inherited from G4VModel
G4String fType
 
G4String fGlobalTag
 
G4String fGlobalDescription
 
G4VisExtent fExtent
 
const G4ModelingParametersfpMP
 

Detailed Description

Definition at line 82 of file G4PhysicalVolumeModel.hh.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
Enumerator
UNLIMITED 

Definition at line 86 of file G4PhysicalVolumeModel.hh.

◆ ClippingMode

Constructor & Destructor Documentation

◆ G4PhysicalVolumeModel()

G4PhysicalVolumeModel::G4PhysicalVolumeModel ( G4VPhysicalVolume pVPV = 0,
G4int  requestedDepth = UNLIMITED,
const G4Transform3D modelTransformation = G4Transform3D(),
const G4ModelingParameters pMP = 0,
G4bool  useFullExtent = false,
const std::vector< G4PhysicalVolumeNodeID > &  baseFullPVPath = std::vector<G4PhysicalVolumeNodeID>() 
)

Definition at line 64 of file G4PhysicalVolumeModel.cc.

71: G4VModel (pMP)
72, fpTopPV (pVPV)
73, fTopPVCopyNo (pVPV? pVPV->GetCopyNo(): 0)
74, fRequestedDepth (requestedDepth)
75, fUseFullExtent (useFullExtent)
76, fTransform (modelTransform)
77, fCurrentDepth (0)
82, fCurrentTransform (modelTransform)
83, fBaseFullPVPath (baseFullPVPath)
84, fAbort (false)
85, fCurtailDescent (false)
88{
89 fType = "G4PhysicalVolumeModel";
90
91 if (!fpTopPV) {
92
93 // In some circumstances creating an "empty" G4PhysicalVolumeModel is
94 // allowed, so I have supressed the G4Exception below. If it proves to
95 // be a problem we might have to re-instate it, but it is unlikley to
96 // be used except by visualisation experts. See, for example, /vis/list,
97 // where it is used simply to get a list of G4AttDefs.
98 // G4Exception
99 // ("G4PhysicalVolumeModel::G4PhysicalVolumeModel",
100 // "modeling0010", FatalException, "Null G4PhysicalVolumeModel pointer.");
101
102 fTopPVName = "NULL";
103 fGlobalTag = "Empty";
104 fGlobalDescription = "G4PhysicalVolumeModel " + fGlobalTag;
105
106 } else {
107
108 fTopPVName = fpTopPV -> GetName ();
109 std::ostringstream oss;
110 oss << fpTopPV->GetName() << ':' << fpTopPV->GetCopyNo()
111 << " BasePath:" << fBaseFullPVPath;
112 fGlobalTag = oss.str();
113 fGlobalDescription = "G4PhysicalVolumeModel " + fGlobalTag;
115 }
116}
G4Material * GetMaterial() const
G4VPhysicalVolume * fpTopPV
std::vector< G4PhysicalVolumeNodeID > fBaseFullPVPath
G4VPhysicalVolume * fpCurrentPV
G4String fGlobalDescription
Definition: G4VModel.hh:100
G4String fType
Definition: G4VModel.hh:98
G4String fGlobalTag
Definition: G4VModel.hh:99
G4LogicalVolume * GetLogicalVolume() const
virtual G4int GetCopyNo() const =0
const G4String & GetName() const

◆ ~G4PhysicalVolumeModel()

G4PhysicalVolumeModel::~G4PhysicalVolumeModel ( )
virtual

Definition at line 118 of file G4PhysicalVolumeModel.cc.

119{
120 delete fpClippingSolid;
121}

Member Function Documentation

◆ Abort()

void G4PhysicalVolumeModel::Abort ( ) const
inline

Definition at line 256 of file G4PhysicalVolumeModel.hh.

256{fAbort = true;}

◆ CalculateExtent()

void G4PhysicalVolumeModel::CalculateExtent ( )

Definition at line 135 of file G4PhysicalVolumeModel.cc.

136{
137 // To handle paramaterisations, set copy number and compute dimensions
138 // to get extent right
139 G4VPVParameterisation* pP = fpTopPV -> GetParameterisation ();
140 if (pP) {
141 fpTopPV -> SetCopyNo (fTopPVCopyNo);
142 G4VSolid* solid = pP -> ComputeSolid (fTopPVCopyNo, fpTopPV);
143 solid -> ComputeDimensions (pP, fTopPVCopyNo, fpTopPV);
144 }
145 if (fUseFullExtent) {
146 fExtent = fpTopPV -> GetLogicalVolume () -> GetSolid () -> GetExtent ();
147 } else {
148 // Calculate extent of *drawn* volumes, i.e., ignoring culled, e.g.,
149 // invisible volumes, by traversing the whole geometry hierarchy below
150 // this physical volume.
151 G4BoundingExtentScene beScene(this);
152 const G4int tempRequestedDepth = fRequestedDepth;
153 const G4Transform3D tempTransform = fTransform;
154 const G4ModelingParameters* tempMP = fpMP;
155 fRequestedDepth = -1; // Always search to all depths to define extent.
156 fTransform = G4Transform3D(); // Extent is in local cooridinates
158 (0, // No default vis attributes needed.
159 G4ModelingParameters::wf, // wireframe (not relevant for this).
160 true, // Global culling.
161 true, // Cull invisible volumes.
162 false, // Density culling.
163 0., // Density (not relevant if density culling false).
164 true, // Cull daughters of opaque mothers.
165 24); // No of sides (not relevant for this operation).
166 mParams.SetSpecialMeshRendering(true); // Avoids traversing parameterisations
167 fpMP = &mParams;
168 DescribeYourselfTo (beScene);
169 fExtent = beScene.GetBoundingExtent();
170 fpMP = tempMP;
171 fTransform = tempTransform;
172 fRequestedDepth = tempRequestedDepth;
173 }
175 if (radius < 0.) { // Nothing in the scene - revert to top extent
176 fExtent = fpTopPV -> GetLogicalVolume () -> GetSolid () -> GetExtent ();
177 }
179}
HepGeom::Transform3D G4Transform3D
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
void DescribeYourselfTo(G4VGraphicsScene &)
G4VisExtent fExtent
Definition: G4VModel.hh:101
const G4VisExtent & GetExtent() const
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:102
G4double GetExtentRadius() const
Definition: G4VisExtent.cc:75
G4VisExtent & Transform(const G4Transform3D &)
Definition: G4VisExtent.cc:102

Referenced by G4PhysicalVolumeModel(), and G4VVisCommandGeometrySet::Set().

◆ CreateCurrentAttValues()

std::vector< G4AttValue > * G4PhysicalVolumeModel::CreateCurrentAttValues ( ) const

Definition at line 892 of file G4PhysicalVolumeModel.cc.

893{
894 std::vector<G4AttValue>* values = new std::vector<G4AttValue>;
895
896 if (!fpCurrentLV) {
898 ("G4PhysicalVolumeModel::CreateCurrentAttValues",
899 "modeling0004",
901 "Current logical volume not defined.");
902 return values;
903 }
904
905 std::ostringstream oss; oss << fFullPVPath;
906 values->push_back(G4AttValue("PVPath", oss.str(),""));
907
908 oss.str(""); oss << fBaseFullPVPath;
909 values->push_back(G4AttValue("BasePVPath", oss.str(),""));
910
911 values->push_back(G4AttValue("LVol", fpCurrentLV->GetName(),""));
912 G4VSolid* pSol = fpCurrentLV->GetSolid();
913
914 values->push_back(G4AttValue("Solid", pSol->GetName(),""));
915
916 values->push_back(G4AttValue("EType", pSol->GetEntityType(),""));
917
918 oss.str(""); oss << '\n' << *pSol;
919 values->push_back(G4AttValue("DmpSol", oss.str(),""));
920
921 const G4RotationMatrix localRotation = fpCurrentPV->GetObjectRotationValue();
922 const G4ThreeVector& localTranslation = fpCurrentPV->GetTranslation();
923 oss.str(""); oss << '\n' << G4Transform3D(localRotation,localTranslation);
924 values->push_back(G4AttValue("LocalTrans", oss.str(),""));
925
926 oss.str(""); oss << '\n' << pSol->GetExtent() << std::endl;
927 values->push_back(G4AttValue("LocalExtent", oss.str(),""));
928
929 oss.str(""); oss << '\n' << fCurrentTransform;
930 values->push_back(G4AttValue("GlobalTrans", oss.str(),""));
931
932 oss.str(""); oss << '\n' << (pSol->GetExtent()).Transform(fCurrentTransform) << std::endl;
933 values->push_back(G4AttValue("GlobalExtent", oss.str(),""));
934
935 G4String matName = fpCurrentMaterial? fpCurrentMaterial->GetName(): G4String("No material");
936 values->push_back(G4AttValue("Material", matName,""));
937
939 values->push_back(G4AttValue("Density", G4BestUnit(matDensity,"Volumic Mass"),""));
940
942 oss.str(""); oss << matState;
943 values->push_back(G4AttValue("State", oss.str(),""));
944
946 values->push_back(G4AttValue("Radlen", G4BestUnit(matRadlen,"Length"),""));
947
948 G4Region* region = fpCurrentLV->GetRegion();
949 G4String regionName = region? region->GetName(): G4String("No region");
950 values->push_back(G4AttValue("Region", regionName,""));
951
952 oss.str(""); oss << fpCurrentLV->IsRootRegion();
953 values->push_back(G4AttValue("RootRegion", oss.str(),""));
954
955 return values;
956}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
G4State
Definition: G4Material.hh:110
@ kStateUndefined
Definition: G4Material.hh:110
#define G4BestUnit(a, b)
G4VSolid * GetSolid() const
G4bool IsRootRegion() const
G4Region * GetRegion() const
const G4String & GetName() const
G4double GetDensity() const
Definition: G4Material.hh:175
G4State GetState() const
Definition: G4Material.hh:176
G4double GetRadlen() const
Definition: G4Material.hh:215
const G4String & GetName() const
Definition: G4Material.hh:172
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
const G4String & GetName() const
const G4ThreeVector GetTranslation() const
G4RotationMatrix GetObjectRotationValue() const
G4String GetName() const
virtual G4GeometryType GetEntityType() const =0

Referenced by G4VSceneHandler::LoadAtts(), G4ASCIITreeSceneHandler::RequestPrimitives(), and G4VisCommandsTouchable::SetNewValue().

◆ CurtailDescent()

void G4PhysicalVolumeModel::CurtailDescent ( ) const
inline

Definition at line 259 of file G4PhysicalVolumeModel.hh.

259{fCurtailDescent = true;}

Referenced by G4ASCIITreeSceneHandler::RequestPrimitives().

◆ DescribeAndDescend()

void G4PhysicalVolumeModel::DescribeAndDescend ( G4VPhysicalVolume pVPV,
G4int  requestedDepth,
G4LogicalVolume pLV,
G4VSolid pSol,
G4Material pMaterial,
const G4Transform3D theAT,
G4VGraphicsScene sceneHandler 
)
protected

Definition at line 383 of file G4PhysicalVolumeModel.cc.

391{
392 // Maintain useful data members...
393 fpCurrentPV = pVPV;
394 fCurrentPVCopyNo = pVPV->GetCopyNo();
395 fpCurrentLV = pLV;
396 fpCurrentMaterial = pMaterial;
397
398 // Create a nodeID for use below - note the "drawn" flag is true
399 G4int copyNo = fpCurrentPV->GetCopyNo();
400 auto nodeID = G4PhysicalVolumeNodeID
402
403 // Update full path of physical volumes...
404 fFullPVPath.push_back(nodeID);
405
406 const G4RotationMatrix objectRotation = pVPV -> GetObjectRotationValue ();
407 const G4ThreeVector& translation = pVPV -> GetTranslation ();
408 G4Transform3D theLT (G4Transform3D (objectRotation, translation));
409
410 // Compute the accumulated transformation...
411 // Note that top volume's transformation relative to the world
412 // coordinate system is specified in theAT == startingTransformation
413 // = fTransform (see DescribeYourselfTo), so first time through the
414 // volume's own transformation, which is only relative to its
415 // mother, i.e., not relative to the world coordinate system, should
416 // not be accumulated.
417 G4Transform3D theNewAT (theAT);
418 if (fCurrentDepth != 0) theNewAT = theAT * theLT;
419 fCurrentTransform = theNewAT;
420
421 const G4VisAttributes* pVisAttribs = pLV->GetVisAttributes();
422 // If the volume does not have any vis attributes, create it.
423 G4VisAttributes* tempVisAtts = nullptr;
424 if (!pVisAttribs) {
426 tempVisAtts = new G4VisAttributes(*fpMP->GetDefaultVisAttributes());
427 } else {
428 tempVisAtts = new G4VisAttributes;
429 }
430 // The user may request /vis/viewer/set/colourByDensity.
431 if (fpMP->GetCBDAlgorithmNumber() == 1) {
432 // Algorithm 1: 3 parameters: Simple rainbow mapping.
433 if (fpMP->GetCBDParameters().size() != 3) {
434 G4Exception("G4PhysicalVolumeModelTouchable::DescribeAndDescend",
435 "modeling0014",
437 "Algorithm-parameter mismatch for Colour By Density");
438 } else {
439 const G4double d = pMaterial? pMaterial->GetDensity(): 0.;
440 const G4double d0 = fpMP->GetCBDParameters()[0]; // Invisible d < d0.
441 const G4double d1 = fpMP->GetCBDParameters()[1]; // Rainbow d0->d1->d2.
442 const G4double d2 = fpMP->GetCBDParameters()[2]; // Blue d > d2.
443 if (d < d0) { // Density < d0 is invisible.
444 tempVisAtts->SetVisibility(false);
445 } else { // Intermediate densities are on a spectrum.
446 G4double red, green, blue;
447 if (d < d1) {
448 red = (d1-d)/(d1-d0); green = (d-d0)/(d1-d0); blue = 0.;
449 } else if (d < d2) {
450 red = 0.; green = (d2-d)/(d2-d1); blue = (d-d1)/(d2-d1);
451 } else { // Density >= d2 is blue.
452 red = 0.; green = 0.; blue = 1.;
453 }
454 tempVisAtts->SetColour(G4Colour(red,green,blue));
455 }
456 }
457 } else if (fpMP->GetCBDAlgorithmNumber() == 2) {
458 // Algorithm 2
459 // ...etc.
460 }
461 pVisAttribs = tempVisAtts;
462 }
463 // From here, can assume pVisAttribs is a valid pointer. This is necessary
464 // because PreAddSolid needs a vis attributes object.
465
466 // Check if vis attributes are to be modified by a /vis/touchable/set/ command.
467 const auto& vams = fpMP->GetVisAttributesModifiers();
468 if (vams.size()) {
469 // OK, we have some VAMs (Vis Attributes Modifiers).
470 for (const auto& vam: vams) {
471 const auto& vamPath = vam.GetPVNameCopyNoPath();
472 if (vamPath.size() == fFullPVPath.size()) {
473 // OK, we have a size match.
474 // Check the volume name/copy number path.
475 auto iVAMNameCopyNo = vamPath.begin();
476 auto iPVNodeId = fFullPVPath.begin();
477 for (; iVAMNameCopyNo != vamPath.end(); ++iVAMNameCopyNo, ++iPVNodeId) {
478 if (!(
479 iVAMNameCopyNo->GetName() ==
480 iPVNodeId->GetPhysicalVolume()->GetName() &&
481 iVAMNameCopyNo->GetCopyNo() ==
482 iPVNodeId->GetPhysicalVolume()->GetCopyNo()
483 )) {
484 // This path element does NOT match.
485 break;
486 }
487 }
488 if (iVAMNameCopyNo == vamPath.end()) {
489 // OK, the paths match (the above loop terminated normally).
490 // Create a vis atts object for the modified vis atts.
491 // It is static so that we may return a reliable pointer to it.
492 static G4VisAttributes modifiedVisAtts;
493 // Initialise it with the current vis atts and reset the pointer.
494 modifiedVisAtts = *pVisAttribs;
495 pVisAttribs = &modifiedVisAtts;
496 const G4VisAttributes& transVisAtts = vam.GetVisAttributes();
497 switch (vam.GetVisAttributesSignifier()) {
499 modifiedVisAtts.SetVisibility(transVisAtts.IsVisible());
500 break;
502 modifiedVisAtts.SetDaughtersInvisible
503 (transVisAtts.IsDaughtersInvisible());
504 break;
506 modifiedVisAtts.SetColour(transVisAtts.GetColour());
507 break;
509 modifiedVisAtts.SetLineStyle(transVisAtts.GetLineStyle());
510 break;
512 modifiedVisAtts.SetLineWidth(transVisAtts.GetLineWidth());
513 break;
515 if (transVisAtts.IsForceDrawingStyle()) {
516 if (transVisAtts.GetForcedDrawingStyle() ==
518 modifiedVisAtts.SetForceWireframe(true);
519 }
520 }
521 break;
523 if (transVisAtts.IsForceDrawingStyle()) {
524 if (transVisAtts.GetForcedDrawingStyle() ==
526 modifiedVisAtts.SetForceSolid(true);
527 }
528 }
529 break;
531 if (transVisAtts.IsForceDrawingStyle()) {
532 if (transVisAtts.GetForcedDrawingStyle() ==
534 modifiedVisAtts.SetForceCloud(true);
535 }
536 }
537 break;
539 modifiedVisAtts.SetForceNumberOfCloudPoints
540 (transVisAtts.GetForcedNumberOfCloudPoints());
541 break;
543 if (transVisAtts.IsForceAuxEdgeVisible()) {
544 modifiedVisAtts.SetForceAuxEdgeVisible
545 (transVisAtts.IsForcedAuxEdgeVisible());
546 }
547 break;
549 modifiedVisAtts.SetForceLineSegmentsPerCircle
550 (transVisAtts.GetForcedLineSegmentsPerCircle());
551 break;
552 }
553 }
554 }
555 }
556 }
557
558 // Check for special mesh rendering
560 G4bool potentialG4Mesh = false;
561 if (fpMP->GetSpecialMeshVolumes().empty()) {
562 // No volumes specified - all are potentially possible
563 potentialG4Mesh = true;
564 } else {
565 // Name and (optionally) copy number of container volume is specified
566 for (const auto& pvNameCopyNo: fpMP->GetSpecialMeshVolumes()) {
567 if (pVPV->GetName() == pvNameCopyNo.GetName()) {
568 // We have a name match
569 if (pvNameCopyNo.GetCopyNo() < 0) {
570 // Any copy number is OK
571 potentialG4Mesh = true;
572 } else {
573 if (pVPV->GetCopyNo() == pvNameCopyNo.GetCopyNo()) {
574 // We have a name and copy number match
575 potentialG4Mesh = true;
576 }
577 }
578 }
579 }
580 }
581 if (potentialG4Mesh) {
582 // Create - or at least attempt to create - a mesh. If it cannot be created
583 // out of this pVPV the type will be "invalid".
584 G4Mesh mesh(pVPV,theNewAT);
585 if (mesh.GetMeshType() != G4Mesh::invalid) {
586 // Create "artificial" nodeID to represent the replaced volumes
587 G4int artCopyNo = 0;
588 auto artPV = mesh.GetParameterisedVolume();
589 auto artDepth = fCurrentDepth + 1;
590 auto artNodeID = G4PhysicalVolumeNodeID(artPV,artCopyNo,artDepth);
591 fFullPVPath.push_back(artNodeID);
592 fDrawnPVPath.push_back(artNodeID);
593 sceneHandler.AddCompound(mesh);
594 fFullPVPath.pop_back();
595 fDrawnPVPath.pop_back();
596 delete tempVisAtts; // Needs cleaning up (Coverity warning!!)
597 return; // Mesh found and processed - nothing more to do.
598 } // else continue processing
599 }
600 }
601
602 // Make decision to draw...
603 G4bool thisToBeDrawn = true;
604
605 // There are various reasons why this volume
606 // might not be drawn...
607 G4bool culling = fpMP->IsCulling();
608 G4bool cullingInvisible = fpMP->IsCullingInvisible();
609 G4bool markedVisible
610 = pVisAttribs->IsVisible() && pVisAttribs->GetColour().GetAlpha() > 0;
611 G4bool cullingLowDensity = fpMP->IsDensityCulling();
612 G4double density = pMaterial? pMaterial->GetDensity(): 0;
613 G4double densityCut = fpMP -> GetVisibleDensity ();
614
615 // 1) Global culling is on....
616 if (culling) {
617 // 2) Culling of invisible volumes is on...
618 if (cullingInvisible) {
619 // 3) ...and the volume is marked not visible...
620 if (!markedVisible) thisToBeDrawn = false;
621 }
622 // 4) Or culling of low density volumes is on...
623 if (cullingLowDensity) {
624 // 5) ...and density is less than cut value...
625 if (density < densityCut) thisToBeDrawn = false;
626 }
627 }
628 // 6) The user has asked for all further traversing to be aborted...
629 if (fAbort) thisToBeDrawn = false;
630
631 // Set "drawn" flag (it was true by default) - thisToBeDrawn may be false
632 nodeID.SetDrawn(thisToBeDrawn);
633
634 if (thisToBeDrawn) {
635
636 // Update path of drawn physical volumes...
637 fDrawnPVPath.push_back(nodeID);
638
639 if (fpMP->IsExplode() && fDrawnPVPath.size() == 1) {
640 // For top-level drawn volumes, explode along radius...
642 G4Transform3D centred = centering.inverse() * theNewAT;
643 G4Scale3D oldScale;
644 G4Rotate3D oldRotation;
645 G4Translate3D oldTranslation;
646 centred.getDecomposition(oldScale, oldRotation, oldTranslation);
647 G4double explodeFactor = fpMP->GetExplodeFactor();
648 G4Translate3D newTranslation =
649 G4Translate3D(explodeFactor * oldTranslation.dx(),
650 explodeFactor * oldTranslation.dy(),
651 explodeFactor * oldTranslation.dz());
652 theNewAT = centering * newTranslation * oldRotation * oldScale;
653 }
654
655 volumeCount++;
656 DescribeSolid (theNewAT, pSol, pVisAttribs, sceneHandler);
657
658 }
659
660 // Make decision to draw daughters, if any. There are various
661 // reasons why daughters might not be drawn...
662
663 // First, reasons that do not depend on culling policy...
664 G4int nDaughters = (G4int)pLV->GetNoDaughters();
665 G4bool daughtersToBeDrawn = true;
666 // 1) There are no daughters...
667 if (!nDaughters) daughtersToBeDrawn = false;
668 // 2) We are at the limit if requested depth...
669 else if (requestedDepth == 0) daughtersToBeDrawn = false;
670 // 3) The user has asked for all further traversing to be aborted...
671 else if (fAbort) daughtersToBeDrawn = false;
672 // 4) The user has asked that the descent be curtailed...
673 else if (fCurtailDescent) daughtersToBeDrawn = false;
674
675 // Now, reasons that depend on culling policy...
676 else {
677 G4bool daughtersInvisible = pVisAttribs->IsDaughtersInvisible();
678 // Culling of covered daughters request. This is computed in
679 // G4VSceneHandler::CreateModelingParameters() depending on view
680 // parameters...
681 G4bool cullingCovered = fpMP->IsCullingCovered();
682 G4bool surfaceDrawing =
685 if (pVisAttribs->IsForceDrawingStyle()) {
686 switch (pVisAttribs->GetForcedDrawingStyle()) {
687 default:
688 case G4VisAttributes::wireframe: surfaceDrawing = false; break;
689 case G4VisAttributes::solid: surfaceDrawing = true; break;
690 }
691 }
692 G4bool opaque = pVisAttribs->GetColour().GetAlpha() >= 1.;
693 // 5) Global culling is on....
694 if (culling) {
695 // 6) ..and culling of invisible volumes is on...
696 if (cullingInvisible) {
697 // 7) ...and the mother requests daughters invisible
698 if (daughtersInvisible) daughtersToBeDrawn = false;
699 }
700 // 8) Or culling of covered daughters is requested...
701 if (cullingCovered) {
702 // 9) ...and surface drawing is operating...
703 if (surfaceDrawing) {
704 // 10) ...but only if mother is visible...
705 if (thisToBeDrawn) {
706 // 11) ...and opaque...
707 if (opaque) daughtersToBeDrawn = false;
708 }
709 }
710 }
711 }
712 }
713
714 if (daughtersToBeDrawn) {
715 for (G4int iDaughter = 0; iDaughter < nDaughters; iDaughter++) {
716 // Store daughter pVPV in local variable ready for recursion...
717 G4VPhysicalVolume* pDaughterVPV = pLV -> GetDaughter (iDaughter);
718 // Descend the geometry structure recursively...
721 (pDaughterVPV, requestedDepth - 1, theNewAT, sceneHandler);
723 }
724 }
725
726 // Clean up
727 delete tempVisAtts;
728
729 // Reset for normal descending of next volume at this level...
730 fCurtailDescent = false;
731
732 // Pop item from paths physical volumes...
733 fFullPVPath.pop_back();
734 if (thisToBeDrawn) {
735 fDrawnPVPath.pop_back();
736 }
737}
@ FatalErrorInArgument
HepGeom::Translate3D G4Translate3D
bool G4bool
Definition: G4Types.hh:86
G4double GetAlpha() const
Definition: G4Colour.hh:155
const G4VisAttributes * GetVisAttributes() const
std::size_t GetNoDaughters() const
Definition: G4Mesh.hh:48
@ invalid
Definition: G4Mesh.hh:52
const G4VisAttributes * GetDefaultVisAttributes() const
const std::vector< VisAttributesModifier > & GetVisAttributesModifiers() const
const G4Point3D & GetExplodeCentre() const
G4bool IsCullingInvisible() const
const std::vector< PVNameCopyNo > & GetSpecialMeshVolumes() const
G4bool IsExplode() const
G4bool IsCulling() const
const std::vector< G4double > & GetCBDParameters() const
G4bool IsDensityCulling() const
G4bool IsSpecialMeshRendering() const
DrawingStyle GetDrawingStyle() const
G4double GetExplodeFactor() const
G4bool IsCullingCovered() const
G4int GetCBDAlgorithmNumber() const
void VisitGeometryAndGetVisReps(G4VPhysicalVolume *, G4int requestedDepth, const G4Transform3D &, G4VGraphicsScene &)
std::vector< G4PhysicalVolumeNodeID > fDrawnPVPath
virtual void DescribeSolid(const G4Transform3D &theAT, G4VSolid *pSol, const G4VisAttributes *pVisAttribs, G4VGraphicsScene &sceneHandler)
virtual void AddCompound(const G4VTrajectory &)=0
G4int GetForcedNumberOfCloudPoints() const
G4double GetLineWidth() const
G4bool IsDaughtersInvisible() const
void SetColour(const G4Colour &)
void SetVisibility(G4bool=true)
void SetForceAuxEdgeVisible(G4bool=true)
void SetForceCloud(G4bool=true)
G4int GetForcedLineSegmentsPerCircle() const
void SetForceWireframe(G4bool=true)
void SetLineWidth(G4double)
LineStyle GetLineStyle() const
const G4Colour & GetColour() const
G4bool IsVisible() const
G4bool IsForceAuxEdgeVisible() const
G4bool IsForcedAuxEdgeVisible() const
ForcedDrawingStyle GetForcedDrawingStyle() const
void SetForceSolid(G4bool=true)
void SetLineStyle(LineStyle)
void SetForceLineSegmentsPerCircle(G4int nSegments)
void SetDaughtersInvisible(G4bool=true)
void SetForceNumberOfCloudPoints(G4int nPoints)
G4bool IsForceDrawingStyle() const
double dy() const
Definition: Transform3D.h:287
double dz() const
Definition: Transform3D.h:290
double dx() const
Definition: Transform3D.h:284
void getDecomposition(Scale3D &scale, Rotate3D &rotation, Translate3D &translation) const
Definition: Transform3D.cc:173
Transform3D inverse() const
Definition: Transform3D.cc:141

Referenced by VisitGeometryAndGetVisReps().

◆ DescribeSolid()

void G4PhysicalVolumeModel::DescribeSolid ( const G4Transform3D theAT,
G4VSolid pSol,
const G4VisAttributes pVisAttribs,
G4VGraphicsScene sceneHandler 
)
protectedvirtual

Reimplemented in G4LogicalVolumeModel.

Definition at line 739 of file G4PhysicalVolumeModel.cc.

744{
745 G4DisplacedSolid* pSectionSolid = fpMP->GetSectionSolid();
746 G4DisplacedSolid* pCutawaySolid = fpMP->GetCutawaySolid();
747
748 if (!fpClippingSolid && !pSectionSolid && !pCutawaySolid) {
749
750 sceneHandler.PreAddSolid (theAT, *pVisAttribs);
751 pSol -> DescribeYourselfTo (sceneHandler); // Standard treatment.
752 sceneHandler.PostAddSolid ();
753
754 } else {
755
756 G4VSolid* pResultantSolid = nullptr;
757
758 if (fpClippingSolid) {
759 switch (fClippingMode) {
760 default:
761 case subtraction:
762 pResultantSolid = new G4SubtractionSolid
763 ("subtracted_clipped_solid", pSol, fpClippingSolid, theAT.inverse());
764 break;
765 case intersection:
766 pResultantSolid = new G4IntersectionSolid
767 ("intersected_clipped_solid", pSol, fpClippingSolid, theAT.inverse());
768 break;
769 }
770 }
771
772 if (pSectionSolid) {
773 pResultantSolid = new G4IntersectionSolid
774 ("sectioned_solid", pSol, pSectionSolid, theAT.inverse());
775 }
776
777 if (pCutawaySolid) {
778 pResultantSolid = new G4SubtractionSolid
779 ("cutaway_solid", pSol, pCutawaySolid, theAT.inverse());
780 }
781
782 sceneHandler.PreAddSolid (theAT, *pVisAttribs);
783 pResultantSolid -> DescribeYourselfTo (sceneHandler);
784 sceneHandler.PostAddSolid ();
785
786 delete pResultantSolid;
787 }
788}
G4DisplacedSolid * GetSectionSolid() const
G4DisplacedSolid * GetCutawaySolid() const
virtual void PostAddSolid()=0
virtual void PreAddSolid(const G4Transform3D &objectTransformation, const G4VisAttributes &visAttribs)=0

Referenced by DescribeAndDescend().

◆ DescribeYourselfTo()

void G4PhysicalVolumeModel::DescribeYourselfTo ( G4VGraphicsScene sceneHandler)
virtual

Implements G4VModel.

Definition at line 181 of file G4PhysicalVolumeModel.cc.

183{
184 if (!fpTopPV) G4Exception
185 ("G4PhysicalVolumeModel::DescribeYourselfTo",
186 "modeling0012", FatalException, "No model.");
187
188 if (!fpMP) G4Exception
189 ("G4PhysicalVolumeModel::DescribeYourselfTo",
190 "modeling0003", FatalException, "No modeling parameters.");
191
192 G4Transform3D startingTransformation = fTransform;
193
194 volumeCount = 0;
195
197 (fpTopPV,
199 startingTransformation,
200 sceneHandler);
201
202// G4cout
203// << "G4PhysicalVolumeModel::DescribeYourselfTo: volume count: "
204// << volumeCount
205// << G4endl;
206
207 // Reset or clear data...
208 fCurrentDepth = 0;
214 fDrawnPVPath.clear();
215 fAbort = false;
216 fCurtailDescent = false;
217}
@ FatalException

Referenced by G4VtkSceneHandler::AddCompound(), CalculateExtent(), DescribeSolid(), G4LogicalVolumeModel::DescribeYourselfTo(), G4VSceneHandler::Draw3DRectMeshAsDots(), G4VSceneHandler::Draw3DRectMeshAsSurfaces(), G4VisManager::DrawGeometry(), G4VSceneHandler::DrawTetMeshAsDots(), G4VSceneHandler::DrawTetMeshAsSurfaces(), G4ASCIITreeSceneHandler::EndModeling(), G4TouchableUtils::FindTouchableProperties(), G4VisCommandSceneAddLocalAxes::SetNewValue(), G4VisCommandSceneAddVolume::SetNewValue(), G4VisCommandSetTouchable::SetNewValue(), G4VisCommandSetVolumeForField::SetNewValue(), G4VisCommandsTouchable::SetNewValue(), and G4VisCommandViewerCentreOn::SetNewValue().

◆ GetAttDefs()

const std::map< G4String, G4AttDef > * G4PhysicalVolumeModel::GetAttDefs ( ) const

Definition at line 811 of file G4PhysicalVolumeModel.cc.

812{
813 G4bool isNew;
814 std::map<G4String,G4AttDef>* store
815 = G4AttDefStore::GetInstance("G4PhysicalVolumeModel", isNew);
816 if (isNew) {
817 (*store)["PVPath"] =
818 G4AttDef("PVPath","Physical Volume Path","Physics","","G4String");
819 (*store)["BasePVPath"] =
820 G4AttDef("BasePVPath","Base Physical Volume Path","Physics","","G4String");
821 (*store)["LVol"] =
822 G4AttDef("LVol","Logical Volume","Physics","","G4String");
823 (*store)["Solid"] =
824 G4AttDef("Solid","Solid Name","Physics","","G4String");
825 (*store)["EType"] =
826 G4AttDef("EType","Entity Type","Physics","","G4String");
827 (*store)["DmpSol"] =
828 G4AttDef("DmpSol","Dump of Solid properties","Physics","","G4String");
829 (*store)["LocalTrans"] =
830 G4AttDef("LocalTrans","Local transformation of volume","Physics","","G4String");
831 (*store)["LocalExtent"] =
832 G4AttDef("LocalExtent","Local extent of volume","Physics","","G4String");
833 (*store)["GlobalTrans"] =
834 G4AttDef("GlobalTrans","Global transformation of volume","Physics","","G4String");
835 (*store)["GlobalExtent"] =
836 G4AttDef("GlobalExtent","Global extent of volume","Physics","","G4String");
837 (*store)["Material"] =
838 G4AttDef("Material","Material Name","Physics","","G4String");
839 (*store)["Density"] =
840 G4AttDef("Density","Material Density","Physics","G4BestUnit","G4double");
841 (*store)["State"] =
842 G4AttDef("State","Material State (enum undefined,solid,liquid,gas)","Physics","","G4String");
843 (*store)["Radlen"] =
844 G4AttDef("Radlen","Material Radiation Length","Physics","G4BestUnit","G4double");
845 (*store)["Region"] =
846 G4AttDef("Region","Cuts Region","Physics","","G4String");
847 (*store)["RootRegion"] =
848 G4AttDef("RootRegion","Root Region (0/1 = false/true)","Physics","","G4bool");
849 }
850 return store;
851}
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)

Referenced by G4VSceneHandler::LoadAtts(), G4ASCIITreeSceneHandler::RequestPrimitives(), G4VisCommandList::SetNewValue(), and G4VisCommandsTouchable::SetNewValue().

◆ GetBaseFullPVPath()

const std::vector< G4PhysicalVolumeNodeID > & G4PhysicalVolumeModel::GetBaseFullPVPath ( ) const
inline

Definition at line 203 of file G4PhysicalVolumeModel.hh.

204 {return fBaseFullPVPath;}

◆ GetClippingSolid()

const G4VSolid * G4PhysicalVolumeModel::GetClippingSolid ( ) const
inline

Definition at line 179 of file G4PhysicalVolumeModel.hh.

180 {return fpClippingSolid;}

◆ GetCurrentDepth()

G4int G4PhysicalVolumeModel::GetCurrentDepth ( ) const
inline

◆ GetCurrentDescription()

G4String G4PhysicalVolumeModel::GetCurrentDescription ( ) const
virtual

Reimplemented from G4VModel.

Definition at line 231 of file G4PhysicalVolumeModel.cc.

232{
233 return "G4PhysicalVolumeModel " + GetCurrentTag ();
234}

◆ GetCurrentLV()

G4LogicalVolume * G4PhysicalVolumeModel::GetCurrentLV ( ) const
inline

◆ GetCurrentMaterial()

G4Material * G4PhysicalVolumeModel::GetCurrentMaterial ( ) const
inline

◆ GetCurrentPV()

◆ GetCurrentPVCopyNo()

G4int G4PhysicalVolumeModel::GetCurrentPVCopyNo ( ) const
inline

Definition at line 191 of file G4PhysicalVolumeModel.hh.

191{return fCurrentPVCopyNo;}

◆ GetCurrentTag()

G4String G4PhysicalVolumeModel::GetCurrentTag ( ) const
virtual

Reimplemented from G4VModel.

Definition at line 219 of file G4PhysicalVolumeModel.cc.

220{
221 if (fpCurrentPV) {
222 std::ostringstream o;
223 o << fpCurrentPV -> GetCopyNo ();
224 return fpCurrentPV -> GetName () + ":" + o.str();
225 }
226 else {
227 return "WARNING: NO CURRENT VOLUME - global tag is " + fGlobalTag;
228 }
229}

Referenced by GetCurrentDescription().

◆ GetCurrentTransform()

const G4Transform3D & G4PhysicalVolumeModel::GetCurrentTransform ( ) const
inline

◆ GetDrawnPVPath()

const std::vector< G4PhysicalVolumeNodeID > & G4PhysicalVolumeModel::GetDrawnPVPath ( ) const
inline

◆ GetFullPVPath()

const std::vector< G4PhysicalVolumeNodeID > & G4PhysicalVolumeModel::GetFullPVPath ( ) const
inline

◆ GetPVNameCopyNoPath()

G4ModelingParameters::PVNameCopyNoPath G4PhysicalVolumeModel::GetPVNameCopyNoPath ( const std::vector< G4PhysicalVolumeNodeID > &  path)
static

Definition at line 123 of file G4PhysicalVolumeModel.cc.

125{
127 for (const auto& node: path) {
128 PVNameCopyNoPath.push_back
130 (node.GetPhysicalVolume()->GetName(),node.GetCopyNo()));
131 }
132 return PVNameCopyNoPath;
133}
std::vector< PVNameCopyNo > PVNameCopyNoPath

Referenced by G4VViewer::TouchableSetColour(), G4VViewer::TouchableSetVisibility(), and G4VVisCommand::Twinkle().

◆ GetRequestedDepth()

G4int G4PhysicalVolumeModel::GetRequestedDepth ( ) const
inline

Definition at line 177 of file G4PhysicalVolumeModel.hh.

177{return fRequestedDepth;}

Referenced by G4ASCIITreeSceneHandler::EndModeling().

◆ GetTopPhysicalVolume()

G4VPhysicalVolume * G4PhysicalVolumeModel::GetTopPhysicalVolume ( ) const
inline

◆ GetTransformation()

const G4Transform3D & G4PhysicalVolumeModel::GetTransformation ( ) const
inline

Definition at line 185 of file G4PhysicalVolumeModel.hh.

185{return fTransform;}

◆ SetClippingMode()

void G4PhysicalVolumeModel::SetClippingMode ( ClippingMode  mode)
inline

Definition at line 249 of file G4PhysicalVolumeModel.hh.

249 {
250 fClippingMode = mode;
251 }

Referenced by G4VisCommandSceneAddVolume::SetNewValue().

◆ SetClippingSolid()

void G4PhysicalVolumeModel::SetClippingSolid ( G4VSolid pClippingSolid)
inline

Definition at line 245 of file G4PhysicalVolumeModel.hh.

245 {
246 fpClippingSolid = pClippingSolid;
247 }

Referenced by G4VisCommandSceneAddVolume::SetNewValue().

◆ SetRequestedDepth()

void G4PhysicalVolumeModel::SetRequestedDepth ( G4int  requestedDepth)
inline

Definition at line 241 of file G4PhysicalVolumeModel.hh.

241 {
242 fRequestedDepth = requestedDepth;
243 }

◆ Validate()

G4bool G4PhysicalVolumeModel::Validate ( G4bool  warn)
virtual

Reimplemented from G4VModel.

Definition at line 790 of file G4PhysicalVolumeModel.cc.

791{
792// Not easy to see how to validate this sort of model. Previously there was
793// a check that a volume of the same name (fTopPVName) existed somewhere in
794// the geometry tree but under some circumstances this consumed lots of CPU
795// time. Instead, let us simply check that the volume (fpTopPV) exists in the
796// physical volume store.
797 const auto& pvStore = G4PhysicalVolumeStore::GetInstance();
798 auto iterator = find(pvStore->begin(),pvStore->end(),fpTopPV);
799 if (iterator == pvStore->end()) {
800 if (warn) {
802 ed << "Attempt to validate a volume that is no longer in the physical volume store.";
803 G4Exception("G4PhysicalVolumeModel::Validate", "modeling0015", JustWarning, ed);
804 }
805 return false;
806 } else {
807 return true;
808 }
809}
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
static G4PhysicalVolumeStore * GetInstance()

Referenced by G4VisCommandSceneAddVolume::SetNewValue().

◆ VisitGeometryAndGetVisReps()

void G4PhysicalVolumeModel::VisitGeometryAndGetVisReps ( G4VPhysicalVolume pVPV,
G4int  requestedDepth,
const G4Transform3D theAT,
G4VGraphicsScene sceneHandler 
)
protected

Definition at line 236 of file G4PhysicalVolumeModel.cc.

241{
242 // Visits geometry structure to a given depth (requestedDepth), starting
243 // at given physical volume with given starting transformation and
244 // describes volumes to the scene handler.
245 // requestedDepth < 0 (default) implies full visit.
246 // theAT is the Accumulated Transformation.
247
248 // Find corresponding logical volume and (later) solid, storing in
249 // local variables to preserve re-entrancy.
250 G4LogicalVolume* pLV = pVPV -> GetLogicalVolume ();
251 G4VSolid* pSol = nullptr;
252 G4Material* pMaterial = nullptr;
253
254 if (!(pVPV -> IsReplicated ())) {
255 // Non-replicated physical volume.
256 pSol = pLV -> GetSolid ();
257 pMaterial = pLV -> GetMaterial ();
258 DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
259 theAT, sceneHandler);
260 }
261 else {
262 // Replicated or parametrised physical volume.
263 EAxis axis;
264 G4int nReplicas;
265 G4double width;
266 G4double offset;
267 G4bool consuming;
268 pVPV -> GetReplicationData (axis, nReplicas, width, offset, consuming);
269 G4int nBegin = 0;
270 G4int nEnd = nReplicas;
271 if (fCurrentDepth == 0) { // i.e., top volume
272 nBegin = fTopPVCopyNo; // Describe only one volume, namely the one
273 nEnd = nBegin + 1; // specified by the given copy number.
274 }
275 G4VPVParameterisation* pP = pVPV -> GetParameterisation ();
276 if (pP) { // Parametrised volume.
277 for (int n = nBegin; n < nEnd; n++) {
278 pSol = pP -> ComputeSolid (n, pVPV);
279 pP -> ComputeTransformation (n, pVPV);
280 pSol -> ComputeDimensions (pP, n, pVPV);
281 pVPV -> SetCopyNo (n);
283 // Create a touchable of current parent for ComputeMaterial.
284 // fFullPVPath has not been updated yet so at this point it
285 // corresponds to the parent.
286 G4PhysicalVolumeModelTouchable parentTouchable(fFullPVPath);
287 pMaterial = pP -> ComputeMaterial (n, pVPV, &parentTouchable);
288 DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
289 theAT, sceneHandler);
290 }
291 }
292 else { // Plain replicated volume. From geometry_guide.txt...
293 // The replica's positions are claculated by means of a linear formula.
294 // Replication may occur along:
295 //
296 // o Cartesian axes (kXAxis,kYAxis,kZAxis)
297 //
298 // The replications, of specified width have coordinates of
299 // form (-width*(nReplicas-1)*0.5+n*width,0,0) where n=0.. nReplicas-1
300 // for the case of kXAxis, and are unrotated.
301 //
302 // o Radial axis (cylindrical polar) (kRho)
303 //
304 // The replications are cons/tubs sections, centred on the origin
305 // and are unrotated.
306 // They have radii of width*n+offset to width*(n+1)+offset
307 // where n=0..nReplicas-1
308 //
309 // o Phi axis (cylindrical polar) (kPhi)
310 // The replications are `phi sections' or wedges, and of cons/tubs form
311 // They have phi of offset+n*width to offset+(n+1)*width where
312 // n=0..nReplicas-1
313 //
314 pSol = pLV -> GetSolid ();
315 pMaterial = pLV -> GetMaterial ();
316 G4ThreeVector originalTranslation = pVPV -> GetTranslation ();
317 G4RotationMatrix* pOriginalRotation = pVPV -> GetRotation ();
318 G4double originalRMin = 0., originalRMax = 0.;
319 if (axis == kRho && pSol->GetEntityType() == "G4Tubs") {
320 originalRMin = ((G4Tubs*)pSol)->GetInnerRadius();
321 originalRMax = ((G4Tubs*)pSol)->GetOuterRadius();
322 }
323 G4bool visualisable = true;
324 for (int n = nBegin; n < nEnd; n++) {
325 G4ThreeVector translation; // Identity.
326 G4RotationMatrix rotation; // Identity - life enough for visualizing.
327 G4RotationMatrix* pRotation = 0;
328 switch (axis) {
329 default:
330 case kXAxis:
331 translation = G4ThreeVector (-width*(nReplicas-1)*0.5+n*width,0,0);
332 break;
333 case kYAxis:
334 translation = G4ThreeVector (0,-width*(nReplicas-1)*0.5+n*width,0);
335 break;
336 case kZAxis:
337 translation = G4ThreeVector (0,0,-width*(nReplicas-1)*0.5+n*width);
338 break;
339 case kRho:
340 if (pSol->GetEntityType() == "G4Tubs") {
341 ((G4Tubs*)pSol)->SetInnerRadius(width*n+offset);
342 ((G4Tubs*)pSol)->SetOuterRadius(width*(n+1)+offset);
343 } else {
344 if (fpMP->IsWarning())
345 G4warn <<
346 "G4PhysicalVolumeModel::VisitGeometryAndGetVisReps: WARNING:"
347 "\n built-in replicated volumes replicated in radius for "
348 << pSol->GetEntityType() <<
349 "-type\n solids (your solid \""
350 << pSol->GetName() <<
351 "\") are not visualisable."
352 << G4endl;
353 visualisable = false;
354 }
355 break;
356 case kPhi:
357 rotation.rotateZ (-(offset+(n+0.5)*width));
358 // Minus Sign because for the physical volume we need the
359 // coordinate system rotation.
360 pRotation = &rotation;
361 break;
362 }
363 pVPV -> SetTranslation (translation);
364 pVPV -> SetRotation (pRotation);
365 pVPV -> SetCopyNo (n);
367 if (visualisable) {
368 DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
369 theAT, sceneHandler);
370 }
371 }
372 // Restore originals...
373 pVPV -> SetTranslation (originalTranslation);
374 pVPV -> SetRotation (pOriginalRotation);
375 if (axis == kRho && pSol->GetEntityType() == "G4Tubs") {
376 ((G4Tubs*)pSol)->SetInnerRadius(originalRMin);
377 ((G4Tubs*)pSol)->SetOuterRadius(originalRMax);
378 }
379 }
380 }
381}
#define G4warn
Definition: G4Scene.cc:41
CLHEP::Hep3Vector G4ThreeVector
#define G4endl
Definition: G4ios.hh:57
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:87
G4bool IsWarning() const
void DescribeAndDescend(G4VPhysicalVolume *, G4int requestedDepth, G4LogicalVolume *, G4VSolid *, G4Material *, const G4Transform3D &, G4VGraphicsScene &)
Definition: G4Tubs.hh:75
EAxis
Definition: geomdefs.hh:54
@ kPhi
Definition: geomdefs.hh:60
@ kYAxis
Definition: geomdefs.hh:56
@ kXAxis
Definition: geomdefs.hh:55
@ kZAxis
Definition: geomdefs.hh:57
@ kRho
Definition: geomdefs.hh:58

Referenced by DescribeAndDescend(), and DescribeYourselfTo().

Member Data Documentation

◆ fAbort

G4bool G4PhysicalVolumeModel::fAbort
mutableprotected

Definition at line 303 of file G4PhysicalVolumeModel.hh.

Referenced by Abort(), DescribeAndDescend(), and DescribeYourselfTo().

◆ fBaseFullPVPath

std::vector<G4PhysicalVolumeNodeID> G4PhysicalVolumeModel::fBaseFullPVPath
protected

◆ fClippingMode

ClippingMode G4PhysicalVolumeModel::fClippingMode
protected

Definition at line 306 of file G4PhysicalVolumeModel.hh.

Referenced by DescribeSolid(), and SetClippingMode().

◆ fCurrentDepth

G4int G4PhysicalVolumeModel::fCurrentDepth
protected

◆ fCurrentPVCopyNo

G4int G4PhysicalVolumeModel::fCurrentPVCopyNo
protected

◆ fCurrentTransform

G4Transform3D G4PhysicalVolumeModel::fCurrentTransform
protected

◆ fCurtailDescent

G4bool G4PhysicalVolumeModel::fCurtailDescent
mutableprotected

◆ fDrawnPVPath

std::vector<G4PhysicalVolumeNodeID> G4PhysicalVolumeModel::fDrawnPVPath
protected

◆ fFullPVPath

◆ fpClippingSolid

G4VSolid* G4PhysicalVolumeModel::fpClippingSolid
protected

◆ fpCurrentLV

G4LogicalVolume* G4PhysicalVolumeModel::fpCurrentLV
protected

◆ fpCurrentMaterial

G4Material* G4PhysicalVolumeModel::fpCurrentMaterial
protected

◆ fpCurrentPV

G4VPhysicalVolume* G4PhysicalVolumeModel::fpCurrentPV
protected

◆ fpTopPV

◆ fRequestedDepth

G4int G4PhysicalVolumeModel::fRequestedDepth
protected

◆ fTopPVCopyNo

G4int G4PhysicalVolumeModel::fTopPVCopyNo
protected

Definition at line 289 of file G4PhysicalVolumeModel.hh.

Referenced by CalculateExtent(), and VisitGeometryAndGetVisReps().

◆ fTopPVName

G4String G4PhysicalVolumeModel::fTopPVName
protected

Definition at line 288 of file G4PhysicalVolumeModel.hh.

Referenced by G4PhysicalVolumeModel().

◆ fTransform

G4Transform3D G4PhysicalVolumeModel::fTransform
protected

◆ fUseFullExtent

G4bool G4PhysicalVolumeModel::fUseFullExtent
protected

Definition at line 292 of file G4PhysicalVolumeModel.hh.

Referenced by CalculateExtent().


The documentation for this class was generated from the following files: