On possibility of determining neutrino mass hierarchy by the charged-current and neutral-current events of supernova neutrinos in scintillation detectors

One of the unresolved mysteries in neutrino physics is the neutrino mass hierarchy. We present a new method to determine neutrino mass hierarchy by comparing the events of inverse beta decays (IBD), ν¯e+p→n+e+\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\bar{\nu }}_e + p\rightarrow n + e^+$$\end{document}, and neutral current (NC) interactions, ν(ν¯)+p→ν(ν¯)+p\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\nu ({\overline{\nu }}) + p\rightarrow \nu ({\overline{\nu }}) + p$$\end{document}, of supernova neutrinos from accretion and cooling phases in scintillation detectors. Supernova neutrino flavor conversions depend on the neutrino mass hierarchy. On account of Mikheyev–Smirnov–Wolfenstein effects, the full swap of the ν¯e\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\bar{\nu }}_e$$\end{document} flux with the ν¯x\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\bar{\nu }}_x$$\end{document} (x=μ,τ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$x=\mu ,~\tau $$\end{document}) one occurs in the inverted hierarchy, while such a swap does not occur in the normal hierarchy. In consequence, the ratio of high energy IBD events to NC events for the inverted hierarchy is higher than in the normal hierarchy. Since the luminosity of ν¯e\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\bar{\nu }}_e$$\end{document} is larger than that of νx\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\nu _x$$\end{document} in accretion phase while the luminosity of ν¯e\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\bar{\nu }}_e$$\end{document} becomes smaller than that of νx\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\nu _x$$\end{document} in cooling phase, we calculate this ratio for both accretion and cooling phases. By analyzing the change of this event ratio from accretion phase to cooling phase, one can determine the neutrino mass hierarchy.


Introduction
Supernovae (SNe) are among the most powerful sources of neutrinos in our Universe. During a supernova explosion, 99% of the emitted energy (∼ 10 53 erg) is released by neutrinos and antineutrinos of all favors, with energy ranging from several to a few tens MeV, which play the role of astrophysical messengers, escaping almost unimpeded from the supernova core. The supernova neutrino flux has been extensively studied as a probe of both fundamental neutrino properties and core-collapse physics. Therefore, a e-mail: kcl@mail.cgu.edu.tw observing the supernova neutrino signal would enable a wide range of opportunities, both in astrophysics and in particle physics.
Based on various oscillation experiments with atmospheric, solar, and terrestrial neutrinos [1], a considerable progress has been achieved in constraining the neutrino mixing parameters [2][3][4]. The flavor states ν e , ν μ , and ν τ are now well recognized to be superpositions of the vacuum mass eigenstates ν 1 , ν 2 , and ν 3 [5]. Thanks to many successful experiments, the three neutrino flavor mixing angles, θ 12 , θ 23 , and θ 13 , and two mass-squared differences, 2 21 = m 2 2 − m 2 1 and 2 31 = m 2 3 − m 2 1 are well constrained, whereas the sign of 2 31 , i.e., the neutrino mass hierarchy, is still unknown. To determine the neutrino mass hierarchy, recent efforts include works based on reactor neutrinos [6][7][8][9], different baseline experiments [10], Earth matter effects on supernova neutrino signals [11,12], spectral swap of SN neutrino flavors [13], rise time of SN ν e light curve [14], ν e andν e light curves on the early accretion phase [15], analysis of meteoritic SN material [16], and detection of atmospheric neutrinos in sea water or ice [17]. Among them, those works using supernova neutrinos are particularly interesting because of the interplay between intrinsic properties of massive neutrinos and the mechanism of SN explosions.
Historically, the detection of neutrinos from SN1987A [18,19] has motivated a huge amount of theoretical works in both SN physics and neutrino physics. Therefore, many proposals to identify neutrino mass hierarchy by studying neutrinos from galactic SNe have been proposed. Originating from deep inside the SN core, neutrinos can experience significant flavor transitions on their way to the terrestrial detectors. Neutrino flavor conversions arising from the Mikheyev-Smirnov-Wolfenstein (MSW) effect [20][21][22] are sensitive to neutrino mass hierarchy. Additionally, it has been pointed out that the collective neutrino oscillation [23][24][25][26][27][28][29][30][31][32][33] (see [34] for a review) results from the coherent ν − ν forward scatterings in the deep region of the core where neutrino densities are large and may lead to collective pair flavor conversion ν eνe ↔ ν xνx (x = μ, τ ) over the entire energy range. However, unlike the status of MSW effects, consensus on collective flavor transitions has not yet been reached. To avoid digression to diverse scenarios of the collective effect, we assume that MSW effect dominates the flavor conversions when SN neutrinos propagate outwards.
The interactions of SN neutrinos with atomic nuclei and free protons are utilized to resolve the neutrino mass hierarchy in most of the methods. The inverse beta decay (IBD), ν e + p → n + e + , is the major interaction channel for neutrino detection in the water Cherenkov and liquid scintillation detectors. On the other hand, the liquid argon detector has a good sensitivity to ν e via charged-current interactions. Because the threshold of visible energy in a liquid scintillation detector can be as low as 0.2 MeV [35] by controlling the abundance of 14 C, the neutral-current (NC) interactions, ν + p → ν + p, will give rise to a large number of events in a channel other than IBD in this case and become very important. As a result, the detection of other species of SN neutrinos was proposed by measuring NC interactions [36,37].
Inspired by the capability of detecting thousands of neutrino events from a galactic supernova with next-generation scintillation detectors, we proposed to identify the neutrino mass hierarchy by comparing IBD and NC interactions inside the scintillators [38]. In the previous work, a set of specific mean energies of different flavors and luminosity equipartition between flavors for SN neutrinos are adopted. Instead of a unique scenario, models with different sets of mean energies and partitions of luminosities are explored in this work and the evolution with time of SN neutrinos is also accounted for in calculating IBD and NC events in the scintillation detector. We study how the way that the IBD and NC events change with time during a SN explosion is related to the neutrino mass hierarchy.
The paper is organized as follows. In Sect. 2, we briefly review the flavor transitions of SN neutrinos as they propagate outward from deep inside a SN and traverse the Earth medium to reach the detector and describe the supernova neutrino fluence in our calculation. In Sect. 3, we describe the event calculation for inverse beta decay and neutral current inside liquid scintillation detectors and briefly discuss interactions of SN neutrinos in the detectors. Then, in Sect. 4, we define a ratio R of total IBD events to total NC events and present our calculations for considered parameter space with statistical uncertainties addressed. Finally, in Sect. 5 we summarize our results and conclude.

Primary neutrino fluence
A SN neutrino burst lasts for t ≈ 10s and includes all six flavors of neutrinos. The total gravitational binding energy released in the explosion is E ≈ 10 53 erg. The neutrino flavors ν μ , ν τ and their antiparticles have similar interactions and thus similar average energies and fluences. Thus the total energy is divided as E = E ν e + Eν e + 4E ν x . In this work, the condition of equipartition of energies and luminosities among the primary neutrino flavors, E ν e ≈ Eν e ≈ E ν x and L ν e ≈ Lν e ≈ L ν x , is relaxed. The primary SN neutrino energy spectrum is typically not purely thermal. We adopt the Keil parametrization [39] for the neutrino fluence where α = E α / < E α > is the time-integrated flux, < E α > is the average neutrino energy, and η α denotes the pinching of the spectrum. In our calculation, we take η α = 3 for all flavors. If flavor conversions do not occur during the propagations of neutrinos from the SN core to the Earth, a SN at a distance d thus yields a neutrino fluence with E α in units of 10 52 erg, d in 10 kpc, and energies in MeV. For the numerical evaluations, we take a representative supernova at the Galactic center region with d = 10 kpc, and a total energy output of E = 3 × 10 53 erg. Further, different sets of the average energies, (< E ν e >, < Eν e >, < E ν x >), will be taken in our calculation. In addition to < E ν e >= 12 MeV, < Eν e >= 15 MeV, < E ν x >= 18 MeV, more hierarchical values, < E ν e >= 10 MeV, < Eν e >= 15 MeV, < E ν x >= 24 MeV, and more degenerate values, < E ν e >= 12 MeV, < Eν e >= 14 MeV, < E ν x >= 16 MeV, will also be considered.

Neutrino fluence on earth
As neutrinos propagate outwards from deep inside a SN and finally reaches the Earth, their flavor contents are modified by the MSW effect. The fluxes of ν e andν e arriving at the detector can be written as: for the normal hierarchy, and for the inverted hierarchy [40]. Here P 2e (P 2e ) is the probability that a mass eigenstate ν 2 (ν 2 ) is observed as a ν e (ν e ) since neutrinos arrive at the Earth as mass eigenstates. We do not consider the regeneration factor due to the Earth matter effect and thus take P 2e = sin 2 θ 12 in this work. From Eqs.
(3)- (6), it is shown that, in the normal hierarchy, ν e completely comes from ν 0 x from the source whileν e comes from bothν 0 e andν 0 x . On the other hand, in the inverted hierarchy, ν e comes from both ν 0 e and ν 0 x whileν e completely comes fromν 0 x . For the rest of flavors, the condition of flux conservation gives and for the normal and inverted hierarchies, respectively.

Events of inverse beta decay and neutral current interaction inside scintillation detectors
In scintillation detectors, inverse beta decays (IBD) are the most dominant interactions. IBD events are obtained in scintillation detectors by measuring the positron energy deposit. The observed event spectrum and total number of IBD events are given by where N p is the number of the target protons in the detector and cross section σ IBD (E ν ) is taken from [41]. The minimum neutrino energy for generating IBD interaction is E min = 1.8 MeV.
Inside the scintillation detector, the yield of νp elastic scatterings is also comparable to that of IBD due to the large number of free protons [37]. The observed event spectrum is given as where F tot ≡ F e + Fē + 4F x is the total fluence of the SN neutrinos and T is the recoil kinetic energy of protons which are scattered by SN neutrinos. To produce a proton recoil energy T requires a minimum neutrino energy E ν,min = m p T /2, with m p the proton mass. In other words, a neutrino of energy E ν can produce a proton recoil energy between 0 and T max = 2E 2 ν /m p . These protons are slow, and hence they are detected with quenched energies T < T . The proton recoil energy T is mapped to an electron-equivalent quenched energy T through the quenching function where k B is Birks constant [42]. The number of NC events is then given by where the differential cross section, dσ νp /dT , is taken from [36,37].We point out that not all signals within the energy range of proton recoils are taken into account. Since the scintillator is made of hydrocarbon, a natural isotope of the carbon, 14 C, decays into 14 N, emitting electrons below 0.2 MeV with a high rate. Below this energy, the signal is flooded by very low energy electrons. Therefore, a threshold of T min = 0.2 MeV is set for recording the signal. The threshold of T min = 0.2 MeV is converted to the threshold of proton recoil energy T min , e.g. T min = 0.93 MeV for JUNO detector.
In our previous work, the ratio used to probe neutrino mass hierarchy is the ratio of the total interactions of NC to those of IBD, which require reconstruction of the entire spectrum of SN neutrinos from detected events, d N/dT and d N/d E e + . By exploring SN neutrino physics with a fixed set of parameters, we have shown that the capability of using SN neutrinos to probe neutrino parameters. In this work we would like to expand the parameter space of SN neutrinos by releasing the energy-equipartition condition and taking into account more combinations of mean energies of different flavors. we would also construct more realistic observables directly and explicitly related with detected events, d N/dT and d N/d E e + , instead of the ratio of interactions.

Resolving neutrino mass hierarchy
Energy equipartition between all flavors is commonly assumed during the entire SN neutrino burst. However, neutrino emissions actually evolve with time as the SN explodes. An important feature of the evolution of SN neutrino emissions is that the hierarchy of luminosities in accretion and cooling phases are reversed. For our present understanding [31,45], L ν e ≈ Lν e > L ν x during the accretion phase and L ν e ≈ Lν e < L ν x during the cooling phase. Therefore, while L ν e ≈ Lν e , L ν x /L ν e grows as SN neutrino emissions evolve from the accretion phase into the cooling phase. In addition to scan over plausible ranges for luminosity ratios, we check three specific scenarios corresponding to the accretion phase, the equipartition model, and the cooling phase, respectively, as in Table 1.
The number of neutrino events inside the scintillator depends on < E α >'s and L α 's, the SN neutrino parameters. In stead of assuming a standard SN neutrino emission model as in [38], we investigate the ratio of IBD events to NC events for different SN neutrino emission models by varying the SN neutrino parameters, < E α >'s and L α 's, in calculating neutrino events inside the scintillator. We define R to be the ratio of the total IBD events to the total NC events, For a specific set of < E α >'s, the event number of IBD is proportional to Lν e while that of NC is proportional to L ν e + Lν e + 4L ν x , the total luminosity. The ratio R is then determined by the ratios of luminosities. While the νp elastic scattering cross section is identical for all flavors and accounts for the total neutrino fluence, the NC spectrum should be the same for both neutrino mass hierarchies. Meanwhile, the IBD spectrum shall be different for neutrino flavor conversions inside the supernova are sensitive to neutrino mass hierarchy. As L ν x /L ν e grows from smaller than one in the accretion phase to larger than one in the cooling phase, the fraction of ν x flux to the total neutrino flux increases. Due to the dense matter inside the SN, theν e flux is fully swapped with theν x flux in the inverted hierarchy by the MSW effect. As a result, Eq. 6 implies that R, the fraction of IBD events to NC events shall increase as the SN neutrino emission evolves from the accretion phase to the cooling phase for the inverted hierarchy. On the contrary, the fraction ofν e flux to the total neutrino flux decreases from the accretion phase to the cooling phase. Eq. 4 then implies that R shall decrease from the accretion phase to the cooling phase for the normal hierarchy.
For both normal and inverted hierarchies, we scan R over the luminosity range of 0.5 ≤ L ν x /L ν e ≤ 2.0 and 0.7 < Lν e /L ν e < 1.5 for three scenarios with mean energies presented in Sect. 2.1. The values of R at JUNO [46] and LENA [47] detectors are shown in Figs. 1 and  2, respectively, in which the plots on the upper panel are for the normal hierarchy while those on the lower panel are for the inverted hierarchy. From left to right, the three columns correspond to three different sets of neutrino mean energies of (< E ν e >, < Eν e >, < E ν x >) = Since L ν e ≈ Lν e , we can fix to Lν e /L ν e = 1 and explore how the R value changes as L ν x /L ν e grows. It is clearly seen, especially along the line at which Lν e /L ν e = 1, that the directions of the change of R are opposed to each other for normal and inverted hierarchies. As the SN explosion evolves from the accretion phase to the cooling phase, L ν x /L ν e increases from less than 1 to greater than 1 along the line from bottom up. Meanwhile, R is decreasing for the normal hierarchy on the upper panel and increasing for the inverted hierarchy on the lower panel. Therefore, the neutrino mass hierarchy can be identified by measuring the change of R from the accretion phase to the cooling phase. To illustrate this, we take a two-phase scenario to model the the time evolution of the SN neutrino emission as proposed in [38]. The energy ratio between flavors in each phase are taken to be where A and C denote the accretion and cooling phases and are marked by triangles and squares on the plots in Fig. 1 and 2 , respectively. The corresponding values of Lν e /L ν e and L ν x /L ν e are shown in Table 1. Assuming a SN explosion with a total energy output of E = 3 × 10 53 erg at a distance of 10 kpc, the values of R and related event numbers are presented in Tables 2 and 3 for JUNO and LENA, rspectively. We note that, for the equipartition scenario marked by circles, the values in the table are calculated for the whole duration of the explosion with a total energy of 3 × 10 53 erg. Hence the numbers are larger than those for the accretion and cooling phases. When the mean energy distribution becomes more hierarchical, the event numbers of both NC and IBD become larger because more neutrinos are shifted to the high energy tail and the the cross sections are larger at higher energies. Since IBD event number depends on the the flux of Eqs. (4) or (6) while NC event number depends on the the total flux of F e + Fē + 4F x , the NC event number grows more than the IBD one. Therefore, we find that R is smaller when the mean energy distribution is more hierarchical in spite of the mass hierarchy as shown in Tables 2 and 3 . For the inverted hierarchy, the IBD event number is determined by F 0 x from Eq. (6). Since the energy fraction of the ν x in the accretion phase is smaller than that in the cooling phase, E ν x ,A /(E ν e ,A + Eν e ,A + 4E ν x ,A ) < E ν x ,C /(E ν e ,C + Eν e ,C + 4E ν x ,C ), the value  Table 2 Numbers of IBD and NC events in three specific scenarios listed in Table 1 Table 3 Numbers of IBD and NC events in three specific scenarios listed in Table 1  of R in the accretion phase should be smaller than that in the cooling phase, R A < R C , for the inverted hierarchy. On the contrary, the energy fraction of theν e in the accretion phase is larger than that in the cooling phase, Eν e ,A /(E ν e ,A + Eν e ,A + 4E ν x ,A ) > Eν e ,C /(E ν e ,C + Eν e ,C + 4E ν x ,C ). Meanwhile, one has |E ν x ,A /(E ν e ,A +Eν e ,A +4E ν x ,A )−E ν x ,C /(E ν e ,C +Eν e ,C + 4E ν x ,C )| < |Eν e ,A /(E ν e ,A +Eν e ,A +4E ν x ,A )−Eν e ,C /(E ν e ,C + Eν e ,C + 4E ν x ,C )|. As a result, R A > R C is inferred from Eq. 4 for the normal hierarchy.
To determine whether R A > R C or R A < R C requires clear discrimination between R A and R C . This can be achieved by requiring the ranges of R A and R C do not over- lap. The ranges of R at JUNO and LENA for the three sets of mean energies are shown in Fig. 3. One can see that, with the measurements of JUNO, NMH can be identified if it is normal but it cannot be determined if it is inverted, since the ranges of R A and R C at JUNO overlap for all the three sets of mean energies in the inverted hierarchy. At LENA, R A and R C are clearly discriminated in between, except for the mean energies of (10 MeV, 15 MeV, 24 MeV) in the inverted hierarchy. Therefore, LENA is capable of determining the NMH for less hierarchical mean energy distributions. Moreover, one can also infer that, even in the most hierarchical distribution of mean energies, the ranges of R A and R C can finally be separated as long as the detector mass is large enough to collect enough events. From Tables 2 and 3 , the difference between R A and R C is smaller in the inverted hierarchy than in the normal hierarchy. This is the reason why the normal hierarchy is easier to be identified than the inverted hierarchy. We note that R values for the same model parameters are different at JUNO and LENA. As shown in Table 4, the scintillation materials in the two detectors are different resulting in different Birk's constants. When applying the same cut of 0.2 MeV to the quenched signal T , the corresponding proton recoil T min 's are different for the two detectors. The higher T min for LENA indicates that fewer fraction of proton recoils are picked as NC events such that the R values are larger for LENA than for JUNO while applied to the same mean energies and luminosity ratios.
The flavor composition in SN neutrino models is conventionally taken to ν e ,ν e , and ν x , with the last accounting for ν μ , ν τ , and their antiparticles. In this work, we consider SN neutrino models with the canonical spectral hierarchy E ν x > Eν e > E ν e and with the lowest E ν e = 10 MeV. Some simulations [48,49] on the SN explosion have indicated a possibility of an even lower E ν e = 8 ∼ 9 MeV with a spectral hierarchy E ν x Eν e > E ν e . With a fixed total energy out burst, the lower average energies implies more SN neutrinos reaching the Earth and hence more events detected. Different sets of average energies will result in different R values. However, in our framework, the NMH identification does not rely on the individual R values. Instead, it is based on the binary relation between R A and R C , which depends on  Fig. 3 Expected R values and uncertainties at JUNO and LENA detectors for both mass hierarchies. Each plot corresponds to a set of mean energies of different flavors denoted at the upper right corner in MeV the energy fractions of flavors as discussed above. Since the condition that L ν e ≈ Lν e > L ν x during the accretion phase and L ν e ≈ Lν e < L ν x during the cooling phase remains unchanged in these simulations, the lower average energies will not change the NMH identification.

Summary and conclusions
We have described how to identify the neutrino mass hierarchy by measuring the time variation of the SN neutrino events inside scintillation detectors. IBD events and NC events at scintillation detectors are taken to define the event ratio R, which can be calculated for given mean energies of and luminosity ratios between different flavors with detector parameters. Our knowledge of the time evolution of SN neutrino emissions indicates that, as the SN explosion evolves from the accretion phase to the cooling phase, the ratio R of IBD events to NC events shall decrease for the normal mass hierarchy and increase for the inverted mass hierarchy. We not only calculate R over physically plausible ranges of luminosity ratios between flavors but also evaluate statistical uncertainties of R arising from measurements for an illustrative model of the SN neutrino emissions. This clarify the detector capability for resolving the neutrino mass hierarchy. We have performed our analysis with three different sets of mean energies of flavors in two scintillation detectors, JUNO and LENA, and found that, as L ν x /L ν e grows, R is increasing for the inverted hierarchy and decreasing for the normal hierarchy in spite of the mean energies of different flavors. With a two-phase scenario to model the SN neutrino emission, we have presented the IBD and NC events, R's, and, σ R 's, deviations of R, numerically to illustrated that R does change in opposite directions as SN neutrinos evolves from the accretion phase to the cooling phase. The detector capability for different mass hierarchies has also been checked and compared for the two detectors, and we have found that the neutrino mass hierarchy will be easier to be identified if it is normal.
In preparation of this article, we noticed a paper [50], in which a "hybrid method" was introduced to determine the NMH. Unlike those efforts on resolving the NMH, the hybrid method, which is similar to our method, determines the NMH by simultaneous measurements of responses from different reaction channels of at least two detectors fromν e and ν e fluxes of accretion and cooling phases. In [50], IBD ofν e on free protons of H 2 O/CH 2 is chosen to be theν e channel and the chosen ν e channel is ν e + 208 Pb → e − + 208 Bi, in which the one-neutron (1n) and two-neutron (2n) decays of bismuth (Bi) are the considered responses. The quantity r is defined as the number ratio ofν e induced events (the IBD events) to the ν e induced events, for which (1n), (2n) and, their sum (tot n) are taken to be different responses and investigated separately. The ranges of the three r ratios in accretion and cooling phases for different SN neutrino models in normal and inverted hierarchies demonstrates the promise to constrain the NMH by using the hybrid method.
The r ratios, like our R ratios, are hierarchy-dependent relative quantities, which are independent of the supernova-Earth distance, of the size of targets, and of the ratios of these quantities. Any specific r ratio is unlikely to determine the NMH. As indicated in [50], a number of these ratios of events from different detectors associated to various particle emissions for different supernova phases should be compared simultaneously for better constraint on the NMH. On the contrary, our method does not rely on simultaneous comparison of a set of numbers to identify the NMH. Instead, the binary relation of a single R ratio in different supernova phases is employed in the identification of the NMH. Compared with the complexity in using the hybrid method in [50], our method is simple, straightforward, and determinate. A SN neutrino model can be characterized by the luminosity, L, the mean energy, < E >, and the shape parameter, η for each flavor. Simulations (see [51] for a review) have shown that η's do not vary much between 2 and 3. Therefore, a common η is assumed throughout this work for all three flavors during both phases for simplicity. In principle, all the parameters, L, < E >, and η, vary with time. To incor-porate the time-dependence of these SN parameters, one can use simulation data of SN neutrino emissions to obtain the event rate, d N/dt. Besides of IBD and NC interactions, SN ν e flux can also be measured in the liquid scintillator [52]. In future studies, we shall work on event rates for different channels of neutrino interactions in various detector to explore the SN explosion with all three flavors of neutrinos taken into account.