© ESO, 2011

1. Introduction

Massive stars end their life exploding as core collapse supernovae (SN) and leaving compact remnants, such as neutron stars or black holes. Most of the binding energy is released by means of neutrino emission, though the details of the emission are not yet completely understood. The SN1987A event provided evidence that the neutrino emission exists and is detectable (Hirata et al. 1987; Bionta et al. 1987; Alekseev et al. 1988). The neutrino signal offers us the most promising opportunity to date to probe this unique astrophysical system.

Neutrino emission from core-collapse supernovae that have exploded in the Universe cumulate in the diffuse supernova neutrino background (DSNB). The prospects of detection seem promising; for reviews and references to the original papers, see Ando & Sato (2004) and Beacom (2010). To predict the DSNB signal, two different quantities are needed: the explosion rate of core collapse SNe as a function of the redshift and the average neutrino emission of the individual supernovae. Each of these quantities is known up to an error which implies uncertainty in the prediction. Until both uncertainties are quantified, it is impossible to evaluate the reliability of the expectations.

The rate of core-collapse supernovae can be identified with the rate of massive-star formation as a function of redshift. This quantity, in turn, can be obtained from knowledge of the cosmic history of star formation and the initial distribution of mass. Horiuchi et al. (2009) used a comprehensive compilation of data to evaluate the uncertainty in the rate.

We attempt to complement the study of the uncertainties by evaluating the impact of the neutrino emission for an individual supernova. Following Fukugita & Kawasaki (2003), we base our inferences on SN1987A observations using for the analysis the neutrino emission model of Pagliaroli et al. (2009). We evaluate the expected signal and its uncertainty for the Super-Kamiokande detector. We consider its present configuration, 22.5 kton of fiducial volume with a threshold of E_th = 19.3 MeV, Malek et al. (2003) and consider the possibility of lowering the threshold to 11.3 MeV by employing gadolinium loading as advocated by Beacom & Vagins (2004).

2. The emission from an individual supernova

We determine the spectrum from SN1987A data, which we then use for the analysis of DSNB. This is done by parameterizing the time and energy distribution, then by fitting the data from Kamiokande-II, IMB, and Baksan and finally by integrating over the time distribution. The main motivation behind this procedure is that we still do not have a definitive theory of supernova explosion; thus, despite their paucity, the data from SN1987A play a very important role in guiding our understanding of supernova neutrino manifestations.

The model for supernova emission.

Motivated by the prospect of exploring, through the neutrino signal, the physics and astrophysics of the gravitational collapse, we proposed in Pagliaroli et al. (2009) a parameterization of emission, based on our present understanding of emission processes improving the model by Loredo & Lamb (2001). Our model has two emission phases and the flux is ${\mathrm{\Phi }}_{\mathit{\nu ̅}\mathrm{e}}\mathrm{\left(}\mathit{t,}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{=}{\mathrm{\Phi }}_{\mathrm{a}}\mathrm{\left(}\mathit{t,}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{+}\mathrm{[}\mathrm{1}\mathrm{-}{\mathit{j}}_{\mathit{k}}\mathrm{(}\mathit{t}\mathrm{\left)}\mathrm{\right]}{\mathrm{\Phi }}_{\mathrm{c}}\mathrm{(}\mathit{t}\mathrm{-}{\mathit{\tau }}_{\mathrm{a}}\mathit{,}{\mathit{E}}^{\mathrm{\prime }}\mathrm{)}\mathit{,}$(1)where t is the emission time and E′ is the emitted neutrino energy. The first term Φ_a is the flux generated during the phase of accretion and above the shock by the interactions between the positrons and the target neutrons. It describes volume emission. It uses three free parameters: the initial accreting mass M_a, the time scale τ_a, and the initial temperature of the positrons T_a. The second term Φ_c is the flux coming from the thermal emission of the new-born proto-neutron star (i.e., the cooling phase). It describes radiation from a surface, whose area is proportional to the radius of the neutrino sphere R_c, with a time scale τ_c and the initial temperature of the emitted antineutrinos T_c. Finally, the function j_k(t) links smoothly the two emission phases, delaying the cooling emission by τ_a. Analytic expressions of these three functions, Φ_a, Φ_c, and j_k(t) are given in Eqs. (10), (13) and (18) of Pagliaroli et al. (2009), respectively. We include neutrino oscillations with normal hierarchy as discussed in Sect. C of Pagliaroli et al. (2009).

SN1987A data analysis.

In Pagliaroli et al. (2009), we tested the parameterized model using the small set of events collected in 1987 by Kamiokande-II, IMB, and Baksan, leaving aside an interpretation of the events recorded by LSD (Aglietta et al. 1987). When we compared the fit based only on the cooling phase, the one adopted in usual SN1987A data analyses, with our fit, we found that the two-phase emission model is 50 times more probable (i.e., we got a 2.5σ indication in its favor). Moreover, the best-fit values of the astrophysical parameters, namely R_c = 16km, M_a = 0.22   M_⊙, T_c = 4.6MeV, T_a = 2.4MeV, τ_c = 4.7s, and τ_a = 0.55s, agree well with the general expectations. For instance, the duration of the accretion phase is shorter than one second, the radius of the neutrino sphere is similar in size to the neutron star, and the total radiated energy 2.2 × 10⁵³ erg is similar to the binding energy. Finally, the luminosity curve and the mean energy as functions of time both resemble the results of numerical simulations. The errors in the parameters and their correlations are large, as expected from the limited number of events (see again Pagliaroli et al. 2009).

The previous analyses of SN1987A data of Lunardini (2006) and Yuksel & Beacom (2007) designed to predict DSNB focused only on the energy spectrum, motivated by the opinion that the time distribution of the events is irrelevant. However, we show that our detailed theoretical description of neutrino emission leads to a peculiar integrated spectrum that retains an imprint of the two different emission phases.

Fig. 1 The SN1987A best-fit spectrum (continuous line) and the two contributions from accretion phase (dotted line) and cooling phase (dashed line). In the inset, we plot the percentage difference with the Fermi-Dirac approximation of the spectrum, given in the text.

The emission spectrum.

By substituting in the flux of Eq. (1) the best-fit values reported in the previous paragraph and integrating both in a 30 s time window and on the emitting surface, we obtain the reference spectrum dN/dE′ for a single supernova event from SN1987A shown in Fig. 1. We compare this spectrum with a Fermi-Dirac distribution, i.e., $\frac{\mathit{k}}{{\mathit{T}}^{\mathrm{3}}}\frac{{\mathit{E}}^{\mathrm{\prime }\mathrm{2}}}{\mathrm{1}\mathrm{+}\exp \mathrm{(}{\mathit{E}}^{\mathrm{\prime }}\mathit{/}\mathit{T}\mathrm{-}\mathit{\eta }\mathrm{)}}$. To obtain the same total integral and the same first two momenta between our spectrum and the Fermi-Dirac distribution, we need k = 2.97 × 10⁵⁷, a temperature T = 3.76   MeV, and a pinching factor η = 0.679. However, the percentage difference between the Fermi-Dirac spectrum and our spectrum is not negligible as shown in the inset of Fig. 1. The agreement worsens using the parametrization of Keil et al. (2003). The emitted electron antineutrinos of our best-fit spectrum are of relatively low energy; half of them are emitted at energies below 11.2 MeV. In particular, when we calculate the number of events N_ev expected from a galactic supernova assuming an energy threshold of E_min > 6.5 MeV, as a function of the upper extreme of integration E_max, we find that there are the 25%, 50%, and 75% of the events below E_max = 14, 18, and 24 MeV, respectively¹.

3. Expectations for the diffuse neutrino flux

We estimate the diffuse neutrino flux accumulated by all the past supernovae exploded in the Universe, assuming that the antineutrino flux discussed in the previous Section represents the typical emission of a core-collapse event.

In a Friedmann-Robertson-Walker flat Universe, the expected DSNB flux is $\frac{\mathrm{d}\mathit{\phi }\mathrm{\left(}\mathit{E}\mathrm{\right)}}{\mathrm{d}\mathit{E}}\mathrm{=}\frac{\mathit{c}}{{\mathit{H}}_{\mathrm{0}}}{\mathrm{\int }}_{\mathrm{0}}\mathrm{d}\mathit{z}\frac{{\mathit{R}}_{\mathrm{CCSN}}\mathrm{\left(}\mathit{z}\mathrm{\right)}}{\sqrt{{\mathrm{\Omega }}_{\mathrm{\Lambda }}\mathrm{+}{\mathrm{\Omega }}_{\mathrm{m}}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}{\mathrm{)}}^{\mathrm{3}}}}\frac{\mathrm{d}\mathit{N}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}}{\mathrm{d}{\mathit{E}}^{\mathrm{\prime }}}\mathrm{·}$(2)The last term of the integrand, dN/dE′, is the spectrum of a single SN emission discussed above and calculated for the red-shifted energy E′ = (1 + z)E. The Hubble constant is H₀ = 71km   s^-1   Mpc^-1, c is the light speed, and the values of the cosmological constants are Ω_m = 0.27, Ω_Λ = 0.73 as measured from the WMAP experiment in Jarosik et al. (2011). These assumptions do not introduce significant errors into the predictions.

The key term we have to discuss is R_CCSN(z), namely the rate of core-collapse SN for comoving volume as a function of the redshift z. This function can be obtained from the number of stars formed in a comoving volume, i.e., from the star formation rate as a function of the redshift, and from the fraction of stars in the appropriate mass range for core collapse SNe. Following Horiuchi et al. (2009), we adopt the initial mass function of Salpeter (1955) and three different analytic fits for the rate of core-collapse to take into account the astrophysical uncertainty. In Fig. 2, we show the quantity $\mathit{n}\mathrm{\left(}\mathit{z}\mathrm{\right)}\mathrm{=}\frac{{\mathit{R}}_{\mathrm{CCSN}}\mathrm{\left(}\mathit{z}\mathrm{\right)}}{\sqrt{{\mathrm{\Omega }}_{\mathrm{\Lambda }}\mathrm{+}{\mathrm{\Omega }}_{\mathrm{m}}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}{\mathrm{)}}^{\mathrm{3}}}}\mathit{,}$(3)which we call the “cosmological factor”, namely the first integrand term in Eq. (2). This plot summarizes all the relevant assumptions about the cosmology and the distribution of the sources. Horiuchi et al. (2009) demonstrated that the use of either the Kroupa (2001) or Baldry & Glazebrook (2003) initial mass functions, or the uncertainties in the lowest mass that forms a supernova, introduces only a negligible error in the estimates.

We now use Eq. (2) to calculate the flux of the diffuse SN neutrino background dφ/dE as a function of the neutrino energy at the Earth. We plot it in Fig. 3 along with the upper bound in the energy region E > E_th = 19.3 MeV of  cm^-2 s^-1 obtained by the Super-Kamiokande collaboration (Malek et al. 2003). The total flux of DSNB with the SN1987A best-fit model is φ = 27.2cm^-2s^-1 for the fiducial rate. Only  ~1% of this flux can be actually observed with a neutrino energy threshold of 19.3 MeV; with the use of gadolinium, the accessible spectrum becomes  ~8%, showing that most of the flux falls in the very low energy region.

Event rate in Super-Kamiokande.

Fig. 2 The cosmological factor n of Eq. (3) for various fits of the core-collapse rate. The continuous line is for the fiducial rate, the dotted line is for the upper one, and the dashed line is for the lower rate.

Fig. 3 Expected DSNB flux in a logarithmic scale for the three different fits of core-collapse rate. Also is shown the upper limit of this quantity given by Super-Kamiokande in 2003 in the energy range between 19.3 and 27.3 MeV. The meaning of continuous, dotted and dashed lines as in Fig. 2.

We consider the events expected in a detector such as Super-Kamiokande, with a fiducial mass of M_d = 22.5 kton of water and a detection efficiency ϵ set equal to the 98% above a neutrino energy threshold (E_th = 19.3 MeV in the present configuration or E_th = 11.3 MeV by loading the detector with gadolium). However, we note that a scintillator-based detector with average chemical formula C₉H₂₁ (or C₆H₃(CH₃)₃) and with 16 (or 25) kton of mass has the same number of target protons N_p = 1.5 × 10³³, thus the same number of interactions. The event rate is calculated easily as $\mathit{Ṅ}\mathrm{ev}\mathrm{=}{\mathit{N}}_{\mathrm{p}}{\mathrm{\int }}_{{\mathit{E}}_{\mathrm{th}}}\mathrm{d}\mathit{E\sigma }\mathrm{\left(}\mathit{E}\mathrm{\right)}\frac{\mathrm{d}\mathit{\phi }}{\mathrm{d}\mathit{E}}\mathit{ϵ}\mathrm{\left(}\mathit{E}\mathrm{\right)}\mathit{.}$(4)For the best-fit model adopted, we expect to see in Super-Kamiokande the number of events per year reported in Table 1 for two different energy thresholds and three different fits of the rate. We expect 0.39–0.65 events per year when we consider a threshold of 19.3 MeV; this increases to 1.35–2.35 events per year lowering the energy threshold to 11.3 MeV. Our results agree well with those of Horiuchi et al. (2009).

4. Uncertainty due to the emission model

In Sect. 3, we discussed the expectations for the DSNB flux based on the SN1987A best-fit model and we considered only the uncertainty related to the core-collapse rate. However, the SN1987A data set is small and the procedure of using only the best-fit model is doubtful, until one estimates the associated theoretical error. The main aim of this work is, indeed, to quantify this uncertainty, which is an entirely new result.

We proceed as follows. The likelihood in the SN1987A data analysis used by Pagliaroli et al. (2009) is a function of six unknown astrophysical parameters. The same is true for the DSNB rate Ṅ_ev. For a given energy threshold, this can be written explicitly as $\mathit{Ṅ}{}_{\mathrm{ev}}\mathrm{=}{\mathit{M}}_{\mathrm{a}}\mathit{f}\mathrm{\left(}{\mathit{T}}_{\mathrm{a}}\mathit{,}{\mathit{\tau }}_{\mathrm{a}}\mathrm{\right)}\mathrm{+}{\mathit{R}}_{\mathrm{c}}^{\mathrm{2}}\mathit{g}\mathrm{\left(}{\mathit{T}}_{\mathrm{c}}\mathit{,}{\mathit{\tau }}_{\mathrm{c}}\mathit{,}{\mathit{\tau }}_{\mathrm{a}}\mathrm{\right)}$, where the two contributions correspond to the two terms of Eq. (1) and f and g are known functions. The idea is simply to rewrite the likelihood, substituting the neutrino-sphere radius R_c with the quantity in which we are interested, Ṅ_ev. We can then apply a standard marginalization procedure: for each fixed value of Ṅ_ev and varying the other five parameters, we find the maximum value of the likelihood L. Thus, we obtain Δχ² = 2(lnL_max − lnL) a function of Ṅ_ev only.

Fig. 4 Δχ2 curves as a function of event rate Ṅev in Super-Kamiokande, as obtained from SN1987A data analysis and for two different energy thresholds. The set of thin lines correspond to the present threshold of 19.3 MeV whereas the set with thick lines correspond to Eth = 11.3 MeV. The meaning of the individual lines is as in the previous figures.

Fig. 5 Cumulative distribution of the events rate in Super-Kamiokande, with the same notation as Fig. 2 for the lines. We give the 2σ error bars on that prediction when Eth = 11.3 MeV and  = 19.3 MeV.

The results of the numerical evaluations of this Δχ² are shown in Fig. 4. The thin lines represent the neutrino energy threshold used in Malek et al. (2003), E_th = 19.3   MeV, whereas the thick ones correspond to the threshold of E_th = 11.3   MeV that can be achieved with gadolinium. Three lines of each type are given, for the three different fits of the cosmic rate for core-collapse events. The 1σ and 2σ errors found from the Δχ² curves are reported in Table 1.

We note that this procedure properly takes into account the correlations among the parameters, that are dictated by the likelihood function, without the need to explicitly calculate the covariance matrix.

5. Summary and discussion

We have evaluated the signal expected for diffuse supernova neutrino background in Super-Kamiokande and its associated uncertainties, using SN1987A data to constrain the model of neutrino emission. The results are summarized in Fig. 5, where we show the cumulative rate of DSNB events in Super-Kamiokande as a function of the energy threshold for the three different fits of the core-collapse rate.

Following Horiuchi et al. (2009), we used the three different fits as a “generous” assessment of the cosmological uncertainty, although we are unable to attach to this range a precise statistical meaning. The result is given in Table 1: the cosmological uncertainty is  ~25% for the higher threshold and  ~27% for the lower threshold, in reasonable agreement with Horiuchi, despite the different emission model used.

On the other hand, selecting the fiducial cosmological model for the rate, we evaluated the uncertainty on the DSNB due to the emission model. The 1σ (resp. 2σ) percentage error is  ~33% (resp. 65%) for the higher threshold and becomes  ~23% (resp. 46%) for the lower threshold. The prediction is more precise in the second case since the SN1987A data provide a tighter constraint of the low energy region of the spectrum. On the basis of the results of Pagliaroli et al. (2009), we conservatively estimate a 10% systematic uncertainty in our spectrum, which is intrinsic to the parameterization used in the fit. Our 2σ ranges are curiously similar to those obtained by Ando & Sato (2004), who use three different emission models obtained from numerical simulations.

It is difficult to combine the two errors in a robust way, since only the one related to the emission model has a precise statistical meaning. We thus construct a global range for the expected rate of DSNB events by the following conservative procedure. The upper value of the range is obtained by summing the best-fit value for the upper core collapse SN rate and its upward 1σ (resp., 2σ) statistical error; similarly, for the lower value of the range. In particular for E_th = 19.3 MeV, the expected event rate in Super-Kamiokande ranges between 0.27 (resp., 0.19) events per year and 0.88 (resp., 1.16) events per year, giving a global uncertainty of  ~53% (resp., 95%). For the lower energy threshold E_th = 11.3 MeV, the range becomes 1.06–2.94 (resp., 0.83–3.59) events per year with a percentage total error of 47% (resp., 76%). As quantified by the percentage errors, most of the global uncertainty is due to the emission model uncertainty and to the low number of SN1987A data.

References

Appendix A: An alternative expression for the rate

We introduce an alternative expression for Ṅ_ev that emphasizes the role of the supernova flux at the source providing additional insigh into the results. From the emission spectrum dN/dE′, we define an effective flux $\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{\ast }}}{\mathrm{d}\mathit{E}}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{\equiv }\frac{\mathrm{1}}{\mathrm{4}\mathit{\pi }{\mathit{d}}^{\mathrm{2}}\mathit{T}}\frac{\mathrm{d}\mathit{N}}{\mathrm{d}{\mathit{E}}^{\mathrm{\prime }}}\mathrm{with}\{\begin{array}{c}\\ \mathit{d}\mathrm{\equiv }\sqrt{\frac{\mathit{c}}{\mathrm{4}\mathit{\pi }{\mathit{n}}_{\mathrm{\ast }}{\mathit{H}}_{\mathrm{0}}\mathit{T}}}\mathrm{=}\mathrm{310}\mathrm{kpc}\\ \mathit{T}\mathrm{\equiv }\mathrm{1}\mathrm{yr}\end{array}$(A.1)where we have introduced the typical value for the cosmic density of n _∗  = 2 × 10^-4 SN/(Mpc³ yr) and the observational time T. By using as an integration variable E′ = E(1 + z) (i.e., the neutrino energy at the emission), we rewrite the signal rate in Eq. (4) mimicking closely the signal from a galactic supernova $\mathit{Ṅ}\mathrm{ev}\mathrm{=}{\mathit{N}}_{\mathrm{p}}{\mathrm{\int }}_{{\mathit{E}}_{\mathrm{th}}}^{\mathrm{\infty }}\mathrm{d}{\mathit{E}}^{\mathrm{\prime }}\mathit{\sigma }\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{\ast }}}{\mathrm{d}{\mathit{E}}^{\mathrm{\prime }}}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}{\mathit{ϵ}}_{\mathrm{\ast }}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\mathit{.}$(A.2)The distribution of cosmic supernovae is contained in the function ${\mathit{ϵ}}_{\mathrm{\ast }}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{\equiv }{\mathrm{\int }}_{\mathrm{0}}^{\frac{{\mathit{E}}^{\mathrm{\prime }}}{{\mathit{E}}_{\mathrm{th}}}\mathrm{-}\mathrm{1}}\frac{\mathrm{d}\mathit{z}}{\mathrm{1}\mathrm{+}\mathit{z}}\frac{\mathit{n}\mathrm{\left(}\mathit{z}\mathrm{\right)}\mathit{\sigma }\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}\mathrm{)}\mathrm{\right)}}{{\mathit{n}}_{\mathrm{\ast }}\mathit{\sigma }\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}}\mathit{ϵ}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}\mathrm{)}\mathrm{\right)}\mathit{.}$(A.3)This function leads to a severe cut in the rate of events given by Eq. (A.2). Above the energy threshold, this function can be approximated linearly as ${\mathit{ϵ}}_{\mathrm{\ast }}\mathrm{\left(}{\mathit{E}}^{\mathrm{\prime }}\mathrm{\right)}\mathrm{\propto }{\mathit{E}}^{\mathrm{\prime }}\mathrm{-}{\mathit{E}}_{\mathrm{th}}\mathit{,}$(A.4)which corresponds to DSNB selecting the highest energy tail of the spectrum. This is important in connection with SN1987A. We have only limited information on the highest energy tail of the spectrum, mostly thanks to IMB, as noted by Fukugita et al. (2003). We can constrain far more reliably the spectrum at lower energies, thanks to data from Kamiokande-II and Baksan. Thus, we expect, by lowering the threshold, not only the number of events due to DSNB to increase, but also the accuracy of the predictions based in SN1987A to improve, in agreement with Fig. 5. We also note that the numerical value of d in Eq. (A.1) agrees with the DSNB signal being modest.

All Tables

All Figures

	Fig. 1 The SN1987A best-fit spectrum (continuous line) and the two contributions from accretion phase (dotted line) and cooling phase (dashed line). In the inset, we plot the percentage difference with the Fermi-Dirac approximation of the spectrum, given in the text.
In the text

	Fig. 2 The cosmological factor n of Eq. (3) for various fits of the core-collapse rate. The continuous line is for the fiducial rate, the dotted line is for the upper one, and the dashed line is for the lower rate.
In the text

	Fig. 3 Expected DSNB flux in a logarithmic scale for the three different fits of core-collapse rate. Also is shown the upper limit of this quantity given by Super-Kamiokande in 2003 in the energy range between 19.3 and 27.3 MeV. The meaning of continuous, dotted and dashed lines as in Fig. 2.
In the text

	Fig. 4 Δχ2 curves as a function of event rate Ṅev in Super-Kamiokande, as obtained from SN1987A data analysis and for two different energy thresholds. The set of thin lines correspond to the present threshold of 19.3 MeV whereas the set with thick lines correspond to Eth = 11.3 MeV. The meaning of the individual lines is as in the previous figures.
In the text

	Fig. 5 Cumulative distribution of the events rate in Super-Kamiokande, with the same notation as Fig. 2 for the lines. We give the 2σ error bars on that prediction when Eth = 11.3 MeV and  = 19.3 MeV.
In the text