SPREAD-SPECTRUM
VLF REMOTE SENSING
OF THE IONOSPHERE
A DISSERTATION
SUBMITTED TO THE DEPARTMENT OF ELECTRICAL ENGINEERING
AND THE COMMITTEE ON GRADUATE STUDIES
OF STANFORD UNIVERSITY
IN PARTIAL FULFILLMENT OF THE REQUIREMENTS
FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY
By
David Christopher Shafer
July, 1994
© Copyright by David Christopher Shafer 1994
All Rights Reserved
Abstract
A novel method for the characterization of transient, localized disturbances to
the lower
ionosphere, making use of the modulated Very Low Frequency (VLF) transmissions
used
for communication and navigation purposes, is presented in this work. These
disturbances include those produced by energetic electron precipitation from the
magnetosphere, which have been intensively studied for several years, and those
produced by other effects, such as localized heating above intense lightning
discharges.
Radio signals transmitted for other purposes, such as the stations in the 20 to
30 KHz
range operated by the United States Navy, provide useful diagnostics of these
disturbances. These signals propagate in the earth-ionosphere waveguide to long
distances from the transmitter, and are influenced by both ionospheric and
ground
conditions along the path. A transient disturbance of the ionospheric
conductivity within
a few hundred kilometers of the signal path thus appears as a transient
perturbation of
the received signal. The disturbance can have a different effect on each of the
frequencies
which make up the modulated signal, particularly when the receiver is near a
point along
the signal path where several waveguide modes interfere destructively. Near
these points,
the predicted steady-state signal strength can vary by several dB between
frequencies only
200 Hz apart, and the predicted perturbation due to changing ionospheric
conditions can
vary by as much as 1 dB across the same range. The measurement method exploits
the
finite bandwidth of these signals of opportunity (several hundred Hz in the case
of most
transmitters), employing a descision-directed technique to infer the probing
signal, and
using input/output analysis to estimate the transfer function of the propagation
channel.
With the addition of frequency as an independent variable, a direct comparison
between
experimental results and theoretical models can be made. This previously
unavailable
frequency dimension is explored by applying this new spread-spectrum analysis
method
to a data set collected in the Antarctic; the results suggest that several
effects previously
neglected in modeling lightning-associated ionospheric disturbances must now be
taken
into account, including coupling of waveguide modes within the disturbance, and
the
structure of the received field at the antenna. In addition, the steady-state
propagation
conditions often show variations with frequency which are indicative of the
existence of
several significant waveguide modes, which should further constrain ionospheric
models.
Acknowlegements
As I write this, itís hard to believe that this long, sometimes frustrating,
sometimes
distracted, and oftentimes confusing process is coming to a close. Left to my
own devices,
I would have abandoned this path a long time ago; the fact that this document
exists is
largely the result of the work, support, encouragement, and love of many other
people. I
would like to thank my advisor, Professor Umran Inan, for his belief in me and
in this
work; it is largely his efforts which have kept this work on track over the
years. The same
is true for Professor Robert Helliwell, whose enthusiasm for the field, and
formation of a
world-class Antarctic scientific effort, brought me to Stanford in the first
place, and for
which I will be ever grateful. I was also very fortunate to have Professor
Arogyaswami
Paulraj on my reading committee, and I would like to thank him for taking the
time to
discuss this work, and to evaluate this dissertation. Much of the research
history of these
fields, and the relationships among the subjects that form the foundation of
this work, are
the result of the efforts of these three men. I would also like to thank
Professor Abbas El-
Gamal, for agreeing to chair my oral defense committee in the middle of a hot
summer,
when it would have been easy to bow out, and turn to more pressing matters.
My education, both specific and general, has benefited greatly from discussions
with many of the research faculty, staff, and fellow students in the Space,
Telecommunications and Radioscience Laboratory. In particular, I would like to
thank
Don Carpenter, Ron Bracewell, Tim Bell, Tony Fraser-Smith, and Martin Walt, for
access
to the wealth of information contained in their files and their heads; Vikas
Sonwalker, for
sharing his knowledge about everything from Jovian lightning to where to find
good
Dosas; Bill Trabucco and Jerry Yarbrough, for all the fine lunch conversations,
and
always knowing how everything worked and where everything was in the lab.
My fellow students at Stanford have also been supportive and helpful, as well as
being the best spotters of a specious argument around. When I arrived at
Stanford, I had
the good fortune to have Tom Wolf as a mentor and friend, who showed me the way
around so well. Tom and Bill Burgess did much of the field work at Palmer
Station, and I
would like to thank them both for doing such a fine job of setting up and
calibrating the
VLF equipment there, so that the results of Chapter 5 actually worked out. Also,
I had a
lot of fun with Bill on returning from Palmer in 1990, and many fine Chilean
meals.
Though I never had the pleasure of working in the field with him, Juan Rodriguez
has
been a good friend, and very calming as we prepared for our orals. Much of the
field work
in this research group has been carried out by students, and I have had the
pleasure of
working with Steve Reising, Paul McGill, Alan Brown, and Dave Lauben in
particular. As
my background is in experimental work, I was most fortunate to have access to
Lee
Poulsenís waveguide propagation model, and to Steve Cummer, Sean Lev-Tov, and
the
aforementioned Juan, to explain how it worked, and what to do when it didnít.
In
addition, I have enjoyed the company and wisdom of many of my fellow students;
this
work is the better for the contributions of Tom Mielke, Jane Oh, Wai-Yeung Yip,
Sasha
Draganov, Yuri Taranenko, Jasna Djurovic, Victor Pasko, Curt Carlson, Hoc Ngo
and
Luis Sa, and I thank them for their discussions, comments, and work.
My working style could only be described (charitably) as peripatetic, and I
would
like to thank our strange crew, including Ken Britten, Inka Mimms, Dave Ellison,
Linda
Fletcher, Kim Fletcher, Aleksey Novicov, Steve Schoenberg, Nancy Serpiello,
Alisa
Koreneff, Stewart Carl, Nina Demediuk, Sarah Kelly, Karen Bollinger, Bob Smith,
Linda
Gass, Jamey Graham, Noam Eisen, Andy Butler, and Sue and Ed Smith, for hanging
in
there. This is doubly true for those friends who are also now clients, and have
been
waiting for this work to be done with bated breath.
Finally, special thanks are due to: my parents, William H. and Polly A. Shafer,
for
everything, and above all teaching me not to fear the truth; John Katsufrakis,
the 1983
Siple Station crew, Bill Trabucco, and Ev Paschal, for showing me the good
results from
hard work; Bill Scott, for sharing all these times, both good and bad, and being
there; and
Kara, for the future, and what we are becoming.
DAVID C. SHAFER
Stanford, California
July 25, 1994
The collection of the field data used in this dissertation, and the subsequent
analysis
and development which led to the bulk of the findings reported here, was made
possible
through the sponsorship of the National Science Foundation, Polar Programs
Division,
under grant DPP-9020687. The initial conception of the work was supported by
Stanford
University, through Seed Research Project 172-C-062 of the Center for Integrated
Systems. The author wishes to thank the people of the United States, and the
Board of
Trustees of Stanford University, for making this work possible.
Contents
Abstract v
Acknowlegements vi
Contents vii
List of Tables x
List of Illustrations xi
1. INTRODUCTION 1
1.1. MOTIVATION
2
1.2. PROBLEM
DESCRIPTION 10
1.3. CHAPTER
DESCRIPTIONS AND CONTRIBUTIONS
12
2. THEORETICAL
BASIS 15
2.1. MEASUREMENT
OF AMBIENT ELECTRON DENSITY PROFILES
15
2.2. MEASUREMENT
OF ELECTRON DENSITY AND
TEMPERATURE PERTURBATIONS 23
2.3. USE OF SPREAD-SPECTRUM SIGNALS
AS MULTI-FREQUENCY
PROBES 24
3. PROCESSING
ALGORITHMS 26
3.1. FSK SIGNALS
28
3.2. MSK SIGNALS
29
3.3. PROCESSING
OVERVIEW 34
3.3.1. Filtering
and Clipping of Broadband Data 35
3.3.2. Demodulation
of Broadband Data 38
3.3.3. Detection
of Transmitted Bit Sequence
41
3.3.4. Reconstruction
of the Transmitted Signal
43
3.4. COMPARISON WITH RECEIVED SIGNALS 44
3.5. SMOOTHING
AND DISPLAYING RESULTS
49
4. SIMULATION
OF PROCESSING SYSTEMS
55
4.1. PROCESSING OF SIMULATED MSK SIGNAL 56
4.1.1. Construction
of the Test msk Signal
57
4.1.2. Simulation
of the Slowly-Varying Propagation System
59
4.1.3. Use of
Spread-Spectrum Processing
65
4.2. PERFORMANCE OF EXISTING AMPLITUDE AND
PHASE
RECEIVERS 66
4.2.1. Amplitude Receivers
67
4.2.2. Phase
Receivers 68
5. APPLICATION
TO DATA 73
5.1. SIGNAL PATHS STUDIED 73
5.2. RECEIVING,
RECORDING AND PROCESSING OF DATA
78
5.2.1. Antenna
Efficiency 80
5.2.2. Preamp
Frequency Response and Noise Figure
82
5.2.3. A/D
Converter and Storage System
85
5.2.4. Noise Calibration Measurements and
Dynamic Range 86
5.3. PROCESSING
SYSTEM 87
5.4. MSK SIGNALS ON EXPECTED
SINGLE-MODE PATHS 88
5.4.1. Time-Averaged
Propagation Channel Transfer Functions
88
5.4.2. LEP-Induced
Perturbations to Transfer Functions
94
6. SUMMARY
AND SUGGESTIONS FOR FUTURE WORK 101
6.1. SUMMARY
OF RESULTS AND CONTRIBUTIONS
101
6.2. SUGGESTIONS
FOR FUTURE WORK 102
6.2.1. Impulsive Noise Removal for Broadband
Recordings 103
6.2.2. Implementation of Real-Time Processing 103
6.2.3. Measurements
at Separated Stations
106
6.2.4. Studying
Shorter, Overland Paths 106
6.2.5. Use of
Other Transmitters 106
BIBLIOGRAPHY 108
List of Tables
1. Transmitter
and Arrival Characteristics 74
2. Processing
Rate Requirements for Field Receiver
105
List of Illustrations
1. Example of
Lightning-induced Electron Precipitation (LEP) events 3
2. Scattering
of energetic electrons by a whistler-mode wave launched by
lightning 6
3. Perturbation
of subionospherically-propagating VLF or LF signals by
energetic electrons 7
4. Decomposing
a wave into a sum of normal modes 11
5. System model used in this work
13
6. Nighttime
ambient profiles of electron density vs. altitude
16
7. Predicted
signal amplitude as a function of distance along the path, for a
fixed ambient electron density
17
8. Predicted
signal amplitude as a function of distance along the path, for a
fixed frequency (26.00 KHz), and two ambient electron density profiles.
18
9. Predicted
signal amplitude as a function of frequency, for the NSU-to-
Stanford path 19
10. Predicted
signal phase as a function of frequency, for the NSU-to-
Stanford path 20
11. Predicted
measured signal phase as a function of frequency, for the
NSU-to-Stanford path 21
12. Comparison of measured propagation
channel transfer function with
the prediction, for the NPM-to-Palmer path using electron density
profile II 22
13. Predicted
amplitude of the propagation channel transfer function
perturbation, as a function of frequency, on the NSU-to-Stanford path
near 26.0 KHz 23
14. Common VLF modulation
methods 27
15. Power Spectral Density of MSK-modulated signal
30
16. Operation of an MSK modulator,
implemented by phase-shift keying
and amplitude modulation of two coherent carriers
31
17. MSK waveform, illustrated for a specific bit sequence
32
18. Phase of the MSK-modulated signal, relative to an
unmodulated carrier
at the same center frequency
33
19. Overview
of the processing steps used in this work, for MSK signals.
34
20. Effect of initial
processing steps on the broadband data 35
21. Amplitude
and phase response of the high-pass filter used as the RF filter
36
22. Setting
of the clipping threshold
38
23. Amplitude and phase of baseband filter
used in the processing
39
24. Amplitude and phase output from
baseband stage. 40
25. Reconstruction of the transmitted phase
from the observed phase of the
decimator output 42
26. Comparison
of output to input fourier transforms to yield transfer
function estimate 45
27. Output of the propagation
channel transfer function estimate, for a
sample data set 49
28. Smoothed
sample data, based on the same data set as Figure 27 50
29. Estimate
of the average transfer function of the propagation channel for
a 4 minute data set, over ±300 Hz in frequency
51
30. Spectrogram
of the propagation channel transfer function residual, after
removing the average transfer function estimate from the data
53
31. An expansion of two
overlapping transfer function perturbations from
Figure 30 54
32. Processing
flow for the simulation
56
33. Autocorrelation of a sample received
bit sequence 58
34. Location of simulated VLF transmitter
("NSU"), and the simulated
ionospheric perturbation on the great-circle propagation path to a
receiver at Stanford
59
35. Electron
density profiles used in the simulation
60
36. Mode structure of the
received signal at Stanford, from the fictitious
26.0 KHz transmitter in Nebraska
61
37. Amplitude and phase of transfer
function for propagation at 26 KHz
from "NSU" to Stanford, under the ambient conditions
62
38. Amplitude of propagation from NSU to SU
as a function of frequency
63
39. Phase of propagation from
NSU to SU as a function of frequency
64
40. Spectrum
of the broadband VLF, with the synthetic NSU signal added at
26.0 KHz 65
41. Spectrogram
resulting from processing the simulated NSU signal
66
42. Typical
block diagram of a frequency-doubling MSK phase receiver,
with amplitude output
69
43. Waveform of MSK signal and an FSK
signal synthesized from it by a
square-law frequency doubler
70
44. Spectra
of MSK signal and frequency-doubled FSK signal generated
from it 71
45. Great-circle propagation
paths for the signals studied in this section
74
46. Mode structure of the arriving VLF signals at Palmer
Station 75
47. Predicted propagation from NPM in
Hawaii to Palmer Station
76
48. Predicted propagation from NSS in
Maryland to Palmer Station
76
49. Predicted propagation from NAA in Maine
to Palmer Station 77
50. Predicted propagation from NLK in
Washington to Palmer Station
77
51. Diagram of the receiving and recording
equipment used to acquire the
data analyzed in this chapter. 79
52. Amplitude
response of the Palmer antenna and preamp
83
53. Phase response of the
Palmer antenna and preamp
84
54. Sensitivity
of the Palmer antenna and preamp
85
55. Data recording system noise
and signal detectability 87
56. Average
propagation channel transfer function and 95% confidence
limits for the NPM-to-Palmer path
89
57. Average propagation channel transfer
function and 95% confidence
limits for the NSS-to-Palmer path 89
58. Average
propagation channel transfer function and 95% confidence
limits for the NAA-to-Palmer path
90
59. Average
propagation channel transfer function and 95% confidence
limits for the NLK-to-Palmer path
90
60. Effect of mode structure of incoming wave on received
signal, for a small
loop antenna 92
61. Average
transfer function of the propagation channel from NPM to
Palmer, as observed on the E/W antenna
94
62. Two transfer function perturbations on the NPM-to-Palmer
path, 95
63. An expansion of the time axis in Figure 62, showing a
single LEP event 96
64. Four transfer function perturbations on the NSS-to-Palmer
path 97
65. An expansion of the time axis in Figure 64, showing two
LEP events 98
66. Four transfer function perturbations on the NAA-to-Palmer
path,
coincident with those on NSS
99
67. An expansion of the time
axis in Figure 66, showing two LEP events
100
68. Possible
implementation of the processing described in this dissertation,
in a form suitable for use in the field 104
1.
INTRODUCTION
Our direct experience shows that the earthís atmosphere is far from a steady
environment; near the surface, weather conditions can change dramatically and
rapidly.
These changes are driven by thermal differences between parts of the Earth, and
depend
on the presence of a substantial atmosphere. As one ascends, however, the
atmosphere
thins out, and the thermal gradients which drive the weather near the surface
become
smaller. Weather phenomena familiar to surface observers rarely reach altitudes
greater
than 10 km or so. At the highest altitudes attainable with aircraft and
balloons, about
40 km, "weather", as we know it on the ground, is non-existent, and
one might expect
that the atmosphere would simply continue to thin out, and blend into the vacuum
of
space.
In fact, this does not happen at all, because ionizing radiation from space,
primarily
in the form of the ultraviolet light from the Sun and cosmic rays from deep
space, strips
electrons from their atoms, and renders the upper reaches of the atmosphere
electrically
conductive. The ionosphere and magnetosphere, as these regions are known,
exhibit a
different type of "weather", in which electromagnetic forces play a
central role. The
ionosphere conventionally begins at about 60 km, and extends up to about 500 km.
At
ionospheric altitudes, particularly the lower regions which are the concern of
this work,
the rate of collisions for the free electrons is high, and chemical processes
are important.
At higher altitudes, the collision rate drops along with the density, and the
medium
begins to be dominated by the motion of the charged particles in the Earthís
magnetic
field; this latter region is known as the magnetosphere.
In the upper ionosphere and magnetosphere, it is possible to make in situ
measurements from satellites, because the atmospheric density is low enough to
allow a
satellite to orbit for a reasonable amount of time. Below about 200 km, however,
the
atmospheric drag is high enough that satellite orbits decay rapidly, so that the
only way to
make in situ measurements of this region is with the use of sounding rockets.
These
measurements can be important, but repetitive measurements are expensive, and it
is
difficult to distinguish between spatial and temporal variations.
Fortunately, we can make use of the electrical conductivity of the lower
ionosphere
to make indirect measurements, by utilizing the fact that the medium reflects
low-
frequency radio waves. The frequencies of interest here are in the VLF and LF
range of
the electromagnetic spectrum, covering the ranges of 3-30 KHz (VLF), and 30-300
KHz
(LF). Radio waves at these frequencies are used for communication and navigation
purposes, and propagate over long distances, being efficiently guided in the
earth-
ionosphere waveguide. Thus, they are highly sensitive to the ionospheric
conditions along
their paths, and form the focus of this work.
1.1. MOTIVATION
Interest in radio propagation in the VLF and LF bands dates back many years, and
several
studies were made to characterize the propagation channel [Davies, 1990, pp.
380-399;
Field and Lewinstein, 1978; Turtle et. al., 1989]. These studies were undertaken
to support
communication and navigation uses of the VLF and LF bands, and thus generally
concern the measurement or modeling of the amplitude and phase of propagation,
or the
atmospheric noise, on coarse frequency and time scales. This type of approach
has been
sufficient for the intended uses - the design and evaluation of VLF and LF
communication systems - but is not well suited to the study of localized and
transient
effects, which are of interest here.
The earth-ionosphere waveguide constitutes a dynamic channel for the
propagation of VLF and LF radio waves, since the waves are sensitive to the
boundary
conditions at both the earthís surface and the ionosphere, and the electrical
characteristics of the ionosphere vary widely. The electron density and
temperature in the
ionosphere result from a balance between the ionization of neutral atoms and
molecules
by energy input from a number of sources, and the removal of the free electrons
by
attachment and recombination processes, resulting in substantial variations over
many
time and distance scales. The main variation in energy input comes as a result
of the
change in the flux of ionizing radiation from the sun, which causes strong
diurnal
variations in the lower ionospheric layers. The solar UV energy during the
daytime
increases the equilibrium ionization levels, and, particularly at the lower
elevations where
the collision rate is high, the ionization levels drop soon after sunset. This
diurnal
variation is important for propagation of VLF and LF radio waves, since these
frequencies
are more sensitive to changes in the ionosphere at lower altitudes than are
higher
frequencies.
Figure 1. An example of Lightning-induced Electron Precipitation (LEP) events,
observed
on the phase of a subionospherically-propagating VLF signal originating in
Hawaii, and
received at Palmer Station, Antarctica. The top two panels show the
"carrier phase" of the
VLF signal at 23.3-23.5 KHz at different time scales, while the bottom panel
shows the
spectrogram of the broadband VLF received at the same time, displayed at the
same time
scale as the middle panel. The ~5 advances
in phase, lasting ~1 second, are LEP events,
and are seen to be associated with the arrivals of magnetospheric whistlers. The
~3 dips in
the phase visible in the expanded scale are artifacts due to bit errors made in
pairs by the
MSK decoder.
At night, when
the background ionization at the lower ionospheric altitudes is low,
other energy sources can have a strong effect on the electron density and
temperature at
these altitudes. Many of these sources are transient in time, and localized in
space, and
are results of complex processes occurring in near-earth space. An important
example of
a transient propagation effect, which has recently emerged as an significant
process in the
D region and for its coupling to the magnetosphere [e.g., Burgess and Inan,
1993], was
first observed on VLF signals propagating to the Antarctic [Helliwell et al.,
1973], and is
illustrated in Figure 1. Such transient VLF perturbation events, originally
named "Trimpi
Events" after the Stanford researcher who first observed the coincidence
between these
VLF signal perturbations and magnetospheric whistlers, are also referred to as
LEP
(Lightning-induced Electron Precipitation) events. These transient perturbations
are
believed to be due to secondary ionization enhancements produced by bursts of
precipitating energetic electrons. The electrons are precipitated as a result of
interactions
in the magnetosphere with whistler waves originating in lightning discharges
(Figure 2).
LEP events have received much attention in recent years, both in the form of
experimental work [e.g., Wolf and Inan, 1990; Inan et. al., 1988b; Dowden and
Adams,
1989; Smith and Cotton, 1990], and theoretical studies [e.g., Poulsen et. al.,
1993a,b;
Cotton and Smith, 1991]. Although all the VLF perturbations analyzed in this
work are of
this type, there are several other known causes for sudden changes in VLF signal
amplitude and phase, including excess ionization in the D region due to intense
heating
from lightning flashes [Inan et. al., 1991; 1993; Taranenko et. al., 1992],
direct heating of
the D region by powerful VLF, LF and HF transmitters [Inan, 1990; Rodriguez et.
al.,
1992], and heating and excess ionization following high-altitude nuclear
explosions
[Allcock et. al, 1963; Helliwell and Carpenter, 1963].
An LEP event occurs when high-energy electrons trapped in the magnetosphere
interact with a whistler-mode wave launched by a lightning flash, are thereby
scattered
out of their stably-trapped orbits, and ëprecipitateí into the lower
ionosphere. There, they
collide with the neutral constituents, and deposit their energy in the form of
heat, X-rays,
ionization, and light. The source of this whistler-mode wave most commonly is
the radio
wave radiated from the impulsive current contained in a lightning flash [Uman,
1987,
p. 114], but other sources are possible as well. Much of this energy remains
confined to
the waveguide formed by the conducting surfaces of the Earth and the ionosphere,
but
some of the energy in the VLF frequencies can couple through the ionosphere, and
into
whistler-mode waves in the magnetosphere [Helliwell, 1965].
In the presence of enhancements to the ionization in the magnetosphere, aligned
with the magnetic field lines, the whistler-mode wave can be guided along the
field lines,
and can again couple through the ionosphere to the ground in the geomagnetically
conjugate hemisphere [Helliwell, 1965]. The whistler-mode wave is circularly
polarized,
and the index of refraction of the magnetized plasma for the whistler mode is
large, so the
velocity of propagation of the wave is relatively slow (typically 0.01-0.1 c,
where c is the
speed of light), so that it may take several seconds for the wave to propagate
back to the
Earthís surface. Furthermore, the medium is dispersive, with the lowest
frequencies
having the lowest propagation velocities, and the impulse which entered the
magnetosphere is observed in the conjugate region as a falling tone spread over
one
second or so; this characteristic accounts for the name of the propagation mode,
i.e., the
whistler.
While in the magnetosphere, the slowly-moving wave encounters a population of
trapped electrons, which execute a helical motion in the magnetic field. The
transverse
component of this motion occurs at the gyrofrequency,
, given by:
(1.1)
where e is the electronic charge, B is the local magnetic field, and m is the
electron mass
[Ramo et al., 1984, p. 705]. The motion consists of two velocity components: , which is
parallel to the magnetic field, and ,
which is perpendicular to the external field. The
sense of the electronsí motion is such that
is opposite to the direction of propagation
of the whistler-mode wave, while is
in the same direction of rotation as the
perpendicular components of the wave electric and magnetic field. Thus, the
electrons
can experience a resonant interaction with the wave, called cyclotron resonance
[Dungey,
1963], when:
(1.2)
where:
w is the frequency of
the whistler wave component,
is
the whistler wave number component parallel to the local magnetic field, and
is
the parallel component of the electron velocity.
Near the geomagnetic equator, where the local magnetic field exhibits a minimum,
electrons can remain in resonance with a particular wave frequency for an
extended time,
and the wave-particle interaction is most efficient for fixed-frequency signals
[Helliwell,
1970]. For whistler signals with varying frequency, the strongest interaction
occurs at
points along the field line which are as much as 10-12
away (in latitude) from the
equator [Inan et. al., 1989]. A common result of the interaction is pitch angle
scattering of
the electrons, where the angular momentum of the particles are redirected due to
their
interaction with the wave, and the angle between the electronsí velocity
vectors and the
magnetic field are changed. While individual electrons may have their pitch
angles either
increased or decreased, there is a minimum pitch angle below which the electronsí
mirror
points will be in the lower, denser, regions of the ionosphere. Those electrons
which have
their pitch angles increased will remain trapped in the magnetosphere, but some
of those
which experience a pitch angle decrease will be scattered into the loss cone,
which is the set
of pitch angles that correspond to mirror points in the lower ionosphere. As
these
scattered electrons subsequently travel down the field line, they are set to
mirror at a
lower altitude, and have a greatly increased probability of collision with
neutrals in the
ionosphere [Inan et. al., 1985], as illustrated in Figure 2.
Figure 2. Scattering of energetic electrons by a whistler-mode wave launched by
lightning.
The lightning flash on the Earthís surface releases an intense impulse of
electromagnetic
energy, some of which is in the VLF region of the spectrum and can penetrate the
ionosphere and travel along the magnetic field in the whistler mode. This wave
can
interact via cyclotron resonance with energetic electrons trapped by the
magnetic field,
which can amplify the wave, and lower the mirror points of the electronsí
orbits. In order
to show both the ionospheric and magnetospheric events in this figure, it is not
to scale.
Figure 3. Perturbation of subionospherically-propagating VLF or LF signals by
energetic
electrons. When electronsí mirroring points are sufficiently lowered (for
example, by the
mechanism in Figure 2), collisions with neutral ionospheric constituents cause
the
electronsí energy to be distributed into the ionosphere, and create a region of
secondary
ionization. This appears to the VLF or LF wave as a change in conductivity of
the earth-
ionosphere waveguide boundary, which changes the propagation characteristics
through
the perturbed region. This change in propagation characteristics causes the
signal
observed by a suitably-located fixed receiver to undergo an amplitude and phase
change,
illustrated here as a phase advance and amplitude decrease by the thickness and
position
of the wave crests.
The incoming >100 KEV electrons carry a much larger energy than is required
to
ionize neutral constituents of the ionosphere, and can create a large patch of
increased
ionization at altitudes from 40-90 km [Rees, 1963]. This excess ionization
changes the
conductivity of the affected ionospheric regions, and alters the propagation
characteristics of VLF and LF waves in regions under or near the disturbance [Galejs,
1972; Poulsen, 1991], as illustrated in Figure 3. It is these changes in
amplitude and phase,
observed with an appropriate receiver, that contain information about the
ionospheric
disturbance, and which are examined in this work.
Following the termination of the burst of precipitating electrons, the disturbed
patch of the ionosphere is no longer in equilibrium, and chemical processes
(recombination and attachment) cause the excess ionization to disappear over a
period of
time [Glukhov et. al., 1992]. This loss of the excess free electrons over time
leads to a
recovery of the perturbed signal levels over a period of 10-100 seconds, and
accounts for
the characteristic form of the disturbance, as seen in Figure 1. The details of
this recovery
process depend on the chemistry in the D region [Glukhov et. al., 1992], and the
resulting
signatures are expected to be an important source of information about the D
region
[Pasko and Inan, 1994].
Previous work on LEP events and associated VLF perturbations have concentrated
on establishing the association of VLF signatures with whistlers [Carpenter and
LaBelle,
1982; Inan and Carpenter, 1986] and lightning discharges [Inan et. al., 1988a;
Yip et. al.,
1991], on the determination of their occurrence characteristics [Carpenter and
Inan,
1987], geographic and geomagnetic distribution [Inan et. al., 1988b; Freidel and
Hughes,
1993], and geomagnetic conjugacy [Burgess and Inan, 1990]. Especially in the
light of a
new realization of the possible importance of this phenomenon for radiation belt
loss
processes on a global scale [Burgess and Inan, 1993], there is a need to extract
quantitative
information from the observed VLF signatures of LEP events. This task requires
the use
of an accurate model of the propagation of VLF waves in the earth-ionosphere
waveguide, and their scattering from localized disturbances.
In the most general case, the propagation of radio waves, and their interactions
with the surfaces of the Earth and the ionosphere, is accurately described by
Maxwellís
equations. However, it is impossible to solve the equations in closed form for
realistic
boundary conditions, and we must turn to approximate, computer-aided methods.
One
approach for the formulation of the propagation problem is to consider all the
possible
paths from the transmitter to the receiver, which satisfy the boundary
conditions at all
points. This leads to an infinite sum of ërayí paths, though the higher-order
rays
contribute less to the total solution, because they are attenuated more at each
ionospheric
reflection, and travel over longer paths [Morfitt and Shellman, 1976]. This ray
approach is
most workable for shorter paths, because the angles of incidence for the rays
rapidly
become small for short separations between transmitter and receiver, with
correspondingly higher attenuation, so that the solution may be written in terms
of only
a few rays. However, for longer paths, many significant rays exist, and because
the path
geometry must be recalculated for each transmitter-receiver separation, the ray
approach
becomes cumbersome.
An alternative approach involves the expansion of the propagating wave in an
infinite sum of normal ëmodesí [Budden, 1985]. The modes are the
eigenfunctions of
Maxwellís equations with the appropriate boundary conditions applied, so that
the
propagation of a single mode (i.e., the evolution of the amplitude and phase of
the mode
as it propagates along the path), in a waveguide with homogeneous boundary
conditions,
can be expressed in terms of a single propagation constant. Because Maxwellís
equations
are linear, an arbitrary input excitation can be expanded as an infinite sum of
such
modes, which can then be propagated independently to the receiver, as long as
any
changes to the boundary are slowly-varying, and satisfy the WKB approximation
[Budden, 1985]. This process is illustrated in Figure 4, which shows the
variation with
height of a field component at the input to the system, and its decomposition
into a sum
of two waveguide modes. The modes propagate independently to the receiver, where
they
can be summed to give the structure of the output wave. Since, in general, the
attenuation
(i.e., the amplitude variation) and propagation speed (i.e., the phase
variation) is different
for different waveguide modes, the complete output wave need not look like the
input in
structure. Coupling between modes can occur at discontinuities (i.e., changes
occurring
over distances of under a wavelength) in either the ionosphere or ground
conductivity,
and can be formulated by imposing the electromagnetic boundary conditions at
such
interfaces [Morfitt and Shellman, 1976].
For long ( km) propagation paths
over seawater and at VLF frequencies, the
effect on the VLF signal of a disturbance in electron density,
, can be modeled using a
single dominant waveguide mode [Inan and Carpenter, 1987]. In such cases, the
received
signal perturbations are predicted to be in the form of an amplitude decrease
and a phase
advance, which qualitatively compares well with the experimental results for
such paths
[Wolf and Inan, 1990]. For paths where more than one waveguide mode is
significant, the
situation is more complex, and we must turn to computer models for predictions.
Two-
dimensional models, such as that of Tolstoy and Rosenberg [1983], can reproduce
some
of the observed features, but for studies of the effects of ionization patches
of finite size
and occurring in general on or off the great-circle propagation path, a full
three-
dimensional model must be employed. Such a three-dimensional model of VLF earth-
ionosphere waveguide propagation was recently developed, and is utilized in this
dissertation [Poulsen et. al., 1993a, 1993b]. The model is based on the
two-dimensional
VLF earth-ionosphere waveguide propagation model known as the Long Wave
Propagation Capability (LWPC) [Ferguson et. al., 1989].
The three-dimensional model, and the underlying LWPC code, is in turn based on
two approximations to the full solution of the modal propagation. First, the
propagation
of the waves from the transmitter to the receiver, from the transmitter to the
scattering
center, and from the scatterer to the receiver, are all assumed to satisfy the
WKB
approximation along a finite set of path segments. That is, the paths are broken
into a
finite number of "slabs", and coupling between propagating modes is
ignored, except at
the boundaries between slabs. The slab boundaries allow variations in ground and
ionospheric conductivities (due to changes in electron density, temperature, and
magnetic field) to be taken into account. Second, the change to the refractive
index in the
perturbed region is assumed small, so that the Born approximation would be
applicable.
This involves the use of the unperturbed field inside the disturbance, rather
than a partly-
scattered field, to determine the total field by summing the contributions of
fields
scattered by each small area of the disturbed region. These assumptions are
justified for
the disturbances of geophysical interest in our work, as well as for the ambient
ionospheric conditions considered here [Poulsen et. al., 1993b].
1.2. PROBLEM
DESCRIPTION
To determine the electrical properties of the ionosphere, and their variations
over
time, we have at our disposal several relevant tools. First, the 3-D waveguide
propagation
model code described above can accurately predict, for a given set of
ionospheric
parameters, the received signal amplitude and phase at any frequency. In
addition, there
is available a network of VLF transmitting stations, operating in the 10-50 KHz
range,
with radiated powers up to 1,200 KW. These VLF transmitters exist for other
purposes,
typically navigation and communication, and the strongest signals are useful for
these
purposes worldwide. While the transmitters employ a variety of modulation
schemes, the
three types that will most concern us here are the communication transmitters
employing
FSK, or, more commonly, MSK modulation, and the Omega navigation transmitters
using a pulsed-cw format, with each station transmitting pulses at several
frequencies
sequentially.
Figure 4. Illustration of the process of decomposing a wave into a sum of normal
modes.
The input wave and output wave are shown schematically at the top panel. The
variation
of field strength with altitude is complex, and changes as the wave propagates.
In a modal
decomposition, the input wave is written as a sum of a (possibly infinite) set
of normal
modes, as illustrated in the bottom panel. The normal modes propagate with only
an
amplitude and phase shift if the medium is homogeneous, and linearity allows the
output
wave to be written as the sum of the modified normal modes.
The measurement of the temporal variations in amplitude and phase of waves
which originate at a transmitter and arrive at the receiver over a given
propagation path,
can be viewed as a system identification problem, with the input consisting of a
modulated signal covering several hundred Hz of bandwidth, and a center
frequency in
the 20-30 KHz range. The approach taken here involves treating the radiation,
propagation, and reception of the modulated signal as a "gray box",
and, using the
modulated VLF signals as a probe, attempting to identify parameters of the
resulting
system model which can be directly compared with the propagation model. Assuming
the
source excitation remains steady, the measured temporal variations are a direct
result of
variations in the characteristics of the propagation path. The system model used
for this
is illustrated in Figure 5.
The receiving system is well-characterized, and is both linear and
time-invariant, so
its transfer function can be transposed with the noise component additions,
which
reduces the system to a known input ,
followed by a linear, but slowly time-varying,
system. The system output is
formed by adding to this the noise components,
transformed by the action of the receiver. The effects of the propagation can be
isolated,
because the time scales over which significant variations occur in the
ionosphere is quite
different from that of the other variations in the system. Variations in the
transfer
functions of the transmitter and antenna, and of the receiver, are much slower
than
variations caused by the ionospheric disturbances of interest. The time scales
for these
variations is measured in hours or more, while significant transient changes in
the
propagation path are seen on time scales ranging from less than a second to many
minutes.
By using the VLF probe signals, the transfer function from the modulator output
to
the receiver output may be measured. The propagation model referred to earlier
gives
theoretical information about the propagation effects at selected frequencies in
the range
of interest. These two may be directly compared, to test and refine estimates of
the
parameters which govern the VLF wave propagation, and which are of geophysical
interest.
1.3. CHAPTER
DESCRIPTIONS AND CONTRIBUTIONS
This work is organized into six chapters, which describe the background of the
work, the
techniques developed for analysis, and the results of the application of these
techniques to
the data set. Below are descriptions of the remaining five chapters, and a
summary of the
contributions they describe.
Chapter 2 describes the theoretical background for this work, and presents
results
of simulations of VLF propagation in the earth-ionosphere waveguide. The model
calculations show that, under realistic ionospheric conditions, the propagation
of the
VLF waves is sensitively dependent on the electron density profile in the
ionosphere, and
that this sensitivity may manifest itself as a strong dependence on the wave
frequency of
the signal amplitude and phase at a fixed receiver location.
Figure 5. System model used in this work. The transmitter is represented by the
top two
blocks, with the unknown factors, which must be estimated from the data or
elsewhere,
listed as inputs. The propagation model is the 3-D propagation code, and
estimating the
ionospheric parameters listed as its inputs is the ultimate goal of this type of
work. The
noise sources represent the additive atmospheric noise, and the receiver model
accounts
for effects in the receiver and processing, which can be measured in advance.
Chapter 3 describes a means of using the modulated signals from VLF and LF
transmitters to probe this frequency dependence, by decoding the message
information
contained in the modulated transmission, and feeding this information back to
estimate
the transmitted signal. This estimate is then compared in the frequency domain
to the
received signal, resulting in a measurement of the transfer function of the
propagation
channel from transmitter to receiver, over a finite frequency range, with high
time
resolution. The contributions of the work described in Chapter 3 include:
* Development
of a new VLF remote sensing technique, which fully utilizes the
characteristics of modulated signals to determine substantially more
information about propagation conditions than techniques currently in use.
* Introduction
of frequency as an independent variable in VLF propagation
measurements and calculations, providing additional experimental tests for
theoretical models of geophysical effects.
Chapter 4 presents the results of a simulation of this processing, as applied to
a
short, overland propagation path. The results of this simulation show that:
* The technique
described in Chapter 3, in the presence of a realistic atmospheric
noise environment, can detect multi-mode propagation and mode coupling
when these effects are important.
Chapter 5 uses the new technique to analyze a data set obtained from Palmer
Station, Antarctica, which was recorded for this purpose digitally with high
bandwidth
and dynamic range. A number of LEP perturbations were recorded and analyzed,
with
the result that:
* Currently
available propagation models give results consistent with the largely
single-mode propagation studied.
The final Chapter of this work, Chapter 6, summarizes the above results. In
addition, since the method used to undertake this study, involving
high-resolution,
broadband recording and post-analysis, is difficult to apply where large amounts
of field
data are desired, the final contribution of this work is:
* Consideration
of the design of a field receiver implementing this technique,
which allows reduction in interference and sensitivity to transmitter conditions
compared with current designs.
2.
THEORETICAL BASIS
In order for a frequency-dependent radio diagnostic technique to yield useful
information on geophysical processes, the propagation paths and/or the
scattering of
energy must exhibit discernible differences over the band of frequencies
illuminated by
the transmitted signals. These differences may be due to many conditions,
including the
ambient electron density profile, the ambient electron collision frequency
profile, or the
perturbations from the ambient levels which occur due to electron precipitation
or
heating. In this section, we investigate the frequency dependence of the
propagation
characteristics of VLF signals over the roughly 300 Hz bandwidth illuminated by
the
available transmitters. For this purpose, we utilize the multiple-mode,
three-dimensional
model of VLF propagation in the earth-ionosphere waveguide in the presence of
localized
D-region disturbances [Poulsen et. al., 1993a], hereafter referred to as the 3-D
waveguide
propagation model.
2.1. MEASUREMENT OF AMBIENT ELECTRON DENSITY
PROFILES
An important input to the 3-D waveguide model is the altitude profile of
electron
density and collision frequency, which must be specified for each segment along
the
propagation path. Typical values of electron density, observed at night at
mid-to-low
latitudes, are encompassed by the three different profiles in Figure 6; these
representative
profiles were also used in other recent work [Inan et. al., 1992]. With the
electrical
conductivity of the ionosphere, which is a function of the electron density ( )
and
collision time ( ), as defined, and using a tabulated data base for the
conductivity of the
ground over the globe, the 3-D waveguide propagation model can be used to
determine
the amplitude and phase as a function of frequency and distance along any great
circle
path between a transmitter and receiver. The path chosen for our analysis is the
same as
that used later in Chapter 4, and runs from a fictitious transmitter at 41.5
N, 97.6 W, to
a receiver at Stanford (37.4 N, 122.2 W),
operating at a center frequency of 26.0 KHz.
Figure 6. Three nighttime ambient profiles of electron density vs. altitude,
which have
been used in recent past work [Inan et. al., 1992], and which encompass a
1:10:100 range
over the important altitude range of 75-100 km. Profile I represents a tenuous
ionosphere,
Profile II represents a normal density, and Profile III represents an increase
in density at
the reflection altitude of a factor of about ten.
Figure 8 shows the signal strength at the receiver, as computed using the 3-D
waveguide propagation model, for the case where the density profile II shown in
Figure 6
was assumed to be in effect all along the propagation path. One feature of the
propagation, which is well-illustrated in Figure 8, is the occurrence of sharp
nulls in the
amplitude of the received signal, as a function of distance from the transmitter
along a
great-circle path. These nulls are due to the fact that the wave energy
propagates in
several waveguide modes simultaneously, and each of these modes has a different
propagation constant. If no conversion of energy from one mode to another
occurred,
and the modes all had the same attenuation, the nulls would repeat every time
destructive
interference occurs between two or more modes with similar amplitudes, due to
differences in the propagation constants. For VLF waves propagating in the
earth-
ionosphere waveguide, the attenuation rates of different modes differ [Poulsen
et. al.,
1993a], and the nulls do not appear in a simple, repetitive pattern. Even though
the
attenuation of the waveguide modes differ, there are still places where two
modes have
nearly the same amplitude, since certain modes having higher attenuation rates
are
excited more efficiently by typical transmitting antennas. This can cause sharp
nulls near
the transmitter, with a transition to a smoother characteristic as many of the
modes are
attenuated with distance. The path illustrated in Figure 8 is relatively short,
at about
2,100 km, and shows nulls of various depths all along the propagation path.
Figure 7. Predicted signal amplitude as a function of distance along the path,
for a fixed
ambient electron density. The solid line is for a frequency of 25.85 KHz, and
the dashed
line for 26.15 KHz.
Although the difference between the
predicted amplitude at 150 Hz below and
above the center frequency is only one dB or so, in many locations along this
particular
path, this difference is easily measurable over time scales of a few seconds, as
will be
shown in Chapter 4. In addition, the modal interference produces numerous places
along
the path where there is much greater sensitivity to frequency; for example, note
the shift
in the position ( km) of the null
near 550 km, and the large shift near the end of
the path at 2,100 km. Differences in phase as a function of frequency are
similar, with
deviations of several degrees over a 300 Hz frequency difference, which are also
easily
measurable.
Figure 8. Predicted signal amplitude as a function of distance along the path,
for a fixed
frequency (26.00 KHz), and two ambient electron density profiles. The solid line
is for
Profile II, and the dashed line for Profile I.
Although these differences are interesting, it is the sensitivity of these
frequency
differences to the selected electron density profile which is of interest here.
These ambient
density profiles, and the associated collision frequency profiles, are the
result of a sum of
many complex physical processes. Some of these processes create free electrons
in the
ionosphere (e.g., solar UV and X-ray energy, cosmic rays, and precipitation of
electrons
from the magnetosphere), while others remove them (e.g., chemical changes
affecting
attachment and recombination processes), and measurement of the actual ambient
profiles under a variety of geophysical conditions can be expected to give
information on
the importance of these processes and the parameters governing them. In
addition, the
ambient profiles themselves enter into the calculation of the importance of
other
ionospheric effects, such as the heating of the ionosphere over powerful VLF
transmitters.
Finally, knowledge of the ambient profiles is important for predicting those
aspects of
VLF propagation of interest in communication and navigation applications. As
shown in
Figure 8, for example, two different density profiles lead to quite different
predictions of
the field strength as a function of distance along the path, with frequency held
fixed at
26.00 KHz; differences of as much as 15 dB are seen at some locations, with at
least 5 dB
being common at most locations along the path. While these differences as a
function of
distance along the path provide another way to constrain the ambient profiles,
direct
utilization of this method requires simultaneous measurements at multiple
points. The
method used in this dissertation provides a means for the determination of the
same
constraints from a single-point measurement.
Figure 9. Predicted signal amplitude as a function of frequency, for the NSU-to-Stanford
path. Two profiles are included; the solid line is for Profile I, and the dashed
line for
Profile II.
When frequency is introduced as an independent variable, a plot of predicted
signal
amplitude vs. frequency, for a set of candidate profiles, may be produced. Such
a plot is
shown in Figure 9, where the predicted signal amplitude is shown over a ±150 Hz
range
around the nominal center frequency of 26.00 KHz, for two electron density
profiles. It is
clear from Figure 9 that an additional experimental test might now be applied to
the
ambient profile, since not only does the overall amplitude differ substantially
between the
two profiles, but there is an additional measurable difference in the dependence
of the
amplitude on frequency over the 300 Hz range.
A similar result may be obtained for the phase of the received signal, as shown
in
Figure 10. Here, phase differences of 15 -40
are predicted, but, as will be discussed in
Chapter 3, absolute phase delay and absolute group delay cannot be measured
without
information about the transmitted carrier phase and bit transition times. An
absolute
phase delay results in a simple phase offset, and an absolute group delay in a
linear trend
of phase vs. frequency, so these must be removed from the predictions to match
what is
available in the data. Since recordings were not made at the transmitter in this
study, the
phase data are replotted in Figure 11, with the linear trends removed.
It is clear from Figure 11 that, even without absolute timing, the phase vs.
frequency
plots are sufficiently different across the 300 Hz band of frequency to
discriminate
between the two profiles.
Figure 10. Predicted signal phase as a function of frequency, for the NSU-to-Stanford
path. Two profiles are included; the solid line is for Profile I, and the dashed
line for
Profile II.
The 3-D waveguide model makes predictions of the received E- and H-fields, but
normalizes the result to 1 KW of radiated power, independent of frequency, and
models
the transmitting and receiving antennas as short vertical electric monopoles,
using the
techniques of Ferguson and Snyder [1980] to determine the initial excitation
factors.
Referring back to Figure 5, it is apparent that the technique developed in this
work
measures the transfer function from the MSK modulator output to the recorded
signal
from the receiver, and therefore measures the product of three possibly
significant
transfer function components: the transmitter, the propagation channel, and the
receiver.
The contribution of the receiver transfer function is easy to estimate, since
the receiver is
accessible for direct measurement, and its characteristics do not change over
the time
scales of this experiment.
Figure 11. Predicted measured signal phase as a function of frequency, for the
NSU-to-
Stanford path. Two profiles are included; the solid line is for Profile I, and
the dashed line
for Profile II. This differs from Figure 10 due to suppression of constant- and
linear-phase
terms, which cannot be measured without knowledge of the transmitterís absolute
phase
and bit transition times.
The contribution from the transmitter requires more analysis to estimate, and,
in
the general case, would need to be measured. However, for the present study, we
can
make use of the regular properties of single-mode propagation over seawater to
check the
validity of modeling the transmitter as a vertical electric monopole over a
ground plane,
as is done in the 3-D waveguide model. To do this, we can use available data
(see
Chapter 5) from the NPM-to-Palmer path, which is a long path (12,349 km),
entirely
over seawater. The propagation along this path is well-understood, and fits the
predictions of a single-mode theory quite well [Inan and Carpenter, 1987].
Because there
is very little modal interference, and since the propagation constant does not
change
radically with frequency, the predicted amplitude and phase as a function of
frequency
are quite smooth, so that any variations which show up must be due to antenna
effects.
Figure 12. Comparison of measured propagation channel transfer function with the
prediction, for the NPM-to-Palmer path using electron density profile II. The
predictions
(shown as circles) have been normalized so the total predicted amplitude is in
agreement
with the total received power, and to remove any constant phase or group delay,
as
required by the measurement technique. The measured transfer function was
computed
from data taken on the magnetic east-west (E/W) antenna, and is displayed as a
solid line,
with the 95% confidence error range indicated by dashed lines.
The results of this comparison, in Figure 12, show that the predicted received
signal
amplitude and phase match that actually observed to within the measurement
error. The
nature of the measurement and its error will be discussed at length in Chapter
3, but it is
useful to note here that the measurement agrees with the prediction (for the
chosen
ambient electron density profile), to within ±0.2 dB and ±2 . The assumption
outlined
above, which treats both the transmitting and receiving antennas as short
vertical electric
monopoles, was used in formulating these predictions. The good agreement between
data
and model predictions indicates that this assumption is reasonable. For the
remainder of
this work, we will use this assumption, unless otherwise stated.
2.2. MEASUREMENT
OF ELECTRON DENSITY AND
TEMPERATURE PERTURBATIONS
Figure 13. Predicted amplitude of the propagation channel transfer function
perturbation,
as a function of frequency, on the NSU-to-Stanford path near 26.0 KHz. The
circles
represent the prediction for the perturbation located 500 km from the receiver,
and the
crosses the prediction for the perturbation being 1000 km distant from the
receiver.
Although the measurement of the ambient profile is of scientific interest, it is
made more
difficult by the important role played by the transmitter antenna excitation
factors, which
are difficult to estimate in the general case. However, changes to the ambient
electron
density and collision frequency profiles can be measured without this
uncertainty,
provided that the time scale of the changes to the density or collision
frequency profiles
are sufficiently different from any time scales associated with the transmitting
antenna,
and that the perturbations are not occurring too close to the transmitting or
receiving
antennas. In this case, the amplitudes and phases of the waveguide modes
launched by
the transmitter depend only on conditions which do not vary quickly with time.
The
receiving antenna, which converts incoming waveguide modes into a signal at the
receiver, acts similarly, and the contribution of the antenna excitation
factors, although
unknown, may be assumed to be fixed, and removed from the data record.
The predictions illustrated in Figure 13 were obtained by utilizing the 3-D
waveguide model for two different density perturbations, on the 26.0 KHz
Nebraska-to-
Stanford path considered earlier.
2.3. USE OF SPREAD-SPECTRUM SIGNALS
AS
MULTI-FREQUENCY PROBES
From the above discussion, it is apparent that the existing theory predicts
significant
differences in received amplitude and phase as a function of frequency, as
either the
ambient electron density profile, or certain characteristics of the perturbation
to that
profile (such as its location along the path), are varied. The performance of
existing
receivers, which treat the MSK signal as if it were a CW signal in many
respects, are
discussed in Chapter 4. However, the real advantage of the techniques presented
here lie
in the ability to make measurements across a band of frequencies, and the
comparison of
those measurements with the existing models.
There have been a few attempts to extend the existing single-frequency
technology
to make measurements at two or more frequencies. All of these efforts have used
a two-
frequency approach, which essentially uses two receivers, one tuned to each of
the two
main sidebands generated when an MSK signal is run through a frequency-doubling
operation to generate an FSK signal (see Chapter 4). In this manner, it is
possible to use
MSK signals to measure perturbations in group, as well as phase, delay, provided
that a
model consisting solely of phase and group delay is adequate for representing
the
propagation (see Chapter 4 for a discussion of this approach).
The first use of the two-frequency technique in the literature appears to be in
1981,
when it was applied to magnetospheric studies [Thompson, 1981]. In this study,
the
subionospheric VLF signal was assumed constant, and removed from the analysis,
in
order to focus on the (much weaker) long-delayed signals resulting from
propagation in
the ducted whistler mode. The MSK coding was demodulated only to examine the
residual whistler-mode signal, and frequency spread of the MSK signal was not
directly
addressed. Nevertheless, this technique allowed the accurate determination of
several
important properties of whistler ducts.
Dowden and Adams [1989, 1992] have applied the two-frequency technique to the
measurement of group delays on subionospheric signal perturbations. Although not
noted in the paper, their method only applies to cases where significant
inter-symbol
interference (ISI) is not present. In such cases, the two-frequency method does
indeed
give an accurate measurement of the phase- and group-delay perturbation. The
Dowden
and Adams papers do not indicate whether the amount of ISI was monitored,
although
this information is in principle available from the MSK discriminator. Hence, it
is not
clear whether the data presented are
truly a measurement of the best-fit phase and group
delays.
Finally, the Omega navigation system has been used to measure the perturbations
to the electron collision frequency profile, with low time resolution [Barr et.
al., 1985].
3.
PROCESSING ALGORITHMS
Since the primary purpose of VLF transmitters is to transmit data, the effects
of this
modulation on the received signal must be considered, and used in the analysis,
in order
to measure the transfer function of the VLF propagation channel. To do this, we
must
know the type of modulation employed, any parameters employed in the modulation,
such as bit rate and carrier phase, and the data being transmitted. In analyzing
the
characteristics (e.g., usable bandwidth) of various modulation schemes, some of
the
modulation parameters are often taken as unknowns, and treated as random
variables,
and the data transmitted are treated as samples from a random process. It is
important to
note, however, that in the case of any real measurement, the transmitted signal
constitutes a single realization or sample function of this random process.
By treating the transmitted data bits as samples from a random process with
known
statistics, we can calculate the power spectral density (PSD) of the modulated
signal. For
this calculation, the data bits are assumed to be equiprobable, independent, and
identically distributed, so that the resulting PSD is determined solely by the
characteristics of the modulation method. Of the many digital bandpass
modulation
methods, three are commonly used in the VLF range: continuous wave (CW),
frequency-
shift keying (FSK), and minimum-shift keying (MSK). CW, in the ideal case,
transmits
no data at all, and the CW signal consists simply of a single sine wave,
extending at
constant amplitude and phase from to . In practice, CW is a good
approximation to those transmissions which have very small bandwidths (i.e.,
constant-
frequency pulse with durations much longer than one period), such as those of
the
Omega navigation system. The PSD of a CW process is simply a pair of delta
functions at
, where
is the carrier frequency, as shown (for positive frequencies only) in the
left-
hand panels of Figure 14. If we take the squared magnitude of the transform of
any
sample of such a process, we obtain a so-called periodogram, which converges to
the PSD
as the sample length tends toward infinity, since no randomness is involved.
Figure 14. Common VLF modulation methods. The top panels show the power spectral
density (PSD) of the three common modulation methods, while the bottom panels
illustrate the finite-time transform of a particular realization of the
underlying random
process. The frequency axes are linearly scaled and zero-suppressed, and would
apply,
with appropriate scaling, to any combination of bit rate and carrier frequency.
These plots
cover ±500 Hz from the carrier, and were generated with a carrier frequency of
10 KHz,
and a frequency deviation of ±50 Hz (bit rate of 100 bps for FSK, and 200 bps
for MSK).
For any signal modulated in a non-trivial manner, however, the periodogram does
not converge to the PSD as the length of the data examined increases. This is
illustrated in
the center- and right-hand panels of Figure 14, where the difference between the
PSD and
the finite-time spectrum is apparent. If longer sequences are used in computing
the
periodogram, the frequency resolution becomes finer, but the variance does not
decrease.
The longer periodograms can be smoothed [Oppenheim & Schaefer, pp. 548-549],
to
yield a better estimate of the PSD, but of interest here is the fact that all
the "noise"
apparent in the lower panels of Figure 14 is a consequence of the particular bit
sequence
transmitted. If the modulation method, the bit sequence, and the other relevant
parameters such as carrier phase, are all known, then the transmitted waveform,
and
hence its spectrum, is completely known.
3.1. FSK SIGNALS
When the modulation employed is coherent binary FSK, the PSD may be derived
by writing the transmitted signal as an infinite sum of pulses at two
frequencies, at
:
(3.1)
where:
is the transmitted signal
is the transmitted amplitude,
is the carrier frequency,
is the carrier phase,
is the frequency offset from the carrier,
are
the message bits (±1),
is the bit rate, and
is
a one-bit pulse of unit amplitude.
An expansion of (3.1) into in-phase and quadrature (I and Q) components yields:
(3.2)
The I component is simply a coherent cosine function, which has zero mean, and
the Q
component is a series of half-sine pulses, whose signs depend on the message
bits . As
long as the are equiprobable
±1, the Q component will also have zero mean, and the I
and Q components will be independent, since the sine and cosine functions are
orthogonal over complete half-cycles. Under this assumption, the I component
produces
a pair of impulses in the lowpass-equivalent spectrum, offset by
, and the Q
component provides all the message-modulated component [Carlson, 1986, p. 520]:
(3.3)
This PSD is illustrated in the top-center panel in Figure 14, while the panel
below it
gives the transform of an actual bit sequence generated with this coherent
binary FSK
modulation. A receiver attempting to measure the amplitude and phase of such a
signal
can do so by measuring the amplitude and phase of either or both of the impulse
components, since these are always transmitted with coherent FSK, and do not
depend
on the message bits being sent. Since the phase may be measured at two
frequencies, the
group delay, as well as the phase delay, of the propagation channel may be
measured.
In practice, such a measurement suffers from two drawbacks. First, as a
practical
matter, it is impossible to suppress the information-carrying, modulated,
portion of the
FSK signal, without also affecting the measurement of the amplitude and phase of
the
carrier component, and any tradeoff which narrows the measurement bandwidth to
reduce the message modulation reduces the time resolution in direct proportion.
Second,
when the propagation channel affects the received signal in a way which differs
from
simple phase and group delays, this information cannot be extracted from the
received
signal, since the spread-out Q component is removed in the processing.
Therefore, even
in the case of FSK, where coherent signal components exist, a measurement which
identifies the bit sequence being transmitted and feeds back the result to
estimate the
propagation characteristics can be expected to perform better than simple
filtering of the
narrowband components.
3.2. MSK SIGNALS
For the frequency-dependent remote sensing needs as outlined in Chapter 2, the
VLF
stations transmitting MSK modulated data are nearly ideal. First, they provide
relatively
high signal levels, with several at or near 1 MW input power. Second, the strong
MSK
stations are geographically well situated, so that the propagation paths to
receivers such as
those at Palmer Station and Stanford cover a variety of terrain and expected
ionospheric
conditions. Finally, the MSK signal covers about 300 Hz of bandwidth with
significant
energy, and does so on a fairly continuous basis, provided only that the bit
sequence
being transmitted does not contain long strings of ones or zeros. Although this
bandwidth is small by normal spread-spectrum standards, it is much wider than
that
explicitly used in previous remote-sensing work, which generally treats all VLF
transmitters, no matter what the modulation, as if they were CW signals, by
attempting
to ignore the modulated components.
MSK modulation, however, differs significantly from the CW and FSK discussed
above, in that the MSK signal carries information on both the I and Q
components, and
there are no discrete spectral impulses present. The MSK carrier can be written
in a form
similar to Equation 3.1 above, except that the bit rate is twice that for the
FSK, and in
order to preserve phase continuity at the bit boundaries, the carrier phase must
depend
on the previous transmitted bits:
(3.4)
where:
is the carrier phase,
is the frequency offset from the carrier
Figure 15. Power Spectral Density of MSK-modulated
signal, illustrated here for a bit rate
of 200 bits/sec, which is a rate commonly used by the transmitters in this
study. The
spectrum covers at least 200 Hz with significant energy.
Since the carrier phase is now a function of the transmitted bits, the I and Q
components become:
(3.5)
and the I component no longer contributes fixed impulses to the spectrum, as was
the
case for FSK. Instead, both the I and Q components contain message information,
which
doubles the bit rate for the same bandwidth, thus making more effective use of
available
transmitter power and bandwidth. It can be shown [Carlson, 1986, pp. 522-525]
that, for
the case where the message bits are independent and identically distributed (I.I.D.),
the
lowpass equivalent spectrum of the MSK signal is:
(3.6)
which is illustrated in Figure 15.
Figure 16. Operation of an MSK modulator, implemented by phase-shift keying and
amplitude modulation of two coherent carriers. The BPSK modulators invert the
phase of
each carrier based on the number of message bits which have been sent on the
other
carrier. This preserves phase continuity at the transitions between bits, and
destroys the
coherence of the carriers.
An alternate view of the MSK modulation operation, which is useful for the
analysis
of existing phase receivers [Wolf, 1990], is illustrated in Figure 16, and may
be obtained
by considering the two carrier frequencies separately. In order to satisfy the
requirement
of phase continuity at the bit boundaries, the phase of the "mark" and
"space"
frequencies (at ) must be inverted
whenever an odd number of bits have been
sent on the other frequency. This can be done by passing each of the carriers
though a
binary phase-shift keying (BPSK) modulator, prior to on/off modulation which
selects
one carrier or the other for transmission. Because of the BPSK modulation, there
is no
coherent component in the output spectrum (for I.I.D. message bits), even though
the
amplitude modulation is always positive. This BPSK modulation spreads out the
carrier
energy over a frequency range on the order of the inverse of the bit time, while
ensuring
that the two (BPSK-modulated) carriers are at the same phase at the bit
boundaries.
Figure 17. MSK waveform, illustrated for
a specific bit sequence of ë10100011í, a bit rate of
200 bits per second, on a carrier of 200 Hz. For a real MSK wave, the carrier
would be
higher in frequency by two orders of magnitude, but this illustrates the
relationship
between the transmitted waveforms for each bit. Note, for example, that the
5-msec
waveform corresponding to the second ë1í is out of phase with the continuation
of that for
the first ë1í. These inversions, which depend on the message bits, destroy the
coherent
carrier waves present in the FSK case.
While the above analyses lead to an understanding of the power spectrum of the
MSK signal, the processing system described in this thesis works initially in
the time
domain to decode the bit sequence. Thus, an understanding of the MSK modulation
in
the time domain is needed. In MSK, at each bit interval, either the mark or
space
frequency is transmitted, but the signal may be inverted as needed to maintain
phase
continuity. This is illustrated in Figure 17 for a particular bit sequence,
transmitted on a
carrier which is at an artificially low frequency to illustrate better the
modulation. In each
bit time, the transmitted signal either gains or loses 90 of phase (ë1í or ë0í bit transmitted,
respectively), with respect to an unmodulated sinusoid at the carrier (center)
frequency.
Because of this, both phases of each frequency are needed to generate the MSK-
modulated signal, and, in general, the coherence of the carrier impulses present
in the
FSK spectrum is broken up.
One of the first steps in the signal processing described below is to mix the
incoming broadband signal with a complex exponential at the carrier frequency to
shift
the desired transmitter down to baseband, where a suilowpass filter can extract
the
transmitter signal of interest. The results of this operation are a time-varying
amplitude
and phase, relative to an unmodulated carrier exactly at the center frequency.
For the
MSK signals, the phase of this complex baseband signal is analyzed to demodulate
the
MSK signal, and thus determine the bit sequence transmitted. Since the MSK
signal either
gains or loses 90 of phase relative
to the center frequency in each bit time, a plot of
baseband signal phase vs. time can be described as a walk on a lattice, with the
steps at
each vertex (corresponding to the bit transitions) determined at the starting
edge of the
message bits. Such a plot is illustrated in Figure 18, and is used in the
working of the
demodulator algorithm described later in Section 3.3.3.
Figure 18. Plot of the phase of the MSK-modulated
signal, relative to an unmodulated
carrier at the same center frequency, for the same ë10100011í sequence used in
Figure 17.
The possible phase paths are shown as light lines, and all transitions occur on
the bit
boundaries, which here are 5 msec apart. The dark line is the phase path taken
by the
sample bit sequence. This phase gains 90
when a ë1í is transmitted, and loses 90
when a
ë0í is transmitted.
3.3. PROCESSING
OVERVIEW
Figure 19. Overview of the processing steps used in this work, for MSK signals.
The data of interest for the development and application of the spread-spectrum
remote sensing method put forth in this thesis concerns the measurements of
transmissions from MSK-modulated VLF stations, operating in the 20-30 KHz range.
The processing flow for these data is shown in Figure 19. The input data, which
consist of
16-bit samples of broadband VLF amplitude taken at a rate of 100 KHz, are pre-
processed to remove as much interfering energy as possible, and then
downconverted
directly to a complex baseband signal sampled at 1 KHz. The phase of the
baseband signal
is then processed to determine the transmitter parameters needed to simulate the
transmitted signal, so that the simulated signal can be used to compare with the
received
signal in the frequency domain. Differences between the simulated transmitted
and the
received signals reflect propagation effects, which is what we aim to measure.
The
following sections describe the processing steps in more detail.
3.3.1. Filtering
and Clipping of Broadband Data
Figure 20. Effect of initial processing steps on the broadband data. The top
panels display
the time series and spectrum of the raw data, sampled at 100 KHz, which shows
the strong
contamination from the Omega transmitter in Argentina, and a pair of sferics
near 2 msec.
The middle panels show the same data after the RF filter; the Omega interference
is gone,
but the sferic remains at high amplitude. The lower panels show the clipped
output, which
is used in all subsequent processing.
As described in Chapter 5, the raw data arrive in the form of 16-bit samples of
broadband
signal amplitude taken at 100 KHz, which are contaminated by various noise
sources. The
main noise source in this frequency band is impulsive noise due to lightning
discharges,
known as "sferics," but there are also other interfering signals in
the received data. These
are due to other transmitters, magnetospherically-generated signals such as
auroral hiss
(observed typically at high latitudes), and thermal noise. A representative
segment of the
raw data is shown in the top panel of Figure 20. For data recorded at Palmer
Station,
Antarctica, which constitutes the data set studied here (Chapter 5),
contamination from a
nearby VLF navigation transmitter (Omega Argentina) is particularly strong, with
a
signal intensity at least 10 dB above that of any of the MSK signals to be
measured. In
Figure 20, the 12.9 KHz Omega signal manifests itself as a dominant oscillation
with a
period of 77 µsec. Also seen in Figure 20 is a sferic at about 2 msec, which
has a long, low-
frequency tail. During the initial sferic burst, both the Omega signal and the
desired MSK
signals are overwhelmed, but they can be recovered in between the individual
bursts of
impulsive noise, and during the slow tail.
Figure 21. Amplitude and phase response of the high-pass filter used as the RF
filter. The
response is tailored to allow accurate measurement above 20 KHz, while rejecting
the
Omega stations and magnetospheric signals below 15 KHz.
To remove most of the magnetospheric signals, and some of the interfering
transmitters, the first operation on the raw data is to pass it through a band
selection
filter. This filter is illustrated in Figure 21, and, since this filter is
applied to the
broadband data before resampling, it is relatively short (101 points). The
results of
applying this filter to the previously-illustrated broadband data are shown in
the second
panel of Figure 20. The amplitude of the interfering Omega signal is
considerably
reduced, and the desired MSK signals are visible as a non-uniform signal
envelope in the
figure. This envelope is due to the fact that there are several MSK signals at
comparable
amplitudes, which interfere constructively and destructively over a time equal
to the
reciprocal of their frequency difference. It should be noted here that this
difference is
several KHz, and is unrelated to the bit time of 5 msec.
Once the data have been filtered, considerable energy may remain in the sferics,
which, because they are impulsive, have spectra that extend into the band of
interest. The
sferics produce two effects: they mask the desired signals, by causing the
receiver to
saturate, and they cause ringing in the narrowband stages to follow,
particularly in the
baseband filter. Since the larger sferics are strong enough to cause the
receiver to saturate,
there is no information available about the behavior of the desired MSK signals
while the
sferic is present (typically for 0.5-3 msec), and the only effect we can
mitigate is the
ringing in the baseband filter. To this end, the filtered output is clipped,
with the
threshold of clipping set on the basis of the distribution of the filtered
signal.
To set the clipping threshold, a histogram is first built up of the amplitudes
reached
during the entire data run. This distribution of signal levels in the data
consists of two
parts. The main component is due to the sum of many MSK transmitters in the
band,
and, as the sum of a large number of independent sinusoids, is approximately
gaussian.
The other component is a high-amplitude tail, caused by the rare, but
high-amplitude,
impulsive sferics. Therefore, a simple analysis of the amplitude histogram can
be made,
and an appropriate cutoff found, above which the large majority of the
excursions are
due to interference from the sferics. This procedure is illustrated in Figure
22, which
determines a clipping level of 0.30. The results of clipping the filtered data
to this level are
shown in the bottom panel of Figure 20.
In the course of this work, attempts were made to go beyond simply clipping out
the sferics, by performing the threshold procedure above, and then setting the
contaminated data to a "no data" value, rather than clipping it.
However, such algorithms
did not produce noticeable gains, relative to the dynamic clipping discussed
above. In
addition, suppression of data brings about substantial additional complexity in
the
subsequent processing, since the "no data" values must actively be
handled at each
processing step. Therefore, all the data presented in this work have been
processed by the
filtering and clipping methods discussed in this section.
Figure 22. Setting of the clipping threshold. The algorithm fits the amplitude
histogram
data (open circles), to a gaussian, and uses the fit to adjust the clip levels.
The bottom
panel, on an expanded scale, shows the "tails" of the distribution,
which are largely due to
sferic energy.
3.3.2. Demodulation
of Broadband Data
The first step in processing the filtered and clipped data is to bring it down
to baseband,
where the subsequent processing takes place. For this work, the baseband
bandwidth
needed is about ±300 Hz, and the baseband is chosen to be sampled at 1 KHz.
Since the
processing is digital, frequency translation is implemented by mixing the
broadband
signal with a complex exponential, passing the result through a lowpass filter,
and
decimating the (complex) result down to the baseband, at zero frequency, and a
sample
rate of 1 KHz. The lowpass filter used in this work is shown in Figure 23, and
is flat
enough in both amplitude and phase so that the errors it introduces are not
measurable
compared to the transmitter and receiver unknowns. This filter is a very long
FIR filter
(3001 points), but, since it is a decimation filter, it only needs to be run
once for every
output point. All the data processed through this filter uses a decimation ratio
of 100,
which results in a substantial computation savings.
Figure 23. Amplitude and phase of baseband filter used in the processing. The
filter
illustrated is suifor processing MSK signals, where the first spectral nulls are
at ±150 Hz
from the center frequency. For this frequency range, the filter is flat to
within 0.001 dB,
and 0.001 . The contribution from a similar interfering transmitter 600 Hz away
from the
center frequency is down by 100 dB.
The baseband signal is a complex (I and Q) signal. This is convenient, since all
the
processing done after this step can be done at the baseband, and since the
bandwidth of
the VLF stations is only a few hundred hertz, the sample rate needed to capture
this
information is low. The use of a complex signal at DC, rather than a real signal
centered
about a suifrequency, allows the various filters in the system to be (possibly
complex-
valued) lowpass filters. For the case where the filter is symmetric about zero
frequency,
which is true of most of the filters used in this work, the filter coefficients
are real, and the
processing load is reduced by nearly a factor of two.
In addition, the amplitude of the baseband signal resulting from this process is
simply the integrated amplitude of the MSK signal at the receiver, as would be
read by a
conventional heterodyne receiver, and the phase of the baseband contains the
message
modulation, in the form illustrated in Figure 18. The amplitude and phase of the
baseband signal for a particular example are shown in Figure 24, where the MSK
modulation can clearly be seen in the lower panel. The positive and negative
phase ramps
from the modulation differ from the theoretical value of 90
per bit time, due to the
additive noise (marked by arrows) and filtering. In addition to small
deviations, there are
three significant sferics during this 250-msec interval: one at 60 msec, one at
200 msec,
and one at 240 msec. The effect of the baseband filter on the data is seen in
the rounding
of the transitions between positive- and negative-going phase ramps (i.e.,
between ë1í and
ë0í bits). This filter suppresses the more distant sidebands of the MSK
signal, and thus
rounds the transition edges, but, since the filter is zero-phase, this rounding
is symmetric,
and the location of the bit transition times can be accurately determined from
the filtered
phase.
Figure 24. Amplitude and phase output from baseband stage. The top panel shows
the
amplitude of the output of the
baseband filter, and the bottom panel the phase. The MSK
modulation is clearly visible in the phase.
3.3.3. Detection
of Transmitted Bit Sequence
In order to estimate the transfer function of the propagation channel, the input
signal
must be known. For an MSK transmitter, five quantities must be known in advance,
or
must be reliably estimated from the data: ,
the modulatorís center frequency, ,
the
sequence of transmitted bits, , the
frequency deviation (or, alternatively, the baud
rate), , the carrier phase, and
, the time of the bit transitions. For the VLF signals
used here, the center frequency is precisely known, and controlled to a high
degree of
precision [ITU, 1979] by a frequency standard. Similarly, the baud rate is also
known,
from an inspection of the received data, since in practice only a few baud rates
are used,
with 200 bits/sec being by far the most common. Hence, only three parameters
must be
estimated from the data: , , and .
The time of the bit transitions, ,
is estimated by considering the entire data
record, and looking for the times at which the phase changes direction.
Referring to
Figure 25, those times ideally occur at bit boundaries. The program estimates
the time of
the transition, modulo the bit time, to a small fraction of the bit time, by
comparing the
phases before and after each stationary point on the phase plot. By this
procedure, it is
possible to get to better
than 1% of a bit time, while only observing the phase at a few
points in each bit time. For each transition from ë0í to ë1í bits, or vice
versa, the phase plot
has an inflection point. There are spurious inflection points due to noise and
interference, but these are randomly distributed in time, while those which
define the bit
time boundaries all occur nearly at multiples of the bit time. The estimation
algorithm for
first builds a histogram of the
inflection point time estimates (modulo the bit time),
and then determines where the best-fit peak of this histogram falls, which
results in a
single value of which is then
used in the subsequent stages for the entire data run. The
data runs used in this work are short enough (5 minutes maximum), so that this
procedure works well, though a dynamic estimate of
could be implemented to allow
the processing of large data sets.
Once the bit transition times are known, there are two possible hypotheses for
the
data received in any bit time: either a ë1í was transmitted, or a ë0í was
transmitted. By
examining the data points for a bit interval, the goodness-of-fit with each
possibility is
computed using a quadratic cost function. A provisional identification of the
transmitted
bit is made based on best fit at this stage.
Figure 25. Reconstruction of the transmitted phase, and hence demodulation of
the MSK
signal, from the observed phase of the decimator output. The solid line in the
above graph
is the received phase, as observed at the output of the baseband. There is
noise, and the
effects of two spherics can be seen in the data, as well as the rounded corners
at the bit
transition points, which is caused by the finite bandwidth of the baseband
filter. The
dashed line is the reconstruction of the transmitter phase from this data.
While the above scheme would be sufficient in the absence of noise, there is a
good
chance that some of the bit decisions will be incorrect. Referring to the
observed data in
Figure 25, the sferic which arrives at about 60 msec disrupts the phase
completely for
nearly 5 msec, and could easily lead to a bit error. A single bit error will
result in the
reconstructed phase advancing by 90 when
it should have retarded, or vice versa, and,
after the disturbance has vanished, the reconstructed phase will differ from the
transmitted phase by a full 180 . This is a particularly egregious form of error
for our
purposes, because when the transfer function is estimated with this data, the
presence of a
180 flip in the reconstruction
causes the magnitude of the estimate to drop to zero,
reverse phase, and rise again. This appears as sharp drops, all the way to zero,
in the
amplitude, and coincident 180 shifts
in the phase, and makes the resulting data very hard
to interpret.
In order to prevent such occurrences, the MSK decoder does not accept the
provisional bit decisions without looking ahead; instead, it uses the
goodness-of-fit values
to correct decisions which are likely to be in error. After a provisional bit
decision is
made on two bits, if the phase during the second bit is close to 180
off, and the
goodness-of-fit is better for the second bit than the first, the decoder assumes
that the
first bit is in error, and reverses its earlier decision. Even though this
scheme is not perfect
at correcting bit errors, it forces the decoder to always make the errors in
pairs, and the
reconstructed phase never gets too far from the received phase. Since the
expected
perturbations dealt with in this work are always small, over the time spans
comparable to
the bit rates, this correction mechanism is not likely to be activated by a
propagation or
LEP perturbation. The reconstructed phase, as determined by the algorithm
described
above, is shown as the dashed line in Figure 25.
One additional possibility should be mentioned, and that is that use could be
made
of any redundancy present in the transmitted message to correct bit errors and
determine
the input signal to the propagation medium to a high degree of confidence. For
many
messages transmitted over such a channel, some form of forward error-correcting
code is
most likely applied, and knowledge of the basic code structure would provide an
ideal
way to determine the transmitted signal, extending the usefulness of these
algorithms to
lower signal-to-noise ratio regimes. If the system has been partitioned into a
source
encoder, which might contain a cryptographic component, and a channel coder,
which
serves only to correct channel errors, only knowledge of the latter would be
needed for
measurements of propagation characteristics. The algorithms investigated in this
thesis
need only obtain a good estimate of the transmitted signal; they do not need to
interpret
the message sent in any other way.
3.3.4. Reconstruction
of the Transmitted Signal
In order to estimate the transfer function of the propagation path, which, as
discussed in
Section 1.2, gives information about the electron density profiles and
perturbations to
these profiles, the signal parameters measured above ( ,
, and ) must be fed back
and used to reconstruct the signal launched from the transmitter into the earth-
ionosphere waveguide. In addition, the output signal from the waveguide must be
measured. Both of these processes involve antennas and transmitting or receiving
electronics, and the contribution of these to the system model (see Figure 5)
must be
accounted for.
On the transmitting side, the signal parameters measured above can directly be
used to reconstruct the modulated VLF signal. Referring to Figure 19, this
function is
performed by the box labeled MSK modulator. The inputs to this stage are
, , and
,
and the output consists of an MSK-modulated signal at the system baseband,
which, as
discussed above, is a complex signal centered at zero frequency. If the
estimated
parameters are correct, the output will be an accurate representation of the
signal which
was presented to the power amplifier of the VLF transmitter. The issue of the
relation
between the transmitter input signal and the EM wave launched into the earth-
ionosphere waveguide was discussed in Chapter 2, and the excitation coefficients
for the
various modes are included in the 3-D waveguide propagation model, so the MSK-
modulated signal constructed here should be suifor comparison with the model
results, at
least over the time scales of interest in this work.
Representation of the receiver is somewhat simpler. First, the receiver is
available
for analysis, and, as the characterization data presented in Chapter 5 show, the
receiving
hardware can be well-represented by a combination of a group delay, a phase
delay and a
frequency offset. The data are processed as 32-bit floating-point numbers, which
have
sufficient dynamic range for the contribution of the processing software to be
largely due
to amplitude and phase distortions in the baseband filter. The first two of
these, the phase
and group delays, are absorbed into the unknown transmitter phase
and bit transition
time . The frequency offset is
small enough to be adjusted for in the data display, which
leaves only the response of the baseband filter to be accounted for. This is
done in
Figure 19 by the step labeled receiver filter, and is important only when
bandwidths wider
than about ±300 Hz are examined, or a narrower than normal baseband filter is
used. For
these cases, the reconstructed MSK signal is resampled at the higher frequency
used in the
original data sampling (100 KHz), the baseband filter applied, and the result
decimated
back down to the baseband sample rate. Thus, if the estimates for the signal
parameters
are correct, and the transmitter is well-modeled by the 3-D waveguide
propagation
model, the differences between the signal at the output of this stage and the
received
signal are solely due to noise and propagation effects.
3.4. COMPARISON
WITH RECEIVED SIGNALS
To perform the comparison between the received signal and the reconstruction of
the
transmitted signal, in a form which allows direct comparison with the model
results, both
signals are converted into the frequency domain, and compared as spectra. This
is
illustrated in Figure 26, for a segment of 256 samples. The transform of both
the
transmitted and received signals contains a large amount of message variation,
but, if the
transfer function from transmitter to receiver is nearly constant, most of these
variations
will be common to both signals, and will drop out of the comparison. This
removes the
effects of the message, leaving only the transfer function estimate, as a
function of
frequency, plus noise. The noise level is roughly constant across the band of
frequencies
processed here, and so the accuracy of the transfer function estimate drops off
at
frequencies where the input signal power is less.
Figure 26. Comparison of output to input fourier transforms to yield transfer
function
estimate. The left-hand panels show the transform of the received signal, over a
time
interval of 256 samples (about 50 bits), and the middle panels show the same for
the
reconstructed transmitted signal. The left-hand panels show the ratio of the two
transforms, which is the raw transfer function estimate.
While Figure 26 illustrates the general procedure employed to estimate the
transfer
function, there are several refinements which must be introduced in order to
obtain
useful estimates in the presence of noise. First, the propagation is modeled as
a simple
linear system, with a known input and output, an unknown transfer function, and
a
stochastic noise term which is independent of the input:
(3.7)
where:
is the unknown channel impulse response,
is
the input signal (inferred from data),
is the output signal (measured), and
is the noise.
If the additive noise is
uncorrelated with the input , and
has zero mean,
then the cross-spectrum of the output and input signals will be given by:
(3.8)
where:
is the unknown transfer function,
is the spectrum of the input signal,
is the cross-spectrum of output and input signals.
This representation suggests forming a simple transfer function estimate as the
ratio of the finite-length transform of the output signal to that of the input
signal. This
"raw" transfer function estimate is performed simply by taking the
ratio of the Fourier
transforms of the input and output functions, neglecting the additive noise
term, and is
referred to as the empirical transfer function estimate [Ljung, 1987, p. 146]:
(3.9)
where:
is the FFT of the output signal,
is the FFT of the input signal, and
is the raw estimate at frequency w.
is an asymptotically unbiased
estimate of the transfer function ,
provided that the additive noise has zero mean, and that the input,
, is not
periodic over the length of the sample window N. This estimate is very general,
in that it
assumes only linearity of the propagation channel, but it is also very noisy,
because it
calculates an independent estimate at each frequency point
, and it is possible to show
[Ljung, 1987, pp. 146-151] that the variance of the empirical transfer function
estimate
does not decrease with N. For this problem, as with many practical systems, we
can
obtain a better estimate by smoothing the raw estimates, based on the assumption
that,
for this propagation channel, the transfer function is not independently
determined at
each frequency. If we introduce the assumption that the true transfer function,
,
varies slowly compared to the frequency difference between each of the N output
points,
then each of these output points is an independent estimate of the same transfer
function
point, and has a variance that depends on the observed output signal spectrum at
that
point:
(3.10)
where:
is the spectrum of the noise, and
is the observed transform of the input.
Therefore, a better estimate of is:
(3.11)
where:
is a weighting function centered about
.
Now, we assume that the noise spectrum
is also smooth, and does not change
much over the interval where the weighting function is non-zero. Then, this term
in the
integrals is a constant, and cancels out, since it appears in both the numerator
and
denominator. Thus, we can simplify (3.10) to:
(3.12)
which is the form of the estimator used in this work. In our implementation, a
rectangular weighting function was generally used, and both the width of the
weighting
function and the length of the transform used to compute
and can be chosen
depending on the presentation desired.
3.5. SMOOTHING AND DISPLAYING RESULTS
Figure 27. Output of the propagation channel transfer function estimate, for a
sample data
set. The top panel gives the amplitude, and the bottom the phase, of the
transfer function,
over a ±300 Hz frequency range. Note that the falloff in the input signal power
past about
120 Hz causes the transfer function estimate to be very noisy at these
frequencies.
The results from this computation, for a typical data set, are shown in Figure
27.
For this Figure, the estimator described above was used, with the FFT size
set to 256,
with 100 points between FFT segments, and with a rectangular weighting function
covering 10 Hz. Thus, there are 10 output slices per second, with 61 separate
frequencies
displayed, covering ±300 Hz from the center frequency. The computed transfer
function
is displayed as a spectrogram, in two panels, the top panel being the amplitude
(in dB) of
the transfer function, and the bottom panel being the phase. This format is used
from
here on throughout this dissertation. The phase data have been processed so as
to make
best use of the spectrogramís limited gray-scale range, by setting the zero
point of the
phase and eliminating any drift between the local oscillator and the average of
the
received data. The darker portions of the display indicate higher signal
amplitudes and
phase advances . Visible in Figure 27 are several LEP-associated perturbations,
as well as
variations in the background propagation channel characteristics. The rapid
falloff in the
MSK signal spectrum beyond about 120 Hz causes the transfer function estimate at
these
frequencies to become dominated by the noise, and sets the limit on the
frequency range
over which the transfer function of the propagation channel may be determined by
this
technique.
Figure 28. Smoothed sample data, based
on the same data set as Figure 27. Note the
changed z-axis scales.
The data presented may be further smoothed, depending on the intended use
of the
display, with some further reduction in the apparent noise, at the expense of
time and
frequency resolution. For example, the data in Figure 28 are averaged over 2
seconds, and
30 Hz. The sideband structure of the MSK-modulated signal is easily visible,
since the
noise which adds to the input becomes more significant as the available signal
power goes
down. However, in this example, as with most of the data, only the central ±150
Hz band
will probably be useful for most applications.
Figure 29. Estimate (solid lines) of the average transfer function of the
propagation
channel for the 4 minute data set, over ±300 Hz in frequency, with the 95%
confidence
limits (dashed lines). The falloff in the power of the input MSK signal beyond
±120 Hz is
visible here as an increase in the confidence limits.
These spectrograms contain all the information available about the propagation
channel transfer function, for a given choice of FFT and smoothing parameters;
however,
for some purposes, further processing is useful to highlight certain features.
For example,
the time-averaged transfer function estimate may be derived by averaging the FFT
"slices," and compared directly against results of the 3-D waveguide
propagation model,
as shown in Figure 29. The high signal-to-noise ratio near the center frequency
results in
very good estimates of the transfer function in this range, though beyond about
±150 Hz
of the center frequency, the input signal is below the noise, and the transfer
function
estimate is biased and noisy. For this reason, the data presentations in the
rest of this
dissertation are generally limited to the useful range of ±150 Hz for the MSK-modulated
signals.
When the background propagation channel transfer function is removed from the
data, the residual highlights transient changes to the transfer function, such
as those
caused by the LEP events discussed in Chapter 2. Figure 30 displays this
residual, as the
deviation from the average transfer function (Figure 29). In this figure, the
full 4 minutes
and ±300 Hz range of the data are displayed.
For detailed examination of an individual perturbation to the propagation
channel
transfer function, the time and frequency scales may be expanded, as shown in
Figure 31.
These transfer function perturbations are typical of the LEP-induced events
analyzed in
Chapter 5. They show a sharp onset (~1 second risetime), and a slow decay (~20
seconds
fall time), and are uniform across the ±150 Hz frequency band.
Figure 30. Spectrogram of the propagation channel transfer function residual,
after
removing the average transfer function estimate from the data. This display
highlights the
transient disturbances to wave propagation, of which at least six events
characteristic of
LEP are visible. All these perturbations are amplitude decreases and phase
advances, and
do not show significant variation with frequency.
Figure 31. An expansion of two overlapping transfer function perturbations from
Figure
30, with an expanded frequency scale showing only frequencies within ±150 Hz of
the
carrier, where the input signal is strong enough to give good estimates of the
transfer
function. The first perturbation has a sharp onset at about 16 seconds, and is
followed by a
larger event at about 20 seconds. Both perturbations are amplitude decreases and
phase
advances which are smooth in frequency.
4.
SIMULATION OF PROCESSING SYSTEMS
Since many physical processes influence the propagation of VLF signals, and the
net
effect of these processes is difficult to predict, it is not always possible or
cost-effective to
assess the performance of a new measurement technique using only data obtained
in the
field. While the theoretical basis for this technique has been described above,
it remains
to be shown that it is capable of resolving geophysically significant
propagation
differences in the presence of ambient interfering signals and atmospheric
noise. To
demonstrate this, we turn to simulation, and generate a VLF signal at a
frequency not
utilized by any of the existing transmitters, alter the signal by passing it
through filters
which simulate the characteristics of a propagation path of interest, superpose
the
modified signal on to wideband VLF data recorded in the field, and analyze the
result
using the techniques described in Chapter 3. The objective of this simulation is
to
determine whether the characteristics of the propagation path (as artificially
simulated by
a filter) can be assessed via measurements of signal amplitude and phase as a
function of
frequency. The degree to which the results of the analysis agree with the
filters used in
simulating the propagation channel reflects the utility of the measurement
technique on
real, noisy, data.
Figure 32. Processing flow for the simulation. By simulating a chosen
propagation path
with the 3-D waveguide model, and using the results to generate a channel
transfer
function filter, a signal can be generated which can be processed as described
in Chapter 3
to yield propagation amplitude and phase information. These results can be
compared
with the channel transfer function filter to test the analysis method.
4.1. PROCESSING
OF SIMULATED MSK SIGNAL
Figure 32 diagrams the simulation process. An assumed set of values of the
ionospheric parameters chosen to match those used in recent work [Inan et. al.,
1992],
combined with an arbitrary transmitter and receiver location and carrier
frequency, are
used in the 3-D waveguide propagation model [Poulsen et. al., 1993a] to
calculate the
propagation channel characteristics at several closely-spaced frequencies. The
model
results are then used to design a time-varying filter which simulates both the
phase and
amplitude changes due to the static propagation path, and a time-varying
perturbation,
representing a localized ionospheric disturbance somewhere along the path. The
filter
was designed using the window method [Oppenheim & Schaefer, pp. 239-250],
modified
to accept an arbitrary amplitude and phase in the desired transfer function. A
4-sample
Blackman-Harris window [Harris, 1978] was used. This filter may contain multiple
stages: a first stage which is time-invariant and which accounts for the effects
of the
ambient earth-ionosphere waveguide; and a second stage, which is slowly
time-varying
and which simulates the effects of an LEP-induced ionospheric density change on
the
propagation of the VLF wave. The output of the filter is then mixed in with
recorded
data, with the gain adjusted to give a signal amplitude similar to those
observed on VLF
transmitter signals in an actual field-recorded data set. This data is then
analyzed using
the measurement methods described in Chapter 3, and the resulting time-varying
propagation channel spectrum is compared with the characteristics of the
propagation
channel filter to determine the degree to which the measurement method can
determine
the characteristics of the filter, over time scales of geophysical interest, and
in the
presence of a realistic noise environment.
4.1.1. Construction
of the Test MSK Signal
As discussed in Chapter 3, most of the high-power VLF transmitters currently in
operation employ MSK modulation, since the bandwidth of these systems is
limited, and
the higher data rates obtainable with MSK, as opposed to FSK, modulation justify
the
slightly greater complexity in the modulator and receiver. All of the
transmitters studied
in this work utilize MSK, which is the modulation chosen for our simulated
signal.
There are six variables which fully determine the MSK-modulated signal (See
Chapter 3), at the input to the transmitter: the amplitude ( ), the carrier
frequency ( ),
the carrier phase ( ), the bit rate ( ), the time offset ( ), and the bit
sequence being
transmitted ( ). Of these, is
normalized to unity, since it is simply a multiplier which
determines the signal level (the relative amplitudes of the antenna excitation
factors for
each of the waveguide modes are computed by the 3-D waveguide propagation
model),
and and
are known. The simulation code must assign values to the other three
variables. For the performance tests described in this Chapter,
and are both set to
zero, since the processing described in Chapter 3 cannot estimate these without
absolute
time and phase information, which is not present in the field data. The code
which
generates the simulated MSK waveform can set these variables to arbitrary
values,
however, and this capability was used in developing and testing the analysis
code. Note
that in real data acquired without absolute phase and time information, these
quantities
would similarly be unknown, and would be suppressed in the data analysis, as
discussed
in Chapter 3.
Figure 33. Autocorrelation of a sample received bit sequence. The value at zero
lag is
normalized to unity, and is off the top of the graph. No simple repeating
structures are
seen, out to a lag of 5 seconds (1000 bits).
With the above assumptions, the MSK-modulated waveform is determined solely
by the bit sequence being transmitted, .
Although the actual data being transmitted has
structure, since it presumably is the output of some channel-coding process at
the
transmitter, an examination of the second-order statistics of actual received
data (in this
case, the NPM signal observed at Palmer Station) failed to reveal any
significant
departure from independence between the bits (a sample autocorrelation function
is
shown in Figure 33). Therefore, in the simulation, the data bits were simulated
as
independent and identically distributed (IID) random variables. For such a case,
it is a
simple matter to construct a MSK waveform which is representative of the signal
at the
input to the power amplifier stages of a real VLF transmitter. Once this test
waveform has
been constructed, the bit sequence used in the modulation, as well as the
carrier phase
and time origin, are not used in further processing, as they are generally not
available to
commonly used receivers in the field.
Figure 34. Location of simulated VLF transmitter ("NSU"), and the
simulated ionospheric
perturbation on the great-circle propagation path to a receiver at Stanford. The
disturbance to the ionospheric density profile is assumed gaussian in shape,
with s=75
km.
4.1.2. Simulation
of the Slowly-Varying Propagation System
The novel feature of the propagation channel measurable by the technique
described here
is the variation of the amplitude and phase over the ~300 Hz frequency band
illuminated
by the signal. In order to test the method, therefore, we need to find a set of
conditions
which are not only of geophysical interest, but which can be expected to lead to
variations
in the signal amplitude and phase over the 300 Hz bandwidth occupied by a
typical VLF
MSK-modulated signal, as a result of propagation through the channel. As
discussed in
Chapter 2, these variations are expected to be most significant when the
received signal
consists of many waveguide modes with similar amplitudes, which occurs for the
case of
short, land-based propagation paths.
Figure 35. Electron density profiles used in the simulation. The ambient profile
is
"Profile II" from Figure 6, while the perturbation is related to a
whistler described in the
text.
The path chosen for our simulation is illustrated in Figure 34, from a
hypothetical
transmitter at the location of the commonly-observed 48.5 KHz Air Force
transmitter in
Silver Creek, Nebraska (41 30¥ N,
97 36¥ W), but operating at 26.0
KHz, to a receiver
located at Stanford (37 25¥ N, 122
10¥ W). The 3-D waveguide propagation model was
utilized, for frequencies from 25.80 KHz to 26.20 KHz, in 50 Hz steps, and for
an ambient
nighttime D-region electron density profile in effect along the entire path,
except under
the disturbance, where the altitude profile of electron density was taken to be
one
representative of conditions following an LEP event. The ambient profile used
corresponds to Profile II of Figure 6, while the disturbed profile is that which
would be
expected for a lightning-induced electron precipitation burst of peak flux
, with a duration of 0.2 seconds, as would be produced by a
0.2-
6.0 KHz whistler, propagating at L=2.5 [Inan et. al., 1988b]. These two density
profiles
are illustrated in Figure 35.
Figure 36. Mode structure of the received signal at Stanford, from the
fictitious 26.0 KHz
transmitter in Nebraska, as computed by the 3-D waveguide propagation model.
Many
modes contribute significantly to the total field at the receiver, so small
perturbations to
individual modes may affect the total field strength disproportionately.
When the 3-D waveguide propagation model is used for the ambient conditions,
there are at least six propagating modes which contribute significantly to the
received
signal (i.e., within 6 dB of the total). These modes are shown in a polar
display in
Figure 36. From this figure, it is apparent that the total predicted field at
the receiver is
less than the amplitude of the single largest mode, due to partial cancellation
by other,
slightly weaker modes. Poulsen [1993b] has shown that the strength of the field
scattered
by a typical LEP disturbance depends strongly on the waveguide mode; in
particular, the
quasi-transverse magnetic (QTM) modes are scattered more effectively than quasi-
transverse electric (QTE) modes, and the higher-order modes are scattered more
effectively than the lower-order modes. When the total received signal consists
of
multiple significant modes, the change in the total field at the receiver is
very sensitive to
the phase changes to the different modes. Since the scattering may also be a
function of
frequency, small differences in scattered fields at different frequencies making
up a single
modulated signal may combine to yield a perturbation which is a strong function
of
frequency.
Figure 37. Amplitude and phase of transfer function for propagation at 26 KHz
from
"NSU" to Stanford, under the ambient conditions. The predicted
amplitude varies by
more than 3 dB over the 300 Hz band of frequency; the predicted phase varies by
nearly
15 over the same interval.
The extent of these differential amplitude and phase perturbations may be
predicted by using the 3-D waveguide propagation model for the above conditions,
for a
closely-spaced set of frequencies covering the range over which the simulated
MSK signal
has significant energy. For this simulation, model calculations were carried out
for nine
frequencies spaced 50 Hz apart, ranging from 25.8 KHz to 26.2 KHz. The predicted
variation of the received signal amplitude and phase, for the chosen ambient
electron
density, are shown in Figure 37. This result illustrates that, for this short,
multi-mode
path, the transfer function of the propagation channel exhibits significant
variations over
the frequency band of interest.
When the localized perturbation to the electron density is introduced, the
amplitude and phase of the signal at the receiver changes, with the difference
due to the
perturbation also being a function of frequency. As discussed in Chapter 3,
these
differential perturbations are more easily measured, since they are not
sensitive to
assumptions about the transmitter and antenna system. The detectability of
differential
changes in signal amplitude and phase, due to typical LEP-induced electron
density
changes, is of most interest in this Chapter. These propagation effects can be
simulated by
using a slowly time-varying filter, which normally has an impulse response
appropriate to
the ambient propagation prediction, but changes to an impulse response
appropriate to
the perturbed condition, over a time scale of several seconds, which is 10-100
times the
length of the impulse response of the filter. These specific propagation channel
transfer
functions, and filters which represent the channel transfer function
predictions, are
shown in Figures 38 and 39.
Figure 38. Amplitude of propagation from NSU to SU as a function of frequency,
as
predicted by the 3-D waveguide propagation model, for both ambient and perturbed
conditions (open circles in the top panel). A 4001-point FIR filter is fit to
each of these
predicted propagation characteristics (solid lines in the top panel). The
difference between
the two characteristics and filters is shown in the bottom panel, and is the
predicted
amplitude of the LEP-induced perturbation to the electron density, as a function
of
frequency. The 3-D waveguide propagation model was run only out to ±200 Hz; the
FIR
filter was designed to have constant amplitude beyond these points, since signal
components beyond ±200 Hz are heavily attenuated by the baseband filter during
processing.
Figure 39. Phase of propagation from NSU to SU as a function of frequency, as
predicted
by the 3-D waveguide propagation model, with the simulation filters (top panel).
The
difference in the bottom panel is the predicted phase of the LEP-induced phase
perturbation. As with the amplitude response, the phase of the FIR filter was
selected to be
constant beyond the ±200 Hz range covered by the 3-D waveguide propagation
model
runs.
To simulate the effects of an LEP-induced ionospheric disturbance, the generated
MSK signal is passed through a time-varying filter, as indicated in Figure 32.
This filter
varies from the "ambient" filter to the "disturbed" filter
in one filter length (40 msec),
and then relaxes to the "ambient" filter with a first-order time
constant of 10 seconds.
Such a time variation approximates the typical onset and decay times seen in LEP-
induced ionospheric disturbances generated by mid-latitude whistlers [Inan and
Carpenter, 1987].
Finally, the output from this filter, which simulates the propagation over the
short,
land-based path indicated in Figure 34, is combined with VLF data recorded in
the field.
The amplitude of the synthetic signal is set to a value at the low end of the
range observed
for real VLF transmitters in the Antarctic, which results in a spectrum similar
to that
shown in Figure 40. This composite signal is then used as the input to the
analysis
procedure, without reference to the filters used to generate the perturbation.
Figure 40. Spectrum of the broadband VLF, with the synthetic NSU signal added at
26.0 KHz. The NSU signal is mixed in at a low enough amplitude to simulate a
relatively
weak VLF transmitter.
4.1.3. Use of
Spread-Spectrum Processing
Using the simulated signal as input, the processing steps detailed in Chapter 3
were
performed, with the result displayed as a spectrogram of the transfer function
in
Figure 41. The two perturbations, at 60 and 180 seconds, are clearly visible,
and show
some variation over frequency. The amplitude perturbation is greatest above the
center
frequency (peaks at about +100 Hz); the phase perturbation is positive below the
center,
and negative above it. Both observations agree with the characteristics of the
filter used to
simulate the propagation, indicating that frequency-dependent amplitude and
phase
variations imposed during the propagation of the signal can be recovered upon
analysis.
Figure 41. Spectrogram resulting from processing the simulated NSU signal. Two
simulated perturbations were introduced, at 60 and 180 seconds from the start of
the
record. The nature of the perturbation is clearly visible in the amplitude
record, and, to a
lesser extent, in the phase as well. In addition, the asymmetry of the
perturbation over
frequency is visible, particularly in the phase.
4.2. PERFORMANCE
OF EXISTING AMPLITUDE AND
PHASE RECEIVERS
Most of the receivers used in previous studies in this field measure the
integrated
amplitude of the VLF signal. When phase receivers are used, they track the phase
of a
single coherent component of the modulated VLF signal using phase-locked loop (PLL)
techniques [Wolf, 1990]. The operation of these receivers, and the limitations
they impose
on the data analysis, are described in this section.
4.2.1. Amplitude
Receivers
Amplitude receivers which measure the integrated signal amplitude are easy to
build, and also easy to analyze, since they are, ideally, linear receivers
except for the
envelope detector. The forms of modulation commonly used in the VLF band are
constant-amplitude, and so the resulting measurement, whether on MSK, FSK, or CW
signals, is a true measurement of the overall attenuation of the signal as a
result of the
propagation. Since the signal amplitude is a function of frequency, as are the
propagation
effects and the receiver filters, the detector in such a receiver ideally
measures the
cumulative result of three components: the transmitted signal, the propagation,
and the
receiver filter shape. Specifically, if
if the spectrum of the transmitted signal (for
one particular bit sequence), is the propagation effect to be measured, and
is
the cumulative transfer function of all the filters used in the receiver, then
the
measurement output from an ideal amplitude detector is:
(4.1)
With careful design, may be
made nearly flat over the bandwidth where
is
significant. If we also assume that the effect of the propagation is also flat
in amplitude
and linear in phase over the occupied bandwidth, then:
(4.2)
so the signal into the detector is simply a time-delayed replica of the input,
multiplied by
a complex amplitude factor. The detector will ideally strip all the phase
information, and
thus will simply measure the propagation as:
(4.3)
since the detector bandwidth is much smaller than the carrier frequency of
, and
therefore is equal to unity.
Although the above analysis is strictly valid only when
the receiver and propagation effects are linear and time-independent, it can be
applied to
the common case in the VLF bands where variations on the propagation are much
slower
than the frequencies contained in the transmitted signal. Finally, since the
receiver filter
characteristics are (ideally) stable, they can be measured, and the only
time-varying effect
at the output is the propagation amplitude
.
If the amplitude and phase of the propagation varies significantly over the
bandwidth occupied by the transmitted signal, then an integrated-amplitude
receiver will
measure the integral of the propagation constant, weighted by the spectrum of
the
transmitted signal:
(4.4)
where B is the bandwidth occupied by .
Thus, this technique cannot resolve differences in propagation as a function of
frequency, and only produces a weighted integral of the propagation amplitude.
While
this is sufficient for many useful measurements, it introduces a practical
problem as well:
if the receiver filters are not perfectly flat, the deviations from flatness
will multiply the
short-time transmitted signal spectrum. The output will exhibit variations which
depend
on the short-time spectrum of the transmitted signal, and thus on the message
bits being
sent. This message modulation gets worse as the integration time of the detector
output is
reduced, because the integral of the short-time spectrum of the product varies
more as
the time interval is reduced.
4.2.2. Phase
Receivers
For CW and FSK-modulated signals, an ideal phase receiver can be described which
is
even simpler than the integrated-amplitude receiver, since CW and FSK-modulated
signals both contain unmodulated signal components. A phase receiver for CW and
FSK
signals can simply measure the phase of one or more of these unmodulated
components
(see Figure 14). In a practical FSK receiver, the finite bandwidth of the
phase-locked loop
(PLL) which tracks the phase of one of the unmodulated carrier components would
admit some of the modulated signal as well, and therefore the output would
exhibit some
amount of message modulation. Further, there would be a tradeoff between
measurement time and the variations due to message modulation; a shorter time
constant in the PLL accepting more of the modulated signal, and resulting in
greater
message modulation in the output.
Figure 42. Typical block diagram of a frequency-doubling MSK phase receiver,
with
amplitude output, after Wolf [1990].
A related disadvantage of this type of phase receiver is that, even with FSK
signals,
measurements may be made at only two frequencies, which correspond to the
unmodulated carrier components. Even though there is significant energy at other
frequencies, it cannot be used to measure phase at those frequencies. The main
advantage
of this type of receiver is that it is relatively simple to build, and needs no
information
about the message being sent, other than the spacing of the mark/space
frequencies from
the center frequency.
For MSK signals, which have no unmodulated components, the coherent-sideband
type of receiver can still be used, if the MSK signal is first put through a
frequency-
doubler circuit. This serves to construct an artificial FSK signal from the MSK
input. In
the MSK input, both carrier frequencies are spread out because the requirement
for phase
continuity forces each carrier to be present in two phases, 180
apart. Where FSK may be
constructed from two coherent carriers, by switching between them at
zero-crossings, the
MSK, with double the bit rate, can only maintain phase continuity by biphase-
modulating the carriers appropriately, and this modulation spreads the carriers
out over
a bandwidth comparable to the signaling rate. Figure 42 is a block diagram of
such a
receiver, currently used by Stanford at a number of field sites for continuous
measurements [Wolf, 1990].
Figure 43. Waveform of MSK signal and an FSK signal synthesized from it by a
square-law
frequency doubler. The center frequency of the MSK signal has been lowered to
400 Hz to
exaggerate the effects, and offset in amplitude for clarity. The lower trace is
the original
MSK signal, and the upper the frequency-doubled result. The message bits being
sent here
are ë01011í.
Figure 43 illustrates the frequency-doubling operation in the time domain. The
original MSK signal contains no coherent carrier components, as can be seen by
examination of the phases of the carrier for the ë0í bits. The first and
second ë0í of the
transmitted message are separated by an odd number of ë1í bits, and the so the
carrier is
switched in phase by 180 between
the two bits. However, in the FSK signal synthesized
from the MSK, the frequency-doubling operation has mapped both 0
and 180 phases of
the carrier into the same, 0 , phase. With reference to Figure 16, this
operation negates
the effect of the BPSK modulators on the "mark" and "space"
frequency carriers, and the
resulting modulation is FSK, not MSK. This restores the spectral peaks at the
mark and
space frequencies, which can be seen in Figure 44.
Figure 44. Spectra of MSK signal and a frequency-doubled FSK signal generated
from it.
The carrier frequency of the MSK signal is set to 1 KHz, and the spectrum
computed for a
random sequence of data bits (2000 total bits).
Receivers which use these techniques exhibit the same limitation described above
for FSK signals: as the tracking loop bandwidth is made wider to obtain higher
time
resolution, more of the modulated message energy is passed to the output, and
the
amount of message modulation increases. In addition, there is a more fundamental
problem involved with the frequency doubling: since it is a non-linear
operation, it
performs frequency conversion on the components of the transmitted MSK signal.
The
carrier components which appear in the FSK signal at the output of the frequency
doubler were not transmitted, and, in the general case where the propagation is
a
sensitive function of frequency, the FSK carriers may be synthesized by the
doubling
operation from MSK components which were attenuated and phase-shifted by
differing
amounts in propagation. Thus, measurement of the amplitude and phase of the
synthesized carriers does not, in general, give measurements of the amplitude
and phase
of the channel transfer function at the carrier frequencies. This deficiency may
be
reconciled with the time-domain analysis by noting that large variations in the
amplitude
or phase of the channel transfer function over the frequency range occupied by
the MSK
signal imply that the impulse response of that transfer function is significant
over a time
comparable to the duration of one transmitted bit. Such a circumstance would be
a case
of large amounts of Inter-Symbol Interference (ISI), and the waveform during one
bit
interval will depend on those before it. As discussed in Section 4.1.2, cases
where this is
expected to occur include short, land-based propagation paths, where many
waveguide
modes are present at comparable amplitudes.
5.
APPLICATION TO DATA
While the simulation results described above in Chapter 4 demonstrate that the
frequency-dependent analysis techniques used in the present work are suifor
accurate
measurements of the transfer function of the VLF radio propagation channel, it
is the
application of these techniques to field-recorded data which is of scientific
interest. In
this chapter, we describe a data set digitally recorded for this purpose at
Palmer Station,
Antarctica, and give examples of the application of the new analysis techniques
to these
data.
5.1. SIGNAL PATHS STUDIED
For the purposes of this study, four signal paths were studied, all terminating
at Palmer
Station, Antarctica (64 W, 65 S). 0 shows the parameters for the signals
observed. All of
the signals measured arrive at Palmer over long paths, two of them (NPM, NLK)
almost
entirely over seawater, from transmitters in the northern hemisphere. All the
listed signals
are modulated with 200 bits/sec MSK, and are strong enough at Palmer to allow
accurate
decoding of the bit sequence. The characteristics of these signals as observed
at Palmer
have been studied extensively in previous work [Carpenter and LaBelle, 1982;
Tolstoy et.
al., 1986; Inan and Carpenter, 1986; Wolf and Inan, 1990].
Figure 45. Great-circle propagation paths for the signals studied in this
section. All
originate in the northern hemisphere, and terminate at Palmer Station,
Antarctica.
Table 1. Transmitter and Arrival Characteristics
Call Sign
Frequency
(KHz)
Power
(KW)
Latitude
(degrees)
Longitude
(degrees)
Mag Bearing
(degrees)
Path Length
(km)
NPM
23.4
502
21 25¥ N
158 09¥ W
275 48¥
12,349
NSS
21.4
213
38 59¥ N
76 27¥ W
350 00¥
11,606
NAA
24.0
1,200
44 39¥ N
67 17¥ W
357 30¥
12,187
NLK
24.8
245
48 12¥ N
121 55¥ W
318 24¥
13,526
The amplitudes and relative phases of the significant waveguide modes that
constitute the received signals at Palmer Station, as predicted by the 3-D
waveguide
propagation model, are shown in Figure 46. For all four signals, the received
mode
structure is dominated by a single mode, although the dominant mode is different
for the
different paths. According to Poulsen [1991, p. 46], under typical conditions,
the lowest
attenuation rates are obtained for the low-order QTE modes (QTE1 AND QTE2), and
for
the lowest-order QTM mode (QTM1); however, the QTM2 mode is most efficiently
excited by the vertical monopole transmitting antennas, and can remain dominant
out to
quite large distances, in spite of its somewhat higher attenuation rate.
Figure 46. Mode structure of the arriving VLF signals at Palmer Station, as
predicted by
the 3-D waveguide propagation model, for the center frequency of each
transmitter. In all
four cases, a single mode is dominant (> 10 dB above the other mode
strengths).
With a single dominant mode, the propagation characteristics of the VLF wave
should be a relatively smooth function of frequency, since the ionospheric and
ground
conductivities vary smoothly with frequency so that rapid changes with frequency
are
most likely caused by the interference between modes of nearly equal amplitude.
Figures 47-50 show the amplitude and phase at the receiver as a function of
frequency in
a ±200 Hz range around the center frequency, as predicted by the 3-D waveguide
propagation model for the indicated paths. The ambient electron density profile
used was
Profile II from Figure 6 for all these calculations. In comparison with Figure
37, which
shows the predicted propagation channel transfer function for a short,
land-based, path
with several significant modes, all the transfer functions in Figures 47-50 vary
smoothly
with frequency, with maximum amplitude variations across the ±200 Hz bands
shown of
about 1 dB, and phase variations of about 0.3 . In addition, there are no sharp
features in
the spectrum; sharp features would only be expected in cases where several
significant
waveguide modes coexisted, in which case destructive interference can occur at
some
frequencies.
Figure 47. Predicted propagation from NPM in Hawaii to Palmer Station. The
predicted
phase has been adjusted to remove constant phase- and group-delay terms, as
described in
Chapter 3, since these cannot be measured from the field data without absolute
timing
information at both the transmitter and receiver.
Figure 48. Predicted propagation from NSS in Maryland to Palmer Station. The
predicted
phase has been adjusted to remove constant phase- and group-delay terms, as
described in
Chapter 3, since these cannot be measured from the field data without absolute
timing
information at both the transmitter and receiver.
Figure 49. Predicted propagation from
NAA in Maine to Palmer Station. The predicted
phase has been adjusted to remove constant phase- and group-delay terms, as
described in
Chapter 3, since these cannot be measured from the field data without absolute
timing
information at both the transmitter and receiver.
Figure 50. Predicted propagation from NLK in Washington to Palmer Station. The
predicted phase has been adjusted to remove constant phase- and group-delay
terms, as
described in Chapter 3, since these cannot be measured from the field data
without
absolute timing information at both the transmitter and receiver.
Although the above results illustrate signal characteristics under ambient
(unperturbed) ionospheric conditions, the perturbation to the signal due to the
appearance of localized disturbances along its propagation path (i.e., due to
lightning-
induced electron precipitation) is also expected to be a smooth function of
frequency,
when a single, low-order, mode is dominant. For a single waveguide mode, the
phase
perturbation ( ) induced by a change in ionospheric reflecting height ( ) is
given by
[Inan and Carpenter, 1987]:
(5.1)
where:
n is the order of
the mode,
d is the length of
the path under the perturbation,
l is the
wavelength,
is the radius of the Earth,
f is the wave
frequency, and
c is the speed of
light.
For low-order modes in the 20-30 KHz range, and a typical reflecting height of
90 km, the first term in the brackets dominates the
term, and the phase perturbation
is simply proportional to frequency, representing essentially a group delay. In
such a
model, the amplitude perturbation is simply due to increased absorption at the
lower
altitudes, and is varies insignificantly with frequency over the range occupied
by these
signals.
5.2. RECEIVING,
RECORDING AND PROCESSING OF DATA
The receiving equipment for this experiment consisted of a set of IGY-style
crossed-loop
antennas, which are single loops of #6 cable, supported in an isosceles triangle
9 meters
tall and 18 meters across at the base [Paschal, 1978], a broadband VLF preamp
and line
receiver, and a digital data recording system, as is illustrated in Figure 51.
All parts of this
system, except for the digital recorder, have been in use since 1988 at Palmer
station for a
variety of ionospheric and magnetospheric studies [Wolf, 1990]. In addition to
the digital
wideband recordings made for this experiment, the existing narrowband data
system was
used to monitor activity on the signal paths, in order to determine when to
record
broadband data. In this way, many examples of lightning-associated VLF signal
perturbations were obtained, despite the limited recording capacity of the
digital
wideband recorder (5 minutes per event, and about 5 hours total).
Figure 51. Diagram of the receiving and
recording equipment used to acquire the data
analyzed in this chapter.
To reduce interference from the generators and other equipment at Palmer
Station,
the loop antennas are situated on a glacier in back of the station, about 1 km
from the
main buildings. The signals received by the loops were amplified by the preamps,
and
driven differentially on lines back to the recording hut. There, a line
receiver unit supplies
gain, filtering, and buffering, and allows the injection of calibration signals.
The
conditioned analog signal is then passed through an anti-alias filter, digitized
at a rate of
100 KHz, and stored on optical disk.
5.2.1. Antenna Efficiency
Since the antennas used are electrically small loops, compared to the
wavelengths of
signals in this frequency range, they are sensitive to the B-field of the
incoming waves,
and are most sensitive to plane waves arriving in the plane of the antenna loop;
thus, the
N/S loop is most sensitive to waves arriving from the North or South, and the
E/W loop
to waves from the East or West. The antennas are aligned to within ±5
of geomagnetic
N/S and E/W.
Most of the data used in this work were recorded from the N/S loop antenna,
since
it gives the largest signals from the North American transmitters NSS, NAA, and
NLK.
Although the N/S antenna orientation nearly nulls the reception of NPM, this is
partly
compensated by the fact that the NPM signal to Palmer propagates almost entirely
over
water, with relatively little attenuation. The great-circle propagation paths
for the other
frequencies are over land for significant parts of their length, and the signals
arriving at
Palmer are attenuated. The net result of these effects is that the amplitude of
the listed
signals in the wideband spectrum observed by the N/S antenna are close to equal,
as can
be seen in Figure 40.
The ultimate limit to the performance of a VLF receiver using a loop antenna is
set
by the resistance of the wire used for the antenna, and its relationship to the
area enclosed
by the loop [Paschal, 1988]. The wire resistance contributes thermal noise to
the receiver
input, and, even with a noiseless receiver, limits the weakest detecsignals.
This noise
voltage is given by:
(5.2)
where:
is the noise voltage density, in ,
is Boltzmannís constant,
is the temperature of the wire, and
is the resistance of the wire.
Since an incoming plane wave, arriving from a direction in the plane of the
loop, induces
a voltage on the antenna of:
(5.3)
where:
is the number of turns of wire in the antenna,
is
the area of the antenna bundle,
is the electric field strength of the wave, and
is the speed of light,
the noise introduced by the resistance of the antenna wire is equivalent to an
arriving
wave of magnitude:
(5.4)
where:
is the antenna sensitivity, in .
This sensitivity is a spot noise density [Motchenbacher and Connelly, 1993, p.
44],
and is a function of frequency. Over the bandwidths of typical VLF signals,
however, the
sensitivity only varies by ~1 percent, so that we can use the carrier frequency
in the above
expression, to calculate the minimum detecsignal from the spot noise density and
the
post-detection bandwidth, B, occupied by the signal:
(5.5)
where is the incoming signal
strength, in , which would give 0
dB signal-to-
noise ratio at the output of a noiseless receiver connected to the antenna in
question. For
the type of loop antennas used at Palmer, Paschal [1978], estimates the antenna
noise
component at about 0.016 µV/m, over a bandwidth of 500 Hz, at a center
frequency of
25 KHz. This is a theoretical limit to the detecsignal with the Palmer antenna,
and does
not account for the noise in the receiver, nor for atmospheric sources of noise.
However,
the sensitivity of the IGY-type antenna is typically about 20 dB better than the
atmospheric noise levels normally present in the 20-30 KHz frequency range
[Fraser-
Smith and Helliwell, 1985].
A loop antenna responds to the rate of change of the magnetic field passing
through the enclosed loop area, so that the voltage induced on the loop
increases linearly
with frequency. However, the loop impedance contains an inductive component,
and,
above the frequency where the loopís inductive reactance equals the input
impedance of
the loop and receiver input stage, the increasing antenna reactance offsets the
increasing
induced voltage, and the current flowing into the receiver input becomes
independent of
frequency. The input turnover frequency is given by [Paschal, 1988]:
(5.6)
where:
is the input resistance of the preamplifier,
is
the turns ratio of the input transformer, and
is the inductance of the loop.
The frequency response of the system is limited on the low end by
, which for the
Palmer system is about 300 Hz .
At the higher-frequency end, the magnetic loop antenna performance begins to
deviate from that described above when the size of the antenna is a significant
fraction of
a wavelength of the incoming wave, so that, at any instant of time, the magnetic
field
strength varies across the loop. However, since this occurs at frequencies
outside the band
of interest by several orders of magnitude, the high-frequency response is
limited only by
circuit and system considerations discussed below.
5.2.2. Preamp
Frequency Response and Noise Figure
Since the antenna is located at some distance from the recording equipment, its
output is
fed into a preamplifier, and then driven on a balanced differential analog
transmission
line, down to the recording hut about 1000 meters from the antenna site. The
combination of the antenna and preamp determine most of the characteristics of
the
recording system, including the low-frequency response, flatness of gain, phase
delay, and
noise level. The only major characteristic controlled by the recording system
itself is the
high-frequency response, which is sharply rolled off by the anti-alias filter.
Figure 52. Amplitude response of the
Palmer antenna and preamp. The system is flat to
within 1 dB from 1-100 KHz, as shown in the top panel. The bottom panel is an
expansion of the response over the range of interest in this work, showing that
the
response is flat to within better than 0.01 dB over the typical 500 Hz
bandwidth occupied
by an MSK-modulated VLF transmitter.
The frequency response of the Palmer antenna and preamp are shown in Figures 52
and 53, and from those figures it is apparent that the gain and phase errors are
insignificant over the 20-30 KHz frequency range occupied by the MSK-modulated
VLF
transmitters studied in this work.
The system noise is dominated by the noise in the front end of the receiver,
which
was designed to provide a noise figure of about 4 dB at the lower frequencies,
where the
antenna sensitivity is lower. At higher frequencies, the noise factor increases
with the
square of the frequency, which offsets the linearly increasing antenna
sensitivity to
produce a broad minimum in system sensitivity. The signals studied in this work
are in
the 20-30 KHz range, for which the combination of antenna and preamp has a noise
density equivalent to an incoming plane wave of about 0.02
. Over the
entire 30 KHz bandwidth of the recording system, this translates into a
broadband system
noise of about 3.5 µV/M, RMS, which is much lower than the noise from
atmospheric
sources [Fraser-Smith and Helliwell, 1985].
Figure 53. Phase response of the Palmer
antenna and preamp. In the bottom panel,
covering the 10-30 KHz frequency range, a constant phase and group delay have
been
removed, to show the phase error. Over any 500 Hz range, the phase error is much
less
than 0.1 .
After transmission to the recording hut, the balanced signal from the line is
converted to a single-ended signal, filtered, and amplified to useful levels by
the receiver.
The receiver consists of two components: the line receiver, which performs most
of the
processing, and a time code mixer, which mixes in a 10 KHz modulated time- and
phase-
reference tone, without altering the signal level. As all this processing is
performed with
signals which have already been amplified well above the noise floor of the
processing
equipment, the effect of noise in the later amplification stages on the gain,
phase, and
noise are minimal. Since the calibration data were taken with this equipment,
this noise is
included in the above figures, but is insignificant in comparison with the
effects of the
preamp and antenna.
5.2.3. A/D
Converter and Storage System
The wideband VLF data are recorded in fully digital form, using commercial
equipment
and custom software for the filtering, conversion and storage of the broadband
data. The
output of the receiver first passes through a "brickwall" low-pass
filter, with a cutoff of
over 100 dB per octave beyond 33 KHz, to remove signal components which would
cause
aliasing. The data are then sampled at 100 KHz, using a 16-bit A/D converter
encapsulated in a metal enclosure to reduce the electromagnetic interference
from the
host PC, so that true 16-bit performance can be achieved.
Figure 54. Sensitivity of the Palmer antenna and preamp, expressed as the
electric field of
an incoming wave which gives an output equal to the system noise density. The
noise
figure of the preamp is the difference between the "System" and
"Antenna" curves, and
approaches 4 dB at frequencies below a few KHz. Above 10 KHz, the preamp noise
figure
increases with the square of the frequency, and the system becomes less
sensitive.
Converted data are stored first on the hard disk, in a circular buffer capable
of
retaining up to five minutes of data. The sampling system can run continuously,
saving
the most recent data, which allows the operator to trigger saving of the data
when an
interesting event is observed. While the system is recording, there is no
analysis possible
of the incoming data. However, other instruments available at Palmer give good
indications of activity, including on-line displays of VLF signal perturbations
on the
paths studied in this work. These available displays were used to control the
capture of
data.
Once five minutes of data are stored on the hard disk, they must be copied to
the
optical (WORM) disk before another recording may be made, due to storage
limitations
on the hard disk. This process takes nearly an hour, after which the system is
again ready
to record. For post-processing, the data from the WORM disk are loaded onto the
hard
disk of the analysis computer, using the algorithms described in Chapter 3.
5.2.4. Noise
Calibration Measurements and Dynamic Range
The recording system was checked for noise performance by replacing the antenna
with a
dummy loop, having the same impedance, but no collecting area, effectively
removing the
incoming VLF signal from the output. What is left at the output is largely
preamp noise,
the time mark, and its harmonics. By injecting a low-level signal into the
calibration
input, we can measure the minimum detecsignal in the presence of the preamp
noise.
A spectrum recorded during such a test is shown in Figure 55, which shows that a
test tone at 5.5 KHz, 95 dB below the clipping level, is easily detectable. Also
visible are
the harmonics of the time mark, and the shape of the anti-alias filter cutoff at
33 KHz.
From Figure 55, it can be seen that the preamp noise is about 10 dB above the
noise floor
set by the A/D converter, so that the recording system does not limit the
detection of low-
level signals. The preamp has a very large dynamic range, about 120 dB in a 1 Hz
bandwidth [Paschal, 1988], larger than the recording range of the A/D converter.
Therefore, the recording level is set so that the sum of all the coherent VLF
transmitters
in the 0-33 KHz range is about 10 dB below clipping, so that large-amplitude
signals,
such as intense spherics, are clipped, and most of the dynamic range is
available for
recording the signals of interest.
Figure 55. Data recording system noise and signal detectability. The lower line
shows the
noise floor of the A/D and processing system only, while the upper line shows a
recording
through the entire system, with the time mark at 10 KHz, and a calibration
signal injected
at 5.5 KHz, 95 dB below clipping. The calibration signal is easily detectable,
even at this
low level. Harmonics of the time mark can be seen, as can the rolloff of the
anti-alias filter
above 33 KHz. The vertical scale is calibrated so that a sine wave input which
just reaches
the clipping levels is 0 dB.
5.3. PROCESSING
SYSTEM
The analysis of the data recorded at Palmer was conducted at Stanford using an
IBM PC-
compatible computer equipped with a board containing a WE DSP-32C chip and
associated memory. Since the data were already recorded digitally, this board
contained
no signal conversion circuitry, but accelerated the processing of the data by
several orders
of magnitude. The software to implement the processing algorithms was written in
a
combination of C and DSP-32C assembler, with the C code on the host PC
performing
data I/O and high-level tasks, and the DSP-32C code optimized for running the
large
filters and transforms involved. With minimal attempts at optimization, the
software
runs in nearly real time on this hardware, with most of the processing work
involved with
the downconversion from broadband data to the complex baseband signal.
The results of the processing were displayed with custom software, either on-
screen, or printed on charts for comparison with existing data, or converted
into Matlab
for further processing or image display.
5.4. MSK SIGNALS ON EXPECTED
SINGLE-MODE PATHS
As discussed in Section 5.1,. the signals recorded and analyzed in this work are
all
predominately single-mode, and are expected to yield propagation channel
transfer
function estimates which are smooth functions of frequency. In addition, the
single-
mode theory predicts that, when LEP-associated perturbations to this transfer
function
are seen, the perturbation should also be constant across the band of
frequencies
illuminated by the signal. We now compare representative data from the field
with the
predictions, both for the time-average propagation channel transfer function,
and for the
LEP-induced perturbations to this transfer function.
5.4.1. Time-Averaged Propagation Channel
Transfer Functions
Following the procedure described in Chapter 3, the time-averaged propagation
channel
transfer functions were computed for selected data runs; a typical set of
results are
presented in Figures 56-59. When these measured transfer functions are compared
directly with the predictions in Figures 47-50, the following differences are
seen:
1) The signal received
from NPM on the N/S antenna is higher than the prediction
by about 7 dB, and has more variation with frequency.
2) The NSS measurement
has an amplitude slope of about +1 dB over 200 Hz,
while the predicted slope is about -1 dB over the same range.
3) The NAA measurement
shows no amplitude slope, and is also about 5 dB lower
than the prediction.
4) The NLK measurement is
about 5 dB lower than the prediction.
Figure 56. Average propagation channel transfer function and 95% confidence
limits for
the NPM-to-Palmer path, for a 5-minute data set collected on the N/S loop
antenna. This
antenna has a null nearly in the direction of NPM; the amplitudes have been
corrected for
the expected antenna response.
Figure 57. Average propagation channel transfer function and 95% confidence
limits for
the NSS-to-Palmer path, for a 5-minute data set collected on the N/S loop
antenna.
Figure 58. Average propagation channel transfer function and 95% confidence
limits for
the NAA-to-Palmer path, for a 5-minute data set collected on the N/S loop
antenna.
Figure 59. Average propagation channel transfer function and 95% confidence
limits for
the NLK-to-Palmer path, for a 5-minute data set collected on the N/S loop
antenna.
There are many uncertainties which affect the predictions, and, to a lesser
extent,
the measurements:
1) The electron density
and collision frequency altitude profiles were assumed
constant ("Profile II" from Figure 6) everywhere along the propagation
path; in
reality, the profiles are not known, and can be expected to vary along the path.
For the predominantly single-mode signals studied here, the main effect of
increasing the electron density would be to lower the reflection height, and the
increased collision frequency at the lower altitudes would be expected to
increase the absorption, and reduce the overall amplitude.
2) Although
the NSS, NAA, and NLK paths are well within the nighttime
hemisphere at the time these data were taken, NPM was near the terminator
(2043 local time). Since the ionospheric electron density is very different
between dark and illuminated portions of the ionosphere, the different
propagation characteristics might cause refraction of some of the signal energy,
resulting in an enhancement at the receiving site.
3) The radiated power from the transmitting antenna is
known only approximately
[Rodriguez, 1994]. This should only affect the comparison of the overall signal
levels.
4) The transmitting and
receiving antennas are assumed to be electrically small
dipoles by the 3-D waveguide propagation model, when computing modal
excitation coefficients. This mirrors the actual construction of typical
transmitting antennas [Watt, 1967, p. 253], but the receiving antennas are
actually small loops (magnetic dipoles). Thus, the comparison of the
measurement to the prediction depends on both the ratio of the E- and B-field
components of the incoming wave, and the orientation of the receiving antenna.
Of these factors, (2) is an important subset of (1), since the effect of the
day/night
terminator on VLF radio propagation is due to the very different electron
density profiles
between daytime and nighttime. The third uncertainty is of no consequence except
when
comparing absolute amplitudes, and may be measurable during the daytime, when
propagation conditions are more stable. Therefore, if we can account for the
effects of the
last item in the above list, the processing described in this work becomes a
tool for the
determination of parameters of geophysical interest: the electron density and
collision
frequency profiles along the signal propagation path.
Figure 60. Effect of mode structure of incoming wave on received signal, for a
small loop
antenna. The small loop is sensitive to the time-varying component of the
magnetic field
which is perpendicular to the plane of the antenna. In the top panel, the wave
has no
magnetic field component in the direction of propagation, and the familiar
cosine antenna
pattern results. In the middle, a mode which has a significant magnetic field in
the
direction of propagation will produce an additional output. In the bottom panel,
when
such a wave is received by an antenna which is nulled (perpendicular to the
direction of
arrival of the signal), the output is due solely to the non-transverse component
of the
magnetic field.
Figure 60 illustrates the importance of knowing the modal structure of the
incoming signal, particularly when that signal is being received on an antenna
which is
nearly nulled to the direction of arrival. In the ideal case, where the incoming
wave is a
pure plane wave, or is pure transverse magnetic (TM), the relation of the
component
to the component is simply
the wave impedance, and there is no component.
In
such a case, the measurements may be compared to the predictions of the 3-D
waveguide
propagation model simply by allowing for the response of the antenna.
If the incoming wave is not TM, then we must consider the effect of the
component on the recorded signal. The non-zero
component rotates the total
magnetic field vector toward the axis of propagation, and it is the angle
between the total
H vector and the antenna normal which must be used when calculating the antenna
response. If the mode in question has non-zero
, then a significant error may be
induced by using the simple response,
when the antenna is nearly nulled, to infer
the total H or total E for comparison with the 3-D waveguide propagation model
results.
In addition, the ratio of to varies with the waveguide mode number, and so this
effect may alter the predicted mode structure as seen at the receiver.
In the data recorded for this work, the only signal which is nearly nulled on
the N/S
loop antenna is NPM. To evaluate the importance of this effect, data recorded on
the
E/W antenna were also analyzed, since for this antenna, the direction of arrival
of the
signal from NPM is nearly in the antenna plane, rather than along a null. The
time-
averaged propagation channel transfer function from these data is shown in
Figure 61.
Unfortunately, the recording equipment used at Palmer does not have the capacity
to
record from both the N/S and E/W antennas simultaneously, so the data in Figure
61 is
from a different day. Nevertheless, the variation with frequency is smoother
(compared
to Figure 56), and the amplitude agrees well with the prediction of the model.
Figure 61. Average transfer function of the propagation channel from NPM to
Palmer, as
observed on the E/W antenna, which has a maximum response nearly in the
direction of
the great circle path to NPM.
5.4.2. LEP-Induced
Perturbations to Transfer Functions
In the data recorded at Palmer, there are a number of perturbations to the
time-averaged
propagation channel transfer function which display the sudden onset and slow
recovery
typical of the LEP-induced perturbations described in Chapter 2. The single-mode
theory
discussed in Section 5.1 predicts that these disturbances should result in
changes to the
propagation channel which are only a slowly-varying function of frequency.
Data showing 10 perturbations on three signal paths (some simultaneous on two of
the paths) are presented in Figures 62-67. Within the resolution of these
displays, the
perturbations all are uniform across frequency, as would be expected for a
single
dominant waveguide mode.
Figure 62. Two transfer function perturbations on the NPM-to-Palmer path, one at
~180
seconds and one at ~290 seconds. The data have been smoothed, and the static
portion of
the transfer function suppressed.
Figure 63. An expansion of the time axis in Figure 62, showing the single LEP
event with
an onset at about 182 seconds from the start of the record (07:18:02 UT). The
perturbation is an amplitude decrease and a phase advance, and is smooth in
frequency.
Figure 64. Four transfer function perturbations on the NSS-to-Palmer path, at
~10, ~50,
~75, and ~110 seconds. The data have been smoothed, and the static portion of
the
transfer function suppressed.
Figure 65. An expansion of the time axis in Figure 64, showing two LEP events
with onsets
at ~50 and ~77 seconds from the start of the record (09:00:50 and 09:01:17 UT).
These
perturbation are amplitude increases and a phase advances, and are smooth in
frequency.
Figure 66. Four transfer function perturbations on the NAA-to-Palmer path,
coincident
with those on NSS, from Figure 64. The data have been smoothed, and the static
portion
of the transfer function suppressed.
Figure 67. An expansion of the time axis in Figure 66, showing two LEP events
with onsets
at ~50 and ~77 seconds from the start of the record (09:00:50 and 09:01:17 UT).
These
perturbation are amplitude increases and a phase advances, and are smooth in
frequency,
and are coincident with the perturbations on NSS.
6.
SUMMARY AND SUGGESTIONS FOR
FUTURE WORK
6.1. SUMMARY
OF RESULTS AND CONTRIBUTIONS
In this dissertation, we have developed a novel technique for the measurement of
the
propagation characteristics of VLF radio waves in the earth-ionosphere waveguide,
which
takes advantage of modulated transmissions to characterize the propagation
channel over
a band of frequencies. The technique for measurement of the transfer function of
the
propagation channel does not rely on assumptions about the characteristics of
that
transfer function, as have all previous techniques reported in the literature.
This introduction of frequency as an independent variable in the propagation
measurements allows direct comparison with existing propagation theory, and thus
constrains the ionospheric parameters (electron density and temperature, as a
function of
altitude) which control the propagation of the VLF waves. Knowledge of these
parameters
is important in studies of the physics of the lower ionosphere, and measurements
of them
are difficult to perform by other means. Comparison of results of model
calculations with
measured transfer function data will allow the remote sensing of these
parameters, on a
continuous-time basis, and on a global scale.
In addition, simulation results presented here show that this technique is
capable of
detecting evidence of multiple-mode propagation and mode coupling, when these
effects
are significant. It is important to consider the mode structure of the incoming
wave at the
receiver, since common receivers use small magnetic loop antennas, which are not
directly sensitive to the vertical electric field of the wave. This effect is of
particular
importance when the direction of arrival of the incoming VLF wave is nearly in
the plane
of the loop.
Finally, application of this technique to a data set collected at Palmer
Station,
Antarctica shows that the measured transfer functions of the propagation channel
do not
exhibit sharp variations in frequency. This result is in general agreement with
the
predictions of the 3-D waveguide propagation model, and indicates that, on these
paths,
the mode structure of the received signal is dominated by a single waveguide
mode.
However, there are quantitative differences between the measurements and the
predictions, which indicate that there are significant differences from the
assumed
electron density and/or temperature profiles.
6.2. SUGGESTIONS
FOR FUTURE WORK
Since the focus of this work has been the development of the measurement
technique,
only a very limited data set has been analyzed. There is much room for the
application of
the techniques described here to larger data sets, and to propagation paths
expected to
show effects of more geophysical interest. For example, study of larger data
sets from
Palmer could determine the diurnal variations in the transfer functions of the
four
propagation channels studied in this work, as well as variations with
geomagnetic
activity. These data, in combination with model calculations using assumed
ionospheric
density and temperature profiles, could be used to test the extent to which the
profiles are
influenced by other factors, since the frequency-dependent processing yields
many points
of comparison with the models. In addition, a study of shorter, land-based paths
might
be expected to capture LEP-induced VLF signal perturbations which show
detecdifferences across the available band of frequency; these would yield
direct
measurements of the transfer function perturbation, which could be directly
compared
with model results for the scattering from such perturbations.
Finally, some of the processing techniques investigated as a part of this work
may
hold promise for other analysis of field data in the future. This is
particularly true given
the rapid advance in the capability of inexpensive signal processing hardware.
Examples
of these techniques would include impulsive noise removal from broadband
recordings,
and the use of the phase of the subionospheric MSK signals as a time reference
for phase
studies of magnetospheric signals.
6.2.1. Impulsive
Noise Removal for Broadband Recordings
One of the processing steps described in Chapter 3 is applied to the broadband
data,
before filtering and decimation to a lower sampling rate: the suppression of
impulsive
interference. In this work, we simply followed standard practice, and filtered
and clipped
the broadband data samples. However, some experimentation was done with more
sophisticated techniques of interference suppression, by recognizing the
impulses and
excising them from the data set. This produced little gain in the narrowband
analyses
performed here, but may be of value to broadband studies, where spectral
representations
are used which cover 10 KHz or so, and in which the masking effect of the
impulsive
interference from spherics is significant.
6.2.2. Implementation
of Real-Time Processing
The data in this study were gathered via broadband recording, followed by post-
processing to yield the propagation channel transfer function estimates. While
this is
useful for development of these techniques, this implementation is not suifor
general use
in the field, because limitations on recording rates and capacity do not allow
the analysis
of the large data sets typically required to reveal novel geophysical phenomena.
In
addition, the upper frequency is limited by the recording rate of the sampling
hardware,
and does not allow the analysis of frequencies above about 33 KHz in the current
implementation.
To overcome these limitations, the techniques described in this work must be
implemented in a field-equipment receiver, and must operate in real time.
Fortunately,
the complex processing is largely done after the individual transmitter has been
isolated
from the broadband data, and resampled at a lower data rate. The isolation of
individual
transmitters is a much simpler process, but must be implemented efficiently,
since this
processing must be applied to the entire broadband data stream.
This partition of the processing into a high data rate, low complexity, front
end,
and a low data rate, high complexity, back end, is not atypical of modern signal
processing systems, and may be used to advantage in the design of a field
receiver. A
block diagram of one possible receiver is shown in Figure 68. In this receiver
design, the
entire processing path is digital. There are identical channel units, each of
which
performs the estimation process on one propagation path, and which plug into a
front
end containing the A/D converter and broadband filtering and clipping functions.
Each
of the narrowband cards can analyze a single transmitter in the 10-60 KHz range.
In
addition to the narrowband outputs (one for each channel card), this receiver
provides a
broadband digital output covering the frequencies below about 10 KHz, which
contains
the whistlers and other magnetospheric signals of interest.
Figure 68. Possible implementation of the processing described in this
dissertation, in a
form suitable for use in the field. In addition to the propagation channel
transfer function
estimation, this receiver also provides a broadband output, which contains the
whistler
data. There is a single A/D converter and front-end processor, and a number of
identical
channel cards, each of which performs the transfer function estimation for one
signal
path.
In this receiver, the processing is divided among DSP elements with different
capabilities, which allows us to perform those filtering and conversion
functions which
occur on the broadband data, and therefore at high rates, using fixed-point
processing,
and to perform the more complex estimation procedure using a slower,
floating-point,
DSP unit. Approximate processing loads for the different types of processors
are:
Table 2. Processing Rate Requirements for Field Receiver
Element
Math Type
Rate
Low-Pass Filter, Clip
Fixed, 12 bits
40 MOPS
High-Pass Filter, Clip, Resample
Fixed, 12 bits
40 MOPS
Downconvert, Filter, Resample
Fixed, 12 bits
50 MOPS
Demodulate, Spectral Estimation
Floating-point
2 MFLOPS
Although the processing algorithms used in the final step could be scaled into
fixed-point form, the low computing rate required (due to the fact that all this
processing
is happening at the lowered sample rate of 1 KHz) makes a floating-point
processor
attractive, for ease of implementation.
In such a receiver, the channel cards are identical, and the parameters which
govern
their operation (center frequency, modulation, bit rate, etc.) could be
downloaded at the
start of operation from a host computer, which could also provide the recording
system.
For example, a receiver designed along these lines could interface to the host
computer
via a SCSI interface, and appear as a peripheral to the host system. This would
allow the
entire code run by all the processors in the receiver to be loaded from the host
computer,
and additional types of processing could be accommodated with software changes.
An additional advantage of such an approach is that channel cards may be
programmed to analyze other types of signals. For example, cards to process
LORAN and
Omega signals could be built, and operate alongside those working on the MSK-
modulated signals. Although not shown in the diagram, a data path from the
front-end
lowpass filter output could be provided to the channel cards, and a channel card
could be
used to detect interesting events, such as whistler arrivals; these triggers
could be used to
"freeze" a buffer on the broadband data, and dump it to the recording
system for later,
off-line, analysis.
6.2.3. Measurements
at Separated Stations
One of the possible applications of the receiver described above is to make
measurements
at separated stations, using clocks synchronized with GPS receivers. Clocks may
be
synchronized with this method to better than 10 nsec accuracy, which corresponds
to 0.1
of phase at 30 KHz. This interferometry technique may be able to reconstruct the
field
scattered by a transient disturbance to the earth-ionosphere waveguide, and thus
provide
information on the source and location of the disturbance.
6.2.4. Studying
Shorter, Overland Paths
Although all the data analyzed in this work was collected on long propagation
paths,
largely over seawater, the simulation results presented in Chapter 4 indicate
that the
shorter, overland, paths should provide the largest deviations with frequency,
and, hence,
the largest sensitivity for measurements of the electron density and temperature
profiles
using these techniques. Such data could be gathered, for example, at Stanford,
from NSS
(Annapolis, MD), NAA (Cutler, ME), and particularly the Air Force transmitter
operating at 48.5 KHz (Silver Creek, NE). The latter is a path which often
experiences
transient disturbances, but the higher frequency of this transmitter would
require a
recording system with larger bandwidth than that used in this work, which only
allowed
measurements to about 33 KHz.
6.2.5. Use of
Other Transmitters
The VLF and LF regions of the spectrum are crowded with transmitters, in
addition to
the large MSK transmitters used in this study. The transmitters most suifor use
in these
studies are the LORAN transmitters, which transmit a coded pulse format centered
at
100 KHz, and the Omega transmitters, which transmit long pulses in the 10-14 KHz
frequency range. Both systems are used for navigation: the Omega system is less
accurate,
but has global coverage; LORAN covers mainly coastal and populated areas, but
provides
more accurate position information.
An advantage of the LORAN transmissions is that they consist of accurately-
controlled pulse trains, which are transmitted first by a master transmitter,
and echoed
after a known delay by a series of slave transmitters. The pulse format is
controlled so that
the third zero crossing is accurately placed, and the duration of the pulses is
short enough
so that the spectrum of a LORAN transmitter extends up to 40 KHz or so away from
100 KHz. This broad spectrum may prove useful in multi-frequency measurements in
the
future, though, again, measurements can only be made with hardware built for the
purpose. Since the transmissions are pulsed, the estimation technique used in
this case
may estimate the impulse response of the earth-ionosphere waveguide, rather than
performing the estimation in the frequency domain.
The other widely-used navigation system operating in the VLF portion of the
spectrum is the Omega system, which consists of eight transmitters scattered
around the
globe. These transmitters operate using signals which consist of a series of
pulses in the
10-14 KHz range. An advantage of the Omega system, as with LORAN, is that no
message information is being sent, and so the demodulator need only estimate the
time of
arrival of the signal of interest, rather than estimating a bit sequence as
well. Also, the
Omega system transmissions cover a range of frequencies much wider than that
occupied
by the MSK transmitters. The main drawback to the Omega transmitters is that,
because
of the pulsed-CW format, they cover their entire list of frequencies every 10
seconds. This
provides sufficient resolution for studies of the ambient electron density and
temperature
profiles, but is not sufficient to study LEP-induced perturbations and similar
events.
One final possibility which should be mentioned in connection with this work is
the eventual use of purpose-built transmitters, sending signals coded to
optimize the
estimation process, rather than to communicate data. Such a modulation would
ideally
provide a wide bandwidth (limited only by antenna matching considerations), with
equal
power density in all the frequency bands of interest. This is akin to
traditional spread-
spectrum techniques, because a spreading code is used to widen the transmitted
bandwidth far beyond the baseband spectrum of the data involved.
BIBLIOGRAPHY
Adams, C. D. D., and R. L. Dowden, VLF group delay of lightning-induced electron
precipitation echoes from measurements of phase and amplitude perturbations at
two frequencies, J. Geophys. Res., 95, 2457, 1990.
Allcock, G. McK., R. A. Helliwell, C. K. Branigan, and J. C. Murphy, Whistler
and other
very low frequency phenomena associated with the high-altitude nuclear explosion
on July 9, 1962, J. Geophys. Res., 68, 735, 1963.
Barr, R., M. T. Reitveld, H. Hopka, and P. Stubbe, Effects of a heated patch of
auroral
ionosphere on VLF-radio wave propagation, J. Geophys. Res., 90, 2861, 1985.
Budden, K. G., The Wave-guide Mode Theory of Wave Propagation, Prentice-Hall,
Englewood Cliffs, New Jersey, 1961.
Budden, K. G., The Propagation of Radio Waves, Cambridge University Press,
Cambridge,
1985.
Burgess, W. C., and U. S. Inan, Simultaneous disturbance of conjugate regions in
association with individual lightning flashes, Geophys. Res. Lett., 17, 259,
1990.
Burgess, W. C., and U. S. Inan, The role of ducted whistlers in the
precipitation loss and
equilibrium flux of radiation belt electrons, J. Geophys. Res., 98, 15643, 1993.
Carlson, A. B., Communication Systems, McGraw-Hill, New York, 1986.
Carpenter, D. L., and J. W. LaBelle, A study of whistlers correlated with bursts
of electron
precipitation near L=2, J. Geophys. Res., 87, 4427, 1982.
Carpenter, D. L., and U. S. Inan, Seasonal, latitudinal, and diurnal
distributions of
whistler-induced electron precipitation events, J. Geophys. Res., 92, 3429,
1987.
Cotton, P. D., and A. J. Smith, Signature of burst particle precipitation on VLF
signals
propagating in the Antarctic earth-ionosphere waveguide, J. Geophys. Res., 96,
19,375, 1991.
Crombie, D. D., Further observations of sunrise and sunset fading of VLF
signals, Radio
Science, 1, 47, 1966.
Crombie, D. D., A. H. Allen, and M. Newman, Phase variations of 16 kc/s
transmissions
from Rugby as received in New Zealand, Proc. IEE, 105(B), 301, 1958.
Davies, K., Ionospheric Radio, Peter Peregrinus, Ltd., London, 1990.
Dowden, R. L., and C. D. D. Adams, Modal effects on amplitude perturbations on
subionospheric signals (Trimpis) deduced from two-frequency measurements, J.
Geophys. Res., 94, 1515, 1989.
Dowden, R. L., and C. D. D. Adams, Ionospheric heating by VLF transmitters?, J.
Atmos.
Terr. Phys., 54, 1051, 1992.
Dungey, J. W., Loss of Van Allen electrons due to whistlers, Planet. Space Sci.,
11, 591,
1963.
Ferguson, J. A., and Snyder, F. P., Approximate VLF/LF waveguide mode conversion
model, Tech. Doc. 400, Naval Ocean Systems Center, San Diego, California, 1980.
Ferguson, J. A., F. P. Snyder, D. G. Morfitt, and C. H. Shellman, Long-wave
propagation
capability and documentation, Tech. Doc. 1518, Naval Ocean Systems Center, San
Diego, California, 1989.
Field, E. C., and Lewinstein, Amplitude-probability distribution model for
VLF/ELF
atmospheric noise, IEEE Trans. Commun., COM-26, 83, 1978.
Fraser-Smith, A. C., and R. A. Helliwell, The Stanford University ELF/VLF
radiometer
project: measurement of the global distribution of ELF/VLF atmospheric noise,
Proc. IEEE Internat. Symp. on Electromag. Compatability, Wakefield,
Massachusetts,
1985.
Freidel, R. H. W., and A. R. W. Hughes, Effect of magnetic field distortions on
gyroresonance at low l values, , J. Geophys. Res., 98, 1,669, 1993.
Gallejs, J., Terrestrial Propagation of Long Electromagnetic Waves, Pergamon
Press, Oxford,
1972.
Glukhov, V. S., V. P. Pasko, and U. S. Inan, Relaxation of transient lower
ionospheric
disturbances caused by lightning-whistler induced electron precipitation bursts,
J.
Geophys. Res., 97, 16971, 1992.
Harris, F. J., On the use of windows for harmonic analysis with the Discrete
Fourier
Transform, Proc. IEEE, 66, 51, 1978.
Helliwell, R. A., Whistlers and Related Ionospheric Phenomena, Stanford
University Press,
Stanford, California, 1965.
Helliwell, R. A., Intensity of discrete VLF emissions, in Particles and Fields
in the
Magnetosphere, B. M. McCormac (ed.), p. 292, 1970.
Helliwell, R. A., and D. L. Carpenter, Whistlers excited by nuclear explosions,
J. Geophys.
Res., 68, 4409, 1963.
Helliwell, R. A., J. P. Katsufrakis, and M. L. Trimpi, Whistler-induced
amplitude
perturbation in VLF propagation, J. Geophys. Res., 78, 4679, 1973.
Imhof, W. L., E. E. Gaines, H. D. Voss, J. B. Reagan, D. W. Datlowe, J. Mobilia,
R. A.
Helliwell, U. S. Inan, J. Katsufrakis, and R. G. Joiner, The results from the
SEEP
active space plasma experiment: effects on the ionosphere, Radio Science, 20,
511,
1985.
Inan, U. S., VLF heating of the lower ionosphere, Geophys. Res. Lett., 17, 729,
1990.
Inan, U. S., and D. L. Carpenter, On the correlation of whistlers and associated
subionospheric VLF/LF perturbations, J. Geophys. Res., 91, 3106, 1986.
Inan, U. S., and D. L. Carpenter, Lightning-induced electron precipitation
events
observed at L~2.4 as phase and amplitude perturbations on subionospheric VLF
signals, J. Geophys. Res., 92, 3293, 1987.
Inan, U. S., D. L. Carpenter, R. A. Helliwell, and J. P. Katsufrakis,
Subionospheric VLF/LF
phase perturbations produced by lightning-whistler induced particle
precipitation,
J. Geophys. Res., 90, 7457, 1985.
Inan, U. S., D. C. Shafer, and W. Y. Yip, Subionospheric VLF signatures of
nighttime D-
region perturbations in the vicinity of lightning discharges, J. Geophys. Res.,
93,
11455, 1988a.
Inan, U. S., T. G. Wolf, and D. L. Carpenter, Gepgraphic distribution of
lightning-
induced electron precipitation observed as VLF/LF perturbation events, J.
Geophys.
Res., 93, 9841, 1988b.
Inan, U. S., M. Walt, H. D. Voss, and W. L. Imhof, Energy spectra and pitch
angle
distributions of lightning-induced electron precipitation: analysis of an event
observed on the S81-1 (SEEP) satellite, J. Geophys. Res., 94, 1379, 1989.
Inan, U. S., T. F. Bell, and J. V. Rodriguez, Heating and ionization of the
lower
ionosphere by lightning, Geophys. Res. Lett., 18, 705, 1991.
Inan, U. S., J. V. Rodriguez, S. Lev-Tov, and J. Oh, Ionospheric modification
with a VLF
transmitter, Geophys. Res. Lett., 19, 2071, 1992.
Inan, U. S., J. V. Rodriguez, and V. P. Idone, VLF signatures of
lightning-induced heating
and ionization of the nighttime D region, Geophys. Res. Lett., 20, 2355, 1993.
ITU, International Frequency List, Volume I, International Telecommunications
Union,
1979.
Ljung, Lennart, System Identification Theory for the User, Prentice-Hall,
Englewood Cliffs,
New Jersey, 1987.
Morfitt, D. G., and C. H. Shellman, MODESRCH: an improved computer program for
obtaining ELF/VLF/LF mode constants in an earth-ionosphere waveguide, Interim
Report Number 77T, Propagation Technology Division, Naval Ocean Systems
Center, San Diego, California, 1976.
Motchenbacher, C. D., and Connelly, J. A., Low-Noise Electronic System Design,
John
Wiley and Sons, New York, 1993.
Naval Ocean Systems Center, Fixed VLF/LF Data Base, Naval Ocean Systems Center,
San
Diego, California, 18-OCT-89.
Oppenheim, A. V., and R. W. Schafer, Digital Signal Processing, Prentice-Hall,
Inc.,
Englewood Cliffs, New Jersey, 1975.
Paschal, E. W., Calibration of very-low-frequency receiving systems for Stanford
University (internal report), Stanford Electronics Laboratories, Stanford,
California, 1978.
Paschal, E. W., The design of broad-band VLF receivers with air-core loop
antennas
(internal report), Stanford Electronics Laboratories, Stanford, California,
1988.
Pasko, V. P., and U. S. Inan, Recovery signatures of lightning-associated VLF
perturbations as a measure of the lower ionosphere, J. Geophys. Res. (in press),
1994.
Poulsen, W. L., Modeling of very low frequency wave propagation and scattering
within
the earth-ionosphere waveguide in the presence of lower ionospheric
disturbances,
PH. D. dissertation, Stanford University, Stanford, California, 1991.
Poulsen, W. L., T. F. Bell, and U. S. Inan, Three-dimensional modeling of
subionospheric
VLF propagation in the presence of localized D-region perturbations associated
with lightning, J. Geophys. Res., 95, 2355, 1990.
Poulsen, W. L., U. S. Inan, and T. F. Bell, A multiple-mode, three-dimensional
model of
VLF propagation in the earth-ionosphere waveguide in the presence of localized
D-
region disturbances, J. Geophys. Res., 98, 1705, 1993a.
Poulsen, W. L., T. F. Bell, and U. S. Inan, The scattering of VLF waves by
localized
ionospheric disturbances produced by lightning-induced electron precipitation,
J.
Geophys. Res., 98, 15553, 1993b.
Ramo, S., J. R. Whinnery, and T. Van Duzer, Fields and Waves in Communication
Electronics, John Wiley and Sons, New York, 1984.
Rees, M. H., Auroral ionization and excitation by incident energetic electrons,
Planet.
Space Sci., 11, 1209, 1963.
Rodriguez, J. V., Modification of the earthís ionosphere by very-low-frequency
transmitters, PH. D. dissertation, Stanford University, Stanford, California,
1994.
Rodriguez, J. V., U. S. Inan, and T. F. Bell, D-region disturbances caused by
electromagnetic pulses from lightning, Geophys. Res. Lett., 19, 2067, 1992.
Smith, A. J., and P. D. Cotton, The Trimpi effect in Antarctica: observations
and models,
J. Atmos. Terr. Phys., 52, 341, 1990.
Taranenko, Y. N., U. S. Inan, and T. F. Bell, Optical signatures of
lightning-induced
heating of the D-region, Geophys. Res. Lett., 19, 1815, 1992.
Thompson, N. R., Whistler mode signals: spectrographic group delays, J. Geophys.
Res.,
86, 4795, 1981.
Tolstoy, A., and T. J. Rosenberg, A quasi three-dimensional propagation model
for
subionospherically propagating VLF radio waves, Radio Science, 20, 535, 1985.
Tolstoy, A., T. J. Rosenberg, U. S. Inan, and D. L. Carpenter, Model predictions
of
subionospheric VLF signal perturbations resulting from localized, electron-
precipitation induced, ionization enhancement regions, J. Geophys. Res., 91,
13473,
1986.
Turtle, J. P., E. C. Field, Jr., C. R. Warber, and P. R. McGill, Low-frequency
transverse
electric atmospheric noise: measurement and theory, Radio Science, 24, 325,
1989.
Uman, M. A., The Lightning Discharge, Academic Press, Orlando, Florida, 1987.
Warber, C. R., and E. C. Field, A long wave TE/TM noise prediction model, Proc.
of
NATO-AGARD conference on ELF/VLF/LF radio propagation and systems aspects,
Brussels, Belgium, 1993.
Watt, A. D., VLF Radio Engineering, Pergamon Press, Oxford, 1967.
Wolf, T. G., Remote sensing of ionospheric effects using very low frequency
radio signals,
PH. D. dissertation, Stanford University, Stanford, California, 1990.
Wolf, T. G., and U. S. Inan, Path-dependent properties of subionospheric VLF
amplitude
and phase perturbations associated with lightning, J. Geophys. Res., 95, 20997,
1990.
Yip, W.-Y., U. S. Inan, and R. E. Orville, On the spatial relationship between
lightning
discharges and propagation paths of perturbed subionospheric VLF/LF signals, J.
Geophys. Res., 96, 249, 1991.