From 877f0b6e0d7bffcb5375c2ee4dce77d3e6e90f87 Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Mon, 25 Sep 2023 22:23:01 +0200 Subject: [PATCH 01/10] Add pre-commit configuration file. --- .pre-commit-config.yaml | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) create mode 100644 .pre-commit-config.yaml diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml new file mode 100644 index 000000000..9ec90cc85 --- /dev/null +++ b/.pre-commit-config.yaml @@ -0,0 +1,17 @@ +repos: + # https://pycqa.github.io/isort/docs/configuration/black_compatibility.html#integration-with-pre-commit + - repo: https://github.com/pycqa/isort + rev: 5.12.0 + hooks: + - id: isort + args: ["--profile", "black", "--filter-files"] + - repo: https://github.com/psf/black + rev: 22.6.0 + hooks: + - id: black-jupyter + # https://black.readthedocs.io/en/stable/guides/using_black_with_other_tools.html?highlight=other%20tools#flake8 + - repo: https://github.com/PyCQA/flake8 + rev: 4.0.1 + hooks: + - id: flake8 + args: [--max-line-length=88, "--extend-ignore=E203,E712"] From 2a3fe8597edb91a9e445d2d8a035046b8b4ce764 Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Mon, 25 Sep 2023 22:23:50 +0200 Subject: [PATCH 02/10] Changes by running pre-commit. --- magicctapipe/__init__.py | 7 +-- magicctapipe/conftest.py | 10 ++-- magicctapipe/image/__init__.py | 7 +-- magicctapipe/image/cleaning.py | 8 +-- magicctapipe/image/leakage.py | 1 - magicctapipe/image/muons/__init__.py | 4 +- magicctapipe/image/muons/muon_analysis.py | 2 +- magicctapipe/io/gadf.py | 3 +- magicctapipe/io/io.py | 3 +- magicctapipe/io/tests/test_gadf.py | 14 +++-- magicctapipe/io/tests/test_io.py | 20 +++---- magicctapipe/io/tests/test_io_monly.py | 22 +++---- magicctapipe/irfs/__init__.py | 12 ++-- magicctapipe/irfs/utils.py | 12 ++-- magicctapipe/reco/__init__.py | 52 +++++----------- magicctapipe/reco/classifier_utils.py | 10 ++-- magicctapipe/reco/direction_utils.py | 2 +- magicctapipe/reco/energy_utils.py | 3 +- magicctapipe/reco/estimators.py | 1 + magicctapipe/reco/event_processing.py | 33 +++++----- magicctapipe/reco/stereo.py | 3 +- magicctapipe/scripts/__init__.py | 6 +- .../scripts/comparison/cleaning_analysis.py | 29 ++++----- .../scripts/comparison/cleaning_comp.py | 39 +++++------- .../comparison/get_images_and_pixel_info.py | 30 +++++----- .../lst1_magic/lst1_magic_create_irf.py | 3 +- .../lst1_magic_dl1_stereo_to_dl2.py | 1 + .../lst1_magic/lst1_magic_dl2_to_dl3.py | 21 +++---- .../lst1_magic_event_coincidence.py | 5 +- .../lst1_magic/lst1_magic_mc_dl0_to_dl1.py | 3 +- .../lst1_magic/lst1_magic_mc_muon_analysis.py | 16 ++--- .../lst1_magic/lst1_magic_stereo_reco.py | 1 + .../lst1_magic/lst1_magic_train_rfs.py | 1 + .../scripts/lst1_magic/magic_calib_to_dl1.py | 9 ++- .../muon_analysis_LST_or_MAGIC_data.py | 5 +- .../scripts/mars/convert_superstar_to_dl1.py | 20 +++---- .../scripts/mars/effective_area_melibea.py | 6 +- .../scripts/mars/magic_calibrated_to_dl1.py | 32 ++++------ .../scripts/mars/mars_images_to_hdf5.py | 8 +-- .../performance/apply_rfs_MAGIC_LST.py | 15 ++--- .../performance/make_irfs_MAGIC_LST.py | 54 ++++++++--------- .../performance/stereo_reco_MAGIC_LST.py | 30 ++++------ .../performance/train_rfs_MAGIC_LST.py | 58 +++++++++--------- magicctapipe/utils/__init__.py | 60 +++++++------------ magicctapipe/utils/badpixels.py | 2 +- magicctapipe/utils/filedir.py | 9 +-- magicctapipe/utils/gti.py | 2 +- magicctapipe/utils/plot.py | 1 + magicctapipe/utils/tels.py | 5 +- magicctapipe/utils/utils.py | 2 + notebooks/check_event_coincidence.ipynb | 4 +- setup.py | 29 ++++----- 52 files changed, 329 insertions(+), 406 deletions(-) diff --git a/magicctapipe/__init__.py b/magicctapipe/__init__.py index 5cf69e71c..a51ad71a3 100644 --- a/magicctapipe/__init__.py +++ b/magicctapipe/__init__.py @@ -1,10 +1,5 @@ +from . import image, irfs, reco, scripts, utils from .version import __version__ -from . import image -from . import irfs -from . import reco -from . import scripts -from . import utils - __all__ = [ "image", diff --git a/magicctapipe/conftest.py b/magicctapipe/conftest.py index c8c2c51ea..73f0c029c 100644 --- a/magicctapipe/conftest.py +++ b/magicctapipe/conftest.py @@ -1,11 +1,13 @@ -import pytest +import subprocess +from math import trunc + import numpy as np +import pandas as pd +import pytest from astropy.io.misc.hdf5 import write_table_hdf5 from astropy.table import Table -import pandas as pd -import subprocess -from math import trunc from ctapipe.utils.download import download_file_cached + from magicctapipe.utils import resource_file maxjoint = 13000 diff --git a/magicctapipe/image/__init__.py b/magicctapipe/image/__init__.py index 34fa8194d..39d0bcb74 100644 --- a/magicctapipe/image/__init__.py +++ b/magicctapipe/image/__init__.py @@ -1,13 +1,10 @@ from .cleaning import ( MAGICClean, PixelTreatment, - get_num_islands_MAGIC, clean_image_params, + get_num_islands_MAGIC, ) - -from .leakage import ( - get_leakage, -) +from .leakage import get_leakage __all__ = [ "MAGICClean", diff --git a/magicctapipe/image/cleaning.py b/magicctapipe/image/cleaning.py index 004ae139c..833f424f7 100644 --- a/magicctapipe/image/cleaning.py +++ b/magicctapipe/image/cleaning.py @@ -1,14 +1,10 @@ import copy import itertools + import numpy as np +from ctapipe.image import hillas_parameters, leakage_parameters, timing_parameters from scipy.sparse.csgraph import connected_components -from ctapipe.image import ( - hillas_parameters, - timing_parameters, - leakage_parameters, -) - __all__ = [ "MAGICClean", "PixelTreatment", diff --git a/magicctapipe/image/leakage.py b/magicctapipe/image/leakage.py index bf92bea1e..2484adcf0 100644 --- a/magicctapipe/image/leakage.py +++ b/magicctapipe/image/leakage.py @@ -1,5 +1,4 @@ import numpy as np - from ctapipe.containers import LeakageContainer __all__ = [ diff --git a/magicctapipe/image/muons/__init__.py b/magicctapipe/image/muons/__init__.py index c72049c19..18a9f4402 100644 --- a/magicctapipe/image/muons/__init__.py +++ b/magicctapipe/image/muons/__init__.py @@ -1,6 +1,4 @@ -from .muon_analysis import ( - perform_muon_analysis, -) +from .muon_analysis import perform_muon_analysis __all__ = [ "perform_muon_analysis", diff --git a/magicctapipe/image/muons/muon_analysis.py b/magicctapipe/image/muons/muon_analysis.py index 26c3eae14..9f0910954 100644 --- a/magicctapipe/image/muons/muon_analysis.py +++ b/magicctapipe/image/muons/muon_analysis.py @@ -1,5 +1,5 @@ import numpy as np -from lstchain.image.muon import tag_pix_thr, fill_muon_event, analyze_muon_event +from lstchain.image.muon import analyze_muon_event, fill_muon_event, tag_pix_thr __all__ = [ "perform_muon_analysis", diff --git a/magicctapipe/io/gadf.py b/magicctapipe/io/gadf.py index f51b8807d..e27c94da7 100644 --- a/magicctapipe/io/gadf.py +++ b/magicctapipe/io/gadf.py @@ -9,10 +9,11 @@ from astropy.io import fits from astropy.table import QTable from astropy.time import Time +from pyirf.binning import split_bin_lo_hi + from magicctapipe import __version__ from magicctapipe.io.io import TEL_COMBINATIONS from magicctapipe.utils.functions import HEIGHT_ORM, LAT_ORM, LON_ORM -from pyirf.binning import split_bin_lo_hi __all__ = [ "create_gh_cuts_hdu", diff --git a/magicctapipe/io/io.py b/magicctapipe/io/io.py index 72d496795..b8c6dd11f 100644 --- a/magicctapipe/io/io.py +++ b/magicctapipe/io/io.py @@ -17,11 +17,12 @@ from ctapipe.coordinates import CameraFrame from ctapipe.instrument import SubarrayDescription from lstchain.reco.utils import add_delta_t_key -from magicctapipe.utils import calculate_mean_direction, transform_altaz_to_radec from pyirf.binning import join_bin_lo_hi from pyirf.simulations import SimulatedEventsInfo from pyirf.utils import calculate_source_fov_offset, calculate_theta +from magicctapipe.utils import calculate_mean_direction, transform_altaz_to_radec + __all__ = [ "format_object", "get_stereo_events", diff --git a/magicctapipe/io/tests/test_gadf.py b/magicctapipe/io/tests/test_gadf.py index f20045368..5853e08a6 100644 --- a/magicctapipe/io/tests/test_gadf.py +++ b/magicctapipe/io/tests/test_gadf.py @@ -1,15 +1,17 @@ +from math import isclose + +import astropy.units as u +import numpy as np +import pytest +from astropy.table import QTable + +from magicctapipe import __version__ from magicctapipe.io.gadf import ( create_event_hdu, create_gh_cuts_hdu, create_gti_hdu, create_pointing_hdu, ) -from magicctapipe import __version__ -import numpy as np -import pytest -import astropy.units as u -from astropy.table import QTable -from math import isclose @pytest.fixture(scope="session") diff --git a/magicctapipe/io/tests/test_io.py b/magicctapipe/io/tests/test_io.py index e6cdfe9a3..3ff6b22f0 100644 --- a/magicctapipe/io/tests/test_io.py +++ b/magicctapipe/io/tests/test_io.py @@ -1,20 +1,20 @@ +import numpy as np +import pandas as pd +import pytest + from magicctapipe.io.io import ( format_object, get_dl2_mean, get_stereo_events, - load_train_data_files, - load_mc_dl2_data_file, + load_dl2_data_file, load_irf_files, - save_pandas_data_in_table, - load_magic_dl1_data_files, load_lst_dl1_data_file, - load_dl2_data_file, + load_magic_dl1_data_files, + load_mc_dl2_data_file, + load_train_data_files, + save_pandas_data_in_table, ) -import pytest -import numpy as np -import pandas as pd - def test_format_object(): """ @@ -163,7 +163,6 @@ def test_load_mc_dl2_data_file_opt(p_dl2, gamma_dl2): assert np.all(data_s["combo_type"] > 0) - def test_load_mc_dl2_data_file_exc(p_dl2, gamma_dl2): """ Check on event_type exceptions @@ -382,7 +381,6 @@ def test_load_dl2_data_file_opt(real_dl2): assert np.all(data_s["combo_type"] > 0) - def test_load_dl2_data_file_exc(real_dl2): """ Check on event_type exceptions diff --git a/magicctapipe/io/tests/test_io_monly.py b/magicctapipe/io/tests/test_io_monly.py index c498b87c4..a8f8d8504 100644 --- a/magicctapipe/io/tests/test_io_monly.py +++ b/magicctapipe/io/tests/test_io_monly.py @@ -1,20 +1,20 @@ +import numpy as np +import pandas as pd +import pytest + from magicctapipe.io.io import ( format_object, get_dl2_mean, get_stereo_events, - load_train_data_files, - load_mc_dl2_data_file, + load_dl2_data_file, load_irf_files, - save_pandas_data_in_table, - load_magic_dl1_data_files, load_lst_dl1_data_file, - load_dl2_data_file, + load_magic_dl1_data_files, + load_mc_dl2_data_file, + load_train_data_files, + save_pandas_data_in_table, ) -import pytest -import numpy as np -import pandas as pd - def test_format_object(): """ @@ -356,7 +356,9 @@ def test_load_dl2_data_file(real_dl2_monly): Checks on default loading """ for file in real_dl2_monly.glob("*"): - data, on, dead = load_dl2_data_file(str(file), "width>0", "magic_only", "simple") + data, on, dead = load_dl2_data_file( + str(file), "width>0", "magic_only", "simple" + ) assert "pointing_alt" in data.colnames assert "timestamp" in data.colnames assert data["reco_energy"].unit == "TeV" diff --git a/magicctapipe/irfs/__init__.py b/magicctapipe/irfs/__init__.py index a392eb1a6..02b5f509d 100644 --- a/magicctapipe/irfs/__init__.py +++ b/magicctapipe/irfs/__init__.py @@ -1,16 +1,16 @@ from .utils import ( - read_simu_info_mcp_sum_num_showers, convert_simu_info_mcp_to_pyirf, - read_dl2_mcp_to_pyirf_MAGIC_LST_list, - plot_sensitivity, - plot_en_res_bias, - plot_en_res_resolution, plot_ang_res, plot_effective_area, + plot_en_res_bias, + plot_en_res_resolution, plot_gamma_eff_gh, plot_irfs_MAGIC_LST, - plot_MARS_sensitivity, plot_MAGIC_reference_sensitivity, + plot_MARS_sensitivity, + plot_sensitivity, + read_dl2_mcp_to_pyirf_MAGIC_LST_list, + read_simu_info_mcp_sum_num_showers, ) __all__ = [ diff --git a/magicctapipe/irfs/utils.py b/magicctapipe/irfs/utils.py index 6c16b1c6e..f44c54796 100644 --- a/magicctapipe/irfs/utils.py +++ b/magicctapipe/irfs/utils.py @@ -1,21 +1,19 @@ -import os import glob +import os +import astropy.units as u +import matplotlib.pylab as plt import numpy as np import pandas as pd from astropy import table -import astropy.units as u from astropy.io import fits from astropy.table import QTable - from pyirf.simulations import SimulatedEventsInfo from magicctapipe.utils.filedir import read_mc_header -from magicctapipe.utils.plot import save_plt, load_default_plot_settings +from magicctapipe.utils.plot import load_default_plot_settings, save_plt from magicctapipe.utils.utils import print_title -import matplotlib.pylab as plt - __all__ = [ "read_simu_info_mcp_sum_num_showers", "convert_simu_info_mcp_to_pyirf", @@ -224,7 +222,7 @@ def plot_sensitivity(data, unit, label, ax=None, **kwargs): e = data["reco_energy_center"] s_mc = e**2 * data["flux_sensitivity"] e_low, e_high = data["reco_energy_low"], data["reco_energy_high"] - ax_ = ax if ax != None else plt + ax_ = ax if ax is not None else plt plt_ = ax_.errorbar( e.to_value(u.GeV), s_mc.to_value(unit), diff --git a/magicctapipe/reco/__init__.py b/magicctapipe/reco/__init__.py index 09d9f740e..a906c21a2 100644 --- a/magicctapipe/reco/__init__.py +++ b/magicctapipe/reco/__init__.py @@ -1,54 +1,32 @@ from .classifier_utils import ( GetHist_classifier, + check_train_test_intersections_classifier, evaluate_performance_classifier, get_weights_classifier, - print_par_imp_classifier, load_init_data_classifier, - check_train_test_intersections_classifier, -) - -from .direction_utils import ( - compute_separation_angle_direction, -) - -from .energy_utils import ( - GetHist2D_energy, - evaluate_performance_energy, - plot_migmatrix, + print_par_imp_classifier, ) - -from .event_processing import ( - RegressorClassifierBase, - # EnergyRegressor, - HillasFeatureSelector, - EventFeatureSelector, - EventFeatureTargetSelector, - EventProcessor, +from .direction_utils import compute_separation_angle_direction +from .energy_utils import GetHist2D_energy, evaluate_performance_energy, plot_migmatrix +from .estimators import DispRegressor, EnergyRegressor, EventClassifier +from .event_processing import ( # EnergyRegressor, + DirectionEstimatorPandas, + DirectionStereoEstimatorPandas, EnergyEstimator, EnergyEstimatorPandas, - DirectionEstimatorPandas, EventClassifierPandas, - DirectionStereoEstimatorPandas, + EventFeatureSelector, + EventFeatureTargetSelector, + EventProcessor, + HillasFeatureSelector, + RegressorClassifierBase, ) - from .global_utils import ( + check_train_test_intersections, compute_event_weights, get_weights_mc_dir_class, - check_train_test_intersections, -) - -from .estimators import ( - DispRegressor, - EnergyRegressor, - EventClassifier, -) - -from .stereo import ( - write_hillas, - check_write_stereo, - check_stereo, - write_stereo, ) +from .stereo import check_stereo, check_write_stereo, write_hillas, write_stereo __all__ = [ "GetHist_classifier", diff --git a/magicctapipe/reco/classifier_utils.py b/magicctapipe/reco/classifier_utils.py index 6f6db717d..9d143bf11 100644 --- a/magicctapipe/reco/classifier_utils.py +++ b/magicctapipe/reco/classifier_utils.py @@ -1,16 +1,16 @@ -import pandas as pd -import numpy as np import glob import os -import sklearn.metrics +import numpy as np +import pandas as pd +import sklearn.metrics -from magicctapipe.utils.utils import info_message +from magicctapipe.reco.global_utils import check_train_test_intersections from magicctapipe.utils.filedir import ( load_dl1_data_stereo_list, load_dl1_data_stereo_list_selected, ) -from magicctapipe.reco.global_utils import check_train_test_intersections +from magicctapipe.utils.utils import info_message __all__ = [ "GetHist_classifier", diff --git a/magicctapipe/reco/direction_utils.py b/magicctapipe/reco/direction_utils.py index c6899e1a0..9aff275a7 100644 --- a/magicctapipe/reco/direction_utils.py +++ b/magicctapipe/reco/direction_utils.py @@ -1,6 +1,6 @@ import pandas as pd from astropy import units as u -from astropy.coordinates import SkyCoord, AltAz +from astropy.coordinates import AltAz, SkyCoord from magicctapipe.utils.tels import get_tel_ids_dl1 diff --git a/magicctapipe/reco/energy_utils.py b/magicctapipe/reco/energy_utils.py index 1c8b389a9..8cebc7b24 100644 --- a/magicctapipe/reco/energy_utils.py +++ b/magicctapipe/reco/energy_utils.py @@ -1,8 +1,7 @@ -import numpy as np import matplotlib.pyplot as plt +import numpy as np from matplotlib import colors - __all__ = [ "GetHist2D_energy", "evaluate_performance_energy", diff --git a/magicctapipe/reco/estimators.py b/magicctapipe/reco/estimators.py index 5bd02a4c3..91d993b62 100644 --- a/magicctapipe/reco/estimators.py +++ b/magicctapipe/reco/estimators.py @@ -7,6 +7,7 @@ import numpy as np import pandas as pd import sklearn.ensemble + from magicctapipe.io.io import TEL_NAMES __all__ = ["EnergyRegressor", "DispRegressor", "EventClassifier"] diff --git a/magicctapipe/reco/event_processing.py b/magicctapipe/reco/event_processing.py index 397441114..b6549bd99 100644 --- a/magicctapipe/reco/event_processing.py +++ b/magicctapipe/reco/event_processing.py @@ -1,26 +1,23 @@ -import re import itertools -import numpy as np -from copy import deepcopy - +import re from abc import ABC, abstractmethod +from copy import deepcopy +import joblib +import numpy as np +import pandas as pd import sklearn import sklearn.ensemble -from sklearn.ensemble import RandomForestRegressor -from sklearn.preprocessing import StandardScaler -import pandas as pd -import joblib - -from ctapipe.image import tailcuts_clean, hillas_parameters -from ctapipe.io import EventSource -from ctapipe.containers import ReconstructedContainer, ReconstructedEnergyContainer -from ctapipe.reco.reco_algorithms import TooFewTelescopesException -from ctapipe.coordinates import CameraFrame, TelescopeFrame - -from astropy.coordinates import AltAz, SkyCoord from astropy import units as u +from astropy.coordinates import AltAz, SkyCoord from astropy.coordinates.angle_utilities import angular_separation +from ctapipe.containers import ReconstructedContainer, ReconstructedEnergyContainer +from ctapipe.coordinates import CameraFrame, TelescopeFrame +from ctapipe.image import hillas_parameters, tailcuts_clean +from ctapipe.io import EventSource +from ctapipe.reco.reco_algorithms import TooFewTelescopesException +from sklearn.ensemble import RandomForestRegressor +from sklearn.preprocessing import StandardScaler __all__ = [ "RegressorClassifierBase", @@ -420,7 +417,7 @@ def show_importances(self): plt.title(cam_id) try: importances = model.feature_importances_ - except: + except Exception: plt.gca().axis("off") continue bins = range(importances.shape[0]) @@ -895,7 +892,7 @@ def _dl1_process(self, event): event.dl1.tel[tel_id].hillas_params = hillas_parameters( camera, event_image_cleaned ) - except: + except Exception: print("Failed") pass diff --git a/magicctapipe/reco/stereo.py b/magicctapipe/reco/stereo.py index f2e841216..422cddee7 100644 --- a/magicctapipe/reco/stereo.py +++ b/magicctapipe/reco/stereo.py @@ -1,6 +1,5 @@ -import numpy as np import astropy.units as u - +import numpy as np from ctapipe.core.container import Container, Field from magicctapipe.utils.tels import tel_ids_2_num diff --git a/magicctapipe/scripts/__init__.py b/magicctapipe/scripts/__init__.py index a645f347f..8fac2f117 100644 --- a/magicctapipe/scripts/__init__.py +++ b/magicctapipe/scripts/__init__.py @@ -4,13 +4,13 @@ dl1_stereo_to_dl2, dl2_to_dl3, event_coincidence, + magic_calib_to_dl1, mc_dl0_to_dl1, + merge_hdf_files, stereo_reconstruction, - train_energy_regressor, train_disp_regressor, + train_energy_regressor, train_event_classifier, - magic_calib_to_dl1, - merge_hdf_files, ) __all__ = [ diff --git a/magicctapipe/scripts/comparison/cleaning_analysis.py b/magicctapipe/scripts/comparison/cleaning_analysis.py index 7ce5fbca4..b5c71360d 100644 --- a/magicctapipe/scripts/comparison/cleaning_analysis.py +++ b/magicctapipe/scripts/comparison/cleaning_analysis.py @@ -1,22 +1,17 @@ -import glob, sys +import argparse +import glob import logging +import sys import numpy as np - -np.set_printoptions(threshold=sys.maxsize) import scipy - -from ctapipe_io_magic import MAGICEventSource - -from ctapipe.instrument import CameraGeometry from ctapipe.image import hillas_parameters +from ctapipe.instrument import CameraGeometry +from ctapipe_io_magic import MAGICEventSource -import argparse - -from magicctapipe.utils import MAGIC_Badpixels +from magicctapipe.utils import MAGIC_Badpixels, MAGIC_Cleaning -# from magicctapipe.utils import bad_pixel_treatment -from magicctapipe.utils import MAGIC_Cleaning +np.set_printoptions(threshold=sys.maxsize) logger = logging.getLogger(__name__) @@ -64,21 +59,21 @@ def main(): all_events = [] if args.mc: - """MC FILES""" + # MC FILES mars_file_mask = "/remote/ceph/group/magic/MAGIC-LST/MCs/MAGIC/ST.03.07/za05to35/Train_sample/1.Calibrated/GA_M1*root" mars_file_list = glob.glob( mars_file_mask - ) # Here makes array which contains files matching the input condition (GA_M1*root). + ) # Here makes array which contains files matching the input condition (GA_M1*root). file_list = list(filter(lambda name: "123" in name, mars_file_list)) print("Number of files : %s" % len(file_list)) else: - """DATA FILES""" + # DATA FILES mars_file_mask = ( "/Users/moritz/Documents/MAGIC-data/Crab/Data/Sorcerer/2020*.root" ) mars_file_list = glob.glob( mars_file_mask - ) # Here makes array which contains files matching the input condition (GA_M1*root). + ) # Here makes array which contains files matching the input condition (GA_M1*root). file_list = list(filter(lambda name: "0" in name, mars_file_list)) print("Number of files : %s" % len(file_list)) @@ -115,7 +110,7 @@ def main(): if scipy.any(event_image_cleaned): try: hillas_params = hillas_parameters(mars_camera, event_image_cleaned) - except: + except Exception: continue logger.debug( diff --git a/magicctapipe/scripts/comparison/cleaning_comp.py b/magicctapipe/scripts/comparison/cleaning_comp.py index dedc0bbe1..43bd32212 100644 --- a/magicctapipe/scripts/comparison/cleaning_comp.py +++ b/magicctapipe/scripts/comparison/cleaning_comp.py @@ -1,28 +1,21 @@ -import glob, sys -import uproot +import argparse import copy -from texttable import Texttable +import csv +import glob +import sys import numpy as np - -np.set_printoptions(threshold=sys.maxsize) - -from ctapipe_io_magic import MAGICEventSource - +import pylab as plt +import uproot from astropy import units as u - -from ctapipe.instrument import CameraGeometry from ctapipe.image import hillas_parameters +from ctapipe.instrument import CameraGeometry +from ctapipe_io_magic import MAGICEventSource +from texttable import Texttable -import pylab as plt - -import argparse -import csv - -from magicctapipe.utils import MAGIC_Badpixels +from magicctapipe.utils import MAGIC_Badpixels, MAGIC_Cleaning -# from magicctapipe.utils import bad_pixel_treatment -from magicctapipe.utils import MAGIC_Cleaning +np.set_printoptions(threshold=sys.maxsize) class bcolors: @@ -187,21 +180,21 @@ def main(): AberrationCorrection = 1.0713 if args.mc: - """MC FILES""" + # MC FILES mars_file_mask = "/remote/ceph/group/magic/MAGIC-LST/MCs/MAGIC/ST.03.07/za05to35/Train_sample/1.Calibrated/GA_M1*root" mars_file_list = glob.glob( mars_file_mask - ) # Here makes array which contains files matching the input condition (GA_M1*root). + ) # Here makes array which contains files matching the input condition (GA_M1*root). file_list = list(filter(lambda name: "123" in name, mars_file_list)) print("Number of files : %s" % len(file_list)) else: - """DATA FILES""" + # DATA FILES mars_file_mask = ( "/home/iwsatlas1/damgreen/CTA/MAGIC_Cleaning/datafile/Calib/*.root" ) mars_file_list = glob.glob( mars_file_mask - ) # Here makes array which contains files matching the input condition (GA_M1*root). + ) # Here makes array which contains files matching the input condition (GA_M1*root). file_list = list(filter(lambda name: "0" in name, mars_file_list)) print("Number of files : %s" % len(file_list)) @@ -242,7 +235,7 @@ def main(): try: hillas_params = hillas_parameters(mars_camera, event_image_cleaned) - except: + except Exception: continue # All of the comparisons between star and ctapipe diff --git a/magicctapipe/scripts/comparison/get_images_and_pixel_info.py b/magicctapipe/scripts/comparison/get_images_and_pixel_info.py index 000cf2495..bb749e99b 100644 --- a/magicctapipe/scripts/comparison/get_images_and_pixel_info.py +++ b/magicctapipe/scripts/comparison/get_images_and_pixel_info.py @@ -2,23 +2,21 @@ # coding: utf-8 import argparse -import yaml + +import matplotlib +import matplotlib.pyplot as plt +import numpy as np import pandas as pd import uproot -from ctapipe_io_magic import MAGICEventSource -from ctapipe.visualization import CameraDisplay -from ctapipe.instrument import CameraGeometry +import yaml from ctapipe.image import hillas_parameters -from ctapipe.image.timing import timing_parameters -from ctapipe.instrument import CameraDescription -from utils import MAGIC_Badpixels -from utils import MAGIC_Cleaning -import matplotlib.pyplot as plt -import matplotlib -import numpy as np - from ctapipe.image.morphology import number_of_islands +from ctapipe.image.timing import timing_parameters +from ctapipe.instrument import CameraDescription, CameraGeometry +from ctapipe.visualization import CameraDisplay +from ctapipe_io_magic import MAGICEventSource from mars_images_to_hdf5 import read_images +from utils import MAGIC_Badpixels, MAGIC_Cleaning # ---------------------------------------------------------------- @@ -28,7 +26,7 @@ # add argparser arguments arg_parser = argparse.ArgumentParser( description=""" -This tool compares the images of events processed by MARS and the magic-cta-pipeline. +This tool compares the images of events processed by MARS and the magic-cta-pipeline. The output is a png file with the camera images and a hdf5 file that contains the pixel information. """ ) @@ -256,7 +254,7 @@ def new_camera_geometry(camera_geom): event_image = event.dl1.tel[tel_id].image event_pulse_time = event.dl1.tel[tel_id].peak_time - """ + """ find different catgories of pixels ------ bad: pixel that is unsuitable, defined subrunwise @@ -265,8 +263,8 @@ def new_camera_geometry(camera_geom): get badrms pixels via badrmspixel_mask get badrms pixel indices via bad_pixel_indices dead: this pixels is not working at all - get dead pixels via deadpixel_mask - get dead pixel indices via dead_pixel_indices + get dead pixels via deadpixel_mask + get dead pixel indices via dead_pixel_indices unsuitable: union of badrms and dead pixels get unsuitable pixels via unsuitable_mask --------- diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_create_irf.py b/magicctapipe/scripts/lst1_magic/lst1_magic_create_irf.py index 58be8c743..8c7d21930 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_create_irf.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_create_irf.py @@ -52,7 +52,6 @@ from astropy import units as u from astropy.io import fits from astropy.table import QTable, vstack -from magicctapipe.io import create_gh_cuts_hdu, format_object, load_mc_dl2_data_file from pyirf.cuts import calculate_percentile_cut, evaluate_binned_cut from pyirf.io.gadf import ( create_aeff2d_hdu, @@ -75,6 +74,8 @@ calculate_event_weights, ) +from magicctapipe.io import create_gh_cuts_hdu, format_object, load_mc_dl2_data_file + __all__ = ["create_irf"] logger = logging.getLogger(__name__) diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py b/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py index d3d129ecb..751111279 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py @@ -29,6 +29,7 @@ from astropy.coordinates import AltAz, SkyCoord, angular_separation from ctapipe.coordinates import TelescopeFrame from ctapipe.instrument import SubarrayDescription + from magicctapipe.io import get_stereo_events, save_pandas_data_in_table from magicctapipe.io.io import TEL_COMBINATIONS from magicctapipe.reco import DispRegressor, EnergyRegressor, EventClassifier diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_dl2_to_dl3.py b/magicctapipe/scripts/lst1_magic/lst1_magic_dl2_to_dl3.py index 3007927e7..07ed050a4 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_dl2_to_dl3.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_dl2_to_dl3.py @@ -33,15 +33,6 @@ from astropy.coordinates.angle_utilities import angular_separation from astropy.io import fits from astropy.table import QTable -from magicctapipe.io import ( - create_event_hdu, - create_gh_cuts_hdu, - create_gti_hdu, - create_pointing_hdu, - format_object, - load_dl2_data_file, - load_irf_files, -) from pyirf.cuts import evaluate_binned_cut from pyirf.interpolation import interpolate_effective_area_per_energy_and_fov from pyirf.io import ( @@ -54,6 +45,16 @@ from pyirf.utils import cone_solid_angle from scipy.interpolate import griddata +from magicctapipe.io import ( + create_event_hdu, + create_gh_cuts_hdu, + create_gti_hdu, + create_pointing_hdu, + format_object, + load_dl2_data_file, + load_irf_files, +) + __all__ = ["dl2_to_dl3"] logger = logging.getLogger(__name__) @@ -142,7 +143,7 @@ def dl2_to_dl3(input_file_dl2, input_dir_irf, output_dir, config): f"\nTarget point: {target_point.round(5).tolist()} with scheme: {scheme}" ) - if 'max_distance' in config_dl3: + if "max_distance" in config_dl3: max_distance = u.Quantity(config_dl3.pop("max_distance")) logger.info(f"selecting only nodes up to {max_distance} from the data") diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py b/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py index 8f3fbf771..bf63c7f32 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py @@ -27,7 +27,7 @@ accidental coincidence rate as much as possible by keeping the number of actual coincident events. -Please note that the time offset depends on the date of observations +Please note that the time offset depends on the date of observations as summarized below: * before June 12 2021: -3.1 us * June 13 2021 to Feb 28 2023: -6.5 us @@ -35,7 +35,7 @@ * April 13 2023 to August 2023: -25.1 us * after Sep 11 2023 : -6.2 us By default, pre offset search is performed using large shower events. -The possible time offset is found among all possible combinations of +The possible time offset is found among all possible combinations of time offsets using those events. Finally, the time offset scan is performed around the possible offset found by the pre offset search. Instead of that, you can also define the offset scan range in the configuration file. @@ -60,6 +60,7 @@ import yaml from astropy import units as u from ctapipe.instrument import SubarrayDescription + from magicctapipe.io import ( get_stereo_events, load_lst_dl1_data_file, diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py b/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py index 41cf519ee..cc0a7582f 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py @@ -53,10 +53,11 @@ random_psf_smearer, set_numba_seed, ) +from traitlets.config import Config + from magicctapipe.image import MAGICClean from magicctapipe.io import SimEventInfoContainer, format_object from magicctapipe.utils import calculate_disp, calculate_impact -from traitlets.config import Config __all__ = ["mc_dl0_to_dl1"] diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_mc_muon_analysis.py b/magicctapipe/scripts/lst1_magic/lst1_magic_mc_muon_analysis.py index 72c5c5be9..734badd6b 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_mc_muon_analysis.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_mc_muon_analysis.py @@ -10,25 +10,27 @@ --config-file ./config.yaml """ +import argparse +import logging import re import time -import yaml -import logging -import argparse -import numpy as np from pathlib import Path + +import numpy as np +import yaml from astropy.table import Table -from traitlets.config import Config -from ctapipe.io import EventSource from ctapipe.calib import CameraCalibrator from ctapipe.image import ( ImageExtractor, - tailcuts_clean, apply_time_delta_cleaning, number_of_islands, + tailcuts_clean, ) +from ctapipe.io import EventSource from lstchain.image.cleaning import apply_dynamic_cleaning from lstchain.image.muon import create_muon_table +from traitlets.config import Config + from magicctapipe.image import MAGICClean from magicctapipe.image.muons import perform_muon_analysis diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py b/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py index accec1cea..c0dde1b16 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py @@ -44,6 +44,7 @@ ) from ctapipe.instrument import SubarrayDescription from ctapipe.reco import HillasReconstructor + from magicctapipe.io import format_object, get_stereo_events, save_pandas_data_in_table from magicctapipe.utils import calculate_impact, calculate_mean_direction diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py b/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py index 8afaaf3b2..b43494fd8 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py @@ -40,6 +40,7 @@ import numpy as np import pandas as pd import yaml + from magicctapipe.io import format_object, load_train_data_files from magicctapipe.io.io import GROUP_INDEX_TRAIN, TEL_NAMES from magicctapipe.reco import DispRegressor, EnergyRegressor, EventClassifier diff --git a/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py b/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py index 645c10053..120058e8d 100644 --- a/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py +++ b/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py @@ -55,6 +55,7 @@ from ctapipe.instrument import SubarrayDescription from ctapipe.io import HDF5TableWriter from ctapipe_io_magic import MAGICEventSource + from magicctapipe.image import MAGICClean from magicctapipe.io import RealEventInfoContainer, SimEventInfoContainer, format_object from magicctapipe.utils import calculate_disp, calculate_impact @@ -94,7 +95,9 @@ def magic_calib_to_dl1(input_file, output_dir, config, max_events, process_run=F # Load the input file logger.info(f"\nInput file: {input_file}") - event_source = MAGICEventSource(input_file, process_run=process_run, max_events=max_events) + event_source = MAGICEventSource( + input_file, process_run=process_run, max_events=max_events + ) is_simulation = event_source.is_simulation logger.info(f"\nIs simulation: {is_simulation}") @@ -393,7 +396,9 @@ def main(): config = yaml.safe_load(f) # Process the input data - magic_calib_to_dl1(args.input_file, args.output_dir, config, args.max_events, args.process_run) + magic_calib_to_dl1( + args.input_file, args.output_dir, config, args.max_events, args.process_run + ) logger.info("\nDone.") diff --git a/magicctapipe/scripts/lst1_magic/muon_analysis_LST_or_MAGIC_data.py b/magicctapipe/scripts/lst1_magic/muon_analysis_LST_or_MAGIC_data.py index 924e3e1bb..cfbd067f9 100644 --- a/magicctapipe/scripts/lst1_magic/muon_analysis_LST_or_MAGIC_data.py +++ b/magicctapipe/scripts/lst1_magic/muon_analysis_LST_or_MAGIC_data.py @@ -14,14 +14,15 @@ --config-file ./config.yaml """ import argparse -import yaml import logging -import numpy as np from pathlib import Path +import numpy as np +import yaml from astropy.table import Table from ctapipe.io import EventSource from lstchain.image.muon import create_muon_table + from magicctapipe.image import MAGICClean from magicctapipe.image.muons import perform_muon_analysis diff --git a/magicctapipe/scripts/mars/convert_superstar_to_dl1.py b/magicctapipe/scripts/mars/convert_superstar_to_dl1.py index 6d2e96fb3..3c5052eb5 100644 --- a/magicctapipe/scripts/mars/convert_superstar_to_dl1.py +++ b/magicctapipe/scripts/mars/convert_superstar_to_dl1.py @@ -1,32 +1,30 @@ #!/usr/bin/env python +import argparse import re import sys -import argparse from pathlib import Path -import uproot +import astropy.units as u import numpy as np - +import uproot from astropy.table import QTable, vstack -import astropy.units as u - from ctapipe.containers import ( CameraHillasParametersContainer, - LeakageContainer, CameraTimingParametersContainer, + LeakageContainer, ReconstructedGeometryContainer, ) -from ctapipe.io import HDF5TableWriter -from ctapipe.core.container import Container, Field from ctapipe.coordinates import CameraFrame +from ctapipe.core.container import Container, Field from ctapipe.instrument import ( - TelescopeDescription, - SubarrayDescription, - OpticsDescription, CameraDescription, CameraReadout, + OpticsDescription, + SubarrayDescription, + TelescopeDescription, ) +from ctapipe.io import HDF5TableWriter magic_optics = OpticsDescription( "MAGIC", diff --git a/magicctapipe/scripts/mars/effective_area_melibea.py b/magicctapipe/scripts/mars/effective_area_melibea.py index a7cd29009..74f2c2a4d 100644 --- a/magicctapipe/scripts/mars/effective_area_melibea.py +++ b/magicctapipe/scripts/mars/effective_area_melibea.py @@ -1,9 +1,7 @@ -import numpy as np - -from astropy.table import QTable, vstack import astropy.units as u - +import numpy as np import uproot +from astropy.table import QTable, vstack from pyirf.simulations import SimulatedEventsInfo melibea_columns = { diff --git a/magicctapipe/scripts/mars/magic_calibrated_to_dl1.py b/magicctapipe/scripts/mars/magic_calibrated_to_dl1.py index 658a9528f..9425e627b 100644 --- a/magicctapipe/scripts/mars/magic_calibrated_to_dl1.py +++ b/magicctapipe/scripts/mars/magic_calibrated_to_dl1.py @@ -2,27 +2,20 @@ # coding: utf-8 import argparse +import copy import glob import re -import yaml -import copy from pathlib import Path import numpy as np - -from ctapipe_io_magic import MAGICEventSource - +import yaml from ctapipe.containers import ( - IntensityStatisticsContainer, ImageParametersContainer, - TimingParametersContainer, + IntensityStatisticsContainer, PeakTimeStatisticsContainer, + TimingParametersContainer, ) - from ctapipe.coordinates import TelescopeFrame - -from ctapipe.io import DataWriter - from ctapipe.image import ( concentration_parameters, descriptive_statistics, @@ -30,16 +23,11 @@ morphology_parameters, timing_parameters, ) +from ctapipe.io import DataWriter +from ctapipe_io_magic import MAGICEventSource -from magicctapipe.image import ( - MAGICClean, - get_leakage, -) - -from magicctapipe.utils import ( - info_message, -) - +from magicctapipe.image import MAGICClean, get_leakage +from magicctapipe.utils import info_message DEFAULT_IMAGE_PARAMETERS = ImageParametersContainer() DEFAULT_TRUE_IMAGE_PARAMETERS = ImageParametersContainer() @@ -318,8 +306,8 @@ def magic_calibrated_to_dl1(input_mask, cleaning_config): ) try: - telescope_type = re.findall("(.*)[_\d]+", telescope)[0] - except: + telescope_type = re.findall("(.*)[_\\d]+", telescope)[0] + except Exception: ValueError( f'Can not recognize the telescope type from name "{telescope}"' ) diff --git a/magicctapipe/scripts/mars/mars_images_to_hdf5.py b/magicctapipe/scripts/mars/mars_images_to_hdf5.py index 745aee39f..daeea1ab1 100644 --- a/magicctapipe/scripts/mars/mars_images_to_hdf5.py +++ b/magicctapipe/scripts/mars/mars_images_to_hdf5.py @@ -7,18 +7,16 @@ by magic-cta-pipe. """ +import argparse import re import sys -import argparse from pathlib import Path import numpy as np -import uproot import tables - +import uproot from ctapipe.core.container import Container, Field -from ctapipe.io import HDF5TableWriter, HDF5TableReader - +from ctapipe.io import HDF5TableReader, HDF5TableWriter from ctapipe_io_magic import MARSDataLevel diff --git a/magicctapipe/scripts/performance/apply_rfs_MAGIC_LST.py b/magicctapipe/scripts/performance/apply_rfs_MAGIC_LST.py index eb2996cce..509b78778 100644 --- a/magicctapipe/scripts/performance/apply_rfs_MAGIC_LST.py +++ b/magicctapipe/scripts/performance/apply_rfs_MAGIC_LST.py @@ -1,24 +1,25 @@ # coding: utf-8 -import time import argparse import glob import os +import time + import pandas as pd from magicctapipe.reco.event_processing import ( - EnergyEstimatorPandas, DirectionEstimatorPandas, + EnergyEstimatorPandas, EventClassifierPandas, ) -from magicctapipe.utils.tels import check_tel_ids, get_array_tel_descriptions -from magicctapipe.utils.utils import print_title, info_message, print_elapsed_time from magicctapipe.utils.filedir import ( - load_cfg_file, - out_file_h5_reco, check_folder, + load_cfg_file, load_dl1_data_stereo, + out_file_h5_reco, ) +from magicctapipe.utils.tels import check_tel_ids, get_array_tel_descriptions +from magicctapipe.utils.utils import info_message, print_elapsed_time, print_title PARSER = argparse.ArgumentParser( description="Apply random forests. For stereo data.", @@ -169,7 +170,7 @@ def apply_rfs_stereo(config_file, only_mc_test, only_data_test): try: mc_ = pd.read_hdf(file, key="dl1/mc_header") mc_.to_hdf(out_file, key="dl2/mc_header") - except: + except Exception: mc_ = pd.read_hdf(file, key="dl2/mc_header") mc_.to_hdf(out_file, key="dl2/mc_header") except Exception as e: diff --git a/magicctapipe/scripts/performance/make_irfs_MAGIC_LST.py b/magicctapipe/scripts/performance/make_irfs_MAGIC_LST.py index 23aa7c243..d75396d3d 100644 --- a/magicctapipe/scripts/performance/make_irfs_MAGIC_LST.py +++ b/magicctapipe/scripts/performance/make_irfs_MAGIC_LST.py @@ -1,54 +1,50 @@ -import os -import time -import shutil +import argparse import logging import operator -import argparse +import os +import shutil +import time +import astropy.units as u import numpy as np from astropy import table -import astropy.units as u from astropy.io import fits - +from pyirf.benchmarks import angular_resolution, energy_bias_resolution from pyirf.binning import ( - create_bins_per_decade, add_overflow_bins, + create_bins_per_decade, create_histogram_table, ) +from pyirf.cut_optimization import optimize_gh_cut from pyirf.cuts import calculate_percentile_cut, evaluate_binned_cut -from pyirf.sensitivity import calculate_sensitivity, estimate_background -from pyirf.utils import calculate_theta, calculate_source_fov_offset -from pyirf.benchmarks import energy_bias_resolution, angular_resolution - -from pyirf.spectral import ( - calculate_event_weights, - PowerLaw, - CRAB_HEGRA, - IRFDOC_PROTON_SPECTRUM, - IRFDOC_ELECTRON_SPECTRUM, +from pyirf.io import ( + create_aeff2d_hdu, + create_background_2d_hdu, + create_energy_dispersion_hdu, + create_psf_table_hdu, + create_rad_max_hdu, ) -from pyirf.cut_optimization import optimize_gh_cut - from pyirf.irf import ( + background_2d, effective_area_per_energy, energy_dispersion, psf_table, - background_2d, ) - -from pyirf.io import ( - create_aeff2d_hdu, - create_psf_table_hdu, - create_energy_dispersion_hdu, - create_rad_max_hdu, - create_background_2d_hdu, +from pyirf.sensitivity import calculate_sensitivity, estimate_background +from pyirf.spectral import ( + CRAB_HEGRA, + IRFDOC_ELECTRON_SPECTRUM, + IRFDOC_PROTON_SPECTRUM, + PowerLaw, + calculate_event_weights, ) +from pyirf.utils import calculate_source_fov_offset, calculate_theta -from magicctapipe.utils.filedir import check_folder, load_cfg_file from magicctapipe.irfs.utils import ( plot_irfs_MAGIC_LST, read_dl2_mcp_to_pyirf_MAGIC_LST_list, ) +from magicctapipe.utils.filedir import check_folder, load_cfg_file from magicctapipe.utils.utils import print_elapsed_time, print_title PARSER = argparse.ArgumentParser( @@ -409,7 +405,7 @@ def make_irfs_MAGIC_LST(config_file): len(gammas_sel_en[gammas_sel_en["selected"] == True]) / len(gammas_sel_en) ) - except: + except Exception: gamma_efficiency["eff_gh"][i] = -1 gamma_efficiency["eff"][i] = -1 hdus.append(fits.BinTableHDU(gamma_efficiency, name="GAMMA_EFFICIENCY")) diff --git a/magicctapipe/scripts/performance/stereo_reco_MAGIC_LST.py b/magicctapipe/scripts/performance/stereo_reco_MAGIC_LST.py index 35a8558e8..b037cdb35 100644 --- a/magicctapipe/scripts/performance/stereo_reco_MAGIC_LST.py +++ b/magicctapipe/scripts/performance/stereo_reco_MAGIC_LST.py @@ -1,42 +1,38 @@ # coding: utf-8 -import os -import time +import argparse import copy import glob -import scipy +import os import select import sys -import argparse -import matplotlib.pyplot as plt +import time import astropy.units as u -from astropy.coordinates import SkyCoord, AltAz, EarthLocation +import matplotlib.pyplot as plt +import scipy +from astropy.coordinates import AltAz, EarthLocation, SkyCoord from astropy.time import Time - -from ctapipe.io import SimTelEventSource -from ctapipe.io import HDF5TableWriter from ctapipe.calib import CameraCalibrator +from ctapipe.containers import ImageParametersContainer from ctapipe.coordinates import TelescopeFrame from ctapipe.image.cleaning import tailcuts_clean from ctapipe.image.morphology import number_of_islands - +from ctapipe.io import HDF5TableWriter, SimTelEventSource from ctapipe.reco import HillasReconstructor -from ctapipe.containers import ImageParametersContainer from ctapipe.visualization import ArrayDisplay from magicctapipe.image import MAGICClean, clean_image_params, get_num_islands_MAGIC -from magicctapipe.utils import calculate_impact, check_folder -from magicctapipe.utils.filedir import load_cfg_file, out_file_h5 -from magicctapipe.utils.utils import print_title, print_elapsed_time -from magicctapipe.utils.tels import check_tel_ids from magicctapipe.reco.stereo import ( + StereoInfoContainer, check_stereo, write_hillas, write_stereo, - StereoInfoContainer, ) - +from magicctapipe.utils import calculate_impact, check_folder +from magicctapipe.utils.filedir import load_cfg_file, out_file_h5 +from magicctapipe.utils.tels import check_tel_ids +from magicctapipe.utils.utils import print_elapsed_time, print_title PARSER = argparse.ArgumentParser( description="Stereo Reconstruction MAGIC + LST", diff --git a/magicctapipe/scripts/performance/train_rfs_MAGIC_LST.py b/magicctapipe/scripts/performance/train_rfs_MAGIC_LST.py index a088b27b5..270d14ec1 100644 --- a/magicctapipe/scripts/performance/train_rfs_MAGIC_LST.py +++ b/magicctapipe/scripts/performance/train_rfs_MAGIC_LST.py @@ -1,49 +1,51 @@ # coding: utf-8 -import time -import os -import glob import argparse -import pandas as pd -import numpy as np +import glob +import os +import time + import matplotlib.pyplot as plt +import numpy as np +import pandas as pd from astropy.coordinates.angle_utilities import angular_separation -from magicctapipe.utils.plot import ( - save_plt, - load_default_plot_settings, - load_default_plot_settings_02, -) -from magicctapipe.utils.utils import print_title, info_message, print_elapsed_time -from magicctapipe.reco.energy_utils import plot_migmatrix, evaluate_performance_energy -from magicctapipe.reco.direction_utils import compute_separation_angle_direction from magicctapipe.reco.classifier_utils import ( - get_weights_classifier, check_train_test_intersections_classifier, + evaluate_performance_classifier, + get_weights_classifier, load_init_data_classifier, print_par_imp_classifier, - evaluate_performance_classifier, +) +from magicctapipe.reco.direction_utils import compute_separation_angle_direction +from magicctapipe.reco.energy_utils import evaluate_performance_energy, plot_migmatrix +from magicctapipe.reco.event_processing import ( + DirectionEstimatorPandas, + EnergyEstimatorPandas, + EventClassifierPandas, ) from magicctapipe.reco.global_utils import ( - get_weights_mc_dir_class, - compute_event_weights, check_train_test_intersections, -) -from magicctapipe.utils.tels import ( - check_tel_ids, - get_tel_name, - get_array_tel_descriptions, + compute_event_weights, + get_weights_mc_dir_class, ) from magicctapipe.utils.filedir import ( - load_cfg_file_check, check_folder, - load_dl1_data_stereo_list_selected, + load_cfg_file_check, load_dl1_data_stereo_list, + load_dl1_data_stereo_list_selected, ) - -from magicctapipe.reco.event_processing import EventClassifierPandas -from magicctapipe.reco.event_processing import EnergyEstimatorPandas -from magicctapipe.reco.event_processing import DirectionEstimatorPandas +from magicctapipe.utils.plot import ( + load_default_plot_settings, + load_default_plot_settings_02, + save_plt, +) +from magicctapipe.utils.tels import ( + check_tel_ids, + get_array_tel_descriptions, + get_tel_name, +) +from magicctapipe.utils.utils import info_message, print_elapsed_time, print_title PARSER = argparse.ArgumentParser( description="Trains random forests for stereo data", diff --git a/magicctapipe/utils/__init__.py b/magicctapipe/utils/__init__.py index 955927f23..4f56846b0 100644 --- a/magicctapipe/utils/__init__.py +++ b/magicctapipe/utils/__init__.py @@ -1,36 +1,22 @@ -from .badpixels import ( - MAGICBadPixelsCalc, -) - -from .camera_geometry import ( - scale_camera_geometry, - reflected_camera_geometry, -) - +from .badpixels import MAGICBadPixelsCalc +from .camera_geometry import reflected_camera_geometry, scale_camera_geometry from .filedir import ( + check_common_keys, + check_folder, + convert_np_list_dict, + drop_keys, load_cfg_file, load_cfg_file_check, - check_folder, - load_dl1_data_stereo_list_selected, - load_dl1_data_stereo_list, - load_dl1_data_stereo, load_dl1_data_mono, - drop_keys, - check_common_keys, - out_file_h5_no_run, + load_dl1_data_stereo, + load_dl1_data_stereo_list, + load_dl1_data_stereo_list_selected, out_file_h5, + out_file_h5_no_run, out_file_h5_reco, read_mc_header, save_yaml_np, - convert_np_list_dict, -) - -from .gti import ( - identify_time_edges, - intersect_time_intervals, - GTIGenerator, ) - from .functions import ( calculate_disp, calculate_impact, @@ -38,31 +24,25 @@ calculate_off_coordinates, transform_altaz_to_radec, ) - -from .plot import ( - save_plt, - load_default_plot_settings, - load_default_plot_settings_02, -) - +from .gti import GTIGenerator, identify_time_edges, intersect_time_intervals +from .plot import load_default_plot_settings, load_default_plot_settings_02, save_plt from .tels import ( - tel_ids_2_num, - num_2_tel_ids, - get_tel_descriptions, + check_tel_ids, + convert_positions_dict, get_array_tel_descriptions, + get_tel_descriptions, get_tel_ids_dl1, - convert_positions_dict, - check_tel_ids, - intersec_tel_ids, get_tel_name, + intersec_tel_ids, + num_2_tel_ids, + tel_ids_2_num, ) - from .utils import ( info_message, - print_elapsed_time, make_elapsed_time_str, - print_title, make_title_str, + print_elapsed_time, + print_title, resource_file, ) diff --git a/magicctapipe/utils/badpixels.py b/magicctapipe/utils/badpixels.py index 3654df6cd..e3a1dd1e3 100644 --- a/magicctapipe/utils/badpixels.py +++ b/magicctapipe/utils/badpixels.py @@ -1,5 +1,5 @@ -from astropy.time import Time import numpy as np +from astropy.time import Time from ctapipe.instrument import CameraGeometry __all__ = [ diff --git a/magicctapipe/utils/filedir.py b/magicctapipe/utils/filedir.py index 543a6d5af..c1ba12819 100644 --- a/magicctapipe/utils/filedir.py +++ b/magicctapipe/utils/filedir.py @@ -1,9 +1,10 @@ +import copy import os import sys -import yaml -import copy -import pandas as pd + import numpy as np +import pandas as pd +import yaml __all__ = [ "load_cfg_file", @@ -216,7 +217,7 @@ def load_dl1_data_stereo(file, drop=False, slope_abs=False): data = data_hillas.merge(data_stereo, on=common_keys) # Index data.set_index(["obs_id", "event_id", "tel_id"], inplace=True) - except: + except Exception: data = pd.read_hdf(file, key="dl2/reco") data.sort_index(inplace=True) # Get absolute value of slope, if needed diff --git a/magicctapipe/utils/gti.py b/magicctapipe/utils/gti.py index 89e3641b7..d1b967a21 100644 --- a/magicctapipe/utils/gti.py +++ b/magicctapipe/utils/gti.py @@ -1,8 +1,8 @@ # coding: utf-8 +import pandas import scipy import uproot -import pandas from .utils import info_message diff --git a/magicctapipe/utils/plot.py b/magicctapipe/utils/plot.py index 64fca9ae9..762acae8c 100644 --- a/magicctapipe/utils/plot.py +++ b/magicctapipe/utils/plot.py @@ -1,4 +1,5 @@ import os + import matplotlib.pylab as plt __all__ = [ diff --git a/magicctapipe/utils/tels.py b/magicctapipe/utils/tels.py index b48ec692e..e934187f7 100644 --- a/magicctapipe/utils/tels.py +++ b/magicctapipe/utils/tels.py @@ -1,7 +1,6 @@ -from astropy import units as u import numpy as np - -from ctapipe.instrument import CameraGeometry, TelescopeDescription, OpticsDescription +from astropy import units as u +from ctapipe.instrument import CameraGeometry, OpticsDescription, TelescopeDescription __all__ = [ "tel_ids_2_num", diff --git a/magicctapipe/utils/utils.py b/magicctapipe/utils/utils.py index 7c90829d0..a291e0e27 100644 --- a/magicctapipe/utils/utils.py +++ b/magicctapipe/utils/utils.py @@ -1,4 +1,5 @@ import datetime + import numpy as np try: @@ -140,6 +141,7 @@ def make_title_str(title, style_char="=", in_space=3, width_char=80): s = f"{s1}{s2}{s3}" return s + def resource_file(filename): """Get the absoulte path of magicctapipe resource files.""" return files("magicctapipe").joinpath("resources", filename) diff --git a/notebooks/check_event_coincidence.ipynb b/notebooks/check_event_coincidence.ipynb index c9c23f245..d626053f9 100644 --- a/notebooks/check_event_coincidence.ipynb +++ b/notebooks/check_event_coincidence.ipynb @@ -64,9 +64,7 @@ "# === Settings ===\n", "# ================\n", "\n", - "input_file_mask = (\n", - " \"/home/evisentin/dataout_magic_lst/dl1/data/coincidence/lst_single_run_stereo/run_03265/\"\n", - ")\n", + "input_file_mask = \"/home/evisentin/dataout_magic_lst/dl1/data/coincidence/lst_single_run_stereo/run_03265/\"\n", "\n", "# ============\n", "# === Main ===\n", diff --git a/setup.py b/setup.py index f8c4a145e..ba50382bf 100644 --- a/setup.py +++ b/setup.py @@ -2,7 +2,8 @@ # Licensed under a 3-clause BSD style license - see LICENSE.rst import os -from setuptools import setup, find_packages + +from setuptools import find_packages, setup entry_points = {} entry_points["console_scripts"] = [ @@ -33,19 +34,19 @@ use_scm_version={"write_to": os.path.join("magicctapipe", "_version.py")}, packages=find_packages(), install_requires=[ - 'lstchain~=0.9.6', - 'ctapipe~=0.12.0', - 'ctapipe_io_magic~=0.4.7', - 'ctaplot~=0.5.5', - 'gammapy~=0.19.0', - 'uproot~=4.1', - 'numpy<1.22.0a0', - 'joblib', - 'pandas', - 'pyirf~=0.6.0', - 'seaborn', - 'scikit-learn', - 'setuptools_scm', + "lstchain~=0.9.6", + "ctapipe~=0.12.0", + "ctapipe_io_magic~=0.4.7", + "ctaplot~=0.5.5", + "gammapy~=0.19.0", + "uproot~=4.1", + "numpy<1.22.0a0", + "joblib", + "pandas", + "pyirf~=0.6.0", + "seaborn", + "scikit-learn", + "setuptools_scm", ], entry_points=entry_points, extras_require={ From 8e9d9cf79dd34fdda7283372a6183259bdbde42e Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Mon, 25 Sep 2023 22:26:04 +0200 Subject: [PATCH 03/10] Add pre-commit in CI. --- .github/workflows/ci.yml | 13 +++++++++++++ 1 file changed, 13 insertions(+) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 90ea27050..a9a3d9894 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -8,6 +8,19 @@ env: PYTEST_ADDOPTS: --color=yes jobs: + lint: + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v4 + with: + fetch-depth: 0 + - uses: actions/setup-python@v4 + with: + python-version: 3.8 + - uses: pre-commit/action@v3.0.0 + with: + extra_args: --files $(git diff origin/main --name-only) + pyflakes: runs-on: ubuntu-latest steps: From 973b25e79f9f8bfe57dbac9ab555061a05062535 Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Tue, 26 Sep 2023 11:14:44 +0200 Subject: [PATCH 04/10] Import only version. --- magicctapipe/__init__.py | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/magicctapipe/__init__.py b/magicctapipe/__init__.py index a51ad71a3..2d735710f 100644 --- a/magicctapipe/__init__.py +++ b/magicctapipe/__init__.py @@ -1,11 +1,3 @@ -from . import image, irfs, reco, scripts, utils from .version import __version__ -__all__ = [ - "image", - "irfs", - "reco", - "scripts", - "utils", - "__version__", -] +__all__ = ["__version__"] From e6f8d799aa3ffd86a27b3e4a03abc82dc888f843 Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Tue, 26 Sep 2023 11:14:59 +0200 Subject: [PATCH 05/10] Fix imports. --- magicctapipe/io/__init__.py | 4 ++++ magicctapipe/io/gadf.py | 6 +++--- magicctapipe/io/io.py | 4 +++- magicctapipe/irfs/utils.py | 4 +--- magicctapipe/reco/classifier_utils.py | 6 +++--- magicctapipe/reco/direction_utils.py | 2 +- magicctapipe/reco/estimators.py | 2 +- magicctapipe/reco/stereo.py | 2 +- magicctapipe/utils/__init__.py | 6 ++++++ magicctapipe/utils/functions.py | 3 +++ 10 files changed, 26 insertions(+), 13 deletions(-) diff --git a/magicctapipe/io/__init__.py b/magicctapipe/io/__init__.py index fa5fb40f7..2dc34a30a 100644 --- a/magicctapipe/io/__init__.py +++ b/magicctapipe/io/__init__.py @@ -13,6 +13,8 @@ create_pointing_hdu, ) from .io import ( + TEL_COMBINATIONS, + TEL_NAMES, format_object, get_dl2_mean, get_stereo_events, @@ -43,4 +45,6 @@ "load_mc_dl2_data_file", "load_train_data_files", "save_pandas_data_in_table", + "TEL_COMBINATIONS", + "TEL_NAMES", ] diff --git a/magicctapipe/io/gadf.py b/magicctapipe/io/gadf.py index e27c94da7..fdfe4b13e 100644 --- a/magicctapipe/io/gadf.py +++ b/magicctapipe/io/gadf.py @@ -11,9 +11,9 @@ from astropy.time import Time from pyirf.binning import split_bin_lo_hi -from magicctapipe import __version__ -from magicctapipe.io.io import TEL_COMBINATIONS -from magicctapipe.utils.functions import HEIGHT_ORM, LAT_ORM, LON_ORM +from .. import __version__ +from ..utils import HEIGHT_ORM, LAT_ORM, LON_ORM +from .io import TEL_COMBINATIONS __all__ = [ "create_gh_cuts_hdu", diff --git a/magicctapipe/io/io.py b/magicctapipe/io/io.py index b8c6dd11f..70b21c777 100644 --- a/magicctapipe/io/io.py +++ b/magicctapipe/io/io.py @@ -21,7 +21,7 @@ from pyirf.simulations import SimulatedEventsInfo from pyirf.utils import calculate_source_fov_offset, calculate_theta -from magicctapipe.utils import calculate_mean_direction, transform_altaz_to_radec +from ..utils import calculate_mean_direction, transform_altaz_to_radec __all__ = [ "format_object", @@ -34,6 +34,8 @@ "load_dl2_data_file", "load_irf_files", "save_pandas_data_in_table", + "TEL_COMBINATIONS", + "TEL_NAMES", ] logger = logging.getLogger(__name__) diff --git a/magicctapipe/irfs/utils.py b/magicctapipe/irfs/utils.py index f44c54796..e8cbbdcbd 100644 --- a/magicctapipe/irfs/utils.py +++ b/magicctapipe/irfs/utils.py @@ -10,9 +10,7 @@ from astropy.table import QTable from pyirf.simulations import SimulatedEventsInfo -from magicctapipe.utils.filedir import read_mc_header -from magicctapipe.utils.plot import load_default_plot_settings, save_plt -from magicctapipe.utils.utils import print_title +from ..utils import load_default_plot_settings, print_title, read_mc_header, save_plt __all__ = [ "read_simu_info_mcp_sum_num_showers", diff --git a/magicctapipe/reco/classifier_utils.py b/magicctapipe/reco/classifier_utils.py index 9d143bf11..341f3afc2 100644 --- a/magicctapipe/reco/classifier_utils.py +++ b/magicctapipe/reco/classifier_utils.py @@ -5,12 +5,12 @@ import pandas as pd import sklearn.metrics -from magicctapipe.reco.global_utils import check_train_test_intersections -from magicctapipe.utils.filedir import ( +from ..utils import ( + info_message, load_dl1_data_stereo_list, load_dl1_data_stereo_list_selected, ) -from magicctapipe.utils.utils import info_message +from .global_utils import check_train_test_intersections __all__ = [ "GetHist_classifier", diff --git a/magicctapipe/reco/direction_utils.py b/magicctapipe/reco/direction_utils.py index 9aff275a7..6d65dcb5a 100644 --- a/magicctapipe/reco/direction_utils.py +++ b/magicctapipe/reco/direction_utils.py @@ -2,7 +2,7 @@ from astropy import units as u from astropy.coordinates import AltAz, SkyCoord -from magicctapipe.utils.tels import get_tel_ids_dl1 +from ..utils import get_tel_ids_dl1 __all__ = ["compute_separation_angle_direction"] diff --git a/magicctapipe/reco/estimators.py b/magicctapipe/reco/estimators.py index 91d993b62..53ddf3100 100644 --- a/magicctapipe/reco/estimators.py +++ b/magicctapipe/reco/estimators.py @@ -8,7 +8,7 @@ import pandas as pd import sklearn.ensemble -from magicctapipe.io.io import TEL_NAMES +from ..io import TEL_NAMES __all__ = ["EnergyRegressor", "DispRegressor", "EventClassifier"] diff --git a/magicctapipe/reco/stereo.py b/magicctapipe/reco/stereo.py index 422cddee7..5fe70a63f 100644 --- a/magicctapipe/reco/stereo.py +++ b/magicctapipe/reco/stereo.py @@ -2,7 +2,7 @@ import numpy as np from ctapipe.core.container import Container, Field -from magicctapipe.utils.tels import tel_ids_2_num +from ..utils import tel_ids_2_num __all__ = ["write_hillas", "check_write_stereo", "check_stereo", "write_stereo"] diff --git a/magicctapipe/utils/__init__.py b/magicctapipe/utils/__init__.py index 4f56846b0..d9c3f84ac 100644 --- a/magicctapipe/utils/__init__.py +++ b/magicctapipe/utils/__init__.py @@ -18,6 +18,9 @@ save_yaml_np, ) from .functions import ( + HEIGHT_ORM, + LAT_ORM, + LON_ORM, calculate_disp, calculate_impact, calculate_mean_direction, @@ -92,4 +95,7 @@ "make_title_str", "resource_file", "get_key_if_exists", + "LON_ORM", + "LAT_ORM", + "HEIGHT_ORM", ] diff --git a/magicctapipe/utils/functions.py b/magicctapipe/utils/functions.py index aebc363a9..49fe50995 100644 --- a/magicctapipe/utils/functions.py +++ b/magicctapipe/utils/functions.py @@ -20,6 +20,9 @@ "calculate_mean_direction", "calculate_off_coordinates", "transform_altaz_to_radec", + "LON_ORM", + "LAT_ORM", + "HEIGHT_ORM", ] # The geographic coordinate of ORM From 5b8935438291b06d0c75841359ce5b632f7d1254 Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Sun, 29 Oct 2023 13:04:46 +0100 Subject: [PATCH 06/10] Run pre-commit after merge with master. --- docs/conf.py | 26 +- magicctapipe/conftest.py | 12 +- magicctapipe/image/__init__.py | 3 +- magicctapipe/image/calib.py | 81 +++--- magicctapipe/image/tests/test_calib.py | 275 ++++++++---------- magicctapipe/io/__init__.py | 13 +- magicctapipe/io/gadf.py | 17 +- magicctapipe/io/io.py | 147 ++++++---- magicctapipe/io/tests/test_io.py | 215 ++++++++------ magicctapipe/io/tests/test_io_monly.py | 19 +- magicctapipe/reco/estimators.py | 9 +- .../lst1_magic_dl1_stereo_to_dl2.py | 4 +- .../lst1_magic_event_coincidence.py | 74 ++--- .../lst1_magic/lst1_magic_mc_dl0_to_dl1.py | 132 +++++---- .../lst1_magic/lst1_magic_stereo_reco.py | 58 ++-- .../lst1_magic/lst1_magic_train_rfs.py | 6 +- .../scripts/lst1_magic/magic_calib_to_dl1.py | 11 +- .../scripts/lst1_magic/merge_hdf_files.py | 2 +- setup.py | 28 +- 19 files changed, 631 insertions(+), 501 deletions(-) diff --git a/docs/conf.py b/docs/conf.py index 44f4ac00d..dd72b0f51 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -3,11 +3,11 @@ import datetime import os -import magicctapipe - # Get configuration information from setup.cfg from configparser import ConfigParser +import magicctapipe + setup_cfg = ConfigParser() setup_cfg.read([os.path.join(os.path.dirname(__file__), "..", "setup.cfg")]) setup_metadata = dict(setup_cfg.items("metadata")) @@ -104,28 +104,28 @@ # -- General configuration extensions = [ - 'sphinx.ext.duration', - 'sphinx.ext.doctest', - 'sphinx.ext.autodoc', - 'sphinx.ext.autosummary', - 'sphinx.ext.intersphinx', + "sphinx.ext.duration", + "sphinx.ext.doctest", + "sphinx.ext.autodoc", + "sphinx.ext.autosummary", + "sphinx.ext.intersphinx", ] intersphinx_mapping = { - 'python': ('https://docs.python.org/3/', None), - 'sphinx': ('https://www.sphinx-doc.org/en/master/', None), + "python": ("https://docs.python.org/3/", None), + "sphinx": ("https://www.sphinx-doc.org/en/master/", None), } -intersphinx_disabled_domains = ['std'] +intersphinx_disabled_domains = ["std"] -templates_path = ['_templates'] +templates_path = ["_templates"] # -- Options for HTML output -html_theme = 'pydata_sphinx_theme' +html_theme = "pydata_sphinx_theme" # The name for this set of Sphinx documents. If None, it defaults to # " v documentation". html_title = f"{project} v{release}" # -- Options for EPUB output -epub_show_urls = 'footnote' +epub_show_urls = "footnote" diff --git a/magicctapipe/conftest.py b/magicctapipe/conftest.py index f22998c3f..0fcce33e8 100644 --- a/magicctapipe/conftest.py +++ b/magicctapipe/conftest.py @@ -4,12 +4,12 @@ import numpy as np import pandas as pd import pytest +import yaml from astropy.io.misc.hdf5 import write_table_hdf5 from astropy.table import Table from ctapipe.utils.download import download_file_cached from magicctapipe.utils import resource_file -import yaml maxjoint = 13000 maxmonly = 500 @@ -50,6 +50,7 @@ Temporary paths """ + @pytest.fixture(scope="session") def temp_DL1_gamma(tmp_path_factory): return tmp_path_factory.mktemp("DL1_gammas") @@ -234,10 +235,12 @@ def temp_DL2_real_monly(tmp_path_factory): def temp_DL3_monly(tmp_path_factory): return tmp_path_factory.mktemp("DL3_monly") + """ Custom data """ + @pytest.fixture(scope="session") def dl2_test(temp_DL2_test): """ @@ -271,10 +274,12 @@ def pd_test(): df = pd.DataFrame(np.array([[1, 2], [3, 4], [5, 6]]), columns=["a", "b"]) return df + """ Remote paths (to download test files) """ + @pytest.fixture(scope="session") def base_url(): return "http://www.magic.iac.es/mcp-testdata" @@ -285,10 +290,12 @@ def env_prefix(): # ENVIRONMENT VARIABLES TO BE CREATED return "MAGIC_CTA_DATA_" + """ Downloads: files """ + @pytest.fixture(scope="session") def dl0_gamma(base_url, env_prefix): gamma_dl0 = [] @@ -379,7 +386,7 @@ def config(): @pytest.fixture(scope="session") def config_monly(): config_path = resource_file("test_config_monly.yaml") - return config_path + return config_path @pytest.fixture(scope="session") @@ -418,6 +425,7 @@ def config_check(): Data processing """ + @pytest.fixture(scope="session") def gamma_l1(temp_DL1_gamma, dl0_gamma, config): """ diff --git a/magicctapipe/image/__init__.py b/magicctapipe/image/__init__.py index 07a7060f8..c9856823a 100644 --- a/magicctapipe/image/__init__.py +++ b/magicctapipe/image/__init__.py @@ -1,3 +1,4 @@ +from .calib import calibrate from .cleaning import ( MAGICClean, PixelTreatment, @@ -5,8 +6,6 @@ get_num_islands_MAGIC, ) from .leakage import get_leakage -from .calib import calibrate - __all__ = [ "MAGICClean", diff --git a/magicctapipe/image/calib.py b/magicctapipe/image/calib.py index 95ae4410d..c11b92480 100644 --- a/magicctapipe/image/calib.py +++ b/magicctapipe/image/calib.py @@ -1,37 +1,40 @@ import numpy as np -from ctapipe.image import ( - apply_time_delta_cleaning, - number_of_islands, - tailcuts_clean, -) +from ctapipe.image import apply_time_delta_cleaning, number_of_islands, tailcuts_clean from ctapipe.instrument import CameraGeometry from lstchain.image.cleaning import apply_dynamic_cleaning from lstchain.image.modifier import ( add_noise_in_pixels, - random_psf_smearer, - set_numba_seed + random_psf_smearer, + set_numba_seed, ) -from magicctapipe.image import MAGICClean +from magicctapipe.image import MAGICClean -__all__ = [ - "calibrate" -] +__all__ = ["calibrate"] -def calibrate(event, tel_id, config, calibrator, is_lst, obs_id=None, camera_geoms=None, magic_clean=None): +def calibrate( + event, + tel_id, + config, + calibrator, + is_lst, + obs_id=None, + camera_geoms=None, + magic_clean=None, +): """ - This function calibrates the camera image for a single event of a telescope + This function calibrates the camera image for a single event of a telescope Parameters ---------- - event: event + event: event From an EventSource tel_id: int - Telescope ID + Telescope ID config: dictionary - Parameters for image extraction and calibration + Parameters for image extraction and calibration calibrator: CameraCalibrator (ctapipe.calib) ctapipe object needed to calibrate the camera is_lst: bool @@ -42,7 +45,7 @@ def calibrate(event, tel_id, config, calibrator, is_lst, obs_id=None, camera_geo Camera geometry. Used in case of LST telescope magic_clean: dictionary (1 entry per MAGIC telescope) Each entry is a MAGICClean object using the telescope camera geometry. Used in case of MAGIC telescope - + Returns ------- @@ -54,23 +57,33 @@ def calibrate(event, tel_id, config, calibrator, is_lst, obs_id=None, camera_geo Array of the signal peak time in the camera pixels """ - if (not is_lst) and (magic_clean==None): - raise ValueError("Check the provided parameters and the telescope type; MAGIC calibration not possible if magic_clean not provided") - if (is_lst) and (obs_id==None): - raise ValueError("Check the provided parameters and the telescope type; LST calibration not possible if obs_id not provided") - if (is_lst) and (camera_geoms==None): - raise ValueError("Check the provided parameters and the telescope type; LST calibration not possible if gamera_geoms not provided") - if (not is_lst) and (type(magic_clean[tel_id])!=MAGICClean): - raise ValueError("Check the provided magic_clean parameter; MAGIC calibration not possible if magic_clean not a dictionary of MAGICClean objects") - if (is_lst) and (type(camera_geoms[tel_id])!=CameraGeometry): - raise ValueError("Check the provided camera_geoms parameter; LST calibration not possible if camera_geoms not a dictionary of CameraGeometry objects") + if (not is_lst) and (magic_clean is None): + raise ValueError( + "Check the provided parameters and the telescope type; MAGIC calibration not possible if magic_clean not provided" + ) + if (is_lst) and (obs_id is None): + raise ValueError( + "Check the provided parameters and the telescope type; LST calibration not possible if obs_id not provided" + ) + if (is_lst) and (camera_geoms is None): + raise ValueError( + "Check the provided parameters and the telescope type; LST calibration not possible if gamera_geoms not provided" + ) + if (not is_lst) and (type(magic_clean[tel_id]) != MAGICClean): + raise ValueError( + "Check the provided magic_clean parameter; MAGIC calibration not possible if magic_clean not a dictionary of MAGICClean objects" + ) + if (is_lst) and (type(camera_geoms[tel_id]) != CameraGeometry): + raise ValueError( + "Check the provided camera_geoms parameter; LST calibration not possible if camera_geoms not a dictionary of CameraGeometry objects" + ) calibrator._calibrate_dl0(event, tel_id) calibrator._calibrate_dl1(event, tel_id) image = event.dl1.tel[tel_id].image.astype(np.float64) peak_time = event.dl1.tel[tel_id].peak_time.astype(np.float64) - + if not is_lst: use_charge_correction = config["charge_correction"]["use"] @@ -82,14 +95,14 @@ def calibrate(event, tel_id, config, calibrator, is_lst, obs_id=None, camera_geo signal_pixels, image, peak_time = magic_clean[tel_id].clean_image( event_image=image, event_pulse_time=peak_time ) - + else: increase_nsb = config["increase_nsb"].pop("use") increase_psf = config["increase_psf"]["use"] use_time_delta_cleaning = config["time_delta_cleaning"].pop("use") use_dynamic_cleaning = config["dynamic_cleaning"].pop("use") use_only_main_island = config["use_only_main_island"] - + if increase_nsb: rng = np.random.default_rng(obs_id) # Add extra noise in pixels @@ -133,8 +146,8 @@ def calibrate(event, tel_id, config, calibrator, is_lst, obs_id=None, camera_geo max_island_label = np.argmax(n_pixels_on_island) signal_pixels[island_labels != max_island_label] = False - config["increase_nsb"]["use"]=increase_nsb - config["time_delta_cleaning"]["use"]=use_time_delta_cleaning - config["dynamic_cleaning"]["use"]=use_dynamic_cleaning - + config["increase_nsb"]["use"] = increase_nsb + config["time_delta_cleaning"]["use"] = use_time_delta_cleaning + config["dynamic_cleaning"]["use"] = use_dynamic_cleaning + return signal_pixels, image, peak_time diff --git a/magicctapipe/image/tests/test_calib.py b/magicctapipe/image/tests/test_calib.py index 7229686e8..6426e0df5 100644 --- a/magicctapipe/image/tests/test_calib.py +++ b/magicctapipe/image/tests/test_calib.py @@ -1,14 +1,10 @@ -from magicctapipe.image.calib import calibrate import pytest - - - from ctapipe.calib import CameraCalibrator - from ctapipe.io import EventSource +from traitlets.config import Config from magicctapipe.image import MAGICClean -from traitlets.config import Config +from magicctapipe.image.calib import calibrate @pytest.fixture(scope="session") @@ -22,26 +18,23 @@ def tel_id_MAGIC(): def test_calibrate_LST(dl0_gamma, config_calib, tel_id_LST): - - assigned_tel_ids = [1,2,3] + + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) - obs_id = event_source.obs_ids[0] - subarray = event_source.subarray + subarray = event_source.subarray - tel_descriptions = subarray.tel - camera_geoms={} + tel_descriptions = subarray.tel + camera_geoms = {} for tel_id, telescope in tel_descriptions.items(): - camera_geoms[tel_id]= telescope.camera.geometry - + camera_geoms[tel_id] = telescope.camera.geometry + config_lst = config_calib["LST"] extractor_type_lst = config_lst["image_extractor"].pop("type") @@ -52,92 +45,86 @@ def test_calibrate_LST(dl0_gamma, config_calib, tel_id_LST): config=Config(config_extractor_lst), subarray=subarray, ) - + for event in event_source: - if (event.count <200) and (tel_id_LST in event.trigger.tels_with_trigger): + if (event.count < 200) and (tel_id_LST in event.trigger.tels_with_trigger): signal_pixels, image, peak_time = calibrate( - event=event, - tel_id=tel_id_LST, - obs_id=obs_id, - config=config_lst, - camera_geoms=camera_geoms, - calibrator=calibrator_lst, + event=event, + tel_id=tel_id_LST, + obs_id=obs_id, + config=config_lst, + camera_geoms=camera_geoms, + calibrator=calibrator_lst, is_lst=True, ) - assert len(signal_pixels)==1855 - assert signal_pixels.dtype==bool - assert len(image)==1855 - assert len(peak_time)==1855 + assert len(signal_pixels) == 1855 + assert signal_pixels.dtype == bool + assert len(image) == 1855 + assert len(peak_time) == 1855 - config_lst["image_extractor"]["type"]=extractor_type_lst + config_lst["image_extractor"]["type"] = extractor_type_lst def test_calibrate_MAGIC(dl0_gamma, config_calib, tel_id_MAGIC): - - assigned_tel_ids = [1,2,3] + + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) + subarray = event_source.subarray - - subarray = event_source.subarray - - - tel_descriptions = subarray.tel - camera_geoms={} + tel_descriptions = subarray.tel + camera_geoms = {} for tel_id, telescope in tel_descriptions.items(): - camera_geoms[tel_id]= telescope.camera.geometry - - + camera_geoms[tel_id] = telescope.camera.geometry + config_magic = config_calib["MAGIC"] config_magic["magic_clean"].update({"find_hotpixels": False}) extractor_type_magic = config_magic["image_extractor"].pop("type") config_extractor_magic = {extractor_type_magic: config_magic["image_extractor"]} magic_clean = {} - for k in [1,2]: - + for k in [1, 2]: + magic_clean[k] = MAGICClean(camera_geoms[k], config_magic["magic_clean"]) calibrator_magic = CameraCalibrator( image_extractor_type=extractor_type_magic, config=Config(config_extractor_magic), subarray=subarray, ) - + for event in event_source: - if (event.count <200) and (tel_id_MAGIC in event.trigger.tels_with_trigger): + if (event.count < 200) and ( + tel_id_MAGIC in event.trigger.tels_with_trigger + ): signal_pixels, image, peak_time = calibrate( - event=event, - tel_id=tel_id_MAGIC, - config=config_magic, - magic_clean=magic_clean, - calibrator=calibrator_magic, + event=event, + tel_id=tel_id_MAGIC, + config=config_magic, + magic_clean=magic_clean, + calibrator=calibrator_magic, is_lst=False, ) - assert len(signal_pixels)==1039 - assert signal_pixels.dtype==bool - assert len(image)==1039 - assert len(peak_time)==1039 + assert len(signal_pixels) == 1039 + assert signal_pixels.dtype == bool + assert len(image) == 1039 + assert len(peak_time) == 1039 - config_magic["image_extractor"]["type"]=extractor_type_magic + config_magic["image_extractor"]["type"] = extractor_type_magic def test_calibrate_exc_1(dl0_gamma, config_calib, tel_id_MAGIC): - assigned_tel_ids = [1,2,3] + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) - subarray = event_source.subarray + subarray = event_source.subarray config_magic = config_calib["MAGIC"] config_magic["magic_clean"].update({"find_hotpixels": False}) extractor_type_magic = config_magic["image_extractor"].pop("type") @@ -147,43 +134,40 @@ def test_calibrate_exc_1(dl0_gamma, config_calib, tel_id_MAGIC): config=Config(config_extractor_magic), subarray=subarray, ) - + for event in event_source: - if (event.count <200) and (tel_id_MAGIC in event.trigger.tels_with_trigger): + if (event.count < 200) and ( + tel_id_MAGIC in event.trigger.tels_with_trigger + ): with pytest.raises( ValueError, match="Check the provided parameters and the telescope type; MAGIC calibration not possible if magic_clean not provided", ): - _,_,_ = calibrate( - event=event, - tel_id=tel_id_MAGIC, - config=config_magic, - calibrator=calibrator_magic, + _, _, _ = calibrate( + event=event, + tel_id=tel_id_MAGIC, + config=config_magic, + calibrator=calibrator_magic, is_lst=False, ) - config_magic["image_extractor"]["type"]=extractor_type_magic + config_magic["image_extractor"]["type"] = extractor_type_magic def test_calibrate_exc_2(dl0_gamma, config_calib, tel_id_LST): - assigned_tel_ids = [1,2,3] + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) + subarray = event_source.subarray - - - subarray = event_source.subarray - - tel_descriptions = subarray.tel - camera_geoms={} + tel_descriptions = subarray.tel + camera_geoms = {} for tel_id, telescope in tel_descriptions.items(): - camera_geoms[tel_id]= telescope.camera.geometry - + camera_geoms[tel_id] = telescope.camera.geometry + config_lst = config_calib["LST"] extractor_type_lst = config_lst["image_extractor"].pop("type") @@ -194,40 +178,35 @@ def test_calibrate_exc_2(dl0_gamma, config_calib, tel_id_LST): config=Config(config_extractor_lst), subarray=subarray, ) - + for event in event_source: - if (event.count <200) and (tel_id_LST in event.trigger.tels_with_trigger): + if (event.count < 200) and (tel_id_LST in event.trigger.tels_with_trigger): with pytest.raises( ValueError, match="Check the provided parameters and the telescope type; LST calibration not possible if obs_id not provided", - ): - _,_,_ = calibrate( - event=event, - tel_id=tel_id_LST, - config=config_lst, - camera_geoms=camera_geoms, - calibrator=calibrator_lst, + ): + _, _, _ = calibrate( + event=event, + tel_id=tel_id_LST, + config=config_lst, + camera_geoms=camera_geoms, + calibrator=calibrator_lst, is_lst=True, ) - config_lst["image_extractor"]["type"]=extractor_type_lst + config_lst["image_extractor"]["type"] = extractor_type_lst def test_calibrate_exc_3(dl0_gamma, config_calib, tel_id_LST): - assigned_tel_ids = [1,2,3] + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) - obs_id = event_source.obs_ids[0] - subarray = event_source.subarray + subarray = event_source.subarray - - config_lst = config_calib["LST"] extractor_type_lst = config_lst["image_extractor"].pop("type") @@ -238,38 +217,36 @@ def test_calibrate_exc_3(dl0_gamma, config_calib, tel_id_LST): config=Config(config_extractor_lst), subarray=subarray, ) - + for event in event_source: - if (event.count <200) and (tel_id_LST in event.trigger.tels_with_trigger): + if (event.count < 200) and (tel_id_LST in event.trigger.tels_with_trigger): with pytest.raises( ValueError, match="Check the provided parameters and the telescope type; LST calibration not possible if gamera_geoms not provided", - ): + ): signal_pixels, image, peak_time = calibrate( - event=event, - tel_id=tel_id_LST, - obs_id=obs_id, - config=config_lst, - calibrator=calibrator_lst, + event=event, + tel_id=tel_id_LST, + obs_id=obs_id, + config=config_lst, + calibrator=calibrator_lst, is_lst=True, ) - config_lst["image_extractor"]["type"]=extractor_type_lst + config_lst["image_extractor"]["type"] = extractor_type_lst def test_calibrate_exc_4(dl0_gamma, config_calib, tel_id_MAGIC): - assigned_tel_ids = [1,2,3] + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) - subarray = event_source.subarray - tel_descriptions = subarray.tel - magic_clean={} + subarray = event_source.subarray + tel_descriptions = subarray.tel + magic_clean = {} for tel_id in range(len(tel_descriptions.items())): - magic_clean[tel_id]= f"camera {tel_id}" + magic_clean[tel_id] = f"camera {tel_id}" config_magic = config_calib["MAGIC"] config_magic["magic_clean"].update({"find_hotpixels": False}) extractor_type_magic = config_magic["image_extractor"].pop("type") @@ -278,45 +255,44 @@ def test_calibrate_exc_4(dl0_gamma, config_calib, tel_id_MAGIC): image_extractor_type=extractor_type_magic, config=Config(config_extractor_magic), subarray=subarray, - ) - + ) + for event in event_source: - if (event.count <200) and (tel_id_MAGIC in event.trigger.tels_with_trigger): + if (event.count < 200) and ( + tel_id_MAGIC in event.trigger.tels_with_trigger + ): with pytest.raises( ValueError, match="Check the provided magic_clean parameter; MAGIC calibration not possible if magic_clean not a dictionary of MAGICClean objects", ): - _,_,_ = calibrate( - event=event, - tel_id=tel_id_MAGIC, - config=config_magic, - calibrator=calibrator_magic, + _, _, _ = calibrate( + event=event, + tel_id=tel_id_MAGIC, + config=config_magic, + calibrator=calibrator_magic, magic_clean=magic_clean, is_lst=False, ) - config_magic["image_extractor"]["type"]=extractor_type_magic + config_magic["image_extractor"]["type"] = extractor_type_magic def test_calibrate_exc_5(dl0_gamma, config_calib, tel_id_LST): - assigned_tel_ids = [1,2,3] + assigned_tel_ids = [1, 2, 3] for input_file in dl0_gamma: event_source = EventSource( - input_file, - allowed_tels=assigned_tel_ids, - focal_length_choice='effective' + input_file, allowed_tels=assigned_tel_ids, focal_length_choice="effective" ) + obs_id = event_source.obs_ids[0] - obs_id = event_source.obs_ids[0] - - subarray = event_source.subarray + subarray = event_source.subarray - tel_descriptions = subarray.tel - camera_geoms={} + tel_descriptions = subarray.tel + camera_geoms = {} for tel_id in range(len(tel_descriptions.items())): - camera_geoms[tel_id]= f"camera {tel_id}" - + camera_geoms[tel_id] = f"camera {tel_id}" + config_lst = config_calib["LST"] extractor_type_lst = config_lst["image_extractor"].pop("type") @@ -327,21 +303,20 @@ def test_calibrate_exc_5(dl0_gamma, config_calib, tel_id_LST): config=Config(config_extractor_lst), subarray=subarray, ) - + for event in event_source: - if (event.count <200) and (tel_id_LST in event.trigger.tels_with_trigger): + if (event.count < 200) and (tel_id_LST in event.trigger.tels_with_trigger): with pytest.raises( ValueError, match="Check the provided camera_geoms parameter; LST calibration not possible if camera_geoms not a dictionary of CameraGeometry objects", - ): - _,_,_ = calibrate( - event=event, - tel_id=tel_id_LST, + ): + _, _, _ = calibrate( + event=event, + tel_id=tel_id_LST, obs_id=obs_id, - config=config_lst, - camera_geoms=camera_geoms, - calibrator=calibrator_lst, - is_lst=True, + config=config_lst, + camera_geoms=camera_geoms, + calibrator=calibrator_lst, + is_lst=True, ) - config_lst["image_extractor"]["type"]=extractor_type_lst - + config_lst["image_extractor"]["type"] = extractor_type_lst diff --git a/magicctapipe/io/__init__.py b/magicctapipe/io/__init__.py index 4931f1717..dfa139b5a 100644 --- a/magicctapipe/io/__init__.py +++ b/magicctapipe/io/__init__.py @@ -6,6 +6,12 @@ RealEventInfoContainer, SimEventInfoContainer, ) +from .gadf import ( + create_event_hdu, + create_gh_cuts_hdu, + create_gti_hdu, + create_pointing_hdu, +) from .io import ( check_input_list, format_object, @@ -22,13 +28,6 @@ save_pandas_data_in_table, telescope_combinations, ) -from .gadf import ( - create_event_hdu, - create_gh_cuts_hdu, - create_gti_hdu, - create_pointing_hdu, -) - __all__ = [ "BaseEventInfoContainer", diff --git a/magicctapipe/io/gadf.py b/magicctapipe/io/gadf.py index 77b9664eb..ac7e58863 100644 --- a/magicctapipe/io/gadf.py +++ b/magicctapipe/io/gadf.py @@ -2,16 +2,19 @@ # coding: utf-8 import logging + import numpy as np from astropy import units as u from astropy.coordinates import SkyCoord from astropy.io import fits from astropy.table import QTable from astropy.time import Time +from pyirf.binning import split_bin_lo_hi + from magicctapipe import __version__ from magicctapipe.utils.functions import HEIGHT_ORM, LAT_ORM, LON_ORM -from pyirf.binning import split_bin_lo_hi -#from .io import telescope_combinations + +# from .io import telescope_combinations __all__ = [ "create_gh_cuts_hdu", @@ -125,11 +128,11 @@ def create_event_hdu( source RA/Dec coordinate is set to None """ TEL_COMBINATIONS = { - "M1_M2": [2, 3], # combo_type = 0 - "LST1_M1": [1, 2], # combo_type = 1 - "LST1_M2": [1, 3], # combo_type = 2 - "LST1_M1_M2": [1, 2, 3], # combo_type = 3 - } #TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) + "M1_M2": [2, 3], # combo_type = 0 + "LST1_M1": [1, 2], # combo_type = 1 + "LST1_M2": [1, 3], # combo_type = 2 + "LST1_M1_M2": [1, 2, 3], # combo_type = 3 + } # TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) mjdreff, mjdrefi = np.modf(MJDREF.mjd) diff --git a/magicctapipe/io/io.py b/magicctapipe/io/io.py index 3bef48fea..8cb27c354 100644 --- a/magicctapipe/io/io.py +++ b/magicctapipe/io/io.py @@ -26,7 +26,7 @@ __all__ = [ "check_input_list", "format_object", - "get_dl2_mean", + "get_dl2_mean", "get_stereo_events", "get_stereo_events_old", "load_dl2_data_file", @@ -60,13 +60,14 @@ DEAD_TIME_LST = 7.6 * u.us DEAD_TIME_MAGIC = 26 * u.us + def check_input_list(config): """ This function checks if the input telescope list is organized as follows: 1) All 4 LSTs and 2 MAGICs must be listed 2) All 4 LSTs must come before the MAGICs And it raises an exception in case these rules are not satisfied. - + Below we give two examples of valid lists: i) mc_tel_ids: @@ -103,24 +104,30 @@ def check_input_list(config): ------- This function will raise an exception if the input list is not properly organized. """ - + list_of_tel_names = list(config["mc_tel_ids"].keys()) standard_list_of_tels = ["LST-1", "LST-2", "LST-3", "LST-4", "MAGIC-I", "MAGIC-II"] if len(list_of_tel_names) != 6: - raise Exception(f"Number of telescopes found in the configuration file is {len(list_of_tel_names)}. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.") + raise Exception( + f"Number of telescopes found in the configuration file is {len(list_of_tel_names)}. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II." + ) else: for tel_name in list_of_tel_names[0:4]: if tel_name in standard_list_of_tels[0:4]: pass else: - raise Exception(f"Entry '{tel_name}' not accepted as an LST. Please make sure that the first four telescopes are LSTs, e.g.: 'LST-1', 'LST-2', 'LST-3', and 'LST-4'") + raise Exception( + f"Entry '{tel_name}' not accepted as an LST. Please make sure that the first four telescopes are LSTs, e.g.: 'LST-1', 'LST-2', 'LST-3', and 'LST-4'" + ) for tel_name in list_of_tel_names[4:6]: if tel_name in standard_list_of_tels[4:6]: pass else: - raise Exception(f"Entry '{tel_name}' not accepted as a MAGIC. Please make sure that the last two telescopes are MAGICs, e.g.: 'MAGIC-I', and 'MAGIC-II'") + raise Exception( + f"Entry '{tel_name}' not accepted as a MAGIC. Please make sure that the last two telescopes are MAGICs, e.g.: 'MAGIC-I', and 'MAGIC-II'" + ) return @@ -140,39 +147,49 @@ def telescope_combinations(config): TEL_COMBINATIONS: dict Dictionary with all telescope combinations with no repetions. """ - - + TEL_NAMES = {} - for k, v in config["mc_tel_ids"].items(): # Here we swap the dictionary keys and values just for convenience. + for k, v in config[ + "mc_tel_ids" + ].items(): # Here we swap the dictionary keys and values just for convenience. if v > 0: - TEL_NAMES[v] = k - + TEL_NAMES[v] = k + TEL_COMBINATIONS = {} keys = list(TEL_NAMES.keys()) - + def recursive_solution(current_tel, current_comb): - - if current_tel == len(keys): # The function stops once we reach the last telescope - return - - current_comb_name = current_comb[0] + '_' + TEL_NAMES[keys[current_tel]] # Name of the combo (at this point it can even be a single telescope) - current_comb_list = current_comb[1] + [keys[current_tel]] # List of telescopes (including individual telescopes) - - if len(current_comb_list) > 1: # We save them in the new dictionary excluding the single-telescope values - TEL_COMBINATIONS[current_comb_name[1:]] = current_comb_list; - - current_comb = [current_comb_name, current_comb_list] # We save the current results in this varible to recal the function recursively ("for" loop below) - - for i in range(1, len(keys)-current_tel): - recursive_solution(current_tel+i, current_comb) - - + + if current_tel == len( + keys + ): # The function stops once we reach the last telescope + return + + current_comb_name = ( + current_comb[0] + "_" + TEL_NAMES[keys[current_tel]] + ) # Name of the combo (at this point it can even be a single telescope) + current_comb_list = current_comb[1] + [ + keys[current_tel] + ] # List of telescopes (including individual telescopes) + + if ( + len(current_comb_list) > 1 + ): # We save them in the new dictionary excluding the single-telescope values + TEL_COMBINATIONS[current_comb_name[1:]] = current_comb_list + + current_comb = [ + current_comb_name, + current_comb_list, + ] # We save the current results in this varible to recal the function recursively ("for" loop below) + + for i in range(1, len(keys) - current_tel): + recursive_solution(current_tel + i, current_comb) + for key in range(len(keys)): - recursive_solution(key, ['',[]]) + recursive_solution(key, ["", []]) - return TEL_NAMES, TEL_COMBINATIONS - + def format_object(input_object): """ @@ -199,7 +216,7 @@ def format_object(input_object): string = string.replace("'", "").replace(",", "") return string - + def get_stereo_events_old( event_data, quality_cuts=None, group_index=["obs_id", "event_id"] @@ -225,11 +242,11 @@ def get_stereo_events_old( Data frame of the stereo events surviving the quality cuts """ TEL_COMBINATIONS = { - "M1_M2": [2, 3], # combo_type = 0 - "LST1_M1": [1, 2], # combo_type = 1 - "LST1_M2": [1, 3], # combo_type = 2 - "LST1_M1_M2": [1, 2, 3], # combo_type = 3 - } #TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) + "M1_M2": [2, 3], # combo_type = 0 + "LST1_M1": [1, 2], # combo_type = 1 + "LST1_M2": [1, 3], # combo_type = 2 + "LST1_M1_M2": [1, 2, 3], # combo_type = 3 + } # TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) event_data_stereo = event_data.copy() # Apply the quality cuts @@ -284,9 +301,12 @@ def get_stereo_events_old( return event_data_stereo - def get_stereo_events( - event_data, config, quality_cuts=None, group_index=["obs_id", "event_id"], eval_multi_combo=True + event_data, + config, + quality_cuts=None, + group_index=["obs_id", "event_id"], + eval_multi_combo=True, ): """ Gets the stereo events surviving specified quality cuts. @@ -298,7 +318,7 @@ def get_stereo_events( ---------- event_data: pandas.core.frame.DataFrame Data frame of shower events - config: dict + config: dict Read from the yaml file with information about the telescope IDs. quality_cuts: str Quality cuts applied to the input data @@ -306,26 +326,26 @@ def get_stereo_events( Index to group telescope events eval_multi_combo: bool If True, multiplicity is recalculated, combination type is assigned to each event and the fraction of events per combination type is shown - + Returns ------- event_data_stereo: pandas.core.frame.DataFrame Data frame of the stereo events surviving the quality cuts """ - + TEL_NAMES, TEL_COMBINATIONS = telescope_combinations(config) - + event_data_stereo = event_data.copy() # Apply the quality cuts if quality_cuts is not None: event_data_stereo.query(quality_cuts, inplace=True) - + # Extract stereo events event_data_stereo["multiplicity"] = event_data_stereo.groupby(group_index).size() event_data_stereo.query("multiplicity > 1", inplace=True) - if eval_multi_combo==True: + if eval_multi_combo == True: # Check the total number of events n_events_total = len(event_data_stereo.groupby(group_index).size()) logger.info(f"\nIn total {n_events_total} stereo events are found:") @@ -579,8 +599,8 @@ def load_magic_dl1_data_files(input_dir, config): ---------- input_dir: str Path to a directory where input MAGIC DL1 data files are stored - config: dict - yaml file with information about the telescope IDs. + config: dict + yaml file with information about the telescope IDs. Returns ------- @@ -594,9 +614,9 @@ def load_magic_dl1_data_files(input_dir, config): FileNotFoundError If any DL1 data files are not found in the input directory """ - + TEL_NAMES, _ = telescope_combinations(config) - + # Find the input files file_mask = f"{input_dir}/dl1_*.h5" @@ -680,11 +700,11 @@ def load_train_data_files( directory """ TEL_COMBINATIONS = { - "M1_M2": [2, 3], # combo_type = 0 - "LST1_M1": [1, 2], # combo_type = 1 - "LST1_M2": [1, 3], # combo_type = 2 - "LST1_M1_M2": [1, 2, 3], # combo_type = 3 - } #TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) + "M1_M2": [2, 3], # combo_type = 0 + "LST1_M1": [1, 2], # combo_type = 1 + "LST1_M2": [1, 3], # combo_type = 2 + "LST1_M1_M2": [1, 2, 3], # combo_type = 3 + } # TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) # Find the input files file_mask = f"{input_dir}/dl1_stereo_*.h5" @@ -736,7 +756,10 @@ def load_train_data_files( return data_train -def load_train_data_files_tel(input_dir, config, offaxis_min=None, offaxis_max=None, true_event_class=None): + +def load_train_data_files_tel( + input_dir, config, offaxis_min=None, offaxis_max=None, true_event_class=None +): """ Loads DL1-stereo data files and separates the shower events per telescope combination type for training RFs. @@ -745,8 +768,8 @@ def load_train_data_files_tel(input_dir, config, offaxis_min=None, offaxis_max=N ---------- input_dir: str Path to a directory where input DL1-stereo files are stored - config: dict - yaml file with information about the telescope IDs. + config: dict + yaml file with information about the telescope IDs. offaxis_min: str Minimum shower off-axis angle allowed, whose format should be acceptable by `astropy.units.quantity.Quantity` @@ -755,7 +778,7 @@ def load_train_data_files_tel(input_dir, config, offaxis_min=None, offaxis_max=N acceptable by `astropy.units.quantity.Quantity` true_event_class: int True event class of the input events - + Returns ------- @@ -769,9 +792,9 @@ def load_train_data_files_tel(input_dir, config, offaxis_min=None, offaxis_max=N If any DL1-stereo data files are not found in the input directory """ - + TEL_NAMES, _ = telescope_combinations(config) - + # Find the input files file_mask = f"{input_dir}/dl1_stereo_*.h5" @@ -857,7 +880,7 @@ def load_mc_dl2_data_file(input_file, quality_cuts, event_type, weight_type_dl2) ValueError If the input event type is not known """ - + # Load the input file df_events = pd.read_hdf(input_file, key="events/parameters") df_events.set_index(["obs_id", "event_id", "tel_id"], inplace=True) @@ -981,7 +1004,6 @@ def load_dl2_data_file(input_file, quality_cuts, event_type, weight_type_dl2): ValueError If the input event type is not known """ - # Load the input file event_data = pd.read_hdf(input_file, key="events/parameters") @@ -1078,6 +1100,7 @@ def load_dl2_data_file(input_file, quality_cuts, event_type, weight_type_dl2): return event_table, on_time, deadc + def load_irf_files(input_dir_irf): """ Loads input IRF data files for the IRF interpolation and checks the diff --git a/magicctapipe/io/tests/test_io.py b/magicctapipe/io/tests/test_io.py index c0e113a5f..59407e740 100644 --- a/magicctapipe/io/tests/test_io.py +++ b/magicctapipe/io/tests/test_io.py @@ -14,91 +14,55 @@ load_mc_dl2_data_file, load_train_data_files, load_train_data_files_tel, - load_mc_dl2_data_file, - load_irf_files, save_pandas_data_in_table, - load_magic_dl1_data_files, - load_lst_dl1_data_file, - load_dl2_data_file, telescope_combinations, ) -import pytest -import numpy as np -import pandas as pd - def test_check_input_list(config_check): """ Test on different dictionaries """ - - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':4, 'MAGIC-I':5, 'MAGIC-II':6}}) - - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':3, 'LST-3':0, 'LST-4':0, 'MAGIC-I':2, 'MAGIC-II':6}}) - - check_input_list({'mc_tel_ids':{'LST-2':1, 'LST-1':3, 'LST-4':0, 'LST-3':0, 'MAGIC-II':2, 'MAGIC-I':6}}) - - with pytest.raises( - Exception, - match="Number of telescopes found in the configuration file is 5. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.", - ): - check_input_list(config_check) - - with pytest.raises( - Exception, - match="Number of telescopes found in the configuration file is 5. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.", - ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'MAGIC-I':4, 'MAGIC-II':5}}) - - with pytest.raises( - Exception, - match="Number of telescopes found in the configuration file is 7. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.", - ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':6, 'LST-5':7, 'MAGIC-I':4, 'MAGIC-II':5}}) - - with pytest.raises( - Exception, - match="Entry 'LSTT-1' not accepted as an LST. Please make sure that the first four telescopes are LSTs, e.g.: 'LST-1', 'LST-2', 'LST-3', and 'LST-4'", - ): - check_input_list({'mc_tel_ids':{'LSTT-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':6, 'MAGIC-I':4, 'MAGIC-II':5}}) - - with pytest.raises( - Exception, - match="Entry 'MAGIC-III' not accepted as a MAGIC. Please make sure that the last two telescopes are MAGICs, e.g.: 'MAGIC-I', and 'MAGIC-II'", - ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':6, 'MAGIC-I':4, 'MAGIC-III':5}}) - - with pytest.raises( - Exception, - match="Entry 'MAGIC-I' not accepted as an LST. Please make sure that the first four telescopes are LSTs, e.g.: 'LST-1', 'LST-2', 'LST-3', and 'LST-4'", - ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'MAGIC-I':4, 'LST-3':3, 'LST-4':6, 'MAGIC-II':5}}) + check_input_list( + { + "mc_tel_ids": { + "LST-1": 1, + "LST-2": 2, + "LST-3": 3, + "LST-4": 4, + "MAGIC-I": 5, + "MAGIC-II": 6, + } + } + ) -def test_telescope_combinations(config_gen, config_gen_4lst): - """ - Simple check on telescope combinations - """ - M_LST, M_LST_comb = telescope_combinations(config_gen) - LSTs, LSTs_comb = telescope_combinations(config_gen_4lst) - assert M_LST == {1: 'LST-1', 2: 'MAGIC-I', 3: 'MAGIC-II'} - assert M_LST_comb == {'LST-1_MAGIC-I': [1, 2], 'LST-1_MAGIC-I_MAGIC-II': [1, 2, 3], 'LST-1_MAGIC-II': [1, 3], 'MAGIC-I_MAGIC-II': [2, 3]} - assert LSTs == {1: 'LST-1', 3: 'LST-2', 2: 'LST-3', 5: 'LST-4'} - assert LSTs_comb == {'LST-1_LST-2': [1, 3], 'LST-1_LST-2_LST-3': [1, 3, 2], 'LST-1_LST-2_LST-3_LST-4': [1, 3, 2, 5], 'LST-1_LST-2_LST-4': [1, 3, 5], 'LST-1_LST-3': [1, 2], 'LST-1_LST-3_LST-4': [1, 2, 5], 'LST-1_LST-4': [1, 5], 'LST-2_LST-3': [3, 2], 'LST-2_LST-3_LST-4': [3, 2, 5], 'LST-2_LST-4': [3, 5], 'LST-3_LST-4': [2, 5]} + check_input_list( + { + "mc_tel_ids": { + "LST-1": 1, + "LST-2": 3, + "LST-3": 0, + "LST-4": 0, + "MAGIC-I": 2, + "MAGIC-II": 6, + } + } + ) + check_input_list( + { + "mc_tel_ids": { + "LST-2": 1, + "LST-1": 3, + "LST-4": 0, + "LST-3": 0, + "MAGIC-II": 2, + "MAGIC-I": 6, + } + } + ) -def test_check_input_list(config_check): - """ - Test on different dictionaries - """ - - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':4, 'MAGIC-I':5, 'MAGIC-II':6}}) - - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':3, 'LST-3':0, 'LST-4':0, 'MAGIC-I':2, 'MAGIC-II':6}}) - - check_input_list({'mc_tel_ids':{'LST-2':1, 'LST-1':3, 'LST-4':0, 'LST-3':0, 'MAGIC-II':2, 'MAGIC-I':6}}) - with pytest.raises( Exception, match="Number of telescopes found in the configuration file is 5. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.", @@ -109,31 +73,86 @@ def test_check_input_list(config_check): Exception, match="Number of telescopes found in the configuration file is 5. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.", ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'MAGIC-I':4, 'MAGIC-II':5}}) + check_input_list( + { + "mc_tel_ids": { + "LST-1": 1, + "LST-2": 2, + "LST-3": 3, + "MAGIC-I": 4, + "MAGIC-II": 5, + } + } + ) with pytest.raises( Exception, match="Number of telescopes found in the configuration file is 7. It must be 6, i.e.: LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.", ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':6, 'LST-5':7, 'MAGIC-I':4, 'MAGIC-II':5}}) + check_input_list( + { + "mc_tel_ids": { + "LST-1": 1, + "LST-2": 2, + "LST-3": 3, + "LST-4": 6, + "LST-5": 7, + "MAGIC-I": 4, + "MAGIC-II": 5, + } + } + ) with pytest.raises( Exception, match="Entry 'LSTT-1' not accepted as an LST. Please make sure that the first four telescopes are LSTs, e.g.: 'LST-1', 'LST-2', 'LST-3', and 'LST-4'", ): - check_input_list({'mc_tel_ids':{'LSTT-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':6, 'MAGIC-I':4, 'MAGIC-II':5}}) + check_input_list( + { + "mc_tel_ids": { + "LSTT-1": 1, + "LST-2": 2, + "LST-3": 3, + "LST-4": 6, + "MAGIC-I": 4, + "MAGIC-II": 5, + } + } + ) with pytest.raises( Exception, match="Entry 'MAGIC-III' not accepted as a MAGIC. Please make sure that the last two telescopes are MAGICs, e.g.: 'MAGIC-I', and 'MAGIC-II'", ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'LST-3':3, 'LST-4':6, 'MAGIC-I':4, 'MAGIC-III':5}}) + check_input_list( + { + "mc_tel_ids": { + "LST-1": 1, + "LST-2": 2, + "LST-3": 3, + "LST-4": 6, + "MAGIC-I": 4, + "MAGIC-III": 5, + } + } + ) with pytest.raises( Exception, match="Entry 'MAGIC-I' not accepted as an LST. Please make sure that the first four telescopes are LSTs, e.g.: 'LST-1', 'LST-2', 'LST-3', and 'LST-4'", ): - check_input_list({'mc_tel_ids':{'LST-1':1, 'LST-2':2, 'MAGIC-I':4, 'LST-3':3, 'LST-4':6, 'MAGIC-II':5}}) + check_input_list( + { + "mc_tel_ids": { + "LST-1": 1, + "LST-2": 2, + "MAGIC-I": 4, + "LST-3": 3, + "LST-4": 6, + "MAGIC-II": 5, + } + } + ) def test_telescope_combinations(config_gen, config_gen_4lst): @@ -142,10 +161,27 @@ def test_telescope_combinations(config_gen, config_gen_4lst): """ M_LST, M_LST_comb = telescope_combinations(config_gen) LSTs, LSTs_comb = telescope_combinations(config_gen_4lst) - assert M_LST == {1: 'LST-1', 2: 'MAGIC-I', 3: 'MAGIC-II'} - assert M_LST_comb == {'LST-1_MAGIC-I': [1, 2], 'LST-1_MAGIC-I_MAGIC-II': [1, 2, 3], 'LST-1_MAGIC-II': [1, 3], 'MAGIC-I_MAGIC-II': [2, 3]} - assert LSTs == {1: 'LST-1', 3: 'LST-2', 2: 'LST-3', 5: 'LST-4'} - assert LSTs_comb == {'LST-1_LST-2': [1, 3], 'LST-1_LST-2_LST-3': [1, 3, 2], 'LST-1_LST-2_LST-3_LST-4': [1, 3, 2, 5], 'LST-1_LST-2_LST-4': [1, 3, 5], 'LST-1_LST-3': [1, 2], 'LST-1_LST-3_LST-4': [1, 2, 5], 'LST-1_LST-4': [1, 5], 'LST-2_LST-3': [3, 2], 'LST-2_LST-3_LST-4': [3, 2, 5], 'LST-2_LST-4': [3, 5], 'LST-3_LST-4': [2, 5]} + assert M_LST == {1: "LST-1", 2: "MAGIC-I", 3: "MAGIC-II"} + assert M_LST_comb == { + "LST-1_MAGIC-I": [1, 2], + "LST-1_MAGIC-I_MAGIC-II": [1, 2, 3], + "LST-1_MAGIC-II": [1, 3], + "MAGIC-I_MAGIC-II": [2, 3], + } + assert LSTs == {1: "LST-1", 3: "LST-2", 2: "LST-3", 5: "LST-4"} + assert LSTs_comb == { + "LST-1_LST-2": [1, 3], + "LST-1_LST-2_LST-3": [1, 3, 2], + "LST-1_LST-2_LST-3_LST-4": [1, 3, 2, 5], + "LST-1_LST-2_LST-4": [1, 3, 5], + "LST-1_LST-3": [1, 2], + "LST-1_LST-3_LST-4": [1, 2, 5], + "LST-1_LST-4": [1, 5], + "LST-2_LST-3": [3, 2], + "LST-2_LST-3_LST-4": [3, 2, 5], + "LST-2_LST-4": [3, 5], + "LST-3_LST-4": [2, 5], + } def test_format_object(): @@ -261,9 +297,9 @@ def test_load_train_data_files_tel_p(p_stereo, config_gen): Check dictionary """ - events = load_train_data_files_tel(str(p_stereo[0]),config_gen) - assert list(events.keys()) == [1,2,3] - data = events[2] + events = load_train_data_files_tel(str(p_stereo[0]), config_gen) + assert list(events.keys()) == [1, 2, 3] + data = events[2] assert "off_axis" in data.columns assert "true_event_class" not in data.columns @@ -274,8 +310,8 @@ def test_load_train_data_files_tel_g(gamma_stereo, config_gen): """ events = load_train_data_files_tel(str(gamma_stereo[0]), config_gen) - assert list(events.keys()) == [1,2,3] - data = events[1] + assert list(events.keys()) == [1, 2, 3] + data = events[1] assert "off_axis" in data.columns assert "true_event_class" not in data.columns @@ -285,7 +321,10 @@ def test_load_train_data_files_tel_off(gamma_stereo, config_gen): Check off-axis cut """ events = load_train_data_files_tel( - str(gamma_stereo[0]), config=config_gen, offaxis_min="0.2 deg", offaxis_max="0.5 deg" + str(gamma_stereo[0]), + config=config_gen, + offaxis_min="0.2 deg", + offaxis_max="0.5 deg", ) data = events[1] assert np.all(data["off_axis"] >= 0.2) @@ -399,7 +438,7 @@ def test_load_irf_files(IRF): """ Check on IRF dictionaries """ - + irf, header = load_irf_files(str(IRF)) assert set(list(irf.keys())).issubset( set( diff --git a/magicctapipe/io/tests/test_io_monly.py b/magicctapipe/io/tests/test_io_monly.py index e7d9b8e02..9ffbf298f 100644 --- a/magicctapipe/io/tests/test_io_monly.py +++ b/magicctapipe/io/tests/test_io_monly.py @@ -13,8 +13,6 @@ load_mc_dl2_data_file, load_train_data_files, load_train_data_files_tel, - load_mc_dl2_data_file, - load_irf_files, save_pandas_data_in_table, ) @@ -132,9 +130,9 @@ def test_load_train_data_files_tel_p(p_stereo_monly, config_gen): Check dictionary """ - events = load_train_data_files_tel(str(p_stereo_monly[0]),config_gen) - assert list(events.keys()) == [2,3] - data = events[2] + events = load_train_data_files_tel(str(p_stereo_monly[0]), config_gen) + assert list(events.keys()) == [2, 3] + data = events[2] assert "off_axis" in data.columns assert "true_event_class" not in data.columns @@ -145,8 +143,8 @@ def test_load_train_data_files_tel_g(gamma_stereo_monly, config_gen): """ events = load_train_data_files_tel(str(gamma_stereo_monly[0]), config_gen) - assert list(events.keys()) == [2,3] - data = events[3] + assert list(events.keys()) == [2, 3] + data = events[3] assert "off_axis" in data.columns assert "true_event_class" not in data.columns @@ -156,7 +154,10 @@ def test_load_train_data_files_tel_off(gamma_stereo_monly, config_gen): Check off-axis cut """ events = load_train_data_files_tel( - str(gamma_stereo_monly[0]), config=config_gen, offaxis_min="0.2 deg", offaxis_max="0.5 deg" + str(gamma_stereo_monly[0]), + config=config_gen, + offaxis_min="0.2 deg", + offaxis_max="0.5 deg", ) data = events[2] assert np.all(data["off_axis"] >= 0.2) @@ -384,7 +385,7 @@ def test_get_stereo_events_data(stereo_monly, config_gen): event_data.set_index(["obs_id", "event_id", "tel_id"], inplace=True) event_data.sort_index(inplace=True) data = get_stereo_events(event_data, config_gen) - assert np.all(data["multiplicity"] == 2) + assert np.all(data["multiplicity"] == 2) assert np.all(data["combo_type"] == 3) diff --git a/magicctapipe/reco/estimators.py b/magicctapipe/reco/estimators.py index 97caed90b..b4dd28b2f 100644 --- a/magicctapipe/reco/estimators.py +++ b/magicctapipe/reco/estimators.py @@ -8,13 +8,16 @@ import pandas as pd import sklearn.ensemble - __all__ = ["EnergyRegressor", "DispRegressor", "EventClassifier"] logger = logging.getLogger(__name__) logger.addHandler(logging.StreamHandler()) logger.setLevel(logging.INFO) -TEL_NAMES = {1: "LST-1", 2: "MAGIC-I", 3: "MAGIC-II"} #TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) +TEL_NAMES = { + 1: "LST-1", + 2: "MAGIC-I", + 3: "MAGIC-II", +} # TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) class EnergyRegressor: @@ -61,7 +64,7 @@ def fit(self, event_data): event_data: pandas.core.frame.DataFrame Data frame of shower events """ - + self.telescope_rfs.clear() # Loop over every telescope diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py b/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py index 70d21c2ff..4532fca1a 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_dl1_stereo_to_dl2.py @@ -29,6 +29,7 @@ from astropy.coordinates import AltAz, SkyCoord, angular_separation from ctapipe.coordinates import TelescopeFrame from ctapipe.instrument import SubarrayDescription + from magicctapipe.io import get_stereo_events_old, save_pandas_data_in_table from magicctapipe.reco import DispRegressor, EnergyRegressor, EventClassifier @@ -42,7 +43,8 @@ "LST1_M1": [1, 2], # combo_type = 1 "LST1_M2": [1, 3], # combo_type = 2 "LST1_M1_M2": [1, 2, 3], # combo_type = 3 -} #TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) +} # TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) + def apply_rfs(event_data, estimator): """ diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py b/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py index d57e2fd44..0280c289f 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_event_coincidence.py @@ -58,21 +58,23 @@ import time from decimal import Decimal from pathlib import Path + import numpy as np import pandas as pd import yaml from astropy import units as u from ctapipe.instrument import SubarrayDescription -from magicctapipe.io import ( + +from magicctapipe.io import ( + check_input_list, get_stereo_events, load_lst_dl1_data_file, load_magic_dl1_data_files, save_pandas_data_in_table, telescope_combinations, - check_input_list, ) -__all__ = ["event_coincidence","telescope_positions"] +__all__ = ["event_coincidence", "telescope_positions"] logger = logging.getLogger(__name__) logger.addHandler(logging.StreamHandler()) @@ -87,32 +89,32 @@ def telescope_positions(config): """ - This function computes the telescope positions with respect to the array baricenter. + This function computes the telescope positions with respect to the array baricenter. The array can have any configuration, e.g.: M1+M2+LST1+LST2, the full MAGIC+LST array, etc. - + Parameters ---------- config: dict dictionary generated from an yaml file with information about the telescope IDs. - + Returns ------- TEL_POSITIONS: dict Dictionary with telescope positions in the baricenter reference frame of the adopted array. """ - - #Telescope positions in meters in a generic reference frame: + + # Telescope positions in meters in a generic reference frame: RELATIVE_POSITIONS = { - "LST-1" : [-70.930, -52.070, 43.00], - "LST-2" : [-35.270, 66.140, 32.00], - "LST-3" : [75.280 , 50.490, 28.70], - "LST-4" : [30.910 , -64.540, 32.00], - "MAGIC-I" : [-23.540, -191.750, 41.25], - "MAGIC-II" : [-94.05, -143.770, 42.42] + "LST-1": [-70.930, -52.070, 43.00], + "LST-2": [-35.270, 66.140, 32.00], + "LST-3": [75.280, 50.490, 28.70], + "LST-4": [30.910, -64.540, 32.00], + "MAGIC-I": [-23.540, -191.750, 41.25], + "MAGIC-II": [-94.05, -143.770, 42.42], } telescopes_in_use = {} - tel_cp=config["mc_tel_ids"].copy() + tel_cp = config["mc_tel_ids"].copy() for k, v in tel_cp.copy().items(): if v <= 0: # Here we remove the telescopes with ID (i.e. "v") <= 0 from the dictionary: @@ -122,13 +124,15 @@ def telescope_positions(config): if k in RELATIVE_POSITIONS.keys(): telescopes_in_use[v] = RELATIVE_POSITIONS[k] else: - raise Exception(f"Telescope {k} not allowed in analysis. The telescopes allowed are LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II.") - - average_xyz=np.array([RELATIVE_POSITIONS[k] for k in tel_cp.keys()]).mean(axis=0) + raise Exception( + f"Telescope {k} not allowed in analysis. The telescopes allowed are LST-1, LST-2, LST-3, LST-4, MAGIC-I, and MAGIC-II." + ) + + average_xyz = np.array([RELATIVE_POSITIONS[k] for k in tel_cp.keys()]).mean(axis=0) TEL_POSITIONS = {} for k, v in telescopes_in_use.items(): - TEL_POSITIONS[k] = list(np.round(np.asarray(v)-average_xyz,2)) * u.m - return TEL_POSITIONS + TEL_POSITIONS[k] = list(np.round(np.asarray(v) - average_xyz, 2)) * u.m + return TEL_POSITIONS def event_coincidence(input_file_lst, input_dir_magic, output_dir, config): @@ -151,7 +155,7 @@ def event_coincidence(input_file_lst, input_dir_magic, output_dir, config): config_coinc = config["event_coincidence"] TEL_NAMES, _ = telescope_combinations(config) - + TEL_POSITIONS = telescope_positions(config) # Load the input LST DL1 data file logger.info(f"\nInput LST DL1 data file: {input_file_lst}") @@ -161,7 +165,9 @@ def event_coincidence(input_file_lst, input_dir_magic, output_dir, config): # Load the input MAGIC DL1 data files logger.info(f"\nInput MAGIC directory: {input_dir_magic}") - event_data_magic, subarray_magic = load_magic_dl1_data_files(input_dir_magic, config) + event_data_magic, subarray_magic = load_magic_dl1_data_files( + input_dir_magic, config + ) # Exclude the parameters non-common to LST and MAGIC data timestamp_type_lst = config_coinc["timestamp_type_lst"] @@ -183,20 +189,24 @@ def event_coincidence(input_file_lst, input_dir_magic, output_dir, config): window_half_width = u.Quantity(window_half_width).to("ns") window_half_width = u.Quantity(window_half_width.round(), dtype=int) pre_offset_search = False - if "pre_offset_search" in config_coinc: #looking for the boolean value of pre_offset_search in the configuration file + if ( + "pre_offset_search" in config_coinc + ): # looking for the boolean value of pre_offset_search in the configuration file pre_offset_search = config_coinc["pre_offset_search"] if pre_offset_search: logger.info("\nPre offset search will be performed.") - n_pre_offset_search_events = config_coinc["n_pre_offset_search_events"] #n_pre_offset_search_events is the number of events used to estimate the time offset. Around 100 events may be enough + n_pre_offset_search_events = config_coinc[ + "n_pre_offset_search_events" + ] # n_pre_offset_search_events is the number of events used to estimate the time offset. Around 100 events may be enough else: logger.info("\noffset scan range defined in the config file will be used.") offset_start = u.Quantity(config_coinc["time_offset"]["start"]) offset_stop = u.Quantity(config_coinc["time_offset"]["stop"]) - + event_data = pd.DataFrame() features = pd.DataFrame() - + profiles = pd.DataFrame(data={"time_offset": []}) # Arrange the LST timestamps. They are stored in the UNIX format in @@ -394,16 +404,12 @@ def event_coincidence(input_file_lst, input_dir_magic, output_dir, config): ) n_coincidences.append(n_coincidence) - - - if not any(n_coincidences): logger.info("\nNo coincident events are found. Skipping...") continue n_coincidences = np.array(n_coincidences) - # Sometimes there are more than one time offset maximizing the # number of coincidences, so here we calculate the mean of them @@ -550,14 +556,16 @@ def event_coincidence(input_file_lst, input_dir_magic, output_dir, config): # Create the subarray description with the telescope coordinates # relative to the center of the LST and MAGIC positions - tel_descriptions = {} - for k, v in TEL_NAMES.items(): + tel_descriptions = {} + for k, v in TEL_NAMES.items(): if v[:3] == "LST": tel_descriptions[k] = subarray_lst.tel[k] elif v[:5] == "MAGIC": tel_descriptions[k] = subarray_magic.tel[k] else: - raise Exception(f"{v} is not a valid telescope name (check the config file). Only MAGIC and LST telescopes can be analyzed --> Valid telescope names are LST-[1-4] and MAGIC-[I-II] ") + raise Exception( + f"{v} is not a valid telescope name (check the config file). Only MAGIC and LST telescopes can be analyzed --> Valid telescope names are LST-[1-4] and MAGIC-[I-II] " + ) subarray_lst_magic = SubarrayDescription( "LST-MAGIC-Array", TEL_POSITIONS, tel_descriptions diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py b/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py index 44d282328..51fba65bb 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_mc_dl0_to_dl1.py @@ -21,8 +21,8 @@ (--config-file config_step1.yaml) """ -import argparse #Parser for command-line options, arguments etc -import logging #Used to manage the log file +import argparse # Parser for command-line options, arguments etc +import logging # Used to manage the log file import re import time from pathlib import Path @@ -40,10 +40,11 @@ ) from ctapipe.instrument import SubarrayDescription from ctapipe.io import EventSource, HDF5TableWriter +from traitlets.config import Config from magicctapipe.image import MAGICClean from magicctapipe.image.calib import calibrate -from magicctapipe.io import SimEventInfoContainer, format_object, check_input_list +from magicctapipe.io import SimEventInfoContainer, check_input_list, format_object from magicctapipe.utils import calculate_disp, calculate_impact __all__ = ["mc_dl0_to_dl1"] @@ -56,9 +57,6 @@ PARTICLE_TYPES = {1: "gamma", 3: "electron", 14: "proton", 402: "helium"} - - - def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): """ Processes LST-1 and MAGIC events of simtel MC DL0 data and computes @@ -74,9 +72,13 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): Configuration for the LST-1 + MAGIC analysis """ - assigned_tel_ids = config["mc_tel_ids"] #This variable becomes the dictionary {'LST-1': 1, 'MAGIC-I': 2, 'MAGIC-II': 3} + assigned_tel_ids = config[ + "mc_tel_ids" + ] # This variable becomes the dictionary {'LST-1': 1, 'MAGIC-I': 2, 'MAGIC-II': 3} - logger.info("\nAssigned telescope IDs:") #Here we are just adding infos to the log file + logger.info( + "\nAssigned telescope IDs:" + ) # Here we are just adding infos to the log file logger.info(format_object(assigned_tel_ids)) # Load the input file @@ -84,7 +86,9 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): event_source = EventSource( input_file, - allowed_tels=list(filter(lambda check_id: check_id > 0, assigned_tel_ids.values())), #Here we load the events for all telescopes with ID > 0. + allowed_tels=list( + filter(lambda check_id: check_id > 0, assigned_tel_ids.values()) + ), # Here we load the events for all telescopes with ID > 0. focal_length_choice=focal_length, ) @@ -122,8 +126,6 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): logger.info("\nLST PSF modifier:") logger.info(format_object(config_lst["increase_psf"])) - - logger.info("\nLST tailcuts cleaning:") logger.info(format_object(config_lst["tailcuts_clean"])) @@ -133,8 +135,6 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): logger.info("\nLST dynamic cleaning:") logger.info(format_object(config_lst["dynamic_cleaning"])) - - use_only_main_island = config_lst["use_only_main_island"] logger.info(f"\nLST use only main island: {use_only_main_island}") @@ -181,29 +181,33 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): azimuth = Angle(sim_config["max_az"]).wrap_at("360 deg").degree logger.info(np.asarray(list(assigned_tel_ids.values()))) LSTs_IDs = np.asarray(list(assigned_tel_ids.values())[0:4]) - LSTs_in_use = np.where(LSTs_IDs > 0)[0] + 1 #Here we select which LSTs are/is in use - + LSTs_in_use = ( + np.where(LSTs_IDs > 0)[0] + 1 + ) # Here we select which LSTs are/is in use + if len(LSTs_in_use) == 0: - LSTs_in_use = ''.join(str(k) for k in LSTs_in_use) + LSTs_in_use = "".join(str(k) for k in LSTs_in_use) elif len(LSTs_in_use) > 0: - LSTs_in_use = 'LST'+'_LST'.join(str(k) for k in LSTs_in_use) - + LSTs_in_use = "LST" + "_LST".join(str(k) for k in LSTs_in_use) + MAGICs_IDs = np.asarray(list(assigned_tel_ids.values())[4:6]) - MAGICs_in_use = np.where(MAGICs_IDs > 0)[0] + 1 #Here we select which MAGICs are/is in use - + MAGICs_in_use = ( + np.where(MAGICs_IDs > 0)[0] + 1 + ) # Here we select which MAGICs are/is in use + if len(MAGICs_in_use) == 0: - MAGICs_in_use = ''.join(str(k) for k in MAGICs_in_use) + MAGICs_in_use = "".join(str(k) for k in MAGICs_in_use) elif len(MAGICs_in_use) > 0: - MAGICs_in_use = 'MAGIC'+'_MAGIC'.join(str(k) for k in MAGICs_in_use) + MAGICs_in_use = "MAGIC" + "_MAGIC".join(str(k) for k in MAGICs_in_use) magic_clean = {} for k in MAGICs_IDs: - if k > 0: + if k > 0: magic_clean[k] = MAGICClean(camera_geoms[k], config_magic["magic_clean"]) - + output_file = ( f"{output_dir}/dl1_{particle_type}_zd_{zenith.round(3)}deg_" f"az_{azimuth.round(3)}deg_{LSTs_in_use}_{MAGICs_in_use}_run{obs_id}.h5" - ) #The files are saved with the names of all telescopes involved + ) # The files are saved with the names of all telescopes involved # Loop over every shower event logger.info("\nProcessing the events...") @@ -216,21 +220,42 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): tels_with_trigger = event.trigger.tels_with_trigger # Check if the event triggers both M1 and M2 or not - magic_stereo=(set(MAGICs_IDs).issubset(set(tels_with_trigger))) and (MAGICs_in_use=="MAGIC1_MAGIC2") #If both have trigger, then magic_stereo = True - - for tel_id in tels_with_trigger: + magic_stereo = (set(MAGICs_IDs).issubset(set(tels_with_trigger))) and ( + MAGICs_in_use == "MAGIC1_MAGIC2" + ) # If both have trigger, then magic_stereo = True + + for tel_id in tels_with_trigger: - if tel_id in LSTs_IDs: ##If the ID is in the LST list, we call calibrate on the LST() + if ( + tel_id in LSTs_IDs + ): # If the ID is in the LST list, we call calibrate on the LST() # Calibrate the LST-1 event - signal_pixels, image, peak_time = calibrate(event=event, tel_id=tel_id, obs_id=obs_id, config=config_lst, camera_geoms=camera_geoms, calibrator=calibrator_lst, is_lst=True) + signal_pixels, image, peak_time = calibrate( + event=event, + tel_id=tel_id, + obs_id=obs_id, + config=config_lst, + camera_geoms=camera_geoms, + calibrator=calibrator_lst, + is_lst=True, + ) elif tel_id in MAGICs_IDs: # Calibrate the MAGIC event - signal_pixels, image, peak_time = calibrate(event=event, tel_id=tel_id, config=config_magic, magic_clean=magic_clean, calibrator=calibrator_magic, is_lst=False) + signal_pixels, image, peak_time = calibrate( + event=event, + tel_id=tel_id, + config=config_magic, + magic_clean=magic_clean, + calibrator=calibrator_magic, + is_lst=False, + ) else: logger.info( f"--> Telescope ID {tel_id} not in LST list or MAGIC list. Please check if the IDs are OK in the configuration file" ) - if not any(signal_pixels): #So: if there is no event, we skip it and go back to the loop in the next event + if not any( + signal_pixels + ): # So: if there is no event, we skip it and go back to the loop in the next event logger.info( f"--> {event.count} event (event ID: {event.index.event_id}, " f"telescope {tel_id}) could not survive the image cleaning. " @@ -256,7 +281,7 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): # Parametrize the image hillas_params = hillas_parameters(camera_geom_masked, image_masked) - # + # if any(np.isnan(value) for value in hillas_params.values()): logger.info( f"--> {event.count} event (event ID: {event.index.event_id}, " @@ -264,7 +289,6 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): ) continue - timing_params = timing_parameters( camera_geom_masked, image_masked, peak_time_masked, hillas_params ) @@ -327,10 +351,10 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): n_islands=n_islands, magic_stereo=magic_stereo, ) - + # Reset the telescope IDs event_info.tel_id = tel_id - + # Save the parameters to the output file writer.write( "parameters", @@ -351,7 +375,6 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): for k in IDs_in_use: tel_positions_lst_magic[k] = tel_positions[k] - position_mean tel_descriptions_lst_magic[k] = tel_descriptions[k] - subarray_lst_magic = SubarrayDescription( "LST-MAGIC-Array", tel_positions_lst_magic, tel_descriptions_lst_magic @@ -360,10 +383,10 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): logger.info("\nTelescope positions:") logger.info(format_object(tel_positions)) - + # Save the subarray description subarray_lst_magic.to_hdf(output_file) - + # Save the simulation configuration with HDF5TableWriter(output_file, group_name="simulation", mode="a") as writer: writer.write("config", sim_config) @@ -373,14 +396,13 @@ def mc_dl0_to_dl1(input_file, output_dir, config, focal_length): def main(): - """ Here we collect the input parameters from the command line, load the configuration file and run mc_dl0_to_dl1()""" - + """Here we collect the input parameters from the command line, load the configuration file and run mc_dl0_to_dl1()""" + start_time = time.time() parser = argparse.ArgumentParser() - - - #Here we are simply collecting the parameters from the command line, as input file, output directory, and configuration file + + # Here we are simply collecting the parameters from the command line, as input file, output directory, and configuration file parser.add_argument( "--input-file", "-i", @@ -407,26 +429,32 @@ def main(): default="./config.yaml", help="Path to a configuration file", ) - + parser.add_argument( "--focal_length_choice", - "-f", + "-f", dest="focal_length_choice", type=str, - choices = ["nominal", "effective"], + choices=["nominal", "effective"], default="effective", help='Choice of focal length, either "effective" or "nominal". The default (and standard) value is "effective"', ) - args = parser.parse_args() #Here we select all 3 parameters collected above + args = parser.parse_args() # Here we select all 3 parameters collected above - with open(args.config_file, "rb") as f: # "rb" mode opens the file in binary format for reading - config = yaml.safe_load(f) #Here we collect the inputs from the configuration file + with open( + args.config_file, "rb" + ) as f: # "rb" mode opens the file in binary format for reading + config = yaml.safe_load( + f + ) # Here we collect the inputs from the configuration file # Checking if the input telescope list is properly organized: check_input_list(config) - config['mc_tel_ids']=dict(sorted(config['mc_tel_ids'].items())) #Sorting needed to correctly name the output file + config["mc_tel_ids"] = dict( + sorted(config["mc_tel_ids"].items()) + ) # Sorting needed to correctly name the output file # Process the input data mc_dl0_to_dl1(args.input_file, args.output_dir, config, args.focal_length_choice) @@ -438,4 +466,4 @@ def main(): if __name__ == "__main__": - main() \ No newline at end of file + main() diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py b/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py index 277de8d43..1d47010e1 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_stereo_reco.py @@ -44,7 +44,13 @@ ) from ctapipe.instrument import SubarrayDescription from ctapipe.reco import HillasReconstructor -from magicctapipe.io import format_object, get_stereo_events, save_pandas_data_in_table, check_input_list + +from magicctapipe.io import ( + check_input_list, + format_object, + get_stereo_events, + save_pandas_data_in_table, +) from magicctapipe.utils import calculate_impact, calculate_mean_direction __all__ = ["calculate_pointing_separation", "stereo_reconstruction"] @@ -71,12 +77,14 @@ def calculate_pointing_separation(event_data, config): Angular distance of the LST array and MAGIC pointing directions in units of degree """ - - assigned_tel_ids = config["mc_tel_ids"] #This variable becomes a dictionary, e.g.: {'LST-1': 1, 'LST-2': 0, 'LST-3': 0, 'LST-4': 0, 'MAGIC-I': 2, 'MAGIC-II': 3} + + assigned_tel_ids = config[ + "mc_tel_ids" + ] # This variable becomes a dictionary, e.g.: {'LST-1': 1, 'LST-2': 0, 'LST-3': 0, 'LST-4': 0, 'MAGIC-I': 2, 'MAGIC-II': 3} LSTs_IDs = np.asarray(list(assigned_tel_ids.values())[0:4]) - LSTs_IDs = list(LSTs_IDs[LSTs_IDs > 0]) #Here we list only the LSTs in use + LSTs_IDs = list(LSTs_IDs[LSTs_IDs > 0]) # Here we list only the LSTs in use MAGICs_IDs = np.asarray(list(assigned_tel_ids.values())[4:6]) - MAGICs_IDs = list(MAGICs_IDs[MAGICs_IDs > 0]) #Here we list only the MAGICs in use + MAGICs_IDs = list(MAGICs_IDs[MAGICs_IDs > 0]) # Here we list only the MAGICs in use # Extract LST events df_lst = event_data.query(f"tel_id == {LSTs_IDs}") @@ -86,8 +94,12 @@ def calculate_pointing_separation(event_data, config): df_magic = df_magic.loc[df_lst.index] # Calculate the mean of the LSTs, and also of the M1 and M2 pointing directions - pnt_az_LST, pnt_alt_LST = calculate_mean_direction(lon=df_lst["pointing_az"], lat=df_lst["pointing_alt"], unit="rad") - pnt_az_magic, pnt_alt_magic = calculate_mean_direction(lon=df_magic["pointing_az"], lat=df_magic["pointing_alt"], unit="rad") + pnt_az_LST, pnt_alt_LST = calculate_mean_direction( + lon=df_lst["pointing_az"], lat=df_lst["pointing_alt"], unit="rad" + ) + pnt_az_magic, pnt_alt_magic = calculate_mean_direction( + lon=df_magic["pointing_az"], lat=df_magic["pointing_alt"], unit="rad" + ) # Calculate the angular distance of their pointing directions theta = angular_separation( @@ -121,8 +133,10 @@ def stereo_reconstruction(input_file, output_dir, config, magic_only_analysis=Fa """ config_stereo = config["stereo_reco"] - assigned_tel_ids = config["mc_tel_ids"] #This variable becomes a dictionary, e.g.: {'LST-1': 1, 'LST-2': 0, 'LST-3': 0, 'LST-4': 0, 'MAGIC-I': 2, 'MAGIC-II': 3} - + assigned_tel_ids = config[ + "mc_tel_ids" + ] # This variable becomes a dictionary, e.g.: {'LST-1': 1, 'LST-2': 0, 'LST-3': 0, 'LST-4': 0, 'MAGIC-I': 2, 'MAGIC-II': 3} + # Load the input file logger.info(f"\nInput file: {input_file}") @@ -152,22 +166,26 @@ def stereo_reconstruction(input_file, output_dir, config, magic_only_analysis=Fa LSTs_IDs = np.asarray(list(assigned_tel_ids.values())[0:4]) if magic_only_analysis: - tel_id=np.asarray(list(assigned_tel_ids.values())[:]) - used_id=tel_id[tel_id!=0] - magic_ids=[item for item in used_id if item not in LSTs_IDs] - event_data.query(f"tel_id in {magic_ids}", inplace=True) # Here we select only the events with the MAGIC tel_ids + tel_id = np.asarray(list(assigned_tel_ids.values())[:]) + used_id = tel_id[tel_id != 0] + magic_ids = [item for item in used_id if item not in LSTs_IDs] + event_data.query( + f"tel_id in {magic_ids}", inplace=True + ) # Here we select only the events with the MAGIC tel_ids logger.info(f"\nQuality cuts: {config_stereo['quality_cuts']}") - event_data = get_stereo_events(event_data, config=config, quality_cuts=config_stereo["quality_cuts"]) + event_data = get_stereo_events( + event_data, config=config, quality_cuts=config_stereo["quality_cuts"] + ) # Check the angular distance of the LST and MAGIC pointing directions tel_ids = np.unique(event_data.index.get_level_values("tel_id")).tolist() - Number_of_LSTs_in_use = len(LSTs_IDs[LSTs_IDs > 0]) + Number_of_LSTs_in_use = len(LSTs_IDs[LSTs_IDs > 0]) MAGICs_IDs = np.asarray(list(assigned_tel_ids.values())[4:6]) Number_of_MAGICs_in_use = len(MAGICs_IDs[MAGICs_IDs > 0]) - Two_arrays_are_used = (Number_of_LSTs_in_use*Number_of_MAGICs_in_use > 0) - + Two_arrays_are_used = Number_of_LSTs_in_use * Number_of_MAGICs_in_use > 0 + if (not is_simulation) and (Two_arrays_are_used): logger.info( @@ -384,10 +402,10 @@ def main(): with open(args.config_file, "rb") as f: config = yaml.safe_load(f) - + # Checking if the input telescope list is properly organized: check_input_list(config) - + # Process the input data stereo_reconstruction(args.input_file, args.output_dir, config, args.magic_only) @@ -398,4 +416,4 @@ def main(): if __name__ == "__main__": - main() \ No newline at end of file + main() diff --git a/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py b/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py index 9a6b70171..b95a98ccc 100644 --- a/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py +++ b/magicctapipe/scripts/lst1_magic/lst1_magic_train_rfs.py @@ -55,7 +55,11 @@ logger = logging.getLogger(__name__) logger.addHandler(logging.StreamHandler()) logger.setLevel(logging.INFO) -TEL_NAMES = {1: "LST-1", 2: "MAGIC-I", 3: "MAGIC-II"} #TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) +TEL_NAMES = { + 1: "LST-1", + 2: "MAGIC-I", + 3: "MAGIC-II", +} # TODO: REMOVE WHEN SWITCHING TO THE NEW RFs IMPLEMENTTATION (1 RF PER TELESCOPE) # True event class of gamma and proton MCs diff --git a/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py b/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py index 707579fbe..971f67041 100644 --- a/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py +++ b/magicctapipe/scripts/lst1_magic/magic_calib_to_dl1.py @@ -52,7 +52,12 @@ from ctapipe_io_magic import MAGICEventSource from magicctapipe.image import MAGICClean -from magicctapipe.io import RealEventInfoContainer, SimEventInfoContainer, format_object, check_input_list +from magicctapipe.io import ( + RealEventInfoContainer, + SimEventInfoContainer, + check_input_list, + format_object, +) from magicctapipe.utils import calculate_disp, calculate_impact __all__ = ["magic_calib_to_dl1"] @@ -394,7 +399,9 @@ def main(): check_input_list(config) # Process the input data - magic_calib_to_dl1(args.input_file, args.output_dir, config, args.max_events, args.process_run) + magic_calib_to_dl1( + args.input_file, args.output_dir, config, args.max_events, args.process_run + ) logger.info("\nDone.") process_time = time.time() - start_time diff --git a/magicctapipe/scripts/lst1_magic/merge_hdf_files.py b/magicctapipe/scripts/lst1_magic/merge_hdf_files.py index 7b747935d..56376382e 100644 --- a/magicctapipe/scripts/lst1_magic/merge_hdf_files.py +++ b/magicctapipe/scripts/lst1_magic/merge_hdf_files.py @@ -38,7 +38,7 @@ import tables from ctapipe.instrument import SubarrayDescription -__all__ = ["merge_hdf_files","write_data_to_table"] +__all__ = ["merge_hdf_files", "write_data_to_table"] logger = logging.getLogger(__name__) logger.addHandler(logging.StreamHandler()) diff --git a/setup.py b/setup.py index 1ae1c4cf1..b283ec672 100644 --- a/setup.py +++ b/setup.py @@ -35,20 +35,20 @@ use_scm_version={"write_to": os.path.join("magicctapipe", "_version.py")}, packages=find_packages(), install_requires=[ - 'lstchain~=0.9.6', - 'ctapipe~=0.12.0', - 'ctapipe_io_magic~=0.4.7', - 'ctaplot~=0.5.5', - 'gammapy~=0.19.0', - 'uproot~=4.1', - 'numba<=0.56', - 'numpy<1.22.0a0', - 'joblib', - 'pandas', - 'pyirf~=0.6.0', - 'seaborn', - 'scikit-learn', - 'setuptools_scm', + "lstchain~=0.9.6", + "ctapipe~=0.12.0", + "ctapipe_io_magic~=0.4.7", + "ctaplot~=0.5.5", + "gammapy~=0.19.0", + "uproot~=4.1", + "numba<=0.56", + "numpy<1.22.0a0", + "joblib", + "pandas", + "pyirf~=0.6.0", + "seaborn", + "scikit-learn", + "setuptools_scm", ], entry_points=entry_points, extras_require={ From 1c1e2d9e707baff579c3fa463e55093c24c74005 Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Sun, 29 Oct 2023 13:32:45 +0100 Subject: [PATCH 07/10] Add link to documentation. --- README.md | 1 + 1 file changed, 1 insertion(+) diff --git a/README.md b/README.md index ae223d5a7..2941a6729 100644 --- a/README.md +++ b/README.md @@ -3,6 +3,7 @@ Repository for the analysis of MAGIC and MAGIC+LST1 data, based on [*ctapipe*](https://github.com/cta-observatory/ctapipe). * Code: https://github.com/cta-observatory/magic-cta-pipe +* Docs (preliminary): https://magic-cta-pipe.readthedocs.io/ v0.3.1 of *magic-cta-pipe* provides all the functionalities to perform a MAGIC+LST-1 or a MAGIC-only analysis. Both types of analyses can be performed using the scripts within the *lst1_magic* folder. See the [README](https://github.com/cta-observatory/magic-cta-pipe/blob/master/magicctapipe/scripts/lst1_magic/README.md) for more details on how to run the analysis. From acff21639336647b088e8a452b1d49df8a081d7f Mon Sep 17 00:00:00 2001 From: Alessio Berti Date: Sun, 29 Oct 2023 13:32:59 +0100 Subject: [PATCH 08/10] Add info about pre-commit hook. --- README.md | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/README.md b/README.md index 2941a6729..2af8d8476 100644 --- a/README.md +++ b/README.md @@ -25,6 +25,8 @@ The following command will set up a conda virtual environment, add the necessary conda activate magic-lst1 pip install . +In general, *magic-cta-pipe* is still in heavy development phase, so expect large changes between different releases. + # Instructions for developers People who would like to join the development of *magic-cta-pipe*, please contact Alessio Berti () to get write access to the repository. @@ -41,6 +43,14 @@ flake8 magicctapipe # checks style and code errors black filename.py # reformats filename.py with black ``` +The *black* and *isort* auto-formatters are used for automatic adherence to the code style. To enforce running these tools whenever you make a commit, setup the [pre-commit hook](https://pre-commit.com/): + +```bash +$ pre-commit install +``` + +The pre-commit hook will then execute the tools with the same settings as when the a pull request is checked on github, and if any problems are reported the commit will be rejected. You then have to fix the reported issues before tying to commit again. + In general, if you want to add a new feature or fix a bug, please open a new issue, and then create a new branch to develop the new feature or code the bug fix. You can create an early pull request even if it is not complete yet, you can tag it as "Draft" so that it will not be merged, and other developers can already check it and provide comments. When the code is ready, remove the tag "Draft" and select two people to review the pull request (at the moment the merge is not blocked if no review is performed, but that may change in the future). When the review is complete, the branch will be merged into the main branch.