Incoherent Scatter Analysis Techniques

Document Sample
Incoherent Scatter Analysis Techniques Powered By Docstoc
					       Incoherent Scatter Analysis
              Techniques


               EISCAT Radar School
                  24 August, 2003
                     John Holt




        Incoherent Scatter Analysis
               Techniques
I.     Introduction
II.    Data and Model Parameters
III.   Inverse Theory
IV.    Incoherent Scatter Spectrum as a Function of Basic
       Ionospheric Parameters
V.     Incoherent Scatter Measurements and Analysis
VI.    Future Trends




                                                            1
Introduction




               I. Data and Model Parameters
      A.       Some History
      B.       Voltage
      C.       Power (Lag Profiles, Spectra, ACFs)
      D.       Basic Derived Parameters
      E.       Inferred Parameters




                                                     2
J. V. Evans, TR 274, 1962




J. V. Evans, TR 274, 1962




                            3
                           J. V. Evans, TR 274, 1962




                          J. V. Evans, TR 274, 1962




1.   Above 300 km assume Ti is constant with height. The increase in the bandwidth of
     the signals is caused solely by an increase in Te/Ti.
2.   Above 300 km assume Te/Ti is constant with height. The increase in the bandwidth
     of the signals is caused by an increase in Ti.

       Option 2 was preferred based on the large diurnal variation in temperature
       inferred from satellite drag measurements.




                                                                                        4
             Perkins and Wand, 1964




     I. Data and Model Parameters
A.   Some History
B.   Voltage
C.   Power (Lag Profiles, Spectra, ACFs)
D.   Basic Derived Parameters
E.   Inferred Parameters




                                           5
                       A/D Data




     I. Data and Model Parameters
A.   Some History
B.   Voltage
C.   Power (Lag Profiles, Spectra, ACFs)
D.   Basic Derived Parameters
E.   Inferred Parameters




                                           6
                                     Power Measurements
Range




                                Delay




                    Incoherent scatter power spectrum
                                      2           2                 4
                   λ            1                     λ 
                           2

                               ∑  D  Fi (ω) N e0 (ω) +  4πD  Fe (ω) ∑ Ni0 (ω)
                                                      2                           2
               1+  
                                                                       2

                   4π        i   i                       e         i
   W (ω) =                                                                    2
                                                1  2         1 
                                                                      2
                                                                          
                                   λ 
                                           2
                                                                         
                               1+              Fe (ω) + ∑   Fi (ω) 
                                   4π         De        i  Di      
                                                                          
                                               

        ω = angular frequency displacement from transmitter frequency
        λ = radar wavelength
        De= electron Debye length = (ε0KTe/4πNe2)1/2
        Di= ion Debye length = (ε0KTi/4πNi2)1/2
                      ∞
                             16π 2 KT                  ∞
                                                                16 π2 KT 
        Fe (ω) = 1 − ω∫ exp  − 2 e τ 2  sin(ωτ)d τ − jω∫ exp  − 2 e τ2  cos(ωτ)d τ
                      0        λ me                    0        λ me   
                  ∞
                         16π2 KT                  ∞
                                                           16 π2 KT 
    Fi (ω) = 1 − ω∫ exp  − 2 i τ2  sin(ωτ)d τ − jω∫ exp  − 2 i τ2  cos(ωτ)d τ
                  0        λ mi                   0        λ mi   
                         ∞
                                16 π2 KT 
         N e0 (ω) = 2N e ∫ exp  − 2 e τ2  cos(ωτ)d τ
                 2


                         0        λ me   




                                                                                         7
         Basic Derived Parameters
 • Electron density
 • Electron temperature
 • Relative ion concentrations (N2+, NO+, N2+,
   O+, He+, H+)
 • Ion temperatures (N2+, NO+, N2+, O+, He+,
   H+)
 • Ion velocities (N2+, NO+, N2+, O+, He+, H+)




Errors in incoherent scatter power measurements

The fundamental result is
     ∆Ps    1  Pn 
         =       1 + 
     Ps    Nsamp  Ps 


Where Ps is the error in the signal power, Ps is the signal power,
Pn is the noise power and Nsamp is the number of samples.

The variances and covariances of the autocorrelation function
have the same dependence on the signal to noise ratio and
number of samples.

Later, we will look at how this dependence translates to the
errors in the ionospheric parameters deduced from the ACF.




                                                                     8
        V. Inferred Parameters
A. Ion production and loss (ion continuity
   equation)
B. Neutral wind (ion momentum equation)
C. Neutral temperature (ion energy equation)
D. Electric field (Maxwells equations)




           II. Inverse Theory
A.   Forward Problem
B.   Inverse Problem
C.   Bayes’ Theorem
D.   Least Squares




                                               9
                 Inverse problem theory

Humans were naked worms; yet they had an internal model of
the world. In the course of time up to the present, this model
has been updated many times, following the development of
new experimental possibilities (i.e., the development of their
senses) or the development of their intellect. Sometimes the
updating has been only quantitative, sometimes it has been
qualitative. Inverse problem theory tries to describe the rules
human beings should use for quantitative updatings.



                                    Albert Tarantola
                                    Inverse Problem Theory




Definition of the Forward and Inverse Problems

• Let S represent a physical system (e.g. the ionosphere)
• Assume that we are able to define a set of model
  parameters which completely describe S (e.g. ion
  composition, temperatures, drifts, density, electron
  temperature…).
• Operationally define some observable parameters whose
  actual values depend on the values of the observable
  parameters, given arbitrary values of the model parameters
  (e.g. I. S. radar spectra or ACFs).

 Forward problem: predict the values of the observable
 parameters, given arbitrary values of the model parameters.

 Inverse Problem: Infer the values of the model parameters
 from given values of the observable parameters.




                                                                  10
 The Bayesian approach to inverse problems

The proper formulation of inverse problems is obtained by
using the language of probability calculus, and by using a
Bayesian interpretation of probability. Inverse problem theory
has to be developed from the consideration of uncertainties
(either experimental or in physical laws), and the right (well-
posed) question to set is: given a certain amount of a priori
information on some model parameters, and given an
uncertain physical law relating some observable parameters to
the model parameters, in what sense should I modify the a
priori information, given the uncertain results of some
experiments.

Tarantola




                          BUT




                                                                  11
          There is another point of view

Courts have consistently held that academic license does not
extend to shouting “Bayesian” in a crowded lecture hall.
Press et al.
Numerical Recipes




     A narrow definition of inverse theory

Consider two positive functionals, A and B. A measures
something like the agreement of a model to the data (e.g. χ2)
and B measures something like the “smoothness” of the
solution. Then the central idea in inverse theory is the
prescription:

     Minimize: A + λB

For various values of 0 < λ < ∞ and then settle on a “best”
value of λ by one or another criterion.




                                                                12
       Bayesian solution to inverse problems

•   Let f(d,m) be a joint probability density on the parameters (d,m)
•   Define the marginal probability densities
       f M (m ) = ∫ d d f ( d , m )              f D (d ) = ∫ d m f ( d , m )
                      M                                         M


•   Define the conditional probability densities
                       f (d 0, m)                                f ( d , m 0)
       f M |D =                                  f D| M =
                  ∫
                  M
                      dm f (d 0, m)                         ∫
                                                            M
                                                                dd f ( d , m 0 )


•   Then
                                 f D| M ( d 0 | m ) f M ( m )
       f M | D ( m | d 0) =                                              Bayes Theorem
                              ∫ dm f   D|M   ( d 0 | m ) f M (m)




                                                                                         Monte Carlo Simulation




                                                                                                                  13
                                      Least squares
If the measurements and a priori information have Gaussian
uncertainties, then the maximum likelihood point in the probability
density representing the posterior information in the model space is
the minimum of the objective function:
               1
     S(m) =      (g(m) - dobs )t C-1 (g(m) - dobs ) + (m - m prior ) t C-1 (m - m prior )
               2                                                                         
                                   D                                     M




Where g(m) is the linear model, d is the data vector, m is the model vector, d=g(m),
The Gaussian uncertainties in the observed data are described by
the covariance operator CD and the Gaussian uncertainties in
the a priori information are described by Cm.

For uncorrelated errors:


               1  nd (g i (m) - d iobs ) 2 n m (m - mprior ) 
                                                  α    α     2

     S(m ) =     ∑                        +∑        α 2       
               2  i=1
                           (σ D )
                              i 2
                                            α =1  (σ M )       
                                                               




                Solution of least squares problem

 •       Normal Equations - possible roundoff errors
        – Gauss Jordan elimination - most common
        – Cholesky decomposition - efficient
        – LU decomposition - no covariance matrix
 •       QR decomposition - fewer roundoff problems than
         normal equations
        – Householder transforms - reflections
        – Givens transforms - rotation. Good for sequential
            addition of new data
 •       Singular value decomposition - Essential if fit
         parameters are not independent




                                                                                              14
               Nonlinear least squares
• Iterative linear least squares. Start with trial values for
  parameters and move to successively better solutions
• Requires partial derivatives of objective function with
  respect to fit parameters. May be computed numerically, or
  use analytic partial derivatives.
• If the fit is not linear enough, successive linear least
  squares solutions may not converge. Levenberg-Marquardt
  method alleviates this problem by including a term in the
  steepest descent direction.
• Multiple local minima may still be a problem. This is not
  usually a problem when fitting the incoherent scatter ion
  line.




 III. Incoherent Scatter Spectrum
        as a Function of Basic
       Ionospheric Parameters

A. Forward Problem
B. Monte Carlo Simulations
C. Inverse Problem




                                                                15
III. Incoherent Scatter Spectrum
       as a Function of Basic
      Ionospheric Parameters

A. Forward Problem
B. Monte Carlo Simulations
C. Inverse Problem




   Incoherent scatter spectrum




                                   16
17
18
19
III. Incoherent Scatter Spectrum
       as a Function of Basic
      Ionospheric Parameters

A. Forward Problem
B. Monte Carlo Simulations
C. Inverse Problem




                                   20
21
22
23
24
25
IV. Incoherent Scatter Measurements
 A.   Lag Profiles
 B.   Traditional Single-Altitude Analysis
 C.   Full-Profile Analysis
 D.   Other Methods




 Basic features of monostatic radar signal return
 • Return includes signals from a volume which extends cT/2 in the radar
   los direction, where T is the length of the transmitted waveform.
 • To achieve good height resolution, must use as short a pulse as
   possible. If the pulse is too long, the measurements are smeared over a
   range extent larger than the spatial scale of the model parameters.
 • But pulse must be long enough to provide sufficient lag extent
   (spectral resolution) to allow model parameters to be extracted with
   sufficient accuracy.
 • A further complication arises from the finite impulse response of the
   receiver, which causes the measured signal to be a time average of the
   signal returning to the radar. The impulse response can be shortened by
   increasing the bandwidth of the receiver, but this introduces additional
   noise into the measurement.
 • Modulated waveform, e.g. alternating codes, can simultaneously
   achieve good spatial and frequency resolution, but only at the expense
   of noise in the form of uncorrelated clutter from unwanted ranges.




                                                                              26
                                         Lag Profiles
 • Cross products are computed from a
   set of N receiver samples z(i) from one
   IPP
 • LPM = z(m)⋅z(n) formed from set
   z(1),z(2),…,z(N) up to maximum lag
   lmax = m-n
 • First diagonal = lag 0,
   second diagonal = lag 1, etc.
 • Products are accumulated element by
   element for a specified number of
   IPPs.
 • The integrated LPM is analyzed to
   determine ionospheric parameters
      – The products are combined to form
        ACFs or spectra which are analyzed
        individually
            • Summation rule
            • Deconvolution techniques
      – The entire matrix is analyzed
        simultaneously to produce profiles of
        ionospheric parameters vs range




         Formation of single-pulse ACFs from lag profile matrices

Triangular                                              Triangular
summation rule.                                         summation rule.
Gated samples                                           Best range resolution




Trapezoidal                                              Lag products
summation rule                                           treated as integrals
                                                         of plasma lag-
                                                         profile matrix




                                                                                27
    Another approach – full-profile analysis

  Another approach is to leave the lag-profile matrix
  matrix unmodified and use our theoretical
  knowledge of the incoherent scatter correlation
  function with lag to represent the lag variation of
  the lag-profile matrix. This brings our (very good)
  a priori knowledge about the lag-profile matrix
  into the analysis as soon as possible. Further, it
  facilitates including any additional a priori
  information about the ionospheric parameter
  profiles to contribute to the analysis.




            OASIS full-profile analysis

• We will discuss the MIT OASIS program here. GUISDAP
  also has a full-profile capability.
• We assume that:
   – The data are Gaussian random variables
   – A priori constraints can be expressed as Gaussian
     random variables
   – The model is linear in the vicinity of the solution
• Then, the problem reduces to minimizing the L2 norm of a
  vector – least squares.




                                                             28
             Components of OASIS model

1.   Model for height variation of ionospheric parameters
2.   Model for the incoherent scatter power spectrum
3.   Model for the radar measurement process, which is
     primarily a function of the transmitted waveform and
     receiver impulse response
4.   A priori constraints on the model parameters




                          1. Parameter profiles

 B-splines are a representation of piecewise polynomials

       pp( x) = p ( x )
                  i             x ∈ [ti , ti+1 ],         (i = 1, 2,..., m)
 Splines are piecewise polynomials which have n-1 continuous derivatives
 at the knots., e.g. step-functions (first-order splines), broken lines (second-
 order splines) and cubic-splines (fourth-order splines).

      1.    Nonzero only on k intervals

                 B k ( x) = 0
                   j             outside            x = [t j , t j +k ]
      2.    Always positive within their range
                 B k ( x) > 0
                   j             inside         x = [t j , t j + k ]
      3.    Well-defined normalization.




                                                                                   29
         Collision Frequency Profile

In the E-region, where ion-neutral collisions are important
we can do better. The collision frequency will have a
diffusive equilibrium height variation determined by the
masses of the dominant molecular species. In OASIS, we
actually take the height variation from the MSIS model and
scale that. The collision frequency determination is then
reduced to finding a single scalar parameter instead of a
collision frequency at each altitude.




      OASIS ionospheric parameter profile model
                                n
                   Si = ∑ aij B k ( r )
                                j
                                j =1

      where
         S = N (r ), T (r ),T ( r ), V (r ), P (r )
          ij       e        i          e   los   O


         a = B-spline coefficients
          ij

               k
         B (r )
              j


         r = range
         N (r ) = electron density
               e


         T (r ) = ion temperature
          i


         T (r ) = electron temperature
          e


         V = ion line-of-sight drift
          los


         P (r ) = %O
          O
                        +




                                                              30
           Components of OASIS model

1.   Model for height variation of ionospheric parameters
2.   Model for the incoherent scatter power spectrum
3.   Model for the radar measurement process, which is
     primarily a function of the transmitted waveform and
     receiver impulse response
4.   A priori constraints on the model parameters




            2. Plasma lag-profile matrix

• Given parameter profiles, we can calculate the
  corresponding incoherent scatter ACF at any point within
  the range of the profiles
• Doing this at a specified set of ranges yields a matrix
  whose rows are ACFs
• The ACF matrix is then normalized to include the effects
  of the single-electron cross section, electron density and
  range from the radar.
• This matrix will be denoted by σ(τ,r) in the following




                                                               31
             Components of OASIS model

1.    Model for height variation of ionospheric parameters
2.    Model for the incoherent scatter power spectrum
3.    Model for the radar measurement process, which is
      primarily a function of the transmitted waveform and
      receiver impulse response
4.    A priori constraints on the model parameters




           Incoherent scatter lag-product data
• The expectation value of the lagged product obtained by multiplying
  complex receiver samples z(t) measured by a monostatic radar at times
  t and t ′ is:
                                 ∞     ∞
           z (t ) ( z (t ') ) = A∫ dr ∫ d τ Wt ,t ' ( τ, r ) σ(τ, t )
                                 0     −∞

• The plasma lag profile matrix σ(t,τ) includes the effect of the plasma
  correlation function, single-electron scattering cross section, electron
  density and range from the radar.

• The radar ambiguity function Wt,t′(τ,r) includes the effects of
  transmitter modulation and receiver filtering.

• The radar normalization factor A depends on transmitter power,
  antenna area, system losses, etc.




                                                                             32
                    Radar ambiguity function
                   ∞
  Wt ,t ' (τ, r ) = ∫ d νWt A (ν, r )Wt 'A (ν − τ, r ),
                  −∞

  where Wt A is
                                             2r
  Wt A (τ, r ) = h(t − τ) P (τ −                )
                                             c
 and where h is the receiver impulse response and P is the
 complex envelope of the transmitted waveform. The
 function WtA(τ,r) is known as the amplitude ambiguity
 function for sampling time t and indicates how the signal
 sampled at time t is composed of signals scattered from
 ranges r and received at times t.




       Discrete form of model prediction for
              measured lag products
Replacing the integrals by an integration rule we get the
model prediction for the measured lag products
                        m2 lmax
             M ij =    ∑ ∑w
                       k = m1 l = 0
                                      jkl   σ i + k ,l (ri + k )

If we take this to be a two-dimensional Riemann sum, then
Wjkl is a matrix representation of the ambiguity function for
lag j. The range summation index k is relative to the range
index i because the shape of the ambiguity function is not a
function of range. The lag summation is over all lags for
which the lag products are not zero, that is lags less than the
length of the transmitted waveform.




                                                                   33
     Two-dimensional radar ambiguity functions for a four-pulse
            modulation and boxcar impulse response




   Two-dimensional radar ambiguity functions for a single-pulse
          modulation and wideband Butterworth filter




Note triangular
Weighting of
single-pulse|
lag products




                                                                  34
Sample two-dimensional radar ambiguity function for a 100-
              µs single-pulse modulation

                                                Waveform and filter
                                                Used to collect data
                                                described below




          Components of OASIS model

1.   Model for height variation of ionospheric parameters
2.   Model for the incoherent scatter power spectrum
3.   Model for the radar measurement process, which is
     primarily a function of the transmitted waveform and
     receiver impulse response
4.   A priori constraints on the model parameters




                                                                       35
             Constraints and other data

• Two types of constraint may be included, the radar
  calibration constant and a linear combination of splines
  and their derivative.

                                                Basically, a regularization term,
                            ∂ j Pi
         C=   ∑a
               j =0,3
                         ij
                            ∂rj
                                                But perhaps easier to relate to physics.
                                                For example, a20=1, a30=-1, else aij=0
              i =1,5
                                                yields Ti=Te


• In addition to incoherent scatter lag product measurements,
  measured values of foF2 can be included




                  OASIS least squares fit

                        ( M ij − mij ) 2        ( Di − di )2     (C − c ) 2
          F =∑                             +∑                + ∑ i 2i
                i, j         ∆mij
                                2
                                            i      ∆di2        i   ∆ci


      The first term is the weighted squared residual
      of the lagged products.

      The second term is the weighted squared residual
      of any other data.

      The third term is the weighted squared residual
      of the a priori constraits.




                                                                                           36
Fit (lines) to smooth simulated data (points) when the electron
 temperature contains a Gaussian feature narrower than cT/2


                                       300 µs pulse → 45 km resolution




                                         8.33 km wide Gaussian




 Experimental lag-profile matrix and Oasis fit to the matrix
                       300 km pulse




                                                                         37
OASIS power profile derived from measurements made with
                  640-µs single pulse




                                                Deconvolved 640 µs
                                                measurement is good
                                                match to 100 µs
                                                measurement




OASIS power profile derived from a simultaneous fit to 100-,
            300- and 640-µs measurements




                                                                      38
Ion and electron temperature profiles derived from fit to 100-,
                    300- and 640-µs data




Ion and electron temperature a posteriori probability densities




                                                                  39
Power profile derived from 100 and 300-µs observations of
      the ERIC-2 ionospheric depletion experiment




            VI. Future Trends
A. Adaptive Integration
B. Analysis on Demand
C. Direct Determination of Inferred
   Parameters from Measurements
D. Incorporation of Measurements Directly
   into Assimilative Models




                                                            40
41
42
43
          VI. Future Trends
A. Adaptive Integration
B. Analysis on Demand
C. Direct Determination of Inferred
   Parameters from Measurements
D. Incorporation of Measurements Directly
   into Assimilative Models




                                            44