DumpGeo.cxx

Go to the documentation of this file.

/*
  Copyright (C) 2002-2026 CERN for the benefit of the ATLAS collaboration
*/
 
//  Initial version:
//  - 2017, Sep -- Riccardo Maria BIANCHI
//                 <riccardo.maria.bianchi@cern.ch>
//  Main updates:
//  - 2024, Feb -- Riccardo Maria BIANCHI
//                 <riccardo.maria.bianchi@cern.ch>
//                 Moved configuration to ComponentAccumulator (CA)
//  - 2024, Mar -- Riccardo Maria BIANCHI
//                 <riccardo.maria.bianchi@cern.ch>
//                 Removed GeoExporter, moved all to DumpGeo
//  - 2024, Apr -- Added missing sub-systems: Cavern, HGTD, FWD detectors, ...
 
#include "DumpGeo/DumpGeo.h"
 
// Athena includes
#include "AthenaBaseComps/AthMsgStreamMacros.h"
#include "GeoModelInterfaces/IGeoDbTagSvc.h"
#include "GeoModelUtilities/GeoModelExperiment.h"
 
// From GeoModel in Externals
#include "GeoModelDBManager/GMDBManager.h"
#include "GeoModelKernel/GeoBox.h"
#include "GeoModelKernel/GeoVDetectorManager.h"
#include "GeoModelKernel/GeoVolumeCursor.h"
#include "GeoModelWrite/WriteGeoModel.h"
// #include "GeoModelHelpers/defineWorld.h" // not available in 24.0...
 
// C++ includes
#include <cstdlib>  //For setenv
#include <filesystem>
#include <iomanip>
#include <iostream>
#include <stdexcept>
#include <system_error>
#include <vector>
 
// Units
#include "GeoModelKernel/Units.h"
#define UNITS GeoModelKernelUnits  // so we can use, e.g., 'UNITS::cm'
 
// TODO: replace this with GeoModelHelpers/defineWorld.h
//_____________________________________________________________________________________
GeoPhysVol* createTheWorld() {
// Define the units
#define gr UNITS::gram
#define mole UNITS::mole
#define cm3 UNITS::cm3
 
  // Define the chemical elements
  GeoElement* Nitrogen =
      new GeoElement("Nitrogen", "N", 7.0, 14.0067 * gr / mole);
  GeoElement* Oxygen = new GeoElement("Oxygen", "O", 8.0, 15.9995 * gr / mole);
  GeoElement* Argon = new GeoElement("Argon", "Ar", 18.0, 39.948 * gr / mole);
  GeoElement* Hydrogen =
      new GeoElement("Hydrogen", "H", 1.0, 1.00797 * gr / mole);
 
  // Define the materials
  double densityOfAir = 0.001214 * gr / cm3;
  GeoMaterial* air = new GeoMaterial("Air", densityOfAir);
  air->add(Nitrogen, 0.7494);
  air->add(Oxygen, 0.2369);
  air->add(Argon, 0.0129);
  air->add(Hydrogen, 0.0008);
  air->lock();
 
  //-----------------------------------------------------------------------------------//
  // create the world volume container and
  // get the 'world' volume, i.e. the root volume of the GeoModel tree
  const GeoBox* worldBox =
      new GeoBox(1000 * UNITS::cm, 1000 * UNITS::cm, 1000 * UNITS::cm);
  const GeoLogVol* worldLog = new GeoLogVol("WorldLog", worldBox, air);
  GeoPhysVol* world = new GeoPhysVol(worldLog);
  return world;
}
 
//____________________________________________________________________
DumpGeo::DumpGeo(const std::string& name, ISvcLocator* svcLocator)
    : AthAlgorithm(name, svcLocator) {}
 
//____________________________________________________________________
StatusCode DumpGeo::initialize() {
  ATH_MSG_DEBUG("in initialize()");
 
  ATH_MSG_INFO("===================================================");
  ATH_MSG_INFO("\t\tLaunching DumpGeo");
  ATH_MSG_INFO("===================================================");
 
  ATH_MSG_INFO("Accessing the ATLAS geometry...");
  const GeoModelExperiment* theExpt = nullptr;
  ATH_CHECK(detStore()->retrieve(theExpt, "ATLAS"));
 
  // Get the ATLAS GEOMETRY 'World' volume
  PVConstLink world(theExpt->getPhysVol());
 
  ATH_MSG_INFO(
      "User filter DetectorManagers: " << m_user_filterDetManagersList);
 
  ServiceHandle<IGeoDbTagSvc> geoDbTag("GeoDbTagSvc", name());
  ATH_CHECK(geoDbTag.retrieve());
  ATH_MSG_INFO(
      "This is the Geometry TAG we are dumping: " << geoDbTag->atlasVersion());
 
  GeoPhysVol* volTop = createTheWorld();
 
  // Get list of managers
  // We fill a set from the output vector,
  // so we can use its built-in 'count' method later,
  // to search for DetManagers
  ATH_MSG_INFO("Dumping the GeoModel tree...");
  std::set<std::string> managersList{};
  {
    std::vector<std::string> blub = theExpt->getListOfManagers();
    managersList.insert(blub.begin(), blub.end());
  }
  if (msgLvl(MSG::INFO)) {
    ATH_MSG_INFO("List of the GeoModel Detector Managers that are being dumped: ");
    for (auto const& man : managersList) {
      // get the DetectorManager
      const GeoVDetectorManager* manager = theExpt->getManager(man);
      unsigned int nTreetops = manager->getNumTreeTops();
      std::cout << "\t" << man << " [contains " << nTreetops << " treetops]"
                << std::endl;
      for (unsigned int i = 0; i < nTreetops; ++i) {
        PVConstLink treetop(manager->getTreeTop(i));
        // get treetop's volume
        const GeoVPhysVol* vol = treetop;
        // get treetop's logvol's name
        std::string volName = vol->getLogVol()->getName();
        std::cout << "\t\t treetop: " << volName << std::endl;
      }
    }
  }
 
  if (!(m_user_filterDetManagersList.empty())) {
    // Convert the list of det managers passed by the user into a set
    std::set<std::string> user_managersList{};
    {
      user_managersList.insert(m_user_filterDetManagersList.begin(),
                               m_user_filterDetManagersList.end());
    }
    // safety check:
    // check that all DetManagers requested by the user are in the list
    // If not, print an error message to warn the user and return
    for (auto& userDet : user_managersList) {
      if (!managersList.count(userDet)) {
        ATH_MSG_FATAL(
            "This Detector Manager you requested to dump is not in the list of "
            "DetectorManagers for the geometry tag you are using: "
            << userDet);
        throw GaudiException(
            "The Detector Manager you requested to dump is not in the list of "
            "DetectorManagers.",
            "DumpGeo", StatusCode::FAILURE);
      }
    }
 
    if (!(managersList.empty())) {
      for (auto const& mm : managersList) {
        // get the DetectorManager
        const GeoVDetectorManager* manager = theExpt->getManager(mm);
 
        // get the name of the DetectorManager
        std::string detManName = manager->getName();
        ATH_MSG_INFO("\tDetectorManager: " << detManName);
 
        // get the DetManager's TreeTops
        unsigned int nTreetops = manager->getNumTreeTops();
        ATH_MSG_INFO("\t" << mm << " - # TreeTops: " << nTreetops);
 
        if ((nTreetops > 0) && user_managersList.count(detManName)) {
 
          for (unsigned int i = 0; i < nTreetops; ++i) {
 
            PVConstLink treetop(manager->getTreeTop(i));
 
            // get treetop's volume
            const GeoVPhysVol* vol = treetop;
 
            // get volume's transform
            // NOTE: we use getDefX() to get the transform without any alignment
            GeoTransform* volXf = new GeoTransform(vol->getDefX());
 
            // get volume's logvol's name
            std::string volName = vol->getLogVol()->getName();
            ATH_MSG_INFO("\t* treetop " << i << ": " << volName);
 
            // Add to the main volume a GeoNameTag with the name of the
            // DetectorManager
            volTop->add(new GeoNameTag(detManName));
            // add Transform and Volume to the main PhysVol
            volTop->add(volXf);
            volTop->add(const_cast<GeoVPhysVol*>(vol));
 
            if (m_showTreetopContent) {
                ATH_MSG_INFO("You enabled the option 'showTreetopContent', so we now print the content of the first layers of child volumes of the Treetops...");
                unsigned v1{0};
                GeoVolumeCursor av(treetop);
                while (!av.atEnd()) {
                    ++v1;
                    if(1==v1) ATH_MSG_INFO("first level of child volumes:");
                    ATH_MSG_INFO("\t\t- child's name: " << av.getName());
                    ATH_MSG_INFO("\t\t- child's n. sub-nodes: " << av.getVolume()->getNChildNodes());
                    ATH_MSG_INFO("\t\t- child's n. sub-volumes: " << av.getVolume()->getNChildVols());
 
 
                    unsigned v2{0};
                    GeoVolumeCursor av2(av.getVolume());
                    while (!av2.atEnd()) {
                        ++v2;
                        if(1==v2) ATH_MSG_INFO("second level of child volumes:");
                        ATH_MSG_INFO("\t\t\t- child's logvol's name: " << av2.getVolume()->getLogVol()->getName());
 
                        unsigned v3{0};
                        GeoVolumeCursor av3(av2.getVolume());
                        while (!av3.atEnd()) {
                            ++v3;
                            if(1==v3) ATH_MSG_INFO("third level of child volumes:");
                            ATH_MSG_INFO("\t\t\t\t- child's logvol's name: " << av3.getVolume()->getLogVol()->getName());
                            av3.next();  // increment volume cursor.
                        }             // end while
 
                        av2.next();  // increment volume cursor.
                    }             // end while
 
                    av.next();  // increment volume cursor.
                }             // end while
            }
          }
        }
      }
    }
  }
 
  // DEBUG inspection
  if (msgLvl(MSG::DEBUG)) {
    ATH_MSG_DEBUG(
        "Looping over top volumes in the GeoModel tree (children of the "
        "'World' volume)...");
    GeoVolumeCursor av(world);
    while (!av.atEnd()) {
      std::string volname = av.getName();
      ATH_MSG_DEBUG("\t* relevant NameTag:" << volname);
      av.next();  // increment volume cursor.
    }
  }
 
  ATH_MSG_INFO("Creating the SQLite DB file...");
  if (m_outFileName.empty()) {
    ATH_MSG_FATAL("The name of the output SQLite file is not set!");
    throw GaudiException("The name of the output SQLite file is not set!",
                         "DumpGeo", StatusCode::FAILURE);
  }
  ATH_MSG_INFO("Output file name: " << m_outFileName);
 
  // Delay destructive overwrite handling until the algorithm is actually
  // initialized, immediately before opening the output database. Use one
  // removal attempt without a preceding existence check: a missing file is
  // harmless, while any real filesystem error must stop the job.
  if (m_forceOverwrite) {
    std::error_code error;
    const bool removed = std::filesystem::remove(m_outFileName.value(), error);
    if (error) {
      ATH_MSG_FATAL("Could not remove existing output file '"
                    << m_outFileName << "': " << error.message());
      return StatusCode::FAILURE;
    }
    if (removed) {
      ATH_MSG_INFO("Removed existing output file '" << m_outFileName << "'.");
    }
  }
 
  // open the DB connection
  GMDBManager db(m_outFileName);
 
  // check the DB connection
  if (db.checkIsDBOpen())
    ATH_MSG_INFO("OK! Database is open!");
  else {
    ATH_MSG_ERROR(" ***** Database ERROR!! Exiting...");
    throw GaudiException(
        "The GeoModel SQLite .db file could not be opened successfully.",
        "DumpGeo", StatusCode::FAILURE);
  }
 
  ATH_MSG_INFO("Traversing the GeoModel tree...");
  // Dump the tree volumes into a DB
  // 1. init the GeoModel node action
  GeoModelIO::WriteGeoModel dumpGeoModelGraph(db);  
  // 2. visit all GeoModel nodes
  if (!(m_user_filterDetManagersList.empty())) {
    volTop->exec(&dumpGeoModelGraph);
  } else {
    world->exec(&dumpGeoModelGraph);
  }
  ATH_MSG_INFO("Saving the GeoModel tree to the DB...");
  // 3. save it to an SQLite file
  dumpGeoModelGraph.saveToDB(); 
  ATH_MSG_ALWAYS("DONE. Geometry saved to " << m_outFileName);
 
  // Quick test if DEBUG
  if (msgLvl(MSG::DEBUG)) {
    ATH_MSG_DEBUG(
        "Test - list of all the GeoMaterial nodes in the persistified "
        "geometry:");
    db.printAllMaterials();
    ATH_MSG_DEBUG(
        "Test - list of all the GeoElement nodes in the persistified "
        "geometry:");
    db.printAllElements();
  }
 
  ATH_MSG_DEBUG("End of DumpGeo::init().");
  return StatusCode::SUCCESS;
}