Single-channel pilot manual control output in closed-tracking tasks is modeled in terms of linear discrete transfer functions which are parsimonious and guaranteed stable. The transfer functions are found by applying a modified super-position time series generation technique. A Levinson-Durbin algorithm is used to determine the filter which prewhitens the input and a projective (least squares) fit of pulse response estimates is used to guarantee identified model stability. Results from two case studies are compared to previous findings, where the source of data are relatively short data records, approximately 25 seconds long. Time delay effects and pilot seasonalities are discussed and analyzed. It is concluded that single-channel time series controller modeling is feasible on short records, and that it is important for the analyst to determine a criterion for best time domain fit which allows association of model parameter values, such as pure time delay, with actual physical and physiological constraints. The purpose of the modeling is thus paramount

Biezad, D. J.

Leban, F.

Mashiko, S.

Schmidt, D. K.

English

NASA Technical Reports Server

dified Superposition: A Simp 
ries Approach to Closed-Loop 
Controller Identification 
Biezad, Lt. Col,  USAP 
fessor David K. Schmidt 
r, Frank Leban 
shiko, Cap+-. USAF 
Purdue University 
Institute of Technology 
trical and Comput Engineering 
terson APB, Ohio 3 3  
ABSTRACT 
Sing le-channel manual control output in closed- 
tracking tasks i modeled in terms of linear discrete transfer 
functions which are parsimonious and guaran ed stable. The 
transfer functions are found by applying modified super- 
position time series generation technique. Levinson-Durbin 
algorithm is used to determine the filter h prewhi’tens the 
input and a projective (least squares) fi pulse response 
estimates is used to guarantee in entified model stability. 
Results from t case studies d to previous 
findings, e source of d ely short data 
records, approximately 25 seconds long. Ti elay effects and 
pilot seasonalities are discussed and It is 
concluded that single-channel time r modeling 
is feesible on short records, and that it is important for the 
analyst to determine a criterion for “best time domain fit” 
which allo s association cf model parameter values, such as 
pure time d actual physical and physiological 
constraints. e” of the modeling is thus paramount. 
33.1 
https://ntrs.nasa.gov/search.jsp?R=19860023537 2020-03-20T13:13:54+00:00Z
: AUTOREGRESSIVE PILOT 
pilot modeling 
utoregressi 
closed loop tern identific 
prewhitening 
superposition 
single input, ingle output 
nual control 
(t) 
T 
v 
w 
tor discrete polynomial in z 
coefficient of z - ~  in b (  
pilot model pul 
splayed to pilot 
coefficient of z-k in g ( z )  
pilot transfer function as a ratio of polynomials 
independent, identically distributed 
lag implying "kA" seconds 
t gain expressed in degrees p 
1 points available 
t input uncorrel ith y(t) in degrees at 
hite noise sequence (i.i.d.) at instant t 
controlled element output signal in degrees pitch 
angle at instant t 
sample interval (seconds) 
pilot output in degrees of elevator deflection at 
instant t 
number of sample times in pure time delay 
sformed frequency 
ing filter in 
33.2 
1. ~ N T R O D U C T ~ O N  
ey question of ho the human bein will be inser- 
control loop of remains an issue 
throughout our t nowhere is it 
more urgent tha ign and analysis 
that a pilot of a modern aircraft is 
systems monitor (Rouse, 1 9 8 3 )  in no 
a controll r (Rouse, 1980;  Sheridan, 
assumption in thi ork is that the 
and machine should under stood much 
almer, 1 9 8 3 ) ,  
Although describing cRuer, 1 9 6 5 )  and optimal 
control (Kle.inman, 1969-1 have been inge- 
niously used to provide insight into piloting strategy 
(Schmidt, 1999;  Bacon, 1983;  Hess, 19791,  they are now supple- 
mented with pilot models d rived from the emerging field of 
time series analysis. Time series modeling of pilot behavior 
offers tremendous potential for discerning ey system charac- 
teristics and relationships, such as the actual effect of 
instabilities (Goto, 19 ) ,  pilot stress (Shinners, 1 9 7 4 ) ,  or 
task effects (Agar 
The key questions in time series models involve not only 
the parsimony of parameters, well established y Breddermann et 
a1 ( 1 9 9 8 1 ,  but of identified model ty and the' model's 
practical application in analysis Shinner 
1 seriously discussed the c ntification 
problem, but the manipulation of ons in his 
fitting procedure contains no gua del sltabi- 
lity. The primary purpose of this ent a theore- 
tically sound relatively loop fitting 
procedure, stil sed firmly in methods of 
ox and Jenki 1976 1 ,  which stability 
ithout sacrific model accuracy. 
2 .  MODEL 
The linear discrete closed-loop model structure is shown 
in Figure 1 ,  Each block represents a discrete pulse response 
en convolved discrete input 
discrete outp . Stable pulse 
sequences, even though infinite in durations, eventually must 
decay for a stable system. When the pulse sequence is ex- 
pressed as a ratio of polynomials, stability is guaranteed if 
the denominator root are less in magnitude th n one. The goal 
33.3  
is to identify the pulse response sequence gp(z) and approxi- 
mate its discrete (z domain) transfer function from actual data 
sets { 6 ( t ) } ,  {y(t)], and {W(t)} which are equispaced in time 
ith their means removed. 
The assumptions are model linearity, time invariance, 
causality, uncorrelated inputs W(t) and R(t), and prewhitenable 
input W(t); that is, W(t) is a linear function of previous 
values plus a white noise "shock." Previous values are mathe- 
matically linked by the backward shift operator z ( - l  I. 
ODIFIED SUPERPOSITION TECHNIQUE 
First, every signal in Figure 1 is decomposed conceptually 
into a part linearly correlated with command disturbance W(t), 
the remainder uncorrelated with W(t). For example, output 
y(t) is the sum of Y L ( ~ ) ,  which is correlated with W(t), and o f  
YR(t), considered the effect of an additional unknown input 
R(t), termed "Remnant," uncorrelated with W(t). The pulse 
response to be found relates, for constant sampling interval 
11 a 11 seconds, the linearly correlated pilot output bL(t) to the 
correlated error signal eL(t); that is, 
This pulse response may be expressed as an infinite sequence or as 
a ratio of polynomials: 
m 
where 
and s 11 imposed constraint 
33.4  
If the integer "k" in Equation (2) is allowed alL values 
( -OD < k < + ~ ) ,  then Equation ( 2 )  defines the discrete transfer 
function relating the z-transform of input sequence eL(t) to 
the z-transform of output sequence d~(t) (Franklin and 
e11,1980, p.15), 
the 'L(t) and eL(t) are not directly 
available, they must be "generated" i f  loop closure effects are 
properly taken into account. T o  do this apply superposition to 
signals Y(t) and (t) of Figure 1 :  
Since W(t) is prewhitenable (defined above) and uncorrelat- 
ed with R(t) the cross correlation identification technique of 
Box and Jenkins (details in Appendix) may be applied to Eind 
an estimate of the initial portion of the pulse response sequence 
g i ( z ) ,  between y(t) and W(t). Then al(z) and bl(z) may be 
determined as shown in the section on model stability, such 
that 
The essence of modified superposition is now to qenerate the time 
series YL(t) using the autoregressive relation 
bl(z)YL(t) = al(z)W(t) ( 1 0 )  
Where al(z) and bl(z) are numerator and denominator polyno- 
mials, 
respectively, with the structure of Equations ( 3 )  and (4). The 
linearly correlated ,signal eL(t) is then generated from 
eL(t) =W(t) - YL(t) (11) 
The above process is then repeated by reapplying super- 
position to obtain the following relation between 6(t) and W(t): 
33.5 
The cross correlation identification (Appendix) applied to,  the 
sequence 6 ( t )  and W(t) yields the initial segment of pulse 
response sequence g 3 ( z ) ,  and the polynomials a3(2) and b 3 ( z )  may 
be determined (see next section) such that 
Pilot output linearly correlated with W(t) is generated from the 
autoregressive relation 
Finally, the cross correlation technique (Appendix) is 
applied to dL(t) and eL(t) to find the initial segment of gp(k), 
defined by the' coefficient set { g p ~ ,  05 k <N], of the pilot model 
p u l s e  response .  Numerator and denominator polynomials a r e  then 
found (see next section) which yields 
~ ~ ( 2 )  = 6L(z)/eL(z) ( 1 6 )  
No multiplication or divisions of transfer functions o c c u r s  
throughtout the above procedure. 
4 .  MODEL STABILITY 
A s  mentioned above, the pulse response sequence idcntif i- 
ed [ g ~ ( z ) ,  93(z) and gp(z)] will be truncated at some finite 
l a g  "k" final task is to find a parsimonious numerator polyno-  
mial stable denominator polynomial which together ace equivalent 
mathematically to the identified pulse response. These polyno- 
mials are chosen to have the structure shown in Equation ( 2 1 ,  
which is re-arranged into the following form: 
> kmax > s  -k 
33.6 
S i n c e  t h e  p u l s e  r e s p o n s e  3k  i s  known f o r  0 < k  < N ,  by e q u a t -  
i n g  c o e f f i c i e n t s  f o r  t h e  o p e r a t o r  'z" a t  each e x p o n e n t i a l  
power"&", o r e l a t i o n s h i p s  may be f o u n d  b e t w e e n  n u m e r a t o r  a n d  d e n o -  
m i a t o r  c o e f f i c i e n t s  a k  a n d  bk .  M o r e o v e r ,  by  e q u a t i n g  c o e E E i c i e n t s  
f o r  t h e  o p e r a t o r  z a b o v e  power f c s " ,  f o r  w h i c h  t h e  r i g h t  s i d e  of  
e q u a t i o n  ( 1 5 )  v a n i s h e s ,  o n e  o b t a i n s  f o r  e v e r y  j > 0 
T h e  a b o v e  r e l a t i o n  e x i s t s  f o r  a f i n i t e  b u t  l a r g e  n u m b e r  of  
'J > O B ' ,  so projec t ion< t h e o r y  ( l e a s t  s q u a r e s )  may be used  t o  coef -  
f i c i e n t s  b k  ( 0  < k  - S I .  B r i n g i n g  term n g ( S + j ) l '  t o  t he  o t h e r  
s i d e  O f  e q u a t i o n  ( 1 8 )  a n d  d i v i d e d  by t ' g ( S + j ) n  o n e  may write 
T T 
A [ b l r  b2, ..., b s ]  = [ - 1 ,  - 1 1  
a n d  t h e  "jWth row of A is g i v e n  b y  
I g s + j - 1  , g.s+j-2 , ..., g s + j - s  I j > O  
g s + j  gs+ j g s + j  
1 by  s 
The s o l u t i o n  f r o m  l i n e a r  a l g e b r a  is 
( 2 0 )  
To p r o v i d e  a p a r s i m o n i o u s  d e n o m i n a t o r ,  t h e  s o l u t i o n  of 
E q u a t i o n  ( 2 1 )  is accepted f o r  t h e  l o w e s t  o r d e r  " s "  w h i c h  has 
b o t h  a s table  c h a r a c t e r i s t i c  e q u a t i o n  (1.e.  r.oots l e s s  t h a n '  1 . 0  
i n  m a g n i t u d e )  a n d  which y i e l d s  a model p u l s e  r e s p o n s e  s i m i l a r  
i n  s h a p e  t o  t he  t r u n c a t e d  p u l s e  r e s p o n s e  i d e n t i f i e d  from t h e  
d a t a .  O n c e  a s t a b l e  d e n o m i n a t o r  i s  f o u n d  t h e  n u m e r a t o r  a ( z )  
a n d  t h e  g a i n  K may be d e t e r m i n e d  b y  o n c e  a g a i n  m a t c h i n g  coef -  
f i c i e n t s  i n  E q u a t i o n  ( 1 7 ) ;  
K = 9, 
33.7 
sidual to be the actual output time 
model output series at each sample 
adjusted by a suitable minimization 
the error residual variance. 
adjusted to provide a steady state 
input to the transfer function is a 
int recommended by Agarwal ( 1 9 8 0 ) .  
is to be included, the final form 
tion ( 2 1 ,  and the indices for 
( 1 7 ) - ( 2 3 )  should be incre- 
ntification (for example 
g, identified from the 
also be applied to the model. There 
ceptability and statistical signifi- 
are common sense checks which com- 
eries verses actual autocorrelation 
tocorrelation of residuals €or whi- 
cks for neg 1 ig ible cross-co r r e lat ion 
gnificance tests may be performed after 
indicate the model is reasonable. Chi- 
e available Erom the (t) and v(t) prewhi- 
sed in the Appendix and shown in Figure 
n safely neglect correlations beyond a lag 
to be computed are, for 
nd, for uncorr 
ould pass the chi-squared significance test for degrees 
(20-1-s-1)  respectively ( B o x  and Jenkins, 
of either significance test is evidence 
or a modeling inadequacy. 
33.8 
0 odified su erposition 
equence relating y(t) and W ( t )  using 
on (Appendi 
unction G1(I) 
e least squares 
the data g 1 ( z )  
(t)], (eL(t)] using 
ng 6(t) and 
etermine a 
e response 
nse g , ( z ) ,  from 
identification, 
e QT s ION 
the Flight 
re accomplished 
using a single 
1 controlled 
disturbance 
nd a standard . Data s e t s  
500 points were 
ly 25 seconds 
ent many low- 
" and the 
ect iden- 
identification 
pilot injected 
reasonable 
33.9 
p le"  c o n t r o l l e d  element s u c h  as K / s .  An a 
p r i o r i  selecte d e l a y  of 0 . 2  s e c o n d s  y i e l d e d  t h e  l o w e s t  
r e v i o u s  r e s u l t s  
r t h a t  a d e l a y  
s p i l o t  beha-  
l i n g  r e s u l t s  
r e t a t i o n  and  
, a t r a n s f o r -  
c c o m p l i s h e d  u s i n g  
t h e  t r a n s f o r m e d  f r e q u e n c y  (u) 
1 Bode i n t e r p r e -  
5 w i t h  F i g u r e s  2 
one  can n c e  be tween  t h e  
( 0  w 2 5  r p s ) .  
Only  t h e  f i r s t  
model ,  *and t h e  
r e a s o n a b l y  ac beyond t h i s  time. T h i s  
e r f i t  (Kashyap ,  1 9 7 6 ) ,  
i n c r e a s e d  e r r o r  r e s i d u a l  when t h e  
i n d e p e n d e n t  of  model d e r i v a t i o n i i n  
F o r  t h e  d o u b l e - i n t e ~ r a t o r  c o n t r o l l e d  e l e m e n t  a more 
as i d e n t i f i e d  and  is shown i n  T a b l e  
riori a e l e ~ ~ e ~  t i m  y ( 0 . 0 5  seconds 
y r e s p o n s e  p l o t  i n  F i g u r e  7 t h e r e  is  some 
he p h a s e  p l  a re  shown i n  F i g u r e s  8 
( 0 . 2  and 0 . 0 5  
h e  t r a n s f o r  ' domain y i e l d s  
f rom t h e s e  nd t h e y  a r e  n o t  
33.10 
In contrast to control of a "simple" K/s, the "best fit 
(minimizing residual error variance) was obtained when time 
delay was set to 0 . 0 5  control for of K/s2. The phase contribu- 
tion (from the poles and zeros) of the discrete transfer func- 
tion is apparent, as time delay changes between 0.2 seconds and 
0 .05  seconds, as may be seen by the phase plots of Figures 10 
and 11 in which the pure time delay has been removed from the 
discrete- transfer function. Selecting the larger pure time 
delay for the model exposes the considerable lead generation 
from the transfer function poles and zeros. This lead genera- 
tion is not a s  apparent when pure time delay is reduced for the 
"best time domain" fit, but the resulting 0.05 seconds might be 
judged too fast to associate with a lumped physiological delay 
for a human operator. A possible explanation is unmodeled 
pilot anticipation; that is, a possible anticipatory loop clo- 
sure not accounted for in Figure 1. 
Further evidence of this is provided in the time history 
for the best fitting model in Figure 12. Note that a seasonal 
pilot residual (where pilot output "leads" model output) occurs 
during some of the longer intervals of large slope. This could 
be caused by momentary anticipatory behavior arising from the 
"pursuit" display including commanded input, a factor not 
accounted for in a time invariant model. Thus in determining 
the "best" model usinq time series analysis, the purpose of the 
model must be given as much consideration as tests for "best 
fit." 
In summary for the K/s2 controlled element, an a priori 
time delay in series with a rate sensitive gain describes "pilot" 
behavior over his usable bandwidth, in agreement with classical 
results (McRuer, 1974). When pure time delay is not set a 
priori but allowed to vary in obtaining the "best time domain 
fit," the minimization of an error variance criterion results 
in a math model where the time delay is perhaps too small to be 
associated with physiological operator delays. This case is 
associated with a pursuit task in which the command as well as 
the plant output is displayed. 
33.11 
6 .  CONCLUSIO 
osition technique as described for 
discrete transfer function, 
model valid tion. Results 
t the time series technique appears 
on "short" data records. The analyst 
etermine the criterion for a "best time 
association of parameter values, such 
th actual physical and physiological 
pilot residual, possibly caused 
re observed as first noted by 
11 modebed time invariate 
ould concentrate on the full potential of 
dels for analyses, especially their ability 
d accurate power spectral densities, and on 
multi-channed closed-loop pilot model-ing. 
7 e ACKNOWLEDG 
uthors express their gratitude to r. Yuan Pin-Jar, 
ible for the laboratory computer programs and 
used in the simulations, and to Mr. Chuck 
tance in computer operation, maintenance, and 
33.12 
8 .  APPEND1 X :  Cross Correlation IdentiEication 
( B o x  and Jenkins, 1976) 
Given the situation in Figure 1 3 ,  the goal is to find the- 
pulse response relating Y(t) and W(t), which i s  prewhitenable 
by u w ( z ) .  The prewhitening is accomplished by applying the  
Levinson-Durbin algorithm as given by Kay and Marple (1981, pp. 
1 3 8 8 - 1 3 8 9 ) .  By reversing the order of the blocks in the for- 
ward path of Figure t and multiplying each signal at the 
summer by 
1, the €0110 ing equation results: 
(t) + $Iw -'(z)V(t) = s ( t )  ( 2 8 )  
d ( L )  = iiw --hz)y(t) ( 2 9 )  
Now multiply Equation ( 2 8 )  by w(t-k) and take the expecta- 
tion, recalling that w(t) is uncorrelated by assumption 
with v(t): 
By expanding G(z) using shift properties of z one obtains 
Since w(t) is an independent, identically distributed sequence of 
random numbers with variance ow2 , one obtains for every lag k 
k > O  
Conventional estimation relations may now be used to estimate 
the terms in Equation ( 3 2 )  and solve for gk; for example, from B o x  
and Jenkins ( 1 9 7 6 ,  pp. 3 2 - 3 3 )  one obtains 
N N 1 
N t=l Id-k t=k 
( 3 3 )  ik { -  (t)} = i2 1 B(t)w(t-k)} 
which detemines the pulse response sequence estimate gk. 
33.13 
REFERENCES 
1. trol Approach to 
riteria, ' A I A A  
8 VOl. 6, NO. 5, 
d Kleinman, Closed Loop Models for 
e ~ ~ ~ r e m e n ~  Simulators," NASA 
e Identifiability 
tern," Identifica- 
an ed., North 
On the Dynamics of Human Pilots in Marginally - 
-Journal of Aircraft, V O ~ .  12, NO. 
7, ling Quabiti of Plight V hicles, "*  internal 
t o  CALSPAM A vanced Technology 
YOKk. 
r@dlCtiQR Of pilot 0 inion Rating Using an Optimal 
Human Factors, V 0. 5, Oct. 1977, 
1 1 ,  
d Levinson 
tornatica, 
33.14 
CRuer, Graham, Kr 1, and Reisener, Human Dynamics in 
ompensatory ~ ~ F ~ L - T ~ - 6 5 - 1 5 ,  ight-Patterson AFB, 
Ohio, July 1 
Human Pilot 
15, Palmer, Cybernetics , Vol. 
, Rosenbrock, H e ,  r Use of Hu bility:: A 
t e r ,  No. 5, Oct. 
terns Engineering achine Inter- 
ction, North Holland, 1980 ,  
. Schmidt, D. ,  Optimal Flight Control Syntheses Via Pilot 
IT PKesS,  1 9 7 4 .  19, Sheridan and Ferrell, 
0. Shinners, S., ling of Human 
I 
c -  33.15 
Kp -'I ( l + a l z - l )  
( l+blzp-l) 
seconds T = 4 (0.2 seconds) 
K* 
P 
K* 
P 
1 
1 
0.64 
0.79 
0.71 
0 .32  
-- 
D i  rec t 
Identification 
- 
0.69 
0.72 
0.69 
0.37 
i ch  minimizes error residual variance 
elds steady state step response o f  u n i t y  
33.16 
dified Superposition Identification Resul ts  
for Controlled El n t  K / S ~  ( K  = 1) 
Kz -' ( 1 +a 2 -  ) 
- 3- - +b2 E -  2*b3 z ) 
I S ~ r ~ c t ~ r ~  G ( z )  = 
P 
Signal t o  noise ratio : 30 
5Q0 A .= 0.05 5 conds T = 4 (0.2 seconds) 
KP 
Up** 
a 1  
bl 
b2 
b3 
Roots 
0.03 
0.033 
10.9 
0.91 
- 0.1 
0.14 
0.64 tj 0.57 
1 = 0.2 secorids 
0.89 
1 . 2 2  
-0.67 
-1 .41  
0.88 
-0.06 
0.08 
0.67 &J  0.57 
ch minimizes error residual var iance  
I d s  s ~ ~ a ~ y  sta nse of u n i t y  
33.17 
c 
& 
33.18 
20 .oo 
10 .oo 
0 .oo ca 
cf 
-10 .oo 
-20 .oo 
Figure 2 Manuel Controller Frequency Rerpon 
Controlled Element K/r 
10 -I 10 * 10 10 * 10 
FREQUENCY-RQD PER SEC 
33.19 
Figure 3 Manual Controller Frequency Rerpon 
Controlled Element K/r 
-90 .o 
-180.0 
w 
W 
E-270 .O 
3 
33 .20  
20 .oo 
10 .00 
0 .oo c3 
(I 
Figure 4 
10 
TRR 
-1 
33.21 
.0 
-90 .o 
A 
c3 
w 
.Q 
Figure 5 W' Response Pha 
10 -' 10 O 10 10 10 
T R ~ N S F ~ R M E ~  FRE - RFID PERS 
33.22 
0 
0 
0 
L n  
N 
0 
0 In 
0 
0 In . 
0 
0 
U3 
0 
0 
L". 
0 
0 
L n  
Lu 
I 
33.23 
-10 .oo 
-20.00 
Figure 7 Msnual Controller Frequency Response M 
K / S ~ .  I = 0.2 seconds 
20 .oo 
10 .oo 
0 .oo c3 
r: a 
10 -' 10 O 10 10 10 
FREQUENCY-RRD PER SEC 
33.24 
spon 
3 
33.25 
ontroller Frequency Rerponre Ph 
0.0 
-90.0 
in 
w 
w 
w I 
0 -I 10 
33.26 
igure 10 Manual Controller Frequency Pole-zero R e s  
K / s ~ ,  T - 0.2 seconds 
33.27 
90.0 - 
0.0 - 
-90 .0 - 
W 
E-180.0 - 
W 
n 
-270 .O - 
Figure 11 Manual Controller Frequency Pole-zero Rerpon 
K h 2 ,  r = 0.5 seconds 
33.28  
Ir! 
0 
U 
b- 
'y, --. 
Y .. 
e-. 
Y 
Q v 
3 
0 
w 
3 
Q 
Y 
w 
0 - 
0, 
U 
0 
I 
33.29 
33.30 


Modified superposition: A simple time series approach to closed-loop manual controller identification

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server