Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Data Structures | Public Types | Public Member Functions | Protected Member Functions | Protected Attributes
G4PhysicalVolumeModel Class Reference

#include <G4PhysicalVolumeModel.hh>

Inheritance diagram for G4PhysicalVolumeModel:
G4VModel G4LogicalVolumeModel

Data Structures

class  G4PhysicalVolumeModelTouchable
 
class  G4PhysicalVolumeNodeID
 

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)
 
virtual ~G4PhysicalVolumeModel ()
 
void DescribeYourselfTo (G4VGraphicsScene &)
 
G4String GetCurrentDescription () const
 
G4String GetCurrentTag () const
 
G4VPhysicalVolumeGetTopPhysicalVolume () const
 
G4int GetRequestedDepth () const
 
const G4VSolidGetClippingSolid () const
 
G4int GetCurrentDepth () const
 
G4VPhysicalVolumeGetCurrentPV () const
 
G4LogicalVolumeGetCurrentLV () const
 
G4MaterialGetCurrentMaterial () 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 SetBaseFullPVPath (const std::vector< G4PhysicalVolumeNodeID > &baseFullPVPath)
 
void SetRequestedDepth (G4int requestedDepth)
 
void SetClippingSolid (G4VSolid *pClippingSolid)
 
void SetClippingMode (ClippingMode mode)
 
G4bool Validate (G4bool warn)
 
void CurtailDescent ()
 
- Public Member Functions inherited from G4VModel
 G4VModel (const G4Transform3D &modelTransformation=G4Transform3D(), const G4ModelingParameters *=0)
 
virtual ~G4VModel ()
 
const G4ModelingParametersGetModelingParameters () const
 
const G4StringGetType () const
 
const G4VisExtentGetExtent () const
 
const G4StringGetGlobalDescription () const
 
const G4StringGetGlobalTag () const
 
const G4Transform3DGetTransformation () const
 
void SetModelingParameters (const G4ModelingParameters *)
 
void SetExtent (const G4VisExtent &)
 
void SetType (const G4String &)
 
void SetGlobalDescription (const G4String &)
 
void SetGlobalTag (const G4String &)
 
void SetTransformation (const G4Transform3D &)
 

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)
 
void CalculateExtent ()
 

Protected Attributes

G4VPhysicalVolumefpTopPV
 
G4String fTopPVName
 
G4int fTopPVCopyNo
 
G4int fRequestedDepth
 
G4bool fUseFullExtent
 
G4int fCurrentDepth
 
G4VPhysicalVolumefpCurrentPV
 
G4LogicalVolumefpCurrentLV
 
G4MaterialfpCurrentMaterial
 
G4Transform3DfpCurrentTransform
 
std::vector
< G4PhysicalVolumeNodeID
fBaseFullPVPath
 
std::vector
< G4PhysicalVolumeNodeID
fFullPVPath
 
std::vector
< G4PhysicalVolumeNodeID
fDrawnPVPath
 
G4bool fCurtailDescent
 
G4VSolidfpClippingSolid
 
ClippingMode fClippingMode
 
- Protected Attributes inherited from G4VModel
G4String fType
 
G4String fGlobalTag
 
G4String fGlobalDescription
 
G4VisExtent fExtent
 
G4Transform3D fTransform
 
const G4ModelingParametersfpMP
 

Detailed Description

Definition at line 82 of file G4PhysicalVolumeModel.hh.

Member Enumeration Documentation

anonymous enum
Enumerator
UNLIMITED 

Definition at line 86 of file G4PhysicalVolumeModel.hh.

Constructor & Destructor Documentation

G4PhysicalVolumeModel::G4PhysicalVolumeModel ( G4VPhysicalVolume pVPV = 0,
G4int  requestedDepth = UNLIMITED,
const G4Transform3D modelTransformation = G4Transform3D(),
const G4ModelingParameters pMP = 0,
G4bool  useFullExtent = false 
)

Definition at line 59 of file G4PhysicalVolumeModel.cc.

63  :
64  G4VModel (modelTransformation, pMP),
65  fpTopPV (pVPV),
66  fTopPVCopyNo (0),
67  fRequestedDepth (requestedDepth),
68  fUseFullExtent (useFullExtent),
69  fCurrentDepth (0),
70  fpCurrentPV (0),
71  fpCurrentLV (0),
74  fCurtailDescent (false),
75  fpClippingSolid (0),
77 {
78  if (fpTopPV) {
79 
80  fType = "G4PhysicalVolumeModel";
81  fTopPVName = fpTopPV -> GetName ();
82  fTopPVCopyNo = fpTopPV -> GetCopyNo ();
83  std::ostringstream o;
84  o << fpTopPV -> GetCopyNo ();
85  fGlobalTag = fpTopPV -> GetName () + "." + o.str();
86  fGlobalDescription = "G4PhysicalVolumeModel " + fGlobalTag;
87 
88  CalculateExtent ();
89  }
90 }
G4String fType
Definition: G4VModel.hh:108
G4Transform3D * fpCurrentTransform
G4VPhysicalVolume * fpCurrentPV
G4String fGlobalTag
Definition: G4VModel.hh:109
G4VModel(const G4Transform3D &modelTransformation=G4Transform3D(), const G4ModelingParameters *=0)
Definition: G4VModel.cc:38
G4String fGlobalDescription
Definition: G4VModel.hh:110
G4VPhysicalVolume * fpTopPV
G4PhysicalVolumeModel::~G4PhysicalVolumeModel ( )
virtual

Definition at line 92 of file G4PhysicalVolumeModel.cc.

References fpClippingSolid.

93 {
94  delete fpClippingSolid;
95 }

Member Function Documentation

void G4PhysicalVolumeModel::CalculateExtent ( )
protected

Definition at line 97 of file G4PhysicalVolumeModel.cc.

References DescribeYourselfTo(), G4VModel::fExtent, G4VModel::fpMP, fpTopPV, fRequestedDepth, G4VModel::fTransform, fUseFullExtent, G4BoundingSphereScene::GetCentre(), G4VModel::GetExtent(), G4BoundingSphereScene::GetRadius(), HepGeom::Transform3D::inverse(), and G4ModelingParameters::wf.

Referenced by Validate().

98 {
99  if (fUseFullExtent) {
100  fExtent = fpTopPV -> GetLogicalVolume () -> GetSolid () -> GetExtent ();
101  }
102  else {
103  G4BoundingSphereScene bsScene(this);
104  const G4int tempRequestedDepth = fRequestedDepth;
105  fRequestedDepth = -1; // Always search to all depths to define extent.
106  const G4ModelingParameters* tempMP = fpMP;
107  G4ModelingParameters mParams
108  (0, // No default vis attributes needed.
109  G4ModelingParameters::wf, // wireframe (not relevant for this).
110  true, // Global culling.
111  true, // Cull invisible volumes.
112  false, // Density culling.
113  0., // Density (not relevant if density culling false).
114  true, // Cull daughters of opaque mothers.
115  24); // No of sides (not relevant for this operation).
116  fpMP = &mParams;
117  DescribeYourselfTo (bsScene);
118  G4double radius = bsScene.GetRadius();
119  if (radius < 0.) { // Nothing in the scene.
120  fExtent = fpTopPV -> GetLogicalVolume () -> GetSolid () -> GetExtent ();
121  } else {
122  // Transform back to coordinates relative to the top
123  // transformation, which is in G4VModel::fTransform. This makes
124  // it conform to all models, which are defined by a
125  // transformation and an extent relative to that
126  // transformation...
127  G4Point3D centre = bsScene.GetCentre();
128  centre.transform(fTransform.inverse());
129  fExtent = G4VisExtent(centre, radius);
130  }
131  fpMP = tempMP;
132  fRequestedDepth = tempRequestedDepth;
133  }
134 }
G4Transform3D fTransform
Definition: G4VModel.hh:112
Transform3D inverse() const
Definition: Transform3D.cc:142
const G4VisExtent & GetExtent() const
int G4int
Definition: G4Types.hh:78
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
G4VPhysicalVolume * fpTopPV
G4VisExtent fExtent
Definition: G4VModel.hh:111
double G4double
Definition: G4Types.hh:76
void DescribeYourselfTo(G4VGraphicsScene &)
std::vector< G4AttValue > * G4PhysicalVolumeModel::CreateCurrentAttValues ( ) const

Definition at line 832 of file G4PhysicalVolumeModel.cc.

References fFullPVPath, fpCurrentLV, fpCurrentMaterial, fpCurrentTransform, G4BestUnit, G4Exception(), G4Material::GetDensity(), G4VSolid::GetEntityType(), G4VSolid::GetName(), G4Region::GetName(), G4Material::GetName(), G4LogicalVolume::GetName(), G4Material::GetRadlen(), G4LogicalVolume::GetRegion(), G4LogicalVolume::GetSolid(), G4Material::GetState(), G4LogicalVolume::IsRootRegion(), JustWarning, and kStateUndefined.

Referenced by G4VSceneHandler::LoadAtts(), and G4XXXStoredSceneHandler::PreAddSolid().

833 {
834  std::vector<G4AttValue>* values = new std::vector<G4AttValue>;
835  std::ostringstream oss;
836  for (size_t i = 0; i < fFullPVPath.size(); ++i) {
837  oss << fFullPVPath[i].GetPhysicalVolume()->GetName()
838  << ':' << fFullPVPath[i].GetCopyNo();
839  if (i != fFullPVPath.size() - 1) oss << '/';
840  }
841 
842  if (!fpCurrentLV) {
844  ("G4PhysicalVolumeModel::CreateCurrentAttValues",
845  "modeling0004",
846  JustWarning,
847  "Current logical volume not defined.");
848  return values;
849  }
850 
851  values->push_back(G4AttValue("PVPath", oss.str(),""));
852  values->push_back(G4AttValue("LVol", fpCurrentLV->GetName(),""));
853  G4VSolid* pSol = fpCurrentLV->GetSolid();
854  values->push_back(G4AttValue("Solid", pSol->GetName(),""));
855  values->push_back(G4AttValue("EType", pSol->GetEntityType(),""));
856  oss.str(""); oss << '\n' << *pSol;
857  values->push_back(G4AttValue("DmpSol", oss.str(),""));
858  oss.str(""); oss << '\n' << *fpCurrentTransform;
859  values->push_back(G4AttValue("Trans", oss.str(),""));
860  G4String matName = fpCurrentMaterial? fpCurrentMaterial->GetName(): G4String("No material");
861  values->push_back(G4AttValue("Material", matName,""));
863  values->push_back(G4AttValue("Density", G4BestUnit(matDensity,"Volumic Mass"),""));
865  oss.str(""); oss << matState;
866  values->push_back(G4AttValue("State", oss.str(),""));
868  values->push_back(G4AttValue("Radlen", G4BestUnit(matRadlen,"Length"),""));
869  G4Region* region = fpCurrentLV->GetRegion();
870  G4String regionName = region? region->GetName(): G4String("No region");
871  values->push_back(G4AttValue("Region", regionName,""));
872  oss.str(""); oss << fpCurrentLV->IsRootRegion();
873  values->push_back(G4AttValue("RootRegion", oss.str(),""));
874  return values;
875 }
G4String GetName() const
const G4String & GetName() const
G4String GetName() const
G4State
Definition: G4Material.hh:114
const G4String & GetName() const
Definition: G4Material.hh:176
G4double GetDensity() const
Definition: G4Material.hh:178
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
G4Region * GetRegion() const
G4Transform3D * fpCurrentTransform
virtual G4GeometryType GetEntityType() const =0
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
G4bool IsRootRegion() const
G4double GetRadlen() const
Definition: G4Material.hh:218
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4State GetState() const
Definition: G4Material.hh:179
double G4double
Definition: G4Types.hh:76
G4VSolid * GetSolid() const
void G4PhysicalVolumeModel::CurtailDescent ( )
inline
void G4PhysicalVolumeModel::DescribeAndDescend ( G4VPhysicalVolume pVPV,
G4int  requestedDepth,
G4LogicalVolume pLV,
G4VSolid pSol,
G4Material pMaterial,
const G4Transform3D theAT,
G4VGraphicsScene sceneHandler 
)
protected

Definition at line 336 of file G4PhysicalVolumeModel.cc.

References density, HepGeom::Transform3D::dx(), HepGeom::Transform3D::dy(), HepGeom::Transform3D::dz(), G4Colour::GetAlpha(), G4VisAttributes::GetColour(), HepGeom::Transform3D::getDecomposition(), G4Material::GetDensity(), G4VisAttributes::GetForcedDrawingStyle(), G4VisAttributes::GetForcedLineSegmentsPerCircle(), G4VisAttributes::GetLineStyle(), G4VisAttributes::GetLineWidth(), G4LogicalVolume::GetNoDaughters(), G4LogicalVolume::GetVisAttributes(), G4ModelingParameters::hlhsr, G4ModelingParameters::hsr, HepGeom::Transform3D::inverse(), G4VisAttributes::IsDaughtersInvisible(), G4VisAttributes::IsForceAuxEdgeVisible(), G4VisAttributes::IsForceDrawingStyle(), G4VisAttributes::IsVisible(), G4VisAttributes::SetColour(), G4VisAttributes::SetDaughtersInvisible(), G4VisAttributes::SetForceAuxEdgeVisible(), G4VisAttributes::SetForceLineSegmentsPerCircle(), G4VisAttributes::SetForceSolid(), G4VisAttributes::SetForceWireframe(), G4VisAttributes::SetLineStyle(), G4VisAttributes::SetLineWidth(), G4VisAttributes::SetVisibility(), G4VisAttributes::solid, G4ModelingParameters::VASColour, G4ModelingParameters::VASDaughtersInvisible, G4ModelingParameters::VASForceAuxEdgeVisible, G4ModelingParameters::VASForceLineSegmentsPerCircle, G4ModelingParameters::VASForceSolid, G4ModelingParameters::VASForceWireframe, G4ModelingParameters::VASLineStyle, G4ModelingParameters::VASLineWidth, G4ModelingParameters::VASVisibility, and G4VisAttributes::wireframe.

343 {
344  // Maintain useful data members...
345  fpCurrentPV = pVPV;
346  fpCurrentLV = pLV;
347  fpCurrentMaterial = pMaterial;
348 
349  const G4RotationMatrix objectRotation = pVPV -> GetObjectRotationValue ();
350  const G4ThreeVector& translation = pVPV -> GetTranslation ();
351  G4Transform3D theLT (G4Transform3D (objectRotation, translation));
352 
353  // Compute the accumulated transformation...
354  // Note that top volume's transformation relative to the world
355  // coordinate system is specified in theAT == startingTransformation
356  // = fTransform (see DescribeYourselfTo), so first time through the
357  // volume's own transformation, which is only relative to its
358  // mother, i.e., not relative to the world coordinate system, should
359  // not be accumulated.
360  G4Transform3D theNewAT (theAT);
361  if (fCurrentDepth != 0) theNewAT = theAT * theLT;
362  fpCurrentTransform = &theNewAT;
363 
364  const G4VisAttributes* pVisAttribs = pLV->GetVisAttributes();
365  if (!pVisAttribs) pVisAttribs = fpMP->GetDefaultVisAttributes();
366  // Beware - pVisAttribs might still be zero - create a temporary default one...
367  G4bool visAttsCreated = false;
368  if (!pVisAttribs) {
369  pVisAttribs = new G4VisAttributes;
370  visAttsCreated = true;
371  }
372  // From here, can assume pVisAttribs is a valid pointer.
373 
374  // Make decision to draw...
375  G4bool thisToBeDrawn = true;
376 
377  // Update full path of physical volumes...
378  G4int copyNo = fpCurrentPV->GetCopyNo();
379  fFullPVPath.push_back
380  (G4PhysicalVolumeNodeID
382 
383  // In case we need to copy the vis atts for modification...
384  G4bool copyForVAM = false;
385  const G4VisAttributes* pUnmodifiedVisAtts = pVisAttribs;
386  G4VisAttributes* pModifiedVisAtts = 0;
387 
388  // Check if vis attributes are to be modified by a /vis/touchable/set/ command.
389  const std::vector<G4ModelingParameters::VisAttributesModifier>& vams =
391  std::vector<G4ModelingParameters::VisAttributesModifier>::const_iterator
392  iModifier;
393  for (iModifier = vams.begin();
394  iModifier != vams.end();
395  ++iModifier) {
397  iModifier->GetPVNameCopyNoPath();
398  if (vamPath.size() == fFullPVPath.size()) {
399  // OK - there's a size match. Check it out.
400 // G4cout << "Size match" << G4endl;
402  std::vector<G4PhysicalVolumeNodeID>::const_iterator iPVNodeId;
403  for (iVAMNameCopyNo = vamPath.begin(), iPVNodeId = fFullPVPath.begin();
404  iVAMNameCopyNo != vamPath.end();
405  ++iVAMNameCopyNo, ++iPVNodeId) {
406 // G4cout
407 // << iVAMNameCopyNo->fName
408 // << ',' << iVAMNameCopyNo->fCopyNo
409 // << "; " << iPVNodeId->GetPhysicalVolume()->GetName()
410 // << ',' << iPVNodeId->GetPhysicalVolume()->GetCopyNo()
411 // << G4endl;
412  if (!(
413  iVAMNameCopyNo->GetName() ==
414  iPVNodeId->GetPhysicalVolume()->GetName() &&
415  iVAMNameCopyNo->GetCopyNo() ==
416  iPVNodeId->GetPhysicalVolume()->GetCopyNo()
417  )) {
418  break;
419  }
420  }
421  if (iVAMNameCopyNo == vamPath.end()) {
422 // G4cout << "Match found" << G4endl;
423  if (!copyForVAM) {
424  pModifiedVisAtts = new G4VisAttributes(*pUnmodifiedVisAtts);
425  pVisAttribs = pModifiedVisAtts;
426  copyForVAM = true;
427  }
428  const G4VisAttributes& transVisAtts = iModifier->GetVisAttributes();
429  switch (iModifier->GetVisAttributesSignifier()) {
431  pModifiedVisAtts->SetVisibility(transVisAtts.IsVisible());
432  break;
434  pModifiedVisAtts->SetDaughtersInvisible
435  (transVisAtts.IsDaughtersInvisible());
436  break;
438  pModifiedVisAtts->SetColour(transVisAtts.GetColour());
439  break;
441  pModifiedVisAtts->SetLineStyle(transVisAtts.GetLineStyle());
442  break;
444  pModifiedVisAtts->SetLineWidth(transVisAtts.GetLineWidth());
445  break;
447  if (transVisAtts.GetForcedDrawingStyle() ==
449  pModifiedVisAtts->SetForceWireframe
450  (transVisAtts.IsForceDrawingStyle());
451  }
452  break;
454  if (transVisAtts.GetForcedDrawingStyle() ==
456  pModifiedVisAtts->SetForceSolid
457  (transVisAtts.IsForceDrawingStyle());
458  }
459  break;
461  pModifiedVisAtts->SetForceAuxEdgeVisible
462  (transVisAtts.IsForceAuxEdgeVisible());
463  break;
465  pModifiedVisAtts->SetForceLineSegmentsPerCircle
466  (transVisAtts.GetForcedLineSegmentsPerCircle());
467  break;
468  }
469  }
470  }
471  }
472 
473  // There are various reasons why this volume
474  // might not be drawn...
475  G4bool culling = fpMP->IsCulling();
476  G4bool cullingInvisible = fpMP->IsCullingInvisible();
477  G4bool markedVisible = pVisAttribs->IsVisible();
478  G4bool cullingLowDensity = fpMP->IsDensityCulling();
479  G4double density = pMaterial? pMaterial->GetDensity(): 0;
480  G4double densityCut = fpMP -> GetVisibleDensity ();
481 
482  // 1) Global culling is on....
483  if (culling) {
484  // 2) Culling of invisible volumes is on...
485  if (cullingInvisible) {
486  // 3) ...and the volume is marked not visible...
487  if (!markedVisible) thisToBeDrawn = false;
488  }
489  // 4) Or culling of low density volumes is on...
490  if (cullingLowDensity) {
491  // 5) ...and density is less than cut value...
492  if (density < densityCut) thisToBeDrawn = false;
493  }
494  }
495 
496  // Record thisToBeDrawn in path...
497  fFullPVPath.back().SetDrawn(thisToBeDrawn);
498 
499  if (thisToBeDrawn) {
500 
501  // Update path of drawn physical volumes...
502  fDrawnPVPath.push_back
503  (G4PhysicalVolumeNodeID
504  (fpCurrentPV,copyNo,fCurrentDepth,*fpCurrentTransform,thisToBeDrawn));
505 
506  if (fpMP->IsExplode() && fDrawnPVPath.size() == 1) {
507  // For top-level drawn volumes, explode along radius...
509  G4Transform3D centred = centering.inverse() * theNewAT;
510  G4Scale3D oldScale;
511  G4Rotate3D oldRotation;
512  G4Translate3D oldTranslation;
513  centred.getDecomposition(oldScale, oldRotation, oldTranslation);
514  G4double explodeFactor = fpMP->GetExplodeFactor();
515  G4Translate3D newTranslation =
516  G4Translate3D(explodeFactor * oldTranslation.dx(),
517  explodeFactor * oldTranslation.dy(),
518  explodeFactor * oldTranslation.dz());
519  theNewAT = centering * newTranslation * oldRotation * oldScale;
520  }
521 
522  DescribeSolid (theNewAT, pSol, pVisAttribs, sceneHandler);
523 
524  }
525 
526  // Make decision to draw daughters, if any. There are various
527  // reasons why daughters might not be drawn...
528 
529  // First, reasons that do not depend on culling policy...
530  G4int nDaughters = pLV->GetNoDaughters();
531  G4bool daughtersToBeDrawn = true;
532  // 1) There are no daughters...
533  if (!nDaughters) daughtersToBeDrawn = false;
534  // 2) We are at the limit if requested depth...
535  else if (requestedDepth == 0) daughtersToBeDrawn = false;
536  // 3) The user has asked that the descent be curtailed...
537  else if (fCurtailDescent) daughtersToBeDrawn = false;
538 
539  // Now, reasons that depend on culling policy...
540  else {
541  G4bool daughtersInvisible = pVisAttribs->IsDaughtersInvisible();
542  // Culling of covered daughters request. This is computed in
543  // G4VSceneHandler::CreateModelingParameters() depending on view
544  // parameters...
545  G4bool cullingCovered = fpMP->IsCullingCovered();
546  G4bool surfaceDrawing =
549  if (pVisAttribs->IsForceDrawingStyle()) {
550  switch (pVisAttribs->GetForcedDrawingStyle()) {
551  default:
552  case G4VisAttributes::wireframe: surfaceDrawing = false; break;
553  case G4VisAttributes::solid: surfaceDrawing = true; break;
554  }
555  }
556  G4bool opaque = pVisAttribs->GetColour().GetAlpha() >= 1.;
557  // 4) Global culling is on....
558  if (culling) {
559  // 5) ..and culling of invisible volumes is on...
560  if (cullingInvisible) {
561  // 6) ...and the mother requests daughters invisible
562  if (daughtersInvisible) daughtersToBeDrawn = false;
563  }
564  // 7) Or culling of covered daughters is requested...
565  if (cullingCovered) {
566  // 8) ...and surface drawing is operating...
567  if (surfaceDrawing) {
568  // 9) ...but only if mother is visible...
569  if (thisToBeDrawn) {
570  // 10) ...and opaque...
571  if (opaque) daughtersToBeDrawn = false;
572  }
573  }
574  }
575  }
576  }
577 
578  // Delete modified vis atts if created...
579  if (copyForVAM) {
580  delete pModifiedVisAtts;
581  pVisAttribs = pUnmodifiedVisAtts;
582  copyForVAM = false;
583  }
584 
585  // Vis atts for this volume no longer needed if created...
586  if (visAttsCreated) delete pVisAttribs;
587 
588  if (daughtersToBeDrawn) {
589  for (G4int iDaughter = 0; iDaughter < nDaughters; iDaughter++) {
590  // Store daughter pVPV in local variable ready for recursion...
591  G4VPhysicalVolume* pDaughterVPV = pLV -> GetDaughter (iDaughter);
592  // Descend the geometry structure recursively...
593  fCurrentDepth++;
595  (pDaughterVPV, requestedDepth - 1, theNewAT, sceneHandler);
596  fCurrentDepth--;
597  }
598  }
599 
600  // Reset for normal descending of next volume at this level...
601  fCurtailDescent = false;
602 
603  // Pop item from paths physical volumes...
604  fFullPVPath.pop_back();
605  if (thisToBeDrawn) {
606  fDrawnPVPath.pop_back();
607  }
608 }
G4bool IsForceAuxEdgeVisible() const
void SetColour(const G4Colour &)
G4double GetAlpha() const
Definition: G4Colour.hh:142
void SetForceWireframe(G4bool)
const G4Point3D & GetExplodeCentre() const
const G4VisAttributes * GetDefaultVisAttributes() const
G4double GetLineWidth() const
void getDecomposition(Scale3D &scale, Rotate3D &rotation, Translate3D &translation) const
Definition: Transform3D.cc:174
void SetLineStyle(LineStyle)
void SetLineWidth(G4double)
void SetVisibility(G4bool)
G4double GetDensity() const
Definition: G4Material.hh:178
G4bool IsVisible() const
const G4Colour & GetColour() const
Transform3D inverse() const
Definition: Transform3D.cc:142
double dx() const
Definition: Transform3D.h:279
G4bool IsExplode() const
G4bool IsDensityCulling() const
G4Transform3D * fpCurrentTransform
void SetForceSolid(G4bool)
int G4int
Definition: G4Types.hh:78
G4VPhysicalVolume * fpCurrentPV
LineStyle GetLineStyle() const
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
G4double density
Definition: TRTMaterials.hh:39
G4bool IsDaughtersInvisible() const
bool G4bool
Definition: G4Types.hh:79
G4bool IsCullingCovered() const
const G4VisAttributes * GetVisAttributes() const
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
virtual void DescribeSolid(const G4Transform3D &theAT, G4VSolid *pSol, const G4VisAttributes *pVisAttribs, G4VGraphicsScene &sceneHandler)
G4bool IsCullingInvisible() const
double dy() const
Definition: Transform3D.h:282
G4int GetNoDaughters() const
void SetForceAuxEdgeVisible(G4bool)
G4int GetForcedLineSegmentsPerCircle() const
std::vector< PVNameCopyNo > PVNameCopyNoPath
double dz() const
Definition: Transform3D.h:285
G4double GetExplodeFactor() const
G4bool IsCulling() const
G4bool IsForceDrawingStyle() const
virtual G4int GetCopyNo() const =0
HepGeom::Translate3D G4Translate3D
std::vector< G4PhysicalVolumeNodeID > fDrawnPVPath
double G4double
Definition: G4Types.hh:76
void VisitGeometryAndGetVisReps(G4VPhysicalVolume *, G4int requestedDepth, const G4Transform3D &, G4VGraphicsScene &)
void SetDaughtersInvisible(G4bool)
ForcedDrawingStyle GetForcedDrawingStyle() const
PVNameCopyNoPath::const_iterator PVNameCopyNoPathConstIterator
DrawingStyle GetDrawingStyle() const
void SetForceLineSegmentsPerCircle(G4int nSegments)
const std::vector< VisAttributesModifier > & GetVisAttributesModifiers() const
void G4PhysicalVolumeModel::DescribeSolid ( const G4Transform3D theAT,
G4VSolid pSol,
const G4VisAttributes pVisAttribs,
G4VGraphicsScene sceneHandler 
)
protectedvirtual

Reimplemented in G4LogicalVolumeModel.

Definition at line 611 of file G4PhysicalVolumeModel.cc.

References G4VGraphicsScene::AddPrimitive(), G4VGraphicsScene::BeginPrimitives(), G4VGraphicsScene::EndPrimitives(), G4cout, G4endl, G4VisAttributes::GetForcedLineSegmentsPerCircle(), G4VSolid::GetName(), G4VSolid::GetPolyhedron(), HepGeom::Transform3D::inverse(), G4VisAttributes::IsForceLineSegmentsPerCircle(), G4VGraphicsScene::PostAddSolid(), G4VGraphicsScene::PreAddSolid(), G4Colour::Red(), HepPolyhedron::ResetNumberOfRotationSteps(), G4VisAttributes::SetColour(), HepPolyhedron::SetNumberOfRotationSteps(), and G4Visible::SetVisAttributes().

615 {
616  sceneHandler.PreAddSolid (theAT, *pVisAttribs);
617 
618  G4VSolid* pSectionSolid = fpMP->GetSectionSolid();
619  G4VSolid* pCutawaySolid = fpMP->GetCutawaySolid();
620 
621  if (!fpClippingSolid && !pSectionSolid && !pCutawaySolid) {
622 
623  pSol -> DescribeYourselfTo (sceneHandler); // Standard treatment.
624 
625  } else {
626 
627  // Clipping, etc., performed by Boolean operations.
628 
629  // First, get polyhedron for current solid...
630  if (pVisAttribs->IsForceLineSegmentsPerCircle())
632  (pVisAttribs->GetForcedLineSegmentsPerCircle());
633  else
635  const G4Polyhedron* pOriginal = pSol->GetPolyhedron();
637 
638  if (!pOriginal) {
639 
640  if (fpMP->IsWarning())
641  G4cout <<
642  "WARNING: G4PhysicalVolumeModel::DescribeSolid: solid\n \""
643  << pSol->GetName() <<
644  "\" has no polyhedron. Cannot by clipped."
645  << G4endl;
646  pSol -> DescribeYourselfTo (sceneHandler); // Standard treatment.
647 
648  } else {
649 
650  G4Polyhedron resultant(*pOriginal);
651  G4VisAttributes resultantVisAttribs(*pVisAttribs);
652  G4VSolid* resultantSolid = 0;
653 
654  if (fpClippingSolid) {
655  switch (fClippingMode) {
656  default:
657  case subtraction:
658  resultantSolid = new G4SubtractionSolid
659  ("resultant_solid", pSol, fpClippingSolid, theAT.inverse());
660  break;
661  case intersection:
662  resultantSolid = new G4IntersectionSolid
663  ("resultant_solid", pSol, fpClippingSolid, theAT.inverse());
664  break;
665  }
666  }
667 
668  if (pSectionSolid) {
669  resultantSolid = new G4IntersectionSolid
670  ("sectioned_solid", pSol, pSectionSolid, theAT.inverse());
671  }
672 
673  if (pCutawaySolid) {
674  resultantSolid = new G4SubtractionSolid
675  ("cutaway_solid", pSol, pCutawaySolid, theAT.inverse());
676  }
677 
678  G4Polyhedron* tmpResultant = resultantSolid->GetPolyhedron();
679  if (tmpResultant) resultant = *tmpResultant;
680  else {
681  if (fpMP->IsWarning())
682  G4cout <<
683  "WARNING: G4PhysicalVolumeModel::DescribeSolid: resultant polyhedron for"
684  "\n solid \"" << pSol->GetName() <<
685  "\" not defined due to error during Boolean processing."
686  "\n Original will be drawn in red."
687  << G4endl;
688  resultantVisAttribs.SetColour(G4Colour::Red());
689  }
690 
691  delete resultantSolid;
692 
693  // Finally, force polyhedron drawing...
694  resultant.SetVisAttributes(resultantVisAttribs);
695  sceneHandler.BeginPrimitives(theAT);
696  sceneHandler.AddPrimitive(resultant);
697  sceneHandler.EndPrimitives();
698  }
699  }
700  sceneHandler.PostAddSolid ();
701 }
virtual void PostAddSolid()=0
virtual G4Polyhedron * GetPolyhedron() const
Definition: G4VSolid.cc:644
G4String GetName() const
G4VSolid * GetSectionSolid() const
virtual void BeginPrimitives(const G4Transform3D &objectTransformation=G4Transform3D())=0
Transform3D inverse() const
Definition: Transform3D.cc:142
G4VSolid * GetCutawaySolid() const
virtual void AddPrimitive(const G4Polyline &)=0
G4GLOB_DLL std::ostream G4cout
G4int GetNoOfSides() const
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
G4int GetForcedLineSegmentsPerCircle() const
virtual void PreAddSolid(const G4Transform3D &objectTransformation, const G4VisAttributes &visAttribs)=0
#define G4endl
Definition: G4ios.hh:61
static G4Colour Red()
Definition: G4Colour.hh:148
virtual void EndPrimitives()=0
static void SetNumberOfRotationSteps(G4int n)
G4bool IsForceLineSegmentsPerCircle() const
G4bool IsWarning() const
static void ResetNumberOfRotationSteps()
void DescribeYourselfTo(G4VGraphicsScene &)
void G4PhysicalVolumeModel::DescribeYourselfTo ( G4VGraphicsScene sceneHandler)
virtual

Implements G4VModel.

Definition at line 137 of file G4PhysicalVolumeModel.cc.

References FatalException, and G4Exception().

Referenced by CalculateExtent(), G4LogicalVolumeModel::DescribeYourselfTo(), G4ASCIITreeSceneHandler::EndModeling(), G4VisCommandSceneAddVolume::SetNewValue(), and Validate().

138 {
139  if (!fpTopPV) G4Exception
140  ("G4PhysicalVolumeModel::DescribeYourselfTo",
141  "modeling0012", FatalException, "No model.");
142 
143  if (!fpMP) G4Exception
144  ("G4PhysicalVolumeModel::DescribeYourselfTo",
145  "modeling0003", FatalException, "No modeling parameters.");
146 
147  // For safety...
148  fCurrentDepth = 0;
149 
150  G4Transform3D startingTransformation = fTransform;
151 
153  (fpTopPV,
155  startingTransformation,
156  sceneHandler);
157 
158  // Clear data...
159  fCurrentDepth = 0;
160  fpCurrentPV = 0;
161  fpCurrentLV = 0;
162  fpCurrentMaterial = 0;
163  if (fFullPVPath.size() != fBaseFullPVPath.size()) {
164  // They should be equal if pushing and popping is happening properly.
166  ("G4PhysicalVolumeModel::DescribeYourselfTo",
167  "modeling0013",
169  "Path at start of modeling not equal to base path. Something badly"
170  "\nwrong. Please contact visualisation coordinator.");
171  }
172  fDrawnPVPath.clear();
173 }
G4Transform3D fTransform
Definition: G4VModel.hh:112
G4VPhysicalVolume * fpCurrentPV
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
std::vector< G4PhysicalVolumeNodeID > fBaseFullPVPath
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4VPhysicalVolume * fpTopPV
std::vector< G4PhysicalVolumeNodeID > fDrawnPVPath
void VisitGeometryAndGetVisReps(G4VPhysicalVolume *, G4int requestedDepth, const G4Transform3D &, G4VGraphicsScene &)
const std::map< G4String, G4AttDef > * G4PhysicalVolumeModel::GetAttDefs ( ) const

Definition at line 757 of file G4PhysicalVolumeModel.cc.

References G4AttDefStore::GetInstance().

Referenced by G4VSceneHandler::LoadAtts(), G4XXXStoredSceneHandler::PreAddSolid(), and G4VisCommandList::SetNewValue().

758 {
759  G4bool isNew;
760  std::map<G4String,G4AttDef>* store
761  = G4AttDefStore::GetInstance("G4PhysicalVolumeModel", isNew);
762  if (isNew) {
763  (*store)["PVPath"] =
764  G4AttDef("PVPath","Physical Volume Path","Physics","","G4String");
765  (*store)["LVol"] =
766  G4AttDef("LVol","Logical Volume","Physics","","G4String");
767  (*store)["Solid"] =
768  G4AttDef("Solid","Solid Name","Physics","","G4String");
769  (*store)["EType"] =
770  G4AttDef("EType","Entity Type","Physics","","G4String");
771  (*store)["DmpSol"] =
772  G4AttDef("DmpSol","Dump of Solid properties","Physics","","G4String");
773  (*store)["Trans"] =
774  G4AttDef("Trans","Transformation of volume","Physics","","G4String");
775  (*store)["Material"] =
776  G4AttDef("Material","Material Name","Physics","","G4String");
777  (*store)["Density"] =
778  G4AttDef("Density","Material Density","Physics","G4BestUnit","G4double");
779  (*store)["State"] =
780  G4AttDef("State","Material State (enum undefined,solid,liquid,gas)","Physics","","G4String");
781  (*store)["Radlen"] =
782  G4AttDef("Radlen","Material Radiation Length","Physics","G4BestUnit","G4double");
783  (*store)["Region"] =
784  G4AttDef("Region","Cuts Region","Physics","","G4String");
785  (*store)["RootRegion"] =
786  G4AttDef("RootRegion","Root Region (0/1 = false/true)","Physics","","G4bool");
787  }
788  return store;
789 }
bool G4bool
Definition: G4Types.hh:79
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)
const G4VSolid* G4PhysicalVolumeModel::GetClippingSolid ( ) const
inline

Definition at line 158 of file G4PhysicalVolumeModel.hh.

References fpClippingSolid.

159  {return fpClippingSolid;}
G4int G4PhysicalVolumeModel::GetCurrentDepth ( ) const
inline
G4String G4PhysicalVolumeModel::GetCurrentDescription ( ) const
virtual

Reimplemented from G4VModel.

Definition at line 187 of file G4PhysicalVolumeModel.cc.

References GetCurrentTag().

188 {
189  return "G4PhysicalVolumeModel " + GetCurrentTag ();
190 }
G4String GetCurrentTag() const
G4LogicalVolume* G4PhysicalVolumeModel::GetCurrentLV ( ) const
inline
G4Material* G4PhysicalVolumeModel::GetCurrentMaterial ( ) const
inline
G4VPhysicalVolume* G4PhysicalVolumeModel::GetCurrentPV ( ) const
inline
G4String G4PhysicalVolumeModel::GetCurrentTag ( ) const
virtual

Reimplemented from G4VModel.

Definition at line 175 of file G4PhysicalVolumeModel.cc.

References G4VModel::fGlobalTag, and fpCurrentPV.

Referenced by GetCurrentDescription().

176 {
177  if (fpCurrentPV) {
178  std::ostringstream o;
179  o << fpCurrentPV -> GetCopyNo ();
180  return fpCurrentPV -> GetName () + "." + o.str();
181  }
182  else {
183  return "WARNING: NO CURRENT VOLUME - global tag is " + fGlobalTag;
184  }
185 }
G4VPhysicalVolume * fpCurrentPV
G4String fGlobalTag
Definition: G4VModel.hh:109
const std::vector<G4PhysicalVolumeNodeID>& G4PhysicalVolumeModel::GetDrawnPVPath ( ) const
inline
const std::vector<G4PhysicalVolumeNodeID>& G4PhysicalVolumeModel::GetFullPVPath ( ) const
inline

Definition at line 173 of file G4PhysicalVolumeModel.hh.

References fFullPVPath.

174  {return fFullPVPath;}
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
G4int G4PhysicalVolumeModel::GetRequestedDepth ( ) const
inline

Definition at line 156 of file G4PhysicalVolumeModel.hh.

References fRequestedDepth.

Referenced by G4ASCIITreeSceneHandler::EndModeling().

G4VPhysicalVolume* G4PhysicalVolumeModel::GetTopPhysicalVolume ( ) const
inline

Definition at line 154 of file G4PhysicalVolumeModel.hh.

References fpTopPV.

Referenced by G4GMocrenFileSceneHandler::AddSolid(), and G4ASCIITreeSceneHandler::EndModeling().

154 {return fpTopPV;}
G4VPhysicalVolume * fpTopPV
void G4PhysicalVolumeModel::SetBaseFullPVPath ( const std::vector< G4PhysicalVolumeNodeID > &  baseFullPVPath)
inline

Definition at line 201 of file G4PhysicalVolumeModel.hh.

References fBaseFullPVPath, and fFullPVPath.

Referenced by G4VisCommandSceneAddVolume::SetNewValue().

202  {
203  fBaseFullPVPath = baseFullPVPath;
204  fFullPVPath = baseFullPVPath;
205  }
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
std::vector< G4PhysicalVolumeNodeID > fBaseFullPVPath
void G4PhysicalVolumeModel::SetClippingMode ( ClippingMode  mode)
inline

Definition at line 215 of file G4PhysicalVolumeModel.hh.

References fClippingMode.

215  {
216  fClippingMode = mode;
217  }
void G4PhysicalVolumeModel::SetClippingSolid ( G4VSolid pClippingSolid)
inline

Definition at line 211 of file G4PhysicalVolumeModel.hh.

References fpClippingSolid.

211  {
212  fpClippingSolid = pClippingSolid;
213  }
void G4PhysicalVolumeModel::SetRequestedDepth ( G4int  requestedDepth)
inline

Definition at line 207 of file G4PhysicalVolumeModel.hh.

References fRequestedDepth.

207  {
208  fRequestedDepth = requestedDepth;
209  }
G4bool G4PhysicalVolumeModel::Validate ( G4bool  warn)
virtual

Reimplemented from G4VModel.

Definition at line 703 of file G4PhysicalVolumeModel.cc.

References CalculateExtent(), DescribeYourselfTo(), G4VModel::fpMP, fpTopPV, fTopPVCopyNo, fTopPVName, G4cout, G4endl, G4ModelingParameters::GetDefaultVisAttributes(), G4PhysicalVolumeSearchScene::GetFoundVolume(), G4TransportationManager::GetNoWorlds(), G4TransportationManager::GetTransportationManager(), G4TransportationManager::GetWorldsIterator(), G4ModelingParameters::SetDefaultVisAttributes(), and G4VModel::SetModelingParameters().

704 {
705  G4TransportationManager* transportationManager =
707 
708  size_t nWorlds = transportationManager->GetNoWorlds();
709 
710  G4bool found = false;
711 
712  std::vector<G4VPhysicalVolume*>::iterator iterWorld =
713  transportationManager->GetWorldsIterator();
714  for (size_t i = 0; i < nWorlds; ++i, ++iterWorld) {
715  G4VPhysicalVolume* world = (*iterWorld);
716  // The idea now is to seek a PV with the same name and copy no
717  // in the hope it's the same one!!
718  G4PhysicalVolumeModel searchModel (world);
719  G4int verbosity = 0; // Suppress messages from G4PhysicalVolumeSearchScene.
720  G4PhysicalVolumeSearchScene searchScene
721  (&searchModel, fTopPVName, fTopPVCopyNo, verbosity);
722  G4ModelingParameters mp; // Default modeling parameters for this search.
724  searchModel.SetModelingParameters (&mp);
725  searchModel.DescribeYourselfTo (searchScene);
726  G4VPhysicalVolume* foundVolume = searchScene.GetFoundVolume ();
727  if (foundVolume) {
728  if (foundVolume != fpTopPV && warn) {
729  G4cout <<
730  "G4PhysicalVolumeModel::Validate(): A volume of the same name and"
731  "\n copy number (\""
732  << fTopPVName << "\", copy " << fTopPVCopyNo
733  << ") still exists and is being used."
734  "\n But it is not the same volume you originally specified"
735  "\n in /vis/scene/add/."
736  << G4endl;
737  }
738  fpTopPV = foundVolume;
739  CalculateExtent ();
740  found = true;
741  }
742  }
743  if (found) return true;
744  else {
745  if (warn) {
746  G4cout <<
747  "G4PhysicalVolumeModel::Validate(): No volume of name and"
748  "\n copy number (\""
749  << fTopPVName << "\", copy " << fTopPVCopyNo
750  << ") exists."
751  << G4endl;
752  }
753  return false;
754  }
755 }
const G4VisAttributes * GetDefaultVisAttributes() const
std::vector< G4VPhysicalVolume * >::iterator GetWorldsIterator()
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
static G4TransportationManager * GetTransportationManager()
G4VPhysicalVolume * fpTopPV
size_t GetNoWorlds() const
#define G4endl
Definition: G4ios.hh:61
void SetDefaultVisAttributes(const G4VisAttributes *pDefaultVisAttributes)
void G4PhysicalVolumeModel::VisitGeometryAndGetVisReps ( G4VPhysicalVolume pVPV,
G4int  requestedDepth,
const G4Transform3D theAT,
G4VGraphicsScene sceneHandler 
)
protected

Definition at line 193 of file G4PhysicalVolumeModel.cc.

References G4cout, G4endl, G4VSolid::GetEntityType(), G4VSolid::GetName(), pyG4VTouchable::GetRotation, kPhi, kRho, kXAxis, kYAxis, kZAxis, n, CLHEP::HepRotation::rotateZ(), and width.

197 {
198  // Visits geometry structure to a given depth (requestedDepth), starting
199  // at given physical volume with given starting transformation and
200  // describes volumes to the scene handler.
201  // requestedDepth < 0 (default) implies full visit.
202  // theAT is the Accumulated Transformation.
203 
204  // Find corresponding logical volume and (later) solid, storing in
205  // local variables to preserve re-entrancy.
206  G4LogicalVolume* pLV = pVPV -> GetLogicalVolume ();
207 
208  G4VSolid* pSol;
209  G4Material* pMaterial;
210 
211  if (!(pVPV -> IsReplicated ())) {
212  // Non-replicated physical volume.
213  pSol = pLV -> GetSolid ();
214  pMaterial = pLV -> GetMaterial ();
215  DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
216  theAT, sceneHandler);
217  }
218  else {
219  // Replicated or parametrised physical volume.
220  EAxis axis;
221  G4int nReplicas;
222  G4double width;
223  G4double offset;
224  G4bool consuming;
225  pVPV -> GetReplicationData (axis, nReplicas, width, offset, consuming);
226  if (fCurrentDepth == 0) nReplicas = 1; // Just draw first
227  G4VPVParameterisation* pP = pVPV -> GetParameterisation ();
228  if (pP) { // Parametrised volume.
229  for (int n = 0; n < nReplicas; n++) {
230  pSol = pP -> ComputeSolid (n, pVPV);
231  pP -> ComputeTransformation (n, pVPV);
232  pSol -> ComputeDimensions (pP, n, pVPV);
233  pVPV -> SetCopyNo (n);
234  // Create a touchable of current parent for ComputeMaterial.
235  // fFullPVPath has not been updated yet so at this point it
236  // corresponds to the parent.
237  G4PhysicalVolumeModelTouchable parentTouchable(fFullPVPath);
238  pMaterial = pP -> ComputeMaterial (n, pVPV, &parentTouchable);
239  DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
240  theAT, sceneHandler);
241  }
242  }
243  else { // Plain replicated volume. From geometry_guide.txt...
244  // The replica's positions are claculated by means of a linear formula.
245  // Replication may occur along:
246  //
247  // o Cartesian axes (kXAxis,kYAxis,kZAxis)
248  //
249  // The replications, of specified width have coordinates of
250  // form (-width*(nReplicas-1)*0.5+n*width,0,0) where n=0.. nReplicas-1
251  // for the case of kXAxis, and are unrotated.
252  //
253  // o Radial axis (cylindrical polar) (kRho)
254  //
255  // The replications are cons/tubs sections, centred on the origin
256  // and are unrotated.
257  // They have radii of width*n+offset to width*(n+1)+offset
258  // where n=0..nReplicas-1
259  //
260  // o Phi axis (cylindrical polar) (kPhi)
261  // The replications are `phi sections' or wedges, and of cons/tubs form
262  // They have phi of offset+n*width to offset+(n+1)*width where
263  // n=0..nReplicas-1
264  //
265  pSol = pLV -> GetSolid ();
266  pMaterial = pLV -> GetMaterial ();
267  G4ThreeVector originalTranslation = pVPV -> GetTranslation ();
268  G4RotationMatrix* pOriginalRotation = pVPV -> GetRotation ();
269  G4double originalRMin = 0., originalRMax = 0.;
270  if (axis == kRho && pSol->GetEntityType() == "G4Tubs") {
271  originalRMin = ((G4Tubs*)pSol)->GetInnerRadius();
272  originalRMax = ((G4Tubs*)pSol)->GetOuterRadius();
273  }
274  G4bool visualisable = true;
275  for (int n = 0; n < nReplicas; n++) {
276  G4ThreeVector translation; // Null.
277  G4RotationMatrix rotation; // Null - life long enough for visualizing.
278  G4RotationMatrix* pRotation = 0;
279  switch (axis) {
280  default:
281  case kXAxis:
282  translation = G4ThreeVector (-width*(nReplicas-1)*0.5+n*width,0,0);
283  break;
284  case kYAxis:
285  translation = G4ThreeVector (0,-width*(nReplicas-1)*0.5+n*width,0);
286  break;
287  case kZAxis:
288  translation = G4ThreeVector (0,0,-width*(nReplicas-1)*0.5+n*width);
289  break;
290  case kRho:
291  if (pSol->GetEntityType() == "G4Tubs") {
292  ((G4Tubs*)pSol)->SetInnerRadius(width*n+offset);
293  ((G4Tubs*)pSol)->SetOuterRadius(width*(n+1)+offset);
294  } else {
295  if (fpMP->IsWarning())
296  G4cout <<
297  "G4PhysicalVolumeModel::VisitGeometryAndGetVisReps: WARNING:"
298  "\n built-in replicated volumes replicated in radius for "
299  << pSol->GetEntityType() <<
300  "-type\n solids (your solid \""
301  << pSol->GetName() <<
302  "\") are not visualisable."
303  << G4endl;
304  visualisable = false;
305  }
306  break;
307  case kPhi:
308  rotation.rotateZ (-(offset+(n+0.5)*width));
309  // Minus Sign because for the physical volume we need the
310  // coordinate system rotation.
311  pRotation = &rotation;
312  break;
313  }
314  pVPV -> SetTranslation (translation);
315  pVPV -> SetRotation (pRotation);
316  pVPV -> SetCopyNo (n);
317  if (visualisable) {
318  DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
319  theAT, sceneHandler);
320  }
321  }
322  // Restore originals...
323  pVPV -> SetTranslation (originalTranslation);
324  pVPV -> SetRotation (pOriginalRotation);
325  if (axis == kRho && pSol->GetEntityType() == "G4Tubs") {
326  ((G4Tubs*)pSol)->SetInnerRadius(originalRMin);
327  ((G4Tubs*)pSol)->SetOuterRadius(originalRMax);
328  }
329  }
330  }
331 
332  return;
333 }
Definition: geomdefs.hh:54
G4String GetName() const
CLHEP::Hep3Vector G4ThreeVector
Definition: G4Tubs.hh:84
#define width
virtual G4GeometryType GetEntityType() const =0
int G4int
Definition: G4Types.hh:78
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
const G4int n
EAxis
Definition: geomdefs.hh:54
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:92
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
Definition: geomdefs.hh:54
void DescribeAndDescend(G4VPhysicalVolume *, G4int requestedDepth, G4LogicalVolume *, G4VSolid *, G4Material *, const G4Transform3D &, G4VGraphicsScene &)
G4bool IsWarning() const

Field Documentation

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

Definition at line 260 of file G4PhysicalVolumeModel.hh.

Referenced by SetBaseFullPVPath().

ClippingMode G4PhysicalVolumeModel::fClippingMode
protected

Definition at line 265 of file G4PhysicalVolumeModel.hh.

Referenced by SetClippingMode().

G4int G4PhysicalVolumeModel::fCurrentDepth
protected

Definition at line 255 of file G4PhysicalVolumeModel.hh.

Referenced by GetCurrentDepth().

G4bool G4PhysicalVolumeModel::fCurtailDescent
protected

Definition at line 263 of file G4PhysicalVolumeModel.hh.

Referenced by CurtailDescent().

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

Definition at line 262 of file G4PhysicalVolumeModel.hh.

Referenced by GetDrawnPVPath().

std::vector<G4PhysicalVolumeNodeID> G4PhysicalVolumeModel::fFullPVPath
protected
G4VSolid* G4PhysicalVolumeModel::fpClippingSolid
protected
G4LogicalVolume* G4PhysicalVolumeModel::fpCurrentLV
protected

Definition at line 257 of file G4PhysicalVolumeModel.hh.

Referenced by CreateCurrentAttValues(), and GetCurrentLV().

G4Material* G4PhysicalVolumeModel::fpCurrentMaterial
protected

Definition at line 258 of file G4PhysicalVolumeModel.hh.

Referenced by CreateCurrentAttValues(), and GetCurrentMaterial().

G4VPhysicalVolume* G4PhysicalVolumeModel::fpCurrentPV
protected

Definition at line 256 of file G4PhysicalVolumeModel.hh.

Referenced by GetCurrentPV(), and GetCurrentTag().

G4Transform3D* G4PhysicalVolumeModel::fpCurrentTransform
protected

Definition at line 259 of file G4PhysicalVolumeModel.hh.

Referenced by CreateCurrentAttValues().

G4VPhysicalVolume* G4PhysicalVolumeModel::fpTopPV
protected

Definition at line 249 of file G4PhysicalVolumeModel.hh.

Referenced by CalculateExtent(), GetTopPhysicalVolume(), and Validate().

G4int G4PhysicalVolumeModel::fRequestedDepth
protected

Definition at line 252 of file G4PhysicalVolumeModel.hh.

Referenced by CalculateExtent(), GetRequestedDepth(), and SetRequestedDepth().

G4int G4PhysicalVolumeModel::fTopPVCopyNo
protected

Definition at line 251 of file G4PhysicalVolumeModel.hh.

Referenced by Validate().

G4String G4PhysicalVolumeModel::fTopPVName
protected

Definition at line 250 of file G4PhysicalVolumeModel.hh.

Referenced by Validate().

G4bool G4PhysicalVolumeModel::fUseFullExtent
protected

Definition at line 254 of file G4PhysicalVolumeModel.hh.

Referenced by CalculateExtent().


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