add stacking scheme for bulk electron tracks

include/RMGHardware.hh

@@ -141,6 +141,9 @@ class RMGHardware : public G4VUserDetectorConstruction {
     /** @brief Get the instance of the world volume, after @ref Construct had been called once. */
     [[nodiscard]] const G4VPhysicalVolume* GetDefinedWorldVolume() const { return fWorld; }
 
+    /** @brief Get the instance of the world volume, after @ref Construct had been called once. */
+    G4VPhysicalVolume* GetDefinedWorldVolume() { return fWorld; }
+
     /** @brief Set the maximum step size.
      *
      * @details This is used as a @c G4UserLimit to limit step sizes to being no larger than the

include/RMGUserInit.hh

@@ -33,6 +33,7 @@
 #include "RMGTrackOutputScheme.hh"
 #include "RMGVGenerator.hh"
 #include "RMGVOutputScheme.hh"
+#include "RMGVolumeDistanceStacker.hh"
 
 /**
  * @brief User initialization class.
@@ -155,6 +156,7 @@ class RMGUserInit {
       AddOptionalOutputScheme<RMGParticleFilterScheme>("ParticleFilter");
       AddOptionalOutputScheme<RMGTrackOutputScheme>("Track");
       AddOptionalOutputScheme<RMGGeometryCheckOutputScheme>("GeometryCheck");
+      AddOptionalOutputScheme<RMGVolumeDistanceStacker>("VolumeStacker");
     }
 
     /**

include/RMGVolumeDistanceStacker.hh

@@ -0,0 +1,56 @@
+// Copyright (C) 2025 Manuel Huber <https://orcid.org/0009-0000-5212-2999>
+//
+// This program is free software: you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License as published by the Free
+// Software Foundation, either version 3 of the License, or (at your option) any
+// later version.
+//
+// This program is distributed in the hope that it will be useful, but WITHOUT
+// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
+// details.
+//
+// You should have received a copy of the GNU Lesser General Public License
+// along with this program.  If not, see <https://www.gnu.org/licenses/>.
+
+#ifndef _RMG_VOLUME_DISTANCE_STACKER_HH_
+#define _RMG_VOLUME_DISTANCE_STACKER_HH_
+
+#include <memory>
+#include <optional>
+#include <string>
+
+#include "G4GenericMessenger.hh"
+
+#include "RMGVOutputScheme.hh"
+
+/** @brief Special scheme to stack electron/positron tracks created in a specific volume. */
+class RMGVolumeDistanceStacker : public RMGVOutputScheme {
+
+  public:
+
+    RMGVolumeDistanceStacker();
+
+    /** @brief Wraps @c G4UserStackingAction::StackingActionClassify
+     *  @details This is used to classify all e-/e+ tracks as @c fWaiting if the conditions are met.
+     */
+    std::optional<G4ClassificationOfNewTrack> StackingActionClassify(const G4Track*, int) override;
+
+    /** @brief Set the minimum distance to any other volume for this track to be stacked. */
+    void SetVolumeSafety(double safety) { fVolumeSafety = safety; }
+
+    /** @brief Set the volume name in which to stack e-/e+ tracks. */
+    void SetVolumeName(std::string volume) { fVolumeName = volume; }
+
+  private:
+
+    std::unique_ptr<G4GenericMessenger> fMessenger;
+    void DefineCommands();
+
+    double fVolumeSafety = -1;
+    std::string fVolumeName;
+};
+
+#endif
+
+// vim: tabstop=2 shiftwidth=2 expandtab

src/CMakeLists.txt

@@ -49,6 +49,7 @@ set(PROJECT_PUBLIC_HEADERS
     ${_root}/include/RMGVGenerator.hh
     ${_root}/include/RMGVOutputScheme.hh
     ${_root}/include/RMGVVertexGenerator.hh
+    ${_root}/include/RMGVolumeDistanceStacker.hh
     ${_root}/include/RMGVertexOutputScheme.hh
     ${_root}/include/RMGWorkerInitialization.hh)
 
@@ -91,7 +92,8 @@ set(PROJECT_SOURCES
     ${_root}/src/RMGUserAction.cc
     ${_root}/src/RMGVertexConfinement.cc
     ${_root}/src/RMGVertexFromFile.cc
-    ${_root}/src/RMGVertexOutputScheme.cc)
+    ${_root}/src/RMGVertexOutputScheme.cc
+    ${_root}/src/RMGVolumeDistanceStacker.cc)
 
 # Write RMGConfig.hh
 # no need to install, it is included in the header list above

src/RMGOutputTools.cc

@@ -21,6 +21,7 @@
 #include "G4Electron.hh"
 #include "G4Gamma.hh"
 #include "G4LogicalVolume.hh"
+#include "G4MultiUnion.hh"
 #include "G4TransportationManager.hh"
 #include "G4VSolid.hh"
 
@@ -423,8 +424,15 @@ double RMGOutputTools::distance_to_surface(const G4VPhysicalVolume* pv, const G4
     sample_tf.Invert();
     const auto sample_point = sample_tf.TransformPoint(position);
     const auto sample_solid = sample_physical->GetLogicalVolume()->GetSolid();
+
+    // MultiUnions have a flag to calculate accurate safeties, instead of just using the outer bounding box.
+    auto mu = dynamic_cast<G4MultiUnion*>(sample_solid);
+    if (mu != nullptr) mu->SetAccurateSafety(true);
+
     const double sample_dist = sample_solid->DistanceToIn(sample_point);
     if (sample_dist < dist) { dist = sample_dist; }
+
+    if (mu != nullptr) mu->SetAccurateSafety(false);
   }
   return dist;
 }

src/RMGVolumeDistanceStacker.cc

@@ -0,0 +1,81 @@
+// Copyright (C) 2025 Manuel Huber <https://orcid.org/0009-0000-5212-2999>
+//
+// This program is free software: you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License as published by the Free
+// Software Foundation, either version 3 of the License, or (at your option) any
+// later version.
+//
+// This program is distributed in the hope that it will be useful, but WITHOUT
+// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
+// details.
+//
+// You should have received a copy of the GNU Lesser General Public License
+// along with this program.  If not, see <https://www.gnu.org/licenses/>.
+
+#include "RMGVolumeDistanceStacker.hh"
+
+#include "G4Electron.hh"
+#include "G4Positron.hh"
+
+#include "RMGHardware.hh"
+#include "RMGManager.hh"
+#include "RMGOutputTools.hh"
+
+namespace u = CLHEP;
+
+RMGVolumeDistanceStacker::RMGVolumeDistanceStacker() { this->DefineCommands(); }
+
+std::optional<G4ClassificationOfNewTrack> RMGVolumeDistanceStacker::StackingActionClassify(
+    const G4Track* aTrack,
+    int stage
+) {
+  // we are only interested in stacking tracks into stage 1 after stage 0 finished.
+  if (stage != 0) return std::nullopt;
+
+  // do not touch the initial track of an event.
+  if (aTrack->GetTrackID() == 0) return std::nullopt;
+
+  // stop if not configured.
+  if (fVolumeName.empty() || fVolumeSafety < 0) return std::nullopt;
+
+  // only defer electron/positron tracks.
+  if (aTrack->GetDefinition() != G4Electron::Definition() &&
+      aTrack->GetDefinition() != G4Positron::Definition())
+    return std::nullopt;
+
+  // note: aTrack->GetLogicalVolumeAtVertex() and aTrack->GetVertexPosition() might not be correctly
+  // set at this point - do not trust it.
+
+  // only defer tracks in the specified volume.
+  const auto vol_name = aTrack->GetVolume()->GetLogicalVolume()->GetName();
+  if (vol_name != fVolumeName) return std::nullopt;
+
+  // only defer tracks that have a minimum distance to other volumes.
+  auto distance = RMGOutputTools::distance_to_surface(aTrack->GetVolume(), aTrack->GetPosition());
+  if (distance < fVolumeSafety) return std::nullopt;
+
+  return fWaiting;
+}
+
+
+void RMGVolumeDistanceStacker::DefineCommands() {
+
+  fMessenger = std::make_unique<G4GenericMessenger>(
+      this,
+      "/RMG/Output/VolumeStacker/",
+      "Commands for controlling stacking tracks in the bulk of a volume."
+  );
+
+  fMessenger->DeclareMethodWithUnit("VolumeSafety", "cm", &RMGVolumeDistanceStacker::SetVolumeSafety)
+      .SetGuidance("Set the minimum distance to any other volume for this track to be stacked.")
+      .SetParameterName("safety", false)
+      .SetStates(G4State_Idle);
+
+  fMessenger->DeclareMethod("VolumeName", &RMGVolumeDistanceStacker::SetVolumeName)
+      .SetGuidance("Set the volume name in which to stack e-/e+ tracks.")
+      .SetParameterName("volume", false)
+      .SetStates(G4State_Idle);
+}
+
+// vim: tabstop=2 shiftwidth=2 expandtab