MEASUREMENT OF MULTIDIMENSIONAL ION VELOCITY DISTRIBUTIONS BY OPTICAL TOMOGRAPHY

We report the development of a new diagnostic capable of measuring plasma ion distributions as a function of all three velocity-space coordinates. The diagnostic makes use of laser-induced fluorescence (LIF) and computer-assisted image reconstruction techniques. LIF yields high resolution, non perturbing measurements of one-dimensional distributions that are integrated in two directions through three-dimensional velocity space. Computer tomography allows for the unambiguous determination of the complete ion velocity distribution. In addition to a description of the diagnostic, examples of recovered distributions obtained from experiments are given, and the effects of the major steps in the data processing are discussed.


INTRODUCTION
The experimental measurement of the full three-dimensional ion velocity distribution f (x, v ,t) in a plasma provides useful information in a number of experiments in both space l -3 and the laboratory4 but has always been difficult to achieve. Ion energy analyzers typically have a wide acceptance angle and consequently poor angular resolution. Since analyzers usually measure the total flux of particles with energies exceeding a threshold value, the distribution function must be inferred from a derivative of the data,5 a procedure which may be a source of error. Additionally, the presence of the analyzer immersed in a laboratory plasma is a perturbation to the system, which may significantly affect the measurement. Although ingenious construction can produce a narrow acceptance angle,4 the other problems cannot be eliminated so easily. Nevertheless, analyzer data have proven to be very useful in obtaining information from both space b . 7 and laboratory8 plasmas. Optical tomography offers a completely different approach to the problem of measuring ion velocity distributions.
Laser induced fluorescence9-1J (LIF) has been used at UC Irvine and other laboratories in recent years to measure ion velocity distributions both parallel and perpendicular to the externally imposed magnetic field. It is non perturbing to the plasma, has good spatial resolution ( 1 mm:!), good speed resolution (3 X 10 3 cm/s), and does not require that a derivative be taken. A single frequency laser beam (WL ,kL ) is used to excite optical electronic transitions in a barium plasma. Some other ion species also may be used. 10.1 I Ions at velocity v; having an electronic transition frequency Wo may be driven to an excited state provided they satisfy the required Doppler shift (1) The emitted photon flux is proportional to the number density of ions satisfying Eq. (1). Scanning the laser frequency yields the one-dimensional distribution function. As indicated by the scalar product in Eq. (1), ion velocity components along the two axes perpendicular to kL are not preferential1y selected by the laser. Thus the measured distribution is (2) where in this case we have taken the laser beam to be along the x direction.
Due to symmetry about the z axis caused by the constraint condition that the ions move in Larmor orbits, a set of LIF velocity scans spanning all angles in the x-z plane (see Fig. 1) contains all information necessary to reconstruct the complete three-dimensional velocity distribution for the experiments described here. Each scan is a different projection of the distribution. The full set of one-dimensional distributions at a continuum of angles in the plane is known as the Radon transform,14 and may be inverted uniquelyl4.15 to yield the actual two-dimensional distribution, from which the three-dimensional dependence can be recovered easily because ofthe aforementioned symmetry. In practice, a discrete number of scans approximates this ideal, while resolution improves if distributions at a larger number of angles are collected. How the LIF scans are collected is described in Sec. I. The details of the reconstruction process are presented in Sec. II. Results and discussion follow in Sec. III.

I. EXPERIMENT At SETUP
Experiments were performed in the VCI Q-machine 1b (see Fig. 1) with main magnetic field Bo = BJ, where Bo = 6.0 kG, density n = 1 X 1010 cm -3 of nearly completely ionized barium, temperature T; = Te ~0.2 eV, plasma length 1.0 m, and diameter 5 cm. Plasma density was inferred by the use of a Langmuir probe, calibrated by comparison to the angle of propagation of a lower hybrid test wave with respect to Bo({1;:::::w1wJ'<') ' The relevant transitions in the Ba + ion are presented in Fig. 2. Barium ions are produced at the hot plate primarily in the ground state; however, a fraction of about 12% wiB be produced in the long-lived ( ~ 10 s) metastable 5 2 D3/2 state. A single-mode tunable scanning dye laser (Coherent CR599-21, rhodamine 6G dye, shown schematical1y in resonant yellow 5854 A wavelength of the 5 D312 to -312 transition. The dye laser is tunable over a range of 30 GHz with a Iinewidth ofless than 2 MHz. This Iinewidth is much narrower than the thermally broadened absorption Iinewidth ( ;(; = 500 MHz). Hence, ion velocity selection occurs, and is described by Eq. (1). Ions excited to the 6 2 P3/2 state decay to the ground st~te almos~ immediately (6.3 ns) emitting violet light (4554 A). The 6 -P 1I2 state can also be used. 13 The laser light is caused to penetrate the plasma at any desired angle relative to Do by means of an optical fiber (diam = O.6mm) carried in a hook-shaped tube as shown in Fig. 4. The light is refocused into a beam at the end of the fiber by a converging lens (diam = 7.0mm, fL = 1O.Omm).
The collection optics are the same as described by Hill et al. 13 The data acquisition and control are carried out with the assistance of a DEC PDP 11/23 minicomputer in con-  junction with ADC and DAC modules in a CAMAC crate. A block diagram of this setup is shown in Fig. 5. Each onedimensional distribution is stored on floppy disk with minimal initial processing to allow for recovery of the original data later, if desired, regardless of what algorithms have been applied.

II. RECONSTRUCTIVE TOMOGRAPHY
The demonstration that a function in two dimensions could be represented by its line integrals (at all angles) was carried out by Radon l4 in 1917, but no practical inversion techniques were available until the work of Cormack 15 in 1963 and those who followed l7 • 1R in recent years. Cormack's inversion of the Radon transform by expansion of functions in circular harmonics is not in general use, partially due to difficulties with singularities similar to those occurring in the better known Abel inversion, 19 which is used for recovering a cylindrically symmetric function from a single projection. PopUlar methods in use today include such techniques as ART (algebraic reconstruction technique) and SIR T (simultaneous iterative reconstruction technique). The technique probably most commonly found in commercial CAT scanners IK is known as the convolution back-projection method, which is also the technique that we use in our experiments. The widest use of reconstructive tomography today is in the medical field; however, the first practical applications were in radio astronomy. 17  have also been employed in plasma physics to determine physical-space distributions in controlled fusion research devices. 20

A. General principles
Consider a function F(r,c/J) defined in the Vx -V z plane, that is finite and continuous inside a circular region of radius R, and 0 outside. The set ofline integrals along parallel rays across the plane may be expressed as another function of two coordinates, G ( p,e) as illustrated in Fig. 6. We may write this dependence as where R denotes the Radon transform. It may be noted that G is not like functions of polar coordinates that are usually encountered, since its value at p = 0 generally depends on e.
This does not, however, introduce any problems in the method. Note also that a complete set of scans spans only 1T radians, since two scans in opposite directions contain the same information. To find F, the transform is inverted by the convolution back-projection method we now describe. Each individual scan G( p,e = const) is replaced by a "filtered" scan which is computed by convolution with a filter function Q H( p,e) = f~ " " G(p',e) where and FA (u) is the "window." A is referred to as the bandwidth. The choice of window and bandwidth depend on the nature of the data and will be discussed later.
Each filtered scan is "back projected" to form the two- dimensional distribution, as illustrated in Fig. 7. Note that the space has been divided into an array of square patches. Each patch is assigned a value equal to the sum of the projecting function values at arguments determined by the intersection of the p axis for each projection and the normal to that axis which passes through the patch. The equivalence l8 of the convolution back-projection method to the inverse Radon transform will not be demonstrated here. Qualitatively, one may gain an intuitive understanding by first considering the operation of back-projection alone. Small structures (compared to the scan width) will be smeared out in a starlike pattern as the contribution from each projection is added. Larger structures will be smeared less. Since any scan may be represented by a Fourier series, this illustrates the fact that back projection suppresses higher harmonic components of a distribution. The convolution, however, enhances higher harmonic components, yielding a cancellation of these effects which results in the correct image. The choice of A in Eq. (5) determines a cutoff in the high harmonics enhancement. The window function FA will determine how this cutoff is approached. The simplest choice is to set FA = 1, known as the "band-limiting" window. The bandwidth A should be selected to correspond to the highest frequency of meaningful harmonics present in the scans. If A is selected at too high a value, noise will seriously affect the reconstruction. If A is too low, information is lost.
In the situation of our interest, the two-dimensional distribution obtained by the procedure just described is a plane projection along the Vv direction of the three-dimensional distribution. We may partition this two-dimensional distribution into what we shall call strips. A strip is a one-dimensional projection of a plane slice through t.he full distribution at constant V z • Note that the symmetry about Bo ensures that the same strip is obtained at any other direction of projection across the slice, i.e., an the projections of a particular slice are the same. Hence for every slice we have a complete set of projections and may reconstruct the slice by means of either Abel inversion or the convolution back-projection method. The set of all these slices is the three-dimensional distribution.

B. Resolution
It may be shown that if 8u is the smallest resolvable width in a scan width D.u, while m = the number of angles scanned, increasing m will fail to improve the resolution when 8ul D.v significantly exceeds 1T/2m. (In our experiment this does not occur due to the rather high resolution of the laser diagnostic.) Alternatively, the resolution of this technique is limited by the number of angles scanned, and may be expressed as the ratio of the resolvable patch diameter d to the reconstructed velocity space diameter D( = 2R) It is important to note that structure smaller than d is detectable, but will be affected by the smearing mentioned earlier.

c. Special considerations
There are some special problems that can arise in processing data collected by optical tomography that were not mentioned in the earlier general discussion. Probably the most significant of these is the presence of other transitions (due to Zeeman splitting) excited by the laser. This results in the collection of a superposition of identical distribution functions, with each extra Zeeman component having a relative intensity and relative velocity shift which can be predicted theoreticall/ I or measured at sufficiently high magnetic field (splitting is proportional to Bo) and low ion temperature. Before the scans may be processed into a two-dimensional image, data from all extra transitions must be systematically removed. We present a method for doing this. Figure 8 shows a typical example of an actual scan in the various processing stages. In (a) is shown the distribution collected and stored in "raw" form. There are two unwanted Zeeman lines present (Bo = 6.0kG). In (b) the distribution has been smoothed by an iterative routine that averages between neighboring points if a threshold value of the local. second derivative is exceeded. This step is necessary, since the presence of noise introduces errors into the Zeeman removal process. A new horizontal axis (/ = 0 line) has also been generated. Figure 8  distribution is normalized (not shown here) to an area that will be identical to that of all the other scans, required by conservation of spatial particle density. In (d) is shown the distribution after convolution with a filter function, as described in Sec. II A. We now describe the details of the Zeeman removal algorithm. One may write i= 1 so that where F(v) is the measured (superposition) distribution, I (v) is the distribution we seek, N is the number of extra transitions present, A j is the relative intensity of the ith extra transition, and V j is its relative velocity shift. To apply the lowest order correction one could replace the expression for I in the second term on the RHS of Eq. (8) by the measured form of F. In fact, this has been done with LIF scans in the past to eliminate the hyperfine structure. 22 To achieve acceptable results with optical tomography one must carry this expansion further. This may be derived by iterating Eq. (8) to yield Even with the assistance of a minicomputer, the number of operations required becomes unwieldy if expanding to more than a few terms is necessary. One may improve the situation dramatically by recognizing that the quantity inside the nested sums of each term is independent of permutations among the values of the indices. Each permutation need be calculated only once, and then be multiplied by the  ing Eq. (13) at some finite value of M, which we shall refer to as the order of the expansion. Note also that: Using Eq. (13), we have computed corrected distributions to as high as 20th order. Going beyond six orders is seldom necessary. It has been found that this alternating series of nested sums will converge wen if either the intensities of the extra Zeeman components are small or the velocity shifts Vi are all the same sign. The former is due to the multiplication of the intensities together [Eq. (11)] and the latter due to the fact that F( v) tends to zero as v tends to + or -infinity. It should be noted that one has the option of preserving any particular Zeeman component and eliminating the others. Convergence and accuracy will vary however. Generally, to achieve the highest accuracy, one should attempt to preserve the most intense peak, as shown in Fig. 8. If convergence becomes a problem it may be necessary to recover a less intense component located on the side of the distribution and translate the result by an appropriate amount. This is usually necessary if the various intensities are of approximately the same amplitude.
An alternative approach to the algorithm described above makes use of the convenient translation properties of the Fourier transform. If one takes the Fourier transform of Eq. (9), a geometric series results, which may be explicitly summed. This leads to a closed form solution for I (v) where k is the Fourier transform variable. Note the possibility of singularities in the denominator whenever the sum of the relative intensities of the extra components exceeds l.
We have not attempted to implement this method. Another important point found in application of optical tomography to velocity space which has no counterpart in physical space tomography is the possible failure to establish a unique origin. We shall explain this by comparison to a physical space example. First, consider the collection of data in recovering a brain image using a common type of medical CT scanner. At each angle, a set of parallel x rays is made to pass through the patient's head and the attenuated rays at the opposite side are collected. These rays yield the set ofline integrals analogous to what we have called one scan. Although the head is positioned carefully, there is no need for the axis of rotating planes to coincide exactly with the "axis" of the patient's head [see Fig. 9 (a) ]. This is because as the angle is varied, the projection of structures in the brain will be transformed properly, yielding an image that is merely a translation in physical space compared to the one obtained if the offset between the axes is zero.
Unlike images of human head, the velocity space of ions in a laboratory plasma has a very well-defined origin, established by the frame of reference of the plasma device. Although the laser is adjusted to scan over the range offrequencies (velocities) inducing fluorescence in the ions, there is no a priori reason to assume that the center of this scan will correspond to zero for the component of ion velocity relative to the laser beam source, unless this laser adjustment is made deliberately and carefully. Failure to do so means that there exists an offset between the true velocity space origin and the center of the laser velocity scan. However, unlike the effect of an offset in a physical space measurement, the center of the laser scan will not remain in correspondence to a fixed point in velocity space as the laser angle is changed, but will trace out a circular arc as shown in Fig. 9 (b). As a result, the reconstruction process will smear the correct image around Optical tomography a semicircle. Thus, it is important to find the true v = 0 (with respect to the laser) in the laser frequency scan.
The necessary calibration may be accomplished easily by examining a scan along a path perpendicular to Bo and through the quiescent plasma column center. Under these circumstances there is no net motion of the ion distribution along the beam direction, which may be verified by rotation of the laser beam in the x-y plane. The frequency offset of the laser is given by the measured difference between the frequency of the peak and the center of the scan, and may be corrected by readjustment of the laser scan range (hardware) or by introducing a shift in the scans artificially after they are conected (software).

m. RtESUl TS AND DISCUSSION
We have applied the optical tomography diagnostic to a laboratory magnetospheric simulation experiment. Upward flowing ion distribution above regions of auroral activity are often found to have conical distributions (in velocity space) about field lines. Ion conics have been observed through the use of energy analyzers on satellites. 3 They are often seen in conjunction with double layers,23.24 electrostatic ion cyclotron waves,25.26 and lower hybrid waves. n . 2H One suggested mechanism of ion conic formation is perpendicular (to the geomagnetic field) ion heating due to waves, followed by /-lV B forces folding the distribution into a conic.
We present reconstructions of ion velocity distributions found in the Q-machine under three distinct conditions: has an axial drift of about 6.4 X 10 4 cm/s. The perpendicular ion temperature is about 0.17 eV and the parallel ion temperature is about 0.15 e V. The drift speed, Til and T,'II are consistent with plasma sheath acceleration effects in front of the hot plate. The distribution shown in this figure is derived from ten scans in the Vx -V z plane. The integration over Vy has not been removed from the reconstructed image. tive of the combined effects of noise, imperfect Zeeman removal (hyperfine structure was ignored), and using a finite number of scans. Therefore, it is appropriate to symmetrize this result. Function values at opposite points about Vx = 0 may be averaged together, yielding the result shown in Fig.  11 (a). Taking strips at constant V z ofthe symmetrized function and once again applying the convolution back-projection method we present in Fig. 11 (b) an example of a parti-First, the plasma was allowed to drift undisturbed from its formation at the hot plate to the opposite end ofthe machine where it was lost. Second, a large amplitude electrostatic ion cyclotron instability29.JO was produced ll • 32 by drawing an electron current in a channel down the axis of the machine. Third, broadband lower hybrid waves were launched from a 12-cm cylindrical antenna 13 which was coaxial with the plasma column. The magnetic field of the Q machine was uniform. Experiments are pl.anned in which Ho will be in the form of a dipole magnetic field so that the effect of /-lVB forces may be studied. cular plane, f ( Vx ,Vy = 0, V z ), embedded in three-dimensional velocity space.

B. Electrostatic ion cyclotron (EIC) instability
Large amplitude (ed!IT~ 1) EIC waves were excited by an axial electron current produced by positively biasing a 6mm button in the center of the plasma column near the end of the machine. Such a configuration generates electron flow and dc potential profiles 34 similar to those observed in the suprauroral region of the earth's magnetosphere. 23 The ion distribution function in the current channel was significantly affected, as shown in Fig. 12, a Vy -integrated reconstruction made from 16 scans. This distribution bears considerable similarity to some magnetospheric observations 3 of ion conics. We conclude that a f.1 VB force may not be needed to generate ion conics under all conditions. Figure 13 (a) shows the symmetrized form of the distribution, while Fig. 13 (b) is a particular plane, f ( Vx ,Vy = 0, V z ), embedded in threedimensional velocity space.

C. Lower hybrid waves
Large amplitude (e¢>IT~ 1) waves were launched broadband with the center of the range at about twice the ion plasma frequency. We have reported a study ofperpendicular ion heating in the presence of such waves. 35 Figure 14(a) shows a Vy integrated distribution obtained from only eight scans. The noise level was such that the generation of a slice in three-dimensional velocity space is not appropriate in this case. However, a symmetrized form is presented in Fig.  14(b). Note that Til IT;II is approximately 2. Experiments with lower hybrid waves will be pursued later when a mirror magnetic geometry is available. Recall that the recovery of the three-dimensional distribution is the result of several operations applied to the input data, including smoothing out noise, elimination of extra Zeeman components, normalization, and procedures to invert the radon transform. Errors may be expected to accumulate, degrading accuracy of the final result. Significant improvements are anticipated with the development of methods to coHect a potentially larger numbers of scans in a smaller amount of time, which will increase tomographic resolution in accordance with Eq. (6). A greater signal to noise ratio is possible upon use of a different Ba + transition. '3 Some other hardware and software modifications are also being considered.