Iterative inversion of global magnetospheric ion distributions using energetic neutral atom (ENA) images recorded by the NUADU/TC2 instrument

. A method has been developed for extracting magnetospheric ion distributions from Energetic Neutral Atom (ENA) measurements made by the NUADU instrument on the TC-2 spacecraft. Based on a constrained linear inversion, this iterative technique is suitable for use in the case of an ENA image measurement, featuring a sharply peaked spatial distribution. The method allows for magnetospheric ion distributions to be extracted from a low-count ENA image recorded over a short integration time (5 min). The technique is demonstrated through its application to a set of representative ENA images recorded in energy Channel 2 (hydrogen: 50–81 keV, oxygen: 138–185 keV) of the NUADU instrument during a geomagnetic storm. It is demonstrated that this inversion method provides a useful tool for extracting ion distribution information from ENA data that are characterized by high temporal and spatial resolution. The recovered ENA images obtained from inverted ion ﬂuxes match most effectively the measurements made at maximum ENA intensity.


Introduction
Energetic Neutral Atom (ENA) imaging is currently established as a powerful method for providing a global characterization of the Earth's ring current. NUADU is a dedicated ENA instrument flown on the TC-2 spacecraft, which is in an elliptical polar-orbit (700×39 000 km) around the Earth at an inclination of 90 • . TC-2 spins about its axis, which is orientated toward the ecliptic pole, thereby providing NUADU with the capability to image ENAs over a 4π solid angle. Each image provides global information concerning ENAs emitted through the charge exchange process. Conventional charged particle measurements are, in contrast, obtained through in-situ measurements.
Each ENA image contains a great deal of information about the magnetosphere, in particular, about the sources of ENAs generated (via charge exchange) between energetic ions spiraling around the geomagnetic field lines and ambient, cold, neutral particles of the exosphere. The information contained in such images is, however, difficult to extract directly and quantitative techniques for processing and presenting this information are required in order to make such images truly useful.
Over the past several years, significant progress has been made in extracting magnetospheric ion distributions from ENA images. In this regard Roelof and Skinner (2000) developed a nonlinear, parametric, retrieval algorithm. Also, Perez et al. (2000) and Perez et al. (2001) developed a method of constrained linear retrieval using a spline function formalism and DeMajistre et al. (2004) used linear formalism to solve  : TH2-TH3  50-81  138-185  31  47  E3: TH3-TH4  81-158 185-300  77  115 the inversion problem. Validation of this latter method was achieved through comparing ion fluxes retrieved by inverting ENA images with those measured in situ by the Cluster constellation (Vallat et al., 2004). Many papers concerning the analysis of extracted magnetospheric ion distributions based on HENA/IMAGE measurements have already been published (e.g. Brandt et al., 2002a and b;Perez et al., 2004;Zhang et al., 2005). In these retrievals, the ENA images were integrated for at least 10 minutes to achieve adequate counting statistics.
The NUADU instrument has, since the launch of TC-2 in 2004, been monitoring, with high temporal and spatial resolution, the magnetospheric ENA emissions that accompany geomagnetic storms (a one frame ENA image covering a full 4π solid angle is recorded every four seconds with a spatial resolution 11.5 • ×2.5 • ). Generally, typical ENA images are recorded for about an hour outside the radiation belts near apogee during magnetic storms. These ENA image data (typical integration time 16 s) display inherent fluctuations due to the limited statistics of ENA sampling in the short integration time, thus accruing to each pixel. Nevertheless, in these images we can identify an obvious evolution sequence in the area of maximum ENA emission. After a long-time integration, an ENA image is obtained that displays a sharp ENA maximum distribution within a relatively smooth and low ENA emission background. Previous inversion methods are not suitable to process such image data due to the strong coefficient of constraint associatively employed that may smooth the inverted ion flux distribution unduly and, thereby, inhibit the proper recovery of the ENA image within its area of maximum emission. The iterative technique applied and described in the present paper extracts and accumulates characteristic ENA information from step solutions, and the iterative process associatively employed is monitored through consideration of the relative errors pertaining between the recovered ENA image and the image obtained through actual measurement. Also, the iterative technique is configured to be tolerant of the inherent noise background and this enables us to investigate the evolution of magnetospheric ion flux distributions through investigating a set of ENA images recorded, while employing a relatively short integration time, during a magnetic storm.
A brief description of NUADU measurements is presented in Sect. 2. Section 3 provides the basic equations governing these measurements. Simulated images are shown in Sect. 4. The retrieval method is applied, discussed and evaluated in Sects. 5 and 6 contains a summary and general conclusions.

NUADU measurements
The Energetic NeUtral Atom Detector Unit (NUADU) aboard TC-2 features 16 solid-state detectors, each with an equal field of view (11.5 • ×2.5 • FWHM), distributed uniformly over a 180 • angle in the elevation plane. Spacecraft spin allows the azimuthal plane to be divided into 128 equal sectors, by counting timing pulses provided by the spacecraft, so that NUADU completes a full 4π image of the ambient particle population by the end of each spacecraft spin (4 s). It is possible to integrate data obtained in N spins onboard within the range N =2−32 to achieve better statistics. The ENAs measured during individual spins are separated into 3 energy channels (Table 1). (See McKenna-Lawlor et al. (2004) for a more detailed description of the NUADU instrument.) NUADU counts individual ENAs without mass discrimination and, under most circumstances, the instrument has a low noise background. Under disturbed conditions, operation of the instrument in toggle mode (when the output voltage is toggled between 0 and +/−5000 V) , allows the instrument to measure ENAs when the high voltage is ON and ENAs plus charged particles when the high voltage is OFF. The simulated cutoff energy of the NUADU deflection system with nominal deflection voltage (10 kV) is 320 keV for singly charged particles. The measured ENA fluxes include hydrogen and other heavy neutrals (mainly oxygen).
The NUADU image, Fig. 1, used herein to test the inversion technique, was obtained by integrating ENA data recorded during 452 spacecraft spins (30 min) in energy Channel 2 of NUADU during the main-phase of a geomagnetic storm which attained a maximum D st =−113 nT ( Fig. 2) on 22 January 2005. This channel measures 50-81 keV hydrogen atoms and 138-185 keV oxygen atoms, but without mass discrimination. During the storm, ENA signatures could be identified in data recorded during each spacecraft spin. However, the accumulation time for each pixel in the left panel of Fig. 1 was 14.125 s. It is noted that the ENA image shown in the left panel of Fig. 1 was cleaned so as to: remove effects due to the Sun pulse; subtract the noise background of each detector according to noise statistics and correct polar overlaps in the field of view of certain individual detectors. ENA-count distribution curves in azimuth and elevation surrounding the ENA maximum pixel are plotted to the right in Fig. 1 (top-right panel (azimuth); bottom-right panel (elevation)), where the maximum count was 12967. The mean noise background (azimuth region 60 • -140 • ) and elevation region (0 • -100 • ) was lower than 100. The average count for all pixels over 4 was 240. After cleaning the data the noise background was found to be reduced to less than 1% of the ENA maximum (which is an   Fig. 1 in units of R E and in GSM coordinates. The ENA count maximum at the nightside was distributed in a narrow area with a rather sharp boundary close to the profile of the Earth (the solid white ring in the left panel of Fig. 1). The shadow of the Earth is represented within the solid white ring by brown-purple and by black colors. The bright region, seen mainly around the outside of the profile of the Earth, is interpreted to be due to ENA emissions from the magnetic L-shells external to L=4.

Simulation equations
Simulation bins for the ENA source were established by dividing the magnetic L-shells (from L=2 to L=8) into segments of equal longitude and latitude in the GSM coordinate system. Each volume element was assigned coordinates described by its L-shell longitude ϕ and latitude θ. The space outside the simulation bins was assumed to be an optically thin medium with respect to ENA flux, so that ENAs could travel along the line-of-sight (LOS) without undergoing charge-exchange.  The counts recorded in each pixel of an ENA image, C δ,ε (pixel with elevation δ and azimuth ε), is represented in the simulation system by the equation where dV is the volume element of the LOS of each pixel; E is the energy range within which neutral atoms are detected; T is the time period during which neutral atoms are 3 +X -X   detected; is the solid angle of the volume element pointing to the δ, ε pixel; A(δ, ε) is the response function of a detector at that pixel with elevation δ and azimuth ε, j ionH and j ionO are the unknown ion fluxes of hydrogen and oxygen; σ H (E) and σ O (E) are the respective charge exchange cross sections between energetic hydrogen and energetic oxygen ions and cold exospheric neutral atoms and n is the exospheric neutral atom density near Earth. Integration is along the LOS of the instrument (note that in the analysis of our data we only consider charge exchange between hydrogen atoms, which comprise the dominant exospheric species). Figure 3 describes the geometry of the NUADU detectors. In this system, the ENA flux was emitted from an Earth-based volume element (at longitude ϕ, magnetic latitude λ, and L) with elevation δ * and azimuth ε * pointing to the volume element from the sensor. The simplest response coefficient of the pixel with normal direction δ and ε at the center of the pixel is given by

Magnetospheric ion flux
On combining variances of energy, pitch angle and magnetic L index, an empirical model of the ion flux in the ring current region deduced from a Kappa distribution (Christon et al., 1991) may be expressed in the form (Shen, 2002) where e=(1+1/κ) κ+1 ≈2.962 (κ=5.5) (Christon et al., 1991), λ is the magnetic latitude, E 0 (7 keV for protons, 16 keV for oxygen ions) (Shen, 2002) is the typical ionenergy of the maximum ion-flux, and L 0 =7.3 is the outer boundary of the ion injection region of the ring current (RC). The pitch angle (α) of the particles within the volume element is expressed in terms of λ (magnetic latitude) and α eq (equatorial pitch angle), such that Since in the absence of mass analysis the proton flux distribution cannot readily be separated from the oxygen ion flux, both are combined in the inversion Eq.
(1) to reflect a mixed ion component, j * (L, ϕ, λ, E, α), as: where j * 0 (L, ϕ) represents the dimensionless unknown part of the ion-flux introduced to appropriately modify the distribution, and j H0 and j O0 are the initial inputs of ion fluxes for protons and oxygen ions with empirical constants.

Neutral density of the exosphere
In general, neutral hydrogen is dominant in the exosphere and we have only considered charge exchange with exospheric hydrogen atoms. Its density, according to Rairden et al. (1986) and Brandt et al. (2002a), in an adjusted Chamberlain model is represented by where H (ϕ, θ ) = 1.46 (1 − 0.3 sin θ cos ϕ) .
Here, r is the geocentric distance, ϕ is the local time angle from noon and θ is the polar angle from the z-axis in the GSM coordinate system. The coefficient of 0.3 in Eq. (9) allows one to represent the geotail in the approximated form of an axis-symmetric (around the Sun-Earth line) exosphere. All numeric coefficients in the above expressions have been obtained from the fit to the Chamberlain model by Rairden et al. (1986) and Brandt et al. (2002a). We would like to stress that although the existence of a geotail is confirmed, quantitative knowledge concerning it is presently sparse.

Numerical quadrature format
The integral in Eq.
(1) can be calculated numerically along the line-of-sight so that where i represents the ϕ coordinate, j the λ coordinate and k the L coordinate of the volume element. The summation operation yields, via j , a set of linear algebraic equations which can be written in (matrix) form where K is the coefficient matrix with elements of the quadrature weights computed from Eq. (10), and J * 0 (L, ϕ) is a column matrix with elements j * 0 (L, ϕ), as defined in Eq. (7), Sect. 3.2. In Eq. (11) the subscripts to the parentheses indicate the dimensions of the matrices, where p<δ×ε and q≤i×k depend on the number of volume elements and on the number of pixels related to the source of ENAs generated through charge-exchange.

ENA image simulation from empirical models
We now use the equations developed in the previous section to simulate NUADU measurements. During a storm, the oxygen ion flux is usually enhanced so that it exceeds the level of proton flux by about 10% (Lennartsson and Sharp, 1982;Mitchell et al., 2003). In-situ measurements made by RAPID/Cluster after the Cluster spacecraft entered the magnetosphere at 13:00 UT on 22 January (that is on the same day, but after the ENA measurements presented in Fig. 1 were recorded) are shown in Fig. 4. Before 13:00 UT, the Cluster spacecraft operated in the solar wind on the dayside. Figure 4 shows that the mean integral heavy ion flux averaged over 4π (E>90 keV), mainly in oxygen ions, was about 5.7×10 4 s −1 sr −1 cm −2 while the mean integral proton flux averaged over 4π (E>27 keV) was 4×10 5 s −1 sr −1 cm −2 . These species maintained this ratio during the evolution in a  composition that occurred during the remainder of the day. Further, as the Cluster spacecraft (C#3) crossed the equatorial plane at L=4.2 dusk-night sector at 20:52 UT (while D st =−69 nT in the recovery phase of the storm), the integral heavy ion flux averaged over 4π (E>90 keV) was 4.1×10 5 s −1 sr −1 cm −2 , and the integral proton flux averaged over 4π (E>27 keV) was 2.9×10 6 s −1 sr −1 cm −2 (i.e. it was about 7 times higher than the former value). Thus, the preliminary parameters' input to the simulation was set at 4×10 5 and 5.7×10 4 s −1 sr −1 cm −2 keV −1 , respectively, for the proton and oxygen ion fluxes. That is, where h(ϕ)= exp [−0.75 (1 − cos(ϕ − ϕ s ))], ϕ s is the azimuth angle shifted to the location where the intensity of the ring current (RC) particles showed its maximum value; the relevant charge exchange cross sections are σ H =2×10 −17 cm 2 and σ O =8×10 −16 cm 2 (Smith and Bewtra, 1987). Figure 5 shows the results of a simulation of NUADU measurements made for 22 January 2005. The azimuth angle ϕ s was chosen to be 160 • , so that the ENA maximum was in the dusk-night region, which is similar to what was found in the actual measurements (see the left panel of Fig. 1). The 5 igure 5. Although it is difficult (Sect. 1) to derive information, even qualitatively, directly from an ENA image, a recovery technique is presented below which is designed to retrieve from the NUADU data within the prevailing constraints, the ion distributions corresponding to the measured ENA images. In the absence (see above) of mass analysis we have to assume the initial ion flux ratio between protons and oxygen ions and extract a mixed distribution based on the empirical model and this assumption.
In the mixed ion component inversion method described below, the maximum ion flux depends on the total percentage of oxygen ions present. Since the charge exchange cross section of oxygen ions is one order of magnitude higher than that pertaining in the case of protons, oxygen ions are more effective in generating ENA emissions than are protons. Thus, the higher the percentage of oxygen ions present, the lower the maximum mixed ion flux is that will be inverted.

Constrained linear inversion by iteration
In a strict sense, Eq. (11) can only be solved uniquely if the matrix K is at once square and non-singular. In practice these two criteria cannot be met simultaneously. Instead of seeking a direct solution, we pursue a least-squares method under constraints (see also DeMajistre et al., 2004); i.e. we determine J * 0 through minimizing the quantities where σ −2 C is the inverse of the measurement covariance matrix, γ is the constraint strength and H is a constraint matrix (here, H=H 0 =I, is adopted as the first order constraint matrix in an iterative solution). For NUADU data, σ −2 C is a diagonal matrix whose elements are 1/σ 2 i , where σ i is the uncertainty corresponding to each pixel. Minimizing Eq. (14) yields When applied to NUADU measurements, the results of a linear least-squares solution were found to oscillate, and even to present negative solutions. To solve this problem we increased the constraint strength (this approach will be discussed further in Section 5.1), and also limited the oscillatory solution region through introducing a mapping-function written as where x indicates a set of linear least-square solutions of Eq. (15) and x max is the maximum, x min the minimum and x mean the mean of this set of step oscillatory solutions, ξ is a step-length coefficient, and x * is a temporary solution which has a limited solution region 0<x * <2 around 1. The ion flux distribution J * 0 could then be approached through an iterative process (which will be described in Sect. 5.2) using the dimensionless temporary solution x * as J x * , where superscript n indicates the iterative sequence, J H0 and J O0 indicate the background inputs of ion fluxes for protons and oxygen ions, which are equivalent to the initial inputs of j H0 and j O0 in Eq. (7).
Instead of comparing the inversion results with the measurement covariance (see DeMajistre et al., 2004), we opted instead to compare the inversion results directly with the measurements. The effect of each converging step (constraint strength variances) could then be justified and monitored by adopting a relative error To better represent properties of the maximum ENA distribution, we chose to include in the statistics only those pixels having measurements above the average count of 240. The relative error minimum then attained among the iterations (constraint strength variances) represents the optimal ion distribution solution. The value of the relative error depends on the number of pixels taken into account in the statistics.  was obtained (Fig. 6). It is noted that relatively low constraint strengths resulted in a high ion flux solution, whereas relatively high constraint strengths resulted in low ion-flux solutions. The higher the empirical ion fluxes input, the higher the constraint strengths were, resulting in an optimal solution.

Constrained linear inversion
The extracted equatorial ion flux with pitch angle 90 • at the optimal constraint strength was found to be distributed rather smoothly with a maximum intensity of 1.56×10 6 s −1 sr −1 cm −2 keV −1 (Fig. 7). The most intense patch in the ion distribution was located in the post-dusk sector (L=4) at about 19-20 h local time and it extended toward both dawn and noon with a gap at dusk (which corresponded with measurements of the high polar angle detectors). This result showed a match with the measurements (see the gap at 18:00 MLT in the left panel of Fig. 1), see also the corresponding recovered ENA image (Fig. 8). The maximum extended outward widely and smoothly in the dusk-night sector, which is also consistent with the extent of the ENA maximum distributions shown in the left panel of Fig. 1 and in Fig. 8. On comparing the recovered image with the left panel of Fig. 1, the maximum ENA count area in the recovered image (Fig. 8) was found to be shifted duskward and to be spread more widely, increasing as a result the relative error. It is stressed that the extracted flux was found by varying the constraint strength γ in Eq. (15) under circumstances where the maximum ENA measured counts were located on the dusk side and the linear constrained method smoothed the inverted ion fluxes around the ENA intensity maximum area.

Iterative inversion
In order to better represent the properties of the intensity patch in the ENA measurements (left panel, Fig. 1), we next performed an iteration to retrieve the ion distribution through employing a constant constraint strength (γ =0.04) that was smaller than that associated with the optimal solution of the previous section. As already indicated above, the ion distri-    O0 x * , where n indicates the iterative sequence. J H0 and J O0 represent the background inputs of ion fluxes for protons and oxygen ions, which are equivalent to the ion values in Eq. (7). After each iterative step the numeric field of the high limit, dimensionless, ion-flux solution was enhanced, as a consequence of the ever increasing input made by the fluxes that characterized each successive adjustment.   The mapping function (Eq. 16) which was associatively introduced, enabled us, through employing a linear, regionlimited, artificial solution, to utilize the enhancing ion flux modifications to guide the step-by-step solutions, so that, gradually, they approached the actual measurements of ENA intensity.
The success of the iterative process in improving, after each step, the input of the ion flux background was determined using the relative error Eq. (17) and the same statistical pixels considered in Sect. 5.1, until a minimum value (0.24, Fig. 9) was obtained. We can follow in this figure how the solution was adjusted step by step. After the relative error minimum, the error profile was found to rise in approximately the same way in which it fell off (even somewhat faster, since the pixels with large ENA counts were still enhanced through the mapping function (Eq. 16)).
The maximum intensity of the extracted equatorial ion flux with pitch angle 90 • at the minimum relative error (0.24) was  1.43×10 6 s −1 sr −1 cm −2 keV −1 (Fig. 10), which is of about the same order of magnitude as that shown in Fig. 7. The intense patches of the ion-flux distribution were shown to have a complex construction featuring relatively sharp boundaries that were consistent with both the measured ENA flux maximum distribution (left panel, Fig. 1) and the recovered image (Fig. 11). The intense ion-flux configuration extended outward on the dusk side but was constricted inward on the night side, which may represent an RC ion distribution region that was compressed earthward during the main phase of the storm (Fig. 2).
The recovered image, characterized by a minimum relative error, is shown in Fig. 11 and it is, in the main, more closely matched to the measurement shown in the left panel of Fig. 1 than is the simulation presented in Fig. 8. To better compare the inversion results obtained using the two methods with the NUADU measurements, we plotted the ENA count distributions versus azimuth and elevation (left and right panels, Fig. 12) for the location of the image maximum, where the black curves represent the measurement in Fig. 1, the red and green curves represent the corresponding recovered images of inverted ion fluxes obtained using, respectively, the iterative scheme (Fig. 11) and the linear constraint scheme (Fig. 8). The blue lines represent the mean ENA count in Fig. 1, which was set as the low limit for relative error calculations. Figure 12 indicates that, above the mean ENA count of the measurements, the iterative scheme results in a closer approach to the observations, especially around the ENA peak region (azimuth 140 • -250 • ; elevation 130 • -180 • ). In regions of lower ENA intensity (elevation 115 • -130 • , right panel Fig. 12), where the ENA count was lower than 110, the   Fig. 1 at the ENA intensity maximum; blue lines represent the mean ENA count in Fig. 1; red and green curves represent recovered images from inverted ion fluxes derived, respectively, by iteration (Fig. 11) and by linear constraint (Fig. 8). For further details see the text.
green curve (obtained using the linear constraint scheme) displays a close approach to the black curve (obtained through measurement).
From the viewpoint of relative error statistics with respect to the high intensity area, the precision of the iterative scheme (with relative error 0.244) is higher by a factor of about two than that derived when using linear constraint inversion (relative error 0.432). It is thus inferred that the iterative technique is the more suitable choice for ion flux extraction from an ENA image characterized by a sharp ENA peak.

Iterative inversions for a short integration time
Since the mapping function Eq. (16) enables us to better cope with oscillatory solutions, the iterative technique can be applied when studying relatively short-integrations of ENA flux measurements made by NUADU during geomagnetic storms. To demonstrate this capability, the method described above is next applied in the case of short integration times. First, we divide the measurements already discussed into six segments and adopt five-minute integration times for each segment (see Fig. 13, which shows the temporal evolution of the ENA measurements concerned). For the same statistical pixels as those already considered above, we derived a set of ion-flux distributions with an average relative error 0.2 during the evolution sequence (Fig. 14) and recovered ENA images (Fig. 15) from each extracted ion-flux distribution.
The short integration measurements in Fig. 13 look rather patchy as compared with those shown in the left panel of Fig. 1, but it is still possible to identify structures within the ENA flux maximum (red color) and in the high latitude ENA distributions (bright-blue color). Figure 13 shows that during the evolution process, the five-minute ENA distributions were, especially in the maximum ENA area, initially enhanced in the first three segments and they then decreased and shifted duskward in the following three segments.
Generally, each ion flux color map in Fig. 14 shared about the same distribution configuration (similar to that seen in Fig. 10), but viewed in detail, the intense ion-flux patches exhibited magnitude changes and shifting locations. The ionflux distribution in the RC region seems to have displayed a configuration that was initially compressed and then recovered during the main phase of the storm (Fig. 2). Figure 15 shows that each recovered ENA image matched well those measurements made (Fig. 13) in the ENA maximum area.

Summary and conclusions
In the present paper we developed a technique for the retrieval of global magnetospheric ion-flux distributions using an iterative procedure. Its capability was demonstrated through being applied to representative ENA images recorded by NUADU aboard TC-2 during a significant magnetic storm. The basic equations adopted in the retrieval process, including the NUADU response function, an ionflux empirical model and the exospheric neutral distribution, were each presented and described. Simulations based on the equations derived were displayed and their important features described. It is noted with regard to the technique adopted that ion-distribution retrieval constitutes an ill-posed problem without a unique solution (noisy measurement data makes least-square solutions oscillatory while the introduction of high-constraint strength smoothes the inversion result too much relative to the measurements). As a result, a mapping function was introduced to select key characteristics from the step solutions and improve the fidelity of the extraction through the application of an iterative process. The latter process was monitored through determining the variance of the relative error between simulation and measurement.      By applying the technique to a representative ENA image recorded during a storm, a recovered ENA image was obtained that showed a close approach to measured reality in the brightest ENA area.
In brief, a technique was developed to extract magnetospheric ion distributions contained in ENA data, together with a means to quantify their quality. This method, which was custom developed for NUADU data recorded over a 4π solid angle, was demonstrated to be applicable in cases of both high and low sensor count rates. It thus constitutes a useful tool for the inversion of global magnetospheric information from ENA images obtained at high temporal and spatial resolution. This method can also support investigations of the evolution of ion-flux distributions extracted from ENA measurements made with short integration times during the development of geomagnetic storms.
The iterative inversion technique described potentially supports the retrieval of magnetospheric ion distributions from multiple ENA images. If, for instance, two ENA images of the ring current were recorded simultaneously from different positions, the ion-flux distribution in the ring current region could be alternately adjusted between the two measurements so that, ultimately, a unique ion-flux distribution solution matching the two simultaneous measurements could be achieved. The NAIK instruments in the straw man payload of the twin spacecraft of China's planned Kuafu-B mission can potentially be utilized to accomplish this goal.