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