Below are the first 10 and last 10 pages of uncorrected machine-read text (when available) of this chapter, followed by the top 30 algorithmically extracted key phrases from the chapter as a whole.
Intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text on the opening pages of each chapter.
Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.
Do not use for reproduction, copying, pasting, or reading; exclusively for search engines.
OCR for page 37
3
Oceanographic and
Atmospheric Applications of
Spatial Statistics and Digital
Image Analysis
James ]. Simpson
Scripps Institution of (oceanography
3.l Introduction
Historically, remote sensing of the environment has helped process-oriented
studies examine individual aspects of the physics, chemistry, and biology of
the earth-ocean-atmosphere system. From such studies, several problems of
global significance have emerged that are cross-disciplinary in nature. Ex-
amples include acid rain, the increase in atmospheric carbon dioxide, antici-
pated depletion of the ozone layer, ET Nino-related mollifications in weather
and ocean circulation with their resultant effects in agriculture en c! fisheries,
and tropical rain forest destruction by fires of human origin. It has been rec-
ognized that the key to progress on these and other cross-clisciplinary issues
in earth science during the decade of the 1990s probably will be addressing
those questions that concern the integrated functioning of the earth as a
system (EOS, 1984~. The hyclrologic cycle, the biogeochemical cycle, and
climate processes are the fundamental processes that integrate the earth as
a system, and thus each of these cycles must be examined in detail and on
a global scale if meaningful progress is to be made on the problems cited
above.
I,arge-scaTe, synoptic observations of a winkle variety of phenomena (e.g.,
sea surface temperature, sea level, winct stress, ozone concentration, radi-
37
OCR for page 38
38
ation heat balance, sea ice, vegetation, and near-surface current velocity)
are required to study these cycles properly, and remote sensing provides the
only practical way to collect synoptically the necessary data on the global
scale. The magnitude of the data set, however, may preclude understanding
unless it can be distilled and synthesized into organized patterns that can
be related meaningfully to the underlying governing physics, chemistry, or
biology of these global cycles and processes. Spatial statistics and math-
ematical methods of digital image analysis provide mechanisms for such a
synthesis.
This chapter discusses a few areas of spatial statistics and data represen-
tation useful in digital image analysis that have immediate application in the
area of remote sensing of the earth. Some special considerations needed for
the correct analysis of remotely-sensed data are presented and a few critical
research areas are identified.
3.2 Selected Analysis Areas
In this section three mathematical methods useful in digital analysis of se-
quences of remotely sensed images are presented. These methods were cho-
sen because they have a broad range of applications in earth science. These
methods are non-Bayesian in nature; a presentation on Bayesian methods
used in image analysis is given in chapter 2 of this report.
Several abbreviations commonly used in remote sensing and the earth
sciences appear in this chapter. Each is defined at the place it first occurs in
the text. For easy reference, there is an appendix containing abbreviations
and their definitions at the end of this chapter.
3.2.] Principal Component Analysis
General Concepts
Principal component analysis (PCA) is a multivariate statistical technique
that can be applied to all forms of multispectral or multi-temporal image
data and is most commonly applied in the general context of arbitrary mul-
tivariate data. Forms of PCA have been used to study pattern classification
(Geladi e! al., 1989), sequential segmentation (Esbensen and Geladi, 1989),
spatial patterns of variability in sea surface temperature and phytoplankton
pigment observed in satellite data (LagerIoef, 1986), and in algorithms for
cloud removal from satellite data (e.g., GalIaudet and Simpson, 1991a).
OCR for page 39
39
The principal component transformation (PCT) is the basis for all these
analyses. The PCT is a linear transformation that isolates uncorrelated
linear combinations of a given set of variables in such a way that each element
in this combination represents a decreasing amount of variance in the original
variables (Ingerbritsen and Lyon, 1985~. This linear transformation defines
a new set of coordinate axes for the data such that (13 the transformed
origin is at the mean of the data distribution (LiDesand and Kiefer, 1987),
(2) the transformed coordinate axes are mutually orthogonal (Jensen, 1986),
and (3) the transformed coordinate axes are in the directions of maximum
variance (Jensen, 1986~. Below, two examples of PCA are developed, one
using multispectral image data and the other using multi-temporal image
data. An extensive treatment of the use of PCA techniques in atmospheric
sciences and oceanography is given by Preisendorfer (1988~.
Multispectral Data Application
The ability to accurately and automatically segment clouds in remotely
sensed imagery is critically important to a broad range of disciplines in earth
science. Clouds significantly affect the net heating of the atmosphere and the
underlying ocean-land surface by modifying solar and terrestrial radiation
(Ohring and CIapp, 1980~. This net radiative heating governs the thermody-
namics and dynamics of the atmosphere, which in turn influence the forma-
tion and dissipation of clouds (e.g., Matveev, 1984~. The potential feedback
effects associated with this cloud-radiation interaction are among the great-
est sources of uncertainty in determining the relation between changes in
external conditions such as solar radiation and atmospheric carbon dioxide
concentration and changes in climate (e.g., Henderson-SeDers, 1982; Ra-
manathan, 1987~. Clouds also affect our ability to remotely sense the prop-
erties of the atmosphere, ocean, and land; such observations are needed, for
example, in weather prediction (e.g., Pailleux, 1986), oceanography (e.g.,
Eckstein and Simpson, 1990a,b), and agriculture (e.g., rainfall, Browning,
1986~. The PCT can be used with multispectral image data to robustly
segment clouds from natural images. An example of such an application is
given below.
The Principal Component Transformation Split-ancI-Merge Clus-
tering (PCTSMC) Algorithm. Here, the development of the PCTSMC
algorithm is presented in abbreviated form. Mathematical details of the al-
gorithm are given in Gallaudet and Simpson (1991a). An AVHRR infrared
OCR for page 40
40
image, designated T. consisting of brightness temperatures calibrated to de-
grees Celsius in bands 3, 4, and 5 (hereinafter referred to as T3, T4, Ts) is
first differenced to construct a 2-banded differenced image D (i.e., T3—T4,
T3—Thy. This differenced image now contains aD the information of the
original infrared image needed for cloud detection, but only in two linearly
independent bands. Next, this differenced image is transformed using the
PCT. The result is a 2-banded transformed differenced image in which all
interband correlation is destroyed (Masher, 1987), thus making it a Togi-
cal preprocessor to the segmentation operation (step 3) of the PCTSMC
algorithm.
Step 3 of the PCTSMC algorithm performs image segmentation using a
split-and-merge clustering procedure (e.g., PavIidis, 1977; Seddon and Hunt,
1985; Richards, 1986) on the PC transform of the di~erencecl image. This
results in a segmented image in which the natural spectra] classes in the
original image are separated into distinct groups (i.e., land versus ocean
versus cloud). The method of clustering that is employed in the PCTSMC
algorithm combines both the partitional and hierarchical approaches. It
consists of a partitional clustering algorithm augmented by a splitting-and-
merging step at each iteration. Combining a partitional with a hierarchical
method has several advantages over the use of either method alone: (1)
pure hierarchical methods are not appropriate for complex data (MuerIe
and Allen, 1968; Fukada, 1980; Jain and Dubes, 1988~; (2) pure hierar-
chical methods are more appropriate for data that is to be partitioned on
the basis of both local and global information, rather than global informa-
tion only plain and Dubes, 1988~; (3) pure hierarchical methods impose a
taxonomic structure on the data (Anderberg, 1973), which is not character-
istic of cloud-containing AVHRR imagery; (4) pure hierarchical methods are
order dependent i.e., the resulting segmentation will vary ~lepending upon
the order in which the regions are split and merged (Cheevasuvit et al.,
1986~; this often results in less than optimal segmentations of the data;
(5) pure partitional algorithms often converge to local minima of the clus-
tering criterion function (Pairman and Kittler, 1986; Jain and Dubes, 19SS);
(6) the combined approach is more efficient than pure merging or pure
splitting methods of region detection (Paviidis, 1977; Richards, 1986~; and
(7) the combined approach is less dependent on the initial segmentation,
and therefore is more capable of recovering all of the natural clusters in the
data (Seddon and Hunt, 1985; Jain and Dubes, 1988~; this is because the
number of clusters in the initial segmentation need not be the same as those
O ~ ~
is
OCR for page 41
41
that actuaby exist in the given data (Seddon and Hunt, 1985; Pairman en c!
Kittler, 1986).
In the final step of the PCTSMC cloud screening procedure, the data
are retransformed back to the feature space, and each of the clusters is
labeled as either a cloud or non-cloud. In the first three steps discussed
above, the operations performed on the data were entirely unsupervised-
i.e., no a priori knowledge was required. Hence, they apply to an arbitrary
image segmentation. In this fourth step, expert knowledge is introduced to
perform a boolean classification. Rules appropriate for land versus cloud
versus ocean separation in AVHRR image data are given in Simpson and
Humphrey (1990) and GaHaudet and Simpson (1991a).
Plate 3.1a shows AVHRR Band 2 data; clouds appear as white or gray
tones. The coastline is white in this panel. Plate 3.1b shows AVHRR Band 4
infrared temperature; coldest temperature is white and warmest tempera-
ture is black. In this panel, the coastline is black. The segmented image
produced by the PCTSMC algorithm is shown in Plate 3.1c and the final
cloud-masked sea surface temperature is shown in Plate 3.16. In this final
panel, the warmest ocean temperature is white, cooler ocean temperatures
appear as shades of gray, and cloud contaminated pixels and land are shown
as black. Land was masked from the Plate 3.1 images using a recursive
polygon fill algorithm (Simpson, 1991) and is rendered either white or black
depending on the gray scale mapping used in the individual panel.
Multi-Temporal Data Application
A major objective of earth science studies is to identify spatial patterns of
variance in temporal sequences of images. Examples include the analysis of
variability in sea surface temperature structure in oceanic current systems
(e.g., Lagerioef, 1986) and in seasonal and interannual variation in phyto-
plankton abundance (e.g., Strub et al., 1990~. The form of PCA used in
such studies is generally referred to as empirical orthogonal function (EOF)
analysis: "empirical" because the functions arise from the data themselves
and "orthogonal" because they are uncorrelated. Note that closed-form
mathematical functions generally cannot be used as the basis functions rep-
resenting the complex images observed in nature.
Empirical orthogonal functions are useful in work with large data sets.
The method separates a data set D(x,t) into spatial components Fi(X) and
OCR for page 42
42
temporal components Alit) such that
NT
D(x, t) = ~ Fi(X) · Ai(t)
i=1
(3.1)
where ~ is the number of EOFs computed. Note that x and t repre-
sent generalized spatial and temporal coordinates, i.e., x represents the set
(x~,...,xArx) and ~ represents the set {ti,...,t~. In matrix notation,
D(x,t) is an NX x IN matrix, with NX representing the number of posi-
tions x and ~ the number of time steps t. Each column of the F matrix is
an EOF. Each EOF has NX values, and there are JVT EOFs. Thus, F is an
NX x JVT matrix. The time series matrix A has ~ rows and JVT columns.
\
Forms of Normalization. The computations may be implemented with
two different normalization schemes because the EOF representation decom-
poses the space-time series, using separation of variables, into a sum of prod-
ucts of temporal amplitudes, Ain't), modulating spatial patterns of variance,
Fi(X). Note that, in physics, separation of variables occurs widely because
of the form of the underlying differential equations. Here, the motivation is
to provide a compact statistical representation of the data in which spatial
patterns in variance can be distinguished from temporal patterns. Either
the Apt) or the Fi(X) can retain the same physical units as the original data.
Historically, a normalization scheme (method 2) has been used in which the
temporal components of the decomposition retained the same units of the
data. More recent studies in oceanography and atmospheric science (e.g.,
Barnett and Patzert, 1980) prefer to use EOFs with the same units as the
original data for easier interpretation of the spatial patterns of variance re-
sulting from the EOF analysis (method I). Both methods of normalization
are equivalent and both are given here for completeness. (The summation
notation is used for consistency with the vast majority of published studies.)
In the first method
Fi(X) · Fjord)
NX
~ Fi~xk)Fj (Xk) = Aitij, (3.2,
k=1
where Fi and Fj are columns of the F matrix and {ij is the Kronecker
delta function. The coefficient Ai is an eigenvalue of the covariance matrix
C = D D , where DT is the transpose of D. The rows of the time series
OCR for page 43
43
matrix are orthonormal, i.e.,
Ai(t) Aj(t)
NT
=
JVT
NT ~ Ai(tk)Aj(tk) = ttj ~
(3.3)
where Ai and Aj are rows of A. The EOFs then give the variance in the
same units as the data.
The second normalization method forces the EOFs to be orthonormal,
ye.,
Hence'
Fi(x) Fj(x) = bij.
TF F—~
_ ,
(3.4)
(3.5)
where ~ is the identity matrix. The rows of the time series matrix are
orthogonal:
Ai(t) Aj(t) = I. (3.6)
With this normalization, the EOFs no longer have the same units as the
data.
Theoretical Basis. The objective of EOF analysis is to represent a given
matrix of data D by the product F A. In the discussion that follows, the
second method of normalization is used, and thus the matrix F satisfies
equation (3.5~.
The equations governing EOF theory can be derived from the eigen-
equation of the covariance matrix,
C · F = F · A ,
(3.7)
where C is the NX x NX covariance matrix = D D , A is the diagonal
matrix of eigenvalues, and F is the matrix of EOFs. Because C is a real
symmetric positive definite matrix, the elements of A are positive. By (lefi-
nition of the covariance matrix C, (3.7) can be rewritten as
{ ~ }
(3.8)
The eigenvectors in the matrix F are now used to define the matrix A =
FT · D. This is the principal component transform. Thus, F · FT = ~ and
necessarily
D = F.A.
(3.9)
OCR for page 44
44
At this point the number of EOFs is equal to JVT. The purpose of EOF
analysis, however, is to produce a representation of the data that is more
compact than the original data set. For this goal to be met, it is necessary
for most of the variance (65 to 80~) in the original data set to be contained
in the first few EOFs. (By convention, the eigenvalues are ordered with the
largest being first.) If this is true, then F may be truncated to produce a
new matrix F having just NT columns, where hT ~ Ad. Then D is the
approximation to D computed from the EOF decomposition:
D~D=F A.
(3.10)
The decomposition is the most efficient representation of D with regard to
a mean square error criterion (Davis, 1976~.
Finally, it should be noted that if D contains a noninteger number of
cycles of a sinusoidally varying variance mode, then the variance represented
by the associated eigenvalue may not agree with the actual variance. For
example, the variance of a sine wave differs for one-half and one full cycle,
but the covariance matrix C = D D is the same. Hence, the EOF will
return the variance for an entire cycle when the data represent one-half
cycle.
A Computational Method. The empirical orthogonal functions can be
computed from a singular value decomposition (SVD) of the data D. When
NX > By, then (Press et al., 1986)
D = U . W . VT
(3.11)
where U is an NX x NX orthogonal matrix with only NT linearly indepen-
(lent columns, W is an NX x IVY matrix with an 1\7 x HI diagonal upper
portion having positive or zero elements, and V is an NT x ~ orthogonal
matrix. Note that the diagonal elements in the IVI x NT upper portion of
W contain the singular values of D. Cad this portion of W the matrix W'.
The U and V matrices are orthogonal in the sense that their columns are
orthonormal, that is,
TU U = I
art
Vl V = I. (3.12)
Using this decomposition, the covariance eigen-equation can be written as
C U=
I NT } U = ~ MU · W · VT) {V · WT · UT) · U. (3.13)
OCR for page 45
45
which simplifies to
C.U= ~ u.W2
NT
.
(3.14)
This is equivalent to the eigen-equation (3.7) with F—the highest order
columns of U and ~—W2/~, using the second method of normalization.
Hence, the diagonal elements Wi of W' relate to the eigenvalues Ai of the
covariance matrix C via the equation
~ = A.
(3.15)
Thus, the SVD returns the EOFs as the highest order columns of U. and
A = W'~VT
Generalizations to Two Dimensions. The EOF decomposition is not
confined to space-time series where space is one-dimensional. If the posi-
tions are actuary (x, y) coordinates, such as latitudes and longitudes or lines
and samples from images, then if each time step covers the same locations,
D(x,y,t) can be reduced to D(x',t). The sequence of positions x' is con-
structed by setting x~ = (at, yip, X2 = (X2, yI), · - -, XNX = (XNX, YI ), · ~ ~ .
XNX.NY = (XNX, ANY ). If one considers x and y to be rows and columns of an
individual time step matrix, then x' is equivalent to concatenating the rows
of D together. Whether one concatenates the rows or the columns is imma-
terial. (Note, EOF analysis with two-dimensional images is actually three
dimensional: x, y, and t. To render the problem tractable, the spatial di-
mensions must be concatenated so that one can work with a two-climensional
space-time matrix.) Then the individual EOFs Fi are also vectors in x' and,
in order to map the patterns of variation associated with the EOFs, this
vector must now be parsed back into its rows and columns.
Example. The purpose of this example is to show that EOF analysis can
be useful for determining the dominant patterns of spatial variance in a
sequence of images. (The normalization scheme used in this example is that
of method 1.) For this purpose, a sequence of eight images was constructed
by superimposing an image of an inclined plane with that of a disk that
is out of phase with the inclined plane. Note also that the range of data
in the image of the inclined plane is about twice that of the image of the
circle. One-half cycle of the image sequence so constructed is shown in
Plate 3.2a. The EOF decomposition of this sequence (Plate 3.2b) identifies
OCR for page 46
46
two dominant patterns of spatial variance in the test sequence. EOF #1 is an
inclined plane that accounts for 67.7% of the total variance in the sequence.
EOF #2 is a circle that accounts for 32.2~o of the total variance. The
corresponding time series for EOF #1 and EOF #2 (Plate 3.2b) correctly
establish the phase relationship between the two patterns.
Additional Considerations. In remote-sensing applications, multiple ob-
servations are generally spread over different times. Moreover, data are gen-
erally not evenly distributed either in time (e.g., due to variations in orbit)
or in space (e.g., cloud cover may obscure some pixels in the scene). These
circumstances can bias the results of an EOF analysis because the EOF
analysis is predicated on evenly distributed data in both space and time.
If the NT images in the data set D (equation (3.1~) are not equispaced in
time, the resultant EOF decomposition may be biased toward specific time
periods. A practical way to minimize such temporal bias is to construct
a weighting scheme for the images in the data set. Generally, temporal
weights are constructed in such a way as to preserve the total mean and to-
tal variance of the data while simultaneously minimizing the temporal bias
(e.g., Kelly, 1985~. Spatial data (i.e., one or more of the NX pixels in a
given image) are often replaced by compositing data from other images that
surround the given image in time. Care must be taken to ensure that the
composite time scale is very much less than the time step between images in
the sequence undergoing EOF analysis. Other interpolation schemes (e.g.,
kriging) can also be used to minimize spatial and temporal bias in EOF
analyses resulting from imperfectly sampled data.
There are three ways to compute EOFs: (1) large covariance matrix
approach, (2) small covariance matrix approach, and (3) singular value de-
composition. The first two methods involve a direct solution of the ei~en-
-- r ~ o-
~ . · ~ . ~ · . · 4. . ~ ~ . . ~ rams
value equation of the covarlance matrix of the data set L'. Mere are two
ways to do this because the data in the image sequence can be stored in
the two-dimensional data array in two different ways. If the data are stored
such that there are NX rows and NT columns, the covariance matrix will
be NX x NX. If, however, the data are stored such that there are NT rows
and NX columns, then the covariance matrix will be NT x NT. For most
remote sensing applications, the number of spatial points NX in a given
image in the sequence is usually very much greater than the number of im-
ages, NT, in the sequence. Thus the method of data storage giving rise
to an NX x NX covariance matrix is often caned the large covariance ma-
OCR for page 47
47
trix approach and generally cannot be used in remote sensing applications.
Likewise, the method of storage resulting in an NT x JVT covariance matrix
is generally caped the small covariance matrix approach. It has been used
often in oceanographic and atmospheric applications (e.g., Preisendorfer,
1988~. Solution of the eigenvalue equation of the covariance matrix by sin-
gular value decomposition was discussed earlier in this chapter. A detailed
discussion of the various methods of solution of the eigenvalue equation of
the covariance matrix, in the context of EOF analysis, is given by GalIaudet
and Simpson (199lb).
3.2.2 Velocity Estimates from Image Sequences
General Comments
AD methods of motion estimation based on image sequence analysis depend
in some way on the detection of image brightness gradients. These gradients
are defined as normal in direction to contours of constant brightness (or sea
surface temperature [SST] for the case of the Advanced Very High Resolution
Radiometer [AVHRR] flying on the NOAA series of operation satellites).
The total velocity v at any point on a contour can be written as
v = zJn+rt,
(3.16)
where ~ and T are the magnitudes of the normal and tangential velocity
components, respectively, and n and t are the unit vectors in the directions
normal and tangential to the contour, respectively. The total velocity vector
can also be decomposed into ordinary Cartesian coordinates
v= u'+v',
(3.17)
where u and v are the magnitudes of the x and y velocity components,
respectively, and ~ and ~ are the unit vectors in the x and y directions,
respectively.
Below, objective methods for computing the total, normal,
and tangential components of near-surface oceanic flow from sequences of
AVHRR data are presented. The formulation follows closely that of Wahl
and Simpson (199Oa,b). Note that atmospheric motions can be computed
using these same methods from sequences of either AVHRR or Geostationary
Operational Environmental Satellite (GOES) data.
OCR for page 60
60
the more familiar two-dimensional satellite images are subsequently con-
structed from this telemetry stream using some set of mathematical trans-
formations.) These data generally are called "level 1" data and are mandated
for primary archives because higher-level data (e.g., the two-dimensional im-
ages) often are produced by irreversible transformations. Moreover, level 1
data usually consist of 10-bit (e.g., AVHRR, CZCS) or less (e.g., DMSP)
data strings packed into zero-filled 16-bit integer format for convenience.
Data compression techniques and optical disc storage technology clearly
are required if the data storage issue is to be adequately addressed. (Data
compression is the process of reducing the number of bits required to store
a given amount of information without loss of information.) For example,
tests with a recent compression algorithm using Lempel-Ziv-Welch (LZW)
coding (Welch, 1984) conservatively show that the average SO-Mb scene can
be compressed to 32 Mb. This represents a 60~o reduction in size from the
original data set. Some scenes may be compressed by as much as 75~o with
EZW coding.
The EZW compression algorithm has two main competitors currently
in common use: Huffman coding and run length coding (R£C). The Hu~-
man coding algorithm is not well suited to satellite data: preliminary tests
indicate that only a 14-20% reduction in size is achievable (Jain, 1989~. Fur-
thermore, it is much slower than the EZW algorithm. The REC algorithm
was originally designed to vectorize bit maps. It is designed to work on
strings of bits which are all Is or Os (Jain, 1989~. This makes it impractical
for satellite data, which tends to vary too much (i.e., the strings of uniform
Is or Os are too short). While the upper 6 bits of each 16-bit sample can be
coded efficiently with the REC algorithm, the remaining 10 bits pose a se-
rious problem for REC methods. Preliminary tests show that the estimated
size reduction obtained from R£C algorithms wiB be in the 10-25~o range.
The speed of the REC is comparable to that of EZW. These considerations
show that the EZW compression algorithm is best suited to the proposed
task.
Preliminary tests also indicate that decompressing a typical satellite pass
so that it can be used for analysis can take as much as 20 minutes per pass.
This time factor depends on the speed of the disc and processing unit. No
compression algorithm can get around this problem; increased access time
is the trade-off for decreased space usage.
Data compression techniques are also needed for the two-climensional
higher level satellite images. Unfortunately, most two-dimensional algo-
rithms either achieve speed by creating distortion in the data or achieve lack
OCR for page 61
61
of distortion by requiring excessively long execution times. Considerable re-
search is needed to develop efficient two-dimensional compression/decompres-
sion algorithms which do not distort the data.
3.3.2 Image Representation
The data structure selected to represent the spatial data in an image wiD
have a critical effect on the implementation and final performance of the
analysis algorithm. The quadtree and octree are hierarchical data struc-
tures often used to represent spatial data. The term quadiree is user! to
describe a class of hierarchical data structures whose common property is
that they are based on the principle of recursive decomposition of space
(Samet, 1989~. They can be differentiated on the following bases: (1) the
principle used to determine the decomposition process, (2) the resolution
(variable or constant), and (3) the type of data they are used to represent.
The prime motivation for the development of the quadtree is the need to
reduce the amount of space necessary to store data through the use of ag-
gregation of homogeneous blocks (Samet, 1989~. An important by-product
of this aggregation is the reduction in execution time of an analysis process.
Quadtrees have proved to be useful data structures for dithering algorithms,
computing geometric properties of images, implementation of linear image
transformations, development of hierarchical hidden-surface algorithms, and
ray tracing. The quadtree is only one of several digital data structures useful
in spatial statistics and digital image analysis.
The constraints on and the need for the large-scale use of remotely sensed
images in studies of global change is clear. Research is required in areas of
both data storage and image representation if optimal use of remotely sensed
data by the earth sciences community is to be achieved.
3.4 Special Considerations
Remote sensing of the environment with earth-observing satellites poses
some additional considerations beyond those normally encountered in labora-
tory-based applications of digital image analysis. In the laboratory, both il-
lumination and viewing geometry can be controlled. Moreover, the imaging
detector is close to the object being cletectecI, and any interfering influence
between imaging detector and object can be minimized. Finally, one gener-
ally has a good notion of what constitutes the detected object. In contrast,
earth observing satellites typically fly 800-900 km above the surface of the
OCR for page 62
62
earth. Viewing geometry and illumination are not controlled and can vary
greatly from orbit to orbit. The SOO-900 km layer of atmosphere between
target and detector acts as a filter that varies spatially and temporally, often
partially corrupting image quality. there is the need to accurately project
the image taken 900 km above the Earth onto a flat map projection suit-
able to the particular application under study. Thus atmospheric correction
algorithms (e.g., Curran and Dungan, 1989; Kaufman, 1988; Gratzki en cl
GerstI, 1989; Simpson and Humphrey, 1990), sensor calibration algorithms
(e.g., Gordon et al., 1983; Eckstein and Simpson, 199Oa), and earth location
algorithms (e.g., Brush, 1985; Goshtasby et al., 1986; Jezching and Jain,
1989) are all pre-processing steps essential prior to meaningful mathemati-
cal analysis.
. .
3.5 Summary
Remote sensing provides the only practical way to obtain the large-scale
synoptic data sets necessary to address major problems of global signifi-
cance in earth science that are fundamentally cross-disciplinary in nature.
The magnitude of the data set, however, may preclude meaningful under-
standing unless it can be distilled and synthesized into organized patterns
of variance that can be meaningfully related to the underlying governing
physics, chemistry, and biology of global-scale cycles and processes. Spa-
tial statistics and mathematical methods of digital image analysis provicle
mechanisms for such a synthesis. The examples cited herein included prin-
cinal component analyses, which are useful for image segmentation and for
determining spatial patterns of variance in large data sets; edge (letection;
pattern matching; optical flow methods, which are useful for determining
fields of motion from sequences of image data; and principal curves, which
are useful for determining the spatial distribution, size, and shape of ice floes
observed from spacecraft data. Atmospheric correction algorithms, sensor
calibration algorithms, and earth location algorithms generally are required
as pre-processes to digital image analysis of remotely-sensed images. Each
of these pre-processing areas contains challenging mathematical problems
which will have to be solved before earth sciences can benefit from the full
potential of remote-sensing technology.
O ~
OCR for page 63
63
Acknowledgments
This work was sponsored in part by the Marine Life Research Group of
the Scripps Institution of Oceanography and by a grant from the Office of
Naval Research. AVHRR data (Plates 3.1, 3.3, and 3.4) were provided by
the Scripps Satellite Oceanography Center. Landsat data (Plate 3.5) and
permission to use material from the Banfield and Raftery ice floe study were
generously provided by Professor Raftery of the University of Washington.
Dr. Barbara Eckstein produced the sample EOF analysis (Plate 3.2) while
working as an ONR-sponsored postdoctoral research assistant in the author's
research group at Scripps. Several other present or former members of this
group (I`. Al-Rawi, D. Atkinson, J. Bloom, T. Gallau(let, C. Humphrey,
J. Toman, and D. WahI) contributed in one way or another to the completion
of this chapter. S. McBride typed initial drafts of the manuscript, and
F. Crowe and G. Tupper assisted with figure preparation.
Bibliography
[1] Aggarwal, J. K., and N. Nandhakuymar, On the computation of motion
from sequences of images A review, Proceedings of ·EEE 76 (1988),
917-935.
[2] Anderberg, M. R., Cluster Analysis for Applications, Academic Press,
New York, 1973.
t3] Banfield, J. D., and A. E. Raft ery, Ice Floe Identification in Satellite
Images using Mathematical Morphology and Clustering about Principal
Curves, Technical Report No. 172, Department of Statistics, University
of Washington, Seattle, 1989.
OCR for page 64
64
[4] Barnett, T. P., and W. C. Patzert, Scales of temperature variability in
the tropical Pacific, J. Phys. Oceanogre 10 (1980), 529-540.
[5] Browning, K. A., Use of radar and satellite imagery for the mea-
surement and short-term prediction of rainfall in the United King-
dom, in Remote Sensing Applications in Meteorology and Climatology,
R. A. Vaughn, ea., NATO AST Series, vol. 201, 1986.
[6] Brush, R. J. H., A method for read time navigation of AVHRR imagery,
IEEE Trans. Geosci. and Remote Sensing 23 (1985), 876-~87.
[7] Campbell, J. D., Edge Structure and the Representation of Pictures,
Ph.D. dissertation, Department of Electrical Engineering, University of
Missouri, Columbia, 1969.
[~] Cheevasuvit, F., H. Maitre, and D. Vidal-Madjar, A robust method for
picture segmentation based on a split and merge procedure, Comput.
Vision, Graph. Image Proc. 34 (1986), 268-281.
[9] Curran, P. J., and J. A. Dungan, Estimation of signal-to-noise: A new
procedure applied toAVIRIS data, IEEE Trans. on Geosci. and Remote
Sensing 27 (1989), 620-627.
[10] Davis, R. E., Predictability of sea surface temperature and sea level
pressure anomalies over the North Pacific Ocean, '7. Phys. Oceanogr. 6
(1976), 249-266.
[11] Davis, R. E., Techniques for statistical analysis and prediction of geo-
physical fluid systems, Geophys. Astrophys. Fluid Dyn. 8 (1977), 245-
277.
[12] Deschamps, P. Y., M. Herman, and D. Tanre, Modeling of the atmo-
spheric effects and its application to remote sensing of ocean color,
Appl. Opt. 22 (1983), 3751-3778.
[13] Dudgeon, D. E., and R. M. Mersereau, Multidlimensional Digital Signal
Processing, Prentice HaD, Englewood Cliffs, Hi., 1984.
[14] Eckstein, B. A., and J. J. Simpson, Aerosol and RayTeigh radiance con-
tributions to Coastal Zone Color Zone Scanner images, Int. ]. Remote
Sensing, in press (199Oa).
[15] Eckstein, B. A., and J. J. Simpson, Cloud screening Coastal Zone Color
tone Scanner Images using channel 5, Int. ]. Remote Sensing, in press
(199Ob).
OCR for page 65
65
[16] EOS, Science and Mission Requirements Working Group Report, Vol-
ume ~ and Appendix, National Aeronautics and Space Administration'
Goddard Space Flight Center, Greenbelt, Md., 1984.
[17] Esbensen, K., and P. Geladi, Strategy of multivariate image analysis
(MIA), Chemometrics and Intelligent Laboratory Systems 7 (1989), 67-
86.
[18] Fukada, Y., Spatial clustering procedures for region analysis, Pattern
Recognition 12 (1980), 395-403.
[19] GalIaudet, T. C., and J. J. Simpson, Automated cloud screening of
AVHRR imagery using split-and-merge clustering, Remote Sensing En-
viron., submitted (1991a).
[20] Galiaudet, T. C., and J. J. Simpson, Seasonal and non-seasonal vari-
ability in sea surface temperature off Punta Eugenia, to be submitted
to Remote Sensing Environ. ~ 199lb).
[21] Geladi, P., H. Isaksson, L. Lin~quist, S. Vold, and K. Esbensen, Prin-
cipal component analysis of multivariate images, Chemometrics and
Intelligent Laboratory Systems 5 (1989), 209-220.
[22] Gordon, H. R., J. W. Brown, O. B. Brown, R. H. Evans, and
D. K. Clark, Nimbus 7 CZCS: Reduction of its radiometric sensitiv-
ity with time, Appl. Opt. 22 (1983), 3929-3931.
[23] Goshtasby, A., G. C. Stockman, and C. V. Page, A regional-based
approach to digital image registration with subpixe! accuracy, IEEE
Trans. Geosci. and Remote Sensing 24 (1986), 390-399.
t24] Gratzki, A., and S. A. W. GerstI, Sensitivity of an atmospheric correc-
tion algorithm for non-Lambertian vegetation surfaces to atmospheric
parameters, IEEE Trans. Geosci. and Remote Sensing 27 (1989), 326-
331.
[25] Hastie, T., and W. X. Stuetzle, Principal curves, A. Am. Stat. Assoc.
84 (1989), 502-516.
t26] Henderson-Sellers, A., Defogging cloud determination algorithms, Na-
t~re 298 (1982), 419-420.
[27] Hildreth, E. C., The Measure of Visual Motion, MTT Press, Cambridge,
Mass., 1983.
OCR for page 66
66
[28] Horn, B. K. P., and B. G. Schunck, Determining optical flow, Artif.
Intell. 17 (1981), 185-203.
[29] Ingerbritsen, S. E., and J. P. Lyon, Principal component analysis of
multi-temporal image pairs, Int. J. Remote Sensing 6 (1985), 687-696.
[30] Jain, A. K., Fundamentals of Digital Image Processing, Prentice Had,
Englewood Cliffs, N.~., 1989.
[31] Jain, A. K., and R. C. Dubes, Algorithms for Clustering Data, Prentice
Hall, En~ewood Cliffs, N.~., 1988.
[32] Jain, J. R., and A. K. Jain, Displacement measurement andits appli-
cation to interframe image coding, IEEE Trans. Common. 29 (1981),
1789-1808.
[33] Jensen, J. R., Introductory Digital Image Processing, Prentice Hall,
Englewood Cliffs, N.J., 1986.
[34] Jezching, T., and A. K. Jain, Registering I'andsat images by point
matching, IEEE Trans. Geosci. and Remote Sensing 27 (1989), 642-
651.
[35] Kaufman, Y. J., Atmospheric effects on spectral signature
Measurements and corrections, IEEE Trans. Geosci. and Remote Sens-
ing 26 (1988), 441-449.
[36] Kelly, K. A., The influence of winds and topography on surface tem-
perature patterns over the northern California slope, if. Geophys. Res.
90 (1985), 11,783-11,793.
[37] Lagerioef, G. S. E., EOF analysis of AVHRR and CZCS data, Third
Conference on Satellite Meteorology and Oceanography, 1986.
[38] I.illesand, T. M., and R. W. Kiefer, Remote Sensing and Image Inter-
pretation, 2nd ea., John Wiley and Sons, New York, 1987.
[39] lLuenberger, D. G. Optimization by Vector Space Methods, John Wiley
and Sons, New York, 1969.
Marr, D., and E. C. Hildreth, Theory of edge detection, Proc. R. Soc.
London, B 207 (1980), 187-217.
[41] Marr, D., and S. Uliman, Directional sensitivity and its use in early
visual primitives, Proc. R. Soc. London, B 211 (1981), 151-180.
OCR for page 67
67
[42] Mather, P. M., Computer Processing of Remotely-Sensed Images: An
Introduction, John Wiley and Sons, New York, 1987.
[43] Matveev, L. T., Computer Processing of Remotely Sensed Images: An
Introduction, John Wiley and Sons, New York, 1984.
[44] MuerIe, J. L., and D. C., Allen, Experimental evaluation of techniques
for automatic segmentation of objects in a complex scene, pp. 3-13 in
Pictorial Pattern Recognition, G. Cheng et al., eds., Thompson, Wash-
ington, D.C., 1968.
[45] Ohring, G., and P. F. CIapp, The effects of changes in cloud amount
on the net radiation at the top of the atmosphere, A. Atmos. Sci. 37
(1980), 447-454.
[46] Pailleux, F. S., The impact of satellite data on global numerical weather
prediction, in Remote Sensing Applications in Meteorology and Clima-
tology, R. A. Vaughn, ea., NATO AST Series, vol. 201, 1986.
[47] Barman, D., and J. Kittler, Clustering algorithms for use with images
of clouds, Int. ~7. Remote Sensing 7 (1986), 855-866.
[48] Pavlidis, T., Structural Pattern Recognition, Springer-VerIag, Berlin,
1977.
[49] Preisendorfer, R. W., Principal Component Analysis in Meteorology
and Oceanography, C. D. Mobley, compiler ea., Developments in Atmo-
spheric Science 17, Elsevier, Amsterdam, 1988.
[50] Press, W. H., B. P. Flannery, S. A. Teukolsky, and W. T. Vetterling,
Numerical Recipes, Cambridge University Press, New York, 1986.
[51] Ramanathan, V., The role of earth radiation budget studies in climate
and global circulation research, '7. Geophys. Res. 92 (1987), 4075-4095.
[52] Richards, J. A., Remote Sensing Digital Image Analysis: An Introduc-
tion, Springer-VerIag, New York, 1986.
[53] ~ --
gamut, tit' Applications of Spatial Data Structures, Computer Graphics,
Image Processing and GIS, Addison-Wesley Publishing Corporation,
New York, 1989.
t54] Seddon, A. M., and G. E. Hunt, Segmentation of clouds using cluster
analysis, Int. '7. Remote Sensing 6 (1985), 717-731.
OCR for page 68
68
[55] Serra, J., Image Analysis and Mathematical Morphology, Academic
Press, New York, 1982.
[56] Simpson, J. J., Image segmentation using recursive polygon fill opera-
tions, Remote Sensing Environ., submitted (1991~.
[57] Simpson, J. J., and J. Bloom, Objective estimates of velocity from
sequences of satellite data using log search and minimum distortion
methods, if. Geophys. Res., submitted (1990~.
[58] Simpson, J. J., and C. Humphrey, An automated cloud screening al-
gorithm for daytime AVHRR imagery, J. Geophys. Res. 95 (1990),
13,459-13,481.
[59] Singh, A., and A. Harrison, Standardized Principal Components, Int.
]. Remote Sensing 6 (1985), 883-896.
[60] Strub, P. T., C. James, A. C. Thomas, and M. A. Abbott, Seasonal and
non-seasonal variability of sateDite-derived surface pigment concentra-
tion in the California Current, '7. Geophys. Res. 95 (1990), 11,501-
11,530.
[61] Verri, A., and T. Poggio, Motion field and optical flow: Qualitative
properties, IEEE Trans. Pattern Anal. Machine Intell. 11 (1989), 490-
498.
[62] WahI, D. D., and J. J. Simpson, Physical processes affecting the objec-
tive determination of near-surface velocity from satellite data, '7. Geo-
phys. Res. 95 (199Oa), 13,511-13,528.
[63] WahI, D. D., and J. J. Simpson, SateDite-derived estimates of the nor-
mal and tangential components of near-surface flow, Int. ]. Remote
Sensing, in press (199Ob).
[64] Welch, Terry A., A technique for high performance data compression,
IEEE Computer 17 (June 1984), 8-19.
OCR for page 69
69
Appendix to Chapter 3
Listed below are definitions and explanations for abbreviations commonly
used in remote sensing.
AVHRR The Advanced Very High Resolution Radiometer on the Na-
tional Oceanic and Atmospheric Administration's polar-orbiting weather
satellite. It measures cloud cover and infrared sea surface temperature.
DFT—Discrete Fourier transform.
EOF Empirical orthogonal function.
EOS Earth Observing System. A proposed National Aeronautics and
Space Administration program for earth observing systems to be launched
between 1997 and 2007.
FFT—Fast Fourier transform.
GOES Geostationary Operational Environmental Satellite.
tional weather satellite used to measure cloud cover.
solar radiation often are computed from GOES data.
LOG Laplacian of the Gaussian operator.
An opera-
Estimates of
EZW Lempel-Ziv-Welch coding used in data compression algorithms.
MN Minimum norm solution.
MU—The Marr-UlIman solution for the normal component of velocity.
OF Optical flow method of computation.
REC Run length coding used in data compression algorithms.
SST Sea surface temperature.
OCR for page 70
Representative terms from entire chapter:
eof analysis
