GeometryCore.h

Go to the documentation of this file.

#ifndef LARCOREALG_GEOMETRY_GEOMETRYCORE_H
#define LARCOREALG_GEOMETRY_GEOMETRYCORE_H

// LArSoft libraries
#include "larcorealg/Geometry/AuxDetGeo.h"
#include "larcorealg/Geometry/AuxDetSensitiveGeo.h"
#include "larcorealg/Geometry/BoxBoundedGeo.h"
#include "larcorealg/Geometry/CryostatGeo.h"
#include "larcorealg/Geometry/GeoObjectSorter.h"
#include "larcorealg/Geometry/GeometryBuilder.h"
#include "larcorealg/Geometry/Iterable.h"
#include "larcorealg/Geometry/OpDetGeo.h"
#include "larcorealg/Geometry/TPCGeo.h"
#include "larcorealg/Geometry/details/GeometryIterationPolicy.h"
#include "larcorealg/Geometry/details/ToGeometryElement.h"
#include "larcorealg/Geometry/details/ZeroIDs.h"
#include "larcorealg/Geometry/fwd.h"
#include "larcorealg/Geometry/geo_vectors_utils.h"       // geo::vect namespace
#include "larcoreobj/SimpleTypesAndConstants/RawTypes.h" // raw::ChannelID_t
#include "larcoreobj/SimpleTypesAndConstants/geo_types.h"
#include "larcoreobj/SimpleTypesAndConstants/geo_vectors.h"

// Framework and infrastructure libraries
#include "fhiclcpp/fwd.h"

// C/C++ standard libraries
#include <cstddef> // size_t
#include <memory>  // std::shared_ptr<>
#include <set>
#include <string>
#include <utility>
#include <vector>

// ROOT class prototypes
class TGeoManager;
class TGeoNode;
class TGeoVolume;
class TGeoMaterial;

namespace geo {

  // BEGIN Geometry group ------------------------------------------------------

  class GeometryCore : Iterable<details::GeometryIterationPolicy, details::ToGeometryElement> {
    using Iteration = Iterable<details::GeometryIterationPolicy, details::ToGeometryElement>;

  public:
    using CryostatList_t = std::vector<CryostatGeo>;

    GeometryCore(fhicl::ParameterSet const& pset,
                 std::unique_ptr<GeometryBuilder> builder,
                 std::unique_ptr<GeoObjectSorter> sorter);
    GeometryCore(GeometryCore const&) = delete;
    GeometryCore(GeometryCore&&) = delete;
    GeometryCore& operator=(GeometryCore const&) = delete;
    GeometryCore& operator=(GeometryCore&&) = delete;

    double DefaultWiggle() const { return fPositionWiggle; }

    std::string const& GDMLFile() const { return fGDMLfile; }

    // BEGIN Detector information

    //
    // global features
    //
    std::string const& DetectorName() const { return fDetectorName; }

    //
    // position
    //

    TGeoVolume const* WorldVolume() const;

    void WorldBox(double* xlo, double* xhi, double* ylo, double* yhi, double* zlo, double* zhi)
      const;

    BoxBoundedGeo WorldBox() const;

    //
    Length_t SurfaceY() const { return fSurfaceY; }

    //
    // object description and information
    //

    TGeoManager* ROOTGeoManager() const;

    std::string const& GetWorldVolumeName() const;

    BoxBoundedGeo DetectorEnclosureBox(std::string const& name = "volDetEnclosure") const;


    std::string VolumeName(Point_t const& point) const;

    std::vector<TGeoNode const*> FindAllVolumes(std::set<std::string> const& vol_names) const;

    std::vector<std::vector<TGeoNode const*>> FindAllVolumePaths(
      std::set<std::string> const& vol_names) const;

    TGeoMaterial const* Material(Point_t const& point) const;

    std::string MaterialName(Point_t const& point) const;

    double TotalMass() const { return TotalMass(GetWorldVolumeName()); }
    double TotalMass(std::string vol) const;


    double MassBetweenPoints(Point_t const& p1, Point_t const& p2) const;

    template <typename Stream>
    void Print(Stream&& out, std::string indent = "  ") const;

    std::string Info(std::string indent = "  ") const;

    // END Detector information


    //
    // group features
    //


    unsigned int Ncryostats() const { return fCryostats.size(); }
    unsigned int NElements() const { return Ncryostats(); }
    unsigned int NSiblingElements(CryostatID const&) const { return Ncryostats(); }

    //
    // access
    //


    bool HasCryostat(CryostatID const& cryoid) const { return cryoid.Cryostat < Ncryostats(); }
    bool HasElement(CryostatID const& cryoid) const { return HasCryostat(cryoid); }


    CryostatGeo const& Cryostat(CryostatID const& cryoid = details::cryostat_zero) const;
    CryostatGeo const& GetElement(CryostatID const& cryoid) const { return Cryostat(cryoid); }


    CryostatGeo const* CryostatPtr(CryostatID const& cryoid) const
    {
      return HasCryostat(cryoid) ? &fCryostats[cryoid.Cryostat] : nullptr;
    }
    CryostatGeo const* GetElementPtr(CryostatID const& cryoid) const { return CryostatPtr(cryoid); }

    CryostatGeo const* PositionToCryostatPtr(Point_t const& point) const;

    CryostatID PositionToCryostatID(Point_t const& point) const;


    CryostatGeo const& PositionToCryostat(Point_t const& point) const;

    //
    // iterators
    //
    using Iteration::begin;
    using Iteration::end;
    using Iteration::Iterate;

    //
    // single object features
    //


    //
    // group features
    //

    unsigned int MaxTPCs() const;

    unsigned int TotalNTPC() const;


    unsigned int NTPC(CryostatID const& cryoid = details::cryostat_zero) const
    {
      CryostatGeo const* pCryo = GetElementPtr(cryoid);
      return pCryo ? pCryo->NElements() : 0;
    }
    unsigned int NElements(CryostatID const& cryoid) const { return NTPC(cryoid); }
    unsigned int NSiblingElements(TPCID const& tpcid) const { return NTPC(tpcid); }

    //
    // access
    //
    bool HasTPC(TPCID const& tpcid) const
    {
      CryostatGeo const* pCryo = CryostatPtr(tpcid);
      return pCryo ? pCryo->HasTPC(tpcid) : false;
    }

    bool HasElement(TPCID const& tpcid) const { return HasTPC(tpcid); }


    TPCGeo const& TPC(TPCID const& tpcid = details::tpc_zero) const
    {
      return Cryostat(tpcid).TPC(tpcid);
    }
    TPCGeo const& GetElement(TPCID const& tpcid) const { return TPC(tpcid); }


    TPCGeo const* TPCPtr(TPCID const& tpcid) const
    {
      CryostatGeo const* pCryo = CryostatPtr(tpcid);
      return pCryo ? pCryo->TPCPtr(tpcid) : nullptr;
    }
    TPCGeo const* GetElementPtr(TPCID const& tpcid) const { return TPCPtr(tpcid); }


    TPCID FindTPCAtPosition(Point_t const& point) const;

    TPCGeo const* PositionToTPCptr(Point_t const& point) const;


    TPCGeo const& PositionToTPC(Point_t const& point) const;

    TPCID PositionToTPCID(Point_t const& point) const;

    //
    // object description
    //


    //
    // group features
    //

    unsigned int NOpDets() const;

    //
    // access
    //
    OpDetGeo const& OpDetGeoFromOpDet(unsigned int OpDet) const;


    unsigned int GetClosestOpDet(Point_t const& point) const;

    //
    // object description
    //

    std::string OpDetGeoName(CryostatID const& cid = details::cryostat_zero) const;


    //
    // group features
    //

    //
    // group features
    //

    unsigned int OpDetFromCryo(unsigned int o, unsigned int c) const;


    //
    // unsorted methods
    //

    CryostatList_t const& Cryostats() const noexcept { return fCryostats; }

  private:
    void LoadGeometryFile();
    void SortGeometry();

    CryostatList_t fCryostats{};

    std::unique_ptr<GeometryBuilder> fBuilder;
    std::unique_ptr<GeoObjectSorter> fSorter;

    TGeoManager* fManager{nullptr};
    std::string fGDMLfile;     
    std::string fDetectorName; 
    double fSurfaceY;          
    double fPositionWiggle;    

    std::vector<GeoNodePathEntry> FindDetectorEnclosure(
      std::string const& name = "volDetEnclosure") const;

    bool FindFirstVolume(std::string const& name, std::vector<GeoNodePathEntry>& path) const;

    void BuildGeometry();

  }; // class GeometryCore

  // END Geometry group --------------------------------------------------------

} // namespace geo

//------------------------------------------------------------------------------
template <typename Stream>
void geo::GeometryCore::Print(Stream&& out, std::string indent /* = "  " */) const
{

  out << "Detector " << DetectorName() << " has " << Ncryostats() << " cryostats:";

  auto const& detEnclosureBox = DetectorEnclosureBox();
  out << "\n"
      << indent << "Detector enclosure: " << detEnclosureBox.Min() << " -- "
      << detEnclosureBox.Max() << " cm => ( " << detEnclosureBox.SizeX() << " x "
      << detEnclosureBox.SizeY() << " x " << detEnclosureBox.SizeZ() << " ) cm^3";

  for (auto const& cryostat : Iterate<CryostatGeo>()) {
    out << "\n" << indent;
    cryostat.PrintCryostatInfo(std::forward<Stream>(out), indent + "  ", cryostat.MaxVerbosity);

    unsigned const int nTPCs = cryostat.NTPC();
    for (unsigned int t = 0; t < nTPCs; ++t) {
      TPCGeo const& tpc = cryostat.TPC(t);

      out << "\n" << indent << "  ";
      tpc.PrintTPCInfo(std::forward<Stream>(out), indent + "    ", tpc.MaxVerbosity);
    } // for TPC

    unsigned int nOpDets = cryostat.NOpDet();
    for (unsigned int iOpDet = 0; iOpDet < nOpDets; ++iOpDet) {
      OpDetGeo const& opDet = cryostat.OpDet(iOpDet);
      out << "\n" << indent << "  [OpDet #" << iOpDet << "] ";
      opDet.PrintOpDetInfo(std::forward<Stream>(out), indent + "  ", opDet.MaxVerbosity);
    } // for
  }   // for cryostat

  out << '\n';

} // geo::GeometryCore::Print()

#endif // LARCOREALG_GEOMETRY_GEOMETRYCORE_H