Double-beta decay of $$^{130}\hbox {Te}$$ to the first $$0^+$$ excited state of $$^{130}\hbox {Xe}$$ with CUORE-0

We report on a search for double beta decay of 130Te to the first 0+ excited state of 130Xe using a 9.8 kg · yr exposure of 130Te collected with the CUORE-0 experiment. In this work we exploit different topologies of coincident events to search for both the neutrinoless and two-neutrino double beta decay modes. We find no evidence for either mode and place lower bounds on the half-lives: T 0ν 0+ 1 > 7.9 · 1023yr and T 2ν 0+ 1 > 2.4 · 1023yr (90% CL). Combining our results with those obtained by the CUORICINO experiment, we achieve the most stringent constraints available for these processes: T 0ν 0+ 1 > 1.4 · 1024yr and T 2ν 0+ 1 > 2.5 · 1023yr


Introduction
Two-neutrino (2νββ) [1] and neutrinoless (0νββ) [2] double beta decay are among the rarest decay processes studied. While the former is allowed by the Standard Model and has been experimentally detected in a number of isotopes [3], the latter has never been observed; its discovery would imply that lepton number is not conserved and that neutrinos are in fact Majorana particles [4,5].
The CUORE experiment (Cryogenic Underground Observatory for Rare Events) [6][7][8][9], which is currently running at Laboratori Nazionali del Gran Sasso (LNGS), is designed to perform a high-sensitivity search for 0νββ decay of 130 Te to the ground state of 130 Xe [10]. The active isotope is contained in TeO 2 crystals, which are operated as thermal detectors in a cryostat capable of reaching temperatures below 10 mK. At this temperature, the crystal heat capacity becomes very small and consequently a release of energy within a crystal results in a detectable increase of its temperature. The soughtafter experimental signature of 0νββ decay is a monochromatic peak in the summed energy spectrum of the final state electrons at 2527.518 ± 0.013 keV [11][12][13], which is the transition energy of the decay. To maximize the sensitivity S. J. Freedman: Deceased. a e-mail: cuore-spokesperson@lngs.infn.it of the search, the radioactive background at the transition energy must be kept as low as possible.
The first tower assembled in the CUORE assembly line was operated as a standalone experiment, named CUORE-0 [14], from 2013 to 2015. CUORE-0 was designed to validate several key aspects of CUORE, including detector construction, data acquisition and the analysis framework. In addition to this, CUORE-0 provided a sensitive probe of several rare decays, including 0νββ [15,16] and 2νββ [17] decay of 130 Te to the ground state of 130 Xe, and the β + EC decay of 120 Te [18]. In this work we focus on a search for double-beta decay of 130 to the first 0 + excited state of 130 Xe (ββ 0 + 1 ) with CUORE-0. As shown in Fig. 1 this decay emits two electrons, which share a maximum energy of 734 keV, followed by a gamma cascade to the ground state of 130 Xe. The most probable de-excitation pattern, which has a branching ratio of 86.0%, involves the emission of two gamma rays with energies of 1257.5 (3) keV and 536.1 keV. Two more patterns are also possible, namely the emission of three gamma rays with energies of 536.1 keV, 586.0 keV and 671.3 keV (branching ratio of 12.2%) and the emission of two gamma rays with energies of 671.3 keV and 1122.1 keV (branching ratio of 1.8%). These gamma lines result in multi-detector coincidence signatures which we exploit in our analysis to achieve very powerful background rejection. As is the case for the decay to the ground state, the summed energy spectrum of the emitted electrons is distinctly different in the 2νββ 0 + 1 vs. the 0νββ 0 + 1 case. The former is a continuous spectrum (0-734 keV), whereas the latter is a monochromatic peak centered at 734 keV.
Prior to the current work, the most stringent constraints on these decays came from the CUORICINO [19] experiment, a predecessor to CUORE-0, which reported the following limits on the decay half-lives [20]: 23 yr.

Experiment
The CUORE-0 tower, just like all the CUORE towers, contains 52 nat TeO 2 crystals (i.e., 130 Te is present at its natural abundance of 34.2% [21]). The crystals are arranged in a copper frame into 13 floors, with each floor containing four crystals [14]. The mass of each crystal is 750 g, for a total detector mass of TeO 2 of 39 kg or 10.8 kg of 130 Te.
The crystals are operated as thermal detectors at a working temperature of ∼ 10 mK to minimize the heat capacity. In CUORE-0, this temperature is obtained using the same dilution refrigerator infrastructure as CUORICINO [22]. The detector thermal link to the fridge is provided by the polytetrafluoroethylene supports which hold each crystal in the copper frames and by the golden wires used to carry the electrical signal.
The temperature of each crystal is continuously sensed by a neutron transmutation doped (NTD) thermistor [23] glued to the crystal surface. The thermistor converts the thermal signal to a voltage output which is digitized at an acquisition rate of 125 samples/s. Signal pulses corresponding to thermal events are identified through a software trigger, and a 5 s long window -1 s before and 4 s after the trigger -is selected for further analysis. The initial second (pre-trigger) is used to establish the baseline temperature of the crystal just prior to the event. The pulse amplitude is used to determine the deposited energy. In order to monitor and correct for changes in detector gain due to temperature drifts, a silicon resistor (heater) [24] is coupled to each crystal and is used to generate reference thermal signals every 300 s [25,26]. A time-coincidence analysis can be performed to search for events that involve multiple crystals simultaneously. To account for the time response of the detector we use a coincidence window of ±5 ms. As the measured event rate is approximately 1 mHz/crystal, the probability of accidental (i.e., causally unrelated) coincidences is extremely small ( 10 −5 ).
In order to reduce background due to environmental radioactivity, the tower and cryogenic infrastructure are surrounded by several layers of shielding, including an internal low-background Roman lead layer and an external antiradon box. The details of the CUORE-0 detector design, operation and performance are described in [14,16]. The dilution refrigerator, shielding, and other cryostat components are those from the CUORICINO experiment [22,27].

Analysis
In this work we exploit the multi-detector coincidence patterns expected to accompany ββ 0 + 1 decays to maximize our sensitivity. The emitted electrons and gamma rays can interact in multiple crystals, producing a variety of experimental scenarios or signatures depending on the number of detectors involved in the process (event multiplicity). We first identify all the possible signatures that can be detected and rank them  [10] by their expected sensitivity to 0νββ 0 + 1 or 2νββ 0 + 1 decays. Our determination of the most significant signatures to be used in the final analysis makes use of both real and simulated data.

Signature identification
To simplify the analysis we restrict ourselves to signatures in which the energy deposition from electrons is limited to the crystal where the decay took place. We further require that each individual de-excitation gamma is completely absorbed in a single crystal, thus discarding events where these gamma rays scatter but subsequently escape the detector. Scenarios where gamma rays are absorbed in the same crystal where they originated from are considered valid. Since the maximum number of gamma rays emitted in the decay is three ( Fig. 1), these choices imply that the maximum number of crystals involved in an event is four (one crystal contains the two electrons). We considered adding a constraint on the proximity of the involved crystals but found that it does not increase our sensitivity significantly.
For every possible signature satisfying these conditions, both for the 2νββ 0 + 1 and the 0νββ 0 + 1 decay modes, we identify a single monochromatic line in one of the up-to-four crystals which can be used to estimate the decay rate through the fit described in Sect. 3.5. For this reason, at least two crystals must be involved in 2νββ 0 + 1 decay signatures, so that at least one of them records a monochromatic peak -when only one crystal is involved the energy deposited is a continuous distribution spread over a 734 keV-wide range. The analysis therefore ultimately involves searching for the peak associated with each signature. Taking as an example the most likely de-excitation pattern (two gammas with energies of 1257.5 keV and 536.1 keV), the allowed scenarios resulting from this set of rules are reported in Table 1. Table 1 Allowed scenarios for the most likely de-excitation pattern (two gammas with energies of 1257.5 keV and 536.1 keV), resulting from the rules defined in Sect. 3.1. The decay always takes place in detector A. ββ, which represents the energy deposited by electrons, corresponds to 734 keV for the 0νββ 0 + 1 decay and to any energy in the 0 ÷ 734 keV range for the 2νββ 0 + 1 decay. An energy deposition of 1793.5 keV is recorded when both gammas deposit their energy in the same crystal. Scenarios where only one crystal is involved (Multiplicity 1) are valid only for 0νββ 0 +

Sensitivity evaluation
Considering all the possible gamma cascade patterns and the aforementioned constraints, a total of 57 signatures remain. However, only a few of them produce a sizable contribution to the half life sensitivity to double beta decay. In the case of a peak search in the presence of nonzero background the half-life sensitivity has the following dependence on the experimental parameters [28]: where is the total detection efficiency, M ·t is the exposure, b is the background rate per unit energy (background index) in the energy region of interest (ROI), ΔE is the energy resolution near the ROI and the constant k incorporates the isotopic abundance of 130 Te, the molecular mass of the source material and the required significance level. While the exposure (35.2 kg·y of TeO 2 ) and k are always the same, the other three parameters vary significantly depending on the signature and associated peak under consideration. We measure the energy resolution and its energy dependence directly from CUORE-0 data. We find the resolution exhibits a linear energy dependence which we parametrize as ΔE(E) = ΔE(2615) × ( p 0 + p 1 × E) [16]; here ΔE(2615) = 4.9 keV is the average FWHM of the 2615 keV line obtained during the calibration runs, p 0 = 0.49 ± 0.04 and p 1 = (2.22 ± 0.15) × 10 −4 . To avoid biasing our ranking procedure, the background index for each signature is estimated from the CUORE-0 background model described in [17] rather than using the CUORE-0 data directly. The model involves a full reconstruction of the measured spectra with highly detailed Monte Carlo simulations based on the Geant4 package [29]. The selection conditions (cuts) associated with each analysis signature are applied to the simulated data and the background index is evaluated in the relevant energy range. Finally, the efficiency term accounts for the probability that a ββ 0 + 1 event is triggered, produces the multi-detector coincidence signature in question, and the peak is properly reconstructed at the expected amplitude. The quantities involved in the efficiency calculation are described in detail in Sect. 3.4. Some 2νββ 0 + 1 scenarios have similar energy deposition patterns, requiring the same event multiplicity and energy of the associated monochromatic peak, but different energy ranges in the crystal where the decay took place. If there is an overlap between these ranges, we apply an additional procedure to properly evaluate the detection efficiency. For example scenarios #6 and #8 from Table 1 both search for a time coincidence between a peak at 1257.5 keV and an event in a wide energy range (0 ÷734 keV for #6, 536.1 ÷ 1270.1 keV for #8). The range between 536.1 and 734 keV is included in both scenarios. If P 6 and P 8 are the probabilities of detecting a decay from scenario #6 and #8 respectively, the probability of detecting a decay in the overlapping energy region is P 6 + P 8 . Any time we observe a situation like this, we modify the original scenarios to exclude the overlapping energy region, which we use instead to create an additional scenario. Referring again to the previous example, a first scenario will require events in the 0 ÷ 536.1 keV range, a second in the 536.1 ÷ 734 keV range and a third in the 734 ÷ 1270.1 range. For this reason, some of the scenarios we considered for our final 2νββ 0 + 1 analysis do not cover the whole 734 keV-wide range.
We define our total sensitivity to 0νββ 0 + 1 /2νββ 0 + 1 decay as the sum in quadrature of the sensitivities given by each signature. We consider a signature to be relevant if its contribution to the total sensitivity to the process exceeds 1%. Only five such scenarios are identified for both the 0νββ 0 + 1 and the 2νββ 0 + 1 channels, and they are ranked by their sensitivity in Tables 2 and 3 respectively. The selected scenarios cover ∼ 97% and ∼ 99% of the total sensitivity for the 2νββ 0 + 1 and 0νββ 0 + 1 cases respectively.

Data selection
We first remove time periods where the data quality are poor; the effect of this is accounted for in the exposure. We next remove events that are either poorly reconstructed by our analysis or are non-signal-like using the pulse shape methods described in [16]. We then impose cuts based on the deposited energy and on the event multiplicity (i.e., the num- Table 2 The five most relevant scenarios that contribute to the total sensitivity to the neutrinoless decay channel by more than 1%. b is the background index (units: counts/keV/kg/y) and the total detection efficiency. The electrons are assumed to be always fully absorbed in the detector 'A'. The numbers within parentheses show how the energy is split between electrons (734 keV) and gammas in detector 'A'. The final fit is performed at the highest energy monochromatic peak measured in each signature, marked here in bold. The observed spectra for each scenario are shown in Fig. 3. Scenario numbers are unrelated to those in Table 1 Scen. Energy (keV Table 3 The five most relevant scenarios that contribute to the total sensitivity to the two neutrino decay channel by more than 1%. b is the background index (units: counts/keV/kg/y) and the total detection efficiency. The electrons are assumed to be always fully absorbed in the detector 'A'. The final fit is performed at the highest energy monochromatic peak measured in each signature, marked here in bold. The energy ranges in detector 'A' do not always cover the full 734 keV-wide range to properly evaluate the efficiency (see Sect. 3.2). The observed spectra for each scenario are shown in Fig. 2. Scenario numbers are unrelated to those in Table 1 Scen. ber of crystals involved in the event). For the 0νββ 0 + 1 case (Table 2), the energy and multiplicity rules a candidate event must pass are listed below.
-In scenario 1, events must involve exactly three hits in the same coincidence time window. One of the three crystals must contain a signal with energy E 1 in the range (734±5σ 734 ) keV. The notation σ 734 indicates the energy resolution at 734 keV, which is estimated from the resolution function reported in Sect. 3.4. Another of the three crystals must have energy, E 2 , in the (536 ± 5σ 536 ) keV range. No requirement is imposed on the energy deposited in the third crystal. -In scenarios 2-5, events must have exactly two crystal hits. No requirement is set on the individual measured energies of the hits (E 1 and E 2 ), but rather on their sum, Labeling as E A and E B the energies indicated in Table 2 that are expected to be deposited in detector A and B respectively, and as E AB their sum (E A + E B ), then E tot must lie in the range (E AB ±5σ AB ), We apply a similar logic to the 2νββ 0 + 1 decay (Table 3), but the range defined by the continuous electron spectrum replaces the ±5σ requirement.
-In scenario 1 events must involve exactly three hits. The energy of the first hit, E 1 , is bound in the range 0 < E 1 < 734; the energy of the second hit, E 2 , must be in the range (536 ± 5σ 536 ). No requirement is set for the hit on the third crystal. -For scenarios 2-5, two hits are required in the same coincidence window. Again, we don't set a requirement on the individual hit energies, but on their sum, The continuous electron spectrum, contained in detector A (Table 3), defines two energy limits, E min A and E max A . We require that E tot satisfy the condition E min

Efficiency evaluation
The detection efficiency is the probability that a ββ 0 + 1 event is triggered and properly reconstructed, including the production of the required multi-detector coincidence signature. The detection efficiency is a product of several factors: -the probability that an event is triggered, which we estimate from the fraction of heater-induced events that are triggered and correctly reconstructed at the expected amplitude; -the probability to include only physical events whose pulse shape does not differ from the average behavior (pulse shape efficiency); -the probability to correctly measure the number of crystals involved in an event (i.e., the event multiplicity); -the fraction of ββ 0 + 1 events that deposit energy according to a particular scenario.
The trigger and pulse shape efficiencies are derived in the same way described in [16]. The multiplicity term is obtained differently depending on the number of crystals involved. The probability that an event involving a single crystal (multiplicity 1) is actually recorded as one is calculated using the 40 K line at 1461 keV [16]: since it's the only γ line emitted in the decay, the only way for it to be measured in a multiplicity > 1 event is by chance. For events with higher multiplicity we don't have a way to evaluate the efficiency directly, so we resort to a statistical derivation. We take into consideration the two main effects that can alter the multiplicity of one event: accidental coincidences and pile-up.
-Accidental coincidences are events that happen simultaneously by chance and not due to causally correlated signals. This leads to an artificially increased multiplicity. -Pile-up refers to the situation where two (or more) signals happen randomly on the same channel within the same 5 seconds-long signal window. Such events are marked as pile-up and ignored by the coincidence calculation. This artificially reduces the multiplicity; if, for example, two signals are coincident but one of them is removed due to pile-up, the second signal is recorded as a multiplicity one event.
The magnitude of these effects can be estimated with Poisson statistics, considering the measured average event rate of 1 mHz/channel, a 5 s long signal window and a 10 ms coincidence window. The resulting efficiencies for events with multiplicities 1 through 4 are listed in Table 4. While multiplicity 4 scenarios are possible, none of them gives a significant contribution to sensitivity. The final term in the efficiency, related to ββ 0 + 1 decay itself, is evaluated using dedicated Monte Carlo simulations. We use the same simulation software employed for the CUORE-0 background model to reproduce the effects of all decay patterns and compute the fraction of fully-contained ββ 0 + 1 events for each of the 57 scenarios. Due to the high statistics of these simulations, the error associated to this efficiency term is extremely small (< 0.1%).

Fitting technique
The final step of the analysis procedure is to obtain for each scenario the energy spectrum selected for the fit. We choose to fit the spectrum of the crystal that records the monochromatic peak with the highest energy as, due to the shape of our observed spectra, it usually has the lowest background.
To simplify the fits, ranges are chosen to exclude any peak from other γ lines. This selection is based on the CUORE-0 background model [17] rather than the data. The final spectra for the 0νββ 0 + 1 decay search are shown in Fig. 2, while those for 2νββ 0 + 1 decay are shown in Fig. 3. The strong background reduction achieved with energy-related cuts and the excellent agreement between real data and the background model are evident in these figures.
The fit function for each signature is the following: where B const and B lin are parameters describing a linear background, is the detection efficiency, t is the live time,  Tables 2  and 3. We perform a simultaneous unbinned extended maximum likelihood (UEML) fit on the five signatures belonging to each ββ 0 + 1 decay, using the RooFit fitting package [30]. The fit parameters are constrained as follows: -the decay rate Γ 0 + 1 ββ is a common parameter for all five signatures; -the two background components are independent for each signature; -the detection efficiency is fixed at the value reported in Tables 2 and 3 Tables 2 and 3 .
Each fit has 11 free parameters: two background parameters for each of the five signatures (for a total of 10) and the decay rate Γ 0 + 1 ββ , which is common for all signatures.

Results
We find no evidence of a ββ 0 + signal, either for the neutrinoless or the two neutrino decay mode. We set a 90% confidence upper limit on the decay rates for the two processes: Γ 0ν 0 + < 8.8 · 10 −25 y −1 , Γ 2ν 0 + < 2.8 · 10 −24 y −1 . In turn, these correspond to the following lower limits for the half lives:  in black are data with energy cuts. The non-shaded area represents the fit range, while the black arrow points to the location of the expected ββ 0 + 1 peak We estimate the systematic uncertainties with a procedure identical to that applied in the analysis for 0νββ decay to the ground state [15,16]. We consider two components for the systematic uncertainty: one factor which scales with the decay rate (σ scaling ), and one which is independent of the decay rate (σ add ). We generate a large number of simulated spectra with a distribution taken from the best fit of each signature, but with the value of a single nuisance parameter modified by 1σ ; we then fit the simulated spectra with the unmodified parameter. To probe the value of σ scaling , we repeat the analysis including in the simulated spectra a fake signal of variable strength. We regress the resulting best-fit decay rates against the simulated values to determine σ add and σ scaling . This procedure is applied separately to get the systematic contributions from the uncertainty on efficiency, energy resolution and peak position. We also run the simulation without changing any parameter, to check for a possible fit bias. The resulting values for σ add and σ scaling are reported in Table 5 for 0νββ 0 + 1 and in Table 6 for 2νββ 0 + 1 . In both cases, the dominant effect is a small negative bias.
We combine the CUORE-0 likelihood curves with those from CUORICINO [20] (Fig. 4 for 0νββ 0 + Thanks to the improved background and analysis techniques, we achieve similar reach as CUORICINO with less than half of the exposure. The lower limits for the half lives of both the 0νββ 0 + 1 and 2νββ 0 + 1 decays obtained by the combination of the results from CUORE-0 and CUORICINO are the best currently available. CUORE will achieve even higher sensitivity, thanks to the improved background and to the powerful coincidence analysis made possible by its closely-packed 988 crystals. ons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. Funded by SCOAP 3 .