From e58aacd68769460de56f0e48b58259f162f58a56 Mon Sep 17 00:00:00 2001 From: Prakash Date: Wed, 14 Feb 2024 16:35:49 -0500 Subject: [PATCH 1/4] Add a very ugly implementation of zyx rotation --- include/remollDetectorConstruction.hh | 3 +- include/remollGDMLReadStructure.hh | 519 ++++++++++++++++++++++++++ remoll.cc | 4 +- src/remollDetectorConstruction.cc | 4 +- 4 files changed, 527 insertions(+), 3 deletions(-) create mode 100644 include/remollGDMLReadStructure.hh diff --git a/include/remollDetectorConstruction.hh b/include/remollDetectorConstruction.hh index 047a2dcba..c39aea00e 100644 --- a/include/remollDetectorConstruction.hh +++ b/include/remollDetectorConstruction.hh @@ -7,6 +7,7 @@ #include "G4GenericMessenger.hh" #include "G4Types.hh" #include "G4Version.hh" +#include "remollGDMLReadStructure.hh" #include #include @@ -25,7 +26,7 @@ class remollDetectorConstruction : public G4VUserDetectorConstruction { public: - remollDetectorConstruction(const G4String& name, const G4String& gdmlfile); + remollDetectorConstruction(const G4String& name, const G4String& gdmlfile, remollGDMLReadStructure* rs=nullptr); virtual ~remollDetectorConstruction(); private: diff --git a/include/remollGDMLReadStructure.hh b/include/remollGDMLReadStructure.hh new file mode 100644 index 000000000..bbbfe4269 --- /dev/null +++ b/include/remollGDMLReadStructure.hh @@ -0,0 +1,519 @@ +#pragma once + +#include "G4LogicalVolume.hh" +#include "G4VPhysicalVolume.hh" +#include "G4PVPlacement.hh" +#include "G4LogicalVolumeStore.hh" +#include "G4PhysicalVolumeStore.hh" +#include "G4AssemblyVolume.hh" +#include "G4ReflectionFactory.hh" +#include "G4PVDivisionFactory.hh" +#include "G4LogicalBorderSurface.hh" +#include "G4LogicalSkinSurface.hh" +#include "G4GDMLReadStructure.hh" + + +class remollGDMLReadStructure : public G4GDMLReadStructure { + public: + remollGDMLReadStructure():G4GDMLReadStructure(){ } + + G4LogicalVolume* FileReads(const xercesc::DOMElement* const fileElement) + { + G4String name; + G4String volname; + + const xercesc::DOMNamedNodeMap* const attributes = + fileElement->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; + ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + + const xercesc::DOMAttr* const attribute = + dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLReadStructure::FileRead()", "InvalidRead", + FatalException, "No attribute found!"); + return nullptr; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + + if(attName == "name") + { + name = attValue; + } + else if(attName == "volname") + { + volname = attValue; + } + } + + const G4bool isModule = true; + remollGDMLReadStructure structure; + structure.Read(name, validate, isModule); + + // Register existing auxiliar information defined in child module + // + const G4GDMLAuxMapType* aux = structure.GetAuxMap(); + if(!aux->empty()) + { + for(auto pos = aux->cbegin(); pos != aux->cend(); ++pos) + { + auxMap.insert(std::make_pair(pos->first, pos->second)); + } + } + + // Return volume structure from child module + // + if(volname.empty()) + { + return structure.GetVolume(structure.GetSetup("Default")); + //return this->GetVolume(this->GetSetup("Default")); + } + else + { + return structure.GetVolume(structure.GenerateName(volname)); + //return this->GetVolume(this->GenerateName(volname)); + } + } + + + void PhysvolReads( const xercesc::DOMElement* const physvolElement, G4AssemblyVolume* pAssembly=0) + { + G4String name; + G4LogicalVolume* logvol = nullptr; + G4AssemblyVolume* assembly = nullptr; + G4ThreeVector position(0.0, 0.0, 0.0); + G4ThreeVector rotation(0.0, 0.0, 0.0); + G4ThreeVector rotationzyx(0.0, 0.0, 0.0); + bool zyx_requested = false; + G4ThreeVector scale(1.0, 1.0, 1.0); + G4int copynumber = 0; + + const xercesc::DOMNamedNodeMap* const attributes = + physvolElement->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; + ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + + const xercesc::DOMAttr* const attribute = + dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", + FatalException, "No attribute found!"); + return; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + + if(attName == "name") + { + name = attValue; + } + if(attName == "copynumber") + { + copynumber = eval.EvaluateInteger(attValue); + } + } + + for(xercesc::DOMNode* iter = physvolElement->getFirstChild(); iter != nullptr; + iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = + dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", + FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + + if(tag == "volumeref") + { + const G4String& child_name = GenerateName(RefRead(child)); + assembly = GetAssembly(child_name); + if(assembly == nullptr) + { + logvol = GetVolume(child_name); + } + } + else if(tag == "file") + { + logvol = FileReads(child); + } + else if(tag == "position") + { + VectorRead(child, position); + } + else if(tag == "rotation") + { + VectorRead(child, rotation); + } + else if(tag == "rotationzyx") + { + VectorRead(child, rotationzyx); + zyx_requested = true; + } + else if(tag == "scale") + { + VectorRead(child, scale); + } + else if(tag == "positionref") + { + position = GetPosition(GenerateName(RefRead(child))); + } + else if(tag == "rotationref") + { + rotation = GetRotation(GenerateName(RefRead(child))); + } + else if(tag == "scaleref") + { + scale = GetScale(GenerateName(RefRead(child))); + } + else + { + G4String error_msg = "Unknown tag in physvol: " + tag; + G4Exception("G4GDMLReadStructure::PhysvolRead()", "ReadError", + FatalException, error_msg); + return; + } + } + + G4Transform3D transform; + if(!zyx_requested){ + transform = G4Transform3D(GetRotationMatrix(rotation).inverse(), position); + } else { + G4RotationMatrix rot; + + rot.rotateY(rotationzyx.y()); + rot.rotateZ(rotationzyx.z()); + rot.rotateX(rotationzyx.x()); + rot.rectify(); // Rectify matrix from possible roundoff errors + + transform = G4Transform3D(rot.inverse(),position); + } + + transform = transform * G4Scale3D(scale.x(), scale.y(), scale.z()); + + if(pAssembly != nullptr) // Fill assembly structure + { + if(assembly != nullptr) // Case of recursive assemblies + { + pAssembly->AddPlacedAssembly(assembly, transform); + } + if(logvol == nullptr) + { + return; + } + pAssembly->AddPlacedVolume(logvol, transform); + } + else // Generate physical volume tree or do assembly imprint + { + if(assembly != nullptr) + { + assembly->MakeImprint(pMotherLogical, transform, 0, check); + } + else + { + if(logvol == nullptr) + { + return; + } + G4String pv_name = logvol->GetName() + "_PV"; + G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()->Place( + transform, pv_name, logvol, pMotherLogical, false, copynumber, check); + + if(pair.first != nullptr) + { + GeneratePhysvolName(name, pair.first); + } + if(pair.second != nullptr) + { + GeneratePhysvolName(name, pair.second); + } + } + } + } + + void ParametersRead(const xercesc::DOMElement* const element) { + + G4ThreeVector rotation(0.0, 0.0, 0.0); + G4ThreeVector rotationzyx(0.0, 0.0, 0.0); + G4ThreeVector position(0.0, 0.0, 0.0); + bool zyx_requested = false; + + G4GDMLParameterisation::PARAMETER parameter; + parameter.pRot = new G4RotationMatrix(); + + for(xercesc::DOMNode* iter = element->getFirstChild(); iter != nullptr; + iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = + dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadParamvol::ParametersRead()", "InvalidRead", + FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + if(tag == "rotation") + { + VectorRead(child, rotation); + } + else if(tag == "rotationzyx") + { + VectorRead(child, rotationzyx); + zyx_requested = true; + } + else if(tag == "position") + { + VectorRead(child, position); + } + else if(tag == "positionref") + { + position = GetPosition(GenerateName(RefRead(child))); + } + else if(tag == "rotationref") + { + rotation = GetRotation(GenerateName(RefRead(child))); + } + else if(tag == "box_dimensions") + { + Box_dimensionsRead(child, parameter); + } + else if(tag == "trd_dimensions") + { + Trd_dimensionsRead(child, parameter); + } + else if(tag == "trap_dimensions") + { + Trap_dimensionsRead(child, parameter); + } + else if(tag == "tube_dimensions") + { + Tube_dimensionsRead(child, parameter); + } + else if(tag == "cone_dimensions") + { + Cone_dimensionsRead(child, parameter); + } + else if(tag == "sphere_dimensions") + { + Sphere_dimensionsRead(child, parameter); + } + else if(tag == "orb_dimensions") + { + Orb_dimensionsRead(child, parameter); + } + else if(tag == "torus_dimensions") + { + Torus_dimensionsRead(child, parameter); + } + else if(tag == "ellipsoid_dimensions") + { + Ellipsoid_dimensionsRead(child, parameter); + } + else if(tag == "para_dimensions") + { + Para_dimensionsRead(child, parameter); + } + else if(tag == "polycone_dimensions") + { + Polycone_dimensionsRead(child, parameter); + } + else if(tag == "polyhedra_dimensions") + { + Polyhedra_dimensionsRead(child, parameter); + } + else if(tag == "hype_dimensions") + { + Hype_dimensionsRead(child, parameter); + } + else + { + G4String error_msg = "Unknown tag in parameters: " + tag; + G4Exception("G4GDMLReadParamvol::ParametersRead()", "ReadError", + FatalException, error_msg); + } + } + + + + if(!zyx_requested){ + parameter.pRot->rotateX(rotation.x()); + parameter.pRot->rotateY(rotation.y()); + parameter.pRot->rotateZ(rotation.z()); + } else { + parameter.pRot->rotateY(rotationzyx.y()); + parameter.pRot->rotateZ(rotationzyx.z()); + parameter.pRot->rotateX(rotationzyx.x()); + } + + + parameter.position = position; + + parameterisation->AddParameter(parameter); + + } + + void Volume_contentRead( const xercesc::DOMElement* const volumeElement) override + { + for(xercesc::DOMNode* iter = volumeElement->getFirstChild(); iter != nullptr; + iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = + dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadStructure::Volume_contentRead()", "InvalidRead", + FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + + if((tag == "auxiliary") || (tag == "materialref") || (tag == "solidref")) + { + // These are already processed in VolumeRead() + } + else if(tag == "paramvol") + { + ParamvolRead(child, pMotherLogical); + } + else if(tag == "physvol") + { + PhysvolReads(child); + } + else if(tag == "replicavol") + { + G4int number = 1; + const xercesc::DOMNamedNodeMap* const attributes = child->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; + ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + const xercesc::DOMAttr* const attribute = + dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLReadStructure::Volume_contentRead()", + "InvalidRead", FatalException, "No attribute found!"); + return; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + if(attName == "number") + { + number = eval.EvaluateInteger(attValue); + } + } + ReplicavolRead(child, number); + } + else if(tag == "divisionvol") + { + DivisionvolRead(child); + } + else if(tag == "loop") + { + LoopRead(child, &G4GDMLRead::Volume_contentRead); + } + else + { + G4cout << "Treating unknown GDML tag in volume '" << tag + << "' as GDML extension..." << G4endl; + } + } + } + + void StructureRead( const xercesc::DOMElement* const structureElement) override + { + G4cout << "G4GDML: Reading structure..." << G4endl; + + for(xercesc::DOMNode* iter = structureElement->getFirstChild(); + iter != nullptr; iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = + dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadStructure::StructureRead()", "InvalidRead", + FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + + if(tag == "bordersurface") + { + BorderSurfaceRead(child); + } + else if(tag == "skinsurface") + { + SkinSurfaceRead(child); + } + else if(tag == "volume") + { + VolumeRead(child); + } + else if(tag == "assembly") + { + AssemblyRead(child); + } + else if(tag == "loop") + { + LoopRead(child, &G4GDMLRead::StructureRead); + } + else + { + G4String error_msg = "Unknown tag in structure: " + tag; + G4Exception("G4GDMLReadStructure::StructureRead()", "ReadError", + FatalException, error_msg); + } + } + } + + +}; diff --git a/remoll.cc b/remoll.cc index b4def9fe3..34f81c52f 100644 --- a/remoll.cc +++ b/remoll.cc @@ -28,6 +28,7 @@ typedef G4RunManager RunManager; #include "remollActionInitialization.hh" #include "remollDetectorConstruction.hh" #include "remollParallelConstruction.hh" +#include "remollGDMLReadStructure.hh" #include "remollSearchPath.hh" @@ -155,8 +156,9 @@ int main(int argc, char** argv) { remollIO::GetInstance(); // Detector geometry + remollGDMLReadStructure* rs = new remollGDMLReadStructure(); G4String material_name = "material"; - remollDetectorConstruction* detector = new remollDetectorConstruction(material_name, geometry_gdmlfile); + remollDetectorConstruction* detector = new remollDetectorConstruction(material_name, geometry_gdmlfile,rs); // Parallel world geometry G4String parallel_name = "parallel"; // Note: name must correspond with name of G4ParallelWorldPhysics remollParallelConstruction* parallel = new remollParallelConstruction(parallel_name, parallel_gdmlfile); diff --git a/src/remollDetectorConstruction.cc b/src/remollDetectorConstruction.cc index 09caaf0cb..5228c01d8 100644 --- a/src/remollDetectorConstruction.cc +++ b/src/remollDetectorConstruction.cc @@ -28,6 +28,7 @@ // GDML export #include "G4GDMLParser.hh" +#include "remollGDMLReadStructure.hh" //CADMesh #ifdef __USE_CADMESH @@ -53,8 +54,9 @@ G4ThreadLocal remollGlobalField* remollDetectorConstruction::fGlobalField = 0; G4UserLimits* remollDetectorConstruction::fKryptoniteUserLimits = new G4UserLimits(0,0,0,DBL_MAX,DBL_MAX); -remollDetectorConstruction::remollDetectorConstruction(const G4String& name, const G4String& gdmlfile) +remollDetectorConstruction::remollDetectorConstruction(const G4String& name, const G4String& gdmlfile, remollGDMLReadStructure* read) : fVerboseLevel(0), + fGDMLParser(read), fGDMLValidate(false), fGDMLOverlapCheck(true), fGDMLPath(""), From 88072f6747ed4d1319407a6d2a39d9bbb01fea7d Mon Sep 17 00:00:00 2001 From: Prakash Date: Wed, 14 Feb 2024 23:26:02 -0500 Subject: [PATCH 2/4] Implement general rotation order --- include/remollGDMLReadStructure.hh | 222 ++++++++++++++--------------- 1 file changed, 107 insertions(+), 115 deletions(-) diff --git a/include/remollGDMLReadStructure.hh b/include/remollGDMLReadStructure.hh index bbbfe4269..47acb9da1 100644 --- a/include/remollGDMLReadStructure.hh +++ b/include/remollGDMLReadStructure.hh @@ -1,33 +1,98 @@ #pragma once +#include "G4UnitsTable.hh" #include "G4LogicalVolume.hh" #include "G4VPhysicalVolume.hh" -#include "G4PVPlacement.hh" -#include "G4LogicalVolumeStore.hh" -#include "G4PhysicalVolumeStore.hh" #include "G4AssemblyVolume.hh" #include "G4ReflectionFactory.hh" -#include "G4PVDivisionFactory.hh" -#include "G4LogicalBorderSurface.hh" -#include "G4LogicalSkinSurface.hh" #include "G4GDMLReadStructure.hh" +void rotatex(G4RotationMatrix& rot,G4ThreeVector anglevec){ + rot.rotateX(anglevec.x()); +} +void rotatey(G4RotationMatrix& rot,G4ThreeVector anglevec){ + rot.rotateY(anglevec.y()); +} + +void rotatez(G4RotationMatrix& rot,G4ThreeVector anglevec){ + rot.rotateZ(anglevec.z()); +} + +std::map>> rotation_map{ + {"xyz", {rotatex,rotatey,rotatez}}, + {"xzy", {rotatex,rotatez,rotatey}}, + {"yxz", {rotatey,rotatex,rotatez}}, + {"yzx", {rotatey,rotatez,rotatex}}, + {"zxy", {rotatez,rotatex,rotatey}}, + {"zyx", {rotatez,rotatey,rotatex}}, +}; class remollGDMLReadStructure : public G4GDMLReadStructure { public: remollGDMLReadStructure():G4GDMLReadStructure(){ } - G4LogicalVolume* FileReads(const xercesc::DOMElement* const fileElement) + void RotationRead(const xercesc::DOMElement* const vectorElement, G4RotationMatrix& rot){ + G4double unit = 1.0; + G4ThreeVector vec; + std::string order = "xyz"; + + const xercesc::DOMNamedNodeMap* const attributes = vectorElement->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + + const xercesc::DOMAttr* const attribute = dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLRead::VectorRead()", "InvalidRead", FatalException, "No attribute found!"); + return; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + + if(attName == "unit") + { + unit = G4UnitDefinition::GetValueOf(attValue); + } + else if(attName == "x") + { + vec.setX(eval.Evaluate(attValue)); + } + else if(attName == "y") + { + vec.setY(eval.Evaluate(attValue)); + } + else if(attName == "z") + { + vec.setZ(eval.Evaluate(attValue)); + } + else if (attName == "order" ){ + order = attValue; + } + } + + vec *= unit; + + for(auto rotate : rotation_map[order]) rotate(rot,vec); + + } + + G4LogicalVolume* FileRead(const xercesc::DOMElement* const fileElement) { G4String name; G4String volname; - const xercesc::DOMNamedNodeMap* const attributes = - fileElement->getAttributes(); + const xercesc::DOMNamedNodeMap* const attributes = fileElement->getAttributes(); XMLSize_t attributeCount = attributes->getLength(); - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; - ++attribute_index) + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) { xercesc::DOMNode* attribute_node = attributes->item(attribute_index); @@ -87,24 +152,22 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } - void PhysvolReads( const xercesc::DOMElement* const physvolElement, G4AssemblyVolume* pAssembly=0) + void PhysvolRead( const xercesc::DOMElement* const physvolElement, G4AssemblyVolume* pAssembly=0) { G4String name; G4LogicalVolume* logvol = nullptr; G4AssemblyVolume* assembly = nullptr; G4ThreeVector position(0.0, 0.0, 0.0); G4ThreeVector rotation(0.0, 0.0, 0.0); - G4ThreeVector rotationzyx(0.0, 0.0, 0.0); - bool zyx_requested = false; + G4RotationMatrix orotm; + bool orotation_requested = false; G4ThreeVector scale(1.0, 1.0, 1.0); G4int copynumber = 0; - const xercesc::DOMNamedNodeMap* const attributes = - physvolElement->getAttributes(); + const xercesc::DOMNamedNodeMap* const attributes = physvolElement->getAttributes(); XMLSize_t attributeCount = attributes->getLength(); - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; - ++attribute_index) + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) { xercesc::DOMNode* attribute_node = attributes->item(attribute_index); @@ -134,28 +197,24 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } } - for(xercesc::DOMNode* iter = physvolElement->getFirstChild(); iter != nullptr; - iter = iter->getNextSibling()) + for(xercesc::DOMNode* iter = physvolElement->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) { if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) { continue; } - const xercesc::DOMElement* const child = - dynamic_cast(iter); + const xercesc::DOMElement* const child = dynamic_cast(iter); if(child == nullptr) { - G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", - FatalException, "No child found!"); + G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", FatalException, "No child found!"); return; } const G4String tag = Transcode(child->getTagName()); if(tag == "volumeref") { - const G4String& child_name = GenerateName(RefRead(child)); - assembly = GetAssembly(child_name); + const G4String& child_name = GenerateName(RefRead(child)); assembly = GetAssembly(child_name); if(assembly == nullptr) { logvol = GetVolume(child_name); @@ -163,7 +222,7 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } else if(tag == "file") { - logvol = FileReads(child); + logvol = FileRead(child); } else if(tag == "position") { @@ -173,10 +232,10 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { { VectorRead(child, rotation); } - else if(tag == "rotationzyx") + else if(tag == "orotation") { - VectorRead(child, rotationzyx); - zyx_requested = true; + RotationRead(child, orotm); + orotation_requested = true; } else if(tag == "scale") { @@ -204,17 +263,10 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } G4Transform3D transform; - if(!zyx_requested){ + if(!orotation_requested){ transform = G4Transform3D(GetRotationMatrix(rotation).inverse(), position); } else { - G4RotationMatrix rot; - - rot.rotateY(rotationzyx.y()); - rot.rotateZ(rotationzyx.z()); - rot.rotateX(rotationzyx.x()); - rot.rectify(); // Rectify matrix from possible roundoff errors - - transform = G4Transform3D(rot.inverse(),position); + transform = G4Transform3D(orotm.inverse(),position); } transform = transform * G4Scale3D(scale.x(), scale.y(), scale.z()); @@ -244,8 +296,7 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { return; } G4String pv_name = logvol->GetName() + "_PV"; - G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()->Place( - transform, pv_name, logvol, pMotherLogical, false, copynumber, check); + G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()->Place(transform, pv_name, logvol, pMotherLogical, false, copynumber, check); if(pair.first != nullptr) { @@ -259,18 +310,17 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } } - void ParametersRead(const xercesc::DOMElement* const element) { - + void ParametersRead(const xercesc::DOMElement* const element) + { G4ThreeVector rotation(0.0, 0.0, 0.0); - G4ThreeVector rotationzyx(0.0, 0.0, 0.0); G4ThreeVector position(0.0, 0.0, 0.0); - bool zyx_requested = false; + bool orotation_requested = false; + G4RotationMatrix orotm; G4GDMLParameterisation::PARAMETER parameter; parameter.pRot = new G4RotationMatrix(); - for(xercesc::DOMNode* iter = element->getFirstChild(); iter != nullptr; - iter = iter->getNextSibling()) + for(xercesc::DOMNode* iter = element->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) { if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) { @@ -290,10 +340,10 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { { VectorRead(child, rotation); } - else if(tag == "rotationzyx") + else if(tag == "orotation") { - VectorRead(child, rotationzyx); - zyx_requested = true; + RotationRead(child,orotm); + orotation_requested = true; } else if(tag == "position") { @@ -368,15 +418,12 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } - - if(!zyx_requested){ + if(orotation_requested){ + parameter.pRot = std::move(&orotm); + } else { parameter.pRot->rotateX(rotation.x()); parameter.pRot->rotateY(rotation.y()); parameter.pRot->rotateZ(rotation.z()); - } else { - parameter.pRot->rotateY(rotationzyx.y()); - parameter.pRot->rotateZ(rotationzyx.z()); - parameter.pRot->rotateX(rotationzyx.x()); } @@ -388,8 +435,7 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { void Volume_contentRead( const xercesc::DOMElement* const volumeElement) override { - for(xercesc::DOMNode* iter = volumeElement->getFirstChild(); iter != nullptr; - iter = iter->getNextSibling()) + for(xercesc::DOMNode* iter = volumeElement->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) { if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) { @@ -416,15 +462,14 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } else if(tag == "physvol") { - PhysvolReads(child); + PhysvolRead(child); } else if(tag == "replicavol") { G4int number = 1; const xercesc::DOMNamedNodeMap* const attributes = child->getAttributes(); XMLSize_t attributeCount = attributes->getLength(); - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; - ++attribute_index) + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) { xercesc::DOMNode* attribute_node = attributes->item(attribute_index); if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) @@ -458,62 +503,9 @@ class remollGDMLReadStructure : public G4GDMLReadStructure { } else { - G4cout << "Treating unknown GDML tag in volume '" << tag - << "' as GDML extension..." << G4endl; - } - } - } - - void StructureRead( const xercesc::DOMElement* const structureElement) override - { - G4cout << "G4GDML: Reading structure..." << G4endl; - - for(xercesc::DOMNode* iter = structureElement->getFirstChild(); - iter != nullptr; iter = iter->getNextSibling()) - { - if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) - { - continue; - } - - const xercesc::DOMElement* const child = - dynamic_cast(iter); - if(child == nullptr) - { - G4Exception("G4GDMLReadStructure::StructureRead()", "InvalidRead", - FatalException, "No child found!"); - return; - } - const G4String tag = Transcode(child->getTagName()); - - if(tag == "bordersurface") - { - BorderSurfaceRead(child); - } - else if(tag == "skinsurface") - { - SkinSurfaceRead(child); - } - else if(tag == "volume") - { - VolumeRead(child); - } - else if(tag == "assembly") - { - AssemblyRead(child); - } - else if(tag == "loop") - { - LoopRead(child, &G4GDMLRead::StructureRead); - } - else - { - G4String error_msg = "Unknown tag in structure: " + tag; - G4Exception("G4GDMLReadStructure::StructureRead()", "ReadError", - FatalException, error_msg); + G4cout << "Treating unknown GDML tag in volume '" << tag << "' as GDML extension..." << G4endl; } } } - }; From 72174b5259de8e935fbf240da84ccf5e89677e37 Mon Sep 17 00:00:00 2001 From: Prakash Date: Fri, 23 Feb 2024 14:40:03 -0500 Subject: [PATCH 3/4] Split into source and header files and cleanup a bit --- include/remollGDMLReadStructure.hh | 520 +---------------------------- src/remollGDMLReadStructure.cc | 502 ++++++++++++++++++++++++++++ 2 files changed, 519 insertions(+), 503 deletions(-) create mode 100644 src/remollGDMLReadStructure.cc diff --git a/include/remollGDMLReadStructure.hh b/include/remollGDMLReadStructure.hh index 47acb9da1..6d461bc59 100644 --- a/include/remollGDMLReadStructure.hh +++ b/include/remollGDMLReadStructure.hh @@ -1,511 +1,25 @@ -#pragma once - -#include "G4UnitsTable.hh" -#include "G4LogicalVolume.hh" -#include "G4VPhysicalVolume.hh" -#include "G4AssemblyVolume.hh" -#include "G4ReflectionFactory.hh" -#include "G4GDMLReadStructure.hh" +// -*- coding: utf-8 -*- +// vim: ai ts=2 sts=2 et sw=2 ft=cpp +// author : Prakash +// date : 2024-02-18 -void rotatex(G4RotationMatrix& rot,G4ThreeVector anglevec){ - rot.rotateX(anglevec.x()); -} -void rotatey(G4RotationMatrix& rot,G4ThreeVector anglevec){ - rot.rotateY(anglevec.y()); -} +#pragma once -void rotatez(G4RotationMatrix& rot,G4ThreeVector anglevec){ - rot.rotateZ(anglevec.z()); -} +#include +#include -std::map>> rotation_map{ - {"xyz", {rotatex,rotatey,rotatez}}, - {"xzy", {rotatex,rotatez,rotatey}}, - {"yxz", {rotatey,rotatex,rotatez}}, - {"yzx", {rotatey,rotatez,rotatex}}, - {"zxy", {rotatez,rotatex,rotatey}}, - {"zyx", {rotatez,rotatey,rotatex}}, -}; +#include "G4GDMLReadStructure.hh" class remollGDMLReadStructure : public G4GDMLReadStructure { - public: - remollGDMLReadStructure():G4GDMLReadStructure(){ } - - void RotationRead(const xercesc::DOMElement* const vectorElement, G4RotationMatrix& rot){ - G4double unit = 1.0; - G4ThreeVector vec; - std::string order = "xyz"; - - const xercesc::DOMNamedNodeMap* const attributes = vectorElement->getAttributes(); - XMLSize_t attributeCount = attributes->getLength(); - - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) - { - xercesc::DOMNode* attribute_node = attributes->item(attribute_index); - - if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) - { - continue; - } - - const xercesc::DOMAttr* const attribute = dynamic_cast(attribute_node); - if(attribute == nullptr) - { - G4Exception("G4GDMLRead::VectorRead()", "InvalidRead", FatalException, "No attribute found!"); - return; - } - const G4String attName = Transcode(attribute->getName()); - const G4String attValue = Transcode(attribute->getValue()); - - if(attName == "unit") - { - unit = G4UnitDefinition::GetValueOf(attValue); - } - else if(attName == "x") - { - vec.setX(eval.Evaluate(attValue)); - } - else if(attName == "y") - { - vec.setY(eval.Evaluate(attValue)); - } - else if(attName == "z") - { - vec.setZ(eval.Evaluate(attValue)); - } - else if (attName == "order" ){ - order = attValue; - } - } - - vec *= unit; - - for(auto rotate : rotation_map[order]) rotate(rot,vec); - - } - - G4LogicalVolume* FileRead(const xercesc::DOMElement* const fileElement) - { - G4String name; - G4String volname; - - const xercesc::DOMNamedNodeMap* const attributes = fileElement->getAttributes(); - XMLSize_t attributeCount = attributes->getLength(); - - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) - { - xercesc::DOMNode* attribute_node = attributes->item(attribute_index); - - if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) - { - continue; - } - - const xercesc::DOMAttr* const attribute = - dynamic_cast(attribute_node); - if(attribute == nullptr) - { - G4Exception("G4GDMLReadStructure::FileRead()", "InvalidRead", - FatalException, "No attribute found!"); - return nullptr; - } - const G4String attName = Transcode(attribute->getName()); - const G4String attValue = Transcode(attribute->getValue()); - - if(attName == "name") - { - name = attValue; - } - else if(attName == "volname") - { - volname = attValue; - } - } - - const G4bool isModule = true; - remollGDMLReadStructure structure; - structure.Read(name, validate, isModule); - - // Register existing auxiliar information defined in child module - // - const G4GDMLAuxMapType* aux = structure.GetAuxMap(); - if(!aux->empty()) - { - for(auto pos = aux->cbegin(); pos != aux->cend(); ++pos) - { - auxMap.insert(std::make_pair(pos->first, pos->second)); - } - } - - // Return volume structure from child module - // - if(volname.empty()) - { - return structure.GetVolume(structure.GetSetup("Default")); - //return this->GetVolume(this->GetSetup("Default")); - } - else - { - return structure.GetVolume(structure.GenerateName(volname)); - //return this->GetVolume(this->GenerateName(volname)); - } - } - - - void PhysvolRead( const xercesc::DOMElement* const physvolElement, G4AssemblyVolume* pAssembly=0) - { - G4String name; - G4LogicalVolume* logvol = nullptr; - G4AssemblyVolume* assembly = nullptr; - G4ThreeVector position(0.0, 0.0, 0.0); - G4ThreeVector rotation(0.0, 0.0, 0.0); - G4RotationMatrix orotm; - bool orotation_requested = false; - G4ThreeVector scale(1.0, 1.0, 1.0); - G4int copynumber = 0; - - const xercesc::DOMNamedNodeMap* const attributes = physvolElement->getAttributes(); - XMLSize_t attributeCount = attributes->getLength(); - - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) - { - xercesc::DOMNode* attribute_node = attributes->item(attribute_index); - - if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) - { - continue; - } - - const xercesc::DOMAttr* const attribute = - dynamic_cast(attribute_node); - if(attribute == nullptr) - { - G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", - FatalException, "No attribute found!"); - return; - } - const G4String attName = Transcode(attribute->getName()); - const G4String attValue = Transcode(attribute->getValue()); - - if(attName == "name") - { - name = attValue; - } - if(attName == "copynumber") - { - copynumber = eval.EvaluateInteger(attValue); - } - } - - for(xercesc::DOMNode* iter = physvolElement->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) - { - if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) - { - continue; - } - - const xercesc::DOMElement* const child = dynamic_cast(iter); - if(child == nullptr) - { - G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", FatalException, "No child found!"); - return; - } - const G4String tag = Transcode(child->getTagName()); - - if(tag == "volumeref") - { - const G4String& child_name = GenerateName(RefRead(child)); assembly = GetAssembly(child_name); - if(assembly == nullptr) - { - logvol = GetVolume(child_name); - } - } - else if(tag == "file") - { - logvol = FileRead(child); - } - else if(tag == "position") - { - VectorRead(child, position); - } - else if(tag == "rotation") - { - VectorRead(child, rotation); - } - else if(tag == "orotation") - { - RotationRead(child, orotm); - orotation_requested = true; - } - else if(tag == "scale") - { - VectorRead(child, scale); - } - else if(tag == "positionref") - { - position = GetPosition(GenerateName(RefRead(child))); - } - else if(tag == "rotationref") - { - rotation = GetRotation(GenerateName(RefRead(child))); - } - else if(tag == "scaleref") - { - scale = GetScale(GenerateName(RefRead(child))); - } - else - { - G4String error_msg = "Unknown tag in physvol: " + tag; - G4Exception("G4GDMLReadStructure::PhysvolRead()", "ReadError", - FatalException, error_msg); - return; - } - } - - G4Transform3D transform; - if(!orotation_requested){ - transform = G4Transform3D(GetRotationMatrix(rotation).inverse(), position); - } else { - transform = G4Transform3D(orotm.inverse(),position); - } - - transform = transform * G4Scale3D(scale.x(), scale.y(), scale.z()); - - if(pAssembly != nullptr) // Fill assembly structure - { - if(assembly != nullptr) // Case of recursive assemblies - { - pAssembly->AddPlacedAssembly(assembly, transform); - } - if(logvol == nullptr) - { - return; - } - pAssembly->AddPlacedVolume(logvol, transform); - } - else // Generate physical volume tree or do assembly imprint - { - if(assembly != nullptr) - { - assembly->MakeImprint(pMotherLogical, transform, 0, check); - } - else - { - if(logvol == nullptr) - { - return; - } - G4String pv_name = logvol->GetName() + "_PV"; - G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()->Place(transform, pv_name, logvol, pMotherLogical, false, copynumber, check); - - if(pair.first != nullptr) - { - GeneratePhysvolName(name, pair.first); - } - if(pair.second != nullptr) - { - GeneratePhysvolName(name, pair.second); - } - } - } - } - - void ParametersRead(const xercesc::DOMElement* const element) - { - G4ThreeVector rotation(0.0, 0.0, 0.0); - G4ThreeVector position(0.0, 0.0, 0.0); - bool orotation_requested = false; - G4RotationMatrix orotm; - - G4GDMLParameterisation::PARAMETER parameter; - parameter.pRot = new G4RotationMatrix(); - - for(xercesc::DOMNode* iter = element->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) - { - if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) - { - continue; - } - - const xercesc::DOMElement* const child = - dynamic_cast(iter); - if(child == nullptr) - { - G4Exception("G4GDMLReadParamvol::ParametersRead()", "InvalidRead", - FatalException, "No child found!"); - return; - } - const G4String tag = Transcode(child->getTagName()); - if(tag == "rotation") - { - VectorRead(child, rotation); - } - else if(tag == "orotation") - { - RotationRead(child,orotm); - orotation_requested = true; - } - else if(tag == "position") - { - VectorRead(child, position); - } - else if(tag == "positionref") - { - position = GetPosition(GenerateName(RefRead(child))); - } - else if(tag == "rotationref") - { - rotation = GetRotation(GenerateName(RefRead(child))); - } - else if(tag == "box_dimensions") - { - Box_dimensionsRead(child, parameter); - } - else if(tag == "trd_dimensions") - { - Trd_dimensionsRead(child, parameter); - } - else if(tag == "trap_dimensions") - { - Trap_dimensionsRead(child, parameter); - } - else if(tag == "tube_dimensions") - { - Tube_dimensionsRead(child, parameter); - } - else if(tag == "cone_dimensions") - { - Cone_dimensionsRead(child, parameter); - } - else if(tag == "sphere_dimensions") - { - Sphere_dimensionsRead(child, parameter); - } - else if(tag == "orb_dimensions") - { - Orb_dimensionsRead(child, parameter); - } - else if(tag == "torus_dimensions") - { - Torus_dimensionsRead(child, parameter); - } - else if(tag == "ellipsoid_dimensions") - { - Ellipsoid_dimensionsRead(child, parameter); - } - else if(tag == "para_dimensions") - { - Para_dimensionsRead(child, parameter); - } - else if(tag == "polycone_dimensions") - { - Polycone_dimensionsRead(child, parameter); - } - else if(tag == "polyhedra_dimensions") - { - Polyhedra_dimensionsRead(child, parameter); - } - else if(tag == "hype_dimensions") - { - Hype_dimensionsRead(child, parameter); - } - else - { - G4String error_msg = "Unknown tag in parameters: " + tag; - G4Exception("G4GDMLReadParamvol::ParametersRead()", "ReadError", - FatalException, error_msg); - } - } - - - if(orotation_requested){ - parameter.pRot = std::move(&orotm); - } else { - parameter.pRot->rotateX(rotation.x()); - parameter.pRot->rotateY(rotation.y()); - parameter.pRot->rotateZ(rotation.z()); - } - - - parameter.position = position; - - parameterisation->AddParameter(parameter); - - } - - void Volume_contentRead( const xercesc::DOMElement* const volumeElement) override - { - for(xercesc::DOMNode* iter = volumeElement->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) - { - if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) - { - continue; - } - - const xercesc::DOMElement* const child = - dynamic_cast(iter); - if(child == nullptr) - { - G4Exception("G4GDMLReadStructure::Volume_contentRead()", "InvalidRead", - FatalException, "No child found!"); - return; - } - const G4String tag = Transcode(child->getTagName()); + public: + remollGDMLReadStructure(); - if((tag == "auxiliary") || (tag == "materialref") || (tag == "solidref")) - { - // These are already processed in VolumeRead() - } - else if(tag == "paramvol") - { - ParamvolRead(child, pMotherLogical); - } - else if(tag == "physvol") - { - PhysvolRead(child); - } - else if(tag == "replicavol") - { - G4int number = 1; - const xercesc::DOMNamedNodeMap* const attributes = child->getAttributes(); - XMLSize_t attributeCount = attributes->getLength(); - for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) - { - xercesc::DOMNode* attribute_node = attributes->item(attribute_index); - if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) - { - continue; - } - const xercesc::DOMAttr* const attribute = - dynamic_cast(attribute_node); - if(attribute == nullptr) - { - G4Exception("G4GDMLReadStructure::Volume_contentRead()", - "InvalidRead", FatalException, "No attribute found!"); - return; - } - const G4String attName = Transcode(attribute->getName()); - const G4String attValue = Transcode(attribute->getValue()); - if(attName == "number") - { - number = eval.EvaluateInteger(attValue); - } - } - ReplicavolRead(child, number); - } - else if(tag == "divisionvol") - { - DivisionvolRead(child); - } - else if(tag == "loop") - { - LoopRead(child, &G4GDMLRead::Volume_contentRead); - } - else - { - G4cout << "Treating unknown GDML tag in volume '" << tag << "' as GDML extension..." << G4endl; - } - } - } + void RotationRead(const xercesc::DOMElement* const vectorElement, G4RotationMatrix& rot); + G4LogicalVolume* FileRead(const xercesc::DOMElement* const fileElement); + void PhysvolRead( const xercesc::DOMElement* const physvolElement, G4AssemblyVolume* pAssembly=0); + void ParametersRead(const xercesc::DOMElement* const element); + void Volume_contentRead( const xercesc::DOMElement* const volumeElement) override; + private: + std::map>> rotations; }; diff --git a/src/remollGDMLReadStructure.cc b/src/remollGDMLReadStructure.cc new file mode 100644 index 000000000..3a6c6a180 --- /dev/null +++ b/src/remollGDMLReadStructure.cc @@ -0,0 +1,502 @@ +// -*- coding: utf-8 -*- +// vim: ai ts=2 sts=2 et sw=2 ft=cpp +// author : Prakash [pgautam@jlab.org] +// date : 2024-02-18 + +#include "G4UnitsTable.hh" +#include "G4LogicalVolume.hh" +#include "G4VPhysicalVolume.hh" +#include "G4AssemblyVolume.hh" +#include "G4ReflectionFactory.hh" + +#include "remollGDMLReadStructure.hh" + +remollGDMLReadStructure::remollGDMLReadStructure():G4GDMLReadStructure() +{ + auto rotatex=[](G4RotationMatrix& rot,G4ThreeVector angle) { rot.rotateX(angle.x()); }; + auto rotatey=[](G4RotationMatrix& rot,G4ThreeVector angle) { rot.rotateY(angle.y()); }; + auto rotatez=[](G4RotationMatrix& rot,G4ThreeVector angle) { rot.rotateZ(angle.z()); }; + rotations = { + {"xyz", {rotatex,rotatey,rotatez}}, + {"xzy", {rotatex,rotatez,rotatey}}, + {"yxz", {rotatey,rotatex,rotatez}}, + {"yzx", {rotatey,rotatez,rotatex}}, + {"zxy", {rotatez,rotatex,rotatey}}, + {"zyx", {rotatez,rotatey,rotatex}} + }; + +} + +void remollGDMLReadStructure::RotationRead(const xercesc::DOMElement* const vectorElement, G4RotationMatrix& rot) +{ + G4double unit = 1.0; + G4ThreeVector vec; + std::string order = "xyz"; + + const xercesc::DOMNamedNodeMap* const attributes = vectorElement->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + + const xercesc::DOMAttr* const attribute = dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLRead::VectorRead()", "InvalidRead", FatalException, "No attribute found!"); + return; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + + if(attName == "unit") + { + unit = G4UnitDefinition::GetValueOf(attValue); + } + else if(attName == "x") + { + vec.setX(eval.Evaluate(attValue)); + } + else if(attName == "y") + { + vec.setY(eval.Evaluate(attValue)); + } + else if(attName == "z") + { + vec.setZ(eval.Evaluate(attValue)); + } + else if (attName == "order" ){ + order = attValue; + } + } + + vec *= unit; + + for(auto rotate : rotations.at(order)) rotate(rot,vec); + +} + +G4LogicalVolume* remollGDMLReadStructure::FileRead(const xercesc::DOMElement* const fileElement) +{ + G4String name; + G4String volname; + + const xercesc::DOMNamedNodeMap* const attributes = fileElement->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + + const xercesc::DOMAttr* const attribute = + dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLReadStructure::FileRead()", "InvalidRead", + FatalException, "No attribute found!"); + return nullptr; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + + if(attName == "name") + { + name = attValue; + } + else if(attName == "volname") + { + volname = attValue; + } + } + + const G4bool isModule = true; + remollGDMLReadStructure structure; + structure.Read(name, validate, isModule); + + // Register existing auxiliar information defined in child module + // + const G4GDMLAuxMapType* aux = structure.GetAuxMap(); + if(!aux->empty()) + { + for(auto pos = aux->cbegin(); pos != aux->cend(); ++pos) + { + auxMap.insert(std::make_pair(pos->first, pos->second)); + } + } + + if(volname.empty()) + { + return structure.GetVolume(structure.GetSetup("Default")); + } + else + { + return structure.GetVolume(structure.GenerateName(volname)); + } +} + + +void remollGDMLReadStructure::PhysvolRead( const xercesc::DOMElement* const physvolElement, G4AssemblyVolume* pAssembly) +{ + G4String name; + G4LogicalVolume* logvol = nullptr; + G4AssemblyVolume* assembly = nullptr; + G4ThreeVector position(0.0, 0.0, 0.0); + G4ThreeVector rotation(0.0, 0.0, 0.0); + G4RotationMatrix orotm; + bool orotation_requested = false; + G4ThreeVector scale(1.0, 1.0, 1.0); + G4int copynumber = 0; + + const xercesc::DOMNamedNodeMap* const attributes = physvolElement->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + + const xercesc::DOMAttr* const attribute = + dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", + FatalException, "No attribute found!"); + return; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + + if(attName == "name") + { + name = attValue; + } + if(attName == "copynumber") + { + copynumber = eval.EvaluateInteger(attValue); + } + } + + for(xercesc::DOMNode* iter = physvolElement->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadStructure::PhysvolRead()", "InvalidRead", FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + + if(tag == "volumeref") + { + const G4String& child_name = GenerateName(RefRead(child)); assembly = GetAssembly(child_name); + if(assembly == nullptr) + { + logvol = GetVolume(child_name); + } + } + else if(tag == "file") + { + logvol = FileRead(child); + } + else if(tag == "position") + { + VectorRead(child, position); + } + else if(tag == "rotation") + { + VectorRead(child, rotation); + } + else if(tag == "orotation") + { + RotationRead(child, orotm); + orotation_requested = true; + } + else if(tag == "scale") + { + VectorRead(child, scale); + } + else if(tag == "positionref") + { + position = GetPosition(GenerateName(RefRead(child))); + } + else if(tag == "rotationref") + { + rotation = GetRotation(GenerateName(RefRead(child))); + } + else if(tag == "scaleref") + { + scale = GetScale(GenerateName(RefRead(child))); + } + else + { + G4String error_msg = "Unknown tag in physvol: " + tag; + G4Exception("G4GDMLReadStructure::PhysvolRead()", "ReadError", + FatalException, error_msg); + return; + } + } + + G4Transform3D transform; + if(!orotation_requested){ + transform = G4Transform3D(GetRotationMatrix(rotation).inverse(), position); + } else { + transform = G4Transform3D(orotm.inverse(),position); + } + + transform = transform * G4Scale3D(scale.x(), scale.y(), scale.z()); + + if(pAssembly != nullptr) // Fill assembly structure + { + if(assembly != nullptr) // Case of recursive assemblies + { + pAssembly->AddPlacedAssembly(assembly, transform); + } + if(logvol == nullptr) + { + return; + } + pAssembly->AddPlacedVolume(logvol, transform); + } + else // Generate physical volume tree or do assembly imprint + { + if(assembly != nullptr) + { + assembly->MakeImprint(pMotherLogical, transform, 0, check); + } + else + { + if(logvol == nullptr) + { + return; + } + G4String pv_name = logvol->GetName() + "_PV"; + G4PhysicalVolumesPair pair = G4ReflectionFactory::Instance()->Place(transform, pv_name, logvol, pMotherLogical, false, copynumber, check); + + if(pair.first != nullptr) + { + GeneratePhysvolName(name, pair.first); + } + if(pair.second != nullptr) + { + GeneratePhysvolName(name, pair.second); + } + } + } +} + +void remollGDMLReadStructure::ParametersRead(const xercesc::DOMElement* const element) +{ + G4ThreeVector rotation(0.0, 0.0, 0.0); + G4ThreeVector position(0.0, 0.0, 0.0); + bool orotation_requested = false; + G4RotationMatrix orotm; + + G4GDMLParameterisation::PARAMETER parameter; + parameter.pRot = new G4RotationMatrix(); + + for(xercesc::DOMNode* iter = element->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = + dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadParamvol::ParametersRead()", "InvalidRead", + FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + if(tag == "rotation") + { + VectorRead(child, rotation); + } + else if(tag == "orotation") + { + RotationRead(child,orotm); + orotation_requested = true; // makes orotation get priority over normal rotation. + } + else if(tag == "position") + { + VectorRead(child, position); + } + else if(tag == "positionref") + { + position = GetPosition(GenerateName(RefRead(child))); + } + else if(tag == "rotationref") + { + rotation = GetRotation(GenerateName(RefRead(child))); + } + else if(tag == "box_dimensions") + { + Box_dimensionsRead(child, parameter); + } + else if(tag == "trd_dimensions") + { + Trd_dimensionsRead(child, parameter); + } + else if(tag == "trap_dimensions") + { + Trap_dimensionsRead(child, parameter); + } + else if(tag == "tube_dimensions") + { + Tube_dimensionsRead(child, parameter); + } + else if(tag == "cone_dimensions") + { + Cone_dimensionsRead(child, parameter); + } + else if(tag == "sphere_dimensions") + { + Sphere_dimensionsRead(child, parameter); + } + else if(tag == "orb_dimensions") + { + Orb_dimensionsRead(child, parameter); + } + else if(tag == "torus_dimensions") + { + Torus_dimensionsRead(child, parameter); + } + else if(tag == "ellipsoid_dimensions") + { + Ellipsoid_dimensionsRead(child, parameter); + } + else if(tag == "para_dimensions") + { + Para_dimensionsRead(child, parameter); + } + else if(tag == "polycone_dimensions") + { + Polycone_dimensionsRead(child, parameter); + } + else if(tag == "polyhedra_dimensions") + { + Polyhedra_dimensionsRead(child, parameter); + } + else if(tag == "hype_dimensions") + { + Hype_dimensionsRead(child, parameter); + } + else + { + G4String error_msg = "Unknown tag in parameters: " + tag; + G4Exception("G4GDMLReadParamvol::ParametersRead()", "ReadError", + FatalException, error_msg); + } + } + + + if(orotation_requested){ + parameter.pRot = std::move(&orotm); + } else { + parameter.pRot->rotateX(rotation.x()); + parameter.pRot->rotateY(rotation.y()); + parameter.pRot->rotateZ(rotation.z()); + } + + + parameter.position = position; + + parameterisation->AddParameter(parameter); + +} + +void remollGDMLReadStructure::Volume_contentRead( const xercesc::DOMElement* const volumeElement) +{ + for(xercesc::DOMNode* iter = volumeElement->getFirstChild(); iter != nullptr; iter = iter->getNextSibling()) + { + if(iter->getNodeType() != xercesc::DOMNode::ELEMENT_NODE) + { + continue; + } + + const xercesc::DOMElement* const child = + dynamic_cast(iter); + if(child == nullptr) + { + G4Exception("G4GDMLReadStructure::Volume_contentRead()", "InvalidRead", + FatalException, "No child found!"); + return; + } + const G4String tag = Transcode(child->getTagName()); + + if((tag == "auxiliary") || (tag == "materialref") || (tag == "solidref")) + { + // These are already processed in VolumeRead() + } + else if(tag == "paramvol") + { + ParamvolRead(child, pMotherLogical); + } + else if(tag == "physvol") + { + PhysvolRead(child); + } + else if(tag == "replicavol") + { + G4int number = 1; + const xercesc::DOMNamedNodeMap* const attributes = child->getAttributes(); + XMLSize_t attributeCount = attributes->getLength(); + for(XMLSize_t attribute_index = 0; attribute_index < attributeCount; ++attribute_index) + { + xercesc::DOMNode* attribute_node = attributes->item(attribute_index); + if(attribute_node->getNodeType() != xercesc::DOMNode::ATTRIBUTE_NODE) + { + continue; + } + const xercesc::DOMAttr* const attribute = + dynamic_cast(attribute_node); + if(attribute == nullptr) + { + G4Exception("G4GDMLReadStructure::Volume_contentRead()", + "InvalidRead", FatalException, "No attribute found!"); + return; + } + const G4String attName = Transcode(attribute->getName()); + const G4String attValue = Transcode(attribute->getValue()); + if(attName == "number") + { + number = eval.EvaluateInteger(attValue); + } + } + ReplicavolRead(child, number); + } + else if(tag == "divisionvol") + { + DivisionvolRead(child); + } + else if(tag == "loop") + { + LoopRead(child, &G4GDMLRead::Volume_contentRead); + } + else + { + G4cout << "Treating unknown GDML tag in volume '" << tag << "' as GDML extension..." << G4endl; + } + } +} From 5001c82329d0da2719c1a0e82e694454f93d1804 Mon Sep 17 00:00:00 2001 From: Ciprian Gal Date: Tue, 27 Feb 2024 09:13:13 -0500 Subject: [PATCH 4/4] update roation of quartz tiles with the new orotation. Work done by Xiang --- .../DetectorArray/DetectorArray.gdml | 2276 ++++++++--------- 1 file changed, 1138 insertions(+), 1138 deletions(-) diff --git a/geometry/detector/ThinQuartz/DetectorArray/DetectorArray.gdml b/geometry/detector/ThinQuartz/DetectorArray/DetectorArray.gdml index ad25af5d4..92bcad72e 100755 --- a/geometry/detector/ThinQuartz/DetectorArray/DetectorArray.gdml +++ b/geometry/detector/ThinQuartz/DetectorArray/DetectorArray.gdml @@ -1,1142 +1,1142 @@ - - - - - + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +