Geant4_10
G4GDMLWriteDefine.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4GDMLWriteDefine.cc 68053 2013-03-13 14:39:51Z gcosmo $
28 //
29 // class G4GDMLWriteDefine Implementation
30 //
31 // Original author: Zoltan Torzsok, November 2007
32 //
33 // --------------------------------------------------------------------
34 
35 #include "G4GDMLWriteDefine.hh"
36 #include "G4SystemOfUnits.hh"
37 
41 
43  : G4GDMLWrite(), defineElement(0)
44 {
45 }
46 
48 {
49 }
50 
52 {
53  G4double x,y,z;
54 
55  const G4double cosb = std::sqrt(mat.xx()*mat.xx()+mat.yx()*mat.yx());
56 
57  if (cosb > kRelativePrecision)
58  {
59  x = std::atan2(mat.zy(),mat.zz());
60  y = std::atan2(-mat.zx(),cosb);
61  z = std::atan2(mat.yx(),mat.xx());
62  }
63  else
64  {
65  x = std::atan2(-mat.yz(),mat.yy());
66  y = std::atan2(-mat.zx(),cosb);
67  z = 0.0;
68  }
69 
70  return G4ThreeVector(x,y,z);
71 }
72 
74 Scale_vectorWrite(xercesc::DOMElement* element, const G4String& tag,
75  const G4String& name, const G4ThreeVector& scl)
76 {
77  const G4double x = (std::fabs(scl.x()-1.0) < kRelativePrecision)
78  ? 1.0 : scl.x();
79  const G4double y = (std::fabs(scl.y()-1.0) < kRelativePrecision)
80  ? 1.0 : scl.y();
81  const G4double z = (std::fabs(scl.z()-1.0) < kRelativePrecision)
82  ? 1.0 : scl.z();
83 
84  xercesc::DOMElement* scaleElement = NewElement(tag);
85  scaleElement->setAttributeNode(NewAttribute("name",name));
86  scaleElement->setAttributeNode(NewAttribute("x",x));
87  scaleElement->setAttributeNode(NewAttribute("y",y));
88  scaleElement->setAttributeNode(NewAttribute("z",z));
89  element->appendChild(scaleElement);
90 }
91 
93 Rotation_vectorWrite(xercesc::DOMElement* element, const G4String& tag,
94  const G4String& name, const G4ThreeVector& rot)
95 {
96  const G4double x = (std::fabs(rot.x()) < kAngularPrecision) ? 0.0 : rot.x();
97  const G4double y = (std::fabs(rot.y()) < kAngularPrecision) ? 0.0 : rot.y();
98  const G4double z = (std::fabs(rot.z()) < kAngularPrecision) ? 0.0 : rot.z();
99 
100  xercesc::DOMElement* rotationElement = NewElement(tag);
101  rotationElement->setAttributeNode(NewAttribute("name",name));
102  rotationElement->setAttributeNode(NewAttribute("x",x/degree));
103  rotationElement->setAttributeNode(NewAttribute("y",y/degree));
104  rotationElement->setAttributeNode(NewAttribute("z",z/degree));
105  rotationElement->setAttributeNode(NewAttribute("unit","deg"));
106  element->appendChild(rotationElement);
107 }
108 
110 Position_vectorWrite(xercesc::DOMElement* element, const G4String& tag,
111  const G4String& name, const G4ThreeVector& pos)
112 {
113  const G4double x = (std::fabs(pos.x()) < kLinearPrecision) ? 0.0 : pos.x();
114  const G4double y = (std::fabs(pos.y()) < kLinearPrecision) ? 0.0 : pos.y();
115  const G4double z = (std::fabs(pos.z()) < kLinearPrecision) ? 0.0 : pos.z();
116 
117  xercesc::DOMElement* positionElement = NewElement(tag);
118  positionElement->setAttributeNode(NewAttribute("name",name));
119  positionElement->setAttributeNode(NewAttribute("x",x/mm));
120  positionElement->setAttributeNode(NewAttribute("y",y/mm));
121  positionElement->setAttributeNode(NewAttribute("z",z/mm));
122  positionElement->setAttributeNode(NewAttribute("unit","mm"));
123  element->appendChild(positionElement);
124 }
125 
126 void G4GDMLWriteDefine::DefineWrite(xercesc::DOMElement* element)
127 {
128  G4cout << "G4GDML: Writing definitions..." << G4endl;
129 
130  defineElement = NewElement("define");
131  element->appendChild(defineElement);
132 }
static const G4double kAngularPrecision
void Scale_vectorWrite(xercesc::DOMElement *, const G4String &, const G4String &, const G4ThreeVector &)
double xx() const
CLHEP::Hep3Vector G4ThreeVector
double x() const
Definition: xmlparse.cc:179
double yy() const
const XML_Char * name
Definition: expat.h:151
xercesc::DOMElement * defineElement
tuple x
Definition: test.py:50
xercesc::DOMElement * NewElement(const G4String &)
Definition: G4GDMLWrite.cc:127
double zx() const
double z() const
Float_t mat
Definition: plot.C:40
Double_t y
Definition: plot.C:279
void Position_vectorWrite(xercesc::DOMElement *, const G4String &, const G4String &, const G4ThreeVector &)
double yz() const
virtual void DefineWrite(xercesc::DOMElement *)
G4GLOB_DLL std::ostream G4cout
tuple degree
Definition: hepunit.py:69
static const G4double kRelativePrecision
#define DBL_EPSILON
Definition: templates.hh:87
xercesc::DOMAttr * NewAttribute(const G4String &, const G4String &)
Definition: G4GDMLWrite.cc:103
static const G4double kLinearPrecision
double zy() const
G4ThreeVector GetAngles(const G4RotationMatrix &)
double y() const
void Rotation_vectorWrite(xercesc::DOMElement *, const G4String &, const G4String &, const G4ThreeVector &)
tuple z
Definition: test.py:28
#define G4endl
Definition: G4ios.hh:61
double yx() const
double G4double
Definition: G4Types.hh:76
double zz() const
virtual ~G4GDMLWriteDefine()