Magnetometer Array System For Magnetically Noisy Environments

Abstract

A system and method for separating spacecraft generated magnetic noise from geomagnetic field data includes detecting magnetic information using a plurality of magnetometers and outputting a representative detection signal, the magnetic information comprising magnetic noise and geomagnetic field data; receiving the representative detection signal including the magnetic information from the plurality of magnetometers and resolving the magnetic noise data and the geomagnetic field data by analyzing a relative gain and phase of the representative detection signal in a time-frequency (TF) domain, if a noise signal of the representative detection signal is sparse in the TF domain then the relative gain and phase is found using cluster analysis, and separating a geomagnetic field data signal from a magnetic noise data signal using cluster centroids in compressive sensing based on the cluster analysis. The system and method further includes outputting a resultant signal representative of at least the geomagnetic field data signal.

Description

FIELD

The present disclosure relates to magnetometers and, more particularly, relates to a magnetometer array system for magnetically noisy environments.

BACKGROUND & SUMMARY

This section provides background information related to the present disclosure which is not necessarily prior art. This section provides a general summary of the disclosure, and is not a comprehensive disclosure of its full scope or all of its features.

Spacecraft equipped with magnetometers can be used to capture in situ measurements of magnetic phenomena in the geospace environment. These measurements are necessary to answer key questions about the nature of the Earth's magnetosphere and its interaction with interplanetary magnetic fields. Understanding how the heliosphere directs the flow of energy, mass, and momentum between the Sun and Earth is critical for applications such as space weather modeling, space exploration, and climate science. A number of missions use spacecraft equipped with magnetometers to measure magnetic fields. For example, The European Space Agency's SWARM mission uses a constellation of three satellites to provide high fidelity magnetic field measurements used to model the Earth's magnetic field and study the Earth's. Magnetometers provide invaluable data for space science research, however, the quality of the data are often limited by magnetic noise generated by the spacecraft. Electrical systems onboard a spacecraft generate stray magnetic fields that interfere with magnetic field measurements. The strength of magnetic fields in the geospace environment ranges several orders of magnitude with natural phenomena such as the interplanetary magnetic field occurring on the order of 6 nT to the Earth's magnetosphere in low-Earth orbit measuring on the order of 60,000 nT. Spacecraft subsystem magnetic fields may completely eclipse the perturbations in natural magnetic fields which are of interest to understanding waves and currents in the solar wind and magnetosphere. The presence of these stray magnetic fields is a significant obstacle for missions that utilize magnetic field data.

On satellites, stray magnetic fields can be generated by subsystems such as solar panels, reaction wheels, battery currents, and magnetorquers. The magnetometer on the CubeSat, Ex-Alta 1, recorded magnetic field noise generated by a magnetorquer which exceeded 7,500 nT peak-to-peak. Satellite magnetometers are typically fixed at the end of a mechanical boom to reduce the magnitude of noise generated by the spacecraft. For example, the mission SWARM uses two magnetometers mounted on a 4.3 m boom. However, the use of a boom is not always possible in designs such as rovers and CubeSats where gravity and cost are limiting factors. Additionally, the implementation of a boom does not always guarantee the elimination of spacecraft noise from magnetic field measurements. The spacecraft, DMSP, employs a single magnetometer on the end of a 5 m boom, but still faces issues with spacecraft noise.

The use of a single magnetometer on a spacecraft requires a careful magnetic cleanliness design and characterization of the spacecraft's magnetic signature in order to minimize stray magnetic fields. In the case of the spacecraft Cassiope, a software update changed the behavior of the spacecraft's fluxgate magnetometer (MGF). Special spacecraft maneuvers to decrease the spacecraft's noise signature were required in order to recalibrate the MGF. Methods to autonomously identify spacecraft noise would allow Cassiope to do in situ MGF calibration without special spacecraft maneuvers.

In spacecraft with multiple magnetometers, the traditional way to cancel stray magnetic field noise is to perform gradiometry. Gradiometry is a technique which compares magnetometer signals from two spatially separated sensors and calculates the gradient between them. Ness et al. uses the gradient to fit a dipole to the spacecraft noise and formulate a coupling matrix. The coupling matrix is then used to subtract the spacecraft noise from the magnetometer measurements. This method can also be applied to higher order magnetic fields but requires arduous pre-flight characterization of the spacecraft's magnetic signature. Ream et al. uses gradients in the frequency domain to identify and suppress spacecraft noise. However, this method assumes that the spectra of the ambient magnetic field and the spacecraft noise do not overlap. Pope et al. uses the axial gradients and fuzzy logic to identify spacecraft noise according to the subsystem that generates it. The identified noise segments are then corrected in the time domain using information about the noise generated by the subsystem. This method is successful at identifying and removing many different individual noise sources, however, it is not designed to correct multiple concurrent noise sources.

Other noise cancellation methods employ state estimation of the magnetic fields generated by spacecraft subsystems by examining spacecraft housekeeping data. Deshmukh et al. uses a supervised machine learning algorithm to estimate the transfer function of housekeeping currents to stray magnetic fields. Total knowledge of a spacecraft's magnetic signature would allow for perfect interference cancellation; however, housekeeping telemetry provides an incomplete mapping of a spacecraft's current distribution. Additionally, housekeeping data are often sampled at a low cadence and may not have the appropriate bandwidth to identify higher frequency noise. For low-cost applications with a large number of spacecraft, such as CubeSat constellations, it is advantageous to use a method that does not require a boom, rely on prior knowledge of the spacecraft's magnetic signature, or requires human analysis.

Recent progress has been made in magnetometer noise cancelation through the application of blind source separation (BSS) algorithms. BSS is the separation of a mixture of source signals without prior knowledge of the signal type or magnetometer location. Constantinescu et al. use maximum variance analysis (MVA) to clean spacecraft magnetometer data. The MVA method finds an orthogonal set of axes to maximize the variance of the measured signals. These axes represent the principal components which are used to identify and remove noise sources. This application of MVA requires that the variance in the noise sources is larger than the variance in the background magnetic field, and can only identify a limited number of noise signals. Imajo et al. proposed the use of independent component analysis (ICA) to separate geomagnetic field data, captured by the satellite Michibiki-1, from stray magnetic field noise. This method separates signals based on statistical independence, and works well when the number of noise sources are not more than the number of magnetometers. The MVA and ICA methods both separate signals through optimizing statistical quantities; however, they are limited by the number of noise signals they can identify. Sheinker and Moldwin proposed a novel BSS method that uses an analytical formulation to estimate the gain of a single noise source between magnetometers. This method is designed for the case in which a single noise source is present, and does not account for the presence of multiple noise sources. Although, the method may be adapted to remove multiple noise sources by adding more magnetometers.

In accordance with the principles of the present teachings, application of an underdetermined blind source separation (UBSS) method based on the unsupervised machine learning algorithm, Density Based Spatial Clustering of Applications with Noise (DBSCAN), and compressive sensing to separate the ambient magnetic field from spacecraft noise is provided. UBSS is a class of problems in which there are M sensors and N unknown source signals such that M<N. In UBSS, no prior knowledge of the source signals is assumed and the number of source signals that can be separated is not limited by the number of sensors. The system used to model UBSS is defined by the following relationship.

$\begin{array}{cc}B\left(k\right)=\mathrm{KS}\left(k\right)& \left(1\right)\end{array}$

UBSS is a topic that has been thoroughly researched in other fields, such as acoustics and radar signal processing. In the field of acoustics, this problem is famously referred to as the cocktail party problem. In the cocktail party problem, there is a room full of people each having conversations. An array of microphones is placed in the room to record the concurrent conversations. The microphone recordings are then used to separate each individual voice. Guo et al. demonstrate the identification of four human voices using three microphones. He et al. also demonstrate the separation of six flutes recorded by three microphones using the DBSCAN method.

Due to the spatial structure of magnetic fields, the same methods developed to solve the cocktail party problem cannot be directly applied to magnetic noise cancelation. When considering a dipole noise source, the vector magnetic field will have a different magnitude and polarity depending on the magnetic latitude and radial distance of the magnetometer. In the present disclosure, we model the spatial structure of magnetic fields with a phase, although magnetic noise signals mix instantaneously. The structure of the magnetic noise signal is not always dipolar and will change depending on the geometry of the noise source. In magnetic UBSS, the mixing matrix, K, is a complex matrix representing the gain and phase of each signal at each magnetometer. In radar signal processing, Bai et al. apply a similar approach by using complex mixing matrices to model time-delayed radar signals with different directions of arrival. According to the present disclosure, DBSCAN is used to estimate the mixing matrix, K. Once K is known, compressive sensing is used to restore the geomagnetic field signal from the noisy magnetometer data.

Two experiments are presented herein to validate this method. The first experiment separates four computer-simulated noise signals from an ambient magnetic field signal. The second experiment separates the same ambient magnetic field signal using real magnetic field data recorded using an experimental CubeSat apparatus with copper coil-generated signals and three PNI RM3100 magnetometers. The aim of this work is to develop a robust signal processing method to remove spacecraft noise and minimize the need for a mechanical boom or a magnetically clean spacecraft. This work focuses on developing a noise cancellation method for geomagnetic field data, but can also be applied to remove noise in measurements of planetary magnetospheres and interplanetary magnetic fields.

According to the principles of the present teachings, in spacecraft equipped with multiple magnetometers, signal processing systems and methods can be used to compare magnetometer measurements and remove stray magnetic noise signals. The present teachings use density-based cluster analysis to identify spacecraft noise signals and compressive sensing to separate spacecraft noise from geomagnetic field data. This method assumes no prior knowledge of the number, location, or amplitude of noise signals, but assumes that they have minimal overlapping spectral properties. The present teachings demonstrate the validity of this method by separating high latitude magnetic perturbations recorded by the low-Earth orbiting satellite, SWARM, from noise signals in simulation and in a laboratory experiment using a mock CubeSat apparatus. In the case of more noise sources than magnetometers, this problem is an instance of underdetermined blind source separation (UBSS). This work presents a UBSS signal processing system and method to remove spacecraft noise and minimize or eliminate the need for a mechanical boom.

In some embodiments, a system and method for separating spacecraft generated magnetic noise from geomagnetic field data is provided that includes detecting magnetic information using a plurality of magnetometers and outputting a representative detection signal, the magnetic information comprising magnetic noise and geomagnetic field data; receiving the representative detection signal including the magnetic information from the plurality of magnetometers and resolving the magnetic noise data and the geomagnetic field data by analyzing a relative gain and phase of the representative detection signal in a time-frequency (TF) domain, if a noise signal of the representative detection signal is sparse in the TF domain then the relative gain and phase is found using cluster analysis, and separating a geomagnetic field data signal from a magnetic noise data signal using cluster centroids in compressive sensing based on the cluster analysis. The system and method further include outputting a resultant signal representative of at least the geomagnetic field data signal.

Further areas of applicability will become apparent from the description provided herein. The description and specific examples in this summary are intended for purposes of illustration only and are not intended to limit the scope of the present disclosure.

DRAWINGS

The drawings described herein are for illustrative purposes only of selected embodiments and not all possible implementations, and are not intended to limit the scope of the present disclosure.

FIG. 1 illustrates three magnetometer measurements of six computer-simulated sinusoidal noise signals. Each magnetometer signal is transformed into the time-frequency (TF) domain using the STFT. The magnitude of the three resulting TF signals is taken and plotted against each other in a scatter plot. The scattered TF points from each magnetometer form straight lines due to Equation 4. This figure does not include the phase subdomain of the H-domain.

FIG. 2 illustrates the scattered time-frequency mixed signals in FIG. 1 are projected onto a half-unit hypersphere through normalization. The six scattered straight lines collapse into six compact clusters. The centroid of each cluster is proportional to each source signals' mixing vector in the mixing matrix, K, due to Equation 4.

FIG. 3 illustrates noise signals in the full H-domain for a two-magnetometer system. The horizontal axes represent the magnitude of the time-frequency magnetometer signals projected onto a unit hypersphere. The vertical axis represents the relative argument of Sensor 2 in radians as defined by Equation 7. The data points are projected onto a plane at Z=−2.5 to distinguish the difference in magnitudes. The phase and magnitude of each noise signal at each magnetometer is discovered by clustering the data in this format.

FIG. 4 illustrates a flowchart process involved in using cluster analysis to discover noise signals and compressive sensing to separate the ambient magnetic field from noise signals.

FIGS. 5A-5E illustrate 10s of four source signals used to simulate spacecraft noise and one signal to simulate the ambient magnetic field, wherein FIG. 5A illustrates the ambient magnetic field signal using SWARM A data starting from 17 Mar. 2015 at 8:53 UTC, FIG. 5B illustrates a 2 Hz sine wave with amplitude of 50 nT, FIG. 5C illustrates a 3 Hz square wave with a magnitude of 100 nT, FIG. 5D illustrates a sine wave with a frequency of 5 Hz and amplitude of 50 nT and FIG. 5E illustrates a sawtooth wave with an amplitude of 110 nT and frequency of 0.7 Hz.

FIGS. 6A-6C illustrate plots (a, b, and c) for one hundred seconds of three magnetometer signals, b(t), created by mixing the five source signals in FIGS. 5A-5E though the mixing matrix defined in Equation 13.

FIGS. 7A-7E illustrate related data, wherein FIG. 7A illustrates the cleaned magnetometer signal in blue with the ambient magnetic field signal overlayed in orange, FIG. 7B illustrates a spectrogram of the uncleaned signal from magnetometer (a) in FIG. 6A, FIG. 7C illustrates a spectrogram of the reconstructed ambient magnetic field signal, FIG. 7D illustrates a spectrogram of the true ambient magnetic field signal (the spectrograms were created using wavelet analysis. The shaded areas indicate where the wavelet does not produce valid results.), and FIG. 7E illustrates a histogram of the signal reconstruction error, s₁-{circumflex over ( )}s₁.

FIG. 8 illustrates a Mock CubeSat Apparatus with three PNI RM3100 magnetometers and four copper coils driven by signal generators. The magnetometers are placed within the mock CubeSat. The Apparatus is placed inside a mu-metal lined copper room that acts as a large magnetic shield.

FIGS. 9A-9C illustrate 18.5 s of three mixed signals recorded by PNI RM3100 magnetometers' z-axis. The five signals present are two sine waves, two square waves, and the added residual magnetic field data. The noise signals have amplitudes between 50 and 500 nT compared to the ambient magnetic field signal with a max amplitude near 300 nT.

FIGS. 10A-10D illustrate reconstructed Sine and Square wave signals from 50 Hz mixed signals in FIGS. 9A-9C.

FIGS. 11A-11E illustrate related data, wherein FIG. 11A illustrates the cleaned magnetometer signal in blue with the ambient magnetic field signal overlayed in orange, FIG. 11B illustrates a spectrogram of the uncleaned signal from magnetometer (a) in FIG. 9A, FIG. 11C illustrates a spectrogram of the reconstructed ambient magnetic field signal, FIG. 11D illustrates a spectrogram of the true ambient magnetic field signal (the spectrograms were created using wavelet analysis. The shaded areas indicate where the wavelet does not produce valid results.), and FIG. 11E illustrates a histogram of the signal reconstruction error, s₁-{circumflex over ( )}s₁.

FIG. 12 is a table summarizing the results of Experiment 1.

FIG. 13. is a table summarizing the results of Experiment 2.

Corresponding reference numerals indicate corresponding parts throughout the several views of the drawings.

DETAILED DESCRIPTION

Example embodiments will now be described more fully with reference to the accompanying drawings. Example embodiments are provided so that this disclosure will be thorough, and will fully convey the scope to those who are skilled in the art. Numerous specific details are set forth such as examples of specific components, devices, and methods, to provide a thorough understanding of embodiments of the present disclosure. It will be apparent to those skilled in the art that specific details need not be employed, that example embodiments may be embodied in many different forms and that neither should be construed to limit the scope of the disclosure. In some example embodiments, well-known processes, well-known device structures, and well-known technologies are not described in detail.

The terminology used herein is for the purpose of describing particular example embodiments only and is not intended to be limiting. As used herein, the singular forms “a,” “an,” and “the” may be intended to include the plural forms as well, unless the context clearly indicates otherwise. The terms “comprises,” “comprising,” “including,” and “having,” are inclusive and therefore specify the presence of stated features, integers, steps, operations, elements, and/or components, but do not preclude the presence or addition of one or more other features, integers, steps, operations, elements, components, and/or groups thereof. The method steps, processes, and operations described herein are not to be construed as necessarily requiring their performance in the particular order discussed or illustrated, unless specifically identified as an order of performance. It is also to be understood that additional or alternative steps may be employed.

When an element or layer is referred to as being “on,” “engaged to,” “connected to,” or “coupled to” another element or layer, it may be directly on, engaged, connected or coupled to the other element or layer, or intervening elements or layers may be present. In contrast, when an element is referred to as being “directly on,” “directly engaged to,” “directly connected to,” or “directly coupled to” another element or layer, there may be no intervening elements or layers present. Other words used to describe the relationship between elements should be interpreted in a like fashion (e.g., “between” versus “directly between,” “adjacent” versus “directly adjacent,” etc.). As used herein, the term “and/or” includes any and all combinations of one or more of the associated listed items.

Although the terms first, second, third, etc. may be used herein to describe various elements, components, regions, layers and/or sections, these elements, components, regions, layers and/or sections should not be limited by these terms. These terms may be only used to distinguish one element, component, region, layer or section from another region, layer or section. Terms such as “first,” “second,” and other numerical terms when used herein do not imply a sequence or order unless clearly indicated by the context. Thus, a first element, component, region, layer or section discussed below could be termed a second element, component, region, layer or section without departing from the teachings of the example embodiments.

Spatially relative terms, such as “inner,” “outer,” “beneath,” “below,” “lower,” “above,” “upper,” and the like, may be used herein for ease of description to describe one element or feature's relationship to another element(s) or feature(s) as illustrated in the figures. Spatially relative terms may be intended to encompass different orientations of the device in use or operation in addition to the orientation depicted in the figures. For example, if the device in the figures is turned over, elements described as “below” or “beneath” other elements or features would then be oriented “above” the other elements or features. Thus, the example term “below” can encompass both an orientation of above and below. The device may be otherwise oriented (rotated 90 degrees or at other orientations) and the spatially relative descriptors used herein interpreted accordingly.

Generally, the use of magnetometers for space exploration is inhibited by magnetic noise generated by spacecraft electrical systems. Mechanical booms are traditionally used to extend magnetometers away from noise sources. If a spacecraft is equipped with multiple magnetometers, signal processing methods can be used to compare magnetometer measurements and remove stray magnetic noise signals. The present teachings provide the use of density-based cluster analysis to identify spacecraft noise signals and compressive sensing to separate spacecraft noise from geomagnetic field data. This method assumes no prior knowledge of the number, location, or amplitude of noise signals, but assumes that they have minimal overlapping spectral properties. The present teachings demonstrate the validity of this method by separating high latitude magnetic perturbations recorded by the low-Earth orbiting satellite, SWARM, from noise signals in simulation and in a laboratory experiment using a mock CubeSat apparatus. In the case of more noise sources than magnetometers, this problem is an instance of underdetermined blind source separation (UBSS). This work presents a UBSS signal processing method to remove spacecraft noise and minimize the need for a mechanical boom using multiple magnetometers.

More particularly, the present teachings enable research quality magnetic field observations (such as the geomagnetic field) from noisy platforms such as satellites. The system comprises three tri-axial magnetometers with two being surface-mounted to the platform and the third attached to any appendage or short boom. The three sensors all measure the large-scale background field (such as the geomagnetic field) and the noise from the platform (consisting of possibly multiple sources of noise), but the platform noise will have different amplitudes at each of the three sensors. Using machine learning algorithms, the mixed signals from all three sensors can be deconvolved to separate the signal of interest from all of the noise.

Magnetometers are instruments designed to measure magnetic fields. They are used for a variety of purposes such as monitoring the magnetic field of the Earth from spacecraft. Spacecraft systems, such as solar panels and reaction wheels, generate magnetic noise that interferes with magnetometer readings. If the spacecraft has multiple magnetometers, each noise source will have a different magnitude at each magnetometer depending on the location of the noise source. The system which describes the magnitude of each noise source at each magnetometer is called a mixing matrix. We provide the use of unsupervised machine learning to estimate the mixing matrix. Once the mixing matrix is estimated, the Earth's magnetic field can be separated from spacecraft magnetic noise using a method called Compressive Sensing. Spacecraft often have many more noise sources than magnetometers, which complicates noise cancellation. The provided method has the ability to clean noisy magnetometer measurements when there are more noise signals present than magnetometers.

METHODOLOGY

According to the principles of the present teachings, we apply a two-step approach to removing spacecraft noise and reconstructing the ambient magnetic field. The first step is to discover the mixing matrix, K, defined in Equation 1. This is achieved by preprocessing the magnetometer data into a clusterable form and applying a clustering algorithm. The second step is to reconstruct the ambient magnetic field and noise signals using compressive sensing. In this step, the mixing matrix, K, is used to demix the magnetometer signals via convex optimization. This two-step process is designed to be applied to each magnetometer axis separately.

Signal Preprocessing

The separation of magnetic field signals from stray magnetic fields is analogous to a problem thoroughly researched in other fields, such as acoustics, and is called UBSS. This problem has been heavily investigated for microphone and radar arrays, but the unique spatial structure of magnetic fields introduces new complications which have not been well-researched. When considering a dipole noise source, the placement of magnetometers at different magnetic latitudes alters the magnitude and polarity of the noise signal. We model this effect as a phase, despite the noise sources mixing instantaneously. The time-frequency (TF) domain mixing model, B(t,k)=KS(t,k), is defined by the following system:

$\begin{array}{cc}\left[\begin{array}{c}{B}_1\left(t,k\right)\\ B_2\left(t,k\right)\\ ⋮\\ B_m\left(t,k\right)\end{array}\right]=\left[\begin{array}{ccccc}1& k_{12}{\mathrm{\angle \varphi }}_{12}& k_{13}{\mathrm{\angle \varphi }}_{13}& \dots & k_{1n}{\mathrm{\angle \varphi }}_{1n}\\ 1& k_{22}{\mathrm{\angle \varphi }}_{22}& k_{23}{\mathrm{\angle \varphi }}_{23}& \dots & k_{2n}{\mathrm{\angle \varphi }}_{2n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 1& k_{m2}{\mathrm{\angle \varphi }}_{m2}& k_{m3}{\mathrm{\angle \varphi }}_{m3}& \dots & k_{\mathrm{mn}}{\mathrm{\angle \varphi }}_{\mathrm{mn}}\end{array}\right]\left[\begin{array}{c}{S}_1\left(t,k\right)\\ S_2\left(t,k\right)\\ ⋮\\ S_n\left(t,k\right)\end{array}\right]& \left(2\right)\end{array}$

In this mixing system, the signals Si(t,k) are the source signals at time t and frequency k. The ambient magnetic field signal we seek to recover, S1(t,k), is assumed to be identical at each magnetometer and is represented by a column of ones in the mixing matrix. In the geospace environment, this allows us to observe phenomena such as ultralow frequency (ULF) waves which have frequencies less than 5 Hz. The phases, φ_ij={0, π}, in the mixing matrix, K, account for the difference of a signal seen by magnetometers at different magnetic latitudes. The phase, φ_ij, is determined by the spatial structure of the noise signal, which depends on the geometry of the noise source. This model does not require that noise sources be dipolar. The value in the mixing matrix k_ij∠φ_ij represents the complex value k_ije^jφ_ij. This value defines presence of the signal S_j(t,k) at magnetometer B_i(t,k).

Sparsity is a precondition of both mixing matrix estimation and compressive sensing, however, spacecraft noise signals are not often sparse in the time domain. The magnetometer signals, b(t), are transformed into the TF domain using a Fourier transform in order to increase signal sparsity. Typically, the Short-Time Fourier Transform (STFT) is used because signals that are present in multiple time windows will provide more data points to be clustered [id=Revision Two]. As a result, periodic signals are easier to identify and remove than aperiodic signals. However, aperiodic signals can be separated with sufficient TF resolution. In the present disclosure, we use the Non-Stationary Gabor Transform (NSGT) to transform magnetometer signals into the TF domain. NSGT has advantages over the STFT because it allows the user to evolve the window size with respect to frequency. As a result, high and low frequencies are not limited to the same window size, and frequency resolution is constant across the frequency spectrum. In order to apply the NSGT, the user specifies a quality, Q, and the lowest frequency they would like to observe. The parameter, Q, is used to automatically calculate the window size with respect to the desired frequency resolution. In comparison to the STFT, the NSGT provides finer frequency resolution at low frequencies and better time resolution at higher frequencies. We perform the NSGT to obtain the UBSS model B(t,k)=KS(t,k). The mixing system of a sparse TF bin where only the signal, S_j(t,k), is present can be defined by a single mixing vector:

$\begin{array}{cc}\left[\begin{array}{c}❘B_1\left(t,k\right)❘\\ ❘B_2\left(t,k\right)❘\\ ⋮\\ ❘B_m\left(t,k\right)❘\end{array}\right]=\left[\begin{array}{c}{k}_{1j}\\ k_{2j}\\ ⋮\\ k_{\mathrm{mj}}\end{array}\right]❘S_j\left(t,k\right)❘& \left(3\right)\end{array}$

The operator |a+jb| applied to the complex value a+jb returns the magnitude of complex value, √(a²+b²). Equation 3 can be rewritten element-wise as:

$\begin{array}{cc}❘S_j\left(t,k\right)❘=\frac{❘B_1\left(t,k\right)❘}{k_{1j}}=\frac{❘B_2\left(t,k\right)❘}{k_{2j}}=\dots =\frac{❘B_m\left(t,k\right)❘}{k_{\mathrm{mj}}}& \left(4\right)\end{array}$

Equation 4 is equivalent to the symmetric form of a line with slope defined by the mixing vector of the noise signal. In order to find the mixing vector of a noise signal, we define a TF space Hϵ^2m in which each phase and magnitude of the m magnetometer signals are an axis. Sparse TF points will draw straight lines through the origin in the H-domain with a slope proportional to the signal's mixing vector. FIG. 1 shows an example of a scatter plot of three mixed TF signals composed of six source signals. The mixed signals form straight lines with slopes defined by Equation 4.

Mixing Matrix Estimation

The slope of the lines drawn through the H-domain are not easily clusterable in their current form as a collection of scattered data points. We transform the scattered data points in H-domain into a clusterable form by projecting the magnitude subdomain onto a unit hypersphere. The H-domain magnitude data are projected onto a half-unit hypersphere by normalizing the TF magnetometer data via the following equation.

$\begin{array}{cc}{B}^{*}\left(t,k\right)=\frac{❘B\left(t,k\right)❘}{B\left(t,k\right)}& \left(5\right)\end{array}$

When the scattered data points have been normalized, they collapse into compact clusters. This is illustrated by the projection of the scattered data points representing six computer generated signals in FIG. 1 onto a half-unit hypersphere in FIG. 2. The centroid of a cluster is proportional to the mixing vector of a noise signal as defined in Equation 2.

The majority of the frequency space is filled with negligible energy points that will project randomly onto the unit hypersphere. We cleanse the data of these points using a magnitude filter. The filter is applied by finding the average signal magnitude and removing data points below a factor, A, of the average signal magnitude. The magnitude filter is applied by removing data points that do not satisfy the following criterion:

$\begin{array}{cc}❘B\left(t,k\right)❘>\lambda ·\mathrm{avg}\left(❘B\left(t,k\right)❘\right)& \left(6\right)\end{array}$

The projected data points form tightly clustered groups on the unit hypersphere that allow us to discover the relative gain between noise signals at different magnetometers. However, we need to find the relative phases between noise signals of magnetometers at different positions. To account for this we join each projected TF point to its relative argument. The relative argument is defined by the following transformation:

$\begin{array}{cc}\arg B\left(t,k\right)=\left\{\arg B_j\left(t,k\right)-\arg (B_0\left(t,k\right)❘j\in \left[0,m\right]\right\}& \left(7\right)\end{array}$

Using the result of Equation 7, we define a new data format, H(t,k), by concatenating the projected magnitude data with the argument of the TF data.

$\begin{array}{cc}H\left(t,k\right)=(B*\left(t,k\right),\arg \left(B\left(t,k\right)\right)& \left(8\right)\end{array}$

The magnetometer data, H(t,k), are now in a format that can be clustered to discover the gain and phase of each signal described in the mixing matrix, K FIG. 3 shows an example of two magnetometer signals transformed into the H-Domain.

Now that the projected magnitude and relative phases are joined, a variety of clustering algorithms can be applied to find the mixing matrix, K In the present disclosure, we use the Density Based Spatial Clustering for Applications with Noise (DBSCAN) method because it does not require user input to discern the number of clusters present, and it will ignore noise points. As a result, the number of noise signals does not need to be defined prior to the application of DBSCAN. DBSCAN has two essential parameters, eps and minPts, that allow this functionality. The maximum distance for two points to become neighbors is the value, eps. If a point has minPts number of neighbors, it is called a core point. Core points are used to define each cluster. If a point is more than eps distance away from any point in a cluster, it is labeled as noise. We use DBSCAN to cluster H(t,k) and use each cluster's centroid as the noise signal's mixing vector. Once the mixing vector of each noise signal is known, the mixing vectors are joined to form the mixing matrix, K. The mixing matrix is used to separate the noise signals from the ambient magnetic field via compressive sensing.

Signal Reconstruction

Compressive sensing is a method used to reconstruct sparse signals with a sampling rate below two times a signal's bandwidth. Reconstructing a signal of length N from a sampled signal of length M, where M<N, is an analogous problem to UBSS. Ordinarily, the system b=Ks, where K is a wide matrix, has infinitely many solutions because if b=Ks is a solution, b=K(s+s′) is also a solution for any vector s′ in the null space of K. Compressive sensing can exactly recover sparse signals and approximate near-sparse signals through minimizing the L1 norm of S with respect to b−Ks<ε. The algorithm works with O(N³) complexity.

We use CVXPY, A Python-Embedded Modeling Language for Convex Optimization, to reconstruct the signals with the estimated mixing matrix, K. The formulation used to recover the signal, s, from b is:

$\begin{array}{cc}\begin{array}{cc}\mathrm{Minimize}& {\omega }^T❘s❘\\ \mathrm{Subject}\mathrm{to}& \mathrm{Ks}=b\end{array}& \left(9\right)\end{array}$

Traditionally, compressive sensing minimizes the L1 norm of the source signals, ∥s∥₁, with respect to Ks=b in order to recover the source signals. Instead of minimizing the L1 norm, we utilize a weighted L1 norm defined by the weighting vector, ω=[ω₁, 1, 1, . . . , 1]^T, where ω₁≥1. The parameter, ω₁, is multiplied with the ambient magnetic field signal, s₁, in order to deter the attribution of energy from other noise signals to it. In the case that the source signals, s, are not sparse at a TF bin, the additional weight increases the cost of attributing energy from other signals to the ambient magnetic field, s₁. The optimal value of the weight, ω₁, depends on the signature of noise signals. Candés et al. (2008) apply a similar approach by iteratively adjusting the weight of each signal with respect to the magnitude of the signal. In the present disclosure, we found the optimal weight, ω₁, experimentally by comparing the reconstructed signal, ŝ₁, to the true signal, s₁.

This system defined in Equation 9 is solved using the Embedded Conic Solver (ECOS) of Domahidi et al., 2013. The ECOS method is a convex optimization algorithm that transforms the problem defined in Equation 9 into a Second Order Cone Problem (SOCP). SOCP problems are convex optimization problems that minimize linear functions with respect to second order cone constraints. The ECOS method applies an interior point solver to converge on the sparse solution of the problem defined by Equation 9.

In some embodiments, the method here can be varied to transform the data into a clusterable format with a few additional updates. The completely transformed magnetometer data are referred to as H(t,k).

The magnetometer measurements, b(t), contain both intrinsic measurement noise and multi-source points (MSP) where the time-frequency distributions of different source signals overlap. The MSPs interfere with the identification of the source signals. Before transforming the data into a clusterable format, we transform the time-series data, b(t), into the time-frequency domain, B(t,k), using the NSGT transform in order to improve sparsity. We then apply a magnitude filter and a single-source point (SSP) filter to remove the intrinsic noise and MSPs.

The measurement uncertainty of the PNI RM3100 magnetometer at 1 Hz is near 2.7 nT. As a result, time-frequency bins with low energy may project randomly onto the unit hypersphere. We cleanse these low energy data points using a magnitude filter which removes data points below a factor of the magnetometer's measurement uncertainty. This filter is defined in Equation 10.

$\begin{array}{cc}❘B\left(t,k\right)❘>\lambda ·\sigma & \left(10\right)\end{array}$

The parameter, σ, defines the standard deviation of the characterized measurement uncertainty of the magnetometer. The parameter, λ, is the factor that defines the magnitude at which low energy points should be removed with respect to the measurement uncertainty, σ. Only statistically meaningful data points will remain after the magnitude filter has been applied.

To identify and remove MSPs, we use a SSP detection algorithm that compares the complex and real components of the mixed magnetometer signals, B(t,k). If only one source signal, Sj(t,k), contributes to a time-frequency bin, (t, k), then both the real and imaginary components of B(t,k) will be proportional to the same mixing vector, K_j. This means that their cosine similarity will be equal to 1. We calculate the cosine similarity using Equation 11 and remove data points that are above a certain threshold, θ₀.

$\begin{array}{cc}\theta ={\cos }^{-1}\left(\frac{{\left\{\mathrm{Re}\left[B\left(t,k\right)\right]\right\}}^T\mathrm{Im}\left[B\left(t,k\right)\right]}{\sqrt{{\left\{\mathrm{Re}\left[B\left(t,k\right)\right]\right\}}^T\mathrm{Re}\left(B\left(t,k\right)\right)}\sqrt{{\left\{\mathrm{Im}\left[B\left(t,k\right)\right]\right\}}^T\mathrm{Im}\left(B\left(t,k\right)\right)}}\right)& \left(11\right)\end{array}$

Once the original magnetometer data, b(t), has been filtered and transformed into a clusterable format H(t,k), a number of clustering algorithms can be used to reconstruct the mixing matrix, K. As described herein, Density Based Spatial Clustering for Applications with Noise (DBSCAN) algorithm can be used. This algorithm requires two data-dependent parameters, eps and minPts, that must be set by the user. In the present embodiment, we use the Hierarchical Density Based Spatial Clustering for Applications with Noise (HDBSCAN) algorithm because it has similar benefits to DBSCAN but does not require the user to set any parameters. HDBSCAN is an extension of DBSCAN that uses a hierarchical clustering scheme. Similar to DBSCAN, HDBSCAN does not require the number of clusters to be defined beforehand, and it will ignore noisy data points. We use HDBSCAN to cluster H(t,k) and use each cluster's centroid as a source signal mixing vector. We then join these vectors to form the mixing matrix, K, which we use to separate the source signals from ambient magnetic field via CS.

Source Signal Separation

To separate the stray magnetic field signals from the ambient magnetic field, we first estimate the mixing matrix, K, using HDBSCAN. Then, we use CS to reconstruct the ambient magnetic field signal from the mixed signals. CS is an algorithm that can recover sparse or near-sparse signals from under-sampled measurements. In some embodiments, one can apply a weighted Basis Pursuit scheme for CS, which solves the following optimization problem in Equation 12.

$\begin{array}{cc}\begin{array}{cc}\mathrm{Minimize}& {\omega }^T❘s❘\\ \mathrm{Subject}\mathrm{to}& b=\mathrm{Ks}\end{array}& \left(12\right)\end{array}$

However, this scheme assumes that there is no measurement error in b. If the measured signals are very noisy, this assumption may not hold. Therefore, we propose a different scheme based on the Dantzig Selector with an iterative weighting scheme. The Dantzig Selector uses the L-infinity norm to measure the error in the residual vector, which makes it more robust to outliers and noise than the least squares method. It also selects fewer source signals to explain the data, which enhances the sparsity of the solution. The formulation of the Dantzig Selector is shown in Equation 13.

$\begin{array}{cc}\begin{array}{cc}\mathrm{Minimize}& {\omega }^T❘s❘\\ \mathrm{Subject}\mathrm{to}& {❘K^T\left(b-K_S\right)}_{\infty }\le \eta \end{array}& \left(13\right)\end{array}$

The parameter, n, in the problem constraint is a tuning parameter that controls how much error is allowed. This formulation enables us to recover the ambient magnetic field signal even when there is some measurement error, e, in the underdetermined measurement system b=Ks+e.

We use a weighting scheme to reduce noise and error in the signal reconstruction process. The weighting scheme assigns different importance to different elements in the solution vector, s, which represents the separated signals. The first element of s, denoted by s₁, corresponds to the ambient magnetic field signal, which is our main interest. The other elements correspond to the stray magnetic field signals, which are considered as noise. In previous embodiments, we use a constant weighting vector that gives more weight to s₁ than to other elements of s. However, this may not be optimal for different types of mixed signals. Therefore, we propose an adaptive weighting scheme that adjusts the weight of s₁ based on whether the mixed time-frequency signal, B(t,k), is a SSP or a MSP. A SSP means that only one source signal contributes to B(t,k), while an MSP means that multiple source signals contribute to B(t,k). If B(t,k) is a SSP, we use the Candés et al. weighting scheme where ω=1/|s|, which assigns less weight to larger elements of s to induce sparsity. If B(t,k) is an MSP, we iteratively increase the weight of s₁ using the formula w₁=w₁+α(s_ratio−w₁), where α is a learning rate, w₁ is the first element of w corresponding to s₁, and s_ratio is the ratio of the sum of absolute values of noise signals to the absolute value of ambient magnetic field signal. The variable, s_ratio, is defined in the following equation where s is the solution vector and ε is a small constant to prevent division by zero.

$\begin{array}{cc}{s}_{\mathrm{ratio}}=\frac{\sum _{l=2}^n❘s_i❘}{❘s_1❘+e}& \left(14\right)\end{array}$

This scheme aims to enhance the sparsity of s by suppressing noise signals and highlighting the ambient magnetic field signal. We use this adaptive weighting scheme with Dantzig Selector, which is a CS constraint that can handle measurement errors. This system allows us to determine and reconstruct the ambient magnetic field with as little noise as possible from under-sampled measurements.

We use CVXPY, a Python-embedded modeling language for convex optimization problems, to solve the system defined by Equation 13. CVXPY automatically transforms the problem into standard form, calls a solver, and unpacks the results. The solver we use is Embedded Conic Solver (ECOS), which converts the problem into a Second Order Cone Problem (SOCP) and applies an interior point solver to find the sparse solution.

System

Generally, in some embodiments as illustrated in FIG. 8, a system 10 is provided for separating spacecraft generated magnetic noise from geomagnetic field data. System 10 can comprise a plurality of magnetometers 12 each configured to detect magnetic information 1000 and outputting a representative detection signal. The magnetic information 1000 comprising magnetic noise and geomagnetic field data. A processing system 14 is configured to receive the representative detection signal including the magnetic information from the plurality of magnetometers 12 and resolve the magnetic noise data and the geomagnetic field data by: analyzing a relative gain and phase of the representative detection signal in a time-frequency (TF) domain, if a noise signal of the representative detection signal is sparse in the TF domain then the relative gain and phase is found using cluster analysis, and separating a geomagnetic field data signal from a magnetic noise data signal using cluster centroids in compressive sensing based on the cluster analysis. The processing system 14 then outputs a resultant signal representative of at least the geomagnetic field data signal.

Experimental Data and Results

We tested the present method of signal and noise separation through two experiments. The first experiment demonstrates the separation of SWARM magnetic field data from computer-simulated signals using virtual magnetometers. The second experiment demonstrates the separation of SWARM magnetic field data from real magnetic noise signals generated with copper coils. The coil-generated magnetic fields were measured using the PNI RM3100 magnetometer and a mock CubeSat described by Deshmukh et al. (2020).

FIG. 4 details the process of identifying noise signals and reconstructing the ambient magnetic field. First (i), the signal offsets are subtracted to center the signals around 0 nT. Second (ii), the signals are transformed into the TF domain using the NSGT to increase signal sparsity. Third (iii), low energy points are filtered out using Equation 6. Fourth (iv), the signals are transformed into H(t,k) by projecting the magnitude, |B(t,k)| onto the unit hypersphere and concatenating it with the phase, arg B(t,k), via Equations 5, 7 and 8. Fifth (e), the data, H(t,k), are clustered using DBSCAN and the cluster centroids are found. Finally, in the last step (vi), compressive sensing is used to reconstruct the ambient magnetic field. The minimum magnitude, A in step (iii), and the parameters eps and MinPts in step (v) may need to be adjusted depending on the length and magnitude of the signals being analyzed.

We evaluate the separation of noise signals via three metrics. The metrics are calculated point-wise using the reconstructed signal, x and the true signal, y, over N data points. The first metric is the Pearson Correlation Coefficient. This measurement gives the covariance between the normalized input and recovered signals.

$\begin{array}{cc}\rho =\frac{{\sum }_{i=0}^{N-1}\left(x_i-\stackrel{\_}{x}\right)\left(y_i-\stackrel{\_}{y}\right)}{\sqrt{{\sum }_{i=0}^{N-1}{❘\left(x_i-\stackrel{\_}{x}\right)❘}^2{\sum }_{i=0}^{N-1}{❘\left(y_i-\stackrel{\_}{y}\right)❘}^2}}& \left(10\right)\end{array}$

The second metric evaluated is the root mean squared error (RMSE). This metric is proportional to the magnitude of the squared error. As a result, the RMSE is very sensitive to large errors.

$\begin{array}{cc}\mathrm{RMSE}=\sqrt{\frac{{\sum }_{i=0}^{N-1}{\left(x_i-y_i\right)}^2}{N}}& \left(11\right)\end{array}$

The final metric is the normalized RMSE (NRMSE). This metric yields the RMSE as a percentage of the magnitude of the signal being measured. It is used to compare the relative error between signals on different orders of magnitude. We calculate the NRMSE by dividing the RMSE of the signal by the max amplitude of the absolute value of the true, detrended signal, |y−y|_max.

$\begin{array}{cc}\mathrm{NRMSE}=\mathrm{RMSE}/{❘y-\stackrel{\_}{y}❘}_{\max }& \left(12\right)\end{array}$

Experiment 1: Computer Simulation

In this experiment, we use four simulated noise signals, s(t)∪[s₂(t), s₃(t), s₄(t), s₅(t)], and three virtual magnetometers b(t)=Ks(t)=[b₁(t), b₂(t), b₃(t)]. The signal, s₁(t), is residual magnetic field data created by subtracting data generated by the IGRF model from SWARM magnetic field data. This process leaves only magnetic perturbations present in the magnetosphere. The magnetic perturbation data we use were measured by the SWARM A satellite on 17 Mar. 2015 between 8:53 and 8:55 UTC. This part of the orbit passes between the 69th and 76th parallel south and was selected to capture perturbations in the southern auroral zone. The present method detailed in FIG. 4 is tested on 100 s of data, although it may be applied to a signal of any length provided that there are enough data points to cluster. The signals are combined through the complex mixing matrix in Equation 13 with phases given in radians.

$\begin{array}{cc}K=\left[\begin{array}{ccccc}1\mathrm{\angle 0}& 0.99\mathrm{\angle 0}& 0.09\mathrm{\angle 0}& 0.7\mathrm{\angle 0}& 0.02\mathrm{\angle 0}\\ 1\mathrm{\angle 0}& 0.09\mathrm{\angle \pi }& 0.99\mathrm{\angle 0}& 0.7\mathrm{\angle 0}& 0.05\mathrm{\angle \pi }\\ 1\mathrm{\angle 0}& 0.12\mathrm{\angle \pi }& 0.12\mathrm{\angle \pi }& 0.13\mathrm{\angle \pi }& 0.99\mathrm{\angle \pi }\end{array}\right]& \left(13\right)\end{array}$

The values in the first column represent the ambient magnetic field signal which appears identically at every magnetometer. FIGS. 5A-5E shows the five source signals used in this simulation. Two of the noise signals are sine waves with frequencies of 2 and 5 Hz. Sine waves are sparse signals that can be represented by a single point in the frequency domain. This makes them easily identifiable by cluster analysis. The two remaining noise signals used are a sawtooth wave with a frequency of 0.7 Hz, and a square wave with a frequency of 3.0 Hz. These signals inhabit a broad frequency spectrum and diminish the sparsity of the mixed signals.

The signals are combined in the mixing system b(t)=Ks(t) with the mixing matrix K from Equation 13. The resulting signals are sampled by the virtual magnetometers at a rate of 50 samples per second. Different noise signals, such as noise generated by reaction wheels, may have higher frequency components and require a higher sampling rate in order to avoid aliasing. A random normal signal with a standard deviation of 6 nT is added to each virtual magnetometer in order to simulate instrument noise. This noise level corresponds to the rated instrument resolution of the PNI RM3100 magnetometer at 50 Hz used in Experiment 2. FIGS. 6A-6C shows the sampled signals.

Following the procedure in FIG. 4, the signals were detrended and transformed into the TF domain using the NSGT. The NSGT is a type of constant-Q transform, so it requires the parameter Q which specifies window size. In this experiment, we used Q=10 and a lower frequency bound of 30 mHz. In step (iv), low energy points were removed using a λ=0.5. The resulting data were transformed into H(t,k) and clustered by DBSCAN with parameters eps=0.3 and MinPts=4. These parameters were optimized experimentally using trial and error, however it may be possible to automate parameter selection based on the signals being analyzed. With this configuration, DBSCAN discovered the five clusters corresponding to each noise source. The clusters, shown in the columns of {circumflex over ( )}K, closely match the original mixing matrix.

Finally, in step (vi), the mixed signals were separated by compressive sensing using the recovered mixing matrix, {circumflex over ( )}K, in Equation 14. The data, H(t,k), are discarded and the raw Fourier transform of the mixed signals are separated by applying the ECOS algorithm to the problem defined in Equation 9 with a weight of ω₁=1.5. The reconstructed SWARM perturbation signal is shown in FIG. 7, as well as a histogram of the reconstruction error and spectrograms of the noisy, cleaned, true SWARM signal.

The reconstructed ambient magnetic field signal resembles the original signal with some additional error. In order to evaluate the reconstruction noise, the Pearson correlation coefficient, RMSE, and NRMSE of each source signal are calculated. The ambient magnetic field was reconstructed with a RMSE of 2.75 nT. The results for the reconstruction of each source signal are shown in the table of FIG. 12. The experiment was repeated without the addition of the 6 nT instrument noise to evaluate the effect of the random noise on the total reconstruction error.

Experiment 2: Magnetic-Coil-Generated Signal Separation

In this experiment, we demonstrate the utility of the present method on real magnetic field data. We use three PNI RM3100 magnetometers to record copper coil-generated noise signals. Four copper coils are driven by signal generators to create the source signals, s(t)∪[s₂(t), s₃(t), s₄(t), s₅(t)]. The signals are combined in the unknown mixing system, b(t)=Ks(t)=[b₁(t), b₂(t), b₃(t)]. The SWARM residual magnetic field data, which is used in experiment one, is added to each magnetometer recording to generate the ambient magnetic field signal, s₁(t).

The present method detailed in FIG. 4 is tested on 100 s of recorded data. The signals, s₂(t) and s₃(t), are sine waves with frequencies of 0.4 and 0.8 Hz. The signals, s₄(t) and s₅(t), are square waves with frequencies of 1 and 2 Hz. The three PNI RM3100 magnetometers and four copper coils are placed on the CubeSat apparatus as shown in FIG. 8. Due to the location and orientation of the four copper coils and three magnetometers, each noise signal will appear at each magnetometer with a different magnitude and magnetic latitude induced phase. Additionally, this experiment was performed in a copper room lined with mu-metal in order to screen out magnetic fields from the surrounding environment.

The PNI RM3100 is a magneto-inductive magnetometer that measures the magnetic field by counting hysteresis loops with a comparator circuit, called a Schmitt Trigger, in an ASIC. The ASIC records magnetic field measurements by adding to a register every time the Schmitt trigger is saturated. This measurement renders the magnetic field when integrated with respect to time. The ASIC has a cycle count register that controls how many clock cycles pass between integrations. The error of the magnetometer will change with respect to the cycle count. In this experiment, each magnetometer is sampled at a rate of 50 Hz with a cycle count of 200 cycles. The PNI RM3100 is rated to have a resolution of 6 nT in this configuration. The mixed signals recorded by the PNI RM3100 magnetometers are shown in FIGS. 9A-9C.

The present method was run on data from the magnetometers' z-axis following the same steps as in FIG. 4 and described herein. The signals were detrended and transformed into the TF domain using the NSGT with a quality factor of Q=20 and a lower frequency bound of 30 mHz. In step 4, low energy points were removed using a λ=2.5. The resulting data were transformed into H(t,k) and clustered by DBSCAN with parameters eps=0.4 and MinPts=4. DBSCAN discovered the following five clusters shown in the columns of {circumflex over ( )}K.

The PNI RM3100 magnetometer was experimentally found to have a lower noise floor when sampled at a higher rate and decimated to a lower rate versus only being sampled at a lower rate. We evaluated this effect by reconstructing the original 50 Hz data in step 6, then downsampling the reconstructed ambient magnetic field signal to 10 Hz, 1 Hz, and averaging the data with a moving mean (N=10). The magnetometer signals were downsampled by applying an eighth order Chebyshev type I anti-aliasing filter and resampling the resulting signal. The mixed signals were separated via weighted compressive sensing using a weight of ω₁=3. The four noise signals reconstructed from the 50 Hz raw data are shown in FIGS. 10A-10D.

The reconstructed coil-generated signals closely resemble square and sine waves with some additional noise. The recovered residual magnetic field data are shown in the top plot of FIG. 11A. The recovered signal is overlayed with the true residual magnetic field signal. The residual data in FIGS. 11A-11E were reconstructed using the mixed signals sampled at the full 50 Hz cadence. The plots show the reconstructed signal, spectrograms of the noisy, cleaned, and true SWARM signal created using wavelet analysis, and a histogram of the signal reconstruction error.

The reconstructed signal closely follows the true geomagnetic perturbation signal with some high frequency noise present. As a result of the geomagnetic field signal being artificially inserted into the magnetometer readings, we are able to calculate the RMSE and Pearson Correlation Coefficient with respect to the original signal. The results for the original, decimated, and moving-mean signals are shown in the table of FIG. 13. These results are also compared to the uncleaned magnetometer data from magnetometer (a) in FIG. 9A.

Discussion

In the present disclosure, we introduced a signal processing method based on UBSS and demonstrated the separation of magnetic noise from geomagnetic field data. In the first experiment, we separated four simulated noise signals from SWARM residual magnetic field data. The noise signals contained both sparse sine wave signals and wideband sawtooth and square wave signals. The method was able to restore the residual magnetic field signal with a correlation coefficient of ρ=0.9988 and RMSE of 2.75 nT. When the experiment was repeated without artificial instrument noise, the method reconstructed the ambient magnetic field signal with a RMSE of 2.84 nT. In the second experiment, we created four magnetic noise signals using copper coils to generate real magnetic field data and placed PNI RM3100 magnetometers within the bus of a mock CubeSat apparatus. The same SWARM magnetic residual data were artificially inserted into the magnetometer measurements. This experiment mimicked the computer-simulated experiment, with two sparse noise signals and two wideband noise signals. At a sampling rate of 50 Hz, the ambient magnetic field signal was reconstructed with a RMSE of 7.94 nT as opposed to 2.75 nT in simulation. The signal separation method was executed using several additional preprocessing techniques such as decimating the sampling rate and applying a moving mean to the magnetometer data. A RMSE of 7.41 nT was achieved by decimating the sample rate to 1 Hz. At 1 Hz, the PNI RM3100 magnetometer is rated to have a measurement error of 2.7 nT due to instrument noise (Regoli et al., 2018). This result places the reconstruction error near the measurement resolution of the magnetometer. When the noisy magnetometer data were decimated, it reduced the RMSE of the signal measured by magnetometer (a) in FIG. 9 from 328.1 to 30.6 nT. In contrast, the decimation of the ambient magnetic field signal reconstructed from the present method did not significantly improve the RMSE. The reconstructed signal decimated to 1 Hz had an RMSE of 7.41 nT compared to 7.94 nT at 50 Hz, however, the UBSS method was able to improve the RMSE by over 20 nT compared to simple downsampling. These results show that the present UBSS method is effective at removing spacecraft noise from magnetic field data.

In general, it is not feasible to adaptively cancel spacecraft noise when a single magnetometer is used. Adaptive noise cancellation requires the removal of noise signals that are time variable. The use of a single magnetometer requires that spacecraft noise be carefully characterized before launch. Otherwise, a change in spacecraft behavior may require special maneuvers to re-characterize noise signatures in situ. The use of multiple magnetometers allows for the discovery of noise signals through the comparison of magnetometer data. Previous systems, such as Sheinker and Moldwin (2016), Deshmukh et al. (2020), and Imajo et al. (2021), each propose methods for noise cancellation using multiple magnetometers. The method proposed by Sheinker and Moldwin (2016) is effective at removing a single noise signal, but is not designed for multiple noise signals. Imajo et al. (2021) propose the use of ICA which is also limited by how many noise signals it can remove. BSS methods require that the number of source signals be less than or equal to the number of mixed signals. Spacecraft contain many electrical systems that could generate magnetic interference, so this condition is rarely met. For example, Pope et al. (2011) identified seven common types of noise signals on Venus Express, which is equipped with two magnetometers. The advantage of the present UBSS method over Imajo et al. (2021) and Sheinker and Moldwin (2016) is that it can cancel noise signals in an underdetermined system. This means that there are more noise signals present than magnetometers. This property of the method provides the flexibility necessary to be applied to many different spacecraft without prior characterization of spacecraft noise. The method also does not require knowledge of magnetometer location and orientation, except that the axis of each magnetometer are aligned. Finally, Deshmukh et al. (2020) designed a state estimation method to transform housekeeping data to magnetic noise signals. Housekeeping currents provide an incomplete mapping of the distribution of currents within a spacecraft. Additionally, housekeeping data are often sampled at a low cadence and may not have the appropriate bandwidth to identify higher frequency noise. The advantage of the present UBSS method over this approach is that it is a blind signal processing algorithm. It requires no housekeeping data to identify and remove noise signals.

The present method functions on the assumption that the noise signals are sparse, meaning that only one noise signal is present at a given frequency. Multiple noise signals may be active at the same time, however, if a signal is not sparse in the frequency domain, then its mixing vector cannot be accurately estimated by cluster analysis. Compressive sensing also requires sparsity in order to accurately reconstruct the separate signals. Compressive sensing can fully reconstruct sparse signals, and approximately reconstruct near-sparse signals. In the present disclosure, we do not exhaustively explore the minimum sparsity required for accurate reconstruction of the ambient magnetic field.

The present method requires that several parameters be set by the user. In the present disclosure, the parameters were manually selected based on the signals being analyzed, but this process could also be automated. The first parameter is the quality factor, Q. This parameter adjusts the window size used in the NSGT. We experimentally selected it, but it may be chosen based on the length of the signal being processed. The parameter, λ, is used to remove low energy noise signals. Data points that are below a fraction, λ, of the average energy data point are removed before clustering occurs. We selected this parameter by analyzing the data projected onto the half-unit hypersphere in FIG. 2, and visually observing if the signals were clusterable. If λ is too small, then the hypersphere will be completely filled with data points, and the noise signals will not be separable. If λ is too large, then small noise signals may not appear at all. Lastly, DBSCAN requires that two parameters, eps, and MinPts, be selected. The parameter, eps, represents the maximum distance allowed for two data points to be considered neighbors. The parameter, MinPts, represents the number of neighbors required for a data point to be considered a core. MinPts may be selected based on the length of signal being processed. A disadvantage of using NSGT and DBSCAN together is that more data points are created for higher frequency signals because the window size is altered based on frequency. Therefore, MinPts should be selected based on the lower frequency signals.

Most heliophysics missions require magnetic field accuracies of better than 1 nT (e.g., the NASA MMS mission (Russell et al., 2016)). Using the PNI RM3100 magnetometer, the method reconstructed the ambient magnetic field signal with an RMSE of 7.94 nT. This error is near the expected measurement noise for the PNI RM3100 magnetometer at 50 Hz, indicating that the accuracy of the method is limited to the total error budget of the magnetometer. Nevertheless, the experiments performed show the successful reconstruction of magnetic perturbation signals measured from within the bus of a mock CubeSat. These results demonstrate the utility of boomless CubeSats for scientific investigation of magnetic field phenomena in the geospace environment. In turn, the low cost of CubeSats enables the use of large constellations of small satellites to measure the geomagnetic field with high temporal and spatial resolution.

CONCLUSION

Based on experimental data and results set forth, the present teachings provide a method for separating spacecraft generated magnetic noise from geomagnetic field data using multiple magnetometers. The method does not require knowledge of the characteristics (location, orientation, amplitude, or spectral signature) and allows the number of noise sources to exceed the number of magnetometers (n>m). The method identifies signals by looking at the relative gain and phase of the magnetometer data in the TF domain. If a noise signal is sparse in this domain, the relative gain and phase is found using cluster analysis. Following the same assumption of sparsity, the signal can be separated from the noisy data using the cluster centroids in compressive sensing.

The method is designed for underdetermined systems in which there are more noise sources than magnetometers. An advantage of this approach is that the UBSS method can be integrated onto any satellite since no prior characterization of noise signals is required. This design eases the assimilation of magnetometers into spacecraft designs by reducing the need for strict magnetic cleanliness requirements and long mechanical booms.

There are several avenues of future development for this method. The most immediate step to be taken is for the selection of parameters to be automated. Our method automates the noise cancellation process, but some rudimentary analysis is still required to select parameters for clustering and preprocessing. It is disclosed that the selection of parameters could be entirely automated. Another avenue of development is to test the limits of the sparsity assumption. Sparsity is a very strict assumption that may not always be met. Finally, an interesting scenario to investigate is where several magnetometers are mounted within the bus of a spacecraft, but one magnetometer is mounted on a short boom, such as on the spacecraft Dellingr. In this scenario, the measurements of one magnetometer may be more accurate than the others. It would be counterproductive if the reconstructed magnetometer signal had more noise than the signal measured by the magnetometer on the boom. It may be possible to account for this by designing a programmable “trust” parameter at the compressive sensing stage. This parameter would indicate an elevated degree of trust in one magnetometer over the others.

We performed two types of experiments to validate the method. The first experiment separated SWARM magnetic perturbation data from four computer-simulated signals. The method was able to reconstruct the ambient magnetic field signal with an RMSE near 3 nT and a correlation of ρ≈0.9988. The reconstruction errors are less than the 6 nT intrinsic instrument noise that was added to each virtual magnetometer. The second experiment used real magnetic noise signals generated by copper coils, and the same SWARM geomagnetic field data. These experiments were able to separate four noise signals and reconstruct the background magnetic perturbation signal with a RMSE of 7.23 nT and a correlation of ρ=0.9958 at a 10 Hz cadence.

These results demonstrate that signal processing methods can identify and remove magnetic noise from spaceborne magnetometer data. The present method diminishes the need to place a magnetometer on a boom or enables significantly shorter booms. This enables the possibility of low cost, boomless spacecraft to capture high fidelity magnetic field measurements.

The resultant output of the present teachings provides robust data relating to low-amplitude and other magnetic fields that are of key importance for geophysical, magnetospheric, and heliospheric investigations. The output is useful for enhancing understanding of the aforementioned investigations and can be obtained and used to enhance navigations systems, detection systems, processing systems, and the like. For example, the system is useful for making scientific quality magnetic measurements from satellites without a boom. Moreover, in some embodiments, the resultant output and, particularly, the magnetic field data (output) can be used for both scientific purposes and navigation using magnetic compassing in addition to magnetic anomaly detection.

The foregoing description of the embodiments has been provided for purposes of illustration and description. It is not intended to be exhaustive or to limit the disclosure. Individual elements or features of a particular embodiment are generally not limited to that particular embodiment, but, where applicable, are interchangeable and can be used in a selected embodiment, even if not specifically shown or described. The same may also be varied in many ways. Such variations are not to be regarded as a departure from the disclosure, and all such modifications are intended to be included within the scope of the disclosure.

Claims

1. A system for separating spacecraft generated magnetic noise from geomagnetic field data, the system comprising: a plurality of magnetometers each configured to detect magnetic information and outputting a representative detection signal, the magnetic information comprising magnetic noise and geomagnetic field data; anda processing system configured to receive the representative detection signal including the magnetic information from the plurality of magnetometers and resolve the magnetic noise data and the geomagnetic field data by: analyzing a relative gain and phase of the representative detection signal in a time-frequency (TF) domain, if a noise signal of the representative detection signal is sparse in the TF domain then the relative gain and phase is found using cluster analysis, andseparating a geomagnetic field data signal from a magnetic noise data signal using cluster centroids in compressive sensing based on the cluster analysis,

2. The system according to claim 1 wherein the cluster analysis comprises a Density Based Spatial Clustering of Applications with Noise (DBSCAN) system by estimating a mixing matrix, K.

3. The system according to claim 2 wherein the compressive sensing is based on the mixing matrix, K.

4. The system according to claim 3 wherein the compressive sensing applies the mixing matrix, K, to demix the geomagnetic field data signal from the magnetic noise data signal using convex optimization.

5. The system according to claim 1 wherein the processing system is configured to receive the representative detection signal and resolve the magnetic noise data and the geomagnetic field data is performed for each axis of the plurality of magnetometers separately.

6. The system according to claim 1 wherein the processing system comprises a machine learning system.

7. A method for separating spacecraft generated magnetic noise from geomagnetic field data, the method comprising: detecting magnetic information using a plurality of magnetometers and outputting a representative detection signal, the magnetic information comprising magnetic noise and geomagnetic field data;receiving the representative detection signal including the magnetic information from the plurality of magnetometers and resolving the magnetic noise data and the geomagnetic field data by: analyzing a relative gain and phase of the representative detection signal in a time-frequency (TF) domain, if a noise signal of the representative detection signal is sparse in the TF domain then the relative gain and phase is found using cluster analysis, andseparating a geomagnetic field data signal from a magnetic noise data signal using cluster centroids in compressive sensing based on the cluster analysis; andoutputting a resultant signal representative of at least the geomagnetic field data signal.

8. The method according to claim 7 wherein the cluster analysis comprises applying Density Based Spatial Clustering of Applications with Noise (DBSCAN) method by estimating a mixing matrix, K.

9. The method according to claim 8 wherein the compressive sensing is based on the mixing matrix, K.

10. The method according to claim 8 wherein the compressive sensing applies the mixing matrix, K, to demix the geomagnetic field data signal from the magnetic noise data signal using convex optimization.

11. The method according to claim 7 wherein the receiving the representative detection signal including the magnetic information from the plurality of magnetometers and resolving the magnetic noise data and the geomagnetic field data is performed for each axis of the plurality of magnetometers separately.