MODELLING THE MAGNETIC DISTURBANCES DUE TO ROAD-TRAFFIC

Magnetic disturbances due to the traffic are tentatively modelled assuming that the sources are moving dipoles. The influencing section of the road (“useful” portion) should be modelled in 3D. The parameters of the model (time of closest position to the magnetometer, velocity, including its sign, dipole moment) are fairly accurately estimated. The fit is improved with the incorporation of a small induction effect.


INTRODUCTION
Disturbances due to human activities are a major concern for magnetic observatories and, in the past, scientists were obliged to move many of the observatories away from disturbing sources. This problem is enhanced by the worldwide tendency to extend records to higher frequencies. Specific examples are to be found in recent literature, for instance Clark et al. (1994), Hejda and Horacek (2001), and Okawa et al. (2007). The last example is closely connected to the one discussed here, in so far as it deals with disturbances due to traffic close to the observatory. The example discussed herein is taken from a geophysical station, named Welschbruch, located in the Vosges mountains, France, where the Earth's magnetic is nearly permanently recorded according to current standards. In our example, the disturbance due to the traffic on the road that passes nearby is nearly smoothed out on one minute data but is clearly identified on the one second records. The aim of the study reported herein is the modelling of this disturbance, based upon a fairly accurate survey of the section of interest of the road, using GPS measurements, on one hand, and the assumption that the magnetic field produced by the passing vehicles may be approximated by a moving dipole, on the other hand.

WELSCHBRUCH STATION : 3D ROAD MODEL
Welschbruch geophysical station is located in the Vosges mountains, about 40 km south of Strasbourg, France. It is equipped like a standard magnetic observatory and serves mainly to testing instruments and training observers yearly. The acquisition system was recently upgraded to one second sampling rate with the device built up at EOST, Strasbourg (Fotzé et al., 2007). Figure 1 shows the position of the variometer house (V), the absolute house (A), and of the road. The shape of the road has been sampled with GPS measurements, along the downhill side of the road, the spacing ranging from 5 to 20 m, according to the curvature of the path. The corresponding curvilinear abscissae have been measured by means of a tape line deployed along the same downhill side. The origin is the point (P) closest to the variometer, and the orientation is positive downwards. Two sets of points, one for the downhill and one for the uphill lane, have been derived from the sampled points and interpolated with a smoothing cubic spline. Figure 2 displays the shape of the two lanes in the geographical reference frame centered on the magnetometer location.
Data Science Journal, Volume 10, 30 August 2011 The separation between the two lanes is large enough to give rise to two distinguishable signals on the magnetometer which are about 50m apart (see Figures 3 and 4).

MODELLING A MOVING DIPOLE
The magnetic effect of the vehicles has been approximated by a moving dipole M . The disturbance field ) Z , Y , (X = S p p p p is expressed in the geographical reference frame by: (1) Where : (2) (X,Y,Z) give the position of the dipole equivalent to the moving vehicle. X, Y, Z are functions of the curvilinear abscissa s , and . As mentioned above, X(s) , Y(s) , Z(s) are approximated by smoothing cubic splines depending on the uphill or downhill direction. Figure  . Of course, as expected, the magnitude is slightly larger for a downhill moving vehicle, which is a little bit closer to the magnetometer. If the independent variable is now time instead of s , the velocity of the vehicle appears to be a key parameter for the shape of the disturbance. We have assumed a constant velocity along the trajectory, which is clearly an oversimplification: obviously, the vehicles slow down in the curve upwards and generally accelerate in the straight section. Thus, a distortion of the modelled signal is to be expected. Under the assumption of constant velocity, the relation between time and s is merely: where o t is the time of the position closest to the magnetometer. v is the velocity, with positive (resp. negative) sign for downhill (resp. uphill) movement. On Figure 4, "warm" colors are for uphill movements whereas "cold" colors are for downhill movements. d is the distance to the magnetometer. Note that Data Science Journal, Volume 10, 30 August 2011 IAGA141 p X and p Z disturbances are symetric with respect to 0 = s whereas p Y is antisymetric. The clear distinction in Y p signal between positive and negative speeds allows discrimination among vehicles moving downhill and uphill. Note the abrupt drop towards 0 for high speeds. It turns out that for these speeds, s(t) is outside the limits of the sampled portion of the road. We have merely set the signal to 0 for these ranges. This slight distortion is not a serious drawback for the estimation of the parameters of the model.  Figure 5 shows an example of an isolated, clearly identified disturbance, recorded by the triaxial magnetometer as well as by a scalar magnetometer. The duration of the disturbance is about 20s. The sampling rate of the scalar magnetometer being only 10s, this record will not been used further.

DISTURBANCE MODELLING
) Y , (X v v have been recalculated in the geographical reference frame although the angle between v X sensor and the true North is small (around 0.37°). In order to process the data, we have selected a 2 minute long interval approximately centered on the zero crossing of the Y record. A parabolic time function based upon 10 values at each end of the selected section has been substracted from the raw signal in order to extract the disturbance.

Modelling an isolated disturbance
The parameters of the model are the components of the dipole and the velocity v (positive or negative) at the time 0 t of 0 = s . In a preliminary version, we assumed that the dipole was parallel to the ambient field. It turned out that the theoretical signal, especially the Z component, was significantly larger than the observed one so that we quickly abandonned this assumption. The model being non-linear in v and 0 t , the estimation of the model parameters classically requires an initial estimate where τ is the shift between the two signals (see Figure 6). In order to compute the cross-correlation ) R(τ , we have used an adapted version of the function xcorr from the Signal Processing Matlab Toolbox. In our case, 1 S stands for one of the recorded components, and 2 S is the modelled signal, computed as described above. The estimation of the parameters M , 0 t , and v is refined further by a non-linear least-squares method using the Levenberg-Marquardt algorithm (Marquardt, 1963). The misfit to be minimised is written: where the time interval defined by 1 N and 2 N is adjusted with respect to velocity in order to keep the curvilinear abscissae within the range of actual sampling. Note that the method requires the calculation of the first derivatives of p S with respect to the parameters. In particular, in order to calculate the first derivatives with respect to the velocity, we need to know the derivatives of the trajectory with respect to the curvilinear abscissa. This is a by-product of the cubic spline approximation easily derived using the spspline Matlab function of the Spline Toolbox. Figure 7 left shows the fit of the model to the observed signal. We may note, as expected, the spurious effect of the inadequate sampling of the trajectory and may suspect that,  We have tried to improve the model by speculating the slight influence of an induction term due to the dipole moving at the surface of a conductive earth. Indeed, the time variation of g S is given by: where A is the matrix explicited in Eq. (2). Neglecting the possible shift between primary and secondary signal, we have empirically incorporated the induction effect via a symmetric matrix α , which allows a small correlation between the three components of the primary disturbance field. p S now written: where I is the 3x3 identity matrix. The elements of the matrix α are 6 parameters more to estimate.
However, the solution is no longer unique with this additional refinement. We have circumvented the difficulty by the introduction of a minimum induction energy constraint. The misfit function takes now the form:     Another small part of the residual disturbance probably arises from the non-constant actual velocity. Unfortunately, there is no mean with the available data to measure the acceleration of the disturbing vehicles. Finally, we may notice that the residuals have roughly the same shape as the disturbance itself however with a slight shift, which might be accounted for in a more elaborated model of the induction part.

Example with two successive disturbances
The following example is provided by the downwards passage of a bus followed some tens of seconds later by the upwards passage of a lorry. On Figure 9 we may notice that, in neither of the disturbances, the signal is as expected with a dipole roughly oriented in the direction of the ambient field. A several hours watch of the traffic has revealed that only buses and lorries give rise to a measurable signal. Cars and motorbikes are not detected.
Data Science Journal, Volume 10, 30 August 2011 Figure 9. Disturbance produced by a bus, followed some tens of seconds later by a lorry. Figure 10 shows that both disturbances are clearly identified by cross-correlation. However, the correlation is less accurate for the lorry due to the actual orientation of its equivalent dipole. The cross-correlations have been calculated with a dipole estimate based upon the disturbance due to the bus, which is obviously not appropriate for the disturbance due to the lorry. Further parameter estimations are restricted to the disturbance caused by the bus.
Data Science Journal, Volume 10, 30 August 2011   Figure  11 that the R.M.S. (provided by the model including some induction effect) is even better than in the previous example. This is partly due to the disturbance being weaker. On the other hand, the residuals displayed on Figure 11 right are very small, and in this example we are probably close to the limit of the disturbance removal that can be performed with this simple model.

CONCLUSION
It appears in the examples analysed in this work that one second records are a powerful tool for assessing the quality of a magnetic observatory. Indeed, in this case, whereas one minute records do not reveal any spurious influence, one second records show that the magnetic environment is far from perfect. In the present case, fortunately to some extent, the numerous disturbances observed on the one second records may be reasonably attributed to heavy vehicles, such as buses and lorries, passing on the road nearby. These vehicles are assumed to be disturbing field sources due to their magnetization. A simple dipolar approximation seems to work well and be appropriate to identify the main properties of the sources, such as the dipolar moment, the velocity, and direction of movement. There are however two complications: 1. the dipole is not parallel to the ambiant field; 2. the model fits better the recorded disturbance if a small induction contribution is incorporated.
The plausibility of this latter assumption would need further investigations in order to estimate its order of magnitude based upon physical assumptions and, in particular, on some knowledge of the ground conductivity. The automatic removal of the disturbances is another issue which has been left aside, all the more as the station has not the status of a registered observatory. However, regarding this specific disturbance, an automatic processing might be easily imagined, at least in the studies where one second data would be necessary.