Tomography by neutrino pair beam

We consider tomography of the Earth's interior using the neutrino pair beam which has recently been proposed. The beam produces a large amount of neutrino and antineutrino pairs from the circulating partially stripped ions and provides the possibility to measure precisely the energy spectrum of neutrino oscillation probability together with a sufficiently large detector. It is shown that the pair beam gives a better sensitivity to probe the Earth's crust compared with the neutrino sources at present. In addition we present a method to reconstruct a matter density profile by means of the analytic formula of the oscillation probability in which the matter effect is included perturbatively to the second order.


Introduction
Our understanding of neutrino has improved greatly since the end of the last century. The observation of flavor oscillations of neutrinos has shown the presence of non-zero masses of neutrinos contrary to the prediction of the Standard Model. This is a clear signature of new physics beyond the Standard Model. Thanks to the remarkable efforts of various neutrino experiments, the mass squared differences and mixing angles of (active) neutrinos have been measured very precisely at present [1].
The origin of neutrino masses is, however, unknown and then it is important to investigate the fundamental theory of neutrino physics. Moreover, it should be useful to consider seriously the application of neutrino physics to various fields of basic science.
In the Standard Model neutrinos are unique matter particles which possess only the weak interaction (in addition to gravity), and their interaction rates are very suppressed accordingly. It is found that most of them can penetrate the Earth without a scattering when the energy is smaller than O (10 5 ) GeV [2]. This shows that neutrinos can be used to probe the deep interior of the Earth.
The idea of the neutrino tomography has been pointed out in the 1970s [3,4]. By measuring the absorption rates of neutrinos passing through the object from different angles, the image of the object can be reconstructed. This is similar to the computed tomography using x-rays, which enables us to probe inside solids without destruction. This method is called as the the neutrino absorption tomography [3]- [23]. In addition to this there have been proposed two other methods of neutrino tomography. One is the method using neutrino oscillations [24]- [44], and the other is the tomography using neutrino diffraction [45,46].
The former one utilizes the energy spectrum of the neutrino oscillation probability, which is distorted, compared to the vacuum one, by the interaction with matter through which neutrinos pass from the production to the detection point. The distortion pattern depends on the profile of the number density of electron in matter, that can be translated into the matter density profile by assuming the charge neutrality and the equality of neutron and proton numbers in matter. It is then possible to probe the deep interior of the Earth by measuring the oscillation probability at the sufficient accuracy.
In this letter we revisit the possibility to realize the neutrino oscillation tomography. The main difficulties of its feasibility include the lack of the powerful neutrino source and no established method to reconstruct the profile of the Earth's interior compared with the medical computed tomography.
As for the first difficulty we consider the neutrino pair beam proposed in Refs. [47,48]. It has been shown that pairs of neutrino and antineutrino can be produced from the partially stripped ions in circular motion at a larger rate than the current neutrino sources from pion and muon decays.
On the other hand, the second one is inherent in the tomography using the oscillation between flavor neutrinos. It has been shown [49,50,51] that the flavor oscillation probability with the density profile ρ(x) is the same as that with ρ(L − x) where x = 0 or L is the production or detection position, if only two flavors of neutrinos are considered. In general the ν α → ν β oscillation probability P (ν α → ν β ) with ρ(x) is equal to P (ν β → ν α ) with ρ(L − x) and an opposite sign of the Dirac-type CP-violating phase [51]. Because of the unitarity conditions 1 = α P (ν α → ν β ) = β P (ν α → ν β ) and the absence of the Dirac-type CP-violation in the two-flavor case P (ν α → ν β ) is invariant under ρ(x) → ρ(L − x). Even for the realistic three flavor case the invariance holds for the oscillations ν e → ν e andν e →ν e [52,53] since they are independent on the Dirac-type phase. #1 In such cases unambiguous reconstruction of ρ(x) is possible only if the profile has the symmetric property with ρ(x) = ρ(L − x). Otherwise, there exist degenerate solutions of reconstruction. It has been, however, proposed that the difficulty can be avoided by using the transition probability of mass eigenstate to flavor eigenstate, which can be realized for the solar and supernova neutrinos [54,32]. Here we focus on the reconstruction of the symmetric density profile with ρ(x) = ρ(L − x), and provide a useful procedure of its reconstruction. Procedures so far proposed are based on the χ 2 analysis (see, for example, Refs. [27,28]), the inverse Fourier transformation [32], and so on. The advantage of ours is that the reconstruction with a sufficient spatial resolution is possible even with a low numerical cost.
This letter is organized as follows: In section 2 we briefly review the neutrino oscillation in matter and present the analytic formula of the oscillation probability based on the perturbation of the matter effect, which will be used to reconstruct the density profile ρ(x). In section 3 it is discussed the oscillation tomography using the neutrino pair beam. We show the possibility of the tomography under the ideal situation. It is then considered how to reconstruct ρ(x) in section 4. Finally, our results are summarized in section 5.
Here and hereafter we assume that all the neutrinos are ultrarelativistic. The free Hamiltonian in the basis of flavor neutrinos is H F 0 is given by #1 This discussion cannot be applied to the case with sterile neutrinos [51].
where E ν is a neutrino energy, m i (i = 1, 2, 3) is a neutrino mass eigenvalue and U αi is the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) mixing matrix [58,59]. The effective potential in the flavor basis is given by where G F is the Fermi constant and n e (x) is the number density of electrons at the distance x. By taking n p = n n = n e and m p = m n in matter n e can be written as The matter density profile is denoted by ρ(x). The oscillation probability measured by the detector at the distance x = L is given by where the initial condition of the amplitude is A αγ (0) = δ αγ . On the other hand, the amplitude of the The oscillation probability for a given matter density profile can be obtained by solving numerically Eq. (2). On the other hand, when the matter effect is sufficiently small, it is obtained analytically based on the perturbation theory [50,60]. We shall expand the amplitude as where A (n) βα (x) is the n-th order correction of the matter effect. The explicit expressions up to the second order are where The oscillation probabilities at x = L up to the second order are then given by When there are only two flavors of neutrinos (ν e and ν µ ), the mixing matrix is given by which leads to the zeroth and first order probabilities where the mass squared difference ∆m 2 = m 2 2 − m 2 1 and On the other hand, the second order expression is given by where Here we have introduced the functions Note that the probabilities for the antineutrino mode is obtained by replacing v(x) → −v(x). The energy spectrum of the oscillation probability obtained by the perturbation theory will be used to reconstruct the matter density profile, which will be discussed in Section 4. In the present analysis the oscillation with two-flavor neutrinos (ν e and ν µ ) is investigated for the sake of simplicity. The study with the realistic three flavor case will be discussed elsewhere [61]. The mixing angle θ and mass squared difference ∆m 2 in vacuum are taken as which correspond to those associated with the solar neutrino oscillation [1].
We consider for instance the matter density profile given by where the background mass density isρ and it is taken asρ = 2.7 g/cm 3 by considering the continental crust. In addition the lump with a density ρ l and a width D l is located at the center L/2.
We show in Fig. 1 the deviation of the oscillation probability from the vacuum one as a function of neutrino energy when the lump density is ρ l = 8.0 g/cm 3 that is close to the iron one. It is seen that the deviation changes in accordance with the density profile. This correspondence is the essence of the oscillation tomography. Furthermore, the deviation due to the matter effect is not so large and then the precise measurement of the energy spectrum is a key for the realization of the tomography.

Tomography by Neutrino Pair Beam
The strong source of neutrino pairs consisting of ν α andν α (α = e, µ, τ) has been proposed recently [47,48]. The neutrino pairs are emitted from circulating heavy ions which are in a quantum mixed state cos θ c |g 〉 + e −i eg t sin θ c |e〉 where the ground state of the ion |g 〉 mixes with its appropriate excited state |e〉 by the mixing angle θ c . The energy gap of these state is denoted by eg and the coherence (mixture) is quantified by ρ eg = sin θ c 2 . The properties of this neutrino pair beam are as follows: the maximal value of the energy sum of neutrino and antineutrino is E m = 2γ eg where γ is the Lorentz boost factor of the circulated ions, the beam is well focused for the large γ region and the effective angular area is approximately given by 1/γ 2 , and the spectrum of a neutrino (or an antineutrino) with energy E ν is given by [47] where the function f (x) is given by and a represents factors of the transition dipole moment and the number of neutrino generations. In this analysis we consider the vector current contribution of neutrino interaction with ionic electron since it gives the oscillation behavior in the probability of the single neutrino detection (while the other neutrino of the pair is undetected) [62]. N and ρ are the number and orbital radius of circulated ions. Note that, when x 1, the function f (x) is approximately given by As representative values we take here aN |ρ eg | 2 = 10 8 , γ = 10 4 , ρ = 4 km and eg = 50 keV as in Ref. [47].
In this case the maximum energy is E m = 1 GeV and the total rate of the pair production is estimated as 5.26 × 10 16 s −1 .
The high intensityν e beam, called as the beta beam [63], has been proposed. Its flux is much larger than the current accelerator experiments as 2.3 cm −2 s −1 . On the other hand, that of the pair beam is estimated as 4.65 × 10 8 cm −2 s −1 , which clearly shows that the precise measurement of the energy spectrum of the oscillation probability is expected.
In this analysis we consider the detector with liquid argon of fiducial volume 10 5 m 3 . The signal event is given by the number of theν e charged current interactionν e + p → e + + n. The cross section is obtained at low energies where m N is a nucleon mass (m N = m p ), g V = 1 and g A = 1.2695. The number of the signal events N (E ν ) for the energy region [E ν , E ν + ∆E ν ] in the detector located at L from the source is given by where V d and n N is the volume and the nucleon density of the detector #2 and T is the duration of the observation. Here the area of the detector is assumed to be smaller than 4πL 2 /γ 2 . It is then for Here we have used the facts that, although the beam produces the pairs of neutrino and antineutrino with all flavors through the charged and neutral current interactions, the dominant one is the pairs of ν e andν e , and that the detection probability ofν e (while the other ν e of the pair is undetected) is approximately given by the neutrino oscillation probability P (ν e →ν e ; E ν ) [62]. It is seen that a large number ofν e events can be detected even at L = 300 km. #3 In order to show the significance of the use of the neutrino pair beam quantitatively, we perform the χ 2 analysis and show how precise the width and density of the lump can be reconstructed when the density profile is given by Eq. (26). We consider the case when the energy spectrum is measured for 100 bins (N b = 100) from E ν = 2 MeV to 100 MeV, #4 and ∆χ 2 is introduced by where the true values of width and density of the lump are taken as D l = 30 km and ρ l = 8.0 g/cm 3 .
Here we take into account only the statistical error and σ(E i ) = N (E i )| D l ,ρ l .
In Fig. 2 contour plots of ∆χ 2 in the D * -ρ * plane are presented for three different cases. The lines with ∆χ 2 = 2.3 (1σ level) and ∆χ 2 = 11.83 (3σ level) [64] are shown. It is seen that there is an approximate degeneracy between D * and ρ * . This is because the leading contribution to the oscillation probability from the lump is proportional to the combination ∆ρ * D * where ∆ρ * = ρ * −ρ.  It is, therefore, found that the neutrino pair beam provides the powerful tool to realize the oscillation tomography. The density and width of the lump can be reconstructed accurately if one takes into account the statistical error only. In order to obtain the realistic precision of the reconstruction we have to include various systematic errors in production, propagation and detection rates as well as the mass squared difference and mixing matrix of neutrinos. This issue is beyond the scope of this analysis.

Reconstruction of Density Profile
Next, we turn to discuss the reconstruction procedure of ρ(x). One of the serious problems to perform the reconstruction is the degeneracy of the oscillation probabilities. In two-flavor neutrino case the probability between flavor neutrinos with the density profile ρ(x) is the same as that with ρ(L − x).
This means that the measurement of the energy spectrum of the oscillation probability is insufficient for the reconstruction when the density profile is asymmetric, i.e. , ρ(x) = ρ(L−x). The reconstruction without ambiguity is possible only for the symmetric density profile with ρ(x) = ρ(L − x). It has been pointed out [54,32] that this difficulty is absent when one uses the solar and supernova neutrinos to probe the interior of the Earth. This is because the initial neutrinos before entering the Earth are not in the flavor eigenstate but in the mass eigenstate due to the MSW effect.
The tomography by the neutrino pair beam under consideration relies on the oscillation probability ofν e →ν e as explained in the previous section, and hence we are faced with this problem. In this analysis we merely assume the symmetric profile to avoid it, which is the first step towards the analyses of more complicated profiles.
Another difficulty of the reconstruction is the cost of numerical calculations. To reduce it we would like to propose a method based on the perturbation of matter effects. Other methods will be discussed elsewhere [61]. Notice that the matter effects can be treated perturbatively for [60] ∆m 2 Thus, the matter density and the neutrino energy should be sufficiently small for a given ∆m 2 .
Our proposed method for the reconstruction is as follows: (i) First, we discretize the neutrino baseline into the N L segments. We take here N L = 60 as an example. The matter densities at each segment ρ i = ρ(x i ) (i = 1, 2, · · · , N L ) are taken as free parameters, which will be determined by using the χ 2 fit to the energy spectrum of the oscillation probability.
(ii) Second, we also divide the energy range of interest into the N E parts. Here the minimum energy is taken as 2 MeV, which is larger than the threshold energy ofν e for the charged current interaction (see the discussion in the previous section). On the other hand, the maximum energy is taken as 100 MeV which is sufficiently small to justify the perturbative treatment of the matter effects. (See Eq. (35)). Here we divide into the N E = 100 parts of even intervals.
(iii) Finally, we determine the densities ρ i by minimizing the χ 2 function which compares the experimental data N obs (E i ) for a given density profile ρ(x) with the theoretical predictions N th (E i ) from unknown ρ i . The explicit form of the χ 2 function is where σ(E i ) = N obs (E i ). In this analysis we take the number of signal events in Eq. (31) estimated from the oscillation probability including the precise matter effect as the observation one N obs (E i ).
On the other hand, the theoretical one N th (E i ) is estimated from the analytical expression of the oscillation probability obtained by the perturbation. The use of the perturbative expression is crucial in reducing the numerical cost to reconstruct ρ i .
First of all, we show the results for the case with the flat density profile ρ(x) =ρ = 2.7 g/cm 3 .
As shown in Fig. 3, the reconstruction by using the probability with the first order correction P (1) (ν e → ν e ; E ν ) is found to be unsuccessful in the present case. On the other hand, when we include the second order correction P (2) (ν e →ν e ; E ν ), the density profile can be reconstructed within a certain error. It is found that the precision of the reconstruction is O (10) % except for the regions near the production and detection points.
Next, we perform the reconstruction of the density profile with a bump or a dip by using the second order expression. The results are shown in Fig. 4. The density at the center for the bump is taken as 8 g/cm 3 which is close to the value for iron, and that for the dip is 1 g/cm 3 which is for water (with a pressure of one atmosphere). It is seen that the profiles can be reconstructed well except for the regions near the production and detection points.
We have found that our proposed method at the second order perturbation works successfully.
A rather complicated profile can be reconstructed as long as the density profile is symmetric ρ(x) = ρ(L − x). See, for example, Fig. 5. It should be stressed that this method can operate even with the limited numerical cost because of the analytical expression of the oscillation probability.
Finally, we have observed that the minimum energy in the oscillation probability used for the reconstruction is important to obtain a better result. In the present analysis we take 2 MeV because of the threshold energy of theν e detection process. As the minimum energy becomes smaller, the result of the reconstruction becomes better. We have found in such cases that the reconstructed densities even near the production and detection points become close to the true values. In addition, the oscillatory behavior of the reconstructed profile at the flat density region becomes changed so that both amplitude (deviation) and wavelength (interval) become smaller, which leads to the better fit. In other words, the main source of the gap from the true profile might be the fact that we can use the oscillation probability for the limited energy range. Especially, those at low energies give a precise information of the density profile. Moreover, the corrections of matter effect at higher orders are also the source of errors in the reconstruction.

Conclusions
We have investigated the oscillation tomography using the neutrino pair beam. It has been shown that the beam can offer the powerful neutrino source to probe the Earth's interior and the precision of the reconstruction of the density profile becomes improved considerably. In addition, we have proposed the tomography method based on the perturbation of matter effects in the oscillation probabilities. This method works only for the symmetric profile with ρ(x) = ρ(L − x) under the situation we have discussed. It has been demonstrated that the profile can be reconstructed well by including the second order correction. We believe that these two ingredients give considerable progress toward the realization of the neutrino tomography.