# Co-Seismic Magnetic Field Perturbations Detected by Swarm Three-Satellite Constellation
## 1 Introduction
### Methods to Analyze the Magnetic Satellite Data  
A series of MATLAB(r) routines have been developed to analyze the Swarm Magnetic data. The first script reads the earthquake catalogue and selects the events in agreement with the magnitude-depth requirements described in Section 2.1 (all earthquakes with Mw \(\geq\) 5.5 and depth according to a linearly varying threshold from 20 km depth for Mw = 5.5 and Mw = 7.0 at 150 km depth). For each event, the algorithm selects the magnetic satellite data on the same date. The polar regions above a 70\({}^{\circ}\) geographical latitude are cut, as they are not suitable for this study due to the typical perturbation of the auroral activity. The goal of this work was to analyze the anomalous magnetic activity associated with earthquakes; thus, other magnetic perturbations, such as the geomagnetic storms, must first be rejected. After that, the script calculates the spatial distance of each sample from the epicentre and the time difference between the acquisition time and the earthquake origin time. If there are some samples close in space-time to the seismic event (distance \(\leq\) 150 km and time difference from 10 min before until 20 min after the earthquake), the script selects and runs an automatic analysis of the close track.  
Two types of automatic analyses have been performed: the first one is called Method 1, and it is the same proposed in previous works on Swarm data related to earthquakes (e.g., [22; 24]. This technique consists of computing the first derivative of the data along the track approximated as the first difference divided for the time interval between two consecutive samples. After this step, a cubic spline is removed and the residuals of X, Y, Z (geomagnetic components) and F (intensity) are plotted with a fixed maximum range of \(\pm\) 1 nT/s.  
The second method (so-called Method 2) has been developed to avoid a possible over-fitting problem of the first method, i.e., that the spline could remove the eventual perturbation induced by the supposed gravity wave. It consists of the subtraction of the International Geomagnetic Reference Field (IGRF-12) model [29] from the original data and the plotting of the residuals in X, Y, Z and F. The IGRF-12 is a global Earth's magnetic field model, which is expected to be sufficient to highlight the studied disturbances. A five-degree polynomial has been further applied to remove the residual trend of each track.  
The output figures of both methods are presented in the Supplementary Materials. Each figure presents a map of the geographical region crossed by the satellite, the projection of the satellite track on the ground (brown for normal conditions and red for some abnormal conditions on the satellite) and the earthquake epicenter (green star). The title of the figure provides the satellites, the time in UTC and Local Time and the geomagnetic indices Dst and ap during the acquisition time. The residuals of X, Y, Z and F are reported together with the Fast Fourier Transform analysis of the signal that helps to select the best cut-off frequency for the following signal processing (for the second method, the mean of the signal is subtracted before calculating the FFT). The figures are grouped in tables that refer to each seismic event.  
All the satellite data of the selected track are stored in a text file for further analysis. The produced material is visually investigated and the tracks that present some possible features linkable to the seismic events are then re-analyzed and plotted again. For this final analysis (presented in Figures 3-14), only the Y component is considered: after the subtraction of IGRF-12 model, a high pass filter at 0.25 Hzwas applied and both quantities were then plotted. The frequency cut-off was chosen after multiple tests as the best value to extract the anomalous signal possibly related to the seismic event. Taking into account that the sampling rate of the data was 1 Hz, the available band was 0.25 Hz-0.50 Hz (following Nyquist theorem).