From aa2cca99cedea1aed9153fa5366b03cf080be26d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Henning=20Schulze=20Ei=C3=9Fing?= Date: Thu, 21 Mar 2024 09:45:20 -0500 Subject: [PATCH 1/5] Add draft for PMT dark count simulation --- fuse/context.py | 1 + .../detector_physics/s1_photon_propagation.py | 2 +- .../detector_physics/s2_photon_propagation.py | 2 +- fuse/plugins/pmt_and_daq/__init__.py | 3 + fuse/plugins/pmt_and_daq/dark_counts.py | 213 ++++++++++++++++++ fuse/plugins/pmt_and_daq/photon_summary.py | 2 +- fuse/plugins/pmt_and_daq/pmt_afterpulses.py | 2 +- fuse/plugins/truth_information/peak_truth.py | 4 +- .../truth_information/records_truth.py | 5 +- 9 files changed, 228 insertions(+), 6 deletions(-) create mode 100644 fuse/plugins/pmt_and_daq/dark_counts.py diff --git a/fuse/context.py b/fuse/context.py index 905cc4c8..d9e7c403 100644 --- a/fuse/context.py +++ b/fuse/context.py @@ -43,6 +43,7 @@ # Plugins to simulate PMTs and DAQ pmt_and_daq_plugins = [ fuse.pmt_and_daq.PMTAfterPulses, + fuse.pmt_and_daq.DarkCounts, fuse.pmt_and_daq.PhotonSummary, fuse.pmt_and_daq.PulseWindow, fuse.pmt_and_daq.PMTResponseAndDAQ, diff --git a/fuse/plugins/detector_physics/s1_photon_propagation.py b/fuse/plugins/detector_physics/s1_photon_propagation.py index 2b43272c..b6e0231e 100644 --- a/fuse/plugins/detector_physics/s1_photon_propagation.py +++ b/fuse/plugins/detector_physics/s1_photon_propagation.py @@ -45,7 +45,7 @@ class S1PhotonPropagationBase(FuseBasePlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2) or PMT AP (0)", "photon_type"), np.int8), + (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), ] dtype = dtype + strax.time_fields diff --git a/fuse/plugins/detector_physics/s2_photon_propagation.py b/fuse/plugins/detector_physics/s2_photon_propagation.py index 70c79941..73480b11 100644 --- a/fuse/plugins/detector_physics/s2_photon_propagation.py +++ b/fuse/plugins/detector_physics/s2_photon_propagation.py @@ -52,7 +52,7 @@ class S2PhotonPropagationBase(FuseBaseDownChunkingPlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2) or PMT AP (0)", "photon_type"), np.int8), + (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), ] dtype = dtype + strax.time_fields diff --git a/fuse/plugins/pmt_and_daq/__init__.py b/fuse/plugins/pmt_and_daq/__init__.py index 0aa69962..c4e34a64 100644 --- a/fuse/plugins/pmt_and_daq/__init__.py +++ b/fuse/plugins/pmt_and_daq/__init__.py @@ -9,3 +9,6 @@ from . import photon_pulses from .photon_pulses import * + +from . import dark_counts +from .dark_counts import * diff --git a/fuse/plugins/pmt_and_daq/dark_counts.py b/fuse/plugins/pmt_and_daq/dark_counts.py new file mode 100644 index 00000000..56f21b3d --- /dev/null +++ b/fuse/plugins/pmt_and_daq/dark_counts.py @@ -0,0 +1,213 @@ +import strax +import straxen +import logging +import numpy as np + +from ...plugin import FuseBasePlugin +from ...common import pmt_gains, build_photon_propagation_output +from ...common import ( + init_spe_scaling_factor_distributions, + photon_gain_calculation, +) + +export, __all__ = strax.exporter() + +logging.basicConfig(handlers=[logging.StreamHandler()]) +log = logging.getLogger("fuse.pmt_and_daq.dark_counts") + +@export +class DarkCounts(FuseBasePlugin): + """Plugin to simulate dark counts in a window around the physics interaction""" + + __version__ = "0.0.1" + + depends_on = "microphysics_summary" + provides = "dark_count_photons" + data_kind = "dark_count_photons" + + save_when = strax.SaveWhen.TARGET + + dtype = [ + (("PMT channel of the photon", "channel"), np.int16), + (("Photon creates a double photo-electron emission", "dpe"), np.bool_), + (("Sampled PMT gain for the photon", "photon_gain"), np.int32), + (("ID of the cluster creating the photon", "cluster_id"), np.int32), + (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), + ] + dtype = dtype + strax.time_fields + + # Config options + + enable_dark_counts = straxen.URLConfig( + default=False, + type=bool, + track=True, + help="Decide if you want to to enable dark count simulation", + ) + + # Get this value from a database. For now lets just set it to 50 Hz per PMT + dark_count_rate = straxen.URLConfig( + default=50 * 494, + type=(int, float), + track=True, + help="Rate of dark counts in Hz combined for all PMTs", + ) + + dark_count_left_window = straxen.URLConfig( + default=2e6, + type=int, + track=True, + help="Left window of the dark count simulation", + ) + + dark_count_right_window = straxen.URLConfig( + default=2e6, + type=int, + track=True, + help="Right window of the dark count simulation", + ) + + #Add a default pointing to the database + dark_count_probability_per_pmt = straxen.URLConfig( + track=True, + help="Probability of dark counts per PMT", + ) + + # reconsider this one for dark counts.... + p_double_pe_emision = straxen.URLConfig( + default="take://resource://" + "SIMULATION_CONFIG_FILE.json?&fmt=json" + "&take=p_double_pe_emision", + type=(int, float), + cache=True, + help="Probability of double photo-electron emission", + ) + + pmt_circuit_load_resistor = straxen.URLConfig( + default="take://resource://" + "SIMULATION_CONFIG_FILE.json?&fmt=json" + "&take=pmt_circuit_load_resistor", + type=(int, float), + cache=True, + help="PMT circuit load resistor [kg m^2/(s^3 A)]", + ) + + digitizer_bits = straxen.URLConfig( + default="take://resource://SIMULATION_CONFIG_FILE.json?&fmt=json&take=digitizer_bits", + type=(int, float), + cache=True, + help="Number of bits of the digitizer boards", + ) + + digitizer_voltage_range = straxen.URLConfig( + default="take://resource://" + "SIMULATION_CONFIG_FILE.json?&fmt=json" + "&take=digitizer_voltage_range", + type=(int, float), + cache=True, + help="Voltage range of the digitizer boards [V]", + ) + + gain_model_mc = straxen.URLConfig( + default="cmt://to_pe_model?version=ONLINE&run_id=plugin.run_id", + infer_type=False, + help="PMT gain model", + ) + + photon_area_distribution = straxen.URLConfig( + default="simple_load://resource://simulation_config://" + "SIMULATION_CONFIG_FILE.json?" + "&key=photon_area_distribution" + "&fmt=csv", + cache=True, + help="Photon area distribution", + ) + + n_tpc_pmts = straxen.URLConfig( + type=int, + help="Number of PMTs in the TPC", + ) + + def setup(self): + super().setup() + + # self.dark_count_probability_per_pmt = np.ones(494) / 494 # uniform distribution + + self.gains = pmt_gains( + self.gain_model_mc, + digitizer_voltage_range=self.digitizer_voltage_range, + digitizer_bits=self.digitizer_bits, + pmt_circuit_load_resistor=self.pmt_circuit_load_resistor, + ) + + self.pmt_mask = np.array(self.gains) > 0 # Converted from to pe (from cmt by default) + self.turned_off_pmts = np.nonzero(np.array(self.gains) == 0)[0] + + self.spe_scaling_factor_distributions = init_spe_scaling_factor_distributions( + self.photon_area_distribution + ) + + def compute(self, interactions_in_roi, start, end): + if not self.enable_dark_counts or (len(interactions_in_roi) == 0): + return np.zeros(0, dtype=self.dtype) + + single_simulation_window = np.zeros(len(interactions_in_roi), dtype=strax.interval_dtype) + #single_simulation_window["length"] = (self.dark_count_left_window + self.dark_count_right_window) + single_simulation_window["time"] = interactions_in_roi["time"] + single_simulation_window["dt"] = np.ones(len(interactions_in_roi)) #length will give time in ns + + + simulation_windows = strax.concat_overlapping_hits( + single_simulation_window, + extensions=(self.dark_count_left_window, self.dark_count_right_window), + pmt_channels=(0, self.n_tpc_pmts), + start = start, + end = end + ) + + #simulation_windows["time"] -= np.int64(left_window) + + # Get the number of dark counts in the simulation window + expected_dark_counts_in_simulation_window = simulation_windows["length"].astype(np.int64) * self.dark_count_rate/1e9 + dark_counts_in_simulation_window =self.rng.poisson(expected_dark_counts_in_simulation_window) + + dark_count_times = get_random_times(simulation_windows["length"], dark_counts_in_simulation_window, self.rng) + dark_count_times += np.repeat(simulation_windows["time"], dark_counts_in_simulation_window) + + # Do we need to add the transit time spread here?? + # I guess no as we just distribute the dark counts randomly in time and it will not make any difference + + + # distribute the dark counts to the PMTs + dark_count_channels = self.rng.choice(self.n_tpc_pmts, len(dark_count_times), p=self.dark_count_probability_per_pmt) + + + # We for sure need to update the inputs for this step + photon_gains, photon_is_dpe = photon_gain_calculation( + _photon_channels=dark_count_channels, + p_double_pe_emision=self.p_double_pe_emision, + gains=self.gains, + spe_scaling_factor_distributions=self.spe_scaling_factor_distributions, + rng=self.rng, + ) + photon_is_dpe = False + + # now build the output + result = build_photon_propagation_output( + dtype=self.dtype, + _photon_timings=dark_count_times, + _photon_channels=dark_count_channels, + _photon_gains=photon_gains, + _photon_is_dpe=photon_is_dpe, + _cluster_id=0, #rethink this part... + photon_type=3, + ) + + return result + +#For sure this can be done more efficiently +def get_random_times(interval_length, number_of_entries, rng): + times = [] + for max_time, n in zip(interval_length, number_of_entries): + times.append(rng.uniform(0, max_time, n)) + return np.concatenate(times) \ No newline at end of file diff --git a/fuse/plugins/pmt_and_daq/photon_summary.py b/fuse/plugins/pmt_and_daq/photon_summary.py index 0f1b5f87..0f1bfc86 100644 --- a/fuse/plugins/pmt_and_daq/photon_summary.py +++ b/fuse/plugins/pmt_and_daq/photon_summary.py @@ -10,7 +10,7 @@ class PhotonSummary(VerticalMergerPlugin): """Plugin that concatenates propagated photons for S1, S2 and PMT afterpulses.""" - depends_on = ("propagated_s2_photons", "propagated_s1_photons", "pmt_afterpulses") + depends_on = ("propagated_s2_photons", "propagated_s1_photons", "pmt_afterpulses", "dark_count_photons") provides = "photon_summary" data_kind = "propagated_photons" diff --git a/fuse/plugins/pmt_and_daq/pmt_afterpulses.py b/fuse/plugins/pmt_and_daq/pmt_afterpulses.py index 4affb54b..f1a00757 100644 --- a/fuse/plugins/pmt_and_daq/pmt_afterpulses.py +++ b/fuse/plugins/pmt_and_daq/pmt_afterpulses.py @@ -35,7 +35,7 @@ class PMTAfterPulses(FuseBasePlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2) or PMT AP (0)", "photon_type"), np.int8), + (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), ] dtype = dtype + strax.time_fields diff --git a/fuse/plugins/truth_information/peak_truth.py b/fuse/plugins/truth_information/peak_truth.py index 3cd4e86e..5d87d9f5 100644 --- a/fuse/plugins/truth_information/peak_truth.py +++ b/fuse/plugins/truth_information/peak_truth.py @@ -9,7 +9,7 @@ @export class PeakTruth(strax.OverlapWindowPlugin): - __version__ = "0.0.3" + __version__ = "0.0.4" depends_on = ( "photon_summary", @@ -26,6 +26,7 @@ class PeakTruth(strax.OverlapWindowPlugin): ("s1_photons_in_peak", np.int32), ("s2_photons_in_peak", np.int32), ("ap_photons_in_peak", np.int32), + ("dark_count_photons_in_peak", np.int32), ("raw_area_truth", np.float32), ("observable_energy_truth", np.float32), ("number_of_contributing_clusters_s1", np.int16), @@ -113,6 +114,7 @@ def compute(self, interactions_in_roi, propagated_photons, peaks): "s1": 1, "s2": 2, "ap": 0, + "dark_count": 3, } for i in range(n_peaks): diff --git a/fuse/plugins/truth_information/records_truth.py b/fuse/plugins/truth_information/records_truth.py index 2f57debd..43075f7c 100644 --- a/fuse/plugins/truth_information/records_truth.py +++ b/fuse/plugins/truth_information/records_truth.py @@ -12,7 +12,7 @@ class RecordsTruth(strax.Plugin): """Plugin that computes the truth information for raw_records.""" - __version__ = "0.0.1" + __version__ = "0.0.2" depends_on = ("raw_records", "photon_summary") provides = "records_truth" @@ -21,6 +21,7 @@ class RecordsTruth(strax.Plugin): (("Number of S1 photons in record", "s1_photons_in_record"), np.int32), (("Number of S2 photons in record", "s2_photons_in_record"), np.int32), (("Number of AP photons in record", "ap_photons_in_record"), np.int32), + (("Number of dark counts in record", "dark_count_photons_in_record"), np.int32), (("Sum of the photon gains", "raw_area"), np.float32), ] @@ -97,6 +98,7 @@ def compute(self, raw_records, propagated_photons): result["s1_photons_in_record"][result_mask] = result_buffer["s1_photons_in_record"] result["s2_photons_in_record"][result_mask] = result_buffer["s2_photons_in_record"] result["ap_photons_in_record"][result_mask] = result_buffer["ap_photons_in_record"] + result["dark_count_photons_in_record"][result_mask] = result_buffer["dark_count_photons_in_record"] return result @@ -109,6 +111,7 @@ def fill_result_buffer(list_input, result_buffer): result_buffer["s1_photons_in_record"][i] = np.sum(photons["photon_type"] == 1) result_buffer["s2_photons_in_record"][i] = np.sum(photons["photon_type"] == 2) result_buffer["ap_photons_in_record"][i] = np.sum(photons["photon_type"] == 0) + result_buffer["dark_count_photons_in_record"][i] = np.sum(photons["photon_type"] == 3) def split_photons_by_channel(propagated_photons): From 0e4166b16d41fdfc143b1a0addea5a92b01f6f3b Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 21 Mar 2024 14:47:07 +0000 Subject: [PATCH 2/5] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../detector_physics/s1_photon_propagation.py | 5 +- .../detector_physics/s2_photon_propagation.py | 5 +- fuse/plugins/pmt_and_daq/dark_counts.py | 57 ++++++++++++------- fuse/plugins/pmt_and_daq/photon_summary.py | 7 ++- fuse/plugins/pmt_and_daq/pmt_afterpulses.py | 5 +- .../truth_information/records_truth.py | 4 +- 6 files changed, 56 insertions(+), 27 deletions(-) diff --git a/fuse/plugins/detector_physics/s1_photon_propagation.py b/fuse/plugins/detector_physics/s1_photon_propagation.py index b6e0231e..56fe8bc6 100644 --- a/fuse/plugins/detector_physics/s1_photon_propagation.py +++ b/fuse/plugins/detector_physics/s1_photon_propagation.py @@ -45,7 +45,10 @@ class S1PhotonPropagationBase(FuseBasePlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), + ( + ("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), + np.int8, + ), ] dtype = dtype + strax.time_fields diff --git a/fuse/plugins/detector_physics/s2_photon_propagation.py b/fuse/plugins/detector_physics/s2_photon_propagation.py index 73480b11..19d0156c 100644 --- a/fuse/plugins/detector_physics/s2_photon_propagation.py +++ b/fuse/plugins/detector_physics/s2_photon_propagation.py @@ -52,7 +52,10 @@ class S2PhotonPropagationBase(FuseBaseDownChunkingPlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), + ( + ("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), + np.int8, + ), ] dtype = dtype + strax.time_fields diff --git a/fuse/plugins/pmt_and_daq/dark_counts.py b/fuse/plugins/pmt_and_daq/dark_counts.py index 56f21b3d..14a59e34 100644 --- a/fuse/plugins/pmt_and_daq/dark_counts.py +++ b/fuse/plugins/pmt_and_daq/dark_counts.py @@ -15,9 +15,11 @@ logging.basicConfig(handlers=[logging.StreamHandler()]) log = logging.getLogger("fuse.pmt_and_daq.dark_counts") + @export class DarkCounts(FuseBasePlugin): - """Plugin to simulate dark counts in a window around the physics interaction""" + """Plugin to simulate dark counts in a window around the physics + interaction.""" __version__ = "0.0.1" @@ -32,7 +34,10 @@ class DarkCounts(FuseBasePlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), + ( + ("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), + np.int8, + ), ] dtype = dtype + strax.time_fields @@ -67,7 +72,7 @@ class DarkCounts(FuseBasePlugin): help="Right window of the dark count simulation", ) - #Add a default pointing to the database + # Add a default pointing to the database dark_count_probability_per_pmt = straxen.URLConfig( track=True, help="Probability of dark counts per PMT", @@ -152,35 +157,42 @@ def compute(self, interactions_in_roi, start, end): return np.zeros(0, dtype=self.dtype) single_simulation_window = np.zeros(len(interactions_in_roi), dtype=strax.interval_dtype) - #single_simulation_window["length"] = (self.dark_count_left_window + self.dark_count_right_window) + # single_simulation_window["length"] = (self.dark_count_left_window + self.dark_count_right_window) single_simulation_window["time"] = interactions_in_roi["time"] - single_simulation_window["dt"] = np.ones(len(interactions_in_roi)) #length will give time in ns - + single_simulation_window["dt"] = np.ones( + len(interactions_in_roi) + ) # length will give time in ns simulation_windows = strax.concat_overlapping_hits( single_simulation_window, - extensions=(self.dark_count_left_window, self.dark_count_right_window), + extensions=(self.dark_count_left_window, self.dark_count_right_window), pmt_channels=(0, self.n_tpc_pmts), - start = start, - end = end - ) + start=start, + end=end, + ) - #simulation_windows["time"] -= np.int64(left_window) + # simulation_windows["time"] -= np.int64(left_window) # Get the number of dark counts in the simulation window - expected_dark_counts_in_simulation_window = simulation_windows["length"].astype(np.int64) * self.dark_count_rate/1e9 - dark_counts_in_simulation_window =self.rng.poisson(expected_dark_counts_in_simulation_window) - - dark_count_times = get_random_times(simulation_windows["length"], dark_counts_in_simulation_window, self.rng) + expected_dark_counts_in_simulation_window = ( + simulation_windows["length"].astype(np.int64) * self.dark_count_rate / 1e9 + ) + dark_counts_in_simulation_window = self.rng.poisson( + expected_dark_counts_in_simulation_window + ) + + dark_count_times = get_random_times( + simulation_windows["length"], dark_counts_in_simulation_window, self.rng + ) dark_count_times += np.repeat(simulation_windows["time"], dark_counts_in_simulation_window) # Do we need to add the transit time spread here?? # I guess no as we just distribute the dark counts randomly in time and it will not make any difference - # distribute the dark counts to the PMTs - dark_count_channels = self.rng.choice(self.n_tpc_pmts, len(dark_count_times), p=self.dark_count_probability_per_pmt) - + dark_count_channels = self.rng.choice( + self.n_tpc_pmts, len(dark_count_times), p=self.dark_count_probability_per_pmt + ) # We for sure need to update the inputs for this step photon_gains, photon_is_dpe = photon_gain_calculation( @@ -189,7 +201,7 @@ def compute(self, interactions_in_roi, start, end): gains=self.gains, spe_scaling_factor_distributions=self.spe_scaling_factor_distributions, rng=self.rng, - ) + ) photon_is_dpe = False # now build the output @@ -199,15 +211,16 @@ def compute(self, interactions_in_roi, start, end): _photon_channels=dark_count_channels, _photon_gains=photon_gains, _photon_is_dpe=photon_is_dpe, - _cluster_id=0, #rethink this part... + _cluster_id=0, # rethink this part... photon_type=3, ) return result -#For sure this can be done more efficiently + +# For sure this can be done more efficiently def get_random_times(interval_length, number_of_entries, rng): times = [] for max_time, n in zip(interval_length, number_of_entries): times.append(rng.uniform(0, max_time, n)) - return np.concatenate(times) \ No newline at end of file + return np.concatenate(times) diff --git a/fuse/plugins/pmt_and_daq/photon_summary.py b/fuse/plugins/pmt_and_daq/photon_summary.py index 0f1bfc86..1bf83587 100644 --- a/fuse/plugins/pmt_and_daq/photon_summary.py +++ b/fuse/plugins/pmt_and_daq/photon_summary.py @@ -10,7 +10,12 @@ class PhotonSummary(VerticalMergerPlugin): """Plugin that concatenates propagated photons for S1, S2 and PMT afterpulses.""" - depends_on = ("propagated_s2_photons", "propagated_s1_photons", "pmt_afterpulses", "dark_count_photons") + depends_on = ( + "propagated_s2_photons", + "propagated_s1_photons", + "pmt_afterpulses", + "dark_count_photons", + ) provides = "photon_summary" data_kind = "propagated_photons" diff --git a/fuse/plugins/pmt_and_daq/pmt_afterpulses.py b/fuse/plugins/pmt_and_daq/pmt_afterpulses.py index f1a00757..a2c3d633 100644 --- a/fuse/plugins/pmt_and_daq/pmt_afterpulses.py +++ b/fuse/plugins/pmt_and_daq/pmt_afterpulses.py @@ -35,7 +35,10 @@ class PMTAfterPulses(FuseBasePlugin): (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), + ( + ("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), + np.int8, + ), ] dtype = dtype + strax.time_fields diff --git a/fuse/plugins/truth_information/records_truth.py b/fuse/plugins/truth_information/records_truth.py index 43075f7c..a78c62be 100644 --- a/fuse/plugins/truth_information/records_truth.py +++ b/fuse/plugins/truth_information/records_truth.py @@ -98,7 +98,9 @@ def compute(self, raw_records, propagated_photons): result["s1_photons_in_record"][result_mask] = result_buffer["s1_photons_in_record"] result["s2_photons_in_record"][result_mask] = result_buffer["s2_photons_in_record"] result["ap_photons_in_record"][result_mask] = result_buffer["ap_photons_in_record"] - result["dark_count_photons_in_record"][result_mask] = result_buffer["dark_count_photons_in_record"] + result["dark_count_photons_in_record"][result_mask] = result_buffer[ + "dark_count_photons_in_record" + ] return result From a2620557367cdf17b8d5e69aecee40d551004d92 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Henning=20Schulze=20Ei=C3=9Fing?= Date: Thu, 21 Mar 2024 16:53:52 +0100 Subject: [PATCH 3/5] Cleanup --- fuse/plugins/pmt_and_daq/dark_counts.py | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/fuse/plugins/pmt_and_daq/dark_counts.py b/fuse/plugins/pmt_and_daq/dark_counts.py index 14a59e34..f2ae4373 100644 --- a/fuse/plugins/pmt_and_daq/dark_counts.py +++ b/fuse/plugins/pmt_and_daq/dark_counts.py @@ -157,11 +157,8 @@ def compute(self, interactions_in_roi, start, end): return np.zeros(0, dtype=self.dtype) single_simulation_window = np.zeros(len(interactions_in_roi), dtype=strax.interval_dtype) - # single_simulation_window["length"] = (self.dark_count_left_window + self.dark_count_right_window) single_simulation_window["time"] = interactions_in_roi["time"] - single_simulation_window["dt"] = np.ones( - len(interactions_in_roi) - ) # length will give time in ns + single_simulation_window["dt"] = np.ones(len(interactions_in_roi)) simulation_windows = strax.concat_overlapping_hits( single_simulation_window, @@ -171,8 +168,6 @@ def compute(self, interactions_in_roi, start, end): end=end, ) - # simulation_windows["time"] -= np.int64(left_window) - # Get the number of dark counts in the simulation window expected_dark_counts_in_simulation_window = ( simulation_windows["length"].astype(np.int64) * self.dark_count_rate / 1e9 @@ -186,9 +181,6 @@ def compute(self, interactions_in_roi, start, end): ) dark_count_times += np.repeat(simulation_windows["time"], dark_counts_in_simulation_window) - # Do we need to add the transit time spread here?? - # I guess no as we just distribute the dark counts randomly in time and it will not make any difference - # distribute the dark counts to the PMTs dark_count_channels = self.rng.choice( self.n_tpc_pmts, len(dark_count_times), p=self.dark_count_probability_per_pmt From 2bbcda3595e7c6ce39050d0d55288b9be91de2e1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Henning=20Schulze=20Ei=C3=9Fing?= Date: Mon, 25 Mar 2024 07:49:55 -0500 Subject: [PATCH 4/5] Make dark_count_probability_per_pmt optional --- fuse/plugins/pmt_and_daq/dark_counts.py | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/fuse/plugins/pmt_and_daq/dark_counts.py b/fuse/plugins/pmt_and_daq/dark_counts.py index f2ae4373..ee95ced8 100644 --- a/fuse/plugins/pmt_and_daq/dark_counts.py +++ b/fuse/plugins/pmt_and_daq/dark_counts.py @@ -74,6 +74,7 @@ class DarkCounts(FuseBasePlugin): # Add a default pointing to the database dark_count_probability_per_pmt = straxen.URLConfig( + default=None, track=True, help="Probability of dark counts per PMT", ) @@ -136,8 +137,6 @@ class DarkCounts(FuseBasePlugin): def setup(self): super().setup() - # self.dark_count_probability_per_pmt = np.ones(494) / 494 # uniform distribution - self.gains = pmt_gains( self.gain_model_mc, digitizer_voltage_range=self.digitizer_voltage_range, @@ -182,9 +181,13 @@ def compute(self, interactions_in_roi, start, end): dark_count_times += np.repeat(simulation_windows["time"], dark_counts_in_simulation_window) # distribute the dark counts to the PMTs - dark_count_channels = self.rng.choice( - self.n_tpc_pmts, len(dark_count_times), p=self.dark_count_probability_per_pmt - ) + if self.dark_count_probability_per_pmt is None: + log.warning("No dark count probability per PMT provided. Using uniform distribution.") + dark_count_channels = self.rng.choice(self.n_tpc_pmts, len(dark_count_times)) + else: + dark_count_channels = self.rng.choice( + self.n_tpc_pmts, len(dark_count_times), p=self.dark_count_probability_per_pmt + ) # We for sure need to update the inputs for this step photon_gains, photon_is_dpe = photon_gain_calculation( From d49e69397ea14bcd039766f03b39d2a1a27db502 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Henning=20Schulze=20Ei=C3=9Fing?= Date: Fri, 19 Apr 2024 14:05:42 +0200 Subject: [PATCH 5/5] Update DarkCounts dtype --- fuse/dtypes.py | 2 +- fuse/plugins/pmt_and_daq/dark_counts.py | 13 ++----------- 2 files changed, 3 insertions(+), 12 deletions(-) diff --git a/fuse/dtypes.py b/fuse/dtypes.py index 65acfafd..058bcef4 100644 --- a/fuse/dtypes.py +++ b/fuse/dtypes.py @@ -74,5 +74,5 @@ (("Photon creates a double photo-electron emission", "dpe"), np.bool_), (("Sampled PMT gain for the photon", "photon_gain"), np.int32), (("ID of the cluster creating the photon", "cluster_id"), np.int32), - (("Type of the photon. S1 (1), S2 (2) or PMT AP (0)", "photon_type"), np.int8), + (("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), np.int8), ] diff --git a/fuse/plugins/pmt_and_daq/dark_counts.py b/fuse/plugins/pmt_and_daq/dark_counts.py index ee95ced8..62ae9142 100644 --- a/fuse/plugins/pmt_and_daq/dark_counts.py +++ b/fuse/plugins/pmt_and_daq/dark_counts.py @@ -4,6 +4,7 @@ import numpy as np from ...plugin import FuseBasePlugin +from ...dtypes import propagated_photons_fields from ...common import pmt_gains, build_photon_propagation_output from ...common import ( init_spe_scaling_factor_distributions, @@ -29,17 +30,7 @@ class DarkCounts(FuseBasePlugin): save_when = strax.SaveWhen.TARGET - dtype = [ - (("PMT channel of the photon", "channel"), np.int16), - (("Photon creates a double photo-electron emission", "dpe"), np.bool_), - (("Sampled PMT gain for the photon", "photon_gain"), np.int32), - (("ID of the cluster creating the photon", "cluster_id"), np.int32), - ( - ("Type of the photon. S1 (1), S2 (2), PMT AP (0) or dark count (3)", "photon_type"), - np.int8, - ), - ] - dtype = dtype + strax.time_fields + dtype = propagated_photons_fields + strax.time_fields # Config options