"""Simulation runner for array simulations."""
from simtools.runners.simtel_runner import SimtelRunner
[docs]
class SimulatorArray(SimtelRunner):
"""
SimulatorArray is the interface with sim_telarray to perform array simulations.
Parameters
----------
simulation_parameters : SimulationParameters
Physical run settings.
array_model : ArrayModel
Telescope and site models.
label: str
Instance label.
"""
def __init__(
self,
simulation_parameters,
array_model,
label=None,
context=None,
):
"""Initialize SimulatorArray."""
super().__init__(
label=label, config=simulation_parameters, array_model=array_model, context=context
)
self.simulation_parameters = simulation_parameters
self.array_model = array_model
self.io_handler = self.context.paths
self._log_file = None
[docs]
def make_run_command(self, run_number=None, input_file=None):
"""
Build and return the command to run sim_telarray.
Parameters
----------
input_file: str
Full path of the input CORSIKA file
run_number: int (optional)
run number
Returns
-------
list
Command to run sim_telarray.
"""
self.file_list = self.runner_service.load_files(run_number=run_number)
command = self._common_run_command(run_number)
if self.simulation_parameters.is_calibration_run():
command += self._make_run_command_for_calibration_simulations()
else:
command += self._make_run_command_for_shower_simulations()
# "-C show=all" should be the last option
return [*command, "-C", "show=all", input_file]
def _make_run_command_for_shower_simulations(self):
"""
Build and return the command to run sim_telarray shower simulations.
Returns
-------
str
Command to run sim_telarray.
"""
power_law = SimulatorArray.get_power_law_for_sim_telarray_histograms(
self.simulation_parameters.primary_particle
)
return [
"-C",
f"power_law={power_law}",
]
def _make_run_command_for_calibration_simulations(self):
"""Build sim_telarray command for calibration simulations."""
cfg = self.context.args
altitude = self.array_model.site_model.get_parameter_value_with_unit(
"reference_point_altitude"
).to_value("m")
options = {"Altitude": altitude}
for key in ("nsb_scaling_factor", "stars"):
if key in cfg:
options[key] = cfg[key]
# stars cannot be 'None' (only 'none')
if key == "stars" and options[key] is None:
options[key] = "none"
run_mode = cfg.get("run_mode")
if run_mode in ("pedestals", "pedestals_nsb_only"):
n_events = cfg.get("number_of_pedestal_events", cfg["number_of_events"])
options["pedestal_events"] = n_events
if run_mode == "pedestals_nsb_only":
options.update(self._pedestals_nsb_only_options())
if run_mode == "pedestals_dark":
n_events = cfg.get("number_of_dark_events", cfg["number_of_events"])
options["dark_events"] = n_events
if run_mode == "direct_injection":
n_events = cfg.get("number_of_flasher_events", cfg["number_of_events"])
options["laser_events"] = n_events
cmd = []
for key, value in options.items():
cmd.extend(["-C", f"{key}={value}"])
return cmd
def _common_run_command(self, run_number):
"""Build generic run command for sim_telarray."""
weak_pointing = self._determine_pointing_option()
config_dir = self.array_model.get_config_directory()
self._log_file = self.runner_service.get_file_name(
file_type="sim_telarray_log", run_number=run_number
)
histogram_file = self.runner_service.get_file_name(
file_type="sim_telarray_histogram", run_number=run_number
)
output_file = self.runner_service.get_file_name(
file_type="sim_telarray_output", run_number=run_number
)
self.array_model.export_config_files("sim_telarray")
cmd = [
str(self.context.software.sim_telarray_exe),
"-c",
str(self.array_model.config_file_path),
f"-I{config_dir}",
]
weak_options = {
"telescope_theta": self.simulation_parameters.zenith_angle,
"telescope_phi": self.simulation_parameters.azimuth_angle,
}
options = {
"histogram_file": histogram_file,
"random_state": "none",
"output_file": output_file,
}
try:
options["random_seed"] = self.array_model.sim_telarray_seed.seed_string
except AttributeError as exc:
raise AttributeError("Error setting sim_telarray seed string") from exc
for key, value in options.items():
cmd.extend(["-C", f"{key}={value}"])
for key, value in weak_options.items():
cmd.extend([("-W" if weak_pointing else "-C"), f"{key}={value}"])
return cmd
def _determine_pointing_option(self):
"""
Determine the pointing option for sim_telarray.
Parameters
----------
label: str
Label of the simulation.
Returns
-------
bool:
True if weak pointing is to be used.
"""
if isinstance(self.label, str):
return any(pointing in self.label for pointing in ["divergent", "convergent"])
return False
def _pedestals_nsb_only_options(self):
"""
Generate options to run sim_telarray for nsb-only pedestal simulations.
Returns
-------
dict
Command to run sim_telarray.
"""
options = {}
null_values = [
"fadc_noise",
"fadc_lg_noise",
"qe_variation",
"gain_variation",
"fadc_var_pedestal",
"fadc_err_pedestal",
"fadc_sysvar_pedestal",
"fadc_dev_pedestal",
]
for param in null_values:
options[param] = 0.0
one_values = [
"fadc_lg_var_pedestal",
"fadc_lg_err_pedestal",
"fadc_lg_dev_pedestal",
"fadc_lg_sysvar_pedestal",
]
for param in one_values:
options[param] = -1.0
return options
def _check_run_result(self, run_number=None):
"""
Check if sim_telarray output file exists.
Parameters
----------
run_number: int
Run number.
Returns
-------
bool
True if sim_telarray output file exists.
Raises
------
FileNotFoundError
If sim_telarray output file does not exist.
"""
output_file = self.runner_service.get_file_name(
file_type="sim_telarray_output", run_number=run_number
)
if not output_file.exists():
raise FileNotFoundError(f"sim_telarray output file {output_file} does not exist.")
self._logger.debug(f"sim_telarray output file {output_file} exists.")
return True
[docs]
@staticmethod
def get_power_law_for_sim_telarray_histograms(primary):
"""
Get the power law index for sim_telarray.
Events will be filled in histograms in sim_telarray with a weight according to
the difference between this exponent and the one used for the shower simulations.
Parameters
----------
primary: str
Primary particle.
Returns
-------
float
Power law index.
"""
power_laws = {
"gamma": 2.5,
"electron": 3.3,
}
if primary.name in power_laws:
return power_laws[primary.name]
return 2.68
