From 8098ab932943dfd25b59c5157f2a04f5a75f24ec Mon Sep 17 00:00:00 2001 From: Pim Nelissen Date: Fri, 20 Mar 2026 10:48:44 +0100 Subject: [PATCH] Add background counts for 3x3 inch NaI. Add auto ROI lookup function from the FWHM. Add ROI test. --- src/pg_rad/background/background.py | 51 ++++ src/pg_rad/configs/defaults.py | 7 + src/pg_rad/data/backgrounds/LU_NaI_3inch.csv | 257 +++++++++++++++++++ tests/test_roi.py | 42 +++ 4 files changed, 357 insertions(+) create mode 100644 src/pg_rad/background/background.py create mode 100644 src/pg_rad/data/backgrounds/LU_NaI_3inch.csv create mode 100644 tests/test_roi.py diff --git a/src/pg_rad/background/background.py b/src/pg_rad/background/background.py new file mode 100644 index 0000000..e37e875 --- /dev/null +++ b/src/pg_rad/background/background.py @@ -0,0 +1,51 @@ +from importlib.resources import files +from typing import Tuple +import numpy as np +from pandas import read_csv + +from pg_rad.configs.defaults import FWHM_PARAMS +from pg_rad.detector.detector import Detector + + +def generate_background( + cps_array: np.ndarray, +) -> np.ndarray: + """ + Generate synthetic background cps for a given detector and energy. + """ + pass + + +def fwhm(A: float, B: float, C: float, E: float) -> float: + return np.sqrt(A + B * E + C * E**2) + + +def get_roi_from_fwhm( + detector: Detector, + energy_keV: float +) -> Tuple[float, float]: + """ + Get the region of interest for given primary gamma energy, using + the 3*FWHM rule. + """ + A, B, C = FWHM_PARAMS.get(detector.type) + delta = 3*fwhm(A, B, C, energy_keV) + return (energy_keV-delta, energy_keV+delta) + + +def get_cps_from_roi( + detector: Detector, + roi_lo: float, + roi_hi: float +) -> float: + + csv = files('pg_rad.data.backgrounds').joinpath(detector.name+'.csv') + data = read_csv(csv) + + # get indices of nearest bins + idx_min = (data["Energy"] - roi_lo).abs().idxmin() + idx_max = (data["Energy"] - roi_hi).abs().idxmin() + idx_start, idx_end = sorted([idx_min, idx_max]) + + cps_sum = data.loc[idx_start:idx_end, "cps"].sum() + return cps_sum diff --git a/src/pg_rad/configs/defaults.py b/src/pg_rad/configs/defaults.py index 9730f70..9ad866e 100644 --- a/src/pg_rad/configs/defaults.py +++ b/src/pg_rad/configs/defaults.py @@ -23,3 +23,10 @@ DETECTOR_EFFICIENCIES = { "NaIR": 0.0216, "NaIF": 0.0254 } + +# A, B, C parameters for different detector types +# for formula FWHM(E) = sqrt(A + B*E + C*E^2) +FWHM_PARAMS = { + "NaI": (-70.55, 2.678, 0.000602), + "HPGe": (1.778, 0.001843, 0.000001) +} diff --git a/src/pg_rad/data/backgrounds/LU_NaI_3inch.csv b/src/pg_rad/data/backgrounds/LU_NaI_3inch.csv new file mode 100644 index 0000000..b75af7d --- /dev/null +++ b/src/pg_rad/data/backgrounds/LU_NaI_3inch.csv @@ -0,0 +1,257 @@ +Energy,Data,cps +0,0,0 +11.6866,0,0 +23.3732,1653,1.181558256 +35.0598,3827,2.735525375 +46.7464,5010,3.581129378 +58.433,7880,5.632594711 +70.1196,11068,7.911365261 +81.8062,13414,9.588277341 +93.4928,14536,10.39027877 +105.179,14580,10.42172981 +116.866,14138,10.10578985 +128.553,13559,9.691922802 +140.239,12725,9.095782702 +151.926,11597,8.289492495 +163.612,10696,7.645461044 +175.299,9847,7.038598999 +186.986,9063,6.478198713 +198.672,8192,5.855611151 +210.359,7627,5.451751251 +222.045,6941,4.961401001 +233.732,6450,4.610436026 +245.419,6046,4.321658327 +257.105,5560,3.974267334 +268.792,4845,3.463187991 +280.478,4267,3.05003574 +292.165,4066,2.906361687 +303.852,3711,2.652609006 +315.538,3413,2.439599714 +327.225,3080,2.201572552 +338.911,3035,2.169406719 +350.598,2789,1.993566833 +362.285,2766,1.977126519 +373.971,2432,1.73838456 +385.658,2182,1.55968549 +397.344,2010,1.436740529 +409.031,1911,1.365975697 +420.718,1871,1.337383846 +432.404,1678,1.199428163 +444.091,1601,1.144388849 +455.777,1538,1.099356683 +467.464,1473,1.052894925 +479.151,1390,0.993566833 +490.837,1369,0.978556112 +502.524,1358,0.970693352 +514.21,1356,0.96926376 +525.897,1323,0.945675482 +537.584,1256,0.897784132 +549.27,1171,0.837026447 +560.957,1072,0.766261615 +572.643,1094,0.781987134 +584.33,1157,0.827019299 +596.017,1073,0.766976412 +607.703,1083,0.774124375 +619.39,1042,0.744817727 +631.076,995,0.711222302 +642.763,938,0.670478914 +654.45,770,0.550393138 +666.136,787,0.562544675 +677.823,761,0.543959971 +689.509,703,0.502501787 +701.196,671,0.479628306 +712.883,611,0.436740529 +724.569,599,0.428162974 +736.256,573,0.40957827 +747.942,623,0.445318084 +759.629,636,0.454610436 +771.316,634,0.453180843 +783.002,602,0.430307362 +794.689,525,0.375268049 +806.375,544,0.388849178 +818.062,557,0.39814153 +829.749,527,0.376697641 +841.435,508,0.363116512 +853.122,509,0.363831308 +864.808,478,0.341672623 +876.495,498,0.355968549 +888.182,515,0.368120086 +899.868,521,0.372408863 +911.555,539,0.385275197 +923.241,546,0.390278771 +934.928,545,0.389563974 +946.615,529,0.378127234 +958.301,489,0.349535382 +969.988,490,0.350250179 +981.674,447,0.319513939 +993.361,437,0.312365976 +1005.05,406,0.290207291 +1016.73,380,0.271622588 +1028.42,357,0.255182273 +1040.11,346,0.247319514 +1051.79,348,0.248749107 +1063.48,337,0.240886347 +1075.17,376,0.268763402 +1086.85,376,0.268763402 +1098.54,363,0.259471051 +1110.23,388,0.277340958 +1121.91,330,0.235882773 +1133.6,364,0.260185847 +1145.29,323,0.230879199 +1156.97,359,0.256611866 +1168.66,330,0.235882773 +1180.35,284,0.203002144 +1192.03,256,0.182987848 +1203.72,252,0.180128663 +1215.41,243,0.173695497 +1227.09,264,0.188706219 +1238.78,226,0.16154396 +1250.47,214,0.152966405 +1262.15,184,0.131522516 +1273.84,179,0.127948535 +1285.53,168,0.120085776 +1297.21,185,0.132237312 +1308.9,177,0.126518942 +1320.59,175,0.12508935 +1332.27,171,0.122230164 +1343.96,187,0.133666905 +1355.65,206,0.147248034 +1367.33,255,0.182273052 +1379.02,331,0.23659757 +1390.71,395,0.282344532 +1402.39,540,0.385989993 +1414.08,608,0.43459614 +1425.77,612,0.437455325 +1437.45,575,0.411007863 +1449.14,516,0.368834882 +1460.82,467,0.333809864 +1472.51,369,0.263759828 +1484.2,254,0.181558256 +1495.88,220,0.157255182 +1507.57,125,0.089349535 +1519.26,136,0.097212294 +1530.94,92,0.065761258 +1542.63,74,0.052894925 +1554.32,67,0.047891351 +1566,68,0.048606147 +1577.69,70,0.05003574 +1589.38,58,0.041458184 +1601.06,66,0.047176555 +1612.75,59,0.042172981 +1624.44,47,0.033595425 +1636.12,60,0.042887777 +1647.81,74,0.052894925 +1659.5,74,0.052894925 +1671.18,84,0.060042888 +1682.87,88,0.062902073 +1694.56,88,0.062902073 +1706.24,65,0.046461758 +1717.93,78,0.05575411 +1729.62,84,0.060042888 +1741.3,62,0.04431737 +1752.99,71,0.050750536 +1764.68,54,0.038598999 +1776.36,53,0.037884203 +1788.05,53,0.037884203 +1799.74,34,0.024303074 +1811.42,36,0.025732666 +1823.11,46,0.032880629 +1834.8,45,0.032165833 +1846.48,41,0.029306648 +1858.17,33,0.023588277 +1869.86,33,0.023588277 +1881.54,29,0.020729092 +1893.23,39,0.027877055 +1904.92,37,0.026447462 +1916.6,28,0.020014296 +1928.29,41,0.029306648 +1939.98,39,0.027877055 +1951.66,46,0.032880629 +1963.35,46,0.032880629 +1975.04,36,0.025732666 +1986.72,39,0.027877055 +1998.41,47,0.033595425 +2010.1,59,0.042172981 +2021.78,56,0.040028592 +2033.47,44,0.031451036 +2045.15,57,0.040743388 +2056.84,50,0.035739814 +2068.53,68,0.048606147 +2080.21,56,0.040028592 +2091.9,53,0.037884203 +2103.59,51,0.03645461 +2115.27,30,0.021443888 +2126.96,40,0.028591851 +2138.65,47,0.033595425 +2150.33,42,0.030021444 +2162.02,49,0.035025018 +2173.71,31,0.022158685 +2185.39,29,0.020729092 +2197.08,49,0.035025018 +2208.77,26,0.018584703 +2220.45,32,0.022873481 +2232.14,23,0.016440315 +2243.83,19,0.013581129 +2255.51,31,0.022158685 +2267.2,16,0.011436741 +2278.89,21,0.015010722 +2290.57,19,0.013581129 +2302.26,20,0.014295926 +2313.95,14,0.010007148 +2325.63,18,0.012866333 +2337.32,18,0.012866333 +2349.01,22,0.015725518 +2360.69,27,0.0192995 +2372.38,21,0.015010722 +2384.07,21,0.015010722 +2395.75,35,0.02501787 +2407.44,49,0.035025018 +2419.13,55,0.039313796 +2430.81,47,0.033595425 +2442.5,68,0.048606147 +2454.19,56,0.040028592 +2465.87,57,0.040743388 +2477.56,63,0.045032166 +2489.25,52,0.037169407 +2500.93,61,0.043602573 +2512.62,37,0.026447462 +2524.31,34,0.024303074 +2535.99,34,0.024303074 +2547.68,26,0.018584703 +2559.37,15,0.010721944 +2571.05,16,0.011436741 +2582.74,9,0.006433167 +2594.43,12,0.008577555 +2606.11,8,0.00571837 +2617.8,11,0.007862759 +2629.48,2,0.001429593 +2641.17,7,0.005003574 +2652.86,5,0.003573981 +2664.54,2,0.001429593 +2676.23,2,0.001429593 +2687.92,0,0 +2699.6,6,0.004288778 +2711.29,2,0.001429593 +2722.98,6,0.004288778 +2734.66,5,0.003573981 +2746.35,5,0.003573981 +2758.04,2,0.001429593 +2769.72,5,0.003573981 +2781.41,6,0.004288778 +2793.1,7,0.005003574 +2804.78,2,0.001429593 +2816.47,2,0.001429593 +2828.16,6,0.004288778 +2839.84,5,0.003573981 +2851.53,4,0.002859185 +2863.22,1,0.000714796 +2874.9,4,0.002859185 +2886.59,4,0.002859185 +2898.28,5,0.003573981 +2909.96,2,0.001429593 +2921.65,3,0.002144389 +2933.34,4,0.002859185 +2945.02,5,0.003573981 +2956.71,3,0.002144389 +2968.4,0,0 +2980.08,0,0 diff --git a/tests/test_roi.py b/tests/test_roi.py new file mode 100644 index 0000000..f4a0122 --- /dev/null +++ b/tests/test_roi.py @@ -0,0 +1,42 @@ +import pytest + +from pg_rad.background.background import get_roi_from_fwhm, get_cps_from_roi +from pg_rad.detector.detector import load_detector + +E_gamma = 661.7 # keV + + +@pytest.fixture +def detector(): + return load_detector("LU_NaI_3inch") + + +@pytest.fixture +def roi_Cs137_NaI(detector): + return get_roi_from_fwhm(detector, E_gamma) + + +@pytest.mark.parametrize("ref_roi_lo, ref_roi_hi", [ + (537.58, 792.69) +]) +def test_roi_acquisition_Cs137_NaI(ref_roi_lo, ref_roi_hi, roi_Cs137_NaI): + """Test against ROI used in GammaVision""" + est_roi_lo, est_roi_hi = roi_Cs137_NaI + + assert pytest.approx(est_roi_lo, rel=0.05) == ref_roi_lo + assert pytest.approx(est_roi_hi, rel=0.05) == ref_roi_hi + + +@pytest.mark.parametrize("ref_cps", [ + 13.61 +]) +def test_cps_acquisition_Cs137_NaI( + ref_cps, + roi_Cs137_NaI, + detector +): + """Test against cps from GammaVision""" + est_roi_lo, est_roi_hi = roi_Cs137_NaI + est_cps = get_cps_from_roi(detector, est_roi_lo, est_roi_hi) + + assert pytest.approx(est_cps, rel=0.1) == ref_cps