Neural Comput & Applic (1993)1:207-214 O 1993 Springer-Verlag London Limited
Neural Computing & Applications
Aeromagnetic Compensation using Neural Networks Peter M. Williams School of Cognitive and Computing Sciences, University of Sussex, Falmer, Brighton BN1 9QH, UK
Airborne magnetic surveys in geophysical exploration can be subject to interference effects from the aircraft. Principal sources are the permanent magnetism of various parts of the aircraft, induction effects created by the earth's magnetic field and eddy-current fields produced by the aircraft's manoeuvres. Neural networks can model these effects as functions of roll, pitch, heading and their time derivatives, together with vertical acceleration, charging currents to the generator, etc., without assuming an explicit physical model. Separation of interference effects from background regional and diurnal fields can also be achieved in a satisfactory way.
normally located on a boom at the rear of the aircraft or in a pod on the wingtip. The magnetometer is then close to major components of the aircraft and interference effects cannot be ignored.
Keywords: Neural networks; Airborne compensation; Regularization; Network pruning
The first component f is a function of geographic coordinates x, y and altitude z, and represents the intensity of the earth's magnetic field at the point (x, y, z). This is the quantity of interest in an actual survey, with both g and h being nuisance terms, g is a function of time and represents diurnal variation in the earth's field. This can vary significantly during a survey flight even on a magnetically quiet day. Diurnal variations are normally corrected by establishing a nearby ground station at which the earth's field is recorded simultaneously. These data are then subtracted from the airborne data, either directly or after smoothing. The third component h is the interference field generated by the aircraft and its manoeuvres. Sources of interference commonly include:
1. Introduction Magnetic prospecting can be used to explore for oil, minerals and even archaeological artifacts. It is effective in searching for targets that are not themselves magnetic, but which are associated with other minerals that have magnetic effects. Airborne surveys over land give more rapid ground coverage than ground surveys and can be performed using either fixed-wing or rotary-wing (helicopter) aircraft. For helicopter installations the magnetometer is usually housed in an aerodynamic bird attached to the aircraft by a cable, well away from any magnetic fields due to the aircraft. For fixed-wing installations the magnetometer is
Original manuscript received 4 June 1993 Correspondence and offprint requests to: P.M. Williams, School
of Cognitive and Computing Sciences, University of Sussex, Falmer, Brighton BN1 9QH, UK.
1.1. Basic Model The total field measured by the airborne magnetometer can be modelled as the sum of three components measured = local + diurnal + interference
= f(x, y, z) + g(t) + h(cb, O, 0 , . . . ) .
9 permanent magnetism of various parts of the aircraft 9 induction effects created by the earth's magnetic field 9 eddy-current magnetic fields produced by the aircraft's manoeuvres. Current carrying wires and movement of control surfaces can also generate magnetic effects.
208
P.M. Williams
1.2. Compensation One approach to the interference problem is to attempt passive compensation by inserting magnets at various places in the aircraft. However, this must be done by trial and error which is time-consuming and assumes a permanent installation. Furthermore, it does not compensate for errors due to the motion of the aircraft. Another approach is to build a system of equations that model these effects analytically [1, 2]. The usual model is a linear function of the products of direction cosines of the earth's total magnetic field vector relative to axes fixed in the aircraft, and their time derivatives. The simplest form uses a system of equations containing 18 unknown coefficients. The interference component is then modelled as a sum of the three terms Hperm + Hind + //eddy with: Hperm = al cos X + a2 cos Y + a 3 cos Z Hin d =
a4Hecos2 X + asHe cos X cos Y + a 6 n e cos S cos Z + a7 He cos 2 Y + a s H e cos Y cos Z + a9H e cos2 Z
/-/eddy =
aloHe cos Xc6s X + an He cos XCOS Y + a12He cos Xc6s Z +
a13H e cos
Ycos S
+
aaaHe cos YcOs Y + alsneCoS YcOs Z + a16He cos Z c6s X + a17He cos Z cos Y + a18He cos Z cOs Z
where cos X, cos Y, cos Z are the direction cosines of the earth's magnetic field along each of the aircraft's major axes and He is the intensity of the field. To determine the coefficients a a , . . . , a~8, the pilot flies a 'figure of merit' in a region of low magnetic gradient to yield measurements for various values of roll, pitch and heading. A high pass filter is then applied to the measured field to remove any residual effects from the earth's field, leaving only variations due to aircraft manoeuvres to be corrected for. The resulting equations are commonly solved using a form of ridge-regression [2]. However, this system of equations is normally poorly conditioned and provides only a first order approximation for some of the effects. Furthermore, some other effects are not easy to calculate unless the exact geometry of the aircraft is known. A neural network model, on the other hand, can provide a non-linear model of the phenomenon incorporating effects such as charging currents to the generators, movement of control cables, attitude
of control surfaces, vertical acceleration etc. without the need for an explicit theoretical model. In addition, it can simultaneously model the earth's regional field without the need for high-pass filtering.
2. Neural Model Resolution of the measured field into the three components discussed in Sect. 1.1 can be represented by a neural network of the form shown in Fig. 1. The single output unit is linear so that the total output is a weighted sum of the outputs from the various hidden units. There is only partial connectivity between input units and hidden units so that the network is made up of three parts that are independent up to the final summation. The first takes inputs from position x, y, z , the second from time t and the third from remaining variables such as attitude, etc. It should be emphasized that this structure is used with data obtained from the figure of merit. When solving the compensation problem, the aim is principally to determine the third component of the network. This is then detached to be used in compensating real data obtained in the survey flight, using variables sensed in flight as input. The first two components are discarded, since they relate to the particular region in which the figure of merit was flown and the diurnal variation experienced on that particular day. If it were possible to fly in a region of constant field, the first component would be superfluous. This is not possible in practice since there are non-zero gradients in the earth's regional field, even at altitude, which must either be filtered or modelled. The figure of merit used in this example is shown in Fig. 2.1 This was flown at approximately 3500
m a g n e t i c field o u t p u t = f + g + h
I
hidden /
position
l
hidden /
I ~ hidden /
time
| attitude ]
Fig. 1. Neural network model of the three components of the measured field. 1 Data supplied by US Geological Survey.
Aeromagnetic Compensation using Neural Networks
209
41
as a function of time alone. An important feature of the application was the physically plausible way in which the separation was made.
2.1. Input Representation
~'01
The input variables to the network were represented as follows:
Position: Coordinates x, y, z were rescaled to the interval [-1, 1]. In addition to x and y inputs, the network was also given the product terms x 2, x y and y2, which allow hidden units to define quadratic regions as well as straight lines. The intention was to allow horizontal gradients to be modelled using fewer hidden units. Since theory implies that the earth's regional field should decrease monotonically with altitude, it was not considered necessary to use product terms involving z. This gave six position variables x, y, z, X 2, xy, y2 .
2351, /537 850
276
Time: The time coordinate t was similarly rescaled to the interval [ - 1 , 1] with t = - 1 corresponding to the beginning of the flight and t = 1 corresponding to the end. Again t2 was provided to give two time variables t, t 2 .
Fig. 2. Figure of merit.
feet above ground level and lasted some 15 minutes. 1700 readings were recorded at half-second intervals. Apart from the wing-tip magnetometer readings, the only sensed variables were roll, pitch, heading and position. Time derivatives of roll, pitch and heading were calculated after the flight. The network was trained by supplying values of variables relating to position, time and attitude as inputs and the field measurements as target outputs. The aim was to train the network to learn the functional relation between aircraft manoeuvres, defined by roll, pitch, heading and their time derivatives, and the interference field these generate. Rather than high-pass filter the data to remove the earth's regional field, the network was allowed to build a model of this field using the first component of the network. It was also allowed to build a model of the diurnal variation using the second component (no base station was established for this particular flight). There is no guarantee that the network will separate these three components in an appropriate way. It could choose, for instance, to model the whole observed variation over the flight
Attitude: Measurements of roll r, pitch p, and heading h were made during flight. These jointly define the instantaneous attitude of the aircraft. Since it can be assumed on physical grounds that aircraft interference depends on these angles only through their sines and cosines, these functions and their time derivatives were used as inputs giving 12 attitude variables in all: sin r, cos r, i" sin r, ~ cos r sinp, c o s p , p s i n p , p cosp sin h, cos h, h sin h, h cos h Sines and cosines automatically scale to [-1, 1]. The trigonometric representation also encodes the equivalence 0~ -= 360~ Time derivatives were rescaled to the same interval by choosing suitable measures of angle and time.
2.2. Output Representation The target training outputs were the measured values of the total field, including diurnal and interference components. These were rescaled to
P.M. Williams
210
lie between - 1 and +1 by using a base level of 54795 gammas and a unit of 25 gammas (1 gamma = 10 -5 Oersted). Note that the absolute DC level is unimportant both for solving the compensation problem and for later exploration purposes. i
if
/ 1
11
/ ~
3162 1
3. Training Networks were trained using optimization and regularization algorithms of Williams [3, 4]. The objective function to be minimised is then N Q = ~- log EL, + Wlog ER
(1)
with ED = 2
(Yi -
t~)2
(2)
i=1
where Yi and ti are corresponding actual and target outputs, N is the number of training items and
Fig. 3. Initial architecture.
sating real surveys, a reduced architecture was therefore used in which there were only 12 + 12 + 12 hidden units as shown in Fig. 4. In this case only 15% of the network's weights were set to zero and 11 + 10 + 10 hidden units were left alive. Despite this difference in the number of active hidden units, the outputs of the two networks were within 0.01 gamma of each other over the whole training set and hence well within the margin of error of the data.
w
ER =
j=l
[w;f
(3)
where W is the number of weights in the network (excluding biases).; The regularizer (3) has a strong tendency to set small weights to zero, principally because of the form of the derivative of Ix[ at the origin. This has the effect of pruning the network of unnecessary connections.
3.1. Architecture Both three and four layer architectures were tried. For a given number of free parameters, no significant improvement was found using 4-layer networks over 3-layer networks and results are reported only for 3-layer architectures. The hyperbolic tangent was used as the transfer function for hidden units. Activation of the output unit was a direct linear sum of weighted inputs and bias. The architecture initially chosen is shown in Fig. 3. This has 18 hidden units for each of the position and time components and 36 for interference. As a result of the pruning effect of the regulariser, 78% of the network's weights were set to zero. Of the 18 + 18 + 36 hidden units available to the network, only 10 +12 +9 were left alive in the sense of having some non-zero output weight. With a view to increased efficiency when compen-
4. Results The total field measured during the flight and the field fitted by the network are shown superimposed in Fig. 5. Vertical scales are in gammas and horizontal scales in seconds from the beginning of the flight. A conventional base level has been chosen for the magnetic field, since only variations are significant. Details are shown in Figs. 6 and 7. The neural network calculated output is the smoother of the two traces in both cases. Fig. 7 shows that the smoothing effect of the neural model can be up to 1 gamma in noisy regions. Fig. 6, however, shows that fluctuations of this magnitude are not always treated as noise. The RMS error over the whole training set was approximately 0.35 gammas. This is roughly the expected error for this type of installation.
9 2 Compare Buntine and Weigend [5] and Mackay [6].
Fig. 4. Pruned architecture.
Aeromagnetic Compensation using Neural Networks
211
30
20
2 10
1 0 100
200
300
400
500
600
700
640
800
650
Fig. 5. Graph showing measured and calculated outputs in gammas for training flight.
660
670
680
Fig. 7. Detail of Fig. 5.
i
1I
J
5oi 40
30
20
10
0
40
50
60
70
80
Fig. 6. Detail of Fig. 5.
0
1O0
200
300
400
500
600
700
800
Fig. 8. Graphs showing separation by the network of the calculated part of Fig. 5 into components of interference (below), regional (centre) and diurnal (above).
4.1. Regional-Residual Separation Separation of the network output into its three components is shown in Fig. 8. There is no significance in the absolute DC levels shown. Suitable levels have been chosen to separate the three curves for display. In fact there is no nonarbitrary way of attaching a separate DC level to the three outputs since the bias on the output unit is shared between the three components. This is no disadvantage in this problem, however, where only changes in magnetic field are of interest. The most important feature of Fig. 8 is that the high frequency components have been placed in
the aircraft interference effect, where they belong, rather than in the regional field. The regional field at this altitude is necessarily smooth and of slowly changing gradient. The oscillations in the calculated regional field along the line of flight correspond to successive north-south traverses indicated in Fig. 2. In fact, it is simple to calculate the neural network model of the regional field over the whole area covered by the flight, not just along the flight path. This is shown in Fig. 9 which was obtained by calculating the outputs of the regional component of the network over all pairs of values of x and y
P.M. Williams
212
OorStb~4 ~
~ ~
I
__ -150
__ 100
150
,ng
Fig. 9. Calculated regional field at central altitude. Fig. 11. Simulated horizontal turn with m i n i m u m at 23 ~ west and m a x i m u m at 147 ~ east.
over a regular grid in the interval [ - 1 , 1]. The plot shown relates to the central altitude, z = 0. Corresponding plots for maximum and minimum altitudes, z = ---l, are similar. 3 The upper curve in Fig. 8 shows the calculated diurnal variation. This is displayed using a larger vertical scale in Fig. 10. Both the total variation of roughly 1.5 gammas and the general shape are physically realistic. 4 The lower curve in Fig. 8 shows the calculated aircraft interference. The major interference effect is clearly the heading error. The three upper levels correspond to north-south traverses and the three lower levels to south-north traverses. It is again possible to calculate the heading error for all values of the heading, not just those experienced during the test flight. This is indeed the purpose of compensation. The predicted error for a general
heading is shown in Fig. 11. This was calculated using constant position and time and assuming zero values for roll, pitch and all time derivatives. Again the DC level has no significance. The important feature is the overall 16 gamma variation. 5 The corresponding predictions for interference due to roll and pitch manoeuvres are shown in Figs. 12 and 13. These have been plotted only for ranges of variables that were experienced during the flight and which were intended to cover the ranges likely to be experienced in survey conditions. It is clear that the major effect comes from the heading error. This is well-known from actual surveys in which uncompensated data exhibits corrugations corresponding to parallel lines flown alternately in opposite directions. An example is given in the upper part of Fig. 14. This shows a small part of a real survey. The lower figure shows the
1.5 1.0 0
0.5 0.0[
200
400
600
10
20
30
800
Fig. 10. Calculated diurnal variation.
3 T h e average variation over the 75 metres sampled vertically was approximately 4 g a m m a s , giving a vertical gradient of approximately - 0 . 0 5 g a m m a / m e t r e . 4 Diurnal variations typically have a period of about a day and an amplitude of around 25 g a m m a s . See Breiner [7] or Dobrin and Savit [8, Ch. 16] for examples.
Fig. 12. Simulated roll manoeuvre.
5 The shape of the curve suggests that the aircraft had already been passively compensated for induction effects, which otherwise would have given an additional 1r periodicity.
Aeromagnetic Compensation using Neural Networks
213
pitch -20
-10
Fig. 13. Simulated pitch manoeuvre.
result of removing the interference effects using the methods described in this paper.
5. Discussion Figures 11-13 reflect characteristics of the particular aircraft used. Other aircraft would exhibit different characteristics. In the same way, another installation in the same aircraft would give rise to different interference effects. In that case another figure of merit needs to be flown. In fact, it is wise to recompensate the aircraft before each survey in case there have been changes in its magnetic characteristics. In addition, the use of roll, pitch and heading as inputs means that the local inclination and declination of the earth's magnetic field are implicitly encoded in the network's weights. This limitation can be overcome if a lower resolution vector (three component fluxgate) magnetometer is also installed. This provides a direct measurement of the three direction cosines of the earth's magnetic field, relative to axes fixed in the aircraft, which can be used together with time derivatives as inputs to the network. This has also been tried with the neural network approach and found effective. The present results relate to an installation with a single high resolution magnetometer. If gradients are required several magnetometers can be installed in various parts of the aircraft, e.g. wingtips, nose, tail. The network shown in Fig. 1 then needs two or four outputs.
6. Conclusion The airborne compensation problem appears to be well-conditioned from a neural network point of
Fig. 14. Magnetic data: raw (above) and compensated (below).
view. Several different architectures were used, including four-layer networks, and different scalings were applied to the input variables. In all cases the actual results were practically indistinguishable. This means that there is a consistent and robust neural network solution to the compensation problem which provides a viable alternative to the classical solution.
Acknowledgements. I am much indebted to Misac Nabighian, Bruno Nilsson, Colin Barnett and Perry Eaton of Newmont Exploration for helpful discussions on the subject of this paper.
P.M. Williams
214
References 1. Leliak P. Identification and evaluation of magneticfield sources of magnetic airborne detector equipped aircraft. 1RE Trans Aerospace Navigational Elect. 1961, 8:95-105 2. Leach BW. Aeromagnetic compensation as a linear regression problem. In: Information Linkage between Applied Mathematics and Industry H. London: Academic Press, 1980 3. Williams PM. A Marquardt algorithm for choosing the step-size in backpropagation learning with conjugate gradients. Cognitive Science Research Paper CSRP 229, University of Sussex, February 1991
4. Williams PM. Improved generalization and network pruning using adaptive Laplace regularization. In: Proceedings 3rd lEE International Conference on Artificial Neural Networks. London: IEE, 1993, 76-80 5. Buntine WL, Weigend AS, Bayesian back-propagation. Complex Systems 1991; 5:603-643 6. MacKay DJC. A practical Bayesian framework for backprop networks. Neural Computation 1992; 4(3):448-472 7. Breiner S. Application Manual for Portable Magnetometers. E G & G Geometrics, 1973 8. Dobrin MB, Savit CH. Introduction to Geophysical Prospecting. New York: McGraw-Hill, 4th ed., 1988
Conference Announcements NEURO-NiMES '93: Sixth International Conference on Neural Networks and their Industrial and Cognitive Applications
25-29 October 1993 Nimes, France
Contact: EC2, 269 Rue de la Garenne, 92024 Nanterre Cedex, France. Artificial Neural Networks in Engineering (ANNIE'93)
14th-17th November 1993 St Louis, Missouri, USA 9 Architectures, pattern recognition, neuro-control, neuroengineering systems.
Contact: ANNIE'93,223 Engineering Management Building, University of Missouri-Rolla, Rolla, Missouri 65401-0249, USA. Second International Workshop on Neural Networks and Genetic Algorithms Applied to SAR and QSAR Studies
8-9 December 1993 Lyon, France 9 Structure-activity and structure-property correlations;
9 9 9 9 9 9
Molecular modelling; Protein classification; Pattern recognition; Comparison of methods; Detection of outliers; Cross-validation and bootstrapping.
Contact: Dr J. Devillers, CTIS, 21 rue de la Banniere, 669003 Lyon, France. Fax: +33 78 62 99 12.
Neural Computing Applications Forum: Quarterly Workshop and Seminar
12-13 January 1994 Kings College, London, UK 9 Neural network architectures; 9 Networks behaviour in practical applications; 9 Medical, engineering, financial, scientific, and industrial applications of neural networks - demand forecasting, pattern recognition, etc. This organisation offers help with putting neural networks into practical applications. Members come from industry and academia, and provides a fertile ground for the development of new ideas and projects. Contact: Karen Edwards, Department of Computing & Applied Mathematics, Aston University, Aston Triangle, Birmingham B4 7ET, UK. Tel: +44 (0)21 359 3611 ext. 4243; Fax: +44 (0)21 333 6215. NIPS*93: Neural Information Processing Systems: Natural and Synthetic
29 November-2 December 1993 Post-Conference Workshops
3rd and 4th December 1993 Denver, CO, USA Neuroscience, theory, implementation and simulation, algorithms and architectures, cognitive science and AI, visual processing control, navigation and planning, applications. Contact: NIPS*93 Registration, NIPS Foundation, PO Box 6003, Pasadena CA 91116-6035, USA.
If you have any problems in making contact with any of the above, please inform the Editor who will help you obtain details.
