00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035 #include "G4GDMLWriteDefine.hh"
00036 #include "G4SystemOfUnits.hh"
00037
00038 const G4double G4GDMLWriteDefine::kRelativePrecision = DBL_EPSILON;
00039 const G4double G4GDMLWriteDefine::kAngularPrecision = DBL_EPSILON;
00040 const G4double G4GDMLWriteDefine::kLinearPrecision = DBL_EPSILON;
00041
00042 G4GDMLWriteDefine::G4GDMLWriteDefine()
00043 : G4GDMLWrite(), defineElement(0)
00044 {
00045 }
00046
00047 G4GDMLWriteDefine::~G4GDMLWriteDefine()
00048 {
00049 }
00050
00051 G4ThreeVector G4GDMLWriteDefine::GetAngles(const G4RotationMatrix& mat)
00052 {
00053 G4double x,y,z;
00054
00055 const G4double cosb = std::sqrt(mat.xx()*mat.xx()+mat.yx()*mat.yx());
00056
00057 if (cosb > kRelativePrecision)
00058 {
00059 x = std::atan2(mat.zy(),mat.zz());
00060 y = std::atan2(-mat.zx(),cosb);
00061 z = std::atan2(mat.yx(),mat.xx());
00062 }
00063 else
00064 {
00065 x = std::atan2(-mat.yz(),mat.yy());
00066 y = std::atan2(-mat.zx(),cosb);
00067 z = 0.0;
00068 }
00069
00070 return G4ThreeVector(x,y,z);
00071 }
00072
00073 void G4GDMLWriteDefine::
00074 Scale_vectorWrite(xercesc::DOMElement* element, const G4String& tag,
00075 const G4String& name, const G4ThreeVector& scl)
00076 {
00077 const G4double x = (std::fabs(scl.x()-1.0) < kRelativePrecision)
00078 ? 1.0 : scl.x();
00079 const G4double y = (std::fabs(scl.y()-1.0) < kRelativePrecision)
00080 ? 1.0 : scl.y();
00081 const G4double z = (std::fabs(scl.z()-1.0) < kRelativePrecision)
00082 ? 1.0 : scl.z();
00083
00084 xercesc::DOMElement* scaleElement = NewElement(tag);
00085 scaleElement->setAttributeNode(NewAttribute("name",name));
00086 scaleElement->setAttributeNode(NewAttribute("x",x));
00087 scaleElement->setAttributeNode(NewAttribute("y",y));
00088 scaleElement->setAttributeNode(NewAttribute("z",z));
00089 element->appendChild(scaleElement);
00090 }
00091
00092 void G4GDMLWriteDefine::
00093 Rotation_vectorWrite(xercesc::DOMElement* element, const G4String& tag,
00094 const G4String& name, const G4ThreeVector& rot)
00095 {
00096 const G4double x = (std::fabs(rot.x()) < kAngularPrecision) ? 0.0 : rot.x();
00097 const G4double y = (std::fabs(rot.y()) < kAngularPrecision) ? 0.0 : rot.y();
00098 const G4double z = (std::fabs(rot.z()) < kAngularPrecision) ? 0.0 : rot.z();
00099
00100 xercesc::DOMElement* rotationElement = NewElement(tag);
00101 rotationElement->setAttributeNode(NewAttribute("name",name));
00102 rotationElement->setAttributeNode(NewAttribute("x",x/degree));
00103 rotationElement->setAttributeNode(NewAttribute("y",y/degree));
00104 rotationElement->setAttributeNode(NewAttribute("z",z/degree));
00105 rotationElement->setAttributeNode(NewAttribute("unit","deg"));
00106 element->appendChild(rotationElement);
00107 }
00108
00109 void G4GDMLWriteDefine::
00110 Position_vectorWrite(xercesc::DOMElement* element, const G4String& tag,
00111 const G4String& name, const G4ThreeVector& pos)
00112 {
00113 const G4double x = (std::fabs(pos.x()) < kLinearPrecision) ? 0.0 : pos.x();
00114 const G4double y = (std::fabs(pos.y()) < kLinearPrecision) ? 0.0 : pos.y();
00115 const G4double z = (std::fabs(pos.z()) < kLinearPrecision) ? 0.0 : pos.z();
00116
00117 xercesc::DOMElement* positionElement = NewElement(tag);
00118 positionElement->setAttributeNode(NewAttribute("name",name));
00119 positionElement->setAttributeNode(NewAttribute("x",x/mm));
00120 positionElement->setAttributeNode(NewAttribute("y",y/mm));
00121 positionElement->setAttributeNode(NewAttribute("z",z/mm));
00122 positionElement->setAttributeNode(NewAttribute("unit","mm"));
00123 element->appendChild(positionElement);
00124 }
00125
00126 void G4GDMLWriteDefine::DefineWrite(xercesc::DOMElement* element)
00127 {
00128 G4cout << "G4GDML: Writing definitions..." << G4endl;
00129
00130 defineElement = NewElement("define");
00131 element->appendChild(defineElement);
00132 }