DD4hep - The AIDA detector description toolkit for high energy physics experiments
DD4hep  Rev:Unversioneddirectory
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
GeometryTreeDump.cpp
Go to the documentation of this file.
1 // $Id$
2 //==========================================================================
3 // AIDA Detector description implementation for LCD
4 //--------------------------------------------------------------------------
5 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
6 // All rights reserved.
7 //
8 // For the licensing terms see $DD4hepINSTALL/LICENSE.
9 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
10 //
11 // Author : M.Frank
12 //
13 //==========================================================================
14 
15 // Framework include files
16 #include "DD4hep/LCDD.h"
17 #include "GeometryTreeDump.h"
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 #include <iostream>
39 
40 using namespace DD4hep::Geometry;
41 using namespace DD4hep;
42 using namespace std;
43 
44 namespace {
45  string indent = "";
46 
48  ostream& m_output = 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* GeometryTreeDump::handleVolume(const 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 << "\">" << endl;
79  m_output << "\t\t\t<solidref ref=\"" << shape->GetName() << "\"/>" << endl;
80  m_output << "\t\t\t<materialref ref=\"" << medium->GetName() << "\"/>" << endl;
81  if (vis.isValid()) {
82  m_output << "\t\t\t<visref ref=\"" << vis.name() << "\"/>" << endl;
83  }
84  if (num > 0) {
85  for (int i = 0; i < num; ++i) {
86  //TGeoNode* n = volume->GetNode(i);
87  TGeoNode* n = vol.ptr()->GetNode(vol->GetNode(i)->GetName());
88  TGeoVolume* v = n->GetVolume();
89  TGeoMatrix* m = n->GetMatrix();
90  m_output << "\t\t\t<physvol>" << endl;
91  m_output << "\t\t\t\t<volumeref ref=\"" << v->GetName() << "\"/>" << endl;
92  if (m) {
93  if (m->IsTranslation()) {
94  m_output << "\t\t\t\t<positionref ref=\"" << n->GetName() << "_pos\"/>" << endl;
95  }
96  if (m->IsRotation()) {
97  m_output << "\t\t\t\t<rotationref ref=\"" << n->GetName() << "_rot\"/>" << endl;
98  }
99  }
100  m_output << "\t\t\t</physvol>" << endl;
101  }
102  }
103  m_output << "\t\t</volume>" << endl;
104  return 0;
105 }
106 
108 void* GeometryTreeDump::handleSolid(const string& name, const TGeoShape* shape) const {
109  if (shape) {
110  if (shape->IsA() == TGeoBBox::Class()) {
111  const TGeoBBox* s = (const TGeoBBox*) shape;
112  m_output << "\t\t<box name=\"" << name << "_shape\" x=\"" << s->GetDX() << "\" y=\"" << s->GetDY() << "\" z=\""
113  << s->GetDZ() << "\" lunit=\"cm\"/>" << endl;
114  }
115  else if (shape->IsA() == TGeoTube::Class()) {
116  const TGeoTube* s = (const TGeoTube*) shape;
117  m_output << "\t\t<tube name=\"" << name << "_shape\" rmin=\"" << s->GetRmin() << "\" rmax=\"" << s->GetRmax() << "\" z=\""
118  << s->GetDz() << "\" startphi=\"0.0\" deltaphi=\"360.0\" aunit=\"deg\" lunit=\"cm\"/>" << endl;
119  }
120  else if (shape->IsA() == TGeoTubeSeg::Class()) {
121  const TGeoTubeSeg* s = (const TGeoTubeSeg*) shape;
122  m_output << "\t\t<tube name=\"" << name << "_shape\" rmin=\"" << s->GetRmin() << "\" rmax=\"" << s->GetRmax() << "\" z=\""
123  << s->GetDz() << "\" startphi=\"" << s->GetPhi1() << "\" deltaphi=\"" << s->GetPhi2()
124  << "\" aunit=\"deg\" lunit=\"cm\"/>" << endl;
125  }
126  else if (shape->IsA() == TGeoTrd1::Class()) {
127  const TGeoTrd1* s = (const TGeoTrd1*) shape;
128  m_output << "\t\t<tube name=\"" << name << "_shape\" x1=\"" << s->GetDx1() << "\" x2=\"" << s->GetDx2() << "\" y1=\""
129  << s->GetDy() << "\" y2=\"" << s->GetDy() << "\" z=\"" << s->GetDz() << "\" lunit=\"cm\"/>" << endl;
130  }
131  else if (shape->IsA() == TGeoTrd2::Class()) {
132  const TGeoTrd2* s = (const TGeoTrd2*) shape;
133  m_output << "\t\t<tube name=\"" << name << "_shape\" x1=\"" << s->GetDx1() << "\" x2=\"" << s->GetDx2() << "\" y1=\""
134  << s->GetDy1() << "\" y2=\"" << s->GetDy2() << "\" z=\"" << s->GetDz() << "\" lunit=\"cm\"/>" << endl;
135  }
136  else if (shape->IsA() == TGeoPgon::Class()) {
137  const TGeoPgon* s = (const TGeoPgon*) shape;
138  m_output << "\t\t<polyhedra name=\"" << name << "_shape\" startphi=\"" << s->GetPhi1() << "\" deltaphi=\"" << s->GetDphi()
139  << "\" numsides=\"" << s->GetNedges() << "\" aunit=\"deg\" lunit=\"cm\">" << endl;
140  for (int i = 0; i < s->GetNz(); ++i) {
141  m_output << "\t\t\t<zplane z=\"" << s->GetZ(i) << "\" rmin=\"" << s->GetRmin(i) << "\" rmax=\"" << s->GetRmax(i)
142  << "\" lunit=\"cm\"/>" << endl;
143  }
144  m_output << "\t\t</polyhedra>" << endl;
145  }
146  else if (shape->IsA() == TGeoPcon::Class()) {
147  const TGeoPcon* s = (const TGeoPcon*) shape;
148  m_output << "\t\t<polycone name=\"" << name << "_shape\" startphi=\"" << s->GetPhi1() << "\" deltaphi=\"" << s->GetDphi()
149  << "\" aunit=\"deg\" lunit=\"cm\">" << endl;
150  for (int i = 0; i < s->GetNz(); ++i) {
151  m_output << "\t\t\t<zplane z=\"" << s->GetZ(i) << "\" rmin=\"" << s->GetRmin(i) << "\" rmax=\"" << s->GetRmax(i)
152  << "\" lunit=\"cm\"/>" << endl;
153  }
154  m_output << "\t\t</polycone>" << endl;
155  }
156  else if (shape->IsA() == TGeoCompositeShape::Class()) {
157  string nn = name;
158  const TGeoCompositeShape* s = (const TGeoCompositeShape*) shape;
159  const TGeoBoolNode* boolean = s->GetBoolNode();
160  TGeoBoolNode::EGeoBoolType oper = boolean->GetBooleanOperator();
161 
162  handleSolid(name + "_left", boolean->GetLeftShape());
163  handleSolid(name + "_right", boolean->GetRightShape());
164 
165  if (oper == TGeoBoolNode::kGeoSubtraction)
166  m_output << "\t\t<subtraction name=\"" << nn << "\">" << endl;
167  else if (oper == TGeoBoolNode::kGeoUnion)
168  m_output << "\t\t<union name=\"" << nn << "\">" << endl;
169  else if (oper == TGeoBoolNode::kGeoIntersection)
170  m_output << "\t\t<intersection name=\"" << nn << "\">" << endl;
171 
172  m_output << "\t\t\t<first ref=\"" << nn << "_left\"/>" << endl;
173  m_output << "\t\t\t<second ref=\"" << nn << "_right\"/>" << endl;
174  indent = "\t";
175  handleTransformation("", boolean->GetRightMatrix());
176  indent = "";
177 
178  if (oper == TGeoBoolNode::kGeoSubtraction)
179  m_output << "\t\t</subtraction>" << endl;
180  else if (oper == TGeoBoolNode::kGeoUnion)
181  m_output << "\t\t</union>" << endl;
182  else if (oper == TGeoBoolNode::kGeoIntersection)
183  m_output << "\t\t</intersection>" << endl;
184  }
185  else {
186  cerr << "Failed to handle unknwon solid shape:" << shape->IsA()->GetName() << endl;
187  }
188  }
189  return 0;
190 }
191 
193 void GeometryTreeDump::handleStructure(const VolumeSet& volset) const {
194  m_output << "\t<structure>" << endl;
195  for (VolumeSet::const_iterator i = volset.begin(); i != volset.end(); ++i)
196  handleVolume((*i)->GetName(), *i);
197  m_output << "\t</structure>" << endl;
198 }
199 
201 void* GeometryTreeDump::handleTransformation(const string& name, const TGeoMatrix* m) const {
202  if (m) {
203  if (m->IsTranslation()) {
204  const Double_t* f = m->GetTranslation();
205  m_output << indent << "\t\t<position ";
206  if (!name.empty())
207  m_output << "name=\"" << name << "_pos\" ";
208  m_output << "x=\"" << f[0] << "\" " << "y=\"" << f[1] << "\" " << "z=\"" << f[2] << "\" unit=\"cm\"/>" << endl;
209  }
210  if (m->IsRotation()) {
211  const Double_t* mat = m->GetRotationMatrix();
212  Double_t theta = 0.0, phi = 0.0, psi = 0.0;
213  getAngles(mat, theta, phi, psi);
214  m_output << indent << "\t\t<rotation ";
215  if (!name.empty())
216  m_output << "name=\"" << name << "_rot\" ";
217  m_output << "x=\"" << theta << "\" " << "y=\"" << psi << "\" " << "z=\"" << phi << "\" unit=\"deg\"/>" << endl;
218  }
219  }
220  return 0;
221 }
222 
225  m_output << "\t<define>" << endl;
226  for (TransformSet::const_iterator i = trafos.begin(); i != trafos.end(); ++i)
227  handleTransformation((*i).first, (*i).second);
228  m_output << "\t</define>" << endl;
229 }
230 
232 void GeometryTreeDump::handleSolids(const SolidSet& solids) const {
233  m_output << "\t<solids>" << endl;
234  for (SolidSet::const_iterator i = solids.begin(); i != solids.end(); ++i)
235  handleSolid((*i)->GetName(), *i);
236  m_output << "\t</solids>" << endl;
237 }
238 
241  m_output << "\t<define>" << endl;
242  for (LCDD::HandleMap::const_iterator i = defs.begin(); i != defs.end(); ++i)
243  m_output << "\t\t<constant name=\"" << (*i).second->name << "\" value=\"" << (*i).second->type << "\" />"
244  << endl;
245  m_output << "\t</define>" << endl;
246 }
247 
250 }
251 
252 static string _path = "";
253 static void dumpStructure(PlacedVolume pv, int level) {
254  TGeoNode* current = pv.ptr();
255  TGeoVolume* volume = current->GetVolume();
256  TObjArray* nodes = volume->GetNodes();
257  int num_children = nodes ? nodes->GetEntriesFast() : 0;
258  char fmt[128];
259 
260  _path += "/";
261  _path += current->GetName();
262  ::snprintf(fmt, sizeof(fmt), "%%4d %%%ds %%7s %%s\n", level * 2 + 5);
263  ::printf(fmt, level, "", " ->LV: ", volume->GetName());
264  ::printf(fmt, level, "", " ->PV: ", current->GetName());
265  ::printf(fmt, level, "", " ->path:", _path.c_str());
266  if (num_children > 0) {
267  for (int i = 0; i < num_children; ++i) {
268  TGeoNode* node = static_cast<TGeoNode*>(nodes->At(i));
269  dumpStructure(PlacedVolume(node), level + 1);
270  }
271  }
272  _path = _path.substr(0, _path.length() - 1 - strlen(current->GetName()));
273 }
274 
275 static void dumpDetectors(DetElement parent, int level) {
276  const DetElement::Children& children = parent.children();
277  PlacedVolume pl = parent.placement();
278  char fmt[128];
279 
280  _path += "/";
281  _path += parent.name();
282  ::snprintf(fmt, sizeof(fmt), "%%4d %%%ds %%7s %%s\n", level * 2 + 5);
283  ::printf(fmt, level, "", "->path:", _path.c_str());
284  if (pl.isValid()) {
285  ::printf(fmt, level, "", " ->placement:", parent.placementPath().c_str());
286  ::printf(fmt, level, "", " ->logvol: ", pl->GetVolume()->GetName());
287  ::printf(fmt, level, "", " ->shape: ", pl->GetVolume()->GetShape()->GetName());
288  }
289  for (DetElement::Children::const_iterator i = children.begin(); i != children.end(); ++i) {
290  dumpDetectors((*i).second, level + 1);
291  }
292  _path = _path.substr(0, _path.length() - 1 - strlen(parent.name()));
293 }
294 
296  //PlacedVolume pv = top.placement();
297  dumpDetectors(top, 0);
298  dumpStructure(top.placement(), 0);
299  //GeometryInfo geo;
300  //collect(top,geo);
301  //handleSetup(LCDD::getInstance().header());
302  //handleDefines(LCDD::getInstance().constants());
303  //handleVisualisation(geo.vis);
304  //handleTransformations(geo.trafos);
305  //handleSolids(geo.solids);
306  //handleStructure(geo.volumes);
307 }
Handle class holding a placed volume (also called physical volume)
Definition: Volumes.h:135
std::set< Volume > VolumeSet
Definition: GeoHandler.h:52
const char * name() const
Access the object name (or "" if not supported by the object)
Definition: Handle.inl:36
bool isValid() const
Check the validity of the object held by the handle.
Definition: Handle.h:124
virtual void handleSolids(const SolidSet &solids) const
Dump all solids in GDML format to output stream.
static void dumpDetectors(DetElement parent, int level)
virtual void handleDefines(const LCDD::HandleMap &defs) const
Dump all constants in GDML format to output stream.
virtual void * handleTransformation(const std::string &name, const TGeoMatrix *matrix) const
Dump single volume transformation in GDML format to output stream.
std::set< TGeoShape * > SolidSet
Definition: GeoHandler.h:56
const Children & children() const
Access to the list of children.
Definition: Detector.cpp:207
TGeoShape * s
Definition: Volumes.cpp:294
return e
Definition: Volumes.cpp:297
VisAttr visAttributes() const
Access the visualisation attributes.
Definition: Volumes.cpp:679
static string _path
T * ptr() const
Access to the held object.
Definition: Handle.h:149
std::vector< std::pair< std::string, TGeoMatrix * > > TransformSet
Definition: GeoHandler.h:55
Handle class holding a placed volume (also called physical volume)
Definition: Volumes.h:237
virtual void * handleSolid(const std::string &name, const TGeoShape *volume) const
Dump solid in GDML format to output stream.
virtual void handleTransformations(const TransformSet &trafos) const
Dump Transformations in GDML format to output stream.
virtual void * handleVolume(const std::string &name, Volume volume) const
Dump logical volume in GDML format to output stream.
void create(DetElement top)
Main entry point: create required object(s)
PlacedVolume placement() const
Access to the physical volume of this detector element.
Definition: Detector.cpp:279
Handle class describing visualization attributes.
Definition: Objects.h:341
Handle class describing a detector element.
Definition: Detector.h:172
View * v
Definition: MultiView.cpp:30
void handleVisualisation(const LCDD::HandleMap &vis) const
Dump all visualisation specs in LCDD format to output stream.
virtual void handleStructure(const VolumeSet &volset) const
Dump structure information in GDML format to output stream.
std::map< std::string, DetElement > Children
Definition: Detector.h:202
std::map< std::string, NamedHandle > HandleMap
Definition: LCDD.h:87
static void dumpStructure(PlacedVolume pv, int level)
TGeoShape TGeoMedium * m
Definition: Volumes.cpp:294
#define RAD_2_DEGREE
Definition: Handle.h:26
const std::string & placementPath() const
Access to the full path to the placed object.
Definition: Detector.cpp:86