Many of the observational data in space physics are obtained as time series. These time series contain measured physical quantities, that may be scalars, vectors, tensors, or multidimensional images. Examples of such quantities are temperatures and densities of the plasma, magnetic or electric field vectors, pressure tensors, and auroral images. These quantities can be measured either in space on a moving platform or on a fixed platform such as the surface of the Earth. Thus the variations in the time series may represent true temporal changes in the system or motion through spatial gradients or some combination of the two. Since most data gathered in space physics are initially in the form of time series, their use is widespread. To analyze these time series data, many data processing methods, analysis techniques and computer algorithms have been developed. In this review, we outline a set of principles for data analysis methods, describe a number of well-established data analysis techniques, discuss the uncertainties and limitations of each technique, and suggest procedures and criteria which may reduce the uncertainties of the results for some analyses. Many of the principles presented are derived for the first time.
Most time-series analyses can be divided into three categories: discontinuity analysis, wave analysis and correlation analysis. Discontinuity analysis determines the orientation, thickness and motion of the interface between two different plasma regions or regimes. The methods for discontinuity analysis are well developed. However due to the lack of wide recognition of the underlying assumptions, uncertainties and validity of each method, there are still many problems in this area. We present a comprehensive discussion of these problems. Wave analysis determines the properties and characteristics of a wave that identify which of the several possible wave modes allowed in a plasma a particular observed fluctuation might be. Techniques for wave analysis are under active development. We discuss the principles of these methods. Since discontinuities can be considered as steepened waves, some of the methods for wave analysis can be used in discontinuity analysis. Correlation analysis determines the relationship between observations at two or more spatially separated locations. Correlation studies can be performed in the time domain, for example to determine the time lags between observers. Correlation analysis can also be performed in the frequency domain. We will give a brief introduction to correlation analysis and discuss how to apply it to the analysis of the data from a cluster of satellites.
From the beginning of the space age, scientists have sought means for quick qualitative visual examination of large amounts of data. With plotted traces, one can easily spot a discontinuity, estimate a wave frequency and correlate two traces by overlaying them. Even today these qualitative visual analyses can be powerful ways to check the results of an analysis, i.e., a quantitative result is suspicious if one cannot verify the consistency of the result by visual inspection of the data. As will be stressed in this paper, relying solely on automated computer analysis may lead to completely wrong results. As a data analyst, one has to frequently return to the examination of the original observations to perform what is often referred to as a "sanity check" or a "reality check".
The first quantitative data analysis method for discontinuities was the minimum variance analysis for discontinuity analysis proposed by Sonnerup and Cahill . In the 1970's, the wave analysis techniques using the magnetic field became mature. Along with the improvements of the plasma measurements, in the 1980's and 1990's, techniques incorporating plasma measurements have been under active development. The correlation analysis of signals made at two or more observation sites has become important since the launches of the ISEE satellites and it will play the key role in analyzing the multiple satellites data in the International Solar Terrestrial Physics (ISTP) program. While in theory the quantitative correlation analysis is most straightforward compared with other analyses, in practice, there are many problems associated with the coherence length of a phenomenon, compared with the spacecraft separation, for example. We suspect that many problems have not yet even been recognized. In the next few years, we expect to see major progress in understanding of correlation analyses as a result of multiple-satellite data analyses.
1.1 Frame of Reference.
Observational data are gathered in a frame of reference at rest with the observer. Physical laws are often stated in a frame of reference that is moving with respect to the observer. In particular, the Rankine-Hugoniot relations are given in the frame of reference at rest with the discontinuity, and the wave dispersion relations usually expressed in the plasma rest frame. In the non-relativistic limit, plasma density, pressure, magnetic field, and wavelength are not dependent on the frame of reference, but velocity, electric field and frequency are. While it is always important to work in the appropriate frame of reference when dealing with a flowing plasma, this is critically true in the supersonic solar wind.
The total time derivative of an observed quantity Q, which can be either a scalar or a vector, is
where v is the velocity of the observer relative to the frame of reference in which the phenomena is described. The first term on the right is the temporal variation in the frame of reference of the plasma (say) and the second term on the right is the apparent temporal variation caused by motion through spatial gradients. Both are observed by a moving instrument as temporal variations. For example, when discussing observations of a wave in the plasma frame, v is the measured plasma flow velocity, and the term is nonzero. For a steady phenomenon, the term is zero (e.g. in the plasma frame) and the observed time variation is due to the motion through stationary gradients. A time-domain correlation analysis assumes that, in a particular frame, for example, a wave frame or a shock frame, = 0. Therefore the velocity derived from the timing difference is the sum of the flow (convection) velocity and the propagation (phase) velocity. We recall that the group velocity measures the velocity of the envelope of the variations whereas the variations themselves move at the phase velocity. The energy of the wave packet flows with the group velocity and it is the group velocity that is restricted to velocities equal to or less than the speed of light, not the phase velocity.
For magnetic field measurements, using the frozen-in Faraday's law and divergence free condition, equation (1.1) becomes,
Assume a coordinate system with variations only along the normal or a direction under one-dimensional assumption, for either a discontinuity or a wave, in which and are zero. It is easy to show that the measured time variation of the field is zero in the direction along n. No stationarity assumption is made here and this is true independent of the frame of reference. This particular feature of the magnetic field is actually the basis of most data analysis techniques. Note that in discontinuity analysis, one cannot always assume that the shock frame is known, nor in wave analysis, that the wave is in steady state.
For other quantities, temporal variations in the frame in which a theory is applied, may affect the results of an analysis. For example, the observed electric field variation is
where is the electric field in the frame of reference in which the phenomenon is described. Note here that and we have assumed v to be uniform in space. When v is not spatially uniform, more terms appear. Under a 1-D assumption, neither component is zero in general. But for 1-D steady state, the tangential components of the variation vanish whereas the normal component does not. This property of the electric field variation provides the foundation of the maximum variance analysis to be discussed in section 2.5. Note that this method requires the steady state assumption and thus is difficult to apply to wave analyses. For discontinuity analysis, not only unsteadiness of the discontinuity but also changes in the motion of the discontinuity relative to the observer affect the results. Since in space, boundaries are often in oscillation and not in simple motion, only on a few occasions is the latter effect unimportant. When comparing measurements from two spacecraft at different times, more effects will occur due to the relative motion between the spacecraft frame and the plasma frame (see more discussion in section 4.1.3).
1.2 Coordinate Systems.
1.2.1. Global coordinate systems.
Three global coordinate systems are often used in space physics. A global coordinate system is useful for studying phenomena of global effects.
GSE Coordinate System. In the Geocentric-Solar-Ecliptic coordinate system, the x direction is from the Earth to the Sun. The z direction is along the normal to the ecliptic plane and pointing to the north. The y direction completes the right handed coordinate system and points to the east, opposite planetary motion. The statistical aberration of the solar wind by the Earth's orbital motion is most readily removed in this system. It is useful for problems in which the orientation of the Earth's dipole axis is not important, such as the bow shock and magnetosheath phenomena.
GSM Coordinate System. In the Geocentric-Solar-Magnetospheric coordinate system, the x direction is also from the Earth to the Sun. The z direction lies in the plane containing the Sun-Earth line and the geomagnetic dipole, and is perpendicular to the Sun-Earth line and positive north of the Sun-Earth line. The y direction completes the right-handed coordinate system. The GSM y and z directions lie at an acute angle (around the x direction) to the GSE counterparts. The GSM coordinates are most useful for studies of solar wind-magnetosphere interaction and where the solar wind determines the geometry of the magnetosphere. These include phenomena near and at the magnetopause/boundary layer, in the outer magnetosphere and near earth magnetotail.
SM Coordinate System. In the Solar Magnetic coordinate system, the z direction is along the Earth's magnetic dipole pointing north. The x-z plane contains the solar direction with roughly toward the sun. Its y direction is on the dawn-dusk meridian and points to dusk in the same direction as the GSM y direction. Thus the SM x and z directions differ from their GSM counterparts by an acute angle of rotation related to the tilt of the dipole. The SM coordinates are usually used in studies of the near earth phenomena, where the geomagnetic field controls the processes, such as the ionosphere, inner magnetosphere, and ground observations.
1.2.2. Local coordinate systems.
To study a local phenomenon or a global phenomenon in a local context, a local coordinate system is most convenient. Three local coordinate systems are often used in the literature.
Boundary Normal Coordinate System (LMN) [Russell and Elphic, 1978]. The direction is along the normal of the boundary, which can be the bow shock or the magnetopause, or a current sheet. The direction is usually chosen to be outward from the Earth for both the bow shock and the magnetopause. The direction is along the geomagnetic field direction for the magnetopause and is along the projection of the upstream magnetic field on the boundary for the bow shock. The direction completes the right-handed coordinate system and points to dawn for the magnetopause (see Fig. 1.1). At the subsolar magnetopause, neglecting the aberration caused by the finite velocity (29.5 km/s) of the Earth through the solar wind, the LMN coordinates are coincident with the GSM coordinates with and .
|Fig. 1.1. Boundary normal coordinate system [Russell and Elphic, 1978].|
The LMN coordinate system is useful only locally unless the boundary is planar. In some cases, a satellite may move along a curved boundary, e.g. the magnetopause, for a long period of time and cross it many times as the boundary rocks back and forth. The normal of the boundary at each crossing may be different. When presenting the data in a fixed LMN coordinate system, one needs to be extremely careful in interpreting the results because the direction of the normal of the boundary may be varying. Similarly one should be cautious in using an LMN coordinate system derived from either the magnetopause or bow shock well away from that boundary such as in presenting the data from the bow shock to the magnetopause. One should not interpret, for example, the N component along a direction determined at the magnetopause as the normal component throughout the magnetosheath.
Field-Aligned Coordinate System. The average magnetic field direction (taken from the in situ measurements) is defined as a preferential direction. A second direction is usually defined according to the symmetry of the system. For example, for magnetospheric wave studies, the azimuthal direction parallel to the direction obtained from the cross product of the magnetic field and the radially outward direction is usually used as the second direction. This direction is eastward in the earth's magnetosphere in the direction of electron drift. This direction is useful for displaying the field perturbations associated with sheets of field aligned currents that are roughly along shells of constant L-value or L-shells. Similarly, for magnetosheath studies, one can define a surface in the sheath which belongs to the family curves that are used for empirical models of the magnetopause and bow shock. The second direction can then be defined perpendicular to the field and along the surface.
Geomagnetic Dipole Coordinate System. The Geomagnetic Dipole coordinate system has its z direction along the geomagnetic dipole axis. The other two directions are defined in terms of geomagnetic latitude and longitude by analogy with geographic coordinates. Therefore, this coordinate system is useful in studying phenomena observed in ionosphere and ground stations. In order to study the global effects of these phenomena, an observation is often presented in geomagnetic local time, which is equivalent to that in SM coordinates, and is the angle between the meridian plane containing the sun and that containing the point of observation converted to hours and increased by 12 hours.
1.2.3. Spacecraft coordinate system.
For some vector quantities only two components are measured by our spacecraft. Examples are the electric field, plasma waves and the plasma moments from certain instruments. The two components are in the plane perpendicular to the spacecraft spin axis. Frequently, the spin axis is nearly along one of the GSE axes (for magnetospheric missions) or perpendicular to the orbit plane (cartwheel mode for ionospheric missions). The fraction and the direction of the background magnetic field projected onto this plane are important to observations which deal with the anisotropies with respect to the magnetic field. The spin modulation of the high time resolution signals can provide useful information. Often artifacts in measurements are most easily detected in spacecraft coordinates.
1.3 Measurements and Their Uncertainties.
A measured physical quantity is derived from measurable quantities by a instrument. Therefore, its accuracy is limited by the instrument capability and the data reduction procedures. As a data analyst, it is very important to understand the principles of the instruments from which the data are measured and the schemes of the data reduction that are used. In this section, we will discuss several general issues concerning the measurements and data reduction. We will assume that the measurements are accurate in other sections in the paper. The word "fluctuation" used in this paper refers to deviations from the average of the measured quantity. The word "noise" refers to the unwanted signals that are not described by the idealized variations of a particular analysis method. It is important to point out that the "noise" in one analysis can be a physical phenomenon for another analysis. Fluctuations include both wanted signals and noise. If one chooses a wrong method, the useful signals will be treated as noise.
Time resolution. Temporal resolution, the time between samples, is a major limit in many data analysis methods. For example, to use the Minimum Variance method (to be discussed in section 2.2) requires the resolution to be at least half of the crossing duration of a discontinuity, i.e. one needs a sample in the middle of the crossing not just at either end. Temporal resolution also sets an upper cutoff frequency, the so-called Nyquist frequency, which is half of the sampling frequency, for Fourier analysis. Stated differently, in order to resolve a waveform, the sampling frequency should be at least twice the wave frequency. If frequencies above the Nyquist frequency are present when a signal is digitized it appears to oscillate at the frequency below the Nyquist frequency. Such a signal is called aliased. The timing difference derived from correlation analysis is limited by the accuracy of the clocks controlling the measurements. Synchronization of clocks of observers becomes very important when the timing differences are small. Sometimes components of vectors are not sampled simultaneously. Sometimes the time separation of sequential measurements is not uniform. Techniques exist for accurately analyzing signals under these circumstances but add complexity to the analysis.
If the bandwidth of the measurements is large compared to the frequency range of the phenomenon of interest, i.e. if the temporal resolution is too high, time domain analysis (which accept signals of any frequency) will be affected by ``unwanted signals.'' For example, if within a discontinuity, there are waves which are polarized in the normal direction of the discontinuity, they will affect the result of a Minimum Variance analysis. Although these high frequency phenomena have their own physical significance they may be treated as ``noise'' in analysis of low frequency phenomena. Therefore, in time domain analysis proper preparation of the data by either running average or filtering is important to limit the band of information to the time scale of the phenomenon of interest.
Amplitude resolution. Digitized measurements have finite amplitude resolution as well as finite temporal resolution. The result of the process is to add a square wave, or finite steps to the data. When a Fourier analysis is performed by digitized data, there is a minimum noise level set by the digitization of D2/12FN, where FN is the Nyquist frequency and D is the digital window. This digital noise is spread uniformly over the bandwidth of the signal i.e. 0 to FN and when decreasing the bandwidth of the digitized signal by subsampling or "decimating" as it is sometime called, it is important to low-pass filter the signal, or else the entire digital noise of the original time series will be added to the new narrower bandwidth. If the digital noise is comparable to the measured signal, the noise is most readily reduced by digitizing more finely than increasing bandwidth [Russell, 1972].
1.3.2 Field measurements.
One of the most accurately determined quantities in space is the magnetic field, if attention is paid to considerations of linearity, magnetic cleanliness and bandwidth relative to the Nyquist frequency. For an elliptical Earth orbiter its magnitude can be calibrated each time at perigee where the field is known and the zero levels of each sensor in the spin plane can be calibrated as the spacecraft rotates [e.g., Kepko et al., 1996]. The intercalibration among different spacecraft can be done in current-free regions [e.g., Khurana et al., 1996].
The electric field measurements are affected by spacecraft charge and Debye length. The understanding of the calibration problem has been significantly improved and the validity of the measurements has been tested [Mozer et al., 1979, 1983]. Nevertheless, special caution should be taken in the regions with significant plasma gradients and low plasma density.
High frequency (> 10 Hz) electric field fluctuations are measured with electric dipole antennae. Often only the components in the spin plane are measured. The reliability of the data is generally good. The major limitation is caused by the finite data rate which results in competition between the frequency and temporal resolutions. For example, to infer the plasma density from waves at the plasma frequency, one needs a fine frequency resolution, but to determine the polarization of a wave based on spin modulation, one needs a temporal resolution a few times the spacecraft spin period. The corresponding oscillating magnetic field is usually measured with wire coil antennae that are sensitive to the time derivative of the magnetic field. A search coil magnetometer is a coil antenna with a highly permeable core and is used at ELF and VLF frequencies, say 10 to 104 Hz.
1.3.3 Plasma moments.
A plasma detector measures the number of particles in an energy range between U and that traverse in a time an area within an element of solid angle around the normal to A. The differential direction (flux) intensity is defined as
In principle, is the cross-section of the detector, and are the angular and energy resolutions of the instrument, and is the sampling time which equals the data rate for continuous measurements. To obtain a complete distribution function, the instrument has to scan all look-directions and energy ranges. Often the angular and energy resolutions (and the signal-noise ratio) compete for the limited telemetry.
The plasma distribution function f is related to J by
where m and v are the mass and velocity of a particle being detected. With the distribution function, the plasma moments, density, velocity, pressure and heat flux (for a particular species), can be derived,
where and are the lower and higher cutoff velocities that are determined by the lower and higher cutoff energies, and , and the spacecraft potential , and are
where q is the electric charge of the particle. For a 2-D detector, an assumption is needed in order to extrapolate the values of fluxes in undetected directions. This assumption becomes important when the anisotropy of the plasma is high. Note that the effects of the spacecraft charge depend on the species and that the cutoff velocities are very large for electrons compared with ions because of their small mass.
The effects of a finite sampling energy range on the moment measurements depend on the temperature and velocity of the plasma being detected [Song et al., 1997]. If the temperature is many times lower than the higher cutoff energy, the higher cutoff has only minor effects. The lower cutoff has more extended effects. In general, see Fig 1.2, the density and pressure are underestimated, and the velocity and temperature overestimated. In order to reduce the errors, some data analysts interpolate the points below the lower cutoff or fill the hole with best-fit to a convective Maxwellian distribution. The energy resolution becomes an important issue in order to accurately derive the moments of cold plasmas. Without a good energy resolution within the bulk of the distribution, the density cannot be derived accurately.
|Fig 1.2. Ratios of the measured to real moments [Song et al., 1997]. The lower cutoff velocity should be evaluated from the lower cutoff energy of the detector and the spacecraft potential at the time. The temperatures, in eV, for = 20 eV ion detector and = 15 eV electron detector, when the spacecraft is uncharged, are given for corresponding values in the x axis. Different lines are for different bulk velocities v0 normalized by (2 T0)1/2, which are close to the Mach numbers.|
In summary, it is very important for a data analyst to understand the scheme of the algorithm with which plasma moments are derived. Caution should be taken in particular for electrons, cold plasmas, plasmas with large anisotropy using 2-D measurements, and plasmas of multiple populations.
For quantitative data analyses, the calibration of a measurement becomes essential. A calibration factor is often a function of time and plasma conditions. For example, degradation of an instrument with time requires the measurements be calibrated and recalibrated, and the calibration factor may change significantly across a shock as the plasma condition differs [Sckopke et al., 1990; Song et al., 1997].
When data analysis involves more than one instrument, intercalibration among these instruments becomes important. For example, the sonic Mach number is the ratio of two moments measured by the same detector. Even though the absolute value of each moment may not be accurate, one may suspect their ratio to be reasonably accurate. The Alfven Mach number, on the other hand, involves three measured quantities from two different instruments. The chance for error is much greater. Another interesting example is that the nonlinear response of two instrument could lead to significant differences in the same physical quantity [Petrinec and Russell, 1993] and these differences may vary with the measured parameters. One way to intercalibrate the magnetic field and plasma moments measurements is to use the force balance requirement under some known conditions, such as near a stagnation region [Song et al., 1993], see Fig. 1.3 and its caption.
|Fig. 1.3. An example of intercalibration between the magnetic field and plasma measurements. At a stagnation region, the sum of the plasma pressure and magnetic pressure should be constant, or the variations in the two should be anticorrelated with a factor of -1. The raw data penal (a) shows that the two pressure are anticorrelated but with a factor differing from -1. A calibration factor is introduced to the plasma pressure to make the slope -1. The intercalibrated plasma density is compared with the densities measured by other instruments to validate the method [Song et al., 1993].|
Quantitative comparison among measurements from different satellites requires knowledge about the above issues for all satellites involved. Without careful intercalibration, one could draw a wrong conclusion. For example, the difference in calibration for two satellites could make normal fluctuations of two physical quantities into two clusters which lead to a linear relationship. This relationship could be mistakenly interpreted as a dependence between the two physical quantities, see Fig. 1.4 for example.
|Fig. 1.4. Random uncorrelated fluctuations in quantities Q1 and Q2 measured from spacecraft A and B which are not intercalibated can sometimes be misinterpreted as a linear relation between Q1 and Q2.|
Advances in technology and accumulation of experience have made multiple-instrument multiple-satellite studies possible. Plasma density can be intercalibrated by comparing particle measurements with plasma frequency measured by plasma wave experiments or wave propagation experiment [e.g. Harvey et al., 1978]. Plasma velocity can be intercalibrated by comparing particle measurements with the field convection velocity, E xB [e.g. Mozer et al., 1983]. However, the latter velocity has components only perpendicular to the magnetic field.
1.4 Principal Axis Analysis.
Principal Axis Analysis provides the mathematical basis for the Minimum Variance Analysis of discontinuity analyses in section 2.2 and for the covariant matrix analyses of wave analyses in section 3.1. For more introductory readings, one is referred to textbooks of multivariate analysis [e.g., Anderson, 1958]. In space physics data analyses, the multivariates are often the three components of the magnetic field, .We define a so-called covariance matrix
where , and are averages of and , respectively. Similarly, the covariance matrix can also be defined in the frequency domain, or are replaced by .When , gives the cross-correlation between the two involved components of the field, and is the auto-correlation. Principal Axis Analysis provides a tool for a coordinate transformation. In the new coordinate system, the cross-correlation, ,between two components vanishes, or
where and are the transformation matrix and its inverse, and is a diagonal matrix. The magnetic field in the new coordinate system is
Mathematically, to find such a transformation is to find the eigenvector and eigenvalue of , or solve for
Because is a 3x3 matrix, there are three solutions, , and with , , . The three eigenvectors referred to as the principal axes (in rows) form the transformation matrix and the three eigenvalues referred to as the maximum, medium, and minimum eigenvalues, respectively, are the diagonal elements of .As will be discussed below, Minimum Variance Analysis (section 2.2 and 3.1) assumes to be the normal direction of a discontinuity or the propagation direction of a wave and Maximum Variance Analysis (section 2.5) assumes (for a different variable) to be the normal direction.
Loosely speaking, Principal Axis Analysis can be visualized as follows. The tip of the measured (magnetic field) vector draws points around the average field in three dimensional space due to variations. A best-fit ellipsoidal surface centered at the tip of the average field that approximates these points is then obtained. The three axes of the ellipsoidal surface are the three principal axes. The lengths of the principal axes represent the standard deviation of the field fluctuations about the average field in the three directions, and their squares are the eigenvalues. The above picture describes very well the perturbations associated with a wave. In the case of a discontinuity the field rotation across it is usually far less than 360. For example, if the field rotates 180, the field vector will vary on only one side of the maximum variation. In this case, the direction of the maximum eigenvector usually, depending on the distribution of the variations, remains parallel to the direction of the maximum variation, but the maximum eigenvalue will be different from the maximum variation. It is worth mentioning that for linearly polarized perturbations (in contrast to rotational perturbations) only one principal axis is determined and the other two have no definitive direction. This behavior causes uncertainty in data analyses and will be discussed in the corresponding subsections.
In general, in the frequency domain, the covariance matrix is complex. The Principal Axis Analysis is concerned only with the real part of the matrix. The meaning of the imaginary part of the matrix will be discussed in section 3.1.