Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4PhysicalVolumeModel.cc
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4PhysicalVolumeModel.cc 110743 2018-06-12 06:33:37Z gcosmo $
28 //
29 //
30 // John Allison 31st December 1997.
31 // Model for physical volumes.
32 
33 #include "G4PhysicalVolumeModel.hh"
34 
35 #include "G4ModelingParameters.hh"
36 #include "G4VGraphicsScene.hh"
37 #include "G4VPhysicalVolume.hh"
38 #include "G4VPVParameterisation.hh"
39 #include "G4LogicalVolume.hh"
40 #include "G4VSolid.hh"
41 #include "G4SubtractionSolid.hh"
42 #include "G4IntersectionSolid.hh"
43 #include "G4Material.hh"
44 #include "G4VisAttributes.hh"
45 #include "G4BoundingSphereScene.hh"
48 #include "G4Polyhedron.hh"
49 #include "HepPolyhedronProcessor.h"
50 #include "G4AttDefStore.hh"
51 #include "G4AttDef.hh"
52 #include "G4AttValue.hh"
53 #include "G4UnitsTable.hh"
54 #include "G4Vector3D.hh"
55 
56 #include <sstream>
57 
60  , G4int requestedDepth
61  , const G4Transform3D& modelTransformation
62  , const G4ModelingParameters* pMP
63  , G4bool useFullExtent)
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)
72 , fpCurrentMaterial (0)
73 , fpCurrentTransform (0)
74 , fAbort (false)
75 , fCurtailDescent (false)
76 , fpClippingSolid (0)
77 , fClippingMode (subtraction)
78 {
79  fType = "G4PhysicalVolumeModel";
80 
81  if (!fpTopPV) {
82 
83  // In some circumstances creating an "empty" G4PhysicalVolumeModel is
84  // allowed, so I have supressed the G4Exception below. If it proves to
85  // be a problem we might have to re-instate it, but it is unlikley to
86  // be used except by visualisation experts. See, for example, /vis/list,
87  // where it is used simply to get a list of G4AttDefs.
88  // G4Exception
89  // ("G4PhysicalVolumeModel::G4PhysicalVolumeModel",
90  // "modeling0010", FatalException, "Null G4PhysicalVolumeModel pointer.");
91 
92  fTopPVName = "NULL";
93  fGlobalTag = "Empty";
94  fGlobalDescription = "G4PhysicalVolumeModel " + fGlobalTag;
95 
96  } else {
97 
98  fTopPVName = fpTopPV -> GetName ();
99  fTopPVCopyNo = fpTopPV -> GetCopyNo ();
100  std::ostringstream o;
101  o << fpTopPV -> GetCopyNo ();
102  fGlobalTag = fpTopPV -> GetName () + "." + o.str();
103  fGlobalDescription = "G4PhysicalVolumeModel " + fGlobalTag;
104 
105  fpCurrentPV = fpTopPV;
106  if (fpCurrentPV) fpCurrentLV = fpCurrentPV->GetLogicalVolume();
107  if (fpCurrentLV) fpCurrentMaterial = fpCurrentLV->GetMaterial();
108  fpCurrentTransform = const_cast<G4Transform3D*>(&modelTransformation);
109 
110  CalculateExtent ();
111  }
112 }
113 
115 {
116  delete fpClippingSolid;
117 }
118 
120 {
121  if (fUseFullExtent) {
122  fExtent = fpTopPV -> GetLogicalVolume () -> GetSolid () -> GetExtent ();
123  }
124  else {
125  G4BoundingSphereScene bsScene(this);
126  const G4int tempRequestedDepth = fRequestedDepth;
127  fRequestedDepth = -1; // Always search to all depths to define extent.
128  const G4ModelingParameters* tempMP = fpMP;
129  G4ModelingParameters mParams
130  (0, // No default vis attributes needed.
131  G4ModelingParameters::wf, // wireframe (not relevant for this).
132  true, // Global culling.
133  true, // Cull invisible volumes.
134  false, // Density culling.
135  0., // Density (not relevant if density culling false).
136  true, // Cull daughters of opaque mothers.
137  24); // No of sides (not relevant for this operation).
138  fpMP = &mParams;
139  DescribeYourselfTo (bsScene);
140  G4double radius = bsScene.GetRadius();
141  if (radius < 0.) { // Nothing in the scene.
142  fExtent = fpTopPV -> GetLogicalVolume () -> GetSolid () -> GetExtent ();
143  } else {
144  // Transform back to coordinates relative to the top
145  // transformation, which is in G4VModel::fTransform. This makes
146  // it conform to all models, which are defined by a
147  // transformation and an extent relative to that
148  // transformation...
149  G4Point3D centre = bsScene.GetCentre();
150  centre.transform(fTransform.inverse());
151  fExtent = G4VisExtent(centre, radius);
152  }
153  fpMP = tempMP;
154  fRequestedDepth = tempRequestedDepth;
155  }
156 }
157 
159 (G4VGraphicsScene& sceneHandler)
160 {
161  if (!fpTopPV) G4Exception
162  ("G4PhysicalVolumeModel::DescribeYourselfTo",
163  "modeling0012", FatalException, "No model.");
164 
165  if (!fpMP) G4Exception
166  ("G4PhysicalVolumeModel::DescribeYourselfTo",
167  "modeling0003", FatalException, "No modeling parameters.");
168 
169  // For safety...
170  fCurrentDepth = 0;
171 
172  G4Transform3D startingTransformation = fTransform;
173 
174  VisitGeometryAndGetVisReps
175  (fpTopPV,
176  fRequestedDepth,
177  startingTransformation,
178  sceneHandler);
179 
180  // Clear data...
181  fCurrentDepth = 0;
182  fpCurrentPV = 0;
183  fpCurrentLV = 0;
184  fpCurrentMaterial = 0;
185  if (fFullPVPath.size() != fBaseFullPVPath.size()) {
186  // They should be equal if pushing and popping is happening properly.
188  ("G4PhysicalVolumeModel::DescribeYourselfTo",
189  "modeling0013",
191  "Path at start of modeling not equal to base path. Something badly"
192  "\nwrong. Please contact visualisation coordinator.");
193  }
194  fDrawnPVPath.clear();
195  fAbort = false;
196  fCurtailDescent = false;
197 }
198 
200 {
201  if (fpCurrentPV) {
202  std::ostringstream o;
203  o << fpCurrentPV -> GetCopyNo ();
204  return fpCurrentPV -> GetName () + "." + o.str();
205  }
206  else {
207  return "WARNING: NO CURRENT VOLUME - global tag is " + fGlobalTag;
208  }
209 }
210 
212 {
213  return "G4PhysicalVolumeModel " + GetCurrentTag ();
214 }
215 
218  G4int requestedDepth,
219  const G4Transform3D& theAT,
220  G4VGraphicsScene& sceneHandler)
221 {
222  // Visits geometry structure to a given depth (requestedDepth), starting
223  // at given physical volume with given starting transformation and
224  // describes volumes to the scene handler.
225  // requestedDepth < 0 (default) implies full visit.
226  // theAT is the Accumulated Transformation.
227 
228  // Find corresponding logical volume and (later) solid, storing in
229  // local variables to preserve re-entrancy.
230  G4LogicalVolume* pLV = pVPV -> GetLogicalVolume ();
231 
232  G4VSolid* pSol = nullptr;
233  G4Material* pMaterial = nullptr;
234 
235  if (!(pVPV -> IsReplicated ())) {
236  // Non-replicated physical volume.
237  pSol = pLV -> GetSolid ();
238  pMaterial = pLV -> GetMaterial ();
239  DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
240  theAT, sceneHandler);
241  }
242  else {
243  // Replicated or parametrised physical volume.
244  EAxis axis;
245  G4int nReplicas;
246  G4double width;
247  G4double offset;
248  G4bool consuming;
249  pVPV -> GetReplicationData (axis, nReplicas, width, offset, consuming);
250  if (fCurrentDepth == 0) nReplicas = 1; // Just draw first
251  G4VPVParameterisation* pP = pVPV -> GetParameterisation ();
252  if (pP) { // Parametrised volume.
253  for (int n = 0; n < nReplicas; n++) {
254  pSol = pP -> ComputeSolid (n, pVPV);
255  pP -> ComputeTransformation (n, pVPV);
256  pSol -> ComputeDimensions (pP, n, pVPV);
257  pVPV -> SetCopyNo (n);
258  // Create a touchable of current parent for ComputeMaterial.
259  // fFullPVPath has not been updated yet so at this point it
260  // corresponds to the parent.
261  G4PhysicalVolumeModelTouchable parentTouchable(fFullPVPath);
262  pMaterial = pP -> ComputeMaterial (n, pVPV, &parentTouchable);
263  DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
264  theAT, sceneHandler);
265  }
266  }
267  else { // Plain replicated volume. From geometry_guide.txt...
268  // The replica's positions are claculated by means of a linear formula.
269  // Replication may occur along:
270  //
271  // o Cartesian axes (kXAxis,kYAxis,kZAxis)
272  //
273  // The replications, of specified width have coordinates of
274  // form (-width*(nReplicas-1)*0.5+n*width,0,0) where n=0.. nReplicas-1
275  // for the case of kXAxis, and are unrotated.
276  //
277  // o Radial axis (cylindrical polar) (kRho)
278  //
279  // The replications are cons/tubs sections, centred on the origin
280  // and are unrotated.
281  // They have radii of width*n+offset to width*(n+1)+offset
282  // where n=0..nReplicas-1
283  //
284  // o Phi axis (cylindrical polar) (kPhi)
285  // The replications are `phi sections' or wedges, and of cons/tubs form
286  // They have phi of offset+n*width to offset+(n+1)*width where
287  // n=0..nReplicas-1
288  //
289  pSol = pLV -> GetSolid ();
290  pMaterial = pLV -> GetMaterial ();
291  G4ThreeVector originalTranslation = pVPV -> GetTranslation ();
292  G4RotationMatrix* pOriginalRotation = pVPV -> GetRotation ();
293  G4double originalRMin = 0., originalRMax = 0.;
294  if (axis == kRho && pSol->GetEntityType() == "G4Tubs") {
295  originalRMin = ((G4Tubs*)pSol)->GetInnerRadius();
296  originalRMax = ((G4Tubs*)pSol)->GetOuterRadius();
297  }
298  G4bool visualisable = true;
299  for (int n = 0; n < nReplicas; n++) {
300  G4ThreeVector translation; // Identity.
301  G4RotationMatrix rotation; // Identity - life enough for visualizing.
302  G4RotationMatrix* pRotation = 0;
303  switch (axis) {
304  default:
305  case kXAxis:
306  translation = G4ThreeVector (-width*(nReplicas-1)*0.5+n*width,0,0);
307  break;
308  case kYAxis:
309  translation = G4ThreeVector (0,-width*(nReplicas-1)*0.5+n*width,0);
310  break;
311  case kZAxis:
312  translation = G4ThreeVector (0,0,-width*(nReplicas-1)*0.5+n*width);
313  break;
314  case kRho:
315  if (pSol->GetEntityType() == "G4Tubs") {
316  ((G4Tubs*)pSol)->SetInnerRadius(width*n+offset);
317  ((G4Tubs*)pSol)->SetOuterRadius(width*(n+1)+offset);
318  } else {
319  if (fpMP->IsWarning())
320  G4cout <<
321  "G4PhysicalVolumeModel::VisitGeometryAndGetVisReps: WARNING:"
322  "\n built-in replicated volumes replicated in radius for "
323  << pSol->GetEntityType() <<
324  "-type\n solids (your solid \""
325  << pSol->GetName() <<
326  "\") are not visualisable."
327  << G4endl;
328  visualisable = false;
329  }
330  break;
331  case kPhi:
332  rotation.rotateZ (-(offset+(n+0.5)*width));
333  // Minus Sign because for the physical volume we need the
334  // coordinate system rotation.
335  pRotation = &rotation;
336  break;
337  }
338  pVPV -> SetTranslation (translation);
339  pVPV -> SetRotation (pRotation);
340  pVPV -> SetCopyNo (n);
341  if (visualisable) {
342  DescribeAndDescend (pVPV, requestedDepth, pLV, pSol, pMaterial,
343  theAT, sceneHandler);
344  }
345  }
346  // Restore originals...
347  pVPV -> SetTranslation (originalTranslation);
348  pVPV -> SetRotation (pOriginalRotation);
349  if (axis == kRho && pSol->GetEntityType() == "G4Tubs") {
350  ((G4Tubs*)pSol)->SetInnerRadius(originalRMin);
351  ((G4Tubs*)pSol)->SetOuterRadius(originalRMax);
352  }
353  }
354  }
355 }
356 
359  G4int requestedDepth,
360  G4LogicalVolume* pLV,
361  G4VSolid* pSol,
362  G4Material* pMaterial,
363  const G4Transform3D& theAT,
364  G4VGraphicsScene& sceneHandler)
365 {
366  // Maintain useful data members...
367  fpCurrentPV = pVPV;
368  fpCurrentLV = pLV;
369  fpCurrentMaterial = pMaterial;
370 
371  const G4RotationMatrix objectRotation = pVPV -> GetObjectRotationValue ();
372  const G4ThreeVector& translation = pVPV -> GetTranslation ();
373  G4Transform3D theLT (G4Transform3D (objectRotation, translation));
374 
375  // Compute the accumulated transformation...
376  // Note that top volume's transformation relative to the world
377  // coordinate system is specified in theAT == startingTransformation
378  // = fTransform (see DescribeYourselfTo), so first time through the
379  // volume's own transformation, which is only relative to its
380  // mother, i.e., not relative to the world coordinate system, should
381  // not be accumulated.
382  G4Transform3D theNewAT (theAT);
383  if (fCurrentDepth != 0) theNewAT = theAT * theLT;
384  fpCurrentTransform = &theNewAT;
385 
386  const G4VisAttributes* pVisAttribs = pLV->GetVisAttributes();
387  // If the volume does not have any vis attributes, create it.
388  G4VisAttributes* tempVisAtts = nullptr;
389  if (!pVisAttribs) {
390  tempVisAtts = new G4VisAttributes; // Default value.
391  // The user may request /vis/viewer/set/colourByDensity.
392  if (fpMP->GetCBDAlgorithmNumber() == 1) {
393  // Algorithm 1: 3 parameters: Simple rainbow mapping.
394  if (fpMP->GetCBDParameters().size() != 3) {
395  G4Exception("G4PhysicalVolumeModelTouchable::DescribeAndDescend",
396  "modeling0014",
398  "Algorithm-parameter mismatch for Colour By Density");
399  } else {
400  const G4double d = pMaterial? pMaterial->GetDensity(): 0.;
401  const G4double d0 = fpMP->GetCBDParameters()[0]; // Invisible d < d0.
402  const G4double d1 = fpMP->GetCBDParameters()[1]; // Rainbow d0->d1->d2.
403  const G4double d2 = fpMP->GetCBDParameters()[2]; // Blue d > d2.
404  if (d < d0) { // Density < d0 is invisible.
405  tempVisAtts->SetVisibility(false);
406  } else { // Intermediate densities are on a spectrum.
407  G4double red, green, blue;
408  if (d < d1) {
409  red = (d1-d)/(d1-d0); green = (d-d0)/(d1-d0); blue = 0.;
410  } else if (d < d2) {
411  red = 0.; green = (d2-d)/(d2-d1); blue = (d-d1)/(d2-d1);
412  } else { // Density >= d2 is blue.
413  red = 0.; green = 0.; blue = 1.;
414  }
415  tempVisAtts->SetColour(G4Colour(red,green,blue));
416  }
417  }
418  } else if (fpMP->GetCBDAlgorithmNumber() == 2) {
419  // Algorithm 2
420  // ...etc.
421  }
422  pVisAttribs = tempVisAtts;
423  }
424  // From here, can assume pVisAttribs is a valid pointer. This is necessary
425  // because PreAddSolid needs a vis attributes object.
426 
427  // Make decision to draw...
428  G4bool thisToBeDrawn = true;
429 
430  // Update full path of physical volumes...
431  G4int copyNo = fpCurrentPV->GetCopyNo();
432  fFullPVPath.push_back
434  (fpCurrentPV,copyNo,fCurrentDepth,*fpCurrentTransform));
435 
436  // Check if vis attributes are to be modified by a /vis/touchable/set/ command.
437  const auto& vams = fpMP->GetVisAttributesModifiers();
438  if (vams.size()) {
439  // OK, we have some VAMs (Vis Attributes Modifiers).
440  for (const auto& vam: vams) {
441  const auto& vamPath = vam.GetPVNameCopyNoPath();
442  if (vamPath.size() == fFullPVPath.size()) {
443  // OK, we have a size match.
444  // Check the volume name/copy number path.
445  auto iVAMNameCopyNo = vamPath.begin();
446  auto iPVNodeId = fFullPVPath.begin();
447  for (; iVAMNameCopyNo != vamPath.end(); ++iVAMNameCopyNo, ++iPVNodeId) {
448  if (!(
449  iVAMNameCopyNo->GetName() ==
450  iPVNodeId->GetPhysicalVolume()->GetName() &&
451  iVAMNameCopyNo->GetCopyNo() ==
452  iPVNodeId->GetPhysicalVolume()->GetCopyNo()
453  )) {
454  // This path element does NOT match.
455  break;
456  }
457  }
458  if (iVAMNameCopyNo == vamPath.end()) {
459  // OK, the paths match (the above loop terminated normally).
460  // Create a vis atts object for the modified vis atts.
461  // It is static so that we may return a reliable pointer to it.
462  static G4VisAttributes modifiedVisAtts;
463  // Initialise it with the current vis atts and reset the pointer.
464  modifiedVisAtts = *pVisAttribs;
465  pVisAttribs = &modifiedVisAtts;
466  const G4VisAttributes& transVisAtts = vam.GetVisAttributes();
467  switch (vam.GetVisAttributesSignifier()) {
469  modifiedVisAtts.SetVisibility(transVisAtts.IsVisible());
470  break;
472  modifiedVisAtts.SetDaughtersInvisible
473  (transVisAtts.IsDaughtersInvisible());
474  break;
476  modifiedVisAtts.SetColour(transVisAtts.GetColour());
477  break;
479  modifiedVisAtts.SetLineStyle(transVisAtts.GetLineStyle());
480  break;
482  modifiedVisAtts.SetLineWidth(transVisAtts.GetLineWidth());
483  break;
485  if (transVisAtts.IsForceDrawingStyle()) {
486  if (transVisAtts.GetForcedDrawingStyle() ==
488  modifiedVisAtts.SetForceWireframe(true);
489  }
490  }
491  break;
493  if (transVisAtts.IsForceDrawingStyle()) {
494  if (transVisAtts.GetForcedDrawingStyle() ==
496  modifiedVisAtts.SetForceSolid(true);
497  }
498  }
499  break;
501  if (transVisAtts.IsForceAuxEdgeVisible()) {
502  modifiedVisAtts.SetForceAuxEdgeVisible
503  (transVisAtts.IsForcedAuxEdgeVisible());
504  }
505  break;
507  modifiedVisAtts.SetForceLineSegmentsPerCircle
508  (transVisAtts.GetForcedLineSegmentsPerCircle());
509  break;
510  }
511  }
512  }
513  }
514  }
515 
516  // There are various reasons why this volume
517  // might not be drawn...
518  G4bool culling = fpMP->IsCulling();
519  G4bool cullingInvisible = fpMP->IsCullingInvisible();
520  G4bool markedVisible = pVisAttribs->IsVisible();
521  G4bool cullingLowDensity = fpMP->IsDensityCulling();
522  G4double density = pMaterial? pMaterial->GetDensity(): 0;
523  G4double densityCut = fpMP -> GetVisibleDensity ();
524 
525  // 1) Global culling is on....
526  if (culling) {
527  // 2) Culling of invisible volumes is on...
528  if (cullingInvisible) {
529  // 3) ...and the volume is marked not visible...
530  if (!markedVisible) thisToBeDrawn = false;
531  }
532  // 4) Or culling of low density volumes is on...
533  if (cullingLowDensity) {
534  // 5) ...and density is less than cut value...
535  if (density < densityCut) thisToBeDrawn = false;
536  }
537  }
538  // 6) The user has asked for all further traversing to be aborted...
539  if (fAbort) thisToBeDrawn = false;
540 
541  // Record thisToBeDrawn in path...
542  fFullPVPath.back().SetDrawn(thisToBeDrawn);
543 
544  if (thisToBeDrawn) {
545 
546  // Update path of drawn physical volumes...
547  fDrawnPVPath.push_back
549  (fpCurrentPV,copyNo,fCurrentDepth,*fpCurrentTransform,thisToBeDrawn));
550 
551  if (fpMP->IsExplode() && fDrawnPVPath.size() == 1) {
552  // For top-level drawn volumes, explode along radius...
553  G4Transform3D centering = G4Translate3D(fpMP->GetExplodeCentre());
554  G4Transform3D centred = centering.inverse() * theNewAT;
555  G4Scale3D oldScale;
556  G4Rotate3D oldRotation;
557  G4Translate3D oldTranslation;
558  centred.getDecomposition(oldScale, oldRotation, oldTranslation);
559  G4double explodeFactor = fpMP->GetExplodeFactor();
560  G4Translate3D newTranslation =
561  G4Translate3D(explodeFactor * oldTranslation.dx(),
562  explodeFactor * oldTranslation.dy(),
563  explodeFactor * oldTranslation.dz());
564  theNewAT = centering * newTranslation * oldRotation * oldScale;
565  }
566 
567  DescribeSolid (theNewAT, pSol, pVisAttribs, sceneHandler);
568 
569  }
570 
571  // Make decision to draw daughters, if any. There are various
572  // reasons why daughters might not be drawn...
573 
574  // First, reasons that do not depend on culling policy...
575  G4int nDaughters = pLV->GetNoDaughters();
576  G4bool daughtersToBeDrawn = true;
577  // 1) There are no daughters...
578  if (!nDaughters) daughtersToBeDrawn = false;
579  // 2) We are at the limit if requested depth...
580  else if (requestedDepth == 0) daughtersToBeDrawn = false;
581  // 3) The user has asked for all further traversing to be aborted...
582  else if (fAbort) daughtersToBeDrawn = false;
583  // 4) The user has asked that the descent be curtailed...
584  else if (fCurtailDescent) daughtersToBeDrawn = false;
585 
586  // Now, reasons that depend on culling policy...
587  else {
588  G4bool daughtersInvisible = pVisAttribs->IsDaughtersInvisible();
589  // Culling of covered daughters request. This is computed in
590  // G4VSceneHandler::CreateModelingParameters() depending on view
591  // parameters...
592  G4bool cullingCovered = fpMP->IsCullingCovered();
593  G4bool surfaceDrawing =
594  fpMP->GetDrawingStyle() == G4ModelingParameters::hsr ||
595  fpMP->GetDrawingStyle() == G4ModelingParameters::hlhsr;
596  if (pVisAttribs->IsForceDrawingStyle()) {
597  switch (pVisAttribs->GetForcedDrawingStyle()) {
598  default:
599  case G4VisAttributes::wireframe: surfaceDrawing = false; break;
600  case G4VisAttributes::solid: surfaceDrawing = true; break;
601  }
602  }
603  G4bool opaque = pVisAttribs->GetColour().GetAlpha() >= 1.;
604  // 5) Global culling is on....
605  if (culling) {
606  // 6) ..and culling of invisible volumes is on...
607  if (cullingInvisible) {
608  // 7) ...and the mother requests daughters invisible
609  if (daughtersInvisible) daughtersToBeDrawn = false;
610  }
611  // 8) Or culling of covered daughters is requested...
612  if (cullingCovered) {
613  // 9) ...and surface drawing is operating...
614  if (surfaceDrawing) {
615  // 10) ...but only if mother is visible...
616  if (thisToBeDrawn) {
617  // 11) ...and opaque...
618  if (opaque) daughtersToBeDrawn = false;
619  }
620  }
621  }
622  }
623  }
624 
625  if (daughtersToBeDrawn) {
626  for (G4int iDaughter = 0; iDaughter < nDaughters; iDaughter++) {
627  // Store daughter pVPV in local variable ready for recursion...
628  G4VPhysicalVolume* pDaughterVPV = pLV -> GetDaughter (iDaughter);
629  // Descend the geometry structure recursively...
630  fCurrentDepth++;
631  VisitGeometryAndGetVisReps
632  (pDaughterVPV, requestedDepth - 1, theNewAT, sceneHandler);
633  fCurrentDepth--;
634  }
635  }
636 
637  delete tempVisAtts;
638 
639  // Reset for normal descending of next volume at this level...
640  fCurtailDescent = false;
641 
642  // Pop item from paths physical volumes...
643  fFullPVPath.pop_back();
644  if (thisToBeDrawn) {
645  fDrawnPVPath.pop_back();
646  }
647 }
648 
650 (const G4Transform3D& theAT,
651  G4VSolid* pSol,
652  const G4VisAttributes* pVisAttribs,
653  G4VGraphicsScene& sceneHandler)
654 {
655  sceneHandler.PreAddSolid (theAT, *pVisAttribs);
656 
657  G4VSolid* pSectionSolid = fpMP->GetSectionSolid();
658  G4VSolid* pCutawaySolid = fpMP->GetCutawaySolid();
659 
660  if (!fpClippingSolid && !pSectionSolid && !pCutawaySolid) {
661 
662  pSol -> DescribeYourselfTo (sceneHandler); // Standard treatment.
663 
664  } else {
665 
666  // Clipping, etc., performed by Boolean operations.
667 
668  // First, get polyhedron for current solid...
669  if (pVisAttribs->IsForceLineSegmentsPerCircle())
671  (pVisAttribs->GetForcedLineSegmentsPerCircle());
672  else
673  G4Polyhedron::SetNumberOfRotationSteps(fpMP->GetNoOfSides());
674  const G4Polyhedron* pOriginal = pSol->GetPolyhedron();
676 
677  if (!pOriginal) {
678 
679  if (fpMP->IsWarning())
680  G4cout <<
681  "WARNING: G4PhysicalVolumeModel::DescribeSolid: solid\n \""
682  << pSol->GetName() <<
683  "\" has no polyhedron. Cannot by clipped."
684  << G4endl;
685  pSol -> DescribeYourselfTo (sceneHandler); // Standard treatment.
686 
687  } else {
688 
689  G4Polyhedron resultant(*pOriginal);
690  G4VisAttributes resultantVisAttribs(*pVisAttribs);
691  G4VSolid* resultantSolid = 0;
692 
693  if (fpClippingSolid) {
694  switch (fClippingMode) {
695  default:
696  case subtraction:
697  resultantSolid = new G4SubtractionSolid
698  ("resultant_solid", pSol, fpClippingSolid, theAT.inverse());
699  break;
700  case intersection:
701  resultantSolid = new G4IntersectionSolid
702  ("resultant_solid", pSol, fpClippingSolid, theAT.inverse());
703  break;
704  }
705  }
706 
707  if (pSectionSolid) {
708  resultantSolid = new G4IntersectionSolid
709  ("sectioned_solid", pSol, pSectionSolid, theAT.inverse());
710  }
711 
712  if (pCutawaySolid) {
713  resultantSolid = new G4SubtractionSolid
714  ("cutaway_solid", pSol, pCutawaySolid, theAT.inverse());
715  }
716 
717  G4Polyhedron* tmpResultant = resultantSolid->GetPolyhedron();
718  if (tmpResultant) resultant = *tmpResultant;
719  else {
720  if (fpMP->IsWarning())
721  G4cout <<
722  "WARNING: G4PhysicalVolumeModel::DescribeSolid: resultant polyhedron for"
723  "\n solid \"" << pSol->GetName() <<
724  "\" not defined due to error during Boolean processing."
725  "\n Original will be drawn in red."
726  << G4endl;
727  resultantVisAttribs.SetColour(G4Colour::Red());
728  }
729 
730  delete resultantSolid;
731 
732  // Finally, force polyhedron drawing...
733  resultant.SetVisAttributes(resultantVisAttribs);
734  sceneHandler.BeginPrimitives(theAT);
735  sceneHandler.AddPrimitive(resultant);
736  sceneHandler.EndPrimitives();
737  }
738  }
739  sceneHandler.PostAddSolid ();
740 }
741 
743 {
744  G4TransportationManager* transportationManager =
746 
747  size_t nWorlds = transportationManager->GetNoWorlds();
748 
749  G4bool found = false;
750 
751  std::vector<G4VPhysicalVolume*>::iterator iterWorld =
752  transportationManager->GetWorldsIterator();
753  for (size_t i = 0; i < nWorlds; ++i, ++iterWorld) {
754  G4VPhysicalVolume* world = (*iterWorld);
755  if (!world) break; // This can happen if geometry has been cleared/destroyed.
756  // The idea now is to seek a PV with the same name and copy no
757  // in the hope it's the same one!!
758  G4PhysicalVolumeModel searchModel (world);
759  G4int verbosity = 0; // Suppress messages from G4PhysicalVolumeSearchScene.
760  G4PhysicalVolumeSearchScene searchScene
761  (&searchModel, fTopPVName, fTopPVCopyNo, verbosity);
762  G4ModelingParameters mp; // Default modeling parameters for this search.
764  searchModel.SetModelingParameters (&mp);
765  searchModel.DescribeYourselfTo (searchScene);
766  G4VPhysicalVolume* foundVolume = searchScene.GetFoundVolume ();
767  if (foundVolume) {
768  if (foundVolume != fpTopPV && warn) {
769  G4cout <<
770  "G4PhysicalVolumeModel::Validate(): A volume of the same name and"
771  "\n copy number (\""
772  << fTopPVName << "\", copy " << fTopPVCopyNo
773  << ") still exists and is being used."
774  "\n But it is not the same volume you originally specified"
775  "\n in /vis/scene/add/."
776  << G4endl;
777  }
778  fpTopPV = foundVolume;
779  CalculateExtent ();
780  found = true;
781  }
782  }
783  if (found) return true;
784  else {
785  if (warn) {
786  G4cout <<
787  "G4PhysicalVolumeModel::Validate(): No volume of name and"
788  "\n copy number (\""
789  << fTopPVName << "\", copy " << fTopPVCopyNo
790  << ") exists."
791  << G4endl;
792  }
793  return false;
794  }
795 }
796 
797 const std::map<G4String,G4AttDef>* G4PhysicalVolumeModel::GetAttDefs() const
798 {
799  G4bool isNew;
800  std::map<G4String,G4AttDef>* store
801  = G4AttDefStore::GetInstance("G4PhysicalVolumeModel", isNew);
802  if (isNew) {
803  (*store)["PVPath"] =
804  G4AttDef("PVPath","Physical Volume Path","Physics","","G4String");
805  (*store)["LVol"] =
806  G4AttDef("LVol","Logical Volume","Physics","","G4String");
807  (*store)["Solid"] =
808  G4AttDef("Solid","Solid Name","Physics","","G4String");
809  (*store)["EType"] =
810  G4AttDef("EType","Entity Type","Physics","","G4String");
811  (*store)["DmpSol"] =
812  G4AttDef("DmpSol","Dump of Solid properties","Physics","","G4String");
813  (*store)["LocalTrans"] =
814  G4AttDef("LocalTrans","Local transformation of volume","Physics","","G4String");
815  (*store)["GlobalTrans"] =
816  G4AttDef("GlobalTrans","Global transformation of volume","Physics","","G4String");
817  (*store)["Material"] =
818  G4AttDef("Material","Material Name","Physics","","G4String");
819  (*store)["Density"] =
820  G4AttDef("Density","Material Density","Physics","G4BestUnit","G4double");
821  (*store)["State"] =
822  G4AttDef("State","Material State (enum undefined,solid,liquid,gas)","Physics","","G4String");
823  (*store)["Radlen"] =
824  G4AttDef("Radlen","Material Radiation Length","Physics","G4BestUnit","G4double");
825  (*store)["Region"] =
826  G4AttDef("Region","Cuts Region","Physics","","G4String");
827  (*store)["RootRegion"] =
828  G4AttDef("RootRegion","Root Region (0/1 = false/true)","Physics","","G4bool");
829  }
830  return store;
831 }
832 
833 #include <iomanip>
834 
835 static std::ostream& operator<< (std::ostream& o, const G4Transform3D t)
836 {
837  using namespace std;
838 
839  G4Scale3D sc;
840  G4Rotate3D r;
841  G4Translate3D tl;
842  t.getDecomposition(sc, r, tl);
843 
844  const int w = 10;
845 
846  // Transformation itself
847  o << setw(w) << t.xx() << setw(w) << t.xy() << setw(w) << t.xz() << setw(w) << t.dx() << endl;
848  o << setw(w) << t.yx() << setw(w) << t.yy() << setw(w) << t.yz() << setw(w) << t.dy() << endl;
849  o << setw(w) << t.zx() << setw(w) << t.zy() << setw(w) << t.zz() << setw(w) << t.dz() << endl;
850 
851  // Translation
852  o << "= translation:" << endl;
853  o << setw(w) << tl.dx() << setw(w) << tl.dy() << setw(w) << tl.dz() << endl;
854 
855  // Rotation
856  o << "* rotation:" << endl;
857  o << setw(w) << r.xx() << setw(w) << r.xy() << setw(w) << r.xz() << endl;
858  o << setw(w) << r.yx() << setw(w) << r.yy() << setw(w) << r.yz() << endl;
859  o << setw(w) << r.zx() << setw(w) << r.zy() << setw(w) << r.zz() << endl;
860 
861  // Scale
862  o << "* scale:" << endl;
863  o << setw(w) << sc.xx() << setw(w) << sc.yy() << setw(w) << sc.zz() << endl;
864 
865  // Transformed axes
866  o << "Transformed axes:" << endl;
867  o << "x': " << r * G4Vector3D(1., 0., 0.) << endl;
868  o << "y': " << r * G4Vector3D(0., 1., 0.) << endl;
869  o << "z': " << r * G4Vector3D(0., 0., 1.) << endl;
870 
871  return o;
872 }
873 
874 std::vector<G4AttValue>* G4PhysicalVolumeModel::CreateCurrentAttValues() const
875 {
876  std::vector<G4AttValue>* values = new std::vector<G4AttValue>;
877  std::ostringstream oss;
878  for (size_t i = 0; i < fFullPVPath.size(); ++i) {
879  oss << fFullPVPath[i].GetPhysicalVolume()->GetName()
880  << ':' << fFullPVPath[i].GetCopyNo();
881  if (i != fFullPVPath.size() - 1) oss << '/';
882  }
883 
884  if (!fpCurrentLV) {
886  ("G4PhysicalVolumeModel::CreateCurrentAttValues",
887  "modeling0004",
888  JustWarning,
889  "Current logical volume not defined.");
890  return values;
891  }
892 
893  values->push_back(G4AttValue("PVPath", oss.str(),""));
894  values->push_back(G4AttValue("LVol", fpCurrentLV->GetName(),""));
895  G4VSolid* pSol = fpCurrentLV->GetSolid();
896  values->push_back(G4AttValue("Solid", pSol->GetName(),""));
897  values->push_back(G4AttValue("EType", pSol->GetEntityType(),""));
898  oss.str(""); oss << '\n' << *pSol;
899  values->push_back(G4AttValue("DmpSol", oss.str(),""));
900  const G4RotationMatrix localRotation = fpCurrentPV->GetObjectRotationValue();
901  const G4ThreeVector& localTranslation = fpCurrentPV->GetTranslation();
902  oss.str(""); oss << '\n' << G4Transform3D(localRotation,localTranslation);
903  values->push_back(G4AttValue("LocalTrans", oss.str(),""));
904  oss.str(""); oss << '\n' << *fpCurrentTransform;
905  values->push_back(G4AttValue("GlobalTrans", oss.str(),""));
906  G4String matName = fpCurrentMaterial? fpCurrentMaterial->GetName(): G4String("No material");
907  values->push_back(G4AttValue("Material", matName,""));
909  values->push_back(G4AttValue("Density", G4BestUnit(matDensity,"Volumic Mass"),""));
911  oss.str(""); oss << matState;
912  values->push_back(G4AttValue("State", oss.str(),""));
914  values->push_back(G4AttValue("Radlen", G4BestUnit(matRadlen,"Length"),""));
915  G4Region* region = fpCurrentLV->GetRegion();
916  G4String regionName = region? region->GetName(): G4String("No region");
917  values->push_back(G4AttValue("Region", regionName,""));
918  oss.str(""); oss << fpCurrentLV->IsRootRegion();
919  values->push_back(G4AttValue("RootRegion", oss.str(),""));
920  return values;
921 }
922 
923 G4bool G4PhysicalVolumeModel::G4PhysicalVolumeNodeID::operator<
925 {
926  if (fpPV < right.fpPV) return true;
927  if (fpPV == right.fpPV) {
928  if (fCopyNo < right.fCopyNo) return true;
929  if (fCopyNo == right.fCopyNo)
930  return fNonCulledDepth < right.fNonCulledDepth;
931  }
932  return false;
933 }
934 
935 std::ostream& operator<<
936  (std::ostream& os, const G4PhysicalVolumeModel::G4PhysicalVolumeNodeID& node)
937 {
938  G4VPhysicalVolume* pPV = node.GetPhysicalVolume();
939  if (pPV) {
940  os << pPV->GetName()
941  << ' ' << node.GetCopyNo()
942 // << '[' << node.GetNonCulledDepth() << ']'
943 // << ':' << node.GetTransform()
944  ;
945 // os << " (";
946 // if (!node.GetDrawn()) os << "not ";
947 // os << "drawn)";
948  } else {
949  os << " (Null node)";
950  }
951  return os;
952 }
953 
954 std::ostream& operator<<
955 (std::ostream& os, const std::vector<G4PhysicalVolumeModel::G4PhysicalVolumeNodeID>& path)
956 {
957  for (const auto& nodeID: path) {
958  os << nodeID << ' ';
959  }
960  return os;
961 }
962 
964 (const std::vector<G4PhysicalVolumeNodeID>& fullPVPath):
965  fFullPVPath(fullPVPath) {}
966 
968 {
969  size_t i = fFullPVPath.size() - depth - 1;
970  if (i >= fFullPVPath.size()) {
971  G4Exception("G4PhysicalVolumeModelTouchable::GetTranslation",
972  "modeling0005",
974  "Index out of range. Asking for non-existent depth");
975  }
976  static G4ThreeVector tempTranslation;
977  tempTranslation = fFullPVPath[i].GetTransform().getTranslation();
978  return tempTranslation;
979 }
980 
982 {
983  size_t i = fFullPVPath.size() - depth - 1;
984  if (i >= fFullPVPath.size()) {
985  G4Exception("G4PhysicalVolumeModelTouchable::GetRotation",
986  "modeling0006",
988  "Index out of range. Asking for non-existent depth");
989  }
990  static G4RotationMatrix tempRotation;
991  tempRotation = fFullPVPath[i].GetTransform().getRotation();
992  return &tempRotation;
993 }
994 
996 {
997  size_t i = fFullPVPath.size() - depth - 1;
998  if (i >= fFullPVPath.size()) {
999  G4Exception("G4PhysicalVolumeModelTouchable::GetVolume",
1000  "modeling0007",
1002  "Index out of range. Asking for non-existent depth");
1003  }
1004  return fFullPVPath[i].GetPhysicalVolume();
1005 }
1006 
1008 {
1009  size_t i = fFullPVPath.size() - depth - 1;
1010  if (i >= fFullPVPath.size()) {
1011  G4Exception("G4PhysicalVolumeModelTouchable::GetSolid",
1012  "modeling0008",
1014  "Index out of range. Asking for non-existent depth");
1015  }
1016  return fFullPVPath[i].GetPhysicalVolume()->GetLogicalVolume()->GetSolid();
1017 }
1018 
1020 {
1021  size_t i = fFullPVPath.size() - depth - 1;
1022  if (i >= fFullPVPath.size()) {
1023  G4Exception("G4PhysicalVolumeModelTouchable::GetReplicaNumber",
1024  "modeling0009",
1026  "Index out of range. Asking for non-existent depth");
1027  }
1028  return fFullPVPath[i].GetCopyNo();
1029 }
G4Transform3D * fpCurrentTransform
Definition: geomdefs.hh:54
void SetDefaultVisAttributes(const G4VisAttributes *pDefaultVisAttributes)
CLHEP::Hep3Vector G4ThreeVector
void SetForceAuxEdgeVisible(G4bool=true)
const G4RotationMatrix * GetRotation(G4int depth) const
G4VPhysicalVolume * fpTopPV
void DescribeYourselfTo(G4VGraphicsScene &)
Definition: test07.cc:36
std::ostream & operator<<(std::ostream &, const BasicVector3D< float > &)
HepGeom::Transform3D G4Transform3D
HepGeom::Translate3D G4Translate3D
virtual G4Polyhedron * GetPolyhedron() const
Definition: G4VSolid.cc:665
void SetLineStyle(LineStyle)
double xy() const
Definition: Transform3D.h:255
Definition: G4Tubs.hh:85
G4VPhysicalVolume * fpCurrentPV
static void ResetNumberOfRotationSteps()
const G4String & GetName() const
double yz() const
Definition: Transform3D.h:267
#define G4endl
Definition: G4ios.hh:61
Transform3D inverse() const
Definition: Transform3D.cc:142
virtual void EndPrimitives()=0
G4double GetRadlen() const
Definition: G4Material.hh:221
ForcedDrawingStyle GetForcedDrawingStyle() const
void VisitGeometryAndGetVisReps(G4VPhysicalVolume *, G4int requestedDepth, const G4Transform3D &, G4VGraphicsScene &)
G4String fGlobalTag
Definition: G4VModel.hh:109
G4VPhysicalVolume * GetFoundVolume() const
G4int GetNoDaughters() const
double dx() const
Definition: Transform3D.h:279
std::vector< G4VPhysicalVolume * >::iterator GetWorldsIterator()
G4VisExtent fExtent
Definition: G4VModel.hh:111
G4bool IsVisible() const
double zy() const
Definition: Transform3D.h:273
const std::vector< G4PhysicalVolumeNodeID > & fFullPVPath
const G4ThreeVector & GetTranslation() const
Definition: test07.cc:36
G4bool IsForceLineSegmentsPerCircle() const
static void SetNumberOfRotationSteps(G4int n)
G4bool IsForceAuxEdgeVisible() const
virtual void BeginPrimitives(const G4Transform3D &objectTransformation=G4Transform3D())=0
G4String GetCurrentDescription() const
void SetModelingParameters(const G4ModelingParameters *)
const G4String & GetName() const
Definition: G4Material.hh:179
double yy() const
Definition: Transform3D.h:264
G4bool IsForcedAuxEdgeVisible() const
void SetForceLineSegmentsPerCircle(G4int nSegments)
double zz() const
Definition: Transform3D.h:276
double yx() const
Definition: Transform3D.h:261
const G4VisExtent & GetExtent() const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
#define width
G4bool Validate(G4bool warn)
G4State GetState() const
Definition: G4Material.hh:182
static const G4double d2
const G4Point3D & GetCentre() const
double xx() const
Definition: Transform3D.h:252
G4String GetCurrentTag() const
static const G4double d1
G4double GetLineWidth() const
const G4ModelingParameters * fpMP
Definition: G4VModel.hh:113
Float_t d
G4RotationMatrix GetObjectRotationValue() const
const std::map< G4String, G4AttDef > * GetAttDefs() const
Double_t radius
const G4VisAttributes * GetVisAttributes() const
size_t GetNoWorlds() const
virtual void AddPrimitive(const G4Polyline &)=0
G4bool IsRootRegion() const
double dz() const
Definition: Transform3D.h:285
HepGeom::Vector3D< G4double > G4Vector3D
Definition: G4Vector3D.hh:35
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
const G4VisAttributes * GetDefaultVisAttributes() const
static G4TransportationManager * GetTransportationManager()
double xz() const
Definition: Transform3D.h:258
void SetForceSolid(G4bool=true)
const G4ThreeVector & GetTranslation(G4int depth) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
void SetVisibility(G4bool=true)
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)
G4PhysicalVolumeModel(G4VPhysicalVolume *=0, G4int requestedDepth=UNLIMITED, const G4Transform3D &modelTransformation=G4Transform3D(), const G4ModelingParameters *=0, G4bool useFullExtent=false)
int G4int
Definition: G4Types.hh:78
G4bool IsForceDrawingStyle() const
EAxis
Definition: geomdefs.hh:54
static G4Colour Red()
Definition: G4Colour.hh:160
G4String GetName() const
virtual void PostAddSolid()=0
void SetDaughtersInvisible(G4bool=true)
G4double GetAlpha() const
Definition: G4Colour.hh:154
void DescribeAndDescend(G4VPhysicalVolume *, G4int requestedDepth, G4LogicalVolume *, G4VSolid *, G4Material *, const G4Transform3D &, G4VGraphicsScene &)
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:92
G4GLOB_DLL std::ostream G4cout
std::vector< G4PhysicalVolumeNodeID > fFullPVPath
double dy() const
Definition: Transform3D.h:282
void SetLineWidth(G4double)
G4VSolid * GetSolid() const
void SetColour(const G4Colour &)
void SetForceWireframe(G4bool=true)
Char_t n[5]
void getDecomposition(Scale3D &scale, Rotate3D &rotation, Translate3D &translation) const
Definition: Transform3D.cc:174
G4State
Definition: G4Material.hh:115
double zx() const
Definition: Transform3D.h:270
Definition: geomdefs.hh:54
virtual G4GeometryType GetEntityType() const =0
virtual void DescribeSolid(const G4Transform3D &theAT, G4VSolid *pSol, const G4VisAttributes *pVisAttribs, G4VGraphicsScene &sceneHandler)
LineStyle GetLineStyle() const
G4int GetForcedLineSegmentsPerCircle() const
G4PhysicalVolumeModelTouchable(const std::vector< G4PhysicalVolumeNodeID > &fullPVPath)
G4Region * GetRegion() const
std::vector< G4AttValue > * CreateCurrentAttValues() const
void SetVisAttributes(const G4VisAttributes *)
Definition: G4Visible.cc:80
G4Transform3D fTransform
Definition: G4VModel.hh:112
const G4Colour & GetColour() const
const G4String & GetName() const
const G4String & GetName() const
G4bool IsDaughtersInvisible() const
virtual void PreAddSolid(const G4Transform3D &objectTransformation, const G4VisAttributes &visAttribs)=0
G4double GetDensity() const
Definition: G4Material.hh:181