
Citation 
 Permanent Link:
 http://ufdc.ufl.edu/AA00013547/00001
Material Information
 Title:
 System modeling using generalized linear feedforward networks
 Creator:
 Jones, Douglas G., 1963
 Publication Date:
 1999
 Language:
 English
 Physical Description:
 vi, 140 leaves : ill. ; 29 cm.
Subjects
 Subjects / Keywords:
 Analytics ( jstor )
Approximation ( jstor ) Error rates ( jstor ) Mathematical vectors ( jstor ) Modeling ( jstor ) Parametric models ( jstor ) Signal processing ( jstor ) Signals ( jstor ) Transfer functions ( jstor ) Weapons ( jstor ) Dissertations, Academic  Electrical and Computer Engineering  UF ( lcsh ) Electrical and Computer Engineering thesis, Ph.D ( lcsh ) Linear control systems ( lcsh ) Mathematical optimization ( lcsh ) Network analysis (Planning) ( lcsh )
 Genre:
 bibliography ( marcgt )
nonfiction ( marcgt )
Notes
 Thesis:
 Thesis (Ph.D.)University of Florida, 1999.
 Bibliography:
 Includes bibliographical references (leaves 132139).
 General Note:
 Typescript.
 General Note:
 Vita.
 Statement of Responsibility:
 by Douglas G. Jones.
Record Information
 Source Institution:
 University of Florida
 Holding Location:
 University of Florida
 Rights Management:
 Copyright [name of dissertation author]. Permission granted to the University of Florida to digitize, archive and distribute this item for nonprofit research and educational purposes. Any reuse of this item in excess of fair use or other copyright exemptions requires permission of the copyright holder.
 Resource Identifier:
 030475416 ( ALEPH )
43164426 ( OCLC )

Downloads 
This item has the following downloads:

Full Text 
SYSTEM MODELING USING
GENERALIZED LINEAR FEEDFORWARD NETWORKS
By
DOUGLAS G. JONES
A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY
UNIVERSITY OF FLORIDA
1999
ACKNOWLEDGE MENTS
There are many people who made the completion of this work possible. First and
foremost I would like to thank my research professor, Dr. Jose C. Principe. He has
taught me the valuable skills of rigor and completeness. Also, his patience with my
progress while I juggled both research and a busy career was greatly appreciated. I
would like to thank my supervisory committee members Dr. Yunmei Chen, Dr.
William W. Edmonson, Dr. John G. Harris, and Dr. Jian Li for their encouragement
and direction.
Many thanks go to the management within the 46 Test Wing, Eglin AFB, FL, who
made it possible for me to attend the University of Florida through the longterm full
time training program. This organization fosters an atmosphere that encourages
academic achievement, and it has continually provided me with the opportunity to
further my education.
There are no words to describe the wisdom, guidance, and support I have received
from my parents. Their love to this day is endless.
Above all I would like to give the highest recognition to my loving wife, Donna,
for her sacrificing devotion throughout my education. To my precious children,
Lindsey and Mark, you are the sunshine in my life.
Lastly, I dedicate this work to my Lord and Savior Jesus Christ. I can do all
things through him who gives me strength.
TABLE OF CONTENTS
page
ACKNOW LEDGM ENTS ......................................................................................... .. ii
ABSTRACT .................................................................................................................. v
CHAPTERS
1 PROBLEM FORM ULATION ................................................................................. 1
1.1 Problem D efinition......................................................................................... 2
1.2 M easurement Objective ......... ............... ........................... ..... ................... 3
1.3 Model Classifications................................................................................... 4
1.4 H historical R eview ............................................................................................... 5
1.5 Modeling Methodology................................................................................ 7
1.6 Research Contribution.................................................................................. 8
2 SYSTEM ASSUM PTIONS...................................................................... .... ........ 10
2.1 Linear, TimeInvariant, Causal Systems............................................................ 10
2.2 R ealizability ................................................................................................. 14
2.3 L2 Signals and Systems ...................................................................................... 14
2.4 System Synthesis.......................................................................................... 15
2.5 Continuous/DiscreteTime Transformation................................... ............ 16
2.6 E xam ples...................... ................................. .............................................. 18
2.7 Conclusion (System Assumptions)............................. ............ ............... 21
3 MODEL FORMULATION....................................................................................... 22
3.1 General Model Form M(p)......................................................................... 22
3.2 D iscreteTim e M odels................................................................................... 23
3.2.1 FIR M odel............................................................................................... .. 23
3.2.2 R M odel............................................................................................... 24
3.2.3 G LFF M odel ........................................................................................... 24
3.2.4 FIR, IIR, GLFF Model Block Diagrams ................................................. 25
3.3 GLFF M odel Kernels ............ .... ................... ................................................ 27
3.3.1 Popular DiscreteTime GLFF Kernels.................................................... 27
3.3.2 Popular ContinuousTime GLFF Kernels............................................... 27
3.4 GLFF Kernel Derivations and Properties..................................... ............ 30
3.4.1 ContinuousTime Laguerre Kernel.............................................................. 30
3.4.2 ContinuousTime Kautz Kernel (Complex Time Scale Vector) ...............30
3.4.3 ContinuousTime Kautz Kernel (Real Time Scale Vector).........................31
3.4.4 ContinuousTime Legendre and Jacobi Kernels..................................... 32
3.5 Impulse and Step Response Modeling.......................................................... 32
iii
3.6 Conclusion (Model Formulation)............................................................... 33
4 OPTIMAL PARAMETER SELECTION............................................................ 35
4.1 M odel D im ension.......................................................................................... 35
4.1.1 C om pleteness ............................................................................... .......... 36
4.1.2 O rthogonal Projection............................................................................... 37
4.1.3 "O ptim al" Selection of N ............................................................................ 38
4.2 W eight V sector ............................................................................................ ...... 39
4.2.1 Optimal Selection of w ........................................................................ 39
4.2.2 Relationship to WienerHopf Equations...................................................41
4.2.3 Minimum Squared Error Performance Function ......................... .............. 42
4.2.4 Computing w Using Step Response Data..................................................43
4.2.5 Adaptation of w .......... ............................................................................46
4.3 T im e Scale V ector............................................................................................ 46
4.3.1 Optimal Selection of A (Historical Review).................................... ..48
4.3.2 Optimal Selection of A (Complex Error Analysis).................................... 56
4.3.3 Further Analysis of the Complex Performance Function ....................... 66
4.3.4 Time Domain Synthesis Using GLFF Networks....................................... 76
4.3.5 Conclusion of the Complex Performance Function.................................. 78
4.4 Conclusion (Optimal Parameter Selection)................... ........... .......... .... 80
5 A PPLIC ATIO N S ........................................... ......................................................... 81
5.1 Examples of GLFF Networks Using Simulated Systems ................................83
5.2 Laboratory Testing of Guided Weapon Systems ............................................. 89
5.2.1 O objectives .............................................................................................. 91
5.2.2 Functional Breakdown of a Missile System ........................................... 92
5.2.3 HITL ClosedLoop Missile Simulations....................................................94
5.3 Modeling the Step Response of a Carco 5Axis FMS ....................................95
5.3.1 DiscreteTime GLFF Models of the FMS .................................................98
5.3.2 Step Response Modeling Procedure ........................................................ 100
5.3.3 Results (FMS Step Response Modeling)................................................ 103
5.3.4 ContinuousTime GLFF Models of the FMS .......................................... 118
5.3.5 Conclusion (Step Reponse Modeling)..................................................... 123
5.4 O their A applications ........................................................................................... 123
5.5 Conclusion (Applications) ............................................................................... 124
APPENDICES
A: ORTHOGONAL SIGNAL SETS ......................................................................... 125
B: HILBERT SPACE CONCEPTS ........ ................................................................. 126
C: SDOMAIN DERIVATION OF MINIMUM J................................................... 128
D: SQUARED ERROR STATIONARITY USING THE LAGUERRE MODEL .... 130
R EFER E N C E S ........................................................................................................ 132
BIOGRAPHICAL SKETCH............................................................................... 140
Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy
SYSTEM MODELING USING
GENERALIZED LINEAR FEEDFORWARD NETWORKS
By
Douglas G. Jones
August 1999
Chairman: Jose C. Principe
Major Department: Electrical and Computer Engineering
The focus of this work is the derivation and development of a class of models
called Generalized Linear FeedForward (GLFF) net%%orks. Included is a description
of the a priori assumptions concerning the systems under study, a development and
formulation of the models (networks) to be used. and a definition of a measurement
criteria used to judge goodnessoffit. The models are constructed from a weighted
sum of parametrically dependent kernels. These kernels represent basis vectors that
span a vector space with dimension equal to the order of the network realization.
Optimization of the model parameters is discussed and a new approach is given; one
that extends the conventional squared error analysis to the complex domain.
After presentation of GLFF model theory and parametric optimization
methodologies, network realizations are constructed to address the system modeling
problem. GLFF networks can be formulated to model discretetime, continuous
time, and sampleddata systems. In this study, it is shown how a flight motion
simulator, used in the laboratory test and evaluation of guided weapon systems, can
be accurately modeled using GLFF networks comprised of Gamma, Laguerre, and
Legendre kernels.
Advantages of GLFF networks are their ability to accurately represent linear time
invariant systems coupled with the simplicity of the procedures employed to optimally
determine their parameters. The network structures have efficient and practical
implementations and applications are offered as proof.
CHAPTER 1
PROBLEM FORMULATION
In this chapter an overview of the entire content of the research is gi\en.
beginning with a definition of the problem to be solved; namely, modeling the
impulse and step responses of causal, linear, timeinvariant systems. The
system/model configuration is illustrated and the primary measurement objective
stated. The model architecture is briefly described and prefaced with a review of the
major contributors leading to the proposed formulation. Finally, the contribution of
this research is summarized. It will be developed throughout the remaining chapters.
It will be shown that the models developed arc capable of approximating systems
with real and/or complex conjugate pair poles, even those \which have slowly
decaying impulse responses (large time constants) and/or oscillatory modes (light
damping). This is difficult to achieve with finite impulse response (FIR) structures.
The models studied are considered parametric since they are functions of parameter
vectors that must be determined. These parameter vectors are fairly easy to optimize
and unlike infinite impulse response (IIR) structures stability can always be
guaranteed.
The networks under consideration in this research are comprised of linear
combinations of special basis filctions. These basis functions are not arbitrarily
chosen but are derived from assumptions of the systems being modeled. The
2
networks therefore become good approximates since they obtain the natural
information of the unknown system. This information resides in a parameter vector
built into the model structure. The parameters are optimally chosen so that the
models approximate the system in a least squares sense. Hence, the proposed theory
is an implementation of parametric least squares error modeling [Cad90].
The modeling methodology can also be described using geometric concepts.
Considering signals as vectors in vector spaces, the models are linear combinations of
parametric basis vectors. They span a subspace of the space spanned by the unknown
system. The model output is a projection of the system vector onto the model
subspace. This space has an orientation that is parametrically established. The error
is a vector that lies outside the model subspace. Projecting orthogonally the system
onto the model space minimizes the error.
1.1 Problem Definition
The objective is to approximate a causal, linear timeinvariant (LTI) system, H,
using a predefined model, M(p), such that a measure of the quality of the
approximation is minimized. Models developed in this study will be called
generalized linear feedforward (GLFF) networks and will be defined in Chapter 3.
The basic system modeling configuration is shown in Figure 11 below. Arguments
are not included with the functions {x, d, y, e, H} to illustrate that the formulation will
apply in both a continuoustime setting {x(t), d(t), y(t), e(t), H(s)} and a sampleddata
or discretetime setting {x(k), d(k), y(k), e(k), H(z)}. Quite often in this thesis we
switch backandforth between continuoustime and discretetime concepts. This is
intentionally done to prove the broad applicability of this work. Most of the theory
however will be derived in continuoustime. It would have been just as appropriate to
develop the concepts in discretetime.
d
H
x e
M(p)
Figure 11: System modeling configuration.
1.2 Measurement Objecti\e
How is it determined whether a model has achieved a good representation of the
system being approximated? For this a measurement of the error or a test for
goodnessoffit of the model is needed. In this thesis the integrated squared error
criteria will be employed; namely, for a given model M(p), the parameter set p is
found such that the energy in the error signal e is minimized. This is one of the
simplest error criteria to work with since the measurement is basically quadratic in
form and its derivatives are linear in certain elements of the parameter set. So the
objective is to minimize
J = e(t)dt= (d(t) y(t)) dt (continuoustime)
or (1.1)
J = e2(k) = J(d(k) y(k))2 (discretetime)
k=O k=O
1.3 Model Classifications
To provide motivation for the GLFF networks (a complete formulation will be
given in Chapter 3), consider the discretetime setting. In the past discretetime
system models have been largely divided into two classes: FIR and IIR. A new class
of models (filters) known as the generalized feedforward (GFF) filters was
popularized with the development of the Gamma filter [Pri93]. This filter has a
structure similar to the tapdelay line (FIR filter) with the delays replaced by a "more
general" transfer function known as the Gamma kernel. Extending this development
can be accomplished by considering other popular kernels, most of which are
orthogonal, such as the Laguerre, Legendre, Jacobi, and Kautz [Lee67], [Arm59],
[Kau54]. The extended formulation with a more comprehensive set of kernels is
called the GLFF model class because its members are used in the context of modeling
linear systems and they obey the superposition principle. How does the GLFF model
class compare to FIR and IIR models? Figure 12 relates these three classes based
upon important issues that are usually considered when deciding which modeling tool
to pull from the signal processing toolbox. GLFF models can be said to fall
somewhere between FIR and IIR models in capability, complexity, stability and
computability.
FIR GLFF IIR
low complexity high
low * { } * high
l capability J
f stability 1
trivial < stability difficult
FcomputabilityI
Figure 12: FIR/GLFF/IIR comparison.
This thesis exploits the strengths of GLFF models. They offer an alternative to
the two extremes in the modeling process. Their advantages are that they have a
reduced complexity compared to IIR models while offering greater capability over
FIR models. Unlike IIR modeling techniques, stability can be easily guaranteed and
their optimal or nearoptimal solutions can be computed with much less difficulty.
The most difficult task is in selection of the time scale r sector A. This is one element
of the parameter set p. The complete parameter set will be defined in Chapter 3.
Optimization of its elements comprises all of Chapter 4.
1.4 Historical Review
Much of the theory on which the GLFF model class is founded began with the
work of Lee [Lee33], [Lee67]. Lee showed how to use orthogonal functions to
synthesize transient networks and how this concept is employed to build electrical
devices that can realize physical systems. Another major contributor is Kautz
[Kau54]. He provided a convenient way to develop the time domain response of a
system using a most general orthogonal signal set. Young and Huggins extended
Kautz' procedure to discretetime [You62a] as did Broome [Bro65]. Nurges and
Yaaksoo [Nur81] formulated state space equations that represent discretetime linear
systems using the Laguerre model. Roberts and Mullis provide a generalized state
variable formulation of orthogonal allpass network structures [Rob87]. These
structures keep data paths local. The fundamental processing elements are alike
except perhaps for local scaling parameters. This simplifies the design and allows for
distributing the computations in a connected array. These orthogonal structures also
have the advantage that higher order models can be constructed by simply adding new
blocks or terms without effecting or forcing a recalculation of the existing lower order
terms. The concept of a generalized feedforward filter was proposed by Principe
[Pri93] with the Gamma model. This idea of creating a more general transversal filter
using the earlier work with orthogonal signal sets leads to the modeling methodology
here described; namely, the generalized linear feedforward networks.
What is the basic objective? Given a desired impulse response, h(t), of a causal
LTI system along with a prescribed tolerance on the allowable error, find an
approximating impulse response, h(t), that is simple in structure and easy to compute.
Let the system and model transfer functions be H(s) and H(s), respectively.
Construct H(s) with the form
m
(s) (S ,) ( A,
H(s) A '= (1.2)
x(s) (s ) = ( p)
i=1
Assuming no coincident poles the impulse response is
n
h(t) = '(H(s)) = A eP (1.3)
i1
where Re(pi) < 0, V i=1, ..., n. It is desired that h(t) and h(t) be as close as possible
under prescribed tolerances.
There are two wellknown approaches to performing this approximation. The first
is to transform h(t) into H(s) and then approximate H(s) in the sdomain by using a
ratio of polynomials. This approach is generally less difficult and several techniques
for accomplishing this approximation have been presented in the literature. These
include Prony's method, Pad6 approximants, continuedfraction expansions, and even
analog filter design [Su71], [Bak81], [Opp89]. The problem with this approach is
that the error in the time domain is difficult to control. The second approach is to
approximate h(t) with a finite sum of exponentially damped sinusoids and damped
exponentials (as described above) and then transform the approximation to find H(s).
Here the modeling error may be controlled directly in the time domain. One way of
accomplishing this is to implement the procedure described by Kautz [Kau54] which
involves approximating h(t) using a weighted sum of special basis functions. These
basis functions are themselves weighted sums of exponentials carefully constructed so
that they are orthonormal. This provides for great simplification when searching for
the optimal parameter set p and gives direct control over the associated error.
1.5 Modeling Methodology
In this thesis the method used to implement GLFF network approximations of
linear systems is based on the work done by Kautz. This is carried out by considering
a construction of h(t) using an infinite sum of weighted time domain basis functions
{(,,,(t)} 0". Since the infinite sum must be truncated to a finite number of terms an
approximation results
N
h(t) =h (t)= w,4,,(tA ) (1.4)
n=0
The weights w,,, w, ... ,WN and basis functions {0,(t,)}I are chosen with the
following goals in mind:
1. For a fixed number of terms, N, make the convergence as rapid as possible.
This is controlled using a time scale vector A.
2. Make the optimal weights independent of the number of taps, N+1. This is
accomplished by using orthogonal basis functions.
3. For a fixed number of taps, N+1, minimize the energy in the error signal. This
requires determining optimal values for all of the elements of the parameter
set p in such a way that J (the integrated squared error) is minimized.
The above procedure is applicable when modeling the impulse response of a
desired system. It may also be of interest to approximate the step response. In fact
the applications given in this thesis involve finding system models using step
response data. If the network input is a step function, this method is employed by
approximating not the desired system output but its derivative or by approximating
the desired output and taking the derivative of the result. This will give an analytic
expression for the impulse response. When only a numerical or tabulated form of the
impulse response is desired, use the step response data to generate a system model.
Once the model is established, use an impulse (or white noise) input to generate the
impulse response data.
1.6 Research Contribution
This study extends the above concepts and contributes to the system modeling
theory in several ways.
First, a unified system modeling framework centered around the GLFF model
class is provided. As will be seen in Chapter 3 this formulation not only includes
models that employ orthogonal signal sets like the Laguerre and Kautz but non
orthogonal ones as well, such as the Gamma. In addition to these generalized models,
the basic framework can also be used to represent FIR and IIR structures.
Optimal parameter selection is a big part of the modeling procedure. The most
difficult task is locating or approximating the optimal time scale vector. A new
approach that extends traditional squared error analysis to the complex domain is
developed. It will be shown how the information in this new space can be used to
locate the optimal time scale parameters.
Finally an application of the theory is given by modeling the step response of a
Carco 5axis flight motion simulator (FMS). This system is used in hardwareinthe
loop (HITL) testing of guided weapon systems. GLFF modeling of other systems and
subsystems used in HITL testing is also considered. It is demonstrated how this class
of models can be employed to accurately represent the dynamics of physical systems
and add value to the task of weapon system test and evaluation.
CHAPTER 2
SYSTEM ASSUMPTIONS
In the previous chapter, a brief discussion was given concerning the assumed
systems to be modeled. A more complete analysis of all the system assumptions is
now offered. Mathematical and physical arguments are used to validate their
necessity. The basic assumptions are that input/output signals have finite energy and
systems are stable, causal, linear, timeinvariant, and realizable. After stating the
assumptions, the technique employed to synthesis these systems by simulation is
given. Higher order systems will be simulated by cascading together 1st and 2nd order
subsystems. This is done in a continuoustime framework. The transformation used
when a digital system is desired is also discussed. Finally some simulation examples
of the kinds of real systems to be modeled are offered.
2.1 Linear, TimeInvariant, Causal Systems
There are many ways to classify systems; linear vs. nonlinear, timeinvariant vs.
time varying, causal vs. noncausal, stochastic vs. deterministic, continuous vs.
discrete, etc. One important step in the study and analysis of systems theory is to
derive a structural or mathematical description of processes and systems. This
assumed structure is essential if plausible analytical models are sought to approximate
them. One could argue that if we had a mathematical description of the process why
would we need additional models. The primary reason for the model is to
approximate input/output behavior of the process. The assumed mathematical
description of the process could be used to give very specific information about its
internal composition (time constants, damping ratios, etc.); however, a component
level specification may be unavailable. In fact even if these system parameters are
known they may not be operating within their stated specification. In this research an
understanding of the internal composition is not the objective, only a model of the
input/output behavior when viewing the system as a black box. The theory presented
in this section does however give the foundation used to formulate and validate this
new class of models and serves to illustrate the derivation of these concepts.
Many systems studied in engineering can be well represented by mathematical
relations that describe their physical makeup or behavior. Quite often systems are
very complex and a realistic characterization may require nonlinear and/or time
varying equations that are difficult to solve. For simplicity, practicality, and in order
to establish an applicable set of tools which can be used in the design and analysis of
systems, approximations and assumptions are made about their physical nature. This
is done primarily to allow for the use of linear systems theory. This reasoning is often
justified in one of two ways:
1. The system is assumed to be linear by nature or it at least operates in a linear
region of its capabilities so that the linearity conditions are satisfied.
2. The system is assumed to be nonlinear by nature and in order to apply linear
systems theory we conceptualize a linearization around a nominal operating
point. With this approach precautions must be taken to ensure the domain of
the linear assumption is not exceeded.
Why is the use of linear systems theory desirable? In addition to the availability
of a wellstudied set of tools, most often a description of a system via the impulse
response and transfer function is desired. These functional system descriptions are
established by employing the theory of linear systems [DeC89].
Consider a linear timeinvariant (LTI) system having input x(t) and output d(t).
The impulse response, h(t), is defined as the output when the input is a unit impulse
function 8(t). The transfer function, H(s), is defined as the Laplace transform of h(t)
with zero initial conditions. This is equivalent to assuming h(t) = 0 for t < 0 (causal).
The above definition of the transfer function is in terms of the impulse response.
However a description of the transfer function taking a different, more mathematical
approach is also available. In practice, a large class of signals and systems in
engineering can be formulated through the use of differential equations. The LTI,
causal assumptions equivalently translate into the requirement that systems be
representable by linear, constant coefficient, ordinary differential equations. Namely,
an /th order system can be written as
d (t) + an_d "l)(t) + + ad(' (t) + a0d(t) =
(2.1)
b,,x(m(t)+ b,_lx('(t)+. .+ bx("(t) + box(t)
where ao, a,, ..., a,_, bo, b\, ..., b,, are all real constants with n 2 m. Systems
described by ordinary differential equations of this form are said to contain lumped
parameters rather than distributed parameters [Kuo87]. Partial differential equations
such as the elliptic equations, hyperbolic equations, and parabolic equations describe
phenomena such as fluid flow, wave propagation, and heat conduction all of which
are considered to be distributed parameter systems [Pow79].
To obtain the transfer function of an LTI system described by the differential
equation above assume zero initial conditions (causality) and take the Laplace
transform
H(s) = D(s) b,,s" + b,, s' +b s+ (2b2)
X(s) s" + a,_,s"'+ +ags + ao
Re%\ writing this in factored form gives
(71sz)
H(s)= A ' (2.3)
S(s p,)
i=1
and in partial fraction form (assuming no repeated poles)
H(s) = I (2.4)
I=, (S p,)
To obtain the impulse response, take the inverse Laplace transform. Assuming no
coincident poles this yields
h(t)= AI'l"t (2.5)
So the LTI, causal assumption leads to the following result: LTI causal systems
can be represented by constant coefficient ordinary differential equations whose
solutions (impulse responses) are comprised of a sum of exponential functions that
are real (when pi is real) and/or complex (when Pi is complex). If some of the poles
are repeated then there are terms with multiplicative factors t, t2, t3, ..., tk1 for a kth
order pole. This representation is the motivation for the modeling structure
considered in this thesis.
2.2 Realizability
It was shown that the LTI, causal system assumption leads to a transfer function
representation of the form
H(s)= A '=' (2.6)
H(s p,)
i=1
An answer is now given to the following question "What are the constraints on
this form that ensures it is a description of a physically realizable system?" The
realizability requirement of H(s) is satisfied with the following constraints:
1. The constant A must be a real.
2. The number of poles, n, must be greater than the number of zeros, m.
3. Complex poles and zeros must occur in conjugate form. That is if pi is
complex then there must be a pj, for some j, such that pj = pi where denotes
complex conjugation.
When H(s) is realizable, it is easy to show that the impulse response is represented
by a weighted sum of exponentials and exponentially decaying sinusoids.
2.3 L2 Signals and Systems
The final constraint imposed is that the signals must belong to the vector space
L2(R), the squareintegrable functions. Here integrable is meant in the Lebesgue
sense [Deb90]. Assume the signal xe L2(R) and x: R R. then
x2 (t)dt < (2.7)
Hence it is said that x has finite energy. Let h(t) be the impulse response of a system.
Since h(t) = 0 for t < 0, and he L2((R) then
Jh'(t)dt
0
In this case it is said that h is a stable system. So by operating in L (R) a restriction is
imposed to work only with finite energy signals and stable systems.
2.4 System Synthesis
Recall the system realizability requirement results in an impulse response
representation that is a weighted sum of exponentials and exponentially decaying
sinusoids. From the transfer function perspective then a system can be thought of as
consisting of a cascade of 1st and 2nd order subsystems. The most general forms of
these subsystem structures are
D, ( s) b
Hi (s) = D =(s) b__
X(s) s+a (2.9)
D, (s) b~s + bo
Hz(s) = =
X(s) s +a,,s +a2o
resulting in differential equations
d (t)+ ald,1() = box(t)
(2.10)
i, (t) + a,,d(t) + ad(t) = b,,x(t) + bzx(t)
During this thesis, examples will be given whereby simulated systems are created.
For such examples, the synthesized systems will consist of these 1st and 2nd order
subsystem blocks. Most often slightly less general forms will be used to keep the
synthesis more in line with engineering (rather than mathematical) concepts. These
special forms convey the engineering information of usual interest: lime constants,
damping ratios, and undamped natural frequencies. These subsystems have the
following transfer functions
H,(s)D, (s) 1 =_,
X(s) rs+l s+Y
2 (2.11)
H, (s) co,
X(s) s2 + 2 ws + Wo
and corresponding differential equations
rd (t)+ d,(t)= x(t)
(2.12)
d,(t) + 2 (),d2(t) + A d(t)= =ox(t)
where
r = time constant
S= damping ratio
t, = undamped natural frequency
2.5 Continuous/DiscreteTime Transformation
The goal of this thesis is to define a model class useful for approximating discrete
time, continuoustime, and sampleddata systems. The above system formulation was
carried out in the continuoustime domain. The same principles however could have
been derived in the discretetime domain. All the concepts discussed have equivalent
discretetime representations. Sometimes a continuoustime system in analytic form
is given when a discretetime representation is desired. This can be resolved by
sampling. In this case, the discretetime system is called a sampleddata system.
Hence a sampling time, T, will necessarily be specified.
Having a continuoustime transfer function, there are several methods used to map
to the discretetime domain. The most popular are the bilinear, impulse or step
invariant, and matchedZ transform [Str88]. There are tradeoffs among these
mappings and some work well where others do not. The most widely used is the
bilinear transform. The transformation (or mapping) is carried out as follows
H(z)=H(SIs 2z1 (2.13)
T z+l
where T is the sampling interval. This is also known as Tustin's method in control
systems circles and is derived by considering the problem of numerical integration of
differential equations [Fra90]. A derivation of the bilinear transform is as follows.
Beginning with the transfer function H(s), determine its differential equation
representation. A difference equation is then derived whose solution matches the
solution of the differential equation at the sample points. The specific transformation
depends on the method with which the differential equation solution is computed.
Consider the 1s order system
D,(s) _
H,(s) = () (2.14)
X (s) s +
or associated differential equation
d,(t)+ d,(t)= x(t) (2.15)
The solution, d\(t), can be written in integral form as
I
d,(t)= J[d,(u)+x(v)]du (2.16)
0
For the sampling interval of time Tthis can be rewritten as
nTT nT
d,(n)= d,(nT)= j[d,(v)+x(u)]dv+ j[d,(v)+x(v)]dv
0 Trr (2.17)
nT
= d(n1)+r I[d,(v)+x(u)]dv
nTT
The method in which the integral is computed over the region [nTT, nT]
determines the discretetime representation, d(n), of d(t) and hence the transformation
used to map H(s) to H(z). Two of the most common integration schemes are Euler's
rule (backward difference method) and the trapezoidal rule. For these two cases the
following mappings result.
Table 21: Continuoustime to discretetime mappings.
Integration Rule Mapping H(s) I4 H(z)
Euler's 1 fzl
s=\ 
T z
Trapezoidal 2(z1)
T _z+l)
So the bilinear transform is derived by computing the integral over [nTT, nT]
using the trapezoidal rule.
In Steiglitz's work [Ste65] a more theoretical treatment and significance of the
bilinear transform is developed. In fact, it is shown that the bilinear transform is an
isomorphism connecting the continuoustime vector space L2(R) and the discretetime
(or digital) vector space of doubleended sequences of complex numbers {x,, }_ that
are square summable. This connection proves that the signal processing theories with
respect to LTI causal signals and systems are identical in both the continuoustime
domain and discretetime domain. In this thesis, when a continuoustodiscrete
transformation is performed the bilinear transform will always be used.
2.6 Examples
Several examples of the 1st and 2nd order systems discussed in Chapter 2.4 are
shown on the next two pages in Figures 21 to 24. These figures illustrate the
impulse and step response behavior that must be approximated by the class of models
proposed in this thesis.
19
2
18
16
1.4
S0T= 5
1.2
T1 1 =
0.8
06
r=2
04
02
0
0 0.5 1 1.5 2 25 3 3.5 4 45 5
t (sec)
Figure 21: 1st order impulse responses, r0.5, 1, 2, 3.
0.9 T=05
r= 1
0.8
0.7 T =2
'=3
0.6
o 0.5
04
0.3
02
0.1
0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
t (sec)
Figure 22: lst order step responses, r=0.5, 1, 2, 3.
0 1 2 3 4 5 6 7 8 9 10
t (sec)
Figure 23: 2nd order impulse responses with at=l.
a) =0.4, under damped; b) &1.0, critically damped; c) 4=1.5, over damped.
 0.6
0.4
0.2
1 2 3 4 5
t (sec)
6 7 8 9 10
Figure 24: 2nd order step responses with a=l1.
a) =0.4, under damped; b) =1.0, critically damped; c) =1.5, over damped.
0.2
0
2.7 Conclusion (System Assumptions)
In this chapter the assumptions made about the systems to be modeled have been
stated; causality, linearity, timeinvariance, realizability. This discussion was to
provide the mathematical foundation and engineering principles that motivate the
development of the system models to be introduced. For the purposes of simulating
real systems, 1st and 2nd order subsystems are used as the basic building blocks. These
subsystems can be cascaded together to create high order systems. Methods were
discussed that allow for transforming continuoustime system descriptions, H(s), into
discretetime descriptions, H(z).
A good model is one with a structure that can capture the dynamics of a system
(time constants, damping coefficients, natural frequencies) without explicit
knowledge of its parameters. For systems with large time constants (slowly decaying
impulse response) and/or small damping coefficients oscillatoryy modes) a good
approximation is more difficult to achieve. The models presented in this thesis are
good approximates of systems even with these characteristics.
CHAPTER 3
MODEL FORMULATION
In this chapter our first major contribution to the problem of system modeling is
presented; namely, a comprehensive framework given by the GLFF model class.
Although many of the components have been around for quite some time it is the
unifying approach taken here that offers new insight. This new formulation
encompasses not only GLFF networks but FIR and IIR models (for the discretetime
case) as well.
3.1 General Model Form M(p)
Recall the general system modeling configuration, shown in Figure 31 below
where H is the system to be approximated using the model M(p). The specific form
of the model M(p) is now discussed.
d
H +
M(p)
Figure 31: System modeling configuration.
p is defined as a parameter set. Its elements are vectors called the parameter
vectors of the model M(p). The parameter set contains three elements (vectors),
p= {N,w,2 }
where (3.1)
N dimension vector
w weight vector
A time scale vector
The model can be written in the form
M(p)= M({N,w,A )= f(N,w,H,(.,2)) (3.2)
where H,(*, A) is called the model kernel at stage n, 0 < n < Ni, and Ni is a component
of the vector N. So M(p) is a function, f, of the parameter vectors and the model
kernel. The kernel argument is used to imply that the domain of operation (s or z)
is universal. When solving specific problems the notation H,(s,A) is used in the
continuoustime case and H,(z,2) in the discretetime case.
3.2 DiscreteTime Models
Consider the definitions above in the discretetime domain. This unifying
representation can be used to illustrate many specific model structures. In this section
it is shown how to describe the FIR, IIR, and GLFF model classes using this
formulation.
3.2.1 FIR Model
In the case of the FIR models
N
M(p)= M({N,w, })= w,,H,(z,A) (3.3)
n=O
where
2=0
n
=nH,(z)
i=0
The functions H,(z) are given by
i=O
i>O
3.2.2 IIR Model
In the case of the IIR models
M(p)= M({N,w,A )
N,
b,, H,, (z,t )
n=0
N1
l=aH,(z,A)
n=l
where
N = (N1,N2)T
w= al,az,...,al ,bob,,..., b
2=0
H(Z,)= H,(z)= H (z)
i=0
Just as in the FIR case the functions H,(z) are given by
H(z) I
ZI
i=0
i>0
3.2.3 GLFF Model
In the case of the GLFF models
M(p)= M({N,w,A])= wnH,,(z,A)
n=0
where
N=(N,)
(3.4)
(3.5)
(3.6)
(3.7)
(3.8)
(3.9)
H, (z, A) = H, (z)
H,(z)= { ,
z1
N =(N,)
w=(Wo Wi1...I^ WN)T
2 (2...L)T (3.10)
H,(z, )= Hn2(z,' )l A(zA)
i=0
The functions H,,(z,A) and H,,(z,A) could be any transfer function. In this thesis a
restriction is made to consider special cases that
1. are founded on the system model assumptions described in Chapter 2.
2. have kernels that are (typically) orthogonal.
3. have kernels that are cascaded lowpass and allpass subs) stems.
4. will allow for modeling the unknown system with real and/or complex poles
via the time scale vector A.
5. have localized feedback.
6. guarantee stability.
3.2.4 FIR, IR, GLFF Model Block Diagrams
To compare the structures of these models examine their block diagrams in
Figures 32, 33, and 34. As can be seen from the FIR model, its lack of any
feedback component limits its applicability. The IR model offers the greatest
capability but at a significant cost and loss of stability control. The GLFF model
offers a practical alternative that rivals the IIR model in performance while
maintaining simple computational procedures similar to the FIR case.
x(k)
x(k)
x(k) 
Figure 32: FIR model block diagram.
y(k)
Figure 33: IIR model block diagram.
y(k)
Figure 34: GLFF model block diagram.
y(k)
3.3 GLFF Model Kernels
The GLFF kernels can be comprised of almost any transfer function. Some of the
most popular discretetime and continuoustime kernels are given in Tables 31 and 3
3, respectively. These tables include each kernel's common name, symbol, and n'h tap
transfer function, H,,(., A). Most GLFF networks are derived from orthogonal signal
sets (see Appendix A). This orthogonality of a kernel is also indicated in the tables.
An important part of a GLFF network is the time scale parameter vector 2. It allows
for local feedback and is responsible for the popularity and usefulness of this entire
class of models. The dimension of the time scale vector, dim(A), determines the
number of degreesoffreedom a kernel offers. dim(A) is also given in the tables.
When fitting a model to an unknown system, the parameter set must be optimally
selected to achieve the minimum integrated squared error.
3.3.1 Popular DiscreteTime GLFF Kernels
The most popular discretetime GLFF kernels are the delay, Gamma, Laguerre,
Legendre, and Kautz. They are listed in Table 31, and restrictions on the time scale
parameters are given in Table 32. Some of the earliest discussions of these kernels
are found in the work by Young [You62a], Perez [Per88], and Principe [Pri93].
3.3.2 Popular ContinuousTime GLFF Kernels
The most popular continuoustime GLFF kernels are the Gamma, Laguerre,
Legendre, Jacobi, and Kautz. They are listed in Table 33, and restrictions on the time
scale parameters are given in Table 34. Some of first to describe these kernels are
Lee [Lee33], Kautz [Kau54], Armstrong Ar[Arm59], and Principe [Pri93].
Table 31: Discretetime GLFF kernels.
Table 32: Restrictions on discretetime time scale parameters.
Kernel Valid time scale
Delay none
Gamma A real, 0 < < 2
Laguerre A real, 1A <1
Legendre Ak real, Ak =k, 0 < 1
Kautz A2 complex, Ak <1
Kernel n htap H2 (z,A) Hk,(z,) H,(z,A) ortho dim(2)
symbol gonal
Delay H,(z) I z' z no 0
Gamma Gn(z,A) I AZ' A1 no 1
___1(12)z 1(1A)z
Laguerre Ln(z,A) _2 z 22 _ I yes 1
1A2z 1Az1 12 1 2z) __
Legendre P,(z,A) 1 z1 1 z A yes >1
12jz 12z 1 1 1'z
Kautz K,(z,A) 1 /1 1 n 1 yes >1
1 I Z1 1 Z i= 1A iz
Z_ i=O
Table 33: Continuoustime GLFF kernels.
Table 34: Restrictions on continuoustime time scale parameters.
Kernel Valid time scale
Gamma A real, A > 0
Laguerre A real, A> 0
Legendre A, real, A, = 2, A > 0
Jacobi 2, real, A, = k + 1
Kautz Ak complex, A, = 2AR + jA/, ARk > 0
Kernel n'h tap H H(s,A) Hk(s,) H,(s,A) ortho dim(A)
___ symbol__ gonal
Gamma G,(s,A) 1 A 2 "A no 1
s
s+A [s+A _
Laguerre L,(s,A) s A s s A "" yes 1
s+A s+A s+A s+A
Legendre P,(s,A) 2j s Ak n2 fSA yes 1
s+A s+k s+ A, =0 S+
Jacobi J,(s,A) s 2 2_ 22 s A, yes >1
s+A_ s++ +k s+ i0 s+
Kautz K,(s,A) 2(sA) (sAkX(s) 2f(sI)( sXs) yes >1
(s+AlXs+A) (s+AkXs+ ) (s+2Xs+x ) s+2+Xs+4)
3.4 GLFF Kernel Derivations and Properties
The Laguerre, Legendre, Jacobi, and Kautz GLFF networks have the interesting
property that their kernels are orthonormal. Hence, we call these orthonormal
generalized linear feedforward (OGLFF) networks. Next these kernels will be
derived from first principles. This is done in the continuoustime domain. The
derivations will make evident their relationship to causal LTI system. The connection
between the network structures and the mathematical description of causal LTI
systems is the primary reason they are able to approximate these systems so well.
3.4.1 ContinuousTime Laguerre Kernel
Consider the sequence I(At) e" }. Regarding this as the solution to an ordinary
differential equation that represents a causal LTI system, the system would consist of
one pole with multiplicity n+l. Orthogonalizing this sequence yields a new sequence
{l((t, )} given by
1,(tA)=V2 I e,nA> 2>0,2>0, t0 (3.11)
Now taking the Laplace transform gives the Laguerre model kernel
L,(s,)= I {l(t)} S=_ An(3.12)
SA S+A
3.4.2 ContinuousTime Kautz Kernel (Complex Time Scale Vector)
Consider the sequence {et '} where the A,'s are distinct complex scalars (i.e.,
components of the time scale vector 2). Since the poles are complex they must exist
in conjugate form in order to satisfy the realizability requirement. Orthogonalizing
this sequence yields a new sequence {k,(t,A)} given by
k,,(t,A)= Re{A,,e1}, n>0, t> 0, Re{A,}>0, Vi (3.13)
i=0
It is now shown how to obtain the coefficients A,i. Taking the Laplace transform of
the new sequence gives the Kautz model kernel
2A(s+In) (s (sA,)(sA) n1
K2,,(t,A)=I {k, (t, )}= ( + m 0,1,..., 2
Sr2AR (sIA) (s1(A,)(sAX)' n1
Kz,,, (t, )= L(k,,(t,A )}= s_+ +,), ,yn i m=0,1,...,
m+(tA)= km+(tA) (s+A +X)(s+ ) ,U (s+A,)(s+A ) 2
Express these functions in partial fraction form. Grouping each pair of time scale
elements (2, ,2 ) gives the kernel expression
K,,(t,)= s+2A, + n=2mandn=2m+1 (3.15)
Hence the coefficients Ani are twice the residue at the pole A, for the nth tap kernel.
3.4.3 ContinuousTime Kautz Kernel (Real Time Scale Vector)
Consider the sequence {e"')} where the A,'s are distinct real scalars (i.e.,
components of the time scale vector A). Orthogonalizing this sequence results in a
new sequence {k,(t, )} given by
k,(tA) = Aie', n >0, t > 0,A, > ,Vi (3.16)
i=0
It is now shown how to obtain the coefficients Ai. Taking the Laplace transform of
the new sequence gives the Kautz model kernel
K,(t, )=I k,,(2)}= 22 k( sA =, A+, (3.17)
S + n i= S +i i=0 S + Ai
where the last expression is simply the partial fraction expansion of the product form
of K,(t, 2). Hence A,i is the residue of the pole at A, for the Inth tap kernel. It will be
real since the time domain function k,(t, 2) is real.
3.4.4 ContinuousTime Legendre and Jacobi Kernels
Consider again the sequence {(e"'} where the An's are distinct real scalars. This is
the same sequence from which the Kautz kernel with a real time scale vector is
derived. This kernel is described by equation (3.17). Now consider the kernel with a
"special" time scale parameter. Namely, if A = (AA,,2.,2,)T with k = 2A, A > 0
then we obtain the Legendre kernel. Likewise, if Ak = k +1 we obtain the Jacobi
kernel.
3.5 Impulse and Step Response Modeling
The models discussed in this chapter can accurately represent any signal or system
belonging to L 2(R). We will be most interested in modeling causal LTI systems using
sampled step response data. In this case, a step function is applied to a continuous
time system to generate its step response. The step response is sampled to obtain a
discretetime representation of the system. Optimal continuoustime and discrete
time GLFF models are determined using this data (optimization will be discussed in
Chapter 4). For example using an N+I tap Kautz filter the approximating transfer
function is
N
H(z) = H(z)= w,K,(z,A) (3.18)
n=0
Suppose the discretetime impulse response is desired. How can it be computed?
Three simple methods are now given.
1. Using the discretetime model /H(z), pass the impulse x(k)=J(k) through this
model to generate an approximation to h(k).
2. Since x(t) is a unit step then X(s) = 1/s. If g(t) is the step response of the system
then g(t) = h(t) x(t) or G(s) = H(s) / s. So
H(s) = sG(s) = D(s)G(s) (3.19)
or
H(z) = D(s)G(s)I, = D(z)G(z) where D(z) = (3.20)
1 T z+1
Using the sampled step response data g(k), pass this through D(z) to generate h(k).
dg(t)
3. Since h(t)= d use a numerical approximation to the derivative (backward
dt
difference perhaps) to obtain h(k). In this case
dg(t) g(kT) g((k)T) g(k) g(k 1)21)
dt ,kT T T
3.6 Conclusion (Model Formulation)
This chapter has dealt with structures of models useful for approximating causal
LTI systems. The generic model M(p) was stated. It was shown how this formulation
can be used to represent FIR, IIR, and GLFF networks. Several GLFF kernels were
derived, and discussion was provided to illustrate their connection to the assumptions
of the systems under study. Tables of discretetime and continuoustime GLFF model
34
kernels were included. These kernels are the most popular in the literature and as the
derivations of the Laguerre, Kautz, Legendre, and Jacobi revealed they are not
theoretical abstractions, but orthogonalized versions of solutions to differential
equation representations of several important classes of systems. Now that the model
structure M(p) has been defined, optimization of the parameter set, p, becomes the
focus.
CHAPTER 4
OPTIMAL PARAMETER SELECTION
In this chapter the selection of the optimal (or near optimal) parameter set for
GLFF networks is discussed. Recall p = {N, w,A } where N is the dimension vector,
w is the weight vector, and A is the time scale vector. Optimization for these three
vectors is handled separately below. In fact, the ability to separate the parameters and
their optimization techniques is another benefit of the GLFF formulation.
4.1 Model Dimension
Consider a GLFF model consisting of a set of time domain basis functions
{o, (t, )}'o with frequency representation {I,(s,A2)}o The model structure is a
weighted sum of these basis functions (kernels). In the sdomain
Y(s)
M(p) = w,(s,A)= WTc(s,2) (4.1)
X(s) ,,=O
where
w= (WoWI...WN)
(4.2)
D(s, 2 )= (o0(s, 2>, ,(s,2),..., N (s,2))T
In the time domain
y(t) = wA(t, )* x(t) = [w O(t,2 )*x(t) (4.3)
where
0(t,2) = (o0(t,2 ),1 0(t,2),..., ON(t, ))T (4.4)
In this thesis, the theory will focus on modeling the impulse response of a
continuoustime system. Given sampled step response data, continuoustime and
discretetime GLFF networks will be constructed to approximate the desired system.
Hence the input, x(t), will be considered an impulse and the desired signal, d(t), will
be the system impulse response.
The selection of N is now discussed. N is a scalar for the GLFF model class.
Suppose N is a fixed finite positive integer. Then d(t) can be approximated by
N
d(t) =y(t)= w,O,(t, A) (4.5)
n=0
4.1.1 Completeness
To ensure that y(t)  d(t) as N  oa the property of completeness is employed. A
set of functions {O,(t,2)}) is called complete [Su71] on the space L2([a,b]) if there
exists no function d(t)e L2([a,b]) such that
b
d(t)O(t, A)dt = 0, Vn = 0,1,2,... (4.6)
a
An equivalent definition is as follows. Given a piecewise continuous function
d(t)E L2([a,b]) and an E > 0 there exists an integer N > 0 such that
b
J(d(t) y(t))2dt < e (4.7)
a
where y(t) is given by the equation (4.5). All of the GLFF kernel sets are complete on
[0,oo). An example of an incomplete set of orthonormal functions {n0(t,2)}) o
defined on [0,2;r) is
(t,=) = 0,(t)= sin(nt) (4.8)
The function d(t)=cos(mt) is nonzero on [0,2r) however it is orthogonal to
{0,(t,')} 0on this interval. Namely,
21r
Jcos(mt)O,(t,2 )dt = 0, Vn,m (4.9)
o
4.1.2 Orthogonal Projection
When {\,(t,2)}Nois an orthonormal set then
N
y(t)= wOn(t,2A) (4.10)
n=0
spans an N+1dimensional Hilbert Space, S (see Appendix B). The Orthogonal
Projection Theorem states that for every d(t) in a Hilbert space H, it can be
decomposed as d(t) = y(t) + e(t) where y(t) S, e(t)e SI. S is the orthogonal
complement of S and H = S E S So y(t) can be considered as the orthogonal
projection of d(t)E H onto the space S. Namely
y(t)= Projs(d(t)) (4.11)
and {, (t, )}N is an orthonormal basis of S.
n=O
4.1.3 "Optimal" Selection of N
Completeness is important because it guarantees that if the approximation of d(t)
using N+1 weighted terms of {0((t,2 )}n0 is not within a prescribed tolerance then
more terms can be included to eventually reach the tolerance requirement. Increasing
N may not improve the approximation of d(t) when a set of incomplete functions is
used.
Having the completeness property provides the following choices when a GLFF
network realization for a given N is inadequate:
1. Increase N until the prescribed tolerance is reached.
2. Select other kernels from the GLFF set and examine the approximation error.
For a given N, if none of the GLFF kernels chosen provide satisfactory results
then N must be increased. If this is done using several kernels the rate of convergence
for each kernel as N is increased could be monitored. Choose from this set the kernel
that yields the minimum N needed to achieve the required tolerance. Other options
are to either raise the tolerance level or choose a different performance criteria. One
suggestion for setting the tolerance is to choose the acceptable error power to be a
certain percent (5% perhaps) of the system power.
Davila and Chiang [Dav95] suggest another approach to estimating the dimension
vector N = (N1, N2) for a general IIR structure. This method examines the
eigenvectors of the input/output data covariance matrix. When the model orders are
overestimated, zeros will appear in the noise subspace eigenvectors. The number of
zeros found can be used to estimate the model orders (dimension vector components).
This procedure could be applied to GLFF networks.
4.2 Weight Vector
One of the greatest advantages to working with a feedforward model is that the
calculation of the optimal weight vector w (in the squared error sense) is simple. This
is true even for GLFF structures. Although the realizability requirement for GLFF
networks demands the weight vector be real, the optimal solution for w is now
derived assuming w, d(t), and {0,(t,2)} N are generally complex. Realizable optimal
solutions are obtained by setting the imaginary components equal to zero.
4.2.1 Optimal Selection of w
Assume d(t) is the impulse response of the desired system (i.e., d(t) = h(t)) and
approximate it with y(t) given by
N
y(t)= w ,(t)
n=0
(4.12)
The parameter set p= {N,w,A } is to be chosen so the error e(t) = d(t) y(t) is
minimum in the integrated squared error sense. Namely, choose p to minimize
If Wn = a, + jbn then
J= f d(t) y(t) ldt
0
N
= ld(t)ldt dw d'(t) (t, A)dt(4
0 i=0 0
E wf d(t)0*(t,A)dt + w,I2 ji,(t,A)2dt
i=0 0 i=0 0
dJ dJ
set = 0 and = 0 to find w, that minimizes J. Hence
da, db,
S= (t),(t, A)t d(t)o0(t, A)t + 2a, f (t,)1dt =0
da, 0 0 0
(4
.13)
.14)
40
so
jRe[d(t)*(t,A)]dt
a= 0o (4.15)
j On (t, A)2dt
0
When {(O(t,2 )}n is an orthonormal set, the bottom integral equals unity yielding
a, = Re[d(t)O (t,A )]dt (4.16)
0
Likewise b, can be derived
Sj d(t)o((t,2 )dt jld(t) (t,A)dt +2b (t, )dt = 0 (4.17)
n 0 0 0
so
0
Im[d(t)*,(t, )]dt
b. = 0(4.18)
Again, if {(,(t,A)}N is orthonormal then equation (4.18) reduces to
b, = Im[d(t)o*(t, A )]dt (4.19)
0
Hence, combining the above results gives
w, = a + jb = Id(t)o(t, A)dt (4.20)
0
When {(Z,(tr,2 )} =is real then w, is real and can be written as
w, = d(t)O,(t, A)dt (4.21)
0
So
w= d(t)O(t, A)dt (4.22)
0
is the weight vector that minimizes the energy of e(t). From this point forward only
real weight vectors will be considered.
4.2.2 Relationship to WienerHopf Equations
It is easy to show that w is in fact the solution to the WienerHopf (or Normal)
equations. The performance function can be expanded as follows
J= e (t)dt = (d(t) y(t))2dt = d(t) w'(t, A ))dt
o o 0
0 0 0 (4.23)
= d2(t)dt 2wTjd(t)o,(t,A)dt +wT jO(t,2A) (t,A)dt w
0 0 0)
now define the following terms
d,= d(t)dt
o
p jd(t)O(t, )dt (4.24)
0
R O(t, )O(t,A)dt
0
Then the performance function J becomes
J = d 2wrp+w Rw (4.25)
To minimize this function take the gradient of J with respect to the weight vector w
and equate to zero
VJ = 2p + 2Rw = 0 (4.26)
Rw = p (4.27)
These are known as the WienerHopf or Normal equations. Solving for w gives the
optimal solution
w = Rp (4.28)
Now the OGLFF kernels {jO(t,/N)} form an orthonormal set so
0
Hence the weight vector becomes
0
which is the solution derived in the previous section.
4.2.3 Minimum Squared Error Performance Function
Now that the optimal weight vector has been derived, the minimum squared error
can be computed. Recall
J = d, 2wrp + wrRw (4.31)
For the optimal weights it was found that (using OGLFF kernels) w = p and R = I.
Plugging these into the equation for J yields
Jopt = de wTW (4.32)
The above optimal weight vector w and minimum squared error Jopt was derived using
a fixed time scale A. Notice that w and Jopt are functions of 2. Hence Jopt can be
written
(4.33)
Jopt(A ) = de wT(A)w( A)
to illustrate the importance of A on the optimal solution. Appendix C presents a
derivation of equation (4.32) in the sdomain.
4.2.4 Computing w Using Step Response Data
In Chapter 5 applications of the GLFF networks are given. The systems there are
continuoustime systems. The available information of these systems is sampled step
response data. GLFF networks will be used to approximate the input/output step
response behavior. Methods for computing the optimal weight vector w using
sampled step response data are now given.
Least Squares Solution. [Hay91] Consider the model dimension N and time scale
vector A fixed or previously chosen to be optimal. Discussion of model dimension
was given in Chapter 4.1. Discussion of the optimal time scale will be given in
Chapter 4.3.
Since the information is sampled (discretetime) step response data the following
definitions are needed. Let the input, x(t), be the unit step and g(t) be the response of
the system under question. If the step response is sampled at interval T then define
g(k) = g(kT), k=0,1, ..., K as the desired signal. Define also
xi(k) = kth sample at the ith tap of the GLFF network (4.34)
y(k) = kth output of the GLFF network
In this discretetime case, the matrix R and vector p are given by
K
R=[r,] where r = (xi,,x) = x(k)xj(k) (4.35)
k=0
and
SPO K
p= : where p, = (g,x,)= g(k)x,(k) (4.36)
k=0
\PN)
This gives the normal equations Rw = p with solution w = R'p. This is the Wiener
Hopf solution. Notice in this case that R will not be the identity matrix. Only when
the input is an impulse (or white noise) will R = I. Nevertheless, the WienerHopf
solution provides the optimal weight vector (for fixed N and A) regardless of the
input.
Integral Equation Solution. [Clu91] A more direct evaluation of the weight vector
solution given in Chapter 4.2.1 is now considered. Recall
w,=fd(t)0,(t,A)dt, n= 0,1,...,N (4.37)
0
is the optimal weight solution assuming an orthonormal basis. For the system
modeling problem this is accomplished by applying an impulse (or white noise) into
the system. In which case d(t)=h(t), the impulse response. Consider, for example, the
Laguerre kernels given by {,(tA)},0={l,(t',)} N0. Some of its properties can be
utilized to derive a simplified integral equation solution. Again work with g(k), k=0,1
,..., K samples of the continuoustime step response g(t). Since the continuoustime
impulse response is h(t) = dg(t)/dt, the optimal weights can be written as
w, = h(t)l,(t, )dt
0
= l,(t,A)dt (4.38)
0 dt
= (g(t)1n(t,A )I j g(t) dt
00 dt
where integration by parts has been used. For a stable system g(t) it must have
lg(oo) < oo. Also using the property of the Laguerre functions [Par71]
i,(t,A2)=Z ei (t"e"2) (4.39)
n i t"
yields
lim l (t,) = 0 (4.40)
Hence
w =g(t) t dt (4.41)
dt
0 t
Let T., be the steady state response time of g(t); namely, T,. is the time which for all
t > T, the output can be consider constant, g(t) = gs, = constant. In this case the
integral can be broken up as
T', In(tA) ddgll( dt
w. = g(t) dt dt g dt
ST (4.42)
= g(t) (t, ) dt + g,J (T,I A)
0ldt
From the Laguerre function property above it is easy to compute the derivative
relationship
S A i(t, A) 2A 1(t,A) (4.43)
dt i=0
The weight computation becomes
n1 I T ,
w.= 22ZJ g(t),(t, A)dt + A (t)l(t, )dt + gsl (T,I,A) (4.44)
i=0 0 0
Up to now the continuoustime step response g(t) has been assumed. To compute w,
using the sampled step response g(k) an integral approximation technique must be
used such as Euler integration
r,, K
Sg(t)1,(t,A)dt = A g(k)1,(k, A) (4.45)
o k=O
where
T
A =  (integration interval) (4.46)
K
Finally this gives the integral equation weight solution
n1 K K
w, = 2/A _g(k)l((k,A)+AAE8g(k)l,(k,2)+ g,sl,(T,) (4.47)
i=0 k=0 k=0
4.2.5 Adaptation of w
Note that the weights can be computed by solving a set of linear (WienerHopf)
equations. This linear relationship among the weights corresponds to a convex
performance function with one minimum point. So the weight vector can be easily
solved in an adaptive manner using a gradient search technique such as Newton's
method or LMS [Wid85]. Unfortunately, this nice relationship does not hold for the
time scale vector. It is (in general) nonconvex with multiple stationary points. The
number of which depends upon the dimension, N, of the model used.
4.3 Time Scale Vector
Optimization of the parameter set elements N and w have been discussed to this
point. Their solutions are fairly easily found for GLFF networks and are the only
elements of p that are required for FIR and IIR models. As has been previously
stated, the power of GLFF networks is in local feedback. This feedback is controlled
by the time scale vector. Although IIR structures also have feedback, it can not be
separated and optimized independently from the taptotap weights like the GLFF
structures.
Unfortunately, the elements of the time scale vector are not linearly related; hence,
the performance function J is usually nonconvex w.r.t. these parameters. As such,
the bulk of the work in finding an optimal p involves finding an optimal 2. This
section is devoted to this issue. It begins by examining several methods that currently
exist for either determining exactly or approximately the optimal A. Each method
usually addresses a particular GLFF kernel, most often the Laguerre. A few papers
have considered the nonorthogonal Gamma and a few discuss optimization of the
more complex Kautz kernel. Although the exact optimal solution is sought, Kautz
has pointed out that nonoptimal solutions (which may be derived in a simplified
manner) can still produce excellent results. His comments suggest that one should
not search to exhaustion for optimality if it is possible to get close with much less
work.
After the overview of optimization methods found in the literature, we propose a
new technique relying heavily on complex variable theory. This is a method that
extends the concepts of squared error analysis to the complex domain. Under certain
assumptions, the optimal time scale can be fairly easily determined by simply looking
for zero crossings.
4.3.1 Optimal Selection of A (Historical Review)
In this section a summary is given of current methods found in the literature to
either locate or approximate the optimal time scale A. Each method usually addresses
a particular GLFF kernel, most often the Laguerre. A few papers have considered the
nonorthogonal Gamma and a few discuss optimization of the more complex Kautz
kernel.
The methods presented in the literature can be grouped into the following
categories according to the employed strategy or desired objective:
1. Approximation via the impulse response.
2. Asymptotic solutions.
3. Achieving prescribed measurement objectives (other than minimum squared
error).
4. Matched filter approach.
5. Moment matching.
6. Satisfying necessary conditions for stationarity of the squared error function.
7. Iterative and adaptive search techniques.
Approximation via the impulse response. These methods approximate A given a
priori knowledge about the impulse response h(t) of the system. The required
information is different for each approach. Either h(t) is needed in analytic form or
tabular form (samples of h(t) at discrete points, t = tk, k = 0,1,2, ... ). These methods
appeal to the necessary requirement for all stable and realizable linear systems;
namely, its transfer function H(s) must be representable as a rational function with no
righthalf plane poles or multiple jaxis poles. Under this assumption, h(t) can be
represented as a sum of exponentials of the type
(A + At+...+Atm")eA' (4.48)
where
Ak = ARk + k, Rk < 0 (4.49)
A frequency domain method for approximating A is as follows [Su 71]. Suppose
h(t) is known (or approximation to it) and its Laplace transform H(s) is determined.
Note that the form of h(t) that is given may not produce a rational H(s). Now find a
rational approximation Ha(s). For instance use the expression for H(s) and expand it
into a Taylor series about s = 0 to obtain
H(s) = ao + as+a2s +* (4.50)
Then a rational function can be found to approximate H(s). One way to accomplish
this is to use a Pad6 approximant. The poles of this approximation can be used as an
initial guess for the poles of H(s). These pole locations could be incrementally
adjusted to locally optimize.
A time domain approach is to use point matching. The basis of this technique is
attributed to Prony [Ham73]. The strategy is to force an approximate impulse
response h(t) to match h(t) at evenly spaced points. Assume h(t) is known at points
separated by the spacing At. Using the approximation
N
(t)= Ae"' (4.51)
it is possible to create N linear equations (using the known points of h(t)) whose
solution gives coefficients ao, a,, ..., aN. of an Nh order characteristic polynomial
xN + a_,XNI + aN_2x + + axx + a0 = 0 (4.52)
Once the coefficients are known, find the roots of this polynomial x, xz, ..., XN. They
are related to the time scale parameters by
.i = nx, (4.53)
At
Another time domain approach involves the use of an orthogonal signal set
[Kau54]. Assume h(t) and its first N derivatives are known. Using these functions,
generate a set of orthogonal functions {(0(t)}N0 satisfying the relationships
o(t)= h(t)
01(t) = h(t)+ b o0o(t)
(t = (t) + b, ,(t) + bz20o(t) (4.54)
N(t) = hI(N(t)+bN(N N_ l(t)+**+ bNlOl(t)+ bNOo(t)
where the b,,n's can be found such that
0
J^.(t)O (tdt= (4.55)
Since each 0~(t) is a linear combination of the derivatives (from the 0th to the nth) of
h(t) then 0n(t) can be rewritten in terms of h)(t), j = 0,1,...,n1. Replace each h()(t)
by A' to get the characteristic equation
A + an(nl)/ + an(n_2),2) +" + +a + an0 = 0 (4.56)
The roots of these equations give the time scale parameters Ani, i=1,2,...,n for model
order n for each n = 1,2,...,N.
Asymptotic solutions. Several methods are offered that are based on power series
representations of the Laguerre model. Schetzen [Sch70], [Sch71] presents an
analytic method for determining the asymptotic optimal time scale LA where
2 = lim AN (4.57)
N)
This is based on the region of convergence (ROC) of the power series equivalence of
the Laguerre model. Justification for A~. in the place of AN is made on the basis that
most engineering applications require the use of a large number of terms N in the
Laguerre expansion to achieve accurate results. In this case, /L is viewed as a good
approximation to the optimal time scale AN.
Bruni [Bru64] offers another approach involving optimization around the ROC of
a Laguerre model expansion. Assuming h(t) is the impulse response of an LTI
lumped parameter system
M
h(t)= Ae p', p,=a,+jb,,a,>0 (4.58)
i=0
decompose a Laguerre model representation of h(t) into M infinite series
h(t) = wnt,l(t, 2) = W"l,(t,2 ) (4.59)
n=0 n=O i=1
where w,i is the nth coefficient of the expansion of the ith exponential in h(t). With
this formulation the Fourier transform of h(t) can be represented as a sum of M
infinite complex geometric series. Mi gives the magnitude of the term in the ith
infinite series. Mi should be minimized for each i in order to achieve the most rapid
convergence of the series. Since M, is a function of A, there is an optimal A (denoted
Ai) for each Mi. Choose as the optimal A for the entire model, the Ai that minimizes
the largest Mi. Methods similar to this one are also given in [War54] and [Hea55].
Achieving prescribed measurement objectives. In addition to locating A that
minimizes the squared error, approaches that minimize other measurement criteria
have also been considered. One such method employed when using the Laguerre
signal set {l(t,2A)}N0is given in [Par71]. Consider the approximation of d(t) using
an N+1 tap Laguerre model
N
d(t) =y(t)= wl, (t, ) (4.60)
n=0
Now define the following quantities
ftd2(t)dt td2(t)dt
M, o and M2 0 (4.61)
Sd2(t)dt d2(t)dt
0 0
M1 is a measure of how rapidly d(t) decays. M2 is a measure of how smoothly d(t)
decays. Parks show that the time scale A given by A = M1/M2 minimizes the
maximum error over the class of all signals with these measurements. The advantage
to this approach is that a complete knowledge of the signal to be approximated is not
required.
In [Fu93] a similar approach is taken when using the discretetime Laguerre
model. In this case, the optimal A is found that minimizes the measurement
J =tnw2 (4.62)
n=0
This performance measure is chosen because it linearly increases the cost of each
additional term used in the Laguerre expansion. This enforces a fast convergence
rate. The optimal A is obtained using a discretetime formulation of M1 and M2
above. A relationship for J is derived that depends only on M1, M2, and A. The
optimal 2 is then determined to be
S2M= 1M2 (4.63)
2M, 1 + 4M ,M M 2M2
Matched filter approach. In [Kuo94] a matched filter bank is used to estimate the
optimal 2 for the Gamma model. Recall that {g,(t,, )}N0 is a nonorthogonal GLFF
signal set. Let d(t) be approximated by
N
d(t)= y(t)= w,g,(t, A) (4.64)
n=O
A matched filter is constructed based on the optimal A condition (g ,(t,2),d(t) = 0
where g ,(t,A) is the component of gN+i(t, 2) that is orthogonal to the previous basis
components. A bank of such filters is created where each bank uses a different A. To
optimize 2 over the region (0, ,,ax) an Mbank filter structure is created where the mth
filter uses the time scale A = A,,x / M. When the desired signal d(t) is applied to the
filter bank, a change in sign of the M outputs will occur between filter m and m+1 for
a value of A that satisfies the optimality condition. When there are multiple
occurrences, the associated A for each case must be used to determine which one
achieves the minimum squared error.
Moment matching. The method presented in [Cel95] uses a Gamma model to
approximate a signal d(t). The Gamma model is conceptualized as a rotating manifold
in space. The kernels {g,(trA)} 0 are the basis vectors that define the N+I
dimensional space of operation and the time scale parameter A is the feature that
allows for rotation. The optimal A is estimated using a moment matching method.
This method attempts to estimate 2 by equating a finite number of time moments of
y(t) (the output of the Gamma model) to those of d(t). Namely, it is required that
mn(y(t))= m,(d(t)), n = ,...,N where
mt(y(t)) t"y(t)dt (4.65)
0
Using these functions an N+1 order polynomial in A, p(2), is formulated. The roots of
p(A) are estimates of the local minima of the squared error performance function J.
The root that minimizes J is chosen as the optimal A.
Satisfying necessary conditions for stationarity of the squared error. Some of the
most popular methods for optimizing the time scale of the Laguerre model center
around the necessary conditions for stationarity of the squared error performance
function J. For instance, approximating d(t) with an N+1 tap Laguerre network
N
d(t) = w,(A)1,(t, A) (4.66)
n=0
where the weight vector components wn(A) have been given an argument in 2 to
emphasis the dependence of these weights on the time scale. For a fixed 2, the
squared error is minimized when the weights are computed as
w,(2) = d(t)l,(t,A)dt (4.67)
0
The squared error function reduces to
N
J= d2(t)dt w,(A) (4.68)
0 n=0
To minimize J (over A) the evaluation of the X must be maximized. This requires that
dN 2 NA = N d
w,() )2w,(A) w,(A2)= 0 (4.69)
dn=O n=O d
which yields the interesting relation
WN(A)WN+(2A)=0 (4.70)
A proof of this assertion is given in Appendix D. Hence the optimal value of 2 is
either a root of WN(2) = 0 or WN+il() = 0. Solving for A analytically is a tedious task
because WN(A) = 0 is a polynomial in 2 of degree (N+M) if d(t) is a system with M
poles. The first proof of this (using mathematical induction) was given by Clowes
[Clo65] and was restricted to the case where d(t) had only simple poles. A more
concise and less restrictive proof was given in [Kin69]. This paper also derives
expressions for the 2nd derivative
d2N [ 1 )[(N+ 4() NWN+,(A)wN(A)] for wN(A)= 0
dt2w( E W (4.71)
dt n=0 [N (N+2)w2(A)w,(A)(N+l)()] for wN+,(a)=0
These equations may be used to determine whether a stationary point is a maximum
or minimum. In [Mas90] a discretetime formulation of the stationary conditions is
formulated and the authors use numerical techniques to find the roots of the
polynomials WN(A) = 0 and WN+I(A) = 0. [Sil93] derives these conditions for the
Gamma model. The derivation given by Clowes was valid only for systems with
rational transfer functions. [Wan94b] extends the results to approximating any L2(R)
stable linear system. [Sil95] describes another efficient way to compute the
dA
weights WN(A) and WN+I(A) via Taylor series and Pad6 approximants. This procedure
is illustrated in both continuoustime and discretetime when approximating a system
excited by an arbitrary input.
Iterative and adaptive search techniques. When representing a system with a
discretetime GLFF network, the elements of the time scale vector are restricted to
bounded regions to guarantee stability. In particular, the discretetime Laguerre
model requires 2e (1,1). One method of optimization is to employ exhaustive search
over this interval, seeking the 2 that minimizes the squared error performance
function. Since this would be done in discrete steps, the only way to get continually
closer to the optimal A is to increase the number of points in the search region.
[Kin77] and [Mai85] discuss minimizing the squared error by employing a Fibonacci
search technique over 2. This is an interval elimination procedure wherein the region
in which the minimum lies is sequentially reduced using a predetermined number of
steps.
[Nur87] uses an adaptive technique whereby the optimal A is chosen to be the
center of gravity of the presumed poles of the system. The method works as follows.
Initialize 2, say Ao. Using the Laguerre model, with this 2, the parameters of a
difference equation representation of the system are estimated. From the
characteristic equation determine the poles of the estimated system. Choose as Ai,
i=1,2,... the center of gravity of these poles. Repeat this process until A2+ ,i <
where Sis a small number. Once the center of gravity converges to a fixed point stop
the adaptation.
4.3.2 Optimal Selection of A (Complex Error Analysis)
A new approach to estimate the optimal time scale is now offered. The
methodology used here is one of extending the squared error performance function to
the complex domain. This is accomplished by defining a new complex performance
function. It will be shown how this function is derived, how it relates to the
conventional squared error function, and what additional information it offers. The
optimal time scale can be determined by looking for zero crossings of the imaginary
component of the new complex performance function.
For systems consisting of a cascade of 2nd order subsystems containing pairs of
complex conjugate poles, examples will be presented that illustrate how the imaginary
component of the new complex performance function can be used to obtain the
optimal (in the minimum squared error sense) time scale parameter. These examples
will solve the system identification problem with a GLFF model containing the
Laguerre kernel.
4.3.2.1 Derivation of the complex performance function. In Appendix C an s
domain description of the squared error performance function is discussed. It is given
by
J = E(s)E(s)ds (4.72)
C
where C is a contour taken in the CCW direction enclosing the entire lefthand plane
(LHP) in s, E(s)= D(s) WT (s,), and w is a real weight vector. For a fixed time
scale A the optimal weight vector is
w= fD(s)(s, A)ds (4.73)
C
yielding the performance function
J = de w'w (4.74)
where
d, D(s)D(s)ds (4.75)
C
Now consider restricting evaluation of the performance function to the upper LHP
only. Namely, create the complex performance function
J , E(s)E(s)ds (4.76)
c.
Make the following definitions:
due, = D(s)D(s)ds (complex energy)
c, (4.77)
wu = D(s)4(s, A)ds (complex weight vector)
The contour Cu is taken in the CCW direction enclosing the upper LHP in s. Now if
F(s) is a function which has real and/or complex conjugate pair poles then
f F(s)ds= 2Re F(s)ds (4.78)
c c,
hence
w = 2 Re[wf]
w 2 "J(4.79)
d, =2Re[d,,]
By making the above definitions, the problem of locating the optimal A using Ju is
studied. By operating in the complex space there is information that can be used to
offer additional insight into the optimal solution of the time scale parameter 2. This
information is not available using squared error analysis because the imaginary
component of Ju is canceled when calculating J. This is due to the contribution from
the poles of D(s) in the lower LHP in s. The squared error performance function can
be considered as a mapping onto R and the complex performance function as a
mapping onto C.
A simplified expression for J, is derived as follows
Ju = 2 E(s)E(s)ds
c=
7 f [D1s) WTq (sA )][D(s) WTD(s,2 )]ds
C9 ,
(4.sU)
= D(s)D(s)ds w' D(s,(s, A)ds
C. c.
w D(s)((s, A)ds +w ~t(s, A)Q~(s,)dsw
Q c,
Using the definition of the complex energy, due, and complex weight vector, w,, as
well as the following relations
and
WT j (s, A)T(s,A)ds w = W'w
C,
jD(s)D(s,A)ds =w
c.
(4.81)
(4.82)
gives
Ju = 4 ~WTW WTW+ WTW
= de WT WU
Note that
2Re[J,] = 2Re[d.,] w2 Re[w,] = d, WTW = J
which is what is expected since
(4.83)
(4.84)
J = E(s)E(s)ds= 2 Re 
c C
= 2Re[J ]
(4.85)
60
4.3.2.2 Complex performance function examples and experimental results. Three
test cases have been chosen to study the complex performance function, Ju. One 2nd
order system and two 4th order systems are modeled using a GLFF network with a
Laguerre kernel. The transfer function for each example follows:
Example 1: (2nd order system, rapidly decaying impulse response)
1 1
H(s) = + (4.86)
s+2j s+2+j
Example 2: (4th order system, rapidly decaying impulse response with slight
oscillation)
1 1 1 1
H(s) = + + + (4.87)
s+2j s+2+j s+l2j s+1+2j
Example 3: (4th order system, slowly decaying impulse response with oscillation)
2 2 3 3
H(s) = +  + + (4.88)
s+0.4j s+0.4+j s+0.62.8j s+0.6+2.8j
Experimental results for each example are tabulated and graphed on the next few
pages. In Figures 41 and 42 the squared error performance function, J (solid curve),
and twice the imaginary component of the complex performance function, 2Im[Ju]
(dash curve), are plotted for the first four model orders for example 1, Figures 43
and 44 for example 2, and Figures 45 and 46 for example 3. Notice the imaginary
component has zero crossings near the stationary points of the squared error
performance function. Im[Jj] crosses negativetopositive near local minima of J. So
the complex performance function provides an estimate of the optimal time scale
through the imaginary component. In fact Im[J,] can be considered an estimate of
V,J(A) that is best near Aop.t
Numerical results are provided in Tables 41, 42, and 43 to illustrate how close
the nearest zero crossing (nearest ZCu) of 2Im[Ju] is to the optimal time scale as the
model order increases. Also the squared error is computed at the optimal, J(Aopt), and
nearest ZCu, J(Anzcu), time scales. For this analysis, Im[Ju] was scaled by 2 to provide
similar scaling to J since J=2Re[Ju]; however, the zero crossings are unaffected by
this scaling. By examining J(Aopt) and J(Anzcu) it is found that when using a Laguerre
network with sufficient order to approximate the unknown system, no degradation
occurs by choosing A at the nearest ZCu. In fact the small difference in ,,op, and Anzc,
is partly due to numerical approximation error when searching for minima of J and
roots of Im[J,].
10
o 20
=L
a,
30
C
40 1
a)
0
z
60
70
0 1 2 3 4 5
Time Scale, X
Figure 41: J for example 1 (orders = 0 to 3).
0.1 '
0.05 
E 0O.
0.05
0.1
0
3
x10
1 r
(order=0)
1 2 3 4 5
(order=2)
05 =  d
0.5 .
E 0
1 S
0 1 2 3 4 5
0.01
0.005
(order=1)
0 .005 O 4 t.............. .................... ................... .
0.01 S
0 1 2 3 4 5
4
x10
4 r
x
(order=3)
0 1 2 3 4 5
Figure 42: J and 2Im[Ju] for example 1.
0
5
0 10
o
S15
u
a)
2 20
CO
25
30
S35
40
45
(order=0)
2
(order=2)
2 4
S0.2
S0.1
 E0.

0.1 *
0.2
6 0
0.1
0.05
0
0.05
0.1
Figure 44 J and 2Im[J for example 2.
Figure 44: J and 2Im[Ju] for example 2.
(order=1)
2
(order=3)
2 4 6
1 2 3 4
Time Scale, X
Figure 43: J for example 2 (orders = 0 to 7).
0.4
0.2
0
0.2
0.4
0.1 r
0.05
0
0.05
0.1
0
1 2 3 4 5 6
Time Scale, X
Figure 45: J for example 3 (orders 0 to 11).
(order=O)

0 2 4 6
(order=2)
.
i
(order=1)
  
(order=3)
0 2 4
Figure 46: J and 2Im[J,] for example 3.
CO
0
c 10
U
0
15
" 20
z25
30
1
j'_0.5
E
0
0.6
0.4
E 0.2
0.
0
0.2
Table 41: Performance function characteristics for example 1.
N oApr A2nzcu J(Ao,) in dB J(Azcu) in dB
0 2.47 2.48 24.7 24.6
1 1.60 1.58 33.8 33.7
2 2.31 2.31 49.9 49.9
3 1.92 1.92 60.7 60.7
Table 42: Performance function characteristics for example 2.
N Aopt Anzcu J(A,,) in dB J(RzcU) in dB
0 2.99 3.12 11.6 11.6
1 1.50 1.43 15.6 15.5
2 2.58 2.60 20.8 20.8
3 1.81 1.78 24.8 24.7
4 2.45 2.47 29.3 29.3
5 1.93 1.92 33.3 33.3
6 2.39 2.40 37.7 37.7
7 2.00 2.00 41.8 41.8
Table 43: Performance function characteristics for example 3.
N Ao, 2nzcu J( A,) in dB J(,,zcu) in dB
0 0.82 0.34 0.5 0.2
1 2.30 2.46 2.7 2.7
2 1.29 1.21 5.2 5.1
3 2.12 2.12 9.9 9.9
4 2.87 2.89 12.7 12.7
5 3.56 3.59 14.1 14.1
6 4.21 4.22 14.7 14.7
7 4.82 4.76 14.9 14.9
8 2.24 2.24 17.1 17.1
9 2.62 2.64 19.6 19.5
10 2.25 2.24 22.9 22.9
11 2.57 2.56 25.5 25.5
4.3.3 Further Analysis of the Complex Performance Function
In the previous section, it was shown that the squared error and complex
performance functions can be expressed as
J = d w w
Sd (4.89)
Ju = due w
where
d, = D(s)D(s)ds du, = D(s)D(s)ds
c c" (4.90)
w= fD(s)(s,A)ds w, = D(s)(s,A)ds
C C.
c c,
and
D(s) = desired system
cp(s,2) = L(s,/) the Laguerre kernel
C = contour enclosing the entire left hand plane (LHP) in s
Cu = contour enclosing the upper left hand plane (LHP) in s
It was also shown that
d, = 2Re[d,]
d = 2Re[de (4.91)
w= 2Re[w,]
Writing J and Ju in normalized form (relative to the total energy of the system) gives
J = (d, WTW)d (4.92)
(4.92)
Ju = (de wT'w)/d
Decomposing due and Wu into their complex component forms
due = dueR + jdue (4.93)
wu = wR + jW.u
yielding
de = 2Re[due] = 2dueR
w = 2RWu(4.94)
w= 2Re[w,] = 2w,
So we can write J and J, in terms of components of the complex energy signal and
complex weight vector.
J = (2due (2wR)T(2wUR))/(2duR)
= (dR 2w uR) dueR
Ju = (due (2wuR)T w)/(2d,R)
=(i ww)/U(4.96)
Now we are interested in 2Im[J,], twice the imaginary component of the complex
performance function. This has the expression
2Im[J,] = (due, 2wWU)/du, (4.97)
This extension to the complex domain has shed light on new information. The
imaginary component of Ju has been shown to have roots (or zero crossings) near the
stationary points of J, one of which provides an approximation to the optimal time
scale solution.
In an effort to capitalize on this discovery, the analytic performance function, Ja,
will be defined and studied as a method of achieving the same information while
requiring less a priori knowledge of the unknown system D(s). As will be seen, Ja is
composed of an energy signal and a weight vector that is the analytic signal
representation of the associated squared error components.
4.3.3.1 Analytic performance function. Consider the components of Ju. The
complex energy signal, due, can be written as
du = 2 D(s)D(s)ds
o0
= J D(o)D(c)da + f D( jw)D( ji)do (4.98)
0 0
~0
= [,D(jw)D(jw)dw] j D(o)D()d7]
Likewise, the complex weight vector, wu, can be written as
w = fD(s)(s,2 )ds
C.
o (4.99)
= f D(a)((c~,A)do + f D(jw)( jw, )dw
00 0
Now collecting real and imaginary components gives
wU= Re D( jo)W( jo,A)do Re [ f D(a)(c, )da +
0 1 1 (4.100)
j Im [i D( jo)(j, A)d J Im J D(o)O(a, ,)d
Consider the components of Wu along the +jm axis only. This is equivalent to the
spectrum of the analytic signal constructed from w. Namely, for a signal f(t), the
analytic signal is specified by
fa(t) f(t) + j(t) (4.101)
where f(t)= H{f(t)}, H is the Hilbert Transform. The spectrum of f,(t) is
equivalent to the spectrum of f(t) for 0 < < oo and zero for oo < m < 0. Now we
define the following functions:
d,,, [real part of d., along + jw axis] + j[imaginary part of du along + jo axis]
= D(joI)D(j))da) + [0] (4.102)
= f D(jo)D(jwo)dw
0
and
w,, =[real part of w, along + jw axis] + j[imaginary part of w, along + jw axis]
= Re 1 D( jo)0(jo, A)dj + j Im [ D( jw)(jo, l)do] (4.103)
0 0
= JD( j)w)( jo,)do
0
dae is called the analytic energy signal and wa is called the analytic weight vector. As
in the case of the complex performance function we define the analytic performance
function, Ja,
Ja d w'w (4.104)
or in normalized form
J= (da, wT a)/de (4.105)
Now decomposing the components of Ja into their complex form gives
dae = dae + Jde' (4.106)
a = WR + Jal
where
doeR = [ D(j))D(j)d (4.107)
dae =0
and
70
a = Re D(j)(iw4(jo9,A)d]o
(4.108)
w = Im [ D( j)W (jI),)dA)
0
The goal is to rewrite J and Ja in terms of analytic energy signal and analytic weight
vector components only. Recall
o J
w=; D(s)W(s, A)ds
C
I Tr f D(jw4)(jw,A2)dw
0
= D(jo)O( jo,)d(+ + D(jmo))( jo,A)do
0 0
= 2 iD( jo)D(jw )dw) + D( j))(jwj,A))dw (4.109)
0 1 0
=2Re[2 i D(jo)o(jtWA)dw]
=2Re[w]
= 2w,,
Likewise
d,= 2Re [ D(jow)D(jm)dw ]
0
= 2Re[de] (4.110)
= 2daeR
So we can write J and Ja in terms of components of dae and Wa as follows
J = (2de (2WaR)T(2waR))(2daeR)
S((4.111)
= (daeR 2w WaR W)/daeR
S=(dae (2WaRT a)/(2daeR)
(4.112)
We are interested in 2Im[Ja], twice the imaginary component of the analytic
performance function. This has expression
2Im[Ja] = 2wW waT /daeR (4.113)
4.3.3.2 Analytic performance function examples and experimental results. Again
consider the three examples used to study the complex performance function. Figures
47, 48, and 49 compare the complex and analytic performance functions. For each
of the three examples the squared error performance function, J (solid curve), twice
the imaginary component of the complex performance function, 2Im[Ju] (dash curve),
and twice the imaginary component of the analytic performance function, 2Im[Ja]
(dotdash curve), are plotted for the first four model orders. Numerical results are
provided in Tables 45, 46, and 47 to illustrate how close the nearest zero crossing
(nearest ZCu) of 2Im[Ju] and the nearest zero crossing (nearest ZCa) of 2Im[Ja] are to
the optimal time scale as the model order increases. Also the squared error is
computed at the optimal, J(A,,p), nearest ZCu, J(Anzcu), and nearest ZCa, J(nzcaa), time
scales.
It is interesting to compare Ju and Ja for these three examples. The difference in
these functions is that Ju includes the line integral along the negative aoaxis. Consider
the complex energy due for each example. From Table 44 we find that duel is 22% of
Table 44: Complex energy signal characteristics.
Example due= dueR + duei % complex
1 0.45+0.10/ 22
2 1.98+0.63/ 32
3 9.91+0.31/ 3
dueR for example 1, 32% for example 2, and only 3% for example 3. Observe from
Figure 49 that J,(A)) = J,(A), for all A, for example 3. Figures 47 and 48 show that
Ju and Ja are significantly different (but have close zero crossings near Aopt) for
examples 1 and 2. So the line integral along the negative o axis contributes
significantly to the evaluation of J for examples 1 and 2 but not example 3. This
results in the similarity of Ju and Ja for the last case. Clearly, the relative values of
dueR and duc, are application dependent.
We have shown that even though Ju and Ja are different functions, their imaginary
components have roots that are similar near 2op,. Since Ja can be computed directly
from a sampling of d(t), it could be used to provide the optimal time scale estimate in
lieu of Ju.
(order=0)
0.1
_U 0.05 
S 0 ... ., ......' ....... .........  ..... .............
E
0.1
0.1 I  j 
0 1 2 3 4 5
3
x10 (order=2)
1
0.5 ....
\
E i
I
I t i I
0 1 2
3 4 5
(order=1)
0.00
 0.00
0.C
i !
4 ..... ......... i ............ ....... ..................
5 
n5 
SI
0 1 2 3 4
0 1 2 3 4 5
S
Figure 47: J, 2Im[Jy], and 2Im[Ja] for example 1.
Table 45: Performance function characteristics for example 1.
N Aopt 2Azcu Azca J(A,,)p in dB J(Azcu) in dB J(Anzca) in dB
0 2.47 2.48 2.48 24.7 24.6 24.6
1 1.60 1.58 1.49 33.8 33.7 31.1
2 2.31 2.31 2.29 49.9 49.9 49.5
3 1.92 1.92 2.62 60.7 60.7 54.6
0.0
(order=1)
, ...... .............. ................................... ..... .....................
i/
. ...... .. .
........;
0.3
0.2
M0.1
4 0
E 0.1
S0.2
0.4
v 0.2
E
;i0
E
,0.2
0.4
0.1
7 0.05
E
 0
E
, 0.05
0.1
0 2 4
 0.05
E
" 0
E
i 0.05
(order=3)
(order=3)
2 4
Figure 48: J, 2Im[J,], and 2Im[Ja] for example 2.
Table 46: Performance function characteristics for example 2.
N Ap A,,zcu Azc J(Aopt) in dB J(AZzcu) in dB J(Azcn) in dB
0 2.99 3.12 3.08 11.6 11.6 11.6
1 1.50 1.43 1.31 15.6 15.5 14.8
2 2.58 2.60 2.19 20.8 20.8 19.0
3 1.81 1.78 3.16 24.8 24.7 22.5
4 2.45 2.47 2.40 29.3 29.3 29.2
5 1.93 1.92 1.86 33.3 33.3 32.8
6 2.39 2.40 2.14 37.7 37.7 34.8
7 2.00 2.00 1.96 41.8 41.8 41.5
0 2 4
(order=2)
(order=O)
(order=0)
i i
0 1 2 3
(order=2)
2 4
Figure 49: J, 2Im[J.], and
(order=1)
0.8
0.6
E
S0.4
E 0.2
0
0.2
0.8
0.6
_E
 0.4
E 0.2
0
0.2
0.6
" 0.4
E
M 0.2
4 0
0.2
2 4
orderr)
(order=3)
2 4
2Im[J] for example 3.
21m[Ja] for example 3.
Table 47: Performance function characteristics for example 3.
N Aopi 2zcu Anzc, J(2op,) in dB J(A zc,) in dB J(Anzca) in dB
0 0.82 0.34 0.55 0.5 0.2 0.4
1 2.30 2.46 2.43 2.7 2.7 2.7
2 1.29 1.21 1.21 5.2 5.1 5.0
3 2.12 2.12 2.13 9.9 9.9 9.9
4 2.87 2.89 2.89 12.7 12.7 12.7
5 3.56 3.59 3.60 14.1 14.1 14.1
6 4.21 4.22 4.22 14.7 14.7 14.7
7 4.82 4.76 4.79 14.9 14.9 14.9
8 2.24 2.24 2.22 17.1 17.1 17.0
9 2.62 2.64 2.60 19.6 19.5 19.6
10 2.25 2.24 2.28 22.9 22.9 22.8
11 2.57 2.56 2.62 25.5 25.5 25.3
0.6
" 0.4
E
4 0.2
E
 o0
0.2
dl
............... ..................... I ..............  ...........
4.3.4 Time Domain Synthesis Using GLFF Networks
For the discretetime setting, one of the fundamental limitations of an FIR model
is its lack of memory depth; namely, its inability to approximate (with low order) a
system with an infinitely decaying impulse response. GLFF networks do not have
this problem due to the presence of local feedback in their structure.
How well a Laguerre network approximates the impulse responses of the three
example systems given by equations (4.86), (4.87), and (4.88) is now illustrated.
Here the system impulse response, d(t), (h(t) = d(t) since the input is an impulse), the
Laguerre network output when using the optimal time scale yo(t), and the Laguerre
network output when using the nearest ZCu time scale yz(t) are plotted. In example 1
the impulse response decays rapidly and has no significant oscillation. yo(t) and yz(t)
are plotted in Figure 410 using a 1st order Laguerre model. In example 2 the impulse
response is also rapidly decaying but with a slight oscillation. Plots of yo(t) and yz(t)
are given in Figure 411 for a 5th order model. Example 3 is a system having a slowly
decaying impulse response with a strong oscillation. Plots of yo(t) and yz(t) are given
in Figure 412 for a 10th order model. This oscillation is difficult to approximate with
a Laguerre model since only a single real pole (time scale) is used in the
approximation. Increasing the order of the network could further reduce the error.
Another option is to use a Kautz model with a pair of complex conjugate poles.
For all three examples we observe that yo(t) and yz(t) are inseparable. The point
to be made in this section is that when using A obtained from the appropriate zero
Example 1: d(t),yo(t),yz(t)
1 2
3 4
Figure 410: Impulse response and approximations using a 1st order Laguerre model.
Example 2: d(t),yo(t),yz(t)
3 4 5
Figure 411: Impulse response and approximations using a 5th order Laguerre model.
Example 3: d(t),yo(t),yz(t)
5 6 7 8
Figure 412: Impulse response and approximations using a 10th order Laguerre model.
2
1.5
1
0.5
0
%I
crossing of Im[J1u (or Im[Ja]) as the time scale for the model, the same result can be
achieved as when using A obtained by searching J for the global minimum. Also the
squared error performance function can always be computed using J = 2Re[Ju].
4.3.5 Conclusion of the Complex Performance Function
For systems composed of sets of real and/or complex conjugate pair poles, a new
complex performance function, J,, was derived when using a GLFF network as a
system model. The imaginary component of Ju has zero crossings near the stationary
points of the squared error performance function J. Experimentally it is determined
that as the order of the model is increased, the zero crossings of Im[J,] near minima
of J converge to these stationary points, one of which will be the optimal time scale,
Aopt. Rather than searching the nonconvex function J for a global minimum, locate
the roots of Im[Ju(A)], A,, and evaluate 2Re[Ju(Ar)]. One of these will be close to Aopt;
namely,
= arg( min{2 Re[Ju (A,)]}j (4.114)
where
2r= {roots of Im[Ju(X)]}.
In addition, we have defined the new analytic performance function and derived its
closed form solution. It was shown that J could be written in terms of components of
Ju and Ja. Expressions for the imaginary components of Ju and Ja were also given.
The value of these functions in determining the optimal time scale of GLFF kernels
was illustrated using three example systems. The important equations from the above
sections are now summarized:
Squared Error Performance Function
(de wTw)/de
J= (deR 2WTRWuR,)/dueR
(daeR 2WRWaR )/daeR
{Squared error components}
{Complex error components}
{Analytic error components}
Complex Performance Function
Ju = (+due WuRWu)/dueR
Analytic Performance Function
Ja = ( dae Wa)/daeR
where
d, = ;D(s)D(s)ds
d,= D(s)D(s)ds
Cu
o
w= D(s(Ds(s,A)ds
C
w = D(s)4(s, A)ds
wa = D(ji)D( jw,A)dw
0
If D(s) is known analytically or in closed form we can employ the residue theorem
to compute de, due, dae, w, Wu, wa. If these functions must be determined numerically
then the analytic performance function can be easily computed. The only requirement
is a sampled data set from an impulse or step response test. In addition, J can always
be computed from Ju or Ja. This illustrates that the squared error solution is
completely embedded within the complex and analytic performance functions.
(4.115)
(4.116)
(4.117)
(4.118)
4.4 Conclusion (Optimal Parameter Selection)
This chapter has addressed optimization techniques for GLFF models, M(p).
Optimization of each element of the parameter set p = {N,w, A} is handled
separately. Model dimension element optimization generally involves increasing N
until a prescribed error tolerance is achieved. Monitoring the rate of convergence for
several GLFF kernels could also aid in choosing the model yielding the smallest N.
The optimal weight vector minimizes the integrated squared error and is found by
solving a set of WienerHopf equations. A direct integral evaluation method was also
given that takes advantage of properties of the model kernels. The weight vector
solution is a function of the time scale vector 2. Optimizing A is the most difficult
task because the minimum squared error solution is nonlinear w.r.t. these parameters.
A number of optimization approaches exist in the literature and they were discussed.
Finally, a new method was proposed that extends the squared error analysis to the
complex domain. Simply locating zero crossings and choosing the one that
minimizes J gives the optimal A. Examples were provided to demonstrate the theory.
CHAPTER 5
APPLICATIONS
Applications of the theory presented in the previous chapters will now be
discussed. In this thesis, GLFF networks will be used to approximate sampleddata
systems. In particular, the step response of a system used in laboratory test and
evaluation of guided weapon systems will be modeled. Other areas where GLFF
models have proven useful are first surveyed.
The GLFF model class encompasses many structures that have been
independently formulated and studied in the past. This includes the tapdelay line,
Gamma, Laguerre, Legendre, Jacobi, and Kautz models defined in both discretetime
and continuoustime domains. The unification of these structures is one contributory
component of this thesis; however, it is worth reviewing the previous history of these
models to gain insight into the many applications that they have been able to address.
Network synthesis. The earliest work dealt primarily with orthogonal model
structures in an effort to synthesize networks and transient responses of linear systems
[Arm57], [Arm59], [Kau54], [Hea55], [Bru64], [Clo65], [Kin77], [War53]. In
[Lee33] and [Hug56] synthesis of electrical circuits was the focus.
System identification and modeling. Similarly as the tools became more
developed the nomenclature and approach shifted to system identification and system
modeling [Bod94], [Eko94], [Kin79], [Lin93], [Mas90], [Mas91], [Nur87], [Nin94],
[Oli87], [Oli94], [Per91], [Si194], [Sil95], [Van94], [Wah91], [Wah94], [Moh88].
The basic principles are the same as in network synthesis only now the modeling
effort involves more than just the network impulse response. A further refinement
includes modeling noncausal systems [Eko95], nonlinear biological systems
[Mam93], probability density functions, [Pab92], [Tia89]. State space formulations
appeared [Dum86], [Gus93], [Heu90], [Heu95], and modeling using step response
data [Wan94b] and frequency response data [Clu92], [Wan95] have been considered.
Function approximation. A representation of an underlying system is not always
desired. Various GLFF models were also useful in function approximation, signal
representation, time series analysis [Bro65], [Cel95], [Con94], [Cle82], [Hwa82],
[Kor91], [Mar69], [Mai85], [McD68], [Ros64], [Sch71], [Wah93], [You62b],
[Yen62]. Particular examples include representation of seismic data [Dea64] and
electrocardiograms [You63].
Speech. Applications of GLFF models used in speech synthesis/analysis [Man63]
and speech compression [A1J94] have appeared.
Control theory. GLFF models are also useful in control theory. Examples include
adaptive controllers [Dum86], [Dum90], [Zer88a], adaptive predictor controllers
[Els94], [Guo94], optimal control [Hag91], [Nur81], robust control [Clu91] and
robust stability tests [Aga92].
Stochastic processes. Stochastic concepts have been addressed including
approximating cross correlations [den93a], correlations and power density spectra
[Lee67], and probability density functions [Tia89], [Pab92].
(Adaptive) filter theory. GLFF models have been used as filters [Kab94]
including adaptive filters [den93a], [den94], [Per88], [Pri93]. Systematic procedures
to design FIR filters have been introduced [Fri81], [Mas90], [Mas91]. Recursive
identification of time varying signals [Gun90], 2D filtering [Tsi94] and hardware
implementations have also been studied [Ji95].
Other applications. Other examples where GLFF networks have been used
include solutions to boundary value problems [Ira92], separation, enhancement, and
tracking multiple sinusoids [Ko94], and approximation of time delays [Lam94],
[Par91], [Moh95].
Although the list of applications above is fairly extensive, it is by no means
exhaustive. Nevertheless, the applicability of the GLFF modeling theory developed in
this research is clearly illustrated.
5.1 Examples of GLFF Networks Using Simulated Systems
The application focus of this thesis is system modeling; namely, it is desired to
approximate (model) LTI systems using GLFF networks. The objective is to obtain
an approximation to the system impulse response. Quite often this is accomplished
by optimizing the parameters of the model using step response data. In this section
several GLFF network realizations are constructed that approximate simulated LTI
systems. Optimization is performed using the step response of the desired systems.
Real (nonsimulated) systems and their associated step response data will be
considered in Chapters 5.3 and 5.4.
Three simulated systems will be studied: 1st order, 2nd order, and 4th order.
Example 1: (1st order)
H(s)= r=1 (5.1)
s + 1
Example 2: (2nd order)
(02
H(s) = w, = 1, = 0.4 (underdamped) (5.2)
s + 2 os + n.
Example 3: (4th order)
H(s)= H,(s)H2(s) (5.3)
where
H (s) = + 2, ,, =1, I = 0.4 (underdamped)
s2 + 241,00s + t (5
2 (5.4)
H,(s)= n2 2 ,2 = 2, 2 = 0.25 (underdamped)
s2 + 22W2s + on2
In each case, the step response is generated and sampled. The sampled step
response data is modeled with Laguerre, Gamma, and Legendre models. The results
are given below in Figures 51 through 56. In each case, the optimal time scale
parameter was determined and used to compute and plot the step response and
impulse response for 6tap (N = 5) networks. Figures 51, 52, and 53 compare the
system and GLFF network step responses. Figures 54, 55, and 56 show the
resulting impulse responses. The optimal time scale parameters of the 6tap models
and associated minimum squared error, J(Aopt), are provided in Tables 51, 52, 53.
1
0.9
0.8
0.7
0.6
0.5
0.4 System
+ + + + Laguerre model
0.3 f Gamma model
x x x x Legendre model
0.2
0.1 
0
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
t (sec)
Figure 51: 1st order step response and model data.
14
1.2
1
0.8
0.6
04 System
+ + + + Laguerre model
Gamma model
0.2 x x x x Legendre model
0
1 2 3 4 5
t (sec)
6 7 8 9 10
Figure 52: 2nd order step response and model data.
1.6
1.4
1.2
1 
0.8
0.6
0.4
0.2
0.2
0
1 2 3 4 5 6 7 8 9 10
t (sec)
Figure 53: 4h order step response and model data.
1 
0.9
0.8
0.7
0.6
0.5 
0.4
0.3
0.2
0.1
0
0 0.5
System
+ + + + Laguerre model
* Gamma model
x x x x Legendre model
Figure 54: 1st order impulse response and model data.
System
+ + + + Laguerre model
* Gamma model
x x x x Legendre model
1 1.5 2 2.5 3 3.5 4 4.5 5
t (sec)
0.8
0.7
xxxx
0.6  System
+ + + + Laguerre model
0.5 Gamma model
f x x x x Legendre model
0.4
x
0.3 x
X
0.2 x
x
0.1
0.1
0 2 ,i 1 i i I I i1
0 1 2 3 4 5 6 7 8 9 10
t (sec)
Figure 55: 2nd order impulse response and model data.
0.8 x
X+.
0.7 x .
x
0.6 x x System
0.5 + ++ + Laguerre model
0 x\ Gamma model
0.4 x \ x x x x Legendre model
x
0.3
0.2
0.1 +* / ;xxxxxxxxxxx>
0 x x +x
x Ixx
0.1 \.x +
+. x .. (*', x
0.2 \ *+*.x
0 1 2 3 4 5 6 7 8 9 10
t (sec)
Figure 56: 4th order impulse response and model data.
Table 51: 1st order model results.
Model N Aot J(2o) in dB
FIR 101 13.7
Laguerre 5 0.88 57.1
Gamma 5 0.05 oo
Legendre 5 Aop = e0.05 56.0
Table 52: 2nd order model results.
Model N A__ pt J(A opt) in dB
FIR 101 13.0
Laguerre 5 0.95 57.1
Gamma 5 0.08 40.7
Legendre 5 opt,i = e0.03i 37.5
Table 53: 4th order model results.
Model N 2opt J(/p,) in dB
FIR 101 10.9
Laguerre 5 0.89 37.7
Gamma 5 0.09 31.2
Legendre 5 opt,i = e003i 28.8
The three examples considered above are typical of the subsystem makeup of
guided weapon systems. With fairly low order models (N=5) the GLFF networks
have demonstrated (in the simulated examples) they can approximate these systems
well. This includes systems having slowly decaying impulse responses and
oscillatory modes. The most difficult task of the modeling process is locating the
optimal time scale. This was done using a numerical search technique in the above
examples. For a given model, if the resulting error is too large the model order can be
increased. This process can be continued until the prescribed error tolerance is
reached.
5.2 Laboratory Testing of Guided Weapon Systems
Guided weapon systems are as important in a military tactical mission as the
aircraft and ships that carry them. They are the tools with which a pilot can carry out
both offensive and defensive plans. As such, their capability and reliability must be
well understood so that they can be operated within the envelope for which they were
developed. To ensure their success it is necessary that the operational envelope is
understood. This knowledge is acquired through extensive test and evaluation of
these systems. The objective is to produce a smart dependable weapon that operates
as expected when launched from an aircraft or ship.
Test and evaluation (T&E) of guided weapon subsystems can be broken down
into three phases [Eic89]. This breakdown is typical when a new weapon or
modification of an existing weapon is under consideration. These phases are
1. Research/Concept Exploration
a) analytic study
b) openloop laboratory and field tests
2. Simulation
a) software simulation and modeling
b) HITL simulation
3. Flight Test
a) captive carry
b) launch
Research and concept exploration are the very beginning of the T&E effort.
Proposals must be carefully considered to determine which (current) weapons could
possibly meet the objective being considered. Can an existing weapon be modified or
does it call for a new design? All these scenarios must be analytically explored.
Flight tests are the final hurdle. It is here where the weapon system's performance is
measured against its expectation. Before flight testing however, many problems can
be worked out by evaluating a weapon system in a simulated environment created in
the laboratory. There are many methods for evaluating the performance of
components of a weapon system such as parametric tests and openloop tests.
Another class of tests studies the behavior of these systems in a closedloop
configuration and involves a simulation of the realworld environment and the
weapon's flight in this environment. This class of tests can be broken down into two
types: alldigital and hardwareintheloop (HITL). In an alldigital closedloop
flight, the entire weapon system is simulated mathematically. If the weapon
subsystem models are accurate, alldigital simulations can give us a good idea of how
the weapon will perform as well as what its limitations and vulnerabilities are. Once
subsystems have been hardware fabricated, these can replace the models in the
closedloop simulation. Ultimately it is desired to have as much of the weapon
hardware in the loop as possible.
5.2.1 Objectives
Several reasons why laboratory HITL testing is an integral part of a weapons
development test program are now given.
1. Asset Availability
During the development of a weapon system, often answers are needed to key
questions about the system to determine the limits of its capability. Rather
than fabricating an entire system and conducting a flight test, many times
these questions can be answered by using simple prototypes or only subsystem
components. This precludes having to completely fabricate a system and
conduct open air flight tests. Flight tests are expensive and consume lots of
resources (fully developed weapon system, aircraft/ships, openair range
time/space, money). Here is where laboratory testing offers a great advantage.
2. Planning and Analysis Requirements
As a missile system reaches development stage and can be flight tested, valid
test conditions (test matrix) must be formulated. Having laboratory results
(alldigital and/or HITL) as a guide, the test matrix can be formulated with
intelligence and unrealistic scenarios can be avoided. Laboratory tests can not
only help define the launch envelope of a weapon, but they can be used to
predict its performance when operating under valid conditions.
3. Repeatability
Another advantage to the laboratory environment is repeatability testing.
Once a weapon system simulation has been developed, it can be used to
conduct hundreds and thousands of tests. Weapon performance when
launched under varying release conditions and "what if" scenarios can be
examined. The laboratory environment is a major contributor to this kind of
testing because real assets do not have to be expended to conduct many tests.
In "live fire" testing assets can usually be used only once.
4. Security/Safety
Finally there is the security advantage. As long as testing is being conducted
in a laboratory, external security vulnerability is reduced. Range safety costs
and concerns that occur in openair testing are also eliminated.
To accomplish valid and accurate laboratory testing, careful attention must be
paid to the test configuration. For HITL tests it must be ensured that the laboratory
equipment can recreate the dynamics and conditions that will be experienced in the
real world. Often it is easy to switch between configurations involving hardware and
software. This requires accurate software models of the hardware being simulated.
This is where GLFF networks discussed in this work enters the application front.
To gain a better understanding of the systems the GLFF networks must model, a
more detailed breakdown of a typical missile system is now given. The functional
components of the missile are discussed as well as how they are integrated into a
HITL test configuration. The GLFF networks can then be used to approximate some
of these components.
5.2.2 Functional Breakdown of a Missile System
This section provides a discussion of components (subsystems) most often found
in missile systems. Figure 57 below illustrates a highlevel functional break down of
a typical missile system [Gar80].
Figure 57: Missile functional block diagram.
Seeker/Guidance Modules. It can be said that the seeker is the eye and the
guidance module is the brain of a weapon system. It is here the missile trajectory is
commanded and controlled. For missiles that are considered lockon after launch, the
guidance module may begin in a target search mode. In this case it drives the seeker
to search for a target in some preprogrammed pattern. During this phase commands
may be sent to the autopilot to have it execute a midcourse maneuver such as pitching
the weapon up in order to gain altitude for target search and future end game
engagement. If target energy is found the guidance module may transition into an
acquisition mode to narrow the search field of the seeker and attempt to select a single
target. At this point it begins target track and will try to keep the seeker pointed in the
direction of the target. During track mode, guidance (acceleration) commands are
sent to the autopilot in accordance with the guidance law being executed (proportional
or pursuit navigation). All the while this module performs seeker head stabilization
and directional pointing to maintain a continuous track of the target until impact.
Autopilot. The job of the autopilot is to maintain a stable flight and deliver on the
commands received from the guidance module. To perform its task, the autopilot
accepts input from various sources and sensors. The three most common are:
guidance module, accelerometers, gyros. The guidance module tells the autopilot that
a certain amount of acceleration is required to maintain the desired trajectory. Upon
receiving a command for a specific acceleration, the autopilot issues commands to the
actuator to drive the control surfaces (fins or canards) in the proper direction.
Measurements are taken from accelerometers to determine whether the weapon has
achieved the desired acceleration. Gyro rotational rate and angle measurements are
taken to ensure proper damping and to keep the motion of the missile in a stable
condition.
Actuator and Control Surfaces. The purpose of the actuator is to accept
commands from the autopilot and in turn drive the control surfaces to achieve the
dynamics needed. As the fins are controlled into the wind the missile body will rotate
and accelerate hopefully in a stable manner sufficient to maintain a target track.
Track angles, Steering commands Fin commands Fin
rates
rates positions
Target Seeker 4 Guidance Autopilot Actuator
signature
Mode, Stabilization
"> I Flight Motion Simulator (FMS) IAa
T arget .............................................. .......................... ... ............................ A ir ram e
Target A
Forces,
............... ....... K inem atics
Missile body
angles, rates,
Target Missile/Target 4 accelerations
Geometry
Target profile
Figure 58: Closedloop missile simulation.
5.2.3 HITL ClosedLoop Missile Simulations
In this section HITL closedloop missile simulations will be discussed. Figure 58
above illustrates a typical test configuration. Usually the seeker, guidance module,
and autopilot are hardware components and the functionality of these systems is the
subject of the test. Special laboratory equipment is used to generate the target scene.
This is dependent upon the kind of sensor in the seeker head. The remaining blocks
are simulated in software. For a HITL closedloop simulation, the seeker, guidance

Full Text 
xml version 1.0 encoding UTF8
REPORT xmlns http:www.fcla.edudlsmddaitss xmlns:xsi http:www.w3.org2001XMLSchemainstance xsi:schemaLocation http:www.fcla.edudlsmddaitssdaitssReport.xsd
INGEST IEID E6RLTDKRC_DJKTO5 INGEST_TIME 20130214T13:55:24Z PACKAGE AA00013547_00001
AGREEMENT_INFO ACCOUNT UF PROJECT UFDC
FILES
PAGE 1
SYSTEM ODEL G GE RALIZED LINE R FEEDFOR B y DO GL G JO T ORKS DIS ERT TIO PRESE T D TO THE GR D TE HOOL OF THE IVER ITY OF FLORID P RTIAL FILLME T OF THE R QUIREME T FOR THE D ORE OF DO TOR OF PHILOSOPHY IVER ITY OF FLORID 1999
PAGE 2
ACKNOWLEDGMENTS There are many people who made the completion of this work possible Fir t and foremost I would like to thank my research professor Dr. Jose C. Principe He h a taught me the valuable skills of rigor and completeness. Also his patience with my progress while I juggled both research and a busy career was greatly appreciated I would like to thank my supervisory committee members Dr. Yunmei Chen Dr. William W. Edmonson, Dr. John G. Harris, and Dr. Jian Li for their encouragement and direction Many thanks go to the management within the 46 Test Wing Eglin AFB FL who made it possible for me to attend the University of Florida through the longterm full time training program. This organization fosters an atmosphere that encourage academic achievement and it has continually provided me with the opportunity to further my education. There are no words to describe the wisdom guidance and support I ha e rece1 ed from my parents Their love to this day is endless. Above all I would like to give the highest recognition to my lo ing wife Donn a, for her sacrificing devotion throughout my education To my preciou children Lindsey and Mark you are the sun hine in my life. La tly I dedicate thi work to my Lord and Sa v ior J u Chri t. I can do all thin g thr u g h him who g i e me trength. 11
PAGE 3
TABLE OF CO TE CKNOWLEDGME TS ............................................................... .. ..... ...................... ii ABSTRACT ................................................................................................................. CHAPTER l PROBLEM FORMUL TIO .................................................................. ................. 1 1.1 Problem Definition ..................................................... ................. .. ..... ... ............. 1.2 Mea urement Objecti .. .. .................. .. ........................... ... ... .. .......... ....... ........... 3 1 .3 Mod l Cla ification ........................................................................................... 4 1.4 Hi tori cal R vi w ................................................................... ............. ..... ........... 5 1.5 Modeling M th d log ........................................................................................ 7 1.6 Re earch ontributi n .......................................................................................... 8 ................................................................. 10 2.1 Linear Tim In riant au al t m ............................................................ 10 2.2 Realizability .......................................................................... ... .... .. ..... ... ............ 14 2.3 L 2 ignal and y t m ...................................................................................... 14 2.4 Sy tern ynth i ...... ................. ..... ................................................................. 15 2.5 ontinu u /Di er t Tim Tran f rm ti n ............... .... ... .. .... .. .... ...... .... .......... 16 2.6 Ex mpl ............................................................................................................ 18 2.7 Conclu i n ( umpti n ................................. ................. .... ........... ... 21 3 MOD LFOR TIO ...................................................... .. .... .. ..... .. ...... ........... 22 3.1 n ralM d I rmM(p ................................ ..... .......... ......... .... ... ...... ..... ....... 22 3.2Di er t Tim Md I ........................................................................................ 23 3.2.l FIRM d l ... .............. ................................ .............. .. ..... .... ... .... ........ ......... 23 3.2.2 IIR M d 1 ..................................................................................................... 24 3. 2. 3 LFF M d I ................................................................................................ 24 3.3.l P pul r Di er t Tim GLFF K m I .......... ....... .............. ......................... 27 3.3.2 P pular ntinu u ime GLFF Kernel .................................................... 27 3.4 GLFF K m I D ri ati n and Pr pertie ........... ...... ..... .. .......... .... .... ................. 30 3.4.l ntinu u Tim Lagu rre Kem l. ................ ......................... .... .... ............ 30 3.4.2 ntinu u Tim Kautz Kem 1 (Complex Time cale V ct r ................. 30 3.4.3 ontinu u Tim K utz Kernel (Real Tim cale Vect r ......................... 31 3.4.4 ontinuou Time Legendre and Jacobi Kem I ... ............. .. ....................... 32 3.5 lmpul e and t p Re pon Modeling ............................................................... 32 111
PAGE 4
3.6 Conclu ion (Mo del Formulation) .............. .. . ... ...... ..... ............. . ......... ..... .... ..... 33 4 OPTIMAL PARAMETER SELECTION ....... ..... ... ............ ... ..... .. ............ . .............. 35 4. 1 Model Dimension ............................................................................................... 35 4.1.1 Completeness ............................................................................................... 36 4.1.2 Orthogonal Projection ... ............................ ... .. .. .. ...... ..... .......................... .. . 37 4.1.3 '' Optim a l '' Selection of N .. ....... ... .............................................................. 38 4.2 Weight Vector ......................... ... ... ................... .......... .. ........... ... ...................... 39 4.2.1 Optimal Selection of w ..... .. ................. ... .. ...... ...... ....................... .... ..... .... 39 4.2.2 Relationship to WienerHopf Equations .. ......................... . ........... .............. 41 4.2.3 Minimum Squared Error Performance Function ......................................... 42 4.2.4 Computing w Using Step Response Data .................................................... 43 4.2.5 Adaptation of w ... . ........ .. ...... ...... .......... ..... ................................ ............... 46 4.3 Time Scale Vector ........ ..... ... ........ ........ ..... .... ............. ....... ..................... .... .. .... 46 4 .3. l Optimal Selection of l (Historical Review) ................................................ 48 4.3.2 Optimal Selection of l (Complex Error Analysis ) .......... ... .. .. ... . ....... .. .... 56 4.3.3 Further Analysis of the Complex Performance Function .... ....... ................ 66 4 .3 .4 Time Domain Synthesis Using GLFF Networks ......................................... 76 4.3.5 Conclusion of the Complex Performance Function ..................................... 78 4.4 Conclusion (Optimal Parameter Selection) . ..... ............................ ................... 80 5 APPLICATIONS ........ .. ........ .... ..... ....... .. ...... ....... .... ..................................... ..... ... .. .. 81 5 1 Examples of GLFF Networks Using Simulated Systems .................................. 83 5.2 Laboratory Testing of Guided Weapon Systems ................................ ...... ..... .. .. 89 5 .2. l Objectives ....................................... .... ...... ........... . . . ..... ... ........................ 91 5 .2.2 Functional Breakdown of a Missile System .. .............................................. 92 5 .2.3 HITL ClosedLoop Missile Simulations .. ............ . .. ........ .. .. .... ....... ...... .. 94 5 .3 Modeling the Step Response of a Carco 5Axis FMS ............................ .. ...... .. 95 5 3 l DiscreteTime GLFF Models of the FMS ........... . ... . ... ..... ... ... .. .. ...... .... .. 98 5 .3.2 Step Response Modeling Procedure ... .. .......... .................. .................... .. 100 5 .3.3 Results (FMS Step Response Modeling) . ... ................................. ...... .. ... . 103 5 .3 .4 ContinuousTime GLFF Models of the FMS ............................................ 118 5.3.5 Conclusion (Step Reponse Modeling) .. ....... ...... ................... ........ .. ......... 1 23 5.4 Other Applications ......... . ...... ... ..... ..... ..... .......... .. ... .. ...... .. ........................... 12 3 5.5 Conclusion (Applications) ................ ......... . .... ............. .... ......... .................... 124 APPENDICES A: ORTHOGONAL SIGNAL SETS .. ......... .......... .. ................................................. 125 B : Iill.,BERT SPACE CONCEPTS .................... .... .......... ............. ........... . .............. 126 C: SDOMAIN DERIVATION OF MINWUM J .............. .. ......... .. ... ... .................. 128 D : SQUARED ERROR STATIONARITY USING THE LAGUERRE MODEL .... 130 REFERENCES ........ .. .. .. . ...... ... . .... ..... ... ...... ... .......... ... ...................... . ... ................. 132 BIOGRAPIDCAL SKETCH .. .... ... ................ .. ........ ....... ............ .. .. ..... ... .......... ..... 140
PAGE 5
Ab tract of Di ertation Pre ented to th Graduate S hool of the ni er it of Florida in Partial Fulfil Im nt f th Requirement for the Degre of Do tor of Phil ph GE R Chairman: Jo e C Prin tp ODEL G R FEEDFOR B D ugl J n ugu t 1 M j r Departm nt : le tri al nd mput r n 0 in nn 0 The f u f thi rk i th T ORK nt f a la call d Gen raliz d Linear F d In I ud d i a d ri pti n f th a pri ri f rmulati n of th m rk ) t b u d and a d finiti n fa m a urement crit ria u d to judg g ffi t. Th m t d fr m a w ight d um f parametrically d p nd nt k m I The k m repr nt b ector that pan a vect r pac with dimen i n qual t th order f th netw rk r alizati n. Optimizati n f th m d I paramet r i di cu d and a new approach i gi n one that extend the c nventi nal quared err r analy i to the c mplex d main After pre entation f GLFF model theory and parametri ptimizati n methodologie n twork realizations are con tructed to addre the y tern modeling pr blem GLFF network can be formulated t mod I di creteti m ntinu u V
PAGE 6
tim and ampled data y tern In thi tudy it i hown how a flight motion imulator u d in th laboratory te t and evaluation of guided weapon sy terns can b a curately mod led u ing GLFF networks compri ed of Gamma Laguerre and L gendre k rnels. Ad antage of GLFF networks are their ability to accurately repre ent linear time invariant systems coupled with the simplicity of the procedures employed to optimally determine their parameters. The network structures have efficient and practical implementations and applications are offered as proof.
PAGE 7
CHAPTER l PROBLE FORMUL TIO In thi chapter an o er 1e of th entir cont nt of th re earch i g1 n beginning with a definition f th probl m to be ol d; name! m d lino th 0 impul e and tep re pon of au al, lin ar tim in ariant tern. Th sy tern/model confiour ti n I illu tr t d and th pnm r m a ur m nt bj ti tated. The model ar hit tur brief! d rib d and pr f major contributor I din 0 t th thi re ar h i umrnrui z d It inall lh n f ut th r matntn hapt r It i 11 b n that th m del d ith r I and/ r mpl njugate pair e n th hi h hav d a ing impul e large time nd/ r illat r m d (light damping) Thi i di fi ult t a hi with finit impul Th mod I tudied are n idered param e rri ince th ar fun ti n f param t r tor that mu t be determin d. The e param l r ect r ar fairl ptimiz and unlike infinite impul e re pon e (IlR) tru tur guaranteed tabilit can al a b The net orks und r c n id ration in thi re ear h are mpri d f linear combination of special basis function Th e ba i functi n are n t arbitrarily cho en but are derived fr m as umpti n f the y tern b ing m d I d. h
PAGE 8
2 netw rk th refore become good approximates since they obtain the natural information of the unknown system. This information resides in a parameter vector built into the model structure. The parameters are optimally chosen so that the model approximate the system in a least squares sense. Hence, the propo se d theory an implementation of parametric l east squa r es e rror modeling [Cad90]. The modeling methodology can also be described usmg geometric concepts. Considering signals as vectors in vector spaces, the models are linear combination of parametric basis vectors. They span a subspace of the space spanned by the unknown system. The model output is a projection of the system vector onto the model subspace. This space has an orientation that is parametrically established. The error is a vector that lies outside the model subspace. Projecting orthogonally the system onto the model space minimizes the error. 1.1 Problem Definition The objective is to approximate a causal, linear timeinvariant ( LTI) system, H using a predefined model M(p ), such that a measure of the quality of the approximation is minimized. Models developed in this study will be called generalized lin ear feedforward (GLFF) networks and will be defined in Chapter 3. The basic system modeling configuration is shown in Figure 11 belo Argument are not included with the function {x, d, y, e, H} to illu trate that the fo1mulation ill app l y in both a continuoustime setting {x(t) d(t) (t) e(t), H( ) } and a amp l eddata or di er tetim setting {x(k), d(k), y(k), e(k) H( ) } Quite often in thi the i witch backandforth between continuou time and discretetime cone pt Thi e
PAGE 9
3 intentionally done to prove the broad applicability of this work. Most of the theory ho w ever will be deri ed in continuou time. It would ha e been ju ta appropriate to develop the concepts in di cretetime. d H X ~ M(p) V Figure 11: y tern rn d ling How i it determin d wh th r m d I ha a hi d a d r pr nt ti n f th s y tern being appr irnat d. r thi a m a ur m nt th IT r r a t t f r goodne offit f th m d I I n d d In th int rat d quar d IT r crit ria will be mpl ycd ; narn ly f r a g i n m d I M(p) th param t r t p f und uch that th rgy in the IT r nal ? i minimi z d. hi i f th imple t IT r criteria t w rk with in e the mea urem nt i ba i all quadrati in f rm and it derivative are Ii near in e11ai n I m n t f th paramet r t. th objective i t m1mm1z J = fo 00 e\ t)dt = [(d(t)y( t)) 2 dt (co ntinuou tim) or (1.1) 00 00 J = L \ k) = L ( d(k) y( k) )1 (di c rete t ime ) k=O k=O
PAGE 10
4 To provide motivation for the GLFF networks (a complete formulation will be given in Chapter 3) con icier the discretetime etting. In the past discretetime y tern models have been largely divided into two classes: FIR and IIR. A new cla s of models (filters) known as the g e n e rali ze d f ee dforward (GFF) filters was popula1ized with the deve l opment of the Gamma filter [Pri93]. This filter has a structure similar to the tapdelay line (FIR filter) with the delays replaced by a more general transfer function known as the Gamma kernel. Extending this development can be accomplished by considering other popular kernels, most of which are orthogonal such as the Laguerre, Legendre, Jacobi, and Kautz [Lee67], [Arm59], [Kau54]. The extended formu l ation with a more comprehensive set of kernels is called the GLFF model class because its members are used in the context of modeling linear systems and they obey the superposition principle. How does the GLFF model class compare to FIR and IIR models? Figure 12 relates these three classes based upon important issues that are usually considered when deciding which modeling tool to pull from the signal processing toolbox. GLFF models can be said to fall omewhere between FIR and IIR models in capability complexity tability and computability. FIR GLFF IIR low { complexity} high capability tri { tability } difficult computability Figur 12: FIR/GLFF/IIR ompan n
PAGE 11
5 Thi the i e ploit the trength of GLFF model The offer an altemati e to the t o extreme in the modeling proce Their ad antage ar that the ha e a reduced comple it compared to IIR model hile off ring greater apabilit o er FIR models. nlike IIR modeling te hnique tabilit can be ea il guarant ed and their optimal or n aroptimal lution n be computed \' ith mu h le diffi ult The mo t diffi ult ta k i rn el cti n of the tim e cale i e tor X Thi i on of the parameter et p Th complete p ram ter et ill b d fin d in h pt r Optimization of it lement c mpri e II f hapter 4 Much f th the ry n hi h th f und d b gan v ith th work f Lee [Le 33] [ 7] ynth iz tran i nt netw rk and h thi al de i that can r aliz ph i al t m n th r maj r ntri but r 1 Kaut [Kau54] H pr n nt p th ti m d matn r r y tern u mg a m t g ner I rth g nal ignal et. unc, and Huggi n nd d Kautz' pr edure t di creteti me [Y u 2a] a did Br m [Br ]. Yaak ur81] formulat d tat pa e equati n that r pr ent di retetim lin ar y tern u ing the Laguerre m d I. R bert and Mulli pr ide a g n raliz d tat variable f rmulati n f rth g nal all pa net rk tru tur [R b87] h tructure ke p data path I al. The fundamental pr c ing element are alik except p rhap f r local caling param ter Thi implifi the de ign and all w f r di tributing the computation in a conn cted array. The rth g nal tru tur al o
PAGE 12
6 h e the advantage that higher order model can be con tructed by s imply adding new block or term without effecting or forcing a recalculation of the existing lower order t rm The concept of a generalized feedforward filter was propo ed by Principe [Pri93] with the Gamma model. This idea of creating a more general transversal filter u ing the earlier work with orthogonal ignal sets leads to the modeling methodolo gy h ere described ; namely the generalized linear feedforward network s. What is the basic objective? Given a desired impulse response h(t) of a ca u a l LTI system along with a prescribed tolerance on the allowable error, find a n approximating impulse response h(t) that is simple in structure and easy to compute. Let the system and model transfer functions be H(s) and H (s) respecti el Construct H(s) with the form m fl (szJ n H(s)= D(s) = A i;i = 2, A X(s) fl (spJ i= l (spi) ( 1.2 ) 1=1 Assuming no coincident poles the impulse response is n h(t) = r 1 ( H(s)) = I, AeP (1.3) 1 = where Re(p ;) 0 V i=l ... n. It is desired that h(t) and h(t) be as c]o ea po ible under prescribed tolerances. There are two wellknown approache to performing thi appro imation. Th e first to tran form h ( t) into H( ) and then approximate H( ) in the domain b u ing a ratio of polynomials. This approach i generally le difficult and e era! t chnique f r acco mpli hing thi appro im at ion h a e be n pre ented in the literature Th e
PAGE 13
7 include Pron method Pade appro imant continuedfraction e pan 10n and e en analog filter de ign [Su7 l], [Bak81] [Opp89]. The problem ith this approa h th t th error in the time dom in i diffi ult to control. The e ond approach i to approximate h(t) with a finite um of ponentiall damped inu oid and damp d e ponential (a de cribed ab ) and then tran forrn th appro imati n to find H( ). Here th modeling error ma b c ntr II d dir ct! in th tim d main. On a of accompli hing thi i to impl m nt th pr dure d crib db K utz [Kau 4] hich involve appr imating h(t) u ing a ighted um ba i function are th m el eighted um of that they ar orth n rmal. Thi pr id f r gr at impli fi ati n ar hin 0 f r the optimal param t r et p nd gi e di r ntrol r th a iat d rr r. 1. In thi th th meth d u d imp! ment net rk a ppr i mati n f lin ar y l m ba don th n b Kautz hi arri d u l b c n id ri n g a c n tructi n f Ii t u ing an in init um f eighted tim d main ba i functi n {Jt A)f = o inc th infinit um mu t be truncated t a finit number ft rm an appr imati n re ult h I = /;(!) = L w n 11 (t A (1.4) n = O The w ight wo, w 1 , w and ba i function /t A)L = o are ch en ith th following g al in mind: 1 F r a fi ed numb r f terrn N mak th c n rg nee a rapid a p ible. Thi i c ntro ll ed u ing a tim ca l e v ctor X
PAGE 14
8 2. Make the optimal weight independent of the number of taps N+l. Thi i accompli hed by using orthogonal basi functions. 3. For a fixed number of tap N+l minimize the energy in the error signal. Thi requires determining optimal values for all of the elements of the parameter et pin such a way that J (the integrated squared error) is minimized. The above procedure is applicable when modeling the impulse response of a desired system. It may also be of interest to approximate the step response. In fact the applications given in thi thesis involve finding system models using step response data. If the network input is a step function this method is employed by approximating not the desired system output but its derivative or by approximating the desired output and taking the derivative of the result. This will give an analytic expression for the impulse response When only a numerical or tabulated form of the impulse response is desired, use the step response data to generate a system model. Once the model is established use an impulse (or white noise) input to generate the impulse response data. 1.6 Research Contribution This study extends the above concepts and contributes to the system modeling theory in everal ways. Fir t, a unified system modeling framework centered around the GLFF model cla pro ided. As will be seen in Chapter 3 this f01mulation not only include model that employ orthogonal ignal et like the Laguerre and Kautz but non orthogonal one a well uch as the Gamma In addition to the e generalized models th ba ic fram ork an al o be u ed to repre ent FIR and IIR tructures.
PAGE 15
9 Optimal parameter selection is a big part of the modeling procedure. The most difficult task is locating or approximating the optimal time scale vector. A new approach that extend traditional squared error analy i to the comple domain is developed. It wiJl be hown how the information in thi new pace can be used to locate the optimal time cale parameter Finally an application of the theory i gi en by mod ling th tep re pon e f a Carco 5axi flight moti n imulator loop (HITL) testing of guided weapon S Thi u d in hardwareinthet m GLFF mod ling of oth r y t m and sub ystem u ed in HITL te ting i al o n id red. It i d m n trated h thi la of model can be employed to accurate! r pr ent the d namic of phy i al t m and add value to the ta k of weapon y tern t t and e aluati n.
PAGE 16
CHAPTER 2 SYSTEM ASSUMPTIONS In the previous chapter, a brief discussion was g1 ven concemmg the assumed systems to be modeled. A more complete analysis of all the system assumptions is now offered. Mathematical and physical arguments are used to validate their necessity. The basic assumptions are that input/output signals have finite energy and systems are stable, causal linear timeinvariant, and realizable. After stating the assumptions, the technique employed to synthesis these systems by simulation is given. Higher order systems will be simulated by cascading together 1 t and 2 nd order subsystems. This is done in a continuoustime framework. The transformation used when a digital system is desired is also discussed. Finally some simulation example of the kinds of real systems to be modeled are offered. 2.1 Linear, TimeInvariant, Causal Systems There are many ways to classify systems linear vs. nonlinear, timeinvariant v time varying causal vs. noncausal, stochastic v deterministic continuous discrete etc One important step in the tudy and analysis of systems theory is to derive a tructural or mathematical description of proce se and system This a urned structure is es ential if plau ible analytical models are ought to appro imate them On could argue that if we had a mathematical de cription of the process why would w ne d additional mod 1 The primary r a on for the model is to 10
PAGE 17
11 appro imate input/output beha ior of the proce The a urned mathematical description of the proces could be u ed to gi e ery pecific information about it internal composition (time con tant damping ratio et howe er a component le el pecification ma be una ailable. In fact e en if the e tern parameter ar known they may not be operating ithin their tated pecification. In thi re earch an under tanding of the internal cornpo ition i not the objecti e on] a model of th input/output behavior hen ie ing the tern a a bla k b The theor pre ented in thi ection doe ho e er g1 e the foundation u ed t formul te and alidat thi new cla of model and er illu trat th d ri ati n f th pt Many y tern tudi d in ngin nn n b ent d b math mati al relati n that d rib th ir ph i al mak up r b ha i r. Quit ft n tm r and a r ali ti hara n ati n ma r quir n nlin ar and/ r tim ar in quati n that ar di ffi ull t r impli it pra ti alit and in rd r la Ii h an ppli abl hi h an b u d in the d I n and an I r y t m appr imati n and urnpti n ar mad i al natur Thi f r th u f linear t m th r Thi r a mng 1 ft n ju lifi din n oftw way : 1. The y tern i a um d t be linear by nature r it at lea t p rate in a Jin ar re g 1 n fit capabilitie that th Jin aiit nditi n ar ati fied 2 The y tern i a urned t be n nlinear b natur and in rd r t appl linear y tern theory c nceptualize a lineari z ati n ar und a n rninal op rating point. With thi appr ach precauti n mu t b taken t en ure th domain f the linear a umpti n i n t e eeded. Why i the u e of linear y tern the ry de irabl ? In additi n t the a ailability f a wel I tudied et of t m t ften a de cri pti n f a y tern via the impul e
PAGE 18
12 response and tran fer function is de si red The se functional system de scri ption are e tablished by employing the theory of linear systems [DeC89]. Consider a linear timeinvariant (LTI) system having input x(t) a nd output d(t). The impul e response h(t), is defined as the output when the input i s a unit impul se function !x_t). The transfer function H(s ), is defined as the Laplace tran sform of h(t) with zero initial conditions. This is equivalent to assuming h(t) = 0 fort< 0 (causal). The above definition of the transfer function is in terms of the impul se re s pon se. However a description of the transfer function taking a different more mathem a ti cal a pproach is also available. In practice a large class of signals a nd systems in engineering can be formulated through the use of differential equations The LTI causal assumptions equivalently translate into the requirement th at systems be representable by linear constant coefficient ordinary differential equations. Namely th d b a n n or er system can e wntten as d < 11 \ t) + a 11 ,d < n l\t) + + a,d < ) (t) + a 0 d(t) = b x ( m ) (t) + b x (ml ) (t) + + b x ( ) (t) + b x(t) m m l I 0 (2. 1 ) where ao, a 1 , a 11 , bo, b 1 , bm are all real constants with n m. Sy ste m s d e cribed b y ordinary differential equations of this form are said to contain lump ed para,neters rather than distributed parameter [Kuo87] Partial differenti a l equations s u c h a the elliptic equations hyperbolic equations, and parabolic equations de scri b e phenomena s uch a fluid flow wave propagation a nd heat conduction a ll of w hich are co n idered to b e di tributed p ara m ter systems [Pow79].
PAGE 19
13 To obtain the transfer function of an LTI system described by the differential equation abo e assume zero initial conditions (causality) and take the Laplace transform (.2) Rewriting this in factored form gi e m IT (s7; ) H(s) = A' 1 ' ; 1 (. ) rr (sP ) 1=! and in partial fraction form (a urning no repeated p n H( )= I, 1=1 ( P ) (_.4) To obtain the impulse re pon tak th m r e Lapla tran f rm umm 0 n coincid nt poles this yield 11 li(r = L eP 1 th LTI au al a umpti n I ad t th f II tn r ult : LTI au al t m can be repr c n tant olution (impul e re p n e ) ar ffi ient ordinary diff r ntial quation who e mpn d of a um f xp n ntial function that are r al (when P i i real) and/ r comp! ( hen P i i c mple ). If ome of the pole ar repeated then there are term with multiplicative factor t t2 t 3, ... l1 for a k th order pol Thi repre ntation i the motivation for the modeling structure considered in this thesis.
PAGE 20
14 2.2 Re a lizability It was hown that the LTI, cau al system assumption leads to a transfer function representation of the form m TI (sz ,) H(s) = Ai ; 1 TI (sP ) i = l (2.6 ) An answer is now given to the following question What are the constraints on this form that ensures it is a description of a physically r e ali z abl e system? The realizability requirement of H(s) is satisfied with the following constraints: 1. The constant A must be a real. 2. The number of poles, n must be greater than the number of zeros, m 3. Complex poles and zeros must occur in conjugate form That is if P i is complex then there must be a pj for some j such that p 1 = p ; where denotes complex conjugation. When H(s) is realizable it is easy to show that the impulse response is represented by a weighted sum of exponentials and exponentially decaying sinusoids. 2.3 L 2 Signals and Systems The final constraint imposed is that the signals must belong to the vector pace L 2(~ ) the squareintegrable functions. Here int e grabl e is meant in the L e b es gu e s ense [Deb90]. Assume the signal x E L 2(~ ) and x : H then 00 f x\ t)dt < 00 ( 2.7 ) Hence it i aid that x ha finit e e n e r y Let h(t) be the impul e response of a sy tern. Since h(t ) = 0 fort< 0 and hE L 2(~ ) then
PAGE 21
15 00 f h\t)dt < oo (2.8) 0 In thi ca e it i aid that hi a stable stem. Sob operating in L ( ~ )are tri tion i impo ed to work only with finite energ ignal and table tern 2.4 Sy tern Synthe i Recall the ystem realizabilit requir m nt r ult in n impul r pon representation that i a eighted um of ponential and p n ntiall d a rng inu oid From the tr n f r fun tion p r p ti th n a t m an b th ught f a con i ting of a ca cad of l I and ... nd rd r ub t m Th m t g n ral f rm f the e ub y tern tru tur are b lO + a lO H j ) = DX :? ( =b_ ::?l _+_b_ ::?O_ 2 + a 11 + a w (2.9) re ulting in diff rential quation . (2 10) d 2 (t)+a 11 d 1 (t)+a w d(t) = b 21 t t) + b 10 ~ t) During thi th xample ill b gi en For uch xampl th ub y tern bl ck M t ften light! I ill C imulated y t m ar creat d. f the e l t and 2 nd ord r g neral f rm will b u d t keep th ynthesi more in line with engineering (rath r than math m tical) concept The e pecial form conv y the engin ring inf rmation of u ual int r t: time con tants damping ratio and undamp d natural frequencies The e sub ystems have t h e following tran fer function
PAGE 22
16 H1(s)= D1(s) =1= _K_ X(s) is+l s+ H (s) = D 2 (s) = m 2 X(s) s2 +2c;m ns +m and corre ponding differential equations where id 1 (t) + d 1 (t) = x(t) d 2 (t) + 2c;m n d 2 (t) + m d(t) = (L} ~x (t) i= time constant c; = damping ratio {J), 1 = undamped natural frequency 2.5 Continuous/DiscreteTime Transformation (2 .11 ) (2.12) The goal of this thesis is to define a model class useful for approximating discrete time, continuoustime, and sampleddata systems. The above system formulation was carried out in the continuoustime domain. The same principles however could h ave been derived in the discretetime domain. All the concepts discussed have equivalent discretetime representations. Sometimes a continuoustime system in analytic form is given when a discretetime representation is desired This can be resolved by sampling. In this case, the discretetime system is called a ampleddata system Henc a sampling time T will necessarily be specified. Having a continuoustime transfer function, there are several methods used to map to the discretetime domain The most popular are the bilinear impulseor step invariant a nd matched Z tran form [Str88]. There are tradeoffs among these mapping a nd some work well where other do not. The most widely used i the bilinear tran iform. The tran formation (or mapping) i s ca rried out as follows
PAGE 23
17 H )=H(s) l s =2:1 2.13) T : + I here T i the ampling inter al. Thi i al o kno n a Tu tin method in control y tern circle and i deri ed by con idering the problem of numerical integration of differential equation [Fra90]. deri ation of the bilinear tran form i a follo Beginning ith the tran fer function H( determjn it differential equation repre entation. difference equation i then deri ed ho olution match the olution of the differential equation at the ample point The pecifi tran formation depend on the m thod ith whi h the diff r ntial equ tion oluti n I omputed Con ider the l t order y tern D ( ==_,14 x ) + or a ciated diff r ntial quation (_, l 5) The olution d 1 t) can b written in int gral form a I d, t) = 7 f[ d, ( u + x u ]du (2.16) 0 For th ampJing interval of time T thi can be r ritten a nT T nT di(n)=d,(nT)= ~ f[di(u + (u)]du +7 f[ di(u + x (u ]du 0 nT T nT 2.17) = d,(n1)+7 f[ d,(u)+x(v)]du nT T The method in which the integral computed o er th region [nTT nT] determines the di cretetime repre entation, d(n), of d(t) and hence th tran formation used to map H( ) to H( z ). Two of the most common integration sch mes are Euler
PAGE 24
18 rul e (backward difference method) and the trap eza idal rul e For these two ca s e s the following mappings result. Table 21 : Continuoustime to discretetime mappings. Integration Rule Mapping H(s) H H( z ) Euler s s=~( z~ l) Trapezoidal s= :( : :~) So the bilinear transform is derived by computing the integral o v er [nTT n11 u s ing the trapezoidal rule. In Steiglitz s work [Ste65] a more theoretical treatment and significance of the bilinear transform is developed. In fact it is shown that the bilinear transform is an isomorphism connecting the continuoustime vector space L 2(~ ) and the discretetime ( or digital) vector space of doubleended sequences of complex numbers { xn }: == that are square summable. This connection proves that the signal processing theorie s w ith respect to LTI causal signals and systems are identical in both the continuoustime domain and discretetime domain. In this thesis when a continuoustodiscrete t ra n s formation is performed the bilinear transform will always be used 2.6 Examples Several examples of the l s t and 2 nd order systems discussed in Chapter 2 .4 are s hown on the next two pages in Figures 2 1 to 24. These figures illustrate the impul se and st e p r es ponse behavior that must be appro x imated by the clas s of models propo se d in thi s the s is.
PAGE 25
1 8 1 6 1 4 1 2 08 06 0 4 02 09 08 0 7 06 ~05 04 0.3 02 0 1 19 05 1 5 2 25 3 35 4 45 5 t (sec) Figure 21: l I order impul er pon es, z=0 5 l 2 3. r=O 5 0 5 1 5 2 25 3 3 5 4 45 5 I (sec) Figure 22: l I order step re ponses z=0 5 1 2, 3.
PAGE 26
20 0 7 0 6 0 5 0.4 0 3 = ..c 0 2 0 1 0 0 1 0 2 0 2 3 4 5 6 7 8 9 10 t (sec) Figure 23: 2 nd order impulse responses with [4 1 =1. a) ~0.4 under damped ; b) ~1.0 critically damped; c) ~1.5 over damped. 1 2 0 8 0 6 0 4 0 2 0 0 2 0 2 3 4 5 6 7 8 9 10 t (sec) Figure 24: 2 nd order step responses with ~=1. a) ~0.4, under damped ; b) ~1.0 critically damped ; c) ~1.5 over damped.
PAGE 27
21 2.7 Conclu ion (System Assumptions) In thi chapter the as umption made about the systems to be modeled have been stated; causality, linearity, timeinvariance realizability. This discussion wa to provide the mathematical foundation and engineering principle that moti ate the development of the y tern model to be introduced. For the purpo e of imulating real y terns 1 t and 2 nd order ub y tern are u ed a the ba ic building block The e ub ystem can be ca caded together to ere te high order tern Method ere di cussed that allow for tran forming continuou time y tern de cripti n H( ) into di cretetime de cription H( A good m del i one with a tructur th t can captur th dynami of a y t m (time con tant damping co ffici nt natur 1 fr qu nci ) ithout e plicit knowledge of it par meter F r y tern with larg time con tant ( 1 wly decaying impulse re pon e) and/or mall damping c fficient ci 11 tory mod ) a good approximation i more difficult to achie e. The model pre ented in thi the i are good approximates of y tern ev n with the e characteri tic
PAGE 28
CHAPTER 3 MODEL FORMULATION In this chapter our first major contribution to the problem of system modeling i presented; namely, a comprehensive framework given by the GLFF model class. Although many of the components have been around for quite some time it is the unifying approach taken here that offers new insight. This new formulation encompasses not only GLFF networks but FIR and IIR models (for the discretetime case) as well. 3.1 General Model Form M(p) Recall the general system modeling configuration, shown in Figure 31 below where His the system to be approximated using the model M(p). The specific form of the model M(p) is now discussed. d H X l e M(p) y Figure 31: System modeling configuration. 22
PAGE 29
23 p is defined as a parameter set. Its elements are vector called the parameter vecto r of the model M(p). The parameter et contains three element (vectors), p={ w,1} where (3.1) N dimension ector w weight ector A ti me scale ector The model can be written in the form M(p)= M({N ,w 1})=/( ,w, H J ,A)) (3.2) where H,i( A) i called the model kemel at ta en, 0 n ,, and ; 1 a compon nt of the vector N. So M(p) i a function f, of the param ter ctor and the mod I kernel. The kernel argument 1 u d to imp! that the d main of peration ( or z) urnv r al. When olving pecifi pr blem th n tation H,i( A continuou tim ca and H,i( z A) in th di r tetim ca u ed m the on icier th definiti n ab in the d main Thi unifying repre entation can be u ed to illu trate many p cific mod I tructure In thi ection it h wn how to de crib the FIR, IIR and GLFF model clas e u ing thi formulation. 3.2.l FIR Model In the case of the FIR model M(p) = M({N ,w, A }) = I w,,H 11 (z,A) (3.3) 11 = 0 where
PAGE 30
24 /4=0 n HJ z, A) = H Jz) = IT H;( z ) i = O The functions H ,(z) are given by { 1 H;( z ) = z 3 2.2 IIR Model i = 0 i > 0 In the case of the IIR models where N, Lb n H n ( z, A) M(p) = M({N, w,A }) = n "" = ~ 1La nHn(z,A) n= l /4=0 n H n(z, A)= H n(z )= IT H;( z) i = O Just as in the FIR case the functions H;( z) are given by { 1 i=O H z = ;( ) I Q z i> 3.2.3 GLFF Model In the case of the GLFF models M ( p) = M ( { N W, /4 } ) = L W 11 H n ( Z, A) n = O where (3.4) (3.5) (3.6) (3.7) (3.8) (3 .9 )
PAGE 31
25 10 n H n ( z, },,,) = H 2 ( z ,},,,) f1 H ,1 ( A) 1=0 The functions H ,1 ( z ,).,) and H n 2 ( z, ).,) could be an tr n f r fun ti n In thi th 1 a restriction is made to consider special ca e that 1. are founded on the y tern model a umption d tib din h pt r 2 have kernel that are (typically) orthogonal. 3 have kernel that are ca caded I pa and allpa ub t m 4. will allow for modeling the unkn n tern ith r al and/ r ompl p le via the time cale vector X 5 have localized feedback. 6. guarantee tability. 3.2.4 FIR, IIR, GLFF Model Block Di gram To compar the tructur f th Figure 32, 33 and 34 feedback component limit it appli abilit capability but at a ignificant c t and lo amm th ir block diagram in th FIR m d I it lack of any Th IIR m d I offer the gr test f tability ontrol. Th GLFF model offer a practical alternati e that 1i al th IIR m de! in p rformance while maintaining simpl computational procedur imilar to th FIR ca e
PAGE 32
26 x( k ) x( k ) x (k1 ) x( k N) y( k ) Figure 32 : FIR model block diagram x (k) y( k ) y( k2) y (k1 ) Figure 33 : IIR model block diagram. y( k ) Figure 34 : GLFF model block diagram
PAGE 33
27 3.3 GLFF Model Kernels The GLFF kernels can be comprised of almost any transfer function. Some of the most popular discretetime and continuoustime kernels are gi en in Table 31 and 3th 3, respectively. These tables include each kernel's common name ymbol and n tap transfer function Hn( A). Most GLFF networks are deri ed from ortlw anal i nal sets (see Appendix A). This orthogonality of a kernel i al o indicated in the tab! An important part of a GLFF network is the time cale parameter ector X It allo for local feedback and is responsible for the popularit and u fuln of thi ntir class of models. The dimen ion of the time ale tor dim /4 d t 1min th number of degreesoffreedom a kernel ffer dim(A) i al o gi n in the tab! When fitting a model to an unkn n t m th param t r mu t b ptimall selected to achie e the minimum int grat d quar d ff r. Leg ndr and Kaut z. Th ar Ii t d in Tab! 31 and r tri tion on th time cale param t r ar gi n in abl me f th arli t di cu ion of the e kernels ar f und in the rk b Y ung [Y u62a], P rez [Per88], and Principe [Pri93]. he m t p pular continuou time GLFF kernel are the Gamma, Laguene, Leg ndre Jacobi, and Kautz They are listed in Table 3 3, and restriction on the time cale parameter are given in Tab l e 34. Some of first to describe these kernel are Lee [Lee33], Kautz [Kau54], Arm trong Ar[Arm59], and Principe [Pri93].
PAGE 34
28 T a bl e 31 : Di sc r e t etim e GL FF ke rn e l s. Kerne l n th ta p H j2 ( Z, A) H n ( z, A) o r t h o d im ( A ) H k 1 ( z, A) sy mb o l go n a l De l ay H n(z) 1 I n n o 0 z z Gamma G n(z, A) 1 k 1 ( k J n o 1 1 ( 1 A) z1 l( lA ) z1 L ag u e rr e L n(z, A ) I A ( 2 1 A J yes 1 z I l A z1 1A z lA z1 lA z1 Lege nd re P ,,(z, A) 1 zl A k 1 n1 z1 A i yes ;?: 1 1A j z1 l A kz1 1 A n I rrl X I z 1=0 z Ka ut z K n(z, A) I X J1I AJ n Z1 A ~ yes 2::: 1 z k 1 A kZ I l A I f1 l A I 1 A z I nZ i=O iz ) Ta bl e 3 2: R es t ric tion s o n di scretet i m e t i m e sca l e p ara m eters Kernel Valid time scale D e l ay n o n e G am m a A rea l 0< A < 2 L ag u erre A r ea l !Al < 1 Lege ndr e A k rea l A k = A k, 0 < A < 1 Ka u tz A k co mpl ex jA k j < 1
PAGE 35
29 Table 33: Continuoustime G LFF kernels. Kernel n th tap H j 2 (s,A) H kl(s, A) H n(s, A) orthodim(A) symbol gonal Gamma G nCs A) 1 A c:AJ no 1 +A Laguerre L n(s,A) sA ye 1 s+A s+A +A +A Legendre P n(s,A) fiI: sA fii: n1 A y l k IT I s+A 1 s+A k s + A n 1 = 0 +A 1 Jacobi J nCs A) fiI: sA fii: n 1 S A l k IT I +A k s+A 1 + A n 1 = 0 + A 1 Kautz K ,i(s,A) (s~A 1 j) (sAk )(sit") ,pX;: hi a ),Xs'\) l (s+A 1 s+"S) (s+ A" )(s+ it ") (s+~I +J,,) I +1i +,t,) Table 34: Re tricti n n ntinu u tim tim al parnm t r K e rn e l Va lid tim e scale Gamma A real A>O Laguerre A real A>O Leg ndre A/.. real A = 1k+ I A A> Q k 1 Jacobi A k r al, Ak = k + l Kautz A k complex, A = An+ }A, ,.,, AR k >O
PAGE 36
30 3.4 GLFF Kernel Derivations and Properties The Laguerre Legendre, Jacobi, and Kautz GLFF networks have the interesting property that their kernels are orthonormal. Hence, we call these orthonormal ge n e rali ze d lin e ar f ee dforward (OGLFF) networks. Next these kernels will be derived from first principles. This is done in the continuoustime domain The derivations will make evident their relationship to causal LTI system The connection between the network structures and the mathematical description of causal LTI systems is the primary reason they are able to approximate these systems so well. 3.4.1 ContinuousTime Laguerre Kernel Consider the sequence { (At f i.1 }. Regarding this as the solution to an ordinary differential equation that represents a causal LTI system the system would consist of one pole with multiplicity n + 1. Orthogonalizing this sequence yields a new sequence {l n (t A)} given by (3 11) Now taking the Laplace transform gives the Laguerre model kernel (3 12) 3.4.2 ContinuousTime Kautz Kernel (Complex Time Scale Vector) Consider the sequence { e .l" } where the l n' s are distinct compl ex scalars (i.e. components of the time scale vector A) Since the poles are complex they must exist
PAGE 37
31 in conjugate form in order to satisfy the realizability requirement. Orthogon alizing this sequence yields a new sequence { k n (t, A)} given by n k n( t A)= I,Re{ ,\,e ,l 1 } 112:0 t2:0 Re{J,}2:0, Vi (3.13) i=O It is now shown how to obtain the coefficients A ni Taking the Laplace transform of the new sequence gives the Kautz model kernel K ( A)= O {J,. (t A)}= ~(s+IA ,,,l) rr ml (sA )(sA ~ ) n1 2 "' t ;., "2m (s + A m)(s + A ~ ,) ;=0 ( + A )(s +A ~ )' 111 = O l 2K t A = t A}= (sA,)(sA ~ ) m=Ol n1 2111+ 1 ( ) { kzm + I ( ) ( A )( X ) ( A )( + X)' ' ., 2 S + 111 S + m 1=() S + I I Expre ss these functions in partial fraction form Grouping each pair of time element (A, ,A ~ ) gives the kernel ex pre ion 11 = 2m and 11 = _,n + l (3 .14 ) ale (3. 15 ) Hence the coefficients Ani are twice the re idue a t the pole Ai for th n th tap kern l. Con ider the sequence { e J" 1 } where the ,1 11 are di tinct r e al calar (i.e., component of the time cale v ctor A). Orthogonalizing thi equence re ult rn a new se quence { k n (t, A)} given by n k n( t ,A) = L, ,\ ,e A, t' 112: 0 t 2: 0 A 2: 0 V i (3 .16) 1 = 0 It i s now shown how to obtain the coefficients A 11 i Taking the Laplace tran form of the new sequence gives the Kautz model kernel
PAGE 38
32 M1 n1 /4 n L1 K n (t A) = {k n (t A)} = ~""" /1,nn_s i = I' ~_ i s+A n i= O s+A ; i=O s+A ; (3.17) where the last expression is simply the partial fraction expansion of the product form of Kn(t A). Hence A 11 ; is the residue of the pole at A; for the n th tap kernel. It will be real since the time domain function kn(t A) is real. 3.4.4 ContinuousTime Legendre and Jacobi Kernels Consider again the sequence { e A,,r } where the A,/ s are distinct r e al scalars. This is the same sequence from which the Kautz kernel with a real time scale vector is derived. This kernel is described by equation (3.17). Now consider the kernel with a special time scale parameter. Namely if A= (A1 A2 , A N f with A k = 2 k t A, A> 0 then we obtain the Legendre kernel. Likewise if A k = k + l we obtain the Jacobi kernel. 3 5 Impulse and Step Response Modeling The models discussed in this chapter can accurately represent any signal or system belonging to L\ lR?. ). We will be most interested in modeling causal LTI systems using sampled step response data. In this case a step function is applied to a continuous time system to generate its step response. The step response is sampled to obtain a discretetime representation of the system. Optimal continuoustime and discrete time GLFF models are determined using this data (optimization will be discussed in Chapter 4 ). For example using an N+ l tap Kautz filter the approximating transfer function is
PAGE 39
33 N H( ) = H(z) = L w n Kn(z,A) (3.18) n = O Suppose the discretetime impul e response is desired How can it be computed? Three simple methods are now given. ,._ l. Using the discretetime model H(z) pass the impulse x(k) = 8J) through thi model to generate an approximation to h(k) 2. Since x(t) is a unit step then X(s) = 1/ If g(t) is the tep response of the ystem then g(t) = h(t) x(t) or G( ) = H( ) I So or H ( s ) = G ( ) =D ( )G() 2 7 1 H( z ) = D( )G( )I ~ = D( )G( ) where D( z ) = " sT :+ 1 T Z + l (3.19) (3.20) Using the ampled tep re pon e data (k) pa thi through D( ) to gen rat h(k). 3. Since h(t) = d (t) u e a numerical appro imation to the deri ati e (backward dt difference perhap ) to obtain h(k) In thi ca e dg(t) = g(kT)g((k l)T) = g(k)g(k1) = h(k) dt t = kT T T (3 21) 3.6 Conclu ion (Model Formulation) Thi chapter has dealt with structures of model useful for appro imating cau al LTI sy terns. The generic model M(p) was tated. It wa hown how this formulation can be used to represent FIR, IIR and GLFF n tworks. Several GLFF kernels were derived, and discus ion was provided to illu trate their connection to the as umption of the systems under study. Tables of discretetime and continuoustime GLFF mode l
PAGE 40
34 kernels were included. These kernels are the most popular in the literature and as the derivations of the Laguerre Kautz Legendre and Jacobi revealed they are not theoretical abstractions, but orthogonalized versions of solutions to differential equation representations of several important classes of systems. Now that the model structure M(p) has been defined optimization of the parameter set p becomes the focus.
PAGE 41
CHAPTER4 OPTIMAL PARAMETER SELECTIO In this chapter the selection of the optimal (or near optimal) parameter set for GLFF networks is discussed. Recall p = { N w A} where Ni the dimension ector w is the weight vector and A is the time scale vector. Optimization for the e three vectors is handled separately below In fact the ability to separate the parameters and their optimization techniques is another benefit of the GLFF formulation 4 1 Model Dimen ion Consider a GLFF model con isting of a t of tim domain ba i function n (t A )L = o with frequency repre entation { n( A ) L = o Th model tructure i a weighted um of the e ba is functions (kem I ) In the s domain Y ( ) T M(p) = = Li wn 11 ( A) = W ( A) X( s ) n = O (4.1) where W = ( WO, W I ... W N ) T (4.2) In the time domain y (t) = [t, w \l) (t /4)} x (t) = [ w r (1,/4 )] x(t) (4.3) where 35
PAGE 42
36 (4.4) In this thesis, the theory will focus on modeling the impulse response of a continuoustime system. Given sampled step response data continuoustime and discretetime GLFF networks will be constructed to approximate the desired system. Hence the input, x(t), will be considered an impulse and the desired signal, d(t), will be the system impulse response. The selection of N is now discussed. N is a scalar for the GLFF model class Suppose N is a fixed finite positive integer. Then d(t) can be approximated by N d(t) = y( t) = L wnn(t,A) (4.5) n = O 4.1 1 Completeness To ensure that y (t) d(t) as N 00 the property of completeness is employed. A set of functions { n (t' A) r = O is called complete [Su71] on the space L 2 ( [a b]) if there exists no function d(t)E L 2c [a b]) such that b f d(t)Jt,A)dt = 0 Vn = 0,1,2 ... (4. 6 ) (l An equivalent definition is as follows Given a p1ecew1se continuous function d(t)E L 2c [a b]) and an c > 0 there exists an integer N > 0 such that b f ( d ( t) y( t) )2 dt < c (4. 7 ) a
PAGE 43
37 where y( t) is given by the equation (4.5). All of the GLFF kernel sets are complete on [0 00 ). An example of an incomplete set of orthonormal functions n( t,A)f =o defined on [O 2.n) is n( t,A ) = n( t ) = sin (n t ) (4.8) The function d(t) = cos(mt) 1s nonzero on [O 2.n) howe er it 1s orthogonal to { n (t' /4) r = O on this interval. amely, 2n J cos(mt) n (t /4 )dt = 0 Vn m (4 9) 0 4 1.2 Orthogonal Projection y( t) = L wnn(t,/4) (4 10) n = O spans an N+ldimen ional Hilbert Space S (see Appendix B). The Orthogonal Proj ec tion Th eo r em state that for every d(t) in a Hilbert space H it can be decomposed as d(t) = y (t) + e( t) where y( t) E S e(t)E S1.. S1. i th rthogonal _l_ complement of S and H = S EB S So y( t) can b con id red a th rthogonal projection of d(t)E H onto the pace S am ly ( l ) = Pro j ( d ( l ) ) (4.11) and n (t,/4 )r = O is an orthonormal ba i of
PAGE 44
38 4.1.3 Optimal Selection of N Completeness is important because it guarantees that if the approximation of d ( t ) using N+l weighted terms of n (t,A)r = o is not within a prescribed tolerance then more terms can be included to eventually reach the tolerance requirement. Increasing N may not improve the approximation of d(t) when a set of incomplete functions i s used. Having the completeness property provides the following choices when a GLFF network realization for a given N is inadequate: 1 Increase N until the prescribed tolerance is reached. 2. Select other kernels from the GLFF set and examine the approximation error. For a given N if none of the GLFF kernels chosen provide satisfactory results then N must be increased. If this is done using several kernels the rate of convergence for each kernel as N is increased could be monitored. Choose from this set the kernel that yields the minimum N needed to achieve the required tolerance Other options are to either raise the tolerance level or choose a different performance criteria One suggestion for setting the tolerance is to choose the acceptable error power to be a certain percent (5 % perhaps) of the system power. Davila and Chiang [Dav95] suggest another approach to estimating the dimen s ion vector N = (N 1 N 2 / for a general IIR structure. This method examines the eigenvectors of the input/output data covariance matrix. When the model orders are overestimated zeros will appear in the noise subspace eigenvectors. The number of ze ro s found can be used to estimate the model orders (dimension vector components ). This procedure could be applied to GLFF networks.
PAGE 45
39 4.2 Weight Vector One of the greatest advantages to working with a feedforward model is that the calculation of the optimal weight vector w (in the squared error sense) is simple. This is true even for GLFF structures. Although the realizability requirement for GLFF networks demands the weight vector be real the optimal solution for w is now derived assuming w, d(t) and {Jt,A) r =o are generally complex. Realizable optimal solutions are obtained by setting the imaginary components equal to zero. 4.2.1 Optimal Selection of w Assume d(t) is the impulse response of the desired system (i.e d(t) = h(t)) and approximate it with y( t) given by y( t) = L w 11 11 (t A) (4. 12 ) n=O The parameter set p = { N w, A} i to be cho en o th err r e(t) = d(t) (t) minimum in the integrated squared rror ense. amely, hoo e p to minimize 00 J = J ld(t)y( t )i1 dt 0 oo N oo = Jl dCt)l 2 dt2, w, J d (1),(t A )dt (4.13) 0 t = O 0 N oo N oo 2, w ; J d(t)f(t,A)dt + 2, JwJ J J,(t,A)j2 dt i = O Q 1 = 0 o If Wn =an+ jbn then set d J = 0 and d J = 0 to find Wn that minimizes J. Hence Jan Jbn dl 00 00 00 J = J d*(t)n(t,A )dt J d(t) ; (t,A )dt + 2an J l,/t,A )12 dt = 0 an O O 0 (4.14)
PAGE 46
40 so 00 f Re[ d(t)(j) : (t A )]dt a =o _____ n oo ( 4.15 ) f lJt Af dt 0 When 11 (t,A)r = o is an orthononnal set, the bottom integral equals unity yielding 00 a n = f Re[d(t)(j) : (t A )]dt ( 4.16 ) 0 Likewise b 11 can be derived d J = j j d(t)(j) : (t A )dt j j d*(t)(j) n (t,/4 )dt + 2bJl J t A f dt = 0 ( 4 17 ) d bn o o o so 00 f Im[ d(t)(j) : (t A) ]dt b ==o _____ n oo ( 4.18 ) f 11 (t,A f dt 0 Again if n (t A)r = O is orthononnal then equation (4.18) reduces to 00 b n = f1m[d(t)(j) : (t J)]dt ( 4.19 ) 0 Hence combining the above results gives 00 w 11 =a n + jb 11 = f d(t)(j) : (t A )dt ( 4.20 ) 0 When 11 (t /4)r = O is real then W11 is real and can be written as 00 w 11 = f d(t)(j)Jt /4 )dt (4.21 ) 0
PAGE 47
41 So 00 w= f d(t)r/J(t,J)dt (4.22) 0 is the weight ector that minimizes the energy of e (t) From this point forward only real weight vectors will be considered. 4.2.2 Relationship to WienerHopf Equations It is easy to show that w is in fact the solution to the WienerHopf (or Normal) equations The performance function can bee panded a follo 00 00 00 J = J e ( t) dt = J ( d ( t) y ( t ) 2 dt = J ( d ( t) w T t A ) ) 2 dt 0 0 0 = j d 2( t)dt 2 wT j d(t)r/J n( t A ) dt + wT [j r/J(t A r/J T (t,A )dt]w 0 0 0 (4.23) now define the following term 00 d e = J d \ t ) dt 0 00 p = J d(t)r/J(! J)dt (4 24) 0 00 R = J r/J(t A)r/J T (t,A)dt 0 Then the performance function J becomes (4.25) To minimize this function take the gradient of J with respect to the weight vector w and equate to zero '\ \ J =2p+2Rw=O (4.26) or
PAGE 48
42 Rw=p (4.27) These are known as the WienerHopf or Normal equations. Solving for w gives the optimal solution R 1 w= p (4 28) Now the OGLFF kernels n (t,A )r = 0 form an orthonormal set so 00 R = f (j)(t A )rjl (t,A )dt = I (4.29) 0 Hence the weight vector becomes 00 w = R 'p = p = f d(t)(j)(t,A )dt (4 30) 0 which is the solution derived in the previous section. 4.2.3 Minimum Squared Error Performance Function Now that the optimal weight vector has been derived the minimum squared error can be computed. Recall (4 31) For the optimal weights it was found that (using OGLFF kernels) w = p and R = I. Plugging these into the equation for J yields (4.32) The above optimal weight vector wand minimum squared error l opt was derived using a fixed time scale X Notice that w and l o pt are functions of X Hence l o p t can be written ( 4.33)
PAGE 49
43 to illustrate the importance of A on the optimal solution. Appendix C presents a derivation of equation (4.32) in the sdomain. 4.2.4 Computing w Using Step Response Data In Chapter 5 applications of the GLFF networks are gi en The sy terns there are continuoustime systems. The available information of the e y tern i ampled tep response data. GLFF networks will be u ed to appro imate the input/output tep response beha 1or. Method for computing the optimal weight ector w u ing sampled step respon e data are now gi en. Least Squares Solution. [Hay9 l] Con ider the model dimen ion N and time scale vector A fixed or pre iou ly cho en to be optimal. Di u ion of model dimen ion wa g1 ven rn Chapter 4.1. Di cu ion of the optimal time cale will be given in Chapter 4.3. Since the information i ampled ( di cret tim ) tep re pon e data the following definitions are needed Let the input x( t ) be the unit tep and (t) be the re ponse of the system under que tion ampled at interval T then define g(k) = g(k1), k=O,l, ... Ka the desired signal. Define al o xi(k) = fh sample at the /h tap of the GLFF network y(k) = k th output of the GLFF network In this discretetime ca e the matrix Rand ector pare given by K R=h] where~ 1 =(x,,x j )= I x,( k)x/k) k = O and (4.34) ( 4.35)
PAGE 50
44 p=[;:] where K P i = (g x i ) = L,g(k)x i (k) k = O ( 4 36) This gives the normal equations Rw = p with solution w = K 1 p. This is the Wiener Hopf solution. Notice in this case that R will not be the identity matrix Only when the input is an impulse (or white noise) will R = I. Nevertheless, the WienerHopf solution provides the optimal weight vector (for fixed N and J) regardless of the input. Integral Equation Solution [Clu91] A more direct evaluation of the weight vector solution given in Chapter 4.2.1 is now considered. Recall 00 wn = f d(t){/} 11 (t,A)dt, n=O,l, .. N (4 37) 0 is the optimal weight solution assummg an orthonormal basis. For the system modeling problem this is accomplished by applying an impulse (or white noise) into the system. In which case d(t)=h(t) the impulse response. Consider, for example the Laguerre kernels given by { (/J /1 (t A) r =o ={l n (t /4) r = O. Some of its properties can be utilized to derive a simplified integral equation solution. Again work with g(k) k=O l .. K samples of the continuoustime step response g(t). Since the continuoustime impulse response is h(t) = dg(t)/dt, the optimal weights can be written as 00 w 11 = f h(t)l n (t A)dt 0 = f 00 dg(t) l (t A)dt d 11 0 t ( 4 38) = (g(t)lJt,A)I ~ j g(t) dlJt A) dt 0 Jt
PAGE 51
45 where integration by parts has been used. For a stable system g(t) it must have Jg( oo )J < 00 Also using the property of the Laguerre functions [Par7 l] ( 4.39) yields liml J t ,A)= O (4.40) Hence J oo () dl n (t A) d W n = g t ''t 0 dt (4.41) Let T 55 be the steady state response time of g(t); namely T ss i the time which for all t > T ss the output can be consider constant, g(t) = g s = con tant. In this case the integral can be broken up as T I ( )df n (l,A)d J oo df (l,A)d W =gt'tg n l n at Sl dt 0 =f (t/ 1 "i/')dt+g,J ,, (T, A) 0 (4.42) From the Laguerre function property above it is easy to compute the derivative relationship (4.43) The weight computation becomes n l T ,, T ,. Wn = 2A L, f g(t)l;(t,A )dt + A f g(t)ln(t,A )dt + g s Jn(~ s A) (4.44) t = O O O
PAGE 52
46 Up to now the continuoustime step response g(t) has been assumed. To compute w 11 using the sampled step response g(k) an integral approximation technique must be used such as Euler integration T ~ K f g(t)l/t A )dt = /). L g(k)((k A) (4.45 ) O k = O where !'1 = ~ s (integration interval) K (4.46 ) Finally this gives the integral equation weight solution n l K K W 11 = 2/4/'1 LL g(k)((k A) +A/'1 L g(k)l 11 (k A) + g ) 11 (~ s A) (4.47 ) i = O k = O k = O 4.2.5 Adaptation of w Note that the weights can be computed by solving a set of linear (WienerHopf) equations. This linear relationship among the weights corresponds to a convex performance function with one minimum point. So the weight vector can be easily solved in an adaptive manner using a gradient search technique such as Newton s method or LMS [Wid85]. Unfortunately this nice relationship does not hold for the time scale vector. It is (in general) nonconvex with multiple stationary points. The number of which depends upon the dimension N of the model used. 4.3 Time Scale Vector Optimization of the parameter set elements N and w have been discussed to this point. Their solutions are fairly easily found for GLFF networks and are the only
PAGE 53
47 elements of p that are required for FIR and IIR models. As has been previously stated the power of GLFF networks is in local feedback. This feedback is controlled by the time scale vector. Although IIR structures also have feedback it can not be separated and optimized independently from the taptotap weights like the GLFF structures Unfortunately the elements of the time scale vector are not linearly related; hence the performance function J is usually noncon ex w.r.t. these parameters. As such the bulk of the work in finding an optimal p involves finding an optimal X This section is devoted to this is ue It begins by examining everal method that currently exist for either determining exactly or approximately the optimal X Each method usually addresses a particular GLFF kernel most often the Laguerre A few papers have considered the nonorthogonal Gamma and a few di cus optimization of the more complex Kautz kernel. Although the ex a c t optimal olution i ought Kautz has pointed out that nonoptimal olution (which may be derived in a simplified manner) can still produce excellent results. His comment uggest that one should not search to exhaustion for optimality if it i po ible to get clo e with much less work. After the overview of optimization methods found in the literature we propose a new technique relying heavily on complex variable theory. This is a method that extends the concepts of quared error analy is to the complex domain. Under certain assumptions the optimal time scale can be fairly easily determined by simply looking for zero crossings.
PAGE 54
48 4.3.1 Optimal Selection of A (Historical Review) In this section a summary is given of current methods found in the literature to either locate or approximate the optimal time scale X Each method usually a ddresses a particular GLFF kernel most often the Laguerre. A few papers have considered the nonorthogonal Gamma and a few discuss optimization of the more complex Kautz kernel The methods presented in the literature can be grouped into the following categories according to the employed strategy or desired objective: 1 Approximation via the impulse response. 2 Asymptotic solutions. 3. Achieving prescribed measurement objectives (other than minimum squared error). 4. Matched filter approach 5. Moment matching. 6. Satisfying necessary conditions for stationarity of the squared error function. 7. Iterative and adaptive search techniques. Approximation via the impulse response. These methods approximate /4 given a priori knowledge about the impulse response h(t) of the system The required information is different for each approach. Either h(t) is needed in analytic form or tabular form (samples of h(t) at discrete points t = tk, k = 0 1 2 .. ) These methods appeal to the necessary requirement for all stable and realizable line a r systems ; namely its transfer function H(s) must be representable as a rational function with no righthalf plane poles or multiple jaxis poles. Under this assumption h ( t) can be represented as a sum of exponentials of the type ( 4.48 ) where
PAGE 55
49 (4.49) A frequency domain method for approximating /4 is as follows [Su 71]. Suppose h(t) is known (or approximation to it) and its Laplace transform H(s) is determined Note that the form of h(t) that is given may not produce a rational H(s). ow find a rational approximation Ha(s). For instance use the expression for H(s) and expand it into a Taylor series abouts = 0 to obtain (4.50) Then a rational function can be found to approximate H(s). One way to accomplish this is to use a Pade approximant. The poles of thi approximation can be used as an initial guess for the poles of H(s) The e pole location could be incrementally adjusted to locally optimize. A time domain approach i to use point matching The ba i of this technique is attributed to Prony [Ham73]. The trat gy i to force an approximate impulse response h(t) to match h(t) at evenly paced point As ume h(t) is known at points separated by the spacing t..t. U ing the approximation h(t) = L Ae J, 1 (4 51) 1 = 1 it is possible to create N linear equations (using the known points of h(t)) w h ose solution gives coefficients ao a 1 , aN l of an Nh order characteristic polynomial (4 52) Once the coefficients are known, find the roots of this polynomial x 1 x 2 ... XN. They are related to the time scale parameters by
PAGE 56
50 1 X = In x I tit I ( 4.53) Another time domain approach involves the use of an orthogonal signal set [Kau54]. Assume h(t) and its first N derivatives are known. Using these functions generate a set of orthogonal functions { Jt) r = O satisfying the relationships 0 (t) = h(t) 1 (t) = h(t) + b 10 0 (t) 2 (t) = h(t) + b 2 1 Jt) + b 20 0 (t) (4 54) where the b 111 n's can be found such that f /t k (t)dt = 8 jk ( 4 55) 0 Since each n(t) is a linear combination of the derivatives (from the 0 th to the n th ) of h(t) then
PAGE 57
51 This is based on the region of convergence (ROC) of the power series equivalence of the Laguerre model. Justification for ;L, in the place of AN is made on the basis that most engineering applications require the use of a large number of terms N in the Laguerre expansion to achieve accurate results. In this case, ;L, is viewed as a good approximation to the optimal time scale A N Bruni [Bru64] offers another approach involving optimization around the ROC of a Laguerre model expansion. Assuming h(t) is the impulse response of an LTI lumped parameter system M h(t) = L Ae p, r, p 1 =a ; + jb 1 a 1 > 0 (4.58) i=O decompose a Laguerre model representation of h(t) into M infinite series 00 00 M h(t) = L wJn(t,A) =LL w,)n(t,A) ( 4 .5 9) n=O n=O 1=l where Wni is the n th coefficient of the expan ion of the /h exponential in h(t) With this formulation the Fourier transform of h(t) can be repre ented a a um of M infinite complex geo metric series M i give the magnitude of the term in the /h infinite serie M i hould be minimi ze d for each i in order to achieve the most rapid convergence of the serie Since Mi is a function of A there is an optimal A (denoted Ai) for each M i, Choose as the optimal A for the entire model the Ai that minimizes the largest M i. Methods similar to this one are also given in [War54] and [Hea55]. Achieving prescribed measurement objectives. In addition to locating A that minimizes the squared error, approaches that minimize other measurement criteria have also been considered. One such method employed when using the Laguerre
PAGE 58
52 signal set {l n (t A) r = o is given in [Par7 l]. Consider the approximation of d(t) using a n N+ l tap Laguerre model N d(t) = y (t) = L w) n (t /4) (4 60) n = O Now define the following quantities 00 00 f td 2 (t )dt f td 2( t)dt M 1 = and M =0 __ 2 00 (4.61) f d 2 (t )dt f d 2( t)dt 0 0 M 1 is a measure of how rapidly d(t) decays. M 2 is a measure of how smoothly d(t) decays. Parks show that the time scale /4 given by /4 = .J M 1 / M 2 minimizes the maximum error over the class of all signals with these measurements The advantage to this approach is that a complete knowledge of the signal to be approximated is not required In [Fu93] a similar approach is taken when usmg the discretetime Laguerre model. In this case the optimal /4 is found that minimizes the measurement 00 J = Inw ~ ( 4.62 ) n = O This performance measure is chosen because it linearly increases the cost of each additional term used in the Laguerre expansion This enforces a fast convergence rate. The optimal /4 is obtained using a discretetime formulation of M1 and M 2 above. A relationship for J is derived that depends only on M1 M 2, and X The optimal /4 is then determined to be
PAGE 59
53 A= 2MllM 2 2MI 1+.,j4M1M 2 M ~ 2M 2 (4.63) Matched filter approach. In [Kuo94] a matched filter bank is used to estimate the optimal A for the Gamma model. Recall that {gJt,A)r = o is a nonorthogonal GLFF signal set. Let d(t) be approximated by N d(t) = y (t) = L, w n g n (t,A) (4.64) n= O A matched filter is constructed based on the optimal A condition ( g t+ i (t A) d(t)) = 0 where g ~+ i (t, A) is the component of g N+ t (t A) that is orthogonal to the previous basis components. A bank of such filters is created where each bank uses a different /4. To optirruze A over the region (0 Am ax ) an Mbank filter structure is created where the m th filter uses the time scale A= Am ax I M When the desired signal d(t) is applied to the filter bank a change in sign of the M outputs wi II occur between filter m and m+ l for a value of A that satisfies the optimality condition. When there are multiple occurrences the associated A for each case must be used to determine which one achieves the minimum squared error. Moment matching. The method presented in [Cel95] uses a Gamma model to approximate a signal d(t). The Gamma model is conceptualized as a rotating manifold rn space. The kernels {g 11 (t,A )r = o are the basis vectors that define the N+l dimensional space of operation and the time scale parameter A is the feature that allows for rotation. The optimal A is estimated using a moment matching method. This method attempts to estimate A by equating a finite number of time moments of
PAGE 60
54 y( t) (the output of the Gamma model) to those of d(t). Namely it is required that m 11 (y(t)) = m 11 ( d (t)), n = O, ... N where 00 m 11 (y( t)) = J t 11 y (t)dt ( 4.65) 0 Using these functions an N+ 1 order polynomial in A p(A) is formulated. The roots of p(A) are estimates of the local minima of the squared error performance function J The root that minimizes J is chosen as the optimal X Satisfying necessary conditions for stationarity of the squared error. Some of the most popular methods for optimizing the time scale of the Laguerre model center around the necessary conditions for stationarity of the squared error performance function J. For instance, approximating d(t) with an N+ 1 tap Laguerre network N d(t) = L wJA )lJt,A) (4.66) n = O where the weight vector components wn(A) have been given an argument in A to emphasis the dependence of these weights on the time scale. For a fixed A, the squared error is minimized when the weights are computed as 00 wn(A) = J d(t)ln(t,A )dt (4.67) 0 The squared error function reduces to oo N J = J d 2 (t)dt L w:( A) (4.68) O n=O To minimize J (over J) the evaluation of the L must be m(l,Ximi zed This requires that (4.69)
PAGE 61
55 which yields the interesting relation (4.70) A proof of this assertion is given in Appendix D Hence the optimal value of A is either a root of w N (A) = 0 or W N+ i(A) = 0. Solving for A analytically is a tedious task because w N (A) = 0 is a polynomial in A of degree (N+M) if d(t) is a system with M poles. The first proof of this (using mathematical induction) was given by Clowes [Clo65] and was restricted to the case where d(t) had only simple poles. A more concise and less restrictive proof was given in [Kin69]. This paper also derives expressions for the 2 nd derivative (4.71) These equations may be used to determine whether a stationary point is a maximum or minimum. In [Mas90] a discretetime formulation of the stationary conditions is formulated and the authors use numerical techniques to find the roots of the polynomials w N (A) = 0 and W N+ i(A) = 0 [Sil93] derives these conditions for the Gamma model. The derivation given by Clowes was valid only for systems with rational transfer functions. [Wan94b] extends the results to approximating any L2( 1R?. ) stable linear system. [Si195] desc1ibes another efficient way to compute the derivatives :A w n (A) and uses this to compute highorder approximations of the weights w N (A) and WN+J(A) via Taylor series and Pade approximants. This procedure is illustrated in both continuoustime and discretetime when approximating a system excited by an arbitrary input.
PAGE 62
56 Iterative and adaptive search techniques. When representing a system with a discretetime GLFF network the elements of the time scale vector are restricted to bounded regions to guarantee stability. In particular the discretetime Laguerre model requires AE (1 1). One method of optimization is to employ exhaustive search over this interval seeking the A that minimizes the squared error performance function. Since this would be done in discrete steps, the only way to get continually closer to the optimal /4 is to increase the number of points in the search region [Kin77] and [Mai85] discuss minimizing the squared error by employing a Fibonacci search technique over X This is an interval elimination procedure wherein the region in which the minimum lies is sequentially reduced using a predetermined number of steps. [Nur87] uses an adaptive technique whereby the optimal /4 is chosen to be the center of gravity of the presumed poles of the system. The method works as follows. Initialize A say AQ. Using the Laguerre model, with this A the parameters of a difference equation representation of the system are estimated. From the characteristic equation determine the poles of the estimated system Choose as A ;, i=l 2 ... the center of gravity of these poles. Repeat this process until l..1 ;+ i A ; I < 8 where 8 is a small number. Once the center of gravity converges to a fixed point stop the adaptation 4 3.2 Optimal Selection of /4 (Complex Error Analysis) A new approach to estimate the optimal time scale 1s now offered. The methodology used here is one of extending the squared error performance function to
PAGE 63
57 the complex domain. This is accomplished by defining a new complex performance fun c tion It will be shown how this function is derived how it relates to the conventional squared error function and what additional information it offers The optimal time scale can be determined by looking for zero crossings of the imaginary component of the new complex performance function For systems consisting of a cascade of 2 nd order subsystems containing pairs of complex conjugate poles examples will be presented that illustrate how the imaginary component of the new complex performance function can be used to obtain the optimal (in the minimum squared error sense) time scale parameter. These examples will solve the system identification problem with a GLFF model containing the Laguerre kernel. 4.3 2 l Deriv a tion of the complex performance function. In Appendix C an s domain description of the squared error performance function is discussed It is given by J = 2 ~ f E (s )E() ds (4.72) C where C is a contour taken in the CCW direction enclosing the entire lefthand plane (LHP) in s E(s) = D(s) w 7 (s /4) and w is a real weight vector. For a fixed time scale A the optimal weight vector is w = 2 ~ f D(s)(s A )ds (4.73) C yielding the performance function (4.74) where
PAGE 64
58 d e = 2 ~ f D(s)D(s)ds (4.75) C Now consider restricting evaluation of the performance function to the upper LHP only. Namely, create the complex performance function Ju= 2 ~ f E(s)E(s)ds (4.76) c,, Make the following definitions: du e = 2 aj jD(s)D(s)ds (complex energy) c,, Wu= 2 aj f D(s)
PAGE 65
59 be considered as a mapping onto and the complex performance function as a mapping onto C A simplified expression for l u is derived as follows J u = 2 ~ f E(s)E(s)ds cu = 2 ~ f[D( s) wT (s,1) ][D (s)wT (s A) ] ds Cu = 2 ~ f D (s) D ( s)ds wT 2 ~ f D (s) (s, A )ds (4.80) cu cu w T 2 ~ fD (s) ( A)d +wT 2 ~ f (s A) T(s,A)dsw Using the definition of the complex energy, d 11 e a nd complex weight vector, w 11 as well as the following relations w T 2 ~ f ( ,A ) T (s,A )d w = w T w (4.81) cu and 2 ~ f D (s) (s A )ds = w (4.82) Cu gives J = d WT w J_ WT w + J_ WT w u u e u 2 2 ( 4.83) Note that 2 Re[ J u ] = 2 Re[ du e ] WT 2 Re[ Wu ] = d e WT W = J (4.84) which is what is expected since J = 2 ~ f E(s)E(s)ds = 2RJ 2 ~ f E(s)E(s)dsj = 2Re[ JJ C l (4.85)
PAGE 66
60 4.3.2.2 Complex performance function examples and experimental results. Three test cases have been chosen to study the complex performance function, lu. One 2 nd order system and two 4 th order systems are modeled using a GLFF network with a Laguerre kernel. The transfer function for each example follows: Example 1: (2 nd order system, rapidly decaying impulse response) l l H(s)=+s+2j s+2+ j Example 2: (4 th order system, rapidly decaying impulse response with slight oscillation) l l 1 l H(s)=+++s+2j s+2+ j s+l2} s+l+2j (4.86) (4.87) Example 3: (4 th order system, slowly decaying impulse response with oscillation) 2 2 3 3 H(s) = +++ s + 0.4 j s + 0.4 + j s + 0.6 2.8 j s + 0.6 + 2.8 j ( 4.88) Experimental results for each example are tabulated and graphed on the next few pages In Figures 41 and 42 the squared error performance function, J (solid curve), and twice the imaginary component of the complex performance function, 2Im[lu] ( dash curve), are plotted for the first four model orders for example 1, Figures 43 and 44 for example 2, and Figures 45 and 46 for example 3. Notice the imaginary component has zero crossings near the stationary points of the squared error
PAGE 67
61 performance function. Im[l u ] crosses negativeto positive near local minima of J. So the complex performance function provides an estimate of the optimal time scale through the imaginary component. In fact Im[l u ] can be considered an estimate of V /4 J (A ) that is best near A opt Numerical results are provided in Tables 41 42 and 43 to illustrate how close the nearest zero crossing (nearest ZCu) of 2Im[Ju] is to the optimal time scale as the model order increases Also the squared error is computed at the optimal J(A op i) and nearest ZCu J(A nzcu), time scales. For this analysis Im[lu] was scaled by 2 to provide similar scaling to J since 1 = 2Re[l u ]; howe v er the zero crossings are unaffected by this scaling By examining J(A 0 p 1 ) and l ( A 11 zcu) it is found that when using a Laguerre network with sufficient order to approximate the unknown ystem no degradation occurs by choo s ing /4 at the ne a re s t ZC u. In fact the small difference in A o p t and /4 112 u is partly due to numerical appro x imation error when earching for minima of J and roots of Im[l u ].
PAGE 68
10 co C 0 20 "B C :::i LL 30 C E 0 40 a.. "O 50 E 0 z 60 62 70.___ ___ ___. ____ ___._ ____ __,__ ____ ~~ 0 2 3 4 5 Time Scale A Figure 41: J for example 1 (orders = 0 to 3). (order=0) (order=1) 0 1 0.01 0 05 0 005 ~::i 2. 0 f ..., ,, 0 05 ~::i 2. f ..., 0.1 '~~''~ 0 3 x10 2 3 A (order=2) 4 5 0 5 0 0 5 :::'' : \\ 1 ._.___.____.______.._____._,.____,__, 0 2 3 4 5 A '=, 2. f ..., 0 0 005 r 0 01 ~....,.__~'~~ '=, ..., 0 4 x10 2 3 4 A (order=3) 5 2 f O ~=~~+! .. : i ,,,1 ..., ( \ 2 1::, ______ ...., : . I I ; \. 4 '~'~~ 0 2 3 4 5 A Figure 42: J and 2Im[lu] for example 1.
PAGE 69
5 (0 :s10 C 0 g 15 :::J LL a> g 20 n1 E .g 25 a> 0.. g 30 N 35 z 40 63 45 ,__ ___ ..__ ___ ..__ ___ ~~~~ 0 2 3 4 5 6 Time Scale ,),.,, Figure 43 : J for ex a mple 2 ( orders= 0 to 7 ). (order=0) 0 4 ~~~~~ '=' 0 2 \i """""""" 2., ,.. ,,.""'_,,. ..E O ......... '. .._\, ... .. ,_ = ~ !! 7 0 2 r++t 0 4 ,__ __ _.__ __ __._ __ ___. 0 0 1 \ ~::J 0.05 '\" """ E o .. L.. ........ ,<~\ f , \ / 2 ),.,, (order=2) 4 6 .............. .......... __ 0 05 + ;, +i \/ 0 1 ~~~~ 0 2 4 6 ),.,, (order=1) 0.2 \ \ ," 0 1 \ ~ ;.., ..E O .... / .... ......... 7 \ __ ,,/ 0 1 0.2 ..__ __ .......__ __ __._ __ 0 0 1 i I 2 ),.,, (order = 3) 4 6 0 05 t """ r++i r O \ta P ~,~ 0 05 +: .: \} 0 1 ..__ __ .......__ __ __._ __ 0 2 4 6 Figure 4 4: J and 2 Im[l u ] for exa mple 2
PAGE 70
5 Q) 0 C ct! E 15 0 't: Q) 0... g 20 N co E 0 z 25 64 30 ~~~~~~~~ 0 2 3 4 5 Ti me Scale A Figure 45: J for example 3 (or ders Oto 11). (order=0) 1 (order=1) 6 ~o 5 ~ ++l ,~0 5 ~ .+4 0 .... .... . ............. ___ __ .. : : 1 0 0 '\j 2 A (order=2) 4 6 ... 0 2 ~~~~ 0 2 4 6 A ,,,........ ...... .. ... .......... o ~ "ct .. 0 .. ....... _, .... j 2 A (order=3) 4 6 0 6 .....,.~, ::?.. 0.2 0 L +, ... [\ ..... 1 /,,,,\ +~ ~ I :, "\, ,,,' \./ \ .. ,, 1 '\,.. .. .. i . 0 2 ~~~~ 0 2 4 6 Figure 46: J a nd 2 Im[Ju] for example 3
PAGE 71
65 Table 41: Perlormance function characteristics for example 1. N A. opt A nzcu l(A 0 p 1 ) in dB l(A nzcu ) in dB 0 2.47 2.48 24.7 24 6 1 1.60 1.58 33 8 33.7 2 2.31 2.31 49 9 49 9 3 1.92 1.92 60 7 60.7 Table 42: Perlormance function characteristics for example 2. N A om A nzC u l(A. 0 ,)/ ) in dB l( AnZCtJ in dB 0 2 99 3.12 11.6 11.6 l 1.50 1.43 15.6 15.5 2 2.58 2 60 20.8 20 8 3 1.81 1.78 24.8 24.7 4 2.45 2.47 29.3 29.3 5 1.93 1.92 33 3 33.3 6 2.39 2.40 37 7 37 7 7 2.00 2.00 41.8 41.8 Table 43: Performance function characteristics for example 3. N A oo t A nzcu l(A. 001 ) in dB l(An zcu ) in dB 0 0.82 0.34 0.5 0 2 l 2.30 2.46 2.7 2.7 2 1.29 1.21 5.2 5 1 3 2.12 2.12 9.9 9.9 4 2.87 2.89 12.7 12.7 5 3.56 3.59 14 1 14.l 6 4 21 4.22 14 7 14 7 7 4.82 4.76 14.9 14.9 8 2 24 2.24 17.1 17. l 9 2.62 2.64 19.6 19 5 10 2.25 2.24 22.9 22.9 11 2.57 2.56 25.5 25.5
PAGE 72
66 4.3.3 Further Analysis of the Complex Performance Function In the previous section, it was shown that the squared error and complex performance functions can be expressed as J =d e WT W J u = du e WT Wu where and d e = 2 ~ f D(s)D(s)ds C w = 2 ~ f D(s)
PAGE 73
67 So we can write J and lu in terms of components of the complex energy signal and complex weight vector. J u = (d ue (2w u Rf w u )/(2d ue R) =( d ue w~Rw u)/d ue R (4.95) (4.96) Now we are interested in 2Im[ JJ twice the imaginary component of the complex performance function. This has the expression (4.97) This extension to the complex domain has shed light on new information. The imaginary component of ] 11 has been shown to have roots (or zero crossings) near the stationary points of J one of which provides an approximation to the optimal time scale solution. In an effort to capitalize on this discovery the analytic performance function l a, will be defined and studied as a method of achieving the same information while requiring less a priori knowledge of the unknown system D(s) As will be seen l a is composed of an energy signal and a weight vector that is the analytic signal representation of the associated squared error components. 4.3.3.1 Analytic performance function. Consider the components of lu. The complex energy signal, du e can be written as
PAGE 74
68 d u e = 2 ~ f D(s)D(s)ds c,, 0 00 = 2 ~ J D(a)D(a)da+ 2 ~ J D(jm)D(jm)dm (4 98) 0 = [ 2 1 D(jw)D(jw)dw ]j[ ,', J D(a)D(a)da] Likewise, the complex weight vector, Wu, can be written as Wu= 2 ~ f D(s)(s,A)ds Cu 0 00 (4.99) = 2 ~ J D(a)(a,A)da+ 2 ~ J D(jm)(jm,A)dm 0 Now collecting real and imaginary components gives w. = { Re[ 2 ~ D(jw)
PAGE 75
69 d 0 e = [real part of d ue along + jcoaxis] + }[imaginary part of d ue along + jcoaxis] = [ 2 1 [ D(jw)D(jw)dw] + j[O] (4.102) 00 = 2 ~ J D(jco)D(jco)dco 0 and w 0 =[real partofw 11 along+ }maxis]+ J[imaginarypartofw 11 along+ jcoaxis] = { Re[ 2 ~ l D(jw)(jw A.)dw ]} + +m[ ,', [ D (jw) (jw A.)dw ]} (4 103) 00 = 2 ~ J D(j co )( j co,)., )dco 0 d ae is called the analytic energy signal and Wa is called the anal y tic weight vector. As in the case of the complex performance function we define the anal y tic performance function, l a l = d w 7 w a ae a or in normalized form Now decomposing the components of l a into their complex form gives where and d ae = d ae R + j d ae / W a = W a R + ]W a l d .. = [ 2 1 D(jw)D(jw)dw] d ae l = 0 (4 104) (4 105 ) (4 106 ) (4.107)
PAGE 76
70 w , = Re[ 2 ~ I D(jm)(jm,A)dm] w , = rm[ 2 ~ I D(jm)(jm,A )dm] (4.108) The goal is to rewrite land la in terms of analytic energy signal and analytic weight vector components only. Recall Likewise w = 2 ~ f D(s)(jm A )dm r + '.I D(jm)(jm,A )dm = 2Re [ 2 ~ I D(jm)(jm,A)dm] =2Re[w 0 ] d = 2Re[ 2 1 J D(jm)D(jm)dm] = 2Re[ d 0 e ] So we can write land la in terms of components of da e and Waas follows l = ( 2da e R (2waR f (2waR) )/(2da e R) = ( d a e R 2 W ; R W aR) / d a e R 1 0 = (da e (2w a Rf w a )/(2da e R) = ae w;Rw a )/d ae R (4.109) (4.110) (4.111) (4.112)
PAGE 77
71 We are interested in 2Im[1 0 ], twice the 1magmary component of the analytic performance function. This has expression (4.113) 4.3.3.2 Analytic performance function examples and experimental results. Again consider the three examples used to study the complex performance function. Figures 47 48 and 49 compare the complex and analytic performance functions. For each of the three examples the squared error performance function J (solid curve) twice the imaginary component of the complex performance function 2Im[Ju] (dash curve) and twice the imaginary component of the analytic performance function 2Im[la] ( dotdash curve ), are plotted for the first four model orders. Numerical results are provided in Tables 45 46 and 47 to illustrate how close the nearest zero crossing (nearest ZCu) of 21m[lu] and the nearest zero crossing (nearest ZC a ) of 21m[l a ] are to the optimal time scale as the model order increases. Also the squared error is computed at the optimal l(A op 1) nearest ZCu J(A nzc u) and nearest ZC a, ](A nzc a) time scales. It is interesting to compare lu and l a for these three examples. The difference in these functions is that Ju includes the line integral along the negative oaxis. Consider the complex energy d ue for each example From Table 44 we find that d ue l is 22 % of Table 44: Complex energy signal characteristics Example d ue = du e R + d ue l % complex 1 0.45+0 l0j 22 2 1.98+0.63} 32 3 9.91+0 31} 3
PAGE 78
72 d ue R for example 1 32 % for example 2, and only 3 % for example 3. Observe from Figure 49 that l ) A ) ::::: 1 0 ( /4), for all /4 for example 3. Figures 47 and 48 show that l u and l a are significantly different (but have close zero crossings near A op i) for examples l and 2. So the line integral along the negative a axis contributes significantly to the evaluation of J for examples 1 and 2 but not example 3. This results in the similarity of lu and l a for the last case. Clearly the relative values of du e R and du e l are application dependent. We have shown that even though l u and l a are different functions, their imaginary components have roots that are similiar near A op t Since l a can be computed directly from a sampling of d(t) it could be used to provide the optimal time scale estimate in lieu of lu,
PAGE 79
73 (order=0) 0.1 \ ~Cll 0.05 ..... +:, + \ 7 E \. \ ~ o ~ ~\ \ 1_ ,,, l ~~/c:.._; ~jj 0 05 .................. \ ,!.: '. < .... [; ,_/ tt1 : t i t I l i 0 .1 ~~~'~'~ 0 3 X 10 2 3 ).. (order=2) 4 5 ,/ 1 ... .. . .. 2 3 4 ).. (order=3) 5 1 ,.....,,....,....,.,~., 1 ~~~~.., u \ N 0 1 2 3 ; \ ,Cl! 0 5 y .... ...... .. t+rrr==I i _; 0 + ...... ......... ...... '.:: .:.. \ ;,.....i \ e. , ..=!""" \' . I E ! \ \ \ .. \ 1 .......__....___..___.____,~__._. _.__. 0 5 l l \ 1 ....____.~~____,~____, ~____, 0 2 3 ).. 4 5 0 2 Figure 47 : J 2Im[lu], and 2Im[J 0 ] for example 1. 3 4 5 Table 45: Performance function characteristics for example 1. A o p, An zc u Anz c a l (A 0 p,) in dB l (An zc u) in dB l (Anzca) in dB 2.47 2.48 2.48 24.7 24.6 24.6 1.60 1.58 1.49 33.8 33.7 31.1 2.31 2.31 2.29 49 9 49 9 49.5 1.92 1.92 2.62 60.7 60 7 54.6
PAGE 80
74 (order=0) (order=1) 0 4 ~~~~~ ~Cll 0.2 \ ..... 7 ', E \ o ~\~2.. _ : : ::: \ , 0.2 .... :. ..___ __ _ ,.,. 0.3 \ 0 2: \ Cll : 7 0 1 \ t++4i E . ': ~~+ i __ 1 ~; 0 .. \ ~ .. :~:; : t ~~ ~~i ~::: : .. .. ..;; .., ~ : ; : t : ~ :.: : : : : : : : : :: : :: 0 1 P., ___, +ii / 7 0 2 \ .... i ;/ ++1 '. 0.4 '~'' 0 2 4 6 0 2 4 6 0.1 A (order=2) 7 ro 0 05 : t, ;;r::;::;:::: = =1 E o .. \ ...... / ,r~: : : ~ J :: ~ """ : : ::; ;~:~::: ~ ... .. ... 0 05 ... i .... .1~,rt1 7 \ii \ / t ; 0 1 ~~~~ 0.1 I A (order=3) 7 ro 0 05 ...... it1 E ,.. a O \ r \~:J?~ :. : ~ :r llf:::. _ :;:=;:;_;;_ ,! ;.~;,,... ~.. :::_._ ..: ::::~ ~~ c. :.:;,,, : : :.. ::.:: , 0 05 \ }; \ /_i ; 0.1 ._______ i __ __ _.___ __ __, 0 2 4 6 0 2 4 6 A A Figure 48: J 2Im[Ju] and 2Im[l a ] for example 2. Table 46: Performance function characteristics for example 2. N Ao pt A. nz C u A n zC a J ( A. 0 p 1 ) in dB l ( A.n z c u) in dB J (A.n zca ) in dB 0 2.99 3 12 3.08 11.6 11.6 11.6 1 1.50 1.43 1.31 15.6 15.5 14.8 2 2.58 2.60 2 .19 20.8 20.8 19.0 3 1.81 1.78 3 16 24 8 24 7 2 2 .5 4 2.45 2.47 2.40 29.3 29.3 2 9 .2 5 1.93 1.9 2 1.86 33.3 33.3 3 2 8 6 2.39 2 .40 2.14 37.7 37 7 34. 8 7 2 00 2.00 1.96 41.8 41.8 41.5
PAGE 81
75 (order=O) ........_____: 0 8 +t1 ~C'tl l o 6 .....=,~ 0 4 +t1 2. f 0 2 +t1 , 0 2 0 , C'tl 0 4 r+J I 2 A ( order=2 ) 3 4 .....=,~ 0 2 +t1 0 ,(~ \ r,,,<.,,~ V ] \ ~,, : ,_: ..' I 0 2 0 2 4 6 A (order=1 ) ~C'tl l o 6 .....=,~ 0.4 t+1 2. f 0 2 , 0 ... ~ \;: : 17 .,.[ >< "' 0 2 ,..__ __ _.__ __ _._ __ ___, 0 2 4 6 A (order=3) 0 6 ~~~~ , C'tl 0 4 ;t,c;j f __;::, 0 2 2. .f :\ f \ i ' ' ,. , o ryt ...... ; <7 .............. J' .,,. \ ;': 1  l = ~ .; : :. ~ : ; , 0 2 ,..__ __ _.__ __ _._ __ ___, 0 2 4 6 Fi g ur e 49 : J 2 Im[l u ] a nd 2 lrn[l a ] fo r exa mpl e 3. T a bl e 47 : P erfo rm a n ce fun c ti o n c h arac t eris ti cs for exa mpl e 3 N A. o p t A. n z cu A. nzca l (A o pr ) in dB l (A. nzcu ) in dB l (A.n zca ) in dB 0 0.82 0.34 0 55 0 5 0 2 0.4 1 2.30 2.46 2.43 2 7 2 7 2 7 2 1.29 1.2 1 1. 2 1 5 2 5 1 5 .0 3 2 .1 2 2 1 2 2 13 9 9 9. 9 9 9 4 2 8 7 2 89 2 89 1 2. 7 1 2 .7 1 2 .7 5 3 .5 6 3 59 3.60 14 1 1 4 .1 14 1 6 4.2 1 4 22 4 22 1 4 7 1 4. 7 14.7 7 4 .82 4 76 4. 79 14 9 14 9 14 9 8 2 24 2.24 2.22 17 1 17 1 17.0 9 2.62 2 64 2 60 1 9.6 19 5 19. 6 10 2 2 5 2 24 2.28 22.9 22 9 22 8 11 2.5 7 2 56 2. 6 2 2 5.5 2 5 5 2 5 .3
PAGE 82
76 4.3.4 Time Domain Synthesis Using GLFF Networks For the discretetime setting one of the fundamental limitations of an FIR model is its lack of memory depth; namely its inability to approximate (with low order) a system with an infinitely decaying impulse response. GLFF networks do not have this problem due to the presence of local feedback in their structure. How well a Laguerre network approximates the impulse responses of the three example systems given by equations (4.86) (4.87) and (4.88) is now illustrated. Here the system impulse response d(t) (h(t) = d(t) since the input is an impulse), the Laguerre network output when using the optimal time scale y o(t), and the Laguerre network output when using the nearest ZCu time scale yz (t) are plotted. In example 1 the impulse response decays rapidly and has no significant oscillation. y o(t) and yz (t) are plotted in Figure 410 using a 1 s t order Laguerre model. In example 2 the impulse response is also rapidly decaying but with a slight oscillation. Plots of y o(t) and yz (t) are given in Figure 411 for a 5 th order model. Example 3 is a system having a slowly decaying impulse response with a strong oscillation. Plots of y o(t) and yz (t) are given in Figure 412 for a 10 th order model. This oscillation is difficult to approximate with a Laguerre model since only a single real pole (time scale) is used in the approximation. Increasing the order of the network could further reduce the error. Another option is to use a Kautz model with a pair of complex conjugate poles For all three examples we observe that y o(t) and y z (t) are inseparable. The point to be made in this section is that when using A obtained from the appropriate zero
PAGE 83
77 Example 1: d(t),yo(t),yz(t) 2 1.5 1 0.5 Or====t 1 2 3 4 5 Figure 410: Impulse response and approximations using a 1 st order Laguerre model. 4 3 \ 2 1 Example 2: d(t),yo(t),yz(t) Of__:\== ...=== = t 1~ 3 4 5 Figure 411: Impulse re s ponse and approximations using a 5 th order Laguerre model. Example 3: d(t),yo(t) yz(t) 5 4 3 2 1 OH~+~/L.:._======~t 1 7 8 2 3 Figure 412 : Impulse response and approximations using a 10 th order Laguerre model.
PAGE 84
78 crossing of Im[l u ] (or Im[l a ]) as the time scale for the model the same result can be achieved as when using /4 obtained by searching l for the global minimum. Also the squared error performance function can always be computed using l = 2Re[l u ]. 4 3.5 Conclusion of the Complex Performance Function For systems composed of sets of real and/or complex conjugate pair poles a new complex performance function l u, was derived when using a GLFF network as a system model. The imaginary component of l u has zero crossings near the stationary points of the squared error performance function l. Experimentally it is determined th a t as the order of the model is increased the zero crossings of Im[lu] near minima of l converge to these stationary points one of which will be the optimal time scale A op t Rather than searching the nonconvex function l for a global minimum locate the roots of Im[l u( A)] Ar and evaluate 2Re[l u( A r) J. One of these will be close to A 0 p 1 ; namely ( 4.114 ) where A r = {roots of Irn[lu(A)]}. In addition, we have defined the new analytic performance function and derived its closed form solution. It was shown that l could be written in terms of components o f lu and l a Expressions for the imaginary components of lu and l a were also given. The value of these functions in determining the optimal time scale of GLFF kernels was illustrated using three example systems The important equations from the abo v e sections are now summarized:
PAGE 85
79 Squared Error Performance Function l ( d e wT w )/ d e {Squared error components} J = ( du e R 2w ~R wuR) / d ueR { Complex error components} (d ae R 2 w;RwaR)/da e R {Analytic error components} Complex Performance Function Analytic Performance Function where d e = 2 ~ f D(s)D(s)ds w = 2 ~ f D(s)(s /4 )ds C C du e = 2 ~ f D(s)D(s)ds Wu= 2 ~ f D(s)(s,A)ds c. c. 00 00 d ae = 2 ~ J D(jw)D(jw)dw w 0 = 2 n f D(jw)(jw,A)dw 0 0 (4.115) (4.116) ( 4.117) (4.118) If D(s) is known analytically or in closed form we can employ the residue theorem to compute d e du e, d ae, w, Wu, Wa. If these functions must be determined numerically then the analytic performance function can be easily computed. The only requirement is a sampled data set from an impulse or step response test. In addition, J can always be computed from l u or l a. This illustrates that the squared error solution 1s completely embedded within the complex and analytic performance functions.
PAGE 86
80 4.4 Conclusion (Optimal Parameter Selection) This chapter has addressed optimization techniques for GLFF models, M(p ). Optimization of each element of the parameter set p = { N w A} is handled separately. Model dimension element optimization generally involves increasing N until a prescribed error tolerance is achieved. Monitoring the rate of convergence for several GLFF kernels could also aid in choosing the model yielding the smallest N. The optimal weight vector minimizes the integrated squared error and is found by solving a set of WienerHopf equations. A direct integral evaluation method was also given that takes advantage of properties of the model kernels. The weight vector solution is a function of the time scale vector X Optimizing A is the most difficult task because the minimum squared error solution is nonlinear w.r.t. these parameters. A number of optimization approaches exist in the literature and they were discussed. Finally, a new method was proposed that extends the squared error analysis to the complex domain. Simply locating zero crossings and choosing the one that minimizes J gives the optimal X Examples were provided to demonstrate the theory.
PAGE 87
CHAPTER 5 APPLICATIONS Applications of the theory presented in the prev10us chapters will now be discussed In this thesis GLFF networks will be used to approximate sampleddata systems. In particular the step response of a system used in laboratory test and evaluation of guided weapon systems will be modeled. Other areas where GLFF models have proven useful are first surveyed The GLFF model class encompasses many structures that have been independently formulated and studied in the past. This includes the tapdelay line Gamma Laguerre Legendre Jacobi and Kautz models defined in both discretetime and continuoustime domains. The unification of these structures is one contributory component of this thesis; however it is worth reviewing the previous history of these models to gain insight into the many applications that they have been able to address. Network synthesis. The earliest work dealt primarily with orthogonal model structures in an effort to synthesize networks and transient responses of linear systems [Arm57], [Arm59], [Kau54], [Hea55] [Bru64], [Clo65] [Kin77], [War53] In [Lee33] and [Hug56] synthesis of electrical circuits was the focus. System identification and modeling Similarly as the tools became more developed the nomenclature and approach shifted to system identification and system modeling [Bod94] [Eko94] [Kin79] [Lin93], [Mas90] [Mas91] [Nur87], [Nin94]
PAGE 88
82 [Oli87], [Oli94], [Per91] [Sil94] [Sil95], [Van94], [Wah91], [Wah94] [Moh88]. The basic principles are the same as in network synthesis only now the modeling effort involves more than just the network impulse response A further refinement includes modeling noncausal systems [Eko95] nonlinear biological system s [Mam93] probability density functions [Pab92] [Tia89]. State space formulations appeared [Dum86] [Gus93] [Heu90] [Heu95] and modeling using step response data [Wan94b] and frequency response data [Clu92], [Wan95] have been considered. Function approximation. A representation of an underlying system is not always desired. Various GLFF models were also useful in function approximation signal representation time series analysis [Bro65] [Cel95] [Con94], [Cle82] [Hwa82], [Kor91], [Mar69] [Mai85] [McD68], [Ros64] [Sch71], [Wah93], [You62b], [Yen62]. Particular examples include representation of seismic data [Dea64] and electrocardiograms [Y ou63]. Speech. Applications of GLFF models used in speech synthesis/analysis [Man63] and speech compression [AIJ94] have appeared. Control theory. GLFF models are also useful in control theory Examples include adaptive controllers [Dum86], [Dum90], [Zer88a] adaptive predictor controllers [Els94] [Guo94] optimal control [Hag91], [Nur81], robust control [Clu91] and robust stability tests [Aga92]. Stochastic processes. Stochastic concepts have been addressed including approximating cross correlations [den93a], correlations and power density spectra [Lee67], and probability density functions [Tia89], [Pab92]
PAGE 89
83 (Adaptive) filter theory GLFF models have been used as filters [Kab94] including adaptive filters [den93a], [den94], [Per88], [Pri93] Systematic procedures to design FIR filters have been introduced [Fri81], [Mas90], [Mas91]. Recursive identification of time varying signals [Gun90] 2D filtering [Tsi94] and hardware implementations have also been studied [Ji95]. Other applications. Other examples where GLFF networks have been used include solutions to boundary value problems [lra92], separation enhancement, and tracking multiple sinusoids [Ko94] and approximation of time delays [Lam94], [Par91], [Moh95]. Although the list of applications above is fairly extensive, it is by no means exhaustive Nevertheless the applicability of the GLFF modeling theory developed in this research is clearly illustrated. 5.1 Examples of GLFF Networks Using Simulated Systems The application focus of this thesis is system modeling ; namely it is desired to approximate (model) LTI systems using GLFF networks The objective is to obtain an approximation to the system impulse response. Quite often this is accomplished by optimizing the parameters of the model using step response data. In this section several GLFF network realizations are constructed that approximate simulated LTI systems. Optimization is performed using the step response of the desired systems Real (nonsimulated) systems and their associated step response data will be considered in Chapters 5.3 and 5.4.
PAGE 90
84 Three simulated systems will be studied : 1 s t order 2 nd order, and 4 th order. Example 1: (1 s t order) 1 H(s)=, r=l rs+l Example 2: (2 nd order) oi H(s) = 2 ; n 2 m 11 = 1 ; = 0.4 (underdamped) s +2 (J)ns+m n Example 3: (4 th order) where 2 HI ( s) = ? (JJ nl 2 (JJ n l = 1 ;I = 0.4 ( underdamped) S + 2;1(J) n lS + (J)nl (JJ 2 H 2 (s) = 2 n 2 2 mn 2 = 2, ; 2 = 0 .2 5 (underdamped) s +2; 2 mn 2 s+mn 2 (5.1) (5.2) (5.3) (5.4) In each case the step response is generated and sampled. The sampled step response data is modeled with Laguerre Gamma, and Legendre models. The results are given below in Figures 51 through 56. In each case the optimal time scale parameter was determined and used to compute and plot the step response and impulse response for 6tap (N = 5) networks. Figures 51 52, and 53 compare the system and GLFF network step responses. Figures 54, 55, and 56 show the resulting impulse responses. The optimal time scale parameters of the 6tap models and associated minimum squared error, J(A o pi) are provided in Tables 51 52, 53.
PAGE 91
0 9 0 8 0 7 0 6 = .:2 0 5 = er, 0.4 0 3 0 2 0 1 0 0 1.4 1 2 0 8 0 6 = er, 0.4 0 2 0 0 2 0 85 0 5 1 5 2 2 5 I (sec) 3 System + + + + Lag uerre mode l Gamma model x x x x Leg endre model 3 5 4 Figure 51: 1 s t order step response and model data. 2 3 4 5 I (sec) 6 Sys t em + + + + Laguerre model Gamma model x x x x L egendre model 7 8 Figure 52 : 2 nd order step response and mode] data. 4 5 5 g 10
PAGE 92
1 6 1.4 1 2 0 8 = >:= Ol 0 6 0 4 0 2 0 0 2 0 0 9 0 8 0 7 0 6 3o 5 = ..c 0.4 0 3 0 2 0 1 0 0 86 2 3 4 5 t (sec) 6 Sys t em + + + + L ag u e r re model Ga m ma model x x x x L ege nd r e m o d e l 7 8 Figure 53: 4 th order step response and model data. 0 5 1 5 2 2 5 t (sec) 3 Sys t e m + + + + L agu e rre mode l Gam m a mode l x x x x L ege n d r e model 3 5 4 Figure 54: 1 s t order impulse response and model data. 9 10 4 5 5
PAGE 93
0 8 0 7 0 6 0 5 0.4 = >:0 3 = ..::. 0 2 0 1 0 0 1 0 2 0 0 8 0 7 0 6 0 5 0.4 = >:0 3 = ..::. 0 2 X 0 1 0 0 1 X 0.2 X X 0 87 2 3 4 5 t (sec) 6 System + + + + Laguerre model Gamma model x x x x Legendre model 7 8 Figure 55 : 2 nd order impulse response and model data. Sy s tem + + + + Laguerre model x Gamma model x x x x Legendre model X X 2 3 4 5 6 7 8 I (sec) Figure 56: 4 th order impulse response and model data 9 10 9 10
PAGE 94
88 Table 51: 1 s t order model results Model N A ov t l(A ov 1) in dB FIR 101 13.7 Laguerre 5 0.88 57.1 Gamma 5 0.05 00 Legendre 5 /4 = 0 05i op t ,i e 56.0 Table 52: 2 nd order model results. Model N A ov t l(A o ot) in dB FIR 101 13.0 Laguerre 5 0.95 57.1 Gamma 5 0.08 40 7 Legendre 5 /4 = eo 03; o pt 1 37.5 Table 53 : 4 th order model results. Model N A ov t l(A o 0 r) in dB FIR 101 10.9 Laguerre 5 0.89 37 7 Gamma 5 0 09 31.2 Legendre 5 /4 = 0 0 3i o pt ,i e 28.8
PAGE 95
89 The three examples considered above are typical of the subsystem makeup of guided weapon systems. With fairly low order models (N=S) the GLFF networks have demonstrated (in the simulated examples) they can approximate these systems well. This includes systems having slowly decaying impulse responses and oscillatory modes. The most difficult task of the modeling process is locating the optimal time scale. This was done using a numerical search technique in the above examples For a given model if the resulting error is too large the model order can be increased. This process can be continued until the prescribed error tolerance is reached 5.2 Laboratory Testing of Guided Weapon Systems Guided weapon systems are as important in a military tactical rruss10n as the aircraft and ships that carry them. They are the tools with which a pilot can carry out both offensive and defensive plans. As such, their capability and reliability must be well understood so that they can be operated within the envelope for which they were developed. To ensure their success it is necessary that the operational envelope is understood. This knowledge is acquired through extensive test and evaluation of these systems. The objective is to produce a smart dependable weapon that operates as expected when launched from an aircraft or ship. Test and evaluation (T&E) of guided weapon subsystems can be broken down into three phases [Eic89]. This breakdown is typical when a new weapon or modification of an existing weapon is under consideration. These phases are
PAGE 96
90 1. Research/Concept Exploration a) analytic study b) openloop laboratory and field tests 2 Simulation a) software simulation and modeling b ) HITL simulation 3. Flight Test a) captive carry b ) launch Research and concept exploration a re the very beginning of the T &E effort. Proposals must be carefully considered to determine which (current) weapons could possibly meet the objective being considered. Can an existing weapon be modified or does it call for a new design? All these scenarios must be analytically explored. Flight tests are the final hurdle. It is here where the weapon system s performance is measured against its e x pectation Before flight testing however many problems can be worked out by ev a luating a weapon system in a simulated environment created in the laboratory. There are many methods for evaluating the performance of components of a weapon system such as parametric tests and openloop tests. Another class of tests studies the behavior of these systems in a clos e dloop configuration and involves a simulation of the realworld environment and the weapon s flight in this environment. This class of tests can be broken down into two types: alldigital and hardwar e inth e loop (HITL). In an alldigital closedloop flight the entire weapon system is simulated mathematically. If the weapon subsystem models a re accurate alldigital simulations can give us a good idea of how the weapon will perform as well as what its limitations and vulnerabilities are. Once subsystems have been hardware fabricated these can replace the models in the
PAGE 97
91 closedloop simulation Ultimately it is desired to have as much of the weapon hardware in the loop as possible. 5.2.l Objectives Several reasons why laboratory HITL testing is an integral part of a weapons development test program are now given. 1. Asset Availability During the development of a weapon system, often answers are needed to key questions about the system to determine the limits of its capability. Rather than fabricating an entire system and conducting a flight test many times these questions can be answered by using simple prototypes or only subsystem components. This precludes having to completely fabricate a system and conduct open air flight tests. Flight tests are expensive and consume lots of resources (fully developed weapon system aircraft/ships openair range time/space money). Here is where laboratory testing offers a great advantage. 2. Planning and Analysis Requirements As a missile system reaches development stage and can be flight tested valid test conditions (test matrix) must be formulated. Having laboratory results (alldigital and/or HITL) as a guide, the test matrix can be formulated with intelligence and unrealistic scenarios can be avoided. Laboratory tests can not only help define the launch envelope of a weapon but they can be used to predict its performance when operating under valid conditions. 3. Repeatability Another advantage to the laboratory environment is repeatability testing. Once a weapon system simulation has been developed, it can be used to conduct hundreds and thousands of tests. Weapon performance when launched under varying release conditions and "what if' scenarios can be examined. The laboratory environment is a major contributor to this kind of testing because real assets do not have to be expended to conduct many tests In live fire" testing assets can usually be used only once. 4. Security/Safety Finally there is the security advantage. As long as testing is being conducted in a laboratory, external security vulnerability is reduced. Range safety costs and concerns that occur in openair testing are also eliminated.
PAGE 98
92 To accomplish valid and accurate laboratory testing, careful attention must be paid to the test configuration. For HITL tests it must be ensured that the laboratory equipment can recreate the dynamics and conditions that will be experienced in the real world. Often it is easy to switch between configurations involving hardware and software. This requires accurate software models of the hardware being simulated. This is where GLFF networks discussed in this work enters the application front. To gain a better understanding of the systems the GLFF networks must model, a more detailed breakdown of a typical missile system is now given. The functional components of the missile are discussed as well as how they are integrated into a HITL test configuration. The GLFF networks can then be used to approximate some of these components. 5.2.2 Functional Breakdown of a Missile System This section provides a discussion of components (subsystems) most often found in missile systems. Figure 57 below illustrates a highlevel functional break down of a typical missile system [Gar80]. Seeker/ Guidance Autopilot Actuator / Control Surfaces Figure 57 : Missile functional block d iagram.
PAGE 99
93 Seeker/Guidance Modules. It can be said that the seeker is the eye and the guidance module is the brain of a weapon system. It is here the missile trajectory is commanded and controlled. For missiles that are considered lockon after launch the guidance module may begin in a target search mode. In this case it drives the seeker to search for a target in some preprogrammed pattern. During this phase commands may be sent to the autopilot to have it execute a midcourse maneuver such as pitching the weapon up in order to gain altitude for target search and future end game engagement. If target energy is found the guidance module may transition into an acquisition mode to narrow the search field of the seeker and attempt to select a single target. At this point it begins target track and will try to keep the seeker pointed in the direction of the target. During track mode guidance (acceleration) commands are sent to the autopilot in accordance with the guidance law being executed (proportional or pursuit navigation) All the while this module performs seeker head stabilization and directional pointing to maintain a continuous track of the target until impact. Autopilot. The job of the autopilot is to maintain a stable flight and deliver on the commands received from the guidance module To perform its task the autopilot accepts input from various sources and sensors The three most common are: guidance module accelerometers gyros The guidance module tells the autopilot that a certain amount of acceleration is required to maintain the desired trajectory. Upon receiving a command for a specific acceleration, the autopilot issues commands to the actuator to drive the control surfaces (fins or canards) in the proper direction. Measurements are taken from accelerometers to determine whether the weapon has achieved the desired acceleration. Gyro rotational rate and angle measurements are
PAGE 100
94 taken to ensure proper damping and to keep the motion of the missile in a stable condition Actu a tor and Control Surfaces. The purpose of the actuator is to accept commands from the autopilot and in tum drive the control surfaces to achieve the dynamics needed. As the fins are controlled into the wind the missile body will rotate and accelerate hopefully in a stable manner sufficient to maintain a target track. T rack a n g l es S t eeri n g co mm a nd s Fin c omm a nd s Ta, ge t l~Se e k e r ~ Fl Guidance I Autopilot ~ A c tu a t_o r s gn, '" re \~ .... L':'~C.,.~ b il .iC! '?. ~ ...i.... .... .......... . ........ ...... .... '. . ' : Fli g ht Motion Simulator ( FMS ) : Target Scene F i n p os iti o n s For ces m o m e nt s 1 , Kinematics ..__ _____ T_ ar ge _t ____ 1 Missile/Target Geometry T arge t p ro fil e Mi ss il e body a n g l es, rat es, acce l e ration s Figure 58 : Closedloop missile simulation. 5.2.3 HITL ClosedLoop Missile Simulations In this section HITL closedloop missile simulations will be discussed. Figure 58 above illustrates a typical test configuration Usually the seeker guidance module, and autopilot are hardware components and the functionality of these systems is the subject of the test. Special laboratory equipment is used to generate the target scene. This is dependent upon the kind of sensor in the seeker head. The remaining blocks are simulated in software. For a HITL closedloop simulation, the seeker guidance
PAGE 101
95 module and autopilot are mounted on a flight motion simulator (FMS). This FMS is commanded to recreate the body dynamics that the weapon would usually experience in a real flight. Clearly the ability of the FMS to create accurate body motion is critical to a valid HITL simulation. To understand the performance capability of the FMS an accurate software model of this system is desired. This model can be used to study the limitations of the FMS without subjecting the hardware to harsh experimentation. Also the alldigital simulation can be modified to provide for body rate updates either directly computed or filtered by the FMS. In this way effects of the FMS on the accuracy of the weapon simulation can be determined. A discussion is now given of the FMS system that will be modeled. Its functional composition is given along with its expected performance specification. 5 3 Modeling the Step Response of a Carco 5Axis FMS The flight motion simulator (FMS) under study is a CARCO model S458R5E fiveaxis system. This simulator is a precision electrohydraulic positioner consisting of three basic elements: a hydraulically driven flight motion table a control console with interface electronics and a hydraulic power supply. The FMS structure provides five axes of rotational motion, usmg prec1s10n electrohydraulic actuators. The inner, middle and outer gimbals of the threeaxis missilemotion table are designed to roll pitch, and yaw the weapon. The two large outer gimbals are used as target axes and are designed to provide azimuth and elevation maneuvers for the target system.
PAGE 102
96 The electronics console houses the electronic servo amplifiers power supplies and other control electronics needed to regulate the table motion. It provides for both analog and digital closedloop control. A highspeed parallel interface allows remote control and monitor of the FMS from a digital simulation computer. This digital interface allows for realtime DMA transfers at a rate of 1 frame per millisecond (1 KHz digital update ) The hydraulic unit provides the power needed for the high dynamic performance of the fivea x is table. It consists of a motor pump and reservoir. As delivered the FMS system has the following capability using a 150lb load 16 in diameter 30 long mounted on the roll plate. The FMS is illustrated in Figure 59. Figure 59: Carco 5axis FMS.
PAGE 103
97 The performance specifications for each gimbal axis are listed in Table 54. It is this closedloop system that will be modeled with several GLFF networks. Table 54 : Carco FMS performance specifications Roll Yaw Pitch Azimuth Elevation Axis Axis Axis Axis Axis Max position (deg ) Continuous Max velocity (deg/sec) 1 000 300 300 100 10 Max acceleration ( deg/sec 2 ) 20 000 15 000 15 000 12 000 12 000 A typical closed loop system ( plant plus feedback controller ) is shown in Figure 510 When identification of only the plant is of interest model design could be carried out with an openloop experiment. If the controller design is known it i s possible to model the plant in a closedloop experiment [W a n94a]. This is often more desirable than openloop system identification because it causes less disruption to the operation of the system since disabling the feedback is not necessary If the open loop system is unstable or poorly damped then openloop experiments can not be performed without concern for safety and impending system damage n .. . ) u + I: e Controller c Plant V .. ..... .. ...... .......... ... .. .............................................. .............. ........ ........ .............. ...... ............. .... .................. ......... ........... Figure 510 : Typical clo s edloop system In this thesis it is a model of the entire closedloop system that is of interest. This is often needed when the system is to be treated as a subcomponent of a larger
PAGE 104
98 configuration. In our case the system is a hardware device that is to be driven by a digital computer. Since the system is expensive prolonged and abusive use is highly undesirable. If it is possible to model the system in software then preliminary development and debugging of the computer algorithms and hardware interface can be performed using a digital system model. Once the software has been tested the hardware system can then replace the digital model. In this way hardware utilization is kept to a minimum. Hence an accurate model of the closedloop system is desired. 5 3 1 DiscreteTime GLFF Models of the FMS Figure 511 illustrates the general form of a discretetime GLFF network realization. Three transfer function kernels (Gamma Laguerre Legendre) will be considered for this application. x( k ) y (k) Figure 511: Discretetime GLFF network block diagram. The GLFF formulation of the discretetime Gamma Laguerre and Legendre kernels were given in Table 31. These kernels result in the network realizations illustrated in Figures 512 513 and 514. Note that in all cases H 0 1 ( z, A)= 1.
PAGE 105
x( k ) l 1 x( k ) 1 x( k ) l 99 Figure 5 12 : Gamma G LFF Form. Figure 513: Laguerre GLFF Form. F i gure 51 4 : Legendre GLFF Form. 1 I l ,l.,N z 1 xtv( k) y( k ) y( k ) y (k )
PAGE 106
100 5.3.2 Step Response Modeling Procedure As stated above, the inner three FMS axes carry a weapon under test and induce body pitch, yaw, and roll motion. A target may be presented on the outer azimuth and elevation axes. Each axis is independently controlled to produce the proper dynamic behavior that the weapon would sense if it were flying in the real world. Obviously there is much work required in developing a valid HITL weapon simulation. During the HITL simulation development cycle a software model of the FMS can be used in the place of the real system. This will eliminate unnecessary wearandtear on the hardware. Hence, a good digital model of the FMS is the goal. A procedure is now described that will employ GLFF networks to accomplish this task. The model parameters are optimized using measured step response data. This data was produced by subjecting the FMS to an analog step input. The input signal and FMS response were digitized and collected at a 1 KHz sample rate. The following analysis will concentrate on the pitch axis only. The same procedure was applied to the other axes. In addition to the step response data FMS pitch yaw and roll gymbal commands and responses were collected during a typical missile flyout engagement. The gymbal commands were also fed into several GLFF models. The model responses were compared to the FMS signals to characterize the accuracy with which the GLFF networks can approximate the FMS under normal operations. This engagement test was conducted using an AGM650 Maverick mjssile. The AGM650 is a short standoffrange launchandleave airtoground missile. It is a rollstabilized cruciformwing tailcontrolled system with an imaging infrared sensor. The scenario
PAGE 107
101 chosen launches the missile at approximately 200 meters altitude and 2000 meters down range from a stationary ground target. The results of this test are presented with the step response data in section 5.3.3 The approach taken to optimize a discretetime GLFF network parameter set using s tep response data of the FMS system H is illustrated in Figure 515. Since the goal is a network that models the impulse response of H via step response data a derivative is necessary ; hence precautions must be taken to ensure the step response g( t ), has lo w signalto noise ratio. Another approach is to transfer the derivative from the response signal to the network kernels This and other methods were discussed in Chapter 3 5. H g(t) s h ( t ) u ( t ) +1 s ~ t) M(p) Figure 515 : Step response modeling structure. Becau s e of the deri v ative of the step response data the parameter adaptation process must work to reduce the overall noise effects. The process involves selecting the parameter set p = { N w, A} in a twostage fashion. Step 1: Fix the order N and determine the time scale parameter A and weight vector w, that minimizes the integrated squared error J, of the output of the system and model im p ul se responses.
PAGE 108
102 Step 2: Let A adapt to minimize the squared error of the output of the system and model step responses. In this way the optimal weight vector represents the impulse response while the optimal time scale allows for a model match of the step response. In a noise free environment A would in fact be the optimal time scale minimizing the impulse response squared error. This method is now demonstrated. The test case illustrates modeling the FMS pitch axis. The set point (step input) is 1 volt (6 Degrees). The step input and pitch axis step response are shown in Figure 516. Carco FMS Pitch Axis Step Response 8 ,.rr,~rr~, 7 6 5 2 4 .. Ol 3 .. I 0 0 1 0.2 I I I I I I 0.3 0.4 0.5 0 6 0.7 0 8 0 9 t (sec) Figure 516: Pitch axis step response.
PAGE 109
103 5.3.3 Results (FMS Step Response Modeling) GLFF model parameters were optimized using the procedure described above for the first eight orders of Gamma Laguerre, and Legendre networks. Tables 56, 57 and 58 summarize the numerical results. Here, A 1opi is the optimal time scale from step 1. l 1 m i n is the minimum squared error resulting from the WienerHopf solution of the impulse responses, J min = J (A o pt). Likewise, A sapi is the optimal time scale from step 2 J Smin is the minimum squared error resulting from the adaptation of the step responses, l smin = J(A s o pt) As can be seen from the tables all three kernels yield similar capabilities. Because of the completeness of the GLFF structure, increasing the order of the network will always produce either the same or smaller ]/m in' For each (Gamma, Laguerre, Legendre) network squared error performance functions (Figures 521 526 and 531) step responses (Figures 522, 527, and 532), and impulse responses (Figures 523, 528 and 533) are plotted for models of order 1 to 8. For the 8 th order models, the step response and associated error (Figures 524 529, and 534), and flyout response and associated error (Figures 525 530 and 535) are plotted. Finally the 8 th order step responses are scaled to degrees and plotted in Figure 536 along with the scaled errors in Figure 537. Notice that the maximum error is no more than 0.15 for a 6 step (less than 3 percent). Once the time scale was found that produced the minimum squared error for the impulse response 1'min, it was allowed to adapt to search for the time scale yielding the minimum squared error for the step response. These are slightly different due to
PAGE 110
104 the noise in the impulse response Choosing J Smin = 30dB as the error specification a 4 th order Gamma 3 rd order Laguerre and 3 rd order Legendre model would be chosen. Figures 538 and 539 plot the transfer function magnitude and phase of the system and models. They match very well out to 30 Hz. The pitch axis frequency response under the specified load condition is 20 Hz The system response increase from 60 100 Hz is a result of this being a closedloop system. To provide a comparison of the GLFF networks with the conventional tapdelay line Table 55 lists the results of several FIR model performance parameters. Figures 517 through 520 plot the step response impulse response step response error and flyout response using a 201tap delay line. A large order FIR model is required to achieve good results. This is a disadvantage due to increased computational requirement as well as the system delay if the model had to operate in a realtime environment. Table 55: FIR optimization characteristics. Order ]! min ] Smin 51 3.2 10.4 101 12 3 13.1 151 14.4 19.4 201 14.9 21.9
PAGE 111
105 FIR Model Pitch Axis Impulse Response (order=201) 30 25 20 15 er u::: 10 .c ::c 5 0 5 10 0 1 0.2 0.3 0.4 0 5 0.6 0 7 0 8 0 9 t (sec) Figure 517: FIR network impulse response (201 taps). FIR Model Pitch Axis Step Response (order=201) 1 4 ,.........., 1 2 0 8 0 6 u::: en 0 4 0 2 0 0 2 .___ __._ __ _.__ __ .___ __._ __ _.__ __ ..__ _.._ __ _._ ___, 0 1 0 2 0 3 0.4 0 5 0 6 0 7 0 8 0 9 t (sec) Figure 518: FIR network step response (2 01 taps).
PAGE 112
106 FIR Model Pitch Axis Step Response (order=201) 1 5 ~~, ~ ~~~!i ~, ~ '.' ; ~~ ; 4 ; 1 T .. ;05"! ) o~~~~~~~~~ 0 1 0 2 0 3 0.4 0.5 0 6 0.7 0 8 0.9 FIR Model Pitch Axis Step Response Error (order=201) ! ! 0 1 ++;r1 0 05 +++++++1 .$ or ~::::1:== ====+=~,it7,,1 ; 0.05 ~ "w t ++t 0 1 1t t~yi _t= "=,,,..,~:::;;;;;:;;;;;;~:::::::;:::;=:t=;::;::;:::::: :;:I I i i i i 0 1 0 2 0 3 0 4 0 5 0.6 0 7 0 8 0.9 t (sec) Figure 519 : Step response and associated error for 201tap FIR network. FIR Model Pitch Axis Flyout Response (order=201) 5~~~~~~~ rn o1i. ::t+t++~ f 5 ~ !t 10 ++, ~ .. % 15 ... .. . ,.. ==+; 20 ''''"i ___ .,_i ___ .,__i __ ___, 3 4 5 6 7 8 9 10 FIR Model Pitch Axis Flyout Response Error (order=201) 2 ~~,:,,. ., .., ____L~ t+t+::;,, .,; g' ~__,, ~ 0 ~ = = = = = 1""" +t1 .._,. Q) 1 11+ ++++t 2 ''''''...__~ 3 4 5 6 7 8 9 10 t (sec) Figure 520: Flyout response and associated error for 201tap FIR network.
PAGE 113
aJ C 0 u C :::l LL (1) (.) C E 0 't: (1) a.. 0 (1) N co E 0 z 107 OT Gamma Model Pitch Axis Performance Functions (orders=1 to 8) 2 4 6 8 10 12 L..._ _____,1 __ _____.i;=L..L..__._ __ 1.... __ ___J,_ __ ....J 0 0 05 0 1 0 15 0.2 A 0 25 0 3 Figure 521: Gamma network performance functions. 0 35 Table 56: Gamma network optimjzation characteristics. Order A 1op1 ]! min A sopt 1 0.0135 2.9 0 0205 2 0.0380 5 9 0 0530 3 0.06 2 0 7 6 0 0260 4 0.0415 9.4 0 0420 5 0.0590 10.8 0 0405 6 0.0715 11.1 0 0295 7 0.0895 11.7 0 0415 8 0.1050 11.9 0 0500 0.4 ] Smin 19.6 22.5 29.9 33 0 33.0 36 0 40 7 41.7
PAGE 114
108 OT Gamma Model Pitch Axis Step Response (orders=1 to 8) 1 2 en i 0 8 2:0 6 CJ Ol 0.4 Ol 0 2 0 0.2 ,..___ ___._ __ __._ __ ....__ ___. __ __._ __ __,__ __ _.__ ____. __ __, 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 t (sec) Figure 522: Gamm a network step responses. OT Gamma Model Pitch Axis Impulse Response (orders=1 to 8) 30 25 20 15 .. 10 <( CJ 5 ..c 1: 0 5 10 15 20 0 1 0 2 0 3 0.4 0 5 0.6 0.7 0 8 0.9 t (sec) Figure 523: G a mma network impulse responses.
PAGE 115
109 OT Gamma Model Pitch Axis Step Response (order=8) 1 5 .r..r.~,, i 1 e j ~ ~ =~ ~ Io 5 .... _,_/ ____ 1t1t i ) OL...,.L __ ....__ __ __.__ __ __.__ __ _.__ __ __._ __ ___._ __ ____. __ __, 0 1 0 2 0 3 0 4 0 5 0.6 0.7 0 8 0.9 OT Gamma Model Pitch Axis Step Response Error (order=8) ! ! 0 1 1++++t .. 0 05 .$ ..... Q) 0 05 1+++++++t 0 1 1+++ i ~ i +i + i + t 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 t (sec) Fi g ur e 524 : G a mm a order 8 s tep res pon s e a nd as sociated e rro r. OT Gamma Model Pitch Axis Flyout Response (order=8) 5~~~.~~r~ 0 ti.~ ~+++tt, 2.. .............__ 5 ~ t1 f C!J 10 r+++ _.::,,._ =+t+, 'O 2 15 r+++t=::;, =~ =!=== ==t 'O 20'~~_.._ ___ ....__ __ ___._ ___ __.__ __ ___. 3 4 5 6 7 8 9 10 OT Gamma Model Pitch Axis Flyout Response Error (order=8) 2 ..r.~.,., 0) Q) 2.. 0 i+.::::a'+;+;.... ,...,,_ .. ~ +1 ..... Q) 1 +++++2 '''''''' 3 4 5 6 7 8 9 10 t (sec) Fi g ure 52 5: G a mm a ord e r 8 flyout res pon se a nd ass oci a ted e rr o r.
PAGE 116
110 OT Laguerre Model Pitch Axis Performance Functions (orders=1 to 8) 0 2 co :s C 0 4 5 C :::) LL Q.) 6 (.) C E .... 0 't: 8 Q.) a.. "O 10 E .... 0 z 12 1 4 .___ __ .___ ____..__ ____..__ ____. __ ____. __ ____. __ ____. __ ___, 0 6 0 65 0 7 0 75 0 8 A 0 85 0.9 0 95 Figure 526 : Laguerre network performance functions. Table 5 7 : Laguerre network optimization characteristics Order A 1op1 ]! min A sopt 1 0 9625 5.8 0.9475 2 0 9390 7.5 0.9745 3 0 9595 9 3 0.9585 4 0.9415 10.8 0.9605 5 0.9290 11.1 0.9710 6 0.9120 11.6 0.9585 7 0.8950 11.9 0.9500 8 0.9375 12.1 0.9505 ] S m in 22 5 29 9 32 9 33 0 35.7 40.3 41.6 41.7
PAGE 117
1 2 0 8 ......... c:5' 0.6 <{ ...J 0) 2 0.4 0) 0 2 0 111 OT Laguerre Model Pitch Axis Step Response (orders=1 to 8) 0 2 .___ __._ __ __._ __ ....._ __ .__ __.._ __ ___._ __ _ _,__ ___, 0.1 0 2 0.3 0.4 0 5 0.6 0 7 0.8 0.9 t (sec) Figure 527: Laguerre network step responses. OT Laguerre Model Pitch Axis Impulse Response (orders=1 to 8) 30 25 20 15 ......... 10 CJ <{ 5 ...J ......... ::c0 5 10 15 20 0 1 0 2 0 3 0 4 0 5 0.6 0 7 0 8 0 9 t (sec) Figure 528: Laguerre network impulse responses.
PAGE 118
11 2 OT Laguerre Model Pitch Axis Step Response (order=8) 1 5 ,.....rr, I 1 e i . . ..... .. .... .. ....,..._ ... ..... ... =. ............................ . rs~ ; 0 ''"'''....L...''''' 0 1 0 2 0.3 0 4 0 5 0 6 0 .7 0 8 0 .9 OT Laguerre Model Pitch Axis Step Response Error (order=8) ! ! 0 1 ~ ...++++r+~ 0 05 t+++ttrt .$ a, 0 05 ++t+ttr; 0 1 t ++ i + i 1i i ~ 0 1 0 2 0 3 0 .4 0 5 0 6 0 7 0.8 0 9 t (sec) Figure 529: Laguerre order 8 step response and associated error. OT Laguerre Model Pitch Axis Flyout Response (order=8) 5 ,.r.r.,, o~ t,..;~ c:+++++~ 5 ~ +__.o, kt+i +l 0 "::i 10 f++..' :,,,,..... ac+++~ 0 S15 ~ ++++~~= .,,,,. ~====I 0 20 ~~~i ___ ....__ __ __, 0) (l) 3 4 5 6 7 8 9 10 OT Laguerre Model Pitch Axis Flyout Response Error (order=8) 2,.r.~~~ 0 r...i"'<;;. ;;r4t,...1+~\.+I (l) 1 +2 ~~~~ 3 4 5 6 7 8 9 10 t (sec) Figure 530: Laguerre order 8 flyout response and associated error.
PAGE 119
. c:o :E, C 0 D C ::, u. (l) (.) C <1:l E 0 't: (l) a.. "O (l) N cij E 0 z 113 OT Legendre Model Pitch Axis Performance Functions (orders=1 to 8) 2 4 6 8 10 12 .__ ____ 1. _____ __,__ _____ __,__ _____ _, 0 8 0 85 0 9 A 0.95 Figure 531: Legendre network performance functions. Table 58: Legendre network optimization characteristics. Order A 1op1 ]/ min A sopt 1 0.9729 5 6 0.9615 2 0 9642 6 9 0.9863 3 0.9822 8.8 0.9812 4 0 9775 10.3 0 9770 5 0 9846 10.6 0 9819 6 0.9891 11.0 0 9883 7 0.9690 11.4 0.9866 8 0 9854 11.8 0 9848 ] Smin 22 2 29 2 33 l 34 3 35 l 38.8 40.8 42 .2
PAGE 120
1 2 0 8 l ..... 0 0.6 w _J 0) S 0.4 0) 0 2 0 114 OT Legendre Model Pitch Axis Step Response (orders=1 to 8) 0 2 ~~~~~~~~~~ 0 1 0.2 0.3 0 4 0 5 0 6 0 7 0 8 0 9 t (sec) Figure 532: Legendre network step responses. OT Legendre Model Pitch Axis Impulse Response (orders=1 to 8) 30 25 20 15 ..... 10 (.!) w _J 5 .c ..... 1: 0 5 10 15 20 0 1 0 2 0.3 0 4 0.5 0 6 0.7 0 8 0 9 t (sec) Figure 533: Legendre network impulse responses.
PAGE 121
115 OT Legendre Model Pitch Axis Step Response (order=8) 1 5 ~~ ! ~ i,, ~ i,,,',,. ~...., 1/ ;os) o.........._ __ ..__ __ ....__ __ _.__ __ _.__ __ __._ __ ___.__ __ _._ __ _____. 0.1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 OT Legendre Model Pitch Axis Step Response Error (order=8) ! ! 0 1 .0 05 ~ ++41+++.$ ..... Q) 0 05 ~ ++41+++0 1 1~ : : ~ , : +1 0.1 0 2 0.3 0 4 0.5 0.6 0 7 0 8 0 9 t (sec ) Figure 534: Legendre order 8 step re s ponse and associated error. OT Legendre Model Pitch Axis Flyout Response (order=8) 5~~~~.~~~ B 0t~ =++++t; Q) ":;:::5 rt~ ttt1 5 10 1+++' '=++ti "O r_ 215 ~ ++++~=::.::;...~ :::t = = = = ==l "O 20 L__ __._ ___ ..__ __ ___.__ ___ ....__ __ ~ i ~' 3 4 5 6 7 8 9 10 OT Legendre Model Pitch Axis Flyout Response Error (order=8) 2~~~~.~~~ 0) Q) ~o~~ ~~7~~~~ = ~ ~ + ==~==~~== ~ 1 1++++f2 L__ __._ ___ ......___ __ L.. ___ ........_ __ ___._ ___ _.___ __ _____. 3 4 5 6 7 8 9 10 t (sec) Figure 535 : Legendre order 8 flyout response and associ a ted error.
PAGE 122
7 0) 6 Q) 25 0 w _J ~4 <( 0 3 2 ~2 _J 0) ~1 0 116 Laguerre(8) / Gamma(8) / Legendre(8) Model Step Responses 1 ......_ ___._ __ __.__ __ ....____ ____......._ __.__ __ _._ __ ....__ ____.......__~ 0 1 0 2 0 3 0.4 0 5 0 6 0.7 0.8 0 9 t (sec) Figure 536: FMS and 8 th order Gamma Laguerre Legendre network step responses. Laguerre(8) / Gamma(8) / Legendre(8) Model Step Response Errors 0 2 O') Q) 0 0 <( _J Q) 0.2 0 1 0 2 0.3 0 4 0 5 0 6 0 7 0 8 0 9 0 2 O') Q) 0 <( 0 Q) 0 2 0.1 0.2 0.3 0.4 0 5 0.6 0 7 0 8 0 9 0 2 O') Q) 0 c5 w _J Q) 0 2 0 1 0 2 0.3 0 4 0 5 0 6 0 7 0 8 0 9 t (sec) Figure 537 : 8 th order Gamma Laguerre Legendre networks step response errors
PAGE 123
10 m15 :s Q) 20 c: 0) ca 25 30 35 117 Laguerre(8) I Gamma(8) I Legendre(8) Transfer Function Estimates FMS ,,.,,. Le re Gamma 40 ........._ ____ ........._ ___ ~........._~~~........._~ 0 20 40 60 80 1 00 Frequency (Hz) Figure 538: Transfer function magnitude estimates for FMS and 8 th order networks. 20 40 60 en 80 Q) "O ';i;' 100 en ca ..c c.. 120 140 160 180 Laguerre(8) I Gamma(8) I Legendre(8) Transfer Function Estimates 200........._ ___ ~........._~~~~ 0 20 40 60 80 100 Frequency (Hz) Figure 539 : Transfer function phase estimates for FMS and 8 th order networks.
PAGE 124
118 5 3.4 ContinuousTime GLFF Models of the FMS In this section the GLFF modeling capability will be further illustrated by constructing a continuoustime Laguerre model of the PMS using the sampled step respon s e d a t a. The general form of a continuoustime Laguerre network using the Laguerre kernel is given in Figure 540 x(t) 1 1~1 s + .1. s /4 s + /4 s + .1. s /4 s + /4 y( t ) Figure 540: Continuoustime GLFF network block diagram. The optimization procedure employed is that which was derived in Chapter 4. Namely the analytic performance function l a, will be computed and the roots Ar of 2Im[l a ] located. A zero crossing detection algorithm was developed to numerically evaluate 2Im[l a ] for these roots The squared error is calculated for each root and the time scale producing the minimum is selected as the nearest ZC a (5.5) Figure 541 depicts the squared error performance functions for Laguerre networks from order 0 to 8. For each order network Table 59 provides a comparison of the optimal time scale A op t, the nearest ZC a time scale An zca and their respective squared
PAGE 125
119 error values. Figures 542 and 543 show the squared error and imaginary component of the analytic performance function Recall J is embedded within l a and is calculated by ( 5 6 ) Finally Figure 544 graphs the step responses and step response errors for 8 th order models using A o pt (light curve) and A nzca (dark curve) as the time scale. Figure 545 graphs the responses and associated errors of these models when driven by the AGM650 flyout commands
PAGE 126
120 CT Laguerre Model Pitch Axis Performance Function (orders=0 to 8) 0 1 co 2 C 0 3 t5 C ::J LL 4 Q) (.) C ctl 5 0 't: Q) 6 a.. 0 Q) N 7 ct1 E 0 z 8 9 10 0 20 40 60 80 100 120 140 Ti me Scale A Figure 541: Laguerre network performance functions. Table 59: Laguerre network optimization characteristics. Order A opr A nZCa 1 (A opr ) J (AnZCa) 0 13.5 11.0 2.6 2.6 1 38.5 41.0 5.3 5.2 2 63.0 69.0 6 7 6 5 3 42.0 39.0 8.0 7.9 4 60.5 60.0 9.0 9.0 5 74.0 55 0 9.2 9.1 6 92.5 39 5 9.5 9.3 7 111.0 49.0 9.6 9 .5 8 65 0 63.0 9.8 9.7
PAGE 127
121 (order=0) 0 8 I / ;..0 6 ~ :,' t1 C1' ...., 0 4 , 0.2 ! jfi . .. .. .. . . _ : .. .. .. .. ... .. . . .. . . .. .. ... ... . 0 .... ..,. ._ / __ .___ ___ ! \ .. 0 2 .___ ___ ,._ ___ ,.__.. 0 50 100 A (order=2) 0 8 .1 0 6 ~ +t....r C1' 2, 0 4 ,. .. 0 f .. : ., ; ~ ; ; : r = ~ 7 :. ~ ""< . ,+~+t ; .... .. .. ... 0 2 .___ ___ ,._ ___ ,.__.. 0 50 100 A (order=1) ~<1'0.6 ...., 0.4 ...., o.2 =. .. .. .. .. .. .. . ) 0 : i/ _..\ __ .,__._ __ _ _ \ ,,, l ,._ ... .. : 0.2 .___ ___ ...__ ___ _.__~ 0 50 100 A (order=3) 0 8 ++t ,(1' 0 6 \ ~ 0.4 r+\___ ./" 1 , 0 2 ............. ~ = ==~ ,,, 0 /\ .. i ... ~ : ~ : ~ : ~ : ~ . ~ : ; : ;,.... c a 1 ~ 0 .2 .___ ___ ...__ ___ _.__~ 0 50 100 A Figure 542: J and 2Im[l a ] for FMS pitch axis (orders= 0 to 3 ) (order=4) (order=5) 0.4 ~~~~ 0.4 ~~I~~~ 0.3 ~ :;<1' 0 2 .. .... . ....... . ... ++t r , 0 : [\!" ,f ~ cc j t i .l : .. .. .. .. ... .. .. .... .. .. .. .. 0.3 ...... ....,EC1' 0.2 ........ .... . .. .,_ \ 1f1 '.. ____ _____ ~:; 0 1 ~ ~ t1 n 0 f' t/ .. l .. .. .. 0 1 .___ ___ ,._ ___ ,.__.. 0 1 .___ ___ ...__ ___ _._____. 0 50 100 0 50 100 A A (order=6) (order=?) 0.2 .. .... .... \'\:=t+:1 0 15 _,.,.,,(1' ""1"~J'_/" 0.1 i , 0 05 Hf 1: ++t QQ~ j f ~ ~ .. . . .. .. .. . .. .. .. .. .. 0.2 \ ~<1' 0 15 ...... .. ...... ... ~ +t1 2, 0.1 i J 0 05 ..... ; __ __ __ __ _, 0~ : v \ ;P "'::., ., .. . . . . .. . . ,,, 0 50 100 0 50 100 A A Figure 543: J and 2Im[l a ] for FMS pitch axis (orders= 4 to 7)
PAGE 128
122 CT Laguerre Model Pitch Axis Step Response (order=8) i 1 5 r~ .. ; o 5 <( _J O'l 0 I / O'l 0 0 Q) Q) 0 1 0 2 0.3 0.4 0 5 0.6 0 7 0 8 CT Laguerre Model Pitch Axis Step Response Error (order=8) 0 9 0 05 ~......_~ __ ___._ __ _.__ __ .....__ __ ......__~ __ ___._ __ __. 0 0 1 0.2 0.3 0 4 0 5 0 6 0 7 0 8 0 9 t (sec) Figure 544: Laguerre order 8 step r esponse and associated error. CT Laguerre Model Pitch Axis Flyout Response (order=8) ~ 4 5 6 7 8 9 10 CT Laguerre Model Pitch Axis Flyout Response Error (order=8) 0 5 ,"T"".,r,~, .. ~ 0 .... .... ..... Q) Q) 0 5 .__ __ __._ ___ _,__ __ ____. ___ __,_ ___ ...__ __ ___,_ ___ _, 3 4 5 6 7 8 9 10 t (sec) Figure 545: Laguerre order 8 flyout response and associated error.
PAGE 129
123 5.3.5 Conclusion (Step Reponse Modeling) The purpose of this chapter has been to review and illustrate the applicability of GLFF networks in the area of system modeling. Focus was on the following concepts: 1. A methodology for determining GLFF network realizations of the impulse response of a closedloop system using measured step response data is possible. The procedure is simple and can be practically implemented in a laboratory setting. 2. A Carco fiveaxis flight motion simulator used in weapon system hardware intheloop test and evaluation, can be modeled very accurately using low order GLFF networks. This can be quite useful since a digital model of the FMS could be substituted for the hardware when missile flyout models are under development or when FMS response issues are raised and the hardware is unavailable. 5.4 Other Applications When conducting laboratory HITL closedloop simulations some of the weapon s subsystems such as the actuators gyros and accelerometers must be modeled If hardware actuators are available they must be pressurized and loaded to provide accurate fin position. This becomes impractical to do on a simulationbysimulation basis. If a single step response test was performed and the tools described in this research employed, then a GLFF network could be constructed to accurately represent and replace the actuators. Gyros are also difficult to work with in HITL testing and they are usually modeled. Accelerometers can not be exercised because no body translational motion is created in the laboratory. Again if hardware is available and can be statically evaluated to provide step response data GLFF networks could be constructed to model these subsystems. This offers an advantage over classical models because the GLFF networks provide a representation of the real systems as
PAGE 130
124 compared to the theoretical approach where the internal parameters (gains, time constants, damping coefficients, etc.) are considered to be known and assumed to be operating with exact precision. 5 5 Conclusion (Applications) This chapter has discussed various applications of the GLFF model class. A historical overview of the many kinds of problems that have been solved was first presented. Next several examples were given to illustrate how these models can be used to create networks that approximate systems when only step response data is available. Having the model parameters the impulse and step responses were recreated and compared to the true systems Since the "real world" applications of interest are systems and subsystems involved in test and evaluation of guided weapons the T&E process was described. This included the importance of laboratory tests such as alldigital and HITL simulations. A HITL test configuration was given and was used to identify various systems that are good candidates for GLFF networks. Primarily, the FMS system was studied and several GLFF networks were created using sampled step response data
PAGE 131
L APPENDIX A ORTHOGONAL SIGNAL SETS A set of complexvalued functions of the real variable t (A is considered fixed here) n (t J)f =o' n = 0 1 2, .. is called orthogonal [Su71] over the range (a b) if b { Q f m (t J : (t A )dt = nonzero mi= n m=n a The set is called orthonormal over the range (a b) if 125 mi= n m=n (A.1) (A.2)
PAGE 132
APPENDIXB ~BERT SPACE CONCEPTS D e finition : (Inner Product) Let Ebe a comple x vector space. A mapping ( ):E x E ~ c is called an inn e r pr o du c t in E if for any x,y,z E E and any a /3 E C the following conditions are satisfied: A. ( x, Y } = ( y,x)* B. (ax+ /J y, z ) = a( x, z ) + /J( y, z ) C ( x,x ) 0 and ( x,x ) = 0 implies x = 0 D e finition: (Inner Product Space) A vector space with an inner product is called an inn e r produ c t spa ce D e finition : (Norm in an Inner Product Space) The norm in an inner product space Eis defined as the functional !! x ii= (B 1 ) ( B 2 ) ( B.3 )
PAGE 133
BI 127 R e mark : The space L \ [a b]) of all Lebesgue square integrable functions on the inter v al [ a, b] with the inner product b (/ g) = f f ( t)g(t ) dt a is an inner product space D e finition: ( Cauchy Sequence ) A sequence { x 11 } n = 1 2 ... is a Cau c h y s e qu e n ce if for every E > 0 there exists a number M such that ll xm xn ll < E for all n m M. D e finition: ( Complete Inner Product Space ) An inner product space is c ompl e t e if every Cauchy sequence from the space converges to an element in the space. D e finition: (Hilbert Space) A Hilb e rt Spa ce is a complete inner product space. (B.4 )
PAGE 134
APPENDIXC SDOMAIN DERIVATION OF MINIMUM J Using Parseval's theorem we can represent the squared error performance function in the sdomain as 00 J = f e 2 (t)dt = 2 ~ f E(s)E(s)ds (C.l) C where E(s) = D(s)Y(s) and C is a contour taken in a CCW direction enclosing the entire left hand plane (LHP) in s. If Y(s) is the OGLFF model of the impulse response then N E(s) = D(s)L wn n(s,/4) = D(s)WT (s ,/4) n = O where (C.2) So the performance function becomes J = 2 ~ f[ D(s)wT (s ,/4 )][ D(s)wT (s ,/4 )Jis C = 2 ~ fD(s)D(s)dsw T 2 ~fD(s)(s,A)ds (C.3) C C wT 2 ~ fD(s)(s ,A)ds +w r 2 ~ f(s,A)T(s,A)ds w C C Since { n(s,A)r = Ois an orthonormal sequence 2 ~ f n(s,A ) m( s ,A )ds = 8nm (C.4) C
PAGE 135
129 or 2 ~ f (s, /4 )T (s, /4 )ds = I (C.5) C Also for the contour C we have 2 ~ f D(s)(s /4 )ds = 2 ~ f D(s)(s /4 )ds (C.6) C C 00 With the definition d e = J d 2(t)dt = 2 ~ f D(s)D(s)ds which 1s the energy of the signal d(t) then 0 C J = d e WT 2 ~ f D(s)(s /4 )ds + WT w C For a fixed A, to minimize J we must satisfy V w l = 0 Hence V w l = 2 ~ f D(s)(s A )ds + 2w = 0 C which yields w = 2 ~ f D(s)(s A)ds C (C.7) (C.8) (C.9) as the optimal weight vector that minimizes the squared error performance function. Substituting w into the equation for J gives (C.10)
PAGE 136
APPENDIXD SQUARED ERROR ST ATIO ARITY USING THE LAGUERRE MODEL The necessary conditions for stationarity of the quared error performance function when using a Laguerre model are now derived [Wan94b]. The deri ation makes use of the following lemma easily proven using induction. L emma : For the Laguerre kernels L ( ,A)= ,,fiJ, ( Ar, A> O ,n = 0,1 .. n ( +Ar + the following relation is true In Chapter 4.3.1 the minimum squared error performance function for a GLFF model (using the Laguerre kernel) was derived N J = d wr w = d w 2 e e L.,; n (D.2) n = O where d e = 2 ~ f D(s)D(s)ds C w = 2 ~ f D(s)L( s,A )ds C (D.3) W = (WO WI' ... W N ) T L(s,A) = ( La(s, A), Li (s,A ), ... L N (s,A) )r 130
PAGE 137
131 N J achieves a minimum when I, w~ is a maximum. So the stationary points are found n = O by (D.4) Now (D.5) Using the above Lemma and simplifying gives (D.6) N So L w~is a maximum when n = O (D.7) or (D.8) Hence the optimal value of A is either a root of w N (A)= 0 or a root of wN+i (A)= 0.
PAGE 138
REFERENCES [Aga92] Agamennoni 0. Paolini E. Dasages A. On Robust Stability Analysis of a Control System Using Laguerre Series ," Automati ca, vol. 28 no 4 pp. 815818 199 2 [AIJ94] Al Jabri A. Alshebeili S. Laguerre Transform for Speech Compression ," C o n ve r e n ce R eco rd IEEE In s trum en ta t i o n & M e asur e m e nt T ec hnolog y Conf e r e n ce, pp 737740 1994 [Arm57] Armstrong H. On Finding an Orthonormal Basis for Representating Transients ," IRE Tran s o n Cir c uit Th e o ry, vol. CT4 pp. 286287 Sep 1957 [Arm59] Armstrong H ., On the Representation of Transients by Series of Orthogonal Functions ," IRE Tran s. o n Cir c uit Th eory, vol. 6 no 4 pp 351354 Dec 1959 [Bak81] Baker G. GravesMorris P ., Pade Approximants. Part I: Basic Theory ," En cy cl o p e dia of Math e mati cs and its Appli c ations vol. 13 AddisonWesley 1981 [Bod94] Bodin P ., Wahlberg B. Thresholding in High Order Transfer Function Estimation ," Pr oc o f th e 3 3rd c onf e r e n ce o n D e cision & Control pp. 34003405 Dec 1994 [Bro65] Broome P. Discrete Orthonormal Sequences ," Journal of th e Assoc. for Computin g Ma c h i n ery, vol. 12 no. 2 pp. 151168 Apr 1965 [Bru64] Bruni C. Analysis of Approximation of Linear and Timeinvariant Systems Pulse Response by Means of Laguerre Finite Term Expansion ," IEEE Tran s on Automati c Control vol. 9 no 4 pp. 580581 Oct 1964 [Cel95] Celebi S. Principe J. Parametric Least Squares Approximation Using Gamma Bases ," IEEE Trans. on Signal Proc e ssing vol. 43, no. 3 pp. 781784 Mar 1995 [Chu84] Churchill R. Brown J., Compl ex Variabl e s And Applications McGrawHill 1984 [Cle82] Clement P. Laguerre Functions in Signal Analysis and Parameter Identification ," Journal of Franklin In s titut e, vol 313 no 2 pp 8595, 1982 132
PAGE 139
133 [Clo65] Clowes G., "Choice of the Timescaling Factor for Linear System Approximations Using Orthonormal Laguerre Functions ," IEEE Trans. on Automatic Control, vol. 10 no. 4, pp. 487489, Oct 1965 [Clu91] Cluett W ., Wang L., "Modelling and Robust Controller Design Using Step Response Data ," Chemical Engineering Science vol. 46 no. 8, pp. 20652077, 1991 [Clu92] Cluett W., Wang L. "Frequency Smoothing Using Laguerre Model ," IEE ProceedingsD: Control Theory & Applications vol. 139 no. 1 pp 8896, Jan 1992 [Coh92] Cohen H. Mathematics for Sci entists and Engineers, PrenticeHall 1992 [Con94] Constantinides A., Stathaki P., "Complex Interpolation for Rational Orthogonal Signal Approximation with Applications ," 1994 IEEE ICASSP vol. 4 pp 465468 [Dav95] Davila C., Chiang H ., "An Algorithm for PoleZero System Model Order Estimation," IEEE Trans. on Signal Processing vol. 43, no. 4 pp. 10131017 Apr 95 [Dea64] Dean W., Seismological Applications of Laguerre Expansions ," Bull etin of the Seismological Society of America vol. 54 no 1 pp. 365407, Feb 1964 [Deb90] Debnath L., Mikusinski P. Introduction to Hilb ert Spac es with Applications, Academic Press 1990 [DeC89] DeCarlo R., Linear Systems: A Stat e Variable Approach with Numerical Implementation PrenticeHall 1989 [den93a] den Brinker A., "Adaptive modified Laguerre Filters," Signal Processing, vol. 31, pp. 6979 1993 [den93b] den Brinker A., "Calculation of the Local CrossCorrelation Function on the Basis of the Laguerre Transform ," IEEE Trans on Signal Processing vol. 41, no. 5 pp. 19801982, May 1993 [den94] den Brinker A. "Laguerredomain Adaptive Filters," IEEE Trans. on Signal Processing, vol. 42, no 4 pp 953956, Apr 1994 [Dum86] Dumont G ., Zervos C. "Adaptive Control Based on Orthonormal Series Representation," IFAC Workshop Adaptive Systems in Control & Signal Proc essing, pp. 105113 1986 [Dum90] Dumont G. Zervos C., Pageau G., "Laguerrebased Adaptive Control of pH in an Industrial Bleach Plant Extraction Stage," Automatica, vol. 26, no. 4 pp. 781787 1990
PAGE 140
134 [Eic89] Eichblatt Jr. E. T est and Evaluation of th e Tactical Missile, Progress in Astronautics & Aeronautics vol. 119 1989 [Eko94] Ekongolo S ., Loverini C. Ragot J. Identification of LTI Systems using Laguerre Models ," Pro c of the 33rd Conf e r e nc e on D ec ision & Control, pp. 17281732 Dec 1994 [Eko95] Ekongolo S. Maquin D. Ragot J. Time Series Expansion Using Discrete Bilateral Laguerre Models ," Electronics L ette rs, vol. 31 no. 12, pp. 955956 Jun 1995 [Els94] Elshafei A. Dumont G. Elnaggar A. Adaptive GPC Based on Laguerre filters Modeling ," Automatica, vol. 30, no. 12 pp. 19131920, 1994 [Fra90] Franklin G ., Powell J. Workman M. Digital Control of D y nami c S ys t e ms 2 nd ed ., Addison Wesley 1990 [Fri81] Friedman D. On Approximating an FIR Filter Using Discrete Orthonormal Exponentials ," IEEE Tran s. on ASSP, vol. 29, no 4 pp. 923926, Aug 1981 [Fu93] Fu Y. Dumont G. An Optimum Time Scale for Discrete Laguerre Network, IEEE Trans. on Automatic Control vol. 38, no. 6 pp 934938 Jun 1993 [Gar80] Gamell P. Guid ed W ea pon Control S ys tems Brassey's Defence Publishers 1980 [Gun90] Gunnarsson S. Wahlberg B ., Some Asymptotic Results in Recursive Identification Using Laguerre Models ," Pro c of the 29th conference on D ec ision & Control pp. 10681073, Dec 1990 [Guo94] Guoqiang L. Dumont G. SampledData GPC (SDGPC) with Integral Action: The State Space Approach ," Proc. of th e 33rd conference on D ec ision & Control pp. 35673572 Dec 1994 [Gus93] Gustafsson T., Makila P., On System Identification and Model Validation via Linear Programming ," Proc. of the 32nd converence on Decision & Control, pp. 20872092 Dec 1993 [Hag91] Hagglund T., Astrom K., "Industrial Adaptive Controllers Based on Frequency Response Techniques ," Automatica vol. 27, no. 4 pp. 599609 1991 [Ham73] Hamming R., Numerical Methods for Scientists & Engineers 2 nd ed. McGrawHill, 1973 [Hay91] Haykin S. Adaptive Filt e r Th e ory 2 nd ed. Prentice Hall 1991
PAGE 141
135 [Hea55] Head J., "Approximation to Transients by Means of Laguerre Series," Proc. of the Cambridge Philosophical Society vol. 52 pp. 640651 1955 [Heu90] Heuberger P., Bosgra 0. "Approximate System Identification Using System Based Orthonormal Functions ," Proc. of the 29th conference on D ecis ion & Control pp. 10861092, Dec 1990 [Heu95] Heuberger P. Van den Hof P ., Bosgra 0., A Generalized Orthonormal Basis for Linear Dynamical Systems," IEEE Trans. on Automatic Control vol. 40 no. 3, pp 451465 Mar 1995 [Hug56] Huggins W. "Signal Theory ," IRE Trans. on Circuit Theory vol. 3, pp. 210215, 1956 [Hwa82] Hwang C., Shih Y., Parameter Identification via Laguerre Polynomials ," International Journal of S ys tem Science, vol. 13 no 2, pp. 209217, 1982 [Ira92] Iranzo V., Falques A. Some Spectral Approximations for Differential Equations in Unbounded Domains ," Comput e r Methods in Applied M ec hanics and Engin ee ring 98 pp. 105126 1992 [Ji95] Ji Z ., Kozick J. Aburdene M. "Implementation of Discrete Legendre Adaptive IlR Filter Using the DSP 56001 ," Midwest S ym posium on Circuits & S ys t e ms vol. 2, pp. 10751078, 1995 [Kab94] Kaber S. Filtering Nonperiodic Functions ," Computer Methods in Appli ed Mechanics & Engineering vol. 116, pp. 123130 1994 [Kau54] Kautz W. "Transient Synthesis in the Time Domain ," IRE Trans. Cir cuit Theory vol. CT1 pp. 2939, 1954 [Kin69] King J. O Canainn T. Optimum Pole Positions for LaguerreFunction Models ," El ect roni cs L ette r s, vol. 5 no. 23, pp. 601602 Nov 1969 [Kin77] King R. Paraskevopoulos P. Digital Laguerre Filters ," Int ern ational Journal of Circuit Th eory & Applications, vol. 5 pp. 8191 1977 [Kin79] King R. Paraskevopoulos P. "Parameter Identification of DiscreteTime SISO Systems ," International Journal of Control vol. 30, no. 6 pp 10231029 1979 [Ko94] Ko C. Li C., An Adaptive IlR Structure for the Separation, Enhancement and Tracking of Multiple Sinusoids ," IEEE Trans. on Signal Pro cessing, vol. 42, no 10 pp. 23322335, Oct 1994
PAGE 142
136 [Kor91] Korenberg M., Paarmann L., "Orthogonal Approaches to TimeSeries Analysis and System Identification ," IEEE Signal Processing Maga z ine, July 1991 [Kuo87] Kuo B., Automatic Control Syst ems, 2 nd ed., Prentice Hall, 1987 [Kuo94] Kuo J. Celebi S ., "Adaptation of Memory Depth on the Gamma Filter," 1994 IEEE ICASSP vol 4., pp. 373376 [Lam94] Lam J. "A nalysis of the Laguerre Formula for Approximating Delay Systems, IEEE Trans. on Automatic Control, vol. 39, no. 7 pp. 15171521 Jul 1994 [Lee33] Lee Y., Synthesis of Electric Networks by Means of the Fourier Transforms of Laguerre 's Functions ," Journal of Mathematics & Physics vol. XI, pp. 83113, 1933 [Lee67] Lee Y., Statistical Th eory of Communications, Wiley, 1967 [Lin93] Lindskog P. Wahlberg B. "A pplications of Kautz Models in System Identification," Proc eedings of the IFAC World Congress, vol. 5, pp. 309312, Jul 1993 [Mai85] Maione B. Turchiano B. "Laguerre Ztransfer Function Representation of Linear DiscreteTime Systems ," International Journal of Control, vol. 41, no. 1, pp. 245257, 1985 [Man63] Manley H. Analysissynthesis of Connected Speech in Terms of Orthogonalized Exponentially Damped Sinusoids ," The Journal of the Acoustical Soci ety of America, vol. 35, no. 4 pp. 464474 Apr 1963 [Mar69] Marzolla A. "On the Mean Square Approximation of a Function With a Linear Combination of Exponentials ," International Journal of Control vol. 9, no. 1, pp. 1726 1969 [Mar93] Marmarelis V. "Identification of Nonlinear Biological Systems Using Laguerre Expansions of Kernels," Annals of Biomedical Engineering, vol. 21, pp. 573589, 1993 [Mas90] MasnadiShirazi M. Ahmed N., "Laguerre Approximation of Nonrecursive DiscreteTime Systems," 1990 IEEE ICASSP pp. 13091312 Apr 1990 [Mas91] MasnadiShirazi M., Ahmed N ., "Optimum Laguerre Networks for a Class of DiscreteTime Systems," IEEE Trans. on Signal Processing, vol. 39, no. 9, pp. 21042108, Sep 1991
PAGE 143
, e 137 [McD68] McDonough R ., Huggins W. Best LeastSquares Representation of Signals by Exponentials ," IEEE Tran s on Automati c Control vol. 13 no. 4 pp 408412 Aug 1968 [Moh88] Mohan B ., Datta K ., Lumped and Distributed Parameter System Identification via Shifted Legendre Polynomials ," Journal o f D y nami c S y st e m s, M e asur e m e nt & Control vol. 110 pp. 436440 Dec 1988 [Moh95] Mohan B. Datta K. Analysis of Linear T i meInvariant TimeDela y Systems via Orthogonal Functions ," Int e rnation a l Journal of S y st e ms S c i e n ce, vol. 2 6 no. 1 pp. 91111 1995 [Nin94] Ninness B ., Gustafsson G ., A Unifying Construction of Orthonormal Bases for System Identification ," Pr oc. o f th e 33 rd c onf e r e n ce on D ec isi o n & Contr o l pp 33883393 Dec 1994 [Nur81] Nurges U ., Jaaksoo U. Laguerre State Equations for a Multivariable Discrete System ," Aut o mat i on R e mot e Control vol. 42 pp. 16011603 1981 [Nur87] Nurges U ., Laguerre Models in Problems of Approximation and Identification of Discrete Systems ," Automation & R e mot e Control vol. 48 pp 346352 1987 [Oli87] Olivier P. ReducedOrder Models Using Optimal Laguerre Approximations ," El ec tro n i cs L e tt e r s, vol. 23 no 6 pp. 257258 Mar 1987 [Oli94] Olivier P ., Online System Identification using Laguerre Series ," IEE Proc PartD: Control Th e o ry Appli c ations vol. 141 no. 4 pp. 249254 Jul 1994 [Opp89] Oppenheim A. Schaeffer R. Dis c ret e Tim e Signal Proc e s si n g, Prentice Hall 1989 [Pab92] PabonOrtiz C. Artoni M. Laguerre polynomials: novel properties and numerical generation scheme for analysis of a discrete probability distribution ," Comput e r Ph y sics Communications vol. 71 no. 3 pp. 215221 Sep 1992 [Par71] Parks T. Choice of Time Scale in Laguerre Approximations using Signal Measurements ," IEEE Trans. on Automatic Control vol. 16 no. 5 pp. 511513 Oct 1971 [Par91] Partington J. Approximation of Delay Systems by FourierLaguerre Series ," Automati c a vol. 27 no. 3 pp 569572 1991 [Per88] Perez H. Tsujii S. IIR Adaptive Filtering via Discrete Legendre Functions ," El ec tr o ni c s L e tt e rs vol. 24 no 8 pp 450451 Apr 1988
PAGE 144
138 [Per91] Perez H ., Tsujii S ., A System Identification Algorithm Using Orthogonal Functions ," IEEE Trans on Signal Proc ess ing vol. 39, no. 3 pp. 752755, Mar 1991 [Pow79] Powers D. Bounda ry Value Probl ems, 2 nd ed. Academic Press 1979 [Pri93] Principe J ., de Vries B ., de Oliveira P ., The Gamma Filter A New Class of Adaptive IlR Filters with Restricted Feedback, IEEE Trans. on Signal Proc es sing vol. 41 no. 2, Feb 1993 [Rob87] Roberts R. Mullis C. Di gita l Si gna l Proc ess ing, AddisonWesley, 1987 [Ros64] Ross D. Orthonormal Exponentials ," IEEE Trans. on Communication and Electronics, vol. 1 2, no 71 pp. 173176 Mar 1964 [Sch70] Schetzen M. Powerseries Equivalence of Some Functional Series with Applications ," IEEE Trans. on Circuit Th eory, vol. 17 no. 3, pp. 305313, Aug 1970 [Sch71] Schetzen M ., Asymptotic Optimum Laguerre Series ," IEEE Trans. on Cir cuit Th eory, vol. 18 no. 5 pp 493500 Sep 1971 [Sil93] Silva T ., Generalized Feedforward Filters: Some Theoretical Results ," 1993 IEEE ICASSP [Sil94] Silva T. Kautz Filters ," technical r e port Dec 29 1994 [Sil95] Silva T ., Rational Orthonormal Functions on the Unit Circle and on the Imaginary Axis with Applications in System Identification ," t ec hnical report Oct 18, 1995 [Ste65] Steiglitz K. "The Equivalence of Digital and Analog Signal Processing ," Information and Control vol. 8 pp. 455467 1965 [Str88] Strum R ., Kirk D ., Di scre t e Tim e Si g nal Processing Prentice Hall, 1988 [Su7 l] Su K. Tim e Domain S ynthes i s of Lin ea r Networks, Prentice Hall 1971 [Tia89] Tian W. Cumulant Based Probabilistic Power System Simulation Using Laguerre Polynomials," IEEE Trans. on Energ y Conversion vol. 4, no. 4 pp. 567573, Dec 1989 [Tse94] Tseng C. Pei S., Two Dimensional Adaptive Filters Using Laguerre Function ," 1994 IEEE ICASSP vol. 3 pp. 417420 [Wah91] Wahlberg B. System Identification Using Laguerre Models, IEEE Tran s. on Automatic Control, vol. 36 no. 5 pp. 551562 May 1991
PAGE 145
139 [Wah93] Wahlberg B. Hannan E. Parametric Signal Modelling Using L a guerre Filters ," Th e Annals of Appli e d Probabili ty, vol. 3 no. 2 pp 467 496 1993 [Wah94] Wahlberg B. System Identification Using Kautz Models ," IEEE T ra n s on Aut o m a ti c C o ntr o l vol. 39 no. 6 Jun 1994 [Wan94a] Wang L. Cluett W. System Identification based on ClosedLoop Step Response Data ," IEE Proc ee dingsD: Control Th e ory & Appli c ations vol. 141 no. 2 pp 107110, Mar 1994 [Wan94b] Wang L. Cluett W. Optimal Choice of TimeScaling Factor for Linear System Approximations using Laguerre Models ," IEEE Trans on Automati c Control vol. 39 no. 7 pp. 14631467 Jul 1994 [Wan95] Wang L. Cluett W. Building Transfer Function Models from Noisy Step Response Data Using the Laguerre Network ," Chemical Engine e ring Sci e nc e, vol. 50 no. 1, pp 149161 1995 [Wid85] Widrow B ., Steams S. Adapti ve S ig nal Pr oces s i n g, Prentice Hall 1985 [Yen62] Yengst W ., Approximation to a Specified Time Response ," IRE Tran s o n Circuit Th e ory vol. CT9 pp. 152162 Jun 1962 [You62a] Young T. Huggins W. Discrete Orthonormal Exponentials ," Proc. National El ec tr o ni cs C o nf e r e n ce, vol. 18 pp. 1018 Oct 1962 [You62b] Young T. Huggins W ., Complementary Signals and Orthogonalized Exponentials ," IRE Tran s on Cir c uit Th e o ry, vol. 9 pp. 362370 1962 [You63] Young T. Huggins W ., On the Representation of Electrocardiograms ," IEEE Tran s BioM e di c al El ec troni c s vol. BME10 pp. 8695 Jul 1963 [Zer88a] Zervos C ., Dumont G ., Deterministic Adapti v e Control Based on Laguerre Series Representation ," Int e rn a tional Journal of Control vol. 48 no 6 pp. 23332359 1988 [Zer88b] Zervos C. Belanger P. Dumont G ., On PID Controller Tuning Using Orthonormal Series Identification A u tomati c a vol. 24 no 2 pp. 165175 1988
PAGE 146
BIOGRAPHICAL SKETCH Douglas G. Jones was born in Valdosta Georgia on May 18 1963. He received his bachelor of science degree in applied mathematics with highest honor from the Georgia Institute of Technology in 1986. After graduation he worked for RCA Services Corporation in the Seeker Evaluation and Test Simulation Laboratory at Eglin AFB FL He became an employee of the Department of Defense in 1989 and now works as an electronics engineer in the Guided Weapons Evaluation Facility a t Eglin AFB. He received the master of engineering degree from the University of Florida in 1992 and a Ph D in 1999. He plans to continue conducting research in the area of digital signal processing with emphasis on system modeling and simulation design. 140
PAGE 147
I certify that I have read this study and that in my opm10n it conforms to acceptable standards of scholarly presentation and is fully a = d=,.= quality as a dissertation for the degree of Doctor of Phi ophy. I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality as a dissertation for the degree of Doctor of ~~ William W. Edmonson Assistant Professor of Electrical and Computer Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate in scope and quality as a dissertation for the degree of Doctor of Philosophy. (l_iL Q,) ~ ~ arris Assistant Professor of Electrical and Computer Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate in scope and quality as a dissertation for the degree of Doctor of Philosophy. Jian Li Associate Pr and Computer Engineering I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate in scope and quality as a dissertation for the degree of Doctor of Philosophy. nmei Chen rofessor of Mathematics
PAGE 148
This dissertation was submitted to the Graduate Faculty of the College of Engineering and to the Graduate School and was accepted as partial fulfillment of the requirements for the degree of Doctor of Philosophy. August 1999 ~<=42/4 M. J. Ohanian Dean College of Engineering Winfred M Phillips Dean Graduate School
PAGE 149
111 1111111111111~111~1~1 11 llll[lil l lllllillll' il~I I III II I I 3 1262 08554 4228
PAGE 150
.. 1 ,...... 1 1 ... __ ..._, __ r"' 11 ..... ..c UNIVERSITY OF FLORIDA 111 1111111111111111111111111111111111111111111111111111111111111 3 1262 08554 4228

