Magnetic gradiometry using frequency-domain filtering

Accurate measurements of ambient planetary and interplanetary magnetic fields using spacecraft magnetometers typically require accounting for interfering magnetic fields generated by the flight system (FS). The most common method for removing FS-generated time-variable magnetic fields is narrow-band and low-pass filtering of magnetic field data in the frequency domain. However, if fluctuations in the ambient field contain frequencies overlapping those in the FS field, it can be difficult to construct a filter that will not affect both signals. Here we present an alternate method for removing FS time-variable signatures from magnetic field measurements. For spacecraft that make use of a magnetic gradiometer (i.e. with two or more instruments on a boom at different distances from the center of the spacecraft), the dominant frequencies in the FS field can be identified using spectra of the differenced field components. The amplitudes of the FS field at those frequencies can then be suppressed without removing spectral peaks present in the ambient field. We demonstrate the successful application of this method, referred to as gradiometry peak suppression, both to modeled data sets and to 128 Hz Venus Express magnetometer data.

Accurate measurements of ambient planetary and interplanetary magnetic fields using spacecraft magnetometers typically require accounting for interfering magnetic fields generated by the flight system (FS). The most common method for removing FS-generated time-variable magnetic fields is narrow-band and low-pass filtering of magnetic field data in the frequency domain. However, if fluctuations in the ambient field contain frequencies overlapping those in the FS field, it can be difficult to construct a filter that will not affect both signals. Here we present an alternate method for removing FS time-variable signatures from magnetic field measurements. For spacecraft that make use of a magnetic gradiometer (i.e. with two or more instruments on a boom at different distances from the center of the spacecraft), the dominant frequencies in the FS field can be identified using spectra of the differenced field components. The amplitudes of the FS field at those frequencies can then be suppressed without removing spectral peaks present in the ambient field. We demonstrate the successful application of this method, referred to as gradiometry peak suppression, both to modeled data sets and to 128 Hz Venus Express magnetometer data.
Supplementary material for this article is available online Keywords: gradiometry, spectral analysis, magnetic fields (Some figures may appear in colour only in the online journal) * Author to whom any correspondence should be addressed.
Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

Introduction
Spacecraft magnetic field measurements are fundamental for the study of the solar system magnetism including planetary dynamos, planetary induction fields, crustal remanence, the interplanetary magnetic field (IMF) and space weather. However, a key challenge in making such observations is that the field measured by spacecraft magnetometers includes both the ambient field plus magnetic fields generated by the fight system. In many cases the flight system (FS) fields can be the dominant noise if not removed.
Major sources of FS field contamination are the alternating current (AC) fields (i.e. having frequency, f, typically >0.01 Hz) generated by permanently magnetized and/or highcurrent components such as stepper motors, reaction wheels, communication components, and heaters [1]. This contamination can be reduced by mounting the magnetometer sensors on a boom at the end of which the FS field will be attenuated and also by implementing a magnetic cleanliness program during mission design to limit the magnetic moments of individual FS components. FS AC fields can also be mitigated following data acquisition using a variety of computation techniques. If the spacecraft is spinning, the FS fields and the ambient fields can be distinguished during post-processing on the ground since the ambient field will be modulated at the spin frequency but the FS field will remain stationary. For non-spinning spacecraft, separating the FS and ambient field fluctuations is less straightforward. One common method for removing AC fields from such datasets is to filter them in the frequency domain based on a priori constraints on the frequency-content of the FS field signatures. Other methods include time-averaging measurements taken by multiple sensors and time-averaging measurements taken by a single sensor to eliminate high-frequency FS fields [2][3][4][5][6][7][8].
For spacecraft that carry only one magnetometer sensor, filtering can only be accomplished by using spacecraft ancillary data (e.g. the timing of changes of state of FS components) along with a priori knowledge of the frequencies generated by the different subsystems and components on the spacecraft. Information about the expected FS fields can be provided by ground-based testing prior to launch as well as from in-flight calibration activities. The difficulty in accurately identifying FS fields is that they are often composed of multiple signals with different amplitudes and frequencies originating from multiple components at different distances and directions from the sensor. As a result, it is difficult to isolate contributions from the FS and ambient field at a given magnetometer sensor location. This is particularly challenging when the FS fields have frequencies similar to those present in the ambient field, making it difficult to design a filter that targets only the FS fields.
For spacecraft that carry two or more magnetometer sensors mounted at different distances from the spacecraft bus, gradiometry methods can be used to remove spacecraft generated fields from the measurements. Magnetic gradiometry has been a staple in spacecraft magnetometry investigations since it was first developed by Ness (1971) for Mariner 10 [9]. This method uses simultaneous measurements from two or more sensors on a boom to remove the FS field by measuring the field gradient. In this configuration, each sensor measures a combination of the ambient field, ⃗ B amb , and the FS field at the location of the sensor , ⃗ B FS (⃗ r k ), where ⃗ r k is the position of the k th sensor in the spacecraft coordinate frame. Typically, the distance to the ambient field source is much larger than the separation between the two sensors. As a result, the ambient field can be recognized by the fact that it is the same at the sensors while the FS field is weaker at the more distal sensors. It is important that the gradiometer sensors are time synchronized to a fraction of the time between measurements in order to accurately determine the differenced field (e.g. to a typical accuracy of a few microseconds). If there is a temporal offset in the measurements some corrections can be made in post processing, but the efficacy of the gradiometry methods will be reduced.
Here we consider the specific configuration of two sensors mounted along a single boom, with an outboard sensor (k = 2) at the end of the boom and the inboard sensor (k = 1) closer to the spacecraft bus. Each sensor measures the three orthogonal components of the magnetic field ( j =1, 2, 3). The coupling coefficients are defined as α ij = Bi,out B j,in , where B i,out is the ith component of the observation at the outboard sensor and B j,in is the jth component of the observation at the inboard sensor. The coupling coefficients form a diagonal matrix whose entries for i ̸ = j are 0. In the case where there the field has a single multipole term, all entries along the diagonal are identical. However, for fields with more than one multipole term α 11 ̸ = α 22 ̸ = α 33 . The latter includes cases with multiple additive FS sources, as well as cases where the field sources do not lie along the boom axis. The differenced field is given by From (1) it can be seen that the differenced field contains only the field generated by the spacecraft because the ambient field terms cancel out. The spacecraft gradiometry equations developed by Ness [9] give estimates for the FS field at the inboard and outboard sensors, ⃗ B est FS_in and ⃗ B est FS_out , and the ambient field ( ⃗ B est am ) as where α is the coupling matrix discussed above. These equations can be used to correct for direct current (DC) (i.e. time-stationary) fields as well as step changes (i.e. abrupt changes in the FS DC field, resulting from a component being turned on and off) when these fields are strong compared to the intrinsic instrument noise. However, to accurately identify the coupling coefficients, the FS field magnitude must be several times larger than the instrument noise. When the magnitude of the FS field is comparable to the instrument noise, uncertainty in the coupling coefficients can cause significant errors in the removal of the FS field components [5,9]. This becomes especially apparent when using gradiometry to remove a spacecraft-generated AC field. Calculating the spacecraft field on a point-by-point basis leads to large uncertainties around the nodes of the FS AC field oscillations because the FS field contributions are comparable to the instrument noise at those times. However, the frequency content and the timing of the FS contributions with magnitudes comparable to the instrument noise can be identified in the differenced field even if the coupling coefficients, and the true amplitude of the FS field, cannot be accurately determined.
Here we describe a method for cleaning FS fields from magnetometry data by which specific frequencies, referred to as spectral peaks, are identified in the power spectrum of the differenced field and then suppressed in the power spectra of the observed fields at each sensor. This method, which we call frequency-domain gradiometry filtering (or, more briefly, gradiometry peak suppression) enables the identification and suppression of the FS AC field while leaving the ambient field fluctuations largely intact. This method of separating spectral peaks is used in numerous signal processing applications such as removing or enhancing emissions from specific molecules in nuclear magnetic resonance spectroscopy [10][11][12], removing digital artifacts due to transmission errors from digital images [13], and correcting exposure thresholds for specific regions in medical radiography images [14].
This manuscript is organized as follows. In section 2, we describe how to use frequency-domain gradiometry filtering with a magnetic gradiometer to identify and remove AC signals generated by FS components. In section 3, we describe synthetic datasets used to test the gradiometry peak suppression method. In section 4, we discuss the results of application of gradiometry peak suppression to the synthetic datasets. In section 5, we apply the method to Venus Express (VEX) data and compare the results to published measurements. Finally, in section 6, we discuss the applications of peak suppression gradiometry to existing and future missions and data and also discuss its limitations.

Method
The gradiometry peak suppression method identifies one or more discrete frequencies contaminated by FS fields using spectra of ∆B x , ∆B y , and ∆B z and then suppresses the power at those frequencies in the spectra of the inboard and outboard sensor measurements. As such, the first step is to calculate the differenced field (∆ ⃗ B). This can be accomplished using either the individual components, ∆B x , ∆B y , and ∆B z , or the total differenced field given by |∆B| = |B in | − |B out |. Our preferred method is to use the total differenced field since the amplitude of the fluctuations in the total field is more likely to exceed the predetermined threshold for detection than each of the three field components individually for FS sources that do not lie along the boom axis. Once the differenced field is calculated, the next step is to determine when a FS AC field is above the noise floor, b n , which we take as the 3σ instrument noise (see below). To accomplish this, we identify the upper and lower envelopes of the differenced field. The upper envelope (∆B ue ) and lower envelope (∆B le ) are determined by taking a rolling maximum and minimum of the differenced field. Similar to a running average, the maximum and minumum values are determined for a sliding window of a predetermined length, δt. An example timeseries of ∆B, ∆B ue , and ∆B le is shown in figure 1(A). Accurately identifying the FS magnetic field signatures requires, at a minimum, one full period for each AC field contribution, such that the value of δt depends on the frequency range of the FS generated AC fields.
FS AC signals are identified as simultaneous changes in ∆B ue and ∆B le . Interval identification is determined based on two parameters, δb and n. δb is defined as the absolute value of the change in ∆B ue and ∆B le at the beginning of the envelope, and n is a specified number of time steps over which the change in the envelope must occur. Intervals are selected for fluctuation removal when δb is larger than a predefined threshold, δb min , over a specified number of time steps, n, such that δb/n ⩾ δb min /n.
The field magnitude threshold, δb min , indicates the smallest change in the envelopes that can result in a time interval being selected. If a FS AC field has an amplitude comparable to the instrument noise, the signal can be identified in the power spectra, but if the amplitude is lower than the instrument noise, the frequency of the fluctuations will be too distorted to accurately identify. As a result, δb min = b n is a reasonable initial threshold.
The other adjustable parameter, n, depends on how rapidly the amplitude of AC FS signals change in time (e.g. time required for a component to turn on or off and whether there is a ramp-up period of gradually increasing power) ( figure 1(A)). For most components, the change in field will likely be abrupt, occurring over a small fraction of a second, such that a relatively small number of timesteps can be used to identify the edges of the envelope [7]. The next step is to calculate the fast Fourier transform (FFT) of each component of ⃗ B in , ⃗ B out , and ∆ ⃗ B for the interval identified using envelope changes. Prior to calculating the FFT, we apply a Kaiser (β = 10) windowing function [15] to the AC interval; this minimizes the energy in the fluctuations at the edges of the time interval and avoids introducing a heavyside function in the time series that would introduce artificial fluctuations in the spectral analysis. In addition, following standard FFT techniques, the dataset is mirrored at both ends (e.g. the whole interval is temporally reversed and appended to both the beginning and the end of the original interval) to minimize edge effects associated with the finite time interval.
Peaks in the FS field are identified in the spectrum of the ∆B as frequencies with amplitudes greater than the 98th percentile for the interval figure 1(B). This threshold can be modified based on the spectral profile of the FS field. If the FS field has field contributions at only a few dominant frequencies, then the 98th percentile is a reasonable threshold for identifying the spectral peaks. However, if the FS field has significant contributions over a wide range of frequencies, then a lower threshold can be set. In practice, we have found that the threshold should not be set lower than the 95th percentile in order to avoid selecting the minor peaks that are contributed by instrument noise for removal. Additionally, the lowest frequencies (f < 0.1 Hz) are not included in the peak identification step in order to avoid removing broadband noise at the lower frequencies.
After identifying the frequencies of the spectral peaks in the differenced field, the equivalent field contributions at the selected frequencies in ⃗ B in and ⃗ B out are suppressed by setting their amplitude to a small fraction, typically between 1% and 5%, of the observed amplitude. This value can be adjusted based on the amplitude of the interfering signal with respect to the amplitude of the ambient frequencies. The ultimate goal is to reduce the FS generated spectral peaks to noise levels without generating a delta function in the spectrum.
The final step is to perform an inverse FFT to clean the identified FS AC signature from the time series. The distortion caused by the windowing function is corrected by dividing by the same windowing function that was applied to the time series prior to performing the FFT. This returns reconstructed ambient field time series, ⃗ B est amb for the inboard and outboard observations.

Magnetic field models
We validated the gradiometry peak suppression method by applying it to simulated FS magnetic field datasets. In particular, we compared reconstructed fields derived from these datasets to the simulated ambient fields, ⃗ B actual amb , to determine the error in the gradiometry peak suppression calculations given by ε = ⃗ B actual amb − ⃗ B est amb . We considered two field models. Model 1 (supplementary material Model1.csv available online at stacks.iop.org/MST/ 33/015104/mmedia) was designed to test the application of our method on a field model containing multiple AC signatures from different FS components that overlap in time. Model 2 (supplementary material Model2.csv available online at stacks.iop.org/MST/33/01510/mmedia) was used to determine the lower limit on the frequency difference between FS and ambient field frequencies before the gradiometry peak suppresison method can no longer reliably identify and remove the FS-generated fields. Both simulated datasets were developed for a virtual spacecraft with two sensors placed on a boom at 1.4 m and 2.1 m above the top deck of the spacecraft bus. Simulated observations at each of the two sensors include an ambient field, ⃗ B amb , FS-generated AC fields, ⃗ B FS_in and ⃗ B FS_out , and instrument noise, ⃗ B noise_in and ⃗ B noise_out . The total field at each sensor is therefore: The FS and ambient field amplitudes were modeled based on the typical minimum detectable signals for space-based magnetometers (∼1 nT). Sensor noise was included as random Gaussian noise with a 3σ value of 0.22 nT in each sensor, ⃗ B noise_in and ⃗ B noise_out . The models are described in detail in the following sections. The sensor locations, sensor performance and the FS AC fields were chosen to represent those expected for the the upcoming NASA Psyche mission [16,17] which incorporates dual threeaxis fluxgate magnetometers for measuring the magnetic field of asteroid (16) Psyche.

Model 1
We used Model 1 to test the gradiometry peak suppression method on a realistic FS field containing temporally coincident AC signatures that turn on and off and have varying amplitudes and durations. This 30 min time series included 20 different synthetic AC sources, modeled as pure sine waves, superposed on actual measurements of the IMF from the Juno magnetometer [18,19] (figure 2(A)). Each of the AC sources was assigned a randomly-chosen frequency between 0.1 and 16 Hz and a random location on the spacecraft bus. The start time and duration for each source were also randomly assigned across the 30 min interval. To make the modeled FS field as realistic as possible, no restriction was placed on the onset time for each AC signal. The result is that several AC signals can overlap in a single interval. The total field at the inboard and outboard sensors and the differenced are shown in figure 3. The peak-to-peak amplitude of the AC field at the inboard sensor is <3 nT.
One example of overlapping AC signals is the interval shown in the the shaded box in figure 3 (394.75-453.19 s). The time series and spectra for this interval, which includes several different FS AC sources with frequencies of 0.65 Hz, 1.7 Hz, 6.5 Hz and 7.6 Hz, contribute to the observed field. Details of the total field for each sensor and the differenced field for that interval are shown in figures 4(A)-(C) along with the corresponding amplitude spectra (figures 4(D)-(F)).

Model 2
We used Model 2 to determine the minimum frequency separation (δf) required by the gradiometry peak suppression algorithm and the lower limit on the frequencies that can be removed using this method. These relate to frequency resolution in the FFT, which is determined by the sampling frequency and the duration of the intervals selected for FS field removal. With these goals, Model 2 was constructed to contain an ambient field consisting of a single frequency and a single FS source that steps through a range of frequencies around the ambient field frequency. We incorporated FS field sources with frequencies near the lower end of the frequency range where peak suppression is likely to be a viable method for

Model 1 results
To show the details of the peak suppression method described in section 1, we focused on the aforementiond 394.75-453.19 s interval in Model 1 consisting of four FS AC contributions. As noted in section 3.1, the interval is defined by the changes in ∆B ue and ∆B le . The start and stop times for the interval coincide with an AC source switching on or off. Because there are only four discrete frequencies in the modeled FS field, the dominant peaks in the differenced spectrum were identified as those that have amplitudes above the 98th percentile in the differenced field. In this example, we have used |∆B| to determine a single set of peaks for all three components. To suppress the fluctuations in the time series we reduce the amplitude at the selected frequencies in the spectrum to the sensor noise levels by setting it to 1% of the observed value.
After the differenced peaks were suppressed in the observed signal, an inverse FFT was performed on the corrected spectrum and the distortion caused by the windowing function is corrected. The final result is a time series that does not contain significant fluctuations at the suppressed frequencies. A comparison of the original (i.e. raw) data and cleaned (i.e. estimated) data for Model 1 is shown in figure 8. The top three panels show the original ambient field, ⃗ B actual amb , and the reconstructed ambient field, ⃗ B est amb . The lower three panels show their diffrence, ε. The four FS AC contributions noted in section 3 are identified in |∆B| (red arrows in figure 4(F)). The spectra for the outboard sensor components (figures 7(D)-(F)) show that only the fluctuations at 1.7 Hz and 7.6 Hz contribute significantly to the outboard  sensor. The time series (figures 7(A)-(C)) show fluctuations from multiple frequencies for both the x and z components but the y component only shows a minor contribution from the higher frequency fluctuations. Fluctuations from the higher frequency source would have been left in the reconstructed signal if the FS contributions were identified using the components of the differenced field rather than the total differenced field. This is because the power at those frequencies is below the threshold for selection in ∆B y but not in the other two components.  interval down to the instrument noise levels as described in section 3.

Model 2 results
Results from Model 2 are shown in figure 9. The only packets that were not cleaned effectively are the two intervals with FS field frequencies separated by ⩽0.05 Hz from the ambient field frequency of 1.85 Hz (1.8 Hz and 1.9 Hz).
In general, the frequency resolution of the FFT determines the minimum frequency separation between ambient and FS generated fields required to acurately remove only the FS fields. The resolution of the FFT is given by δf = f s /N s where f s is the sampling frequency and N s is the number of samples in the window. Since N s = dt × f s , where dt is the duration of the window in seconds, δf = 1/dt. For example, to obtain a resolution better than 0.05 Hz, the windows should be a minimum of 20 s. This also sets the lower limit on the frequencies that can be removed using peak suppression for a given window, which is given by f min = 2 × δf. The frequency resolution can be increased by mirroring the beginning and end of the data set by effectively increasing the length of the data set. However, it is possible to artificially enhance the power at some frequencies if only a portion of the data set is included in the mirror. The upper limit on the frequencies that can be identified is set by the measurement Nyquist frequency.

Venus Express mission and data
We next apply our method to actual data from the VEX spacecraft mission, which arrived in orbit around Venus in 2006. The goals of the mission included investigation of the plasma environment including lightning on Venus and atmospheric loss mechanisms related to interactions between the Venusian atmosphere and the solar wind [20]. A magnetic gradiometer was included on the spacecraft to help accomplish these goals. The gradiometer consisted of two fluxgate magnetometers on a 1 m boom. One sensor was mounted at the end of the boom and the other was mounted at the base of the boom, 10 cm from the spacecraft bus [21]. The instrument had three measurement modes with different frequencies: 1, 32 and 128 Hz. Given that the sensors were positioned so close to the spacecraft, the FS fields observed by the magnetometer sensors were ∼200 nT at the outboard sensor [22]. Several methods to remove the FS fields from the observations were developed and the mission was able to report the Venusian fields to within an uncertainty of ∼1 nT [7].
During the closest approach to Venus in each orbit, the magnetometer used the 128 Hz burst mode to study electromagnetic waves thought to be generated by lightning on Venus The spectrum for the differenced field. The red arrows indicate the frequencies with spectral peaks above the 98th percentile. [8,23]. To test the gradiometry peak suppression method, we selected several VEX orbits previously analyzed by [24]. Outboard sensor measurements from one of those days, 9 June 2006 143:18-143:45 UT, are shown in figure 10. The average value has been removed from each component so that only the fluctuations are shown.
The spectrum for the total field at the outboard sensor and the differenced field are shown in figure 11. Based on the differenced field spectrum, it is apparent that the prominent fluctuations at 32-40 Hz as well as peaks 48, 50, 56, and 60-61 Hz are FS-generated. These peaks in the FS field are identified in the differenced field spectrum ( figure 11(B)) by the red arrows. The difference field, shown in figure 12(A), also shows a lot of interference from the reaction wheels throughout the time series. FS-generated fields have amplitudes on the order of 5 nT at the outboard sensor. (B) Differenced field for the same interval as (A) after the flight system field has been removed using gradiometry peak suppression.

Method comparison
The cleaning method used by [24], referred to here as the VEX method, and the gradiometry peak suppression method presented here are very similar in the way they identify the frequencies of the FS fields. Both identify frequencies in the the differenced field spectrum based on whether their power is above a specified threshold. There are two subtle differences in this process. First, the VEX method identifies FS frequencies using the power spectral density (PSD) while the peak suppression method uses the amplitude spectral density (ASD). The effect on the processing is negligible because the differenced field does not contain a significant amount of power in the lower frequencies compared to the FS fluctuations. Second, in the VEX method the spectral content of ⃗ B in , ⃗ B out and ∆B, were compared prior to selecting a frequency for removal. Only the frequency bands with power above a predefined threshold in the differenced field and a coupling coefficient above a predefined threshold are selected for removal. In gradiometry peak suppression, the coupling coefficients are not evaluated. All of the frequencies above the 98th percentile in the differenced field are selected and suppressed in the observed spectra for each sensor.
More importantly, the two methods also differ in how the identified FS fields are removed from the observations. In the VEX method, an inverse FFT based on the frequencies selected for removal was used to create a time series of the FS-generated signal. Coupling coefficients were then used to determine the amplitude of the FS signals at the outboard sensor. The time series for the FS-generated field was shifted to match the phase of the fluctuations in the observed time series and then subtracted from the observed field [24].
In comparison, the gradiometry peak suppression method corrects the field in the frequency domain rather than generating a time series for the FS field and using it to make the correction in the time domain. This removes the complication of matching up the phases of the fluctuations in the modeled spacecraft field to the observed field. The VEX method also makes use of the gradiometry equations [9] to determine the amplitude of the FS generated field. This means that the minimum amplitude threshold for identifying a FS field for the VEX method is limited by the need to accurately identify the coupling coefficients. Gradiometry peak suppression only uses the differenced field to identify the frequency of the FS fields rather than calculating the true amplitude of those field contributions. As a result, there is no need to accurately determine the coupling coefficient so the amplitude threshold for FS signals to be identified is reduced compared to more traditional gradiometry methods.

Results
The PSD of the cleaned measurements using two different processing methods are shown in figure 13. The results for the methods show some differences between 30 and 40 Hz. The original VEX results show residual peaks at 33 and 35 Hz, as well as introducing some noise at 31 Hz that was not present in the measured spectrum. The gradiometry peak suppression method completely removes peaks between 30 Hz and 39 Hz, but leaves some residual power at 40 Hz since the power in that frequency is below the 98th percentile for magnitude of the peaks in the differenced field spectrum. The reconstructed time series are very similar for both methods when considering both the wave forms and the noise floor for the results. The red lines show the difference between the results from the VEX method and those from the peak suppression method. During the first 8 s (01:43:18-01:43:26 UT) the difference in the two results is as large as 0.3 nT. For the rest of the time series the differences are <0.15 nT. The differenced field after the FS field has been removed using peak suppression is shown in figures 13(B) and (D). The amplitude fluctuations in the differenced field are reduced from ∼5.0 nT to ∼0.3 nT after gradiometry peak suppression has been applied ( figure 13(D)). There is an improvement of ∼10 pT (14%) in the RMS noise for the 28-42 Hz frequency band when using gradiometry peak suppression compared to the results in [24] ( figure 13(B)).
Results for the whistler waves analysis for 9 June 2006, 01:43:18-01:43:45 are shown in figure 14. Whistler waves are identified using the ellipticity, whether a signal is right or left hand polarized, and the propagation angle of the fluctuations with respect to the background field [25,26]. The coherence is also used to refine the selection of the whistler waves in the data, where the plane with the best coherence between the wave forms in its respective axes is used to identify the times and frequencies where the coherence is above a predefined threshold. The ellipticity is calculated from the ratio of the minimum and maximum eigenvalues of the cospectrum, given by (Re [ f · f * ]), and the propagation angle is determined by calculating the direction of minimum variance for the wave field, ⃗ k, and comparing it to the direction of the background magnetic field, ⃗ B. The whistler-mode waves are signals with ellipticity >0.5, propagation angle <45º, and coherence >0.3 that span at least 1 s in time and 6 Hz in frequency. Times are identified when the parameters are above the predefined thresholds as shown by the thick black lines at the botom of each panel and the results are manually checked for false positives [27].
The whistler wave identified in the data processed using peak suppression are very similar to those identified in the data processed using the VEX method. The most significant differences are during the interval between 01:43:18 and 01:43:26 where some whistlers were falsely identified in the reaction wheel frequency band (32-40 Hz) in the results from the VEX method but not in the peak suppression results. One of the reasons it is so important to accurately remove the reaction wheel noise is that it can mimic whistler waves. Power in the transverse and compressional modes can be used to identify false positives. The power in the transverse waves should be larger than that in the compressional waves for whistler waves. One example is between 01:43:24 and 01:43:25 where the compressional wave power is larger than the transverse wave power indicating that is was a false positive in the VEX method results. The result is that individual packets of whistler waves are identified in the peak suppression results whereas whistler waves were identified continuously from 01:43:19 to 01:43:34 in the VEX results.
Residual reaction wheel noise shows up in the dynamic spectra in figure 14 as a band between 32 Hz and 40 Hz. The differences in the processing results are most distinctive in the ellipticity (figures 14(B) and (G)) and propagation angle (figures 14(C) and (H)) results. The differences between the two methods are shown in figure 15. The time series in panels a-d are the same as the red lines in figure 13. The largest differences are between 01:43:18 and 01:43:26. The absolute difference in the dynamic spectra from the two methods for the three components that are used to identify the whistler signatures are Based on both the time series results and the whistler wave analysis, the gradiometry peak suppression method was able to effectively remove the reaction wheel noise from the the 128 Hz magnetometer data for 9 June 2006. We were able to reproduce the timing, amplitude, ellipticity, propagation angle and wave power for the fluctuations that have been identified as whistler waves generated by lightning on Venus. Peak suppression results also show two distinct packets of whistler waves in agreement with the fluctuations in the time series shown in figure 13.

Discussion
The primary benefit of the gradiometry peak suppression technique is that frequencies of signals generated by the FS can be specifically targeted for removal instead of filtering all fluctuations within a given frequency band or above/below a given threshold in the sensor data. Additionally, as demonstrated in both the results from Model 1, as well as the results from processing the VEX data, multiple FS generated AC signals can be removed simultaneously from a given interval. The interval shown in figure 7 includes four spectral peaks with frequencies between 0.1 and 8 Hz that were identified in the differenced spectrum and removed from the observed spectrum. This demonstrates that the peak suppression method works well even when we have signals from multiple components combined in the total observed signal. AC signals are removed to an uncertainty <0.4 nT, less than twice the instrument noise included in the model. In addition, there are multiple spacecraft AC frequencies contained in the VEX observations. All frequencies were effectively removed, as shown in the final difference field (figure 13).
The limit on the algorithm's ability to accurately remove a FS field is imposed by the difference between the FS field frequencies and the ambient field frequencies. In general, the frequencies contained in the ambient and FS fields must differ by more than 1/dt since this determines the resolution of the frequencies in the FFT. It is important to make sure that the intervals selected are sufficiently long in duration to meet the required frequency resolution for a given application. The duration of the intervals selected can be controlled by the specifications used for the upper and lower envelopes. The duration for the rolling maximum and minimum will determine the shortest interval that can be selected.
In Model 2, each burst of FS-generated fields has a duration of 60 s, which results in a spectral resolution of ∼0.02 Hz. For any spectral peaks in the FS field, there must not be an ambient field peak that falls with ∼0.04 Hz of that signal, or else both the ambient and spacecraft field contributions will be completely suppressed. Only FS signals with a frequency 1.75 ⩽ f ⩽ 1.95 Hz, given a step size of 0.05 Hz in the frequencies tested, were identified and removed from the time series to within 0.6 nT. During the intervals with FS frequencies of 1.8 Hz and 1.9 Hz, the error in field removal was ∼1.2 nT while the error for the rest of the intervals was <0.3 nT.
Another benefit of this method is that it is very flexible. The only requirements are that the gradiometer sensors be time synchronized and spacially separated enough for the FS AC signals to have different amplitudes at each sensor. Having the sensors further apart is ideal so that the differenced field is large. If the sensors are very close together then the differenced field will be at or near noise levels. There is no requirement for the FS fields to be dipolar at the sensors. The VEX spacecraft had sensors separated by 2 m with one very close to the spacecraft. This results in a very large differenced field. It also results in lot of stray fields being observed at the inboard sensor which were completely absent from the outboard sensor measurements. The peak suppression method performs as well in this case as it does with the modeled data shown in section 3 which were built with the assumption that the sensors are placed further from the spacecraft bus (>1.4 m).
Although the method was demonstrated here only for modeled data and for observations collected by VEX, application of the peak suppression method is not limited to spin stabilized spacecraft, nor is it limited to observations collected by fluxgate magnetometers. The only requirement is that the spacecraft carry a gradiometer consisting of two or more time synchronized instruments so that the differenced field can be used to determine the frequencies of the FS generated fields. Peak suppression could similarly be applied to measurements taken by search-coil or vector helium magnetometers if they are used in a gradiometer configuration. If the method is to be used on a spinning spacecraft the data should be de-spun prior to the application of the peak suppression method. This will eliminate the risk of distorting the spin related fluctuations if there are similar frequencies in the FS AC fields.

Conclusions and future work
Our frequency-domain gradiometry filtering can be used to remove FS magnetic fields from magnetic gradiometer data. Both the minimum frequency that can be removed and the minimum frequency resolution that can be achieved are determined by the frequency resolution of the FFT which is, in turn, determined by the duration of the interval being processed.
Using Model 1, we showed that the method is capable of removing concurrent signals from multiple FS sources at once. Figure 7 shows that the FS fluctuations are significantly reduced by the gradiometry peak suppression method with a resulting amplitude of ∼0.25 nT in all components of the reconstructed ambient field (figures 7(J)-(L)). The resulting error in the calculations for the entire 30 min synthetic dataset compared to the original ambient field measurments used to build the model is <0.6 nT.
Model 2 results indicate that the gradiometry peak suppression method is capable of removing FS generated peaks that are separated from the ambient field frequencies by greater than 0.05 Hz. This is driven by the frequency resolution of the FFTs for the individual intervals where FS fields were identified which is, in turn, driven by the number of measurements in those intervals and the sampling frequency. Lower frequencies may be removed by using longer duration time series.
To validate the method on real data, we reprocessed 128 Hz magnetic field data from the VEX mission. We were able to show that the method is effective at identifying and removing the multiple spectral peaks generated by the reaction wheels on the spacecraft. The resulting spectra and wave forms are consistent with the fluctuations that have been identified as whistler waves generated by lightning in the Venusian atmosphere [8,24,27].
Gradiometry peak suppression is a straightforward method that can autonomously identify and remove several FS generated signals in a single step. A list of discrete frequencies is constructed based on the frequencies that have an amplitude above the 98th percentile in the differenced field and each of those frequenices is removed from the observed spectra. There is no need to identify the true amplitude of the FS signals prior to removing them. If the fluctuations are strong enough to cause an observable peak in the field spectrum, they can be removed. The lower limit on the FS field amplitude that can be removed using this method has been identified as δb = b n using modeled data, where b n is the instrument noise. This is much lower than the limits on standard gradiometry methods which are typically several times the instrument noise.
Several possibilities exist for improving the technique. In cases where the ambient field frequency is too close to the spacecraft field frequency, it may be possible to use gradiometry peak suppression with coupling coefficients in frequency space [9] to determine the amplitude of the FS field at the inboard and outboard sensors. This could be used to determine how much of the signals to suppress. The trade off is that this would require accurately identifying the coupling coefficients which would increase the minimum amplitude required for signal identification. Continuing investigations into the best methods for integrating gradiometry peak suppression with coupling coefficients are underway. There are also cases where the FS fields have power over a wide and continuous range of frequencies. While peak suppression works best for well defined FS AC frequencies, it can select any number of frequencies for removal in a given window. However, the challenge in this situation is that, for a wide band of FS frequencies, ambient field frequencies are likely to be removed as well.

Data availability statement
All data that support the findings of this study are included within the article (and any supplementary files).