National Air Traffic Services (NATS) are concerned with ensuring low probabilities of errors in determining aircraft positions. In general, error probabilities depend on the tails of some probability distributions for which there has been no theoretical model. Analysis of radar performance is regularly undertaken by NATS to ensure radar performance is within safety limits, with the maximum range being dependent on the declared separation between aircraft. NATS brought two questions to the Study Group, involving the horizontal (azimuthal) errors in radar data and the vertical errors in altimetry system data. In both cases, NATS asked the Study Group to analyse the data and assess whether the probability distributions that are currently used are good models for the errors

Andre Leger, Andre

Budd, Chris

Hill, Andrew

English

Mathematics in Industry Information Service Eprints Archive

Distribution-independent safety analysis
Problem presented by
Keith Slater
National Air Traffic Services
Problem statement
National Air Traffic Services (NATS) are concerned with ensuring low
probabilities of errors in determining aircraft positions. In general, error
probabilities depend on the tails of some probability distributions for which
there has been no theoretical model. Analysis of radar performance is
regularly undertaken by NATS to ensure radar performance is within safety
limits, with the maximum range being dependent on the declared separation
between aircraft. NATS brought two questions to the Study Group, involving
the horizontal (azimuthal) errors in radar data and the vertical errors in
altimetry system data. In both cases, NATS asked the Study Group to
analyse the data and assess whether the probability distributions that are
currently used are good models for the errors.
Study Group contributors
Chris Budd (University of Bath)
Jeff Dewynne
Andrew Hill (University of Bath)
Sam Howison (University of Oxford)
Andre´ Le´ger (University of Bath)
Henrik Petersen (Maersk Institute, SDU)
Colin Please (University of Southampton)
Jonathan Rougier (University of Durham)
Report prepared by
Chris Budd (University of Bath)
Andrew Hill (University of Bath)
Andre´ Le´ger (University of Bath)
C-1
1 Introduction
NATS brought to the Study Group two problems involving aircraft horizontal and
vertical overlap safety margins. In the first problem, azimuthal error data is currently
modelled by a double exponential distribution, where the central density is captured
by one exponential and the tails by the other. The horizontal overlap probability
(the probability that two aircraft are on the same bearing) is determined from tail-
tail convolution. We have discovered the azimuthal outliers arise through random errors
but also from systematic errors introduced by the way the data is processed. Gibbs
phenomena are seen as a result of interpolation error in plotting aircraft paths and we
have shown they can be removed by post filtering. In the second problem, vertical error
data is currently modelled by a mixed Gaussian exponential distribution, and here we
have shown there is insufficient outlier data to fit any distribution, making any estimate
of the vertical overlap probability invalid.
2 Problem 1: horizontal safety margins
2.1 Introduction
The National Air Traffic Services are responsible for monitoring and ensuring safe ranges
of distance between aircraft. The safety system at NATS requires all radars at its 22
sites to operate at a minimum acceptable level of declared range safety (Target Level
of Safety). A typical declaration may take form ‘let radar X support 5Nm separations
between pairs of aircraft at any range up to 120Nm from the radar’. Radar performance
is inconsistent and varies depending on factors such as icing at the radar, other weather
factors, wind farms and modification of the radar itself and its associated equipment.
NATS apply the following process to determine the maximum safe range:
• Several hours of azimuthal and range data are recorded from many radars
simultaneously.
• The recorded data is post-processed to determine the true position of targets and
hence the individual position azimuthal errors θ for each radar return and for each
radar (see Figure 1).
• The error data are then analysed statistically by another partly automated process
to produce a maximum safe range estimate for a particular radar.
At the last stage, the distribution of azimuthal errors is currently fitted by NATS as a
sum of two symmetric exponentials:
f(θ) =
1− a
2λ
e−
|θ|
λ︸ ︷︷ ︸
Core
+
a
2µ
e−
|θ|
µ
︸ ︷︷ ︸
Tail
(1)
for means λ, µ and weighting coefficient a  1, ensuring the central part of the density
is dominated by the first term and the tails by the second. This mixed exponential
distribution is an empirical choice, and there is no theoretical underpinning for the actual
C-2
Figure 1: Determining azimuthal error
Figure 2: Mixed exponential model
distribution of the data. However, a fitted distribution is required to estimate Horizontal
Overlap Probabilities (HOPs), that is the probability two aircraft are overlapping when
the display in fact shows them have at least the minimum declared separation (e.g. 3, 5
or 10Nm). In terms of a convolution, we have
HOP = P (|Θ1 −Θ2| ≥ θ0) = f ∗ f(θ0) , (2)
where Θ1 and Θ2 are independent samples from the distribution. NATS determine HOP
from the best fit line of the mixed exponential and consequently can find the maximum
safe range. However, the HOP depends on the tail probabilities of the convolution (2),
which are poorly understood, and curve fitting, which is very subjective.
2.2 Mixed exponential distribution of data
Here we test the appropriateness of applying a mixed exponential to azimuthal error
bin data collected from a radar station at Gatwick Airport. All data was exported to
Matlab [1] for analysis, making use of Matlab’s graphical packages. Taking the logarithm
of the error versus the scaled range, we plotted the blue graph in Figure 4, where we
have superimposed in green the logarithm of a mixed exponential distribution. From
the approximate straight line in the region [−0.2, 0.2] we see the core is approximately
exponentially distributed, however, the tail data looks too random to be modelled
by a straight line. The current model assumes tail data is random and identically
C-3
Figure 3: Horizontal overlap convolution
−1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1
−6
−4
−2
0
2
4
6
8
10
12
Log histogram of azimuth error and mixed exponential
Figure 4: Range and azimuthal data for Gatwick radar
independently distributed (iid). However, upon closer examination of the data this
assumption is false.
2.3 Analysis of the azimuthal and range data
Error data was provided for 7 radar stations situated around the UK’s main airports and
across the UK1. As a aircraft comes into range of the radar and is acquired, it is given
an id number. Azimuthal and range data was organised as a function of id, which is
analagous to time. The rotating radar will pick up the aircraft again on its next rotation,
assigning it another id number. Gaps in id number refer to lapses between acquiring
a new aircraft or losing track of the current aircraft. Looking at the Gatwick data, we
plot in Figure 5 the azimuthal data as a function of id and superimpose range data
scaled by maximum range. Data was also manipulated using Matlab to obtain figures
6 and 7. Highlighting a region of ranged data when there is large azimuthal error and
zooming in, we notice two types of radar-aircraft interaction: (1) aircraft are picked up
by the radar at long range and fly-by resulting in a rectangular hyperbola on the range-id
graph; (2) aircraft are suddenly acquired at short range and followed until they are out
of range, seen by a ‘tick’ on Figure 6(a). Looking at Figure 6(b) it is clear that random
error coincides with the start of the acquisition of an aircraft path, resulting in highest
azimuthal error.
1The 7 stations were Burrington, Claxby, Gatwick, Heathrow10, Heathrow23, Jersey and Mount
Gabriel.
C-4
1 2 3 4 5 6 7 8 9 10 11
x 105
−0.5
0
0.5
1
1.5
2
id
Scaled Range
Azimuthal Error
Figure 5: Log histogram of Gatwick data
6.29 6.295 6.3 6.305 6.31 6.315 6.32 6.325 6.33 6.335
x 105
1
1.2
1.4
1.6
1.8
2
id
Gatwick
(a) Typical radar-
aircraft interaction
5.66 5.68 5.7 5.72 5.74 5.76 5.78 5.8 5.82
x 105
−0.5
0
0.5
1
1.5
2
id
e
rr
o
r/s
ca
le
d 
ra
ng
e
Gatwick
Scaled Range
Azimuthal error
(b) Zoom on region of
large error
5.749 5.7495 5.75 5.7505 5.751 5.7515 5.752 5.7525 5.753 5.7535
x 105
−0.6
−0.4
−0.2
0
0.2
0.4
0.6
0.8
1
1.2
(c) Acquired aircraft
azimuthal error
Figure 6: Tracking azimuthal errors
The width of the hyperbolae seen in Figure 6 is determined by the distance of the aircraft
from the radar. Consider a fly-by, where the range from aircraft to radar is denoted by
r, the speed of the aircraft by v, and L is the minimum distance to the radar which we
scale to time t = 0. Then the range is given by the hyperbolic equation r2 = L2 + v2t2.
L
vt
r
Radar
Figure 7: Calculating range on a fly-by
Now by zooming into the azimuthal error when an aircraft is acquired, shown in
Figure 6(c), we see oscillatory errors resembling Gibbs phenomena. This is often related
to interpolating square functions with polynomials (for example, see Burden and Faires
[2], Chapter 3.4). Hence applying Fourier analysis, we take the fast Fourier transform
(FFT) of the azimuth error. The FFT plotted against id is shown in Figure 8. It is very
apparent from the FFT there is a signal peak in this data corresponding to systematic
C-5
Figure 8: Fourier transform of Gatwick error data
oscillatory error with very marked peaks, quite different from what one would expect
from random errors. Given that there is a total of 157, 000 data points (N), we can
determine the frequency of oscillation from the location of a peak. Looking at the peak
occurring at 13, 000, the frequency of oscillation is such that
fpeak =
idpeak
N
=
13000
157000
≈ 1
12
. (3)
Hence, there are 12 points per wave from this systematic error (see Figure 6(c)). The
same Gibbs phenomenon is seen in the data for the other radar stations. We now
investigate the source of error, by looking at the way radars acquire aircraft.
2.4 Systematic source of error
Azimuth and range data were provided for each radar station. In general most radars
tracked aircraft flying at close range (< 180Nm), but a few tracked aircraft at both close
and far range (> 1, 000Nm). We hypothesise that on take-off, the aircraft is tracked
by an airport radar station capturing the aircraft on every rotation. The path of the
aircraft must then be interpolated, leading to large azimuthal errors from close-range
observation. However, the farther the aircraft is from the airport the more likely it will
be picked up by another radar from long range, providing more radar traces and a better
interpolated path.
Let us consider a hypothetical path of an aircraft described by f = 400(1−e−t ln(2)/20).
Radar readings are taken every ∆t seconds and we perform a cubic spline to approximate
the path. After 20 seconds, the aircraft is picked up by a second radar and we
approximate the path again with a cubic spline.
Hence, acquiring aircraft introduces interpolation errors which become recorded as
false azimuth errors in the data.
2.5 Removal of systematic error
We apply a band-gap filter to remove the component of frequency 1/12 in the FFT data
(also seen in the Heathrow radar). An inverse FFT is then applied. In Figure 10 we
C-6
5 10 15 20 25 30 35 40
−0.6
−0.5
−0.4
−0.3
−0.2
−0.1
0
0.1
0.2
0.3
0.4
time
di
st
an
ce
 e
rr
or
(a) Simulated spline
error
(b) Aircraft path tracked by radar
Figure 9: Simulated spline error in acquiring an aircraft path by radar
plot our original noisy azimuthal data (in blue) with the new filtered data (in red). We
8.484 8.485 8.486 8.487 8.488 8.489 8.49 8.491 8.492 8.493 8.494
x 104
−0.5
0
0.5
1
1.5
2
2.5
id
a
zi
m
ut
h 
er
ro
r/s
ca
le
d 
ra
ng
e
Heathrow23cm 
scaled range
unfiltered error
filtered error
(a) Filtered and noisy
azimuthal error
−1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1
−6
−4
−2
0
2
4
6
8
10
12
Log Histogram of azimuth error
Unfiltered
Filtered
(b) Log histogram of
filtered error
Figure 10: Filtered errors
have clearly removed the systematic error in Figure 10(a) and taking the log histogram
of the azimuthal error (Figure 10(b)) shows a dramatic change in the tail distribution.
The clearer data appears to have a core normal distribution, and in Figure 11 we test
the sample data against a standard normal using a quantile-quantile plot. In fact, the
core data reasonably fits a normal distribution, but the tail data does not.
2.6 Conclusion and recommendations
Analysis of the azimuthal error shows that systematic error is responsible. Errors in
interpolating aircraft paths lead to large errors in azimuth data. Applying an FFT,
filtering out these errors, and using an inverse FFT yields clearer data, reducing the
magnitude of the largest errors. At present NATS are using noisy data to calculate the
HOP, and are producing safety margins that are larger than necessary.
Whilst we can remove the error we could also lose information, so we recommend to
NATS that in the calculation of aircraft paths they remove the systematic error by using
C-7
−5 −4 −3 −2 −1 0 1 2 3 4 5
−0.8
−0.6
−0.4
−0.2
0
0.2
0.4
0.6
Standard Normal Quantiles
Qu
an
tile
s o
f I
np
ut
 S
am
ple
QQ Plot of Sample Data versus Standard Normal
Figure 11: Quantile-quantile plot
an alternative to interpolation, such as the use of GPS data. Having removed the error,
maximum likelihood estimation should be used to fit data to the distribution.
3 Problem 2: vertical safety margins
In problem 1, we dealt with azimuthal error, an angular error in determining the position
of an aircraft. However, NATS are also interested in total vertical error, that is the
difference in height between an aircraft’s assigned altitude and its actual altitude. An
aircraft determines its altitude from its altimeter which works out height from outside
pressure. In bad weather conditions an altimeter may display an inaccurate altitude. The
difference between actual altitude and displayed altitude is known as altimetry system
error.
(a) Height-keeping
performance
(b) Vertical overlap
convolution
Figure 12: Vertical height error and overlap
Total vertical error data has been collected from two Height Monitoring Units
(HMUs) located near Strumble and Gander. From this set of data, the altimetry system
error has been calculated and was provided to the Study Group, to be fitted with a pdf.
Determining the distribution will allow NATS to calculate the vertical overlap probability
Pz. According to the so-called Reich Model, the collision risk Naz is a function of vertical
C-8
overlap probability, horizontal overlap frequency Nr and a weighting kinematic factor K
such that Naz = Pz ×Nr ×K. The height data is believed to be distributed by a Mixed
Gaussian Exponential modelled by
f(d) = (1− α)
(
1
σ
√
2π
exp
(
− d
2
2σ2
))
︸ ︷︷ ︸
Gaussian
+α
(
1
2λ
exp
(
−|d|
λ
))
︸ ︷︷ ︸
Exponential
(4)
for weighting parameter α  1 so that the core is dominated by the Gaussian
distribution. The Study Group was asked to fit the data to this distribution and find
the probability of a tail overlap. The data was supplied as 103 data sets, where each set
i has Ni observations of height errors fi for aircraft type i.
Consider two examples: the Boeing 737-400 series and the Boeing 737-500 series, for
which data is shown in Figure 13. The Boeing 400 series has mean error −56.7 and
1 2 3 4 5 6
x 104
−200
−150
−100
−50
0
50
100
150
200
v
e
rti
ca
l h
ei
gh
t e
rro
r
(a) Boeing 737-400
Series
1 2 3 4 5 6
x 104
−200
−150
−100
−50
0
50
100
150
200
250
v
e
rt
ic
al
 h
ei
gh
t e
rro
r
(b) Boeing 737-500
Series
Figure 13: Boeing 737 ASE data
standard deviation 36.2, and the 500 Series has mean −28.0 with standard deviation
41.2. Observing the data, there is little consistency. Other data sets have positive mean
and smaller standard deviation. Unfortunately the data was given as pure height error
alone for a series of aircraft. A trend could appear if the height error had a time stamp
and was given for a particular aircraft. In Figure 14 we consider the aggregate data and
plot a histogram with log histogram superimposed.
Again we see inconsistent outliers which would make fitting a pdf very difficult. To
determine whether the data fits a mixed Gaussian exponential distribution, we calculate
the expectation, defined by
E(d) = E(X − d|X > d) =
∫∞
u
(s− u)fx(s)ds∫∞
u
fx(s)ds
. (5)
The expectation of the Gaussian core is given by E(d) = σ2/d, a hyperbola in d, and an
exponential tail expectation E(d) = 1/λ would give a straight line as a function of d. In
Figure 15 the expectation of the aggregate data is plotted.
The core appears to be Gaussian, but the tail data is far from constant. It is too
random to fit and should not be used to calculate probability overlap.
C-9
−400 −300 −200 −100 0 100 200 300 400
0
1
2
3
4
5
6
7
x 104
error in feet
hi
st
og
ra
m
−6
−4
−2
0
2
4
6
8
10
12
lo
g 
of
 h
ist
og
ra
m
Figure 14: Aggregate data log histogram and histogram plot
Figure 15: Expectation of aggregate data
3.1 Conclusions and recommendations
NATS use a mixed Gaussian-exponential model of ASE data to calculate the vertical
probability overlap of aircrafts. However, analysis of the data suggests that more data
is required for the outliers. There is currently not enough tail data to fit conclusively to
any pdf. The data would be of better use if it were time stamped, since it would then
elucidate whether the altimeters are being maintained to compensate for ASE; if so, this
would suggest the data is no long iid and cannot be fitted to a pdf.
References
[1] D.J. Higham, N.J. Higham Matlab Guide (SIAM, Philadelphia, 2000)
[2] R.L. Burden, J. Douglas Faires. Numerical Analysis (Brooks/Cole Publishing
Company, USA, 1997)
C-10


National Air Traffic Services

Abstract

Similar works

Full text

Available Versions

Mathematics in Industry Information Service Eprints Archive