G4GDMLReadStructure Class Reference

#include <G4GDMLReadStructure.hh>

Inheritance diagram for G4GDMLReadStructure:

G4GDMLReadParamvol G4GDMLReadSetup G4GDMLReadSolids G4GDMLReadMaterials G4GDMLReadDefine G4GDMLRead

Public Member Functions

 G4GDMLReadStructure ()
virtual ~G4GDMLReadStructure ()
G4VPhysicalVolumeGetPhysvol (const G4String &) const
G4LogicalVolumeGetVolume (const G4String &) const
G4AssemblyVolumeGetAssembly (const G4String &) const
G4GDMLAuxListType GetVolumeAuxiliaryInformation (G4LogicalVolume *) const
G4VPhysicalVolumeGetWorldVolume (const G4String &)
const G4GDMLAuxMapTypeGetAuxMap () const
void Clear ()
virtual void VolumeRead (const xercesc::DOMElement *const)
virtual void Volume_contentRead (const xercesc::DOMElement *const)
virtual void StructureRead (const xercesc::DOMElement *const)

Protected Member Functions

G4GDMLAuxPairType AuxiliaryRead (const xercesc::DOMElement *const)
void AssemblyRead (const xercesc::DOMElement *const)
void DivisionvolRead (const xercesc::DOMElement *const)
G4LogicalVolumeFileRead (const xercesc::DOMElement *const)
void PhysvolRead (const xercesc::DOMElement *const, G4AssemblyVolume *assembly=0)
void ReplicavolRead (const xercesc::DOMElement *const, G4int number)
void ReplicaRead (const xercesc::DOMElement *const replicaElement, G4LogicalVolume *logvol, G4int number)
EAxis AxisRead (const xercesc::DOMElement *const axisElement)
G4double QuantityRead (const xercesc::DOMElement *const readElement)
void BorderSurfaceRead (const xercesc::DOMElement *const)
void SkinSurfaceRead (const xercesc::DOMElement *const)

Protected Attributes

G4GDMLAuxMapType auxMap
G4GDMLAssemblyMapType assemblyMap
G4LogicalVolumepMotherLogical
std::map< std::string, G4VPhysicalVolume * > setuptoPV

Detailed Description

Definition at line 62 of file G4GDMLReadStructure.hh.


Constructor & Destructor Documentation

G4GDMLReadStructure::G4GDMLReadStructure (  ) 

Definition at line 48 of file G4GDMLReadStructure.cc.

00049   : G4GDMLReadParamvol(), pMotherLogical(0)
00050 {
00051 }

G4GDMLReadStructure::~G4GDMLReadStructure (  )  [virtual]

Definition at line 53 of file G4GDMLReadStructure.cc.

00054 {
00055 }


Member Function Documentation

void G4GDMLReadStructure::AssemblyRead ( const xercesc::DOMElement *  const  )  [protected]

Definition at line 635 of file G4GDMLReadStructure.cc.

References assemblyMap, FatalException, G4cout, G4endl, G4Exception(), G4GDMLRead::GenerateName(), PhysvolRead(), and G4GDMLRead::Transcode().

Referenced by StructureRead().

00636 {
00637    XMLCh *name_attr = xercesc::XMLString::transcode("name");
00638    const G4String name = Transcode(assemblyElement->getAttribute(name_attr));
00639    xercesc::XMLString::release(&name_attr);
00640 
00641    G4AssemblyVolume* pAssembly = new G4AssemblyVolume();
00642    assemblyMap.insert(std::make_pair(GenerateName(name), pAssembly));
00643 
00644    for (xercesc::DOMNode* iter = assemblyElement->getFirstChild();
00645         iter != 0; iter = iter->getNextSibling())
00646    {
00647       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) { continue; }
00648       const xercesc::DOMElement* const child
00649             = dynamic_cast<xercesc::DOMElement*>(iter);
00650       if (!child)
00651       {
00652         G4Exception("G4GDMLReadStructure::AssemblyRead()",
00653                     "InvalidRead", FatalException, "No child found!");
00654         return;
00655       }
00656       const G4String tag = Transcode(child->getTagName());
00657 
00658       if (tag=="physvol")
00659       {
00660         PhysvolRead(child, pAssembly);
00661       }
00662       else
00663       {
00664         G4cout << "Unsupported GDML tag '" << tag
00665                << "' for Geant4 assembly structure !" << G4endl;
00666       }
00667    }
00668 }

G4GDMLAuxPairType G4GDMLReadStructure::AuxiliaryRead ( const xercesc::DOMElement *  const  )  [protected]

Definition at line 58 of file G4GDMLReadStructure.cc.

References FatalException, G4Exception(), G4GDMLRead::Transcode(), G4GDMLAuxPairType::type, and G4GDMLAuxPairType::value.

Referenced by VolumeRead().

00059 {
00060    G4GDMLAuxPairType auxpair = {"",""};
00061 
00062    const xercesc::DOMNamedNodeMap* const attributes
00063          = auxiliaryElement->getAttributes();
00064    XMLSize_t attributeCount = attributes->getLength();
00065 
00066    for (XMLSize_t attribute_index=0;
00067         attribute_index<attributeCount; attribute_index++)
00068    {
00069       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00070 
00071       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00072         { continue; }
00073 
00074       const xercesc::DOMAttr* const attribute
00075             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00076       if (!attribute)
00077       {
00078         G4Exception("G4GDMLReadStructure::AuxiliaryRead()",
00079                     "InvalidRead", FatalException, "No attribute found!");
00080         return auxpair;
00081       }
00082       const G4String attName = Transcode(attribute->getName());
00083       const G4String attValue = Transcode(attribute->getValue());
00084 
00085       if (attName=="auxtype") { auxpair.type = attValue; } else
00086         //      if (attName=="auxvalue") { auxpair.value = eval.Evaluate(attValue); }
00087       if (attName=="auxvalue") { auxpair.value = attValue; }
00088    }
00089 
00090    return auxpair;
00091 }

EAxis G4GDMLReadStructure::AxisRead ( const xercesc::DOMElement *const   axisElement  )  [protected]

Definition at line 517 of file G4GDMLReadStructure.cc.

References G4GDMLRead::eval, G4GDMLEvaluator::Evaluate(), FatalException, G4Exception(), kPhi, kRho, kUndefined, kXAxis, kYAxis, kZAxis, and G4GDMLRead::Transcode().

Referenced by ReplicaRead().

00518 {
00519    EAxis axis = kUndefined;
00520 
00521    const xercesc::DOMNamedNodeMap* const attributes
00522          = axisElement->getAttributes();
00523    XMLSize_t attributeCount = attributes->getLength();
00524 
00525    for (XMLSize_t attribute_index=0;
00526         attribute_index<attributeCount; attribute_index++)
00527    {
00528       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00529 
00530       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00531         { continue; }
00532 
00533       const xercesc::DOMAttr* const attribute
00534             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00535       if (!attribute)
00536       {
00537         G4Exception("G4GDMLReadStructure::AxisRead()",
00538                     "InvalidRead", FatalException, "No attribute found!");
00539         return axis;
00540       }
00541       const G4String attName = Transcode(attribute->getName());
00542       const G4String attValue = Transcode(attribute->getValue());
00543       if (attName=="x")
00544         { if( eval.Evaluate(attValue)==1.) {axis=kXAxis;} }
00545       else if (attName=="y")
00546         { if( eval.Evaluate(attValue)==1.) {axis=kYAxis;} }
00547       else if (attName=="z")
00548         { if( eval.Evaluate(attValue)==1.) {axis=kZAxis;} }
00549       else if (attName=="rho")
00550         { if( eval.Evaluate(attValue)==1.) {axis=kRho;}   }
00551       else if (attName=="phi")
00552         { if( eval.Evaluate(attValue)==1.) {axis=kPhi;}   }
00553    }
00554 
00555    return axis;
00556 }

void G4GDMLReadStructure::BorderSurfaceRead ( const xercesc::DOMElement *  const  )  [protected]

Definition at line 94 of file G4GDMLReadStructure.cc.

References FatalException, G4Exception(), G4GDMLRead::GenerateName(), GetPhysvol(), G4GDMLReadSolids::GetSurfaceProperty(), G4GDMLReadDefine::RefRead(), G4GDMLRead::Strip(), and G4GDMLRead::Transcode().

Referenced by StructureRead().

00095 {
00096    G4String name;
00097    G4VPhysicalVolume* pv1 = 0;
00098    G4VPhysicalVolume* pv2 = 0;
00099    G4SurfaceProperty* prop = 0;
00100    G4int index = 0;
00101 
00102    const xercesc::DOMNamedNodeMap* const attributes
00103          = bordersurfaceElement->getAttributes();
00104    XMLSize_t attributeCount = attributes->getLength();
00105 
00106    for (XMLSize_t attribute_index=0;
00107         attribute_index<attributeCount; attribute_index++)
00108    {
00109       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00110 
00111       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00112         { continue; }
00113 
00114       const xercesc::DOMAttr* const attribute
00115             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00116       if (!attribute)
00117       {
00118         G4Exception("G4GDMLReadStructure::BorderSurfaceRead()",
00119                     "InvalidRead", FatalException, "No attribute found!");
00120         return;
00121       }
00122       const G4String attName = Transcode(attribute->getName());
00123       const G4String attValue = Transcode(attribute->getValue());
00124 
00125       if (attName=="name")
00126         { name = GenerateName(attValue); } else
00127       if (attName=="surfaceproperty")
00128         { prop = GetSurfaceProperty(GenerateName(attValue)); }
00129    }
00130 
00131    for (xercesc::DOMNode* iter = bordersurfaceElement->getFirstChild();
00132         iter != 0; iter = iter->getNextSibling())
00133    {
00134       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00135 
00136       const xercesc::DOMElement* const child
00137             = dynamic_cast<xercesc::DOMElement*>(iter);
00138       if (!child)
00139       {
00140         G4Exception("G4GDMLReadStructure::BorderSurfaceRead()",
00141                     "InvalidRead", FatalException, "No child found!");
00142         return;
00143       }
00144       const G4String tag = Transcode(child->getTagName());
00145 
00146       if (tag != "physvolref")  { continue; }
00147  
00148       if (index==0)
00149         { pv1 = GetPhysvol(GenerateName(RefRead(child))); index++; } else
00150       if (index==1)
00151         { pv2 = GetPhysvol(GenerateName(RefRead(child))); index++; } else
00152       break;
00153    }
00154 
00155    new G4LogicalBorderSurface(Strip(name),pv1,pv2,prop);
00156 }

void G4GDMLReadStructure::Clear (  ) 

Definition at line 923 of file G4GDMLReadStructure.cc.

References G4GDMLEvaluator::Clear(), G4GDMLRead::eval, and setuptoPV.

Referenced by G4GDMLParser::Clear().

00924 {
00925   eval.Clear();
00926   setuptoPV.clear();
00927 }

void G4GDMLReadStructure::DivisionvolRead ( const xercesc::DOMElement *  const  )  [protected]

Definition at line 159 of file G4GDMLReadStructure.cc.

References G4ReflectionFactory::Divide(), G4GDMLRead::eval, G4GDMLEvaluator::Evaluate(), G4GDMLEvaluator::EvaluateInteger(), FatalException, G4String::first(), G4Exception(), G4GDMLRead::GenerateName(), G4GDMLRead::GeneratePhysvolName(), G4PVDivisionFactory::GetInstance(), G4LogicalVolume::GetName(), GetVolume(), G4ReflectionFactory::Instance(), kPhi, kRho, kUndefined, kXAxis, kYAxis, kZAxis, pMotherLogical, G4GDMLReadDefine::RefRead(), and G4GDMLRead::Transcode().

Referenced by Volume_contentRead().

00160 {
00161    G4String name;
00162    G4double unit = 1.0;
00163    G4double width = 0.0;
00164    G4double offset = 0.0;
00165    G4int number = 0;
00166    EAxis axis = kUndefined;
00167    G4LogicalVolume* logvol = 0;
00168    
00169    const xercesc::DOMNamedNodeMap* const attributes
00170          = divisionvolElement->getAttributes();
00171    XMLSize_t attributeCount = attributes->getLength();
00172 
00173    for (XMLSize_t attribute_index=0;
00174         attribute_index<attributeCount; attribute_index++)
00175    {
00176       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00177 
00178       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00179         { continue; }
00180 
00181       const xercesc::DOMAttr* const attribute
00182             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00183       if (!attribute)
00184       {
00185         G4Exception("G4GDMLReadStructure::DivisionvolRead()",
00186                     "InvalidRead", FatalException, "No attribute found!");
00187         return;
00188       }
00189       const G4String attName = Transcode(attribute->getName());
00190       const G4String attValue = Transcode(attribute->getValue());
00191 
00192       if (attName=="name") { name = attValue; } else
00193       if (attName=="unit") { unit = eval.Evaluate(attValue); } else
00194       if (attName=="width") { width = eval.Evaluate(attValue); } else
00195       if (attName=="offset") { offset = eval.Evaluate(attValue); } else
00196       if (attName=="number") { number = eval.EvaluateInteger(attValue); } else
00197       if (attName=="axis")
00198       {
00199          if (attValue=="kXAxis") { axis = kXAxis; } else
00200          if (attValue=="kYAxis") { axis = kYAxis; } else
00201          if (attValue=="kZAxis") { axis = kZAxis; } else
00202          if (attValue=="kRho") { axis = kRho; } else
00203          if (attValue=="kPhi") { axis = kPhi; }
00204       }
00205    }
00206 
00207    width *= unit;
00208    offset *= unit;
00209 
00210    for (xercesc::DOMNode* iter = divisionvolElement->getFirstChild();
00211         iter != 0;iter = iter->getNextSibling())
00212    {
00213       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) { continue; }
00214 
00215       const xercesc::DOMElement* const child
00216             = dynamic_cast<xercesc::DOMElement*>(iter);
00217       if (!child)
00218       {
00219         G4Exception("G4GDMLReadStructure::DivisionvolRead()",
00220                     "InvalidRead", FatalException, "No child found!");
00221         return;
00222       }
00223       const G4String tag = Transcode(child->getTagName());
00224 
00225       if (tag=="volumeref") { logvol = GetVolume(GenerateName(RefRead(child))); }
00226    }
00227 
00228    if (!logvol)  { return; }
00229 
00230    G4PVDivisionFactory::GetInstance();
00231    G4PhysicalVolumesPair pair;
00232 
00233    G4String pv_name = logvol->GetName() + "_div";
00234    if ((number != 0) && (width == 0.0))
00235    {
00236      pair = G4ReflectionFactory::Instance()
00237             ->Divide(pv_name,logvol,pMotherLogical,axis,number,offset);
00238    }
00239    else if ((number == 0) && (width != 0.0))
00240    {
00241      pair = G4ReflectionFactory::Instance()
00242             ->Divide(pv_name,logvol,pMotherLogical,axis,width,offset);
00243    }
00244    else
00245    {
00246      pair = G4ReflectionFactory::Instance()
00247             ->Divide(pv_name,logvol,pMotherLogical,axis,number,width,offset);
00248    }
00249 
00250    if (pair.first != 0) { GeneratePhysvolName(name,pair.first); }
00251    if (pair.second != 0) { GeneratePhysvolName(name,pair.second); }
00252 }

G4LogicalVolume * G4GDMLReadStructure::FileRead ( const xercesc::DOMElement *  const  )  [protected]

Definition at line 255 of file G4GDMLReadStructure.cc.

References auxMap, FatalException, G4Exception(), G4GDMLRead::GenerateName(), GetAuxMap(), G4GDMLReadSetup::GetSetup(), GetVolume(), G4GDMLRead::Read(), G4GDMLRead::Transcode(), and G4GDMLRead::validate.

Referenced by PhysvolRead().

00256 {
00257    G4String name;
00258    G4String volname;
00259 
00260    const xercesc::DOMNamedNodeMap* const attributes
00261          = fileElement->getAttributes();
00262    XMLSize_t attributeCount = attributes->getLength();
00263 
00264    for (XMLSize_t attribute_index=0;
00265         attribute_index<attributeCount; attribute_index++)
00266    {
00267       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00268 
00269       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00270         { continue; }
00271 
00272       const xercesc::DOMAttr* const attribute
00273             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00274       if (!attribute)
00275       {
00276         G4Exception("G4GDMLReadStructure::FileRead()",
00277                     "InvalidRead", FatalException, "No attribute found!");
00278         return 0;
00279       }
00280       const G4String attName = Transcode(attribute->getName());
00281       const G4String attValue = Transcode(attribute->getValue());
00282 
00283       if (attName=="name") { name = attValue; } else
00284       if (attName=="volname") { volname = attValue; }
00285    }
00286 
00287    const G4bool isModule = true;
00288    G4GDMLReadStructure structure;
00289    structure.Read(name,validate,isModule);
00290 
00291    // Register existing auxiliar information defined in child module
00292    //
00293    const G4GDMLAuxMapType* aux = structure.GetAuxMap();
00294    if (!aux->empty())
00295    {
00296      G4GDMLAuxMapType::const_iterator pos;
00297      for (pos = aux->begin(); pos != aux->end(); ++pos)
00298      {
00299        auxMap.insert(std::make_pair(pos->first,pos->second));
00300      }
00301    }
00302 
00303    // Return volume structure from child module
00304    //
00305    if (volname.empty())
00306    {
00307      return structure.GetVolume(structure.GetSetup("Default"));
00308    }
00309    else
00310    {
00311      return structure.GetVolume(structure.GenerateName(volname));
00312    }
00313 }

G4AssemblyVolume * G4GDMLReadStructure::GetAssembly ( const G4String  )  const

Definition at line 881 of file G4GDMLReadStructure.cc.

References assemblyMap.

Referenced by PhysvolRead().

00882 {
00883    G4GDMLAssemblyMapType::const_iterator pos = assemblyMap.find(ref);
00884    if (pos != assemblyMap.end()) { return pos->second; }
00885    return 0;
00886 }

const G4GDMLAuxMapType * G4GDMLReadStructure::GetAuxMap (  )  const

Definition at line 897 of file G4GDMLReadStructure.cc.

References auxMap.

Referenced by FileRead(), and G4GDMLParser::GetAuxMap().

00898 {
00899    return &auxMap;
00900 }

G4VPhysicalVolume * G4GDMLReadStructure::GetPhysvol ( const G4String  )  const

Definition at line 849 of file G4GDMLReadStructure.cc.

References FatalException, G4Exception(), G4PhysicalVolumeStore::GetInstance(), and G4PhysicalVolumeStore::GetVolume().

Referenced by BorderSurfaceRead().

00850 {
00851    G4VPhysicalVolume* physvolPtr =
00852      G4PhysicalVolumeStore::GetInstance()->GetVolume(ref,false);
00853 
00854    if (!physvolPtr)
00855    {
00856      G4String error_msg = "Referenced physvol '" + ref + "' was not found!";
00857      G4Exception("G4GDMLReadStructure::GetPhysvol()", "ReadError",
00858                  FatalException, error_msg);
00859    }
00860 
00861    return physvolPtr;
00862 }

G4LogicalVolume * G4GDMLReadStructure::GetVolume ( const G4String  )  const [virtual]

Implements G4GDMLRead.

Definition at line 865 of file G4GDMLReadStructure.cc.

References FatalException, G4Exception(), G4LogicalVolumeStore::GetInstance(), and G4LogicalVolumeStore::GetVolume().

Referenced by DivisionvolRead(), FileRead(), G4GDMLParser::GetVolume(), GetWorldVolume(), PhysvolRead(), ReplicavolRead(), and SkinSurfaceRead().

00866 {
00867    G4LogicalVolume *volumePtr
00868    = G4LogicalVolumeStore::GetInstance()->GetVolume(ref,false);
00869 
00870    if (!volumePtr)
00871    {
00872      G4String error_msg = "Referenced volume '" + ref + "' was not found!";
00873      G4Exception("G4GDMLReadStructure::GetVolume()", "ReadError",
00874                  FatalException, error_msg);
00875    }
00876 
00877    return volumePtr;
00878 }

G4GDMLAuxListType G4GDMLReadStructure::GetVolumeAuxiliaryInformation ( G4LogicalVolume  )  const

Definition at line 889 of file G4GDMLReadStructure.cc.

References auxMap.

Referenced by G4GDMLParser::GetVolumeAuxiliaryInformation().

00890 {
00891    G4GDMLAuxMapType::const_iterator pos = auxMap.find(logvol);
00892    if (pos != auxMap.end()) { return pos->second; }
00893    else { return G4GDMLAuxListType(); }
00894 }

G4VPhysicalVolume * G4GDMLReadStructure::GetWorldVolume ( const G4String  ) 

Definition at line 903 of file G4GDMLReadStructure.cc.

References G4LogicalVolume::GetName(), G4GDMLReadSetup::GetSetup(), GetVolume(), G4VisAttributes::Invisible, setuptoPV, G4LogicalVolume::SetVisAttributes(), and G4GDMLRead::Strip().

Referenced by G4GDMLParser::GetWorldVolume().

00904 {    
00905    G4LogicalVolume* volume = GetVolume(Strip(GetSetup(setupName)));
00906    volume->SetVisAttributes(G4VisAttributes::Invisible);
00907 
00908    G4VPhysicalVolume* pvWorld = 0;
00909 
00910    if(setuptoPV[setupName])
00911    {
00912      pvWorld = setuptoPV[setupName];
00913    }
00914    else
00915    {
00916      pvWorld = new G4PVPlacement(0,G4ThreeVector(0,0,0),volume,
00917                                  volume->GetName()+"_PV",0,0,0);
00918      setuptoPV[setupName] = pvWorld;
00919    }
00920    return pvWorld;
00921 }

void G4GDMLReadStructure::PhysvolRead ( const xercesc::DOMElement *  const,
G4AssemblyVolume assembly = 0 
) [protected]

Definition at line 316 of file G4GDMLReadStructure.cc.

References G4AssemblyVolume::AddPlacedVolume(), G4GDMLRead::check, FatalException, FileRead(), G4String::first(), G4Exception(), G4GDMLRead::GenerateName(), G4GDMLRead::GeneratePhysvolName(), GetAssembly(), G4LogicalVolume::GetName(), G4GDMLReadDefine::GetPosition(), G4GDMLReadDefine::GetRotation(), G4GDMLReadDefine::GetRotationMatrix(), G4GDMLReadDefine::GetScale(), GetVolume(), G4ReflectionFactory::Instance(), G4AssemblyVolume::MakeImprint(), G4ReflectionFactory::Place(), pMotherLogical, position, G4GDMLReadDefine::RefRead(), G4GDMLRead::Transcode(), and G4GDMLReadDefine::VectorRead().

Referenced by AssemblyRead(), and Volume_contentRead().

00318 {
00319    G4String name;
00320    G4LogicalVolume* logvol = 0;
00321    G4AssemblyVolume* assembly = 0;
00322    G4ThreeVector position(0.0,0.0,0.0);
00323    G4ThreeVector rotation(0.0,0.0,0.0);
00324    G4ThreeVector scale(1.0,1.0,1.0);
00325 
00326    const xercesc::DOMNamedNodeMap* const attributes
00327          = physvolElement->getAttributes();
00328    XMLSize_t attributeCount = attributes->getLength();
00329 
00330    for (XMLSize_t attribute_index=0;
00331         attribute_index<attributeCount; attribute_index++)
00332    {
00333      xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00334 
00335      if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00336        { continue; }
00337 
00338      const xercesc::DOMAttr* const attribute
00339            = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00340      if (!attribute)
00341      {
00342        G4Exception("G4GDMLReadStructure::PhysvolRead()",
00343                    "InvalidRead", FatalException, "No attribute found!");
00344        return;
00345      }
00346      const G4String attName = Transcode(attribute->getName());
00347      const G4String attValue = Transcode(attribute->getValue());
00348 
00349      if (attName=="name") { name = attValue; }
00350    }
00351 
00352    for (xercesc::DOMNode* iter = physvolElement->getFirstChild();
00353         iter != 0; iter = iter->getNextSibling())
00354    {
00355      if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00356 
00357      const xercesc::DOMElement* const child
00358            = dynamic_cast<xercesc::DOMElement*>(iter);
00359      if (!child)
00360      {
00361        G4Exception("G4GDMLReadStructure::PhysvolRead()",
00362                    "InvalidRead", FatalException, "No child found!");
00363        return;
00364      }
00365      const G4String tag = Transcode(child->getTagName());
00366 
00367      if (tag=="volumeref")
00368        {
00369          const G4String& child_name = GenerateName(RefRead(child));
00370          assembly = GetAssembly(child_name);
00371          if (!assembly) { logvol = GetVolume(child_name); }
00372        }
00373      else if (tag=="file")
00374        { logvol = FileRead(child); }
00375      else if (tag=="position")
00376        { VectorRead(child,position); }
00377      else if (tag=="rotation")
00378        { VectorRead(child,rotation); }
00379      else if (tag=="scale")
00380        { VectorRead(child,scale); }
00381      else if (tag=="positionref")
00382        { position = GetPosition(GenerateName(RefRead(child))); }
00383      else if (tag=="rotationref")
00384         { rotation = GetRotation(GenerateName(RefRead(child))); }
00385      else if (tag=="scaleref")
00386         { scale = GetScale(GenerateName(RefRead(child))); }
00387      else
00388         {
00389           G4String error_msg = "Unknown tag in physvol: " + tag;
00390           G4Exception("G4GDMLReadStructure::PhysvolRead()", "ReadError",
00391                       FatalException, error_msg);
00392           return;
00393         }
00394    }
00395 
00396    G4Transform3D transform(GetRotationMatrix(rotation).inverse(),position);
00397    transform = transform*G4Scale3D(scale.x(),scale.y(),scale.z());
00398 
00399    if (pAssembly)   // Fill assembly structure
00400    {
00401      if (!logvol) { return; }
00402      pAssembly->AddPlacedVolume(logvol, transform);
00403    }
00404    else             // Generate physical volume tree or do assembly imprint
00405    {
00406      if (assembly)
00407      {
00408        assembly->MakeImprint(pMotherLogical, transform, 0, check);
00409      }
00410      else
00411      {
00412        if (!logvol) { return; }
00413        G4String pv_name = logvol->GetName() + "_PV";
00414        G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()
00415          ->Place(transform,pv_name,logvol,pMotherLogical,false,0,check);
00416 
00417        if (pair.first != 0) { GeneratePhysvolName(name,pair.first); }
00418        if (pair.second != 0) { GeneratePhysvolName(name,pair.second); }
00419      }
00420    }
00421 }

G4double G4GDMLReadStructure::QuantityRead ( const xercesc::DOMElement *const   readElement  )  [protected]

Reimplemented from G4GDMLReadDefine.

Definition at line 559 of file G4GDMLReadStructure.cc.

References G4GDMLRead::eval, G4GDMLEvaluator::Evaluate(), FatalException, G4Exception(), and G4GDMLRead::Transcode().

Referenced by ReplicaRead().

00560 {
00561    G4double value = 0.0;
00562    G4double unit = 0.0;
00563    const xercesc::DOMNamedNodeMap* const attributes
00564          = readElement->getAttributes();
00565    XMLSize_t attributeCount = attributes->getLength();
00566 
00567    for (XMLSize_t attribute_index=0;
00568         attribute_index<attributeCount; attribute_index++)
00569    {
00570       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00571 
00572       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00573         { continue; }
00574       const xercesc::DOMAttr* const attribute
00575             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00576       if (!attribute)
00577       {
00578         G4Exception("G4GDMLReadStructure::QuantityRead()",
00579                     "InvalidRead", FatalException, "No attribute found!");
00580         return value;
00581       }
00582       const G4String attName = Transcode(attribute->getName());
00583       const G4String attValue = Transcode(attribute->getValue());
00584 
00585       if (attName=="unit") { unit = eval.Evaluate(attValue); } else
00586       if (attName=="value"){ value= eval.Evaluate(attValue); } 
00587    }
00588 
00589    return value*unit;
00590 }

void G4GDMLReadStructure::ReplicaRead ( const xercesc::DOMElement *const   replicaElement,
G4LogicalVolume logvol,
G4int  number 
) [protected]

Definition at line 460 of file G4GDMLReadStructure.cc.

References AxisRead(), FatalException, G4String::first(), G4Exception(), G4GDMLRead::GenerateName(), G4GDMLRead::GeneratePhysvolName(), G4LogicalVolume::GetName(), G4GDMLReadDefine::GetPosition(), G4GDMLReadDefine::GetRotation(), G4ReflectionFactory::Instance(), kUndefined, pMotherLogical, position, QuantityRead(), G4GDMLReadDefine::RefRead(), G4ReflectionFactory::Replicate(), G4GDMLRead::Transcode(), and G4GDMLReadDefine::VectorRead().

Referenced by ReplicavolRead().

00462 {
00463    G4double width = 0.0;
00464    G4double offset = 0.0;
00465    G4ThreeVector position(0.0,0.0,0.0);
00466    G4ThreeVector rotation(0.0,0.0,0.0);
00467    EAxis axis = kUndefined;
00468    G4String name;
00469  
00470    for (xercesc::DOMNode* iter = replicaElement->getFirstChild();
00471                           iter != 0; iter = iter->getNextSibling())
00472    {
00473       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00474 
00475       const xercesc::DOMElement* const child
00476             = dynamic_cast<xercesc::DOMElement*>(iter);
00477       if (!child)
00478       {
00479         G4Exception("G4GDMLReadStructure::ReplicaRead()",
00480                     "InvalidRead", FatalException, "No child found!");
00481         return;
00482       }
00483       const G4String tag = Transcode(child->getTagName()); 
00484 
00485       if (tag=="position")
00486         { VectorRead(child,position); } else
00487       if (tag=="rotation")
00488         { VectorRead(child,rotation); } else
00489       if (tag=="positionref")
00490         { position = GetPosition(GenerateName(RefRead(child))); } else
00491       if (tag=="rotationref")
00492         { rotation = GetRotation(GenerateName(RefRead(child))); } else
00493       if (tag=="direction")
00494         { axis=AxisRead(child); } else
00495       if (tag=="width")
00496         { width=QuantityRead(child); } else
00497       if (tag=="offset")
00498         { offset=QuantityRead(child); }
00499       else
00500       {
00501         G4String error_msg = "Unknown tag in ReplicaRead: " + tag;
00502         G4Exception("G4GDMLReadStructure::ReplicaRead()", "ReadError",
00503                     FatalException, error_msg);
00504       }
00505    }
00506 
00507    G4String pv_name = logvol->GetName() + "_PV";
00508    G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()
00509      ->Replicate(pv_name,logvol,pMotherLogical,axis,number,width,offset);
00510 
00511    if (pair.first != 0) { GeneratePhysvolName(name,pair.first); }
00512    if (pair.second != 0) { GeneratePhysvolName(name,pair.second); }
00513 
00514 }

void G4GDMLReadStructure::ReplicavolRead ( const xercesc::DOMElement *  const,
G4int  number 
) [protected]

Definition at line 424 of file G4GDMLReadStructure.cc.

References FatalException, G4Exception(), G4GDMLRead::GenerateName(), GetVolume(), G4GDMLReadDefine::RefRead(), ReplicaRead(), and G4GDMLRead::Transcode().

Referenced by Volume_contentRead().

00425 {
00426   G4LogicalVolume* logvol = 0;
00427   for (xercesc::DOMNode* iter = replicavolElement->getFirstChild();
00428                          iter != 0; iter = iter->getNextSibling())
00429   {
00430     if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00431 
00432     const xercesc::DOMElement* const child
00433           = dynamic_cast<xercesc::DOMElement*>(iter);
00434     if (!child)
00435     {
00436       G4Exception("G4GDMLReadStructure::ReplicavolRead()",
00437                   "InvalidRead", FatalException, "No child found!");
00438       return;
00439     }
00440     const G4String tag = Transcode(child->getTagName());
00441 
00442     if (tag=="volumeref")
00443     {
00444       logvol = GetVolume(GenerateName(RefRead(child)));
00445     }
00446     else if (tag=="replicate_along_axis")
00447     {
00448       if (logvol)  { ReplicaRead(child,logvol,number); }
00449     }
00450     else
00451     {
00452       G4String error_msg = "Unknown tag in ReplicavolRead: " + tag;
00453       G4Exception("G4GDMLReadStructure::ReplicavolRead()",
00454                   "ReadError", FatalException, error_msg);
00455     }
00456   }
00457 }

void G4GDMLReadStructure::SkinSurfaceRead ( const xercesc::DOMElement *  const  )  [protected]

Definition at line 671 of file G4GDMLReadStructure.cc.

References FatalException, G4Exception(), G4GDMLRead::GenerateName(), G4GDMLReadSolids::GetSurfaceProperty(), GetVolume(), G4GDMLReadDefine::RefRead(), G4GDMLRead::Strip(), and G4GDMLRead::Transcode().

Referenced by StructureRead().

00672 {
00673    G4String name;
00674    G4LogicalVolume* logvol = 0;
00675    G4SurfaceProperty* prop = 0;
00676 
00677    const xercesc::DOMNamedNodeMap* const attributes
00678          = skinsurfaceElement->getAttributes();
00679    XMLSize_t attributeCount = attributes->getLength();
00680 
00681    for (XMLSize_t attribute_index=0;
00682         attribute_index<attributeCount; attribute_index++)
00683    {
00684       xercesc::DOMNode* attribute_node = attributes->item(attribute_index);
00685 
00686       if (attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE)
00687         { continue; }
00688 
00689       const xercesc::DOMAttr* const attribute
00690             = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00691       if (!attribute)
00692       {
00693         G4Exception("G4GDMLReadStructure::SkinsurfaceRead()",
00694                     "InvalidRead", FatalException, "No attribute found!");
00695         return;
00696       }
00697       const G4String attName = Transcode(attribute->getName());
00698       const G4String attValue = Transcode(attribute->getValue());
00699 
00700       if (attName=="name")
00701         { name = GenerateName(attValue); } else
00702       if (attName=="surfaceproperty")
00703         { prop = GetSurfaceProperty(GenerateName(attValue)); }
00704    }
00705 
00706    for (xercesc::DOMNode* iter = skinsurfaceElement->getFirstChild();
00707         iter != 0; iter = iter->getNextSibling())
00708    {
00709       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00710 
00711       const xercesc::DOMElement* const child
00712             = dynamic_cast<xercesc::DOMElement*>(iter);
00713       if (!child)
00714       {
00715         G4Exception("G4GDMLReadStructure::SkinsurfaceRead()",
00716                     "InvalidRead", FatalException, "No child found!");
00717         return;
00718       }
00719       const G4String tag = Transcode(child->getTagName());
00720 
00721       if (tag=="volumeref")
00722       {
00723         logvol = GetVolume(GenerateName(RefRead(child)));
00724       }
00725       else
00726       {
00727         G4String error_msg = "Unknown tag in skinsurface: " + tag;
00728         G4Exception("G4GDMLReadStructure::SkinsurfaceRead()", "ReadError",
00729                     FatalException, error_msg);
00730       }
00731    }
00732 
00733    new G4LogicalSkinSurface(Strip(name),logvol,prop);
00734 }

void G4GDMLReadStructure::StructureRead ( const xercesc::DOMElement *  const  )  [virtual]

Implements G4GDMLRead.

Definition at line 815 of file G4GDMLReadStructure.cc.

References AssemblyRead(), BorderSurfaceRead(), FatalException, G4cout, G4endl, G4Exception(), G4GDMLRead::LoopRead(), SkinSurfaceRead(), G4GDMLRead::StructureRead(), G4GDMLRead::Transcode(), and VolumeRead().

00816 {
00817    G4cout << "G4GDML: Reading structure..." << G4endl;
00818 
00819    for (xercesc::DOMNode* iter = structureElement->getFirstChild();
00820         iter != 0; iter = iter->getNextSibling())
00821    {
00822       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00823 
00824       const xercesc::DOMElement* const child
00825             = dynamic_cast<xercesc::DOMElement*>(iter);
00826       if (!child)
00827       {
00828         G4Exception("G4GDMLReadStructure::StructureRead()",
00829                     "InvalidRead", FatalException, "No child found!");
00830         return;
00831       }
00832       const G4String tag = Transcode(child->getTagName());
00833 
00834       if (tag=="bordersurface") { BorderSurfaceRead(child); } else
00835       if (tag=="skinsurface") { SkinSurfaceRead(child); } else
00836       if (tag=="volume") { VolumeRead(child); } else
00837       if (tag=="assembly") { AssemblyRead(child); } else
00838       if (tag=="loop") { LoopRead(child,&G4GDMLRead::StructureRead); }
00839       else
00840       {
00841         G4String error_msg = "Unknown tag in structure: " + tag;
00842         G4Exception("G4GDMLReadStructure::StructureRead()",
00843                     "ReadError", FatalException, error_msg);
00844       }
00845    }
00846 }

void G4GDMLReadStructure::Volume_contentRead ( const xercesc::DOMElement *  const  )  [virtual]

Implements G4GDMLRead.

Definition at line 737 of file G4GDMLReadStructure.cc.

References DivisionvolRead(), G4GDMLRead::eval, G4GDMLEvaluator::EvaluateInteger(), FatalException, G4cout, G4endl, G4Exception(), G4GDMLRead::LoopRead(), G4GDMLReadParamvol::ParamvolRead(), PhysvolRead(), pMotherLogical, ReplicavolRead(), G4GDMLRead::Transcode(), and G4GDMLRead::Volume_contentRead().

Referenced by VolumeRead().

00738 {
00739    for (xercesc::DOMNode* iter = volumeElement->getFirstChild();
00740         iter != 0; iter = iter->getNextSibling())
00741    {
00742       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) { continue; }
00743 
00744       const xercesc::DOMElement* const child
00745             = dynamic_cast<xercesc::DOMElement*>(iter);
00746       if (!child)
00747       {
00748         G4Exception("G4GDMLReadStructure::Volume_contentRead()",
00749                     "InvalidRead", FatalException, "No child found!");
00750         return;
00751       }
00752       const G4String tag = Transcode(child->getTagName());
00753 
00754       if ((tag=="auxiliary") || (tag=="materialref") || (tag=="solidref"))
00755       {
00756         // These are already processed in VolumeRead()
00757       }
00758       else if (tag=="paramvol")
00759       {
00760         ParamvolRead(child,pMotherLogical);
00761       }
00762       else if (tag=="physvol")
00763       {
00764         PhysvolRead(child);
00765       }
00766       else if (tag=="replicavol")
00767       {
00768         G4int number = 1;
00769         const xercesc::DOMNamedNodeMap* const attributes
00770               = child->getAttributes();
00771         XMLSize_t attributeCount = attributes->getLength();
00772         for (XMLSize_t attribute_index=0;
00773                        attribute_index<attributeCount; attribute_index++)
00774         {
00775           xercesc::DOMNode* attribute_node
00776                  = attributes->item(attribute_index);
00777           if (attribute_node->getNodeType()!=xercesc::DOMNode::ATTRIBUTE_NODE)
00778           {
00779             continue;
00780           }
00781           const xercesc::DOMAttr* const attribute
00782                 = dynamic_cast<xercesc::DOMAttr*>(attribute_node);   
00783           if (!attribute)
00784           {
00785             G4Exception("G4GDMLReadStructure::Volume_contentRead()",
00786                         "InvalidRead", FatalException, "No attribute found!");
00787             return;
00788           }
00789           const G4String attName = Transcode(attribute->getName());
00790           const G4String attValue = Transcode(attribute->getValue());
00791           if (attName=="number")
00792           {
00793             number = eval.EvaluateInteger(attValue);
00794           }
00795         }
00796         ReplicavolRead(child,number); 
00797       }
00798       else if (tag=="divisionvol")
00799       {
00800         DivisionvolRead(child);
00801       }
00802       else if (tag=="loop")
00803       {
00804         LoopRead(child,&G4GDMLRead::Volume_contentRead);
00805       }
00806       else
00807       {
00808         G4cout << "Treating unknown GDML tag in volume '" << tag
00809                << "' as GDML extension..." << G4endl;
00810       }
00811    }
00812 }

void G4GDMLReadStructure::VolumeRead ( const xercesc::DOMElement *  const  )  [virtual]

Definition at line 593 of file G4GDMLReadStructure.cc.

References AuxiliaryRead(), auxMap, FatalException, G4Exception(), G4GDMLRead::GenerateName(), G4GDMLReadMaterials::GetMaterial(), G4GDMLReadSolids::GetSolid(), pMotherLogical, G4GDMLReadDefine::RefRead(), G4GDMLRead::Transcode(), and Volume_contentRead().

Referenced by StructureRead().

00594 {
00595    G4VSolid* solidPtr = 0;
00596    G4Material* materialPtr = 0;
00597    G4GDMLAuxListType auxList;
00598    
00599    XMLCh *name_attr = xercesc::XMLString::transcode("name");
00600    const G4String name = Transcode(volumeElement->getAttribute(name_attr));
00601    xercesc::XMLString::release(&name_attr);
00602 
00603    for (xercesc::DOMNode* iter = volumeElement->getFirstChild();
00604         iter != 0; iter = iter->getNextSibling())
00605    {
00606       if (iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE)  { continue; }
00607 
00608       const xercesc::DOMElement* const child
00609             = dynamic_cast<xercesc::DOMElement*>(iter);
00610       if (!child)
00611       {
00612         G4Exception("G4GDMLReadStructure::VolumeRead()",
00613                     "InvalidRead", FatalException, "No child found!");
00614         return;
00615       }
00616       const G4String tag = Transcode(child->getTagName());
00617 
00618       if (tag=="auxiliary")
00619         { auxList.push_back(AuxiliaryRead(child)); } else
00620       if (tag=="materialref")
00621         { materialPtr = GetMaterial(GenerateName(RefRead(child),true)); } else
00622       if (tag=="solidref")
00623         { solidPtr = GetSolid(GenerateName(RefRead(child))); }
00624    }
00625 
00626    pMotherLogical = new G4LogicalVolume(solidPtr,materialPtr,
00627                                         GenerateName(name),0,0,0);
00628 
00629    if (!auxList.empty()) { auxMap[pMotherLogical] = auxList; }
00630 
00631    Volume_contentRead(volumeElement);
00632 }


Field Documentation

G4GDMLAssemblyMapType G4GDMLReadStructure::assemblyMap [protected]

Definition at line 101 of file G4GDMLReadStructure.hh.

Referenced by AssemblyRead(), and GetAssembly().

G4GDMLAuxMapType G4GDMLReadStructure::auxMap [protected]

Definition at line 100 of file G4GDMLReadStructure.hh.

Referenced by FileRead(), GetAuxMap(), GetVolumeAuxiliaryInformation(), and VolumeRead().

G4LogicalVolume* G4GDMLReadStructure::pMotherLogical [protected]

Definition at line 102 of file G4GDMLReadStructure.hh.

Referenced by DivisionvolRead(), PhysvolRead(), ReplicaRead(), Volume_contentRead(), and VolumeRead().

std::map<std::string, G4VPhysicalVolume*> G4GDMLReadStructure::setuptoPV [protected]

Definition at line 103 of file G4GDMLReadStructure.hh.

Referenced by Clear(), and GetWorldVolume().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:03 2013 for Geant4 by  doxygen 1.4.7