Fix user limits (step limiter)

docs/source/developer_guide/developer_guide_physics.rst

@@ -0,0 +1,161 @@
+Physics
+=======
+
+This chapter documents implementation details of GATE's physics handling for
+developers. It focuses on the Python-side orchestration around Geant4 physics
+objects, i.e. the code that turns user-facing settings into Geant4 regions,
+physics constructors, processes, and cuts.
+
+Step Limiter
+------------
+
+The step limiter is implemented through Geant4 user limits. On the GATE side,
+the relevant code is split across three places:
+
+* ``opengate.managers.PhysicsManager`` exposes the user-facing API.
+* ``opengate.physics.Region`` stores per-region user limits and builds the
+  corresponding Geant4 region objects.
+* ``opengate.physics.UserLimitsPhysics`` is a Geant4 physics constructor that
+  adds the required Geant4 processes to the selected particles.
+
+User-facing entry points
+~~~~~~~~~~~~~~~~~~~~~~~~
+
+Users set the maximum step size either through the physics manager:
+
+.. code-block:: python
+
+    sim.physics_manager.set_max_step_size(volume.name, 1 * gate.g4_units.mm)
+
+or through the volume convenience method:
+
+.. code-block:: python
+
+    volume.set_max_step_size(1 * gate.g4_units.mm)
+
+Internally, both paths associate the volume with a ``Region`` object and store
+the value in ``region.user_limits["max_step_size"]``. If no region exists yet
+for the volume, ``PhysicsManager.find_or_create_region()`` creates one. The
+world volume is associated with Geant4's ``DefaultRegionForTheWorld``; other
+volumes use a GATE-created region named ``<volume_name>_region`` by default.
+
+The particles to which user limits are applied are selected with:
+
+.. code-block:: python
+
+    sim.physics_manager.user_limits_particles = ["proton", "GenericIon"]
+
+This setting is a list-like user input. A single string is accepted by the
+setter hook and converted to a one-element list. The special values are:
+
+* ``"all"``: apply user limits to all Geant4 particles known to the particle
+  table.
+* ``"all_charged"``: add the step limiter to charged particles, following
+  Geant4's ``G4StepLimiterPhysics`` behaviour.
+
+Particle names are intentionally not limited to the production-cut aliases
+(``gamma``, ``electron``, ``positron``, ``proton``). Any Geant4 particle name is
+accepted, after applying the small GATE-to-Geant4 alias translation implemented
+by ``translate_particle_name_gate_to_geant4()``.
+
+Region initialization
+~~~~~~~~~~~~~~~~~~~~~
+
+During simulation initialization, every ``Region`` runs
+``Region.initialize_g4_user_limits()``. If none of the user-limit fields is set,
+no ``G4UserLimits`` object is created. Otherwise, GATE creates one and fills all
+limits:
+
+* ``max_step_size`` maps to ``G4UserLimits.SetMaxAllowedStep()``.
+* ``max_track_length`` maps to ``SetUserMaxTrackLength()``.
+* ``max_time`` maps to ``SetUserMaxTime()``.
+* ``min_ekine`` maps to ``SetUserMinEkine()``.
+* ``min_range`` maps to ``SetUserMinRange()``.
+
+Unset upper limits are replaced by ``FLOAT_MAX`` and unset lower limits by
+``0``. This lets one ``G4UserLimits`` object represent only the limits that the
+user actually requested.
+
+After that, ``Region.initialize_g4_region()`` creates or finds the Geant4
+``G4Region``, attaches the ``G4UserLimits`` object to it, and adds the root
+logical volumes associated with the GATE region. This is the part that tells
+Geant4 where the user-limit values apply.
+
+Registering the Geant4 processes
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+A ``G4UserLimits`` object alone is not sufficient. Geant4 also needs tracking
+processes attached to particles so that the limits are actually enforced.
+
+The physics engine checks all regions in
+``PhysicsEngine.initialize_user_limits_physics()``:
+
+.. code-block:: python
+
+    if region.need_step_limiter():
+        need_step_limiter = True
+    if region.need_user_special_cut():
+        need_user_special_cut = True
+
+If at least one region needs a maximum step size or another user special cut,
+GATE registers ``UserLimitsPhysics`` on the active Geant4 physics list.
+
+``UserLimitsPhysics.ConstructProcess()`` then iterates over the full Geant4
+particle table and decides, particle by particle, whether to add:
+
+* ``G4StepLimiter("StepLimiter")`` for ``max_step_size``.
+* ``G4UserSpecialCuts("UserSpecialCut")`` for the other user limits.
+
+For explicit particle names and for ``"all"``, GATE adds both processes. For
+``"all_charged"``, GATE only adds ``G4StepLimiter`` to charged particles. This
+mirrors Geant4's default step-limiter constructor behaviour and avoids applying
+the other special cuts unless the user explicitly asked for the particle.
+
+The created Geant4 process objects are stored in
+``UserLimitsPhysics.g4_step_limiter_storage`` and
+``UserLimitsPhysics.g4_special_user_cuts_storage``. This storage is important:
+the objects are created from Python via pybind11, and keeping references avoids
+their garbage collection after ``ConstructProcess()`` returns.
+
+Validation and name translation
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+Before processes are added, ``UserLimitsPhysics`` validates all explicitly named
+particles against ``G4ParticleTable``. Unknown particles are rejected with a
+fatal error. The special selectors ``"all"`` and ``"all_charged"`` are removed
+before validation.
+
+The translation helper ``translate_particle_name_gate_to_geant4()`` maps the
+historical GATE names used for production cuts to Geant4 names:
+
+* ``electron`` -> ``e-``
+* ``positron`` -> ``e+``
+* ``gamma`` -> ``gamma``
+* ``proton`` -> ``proton``
+
+Names not present in this alias table are passed through unchanged. This is what
+allows inputs such as ``"GenericIon"`` or other Geant4 particle names to work
+without extending the production-cut particle list.
+
+Important distinction from production cuts
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+Production cuts and user limits have different particle-name scopes in GATE.
+Production cuts are still constrained to the particle aliases in
+``cut_particle_names`` because that mirrors Geant4 production-cut usage in GATE.
+User limits, including the step limiter, are more general: Geant4 can apply them
+to any particle with a process manager.
+
+This distinction is why ``sim.physics_manager.user_limits_particles = "all"``
+means all Geant4 particles for user limits. Tests that need the historical
+GATE behaviour should explicitly request:
+
+.. code-block:: python
+
+    sim.physics_manager.user_limits_particles = [
+        "proton",
+        "gamma",
+        "electron",
+        "positron",
+    ]
+

docs/source/developer_guide/index.rst

@@ -18,6 +18,12 @@ Developer guide
    developer_guide_dynamic_parametrisations
    developer_guide_auxiliary_attributes
 
+.. toctree::
+   :caption: Physics
+   :maxdepth: 2
+
+   developer_guide_physics
+
 .. toctree::
    :caption: Good to know
    :maxdepth: 2

docs/source/user_guide/user_guide_physics.rst

@@ -131,7 +131,7 @@ Additionally, you need to tell GATE to which particles you want to apply the ste
 
 .. code-block:: python
 
-    sim.physics_manager.set_user_limits_particles(['electron', 'positron'])
+    sim.physics_manager.user_limits_particles = ['electron', 'positron']
 
 There are other user limits like ''maximum track length'' and ''minimium kinetic energy'', that are used in analogy to the ''maximum step size''.
 You can also use Regions if your geometry is complex. Have a look at the section :ref:`user-limits-details-label` in the detailed part of this user guide for more info.

opengate/managers.py

@@ -573,6 +573,13 @@ def _setter_hook_physics_list_name(self, physics_list_name):
     return physics_list_name
 
 
+def _setter_hook_user_limits_particles(self, particle_names):
+    if not isinstance(particle_names, (list, set, tuple)):
+        return list([particle_names])
+    else:
+        return particle_names
+
+
 class PhysicsManager(GateObject):
     """
     Everything related to the physics (lists, cuts, etc.) should be here.
@@ -629,18 +636,11 @@ class PhysicsManager(GateObject):
             },
         ),
         "user_limits_particles": (
-            Box(
-                [
-                    ("all", False),
-                    ("all_charged", True),
-                    ("gamma", False),
-                    ("electron", False),
-                    ("positron", False),
-                    ("proton", False),
-                ]
-            ),
+            ["all_charged"],
             {
-                "doc": "Switch on (True) or off (False) UserLimits, e.g. step limiter, for individual particles. Default: Step limiter is applied to all charged particles (in accordance with G4 default)."
+                "doc": "List of particles to which UserLimits, e.g. step limiter, are applied. Default: Step limiter "
+                "is applied to all charged particles (in accordance with G4 default).",
+                "setter_hook": _setter_hook_user_limits_particles,
             },
         ),
         "em_parameters": (
@@ -879,8 +879,12 @@ def _normalize_volume_name(self, volume):
             fatal(f"Expected a volume name or a volume object, but received: {volume}")
 
     def find_or_create_region(self, volume_name):
+        if volume_name == self.simulation.world.name:
+            region_name = "DefaultRegionForTheWorld"
+        else:
+            region_name = volume_name + "_region"
         if volume_name not in self.volumes_regions_lut:
-            region = self.add_region(volume_name + "_region")
+            region = self.add_region(region_name)
             region.associate_volume(volume_name)
         else:
             region = self.volumes_regions_lut[volume_name]
@@ -981,19 +985,29 @@ def set_track_structure_em_physics_in_region(
         region.track_structure_em_physics = track_structure_em_physics
 
     def set_user_limits_particles(self, particle_names):
-        if not isinstance(particle_names, (list, set, tuple)):
-            particle_names = list([particle_names])
-        for pn in list(particle_names):
-            # try to get current value to check if particle_name is eligible
-            try:
-                _ = self.user_info.user_limits_particles[pn]
-            except KeyError:
-                fatal(
-                    f"Found unknown particle name '{pn}' in set_user_limits_particles(). Eligible names are "
-                    + ", ".join(list(self.user_info.user_limits_particles.keys()))
-                    + "."
-                )
-            self.user_info.user_limits_particles[pn] = True
+        raise GateDeprecationError(
+            "The function set_user_limits_particles has been removed. Set the particle(s) directly via: \n"
+            "sim.physics_manager.user_limits_particles = XXX, \n"
+            "e.g.\n"
+            "sim.physics_manager.user_limits_particles = 'proton'\n"
+            "or \n"
+            "sim.physics_manager.user_limits_particles = ['gamma', 'electron']"
+        )
+        # if not isinstance(particle_names, (list, set, tuple)):
+        #     particle_names = list([particle_names])
+        # self.user_info.user_limits_particles = particle_names
+
+        # for pn in list(particle_names):
+        #     # try to get current value to check if particle_name is eligible
+        #     try:
+        #         _ = self.user_info.user_limits_particles[pn]
+        #     except KeyError:
+        #         fatal(
+        #             f"Found unknown particle name '{pn}' in set_user_limits_particles(). Eligible names are "
+        #             + ", ".join(list(self.user_info.user_limits_particles.keys()))
+        #             + "."
+        #         )
+        #     self.user_info.user_limits_particles[pn] = True
 
     def freeze_config(self):
         # Freeze the Python-side configuration before any SimulationEngine is

opengate/physics.py

@@ -74,42 +74,44 @@ def ConstructProcess(self):
         """
         ui = self.physics_engine.user_info_physics_manager
 
-        particle_keys_to_consider = []
-        # 'all' overrides individual settings
-        if ui.user_limits_particles["all"] is True:
-            particle_keys_to_consider = list(ui.user_limits_particles.keys())
-        else:
-            keys_to_exclude = ("all", "all_charged")
-            particle_keys_to_consider = [
-                p
-                for p, v in ui.user_limits_particles.items()
-                if v is True and p not in keys_to_exclude
-            ]
-
-        if len(particle_keys_to_consider) == 0:
+        if len(ui.user_limits_particles) == 0:
             self.physics_engine.simulation_engine.simulation.warn_user(
-                "user_limits_particles is False for all particles. No tracking cuts will be applied. Use sim.physics_manager.set_user_limits_particles()."
+                "No particles selected to which UserLimits (incl. step limiter) will be applied. No tracking cuts will be applied. Use sim.physics_manager.set_user_limits_particles()."
             )
 
-        # translate to Geant4 particle names
+        # process user input parameters
+        flag_all = "all" in ui.user_limits_particles
+        flag_all_charged = "all_charged" in ui.user_limits_particles
         particles_to_consider = [
-            translate_particle_name_gate_to_geant4(k) for k in particle_keys_to_consider
+            translate_particle_name_gate_to_geant4(p)
+            for p in ui.user_limits_particles
+            if p not in ["all", "all_charged"]
         ]
 
-        for particle in g4.G4ParticleTable.GetParticleTable().GetParticleList():
+        known_g4_particles = list(
+            g4.G4ParticleTable.GetParticleTable().GetParticleList()
+        )
+        known_g4_particle_names = [p.GetParticleName() for p in known_g4_particles]
+
+        # check for consistency of particle names
+        wrong_particle_names = [
+            p for p in particles_to_consider if p not in known_g4_particle_names
+        ]
+        if len(wrong_particle_names) > 0:
+            fatal(
+                f"UserLimits were requested for the following unknown particles: {wrong_particle_names}"
+            )
+
+        for particle in known_g4_particles:
             add_step_limiter = False
             add_user_special_cuts = False
             p_name = str(particle.GetParticleName())
 
-            if p_name in particles_to_consider:
+            if flag_all or p_name in particles_to_consider:
                 add_step_limiter = True
                 add_user_special_cuts = True
-
             # this reproduces the logic of the Geant4's G4StepLimiterPhysics class
-            if (
-                ui.user_limits_particles["all_charged"] is True
-                and particle.GetPDGCharge() != 0
-            ):
+            elif flag_all_charged and particle.GetPDGCharge() != 0:
                 add_step_limiter = True
 
             if add_step_limiter is True or add_user_special_cuts is True:

opengate/tests/src/actors/test023_interaction_counter_auxiliary_attribute.py

@@ -69,7 +69,7 @@ def run_simulation(
     sim.physics_manager.physics_list_name = "G4EmStandardPhysics_option4"
     if max_step_size is not None:
         sim.physics_manager.set_max_step_size(slab.name, max_step_size)
-        sim.physics_manager.set_user_limits_particles("gamma")
+        sim.physics_manager.user_limits_particles = "gamma"
     sim.run(start_new_process=True)
 
     tree = uproot.open(phsp.get_output_path())["phsp"]

opengate/tests/src/actors/test065_dose_from_edep_ct.py

@@ -117,7 +117,12 @@
     # physics
     sim.physics_manager.physics_list_name = "FTFP_INCLXX_EMZ"
     sim.physics_manager.set_production_cut("world", "all", 1000 * km)
-    sim.physics_manager.user_limits_particles.all = True
+    sim.physics_manager.user_limits_particles = [
+        "proton",
+        "gamma",
+        "electron",
+        "positron",
+    ]
 
     # add dose actor
     dose_postprocess = sim.add_actor("DoseActor", "dose_postprocess")

opengate/tests/src/actors/test081_pgtle_stage1.py

@@ -230,7 +230,7 @@ def voxel_to_mat_name(UH, mat_data):  # from UH to material name
     sim.physics_manager.global_production_cuts.proton = 0.1 * mm
 
     sim.physics_manager.set_max_step_size("ct", 0.1 * mm)
-    sim.physics_manager.set_user_limits_particles(["proton"])
+    sim.physics_manager.user_limits_particles = ["proton"]
 
     # source of proton
     # FIXME to replace by a more realistic proton beam, see tests 044

opengate/tests/src/actors/test081_tle_1_geom.py

@@ -46,7 +46,7 @@ def main(argv):
     sim.physics_manager.physics_list_name = "G4EmStandardPhysics_option3"
     sim.physics_manager.global_production_cuts.all = 1 * mm
     sim.physics_manager.set_max_step_size("waterbox", 1 * mm)
-    sim.physics_manager.set_user_limits_particles("gamma")
+    sim.physics_manager.user_limits_particles = "gamma"
 
     s = f"/process/eLoss/CSDARange true"
     sim.g4_commands_before_init.append(s)

opengate/tests/src/actors/test081_tle_4_high_energy.py

@@ -46,7 +46,7 @@
     sim.physics_manager.physics_list_name = "G4EmStandardPhysics_option3"
     sim.physics_manager.global_production_cuts.all = 1 * mm
     sim.physics_manager.set_max_step_size("waterbox", 1 * mm)
-    sim.physics_manager.set_user_limits_particles("gamma")
+    sim.physics_manager.user_limits_particles = "gamma"
 
     s = f"/process/eLoss/CSDARange true"
     sim.g4_commands_before_init.append(s)