IceCube HESE 12-year data release

<p>This data release describes a revisited analysis of the IceCube high-energy starting events (HESE) sample with an additional 4.5 years of data that employs a full-resimulation reconstruction approach, <a href="https://arxiv.org/abs/1309.7010">DirectFit</a>, which exactly incorporates recent updates in glacial ice modeling. These updates include a microscopic description of the observed ice anisotropy due to <a href="http://doi.org/10.5194/tc-2022-174">birefringence</a> and a full description of ice layer undulations in three dimensions. Each event is reconstructed with the Approximate Bayesian Computation method, with the posterior sample points then fitted to an asymmetric directional probability density function.
</p>
<p>For further details, refer to the IceCube publication <a href="https://pos.sissa.it/444/1030/pdf">PoS-ICRC23-1030</a> / <a href="https://arxiv.org/abs/2307.13878">arXiv:2307.13878</a></p>

<h3>Note on interpretation</h3>
<p>
Reconstructed quantities for subsets of these events have been published previously. IceCube utilizes several event reconstructions, each with different strengths and employed in different contexts. DirectFit, the reconstruction in this data release, performs accurate event reconstruction, but is too computationally expensive to run on large scale Monte Carlo (MC). <b>Therefore this data release should not be employed in a forward-folding likelihood fit for MC-based spectral analyses</b>. Instead, as the reconstruction utilizes a complete description of the ice, the dataset is meant to provide accurate directional constraints for HESE events and for use in neutrino source searches. We highlight the differences between the reconstructions listed here, and those available elsewhere for subsets of this sample, below.
</p>

<p>
For overlapping events already present in the previous <a href="https://doi.org/10.1103/PhysRevD.104.022002"> 7.5 year dataset</a> (arxiv:2011.03545), the reconstructions provided here may be compared to those in Appendix J. The reconstructions in Appendix J were obtained with DirectFit, but using an older ice model. Differences primarily arise from the updates to ice optical modeling and/or changes in the reconstruction routine (shower vs track) applied for select events. As the new ice model leads to a generally improved data-MC agreement in both calibration and physics data, the resulting fits are as a result preferred. The per-event differences between the fits may allow to conservatively gauge the impact of systematic detector uncertainties.
</p>

<p>
For overlapping events already present in the previous <a href="https://doi.org/10.1103/PhysRevD.104.022002"> 7.5 year dataset</a> (arxiv:2011.03545), the reconstructions provided here utilize a different algorithm from those in Appendix F and listed as part of that <a href="https://doi.org/10.21234/4EQJ-BB17">data release</a>. Specifically, <b> the energies and directions published here should not be combined with the spectral fitting tools accompanying the 7.5 year dataset</b>. Diffuse spectral analysis rely on likelihood-based reconstructions (arXiv:1311.4767) that leverage tabulated light yield approximations. Such algorithms are substantially more computationally efficient than DirectFit and are fast enough to be applied to large scale Monte Carlo, as required for forward-folding fits, at a cost of reduced accuracy.
</p>

<p>
For overlapping events already present in the <a href="https://arxiv.org/abs/2304.01174">IceCat-1 catalog of alert tracks</a> and released in real time via <a href="https://gcn.nasa.gov/missions/icecube">GCN/AMON</a>, the visible energy reported here should not be compared to arxiv:2304.01174 (Tab. 3) and published in the accompanying <a href="https://doi.org/10.7910/DVN/SCRUCD">data release </a>. The IceCat-1 data release reports the most likely neutrino energy assuming a power-law spectral index of -2.19 at a given declination. The energies reported here correspond to the visible energy detected by IceCube.
</p>

<h3>Note on event times</h3>
The Modified Julian Date (MJD) is provided here for each event. In previously published HESE datasets that include the MJD, events that occurred between 2012-06-30 and 2015-05-17 were one second behind UTC time due to a master clock firmware issue. This offset has been corrected for in this data release; the overlapping events that fall within those dates will have a one second difference when compared to previous data releases. For those cases the MJDs provided here should be preferred.

<h3>Data release</h3>
Provided in this release is <code>data.npy</code>, which is a numpy structured array that contains the 164 reconstructed events that pass the HESE selection. The structured array contains the following fields:
<pre>
<code>
('id', '< i4')
('mjd', '< f8')
('ra', '< f4')
('dec', '< f4')
('params', [('f0', '< f8'), ('f1', '< f8'), ('f2', '< f8'), ('f3', '< f8'), ('f4', '< f8'), ('f5', '< f8'), ('f6', '< f8'), ('f7', '< f8')])
('reconstruction', 'S10')
('energy', '< f4')
('drlogl', '< f4')
</code>
</pre>
For each event, <code>id</code> corresponds to the event-ID consistent with prior HESE data releases, <code>mjd</code> is the Modified Julian Date when the event was observed, <code>ra</code> and <code>dec</code> correspond to the most-probable arrival direction, <code>params</code> corresponds to the eight parameters that describe the directional probability density function (PDF) using a generalized Fisher-Bingham distribution, <code>reconstruction</code> corresponds to the algorithm used by DirectFit, and <code>energy</code> the most-probable visible energy. For more details on the directional PDF construction see this <a href="https://doi.org/10.1007/s00180-020-01023-w">paper</a>. For <code>reconstruction</code> a Shower (Track) hypothesis indicates that the parameters of interest for the event were obtained with the corresponding algorithm in DirectFit. A few event reconstructions have changed from previous iterations of the analysis due to an updated classification scheme. The indicator for this new classifier is given by <code>drlogl</code>, which is the difference in the reduced log-likelihood between the track and shower hypothesis obtained from a fast, approximative model. More negative values in <code>drlogl</code> correspond to a higher degree of confidence that the event is track-like. In reference to the distinction between reconstruction routines, the terminology is the same as used historically.
</p>

<p>In addition, pre-generated full-sky probability density maps in FITS format are included and collected under <code>fits.[01].tar</code>. Each FITS file contains a event skymap discretized from the PDF described in <code>data.npy</code>. It is included to provide additional flexibility and ease of use at the expense of a slight information loss.
</p>

<h3>Example Python usage</h3>
It's easy to generate directional PDFs from <code>data.npy</code>. First get the dependencies with <code>pip install fb8 healpy</code>. Then a minimum working example to produce probability HEALPix maps with nside=256 is:
<pre>
<code>
import numpy as np
import healpy
from sphere.distribution import fb8
import logging
logger = logging.getLogger()
logger.setLevel(level=logging.ERROR)

nside=256
scan_theta, scan_ra = healpy.pix2ang(nside, np.r_[:healpy.nside2npix(nside)])
scan_dec = np.pi/2 - scan_theta

def fb8_to_pdf(f, i):
print('\rprocessing event {} ...'.format(i), end='')
return f.pdf(f.spherical_coordinates_to_nu(scan_theta, scan_ra))

arr = np.load('data.npy')
fb8s = [fb8(*_) for _ in arr['params']]

# healpix probability maps
maps = np.array([fb8_to_pdf(f, i) for (i, f) in enumerate(fb8s)])
</code>
</pre>

An example script to create FITS files and, optionally, plot skymaps for each event is included in <code>create_fits.py</code>. Plotting relies additionally on the matplotlib library, which can be installed with <code>pip install matplotlib</code>. For usage examples <code>python create_fits.py -h</code>
<h3>Contact</h3>
<p>[e-mail]</p>