diff --git a/GeoModelCore/GeoModelHelpers/src/GeoShapeUtils.cxx b/GeoModelCore/GeoModelHelpers/src/GeoShapeUtils.cxx index e026be00a60e1a4c141aebeca1647fed14ed8515..dd0aaa32d4332e923185535e8b7bd557021cbd28 100644 --- a/GeoModelCore/GeoModelHelpers/src/GeoShapeUtils.cxx +++ b/GeoModelCore/GeoModelHelpers/src/GeoShapeUtils.cxx @@ -28,6 +28,10 @@ #include <vector> +namespace{ + constexpr unsigned tubeApproxEdges = 8; +} + std::pair<const GeoShape* , const GeoShape*> getOps(const GeoShape* composed) { if (composed->typeID() == GeoShapeUnion::getClassTypeID()) { const GeoShapeUnion* unionShape = dynamic_cast<const GeoShapeUnion*>(composed); @@ -48,7 +52,7 @@ std::pair<const GeoShape* , const GeoShape*> getOps(const GeoShape* composed) { unsigned int countComposedShapes(const GeoShape* shape) { std::pair<const GeoShape*, const GeoShape*> ops = getOps(shape); return 1 + (ops.first ? countComposedShapes(ops.first) : 0) + - (ops.second ? countComposedShapes(ops.second) : 0); + (ops.second ? countComposedShapes(ops.second) : 0); } GeoIntrusivePtr<const GeoShape> compressShift(const GeoShape* shift) { @@ -229,6 +233,46 @@ std::vector<GeoTrf::Vector3D> getPolyShapeEdges(const GeoShape* shape, } } + } else if (shape->typeID() == GeoTube::getClassTypeID()) { + constexpr double stepLength = 2.*M_PI / tubeApproxEdges; + const GeoTube* tube = static_cast<const GeoTube*>(shape); + for (const double z :{-1., 1.}){ + for (unsigned int e = 0; e < tubeApproxEdges; ++e) { + edgePoints.emplace_back(tube->getRMax()*std::cos(stepLength*e), + tube->getRMax()*std::sin(stepLength*e), + z*tube->getZHalfLength()); + } + } + } else if (shape->typeID() == GeoTubs::getClassTypeID()) { + const GeoTubs* tubs = static_cast<const GeoTubs*>(shape); + const double stepSize = tubs->getDPhi() / tubeApproxEdges; + for (const double z : {-1., 1.}) { + for (unsigned e = 0; e <= tubeApproxEdges; ++e) { + edgePoints.emplace_back(tubs->getRMax()*std::cos(tubs->getSPhi() + stepSize*e), + tubs->getRMax()*std::sin(tubs->getSPhi() + stepSize*e), + z*tubs->getZHalfLength()); + } + } + } else if (shape->typeID() == GeoEllipticalTube::getClassTypeID()) { + const GeoEllipticalTube* tube = static_cast<const GeoEllipticalTube*>(shape); + constexpr double stepSize = 2.*M_PI / tubeApproxEdges; + for (const double z : {-1.,1.}) { + for (unsigned e = 0; e<tubeApproxEdges; ++e){ + edgePoints.emplace_back(tube->getXHalfLength()*std::cos(stepSize*e), + tube->getYHalfLength()*std::sin(stepSize*e), + z* tube->getZHalfLength()); + } + } + } else if (shape->typeID() == GeoTorus::getClassTypeID()) { + constexpr double stepSize = 2.*M_PI / tubeApproxEdges; + const GeoTorus* torus = static_cast<const GeoTorus*>(shape); + for (unsigned int donut =0; donut < tubeApproxEdges; ++donut) { + for (unsigned int slice =0; slice < tubeApproxEdges; ++slice) { + edgePoints.emplace_back((torus->getRTor() + torus->getRMax()* std::cos(slice *stepSize))*std::cos(donut*stepSize), + (torus->getRTor() + torus->getRMax()* std::cos(slice *stepSize))*std::sin(donut*stepSize), + torus->getRMax()* std::sin(slice*stepSize)); + } + } } else { THROW_EXCEPTION("The shape "<<shape->type()<<" is not supported. Please add it to the list"); } diff --git a/GeoModelCore/GeoModelHelpers/src/MaterialManager.cxx b/GeoModelCore/GeoModelHelpers/src/MaterialManager.cxx index f9b996e2364769588d62601bb6a037a5d5817fc0..70e370fc94f5f6bc402cc21db4e4f71447f13c59 100644 --- a/GeoModelCore/GeoModelHelpers/src/MaterialManager.cxx +++ b/GeoModelCore/GeoModelHelpers/src/MaterialManager.cxx @@ -33,15 +33,14 @@ namespace { const GeoMaterial* MaterialManager::MaterialFactory::get() const { return m_material; } MaterialManager::MaterialFactory::~MaterialFactory() { if (!m_material) return; - // std::cout<<"MaterialFactory() -- finalize " <<m_material->getName()<<std::endl; + if (m_components.empty()) { m_material->lock(); return; } const double inv_totalFraction = 1. / m_totFraction; for(const auto& [element, fraction] : m_components) { - m_material->add(element, - fraction * element->getA() * inv_totalFraction); + m_material->add(element, fraction * element->getA() * inv_totalFraction); } m_material->lock(); } @@ -56,8 +55,7 @@ void MaterialManager::MaterialFactory::addComponent(const ConstMaterialPtr& mat, addComponent(elePtr, mat->getFraction(ele) * fraction); } } -void MaterialManager::MaterialFactory::addComponent(const ConstElementPtr& ele, - double fraction) { +void MaterialManager::MaterialFactory::addComponent(const ConstElementPtr& ele, double fraction) { m_components.emplace_back(std::make_pair(ele, fraction)); m_totFraction += fraction; } @@ -135,7 +133,7 @@ void MaterialManager::printAll() const { } void MaterialManager::addElement(const std::string &name, const std::string &symbol, double z, double a) { - GeoIntrusivePtr<GeoElement> newElement{new GeoElement(name,symbol,z,a*atomicDensity)}; + GeoIntrusivePtr<GeoElement> newElement{make_intrusive<GeoElement>(name,symbol,z,a*atomicDensity)}; addElement(newElement); } @@ -158,7 +156,7 @@ void MaterialManager::setMaterialNamespace(const std::string &name) { void MaterialManager::lockMaterial() { m_factory.reset(); } void MaterialManager::addMaterial(const std::string &name, double density) { - MaterialPtr newMat{new GeoMaterial(name, density * volDensity)}; + MaterialPtr newMat{make_intrusive<GeoMaterial>(name, density * volDensity)}; addMaterial(newMat); } @@ -200,14 +198,14 @@ void MaterialManager::addMatComponent(const std::string &name, double fraction) } void MaterialManager::buildSpecialMaterials() { // Ether - GeoIntrusivePtr<GeoElement> ethElement{new GeoElement("Ether","ET",500.0,0.0)}; + GeoIntrusivePtr<GeoElement> ethElement{make_intrusive<GeoElement>("Ether","ET",500.0,0.0)}; m_elements.insert(std::make_pair("Ether",ethElement)); - GeoIntrusivePtr<GeoMaterial> ether{new GeoMaterial("special::Ether",0.0)}; + GeoIntrusivePtr<GeoMaterial> ether{make_intrusive<GeoMaterial>("special::Ether",0.0)}; ether->add(ethElement,1.); ether->lock(); m_materials.insert(std::make_pair("special::Ether", ether)); // HyperUranium - GeoIntrusivePtr<GeoMaterial> hu{new GeoMaterial("special::HyperUranium",0.0)}; + GeoIntrusivePtr<GeoMaterial> hu{make_intrusive<GeoMaterial>("special::HyperUranium",0.0)}; hu->add(ethElement,1.); hu->lock(); m_materials.insert(std::make_pair("special::HyperUranium", hu));