Merge pull request #68 from pim-n/feature-opt-bkg

optional background and fix cps bug

src/pg_rad/background/background.py

@@ -11,13 +11,17 @@ def generate_background(
     cps_array: np.ndarray,
     detector: Detector,
     energy_keV: float,
+    lam_inp: int | None = None
 
 ) -> np.ndarray:
     """
     Generate synthetic background cps for a given detector and energy.
     """
-    ROI_lo, ROI_hi = get_roi_from_fwhm(detector, energy_keV)
-    lam = get_cps_from_roi(detector, ROI_lo, ROI_hi)
+    if not lam_inp:
+        ROI_lo, ROI_hi = get_roi_from_fwhm(detector, energy_keV)
+        lam = get_cps_from_roi(detector, ROI_lo, ROI_hi)
+    else:
+        lam = lam_inp
 
     rng = np.random.default_rng()
     return rng.poisson(lam=lam, size=cps_array.shape)

src/pg_rad/inputparser/parser.py

@@ -104,7 +104,7 @@ class ConfigParser:
     def _parse_options(self) -> SimulationOptionsSpec:
         options = self.config.get("options", {})
 
-        allowed = {"air_density_kg_per_m3", "seed"}
+        allowed = {"air_density_kg_per_m3", "seed", "bkg_cps"}
         self._warn_unknown_keys(
             section="options",
             provided=set(options.keys()),
@@ -116,12 +116,14 @@ class ConfigParser:
             defaults.DEFAULT_AIR_DENSITY
         )
         seed = options.get("seed")
+        bkg_cps = options.get("bkg_cps")
 
         if not isinstance(air_density, float) or air_density <= 0:
             raise InvalidConfigValueError(
                 "options.air_density_kg_per_m3 must be a positive float "
                 "in kg/m^3."
             )
+
         if (
             seed is not None or
             (isinstance(seed, int) and seed <= 0)
@@ -130,9 +132,17 @@ class ConfigParser:
                 "Seed must be a positive integer value."
             )
 
+        if bkg_cps is not None and (
+            not isinstance(bkg_cps, int) or bkg_cps < 0
+        ):
+            raise InvalidConfigValueError(
+                "Background CPS must be an integer >= 0."
+            )
+
         return SimulationOptionsSpec(
             air_density=air_density,
             seed=seed,
+            bkg_cps=bkg_cps
         )
 
     def _parse_path(self) -> PathSpec:

src/pg_rad/inputparser/specs.py

@@ -18,6 +18,7 @@ class RuntimeSpec:
 class SimulationOptionsSpec:
     air_density: float
     seed: int | None = None
+    bkg_cps: int | None = None
 
 
 @dataclass

src/pg_rad/physics/fluence.py

@@ -109,6 +109,7 @@ def calculate_counts_along_path(
     detector: "Detector",
     velocity: float,
     points_per_segment: int = 10,
+    bkg_cps_input: int | None = None
 ) -> Tuple[np.ndarray, np.ndarray]:
     """Compute the counts recorded in each acquisition period in the landscape.
 
@@ -116,6 +117,7 @@ def calculate_counts_along_path(
         landscape (Landscape): _description_
         detector (Detector): _description_
         points_per_segment (int, optional): _description_. Defaults to 100.
+        bkg_cps_input (int | None, optional): Optional background CPS.
 
     Returns:
         Tuple[np.ndarray, np.ndarray]: Array of acquisition points and
@@ -148,17 +150,23 @@ def calculate_counts_along_path(
         landscape, full_positions, detector
     )
 
-    bkg = generate_background(
-        cps, detector, landscape.point_sources[0].isotope.E
-    )
+    if bkg_cps_input is None:
+        bkg = generate_background(
+            cps, detector, landscape.point_sources[0].isotope.E
+        )
+    elif bkg_cps_input == 0:
+        bkg = bkg_cps_input
+    else:
+        bkg = generate_background(
+            cps, detector, landscape.point_sources[0].isotope.E,
+            lam_inp=bkg_cps_input
+        )
 
     cps_with_bg = cps + bkg
     # reshape so each segment is on a row
     cps_per_seg = cps_with_bg.reshape(num_segments, points_per_segment)
 
     du = s[1] - s[0]
-    integrated_counts = np.trapezoid(cps_per_seg, dx=du, axis=1) / velocity
-    int_counts_result = np.zeros(num_points)
-    int_counts_result[1:] = integrated_counts
+    int_counts = np.trapezoid(cps_per_seg, dx=du, axis=1) / velocity
 
-    return original_distances, s, cps_with_bg, int_counts_result, np.mean(bkg)
+    return original_distances[1:], s, cps_with_bg, int_counts, np.mean(bkg)

src/pg_rad/simulator/engine.py

@@ -47,7 +47,8 @@ class SimulationEngine:
             calculate_counts_along_path(
                 self.landscape,
                 self.detector,
-                velocity=self.runtime_spec.speed
+                velocity=self.runtime_spec.speed,
+                bkg_cps_input=self.sim_spec.bkg_cps
             )
         )