Skip to content

Irf gammaness efficiency offset #317

New issue

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

Sign up for GitHub

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

Already on GitHub? Sign in to your account

Open
jsitarek wants to merge 4 commits into
base: master
Choose a base branch
from
Open

Irf gammaness efficiency offset #317

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
1273559
fixed the naming of POINT-LIKE vs FULL-ENCLOSURE IRFs
jsitarek Dec 3, 2025
9c9fa6f
lowered the absolute tolerance for the test since now the computed az…
jsitarek Dec 5, 2025
90a47a1
reverting a change that went into a wrong branch
jsitarek Dec 5, 2025
d09bf43
applying offset cuts when computing cut efficiencies
jsitarek Dec 5, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 6 additions & 2 deletions magicctapipe/io/gadf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,9 @@


@u.quantity_input(reco_energy_bins=u.TeV, fov_offset_bins=u.deg)
def create_gh_cuts_hdu(gh_cuts, reco_energy_bins, fov_offset_bins, **header_cards):
def create_gh_cuts_hdu(
gh_cuts, reco_energy_bins, fov_offset_bins, point_like, **header_cards
):
"""
Creates a fits binary table HDU for dynamic gammaness cuts.

Expand All @@ -44,6 +46,8 @@ def create_gh_cuts_hdu(gh_cuts, reco_energy_bins, fov_offset_bins, **header_card
Bin edges in the reconstructed energy
fov_offset_bins : astropy.units.quantity.Quantity
Bin edges in the field of view offset
point_like : bool
Meaning: true = IRFs are point-like, false = IRFS are full-enclosure
**header_cards : dict
Additional metadata to add to the header

Expand Down Expand Up @@ -73,7 +77,7 @@ def create_gh_cuts_hdu(gh_cuts, reco_energy_bins, fov_offset_bins, **header_card
("CREATOR", f"magicctapipe v{MCP_VERSION}"),
("HDUCLAS1", "RESPONSE"),
("HDUCLAS2", "GH_CUTS"),
("HDUCLAS3", "POINT-LIKE"),
("HDUCLAS3", "POINT-LIKE" if point_like else "FULL-ENCLOSURE"),
("HDUCLAS4", "GH_CUTS_2D"),
("DATE", Time.now().utc.iso),
]
Expand Down
4 changes: 3 additions & 1 deletion magicctapipe/io/tests/test_gadf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -60,7 +60,9 @@ def event():


def test_create_gh_cuts_hdu(gammaness_cut, energy_bins, fov_bins, header):
g_cuts_fits = create_gh_cuts_hdu(gammaness_cut, energy_bins, fov_bins, **header)
g_cuts_fits = create_gh_cuts_hdu(
gammaness_cut, energy_bins, fov_bins, True, **header
)
g_cuts = g_cuts_fits.data
assert np.allclose(g_cuts["ENERG_LO"][0], np.array([0.1, 1.0, 10.0]))
assert np.allclose(g_cuts["ENERG_HI"][0], np.array([1.0, 10.0, 100.0]))
Expand Down
20 changes: 15 additions & 5 deletions magicctapipe/scripts/lst1_magic/lst1_magic_create_irf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -366,10 +366,12 @@ def create_irf(

# Calculate dynamic gammaness cuts
gh_percentile = 100 * (1 - gh_efficiency)

withinfov = (
event_table_gamma["true_source_fov_offset"].to("deg") < fov_offset_bins[-1]
) * (event_table_gamma["true_source_fov_offset"].to("deg") > fov_offset_bins[0])
cut_table_gh = calculate_percentile_cut(
values=event_table_gamma["gammaness"],
bin_values=event_table_gamma["reco_energy"],
values=event_table_gamma[withinfov]["gammaness"],
bin_values=event_table_gamma[withinfov]["reco_energy"],
bins=energy_bins,
fill_value=gh_cut_min,
percentile=gh_percentile,
Expand Down Expand Up @@ -423,6 +425,7 @@ def create_irf(
gh_cuts=gh_cuts,
reco_energy_bins=energy_bins,
fov_offset_bins=fov_offset_bins,
point_like=is_point_like,
**extra_header,
)

Expand Down Expand Up @@ -472,9 +475,16 @@ def create_irf(
# Calculate dynamic theta cuts
theta_percentile = 100 * theta_efficiency

withinfov = (
event_table_gamma["true_source_fov_offset"].to("deg")
< fov_offset_bins[-1]
) * (
event_table_gamma["true_source_fov_offset"].to("deg")
> fov_offset_bins[0]
)
cut_table_theta = calculate_percentile_cut(
values=event_table_gamma["theta"],
bin_values=event_table_gamma["reco_energy"],
values=event_table_gamma[withinfov]["theta"],
bin_values=event_table_gamma[withinfov]["reco_energy"],
bins=energy_bins,
fill_value=theta_cut_max,
percentile=theta_percentile,
Expand Down
6 changes: 4 additions & 2 deletions magicctapipe/scripts/lst1_magic/lst1_magic_dl2_to_dl3.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -206,11 +206,12 @@ def dl2_to_dl3(input_file_dl2, input_dir_irf, output_dir, config):
aeff_interp = irf_data["effective_area"][0]
print("skipping interpolation since only one point is given")

point_like = irf_data["effective_area"].shape[-1] == 1
aeff_hdu = create_aeff2d_hdu(
effective_area=aeff_interp,
true_energy_bins=irf_data["energy_bins"],
fov_offset_bins=irf_data["fov_offset_bins"],
point_like=True,
point_like=point_like,
extname="EFFECTIVE AREA",
**extra_header,
)
Expand Down Expand Up @@ -247,7 +248,7 @@ def dl2_to_dl3(input_file_dl2, input_dir_irf, output_dir, config):
true_energy_bins=irf_data["energy_bins"],
migration_bins=irf_data["migration_bins"],
fov_offset_bins=irf_data["fov_offset_bins"],
point_like=True,
point_like=point_like,
extname="ENERGY DISPERSION",
)

Expand Down Expand Up @@ -328,6 +329,7 @@ def dl2_to_dl3(input_file_dl2, input_dir_irf, output_dir, config):
gh_cuts=gh_cuts_interp,
reco_energy_bins=irf_data["energy_bins"],
fov_offset_bins=irf_data["fov_offset_bins"],
point_like=point_like,
**extra_header,
)

Expand Down