### Setup conda development environment

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/install.rst

Create a dedicated conda environment from the provided configuration file and install the package in editable mode.

```bash
$ conda env create -f environment.yml
$ conda activate pyirf
$ pip install -e '.[all]'
```

--------------------------------

### Run unit tests

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/install.rst

Execute the test suite to verify the installation.

```bash
$ pytest
```

--------------------------------

### Install pyirf for development

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/install.rst

Install the package in editable mode with all optional dependencies for testing and documentation.

```bash
$ pip install -e '.[all]'  # or [docs,tests] to install them separately
```

--------------------------------

### Initialize Extrapolator with MomentMorphNearestSimplexExtrapolator

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Demonstrates initializing a MomentMorphNearestSimplexExtrapolator. This snippet is incomplete and likely part of a larger example.

```python
import numpy as np
from pyirf.interpolation import MomentMorphNearestSimplexExtrapolator
from scipy.stats import norm 

bins = np.linspace(-10, 10, 51)
grid = np.array([[1], [2], [3]])

```

--------------------------------

### Install released pyirf version

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/install.rst

Use pip to install the latest released version of the package.

```bash
$ pip install pyirf
```

--------------------------------

### Download Private Data for pyirf

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/examples.rst

Downloads required private data files for EventDisplay examples. Requires curl and unzip utilities.

```bash
./download_private_data.sh
```

--------------------------------

### Serve documentation

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/contribute.rst

Serve the built documentation locally to preview changes in a browser.

```bash
$ python -m http.server _build/html
```

--------------------------------

### Build documentation

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/contribute.rst

Commands to build the project documentation using Sphinx.

```bash
make html
```

```bash
make clean html
```

--------------------------------

### Import Binning Utilities

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Imports functions for creating energy and field-of-view offset bins, essential for IRF calculations.

```python
from pyirf.binning import create_bins_per_decade, bin_center
```

--------------------------------

### Instantiate and Use GaussianEstimator

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Demonstrates how to create a GaussianEstimator with predefined grid points, bin edges, and pre-calculated Gaussians. The estimator is then used to predict Gaussian values at target points.

```python
gaussians = np.array([np.diff(norm(loc=x, scale=1/x).cdf(bins))/np.diff(bins) for x in grid])

estimator = GaussianEstimatior(
   grid_points = grid,
   bin_edges = bins,
   gaussians = gaussians * u.m,
   interpolator_cls = MomentMorphInterpolator,
   extrapolator_cls = MomentMorphNearestSimplexExtrapolator
)

targets = np.array([[0.9], [1.5]])

results = u.Quantity([estimator(target).squeeze() for target in targets])
```

--------------------------------

### pyirf.io Module Overview

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/io/index.rst

Overview of the I/O module capabilities for handling GADF-compliant data and EventDisplay conversions.

```APIDOC
## pyirf.io Module

### Description
The `pyirf.io` module provides functions to read input data and write IRFs in the GADF (Gamma Astro Data Format) standard. It is designed to facilitate the processing of DL2 data files.

### Supported Formats
- **EventDisplay DL2 FITS**: Supported via conversion from ROOT files using the official EventDisplay conversion scripts.
```

--------------------------------

### Manage Energy Binning

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Creates logarithmic energy bins and handles overflow/underflow. Useful for preparing data for histogramming.

```python
import numpy as np
import astropy.units as u
from pyirf.binning import (
    create_bins_per_decade,
    add_overflow_bins,
    calculate_bin_indices,
    create_histogram_table,
    bin_center,
)
from astropy.table import QTable

# Create bins with 5 bins per decade
energy_bins = create_bins_per_decade(
    e_min=0.01 * u.TeV,
    e_max=100 * u.TeV,
    bins_per_decade=5,
)
print(f"Number of bins: {len(energy_bins) - 1}")

# Add under/overflow bins (0 to inf)
energy_bins_overflow = add_overflow_bins(energy_bins, positive=True)
print(f"First bin edge: {energy_bins_overflow[0]}")  # 0 TeV
print(f"Last bin edge: {energy_bins_overflow[-1]}")  # inf TeV

# Get bin centers (geometric mean for log spacing)
centers = bin_center(energy_bins)

# Calculate bin indices for events
energies = np.array([0.05, 0.5, 5.0, 50.0]) * u.TeV
bin_indices, valid = calculate_bin_indices(energies, energy_bins)
print(f"Bin indices: {bin_indices}")
print(f"Valid (in range): {valid}")

# Create histogram table from events
events = QTable({
    'reco_energy': np.random.uniform(0.1, 10, 1000) * u.TeV,
    'weight': np.ones(1000),
})
hist = create_histogram_table(events, energy_bins, key='reco_energy')
print(f"Histogram columns: {hist.colnames}")
```

--------------------------------

### Display First 10 Events

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Displays the first 10 rows of the loaded gammas QTable to provide a preview of the event data.

```python
gammas[:10]
```

--------------------------------

### Run tests with coverage

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/contribute.rst

Run tests and generate an HTML coverage report in the htmlcov directory.

```bash
$ pytest --cov=pyirf --cov-report=html -v
```

--------------------------------

### Run tests locally

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/contribute.rst

Execute the test suite locally to verify changes before pushing.

```bash
$ pytest -v
```

--------------------------------

### Download test IRFs

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/install.rst

Download required CTA IRF files from Zenodo to the local directory.

```bash
$ python download_irfs.py 
```

--------------------------------

### Serve coverage report

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/contribute.rst

Serve the generated HTML coverage report locally using a Python HTTP server.

```bash
$ python -m http.server -d htmlcov
```

--------------------------------

### GaussianEstimator Usage

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Demonstrates how to create and use a GaussianEstimator to estimate Gaussian distributions at arbitrary target points.

```APIDOC
## GaussianEstimator Usage

### Description
This section shows how to initialize a `GaussianEstimator` with predefined grid points, bin edges, and pre-calculated Gaussian distributions. It then demonstrates how to use the estimator to predict Gaussian values for new target points.

### Method
N/A (Illustrative Code)

### Endpoint
N/A (Local Object Usage)

### Request Body
```json
{
  "grid_points": "numpy.ndarray",
  "bin_edges": "numpy.ndarray",
  "gaussians": "astropy.units.Quantity",
  "interpolator_cls": "class",
  "extrapolator_cls": "class"
}
```

### Request Example
```python
import numpy as np
import astropy.units as u
from scipy.stats import norm
from pyirf.interpolation import GaussianEstimator, MomentMorphInterpolator, MomentMorphNearestSimplexExtrapolator

# Assuming grid and bins are defined elsewhere
grid = np.linspace(0, 10, 100)
bins = np.linspace(0, 10, 50)

gaussians = np.array([np.diff(norm(loc=x, scale=1/x).cdf(bins))/np.diff(bins) for x in grid])

estimator = GaussianEstimatior(
   grid_points = grid,
   bin_edges = bins,
   gaussians = gaussians * u.m, # Example unit
   interpolator_cls = MomentMorphInterpolator,
   extrapolator_cls = MomentMorphNearestSimplexExtrapolator
)

targets = np.array([[0.9], [1.5]])

results = u.Quantity([estimator(target).squeeze() for target in targets])
```

### Response
#### Success Response (200)
- **results** (astropy.units.Quantity) - The estimated Gaussian values for the target points.
```

--------------------------------

### Download and Unpack Data

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Use curl and unzip to download and extract the necessary data files for the analysis.

```bash
$ curl -fL -o data.zip https://nextcloud.e5.physik.tu-dortmund.de/index.php/s/Cstsf8MWZjnz92L/download
$ unzip data.zip
$ mv eventdisplay_dl2 data
```

--------------------------------

### Manage Simulation Metadata with SimulatedEventsInfo

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Use SimulatedEventsInfo to store MC production parameters and calculate expected shower counts per energy or FOV bin.

```python
import astropy.units as u
import numpy as np
from pyirf.simulations import SimulatedEventsInfo

# Create simulation info from MC production parameters
simulation_info = SimulatedEventsInfo(
    n_showers=1_000_000,
    energy_min=0.01 * u.TeV,
    energy_max=100 * u.TeV,
    max_impact=2000 * u.m,
    spectral_index=-2.0,
    viewcone_min=0 * u.deg,
    viewcone_max=5 * u.deg,
)

# Calculate expected number of showers per energy bin
energy_bins = np.geomspace(0.01, 100, 21) * u.TeV
n_showers_per_bin = simulation_info.calculate_n_showers_per_energy(energy_bins)
print(f"Showers in first bin: {n_showers_per_bin[0]:.1f}")

# Calculate showers per energy and FOV offset (for 2D effective area)
fov_bins = np.linspace(0, 5, 6) * u.deg
n_showers_2d = simulation_info.calculate_n_showers_per_energy_and_fov(energy_bins, fov_bins)
print(f"Shape: {n_showers_2d.shape}")  # (n_energy_bins, n_fov_bins)
```

--------------------------------

### pyirf.cut_optimization API

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/cut_optimization.rst

Reference documentation for the cut_optimization module.

```APIDOC
## pyirf.cut_optimization

### Description
This module provides tools for optimizing event selection cuts, which are essential for separating gamma-ray signals from background events in Cherenkov telescope data.

### Reference
For detailed class and function signatures, please refer to the official pyirf documentation generated via automodapi for the `pyirf.cut_optimization` namespace.
```

--------------------------------

### Write a Changelog Message

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/changes/README.md

Write a suitable message for your change. This message will appear in the changelog.

```text
Fixed ``crazy_function`` to be consistent with ``not_so_crazy_function``
```

--------------------------------

### Generate Changelog with Towncrier

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/changes/README.md

Execute this command in the base directory of the project to generate the changelog. Replace <VERSION NUMBER> with the actual version.

```bash
towncrier build --version=<VERSION NUMBER>
```

--------------------------------

### pyirf.binning Module API

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/binning.rst

Reference for the pyirf.binning module, covering binning and histogram utilities.

```APIDOC
## pyirf.binning Module

### Description
Provides utilities for binning data and creating histograms.

### Methods
This module contains various functions for binning and histogram creation. Please refer to the specific function documentation for details on parameters and usage.

### Example Usage
```python
# Example of using a binning function (specific function not detailed in source)
# from pyirf.binning import create_binning
# bins = create_binning(edges)
```

### Available Functions (Refer to source for full details)
- Functions for creating and manipulating bins.
- Functions for generating histograms from data.

### Error Handling
Refer to individual function documentation for specific error handling.
```

--------------------------------

### Create Background and RAD_MAX HDUs

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Generate background rate and RAD_MAX HDUs for inclusion in a GADF-compliant FITS file.

```python
bg_rate = np.random.uniform(1e-6, 1e-4, (n_energy, n_fov)) / (u.MeV * u.s * u.sr)
hdus.append(create_background_2d_hdu(
    bg_rate, reco_energy_bins, fov_offset_bins,
))

# RAD_MAX (theta cut) HDU
rad_max = np.full((n_energy, n_fov), 0.2) * u.deg
hdus.append(create_rad_max_hdu(
    rad_max, reco_energy_bins, fov_offset_bins,
))

# Write FITS file
hdulist = fits.HDUList(hdus)
# hdulist.writeto('irfs.fits.gz', overwrite=True)
print(f"Created {len(hdus)} HDUs for GADF-compliant IRF file")
```

--------------------------------

### Import Necessary Libraries

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Imports core libraries for numerical operations, units, plotting, and subprocess management.

```python
import numpy as np
import astropy.units as u
import matplotlib.pyplot as plt
import subprocess as sp
```

--------------------------------

### Define Simulated Event Information

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Creates a SimulatedEventsInfo object to store parameters of simulated events, including energy range, spectral index, and shower properties. This is crucial for subsequent IRF calculations.

```python
from pyirf.simulations import SimulatedEventsInfo

simulation_info = SimulatedEventsInfo(
    energy_min=200 * u.GeV,
    energy_max=50 * u.TeV,
    spectral_index=-2.7,
    n_showers=12600000,
    max_impact=300 * u.m,
    viewcone_min=0 * u.deg,
    viewcone_max=0 * u.deg,
)
```

--------------------------------

### Load EventDisplay IRF Data

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Opens and reads EventDisplay IRF data from a ROOT file using uproot. Ensure the 'indir' variable points to the correct data directory.

```python
# Path of EventDisplay IRF data in the user's local setup
# Please, empty the indir_EventDisplay variable before pushing to the repo
indir = "../../data/"
irf_file_event_display = "DESY.d20180113.V3.ID0NIM2LST4MST4SST4SCMST4.prod3b-paranal20degs05b-NN.S.3HB9-FD.180000s.root"

irf_eventdisplay = uproot.open(os.path.join(indir, irf_file_event_display))
```

--------------------------------

### Plot Gamma/Hadron Cut Comparison

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Compares the gamma/hadron (G/H) cut between pyirf and EventDisplay as a function of reconstructed energy. Requires the 'GH_CUTS' HDU from the pyirf file.

```python
from astropy.table import QTable


gh_cut = QTable.read(pyirf_file, hdu='GH_CUTS')[1:-1]


plt.errorbar(
    0.5 * (gh_cut['low'] + gh_cut['high']).to_value(u.TeV),
    gh_cut['cut'],
    xerr=0.5 * (gh_cut['high'] - gh_cut['low']).to_value(u.TeV),
    ls='',
    label='pyirf',
)

plt.legend()
plt.ylabel('G/H-cut')
plt.xlabel(r'$E_\mathrm{reco} / \mathrm{TeV}$')
plt.xscale('log')

None # to remove clutter by mpl objects
```

--------------------------------

### Load and Estimate Energy Dispersions

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Loads IRF data and uses EnergyDispersionEstimator to interpolate and extrapolate energy dispersions. Requires IRF data files to be present.

```python
import numpy as np

from gammapy.irf import load_irf_dict_from_file
from glob import glob
from pyirf.interpolation import (
   EnergyDispersionEstimator, 
   MomentMorphInterpolator,
   DiscretePDFNearestNeighborSearcher
)

# Load IRF data, replace path with actual path
PROD5_IRF_PATH = "pyirf/irfs/*.fits.gz"

irfs = [load_irf_dict_from_file(path) for path in sorted(glob(PROD5_IRF_PATH))]

edisps = np.array([irf["edisp"].quantity for irf in irfs])
bin_edges = irfs[0]["edisp"].axes["migra"].edges
# IRFs are for zenith distances of 20, 40 and 60 deg 
zen_pnt = np.array([[20], [40], [60]])

# Create estimator instance
edisp_estimator = EnergyDispersionEstimator(
     grid_points=zen_pnt,
     migra_bins=bin_edges,
     energy_dispersion=edisps,
     interpolator_cls=MomentMorphInterpolator,
     interpolator_kwargs=None,
     extrapolator_cls=DiscretePDFNearestNeighborSearcher,
     extrapolator_kwargs=None,
 )

 # Estimate energy dispersions
interpolated_edisp = edisp_estimator(np.array([[30]]))
extrapolated_edisp = edisp_estimator(np.array([[10]]))

```

--------------------------------

### pyirf.utils Module Reference

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/utils.rst

Reference documentation for the utility functions contained within the pyirf.utils module.

```APIDOC
## pyirf.utils

### Description
The pyirf.utils module contains a collection of utility functions designed to assist with common tasks in the pyirf library, such as data handling and mathematical operations relevant to CTA Observatory data analysis.
```

--------------------------------

### Apply Style Settings for Plots

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Sets plot limits, scales, labels, and legend for astronomical data visualization. Ensures consistent styling across plots.

```python
plt.xlim(1.e-2, 2.e2)
plt.ylim(2.e-2, 1)
plt.xscale("log")
plt.yscale("log")
plt.xlabel("Reconstructed energy / TeV")
plt.ylabel("Angular Resolution / deg")
plt.grid(which="both")
plt.legend(loc="best")

None # to remove clutter by mpl objects
```

--------------------------------

### Calculate EventDisplay IRFs

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/examples.rst

Executes the script to optimize cuts and calculate sensitivity and IRFs for EventDisplay DL2 data.

```bash
python examples/calculate_eventdisplay_irfs.py
```

--------------------------------

### pyirf.spectral Module Reference

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/spectral.rst

Reference documentation for the spectral module in pyirf, covering event weighting and spectrum definitions.

```APIDOC
## pyirf.spectral

### Description
The spectral module provides tools for event weighting and defining spectra within the pyirf framework. It includes functions for calculating weights based on target spectra and instrument response functions.

### Reference
This module is automatically generated from the pyirf.spectral source code. It includes all objects defined within the module, excluding inheritance diagrams.
```

--------------------------------

### pyirf.sensitivity Module Reference

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/sensitivity.rst

Reference documentation for the sensitivity module in pyirf.

```APIDOC
## pyirf.sensitivity

### Description
This module provides tools for calculating the sensitivity of Cherenkov Telescope Array (CTA) observations. It is part of the pyirf library for Instrument Response Function (IRF) analysis.
```

--------------------------------

### Base Classes

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Lists available base classes for estimators and interpolators.

```APIDOC
## Base Classes

### Description
This section lists the base classes provided for creating custom estimators, interpolators, and extrapolators within the `pyirf.interpolation` module.

### Autosummary
- BaseComponentEstimator
- ParametrizedComponentEstimator
- DiscretePDFComponentEstimator
- BaseInterpolator
- ParametrizedInterpolator
- DiscretePDFInterpolator
- BaseExtrapolator
- ParametrizedExtrapolator
- DiscretePDFExtrapolator
- BaseNearestNeighborSearcher
```

--------------------------------

### Plot Background Rate Comparison

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Compares background rate data from EventDisplay and PyIRF. This involves data preprocessing, including calculating solid angles and interpolating theta cuts. Uses errorbar plots with logarithmic x and y scales.

```python
from pyirf.utils import cone_solid_angle

# Data from EventDisplay
y, edges = irf_eventdisplay["BGRate"].to_numpy()
yerr = irf_eventdisplay["BGRate"].errors()
x = bin_center(10**edges)
xerr = np.diff(10**edges) / 2

# pyirf data

bg_rate = QTable.read(pyirf_file, hdu='BACKGROUND')[0]

reco_bins = np.append(bg_rate['ENERG_LO'], bg_rate['ENERG_HI'][-1])

# first fov bin, [0, 1] deg
fov_bin = 0
rate_bin = bg_rate['BKG'].T[:, fov_bin]


# interpolate theta cut for given e reco bin
e_center_bg = 0.5 * (bg_rate['ENERG_LO'] + bg_rate['ENERG_HI'])
e_center_theta = 0.5 * (rad_max['ENERG_LO'] + rad_max['ENERG_HI'])
theta_cut = np.interp(e_center_bg, e_center_theta, rad_max['RAD_MAX'].T[:, 0])

# undo normalization
rate_bin *= cone_solid_angle(theta_cut)
rate_bin *= np.diff(reco_bins)


# Plot function
plt.errorbar(x, y, xerr=xerr, yerr=yerr, ls='', label="EventDisplay")

plt.errorbar(
    0.5 * (bg_rate['ENERG_LO'] + bg_rate['ENERG_HI']).to_value(u.TeV)[1:-1],
    rate_bin.to_value(1 / u.s)[1:-1],
    xerr=np.diff(reco_bins).to_value(u.TeV)[1:-1] / 2,
    ls='',
    label='pyirf',
)

# Style settings
plt.xscale("log")
plt.xlabel(r"$E_\mathrm{Reco} / \mathrm{TeV}$")
plt.ylabel("Background rate / (s⁻¹) ")
plt.grid(which="both")
plt.legend(loc="best")
plt.yscale('log')

None # to remove clutter by mpl objects
```

--------------------------------

### Inter- and Extrapolation Algorithms

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Summary of available classes for inter- and extrapolation, categorized by component type and grid dimensionality support.

```APIDOC
## Inter- and Extrapolation Algorithms

### Description
Classes used to perform the mathematical operations for IRF estimation. These are plugged into the estimator classes.

### Parametrized Components (AEFF, Rad-Max)
- **GridDataInterpolator**: Interpolation, Arbitrary dimensions.
- **ParametrizedNearestSimplexExtrapolator**: Extrapolation, 1D or 2D.
- **ParametrizedVisibleEdgesExtrapolator**: Extrapolation, 1D or 2D.
- **ParametrizedNearestNeighborSearcher**: Nearest Neighbor, Arbitrary dimensions.

### Discretized PDFs (PSF, EDISP)
- **QuantileInterpolator**: Interpolation, Arbitrary dimensions.
- **MomentMorphInterpolator**: Interpolation, 1D or 2D.
- **MomentMorphNearestSimplexExtrapolator**: Extrapolation, 1D or 2D.
- **DiscretePDFNearestNeighborSearcher**: Nearest Neighbor, Arbitrary dimensions.
```

--------------------------------

### Define Energy and Offset Bins

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Defines the binning scheme for true energy and field-of-view offsets. Energy bins are created per decade, and a single offset bin is defined around the FACT wobble offset.

```python
true_energy_bins = create_bins_per_decade(simulation_info.energy_min, simulation_info.energy_max, 5)

# single offset bin around the wobble distance
# since we are dealing with point-like simulations 
wobble_offset = 0.6 * u.deg
fov_offset_bins = [0.59, 0.61] * u.deg
```

--------------------------------

### Export IRFs to GADF FITS Files

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Exports IRF objects (effective area, energy dispersion, PSF) to GADF FITS files using `pyirf.io.gadf` and `astropy.io.fits`. Common metadata like creator, telescope, instrument, and date are included.

```python
from pyirf.io.gadf import create_aeff2d_hdu, create_energy_dispersion_hdu, create_psf_table_hdu
from astropy.io import fits
from astropy.time import Time
from pyirf import __version__

# set some common meta data for all hdus
meta = dict(
    CREATOR='pyirf-v' + __version__,
    TELESCOP='FACT',
    INSTRUME='FACT',
    DATE=Time.now().iso,
)

hdus = []

# every fits file has to have an Image HDU as first HDU.
# GADF only uses Binary Table HDUs, so we need to add an empty HDU in front
hdus.append(fits.PrimaryHDU(header=fits.Header(meta)))

hdus.append(create_aeff2d_hdu(aeff_selected, true_energy_bins, fov_offset_bins, **meta))
hdus.append(create_energy_dispersion_hdu(edisp, true_energy_bins, migration_bins, fov_offset_bins, **meta))
hdus.append(create_psf_table_hdu(psf, true_energy_bins,    source_offset_bins, fov_offset_bins, **meta))

fits.HDUList(hdus).writeto('fact_irf.fits.gz', overwrite=True)
```

--------------------------------

### Helper Classes

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Lists available helper classes within the interpolation module.

```APIDOC
## Helper Classes

### Description
This section lists the helper classes available in the `pyirf.interpolation` module.

### Autosummary
- PDFNormalization
```

--------------------------------

### Optimize Gamma/Hadron Separation Cut

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Optimizes gamma/hadron separation cuts for best sensitivity by scanning over cut efficiencies. Requires signal and background event tables, energy bins, and pre-computed theta cuts.

```python
import numpy as np
import astropy.units as u
from astropy.table import QTable
import operator
from pyirf.cut_optimization import optimize_gh_cut
from pyirf.cuts import calculate_percentile_cut
from pyirf.binning import add_overflow_bins, create_bins_per_decade

# Signal (gamma) events
signal = QTable({
    'reco_energy': np.random.uniform(0.1, 100, 5000) * u.TeV,
    'theta': np.abs(np.random.normal(0, 0.1, 5000)) * u.deg,
    'gh_score': np.random.uniform(0.3, 1.0, 5000),  # Gamma-like
    'weight': np.ones(5000) * 0.5,
})

# Background (hadron) events
background = QTable({
    'reco_energy': np.random.uniform(0.1, 100, 20000) * u.TeV,
    'reco_source_fov_offset': np.random.uniform(0, 3, 20000) * u.deg,
    'gh_score': np.random.uniform(0, 0.7, 20000),  # Hadron-like
    'weight': np.ones(20000) * 0.1,
})

# Energy bins
reco_energy_bins = add_overflow_bins(
    create_bins_per_decade(0.1 * u.TeV, 100 * u.TeV, bins_per_decade=5)
)

# Pre-computed theta cuts (68% containment)
theta_cuts = calculate_percentile_cut(
    signal['theta'], signal['reco_energy'], reco_energy_bins,
    fill_value=0.3 * u.deg, percentile=68,
)

# Optimize G/H cut for best sensitivity
gh_efficiencies = np.arange(0.1, 0.9, 0.1)
best_sensitivity, best_gh_cuts = optimize_gh_cut(
    signal=signal,
    background=background,
    reco_energy_bins=reco_energy_bins,
    gh_cut_efficiencies=gh_efficiencies,
    theta_cuts=theta_cuts,
    op=operator.ge,  # Keep events with score >= cut
    alpha=0.2,
    fov_offset_max=3 * u.deg,
    progress=False,
)

print(f"Best cuts columns: {best_gh_cuts.colnames}")
# ['low', 'center', 'high', 'cut', 'efficiency']
print(f"Sensitivity columns: {best_sensitivity.colnames}")
```

--------------------------------

### Plot Energy Dispersion Matrix

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Visualizes the energy dispersion matrix from PyIRF data. Requires reading data from a specified HDU and uses pcolormesh for display. Sets log scales for both axes and adds a colorbar.

```python
edisp = QTable.read(pyirf_file, hdu='ENERGY_DISPERSION')[0]

e_bins = edisp['ENERG_LO'][1:]
migra_bins = edisp['MIGRA_LO'][1:]

plt.title('pyirf')
plt.pcolormesh(e_bins.to_value(u.TeV), migra_bins, edisp['MATRIX'].T[1:-1, 1:-1, 0].T, cmap='inferno')

plt.xscale('log')
plt.yscale('log')
plt.colorbar(label='PDF Value')

plt.xlabel(r'$E_\mathrm{True} / \mathrm{TeV}$')
plt.ylabel(r'$E_\mathrm{Reco} / E_\mathrm{True}$')

None # to remove clutter by mpl objects
```

--------------------------------

### Perform Spectral Weighting

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Reweights simulated events to match a target spectrum. Uses built-in reference spectra or custom power laws.

```python
import numpy as np
import astropy.units as u
from pyirf.spectral import (
    PowerLaw, LogParabola, calculate_event_weights,
    CRAB_HEGRA, CRAB_MAGIC_JHEAP2015,
    IRFDOC_PROTON_SPECTRUM, IRFDOC_ELECTRON_SPECTRUM,
)
from pyirf.simulations import SimulatedEventsInfo

# Built-in reference spectra
print(f"Crab HEGRA: {CRAB_HEGRA}")
print(f"Crab MAGIC: {CRAB_MAGIC_JHEAP2015}")

# Define custom power law spectrum
custom_spectrum = PowerLaw(
    normalization=3.0e-11 / (u.TeV * u.cm**2 * u.s),
    index=-2.5,
    e_ref=1 * u.TeV,
)

# Calculate flux at specific energy
energy = np.array([0.1, 1.0, 10.0]) * u.TeV
flux = custom_spectrum(energy)
print(f"Flux at 1 TeV: {flux[1]:.2e}")

# Create simulated spectrum from simulation info
simulation_info = SimulatedEventsInfo(
    n_showers=1_000_000,
    energy_min=0.01 * u.TeV,
    energy_max=100 * u.TeV,
    max_impact=1500 * u.m,
    spectral_index=-2.0,
    viewcone_min=0 * u.deg,
    viewcone_max=0 * u.deg,  # Point source
)
t_obs = 50 * u.hour
simulated_spectrum = PowerLaw.from_simulation(simulation_info, t_obs)

# Calculate event weights to reweight to Crab spectrum
true_energies = np.random.uniform(0.1, 10, 1000) * u.TeV
weights = calculate_event_weights(
    true_energies,
    target_spectrum=CRAB_HEGRA,
    simulated_spectrum=simulated_spectrum,
)
print(f"Mean weight: {weights.mean():.3f}")
```

--------------------------------

### Import Data Handling Libraries

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Imports libraries required for reading and handling tabular data with units.

```python
from astropy.table import QTable
import astropy.units as u
import tables
```

--------------------------------

### Plotting Flux Sensitivity

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Plots the flux sensitivity from EventDisplay and pyirf data. Requires matplotlib and astropy units.

```python
y, edges = irf_eventdisplay["DiffSens"].to_numpy()
yerr = irf_eventdisplay["DiffSens"].errors()
bins = 10**edges
x = bin_center(bins)
width = np.diff(bins)

fig, (ax_sens, ax_ratio) = plt.subplots(
    2, 1,
    gridspec_kw={'height_ratios': [4, 1]},
    sharex=True,
)

ax_sens.errorbar(
    x,
    y, 
    xerr=width/2,
    yerr=yerr,
    label="EventDisplay",
    ls=''
)

unit = u.Unit('erg cm-2 s-1')


e = sensitivity['reco_energy_center']
w = (sensitivity['reco_energy_high'] - sensitivity['reco_energy_low'])
s = (e**2 * sensitivity['flux_sensitivity'])


ax_sens.errorbar(
    e.to_value(u.TeV),
    s.to_value(unit),
    xerr=w.to_value(u.TeV) / 2,
    ls='',
    label='pyirf'
)

ax_ratio.errorbar(
    e.to_value(u.TeV), s.to_value(unit) / y,
    xerr=w.to_value(u.TeV)/2,
    ls=''
)
ax_ratio.set_yscale('log')
ax_ratio.set_xlabel("Reconstructed energy / TeV")
ax_ratio.set_ylabel('pyirf / eventdisplay')
ax_ratio.grid()
ax_ratio.yaxis.set_major_formatter(ScalarFormatter())

ax_ratio.set_ylim(0.5, 2.0)
ax_ratio.set_yticks([0.5, 2/3, 1, 3/2, 2])
ax_ratio.set_yticks([], minor=True)

# Style settings
ax_sens.set_title('Minimal Flux Satisfying Requirements for 50 hours')
ax_sens.set_xscale("log")
ax_sens.set_yscale("log")
ax_sens.set_ylabel(rf"$(E^2 \cdot \mathrm{{Flux Sensitivity}}) /$ ({unit.to_string('latex')})")

ax_sens.grid(which="both")
ax_sens.legend()
fig.tight_layout(h_pad=0)


None # to remove clutter by mpl objects
```

--------------------------------

### Plot Energy Resolution Comparison

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Compares energy resolution data from EventDisplay and PyIRF. Uses errorbar plots and sets a logarithmic x-axis scale. Requires data to be pre-processed into x, y, xerr, and yerr format.

```python
# Data from EventDisplay
y, edges = irf_eventdisplay["ERes"].to_numpy()
yerr = irf_eventdisplay["ERes"].errors()
x = bin_center(10**edges)
xerr = np.diff(10**edges) / 2


# Data from pyirf
bias_resolution = QTable.read(pyirf_file, hdu='ENERGY_BIAS_RESOLUTION')[1:-1]

# Plot function
plt.errorbar(x, y, xerr=xerr, yerr=yerr, ls='', label="EventDisplay")
plt.errorbar(
    0.5 * (bias_resolution['reco_energy_low'] + bias_resolution['reco_energy_high']).to_value(u.TeV),
    bias_resolution['resolution'],
    xerr=0.5 * (bias_resolution['reco_energy_high'] - bias_resolution['reco_energy_low']).to_value(u.TeV),
    ls='',
    label='pyirf'
)
plt.xscale('log')

# Style settings
plt.xlabel(r"$E_\mathrm{reco} / \mathrm{TeV}$")
plt.ylabel("Energy resolution")
plt.grid(which="both")
plt.legend(loc="best")

None # to remove clutter by mpl objects
```

--------------------------------

### Visualizing Point Spread Function (PSF) in 2D

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Generates a 2D polar plot of the radially symmetric Point Spread Function (PSF) for different energy bins. Requires astropy QTable and matplotlib.

```python
psf_table = QTable.read(pyirf_file, hdu='PSF')[0]
# select the only fov offset bin
psf = psf_table['RPSF'].T[:, 0, :].to_value(1 / u.sr)

offset_bins = np.append(psf_table['RAD_LO'], psf_table['RAD_HI'][-1])
phi_bins = np.linspace(0, 2 * np.pi, 100)



# Let's make a nice 2d representation of the radially symmetric PSF
r, phi = np.meshgrid(offset_bins.to_value(u.deg), phi_bins)

# look at a single energy bin
# repeat values for each phi bin
center = 0.5 * (psf_table['ENERG_LO'] + psf_table['ENERG_HI'])


fig = plt.figure(figsize=(15, 5))
axs = [fig.add_subplot(1, 3, i, projection='polar') for i in range(1, 4)]


for bin_id, ax in zip([10, 20, 30], axs):
    image = np.tile(psf[bin_id], (len(phi_bins) - 1, 1))
    
    ax.set_title(f'PSF @ {center[bin_id]:.2f} TeV')
    ax.pcolormesh(phi, r, image)
    ax.set_ylim(0, 0.25)
    ax.set_aspect(1)
    
fig.tight_layout()

None # to remove clutter by mpl objects
```

--------------------------------

### Plot Direction Cut Comparison

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Compares the theta cut (direction cut) between EventDisplay and pyirf as a function of reconstructed energy. Requires 'RAD_MAX' HDU from pyirf file and 'ThetaCut;1' from EventDisplay file.

```python
from astropy.table import QTable


rad_max = QTable.read(pyirf_file, hdu='RAD_MAX')[0]


theta_cut_ed, edges = irf_eventdisplay['ThetaCut;1'].to_numpy()



plt.errorbar(
    bin_center(10**edges),
    theta_cut_ed,
    xerr=np.diff(10**edges),
    ls='',
    label='EventDisplay',
)

plt.errorbar(
    0.5 * (rad_max['ENERG_LO'] + rad_max['ENERG_HI'])[1:-1].to_value(u.TeV),
    rad_max['RAD_MAX'].T[1:-1, 0].to_value(u.deg),
    xerr=0.5 * (rad_max['ENERG_HI'] - rad_max['ENERG_LO'])[1:-1].to_value(u.TeV),
    ls='',
    label='pyirf',
)

plt.legend()
plt.ylabel('θ-cut / deg²')
plt.xlabel(r'$E_\mathrm{reco} / \mathrm{TeV}$')
plt.xscale('log')

None # to remove clutter by mpl objects
```

--------------------------------

### Define pyirf Output File Path

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Specifies the path to the pyirf output FITS file, which contains IRF and sensitivity information.

```python
pyirf_file = '../../pyirf_eventdisplay.fits.gz'
```

--------------------------------

### Create Custom Gaussian Estimator Class

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/interpolation.rst

Defines a custom GaussianEstimator by inheriting from DiscretePDFComponentEstimator, adding unit handling for Gaussian parameters. Requires astropy.units.

```python
import astropy.units as u
from pyirf.interpolation import (DiscretePDFComponentEstimator,
                                    MomentMorphInterpolator)

class GaussianEstimatior(DiscretePDFComponentEstimator):
   @u.quantity_input(gaussians=u.m)
   def __init__(
      self,
      grid_points,
      bin_edges,
      gaussians,
      interpolator_cls=MomentMorphInterpolator,
      interpolator_kwargs=None,
      extrapolator_cls=None,
      extrapolator_kwargs=None,
   ):
      if interpolator_kwargs is None:
         interpolator_kwargs = {}

      if extrapolator_kwargs is None:
         extrapolator_kwargs = {}

      self.unit = gaussians.unit

      super().__init__(
         grid_points=grid_points,
         bin_edges=bin_edges,
         binned_pdf=gaussians.to_value(u.m),
         interpolator_cls=interpolator_cls,
         interpolator_kwargs=interpolator_kwargs,
         extrapolator_cls=extrapolator_cls,
         extrapolator_kwargs=extrapolator_kwargs,
      )

   def __call__(self, target_point):
      res = super().__call__(target_point)

      # Return result with correct unit
      return u.Quantity(res, u.m, copy=False).to(self.unit)

```

--------------------------------

### Plotting Point Spread Function (PSF) Profile

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Plots the profile of the Point Spread Function (PSF) for selected energy bins. Requires astropy QTable and matplotlib.

```python
# Profile
center = 0.5 * (offset_bins[1:] + offset_bins[:-1])
xerr = 0.5 * (offset_bins[1:] - offset_bins[:-1])

for bin_id in [10, 20, 30]:
    plt.errorbar(
        center.to_value(u.deg),
        psf[bin_id],
        xerr=xerr.to_value(u.deg),
        ls='',
        label=f'Energy Bin {bin_id}'
    )
    
#plt.yscale('log')
plt.legend()
plt.xlim(0, 0.25)
plt.ylabel('PSF PDF / sr⁻¹')
plt.xlabel('Distance from True Source / deg')

None # to remove clutter by mpl objects
```

--------------------------------

### Interpolate PSF Table

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Interpolate PSF tables between grid points using quantile-based interpolation.

```python
import numpy as np
import astropy.units as u
from pyirf.interpolation import PSFTableEstimator, QuantileInterpolator

# Grid points (e.g., zenith angles)
grid_points = np.array([[20], [40], [60]])
n_energy = 10
n_fov = 2
n_rad = 30

# PSF tables at each grid point (shape: n_points, n_energy, n_fov, n_rad)
psf_tables = np.random.uniform(0, 100, (3, n_energy, n_fov, n_rad)) / u.sr

# Source offset bins
source_offset_bins = np.linspace(0, 0.5, n_rad + 1) * u.deg

# Create PSF estimator
psf_estimator = PSFTableEstimator(
    grid_points=grid_points,
    source_offset_bins=source_offset_bins,
    psf=psf_tables,
    interpolator_cls=QuantileInterpolator,
)

# Interpolate to new observation condition
target = np.array([[35]])
psf_interp = psf_estimator(target)
print(f"Interpolated PSF shape: {psf_interp.shape}")  # (n_energy, n_fov, n_rad)
print(f"Unit: {psf_interp.unit}")
```

--------------------------------

### View Git Contributors

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/AUTHORS.rst

Use this bash command to view a sorted list of git contributors with their email addresses. This helps in identifying individuals who have contributed to the project.

```bash
git shortlog -sne
```

--------------------------------

### Initialize Plot for Sensitivity

Source: https://github.com/cta-observatory/pyirf/blob/main/examples/comparison_with_EventDisplay.ipynb

Initializes a matplotlib figure with a specified size (12x8 inches) for plotting sensitivity data.

```python
plt.figure(figsize=(12,8))
```

--------------------------------

### Gammapy Interoperability Functions

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/gammapy.rst

This section details the functions available for converting pyirf IRFs and binning information into gammapy compatible formats.

```APIDOC
## pyirf.gammapy Module

### Description
This module provides functions to convert the ``pyirf`` quantities for IRFs and the binning to the corresponding ``gammapy`` classes.

### Reference/API

.. automodapi:: pyirf.gammapy
   :no-inheritance-diagram:
```

--------------------------------

### Interpolate Energy Dispersion

Source: https://context7.com/cta-observatory/pyirf/llms.txt

Perform quantile-based interpolation on energy dispersion matrices.

```python
import numpy as np
from pyirf.interpolation import EnergyDispersionEstimator, QuantileInterpolator

# Grid of energy dispersions at different observation conditions
grid_points = np.array([[20], [40], [60]])  # e.g., zenith angles
n_energy = 10
n_migra = 50
n_fov = 1

# Energy dispersion matrices (shape: n_points, n_energy, n_migra, n_fov)
energy_dispersions = np.random.uniform(0, 1, (3, n_energy, n_migra, n_fov))
# Normalize to valid PDFs
energy_dispersions /= energy_dispersions.sum(axis=2, keepdims=True)

# Migration bins (reco/true energy ratio)
migra_bins = np.geomspace(0.2, 5, n_migra + 1)

# Create estimator
edisp_estimator = EnergyDispersionEstimator(
    grid_points=grid_points,
    migra_bins=migra_bins,
    energy_dispersion=energy_dispersions,
    interpolator_cls=QuantileInterpolator,
)

# Interpolate to new point
target = np.array([[30]])
edisp_interp = edisp_estimator(target)
print(f"Interpolated shape: {edisp_interp.shape}")  # (n_energy, n_migra, n_fov)
print(f"Sum check (should be ~1): {edisp_interp.sum(axis=1).mean():.3f}")
```

--------------------------------

### Download Data File

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Downloads the DL3 data file using curl, checking for existence and freshness. Raises an error if the download fails.

```python
path = "gamma_test_dl3.hdf5"
url = f"https://factdata.app.tu-dortmund.de/dl3/FACT-Tools/v1.1.2/{path}"
ret = sp.run(["curl", "-z", path, "-fsSLO", url], stdout=sp.PIPE, stderr=sp.PIPE, encoding='utf-8')
if ret.returncode != 0:
    raise IOError(ret.stderr)
```

--------------------------------

### Statistics Module API

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/statistics.rst

Reference documentation for the pyirf.statistics module, generated using automodapi.

```APIDOC
## Statistics Module API

### Description
This section provides the API reference for the `pyirf.statistics` module. It details the functions and classes available for statistical analysis within the PyIRF library.

### Module
`pyirf.statistics`

### Content
This documentation is auto-generated and reflects the structure and docstrings of the module's components. For detailed usage, please refer to the individual function or class docstrings within the library.
```

--------------------------------

### Simulation Information API

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/simulation.rst

Provides access to simulation-related functionalities within the PyIRF library.

```APIDOC
## Simulation Information API

### Description
This section covers the API endpoints and functionalities related to generating and managing simulations using the `pyirf.simulations` module.

### Method
N/A (This is a module-level documentation)

### Endpoint
N/A (This is a module-level documentation)

### Parameters
N/A

### Request Example
N/A

### Response
N/A
```

--------------------------------

### Convert PyIrf PSF to Gammapy 3D

Source: https://github.com/cta-observatory/pyirf/blob/main/docs/notebooks/fact_example.ipynb

Converts a PyIrf PSF table to a Gammapy `PSF3D` object using `pyirf.gammapy.create_psf_3d`. This function requires the PSF table, true energy bins, source offset bins, and FOV offset bins.

```python
from pyirf.gammapy import create_psf_3d

psf_gammapy = create_psf_3d(psf, true_energy_bins, source_offset_bins,  fov_offset_bins)

plt.figure()
psf_gammapy.plot_psf_vs_rad(offset=[wobble_offset], energy_true=[1., 10.]*u.TeV)
plt.legend(plt.gca().lines, ['1 TeV', '10 TeV'])
```