Skip to content
New issue

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

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

Fix the shape of the interpolated Effective Area #1270

Merged
merged 1 commit into from
Jun 26, 2024
Merged

Fix the shape of the interpolated Effective Area #1270

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
24 changes: 12 additions & 12 deletions lstchain/high_level/interpolate.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -422,7 +422,7 @@ def interpolate_irf(irfs, data_pars, interp_method="linear"):
main_headers = fits.open(irfs[0])[1].header

point_like = main_headers["HDUCLAS3"] == "POINT-LIKE"

# Update headers to be added to the final IRFs
extra_headers = dict()

Expand Down Expand Up @@ -471,19 +471,19 @@ def interpolate_irf(irfs, data_pars, interp_method="linear"):
}

for hdu_name in hdu_names:

if hdu_name == "BACKGROUND":
log.warning("The interpolation of BACKGROUND is not yet supported")
continue

if interp_col_keys[hdu_name] == 'cut':
gadf_irf = False
else:
gadf_irf = True
gadf_irf = True

irf_list = load_irf_grid(
irfs,
extname=hdu_name,
irfs,
extname=hdu_name,
interp_col=interp_col_keys[hdu_name],
gadf_irf=gadf_irf
)
Expand All @@ -492,17 +492,17 @@ def interpolate_irf(irfs, data_pars, interp_method="linear"):
if hdu_name in ['EFFECTIVE AREA', 'ENERGY DISPERSION', 'PSF']:
e_true = np.append(temp_irf["ENERG_LO"][0], temp_irf["ENERG_HI"][0][-1])
fov_off = np.append(temp_irf["THETA_LO"][0], temp_irf["THETA_HI"][0][-1])

if hdu_name == "EFFECTIVE AREA":
aeff_estimator = EffectiveAreaEstimator(
grid_points=irf_pars_sel,
effective_area=irf_list,
interpolator_kwargs={"method": interp_method},
)
aeff_interp = aeff_estimator(interp_pars_sel)

aeff_hdu_interp = create_aeff2d_hdu(
effective_area=aeff_interp.T[0],
effective_area=aeff_interp[0],
true_energy_bins=e_true,
fov_offset_bins=fov_off,
point_like=point_like,
Expand All @@ -519,7 +519,7 @@ def interpolate_irf(irfs, data_pars, interp_method="linear"):
energy_dispersion=irf_list,
)
edisp_interp = edisp_estimator(interp_pars_sel)

edisp_hdu_interp = create_energy_dispersion_hdu(
energy_dispersion=edisp_interp[0],
true_energy_bins=e_true,
Expand Down Expand Up @@ -556,7 +556,7 @@ def interpolate_irf(irfs, data_pars, interp_method="linear"):
cut_header = fits.Header()
cut_header["CREATOR"] = f"lstchain v{__version__}"
cut_header["DATE"] = Time.now().utc.iso

for k, v in extra_headers.items():
cut_header[k] = v

Expand Down
10 changes: 7 additions & 3 deletions lstchain/high_level/tests/test_interpolate.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -257,17 +257,17 @@ def test_interp_irf(
psf_1_meta = fits.open(irf)["PSF"].header
psf_2 = psf_1.copy()
psf_2_meta = psf_1_meta.copy()

psf_1["RPSF"][0] *= factor_czd
psf_2["RPSF"][0] *= factor_czd * factor_sd

psf_1_meta["ZEN_PNT"] = (zen_2 * 180 / np.pi, "deg")
psf_1_meta["AZ_PNT"] = (az_1 * 180 / np.pi, "deg")
psf_1_meta["B_DELTA"] = (del_1 * 180 / np.pi, "deg")
psf_2_meta["ZEN_PNT"] = (zen_2 * 180 / np.pi, "deg")
psf_2_meta["AZ_PNT"] = (az_2 * 180 / np.pi, "deg")
psf_2_meta["B_DELTA"] = (del_2 * 180 / np.pi, "deg")

psf_hdu_1 = fits.BinTableHDU(
psf_1, header=psf_1_meta, name="PSF"
)
Expand Down Expand Up @@ -333,10 +333,14 @@ def test_interp_irf(
hdu_en_srcdep = interpolate_irf(irfs_en_srcdep, data_pars)
hdu_en_srcdep.writeto(irf_file_en_srcdep_final, overwrite=True)

aeff_shape_final = Table.read(irf_file_g_final, hdu=1)["EFFAREA"].shape
aeff_shape_2 = Table.read(irf_file_g_2, hdu=1)["EFFAREA"].shape

assert hdu_g[1].header["ZEN_PNT"] == zen_t
assert irf_file_g_2.exists()
assert irf_file_g_3.exists()
assert irf_file_g_final.exists()
assert aeff_shape_final == aeff_shape_2

assert hdu_en[1].header["ZEN_PNT"] == zen_t
assert irf_file_en_2.exists()
Expand Down