1 Introduction↩︎

The discovery of neutrino oscillations [1], [2] stands as one of the landmark results of modern particle physics, establishing that neutrinos are massive and that lepton flavors mix. Yet the minimal extension of the Standard Model required to accommodate this fact leaves many fundamental questions unanswered: the absolute neutrino mass scale; the Dirac, Pseudo-Dirac, or Majorana nature of the mass; the mass ordering; and, most relevant to the present work, whether the leptonic sector violates CP symmetry. The CP-violating phase \(\delta_{\rm CP}\) enters the Pontecorvo–Maki–Nakagawa–Sakata (PMNS) mixing matrix and governs a matter–antimatter asymmetry in neutrino oscillations that next-generation long-baseline experiments, chief among them DUNE [3] and Hyper-K [4], are specifically designed to measure.

A reliable extraction of \(\delta_{\rm CP}\), however, rests on a tacit assumption that is rarely scrutinized in practice: CPT conservation. Long-baseline oscillation experiments are sensitive to the difference between neutrino and antineutrino oscillation probabilities, \(\Delta P \equiv P(\nu_\mu \to \nu_e) - P(\bar\nu_\mu \to \bar\nu_e)\), which is sourced both by genuine CP violation and by any CPT-violating asymmetry between the mass spectra and mixing parameters of neutrinos and antineutrinos. If CPT is violated at even a modest level, a CPT-conserving analysis will absorb the extra asymmetry into a biased inference of \(\delta_{\rm CP}\), potentially reporting a spurious signal of CP violation or obscuring a real one. This concern is not merely theoretical: it was shown in [5] that the hints of CP violation reported by T2K, a predecessor of Hyper-K, are not robust against a CPT-violating explanation, and that a CPT-violating, CP-conserving scenario is in perfect agreement with the neutrino oscillation data available at the time, even after combining T2K with NO\(\nu\)A and reactor experiments. Testing CPT conservation is therefore not merely an exotic pursuit; it is a necessary prerequisite for the program of discovery of leptonic CP violation.

Beyond its immediate practical consequences for oscillation physics, CPT symmetry occupies a singular position in the architecture of modern theoretical physics. CP violation, however striking, is accommodated within the Standard Model through the Cabibbo–Kobayashi–Maskawa and PMNS phases; its discovery refined the theory but did not overturn its foundations. CPT violation, by contrast, would constitute a far more radical departure. The CPT theorem [6]–[9] guarantees that any theory satisfying three pillars of quantum field theory—Lorentz invariance, locality, and the spin-statistics connection in a unitary, local, relativistic framework—must be invariant under the combined operation of charge conjugation C, parity inversion P, and time reversal T. The theorem does not require that any of C, P, or T be individually conserved; it demands only their product. A confirmed violation of CPT would therefore signal the breakdown of at least one of these foundational principles: it would imply either a departure from Lorentz invariance, a non-local interaction structure, or a failure of the standard spin-statistics relation. Any of these would necessitate physics beyond the framework of local relativistic quantum field theory—the mathematical language underlying every confirmed fundamental theory we possess.

Proposals for CPT violation arise naturally in several beyond-Standard-Model contexts. In string-field theory and non-commutative geometry, spontaneous breaking of Lorentz symmetry generates CPT-odd operators in the low-energy effective Lagrangian [10], [11]. The Standard-Model Extension (SME) of Kostelecký and collaborators provides a systematic effective-field-theory catalog of all Lorentz- and CPT-violating operators consistent with the remaining symmetries of the Standard Model [12]–[14]. A qualitatively different and particularly minimal realization was proposed in [15]: CPT violation in the Dirac mass terms of the three neutrino flavors, which generates independent mass spectra and mixing parameters for neutrinos and antineutrinos while preserving Lorentz invariance. This scenario, strongly motivated by braneworld models with extra dimensions where neutrinos propagate into the bulk, parametrizes the CPT-violating observable most directly as \(\delta\Delta m^2_{31} \equiv \Delta\bar{m}^2_{31} - \Delta m^2_{31} \neq 0\), and is the framework adopted throughout this work. Quantum-gravity scenarios involving space-time foam or decoherence induced by virtual black holes also generically predict CPT non-invariance [16], [17]. The neutrino sector thus provides a low-energy window onto Planck-scale physics.

Neutrinos are, in several respects, the most powerful laboratory for CPT tests among elementary particles. First, unlike the kaon or \(B\)-meson systems, where CPT tests rely on particle–antiparticle mass differences, oscillation experiments are directly sensitive to the interference between mass eigenstates, giving access to the full matrix of mass-squared differences and mixing angles for both neutrinos and antineutrinos separately. Second, the origin of neutrino mass remains unknown. If neutrinos are Majorana particles—their own antiparticles—the usual notion of a CPT-conjugate state requires re-examination; the relevant CPT-odd quantity becomes a difference in the effective Majorana masses or phases. If they are Dirac, neutrino and antineutrino mass matrices are in principle independent, and an asymmetry \(\delta\Delta m^2_{31}\) is a direct, gauge-invariant observable. In either case, the lack of a confirmed mass-generation mechanism means that new physics in the neutrino mass sector is not only possible but perhaps expected, and CPT-violating mass terms are among the most minimal such extensions. Third, cosmological and astrophysical neutrinos traverse macroscopic and even Hubble-scale distances, providing extraordinary sensitivity to Lorentz- and CPT-violating effects that accumulate coherently with baseline [18], [19]. Fourth, the atmospheric neutrino flux delivers a broad-band, isotropic sample of both neutrinos and antineutrinos spanning energies from sub-GeV to hundreds of TeV, and experiments such as IceCube-DeepCore and KM3NeT/ORCA are able to exploit this diversity of \(L/E\) to disentangle CPT-violating mass asymmetries from matter effects and standard-oscillation parameters. In fact, the next-generation atmospheric neutrino experiments will perform some of the most precise measurements of the neutrino oscillation parameters [20], including the CP-violating phase [21]–[24].

In this work, we study the interplay between CPT violation and CP violation in long-baseline neutrino oscillation experiments. We demonstrate, analytically and numerically, that CPT-violating mass-squared asymmetries generate a contribution to \(\Delta P\) that is functionally degenerate with a shift in \(\delta_{\rm CP}\), inducing an effective, energy-dependent phase \(\phi_{\rm eff}(E)\) that causes a CPT-conserving analysis to infer a biased CP phase. Using the DUNE experiment as a concrete example, we show that current experimental bounds on \(\delta\Delta m^2_{31}\) are already sufficient to produce significant biases in the reconstructed \(\delta_{\rm CP}\) and, in CP-conserving scenarios, to fake a near-maximal CP-violation signal. We then show how atmospheric neutrino telescopes—specifically IceCube-DeepCore and KM3NeT/ORCA—can break this degeneracy by providing independent, high-statistics constraints on \(\delta\Delta m^2_{31}\), and we project the combined sensitivity of present and future configurations as a function of time.

The rest of this article is structured as follows. Section 2 develops the analytic framework, decomposing \(\Delta P\) into CP, CPT, and matter contributions and deriving the effective phase shift \(\phi_{\rm eff}(E)\) that quantifies the CP–CPT degeneracy. Section 3 presents the DUNE sensitivity analysis using GLoBES, illustrating the bias on \(\delta_{\rm CP}\) as a function of \(\delta\Delta m^2_{31}\), describes the atmospheric neutrino analysis with the Pynu framework [20], and presents current constraints and future projections from IceCube-DeepCore and ORCA. We conclude in Section 4.

2 Oscillations Phenomenology↩︎

2.1 Decomposition of \(\Delta P\)↩︎

The \(\nu_e\) and \(\bar\nu_e\) appearance probabilities, in the presence of a mass-squared-splitting asymmetry \(\delta \Delta \equiv \Delta_{\nu\bar\nu}(\Delta m^2) = \bar{\Delta} m^2 - \Delta m^2\), are given by the Cervera formula [25] as \[\begin{align} P(\nu_\mu \rightarrow \nu_e) = \alpha^2 + \beta^2 + 2\alpha\beta \cos(\Delta_{31} + \delta_{\textrm{CP}}), \\ P(\bar\nu_\mu \rightarrow \bar\nu_e) = \bar\alpha^2 + \beta^2 + 2\bar\alpha\beta \cos(\bar\Delta_{31} - \delta_{\textrm{CP}}). \end{align}\] Here, \(\Delta_{ij} \equiv \Delta m^2_{ij} L / (4E)\) is the standard oscillation phase. Following similar conventions as used in [26], the amplitudes \(\alpha\), \(\bar\alpha\), and \(\beta\) are defined by \[\begin{align} \alpha = \sqrt{P_{\rm{atm}}} &= s_{23}\sin2 \theta_{13}\frac{\sin((1-\hat{A})\Delta_{31})}{1-\hat{A}} ,\\ \bar \alpha = \sqrt{\bar P^{\cancel{\rm{CPT}}}_{\rm{atm}}} &= s_{23}\sin2 \theta_{13}\frac{\sin((1+ \hat{\bar A})\bar\Delta_{31})}{1+\hat{ \bar A}} ,\\ \beta = \sqrt{P_{\rm{sol}}} &= c_{23}\sin2 \theta_{12}\frac{\sin(\hat{A} \Delta_{21})}{\hat{A}}. \end{align}\] where \(A\) is the matter potential, and we have introduced the dimensionless ratios \(\hat{A} \equiv A/\Delta\) and \(\hat{\bar A} \equiv A/\bar\Delta\), with \(\Delta\) (\(\bar\Delta\)) the mass-squared splitting appropriate to the channel.

Here, the parameters that are related to CPT violation include \(\bar \alpha_{\rm{CPT}}\), \(\hat{\bar A}\), and \(\bar \Delta\). Expanding \(\cos(\Delta \pm \delta_{\rm{CP}})\), we can obtain the difference in neutrino and antineutrino oscillation probabilities \[\label{eq:decomp} \begin{align} \Delta P = &(\alpha^2 - \bar \alpha^2) \\ &+ 2\beta \cos \delta_{\rm{CP}}(\alpha \cos \Delta - \bar\alpha \cos \bar \Delta ) \\ &- 2\beta \sin \delta_{\rm{CP}} (\alpha \sin \Delta + \bar \alpha \sin \bar \Delta ). \end{align}\tag{1}\]

Observe that the second term is CP-even, and the third term is CP-odd: in vacuum, without CPT violation, and with \(\alpha = \bar \alpha\), this term is proportional to the Jarlskog invariant. However, due to interference terms, it is hard to disentangle the effects of CP, CPT, and the matter potential. To disentangle these contributions, denoting the difference in atmospheric oscillation probability difference with and without CPT conservation as \(\bar \alpha_{\rm{CPT}} - \bar \alpha_{\cancel{\rm{CPT}}}\), we can define the CPT contribution \(\Delta P_{\rm{CPT}}\) as \[\begin{align} \Delta P_{\rm{CPT}} = &(\bar \alpha_{\rm{CPT}} - \bar \alpha_{\cancel{\rm{CPT}}}) \\ + &(\sin \delta_{\rm{CP}} + \cos \delta_{\rm{CP}})\\ \times&(2 (\delta \bar \alpha) \beta \cos \bar \Delta + 2 \bar \alpha \beta (\cos \Delta - \cos \bar \Delta)). \end{align}\]

One identifies the matter-CPT interference term as well as the matter-CP-CPT interference terms, and verifies that \(\Delta P_{\rm{CPT}} = 0\) with \(\bar \Delta = \Delta\) and \(\delta \bar \alpha = 0\). With this parameterization, we obtain back the decomposition from Ref. [26] \[\Delta P = \Delta P_{\rm{CP, matter}} + \Delta P_{\rm{CPT}}.\]

For simplicity, we can define \(\Delta P_{\rm{CP}} = \Delta P_{\rm{CP}}^{\rm{vacuum}}\) to be the vacuum genuine CP violation, and let \(\Delta P_{\rm{matter}} = \Delta P - \Delta P_{\rm{CP}}\) absorb the leftover probability difference in CPT-conserving scenario. Another choice, as done in Ref. [26], would be to disentangle matter and CP contributions by defining \(\Delta P_{\rm{matter}} = \Delta P (\delta_{\rm{CP}} = \pi /2) - \Delta P (\delta_{\rm{CP}} = 0)\).

We use the nuSQuIDS package [27] in a constant-matter-density configuration, taking the DUNE baseline and matter density as our reference accelerator setup, and choose a CPT-violating mass-squared splitting of \(\delta \Delta = 1.0 \times 10^{-3}~\text{eV}^2\). Figure 1 (a) shows the resulting decomposition of the total \(\Delta P\) into its three components. Panels (b) and (c) further demonstrate that, at a true \(\delta_{\rm{CP}} = -112^\circ\), the CP and CPT contributions to \(\Delta P\) conspire to produce a degeneracy: an impostor CPT-conserving solution at \(\delta_{\rm{CP}} = -84^\circ\), \(\delta\Delta = 0\) that reproduces the true \(\Delta P\).

This example of an impostor solution leads us to understand qualitatively first the physical implications of the decomposition of \(\Delta P\) in Eq. 1 : CPT violation can effectively mimic CP violation in long-baseline experiments. Two features are crucial:

1. CPT-violating terms interfere with CP-sensitive contributions, producing an energy-dependent effect functionally similar to genuine CP violation in vacuum.

2. Matter effects further modify the interference, enhancing the degeneracy between CP and CPT effects.

This interplay allows for situations in which a true CP phase with small CPT violation produces the same \(\Delta P\) as a CPT-conserving scenario with a shifted CP phase.

Additional considerations to break or mitigate this degeneracy include:

• Energy dependence: The mimicking effect varies with neutrino energy, so combining spectral information across the energy range can help distinguish CP and CPT contributions.

• Channel complementarity: Using multiple oscillation channels, such as \(\nu_\mu\) disappearance or \(\nu_e\) appearance at different baselines, provides independent constraints on CPT-violating parameters.

• Baseline dependence: Different experiments experience different matter effects, so combining long-baseline setups improves sensitivity to genuine CP violation.

In summary, an observed CP-odd asymmetry in \(\nu_e\) appearance does not necessarily indicate genuine CP violation in vacuum. Only by combining multiple channels, baselines, and energy spectra, or by imposing precise constraints on CPT-violating parameters, can this degeneracy be resolved. The decomposition of \(\Delta P\) into CP, matter, and CPT contributions provides a clear framework to quantify the extent to which CPT can mimic CP violation.

2.2 Effective CP Phase Shift from CPT Violation↩︎

Now, we turn to a quantitative study of this degeneracy by further parameterizing the decomposition in Eq. 1 through a coordinate transformation. We start from the CP-dependent part of the probability difference, \[\Delta P_{\rm int}(E) = A(E)\cos\delta_{\rm CP} + B(E)\sin\delta_{\rm CP},\] where \[\begin{align} A(E) &= 2\beta \left(\alpha \cos\Delta - \bar\alpha\cos\bar\Delta \right), \\ B(E) &= -2\beta \left(\alpha \sin\Delta + \bar\alpha\sin\bar\Delta \right), \end{align}\] and recast this expression in the form \[\Delta P_{\rm int}(E) = R(E)\,\sin\big(\delta_{\rm CP} + \phi_{\rm eff}(E)\big),\] with amplitude \(R(E) = \sqrt{A^2(E) + B^2(E)}\). The effective phase shift \(\phi_{\rm eff}(E)\) is defined through \[\tan\phi_{\rm eff}(E) = \frac{A(E)}{B(E)},\] which yields \[\phi_{\rm eff}(E) = \arctan\!\left[ \frac{\alpha \cos\Delta - \bar\alpha\cos\bar\Delta}{-\alpha \sin\Delta - \bar\alpha\sin\bar\Delta} \right]. \label{eq:phieff}\tag{2}\]

This form makes explicit that CPT violation modifies the CP-phase dependence by introducing an effective, energy-dependent phase shift, where we show an illustrative example in Fig. 2. See Section 3.3 for a description of the energy range in which these effects are most significant.

We can further make this more transparent by considering an expansion of the CPT-violating quantities and showing the first-order contributions: we can denote \[\bar\Delta = \Delta + \delta\Delta, \qquad \bar\alpha_{\rm CPT} = \alpha + \delta\bar\alpha,\]

To first order,

\[\begin{align} \phi_{\rm eff}(E) \simeq -\frac{1}{2}\delta \Delta_{31}(E) &+ \frac{1}{2}\frac{\delta\bar\alpha(E)}{\alpha(E)}\,\cot\Delta_{31}(E). \end{align}\]

This result provides a direct physical interpretation of the CP–CPT degeneracy: CPT violation effectively induces a shift \[\delta_{\rm CP} \;\longrightarrow\; \delta_{\rm CP} + \phi_{\rm eff}(E),\] with \(\phi_{\rm eff}(E)\) depending nontrivially on energy through both the oscillation phase and matter effects. The term proportional to \(\delta\bar\alpha/\alpha\) encodes amplitude-level CPT violation and matter–CPT interference. The term proportional to \(\delta\Delta\) represents a genuine phase shift arising from differences in the mass-squared splittings.

Because \(\phi_{\rm eff}(E)\) is energy dependent, a CPT-violating scenario can closely reproduce the oscillation probability of a CPT-conserving theory with a shifted \(\delta_{\rm CP}\) at a fixed energy, while deviations emerge once spectral information is taken into account. Supposing that only monochromatic measurements of the CP phase are taken at two energies \(E_1\) and \(E_2\), then at first order,

\[\Delta \delta_{\rm CP}^{\rm exp} \simeq \phi_{\rm eff}(E_2) - \phi_{\rm eff}(E_1).\] Near the first oscillation maximum (\(\Delta \sim \pi/2\)), this reduces to \[\Delta \delta_{\rm CP}^{\rm exp} \sim - \frac{1}{2}\delta\Delta_{31}(E)\]

On the one hand, percent-level shifts in the inferred CP phase can arise from \(\mathcal{O}(10^{-3}\,\text{eV}^2)\) differences in the atmospheric mass-squared splitting or comparable fractional changes in the oscillation amplitude (while such effects are constrained by existing oscillation data, this estimate demonstrates that CPT violation at a level not yet excluded could induce experiment-dependent shifts in \(\delta_{\rm CP}\)); and on the other hand, \(\phi_{\rm{eff}}(E)\) being an energy-dependent quantity suggests that different long-baseline experiments probing distinct \(L/E\) regimes may infer different effective values of \(\delta_{\rm CP}\).

It is worth noting here that the expansion is only valid when \(|\delta\Delta| \ll 1\) and \(|\delta\bar\alpha| \ll \alpha\); however, the expansion is qualitatively valid in showing the energy dependence of the \(\phi_{\rm{eff}}(E)\) term, and the full prescription (Eq. 2 ) is always valid.

3 Accelerators and Atmospheric Neutrino Telescopes↩︎

3.1 Accounting for tension in NOvA-T2K↩︎

The T2K and NOvA collaborations have performed a joint fit on neutrino oscillation parameters, and presented a mild tension in their measurements of the CP-violating phase [28]. Our parameterization of the CPT violation in oscillations provides a natural mechanism and explanation to account for the difference in best-fit points \(\delta_{\rm CP}\) between T2K and NO\(\nu\)A. Here, we provide viewpoints from two directions as shown in Fig. 3, where we show that the true \(\Delta P_{\nu \bar \nu}\) is equivalence classes forming submanifolds in the \(\delta_{\rm{CP}}^{\rm{true}}\) and \(\delta\Delta^{\rm{true}}\) plane. The physical observable is \[\langle\Delta P\rangle = \langle P(\nu_\mu\!\to\!\nu_e)\rangle_{\Phi_\nu} - \langle P(\bar\nu_\mu\!\to\!\bar\nu_e)\rangle_{\Phi_{\bar\nu}}.\] For simplicity, we also assume that the flux is Gaussian distributed around the characteristic flux of the two experiments with spread mimicking the on- and off-axis effects: T2K assumes \(\mu_E = 0.6\) GeV with \(\sigma_E = 0.1\) GeV, and NOvA assumes \(\mu_E = 1.9\) GeV with \(\sigma_E = 0.45\) GeV.

On the one hand, we assume a CPT-violating but CP-conserving truth, and look at the equivalent points in the submanifold where the physical observable \(\Delta P\) is invariant, and we find that they cross the CPT-conserving axis at different points due to the different baselines. The key observation is that any point along this submanifold is perfectly equivalent to any other point, so there’s nothing special about these labeled CPT-conserving solutions. On the other hand, we ask the preimage question: given the best-fit points of the NOvA and T2K experiments, respectively, where the \(\Delta P\) is computed with the CPT-conserving hypothesis, what are the equivalence manifolds in the \((\delta \Delta, \, \delta_{\rm{CP}})\) plane of the two experiments? We can even further trace the locus on these points to find a best tension-reconciling CPT-non-conserving truth. One can see that the normalized difference in \(\Delta P\) across the flux range reduces by 75% if we allow for a CPT-non-conserving truth. The mere existence of the distinct intersection points constitutes a CPT-based explanation for the T2K–NOvA tension that is invisible to any CPT-conserving analysis, confirmed by the proximity of the two preimage contours.

The two experiments probe the same \(L/E\) regime, but with neutrino beams that span distinct baselines and peak at different energies, implying a varying impact of matter effects in each case. This leads to separate degeneracy manifolds in the \((\delta_{\rm CP}, \delta \Delta m^2_{31})\) plane, whose CPT-conserving intersections are generically displaced from one another. The degeneracy is exact at the level of the flux-averaged asymmetry \(\langle\Delta P\rangle\) but approximate in the spectral sense, since \(\phi_{\rm eff}(E)\) varies across the flux window; this is the handle that future high-statistics experiments can exploit to break it.

Even with the flux information \(\phi(E)\) and reconstruction, a single accelerator experiment with a fixed matter effect profile and some energy spectral information typically still suffers from the CP-CPT degeneracy if a fit for oscillation parameters is performed in both horn currents combined. Supp. Fig. [fig:nova_combines] in the Appendix shows the CP phase measurements for NOvA following the treatment discussed here. In the next subsection, we also show a more detailed discussion for DUNE, where we further include flux information and detector energy reconstruction effects. We find in these plots that the measured \(\delta_{\rm{CP}}\) values follow the \(\Delta P\) isocontours approximately, validating its usage as an effective observable at truth level.

This result motivates the central question of the remainder of the paper: can an independent, external constraint on \(\delta \Delta m^2_{31}\) from atmospheric neutrinos break this degeneracy before DUNE begins its CP measurement?

3.2 CP-CPT degeneracy in DUNE↩︎

We first present the problem of the degeneracy between CP and CPT in accelerators, using the DUNE experiment as an example. We use the GLoBES package with the DUNE configuration [29], [30]. In this fit, we assume 7.5 years of FHC (Forward Horn Current) and RHC (Reverse Horn Current) equally, but we assume that only a combined fit will be performed, assuming that CPT is conserved, where \(\delta\Delta\equiv \Delta\bar{m}^2_{31} - \Delta m^2_{31} = 0\).

CPT violation in the context of DUNE has a significant impact on the measurement not only due to the combination of FHC and RHC in analysis \(R_{\nu_e} = R_{\nu_e}^{\rm{FHC}} + R_{\nu_e}^{\rm{RHC}}\), but also due to the wrong sign contamination of the beam. This is because, taking \(\nu_e\) appearance rate in RHC as an example, the event rate is computed by the sum of \(\nu_\mu\) and \(\bar \nu_\mu\) components \[\begin{align} R^{\rm{RHC}}_{\nu_e, \bar\nu_e} &= R^{\rm{RHC, wrong}}_{\nu_\mu \rightarrow \nu_e} + R^{\rm{RHC}}_{\bar \nu_\mu \rightarrow \bar \nu_e}, \\ &=\phi^{\rm{RHC}}_{\nu_\mu} (E) \cdot P_{\mu e} (E; \Delta_{31}) \cdot \sigma ^{\rm{CC}}_{\nu_e}(E) \cdot \epsilon_{\rm{det}}(E) \\ &+\phi^{\rm{RHC}}_{\bar \nu_\mu} (E) \cdot \bar P_{\bar \mu \bar e} (E; \bar \Delta_{31}) \cdot \sigma ^{\rm{CC}}_{\bar \nu_e}(E) \cdot \epsilon_{\rm{det}}(E), \end{align} \label{eq:DUNE-rate}\tag{3}\] where the first contribution in flux alone is \(15\%\) [31], which can correspond to \(35\%\) in event rate due to oscillation and cross section. With CPT violation, the neutrino and antineutrino oscillation probabilities are different, leading to different distortions of the event distribution that hide in plain sight and can be mistaken for a CP-violation effect.

In GLoBES, the calculation of event rates and therefore likelihood for a set of oscillation parameters is obtained via the usage of “rules": an analysis-level reconstructed spectrum. Conventionally, such as shipped in the DUNE configuration [29], [30], there are four such spectra: FHC appearance, RHC appearance, FHC disappearance, and RHC disappearance. However, each of these spectra, as shown in 3 , contains both the right and the wrong sign contributions, and assigning a single \(\Delta m^2_{31}\) value cannot account for the effect of CPT violation. We therefore use an equivalent 8-rule configuration in which each physical rule is split into a neutrino-only and an antineutrino-only sub-rule. The CPT-violating Asimov prediction is generated by evaluating the neutrino sub-rules with \(\Delta m^2_{31}\) and the antineutrino sub-rules with \(\Delta\bar m^2_{31}\), thereby realizing Eq. 3 at the rate level. It is worth noting that the 8 sub-rules are not treated as eight independently observed samples. After the neutrino and antineutrino components are generated with their respective mass splittings, the corresponding sub-rules are summed back into the 4 physical DUNE spectra in order to enforce shared signal and background normalization nuisance parameters for each physical spectrum.

The right two panels of Fig. 4 show two panels assuming a near-maximal CP violation and zero CP violation Asimov truths, respectively. Varying the Asimov data set true \(\delta \Delta\) on a grid, \(\delta_{CP}^{\textrm{meas}}\) is scanned over \([-\pi, \pi]\) while profiling over \(\theta_{23}\), \(\Delta m^2_{31}\), and \(\theta_{13}\). We show that the band of fitted \(\delta_{CP}^{\textrm{meas}}\) assuming CPT conservation shifts from the true value. In the case of existing CP-violation with \(\delta_{\rm{CP}} = -112 °\), assuming \(\delta \Delta \leq -0.3 \times 10^{-3}\) eV\(^2\), which is well below current bounds of \(|\delta \Delta| \leq 0.8 \times 10^{-3}\) eV\(^2\) [32], the allowed value of \(\delta_{\rm{CP}}\) incorrectly assuming CPT conservation completely misses the true value of CP violation. On the other hand, if the nature were CP-conserving, a CPT violation of \(\delta \Delta \geq 0.8 \times 10^{-3}\) eV\(^2\) could lead to double-minima structure in the \(\chi^2\) profile of the \(\delta_{\rm{CP}}\) fit when the true \(\delta_{\rm{CP}}\) is close to 0, leading to a potential false detection of near-maximal CP violation.

This double minima and basin switch behavior can be understood intuitively as two competing fit minima. For neutrinos and antineutrinos respectively, the \(\delta_{\rm{CP}}\) term enters in the interference probability terms in \(+\) and \(-\) signs, i.e., \[\begin{align} \Delta P ^{\rm{int}}_{\nu}(\delta_{\rm{CP}}) &= 2\alpha \beta \cos(\Delta_{31}+\delta_{\rm{CP}}),\\ \Delta P ^{\rm{int}}_{\bar \nu}(\delta_{\rm{CP}}) &= 2\bar\alpha \bar\beta \cos(\bar\Delta_{31}-\delta_{\rm{CP}}). \end{align}\] Therefore, when we shift \(\delta_{\rm{CP}}\rightarrow \delta_{\rm{CP}}+ \pi\), we flip the sign of both interference terms: \[\Delta P ^{\rm{int}}_{\nu, \bar \nu}(\delta_{\rm{CP}}+ \pi) = -\Delta P ^{\rm{int}}_{\nu, \bar \nu}(\delta_{\rm{CP}}).\] This creates the two competing minima near \(\delta_{\rm{CP}}\approx 0 °\). At small \(\delta\Delta\), the fit stays near the true \(\delta_{\rm{CP}}\) and absorbs the \(\bar \nu\) rate mismatch induced by CPT violation through shifts in \(\theta_{23}^{\rm{meas}}\) and \(\Delta_{31}^{\rm{meas}}\) in the non-interference \(\alpha^2 + \beta^2\) terms, and the interference terms \(\Delta P_{\rm{int}}\) remains close to their true value. However, at larger \(\delta\Delta\), such a spectral difference in the \(\bar \nu\) channel cannot be absorbed as well by the non-interference terms, and as a result, a competing fit method that better accommodates the enhanced \(\bar\nu_e\) appearance rate wins when \(\delta\Delta > 0\) at the cost of worsening the \(\nu_e\) appearance fit. In the bottom right panel of Fig. 5, where the near-truth fit is around \(\delta_{\rm{CP}}^{\rm{meas}} \approx -31 °\) and flips to \(\delta_{\rm{CP}}^{\rm{meas}} \approx +169 °\). One can further observe from the top panel that such basin flip occurs at \(\delta\Delta>0\) for \(\delta_{\rm{CP}}\approx 0 °\) and at \(\delta\Delta<0\) for \(\delta_{\rm{CP}}\approx 180 °\). Furthermore, one can see that in each region across the discontinuity, the \(\delta_{\rm{CP}}^{\rm{meas}}\) roughly follows the \(\Delta P\) isocontours for analytic calculation. Specific event distributions illustrating the two competing fit strategies as favoring the \(\nu_e\) and \(\bar \nu_e\) channels can be found in the Appendix at Supp. Fig. [sppfig:dune_evt_dist].

The sensitivity of DUNE to neutrino oscillation parameters can also be enhanced through measurements of the atmospheric neutrino flux [22], [33], [34]. In contrast to the beam-based analysis discussed above, atmospheric neutrinos probe a broad range of baselines and energies, with the strongest sensitivity to CP-violating effects arising in the sub-GeV region [35]. We investigate how CPT violation modifies the extraction of the CP phase in this complementary channel and find that, when CPT is violated, the value of \(\delta_{\textrm{CP}}\) inferred under the assumption of CPT conservation can deviate significantly from the true value, as illustrated in Fig. 4. Moreover, we observe that the correlation between \(\delta_{CP}^{\rm meas}\) and \(\delta(\Delta m^2_{31})\) differs qualitatively from that obtained in the beam analysis, reflecting the distinct energy dependence of CPT-violating effects and the dominance of lower-energy events in the atmospheric sample (see Section 2).

3.3 Uncovering the impostor CP with IceCube and ORCA↩︎

To extract meaningful \(\delta_{\rm CP}\) measurements from long-baseline accelerator neutrino experiments without fitting for oscillation parameters in FHC and RHC separately, we turn to atmospheric neutrino telescopes to place an independent constraint on \(\delta\Delta m^2_{31}\). Atmospheric neutrinos offer a qualitative advantage that long-baseline beams cannot replicate: a natural decoupling between the CP-sensitive and CPT-sensitive observables. The CP-violating phase \(\delta_{\rm CP}\) enters the oscillation probability through the interference between the solar and atmospheric amplitudes, a term suppressed by \(\alpha_{21} = \Delta m^2_{21}/\Delta m^2_{31}\) that is only significant below \(\sim\!1~\text{GeV}\), where the solar oscillation length is commensurate with trans-Earth baselines. Above this threshold the rapid \(\Delta m^2_{21}\)-driven oscillations average out over any realistic energy resolution, effectively decoupling \(\delta_{\rm CP}\) from the multi-GeV sample. The CPT-violating asymmetry \(\delta\Delta m^2_{31}\), by contrast, acts through the atmospheric channel and is most prominent above \(1~\text{GeV}\), where it shifts the \(\bar\nu_\mu\) disappearance minimum relative to the \(\nu_\mu\) minimum, generating a characteristic dipole pattern in the \((E,\cos\theta_z)\) plane (Supp. Figs. 3–6). The two effects therefore, occupy largely distinct regions of the observable phase space, allowing a combined fit over the full energy range to constrain \(\delta\Delta m^2_{31}\) with little contamination from the unknown \(\delta_{\rm CP}\). An additional handle comes from the MSW resonance: in normal ordering the resonance enhances neutrino conversion at \(E\sim 3\)–\(6~\text{GeV}\), while CPT violation shifts the resonance condition independently in the antineutrino sector, producing a matter-asymmetric signature that is qualitatively distinct from any CP effect. The sensitivity to \(\delta\Delta m^2_{31}\) is therefore concentrated in the antineutrino-driven disappearance sample for normal ordering and in the neutrino sample for inverted ordering, and the inclusive atmospheric flux — containing both neutrinos and antineutrinos across a wide range of energies and baselines — exploits both simultaneously in a single combined fit.

This decoupling is made explicit in Fig. 6, which shows the derivatives \(\partial P_{\alpha\beta}/\partial\Delta m^2_{31}\) and \(\partial P_{\alpha\beta}/\partial\delta_{\rm CP}\) in the \((E,\cos\theta_z)\) plane for the \(\nu_\mu\) disappearance and \(\nu_\mu\to\nu_e\) appearance channels, for both neutrinos and antineutrinos under normal and inverted ordering, using the CHIC code [36], [37]. The \(\partial P_{\mu\mu}/\partial\Delta m^2_{31}\) panel is dominated by a large localized structure at the MSW resonance (\(E\sim6\)–\(8\) GeV, \(\cos\theta_z\lesssim-0.8\)) in normal ordering, while the \(\partial P_{\mu\mu}/\partial\delta_{\rm CP}\) panel is nearly structureless across the same region. The two parameters are therefore, effectively orthogonal in the observable space, most accessible to IceCube and ORCA, providing a clean experimental handle on \(\delta\Delta m^2_{31}\) that is not contaminated by the unknown value of \(\delta_{\rm CP}\).

Using the Pynu framework, introduced in previous sensitivity studies of atmospheric neutrinos to oscillation parameters [27], [38], we perform data fits based on existing data release from the IceCube (IceCube-DeepCore) and ORCA (ORCA-6) experiments [39]–[41].

For the data fit, we are using 8 years of IceCube DeepCore data, and 433 kt-yr of ORCA-6 data; for the sensitivity projection, we are assuming 10 years of 7-string IceCube Upgrade [42], five of which have been recently deployed, along with 5 years of the full ORCA detector, whose simulation we have developed in previous work [27]. In both cases, the combined results profile over \(\sin^2\theta_{23}\) and \(\delta_{\rm{CP}}\) in addition to scanning over \(\Delta m^2_{31}\) and \(\delta \Delta\). The systematic and statistical treatments are partially inherited from previous work [27] but with important changes to account for the changes in data releases, which includes hyperplane systematics and Monte Carlo statistics uncertainties. We do not include DUNE or Hyper-Kamiokande in this analysis, as both detectors are expected to collect limited statistics in the multi-GeV region, where the atmospheric neutrino flux is relatively low.

The Monte Carlo datasets are collections of simulated events \({e}\), each carrying a weight \(w_e\). For the ORCA release each event also ships with an MC statistical variance \(\sigma_e^2\); the IceCube DeepCore release does not provide one, so we adopt the Poisson estimator \(\sigma_e^2 = w_e^2\). The MC variance in analysis bin \(i\) is then \(\sigma_{{\rm MC},i}^2 = \sum_{e\in i}\sigma_e^2\). To mimic what is done in the latest ORCA-6 analysis, our test statistics follow the Barlow–Beeston light prescription [43], [44]: bin \(i\) is assigned an auxiliary scale factor \(\beta_i\) that absorbs its MC statistical uncertainty on the predicted yield. Writing \(O_i\) for the observed count, \(E_i(\boldsymbol{\eta})\) for the total expectation (neutrino plus muon background), and \(\sigma^2_{{\rm MC},i}\) for the corresponding MC variance, the test statistic minimised at each grid point reads

with \(\tau_i \equiv \sigma^2_{{\rm MC},i}/E_i^2(\boldsymbol{\eta})\). The first sum is the saturated Poisson term used in our previous analysis, the second is the penalty term for the departure of \(\beta_i\) from unity, and the third are the penalty terms for Gaussian-modeled nuisance parameters; the mean and spread of the nuisance parameters are given in Supp. Table 1.

The optimal \(\beta_i\) at each bin minimises Eq. 4 in closed form: it is the positive root of a quadratic, \[\label{eq:beta-closed-form} \beta_i \;=\; \tfrac{1}{2}\!\left(-b_i + \sqrt{b_i^{\,2}-4 c_i}\right),\tag{5}\] where \(b_i = E_i\tau_i - 1\) and \(c_i = -O_i\tau_i\). With this minimization \(\beta_i\) is eliminated analytically at every evaluation of \(\chi^2\) rather than minimised numerically.

The total per-bin expectation includes an additive atmospheric-muon background scaled by a per-experiment normalisation nuisance \(\alpha_\mu^{\rm exp}\): \[\label{eq:muon-add} \begin{align} E_i^{\rm exp}(\boldsymbol{\eta}) \; &=\; E_i^{\nu,\,{\rm exp}}(\boldsymbol{\eta}) \;+\; \alpha_\mu^{\rm exp}\, M_i^{\rm exp}, \\ \sigma^2_{{\rm MC},i} \; &=\; \sigma^2_{\nu,i} + (\alpha_\mu^{\rm exp})^2\,\sigma^2_{\mu,i}. \end{align}\tag{6}\]

Using this statistical treatment, we have obtained the confidence regions, shown in Fig. 5. With current data from IceCube and ORCA, we can constrain \(| \delta\Delta m^2_{31}|\leq 0.57\times10^{-3}~\text{eV}^2\) (see Suppl. Fig. 9 for 1-dimensional \(\chi^2\) profile). However, this is with IceCube DeepCore and ORCA6, where these two low-energy (\(\mathcal{O}(1)~\text{GeV}\)) experiments will be superseded by their respective successors IceCube Upgrade and the full ORCA detector. In the next decade, we will be able to constrain the \(\delta \Delta\) CPT violation to be less than \(10^{-4}\) eV\(^2\) at 90% confidence level, providing sufficient support to the upcoming DUNE measurement, ruling out the parameter spaces where an impostor solution might exist. Assuming that the nature is CPT-conserving, Fig. 7 shows the amount of significance we can exclude the current best fit, which is only very slightly biased towards a CPT-violating truth at \(5\times10^{-5}\) eV\(^2\), as a function of time in years. Starting with the 7-string IceCube-Upgrade configuration, we consider the addition of ORCA full detector in 2030 [45], we can rule out this point at 1\(\sigma\) in 10 years, effectively ruling out CPT as a source of degeneracy where it is no longer able to produce impostor solutions that differ significantly in \(\delta_{\rm{CP}}^{\rm{meas}}\).

4 Conclusion↩︎

We have studied the interplay between CPT violation and CP violation in long-baseline neutrino oscillation, with the concrete goal of understanding whether a CPT-conserving analysis of next-generation experiments can be trusted to deliver an unambiguous measurement of \(\delta_{\rm CP}\).

The answer, in the absence of external constraints on \(\delta\Delta m^2_{31}\), is no. CPT violation parametrized as an asymmetry between the neutrino and antineutrino atmospheric mass-squared splittings induces an effective, energy-dependent phase shift \(\phi_{\rm eff}(E)\) in the \(\nu_\mu \to \nu_e\) appearance asymmetry \(\Delta P\). Because \(\phi_{\rm eff}(E)\) varies with energy, it cannot be absorbed into a single redefinition of \(\delta_{\rm CP}\) across the full spectrum; nevertheless, once integrated over a realistic flux window, it produces a flux-averaged \(\langle\Delta P\rangle\) that is indistinguishable from that of a CPT-conserving scenario with a shifted CP phase. The result is a continuous family of \((\delta_{\rm CP}, \delta\Delta m^2_{31})\) truth points — a degeneracy manifold — that are experimentally equivalent at the level of the observable a single long-baseline experiment can access.

This manifold has direct consequences for the current experimental landscape. T2K and NO\(\nu\)A probe different \(L/E\) regimes and therefore inhabit geometrically distinct degeneracy manifolds in the \((\delta_{\rm CP}, \delta\Delta m^2_{31})\) plane. Their CPT-conserving best-fit \(\delta_{\rm CP}\) values correspond to different intersection points of their respective manifolds with the \(\delta\Delta = 0\) axis. The tension between those intersections is therefore not a fundamental disagreement between the experiments — it is the geometric consequence of fitting a single parameter to data that is sensitive to two. We showed explicitly that the observed T2K and NO\(\nu\)A central values are simultaneously explained by a CPT-violating truth at \(\delta_{\rm CP}^{\rm true} \approx 109^\circ\), \(\delta\Delta m^2_{31} \approx 1.48 \times 10^{-3}~\text{eV}^2\), a point well within the reach of current experimental bounds.

For DUNE, the consequences are more severe. Using a full GLoBES implementation with an 8-rule CPT-extended treatment that correctly accounts for both the right-sign and wrong-sign beam contributions in each horn polarity, we showed two qualitatively distinct failure modes of a CPT-conserving analysis. First, for a truth with genuine CP violation (\(\delta_{\rm CP}^{\rm true} = -112^\circ\)), a CPT asymmetry \(\delta\Delta m^2_{31} \lesssim -0.3 \times 10^{-3}~\text{eV}^2\) — well within current bounds — is sufficient to cause the CPT-conserving fit to miss the true CP phase entirely, with the inferred \(\delta_{\rm CP}^{\rm meas}\) drifting by tens of degrees. Second, and more strikingly, for a CP-conserving truth (\(\delta_{\rm CP}^{\rm true} = 0^\circ\)), a CPT asymmetry \(\delta\Delta m^2_{31} \gtrsim +0.76 \times 10^{-3}~\text{eV}^2\) triggers a basin switch in the \(\chi^2\) landscape: the global minimum jumps from near \(\delta_{\rm CP} \approx 0^\circ\) to near \(\delta_{\rm CP} \approx \pm\pi\), producing a false detection of near-maximal CP violation from a CP-conserving truth. In both cases the measured \(\delta_{\rm CP}\) follows the analytic \(\langle\Delta P\rangle\) isocontours on either side of the discontinuity, confirming that the bias is a genuine feature of the degeneracy structure rather than a numerical artifact.

The resolution lies outside the long-baseline program itself. Atmospheric neutrino telescopes are the natural tool: their inclusive flux contains both neutrinos and antineutrinos across a broad range of energies and baselines, giving them direct sensitivity to \(\delta\Delta m^2_{31}\) in a single combined fit. Using existing data from IceCube-DeepCore (8 years) and KM3NeT/ORCA-6 (433 kt-yr), we obtain the current world-leading constraint \(|\delta\Delta m^2_{31}| \leq 0.57 \times 10^{-3}~\text{eV}^2\) at 90% CL. With the IceCube Upgrade (10 years from 2026) and the full KM3NeT/ORCA detector (5 years from 2031), this bound will reach \(10^{-4}~\text{eV}^2\) at 90% CL, shrinking the allowed CPT-violating parameter space to a level where no impostor solution can produce a \(\delta_{\rm CP}^{\rm meas}\) that differs significantly from the true CP phase. The combined IceCube+ORCA significance against the current best-fit CPT-violating point at \(\delta\Delta m^2_{31} = 5 \times 10^{-5}~\text{eV}^2\) exceeds \(1\sigma\) within a decade.

The broader lesson is conceptual as much as practical. The \(\Delta P\) observable accessible to any single long-baseline experiment is not a function of \(\delta_{\rm CP}\) alone: it is a function of the full truth \((\delta_{\rm CP}, \delta\Delta m^2_{31})\), and the two parameters are entangled on a continuous manifold. The CPT-conserving assumption does not extract \(\delta_{\rm CP}\) from that manifold; it selects one particular point on it. Only by combining the spectral information of long-baseline experiments with independent external constraints on \(\delta\Delta m^2_{31}\) from atmospheric telescopes can the degeneracy be broken and the CP measurement be trusted. The upcoming generation of experiments — DUNE, the IceCube Upgrade, and full KM3NeT/ORCA — is precisely the combination needed to achieve this, provided their results are analyzed jointly rather than sequentially under the assumption that CPT is conserved.

Code Availability and AI Usage↩︎

The code used to perform analysis and to produce the plots in this paper are available at:

https://github.com/MiaochenJin/CPT-Hides-in-Plain-Sight.

Large Language Models (LLMs) have been used in the development of this paper with heavy human supervision. All data analysis are end-to-end deterministic code pipeline (up to random seeds) and no text or code content has been created ground-up by LLMs. LLMs are used in structural and efficiency improvement of existing code architecture as well as slight language refining of human-written paper text and styling of existing plots. The code release is also rearranged from development-time code with LLM to maximize clarity and reproducibility.

5 Details on Accelerator Experiments’ Distributions and Measurements↩︎

6 Details on Neutrino Telescopes’ Distributions↩︎

Supp. Table 1 shows the groups of common and experiment-specific systematic parameters used in the atmospheric analysis. Supp. Fig. 8 show for IceCube DeepCore and Upgrade (7-String configuration) the true and reconstructed binned event count differences in the leading track sample between a representative CPT-violating scenario (\(\Delta m^2_{31} = 2.511\times10^{-3}~\text{eV}^2\), \(\Delta\bar{m}^2_{31} = 2.0\times10^{-3}~\text{eV}^2\)) and the standard CPT-conserving case. For IceCube DeepCore and IceCube Upgrade-7 respectively, we show the track samples’ event distributions since the track sample contains more information on the oscillatory behaviors, where patterns still exist after the energy and zenith reconstruction smearing. This oscillatory pattern gives rise to the sensitivity of atmospheric telescopes to \(\delta\Delta m^2_{31}\), where the granularity of the pattern after reconstruction smearing dominates the sensitivity.

Comparing IceCube Upgrade-7 with IceCube DeepCore, the minimum neutrino energy goes down to 1 GeV, enabling the detector to see 2 more oscillation maxima, and the reconstruction improves significantly (Supp. Figs. 5–6). Coupling this improvement to a more finely-binned event distribution histograms gives rise to the improvement to the sensitivity of CPT violation as shown in Fig. 5.

In 1 we also provide more information on the sources of systematic uncertainties used in this atmospheric oscillations analysis. Aside from the separate detector systematics, all cross-section systematics and flux systematics are shared between the experiments, thereby realizing a true synergetic combined fit. This prescription is similar to our previous article  [27]. Specifically, IceCube DeepCore employs hypersurface (HS) parametrization for detector systematics, so we also include in Eq. 7 a brief discussion on how these HS systematics are applied in the analysis framework. \[\label{eq:ic-hypersurface} E_i^{\rm IC}(\boldsymbol{\eta}) = \sum_{c} H_{c i}^{\rm IC}\!\left(\boldsymbol{\eta}_{\rm flux}\right) \left[ a_{c i}\!\left(\Delta m^2_{31}\right) + \sum_{j\in{\rm HS}} s_{c i j}\!\left(\Delta m^2_{31}\right) \left(\eta_j-\eta_j^{\rm nom}\right) \right],\tag{7}\] Here \(i\) labels the reconstructed analysis bin and \(c\) labels the IceCube DeepCore hypersurface category, taken to be \(c\in\{\mathrm{NC}+\nu_e\mathrm{CC},\,\nu_\mu\mathrm{CC},\,\nu_\tau\mathrm{CC}\}\). The quantity \(H_{ci}^{\rm IC}(\boldsymbol{\eta}_{\rm flux})\) is the flux-reweighted Monte Carlo expectation in bin \(i\) for category \(c\), after applying the atmospheric-flux nuisance parameters at the event level. The bracketed term is the hypersurface correction factor for that category and bin. Its intercept \(a_{ci}(\Delta m^2_{31})\) gives the nominal detector-response correction, while \(s_{cij}(\Delta m^2_{31})\) is the linear response slope to the IceCube detector nuisance parameter \(\eta_j\). The subtraction \(\eta_j-\eta_j^{\rm nom}\) ensures that the slope terms vanish at the nominal hypersurface point, leaving only the intercept correction. The final expectation \(E_i^{\rm IC}\) is obtained by applying these category-wise correction factors to the flux-reweighted histograms and summing over the three categories.

Additionally, we report in Suppl. Fig. 9 the 1-dimensional \(\chi^2\) profiles of the atmospheric CPT violation data-fit as well as the exclusion sensitivity projection.

References↩︎