Probing New Physics Models of Neutrinoless Double Beta Decay with SuperNEMO

The possibility to probe new physics scenarios of light Majorana neutrino exchange and right-handed currents at the planned next generation neutrinoless double beta decay experiment SuperNEMO is discussed. Its ability to study different isotopes and track the outgoing electrons provides the means to discriminate different underlying mechanisms for the neutrinoless double beta decay by measuring the decay half-life and the electron angular and energy distributions.


Introduction
Oscillation experiments have convincingly shown that at least two of the three active neutrinos have a finite mass and that flavour is violated in the leptonic sector [1]. Despite this success, oscillation experiments are unable to determine the absolute magnitude of neutrino masses. Upper limits on the effective electron neutrino mass of 2.3 eV [2] and 2.05 eV [3] can be set from the analysis of tritium β decay experiments. Astronomical observations combined with cosmological considerations yield an upper bound to be set on the sum of the three neutrino masses of the order of 0.7 eV [4]. However, the most sensitive probe of the absolute mass scale of Majorana neutrinos is neutrinoless double β decay (0νββ) [5,6,7,8]. In this process, an atomic nucleus with Z protons decays into a nucleus with Z + 2 protons and the same mass number A under the emission of two electrons, This process can be described by the exchange of a light neutrino connecting two V-A weak interactions, see Fig. 1 (a). The process (1) is lepton number violating and, in the standard picture of light neutrino exchange, it is only possible if the neutrino is identical to its own anti-particle, i.e. if neutrinos are Majorana particles. Combined with the fact that neutrino masses are more than five orders of magnitude smaller than the masses of other fermions, such a possibility suggests that the origin of neutrino masses is different from that of charged fermions. Several mechanisms of mass generation have been suggested in the literature, the most prominent example being the seesaw mechanism [9] in which heavy right-handed neutrinos mix with the left-handed neutrinos and generate light Majorana masses for the observed active neutrinos. The Majorana character of the active neutrinos can then be connected to a breaking of lepton number symmetry close to the GUT scale and might be related to the baryon asymmetry of the Universe through the baryogenesis via leptogenesis mechanism [10].
Because of its sensitivity to the nature and magnitude of the neutrino mass, 0νββ decay is a crucial experimental probe for physics beyond the Standard Model and its discovery will be of the utmost importance. It will prove lepton number to be broken, and in most models it will also provide direct evidence that the light active neutrinos are Majorana particles 1 [12]. However, the measurement of 0νββ decay in a single isotope is not sufficient to prove that the standard mechanism of light Majorana neutrino exchange is the dominant source for the decay. There are a host of other models, such as Left-Right symmetry [5], Rparity violating Supersymmetry (SUSY) [13] or Extra Dimensions [11], which can provide alternative mechanisms to trigger 0νββ decay. In some of these models, additional sources of lepton number violation can supplement light neutrino exchange. For example, in Left-Right symmetric models, there are additional contributions from righthanded currents and the exchange of heavy neutrinos. In other models, such as R-parity violating SUSY, 0νββ decay can be mediated by other heavy particles that are not directly related to neutrinos.
There are several methods proposed in the literature to disentangle the many possible contributions or at least to determine the class of models that give rise to the dominant mechanism for 0νββ decay. Results from 0νββ decay can be compared with other neutrino experiments and observations such as tritium decay to determine if they are consistent. At the LHC there could also be signs of new physics exhibiting lepton number violation that is related to 0νββ (see [14] for such an example in R-parity violating Supersymmetry). Such analyses would compare results for 0νββ with other experimental searches, but there are also ways to gain more information within the realm of 0νββ decay and related nuclear processes. Possible techniques include the analysis of angular and energy correlations between the electrons emitted in the 0νββ decay [5,15,16,17,18] or a comparison of results for 0νββ decay in two or more isotopes [19,20,21,22]. These two approaches are studied in this paper. Other proposed methods are the comparative analysis of 0νββ decay to the ground state with either 0νβ + β + or electron capture decay [23] and 0νββ decay to excited states [24].
Currently, the best limit on 0νββ decay comes from the search for 0νββ decay of the isotope 76 Ge giving a half-life of T 1/2 > 1.9 · 10 25 years [25]. This results in an upper bound on the 0νββ Majorana neutrino mass of m ν ≤ 300 − 600 meV. A controversial claim of observation of 0νββ decay in 76 Ge gives a half-life of T 1/2 = (0.8 − 18.3) · 10 25 y [26] and a resulting effective Majorana neutrino mass of m ν = 110 − 560 meV. Next generation experiments such as CUORE, EXO, GERDA, MA-JORANA or SuperNEMO aim to increase the half-life exclusion limit by one order of magnitude and confirm or exclude the claimed observation. The planned experiment SuperNEMO allows the measurement of 0νββ decay in several isotopes ( 82 Se, 150 Nd and 48 Ca are currently being considered) to the ground and excited states, and is able to track the trajectories of the emitted electrons and determine their individual energies. In this respect, the SuperNEMO experiment has a unique potential to disentangle the possible mechanisms for 0νββ decay.
This paper addresses the question of how measurements at SuperNEMO can be used to gain information on the underlying physics mechanism of the 0νββ decay. The sensitivity of SuperNEMO to new physics parameters in two models is determined by performing a detailed simulation of the SuperNEMO experimental set-up. By analysing both the angular and energy distributions in the standard mass mechanism and in a model incorporating right-handed currents, the prospects of discriminating 0νββ decay mechanisms are examined. The two models are specifically chosen to represent all possible mechanisms, as they maximally deviate from each other in their angular and energy distributions. This paper is organised as follows. In Section 2 a short description of the theoretical framework on which the calculations of the 0νββ decay rate and the angular and energy correlations are based is shown. The example physics models are introduced and reviewed. Section 3 gives a brief overview of the SuperNEMO experiment design and in Section 4 a detailed account of the simulation analysis and its results are presented. In Section 5, the expected constraints from SuperNEMO on new physics are shown and the prospects of disentangling 0νββ mechanisms by analysing the angular and energy distributions and by comparing rates in different isotopes are addressed. Our conclusions are presented in Section 6.

Effective Description
Contributions to 0νββ decay can be categorised as either long-range or short-range interactions. In the first case, the corresponding diagram involves two vertices which are both point-like at the Fermi scale, and connected by the exchange of a light neutrino. Such long-range interactions are described by an effective Lagrangian [27,28] where G F is the Fermi coupling constant and the leptonic and hadronic Lorentz currents are defined as j a =ēO a ν and J a =ūO a d, respectively. Here, O a denotes the corresponding transition operator, with a = V − A, V + A, S − P, S + P, T L , T R [27]. In Equation (2), the contribution from V − A currents originating from standard weak couplings has been separated off and the summation runs over all Lorentz invariant and non-vanishing combinations of the leptonic and hadronic currents, except for the case a = b = V − A. The effective coupling strengths for longrange contributions are denoted as ǫ lr ab . For short-ranged contributions, the interactions are represented by a single vertex which is point-like at the Fermi scale, and they are described by the Lagrangian [28,29] Here, m p denotes the proton mass and the leptonic and hadronic currents are given by J a = uO a d and j ′ a = eO a e C , respectively. The transition operators O a are defined as in the long-range case above, and the summation runs over all Lorentz invariant and non-vanishing combinations of the hadronic and leptonic currents. The effective coupling strengths for the short-range contributions are denoted as ǫ sr abc . Described by the first term in Equation (2), the exchange of light left-handed Majorana neutrinos leads to the 0νββ decay rate where m ν is the effective Majorana neutrino mass in which the contributions of the individual neutrino masses m i are weighted by the squared neutrino mixing matrix elements, U 2 ei , m ν = | i U 2 ei m i |. Analogously, other new physics (NP) contributions, of both long-and short-range nature, can in general be expressed as where ǫ N P denotes the corresponding effective coupling strength, i.e. is either given by ǫ lr ab for a long-range mechanism or by ǫ sr abc for a short-range mechanism. In Equations (4) and (5), the nuclear matrix elements for the mass mechanism and alternative new physics contributions are given by M mν and M N P , respectively, and G 01 , G N P denote the phase space integrals of the corresponding nuclear processes. It is assumed that one mechanism dominates the double β decay rate.

Left-Right Symmetry
The focus in this paper is on a subset of the Left-Right symmetric model [5], which incorporates left-handed and right-handed currents under the exchange of light and heavy neutrinos. Left-Right symmetric models generally predict new gauge bosons of the extra right-handed SU(2) gauge symmetry as well as heavy right-handed neutrinos   [30].
which give rise to light observable neutrinos via the seesaw mechanism. The 0νββ decay half-life in the Left-Right symmetric model can be written as a function of the effective parameters µ, η, λ [30], where contributions from the exchange of heavy neutrinos are omitted. The coefficients C mm , C ηη etc. are combinations of phase space factors and nuclear matrix elements. The first three terms give the contributions from the following processes: 1. C mm µ 2 : Fully left-handed current neutrino exchange, see Fig. 1 (a) (mass mechanism); 2. C λλ λ 2 : Right-handed leptonic and right-handed hadronic current neutrino exchange, see Fig. 1 (b); 3. C ηη η 2 : Right-handed leptonic and left-handed hadronic current neutrino exchange.
The remaining terms in Equation (6) describe interference effects between these three processes. The effective parameters µ, η, λ in (6) are given in terms of the underlying physics parameters as with the electron mass m e , the left-and right-handed W boson masses M WL and M WR , respectively, and the mixing angle ζ between the W bosons. The 3 × 3 matrices U 11 and U 12 connect the weak eigenstates (ν e , ν µ , ν τ ) of the light neutrinos with the mass eigenstates of the light neutrinos (ν 1 , ν 2 , ν 3 ), and heavy neutrinos, (N 1 , N 2 , N 3 ), respectively. We assume that the neutrino sector consists of three light neutrino states, m νi ≪ m e , and three heavy neutrino states, M Ni ≫ m p , i = 1, 2, 3. Consequently, the summations in (7,8,9) are only over the light neutrino states. For a simple estimate of the sensitivity of 0νββ decay to the model parameters, we neglect the flavour structure in U 11 and U 12 ; using the assumption that the elements in U 11 are of order unity (almost unitary mixing), and those in U 12 are of order m D /M R ∼ m ν /M R , with the effective magnitude m D of the neutrino Dirac mass matrix, and the light and heavy neutrino mass scales, m ν and M R , leads to the approximate relations: In the following analysis a simplified model incorporating only an admixture of mass mechanism (MM) due to a neutrino mass term µ = m ν /m e and right-handed current due to the λ term (RHC λ ) is considered: As we will see in Section 2.4, these two mechanisms exhibit maximally different angular and energy distributions, and with an admixture between them, to a good approximation any possible angular and energy distribution can be produced. In our numerical calculation we use the values as given in Table 1 for the coefficients C mm , C λλ and C mλ in Equation (13). Furthermore, we assume that the parameter µ is real-valued positive and λ is real-valued.

Nuclear Matrix Elements
As demonstrated in Equations (4) and (5), a calculation of the nuclear matrix elements (NMEs) is required to convert the measured half-life rates or limits into new physics parameters. Exact solutions for the NMEs do not exist, and approximations have to be used. Calculations using the nuclear shell model exist for lighter nuclei such as 76 Ge and 82 Se, though the only reliable results are for 48 Ca. Quasi-particle random phase approximation (QRPA) calculations are applied for most isotopes as a greater number of intermediate states can be included. In this paper, a comparison between two possible SuperNEMO isotopes ( 82 Se and 150 Nd) and the isotope that gives the current best limit ( 76 Ge) is made. Consistent calculations of the NMEs for these three isotopes in both the MM and RHC are rare (only [30] and [31]). All the results are shown using NMEs from [30], displayed in Table 1. Recent work on the calculation of NMEs for the heavy isotope 150 Nd suggests that nuclear deformation must be included, as QRPA calculations usually consider the nucleus to be spherical. To compensate for this a suppression factor of 2.7 is introduced into the NME due to an approximation arising from the BCS overlap factor [32], M( 150 Nd)/2.7. This gives a suppression C mm,λλ,mλ /(2.7) 2 in Table 1. The 82 Se nuclei are assumed to be spherical and no correction is added in this paper.
The NMEs are a significant source of uncertainty in double β decay physics and quantitative results in this paper could change with different calculations (particularly for 150 Nd). For example, more recent studies [33] suggest the NMEs from 150 Nd for the MM are an additional factor 1.3-1.7 smaller. In our analysis we assume a theoretical uncertainty of 30% in the NMEs of all isotopes and mechanisms considered throughout. Even though the choice of NME changes quantitative results for the extracted physics parameters, the conclusions about the advantages of using different kinematic variables will not be affected.

Angular and Energy Distributions in the Left-Right Symmetric Model
For our event simulation, the three-dimensional distribution of the 0νββ decay rate in terms of the kinetic energies t 1,2 of the two emitted electrons and the cosine of the angle between the electrons cos θ 12 is used: The distributions for the MM and for the RHC λ mechanism are given by with the electron momenta p i = t i (t i + 2) and velocities β i = p i /(t i + 1), and the mass difference Q between the mother and daughter nucleus. All energies and momenta are expressed in units of the electron mass and c 1 and c 2 are normalisation constants. The Fermi function F is given by where s = 1 − (αZ) 2 , u = αZ(t + 1)/p, α = 1/137.036, Γ is the Gamma function and c 3 is a normalisation constant. Here, Z is the atomic number of the daughter nucleus. The normalisation of the distributions is irrelevant when discussing energy and angular correlations. Using Equations (15) and (16), the differential decay widths with respect to the cosine of the angle θ 12 , may be determined.
The differential width in Equation (18) can be written as [5,18] dΓ with the total decay width Γ . The distribution shape is linear in cos θ 12 , with the slope determined by the parameter k θ which can range between −1 ≤ k θ ≤ 1, depending on the underlying decay mechanism. Assuming the dominance of one scenario, the shape does not depend on the precise values of new physic parameters (mass scales, coupling constants) but is a model specific signature of the mechanism. For the MM and RHC λ mechanisms, the theoretically predicted k θ is found from Equation (18) and is given by The correlation coefficient k θ is modified when taking into account nuclear physics effects and exhibits only a small dependence on the type of nucleus. The MM and the RHC λ mechanisms give the maximally and minimally possible values for the angular correlation coefficient k θ in a given isotope, respectively. Experimentally, k θ can be determined via the forwardbackward asymmetry of the decay distribution, Here, N + (N − ) counts the number of signal events with the angle θ 12 larger (smaller) than 90 • . Analogously, the MM and RHC λ mechanism also differ in the shape of the electron energy difference distribution, Equation (19). For the isotopes 82 Se and 150 Nd, these distributions are shown in Fig. 2. Again, the shape is largely independent of the isotope under inspection. The following asymmetry in the electron energy-difference distribution is determined, thereby defining an energy correlation coefficient k E , where Q is the energy release of the decay. The rate N + (N − ) counts the number of signal events with an electron energy difference smaller (larger) than Q/2. For the MM and RHC λ mechanism, the theoretical k E parameter may be found from Equation (19) and is given by in the isotopes 82 Se and 150 Nd. As can be seen in the angular distribution, the corresponding energy correlation coefficients in the two mechanisms considered are, to a good approximation, at their upper and lower limits in a given isotope.

SuperNEMO
SuperNEMO is a next generation experiment building on technology used by the currently running NEMO-III experiment [34,35,36,37,38,39,40]. The design of the detector consists of 20 modules each containing approximately 5 kg of enriched and purified double β emitting isotope in the form of a thin foil (with a surface density of 40 mg/cm 2 ). Isotopes under consideration for Super-NEMO are 82 Se, 150 Nd and 48 Ca.
The foil is surrounded by a tracking chamber comprising nine planes of drift cells (44 mm diameter) on each side operating in Geiger mode in a magnetic field of 25 Gauss. The tracking chamber has overall dimensions of 4 m height (parallel to the drift cells), 5 m length and 1 m width (perpendicular to the foil); the foil is centred in this volume with dimensions of 3 m height and 4.5 m length. The tracking allows particle identification (e − , e + , γ, α) and vertex reconstruction to improve background rejection and to allow measurement of double β decay angular correlations.
Calorimetry consisting of 25×25 cm 2 square blocks of 5 cm thickness scintillating material connected to low activity photomultiplier tubes (PMTs) surrounds the detector on four sides. An energy resolution of 7% (FWHM) and time resolution of 250 ps (Gaussian σ) at 1 MeV for the blocks is required. The granularity of the calorimetry allows the energy of individual particles to be measured. Additional γ-veto calorimetry to identify photons from background events of thickness 10 cm surrounds the detector on all sides. The modules are contained in shared background shielding and will be housed in an underground laboratory to reduce the cosmic ray flux. A diagram of the planned SuperNEMO module design is shown in Fig. 3.

Simulation Description
A full simulation of the SuperNEMO detector was performed including realistic digitisation, tracking and event selection. Signals for two mechanisms of 0νββ decay (mass mechanism MM and right-handed current via the λ parameter RHC λ ) and the principal internal backgrounds were generated using DECAY0 [41]. This models the full event kinematics, including angular and energy distributions.
A GEANT-4 Monte Carlo simulation of the detector was constructed. Digitisation of the hits in cells was obtained by assuming a Geiger hit model validated with NEMO-III with a transverse resolution of 0.6 mm and a longitudinal resolution of 0.3 cm. The calorimeter response was simulated assuming a Gaussian energy resolution of 7%/ √ E (FWHM) and timing resolution of 250 ps (Gaussian σ at 1 MeV). Inactive material in front of the γ-veto was partially simulated.
Full tracking was developed consisting of pattern recognition and helical track fitting. The pattern recognition uses a cellular automaton to select adjacent hits in the tracking layers. Helical tracks are fitted to the particles. Tracks are extrapolated into the foil to find an appropriate event origin and into the calorimeter where they may be associated with calorimeter energy deposits. The realistic event selection (validated using NEMO-III) was optimised for double β decay electrons (two electrons with a common vertex in the foil). The selection criteria are: events must include only two negatively charged particles each associated with one calorimeter hit; event vertices must be within the foil and the tracks must have a common vertex of <10 standard deviations between intersection points in the plane of the source foil; the time of flight of the electrons in the detector must be consistent with the hypothesis of the electrons originating in the source foil; the number of Geiger drift cell hits unassociated with a track must be less than 3; the energy deposited in individual calorimeter blocks must be > 50 keV; there are zero calorimeter hits not associated with a track; tracks must have hits in at least one of the first three and one of the last three planes of Geiger drift cells; the number of delayed Geiger drift cell hits due to α particles from 214 Bi-214 Po events must be zero; there are no hits in the γ-veto detectors with energy > 50 keV.
Using these experimental selection criteria the signal efficiency was found to be 28.2% for the MM and 17.0% for the RHC λ in 82 Se and 29.1% for the MM and 17.3% for the RHC λ in 150 Nd. This is higher than the efficiency for MM detection in 100 Mo decays in NEMO-III of 17.4% (in the electron energy sum range 2-3.2 MeV) [40].
The variables reconstructed from the simulation are the energy sum, where a peak is expected at the energy release, Q, of the 0νββ decay, the energy difference and the cosine of the opening angle of the two electrons. Simulations of the angular and energy difference distributions of the two electrons in a signal sample are shown in Fig. 4 for the isotope 82 Se (similar results hold for 150 Nd). The reconstructed distributions, normalised to the theoretical distributions, show detector effects which arise due to multiple scattering in the source foil, compared to the theoretically predicted distributions based on Equations (18) and (19). This influence is particularly strong in the righthanded current as one electron usually has low energy so the shape of the distribution is changed (on average a 30 • deviation from the generated distribution). The reconstruction efficiency is also low for small angular separation between the electrons when they travel through the same drift cells.
The backgrounds were processed by the same detector simulation and reconstruction programs as the signal. The dominant two neutrino double β decay (2νββ) background and the background due to foil contamination were normalised assuming a detector exposure of 500 kg y. Due to the high decay energy Q for 0νββ in 150 Nd, the 214 Bi background may be neglected. The activities were assumed to be 2 µBq/kg for 208 Tl and 10 µBq/kg for 214 Bi. These are the target radioactive background levels in the baseline SuperNEMO design. Reconstructed distributions of the experimental variables including background events for the MM at an example signal half-life of 10 25 y are shown in Figure 5.

Limit Setting
To determine the longest half-life that can be probed with SuperNEMO, exclusion limits at 90% CL on the halflife using the distribution of the sum of electron energies ( Fig. 5 (a)) were set using a Modified Frequentist (CL s ) [42] method. This method uses a log-likelihood ratio (LLR) of the signal-plus-background hypothesis and the background-only hypothesis, where the signal is due to the 0νββ process. The effect of varying the 214 Bi background activities on the expected limit to the MM is shown in Fig. 6. The expected limit is given by the median of the distribution of the LLR and the widths of the bands shown represent one and two standard deviations of the LLR distributions for a given 214 Bi activity. For comparison, the NEMO-III internal 214 Bi background is < 100 µBq/kg in 100 Mo and 530 ± 180 µBq/kg in 82 Se. The NEMO-III internal 208 Tl background is 110 ± 10 µBq/kg in 100 Mo, 340 ± 50 µBq/kg in 82 Se and 9320 ± 320 µBq/kg in 150 Nd [39]. The γ-veto used reduces the number of radioactive background events by 30% for 214 Bi in the electron energy sum window > 2.7 MeV.
All external backgrounds from outside the foil, apart from radon in the tracking volume, are expected to be negligible and were not simulated. The energy distribution of the external radon background is similar to the internal background. Simulations have shown that a contamination of 10 µBq/kg of 214 Bi in the foil is equivalent to 280 µBq/m 3 of 214 Bi in the gas volume and 2 µBq/kg of 208 Tl in the foil is equivalent to 26 µBq/m 3 of 208 Tl in the gas volume. Figure 6 shows that this level of external background would lead to a ∼ 15% reduction in the halflife limit. The dominant 2νββ background is measured by SuperNEMO and statistical uncertainties on its half-life are expected to be negligible. Inclusion of an estimated 7% correlated systematic uncertainty on the signal and background distributions [35] leads to a ∼ 5% reduction in the MM half-life limit. The effects of external background and of systematic uncertainties on the 2νββ background are not included in the results of this paper.
Expected exclusion limits at 90% confidence level on the half-life are shown in Fig. 7. Results are displayed as a function of RHC λ admixture, where the signal distribution is produced by combining weighted combinations per bin of the MM and RHC λ contributions at the event level. An admixture of 0% corresponds to a pure MM contribution, and an admixture of 100% to pure RHC λ . Interference terms are assumed to be small and are neglected in the experimental simulation. The lower efficiency in the case of RHC λ results in a lower limit for larger admixtures. The half-life limit is approximately twice as sensitive in measurements of 82 Se due to the lower mass number and higher 2νββ decay half-life, though this is compensated in 150 Nd by more favourable phase space when converting into physics parameter space. In the case where one mechanism dominates SuperNEMO is expected to be able to exclude 0νββ half-lives up to 1.2 · 10 26 y (MM) and 6.1 · 10 25 y (RHC λ ) for 82 Se, and 5.1 · 10 25 y (MM) and 2.6 · 10 25 y (RHC λ ) for 150 Nd.

Observation
A 0νββ signal rate with significant excess over the background expectation, as for example shown in The angular asymmetry parameter k θ in Equation (23) is experimentally accessible by defining N + as the number of events with measured angle cos θ < 0 and N − as the number of events with cos θ > 0. Similarly, an energy difference asymmetry k E can be obtained where N + is the number of events with energy difference < Q/2 (half the energy of the 0νββ decay) and N − is the number of events with energy difference > Q/2. The electron energy sum is required to be greater than 2.7 MeV for 82 Se and 3.1 MeV for 150 Nd to maximise signal to background ratio. This results in signal efficiencies of 23.2% for the MM and 13.2% for the RHC λ in 82 Se and 19.1% for the MM and 10.4% for the RHC λ in 150 Nd.
Experimentally, the distributions are only available as a sum of signal plus background so the measured values differ from the theoretically expected values due to the background distributions. This generally results in reconstructed correlation factors that are biased towards positive values. The measured values of k θ,E are shown in Fig. 9 for a number of half-lives in the two isotopes. Statistical uncertainties are shown as the width of the bands.
All reconstructed k θ,E values are displayed as a function of the corresponding theoretical k T θ,E parameter, to allow for a model independent generalisation. It can be seen that the energy difference distribution allows stronger model discrimination than the angular distribution.

Model Parameter Constraints
Having performed a detailed experimental analysis including a realistic simulation of the detector setup, the results are interpreted in terms of the expected reach of the SuperNEMO experiment to new physics parameters of the combined MM and RHC λ model of 0νββ decay.
Using Equation (13) for the 0νββ decay half-life together with the coefficients listed in Table 1, the expected 90% CL limit on T 1/2 shown in Fig. 7 can be translated into a constraint on the model parameters m ν and λ. Assuming all other contributions are negligible this is shown in Fig. 10 (a), as a contour in the m ν − λ parameter plane. In the case SuperNEMO does not see a signal these parameters would be constrained to be located within the coloured contour. The odd shape of the coloured contour is a direct consequence of the SuperNEMO 90% CL exclusion limit as a function of the specific admixture between the MM and the RHC λ shown in Fig. 7. The small interference term, though not included in the experimental simulation, is taken into account through Equation (13) in this figure and results in the asymmetry of the distribution with respect to the sign of the parameter λ.
As shown in Section 4, SuperNEMO is expected to be more sensitive to the 0νββ half-life when using the isotope 82 Se, but this is compensated by the larger phase space of 150 Nd. As a result, the constraint on the model parameters is slightly stronger for 150 Nd. Due to the large uncertainty in the NMEs, this might be different for other NME calculations. To demonstrate the improvement over existing experimental bounds, the parameter constraints are shown in Fig. 10 (b) on a logarithmic scale (for positive values of λ), comparing the expected SuperNEMO reach with the current constraints from the 0νββ experiments NEMO-III [40,38] and Heidelberg Moscow [25]. Figure 10 shows that SuperNEMO is expected to constrain model parameters at 90% CL down to m ν =70-73 meV and λ=(1-1.3)·10 −7 . This would be an improvement by a factor 5-6 over the current best limit from the Heidelberg Moscow experiment and more than an order of magnitude compared to the NEMO-III results.

Angular and Energy Correlations
As a more intriguing scenario it is now assumed that SuperNEMO actually observes a 0νββ decay signal in 82 Se or 150 Nd. Because of the tracking abilities described in Section 4 this opens up the additional possibility of measuring the angular and energy distribution of the decays. Depending on the number of signal events detected, this can be crucial in distinguishing between different 0νββ decay mechanisms. In the analysis a reconstruction of the angular and energy correlation coefficients k θ and k E is used to determine the theoretical coefficients, and thereby the admixture between the left-and right-handed currents in the combined MM and RHC λ model.
For the isotope 82 Se, this is shown in Fig. 11 for different RHC λ admixtures. The two blue elliptical contours correspond to the allowed one standard deviation (m ν −λ) parameter space at SuperNEMO when observing a signal at T 1/2 = 10 25 y and T 1/2 = 10 26 y, respectively. This takes into account a nominal theoretical uncertainty on the NME of 30% and a one standard deviation statistical uncertainty on the measurement determined from the simulation (Fig. 8). The blue elliptical error bands therefore give the allowed parameter region when only considering the total 0νββ rate, which does not allow to distinguish between the MM and RHC λ contributions.
When taking into account the information provided by the angular and energy difference distribution shape, the parameter space can be constrained significantly. This is shown using the green contours in Fig. 11 for (a) a pure MM model, (b) a RHC λ admixture of 30%, corresponding to an angular correlation of k θ ≈ 0.4 and (c) a pure RHC λ model. This fixes two specific directions in the m ν − λ plane (one for positive and one for negative λ). The widths of the contours are determined by the uncertainty in determining the theoretical correlation and admixture from the apparent distribution shape, see Fig. 9. The outer (light green) contours in Fig. 11 give the one standard deviation uncertainty on the parameters from reconstructing the angular distribution, while the inner (darker green) contour gives the one standard deviation uncertainty when using the distributions of the electron energy difference. As was outlined in Section 4, the energy difference distribution is expected to be easier to reconstruct and therefore gives a better determination of the RHC λ admixture and a better constraint. While interference between MM and RHC λ is neglected in the simulation, it is taken into account in Equation (13) through the term C mλ µλ resulting in the slightly tilted elliptical contours and the asymmetry for λ ↔ −λ. Finally, the red contours in Fig. 11 show the constraints on the model parameters when combining both the determination of the 0νββ decay rate and the decay energy distribution. This demonstrates that such a successful combination can make it possible to determine the mechanism (i.e. the degree of MM and RHC λ admixture in this case), and provide a better constraint on the model parameters. From Fig. 11 (a), the Majorana mass term can be determined at m ν = 245 +56 −41 meV while the λ parameter is constrained to be −0.87 · 10 −7 < λ < 0.92 · 10 −7 in the case of a measured 0νββ decay half-life of 82 Se of T 1/2 = 10 25 y. For a 82 Se half-life of T 1/2 = 10 26 y, the uncertainty on the decay rate increases as SuperNEMO reaches its exclusion limit for RHC λ admixtures. It is therefore only possible to extract upper limits on the model parameters from Fig. 11 for T 1/2 = 10 26 y. However, the shape information provides additional constraints on the parameter space. In Fig. 12 we show the analogous plots for the isotope 150 Nd assuming a decay half-life of T 1/2 = 10 25 y.  (1) and (2) using decay rate and energy distribution shape reconstruction (red contours). The admixture of the MM and RHC λ contributions is assumed to be: (a) pure MM contribution; (b) 30% RHC λ admixture; and (c) pure RHC λ contribution. NME uncertainties are assumed to be 30% and experimental statistical uncertainties are determined from the simulation.

Rate Comparison of 150 Nd and 82 Se
While reconstruction of the decay distribution can be an ideal way to distinguish between different mechanisms, it might be of little help if 0νββ decay is observed close to the exclusion limit of SuperNEMO, or not at all. This is demonstrated in Fig. 11 where, for a half-life of T 1/2 = 10 26 y, the reconstruction of the energy difference distribution will be problematic due to the low number of events (compare Fig. 9). As an alternative, it is possible to compare the 0νββ rate in different isotopes. This method, which could provide crucial information close to the exclusion limit, is especially relevant for SuperNEMO which could potentially measure 0νββ decay in two (or more) isotopes. Such a comparative analysis was used in [21] to distinguish between several new physics mechanisms. A combined analysis of several isotopes, potentially measured in other experiments, will improve the statistical significance [22].
The possibility of sharing the two isotopes equally in SuperNEMO, each with a total exposure of 250 kg y, is now considered. In the cases where the MM or the RHC λ contributions dominate, the following half-life ratios can be found: These ratios and their uncertainties are determined by the 0νββ decay NMEs and phase spaces. The factor 2.7 is the correction added to the 150 Nd NMEs as described in Section 2.3. It has recently been suggested that uncertainties in NME calculations are highly correlated [43] so measurements made with two or more isotopes could reduce the uncertainty on the physics parameters significantly. Additionally, most experimental systematic uncertainties would cancel if different isotopes are studied in a single experiment such as SuperNEMO. This would not be possible when comparing results with other experiments. The statistical uncertainties are naturally greater than in the single-isotope case, due to the exposure being halved for each isotope.
The results of the combined NME and statistical uncertainties analysis, including a possible correlation of the NMEs, are illustrated in Fig. 13. It shows the 0νββ halflife of 150 Nd as a function of the half-life in 82 Se assuming a pure MM model, with the coloured contours giving the deviation from the hypothesis that the mass mechanism is the single source of 0νββ decay in both isotopes at the 1, 2 and 5 standard deviation level. The statistical uncertainties used in Fig. 13 are derived from our experimental simulation and the standard 30% NME uncertainties are applied. The effect of a possible correlation of the NMEs is shown by assuming the covariance coefficient between the NME uncertainties of 82 Se and 150 Nd to be (a) zero (no correlation), (b) 0.7 and (c) 1.0 (full correlation). The experimental uncertainties and expected sensitivity (90% CL exclusion) limits are calculated for 250 kg y of exposure of each isotope and assume a 50% 82 Se and 50% 150 Nd option for SuperNEMO. The red line shows the relationship for the half-life ratio in the pure RHC λ model (Equation (28)). It can be seen that an exclusion at two standard deviations is possible if the NME errors are perfectly correlated and at the one standard deviation level if the correlation is 70%, which is a more realistic assumption.
Other mechanisms have different half-life ratios [21] so they could be excluded with different CLs at Super-NEMO. One important advantage of this method is that it provides a possibility to falsify the mass mechanism as the sole source for 0νββ. A measurement within the blue contour would indicate that the pure MM model can be excluded at the 5 standard deviation level and new physics is required to explain the measured half-lives.

Combined Energy and Rate Comparison of 150 Nd and 82 Se
In the most favourable case, signal event rates in two isotopes could be high enough (0νββ decay half-lives small enough) that the distribution method and the two isotope rate analysis can be combined to put further constraints on the parameter space. The effect of such a combined analysis on the allowed parameter space is shown in Fig. 14, where the 50% 150 Nd -50% 82 Se two-isotope option (red contours) is compared to the single-isotope options 100% 82 Se (green contours) and 100% 150 Nd (blue contours). The 0νββ decay half-life of 82 Se is assumed to be 10 25 y, and the half-life of 150 Nd is determined by the respective MM-RHC λ admixture, i.e. (a) T N d 1/2 = 10 25 /2.45 y, (b) 10 25 /2.73 y and (c) 10 25 /3.64 y. The NME uncertainties are assumed to be 30% with a 0.7 covariance between the uncertainties of the NMEs of 82 Se and 150 Nd. As can be seen in Fig. 14, the two-isotope option can improve the constraints on the parameter space along the radial direction, e.g. it allows a more accurate determination of the MM neutrino mass m ν in Fig. 14 (a). On the other hand, the accuracy in the lateral direction (the parameter λ in Fig. 14 (a)) becomes worse compared to the best single-isotope option due to the reduced statistics for a given isotope.

Conclusion
The 0νββ decay is a crucial process for physics beyond the Standard Model, and the next generation SuperNEMO experiment is designed to be a sensitive probe of this lepton number violating observable. In addition to being able to measure the 0νββ half-life of one or more isotopes, it also allows the determination of the angular and energy difference distributions of the outgoing electrons.
In this paper we have focussed on the sensitivity of SuperNEMO to new physics and its ability to discriminate between different 0νββ mechanisms. This was achieved by a detailed analysis of two important models, namely the standard mass mechanism via light left-handed Majorana neutrino exchange and a contribution from right-handed current via the effective λ parameter stemming from Left-Right symmetry. The study included a full simulation of the process and the SuperNEMO detector at the event level, allowing a realistic estimation of the experimental 90% CL exclusion limit and statistical uncertainties.
SuperNEMO is expected to exclude 0νββ half-lives up to 1.2 · 10 26 y (MM) and 6.1 · 10 25 y (RHC λ ) for 82 Se and 5.1·10 25 y (MM) and 2.6·10 25 y (RHC λ ) for 150 Nd at 90% CL with a detector exposure of 500 kg y. This corresponds to a Majorana neutrino mass of m ν ≈ 70 meV and a λ parameter of λ ≈ 10 −7 , giving an improvement of more than one order of magnitude compared to the NEMO-III experiment.
It has been shown that the angular and electron energy difference distributions can be used to discriminate new physics scenarios. In the framework of the two mechanisms analysed, it was demonstrated that using this technique the individual new physics model parameters can be determined. For a half-life of T 1/2 = 10 25 y with an exposure of 500 kg y, the Majorana neutrino mass can be determined to be 245 meV with an uncertainty of 30% while the λ parameter can be constrained at the same time to be smaller than |λ| < 0.9 · 10 −7 . Such a decay distribution analysis could be easily extended further to include other new physics scenarios with distinct distributions and the 10 Fig. 13: The 0νββ half-life of 150 Nd as a function of measured half-life in 82 Se for the hypothesis that the MM is the single 0νββ decay source. The contours show the 1, 2 and 5 standard deviation levels assuming statistical uncertainties derived from the experimental simulation and 30% NME errors assumed to have (a) no, (b) 0.7 and (c) perfect correlation. The experimental uncertainties and expected sensitivity (90% CL exclusion) limit are calculated for 250 kg y of exposure (assuming a 50% 82 Se and 50% 150 Nd option). The red line shows the relationship for the RHC λ . The blue contour shows the 5σ exclusion of the MM.  (1) and (2) with an exposure of 250 kg y in 82 Se and 150 Nd (red contour). The admixture of the MM and RHC λ contributions is assumed to be: (a) pure MM contribution; (b) 30% RHC λ admixture; and (c) pure RHC λ contribution. NME uncertainties are assumed to be 30% with a correlation of the uncertainties of 0.7, and experimental statistical uncertainties are determined from the simulation.
results are quoted in terms of a generalised distribution asymmetry parameter to allow new physics scenarios to be compared. As the two example mechanisms considered exhibit maximally different angular and energy distribution shapes, they serve as representative scenarios covering a large spectrum of the model space. For example, the right-handed current contribution due the effective η parameter, also arising in Left-Right symmetrical models, can be distinguished from the mass mechanism and the right-handed current λ contribution by looking at both the angular and energy difference decay distribution. This would allow a determination of all three model parameters m ν , λ and η by looking at the total rate and the angular and energy difference distribution shapes.
Further insight into the mechanism of 0νββ can be gained by using multiple isotopes within the SuperNEMO setup. This possibility was explored by studying the option of having 50% 150 Nd and 50% 82 Se, each with an exposure of 250 kg y. While the statistics per isotope is reduced compared to the individual 100% options, the ability to measure the ratio between the half-lives of the two isotopes can be used as additional information on the underlying physics mechanism responsible for 0νββ decay. As was shown for the isotopes 82 Se and 150 Nd at Super-NEMO, this could be a powerful method to falsify the mass mechanism as the dominant 0νββ mechanism. A correlation between the uncertainties of nuclear matrix elements, which is generally expected on theoretical grounds, has proven to be of importance and its impact on the falsification potential was analysed. Within SuperNEMO such a correlation could also be found between the systematic uncertainties in the measurements of different isotopes.
SuperNEMO also has a number of other possibilities to disentangle the underlying physics. The detection technology is not dependent on one particular isotope and any double β emitting source could be studied in the detector. In this paper 82 Se and 150 Nd have been considered but other isotopes such as 48 Ca or 100 Mo are feasible. The analysis can be extended to cover more than two isotopes thereby achieving a higher significance and a comparison with other experimental results will provide additional information. SuperNEMO is also able to measure a 0νββ decay to an excited state, by measuring two electrons and an accompanying photon. This again could be used to aid the analysis to discriminate between new physics mechanisms.
A combination of the above methods makes Super-NEMO an exciting test of new physics. These methods would prove invaluable in excluding or confirming dominating mechanisms of lepton number violation in the reach of the next generation 0νββ experiments.