Citation
System modeling using generalized linear feedforward networks

Material Information

Title:
System modeling using generalized linear feedforward networks
Creator:
Jones, Douglas G., 1963-
Publication Date:
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 )
non-fiction ( marcgt )

Notes

Thesis:
Thesis (Ph.D.)--University of Florida, 1999.
Bibliography:
Includes bibliographical references (leaves 132-139).
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 non-profit 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 long-term 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, Time-Invariant, Causal Systems............................................................ 10
2.2 R ealizability ................................................................................................. 14
2.3 L2 Signals and Systems ...................................................................................... 14
2.4 System Synthesis.......................................................................................... 15
2.5 Continuous/Discrete-Time 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 iscrete-Tim 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 Discrete-Time GLFF Kernels.................................................... 27
3.3.2 Popular Continuous-Time GLFF Kernels............................................... 27
3.4 GLFF Kernel Derivations and Properties..................................... ............ 30
3.4.1 Continuous-Time Laguerre Kernel.............................................................. 30
3.4.2 Continuous-Time Kautz Kernel (Complex Time Scale Vector) ...............30
3.4.3 Continuous-Time Kautz Kernel (Real Time Scale Vector).........................31
3.4.4 Continuous-Time 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 Wiener-Hopf 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 Closed-Loop Missile Simulations....................................................94
5.3 Modeling the Step Response of a Carco 5-Axis FMS ....................................95
5.3.1 Discrete-Time 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 Continuous-Time 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: S-DOMAIN 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 goodness-of-fit. 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 discrete-time, continuous-







time, and sampled-data 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, time-invariant 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 time-invariant (LTI) system, H,

using a pre-defined 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 1-1 below. Arguments

are not included with the functions {x, d, y, e, H} to illustrate that the formulation will

apply in both a continuous-time setting {x(t), d(t), y(t), e(t), H(s)} and a sampled-data

or discrete-time setting {x(k), d(k), y(k), e(k), H(z)}. Quite often in this thesis we

switch back-and-forth between continuous-time and discrete-time concepts. This is







intentionally done to prove the broad applicability of this work. Most of the theory

however will be derived in continuous-time. It would have been just as appropriate to

develop the concepts in discrete-time.



d
H

x e

M(p)


Figure 1-1: 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

goodness-of-fit 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 (continuous-time)
or (1.1)
J = e2(k) = J(d(k)- y(k))2 (discrete-time)
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 discrete-time setting. In the past discrete-time

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 tap-delay 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 1-2 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 1-2: 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 near-optimal 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 discrete-time [You62a] as did Broome [Bro65]. Nurges and

Yaaksoo [Nur81] formulated state space equations that represent discrete-time linear

systems using the Laguerre model. Roberts and Mullis provide a generalized state

variable formulation of orthogonal all-pass 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)
i--1

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 well-known approaches to performing this approximation. The first

is to transform h(t) into H(s) and then approximate H(s) in the s-domain 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, continued-fraction 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 5-axis flight motion simulator (FMS). This system is used in hardware-in-the-

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, time-invariant, 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 continuous-time 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, Time-Invariant, Causal Systems

There are many ways to classify systems; linear vs. nonlinear, time-invariant vs.

time varying, causal vs. non-causal, 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 well-studied 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 time-invariant (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


(71s-z)
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 square-integrable 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/Discrete-Time Transformation

The goal of this thesis is to define a model class useful for approximating discrete-

time, continuous-time, and sampled-data systems. The above system formulation was

carried out in the continuous-time domain. The same principles however could have

been derived in the discrete-time domain. All the concepts discussed have equivalent

discrete-time representations. Sometimes a continuous-time system in analytic form

is given when a discrete-time representation is desired. This can be resolved by

sampling. In this case, the discrete-time system is called a sampled-data system.

Hence a sampling time, T, will necessarily be specified.

Having a continuous-time transfer function, there are several methods used to map

to the discrete-time domain. The most popular are the bilinear, impulse- or step-

invariant, and matched-Z 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 2z-1 (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

nT-T nT
d,(n)= d,(nT)= j[-d,(v)+x(u)]dv+ j[-d,(v)+x(v)]dv
0 Tr-r (2.17)
nT
= d(n-1)+r I[-d,(v)+x(u)]dv
nT-T

The method in which the integral is computed over the region [nT-T, nT]

determines the discrete-time 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 2-1: Continuous-time to discrete-time mappings.

Integration Rule Mapping H(s) I-4 H(z)
Euler's 1 fz-l
s=-\ --
T z
Trapezoidal 2(z-1)
T _z+l)


So the bilinear transform is derived by computing the integral over [nT-T, 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 continuous-time vector space L2(R) and the discrete-time

(or digital) vector space of double-ended 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 continuous-time

domain and discrete-time domain. In this thesis, when a continuous-to-discrete

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 2-1 to 2-4. 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 2-1: 1st order impulse responses, r-0.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 2-2: lst order step responses, r=0.5, 1, 2, 3.




































0 1 2 3 4 5 6 7 8 9 10
t (sec)

Figure 2-3: 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 2-4: 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, time-invariance, 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 continuous-time system descriptions, H(s), into

discrete-time 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 discrete-time

case) as well.


3.1 General Model Form M(p)

Recall the general system modeling configuration, shown in Figure 3-1 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 3-1: 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

continuous-time case and H,(z,2) in the discrete-time case.


3.2 Discrete-Time Models

Consider the definitions above in the discrete-time 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
Z-I


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)= { ,
z-1







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 3-2, 3-3, and 3-4. 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 3-2: FIR model block diagram.


y(k)


Figure 3-3: IIR model block diagram.


y(k)


Figure 3-4: 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 discrete-time and continuous-time kernels are given in Tables 3-1 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 degrees-of-freedom 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 Discrete-Time GLFF Kernels

The most popular discrete-time GLFF kernels are the delay, Gamma, Laguerre,

Legendre, and Kautz. They are listed in Table 3-1, and restrictions on the time scale

parameters are given in Table 3-2. 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 Continuous-Time GLFF Kernels

The most popular continuous-time GLFF kernels are the Gamma, Laguerre,

Legendre, Jacobi, and Kautz. They are listed in Table 3-3, and restrictions on the time

scale parameters are given in Table 3-4. Some of first to describe these kernels are

Lee [Lee33], Kautz [Kau54], Armstrong Ar[Arm59], and Principe [Pri93].
















Table 3-1: Discrete-time GLFF kernels.


Table 3-2: Restrictions on discrete-time 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-' A-1 no 1
___1-(1-2)z 1-(1-A)z
Laguerre Ln(z,A) _2 z 22 -_ I yes 1
1-A2z- 1-Az-1 1-2 1- 2z-) __
Legendre P,(z,A) 1 z-1 1 z A yes >1
1-2jz- 1-2z- 1 1 1-'z-
Kautz K,(z,A) 1 -/-1 1 n -1 yes >1

1 I Z-1 1- Z- i= 1-A iz-
Z_ i=O














Table 3-3: Continuous-time GLFF kernels.


Table 3-4: Restrictions on continuous-time 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 n-2 fS-A 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) (s-AkX(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 continuous-time 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 Continuous-Time 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)
S-A S+A


3.4.2 Continuous-Time 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,,e-1}, 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

2-A(s+In) (s (s-A,)(s-A) n-1
K2,,(t,A)=I {k, (t, )}= ( + m 0,1,..., 2

Sr2AR (s-IA) (s1(-A,)(s-AX)' n-1
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 Continuous-Time 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( s-A =, 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 Continuous-Time 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

discrete-time representation of the system. Optimal continuous-time 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 discrete-time impulse response is desired. How can it be computed?

Three simple methods are now given.

1. Using the discrete-time 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 discrete-time and continuous-time 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 s-domain

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

continuous-time system. Given sampled step response data, continuous-time and

discrete-time 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+1-dimensional 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 Wiener-Hopf Equations

It is easy to show that w is in fact the solution to the Wiener-Hopf (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 Wiener-Hopf or Normal equations. Solving for w gives the

optimal solution

w = R-p (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 s-domain.


4.2.4 Computing w Using Step Response Data

In Chapter 5 applications of the GLFF networks are given. The systems there are

continuous-time 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 (discrete-time) 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 discrete-time 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 Wiener-Hopf

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 continuous-time step response g(t). Since the continuous-time

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) dd-gll( dt
w. =- g(t) dt dt g dt
ST (4.42)
= g(t) (t, ) dt + g,J (T,I A)
0-ldt

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

n-1 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 continuous-time 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

n-1 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 (Wiener-Hopf)

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) non-convex 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 tap-to-tap weights like the GLFF

structures.

Unfortunately, the elements of the time scale vector are not linearly related; hence,

the performance function J is usually non-convex 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 non-orthogonal Gamma and a few discuss optimization of the

more complex Kautz kernel. Although the exact optimal solution is sought, Kautz

has pointed out that non-optimal 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

non-orthogonal 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

right-half plane poles or multiple j-axis 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_,XN-I + 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,...,n-1. Replace each h()(t)

by A' to get the characteristic equation

A + an(n-l)/ + 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)= A-e- 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 discrete-time 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 discrete-time 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= -1-M2 (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 non-orthogonal 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 M-bank 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 discrete-time 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 continuous-time and discrete-time when approximating a system

excited by an arbitrary input.








Iterative and adaptive search techniques. When representing a system with a

discrete-time GLFF network, the elements of the time scale vector are restricted to

bounded regions to guarantee stability. In particular, the discrete-time 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 left-hand 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+2-j s+2+j



Example 2: (4th order system, rapidly decaying impulse response with slight
oscillation)


1 1 1 1
H(s) = + + +-- (4.87)
s+2-j s+2+j s+l-2j s+1+2j



Example 3: (4th order system, slowly decaying impulse response with oscillation)


2 2 -3 -3
H(s) = --+ -- + -+ -(4.88)
s+0.4-j s+0.4+j s+0.6-2.8j s+0.6+2.8j



Experimental results for each example are tabulated and graphed on the next few

pages. In Figures 4-1 and 4-2 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 4-3

and 4-4 for example 2, and Figures 4-5 and 4-6 for example 3. Notice the imaginary

component has zero crossings near the stationary points of the squared error







performance function. Im[Jj] crosses negative-to-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[J,] can be considered an estimate of

V,J(A) that is best near Aop.t

Numerical results are provided in Tables 4-1, 4-2, and 4-3 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 4-1: 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 4-2: J and 2Im[Ju] for example 1.










0


-5


0 -10
o
S-15
u-
a)
2 -20
CO

-25


-30


S-35

-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 4-4 J and 2Im[J for example 2.
Figure 4-4: J and 2Im[Ju] for example 2.


(order=1)


2

(order=3)


2 4 6


1 2 3 4
Time Scale, X

Figure 4-3: 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 4-5: 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 4-6: J and 2Im[J,] for example 3.


CO

0
c -10
U-










0
-15



" -20



z-25



-30


1




j'_0.5
E


0


0.6


0.4


E 0.2
-0.
0


-0.2













Table 4-1: 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 4-2: 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 4-3: 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

4-7, 4-8, and 4-9 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]

(dot-dash curve), are plotted for the first four model orders. Numerical results are

provided in Tables 4-5, 4-6, and 4-7 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 4-4 we find that duel is 22% of


Table 4-4: 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 4-9 that J,(A)) = J,(A), for all A, for example 3. Figures 4-7 and 4-8 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 4-7: J, 2Im[Jy], and 2Im[Ja] for example 1.


Table 4-5: 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
S-0.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 4-8: J, 2Im[J,], and 2Im[Ja] for example 2.




Table 4-6: 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 4-9: 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 4-7: 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 discrete-time 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 4-10 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 4-11 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 4-12 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 4-10: Impulse response and approximations using a 1st order Laguerre model.


Example 2: d(t),yo(t),yz(t)


3 4 5


Figure 4-11: Impulse response and approximations using a 5th order Laguerre model.


Example 3: d(t),yo(t),yz(t)


5 6 7 8


Figure 4-12: 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 non-convex 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 Wiener-Hopf 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 sampled-data

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 tap-delay line,

Gamma, Laguerre, Legendre, Jacobi, and Kautz models defined in both discrete-time

and continuous-time 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 non-causal 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 (non-simulated) 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 5-1 through 5-6. In each case, the optimal time scale

parameter was determined and used to compute and plot the step response and

impulse response for 6-tap (N = 5) networks. Figures 5-1, 5-2, and 5-3 compare the

system and GLFF network step responses. Figures 5-4, 5-5, and 5-6 show the

resulting impulse responses. The optimal time scale parameters of the 6-tap models

and associated minimum squared error, J(Aopt), are provided in Tables 5-1, 5-2, 5-3.










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 5-1: 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 5-2: 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 5-3: 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 5-4: 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 i-1
0 1 2 3 4 5 6 7 8 9 10
t (sec)


Figure 5-5: 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 5-6: 4th order impulse response and model data.









Table 5-1: 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 5-2: 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 5-3: 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) open-loop 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 open-loop tests.

Another class of tests studies the behavior of these systems in a closed-loop

configuration and involves a simulation of the real-world environment and the

weapon's flight in this environment. This class of tests can be broken down into two

types: all-digital and hardware-in-the-loop (HITL). In an all-digital closed-loop

flight, the entire weapon system is simulated mathematically. If the weapon

subsystem models are accurate, all-digital 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







closed-loop 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, open-air 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
(all-digital 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 open-air 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 5-7 below illustrates a high-level functional break down of

a typical missile system [Gar80].


Figure 5-7: 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 pre-programmed 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 5-8: Closed-loop missile simulation.


5.2.3 HITL Closed-Loop Missile Simulations


In this section HITL closed-loop missile simulations will be discussed. Figure 5-8

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 closed-loop simulation, the seeker, guidance




Full Text
xml version 1.0 encoding UTF-8
REPORT xmlns http:www.fcla.edudlsmddaitss xmlns:xsi http:www.w3.org2001XMLSchema-instance xsi:schemaLocation http:www.fcla.edudlsmddaitssdaitssReport.xsd
INGEST IEID E6RLTDKRC_DJKTO5 INGEST_TIME 2013-02-14T13: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 long-term 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 Wiener-Hopf 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 Closed-Loop Missile Simulations .. ............ . .. ........ .. .. .... ....... ...... .. 94 5 .3 Modeling the Step Response of a Carco 5-Axis FMS ............................ .. ...... .. 95 5 3 l Discrete-Time 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 Continuous-Time 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 : I-ill.,BERT SPACE CONCEPTS .................... .... .......... ............. ........... . .............. 126 C: S-DOMAIN 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 f-fi 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 crete-ti 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 time-invariant ( LTI) system, H using a pre-defined 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 1-1 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 continuous-time setting {x(t) d(t) (t) e(t), H( ) } and a amp l ed-data or di er te-tim setting {x(k), d(k), y(k), e(k) H( ) } Quite often in thi the i witch back-and-forth between continuou -time and discrete-time 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 crete-time. d H X --~ M(p) V Figure 1-1: 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 -of-fit 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 discrete-time etting. In the past discrete-time 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 tap-delay 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 1-2 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 1-2: 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 ar-optimal 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 crete-ti me [Y u 2a] a did Br m [Br ]. Yaak ur81] formulat d tat pa e equati n that r pr ent di rete-tim 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 well-known 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 continued-fraction 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 5-axi flight moti n imulator loop (HITL) testing of guided weapon S Thi u d in hardware-in-thet 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 time-invariant, 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 continuous-time 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, Time-Invariant, Causal Systems There are many ways to classify systems linear vs. nonlinear, time-invariant v time varying causal vs. non-causal, 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 time-invariant (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 square-integrable 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) i-s+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/Discrete-Time Transformation (2 .11 ) (2.12) The goal of this thesis is to define a model class useful for approximating discrete time, continuous-time, and sampled-data systems. The above system formulation was carried out in the continuous-time domain. The same principles however could h ave been derived in the discrete-time domain. All the concepts discussed have equivalent discrete-time representations. Sometimes a continuous-time system in analytic form is given when a discrete-time representation is desired This can be resolved by sampling. In this case, the discrete-time system is called a ampled-data system Henc a sampling time T will necessarily be specified. Having a continuous-time transfer function, there are several methods used to map to the discrete-time 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,(n-1)+7 f[ d,(u)+x(v)]du nT T The method in which the integral computed o er th region [nT-T nT] determines the di crete-time 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 2-1 : Continuous-time to discrete-time 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 [nT-T 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 continuous-time vector space L 2(~ ) and the discrete-time ( or digital) vector space of double-ended 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 continuous-time domain and discrete-time domain. In this thesis when a continuous-to-discrete 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 2-4. 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 2-1: 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 2-2: 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 2-3: 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 2-4: 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, time-invariance 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 crete-time 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 discrete-time case) as well. 3.1 General Model Form M(p) Recall the general system modeling configuration, shown in Figure 3-1 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 3-1: 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 te-tim 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 3-2, 3-3 and 3-4 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 (k-1 ) x( k N) y( k ) Figure 3-2 : FIR model block diagram x (k) y( k ) y( k-2) y (k-1 ) Figure 3-3 : IIR model block diagram. y( k ) Figure 3-4 : 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 discrete-time and continuous-time kernels are gi en in Table 3-1 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 degrees-of-freedom 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! 3-1 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, Laguen-e, 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 3-4. 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( l-A ) z1 L ag u e rr e L n(z, A ) I -A ( 2 1 -A J yes 1 z I l -A z1 1-A 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 -J1-I 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 3-3: Continuous-time 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) s-A ye 1 --s+A s+A +A +A Legendre P n(s,A) fiI: s-A fii: n-1 -A y l k IT I s+A 1 s+A k s + A n 1 = 0 +A 1 Jacobi J nCs A) fiI: s-A 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) (s-Ak )(s-it") ,pX;: hi a -),Xs-'\) l (s+A 1 s+"S) (s+ A" )(s+ it ") (s+~I +J,,) I +-1i +,-t,) Table 3-4: 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 continuous-time 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 Continuous-Time 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 Continuous-Time 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 m-l (s-A )(s-A ~ ) n-1 2 "' t ;., "-2m (s + A m)(s + A ~ ,) ;=0 ( + A )(s +A ~ )' 111 = O l -2K t A = t A}= (s-A,)(s-A ~ ) m=Ol n-1 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 Continuous-Time 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 discrete-time representation of the system. Optimal continuous-time 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 discrete-time impul e response is desired How can it be computed? Three simple methods are now given. ,._ l. Using the discrete-time 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 ) = -" s-T :+ 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(k-1) = 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 discrete-time and continuous-time 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 continuous-time system. Given sampled step response data continuous-time and discrete-time 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+l-dimen 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 S-1.. S-1. 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 Wiener-Hopf Equations It is easy to show that w is in fact the solution to the Wiener-Hopf (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 Wiener-Hopf 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 s-domain. 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 continuous-time 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 discrete-time 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 Wiener-Hopf 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 continuous-time step response g(t). Since the continuous-time 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 continuous-time 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 (Wiener-Hopf) 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) non-convex 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 tap-to-tap weights like the GLFF structures Unfortunately the elements of the time scale vector are not linearly related; hence the performance function J is usually non-con 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 non-orthogonal 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 non-optimal 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 non-orthogonal 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 right-half plane poles or multiple j-axis 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 discrete-time 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 discrete-time 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= 2Ml-l-M 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 non-orthogonal 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 M-bank 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 discrete-time 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 high-order approximations of the weights w N (A) and WN+J(A) via Taylor series and Pade approximants. This procedure is illustrated in both continuous-time and discrete-time when approximating a system excited by an arbitrary input.

PAGE 62

56 Iterative and adaptive search techniques. When representing a system with a discrete-time GLFF network the elements of the time scale vector are restricted to bounded regions to guarantee stability. In particular the discrete-time 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 left-hand 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+l-2} 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 4-1 and 4-2 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 4-3 and 4-4 for example 2, and Figures 4-5 and 4-6 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 negative-to 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 4-1 4-2 and 4-3 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 4-1: 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 4-2: J and 2Im[lu] for example 1.

PAGE 69

-5 (0 :s-10 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 4-3 : 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 \t-a 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 4-5: 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 ~ -"c---t .. 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 4-6: J a nd 2 Im[Ju] for example 3

PAGE 71

65 Table 4-1: 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 4-2: 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 4-3: 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 4-7 4-8 and 4-9 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] ( dot-dash curve ), are plotted for the first four model orders. Numerical results are provided in Tables 4-5 4-6 and 4-7 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 4-4 we find that d ue l is 22 % of Table 4-4: 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 4-9 that l ) A ) ::::: 1 0 ( /4), for all /4 for example 3. Figures 4-7 and 4-8 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:.._; ~-j-j -0 05 .................. \ ,!.: '. < .... --[; ,_/ --t----t---1 : 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--+--r---r---r--==I i _; 0 + ...... ......... ...... '.:: .:.. \ ;,.....i \ e-. , ..-=-!""" \' -. I E ! \ \ \ .. \ -1 .......__....___..___.____,~__._. _.__. -0 5 l l \ -1 ....____.~-~____,~____, ~____, 0 2 3 ).. 4 5 0 2 Figure 4-7 : J 2Im[lu], and 2Im[J 0 ] for example 1. 3 4 5 Table 4-5: 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---+---+---4---i E . ': ~~---+ i __ ---1 ~; 0 .. \ ~ .. :~:; : t ~~ ~~i ~::: : .. .. -..;; .., ~ : ; : t : ~ :.: : : : : : : : : :: : :: -0 1 P., ___, ---+----i----i / 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~--,r--t-----1 7 \ii \ / t ; -0 1 ~--~--~--~ 0.1 I A (order=3) 7 ro 0 05 ...... i-------t-----------1 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 4-8: J 2Im[Ju] and 2Im[l a ] for example 2. Table 4-6: 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 ------+---t------1 ~C'tl l o 6 .....=,~ 0 4 ------+---t------1 2. f 0 2 -------+---t------1 --, -0 2 0 --, C'tl 0 4 r---+J I 2 A ( order=2 ) 3 4 .....=,~ 0 2 ---+-----t------1 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 4-9 : J 2 Im[l u ] a nd 2 lrn[l a ] fo r exa mpl e 3. T a bl e 4-7 : 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 discrete-time 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 4-10 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 4-11 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 4-12 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 4-10: 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 4-11: 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 4-12 : 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 non-convex 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 Wiener-Hopf 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 sampled-data 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 tap-delay line Gamma Laguerre Legendre Jacobi and Kautz models defined in both discrete-time and continuous-time 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 non-causal 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 (non-simulated) 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 5-1 through 5-6. In each case the optimal time scale parameter was determined and used to compute and plot the step response and impulse response for 6-tap (N = 5) networks. Figures 5-1 5-2, and 5-3 compare the system and GLFF network step responses. Figures 5-4, 5-5, and 5-6 show the resulting impulse responses. The optimal time scale parameters of the 6-tap models and associated minimum squared error, J(A o pi) are provided in Tables 5-1 5-2, 5-3.

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 5-1: 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 5-2 : 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 5-3: 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 5-4: 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 5-5 : 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 5-6: 4 th order impulse response and model data 9 10 9 10

PAGE 94

88 Table 5-1: 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 5-2: 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 5-3 : 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) open-loop 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 open-loop tests. Another class of tests studies the behavior of these systems in a clos e d-loop configuration and involves a simulation of the real-world environment and the weapon s flight in this environment. This class of tests can be broken down into two types: all-digital and hardwar e -in-th e -loop (HITL). In an all-digital closed-loop flight the entire weapon system is simulated mathematically. If the weapon subsystem models a re accurate all-digital 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 closed-loop 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 open-air 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 (all-digital 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 open-air 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 5-7 below illustrates a high-level 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 pre-programmed 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 5-8 : Closed-loop missile simulation. 5.2.3 HITL Closed-Loop Missile Simulations In this section HITL closed-loop missile simulations will be discussed. Figure 5-8 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 closed-loop 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 all-digital 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 5-Axis FMS The flight motion simulator (FMS) under study is a CARCO model S-458R-5E five-axis system. This simulator is a precision electro-hydraulic 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 electro-hydraulic actuators. The inner, middle and outer gimbals of the three-axis missile-motion 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 closed-loop control. A high-speed parallel interface allows remote control and monitor of the FMS from a digital simulation computer. This digital interface allows for real-time 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 five-a x is table. It consists of a motor pump and reservoir. As delivered the FMS system has the following capability using a 150-lb load 16 in diameter 30 long mounted on the roll plate. The FMS is illustrated in Figure 5-9. Figure 5-9: Carco 5-axis FMS.

PAGE 103

97 The performance specifications for each gimbal axis are listed in Table 5-4. It is this closed-loop system that will be modeled with several GLFF networks. Table 5-4 : 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 5-10 When identification of only the plant is of interest model design could be carried out with an open-loop experiment. If the controller design is known it i s possible to model the plant in a closed-loop experiment [W a n94a]. This is often more desirable than open-loop 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 open-loop experiments can not be performed without concern for safety and impending system damage ------------n .. . ) u + I: e Controller c Plant V .. ..... .. ...... .......... ... .. .............................................. .............. ........ ........ .............. ...... ............. .... .................. ......... ........... Figure 5-10 : Typical clo s ed-loop system In this thesis it is a model of the entire closed-loop system that is of interest. This is often needed when the system is to be treated as a sub-component 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 closed-loop system is desired. 5 3 1 Discrete-Time GLFF Models of the FMS Figure 5-11 illustrates the general form of a discrete-time GLFF network realization. Three transfer function kernels (Gamma Laguerre Legendre) will be considered for this application. x( k ) y (k) Figure 5-11: Discrete-time GLFF network block diagram. The GLFF formulation of the discrete-time Gamma Laguerre and Legendre kernels were given in Table 3-1. These kernels result in the network realizations illustrated in Figures 5-12 5-13 and 5-14. 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 5-13: 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 wear-and-tear 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 AGM-650 Maverick mjssile. The AGM-650 is a short standoff-range launch-and-leave air-to-ground missile. It is a roll-stabilized cruciform-wing tail-controlled 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 discrete-time GLFF network parameter set using s tep response data of the FMS system H is illustrated in Figure 5-15. 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 signal-to 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 5-15 : 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 two-stage 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 5-16. Carco FMS Pitch Axis Step Response 8 ,-------.------r----r-----,--~----r-----r---~----, 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 5-16: 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 5-6, 57 and 5-8 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 Wiener-Hopf 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 5-21 5-26 and 5-31) step responses (Figures 5-22, 5-27, and 532), and impulse responses (Figures 5-23, 5-28 and 5-33) are plotted for models of order 1 to 8. For the 8 th order models, the step response and associated error (Figures 5-24 5-29, and 5-34), and flyout response and associated error (Figures 5-25 5-30 and 5-35) are plotted. Finally the 8 th order step responses are scaled to degrees and plotted in Figure 5-36 along with the scaled errors in Figure 5-37. 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 5-38 and 5-39 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 closed-loop system. To provide a comparison of the GLFF networks with the conventional tap-delay line Table 5-5 lists the results of several FIR model performance parameters. Figures 5-17 through 5-20 plot the step response impulse response step response error and flyout response using a 201-tap 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 real-time environment. Table 5-5: 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 5-17: 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 5-18: 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 --+----+---------------;-----r----1 0 05 --+-----+----+---+----+---+---+----1 .$ or ---~::::1:== ====+=~-,---i----t--7,--,---1 -; -0.05 ----~ -"w -t ----+------+----t -0 1 1-------t ---t--~-y-i _-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 5-19 : Step response and associated error for 201-tap FIR network. FIR Model Pitch Axis Flyout Response (order=201) 5~------~--~---~--~---~--~ rn o1----i-. ::-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-""" --+---t-----1 .._,. Q) -1 1----1---+ ----+----+--+--+-----t -2 '-------'----'-------'----'-------'----...__--~ 3 4 5 6 7 8 9 10 t (sec) Figure 5-20: Flyout response and associated error for 201-tap 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 5-21: Gamma network performance functions. 0 35 Table 5-6: 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 5-22: 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 .... _,_/ ____ 1-----t----1----------------t 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 t-----i--.~ ~--+----+---+-----t----t---, 2.. .............__ -5 -----~ -----------------t------1 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 5-26 : 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 5-27: 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 5-28: Laguerre network impulse responses.

PAGE 118

11 2 OT Laguerre Model Pitch Axis Step Response (order=8) 1 5 ,-----.------.------.------.------.------r------r----, 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--+---+-----+----t----t------r----t .$ a, -0 05 --+------+---t--+---t---t------r----; -0 1 --t ---+---+ i ----+ i ------1-i --i -~ 0 1 0 2 0 3 0 .4 0 5 0 6 0 7 0.8 0 9 t (sec) Figure 5-29: 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, k---t--+---i ----+----l 0 "--::i -10 f------+----+---------..' -----:,,,,..... ac--+--------+----+--~ -0 ---S-15 ~ --+--+----+----+---~-~= .,,,,. ~====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--"'<;;. ;;r---4-----t-,...----1-----+----~\.-+-------I (l) -1 -----+------------------2 ~--~-----------------~--~ 3 4 5 6 7 8 9 10 t (sec) Figure 5-30: 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 5-31: Legendre network performance functions. Table 5-8: 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 5-32: 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 5-33: 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 ~ --+---+-----4----1-----+-----+-----+-.$ ..... Q) -0 05 ~ --+---+-----4----1---+----+----+---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 5-34: 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 r------t-~ -----t---t--------t-----1 5 -10 1------+----+----+-' '------=--+----+-----t---i "O r----_ 2-15 ~ ---+----+-----+----+--~-=::.::;...~ :::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-----+----+-----+----+-----f-------2 L-__ __._ ___ ......___ __ --L.. ___ ........_ __ ___._ ___ _.___ __ _____. 3 4 5 6 7 8 9 10 t (sec) Figure 5-35 : 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 5-36: 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 5-37 : 8 th order Gamma Laguerre Legendre networks step response errors

PAGE 123

-10 m-15 :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 5-38: 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 5-39 : Transfer function phase estimates for FMS and 8 th order networks.

PAGE 124

118 5 3.4 Continuous-Time GLFF Models of the FMS In this section the GLFF modeling capability will be further illustrated by constructing a continuous-time Laguerre model of the PMS using the sampled step respon s e d a t a. The general form of a continuous-time Laguerre network using the Laguerre kernel is given in Figure 5-40 x(t) 1 1---~---1 s + .-1. s /4 s + /4 s + .-1. s /4 s + /4 y( t ) Figure 5-40: Continuous-time 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 5-41 depicts the squared error performance functions for Laguerre networks from order 0 to 8. For each order network Table 5-9 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 5-42 and 5-43 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 5-44 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 5-45 graphs the responses and associated errors of these models when driven by the AGM-650 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 5-41: Laguerre network performance functions. Table 5-9: 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 ~ :,' ----t----1 C1' ...., 0 4 -, 0.2 -------! -------jf----i . .. .. .. -. . _ : .. .. .. .. ... .. . . .. . . .. .. ... ... . 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 5-42: 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 ........ .... . .. .,_ \ --1f-----1 '-------.. ____ _____ ~--:; 0 1 ~ --~ ---t--1 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 ...... .. ...... ... ~ ---+----t----1 2, 0.1 i J 0 05 ..... ; __ __ __ __ _, 0~ : v \ ;P "'::., ., .. . . . . .. . . ,,, 0 50 100 0 50 100 A A Figure 5-43: 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 5-44: 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 5-45: 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 closed-loop system using measured step response data is possible. The procedure is simple and can be practically implemented in a laboratory setting. 2. A Carco five-axis flight motion simulator used in weapon system hardware in-the-loop 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 closed-loop 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 simulation-by-simulation 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 all-digital 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 complex-valued 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 S-DOMAIN DERIVATION OF MINIMUM J Using Parseval's theorem we can represent the squared error performance function in the s-domain 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)ds-w 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. 815-818 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 737-740 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. CT-4 pp. 286-287 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 351-354 Dec 1959 [Bak81] Baker G. Graves-Morris 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 Addison-Wesley 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. 3400-3405 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. 151-168 Apr 1965 [Bru64] Bruni C. Analysis of Approximation of Linear and Time-invariant Systems Pulse Response by Means of Laguerre Finite Term Expansion ," IEEE Tran s on Automati c Control vol. 9 no 4 pp. 580-581 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. 781-784 Mar 1995 [Chu84] Churchill R. Brown J., Compl ex Variabl e s And Applications McGraw-Hill 1984 [Cle82] Clement P. Laguerre Functions in Signal Analysis and Parameter Identification ," Journal of Franklin In s titut e, vol 313 no 2 pp 85-95, 1982 132

PAGE 139

133 [Clo65] Clowes G., "Choice of the Time-scaling Factor for Linear System Approximations Using Orthonormal Laguerre Functions ," IEEE Trans. on Automatic Control, vol. 10 no. 4, pp. 487-489, Oct 1965 [Clu91] Cluett W ., Wang L., "Modelling and Robust Controller Design Using Step Response Data ," Chemical Engineering Science vol. 46 no. 8, pp. 2065-2077, 1991 [Clu92] Cluett W., Wang L. "Frequency Smoothing Using Laguerre Model ," IEE Proceedings-D: Control Theory & Applications vol. 139 no. 1 pp 88-96, Jan 1992 [Coh92] Cohen H. Mathematics for Sci entists and Engineers, Prentice-Hall 1992 [Con94] Constantinides A., Stathaki P., "Complex Interpolation for Rational Orthogonal Signal Approximation with Applications ," 1994 IEEE ICASSP vol. 4 pp 465-468 [Dav95] Davila C., Chiang H ., "An Algorithm for Pole-Zero System Model Order Estimation," IEEE Trans. on Signal Processing vol. 43, no. 4 pp. 1013-1017 Apr 95 [Dea64] Dean W., Seismological Applications of Laguerre Expansions ," Bull etin of the Seismological Society of America vol. 54 no 1 pp. 365-407, 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 Prentice-Hall 1989 [den93a] den Brinker A., "Adaptive modified Laguerre Filters," Signal Processing, vol. 31, pp. 69-79 1993 [den93b] den Brinker A., "Calculation of the Local Cross-Correlation Function on the Basis of the Laguerre Transform ," IEEE Trans on Signal Processing vol. 41, no. 5 pp. 1980-1982, May 1993 [den94] den Brinker A. "Laguerre-domain Adaptive Filters," IEEE Trans. on Signal Processing, vol. 42, no 4 pp 953-956, Apr 1994 [Dum86] Dumont G ., Zervos C. "Adaptive Control Based on Orthonormal Series Representation," IFAC Workshop Adaptive Systems in Control & Signal Proc essing, pp. 105-113 1986 [Dum90] Dumont G. Zervos C., Pageau G., "Laguerre-based 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. 955-956 Jun 1995 [Els94] Elshafei A. Dumont G. Elnaggar A. Adaptive GPC Based on Laguerre filters Modeling ," Automatica, vol. 30, no. 12 pp. 1913-1920, 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. 923-926, 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 934-938 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. 1068-1073, Dec 1990 [Guo94] Guoqiang L. Dumont G. Sampled-Data GPC (SDGPC) with Integral Action: The State Space Approach ," Proc. of th e 33rd conference on D ec ision & Control pp. 3567-3572 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. 2087-2092 Dec 1993 [Hag91] Hagglund T., Astrom K., "Industrial Adaptive Controllers Based on Frequency Response Techniques ," Automatica vol. 27, no. 4 pp. 599-609 1991 [Ham73] Hamming R., Numerical Methods for Scientists & Engineers 2 nd ed. McGraw-Hill, 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. 640-651 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. 1086-1092, 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 451-465 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. 209-217, 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. 105-126 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. 1075-1078, 1995 [Kab94] Kaber S. Filtering Non-periodic Functions ," Computer Methods in Appli ed Mechanics & Engineering vol. 116, pp. 123-130 1994 [Kau54] Kautz W. "Transient Synthesis in the Time Domain ," IRE Trans. Cir cuit Theory vol. CT-1 pp. 29-39, 1954 [Kin69] King J. O Canainn T. Optimum Pole Positions for Laguerre-Function Models ," El ect roni cs L ette r s, vol. 5 no. 23, pp. 601-602 Nov 1969 [Kin77] King R. Paraskevopoulos P. Digital Laguerre Filters ," Int ern ational Journal of Circuit Th eory & Applications, vol. 5 pp. 81-91 1977 [Kin79] King R. Paraskevopoulos P. "Parameter Identification of Discrete-Time SISO Systems ," International Journal of Control vol. 30, no. 6 pp 1023-1029 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. 2332-2335, Oct 1994

PAGE 142

136 [Kor91] Korenberg M., Paarmann L., "Orthogonal Approaches to Time-Series 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. 373-376 [Lam94] Lam J. "A nalysis of the Laguerre Formula for Approximating Delay Systems, IEEE Trans. on Automatic Control, vol. 39, no. 7 pp. 1517-1521 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. 83-113, 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. 309-312, Jul 1993 [Mai85] Maione B. Turchiano B. "Laguerre Z-transfer Function Representation of Linear Discrete-Time Systems ," International Journal of Control, vol. 41, no. 1, pp. 245-257, 1985 [Man63] Manley H. Analysis-synthesis of Connected Speech in Terms of Orthogonalized Exponentially Damped Sinusoids ," The Journal of the Acoustical Soci ety of America, vol. 35, no. 4 pp. 464-474 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. 17-26 1969 [Mar93] Marmarelis V. "Identification of Nonlinear Biological Systems Using Laguerre Expansions of Kernels," Annals of Biomedical Engineering, vol. 21, pp. 573-589, 1993 [Mas90] Masnadi-Shirazi M. Ahmed N., "Laguerre Approximation of Nonrecursive Discrete-Time Systems," 1990 IEEE ICASSP pp. 1309-1312 Apr 1990 [Mas91] Masnadi-Shirazi M., Ahmed N ., "Optimum Laguerre Networks for a Class of Discrete-Time Systems," IEEE Trans. on Signal Processing, vol. 39, no. 9, pp. 2104-2108, Sep 1991

PAGE 143

, e 137 [McD68] McDonough R ., Huggins W. Best Least-Squares 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. 436-440 Dec 1988 [Moh95] Mohan B. Datta K. Analysis of Linear T i me-Invariant Time-Dela 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. 91-111 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 3388-3393 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. 1601-1603 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. Reduced-Order Models Using Optimal Laguerre Approximations ," El ec tro n i cs L e tt e r s, vol. 23 no 6 pp. 257-258 Mar 1987 [Oli94] Olivier P ., Online System Identification using Laguerre Series ," IEE Proc Part-D: Control Th e o ry Appli c ations vol. 141 no. 4 pp. 249-254 Jul 1994 [Opp89] Oppenheim A. Schaeffer R. Dis c ret e -Tim e Signal Proc e s si n g, Prentice Hall 1989 [Pab92] Pabon-Ortiz 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. 215-221 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. 511-513 Oct 1971 [Par91] Partington J. Approximation of Delay Systems by Fourier-Laguerre Series ," Automati c a vol. 27 no. 3 pp 569-572 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 450-451 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. 752-755, 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, Addison-Wesley, 1987 [Ros64] Ross D. Orthonormal Exponentials ," IEEE Trans. on Communication and Electronics, vol. 1 2, no 71 pp. 173-176 Mar 1964 [Sch70] Schetzen M. Power-series Equivalence of Some Functional Series with Applications ," IEEE Trans. on Circuit Th eory, vol. 17 no. 3, pp. 305-313, Aug 1970 [Sch71] Schetzen M ., Asymptotic Optimum Laguerre Series ," IEEE Trans. on Cir cuit Th eory, vol. 18 no. 5 pp 493-500 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. 455-467 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. 417-420 [Wah91] Wahlberg B. System Identification Using Laguerre Models, IEEE Tran s. on Automatic Control, vol. 36 no. 5 pp. 551-562 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 Closed-Loop Step Response Data ," IEE Proc ee dings-D: Control Th e ory & Appli c ations vol. 141 no. 2 pp 107-110, Mar 1994 [Wan94b] Wang L. Cluett W. Optimal Choice of Time-Scaling Factor for Linear System Approximations using Laguerre Models ," IEEE Trans on Automati c Control vol. 39 no. 7 pp. 1463-1467 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 149-161 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. CT-9 pp. 152-162 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. 10-18 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. 362-370 1962 [You63] Young T. Huggins W ., On the Representation of Electrocardiograms ," IEEE Tran s Bio-M e di c al El ec troni c s vol. BME-10 pp. 86-95 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. 165-175 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