Integrate Kinkal volumes and update Calo geometry

CommonMC/src/MakeSurfaceSteps_module.cc

@@ -60,6 +60,16 @@ namespace mu2e {
       vdmap_[VirtualDetectorId(VirtualDetectorId::TT_Back)] = SurfaceId("TT_Back");
       vdmap_[VirtualDetectorId(VirtualDetectorId::TT_OutSurf)] = SurfaceId("TT_Outer");
       vdmap_[VirtualDetectorId(VirtualDetectorId::TT_InSurf)] = SurfaceId("TT_Inner");
+
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_0_SurfIn)] = SurfaceId("EMC_Disk_0_Front");
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_1_SurfIn)] = SurfaceId("EMC_Disk_1_Front");
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_0_SurfOut)] = SurfaceId("EMC_Disk_0_Back");
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_1_SurfOut)] = SurfaceId("EMC_Disk_1_Back");
+
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_0_EdgeIn)] = SurfaceId("EMC_Disk_0_Inner");
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_1_EdgeIn)] = SurfaceId("EMC_Disk_1_Inner");
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_0_EdgeOut)] = SurfaceId("EMC_Disk_0_Outer");
+      vdmap_[VirtualDetectorId(VirtualDetectorId::EMC_Disk_1_EdgeOut)] = SurfaceId("EMC_Disk_1_Outer");
     }
 
   void MakeSurfaceSteps::produce(art::Event& event) {
@@ -72,6 +82,7 @@ namespace mu2e {
     for(auto const& vdspmc : vdspmccol) {
       // only some VDs are kept
       auto isid = vdmap_.find(vdspmc.virtualDetectorId());
+      if(debug_ > 0) std::cout<<" VID "<<vdspmc.virtualDetectorId().name()<<std::endl;
       if(isid != vdmap_.end())ssc->emplace_back(isid->second,vdspmc,det); // no aggregation of VD hits
     }
     auto nvdsteps = ssc->size();

DataProducts/inc/SurfaceId.hh

@@ -36,10 +36,12 @@ namespace mu2e {
         CRV_C3=262, CRV_C4,                                   // legacy MDC2020/crv_counters_v07-v09; retire once those datasets are no longer used
         CRV_M1=270, CRV_M2, CRV_M3, CRV_M4, CRV_M5, CRV_M6, CRV_M7, CRV_M8, // CRV-Muon-Taggers (Mu2e/Offline PR #1864)
         CRV_StrongBack=280,                                   // CRV module Al strongback (tracker-side support plate)
-        DS_HatchConcrete=300                                  // detector-area hatch concrete block approximation
+        DS_HatchConcrete=300,                                 // detector-area hatch concrete block approximation
+        EMC_Disk_0_Outer=310, EMC_Disk_0_Inner, EMC_Disk_1_Inner, EMC_Disk_1_Outer, // calorimeter disk bounding surfaces
+        EMC_Disk_0_Front, EMC_Disk_1_Front, EMC_Disk_0_Back, EMC_Disk_1_Back
       };
       // Update this counter whenever you add/remove surface IDs from the enum above.
-      static constexpr std::size_t nSurfaceIds = 64;
+      static constexpr std::size_t nSurfaceIds = 72;
       static std::string const& typeName();
       static std::map<enum_type,std::string> const& names();
   };

DataProducts/src/SurfaceId.cc

@@ -80,7 +80,15 @@ namespace mu2e {
     std::make_pair(SurfaceIdEnum::CRV_M7, "CRV_M7"),
     std::make_pair(SurfaceIdEnum::CRV_M8, "CRV_M8"),
     std::make_pair(SurfaceIdEnum::CRV_StrongBack, "CRV_StrongBack"),
-    std::make_pair(SurfaceIdEnum::DS_HatchConcrete, "DS_HatchConcrete")
+    std::make_pair(SurfaceIdEnum::DS_HatchConcrete, "DS_HatchConcrete"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_0_Outer, "EMC_Disk_0_Outer"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_0_Inner, "EMC_Disk_0_Inner"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_1_Inner, "EMC_Disk_1_Inner"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_1_Outer, "EMC_Disk_1_Outer"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_0_Front, "EMC_Disk_0_Front"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_1_Front, "EMC_Disk_1_Front"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_0_Back, "EMC_Disk_0_Back"),
+    std::make_pair(SurfaceIdEnum::EMC_Disk_1_Back, "EMC_Disk_1_Back")
     };
 
     constexpr std::size_t nSurfaceIdNames = sizeof(surfaceIdNames)/sizeof(surfaceIdNames[0]);

GeometryService/inc/KinKalGeomMaker.hh

@@ -18,6 +18,7 @@ namespace mu2e {
       void makeTracker();
       void makeDS();
       void makeTarget();
+      void makeCalo();
       void makeCRV();
       void makePassiveMaterials();
       // append one averaged rectangular concrete passive-material plane (see .cc)

GeometryService/src/KinKalGeomMaker.cc

@@ -16,6 +16,8 @@
 #include "Offline/KinKalGeom/inc/DetectorSolenoid.hh"
 #include "Offline/KinKalGeom/inc/StoppingTarget.hh"
 #include "Offline/KinKalGeom/inc/CRV.hh"
+#include "Offline/KinKalGeom/inc/Calo.hh"
+#include "Offline/CalorimeterGeom/inc/Calorimeter.hh"
 #include "Offline/BeamlineGeom/inc/Beamline.hh"
 #include "Offline/GeometryService/inc/DetectorSystem.hh"
 #include "Offline/DetectorSolenoidGeom/inc/DetectorSolenoid.hh"
@@ -32,7 +34,6 @@
 namespace mu2e {
   using KinKal::VEC3;
   using KinKal::Cylinder;
-  using KinKal::Disk;
   using KinKal::Surface;
   using KinKal::Frustrum;
   using KinKal::Annulus;
@@ -54,9 +55,49 @@ namespace mu2e {
     makeTarget();
     makeCRV();
     makePassiveMaterials();
+    makeCalo();
     return kkg_;
   }
 
+  void KinKalGeomMaker::makeCalo() {
+    // construct calorimeter geometry from GeometryService
+    art::ServiceHandle<GeometryService> geomService;
+    if(!geomService->hasElement<Calorimeter>()) return;
+    auto const& calo_det = *(GeomHandle<mu2e::Calorimeter>());
+    GeomHandle<DetectorSystem> det;
+    auto toDetectorZ = [&det](double zmu2e) {
+      return VEC3(det->toDetector(CLHEP::Hep3Vector(0.0,0.0,zmu2e))).Z();
+    };
+
+    // Extract geometry from the first two disks (D0 and D1)
+    auto const& geom0 = calo_det.disk(0).diskInfo();
+    auto const& geom1 = calo_det.disk(1).diskInfo();
+    double r0_inner = geom0.innerEnvelopeR();
+    double r0_outer = geom0.outerEnvelopeR();
+    double r1_inner = geom1.innerEnvelopeR();
+    double r1_outer = geom1.outerEnvelopeR();
+
+    // disk front/back faces, converted from Mu2e to detector coordinates
+    double z0_front = toDetectorZ(geom0.origin().z() - 0.5*geom0.size().z());
+    double z0_back  = toDetectorZ(geom0.origin().z() + 0.5*geom0.size().z());
+    double z1_front = toDetectorZ(geom1.origin().z() - 0.5*geom1.size().z());
+    double z1_back  = toDetectorZ(geom1.origin().z() + 0.5*geom1.size().z());
+
+    // Construct Calo with geometry in detector coordinates
+    kkg_->calo_ = std::make_unique<KKGeom::Calo>(r0_inner, r0_outer, r1_inner, r1_outer, z0_front, z0_back, z1_front, z1_back);
+    auto const& calo = *kkg_->calo_;
+
+    // Add all calo surfaces to the map
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_0_Outer),std::static_pointer_cast<Surface>(calo.EMC_Disk_0_OuterPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_0_Inner),std::static_pointer_cast<Surface>(calo.EMC_Disk_0_InnerPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_1_Inner),std::static_pointer_cast<Surface>(calo.EMC_Disk_1_InnerPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_1_Outer),std::static_pointer_cast<Surface>(calo.EMC_Disk_1_OuterPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_0_Front),std::static_pointer_cast<Surface>(calo.EMC_Disk_0_FrontPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_1_Front),std::static_pointer_cast<Surface>(calo.EMC_Disk_1_FrontPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_0_Back),std::static_pointer_cast<Surface>(calo.EMC_Disk_0_BackPtr())));
+    kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::EMC_Disk_1_Back),std::static_pointer_cast<Surface>(calo.EMC_Disk_1_BackPtr())));
+  }
+
   // sort by transverse distance
   struct sortCRVSectors {
     bool operator () (KKCRVSector const& sect1, KKCRVSector const& sect2) {
@@ -89,9 +130,9 @@ namespace mu2e {
     auto outer = std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,zMidLocal),orvd,halfLen);
     auto inner = std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,zMidLocal),irvd,halfLen);
     // expand the disk radii to the DS
-    auto front = std::make_shared<Disk>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,zFrontLocal),irds);
-    auto mid = std::make_shared<Disk>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,zMidLocal),irds);
-    auto back = std::make_shared<Disk>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,zBackLocal),irds);
+    auto front = std::make_shared<KinKal::Disk>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,zFrontLocal),irds);
+    auto mid = std::make_shared<KinKal::Disk>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,zMidLocal),irds);
+    auto back = std::make_shared<KinKal::Disk>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,zBackLocal),irds);
     // add all these to the map
     kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::TT_Front),std::static_pointer_cast<Surface>(front)));
     kkg_->map_.emplace(std::make_pair(SurfaceId(SurfaceIdEnum::TT_Mid),std::static_pointer_cast<Surface>(mid)));
@@ -133,8 +174,8 @@ namespace mu2e {
     //
     auto inner= std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,-1482),ds->rIn1(),ds->halfLength());
     auto outer= std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,-1482),ds->rOut2(),ds->halfLength()); // bounding surfaces
-    auto front= std::make_shared<Disk>(outer->frontDisk());
-    auto back= std::make_shared<Disk>(outer->backDisk());
+    auto front= std::make_shared<KinKal::Disk>(outer->frontDisk());
+    auto back= std::make_shared<KinKal::Disk>(outer->backDisk());
     KKGeom::DetectorSolenoid::MaterialCylinderCollection materialCylinders;
     auto toDetectorZ = [&det](double zmu2e) {
       return VEC3(det->toDetector(CLHEP::Hep3Vector(0.0,0.0,zmu2e))).Z();
@@ -203,8 +244,8 @@ namespace mu2e {
 
     // hard-coded for now
     auto ipa= std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,-2770),300.0,500.0);
-    auto ipafront= std::make_shared<Disk>(ipa->frontDisk());
-    auto ipaback= std::make_shared<Disk>(ipa->backDisk());
+    auto ipafront= std::make_shared<KinKal::Disk>(ipa->frontDisk());
+    auto ipaback= std::make_shared<KinKal::Disk>(ipa->backDisk());
     auto opa= std::make_shared<Frustrum>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,-3766),454.0,728.4,2125.0); // inner surface
     auto tsda= std::make_shared<Annulus>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,-5967),235.0,525.0); // back surface
 
@@ -225,8 +266,8 @@ namespace mu2e {
     // currently use hard-coded geometry
     auto outer = std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,-4300),75,400.0);
     auto inner= std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,-4300),21.5,400.0);
-    auto front= std::make_shared<Disk>(outer->frontDisk());
-    auto back= std::make_shared<Disk>(outer->backDisk());
+    auto front= std::make_shared<KinKal::Disk>(outer->frontDisk());
+    auto back= std::make_shared<KinKal::Disk>(outer->backDisk());
     double startz = -4700;
     double endz = -3900;
     double dz = (endz-startz)/36.0;

KinKalGeom/CMakeLists.txt

@@ -1,5 +1,6 @@
 cet_make_library(
     SOURCE
+      src/Calo.cc
       src/KinKalGeom.cc
       src/KKMaterial.cc
       src/KKStrawMaterial.cc

KinKalGeom/inc/Calo.hh

@@ -0,0 +1,48 @@
+//
+//  Define the nominal calorimeter boundary and reference surfaces, used to extrapolate and sample KinKal track fits, and to build
+//  the passive materials in the fit
+//  original author: Sophie Middleton (2025)
+//
+#ifndef KinKalGeom_Calo_hh
+#define KinKalGeom_Calo_hh
+#include "KinKal/Geometry/Cylinder.hh"
+#include "KinKal/Geometry/Annulus.hh"
+#include <memory>
+namespace mu2e {
+  namespace KKGeom {
+    class Calo {
+      public:
+        using CylPtr = std::shared_ptr<KinKal::Cylinder>;
+        using AnnPtr = std::shared_ptr<KinKal::Annulus>;
+        // constructor with geometry parameters from Calorimeter service
+        Calo(double r0_inner, double r0_outer, double r1_inner, double r1_outer,
+             double z0_front, double z0_back, double z1_front, double z1_back);
+        // return by reference
+        auto const& EMC_Disk_0_Outer() const { return *EMC_Disk_0_Outer_;}
+        auto const& EMC_Disk_0_Inner() const { return *EMC_Disk_0_Inner_;}
+        auto const& EMC_Disk_1_Inner() const { return *EMC_Disk_1_Inner_;}
+        auto const& EMC_Disk_1_Outer() const { return *EMC_Disk_1_Outer_;}
+
+        auto const& EMC_Disk_0_Front() const { return *EMC_Disk_0_Front_;}
+        auto const& EMC_Disk_1_Front() const { return *EMC_Disk_1_Front_;}
+        auto const& EMC_Disk_0_Back() const  { return *EMC_Disk_0_Back_;}
+        auto const& EMC_Disk_1_Back() const  { return *EMC_Disk_1_Back_;}
+
+        // return by ptr
+        auto const& EMC_Disk_0_OuterPtr() const { return EMC_Disk_0_Outer_;}
+        auto const& EMC_Disk_0_InnerPtr() const { return EMC_Disk_0_Inner_;}
+        auto const& EMC_Disk_1_InnerPtr() const { return EMC_Disk_1_Inner_;}
+        auto const& EMC_Disk_1_OuterPtr() const { return EMC_Disk_1_Outer_;}
+        auto const& EMC_Disk_0_FrontPtr() const { return EMC_Disk_0_Front_;}
+        auto const& EMC_Disk_1_FrontPtr() const { return EMC_Disk_1_Front_;}
+        auto const& EMC_Disk_0_BackPtr() const  { return EMC_Disk_0_Back_;}
+        auto const& EMC_Disk_1_BackPtr() const  { return EMC_Disk_1_Back_;}
+
+      private:
+        CylPtr EMC_Disk_0_Inner_, EMC_Disk_0_Outer_ , EMC_Disk_1_Inner_, EMC_Disk_1_Outer_;
+        AnnPtr EMC_Disk_0_Front_, EMC_Disk_0_Back_, EMC_Disk_1_Front_,  EMC_Disk_1_Back_; // annular disk faces
+    };
+  }
+}
+
+#endif

KinKalGeom/inc/KinKalGeom.hh

@@ -9,6 +9,7 @@
 #include "Offline/KinKalGeom/inc/StoppingTarget.hh"
 #include "Offline/KinKalGeom/inc/DetectorSolenoid.hh"
 #include "Offline/KinKalGeom/inc/CRV.hh"
+#include "Offline/KinKalGeom/inc/Calo.hh"
 #include "Offline/DataProducts/inc/SurfaceId.hh"
 #include "KinKal/Geometry/Surface.hh"
 #include "KinKal/Geometry/Plane.hh"
@@ -49,13 +50,15 @@ namespace mu2e {
       auto const& DS() const {return ds_; }
       auto const& ST() const {return st_; }
       auto const& tracker() const {return tracker_; }
+      auto const& calo() const {return calo_; }
       auto const& CRV() const {return crv_; }
       auto const& passiveMaterialPlanes() const { return passiveMaterialPlanes_; }
     private:
       // local copy of detector objects; these hold the actual (typed) surface objects
       std::unique_ptr<KKGeom::Tracker> tracker_;
       std::unique_ptr<KKGeom::DetectorSolenoid> ds_;
       std::unique_ptr<KKGeom::StoppingTarget> st_;
+      std::unique_ptr<KKGeom::Calo> calo_;
       std::unique_ptr<KKGeom::CRV> crv_;
       PassiveMaterialPlaneCollection passiveMaterialPlanes_;
       // the map used to find surfaces by Id

KinKalGeom/src/Calo.cc

@@ -0,0 +1,21 @@
+#include "Offline/KinKalGeom/inc/Calo.hh"
+namespace mu2e {
+  namespace KKGeom {
+    using KinKal::VEC3;
+    using KinKal::Cylinder;
+    using KinKal::Annulus;
+
+    Calo::Calo(double r0_inner, double r0_outer, double r1_inner, double r1_outer,
+               double z0_front, double z0_back, double z1_front, double z1_back) :
+      EMC_Disk_0_Inner_ { std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,0.5*(z0_front+z0_back)),r0_inner, 0.5*(z0_back-z0_front))},
+      EMC_Disk_0_Outer_ { std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,0.5*(z0_front+z0_back)),r0_outer, 0.5*(z0_back-z0_front))},
+      EMC_Disk_1_Inner_ { std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,0.5*(z1_front+z1_back)),r1_inner, 0.5*(z1_back-z1_front))},
+      EMC_Disk_1_Outer_ { std::make_shared<Cylinder>(VEC3(0.0,0.0,1.0),VEC3(0.0,0.0,0.5*(z1_front+z1_back)),r1_outer, 0.5*(z1_back-z1_front))},
+
+      EMC_Disk_0_Front_ { std::make_shared<Annulus>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,z0_front),r0_inner,r0_outer)},
+      EMC_Disk_0_Back_ { std::make_shared<Annulus>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,z0_back),r0_inner,r0_outer)},
+      EMC_Disk_1_Front_ { std::make_shared<Annulus>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,z1_front),r1_inner,r1_outer)},
+      EMC_Disk_1_Back_ { std::make_shared<Annulus>(VEC3(0.0,0.0,1.0),VEC3(1.0,0.0,0.0),VEC3(0.0,0.0,z1_back),r1_inner,r1_outer)}
+    {}
+  }
+}

Mu2eKinKal/fcl/prolog.fcl

@@ -390,6 +390,8 @@ Mu2eKinKal : {
     BackToTracker : false
     ToOPA : true
     ToCRV : false
+    ToCaloD0 : true
+    ToCaloD1 : true
   }
 
   LHSEEDXTRAP : {
@@ -403,6 +405,8 @@ Mu2eKinKal : {
     BackToTracker : false
     ToOPA : false
     ToCRV : false
+    ToCaloD0 : false
+    ToCaloD1 : false
   }
 
   KLXTRAP : {
@@ -416,6 +420,8 @@ Mu2eKinKal : {
     BackToTracker : false
     ToOPA : false
     ToCRV : true
+    ToCaloD0 : true
+    ToCaloD1 : true
   }
 
   CHXTRAP : {
@@ -429,6 +435,8 @@ Mu2eKinKal : {
     BackToTracker : false
     ToOPA : false
     ToCRV : true
+    ToCaloD0 : true
+    ToCaloD1 : true
   }
 
   CENTRALHELIX : {