Add functionality to track in a BeamPipeChain (#768)

### Briefly, what does this PR introduce?
An option had been added to include sensitive elements at the start
and/or end of beam pipe vacuum. This allows easy tracking of the shape
and position of beam electrons/protons to compare with accelerator group
and stand alone simulations. The segment in which scattered electrons of
different kinematics can also be studied/verified.

### What kind of change does this PR introduce?
- [ ] Bug fix (issue #__)
- [x] New feature (issue #__)
- [ ] Documentation update
- [ ] Other: __

### Please check if this PR fulfills the following:
- [ ] Tests for the changes have been added
- [ ] Documentation has been added / updated
- [ ] Changes have been communicated to collaborators

### Does this PR introduce breaking changes? What changes might users
need to make to their code?
No
### Does this PR change default behavior?
No

---------

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>

compact/definitions.xml

@@ -314,6 +314,7 @@ The unused IDs below are saved for future use.
     <constant name="LumiSpecCAL_ID"          value="194"/>
     <constant name="LumiDirectPCAL_ID"       value="195"/>
 
+    <constant name="BackwardsBeamline_ID"    value="197"/>
     <constant name="TaggerTracker_ID"        value="198"/>
     <constant name="TaggerCalorimeter_ID"    value="199"/>
 

compact/far_backward/beamline_extension_electron.xml

@@ -25,7 +25,8 @@
       <pipe id="2" name="Pipe_in_Q4eR"
         xcenter="Q4eR_CenterX" zcenter="Q4eR_CenterZ"
         length="Q4eR_Length" theta="Q4eR_Theta"
-        rout1="Q4eR_InnerRadius" rout2="Q4eR_InnerRadius">
+        rout1="Q4eR_InnerRadius" rout2="Q4eR_InnerRadius"
+        limits="kill_limits">
       </pipe>
       <pipe id="3" name="Pipe_Q4eR_to_B3eR"/>
       <pipe id="4" name="Pipe_in_B3eR"

compact/far_backward/beamline_tracking.xml

@@ -0,0 +1,72 @@
+
+<lccdd>
+
+  <detectors>
+    <detector
+    id="BackwardsBeamline_ID"
+    type="BeamPipeTracking"
+    name="Pipe_tracker"
+    readout="BackwardsBeamlineHits">
+
+        <slice
+        pipe_id="0"
+        grandmother="Pipe_Q1eR_to_B2BeR"
+        mother="Pipe_to_Q1eR_vacuum"
+        name="Pipe_to_Q1eR_tracker"/>
+
+        <slice
+        pipe_id="1"
+        grandmother="Pipe_Q1eR_to_B2BeR"
+        mother="Pipe_in_Q1eR_vacuum"
+        name="Pipe_in_Q1eR_tracker"/>
+
+        <slice
+        pipe_id="2"
+        grandmother="Pipe_Q1eR_to_B2BeR"
+        mother="Pipe_in_Q2eR_vacuum"
+        name="Pipe_in_Q2eR_tracker"/>
+
+        <slice
+        pipe_id="3"
+        grandmother="Pipe_Q1eR_to_B2BeR"
+        mother="Pipe_in_B2AeR_vacuum"
+        name="Pipe_in_B2AeR_tracker"/>
+
+        <slice
+        pipe_id="4"
+        grandmother="Pipe_Q1eR_to_B2BeR"
+        mother="Pipe_in_B2BeR_vacuum"
+        name="Pipe_in_B2BeR_tracker"/>
+
+        <slice
+        pipe_id="5"
+        grandmother="Pipe_Q3eR_to_B7eR"
+        mother="Pipe_in_Q3eR_vacuum"
+        name="Pipe_in_Q3eR_tracker_start"
+        end="false"/>
+
+        <slice
+        pipe_id="6"
+        grandmother="Pipe_Q3eR_to_B7eR"
+        mother="Pipe_in_Q3eR_vacuum"
+        name="Pipe_in_Q3eR_tracker_end"/>
+
+        <slice
+        pipe_id="7"
+        grandmother="Pipe_Q3eR_to_B7eR"
+        mother="Pipe_in_Q4eR_vacuum"
+        name="Pipe_in_Q4eR_tracker"/>
+
+    </detector>
+  </detectors>
+
+  <readouts>
+
+    <readout name="BackwardsBeamlineHits">
+      <segmentation type="NoSegmentation"/>
+      <id>system:8,pipe:8,end:2</id>
+    </readout>
+
+  </readouts>
+
+</lccdd>

compact/far_backward/extended.xml

@@ -13,4 +13,6 @@
   <!-- Luminosity detector -->
   <include ref="lumi.xml"/>
 
+  <include ref="beamline_tracking.xml"/>
+
 </lccdd>

src/BeamPipeChain_geo.cpp

@@ -31,19 +31,22 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector /
   double thickness  = getAttrOrDefault<double>(x_det, _Unicode(wall_thickness), 0);
 
   vector<string> names;
+  vector<int> ids;
   vector<double> xCenters;
   vector<double> zCenters;
   vector<double> lengths;
   vector<double> thetas;
   vector<double> rOuters1;
   vector<double> rOuters2;
+  vector<string> limits;
 
   // Grab info for beamline magnets
   for (xml_coll_t pipe_coll(x_det, _Unicode(pipe)); pipe_coll; pipe_coll++) { // pipes
 
     xml_comp_t pipe(pipe_coll);
 
     names.push_back(getAttrOrDefault<string>(pipe, _Unicode(name), ""));
+    ids.push_back(getAttrOrDefault<int>(pipe, _Unicode(id), 0));
 
     // Vectors momentarily filled with zeros for pipes in between magnets
     xCenters.push_back(getAttrOrDefault<double>(pipe, _Unicode(xcenter), 0));
@@ -52,6 +55,7 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector /
     thetas.push_back(getAttrOrDefault<double>(pipe, _Unicode(theta), 0));
     rOuters1.push_back(getAttrOrDefault<double>(pipe, _Unicode(rout1), 0));
     rOuters2.push_back(getAttrOrDefault<double>(pipe, _Unicode(rout2), 0));
+    limits.push_back(getAttrOrDefault<std::string>(pipe, _Unicode(limits), "world_limits"));
   }
 
   // Calculate parameters for connecting pipes in between magnets
@@ -105,11 +109,16 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector /
     Volume v_vacuum("v_vacuum_" + names[pipeN], s_vacuum, m_Vacuum);
 
     v_tube.setVisAttributes(description.visAttributes(vis_name));
+    v_vacuum.setLimitSet(description, limits[pipeN]);
 
     assembly.placeVolume(v_tube, Transform3D(RotationY(thetas[pipeN]),
                                              Position(xCenters[pipeN], 0, zCenters[pipeN])));
-    assembly.placeVolume(v_vacuum, Transform3D(RotationY(thetas[pipeN]),
-                                               Position(xCenters[pipeN], 0, zCenters[pipeN])));
+    auto placed_vacuum =
+        assembly.placeVolume(v_vacuum, Transform3D(RotationY(thetas[pipeN]),
+                                                   Position(xCenters[pipeN], 0, zCenters[pipeN])));
+
+    DetElement vacuum_element(sdet, names[pipeN] + "_vacuum", ids[pipeN]);
+    vacuum_element.setPlacement(placed_vacuum);
   }
 
   // Final placement

src/BeamPipeTracking_geo.cpp

@@ -0,0 +1,85 @@
+// SPDX-License-Identifier: LGPL-3.0-or-later
+// Copyright (C) 2024 Simon Gardner
+
+//==========================================================================
+//
+// Places a small sensitive disk of vacuum at the end of beam pipes
+//
+//==========================================================================
+
+#include "DD4hep/DetFactoryHelper.h"
+#include "DD4hep/Printout.h"
+#include "TMath.h"
+#include <XML/Helper.h>
+
+using namespace std;
+using namespace dd4hep;
+
+static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens) {
+
+  using namespace ROOT::Math;
+  xml_det_t x_det   = e;
+  string det_name   = x_det.nameStr();
+  int det_id        = x_det.id();
+  Material m_Vacuum = description.material("Vacuum");
+  string vis_name   = dd4hep::getAttrOrDefault<std::string>(x_det, _Unicode(vis), "BeamPipeVis");
+
+  sens.setType("tracker");
+
+  DetElement sdet(det_name, det_id);
+  Assembly assembly(det_name + "_assembly");
+
+  // Grab info for beamline magnets
+  for (xml_coll_t slice_coll(x_det, _Unicode(slice)); slice_coll; slice_coll++) { // pipes
+
+    string grandmotherName = slice_coll.attr<string>(_Unicode(grandmother));
+    string motherName      = slice_coll.attr<string>(_Unicode(mother));
+    bool detStart          = getAttrOrDefault<bool>(slice_coll, _Unicode(end), true);
+    int pipe_id            = getAttrOrDefault<int>(slice_coll, _Unicode(pipe_id), 0);
+    string slice_name      = slice_coll.attr<string>(_Unicode(name));
+    DetElement mother      = description.detector(grandmotherName).child(motherName);
+
+    // Get the mother volume
+    Volume mother_vol = mother.volume();
+
+    // Get mother volume shape as cone segment
+    ConeSegment mother_shape = mother_vol.solid();
+
+    // Get the parameters of the mother volume
+    double rOuter1 = mother_shape.rMax1();
+    double rOuter2 = mother_shape.rMax2();
+    double length  = 2 * mother_shape.dZ();
+
+    double sensitive_thickness = 0.1 * mm;
+
+    //Calculate R or cone after sensitive layer
+
+    double rEnd = rOuter2 - (rOuter2 - rOuter1) * sensitive_thickness / length;
+    double zPos = length / 2.0 - sensitive_thickness / 2.0;
+    if (detStart) {
+      rEnd = rOuter1 - (rOuter1 - rOuter2) * sensitive_thickness / length;
+      zPos = -length / 2.0 + sensitive_thickness / 2.0;
+    }
+
+    ConeSegment s_start_disk(sensitive_thickness / 2, 0.0, rOuter2, 0.0, rEnd);
+    Volume v_start_disk("v_start_disk_" + motherName, s_start_disk, m_Vacuum);
+    v_start_disk.setSensitiveDetector(sens);
+
+    auto disk_placement = mother_vol.placeVolume(v_start_disk, Position(0.0, 0.0, zPos));
+    disk_placement.addPhysVolID("end", detStart);
+    disk_placement.addPhysVolID("pipe", pipe_id);
+    disk_placement.addPhysVolID("system", det_id);
+
+    DetElement slice_element(sdet, slice_name, pipe_id);
+
+    slice_element.setPlacement(disk_placement);
+    description.declareParent(slice_name, mother);
+  }
+
+  auto pv_assembly = description.worldVolume().placeVolume(assembly, Position(0.0, 0.0, 0.0));
+  sdet.setPlacement(pv_assembly);
+
+  return sdet;
+}
+
+DECLARE_DETELEMENT(BeamPipeTracking, create_detector)

templates/epic.xml.jinja2

@@ -71,6 +71,9 @@
       <limit name="ekin_min" particles="*" value="0.001" unit="MeV" />
       <limit name="range_min" particles="*" value="0.1" unit="mm" />
     </limitset>
+    <limitset name="kill_limits">
+      <limit name="time_max" particles="*" value="0.0" unit="ns"/>
+    </limitset>
     <limitset name="cal_limits">
       <limit name="step_length_max" particles="*" value="5.0" unit="mm"/>
     </limitset>