Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4ASCIITreeSceneHandler.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: G4ASCIITreeSceneHandler.cc 108859 2018-03-12 07:44:53Z gcosmo $
28 //
29 //
30 // John Allison 5th April 2001
31 // A scene handler to dump geometry hierarchy.
32 // Based on a provisional G4ASCIITreeGraphicsScene (was in modeling).
33 
35 
36 #include "G4ASCIITree.hh"
37 #include "G4ASCIITreeMessenger.hh"
38 #include "G4VSolid.hh"
39 #include "G4PhysicalVolumeModel.hh"
40 #include "G4VPhysicalVolume.hh"
41 #include "G4LogicalVolume.hh"
42 #include "G4VPVParameterisation.hh"
43 #include "G4Polyhedron.hh"
44 #include "G4UnitsTable.hh"
45 #include "G4Material.hh"
46 #include "G4Scene.hh"
47 #include "G4ModelingParameters.hh"
49 #include "G4VSensitiveDetector.hh"
50 #include "G4VReadOutGeometry.hh"
52 #include "G4AttCheck.hh"
53 #include "G4AttValue.hh"
54 
57  const G4String& name):
58  G4VTreeSceneHandler(system, name),
59  fpOutFile(0),
60  fpLastPV(0),
61  fLastCopyNo(-99),
62  fLastNonSequentialCopyNo(-99)
63 {}
64 
66 
68 
69  G4VTreeSceneHandler::BeginModeling (); // To re-use "culling off" code.
70 
71  const G4ASCIITree* pSystem = (G4ASCIITree*)GetGraphicsSystem();
72  const G4String& outFileName = pSystem -> GetOutFileName();
73  if (outFileName == "G4cout") {
74  fpOutFile = &G4cout;
75  } else {
76  fOutFile.open (outFileName);
78  }
79 
80  static G4bool firstTime = true;
81  if (firstTime) {
82  firstTime = false;
83  G4cout << "G4ASCIITreeSceneHandler::BeginModeling: writing to ";
84  if (outFileName == "G4cout") {
85  G4cout << "G4 standard output (G4cout)";
86  } else {
87  G4cout << "file \"" << outFileName << "\"";
88  }
89  G4cout << G4endl;
90 
92  }
93 
94  if (outFileName != "G4cout") {
95  WriteHeader (fOutFile); fOutFile << std::endl;
96  }
97 }
98 
99 void G4ASCIITreeSceneHandler::WriteHeader (std::ostream& os)
100 {
101  const G4ASCIITree* pSystem = (G4ASCIITree*)GetGraphicsSystem();
102  const G4int verbosity = pSystem->GetVerbosity();
103  const G4int detail = verbosity % 10;
104  os << "# Set verbosity with \"/vis/ASCIITree/verbose <verbosity>\":";
105  for (size_t i = 0;
107  os << "\n# " << G4ASCIITreeMessenger::fVerbosityGuidance[i];
108  }
109  os << "\n# Now printing with verbosity " << verbosity;
110  os << "\n# Format is: PV:n";
111  if (detail >= 1) os << " / LV (SD,RO)";
112  if (detail >= 2) os << " / Solid(type)";
113  if (detail >= 3) os << ", volume, density";
114  if (detail >= 5) os << ", daughter-subtracted volume and mass";
115  if (detail >= 6) os << ", physical volume dump";
116  if (detail >= 7) os << ", polyhedron dump";
117  os <<
118  "\n# Abbreviations: PV = Physical Volume, LV = Logical Volume,"
119  "\n# SD = Sensitive Detector, RO = Read Out Geometry.";
120 }
121 
123  const G4ASCIITree* pSystem = (G4ASCIITree*) GetGraphicsSystem();
124  const G4int verbosity = pSystem->GetVerbosity();
125  const G4int detail = verbosity % 10;
126  const G4String& outFileName = pSystem -> GetOutFileName();
127 
128  // Output left over copy number, if any...
130  if (fLastCopyNo == fLastNonSequentialCopyNo + 1) *fpOutFile << ',';
131  else *fpOutFile << '-';
133  }
134  // Output outstanding rest of line, if any...
135  if (!fRestOfLine.str().empty()) *fpOutFile << fRestOfLine.str();
136  fRestOfLine.str("");
137  fpLastPV = 0;
138  fLastPVName.clear();
139  fLastCopyNo = -99;
141 
142  // This detail to G4cout regardless of outFileName...
143  if (detail >= 4) {
144  G4cout << "Calculating mass(es)..." << G4endl;
145  const std::vector<G4Scene::Model>& models = fpScene->GetRunDurationModelList();
146  std::vector<G4Scene::Model>::const_iterator i;
147  for (i = models.begin(); i != models.end(); ++i) {
148  G4PhysicalVolumeModel* pvModel =
149  dynamic_cast<G4PhysicalVolumeModel*>(i->fpModel);
150  if (pvModel) {
151  const G4ModelingParameters* tempMP =
152  pvModel->GetModelingParameters();
153  G4ModelingParameters mp; // Default - no culling.
154  pvModel->SetModelingParameters (&mp);
155  G4PhysicalVolumeMassScene massScene(pvModel);
156  pvModel->DescribeYourselfTo (massScene);
157  G4double volume = massScene.GetVolume();
158  G4double mass = massScene.GetMass();
159 
160  G4cout << "Overall volume of \""
161  << pvModel->GetTopPhysicalVolume()->GetName()
162  << "\":"
163  << pvModel->GetTopPhysicalVolume()->GetCopyNo()
164  << ", is "
165  << G4BestUnit (volume, "Volume")
166  << " and the daughter-included mass";
167  G4int requestedDepth = pvModel->GetRequestedDepth();
168  if (requestedDepth == G4PhysicalVolumeModel::UNLIMITED) {
169  G4cout << " to unlimited depth";
170  } else {
171  G4cout << ", ignoring daughters at depth "
172  << requestedDepth
173  << " and below,";
174  }
175  G4cout << " is " << G4BestUnit (mass, "Mass")
176  << G4endl;
177 
178  pvModel->SetModelingParameters (tempMP);
179  }
180  }
181  }
182 
183  if (outFileName != "G4cout") {
184  fOutFile.close();
185  G4cout << "Output file \"" << outFileName << "\" closed." << G4endl;
186  }
187  fLVSet.clear();
188  fReplicaSet.clear();
189  G4cout << "G4ASCIITreeSceneHandler::EndModeling" << G4endl;
190  G4VTreeSceneHandler::EndModeling (); // To re-use "culling off" code.
191 }
192 
194 
195  G4PhysicalVolumeModel* pPVModel =
196  dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
197  if (!pPVModel) return; // Not from a G4PhysicalVolumeModel.
198 
199  // This call comes from a G4PhysicalVolumeModel. drawnPVPath is
200  // the path of the current drawn (non-culled) volume in terms of
201  // drawn (non-culled) ancesters. Each node is identified by a
202  // PVNodeID object, which is a physical volume and copy number. It
203  // is a vector of PVNodeIDs corresponding to the geometry hierarchy
204  // actually selected, i.e., not culled.
205  // The following typedef's already set in header file...
206  //typedef G4PhysicalVolumeModel::G4PhysicalVolumeNodeID PVNodeID;
207  //typedef std::vector<PVNodeID> PVPath;
208  const PVPath& drawnPVPath = pPVModel->GetDrawnPVPath();
209  //G4int currentDepth = pPVModel->GetCurrentDepth();
210  G4VPhysicalVolume* pCurrentPV = pPVModel->GetCurrentPV();
211  const G4String& currentPVName = pCurrentPV->GetName();
212  const G4int currentCopyNo = pCurrentPV->GetCopyNo();
213  G4LogicalVolume* pCurrentLV = pPVModel->GetCurrentLV();
214  G4Material* pCurrentMaterial = pPVModel->GetCurrentMaterial();
215 
217  G4int verbosity = pSystem->GetVerbosity();
218  G4int detail = verbosity % 10;
219 
220  // If verbosity < 10 suppress unnecessary repeated printing.
221  // Repeated simple replicas can always be suppressed.
222  // Paramaterisations can only be suppressed if verbosity < 3, since their
223  // size, density, etc., are in principle different.
224  const G4bool isParameterised = pCurrentPV->GetParameterisation();
225  const G4bool isSimpleReplica = pCurrentPV->IsReplicated() && !isParameterised;
226  const G4bool isAmenableToSupression =
227  (verbosity < 10 && isSimpleReplica) || (verbosity < 3 && isParameterised);
228  if (isAmenableToSupression) {
229  // See if this has been found before with same mother LV...
230  PVPath::const_reverse_iterator thisID = drawnPVPath.rbegin();
231  PVPath::const_reverse_iterator motherID = ++drawnPVPath.rbegin();
232  G4bool ignore = false;
233  for (ReplicaSetIterator i = fReplicaSet.begin(); i != fReplicaSet.end();
234  ++i) {
235  if (i->back().GetPhysicalVolume()->GetLogicalVolume() ==
236  thisID->GetPhysicalVolume()->GetLogicalVolume()) {
237  // For each one previously found (if more than one, they must
238  // have different mothers)...
239  // To avoid compilation errors on VC++ .Net 7.1...
240  // Previously:
241  // PVNodeID previousMotherID = ++(i->rbegin());
242  // (Should that have been: PVNodeID::const_iterator previousMotherID?)
243  // Replace
244  // previousMotherID == i->rend()
245  // by
246  // i->size() <= 1
247  // Replace
248  // previousMotherID != i->rend()
249  // by
250  // i->size() > 1
251  // Replace
252  // previousMotherID->
253  // by
254  // (*i)[i->size() - 2].
255  if (motherID == drawnPVPath.rend() &&
256  i->size() <= 1)
257  ignore = true; // Both have no mother.
258  if (motherID != drawnPVPath.rend() &&
259  i->size() > 1 &&
260  motherID->GetPhysicalVolume()->GetLogicalVolume() ==
261  (*i)[i->size() - 2].GetPhysicalVolume()->GetLogicalVolume())
262  ignore = true; // Same mother LV...
263  }
264  }
265  if (ignore) {
266  pPVModel->CurtailDescent();
267  return;
268  }
269  }
270 
271  // Now suppress printing for volumes with the same name but different
272  // copy number - but not if they are parameterisations that have not been
273  // taken out above (those that are "amenable to suppression" will have been
274  // taken out).
275  if (verbosity < 10 && !isParameterised &&
276  currentPVName == fLastPVName &&
277  currentCopyNo != fLastCopyNo) {
278  // Check...
279  if (isAmenableToSupression) {
280  G4Exception("G4ASCIITreeSceneHandler::RequestPrimitives",
281  "vistree0001",
282  JustWarning,
283  "Volume amenable to suppressed printing unexpected");
284  }
285  // Check mothers are identical...
286  else if (pCurrentLV == (fpLastPV? fpLastPV->GetLogicalVolume(): 0)) {
287  if (currentCopyNo != fLastCopyNo + 1) {
288  // Non-sequential copy number...
289  *fpOutFile << ',' << currentCopyNo;
290  fLastNonSequentialCopyNo = currentCopyNo;
291  }
292  fLastCopyNo = currentCopyNo;
293  pPVModel->CurtailDescent();
294  return;
295  }
296  }
297  fpLastPV = pCurrentPV;
298 
299  // High verbosity or a new or otherwise non-amenable volume...
300  // Output copy number, if any, from previous invocation...
302  if (fLastCopyNo == fLastNonSequentialCopyNo + 1) *fpOutFile << ',';
303  else *fpOutFile << '-';
305  }
306  // Output rest of line, if any, from previous invocation...
307  if (!fRestOfLine.str().empty()) *fpOutFile << fRestOfLine.str();
308  fRestOfLine.str("");
309  fLastPVName = currentPVName;
310  fLastCopyNo = currentCopyNo;
311  fLastNonSequentialCopyNo = currentCopyNo;
312  // Start next line...
313  // Indent according to drawn path depth...
314  for (size_t i = 0; i < drawnPVPath.size(); i++ ) *fpOutFile << " ";
315  *fpOutFile << "\"" << currentPVName
316  << "\":" << currentCopyNo;
317 
318  if (pCurrentPV->IsReplicated()) {
319  if (verbosity < 10) {
320  // Add printing for replicas (when replicas are ignored)...
321  EAxis axis;
322  G4int nReplicas;
323  G4double width;
324  G4double offset;
325  G4bool consuming;
326  pCurrentPV->GetReplicationData(axis,nReplicas,width,offset,consuming);
327  G4VPVParameterisation* pP = pCurrentPV->GetParameterisation();
328  if (pP) {
329  if (detail < 3) {
330  fReplicaSet.insert(drawnPVPath);
331  if (nReplicas > 2) fRestOfLine << '-';
332  else fRestOfLine << ',';
333  fRestOfLine << nReplicas - 1
334  << " (" << nReplicas << " parametrised volumes)";
335  }
336  }
337  else {
338  fReplicaSet.insert(drawnPVPath);
339  if (nReplicas > 2) fRestOfLine << '-';
340  else fRestOfLine << ',';
341  fRestOfLine << nReplicas - 1
342  << " (" << nReplicas << " replicas)";
343  }
344  }
345  } else {
346  if (fLVSet.find(pCurrentLV) != fLVSet.end()) {
347  if (verbosity < 10) {
348  // Add printing for repeated LV (if it has daughters)...
349  if (pCurrentLV->GetNoDaughters()) fRestOfLine << " (repeated LV)";
350  // ...and curtail descent.
351  pPVModel->CurtailDescent();
352  }
353  }
354  }
355 
356  if (detail >= 1) {
357  fRestOfLine << " / \""
358  << pCurrentLV->GetName() << "\"";
359  G4VSensitiveDetector* sd = pCurrentLV->GetSensitiveDetector();
360  if (sd) {
361  fRestOfLine << " (SD=\"" << sd->GetName() << "\"";
362  G4VReadOutGeometry* roGeom = sd->GetROgeometry();
363  if (roGeom) {
364  fRestOfLine << ",RO=\"" << roGeom->GetName() << "\"";
365  }
366  fRestOfLine << ")";
367  }
368  }
369 
370  if (detail >= 2) {
371  fRestOfLine << " / \""
372  << solid.GetName()
373  << "\"("
374  << solid.GetEntityType() << ")";
375  }
376 
377  if (detail >= 3) {
378  fRestOfLine << ", "
379  << G4BestUnit(((G4VSolid&)solid).GetCubicVolume(),"Volume")
380  << ", ";
381  if (pCurrentMaterial) {
383  << G4BestUnit(pCurrentMaterial->GetDensity(), "Volumic Mass")
384  << " (" << pCurrentMaterial->GetName() << ")";
385  } else {
386  fRestOfLine << "(No material)";
387  }
388  }
389 
390  if (detail >= 5) {
391  if (pCurrentMaterial) {
392  G4Material* pMaterial = const_cast<G4Material*>(pCurrentMaterial);
393  // ...and find daughter-subtracted mass...
394  G4double daughter_subtracted_mass = pCurrentLV->GetMass
395  (pCurrentPV->IsParameterised(), // Force if parametrised.
396  false, // Do not propagate - calculate for this volume minus
397  // volume of daughters.
398  pMaterial);
399  G4double daughter_subtracted_volume =
400  daughter_subtracted_mass / pCurrentMaterial->GetDensity();
401  fRestOfLine << ", "
402  << G4BestUnit(daughter_subtracted_volume,"Volume")
403  << ", "
404  << G4BestUnit(daughter_subtracted_mass,"Mass");
405  }
406  }
407 
408  if (detail >= 6) {
409  std::vector<G4AttValue>* attValues = pPVModel->CreateCurrentAttValues();
410  const std::map<G4String,G4AttDef>* attDefs = pPVModel->GetAttDefs();
411  fRestOfLine << '\n' << G4AttCheck(attValues, attDefs);
412  delete attValues;
413  }
414 
415  if (detail >= 7) {
416  G4Polyhedron* polyhedron = solid.GetPolyhedron();
417  fRestOfLine << "\nLocal polyhedron coordinates:\n" << *polyhedron;
418  G4Transform3D* transform = pPVModel->GetCurrentTransform();
419  polyhedron->Transform(*transform);
420  fRestOfLine << "\nGlobal polyhedron coordinates:\n" << *polyhedron;
421  }
422 
423  if (fLVSet.find(pCurrentLV) == fLVSet.end()) {
424  fLVSet.insert(pCurrentLV); // Record new logical volume.
425  }
426 
427  fRestOfLine << std::endl;
428 
429  return;
430 }
const std::vector< Model > & GetRunDurationModelList() const
const XML_Char * name
Definition: expat.h:151
const G4VPhysicalVolume * fpLastPV
system("rm -rf microbeam.root")
void DescribeYourselfTo(G4VGraphicsScene &)
G4LogicalVolume * GetLogicalVolume() const
virtual G4bool IsParameterised() const =0
G4VSensitiveDetector * GetSensitiveDetector() const
std::set< G4LogicalVolume * > fLVSet
std::vector< PVNodeID > PVPath
virtual G4Polyhedron * GetPolyhedron() const
Definition: G4VSolid.cc:665
virtual G4VPVParameterisation * GetParameterisation() const =0
#define G4endl
Definition: G4ios.hh:61
G4int GetNoDaughters() const
const G4ModelingParameters * GetModelingParameters() const
virtual void RequestPrimitives(const G4VSolid &)
void WriteHeader(std::ostream &)
G4ASCIITreeSceneHandler(G4VGraphicsSystem &system, const G4String &name)
void SetModelingParameters(const G4ModelingParameters *)
const G4String & GetName() const
Definition: G4Material.hh:179
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
#define width
G4VPhysicalVolume * GetTopPhysicalVolume() const
G4VReadOutGeometry * GetROgeometry() const
G4Transform3D * GetCurrentTransform() const
G4String GetName() const
virtual G4int GetCopyNo() const =0
const std::map< G4String, G4AttDef > * GetAttDefs() const
G4LogicalVolume * GetCurrentLV() const
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
G4int GetVerbosity() const
Definition: G4ASCIITree.hh:47
G4Material * GetCurrentMaterial() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
int G4int
Definition: G4Types.hh:78
G4VPhysicalVolume * GetCurrentPV() const
EAxis
Definition: geomdefs.hh:54
G4String GetName() const
G4VGraphicsSystem * GetGraphicsSystem() const
const std::vector< G4PhysicalVolumeNodeID > & GetDrawnPVPath() const
virtual void GetReplicationData(EAxis &axis, G4int &nReplicas, G4double &width, G4double &offset, G4bool &consuming) const =0
G4GLOB_DLL std::ostream G4cout
static std::vector< G4String > fVerbosityGuidance
virtual G4bool IsReplicated() const =0
virtual void EndModeling()
std::set< PVPath >::iterator ReplicaSetIterator
virtual G4GeometryType GetEntityType() const =0
std::vector< G4AttValue > * CreateCurrentAttValues() const
G4double GetMass(G4bool forced=false, G4bool propagate=true, G4Material *parMaterial=0)
const G4String & GetName() const
const G4String & GetName() const
virtual void BeginModeling()
G4double GetDensity() const
Definition: G4Material.hh:181