Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Waveform noise tuner #1282

Merged
merged 21 commits into from
Aug 1, 2024
Merged
Show file tree
Hide file tree
Changes from 12 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion lstchain/data/lstchain_lhfit_config.json
Original file line number Diff line number Diff line change
Expand Up @@ -279,7 +279,7 @@
},
"waveform_nsb_tuning":{
"nsb_tuning": false,
"nsb_tuning_ratio": 0.5,
"nsb_tuning_rate_GHz": 0.15,
"spe_location": null,
"pre_computed_multiplicity": 10
},
Expand Down
2 changes: 1 addition & 1 deletion lstchain/data/lstchain_standard_config.json
Original file line number Diff line number Diff line change
Expand Up @@ -274,7 +274,7 @@
},
"waveform_nsb_tuning":{
"nsb_tuning": false,
"nsb_tuning_ratio": 0.5,
"nsb_tuning_rate_GHz": 0.15,
"spe_location": null,
"pre_computed_multiplicity": 10
},
Expand Down
302 changes: 208 additions & 94 deletions lstchain/image/modifier.py

Large diffs are not rendered by default.

8 changes: 3 additions & 5 deletions lstchain/image/tests/test_modifier.py
Original file line number Diff line number Diff line change
Expand Up @@ -87,7 +87,7 @@ def test_tune_nsb_on_waveform():
[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]]
)
added_nsb_fraction, original_nsb = 1.0, {1: 0.1 * u.GHz}
added_nsb_rate = {1: 0.1 * u.GHz}
subarray = LSTEventSource.create_subarray(1)
amplitude_HG = np.zeros(40)
amplitude_LG = np.zeros(40)
Expand All @@ -106,8 +106,7 @@ def test_tune_nsb_on_waveform():
charge_spe_cumulative_pdf = interp1d(spe_integral, spe[0], kind='cubic',
bounds_error=False, fill_value=0.,
assume_sorted=True)
nsb_tunner = WaveformNsbTunner(added_nsb_fraction,
original_nsb,
nsb_tunner = WaveformNsbTunner(added_nsb_rate,
pulse_templates,
charge_spe_cumulative_pdf,
pre_computed_multiplicity=0)
Expand All @@ -118,8 +117,7 @@ def test_tune_nsb_on_waveform():
[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]]
)
nsb_tunner = WaveformNsbTunner(added_nsb_fraction,
original_nsb,
nsb_tunner = WaveformNsbTunner(added_nsb_rate,
pulse_templates,
charge_spe_cumulative_pdf,
pre_computed_multiplicity=10)
Expand Down
45 changes: 26 additions & 19 deletions lstchain/reco/r0_to_dl1.py
Original file line number Diff line number Diff line change
Expand Up @@ -59,7 +59,7 @@
write_subarray_tables
)

from ..io.io import add_column_table, extract_simulation_nsb, dl1_params_lstcam_key, get_resource_path
from ..io.io import add_column_table, dl1_params_lstcam_key, get_resource_path
from ..io.lstcontainers import ExtraImageInfo, DL1MonitoringEventIndexContainer
from ..paths import parse_r0_filename, run_to_dl1_filename, r0_to_dl1_filename
from ..visualization.plot_reconstructor import plot_debug
Expand Down Expand Up @@ -403,25 +403,27 @@ def r0_to_dl1(
metadata=metadata,
)
nsb_tuning = False
nsb_tunner = None
nsb_tuner = None
if 'waveform_nsb_tuning' in config.keys():
nsb_tuning = config['waveform_nsb_tuning']['nsb_tuning']
if nsb_tuning:
if is_simu:
nsb_original = extract_simulation_nsb(input_filename)
pulse_templates = {tel_id: NormalizedPulseTemplate.load_from_eventsource(
subarray.tel[tel_id].camera.readout, resample=True)
for tel_id in config['source_config']['LSTEventSource']['allowed_tels']}
if 'nsb_tuning_ratio' in config['waveform_nsb_tuning'].keys():
if 'nsb_tuning_rate_GHz' in config[
'waveform_nsb_tuning'].keys():
moralejo marked this conversation as resolved.
Show resolved Hide resolved
# get value from config to possibly extract it beforehand on multiple files for averaging purposes
# or gain time
nsb_tuning_ratio = config['waveform_nsb_tuning']['nsb_tuning_ratio']
nsb_tuning_rate = config['waveform_nsb_tuning'][
'nsb_tuning_rate_GHz']
else:
# extract the pedestal variance difference between the current MC file and the target data
# FIXME? fails for multiple telescopes
nsb_tuning_ratio = calculate_required_additional_nsb(input_filename,
config['waveform_nsb_tuning']['target_data'],
config=config)[0]
nsb_tuning_rate, _, _ = calculate_required_additional_nsb(
input_filename,
config['waveform_nsb_tuning']['target_data'],
config=config)
spe_location = (config['waveform_nsb_tuning']['spe_location']
or get_resource_path("data/SinglePhE_ResponseInPhE_expo2Gaus.dat"))
spe = np.loadtxt(spe_location).T
Expand All @@ -430,15 +432,20 @@ def r0_to_dl1(
bounds_error=False, fill_value=0.,
assume_sorted=True)
pre_computed_multiplicity = config['waveform_nsb_tuning'].get('pre_computed_multiplicity', 10)
logger.info('Tuning NSB on MC waveform from '
+ str(np.asarray(nsb_original))
+ ' to {0:d}%'.format(int(nsb_tuning_ratio * 100 + 100.5))
+ ' for telescopes ids ' + str(config['source_config']['LSTEventSource']['allowed_tels']))
nsb_tunner = WaveformNsbTunner(nsb_tuning_ratio,
nsb_original,
pulse_templates,
charge_spe_cumulative_pdf,
pre_computed_multiplicity)

allowed_tels = config['source_config']['LSTEventSource'][
'allowed_tels']
logger.info('Tuning NSB on MC waveform by adding ')
logger.info(f'{nsb_tuning_rate:.3f} GHz for telescope ids:')
logger.info(f'{allowed_tels}')

nsb_per_tel = {tel_id: nsb_tuning_rate * u.GHz for tel_id in
allowed_tels}

nsb_tuner = WaveformNsbTunner(nsb_per_tel,
moralejo marked this conversation as resolved.
Show resolved Hide resolved
pulse_templates,
charge_spe_cumulative_pdf,
pre_computed_multiplicity)
else:
logger.warning('NSB tuning on waveform active in config but file is real data, option will be ignored')
nsb_tuning = False
Expand All @@ -451,7 +458,7 @@ def r0_to_dl1(
tmp_source = EventSource(input_url=input_filename,
config=Config(config["source_config"]))
if is_simu:
lhfit_fitter.get_ped_from_true_signal_less(tmp_source, nsb_tunner)
lhfit_fitter.get_ped_from_true_signal_less(tmp_source, nsb_tuner)
else:
lhfit_fitter.get_ped_from_interleaved(tmp_source)
del tmp_source
Expand Down Expand Up @@ -570,7 +577,7 @@ def r0_to_dl1(
waveform = event.r1.tel[tel_id].waveform
selected_gains = event.r1.tel[tel_id].selected_gain_channel
mask_high = selected_gains == 0
nsb_tunner.tune_nsb_on_waveform(waveform, tel_id, mask_high, subarray)
nsb_tuner.tune_nsb_on_waveform(waveform, tel_id, mask_high, subarray)

# create image for all events
r1_dl1_calibrator(event)
Expand Down
4 changes: 2 additions & 2 deletions lstchain/scripts/lstchain_tune_nsb_waveform.py
Original file line number Diff line number Diff line change
Expand Up @@ -82,13 +82,13 @@ def main():

config = read_configuration_file(args.config)

nsb_correction_ratio, data_ped_variance, mc_ped_variance = \
extra_nsb_rate, data_ped_variance, mc_ped_variance = \
calculate_required_additional_nsb(args.input_mc, args.input_data,
config=Config(config))

dict_nsb = {
"nsb_tuning": True,
"nsb_tuning_ratio": np.round(nsb_correction_ratio, decimals=2),
"nsb_tuning_rate_GHz": np.round(extra_nsb_rate, decimals=3),
"spe_location": str(get_resource_path("data/SinglePhE_ResponseInPhE_expo2Gaus.dat"))
}

Expand Down
10 changes: 5 additions & 5 deletions lstchain/scripts/tests/test_lstchain_scripts.py
Original file line number Diff line number Diff line change
Expand Up @@ -236,16 +236,16 @@ def test_validity_tune_nsb(tune_nsb):

def test_validity_tune_nsb_waveform(tune_nsb_waveform):
"""
The resulting nsb_tuning_ratio value of -1 expected in this test is meaningless
because the input data do not allow a full test of the functionality.
This test is only a formal check that the script runs.
The resulting nsb_tuning_rate value of -1 expected in this test is
meaningless because the input data do not allow a full test of the
functionality. This test is only a formal check that the script runs.
"""
output_lines = tune_nsb_waveform.stdout.splitlines()
for line in output_lines:
if '"nsb_tuning"' in line:
assert line == ' "nsb_tuning": true,'
if '"nsb_tuning_ratio"' in line:
assert line == ' "nsb_tuning_ratio": -1.0,'
if '"nsb_tuning_rate"' in line:
assert line == ' "nsb_tuning_rate": -1.0,'
if '"spe_location"' in line:
assert line == f' "spe_location": "{get_resource_path("data/SinglePhE_ResponseInPhE_expo2Gaus.dat")}"'

Expand Down
Loading