DD4hep  1.30.0
Detector Description Toolkit for High Energy Physics
GeometryTreeDump.cpp
Go to the documentation of this file.
1 //==========================================================================
2 // AIDA Detector description implementation
3 //--------------------------------------------------------------------------
4 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
5 // All rights reserved.
6 //
7 // For the licensing terms see $DD4hepINSTALL/LICENSE.
8 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
9 //
10 // Author : M.Frank
11 //
12 //==========================================================================
13 
14 // Framework include files
15 #include "GeometryTreeDump.h"
16 #include <DD4hep/Detector.h>
17 
18 // ROOT includes
19 #include <TROOT.h>
20 #include <TColor.h>
21 #include <TGeoShape.h>
22 #include <TGeoCone.h>
23 #include <TGeoParaboloid.h>
24 #include <TGeoPcon.h>
25 #include <TGeoPgon.h>
26 #include <TGeoSphere.h>
27 #include <TGeoTorus.h>
28 #include <TGeoTube.h>
29 #include <TGeoTrd1.h>
30 #include <TGeoTrd2.h>
31 #include <TGeoArb8.h>
32 #include <TGeoMatrix.h>
33 #include <TGeoBoolNode.h>
34 #include <TGeoCompositeShape.h>
35 #include <TClass.h>
36 #include <TGeoManager.h>
37 #include <TMath.h>
38 
39 // C/C++ include files
40 #include <iostream>
41 
42 using namespace dd4hep;
43 
44 namespace {
45  std::string indent = "";
46 
48  std::ostream& m_output = std::cout;
49 
50  void getAngles(const Double_t* m, Double_t &phi, Double_t &theta, Double_t &psi) {
51  // Retreive Euler angles.
52  // Check if theta is 0 or 180.
53  if (TMath::Abs(1. - TMath::Abs(m[8])) < 1.e-9) {
54  theta = TMath::ACos(m[8]) * RAD_2_DEGREE;
55  phi = TMath::ATan2(-m[8] * m[1], m[0]) * RAD_2_DEGREE;
56  psi = 0.; // convention, phi+psi matters
57  return;
58  }
59  // sin(theta) != 0
60  phi = TMath::ATan2(m[2], -m[5]);
61  Double_t sphi = TMath::Sin(phi);
62  if (TMath::Abs(sphi) < 1.e-9)
63  theta = -TMath::ASin(m[5] / TMath::Cos(phi)) * RAD_2_DEGREE;
64  else
65  theta = TMath::ASin(m[2] / sphi) * RAD_2_DEGREE;
66  phi *= RAD_2_DEGREE;
67  psi = TMath::ATan2(m[6], m[7]) * RAD_2_DEGREE;
68  }
69 }
70 
72 void* detail::GeometryTreeDump::handleVolume(const std::string& name, Volume vol) const {
73  VisAttr vis = vol.visAttributes();
74  TGeoShape* shape = vol->GetShape();
75  TGeoMedium* medium = vol->GetMedium();
76  int num = vol->GetNdaughters();
77 
78  m_output << "\t\t<volume name=\"" << name << "\">" << std::endl;
79  m_output << "\t\t\t<solidref ref=\"" << shape->GetName() << "\"/>" << std::endl;
80  m_output << "\t\t\t<materialref ref=\"" << medium->GetName() << "\"/>" << std::endl;
81  if (vis.isValid()) {
82  m_output << "\t\t\t<visref ref=\"" << vis.name() << "\"/>" << std::endl;
83  }
84  if (num > 0) {
85  for (int i = 0; i < num; ++i) {
86  TGeoNode* geo_nod = vol.ptr()->GetNode(vol->GetNode(i)->GetName());
87  TGeoVolume* geo_vol = geo_nod->GetVolume();
88  TGeoMatrix* geo_mat = geo_nod->GetMatrix();
89  m_output << "\t\t\t<physvol>" << std::endl;
90  m_output << "\t\t\t\t<volumeref ref=\"" << geo_vol->GetName() << "\"/>" << std::endl;
91  if (geo_mat) {
92  if (geo_mat->IsTranslation()) {
93  m_output << "\t\t\t\t<positionref ref=\"" << geo_nod->GetName() << "_pos\"/>" << std::endl;
94  }
95  if (geo_mat->IsRotation()) {
96  m_output << "\t\t\t\t<rotationref ref=\"" << geo_nod->GetName() << "_rot\"/>" << std::endl;
97  }
98  }
99  m_output << "\t\t\t</physvol>" << std::endl;
100  }
101  }
102  m_output << "\t\t</volume>" << std::endl;
103  return 0;
104 }
105 
107 void* detail::GeometryTreeDump::handleSolid(const std::string& name, const TGeoShape* shape) const {
108  if (shape) {
109  if (shape->IsA() == TGeoBBox::Class()) {
110  const TGeoBBox* sh = (const TGeoBBox*) shape;
111  m_output << "\t\t<box name=\"" << name << "_shape\" x=\"" << sh->GetDX() << "\" y=\"" << sh->GetDY() << "\" z=\""
112  << sh->GetDZ() << "\" lunit=\"cm\"/>" << std::endl;
113  }
114  else if (shape->IsA() == TGeoTube::Class()) {
115  const TGeoTube* sh = (const TGeoTube*) shape;
116  m_output << "\t\t<tube name=\"" << name << "_shape\" rmin=\"" << sh->GetRmin() << "\" rmax=\"" << sh->GetRmax() << "\" z=\""
117  << sh->GetDz() << "\" startphi=\"0.0\" deltaphi=\"360.0\" aunit=\"deg\" lunit=\"cm\"/>" << std::endl;
118  }
119  else if (shape->IsA() == TGeoTubeSeg::Class()) {
120  const TGeoTubeSeg* sh = (const TGeoTubeSeg*) shape;
121  m_output << "\t\t<tube name=\"" << name << "_shape\" rmin=\"" << sh->GetRmin() << "\" rmax=\"" << sh->GetRmax() << "\" z=\""
122  << sh->GetDz() << "\" startphi=\"" << sh->GetPhi1() << "\" deltaphi=\"" << sh->GetPhi2()
123  << "\" aunit=\"deg\" lunit=\"cm\"/>" << std::endl;
124  }
125  else if (shape->IsA() == TGeoTrd1::Class()) {
126  const TGeoTrd1* sh = (const TGeoTrd1*) shape;
127  m_output << "\t\t<tube name=\"" << name << "_shape\" x1=\"" << sh->GetDx1() << "\" x2=\"" << sh->GetDx2() << "\" y1=\""
128  << sh->GetDy() << "\" y2=\"" << sh->GetDy() << "\" z=\"" << sh->GetDz() << "\" lunit=\"cm\"/>" << std::endl;
129  }
130  else if (shape->IsA() == TGeoTrd2::Class()) {
131  const TGeoTrd2* sh = (const TGeoTrd2*) shape;
132  m_output << "\t\t<tube name=\"" << name << "_shape\" x1=\"" << sh->GetDx1() << "\" x2=\"" << sh->GetDx2() << "\" y1=\""
133  << sh->GetDy1() << "\" y2=\"" << sh->GetDy2() << "\" z=\"" << sh->GetDz() << "\" lunit=\"cm\"/>" << std::endl;
134  }
135  else if (shape->IsA() == TGeoPgon::Class()) {
136  const TGeoPgon* sh = (const TGeoPgon*) shape;
137  m_output << "\t\t<polyhedra name=\"" << name << "_shape\" startphi=\"" << sh->GetPhi1() << "\" deltaphi=\"" << sh->GetDphi()
138  << "\" numsides=\"" << sh->GetNedges() << "\" aunit=\"deg\" lunit=\"cm\">" << std::endl;
139  for (int i = 0; i < sh->GetNz(); ++i) {
140  m_output << "\t\t\t<zplane z=\"" << sh->GetZ(i) << "\" rmin=\"" << sh->GetRmin(i) << "\" rmax=\"" << sh->GetRmax(i)
141  << "\" lunit=\"cm\"/>" << std::endl;
142  }
143  m_output << "\t\t</polyhedra>" << std::endl;
144  }
145  else if (shape->IsA() == TGeoPcon::Class()) {
146  const TGeoPcon* sh = (const TGeoPcon*) shape;
147  m_output << "\t\t<polycone name=\"" << name << "_shape\" startphi=\"" << sh->GetPhi1() << "\" deltaphi=\"" << sh->GetDphi()
148  << "\" aunit=\"deg\" lunit=\"cm\">" << std::endl;
149  for (int i = 0; i < sh->GetNz(); ++i) {
150  m_output << "\t\t\t<zplane z=\"" << sh->GetZ(i) << "\" rmin=\"" << sh->GetRmin(i) << "\" rmax=\"" << sh->GetRmax(i)
151  << "\" lunit=\"cm\"/>" << std::endl;
152  }
153  m_output << "\t\t</polycone>" << std::endl;
154  }
155  else if (shape->IsA() == TGeoCompositeShape::Class()) {
156  std::string nn = name;
157  const TGeoCompositeShape* sh = (const TGeoCompositeShape*) shape;
158  const TGeoBoolNode* boolean = sh->GetBoolNode();
159  TGeoBoolNode::EGeoBoolType oper = boolean->GetBooleanOperator();
160 
161  handleSolid(name + "_left", boolean->GetLeftShape());
162  handleSolid(name + "_right", boolean->GetRightShape());
163 
164  if (oper == TGeoBoolNode::kGeoSubtraction)
165  m_output << "\t\t<subtraction name=\"" << nn << "\">" << std::endl;
166  else if (oper == TGeoBoolNode::kGeoUnion)
167  m_output << "\t\t<union name=\"" << nn << "\">" << std::endl;
168  else if (oper == TGeoBoolNode::kGeoIntersection)
169  m_output << "\t\t<intersection name=\"" << nn << "\">" << std::endl;
170 
171  m_output << "\t\t\t<first ref=\"" << nn << "_left\"/>" << std::endl;
172  m_output << "\t\t\t<second ref=\"" << nn << "_right\"/>" << std::endl;
173  indent = "\t";
174  handleTransformation("", boolean->GetRightMatrix());
175  indent = "";
176 
177  if (oper == TGeoBoolNode::kGeoSubtraction)
178  m_output << "\t\t</subtraction>" << std::endl;
179  else if (oper == TGeoBoolNode::kGeoUnion)
180  m_output << "\t\t</union>" << std::endl;
181  else if (oper == TGeoBoolNode::kGeoIntersection)
182  m_output << "\t\t</intersection>" << std::endl;
183  }
184  else {
185  std::cerr << "Failed to handle unknwon solid shape:" << shape->IsA()->GetName() << std::endl;
186  }
187  }
188  return 0;
189 }
190 
192 void detail::GeometryTreeDump::handleStructure(const std::set<Volume>& volset) const {
193  m_output << "\t<structure>" << std::endl;
194  for (const auto& vol : volset)
195  handleVolume(vol->GetName(), vol);
196  m_output << "\t</structure>" << std::endl;
197 }
198 
200 void* detail::GeometryTreeDump::handleTransformation(const std::string& name, const TGeoMatrix* mat) const {
201  if (mat) {
202  if (mat->IsTranslation()) {
203  const Double_t* f = mat->GetTranslation();
204  m_output << indent << "\t\t<position ";
205  if (!name.empty())
206  m_output << "name=\"" << name << "_pos\" ";
207  m_output << "x=\"" << f[0] << "\" " << "y=\"" << f[1] << "\" " << "z=\"" << f[2] << "\" unit=\"cm\"/>" << std::endl;
208  }
209  if (mat->IsRotation()) {
210  const Double_t* matrix = mat->GetRotationMatrix();
211  Double_t theta = 0.0, phi = 0.0, psi = 0.0;
212  getAngles(matrix, theta, phi, psi);
213  m_output << indent << "\t\t<rotation ";
214  if (!name.empty())
215  m_output << "name=\"" << name << "_rot\" ";
216  m_output << "x=\"" << theta << "\" " << "y=\"" << psi << "\" " << "z=\"" << phi << "\" unit=\"deg\"/>" << std::endl;
217  }
218  }
219  return 0;
220 }
221 
223 void detail::GeometryTreeDump::handleTransformations(const std::vector<std::pair<std::string, TGeoMatrix*> >& trafos) const {
224  m_output << "\t<define>" << std::endl;
225  for ( const auto& t : trafos )
226  handleTransformation(t.first, t.second);
227  m_output << "\t</define>" << std::endl;
228 }
229 
231 void detail::GeometryTreeDump::handleSolids(const std::set<TGeoShape*>& solids) const {
232  m_output << "\t<solids>" << std::endl;
233  for ( const auto sh : solids )
234  handleSolid(sh->GetName(), sh);
235  m_output << "\t</solids>" << std::endl;
236 }
237 
239 void detail::GeometryTreeDump::handleDefines(const Detector::HandleMap& defs) const {
240  m_output << "\t<define>" << std::endl;
241  for ( const auto& def : defs )
242  m_output << "\t\t<constant name=\"" << def.second->name << "\" value=\"" << def.second->type << "\" />"
243  << std::endl;
244  m_output << "\t</define>" << std::endl;
245 }
246 
248 void detail::GeometryTreeDump::handleVisualisation(const Detector::HandleMap&) const {
249 }
250 
251 static std::string _path = "";
252 static void dumpStructure(PlacedVolume pv, int level) {
253  TGeoNode* current = pv.ptr();
254  TGeoVolume* volume = current->GetVolume();
255  TObjArray* nodes = volume->GetNodes();
256  int num_children = nodes ? nodes->GetEntriesFast() : 0;
257  char fmt[128];
258 
259  _path += "/";
260  _path += current->GetName();
261  ::snprintf(fmt, sizeof(fmt), "%%4d %%%ds %%7s %%s\n", level * 2 + 5);
262  ::printf(fmt, level, "", " ->LV: ", volume->GetName());
263  ::printf(fmt, level, "", " ->PV: ", current->GetName());
264  ::printf(fmt, level, "", " ->path:", _path.c_str());
265  if (num_children > 0) {
266  for (int i = 0; i < num_children; ++i) {
267  TGeoNode* node = static_cast<TGeoNode*>(nodes->At(i));
268  dumpStructure(PlacedVolume(node), level + 1);
269  }
270  }
271  _path = _path.substr(0, _path.length() - 1 - strlen(current->GetName()));
272 }
273 
274 static void dumpDetectors(DetElement parent, int level) {
275  const DetElement::Children& children = parent.children();
276  PlacedVolume pl = parent.placement();
277  char fmt[128];
278 
279  _path += "/";
280  _path += parent.name();
281  ::snprintf(fmt, sizeof(fmt), "%%4d %%%ds %%7s %%s\n", level * 2 + 5);
282  ::printf(fmt, level, "", "->path:", _path.c_str());
283  if (pl.isValid()) {
284  ::printf(fmt, level, "", " ->placement:", parent.placementPath().c_str());
285  ::printf(fmt, level, "", " ->logvol: ", pl->GetVolume()->GetName());
286  ::printf(fmt, level, "", " ->shape: ", pl->GetVolume()->GetShape()->GetName());
287  }
288  for (const auto& c : children)
289  dumpDetectors(c.second, level + 1);
290 
291  _path = _path.substr(0, _path.length() - 1 - strlen(parent.name()));
292 }
293 
294 void detail::GeometryTreeDump::create(DetElement top) {
295  dumpDetectors(top, 0);
296  dumpStructure(top.placement(), 0);
297 }
dd4hep::DetElement::children
const Children & children() const
Access to the list of children.
Definition: DetElement.cpp:207
dd4hep::detail::GeometryTreeDump::handleSolid
virtual void * handleSolid(const std::string &name, const TGeoShape *volume) const
Dump solid in GDML format to output stream.
Definition: GeometryTreeDump.cpp:107
GeometryTreeDump.h
Detector.h
dd4hep::detail::GeometryTreeDump::handleTransformation
virtual void * handleTransformation(const std::string &name, const TGeoMatrix *matrix) const
Dump single volume transformation in GDML format to output stream.
Definition: GeometryTreeDump.cpp:200
dd4hep::detail::GeometryTreeDump::handleSolids
virtual void handleSolids(const std::set< TGeoShape * > &solids) const
Dump all solids in GDML format to output stream.
Definition: GeometryTreeDump.cpp:231
dd4hep::detail::GeometryTreeDump::handleVolume
virtual void * handleVolume(const std::string &name, Volume volume) const
Dump logical volume in GDML format to output stream.
Definition: GeometryTreeDump.cpp:72
dd4hep::PlacedVolume
Handle class holding a placed volume (also called physical volume)
Definition: Volumes.h:173
dd4hep::VisAttr
Handle class describing visualization attributes.
Definition: Objects.h:324
dd4hep::DetElement::placement
PlacedVolume placement() const
Access to the physical volume of this detector element.
Definition: DetElement.cpp:321
dd4hep::detail::GeometryTreeDump::handleDefines
virtual void handleDefines(const Detector::HandleMap &defs) const
Dump all constants in GDML format to output stream.
Definition: GeometryTreeDump.cpp:239
dd4hep::Handle::isValid
bool isValid() const
Check the validity of the object held by the handle.
Definition: Handle.h:128
dd4hep::Handle::name
const char * name() const
Access the object name (or "" if not supported by the object)
dd4hep::DetElement
Handle class describing a detector element.
Definition: DetElement.h:188
dd4hep::Volume
Handle class holding a placed volume (also called physical volume)
Definition: Volumes.h:378
dd4hep::detail::GeometryTreeDump::handleTransformations
virtual void handleTransformations(const std::vector< std::pair< std::string, TGeoMatrix * > > &trafos) const
Dump Transformations in GDML format to output stream.
Definition: GeometryTreeDump.cpp:223
dd4hep::DetElement::placementPath
const std::string & placementPath() const
Access to the full path to the placed object.
Definition: DetElement.cpp:85
dd4hep::detail::GeometryTreeDump::create
void create(DetElement top)
Main entry point: create required object(s)
Definition: GeometryTreeDump.cpp:294
dd4hep::Detector::HandleMap
std::map< std::string, Handle< NamedObject > > HandleMap
Type definition of a map of named handles.
Definition: Detector.h:93
dd4hep::DetElement::Children
std::map< std::string, DetElement > Children
Definition: DetElement.h:206
dd4hep::Handle::ptr
T * ptr() const
Access to the held object.
Definition: Handle.h:153
dd4hep
Namespace for the AIDA detector description toolkit.
Definition: AlignmentsCalib.h:28
dd4hep::detail::GeometryTreeDump::handleVisualisation
void handleVisualisation(const Detector::HandleMap &vis) const
Dump all visualisation specs in Detector format to output stream.
Definition: GeometryTreeDump.cpp:248
dd4hep::detail::GeometryTreeDump::handleStructure
virtual void handleStructure(const std::set< Volume > &volset) const
Dump structure information in GDML format to output stream.
Definition: GeometryTreeDump.cpp:192
dd4hep::Volume::visAttributes
VisAttr visAttributes() const
Access the visualisation attributes.
Definition: Volumes.cpp:1210
RAD_2_DEGREE
#define RAD_2_DEGREE
Definition: Handle.h:26
Generated on Thu Nov 7 2024 13:31:55 for DD4hep by   doxygen 1.8.18