Belle II Software  release-05-02-19
VRMLWriterModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2011 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Leo Piilonen *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <geometry/modules/vrmlWriter/VRMLWriterModule.h>
12 #include <geometry/GeometryManager.h>
13 
14 #include "G4PhysicalVolumeStore.hh"
15 #include "G4LogicalVolumeStore.hh"
16 #include "G4SolidStore.hh"
17 #include "G4VPhysicalVolume.hh"
18 #include "G4LogicalVolume.hh"
19 #include "G4VSolid.hh"
20 #include "G4VisAttributes.hh"
21 #include "G4VisExtent.hh"
22 #include "G4Material.hh"
23 #include "G4ThreeVector.hh"
24 #include "G4RotationMatrix.hh"
25 #include "G4Transform3D.hh"
26 #include "G4VPVParameterisation.hh"
27 #include <G4Tubs.hh>
28 #include <G4Box.hh>
29 #include <G4Polyhedron.hh>
30 #include <G4DisplacedSolid.hh>
31 
32 #include <iomanip>
33 
34 using namespace Belle2;
35 
36 REG_MODULE(VRMLWriter);
37 
38 VRMLWriterModule::VRMLWriterModule()
39 {
40  //Set module properties and the description
41  setDescription("Write the detector geometry in a hierarchical VRML format.");
42 
43  //Parameter definition
44  addParam("outputFile", m_Filename, "Output filename", std::string("belle2.wrl"));
45 }
46 
47 void VRMLWriterModule::initialize()
48 {
49  m_First = true;
50 }
51 
52 void VRMLWriterModule::event()
53 {
54  if (!m_First) return;
55  m_First = false;
56  G4VPhysicalVolume* topVol = geometry::GeometryManager::getInstance().getTopVolume();
57  if (!topVol) {
58  B2ERROR("No geometry found: add the Geometry module to the path before the VRMLWriter module.");
59  return;
60  }
61  m_File.open(m_Filename, std::ios_base::trunc);
62  writePreamble();
63 
64  G4PhysicalVolumeStore* pvStore = G4PhysicalVolumeStore::GetInstance();
65  G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
66  G4SolidStore* solidStore = G4SolidStore::GetInstance();
67 
68  m_PVName = new std::vector<std::string>(pvStore->size(), "");
69  m_LVName = new std::vector<std::string>(lvStore->size(), "");
70  m_SolidName = new std::vector<std::string>(solidStore->size(), "");
71 
72  // Assign legal and unique names to each used physical volume, logical volume and solid
73  for (G4VPhysicalVolume* physVol : *pvStore) {
74  int pvIndex = std::find(pvStore->begin(), pvStore->end(), physVol) - pvStore->begin();
75  if ((*m_PVName)[pvIndex].length() == 0)
76  assignName(m_PVName, pvIndex, (*pvStore)[pvIndex]->GetName(), 0);
77  G4LogicalVolume* logVol = physVol->GetLogicalVolume();
78  int lvIndex = std::find(lvStore->begin(), lvStore->end(), logVol) - lvStore->begin();
79  if ((*m_LVName)[lvIndex].length() == 0)
80  assignName(m_LVName, lvIndex, (*lvStore)[lvIndex]->GetName(), 1);
81  G4VSolid* solid = logVol->GetSolid();
82  int solidIndex = std::find(solidStore->begin(), solidStore->end(), solid) - solidStore->begin();
83  if ((*m_SolidName)[solidIndex].length() == 0)
84  assignName(m_SolidName, solidIndex, (*solidStore)[solidIndex]->GetName(), 2);
85  }
86 
87  // Write all explicitly-referenced solids as PROTOs (replicas are written later)
88  // Use implicit prefix "solid_" to avoid name clashes with logical- and physical-volume names
89  m_IsCylinder = new std::vector<bool>(solidStore->size(), false);
90  for (unsigned int solidIndex = 0; solidIndex < solidStore->size(); ++solidIndex) {
91  if ((*m_SolidName)[solidIndex].length() != 0) {
92  if ((*solidStore)[solidIndex]->GetEntityType() == "G4Tubs") {
93  auto* tube = (G4Tubs*)((*solidStore)[solidIndex]);
94  (*m_IsCylinder)[solidIndex] = ((tube->GetInnerRadius() == 0.0) && (tube->GetDeltaPhiAngle() == 2.0 * M_PI));
95  }
96  describeSolid((*solidStore)[solidIndex], (*m_SolidName)[solidIndex], (*m_IsCylinder)[solidIndex]);
97  }
98  }
99 
100  // Recursively write all physical volumes (as DEFs) and logical volumes (as PROTOs).
101  // Deepest volumes are written first; top volume is written last. Use implicit prefix
102  // "lv_" for logical-volume names to avoid name clashes with solid and physical-volume names.
103  m_PVIndex = new std::vector<std::vector<int> >(lvStore->size(), std::vector<int>());
104  m_LVWritten = new std::vector<bool>(lvStore->size(), false);
105  m_PVWritten = new std::vector<bool>(pvStore->size(), false);
106  describePhysicalVolume(topVol);
107 
108  // Now tell VRML to draw the top physical volume (and, recursively, all daughters)
109  int pvIndex = std::find(pvStore->begin(), pvStore->end(), topVol) - pvStore->begin();
110  int lvIndex = std::find(lvStore->begin(), lvStore->end(), topVol->GetLogicalVolume()) - lvStore->begin();
111  m_File << "# PhysicalVolume " << topVol->GetName() << std::endl <<
112  "DEF " << (*m_PVName)[pvIndex] << " Transform {" << std::endl <<
113  " children lv_" << (*m_LVName)[lvIndex] << "{}" << std::endl <<
114  "}" << std::endl;
115  m_File.close();
116  B2INFO("VRML written to " << m_Filename);
117 
118  delete m_PVName;
119  delete m_LVName;
120  delete m_SolidName;
121  delete m_IsCylinder;
122  delete m_PVIndex;
123  delete m_LVWritten;
124  delete m_PVWritten;
125 
126 }
127 
128 void VRMLWriterModule::assignName(std::vector<std::string>* names, unsigned int index, const G4String& originalName, int select)
129 {
130  G4SolidStore* solidStore = G4SolidStore::GetInstance();
131  G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
132  G4PhysicalVolumeStore* pvStore = G4PhysicalVolumeStore::GetInstance();
133 
134  G4String name = originalName;
135  if (name.length() == 0) { name = "anonymous"; }
136  // Replace problematic characters with underscore (there may be more!)
137  for (char c : " .,:;?'\"*+-=|^!/@#$\\%{}[]()<>") std::replace(name.begin(), name.end(), c, '_');
138  // Avoid duplicate names for solids that will be written to VRML file
139  for (int j = (int)index - 1; j >= 0; --j) {
140  if ((*names)[j].length() == 0) continue;
141  int match = 0;
142  switch (select) {
143  case 0:
144  match = (*pvStore)[j]->GetName().compare((*pvStore)[index]->GetName()); break;
145  case 1:
146  match = (*lvStore)[j]->GetName().compare((*lvStore)[index]->GetName()); break;
147  case 2:
148  match = (*solidStore)[j]->GetName().compare((*solidStore)[index]->GetName()); break;
149  }
150  if (match == 0) {
151  if (name.length() == (*names)[j].length()) {
152  (*names)[j].append("_1");
153  }
154  int n = std::stoi((*names)[j].substr(name.length() + 1), nullptr);
155  name.append("_");
156  name.append(std::to_string(n + 1));
157  break;
158  }
159  }
160  (*names)[index] = name;
161 }
162 
163 void VRMLWriterModule::describeSolid(G4VSolid* solid, const std::string& solidName, bool isCylinder)
164 {
165  m_File << "# Solid " << solid->GetName() << " of type " << solid->GetEntityType() << std::endl;
166  if (isCylinder) {
167  auto* tube = (G4Tubs*)solid;
168  // VRML cylinder is along y axis but G4Tubs is along z axis => rotate in logical volume
169  m_File << "PROTO solid_" << solidName << " [ ] {" << std::endl <<
170  " Cylinder {" << std::endl << std::setprecision(10) <<
171  " radius " << tube->GetOuterRadius() << std::endl <<
172  " height " << tube->GetZHalfLength() * 2.0 << std::endl <<
173  " }" << std::endl <<
174  "}" << std::endl;
175  } else if (solid->GetEntityType() == "G4Box") {
176  auto* box = (G4Box*)solid;
177  m_File << "PROTO solid_" << solidName << " [ ] {" << std::endl <<
178  " Box {" << std::endl << std::setprecision(10) <<
179  " size " << box->GetXHalfLength() * 2.0 << " " <<
180  box->GetYHalfLength() * 2.0 << " " <<
181  box->GetZHalfLength() * 2.0 << std::endl <<
182  " }" << std::endl <<
183  "}" << std::endl;
184  } else if ((solid->GetEntityType() == "G4IntersectionSolid") ||
185  (solid->GetEntityType() == "G4UnionSolid") ||
186  (solid->GetEntityType() == "G4SubtractionSolid") ||
187  (solid->GetEntityType() == "G4BooleanSolid")) {
188  HepPolyhedron* polyhedron = getBooleanSolidPolyhedron(solid);
189  auto* g4polyhedron = new G4Polyhedron(*polyhedron);
190  writePolyhedron(g4polyhedron, solidName);
191  delete polyhedron;
192  delete g4polyhedron;
193  } else {
194  writePolyhedron(solid->GetPolyhedron(), solidName);
195  }
196 }
197 
198 void VRMLWriterModule::describeLogicalVolume(G4LogicalVolume* logVol, const std::string& lvName, const std::string& solidName,
199  bool isCylinder)
200 {
201  G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
202  int lvIndex = std::find(lvStore->begin(), lvStore->end(), logVol) - lvStore->begin();
203  if ((*m_LVWritten)[lvIndex]) return; // at most one PROTO in VRML for each logical volume
204  (*m_LVWritten)[lvIndex] = true;
205  G4Color color(0.0, 1.0, 0.0, 0.5); // default is semi-transparent green
206  if ((lvName.compare(0, 23, "eclBarrelCrystalLogical") == 0) ||
207  (lvName.compare(0, 20, "eclFwdCrystalLogical") == 0) ||
208  (lvName.compare(0, 20, "eclBwdCrystalLogical") == 0)) {
209  color = G4Color(1.0, 0.25, 0.0, 0.7); // orange since ECL crystals have no G4VisAttribute :(
210  }
211  std::string visible = "";
212  G4String materialName = logVol->GetMaterial()->GetName();
213  // Hide containers that have vacuum, air or gas
214  if (materialName == "Vacuum") visible = "#";
215  if (materialName == "G4_AIR") visible = "#";
216  if (materialName == "CDCGas") visible = "#";
217  if (materialName == "ColdAir") visible = "#";
218  if (materialName == "STR-DryAir") visible = "#";
219  if (materialName == "TOPAir") visible = "#";
220  if (materialName == "TOPVacuum") visible = "#";
221  const G4VisAttributes* visAttr = logVol->GetVisAttributes();
222  if (visAttr) {
223  color = const_cast<G4Color&>(logVol->GetVisAttributes()->GetColor());
224  if (!(visAttr->IsVisible())) visible = "#";
225  } else {
226  visible = "#";
227  }
228  bool hasDaughters = (visible.length() == 0) || ((*m_PVIndex)[lvIndex].size() > 0);
229  if (logVol->GetSensitiveDetector() != nullptr) visible = "";
230  m_File << "# LogicalVolume " << logVol->GetName() << " containing " << materialName << std::endl <<
231  "PROTO lv_" << lvName << " [ ] {" << std::endl <<
232  " Group {" << std::endl <<
233  " " << (hasDaughters ? "" : "#") << "children [" << std::endl;
234  if (isCylinder) {
235  // VRML cylinder is along y axis but G4Tubs is along z axis => rotate here
236  m_File << " " << visible << "Transform {" << std::endl <<
237  " " << visible << " rotation 1 0 0 " << M_PI_2 << std::endl <<
238  " " << visible << " children [" << std::endl;
239  visible.append(" ");
240  }
241  m_File << " " << visible << "Shape {" << std::endl <<
242  " " << visible << " appearance Appearance {" << std::endl <<
243  " " << visible << " material Material {" << std::endl <<
244  " " << visible << " diffuseColor " << color.GetRed() << " " <<
245  color.GetGreen() << " " <<
246  color.GetBlue() << std::endl <<
247  " " << visible << " " << (color.GetAlpha() == 1.0 ? "#" : "") <<
248  "transparency " << 1.0 - color.GetAlpha() << std::endl <<
249  " " << visible << " }" << std::endl <<
250  " " << visible << " }" << std::endl <<
251  " " << visible << " geometry solid_" << solidName << "{}" << std::endl <<
252  " " << visible << "}" << std::endl;
253  if (isCylinder) {
254  visible.resize(visible.length() - 2);
255  m_File << " " << visible << " ]" << std::endl <<
256  " " << visible << "}" << std::endl;
257  }
258  G4PhysicalVolumeStore* pvStore = G4PhysicalVolumeStore::GetInstance();
259  for (int daughter : (*m_PVIndex)[lvIndex]) {
260  G4VPhysicalVolume* physVol = (*pvStore)[daughter];
261  std::string pvNameDaughter = (*m_PVName)[daughter];
262  int lvDaughter = std::find(lvStore->begin(), lvStore->end(), physVol->GetLogicalVolume()) - lvStore->begin();
263  std::string lvNameDaughter = (*m_LVName)[lvDaughter];
264  if (physVol->IsReplicated()) {
265  G4VSolid* solid = logVol->GetSolid();
266  EAxis axis;
267  G4int nReplicas;
268  G4double width;
269  G4double offset;
270  G4bool consuming;
271  physVol->GetReplicationData(axis, nReplicas, width, offset, consuming);
272  G4VPVParameterisation* physParameterisation = physVol->GetParameterisation();
273  if (physParameterisation) { // parameterised volume
274  for (int replica = 0; replica < nReplicas; ++replica) {
275  std::string pvNameReplica = pvNameDaughter;
276  std::string lvNameReplica = lvNameDaughter;
277  pvNameReplica.append("_");
278  pvNameReplica.append(std::to_string(replica));
279  physVol->SetCopyNo(replica);
280  physParameterisation->ComputeTransformation(replica, physVol);
281  G4VSolid* solidReplica = physParameterisation->ComputeSolid(replica, physVol);
282  if (solidReplica != solid) { // not sure if this works ...
283  lvNameReplica.append("_");
284  lvNameReplica.append(std::to_string(replica));
285  }
286  writePhysicalVolume(physVol, pvNameReplica, lvNameReplica, false);
287  }
288  } else { // plain replicated volume
289  G4ThreeVector originalTranslation = physVol->GetTranslation();
290  G4RotationMatrix* originalRotation = physVol->GetRotation();
291  for (int replica = 0; replica < nReplicas; ++replica) {
292  G4ThreeVector translation; // No translation
293  G4RotationMatrix rotation; // No rotation
294  std::string pvNameReplica = pvNameDaughter;
295  std::string lvNameReplica = lvNameDaughter;
296  pvNameReplica.append("_");
297  pvNameReplica.append(std::to_string(replica));
298  physVol->SetCopyNo(replica);
299  switch (axis) {
300  default:
301  case kXAxis:
302  translation.setX(width * (replica - 0.5 * (nReplicas - 1)));
303  physVol->SetTranslation(translation);
304  break;
305  case kYAxis:
306  translation.setY(width * (replica - 0.5 * (nReplicas - 1)));
307  physVol->SetTranslation(translation);
308  break;
309  case kZAxis:
310  translation.setZ(width * (replica - 0.5 * (nReplicas - 1)));
311  physVol->SetTranslation(translation);
312  break;
313  case kRho:
314  if (solid->GetEntityType() == "G4Tubs") {
315  lvNameReplica.append("_");
316  lvNameReplica.append(std::to_string(replica));
317  } else {
318  B2WARNING("Built-in volumes replicated along radius for " << solid->GetEntityType() <<
319  " (solid " << solid->GetName() << ") are not visualisable.");
320  }
321  break;
322  case kPhi:
323  physVol->SetRotation(&(rotation.rotateZ(-(offset + (replica + 0.5) * width))));
324  break;
325  }
326  writePhysicalVolume(physVol, pvNameReplica, lvNameReplica, false);
327  }
328  // Restore originals...
329  physVol->SetTranslation(originalTranslation);
330  physVol->SetRotation(originalRotation);
331  }
332  } else {
333  writePhysicalVolume(physVol, pvNameDaughter, lvNameDaughter, (*m_PVWritten)[daughter]);
334  }
335  (*m_PVWritten)[daughter] = true;
336  }
337  G4VisExtent extent = logVol->GetSolid()->GetExtent();
338  double xMin = extent.GetXmin();
339  double xMax = extent.GetXmax();
340  double yMin = extent.GetYmin();
341  double yMax = extent.GetYmax();
342  double zMin = extent.GetZmin();
343  double zMax = extent.GetZmax();
344  double xCenter = 0.5 * (xMin + xMax);
345  double yCenter = 0.5 * (yMin + yMax);
346  double zCenter = 0.5 * (zMin + zMax);
347  bool atOrigin = (std::fabs(xCenter) + std::fabs(yCenter) + std::fabs(zCenter) > 1.0e-12);
348  m_File << " " << (hasDaughters ? "" : "#") << "]" << std::endl <<
349  " " << (atOrigin ? "" : "#") <<
350  "bboxCenter " << xCenter << " " << yCenter << " " << zCenter << std::endl <<
351  " bboxSize " << xMax - xMin << " " <<
352  (isCylinder ? zMax - zMin : yMax - yMin) << " " <<
353  (isCylinder ? yMax - yMin : zMax - zMin) << std::endl <<
354  " }" << std::endl <<
355  "} # end of " << lvName << std::endl;
356 }
357 
358 void VRMLWriterModule::describePhysicalVolume(G4VPhysicalVolume* physVol)
359 {
360  G4LogicalVolume* logVol = physVol->GetLogicalVolume();
361  G4VSolid* solid = logVol->GetSolid();
362  G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
363  int lvIndex = std::find(lvStore->begin(), lvStore->end(), logVol) - lvStore->begin();
364  std::string lvName = (*m_LVName)[lvIndex];
365  if (physVol->IsReplicated()) {
366  G4SolidStore* solidStore = G4SolidStore::GetInstance();
367  int solidIndex = std::find(solidStore->begin(), solidStore->end(), solid) - solidStore->begin();
368  std::string solidName = (*m_SolidName)[solidIndex];
369  EAxis axis;
370  G4int nReplicas;
371  G4double width;
372  G4double offset;
373  G4bool consuming;
374  physVol->GetReplicationData(axis, nReplicas, width, offset, consuming);
375  G4VPVParameterisation* physParameterisation = physVol->GetParameterisation();
376  if (physParameterisation) { // user-parameterised volume
377  for (int replica = 0; replica < nReplicas; ++replica) {
378  std::string lvNameReplica = lvName;
379  G4VSolid* solidReplica = physParameterisation->ComputeSolid(replica, physVol);
380  if (solidReplica != solid) { // not sure if this works ...
381  solidReplica->ComputeDimensions(physParameterisation, replica, physVol);
382  std::string solidNameReplica = solidName;
383  solidNameReplica.append("_");
384  solidNameReplica.append(std::to_string(replica));
385  lvNameReplica.append("_");
386  lvNameReplica.append(std::to_string(replica));
387  describeSolid(solidReplica, solidNameReplica, (*m_IsCylinder)[solidIndex]);
388  describeLogicalVolume(logVol, lvNameReplica, solidNameReplica, (*m_IsCylinder)[solidIndex]);
389  }
390  descendAndDescribe(physVol, lvNameReplica, replica);
391  }
392  } else { // plain replicated volume - "divisions"
393  for (int replica = 0; replica < nReplicas; ++replica) {
394  if (axis == kRho) {
395  if (solid->GetEntityType() == "G4Tubs") {
396  double originalRMin = ((G4Tubs*)solid)->GetInnerRadius();
397  double originalRMax = ((G4Tubs*)solid)->GetOuterRadius();
398  ((G4Tubs*)solid)->SetInnerRadius(offset + width * replica);
399  ((G4Tubs*)solid)->SetOuterRadius(offset + width * (replica + 1));
400  std::string solidNameReplica = solidName;
401  solidNameReplica.append("_");
402  solidNameReplica.append(std::to_string(replica));
403  std::string lvNameReplica = lvName;
404  lvNameReplica.append("_");
405  lvNameReplica.append(std::to_string(replica));
406  describeSolid(solid, solidNameReplica, (*m_IsCylinder)[solidIndex]);
407  ((G4Tubs*)solid)->SetInnerRadius(originalRMin);
408  ((G4Tubs*)solid)->SetOuterRadius(originalRMax);
409  describeLogicalVolume(logVol, lvNameReplica, solidNameReplica, (*m_IsCylinder)[solidIndex]);
410  descendAndDescribe(physVol, lvNameReplica, replica);
411  } else {
412  B2WARNING("Built-in volumes replicated along radius for " << solid->GetEntityType() <<
413  " (solid " << solid->GetName() << ") are not visualisable.");
414  }
415  } else {
416  descendAndDescribe(physVol, lvName, replica);
417  }
418  }
419  }
420  } else { // non-replicated volume
421  descendAndDescribe(physVol, lvName, -1);
422  }
423 }
424 
425 void VRMLWriterModule::descendAndDescribe(G4VPhysicalVolume* physVol, const std::string& lvName, int replica)
426 {
427  G4PhysicalVolumeStore* pvStore = G4PhysicalVolumeStore::GetInstance();
428  G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
429  G4SolidStore* solidStore = G4SolidStore::GetInstance();
430 
431  int pvIndex = std::find(pvStore->begin(), pvStore->end(), physVol) - pvStore->begin();
432  if ((*m_PVWritten)[pvIndex]) return;
433 
434  // Descend to the leaves of the tree
435  G4LogicalVolume* logVol = physVol->GetLogicalVolume();
436  int lvIndex = std::find(lvStore->begin(), lvStore->end(), logVol) - lvStore->begin();
437  for (int daughter = 0; daughter < logVol->GetNoDaughters(); ++daughter) {
438  G4VPhysicalVolume* dPhysVol = logVol->GetDaughter(daughter);
439  (*m_PVIndex)[lvIndex].push_back(std::find(pvStore->begin(), pvStore->end(), dPhysVol) - pvStore->begin());
440  describePhysicalVolume(dPhysVol);
441  }
442  // Write out the physical volume and its corresponding logical volume as we ascend the recursive tree
443  int solidIndex = std::find(solidStore->begin(), solidStore->end(), logVol->GetSolid()) - solidStore->begin();
444  if (replica <= 0) describeLogicalVolume(logVol, lvName, (*m_SolidName)[solidIndex], (*m_IsCylinder)[solidIndex]);
445 }
446 
447 void VRMLWriterModule::writePreamble()
448 {
449  m_File << "#VRML V2.0 utf8" << std::endl << std::endl <<
450  "WorldInfo {" << std::endl <<
451  " info [ \"(c) The Belle II Collaboration\" ]" << std::endl <<
452  " title \"The Belle II Detector\"" << std::endl <<
453  "}" << std::endl <<
454  "Viewpoint {" << std::endl <<
455  " position 0 0 15000" << std::endl <<
456  " description \"front\"" << std::endl <<
457  "}" << std::endl <<
458  "Viewpoint {" << std::endl <<
459  " position 0 0 -15000" << std::endl <<
460  " orientation 0 1 0 3.141593" << std::endl <<
461  " description \"back\"" << std::endl <<
462  "}" << std::endl <<
463  "Viewpoint {" << std::endl <<
464  " position 0 15000 0" << std::endl <<
465  " orientation 1 0 0 -1.570796" << std::endl <<
466  " description \"top\"" << std::endl <<
467  "}" << std::endl <<
468  "Viewpoint {" << std::endl <<
469  " position 0 -15000 0" << std::endl <<
470  " orientation 1 0 0 1.570796" << std::endl <<
471  " description \"bottom\"" << std::endl <<
472  "}" << std::endl <<
473  "Viewpoint {" << std::endl <<
474  " position 15000 0 0" << std::endl <<
475  " orientation 0 1 0 1.570796" << std::endl <<
476  " description \"right\"" << std::endl <<
477  "}" << std::endl <<
478  "Viewpoint {" << std::endl <<
479  " position -15000 0 0" << std::endl <<
480  " orientation 0 1 0 -1.570796" << std::endl <<
481  " description \"left\"" << std::endl <<
482  "}" << std::endl;
483 }
484 
485 void VRMLWriterModule::writePolyhedron(const G4Polyhedron* polyhedron, const std::string& name)
486 {
487  if (polyhedron) {
488  m_File << "PROTO solid_" << name << " [ ] {" << std::endl <<
489  " IndexedFaceSet {" << std::endl <<
490  " coord Coordinate {" << std::endl <<
491  " point [" << std::endl;
492  for (int j = 1; j <= polyhedron->GetNoVertices(); ++j) {
493  m_File << " " << polyhedron->GetVertex(j).x() << " " <<
494  polyhedron->GetVertex(j).y() << " " <<
495  polyhedron->GetVertex(j).z() << std::endl;
496  }
497  m_File << " ]" << std::endl <<
498  " }" << std::endl <<
499  " coordIndex [" << std::endl;
500  for (int j = 0; j < polyhedron->GetNoFacets(); ++j) {
501  G4bool notLastEdge = true;
502  G4int ndx = -1, edgeFlag = 1;
503  m_File << " ";
504  do {
505  notLastEdge = polyhedron->GetNextVertexIndex(ndx, edgeFlag);
506  m_File << " " << ndx - 1;
507  } while (notLastEdge);
508  m_File << " -1" << std::endl;
509  }
510  m_File << " ]" << std::endl <<
511  " }" << std::endl <<
512  "}" << std::endl;
513  } else {
514  B2INFO("Polyhedron representation of solid " << name << " cannot be created");
515  }
516 }
517 
518 void VRMLWriterModule::writePhysicalVolume(const G4VPhysicalVolume* physVol, const std::string& pvName, const std::string& lvName,
519  bool written)
520 {
521  if (written) {
522  m_File << " USE " << pvName << std::endl; // this never happens, as it turns out
523  } else {
524  m_File << " # PhysicalVolume " << physVol->GetName() << " copy " << physVol->GetCopyNo() << std::endl <<
525  " DEF " << pvName << " Transform {" << std::endl; // DEF provides object name when importing VRML into Unity
526  G4RotationMatrix* rot = physVol->GetObjectRotation();
527  G4ThreeVector move = physVol->GetObjectTranslation();
528  double yaw = std::atan2(rot->yx(), rot->xx());
529  double pitch = -std::asin(rot->zx());
530  double roll = std::atan2(rot->zy(), rot->zz());
531  std::cout << physVol->GetName() << " translation: " << move.x() << "," << move.y() << "," << move.z() << std::endl;
532  std::cout << physVol->GetName() << " rotation: " << yaw * 180.0 / M_PI << "," << pitch * 180.0 / M_PI << "," << roll * 180.0 / M_PI
533  << std::endl;
534  if (move.mag2() != 0.0) {
535  m_File << " translation " << std::setprecision(10) <<
536  move.x() << " " << move.y() << " " << move.z() << std::endl;
537  }
538  if (!(rot->isIdentity())) {
539  double trace = std::max(-1.0, std::min(3.0, rot->xx() + rot->yy() + rot->zz()));
540  double angle = std::acos(0.5 * (trace - 1.0));
541  G4ThreeVector rotA(rot->zy() - rot->yz(), rot->xz() - rot->zx(), rot->yx() - rot->xy());
542  if ((rotA.x() == 0.0) && (rotA.y() == 0.0) && (rotA.z() == 0.0)) {
543  // The assignment along x, y or z axis in Hep3Vector::axis() is wrong. The proper axis of
544  // rotation is given by kernel(r - I), which is the cross product of two non-proportional rows
545  if (rot->xx() > 0.0) {
546  rotA.setX((rot->yy() - 1.0) * (rot->zz() - 1.0) - rot->yz() * rot->zy());
547  rotA.setY(rot->yz() * rot->zx() - rot->yx() * (rot->zz() - 1.0));
548  rotA.setZ(rot->yx() * rot->zy() - rot->zx() * (rot->yy() - 1.0));
549  } else if (rot->yy() > 0.0) {
550  rotA.setX(rot->xy() * (rot->zz() - 1.0) - rot->xz() * rot->zy());
551  rotA.setY(rot->xz() * rot->zx() - (rot->xx() - 1.0) * (rot->zz() - 1.0));
552  rotA.setZ((rot->xx() - 1.0) * rot->zy() - rot->xy() * rot->zx());
553  } else {
554  rotA.setX(rot->xy() * rot->yz() - rot->xz() * (rot->yy() - 1.0));
555  rotA.setY(rot->xz() * rot->yx() - rot->yz() * (rot->xx() - 1.0));
556  rotA.setZ((rot->xx() - 1.0) * (rot->yy() - 1.0) - rot->xy() * rot->yx());
557  }
558  }
559  rotA.setMag(1.0);
560  m_File << " rotation " << std::setprecision(10) <<
561  rotA.x() << " " << rotA.y() << " " << rotA.z() << " " << angle << std::endl;
562  }
563  m_File << " children lv_" << lvName << "{}" << std::endl <<
564  " }" << std::endl;
565  }
566 }
567 
568 // The code in GEANT4 geometry/solids/Boolean/src/G4BooleanSolid.cc is buggy so avoid it.
569 // Recursively combine the polyhedrons of two solids to make a resulting polyhedron.
570 HepPolyhedron* VRMLWriterModule::getBooleanSolidPolyhedron(G4VSolid* solid)
571 {
572  G4VSolid* solidA = solid->GetConstituentSolid(0);
573  G4VSolid* solidB = solid->GetConstituentSolid(1);
574  HepPolyhedron* polyhedronA = nullptr;
575  if ((solidA->GetEntityType() == "G4IntersectionSolid") ||
576  (solidA->GetEntityType() == "G4UnionSolid") ||
577  (solidA->GetEntityType() == "G4SubtractionSolid") ||
578  (solidA->GetEntityType() == "G4BooleanSolid")) {
579  polyhedronA = getBooleanSolidPolyhedron(solidA);
580  } else {
581  polyhedronA = new HepPolyhedron(*(solidA->GetPolyhedron()));
582  }
583  HepPolyhedron* polyhedronB = nullptr;
584  G4VSolid* solidB2 = solidB;
585  if (solidB->GetEntityType() == "G4DisplacedSolid") {
586  solidB2 = ((G4DisplacedSolid*)solidB)->GetConstituentMovedSolid();
587  }
588  if ((solidB2->GetEntityType() == "G4IntersectionSolid") ||
589  (solidB2->GetEntityType() == "G4UnionSolid") ||
590  (solidB2->GetEntityType() == "G4SubtractionSolid") ||
591  (solidB2->GetEntityType() == "G4BooleanSolid")) {
592  polyhedronB = getBooleanSolidPolyhedron(solidB2);
593  if (solidB != solidB2) { // was solidB a G4DisplacedSolid?
594  polyhedronB->Transform(G4Transform3D(
595  ((G4DisplacedSolid*)solidB)->GetObjectRotation(),
596  ((G4DisplacedSolid*)solidB)->GetObjectTranslation()));
597  }
598  } else {
599  polyhedronB = new HepPolyhedron(*(solidB->GetPolyhedron()));
600  }
601  auto* result = new HepPolyhedron();
602  if (solid->GetEntityType() == "G4UnionSolid") {
603  *result = polyhedronA->add(*polyhedronB);
604  } else if (solid->GetEntityType() == "G4SubtractionSolid") {
605  *result = polyhedronA->subtract(*polyhedronB);
606  } else if (solid->GetEntityType() == "G4IntersectionSolid") {
607  *result = polyhedronA->intersect(*polyhedronB);
608  } else {
609  B2WARNING("getBooleanSolidPolyhedron(): Unrecognized boolean solid " << solid->GetName() <<
610  " of type " << solid->GetEntityType());
611  }
612  delete polyhedronA;
613  delete polyhedronB;
614  return result;
615 }
Belle2::VRMLWriterModule::assignName
void assignName(std::vector< std::string > *, unsigned int, const G4String &, int)
Create unique and legal name for each solid.
Definition: VRMLWriterModule.cc:128
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::VRMLWriterModule::writePreamble
void writePreamble(void)
Emit VRML for the start of the file.
Definition: VRMLWriterModule.cc:447
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::geometry::GeometryManager::getTopVolume
G4VPhysicalVolume * getTopVolume()
Return a pointer to the top volume.
Definition: GeometryManager.h:59
Belle2::VRMLWriterModule::m_First
bool m_First
Once-only flag to write VRML only on the first event.
Definition: VRMLWriterModule.h:87
Belle2::VRMLWriterModule::m_IsCylinder
std::vector< bool > * m_IsCylinder
Flag to indicate that a solid can be rendered as a VMRL cylinder.
Definition: VRMLWriterModule.h:105
Belle2::VRMLWriterModule::m_File
std::ofstream m_File
Output file.
Definition: VRMLWriterModule.h:93
Belle2::geometry::GeometryManager::getInstance
static GeometryManager & getInstance()
Return a reference to the instance.
Definition: GeometryManager.cc:98
Belle2::VRMLWriterModule::m_Filename
std::string m_Filename
User-specified output filename.
Definition: VRMLWriterModule.h:90
Belle2::VRMLWriterModule::VRMLWriterModule
VRMLWriterModule()
Constructor of the module.
Definition: VRMLWriterModule.cc:38
Belle2::VRMLWriterModule::m_SolidName
std::vector< std::string > * m_SolidName
Modified (legal-character and unique) solid name.
Definition: VRMLWriterModule.h:102
Belle2::VRMLWriterModule::m_LVWritten
std::vector< bool > * m_LVWritten
Flag to indicate that the logical volume has already been written.
Definition: VRMLWriterModule.h:111
Belle2::VRMLWriterModule::writePhysicalVolume
void writePhysicalVolume(const G4VPhysicalVolume *, const std::string &, const std::string &, bool)
Emit VRML for each daughter of a logical volume.
Definition: VRMLWriterModule.cc:518
Belle2::VRMLWriterModule::event
void event() override
Called for each event: this runs the VRML writer only for the first event.
Definition: VRMLWriterModule.cc:52
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VRMLWriterModule::describeSolid
void describeSolid(G4VSolid *, const std::string &, bool)
Emit VRML for each solid.
Definition: VRMLWriterModule.cc:163
Belle2::VRMLWriterModule::getBooleanSolidPolyhedron
HepPolyhedron * getBooleanSolidPolyhedron(G4VSolid *)
Create polyhedron for a boolean solid (recursive)
Definition: VRMLWriterModule.cc:570
Belle2::VRMLWriterModule::m_LVName
std::vector< std::string > * m_LVName
Modified (legal-character and unique) logical-volume name.
Definition: VRMLWriterModule.h:99
Belle2::VRMLWriterModule::describeLogicalVolume
void describeLogicalVolume(G4LogicalVolume *, const std::string &, const std::string &, bool)
Emit VRML for each logical volume.
Definition: VRMLWriterModule.cc:198
Belle2::VRMLWriterModule::m_PVName
std::vector< std::string > * m_PVName
Modified (legal-character and unique) physical-volume name.
Definition: VRMLWriterModule.h:96
Belle2::VRMLWriterModule::descendAndDescribe
void descendAndDescribe(G4VPhysicalVolume *, const std::string &, int)
Emit VRML for a physical volume (recursive)
Definition: VRMLWriterModule.cc:425
Belle2::Module::addParam
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition: Module.h:562
Belle2::VRMLWriterModule::initialize
void initialize() override
Initialize at the start of a job.
Definition: VRMLWriterModule.cc:47
Belle2::VRMLWriterModule::describePhysicalVolume
void describePhysicalVolume(G4VPhysicalVolume *)
Access next physical volume in the tree (recursive)
Definition: VRMLWriterModule.cc:358
Belle2::VRMLWriterModule::writePolyhedron
void writePolyhedron(const G4Polyhedron *, const std::string &)
Emit VRML for the solid's polyhedron.
Definition: VRMLWriterModule.cc:485
Belle2::VRMLWriterModule::m_PVWritten
std::vector< bool > * m_PVWritten
Flag to indicate that the physical volume has already been written.
Definition: VRMLWriterModule.h:114
Belle2::VRMLWriterModule::m_PVIndex
std::vector< std::vector< int > > * m_PVIndex
Indices (in G4PhysicalVolumeStore) of the logical volume's physical-volume daughters.
Definition: VRMLWriterModule.h:108