A new method of measuring Forbush decreases

Forbush decreases (FDs) are short-term depressions in the galactic cosmic ray flux and one of the common signatures of coronal mass ejections (CMEs) in the heliosphere. They often show a two-step profile, the second one associated with the CMEs magnetic structure (flux rope, FR), which can be described by the recently developed model ForbMod. The aim of this study is to utilise ForbMod to develop a best-fit procedure to be applied on FR-related FDs as a convenient measurement tool. We develop a best-fit procedure that can be applied to a data series from an arbitrary detector. Thus, the basic procedure facilitates measurement estimation of the magnitude of the FR-related FD, with the possibility of being adapted for the energy response of a specific detector for a more advanced analysis. The non-linear fitting was performed by calculating all possible ForbMod curves constrained within the FR borders to the designated dataset and minimising the mean square error (MSE). In order to evaluate the performance of the ForbMod best-fit procedure, we used synthetic measurements produced by calculating the theoretical ForbMod curve for a specific example CME and then applying various effects to the data to mimic the imperfection of the real measurements. We also tested the ForbMod best-fit function on the real data, measured by detector F of the SOHO-EPHIN instrument on a sample containing 30 events, all of which have a distinct FD corresponding to the CMEs magnetic structure. Overall, we find that the ForbMod best-fit procedure performs similar to the traditional algorithm-based observational method, but with slightly smaller values for the FD amplitude, as it is taking into account the noise in the data. Furthermore, we find that the best-fit procedure has an advantage compared to the traditional method as it can estimate the FD amplitude even when there is a data gap at the onset of the FD.


Introduction
Forbush decreases (FDs) are short-term depressions in the galactic cosmic ray (GCR) flux that result from the interaction of GCRs with coronal mass ejections (CMEs).They are thus one of the common signatures of CMEs in the heliosphere (see overviews by e.g.Lockwood 1971;Cane 2000;Belov 2009).Detecting CMEs at various different locations in the heliosphere is crucial for the analysis of their evolution, as well as the validation of the CME propagation and evolution models.They are typically observed by in situ detectors that measure magnetic field and plasma properties (e.g.Zurbuchen & Richardson 2006;Kilpua et al. 2017); however, these are not available as frequently as particle detectors that can detect FDs.Therefore, a better understanding of FDs can lead to CME analysis where other types of observation are not possible.
FDs caused by interplanetary coronal mass ejections (ICMEs) often show a two-step profile, one associated with the turbulent and compressed shock-sheath region ahead of the CME magnetic structure and the other with the CME magnetic structure itself (Cane 2000).Because of their different physical properties, the mechanism through which these two regions produce the depression is different, and thus should be modelled differently.The current paradigm is that the CME magnetic structure is a flux rope (FR), a structure with field lines helicoidally winding around a central axis (e.g.Zhang et al. 2021).It was first proposed by Morrison (1956) that FDs might be caused by clouds of closed magnetised plasma, whereby they are initially empty of particles and slowly fill as they propagate through the interplanetary space.This approach has since been utilised to explain FR-associated FDs, where the particles are allowed to enter via "cross-field" transport; in other words, perpendicular diffusion (e.g.Cane et al. 1995;Quenby et al. 2008).In addition to the particle diffusion, another important contribution comes from the expansion of the FR.Namely, as they propagate through the IP space, CMEs expand due to the pressure imbalance (e.g.Klein & Burlaga 1982;Démoulin & Dasso 2009), where, consequently, as the size of the magnetic structure increases, its magnetic field weakens (Bothmer & Schwenn 1998;Démoulin et al. 2008;Vršnak et al. 2019).It was first proposed by Laster et al. (1962) r FD [%] 0 A FD,max 0 r/a 1 -1 Fig. 1.Sketch of the ForbMod model for FR-related FDs.The GCRs diffuse into the self-similarly expanding FR (r/a = const.),which is initially empty at its centre (r/a = 0), where a(t) is the FR radius and r(t) the radial distance from the centre of the FR to its borders.As a result, at an arbitrary time after the FR eruption from the Sun (i.e. at an arbitrary heliocentric distance), a distribution of GCRs inside the FR is observed, given by the Bessel function.The relative GCR density (FD [%]) at the borders of the FR is 0%, whereas in the centre of the FR it is given by the amplitude of the depression (A FD,max ).For model details, see the main text.and Singer et al. (1962) that expansion is needed to explain the observational properties of FDs.More recently, the diffusionexpansion approach was applied by Munakata et al. (2006) and Kuwabara et al. (2009) in a numerical model and by Subramanian et al. (2009) and Arunbabu et al. (2013) in an analytical model, considering isotropic self-similar expansion.In the recent FD analytical model ForbMod (Dumbović et al. 2018), different types of expansion are considered for the expansion-diffusion approach, constrained by the observational ranges of the expansion factors (Gulisano et al. 2012).ForbMod was advanced to take into account the energy dependence of the detector when comparing the model results to measurements (Dumbović et al. 2020), since FDs are energy-dependent phenomena (i.e.their amplitudes depend on the energy range and cut-off of the detector which detects them Cane 2000).This facilitates direct comparison of modelled and observed FDs, as was shown by, for example, Freiherr von Forstner et al. (2021).
It should be noted that there are other approaches to FD modelling in the FR, such as perpendicular transport by guiding centre drifts (e.g.Krittinatham & Ruffolo 2009;Tortermpun et al. 2018), full trajectory integration using CME FR-type models (Petukhova et al. 2019), or CME magnetic field reconstructions from in situ measurements (Benella et al. 2020), as well as describing FDs via the change in the single GCR spectrum modulation parameter attributed with a CME (Guo et al. 2020).However, a recent study by Laitinen & Dalla (2021), in which they performed full-orbit particle simulations, with the interface between the external interplanetary magnetic field and the isolated FR field lines modelled analytically, revealed that the transport around the FR border is fast compared to diffusive radial propagation within the FR.As a result, the propagation of GCRs into the FR can be modelled as diffusion into a cylinder.Moreover, the diffusive approach is further supported by most recent observations (see e.g.recent studies by Janvier et al. 2021;Davies et al. 2023).
The aim of this study is to utilise ForbMod to develop a best-fit procedure that can be applied to FR-related FDs as a convenient measurement tool.Our motivation is to develop a best-fit procedure that can be applied to a data series from an arbitrary detector.Thus, the basic procedure would facilitate measurement estimation of the onset, duration, and magnitude of the FR-related FD, with the possibility of being adapted for the energy response of a specific detector for a more advanced analysis.It should be noted that in many cases the onset of the FR-related FD does not necessarily correspond to the onset of the FD itself, as the latter might correspond to the preceding shock-sheath region.In such cases we refer to the onset of the second decrease as the onset of FD (or FR-related FD).In Sect. 2 we describe the mathematical formalism of the best-fit procedure as well as two samples used in the study to test its performance.We first used synthetic measurements to test the performance of the best-fit procedure under controlled conditions.Next we used a real sample of events detected by the F-detector of the Electron Proton Helium INstrument instrument aboard the SOHO spacecraft (SOHO/EPHIN Müller-Mellin et al. 1995).The results are presented and discussed in Sect.3.

ForbMod best-fit procedure
ForbMod is an analytical diffusion-expansion model that assumes that the CME magnetic structure (i.e.FR) can locally be represented by a straight cylinder, which is initially empty at its centre.The FR is assumed to be propagating at a constant speed, expanding self-similarly, and the particles can only enter it via perpendicular diffusion.It is assumed that GCRs are mostly protons.The time-dependent transport equation with these assumptions simplifies to a time-dependent radial diffusion of protons into an infinite straight cylinder.The corresponding solution is a combination of the radial-dependent solution, which is given by the cylindrical Bessel function of the zero-order, and the time-dependent solution, which is given by the exponential function (for details on the model see Dumbović et al. 2018).The amplitude of the FD, defined as normalised to border conditions, depends on the diffusion rate of particles into the FR.The diffusion rate is governed by the behaviour of the magnetic field, which is generally decreasing due to expansion.Thus, the FD amplitude is expected to decrease with diffusion time (i.e.helisopheric distance).The sketch visualising ForbMod is given in Fig. 1.The diffusion rate of particles also depends on the particle energy, as higher-energy particles enter the FR more easily.Therefore, it is expected that the FD amplitude is energydependent and that it should appear differently in detectors with a different energy response (for details see Dumbović et al. 2018Dumbović et al. , 2020)).
We developed the ForbMod best-fit procedure by applying a non-linear fitting to the data.For that purpose we first used A168, page 2 of 13 synthetic measurements, as is described in Sect.2.2.In general, the non-linear fitting was performed by calculating all possible ForbMod curves constrained within FR borders to the designated dataset and minimising the mean square error, MSE: where ŷi is the predicted value of the ith sample, y i is the corresponding true value, and there are n SAMPLES samples.The ForbMod best-fit procedure actually performs non-linear fitting of the Bessel function of the zero order (i.e. the radial part of the ForbMod solution) to the designated time series.The input for the procedure is the designated dataset with defined FR borders (i.e. the dataset translated into the FR coordinate system radial direction, with the centre of the FR at r = 0 and FR borders at r = −1 and r = 1, respectively; r = r/a, where a is the radius of the FR).The output of the procedure is the best-fitted amplitude of the measured FD, A FR,max , and the corresponding goodness of fit, that is, MSE (as defined by Eq. ( 1)).The FD profiles needed to be translated from a timescale into the spatial scale of the FR in order for us to apply the ForbMod best-fit.For that purpose, we assumed that the spacecraft crosses the full diameter of the FR, which has a circular cross-section.Furthermore, we assumed that the spacecraft trajectory is oriented vertically with respect to the FR axis.We note that with this assumption the model might not realistically capture FRs with orientations deviating significantly from low or high inclinations.Nevertheless, with this assumption, due to self-similar expansion, the spacecraft observes deceleration of plasma flow across the FR.Thus, the diameter can be given by the equation for accelerated linear motion, s = 0.5at 2 +v 0 t.The borders of the FR were set at [-1,1], per the definition of the ForbMod best-fit function.

Synthetic Forbush decrease measurements
In order to evaluate the performance of the ForbMod best-fit procedure, we used synthetic measurements.Synthetic measurements were produced by calculating the theoretical ForbMod curve for a specific example CME and then applying various effects to the data to mimic the "imperfection" of the real measurements.
The arbitrary CME example used to produce synthetic measurements (loosely based on event analysed by Freiherr von Forstner et al. 2021) had the following inputs: -detector response function for perfect detector: Y(E) = 1 -detector cutoff: E cutoff = 0.05 GeV -event month and year: April 2020 -ICME transit time in hours (i.e.diffusion time), T T = 123.38h -starting distance from the Sun, R 0 = 9.54 R SUN -final distance from the Sun, R fin = 1.005 au -initial FR radius, a 0 = 1.81 R SUN -expansion factor for FR radius, n a = 1 -in situ magnetic field strength, B = 16.2 nT -expansion factor for mag.field, n B = 2.The effects applied to the data to mimic the "imperfection" of the real measurements were: -Y-noise: -random noise in the y direction -X-noise: -random noise in the x direction -asymmetry: -systematic left-shift of the FD minimum in the x direction (to simulate the asymmetric FD profile) -data gaps: -simulating a data gap by removing certain data points around the selected position within the FD -low data resolution: -reducing the number of data points in the FD profile.The Y-noise value was calculated and added to each data point of the ForbMod solution in the following way: where y i,SYNTH is the synthetic measurement, y i,ForbMod is the ForbMod solution data point, δy i is the noise, and A norm is the randomly selected value from the normal distribution defined by the ForbMod solution amplitude as the mean.A norm was determined using the Python numpy package function random.normal.δ n was an arbitrary percentage applied to the randomly selected A norm .The parameter δ n could be varied in order to obtain different percentages (i.e.different levels) of noise.The X-noise value was calculated and applied in the same manner as the Y-noise, only in the x direction.The Y-noise and X-noise of different percentages are shown in the upper and middle panels of Fig. 2, respectively.The asymmetry was introduced by incrementally shifting positive r values corresponding to data points of the ForbMod solution: where r i is the positive r value corresponding to the ith data point of the ForbMod solution (counted for positive r values only) and δ a is a shift percentage between two iterations, which could be varied in order to analyse different percentages (i.e.levels) of asymmetry.After applying the asymmetry, the dataset was rescaled to be defined on a [-1,1] interval using the sklearn package function minmax_scale.The asymmetry of different percentages is shown in the lower panels of Fig. 2. The data gaps were introduced by removing a certain number of data points from a certain location within the FD profile.Thus, we could vary both the location of the data gap, as well as its size.This is shown in the example in Fig. 3.In the upper panels, a gap 10% in size (i.e.200 data points out of 2000) is consecutively shifted from the start of the FD to its minimum (five different positions).In the middle panels of Fig. 3 various gap sizes are shown, where the gap is centred around the middle of the FD main phase (i.e. the decrease phase).Finally, as the last effect applied to the data to mimic the "imperfection" of the real measurements, we considered data resolution, as is shown in the bottom panels of Fig. 3.

Event sample and measurements
A sample of events was created by associating SOHO/EPHIN detected FDs with in situ detected ICMEs and remotely observed CMEs.We used a list of SOHO/EPHIN detected FDs compiled by Belov et al. (2021) and a catalogue of near-Earth ICMEs compiled by Richardson & Cane (2010) and regularly updated at CALTECH 1 .Since we aimed to use the sample for the followup analysis, where a 3D reconstruction using stereoscopic data to obtain CME remote properties would be performed, we constrained our analysis to the time period when at least one STEREO spacecraft was operational, where a reliable CME-ICME association could be made, and where CME signatures in white light coronagraphs were significant enough to perform a 3D reconstruction.The sample contains 30 events in the time period 2007-2019.In this study, we only used the part of the 1 http://www.srl.caltech.edu/ACE/ASC/DATA/level3/icmetable2.htmA168, page 3 of 13 Dumbović, M., et al.: A&A, 683, A168 (2024 asymmetry 0.25% Fig. 2. Effects applied to the data to mimic the "imperfection" of the real measurements: Y-noise (upper panels), X-noise (middle panels), and asymmetry (lower panels) for three different levels of noise and asymmetry, respectively.The original theoretical curve is given in black, whereas the data with the applied effect are coloured blue.For details, see the main text.
sample containing the interplanetary in situ signatures of ICMEs and the corresponding FD, and thus we only provided a description of the ICME-FD sample in detail, whereas the description of the CME-part of the sample will be given in a future study.
Interplanetary coronal mass ejections are roughly simultaneous with FD events (Cane et al. 1996;Dumbović et al. 2011), and therefore this association was trivial and performed based on the timing provided by Belov et al. (2021) or Richardson & Cane (2010)'s catalogues.To analyse the in situ properties, we used one-minute plasma and magnetic field data in the geocentric solar ecliptic (GSE) system provided by the OMNIWeb 2 database (King & Papitashvili 2005).We identified ICME signatures using criteria described by (e.g.Zurbuchen & Richardson 2006;Kilpua et al. 2017), considering the notification of Rouillard (2011), according to which an ICME consists of the compressed, heated, and turbulent shocksheath region, followed by a magnetic cloud or ejecta.The magnetic cloud is characterised by a low-temperature, low-plasma beta parameter, an enhanced smoothly rotating magnetic field, 2 https://omniweb.gsfc.nasa.gov/html/ow_data.html#norm_pla and an expanding flow speed profile, whereas the term magnetic ejecta is used when these criteria are only partially met.We note that both magnetic cloud and ejecta are nowadays broadly considered as magnetic obstacles (MO, Jian et al. 2006;Nieves-Chinchilla et al. 2018), and therefore we use this more general term.
Example events are given in Fig. 4. For each event we analysed: the magnetic field strength, B, and its fluctuations, dBrms 3 ; the magnetic field components, B x , B y , and B z ; the plasma density, temperature, and expected temperature 4 ; the plasma flow speed and plasma beta parameter; and the SOHO/EPHIN F-detector particle count, CR count 5 .We used hourly averaged detector F data of the SOHO/EPHIN, as they were shown to be suitable for a FD analysis, especially the MO part (Heber et al. 2015;Dumbović et al. 2020;Belov et al. 2021).The SOHO/EPHIN data were time-shifted to the Earth's bow 3 The root mean square variation in the vector of the interplanetary magnetic field. 4Empirical formula based on the speed of the solar wind calculated according to Lopez (1987) and Richardson & Cane (1995). 5Normalised to the pre-decrease value for each event.
A168, page 4 of 13 Dumbović, M., et al.: A&A, 683, A168 (2024 20 data points Fig. 3. Effects applied to the data to mimic the "imperfection" of the real measurements: data gaps at certain positions (upper panels), data gaps of a certain size (middle panels), and different data resolutions (lower panels).The original theoretical curve is given in black, whereas the data with the applied effect are coloured blue.For details, see the main text.
shock using the OMNIWeb formula for time-shift: where X and Y are GSE X and Y components of the spacecraft position vector (we used daily positions of SOHO), V is the observed solar wind speed (assumed to be radial), V e is the speed of the Earth's orbital motion (30 km s −1 ), and W = tan(0.5atan(V/428)) is the geometry factor (for details, see the OMNIWeb webpage 6 .In order to achieve better visibility of smaller FD variations, SEP events were removed from the data using the <600 counts threshold on the B-detector, which is sensitive to lower energy particles, that is, sharp increases due to SEPs.We visually inspected each event and discarded events where either no MO or corresponding FD was observed.In addition, we classified events according to the quality index (QI), subjectively determined by the observer, in the range 1-4, with 1 being the worst quality and 4 being the best quality (see Fig. 4).The QI is based 6 https://omniweb.gsfc.nasa.gov/html/ow_data.html on criteria regarding the clarity of MO signatures and borders, as well as the clarity of FD signatures and borders.Both are less clear for events designated as having a lower QI.We note that in QI determination, some more weight was given to the clarity of an onset corresponding to both MO and FD, than to one corresponding to other MO/FD properties, due to its importance for the best-fit application.We also note that although we analysed the clarity of MO borders visually, we did not actually determine or measure MO borders, but FD borders.
It can be seen in Fig. 4 that the MO signatures are weakest for the QI1 event (bottom right).The magnetic field is increased and the B z is rotating, although the magnetic field is not very smooth.Plasma beta is low, although the temperature is only partially lower than expected and the expansing profile of the flow speed is missing.The onset of the FD is clear; however, the recovery of the FD seems to be missing.This is likely related to another structure that can be observed trailing the MO.In the QI2 event (bottom left) the MO signatures are clearer.The magnetic field is again increased and the B y is rotating smoothly.Plasma beta is low, but the temperature is again only partially lower than expected.The flow speed seems to show a globally expanding A168, page 5 of 13 Fig. 4. In situ measurements of four different events assigned with four different quality indices, QI4-QI1, from top left to bottom right.Top-tobottom panels for each event show: (1) magnetic field strength (black, in nT) and its fluctuations (grey, in nT); (2) magnetic field GSE components, x, y, and z, coloured red, blue, and green, respectively (in nT); (3) plasma density (black, in cm −3 ), temperature (red, in 1e5 K), and expected temperature (blue, in 1e5 K); (4) plasma flow speed (black, in km s −1 ) and beta (grey, non-dimensional); (5) the SOHO/EPHIN F-detector particle count (in %).The shaded region represents the magnetic obstacle, where the borders were selected according to the low-beta condition.
profile, but with local irregularities.The corresponding FD has a clear onset; however, the onset is slightly shifted in time with respect to the MO onset (according to the low-beta condition).Moreover, the FD shows a double depression without the full recovery to the pre-decrease level.We note that the substructuring of the FD corresponds to the substructuring of the MO.The first depression seems to correspond to the rotation of the Bz component from positive to negative values and a small density increase, whereas the second depression seems to correspond to the rotation of the By component from negative to positive values and the low density region.This substructuring might indicate a complex inner structure of the MO, or the passage of the spacecraft under some oblique angle.However, although they both show substructuring, it appears the global structure may correspond to that of a single MO, as was found for example by Nieves-Chinchilla et al. (2018).The QI3 event shows clear MO signatures as, in addition to other properties that the QI2 event shows, the flow speed in the QI3 event shows a clearly expanding profile.Although a depression can be seen to start with the shock-sheath, the second depression of the FD has a clear onset that corresponds to the MO with a clear tendency to return to the pre-decrease level, although it does not (most likely due to the high-speed stream that follows).Moreover, the FD seems to have different recovery rate at the back (which again corresponds to the substructuring of the MO, similar to the QI2 event, but much less pronounced).Finally, the QI4 event shows clear MO A168, page 6 of 13 signatures and a clear, symmetric FD with borders that roughly match the FR borders.However, we do note that although globally this event shows clear MO signatures, again some substructuring can be observed.To summarise, in general the best-quality events have nicely observed FD and MO signatures and clear borders.The worst-quality events have unclear borders, pronounced substructures, and/or an unclear recovery.Out of the 30 events in the sample, 6, 12, 8, and 4 have QIs 1, 2, 3, and 4, respectively.The full sample of events is given in Table A.1.

Extraction of FD profiles: The border selection problem
In order to apply the best-fit procedure, the FD profiles needed to be extracted, that is, the borders of the FD needed to be selected.The onset of FD is usually defined as the point where the GCR count starts to drop, whereas the end of FD, in a textbook example, would be the point where it returns to its pre-decrease level (e.g.Dumbović et al. 2011).In textbook examples, the recovery of two-step FD is prolonged after the passage of the ICME due to the "shading effect" of the shock-sheath region, whereas for a FR, FD (which can be self-standing or a second decrease in a two-step FD), the end of FD corresponds to the end of the FR (Dumbović et al. 2020).However, in many non-textbook cases, the determination of the onset and end of FD is not a trivial procedure that can be defined by a simple algorithm applicable evenly to all events, as can clearly be seen from Fig. 4. Therefore, the extraction of FD profiles (from the onset to the end) was performed manually by an observer.This may result in a certain subjectivity when determining the borders, especially the end of FD.Therefore, we applied two different versions of border selection, as is shown in Fig. 5: -extended ICME borders -these encompass typical ICME signatures such as, for example, low temperature or low plasma beta, but with borders set based on the GCR depression onset or end and/or based on the total magnetic field increase onset or end.They thus overlap with typical MO signatures only in a subset of events and they often contain external disturbed parts of the FR.
-inner structure borders -we find that many events show clearer MO signatures in a sub-interval of what would normally be selected as the global MO structure according to, for example, the low-beta condition.Therefore, first we selected borders to encompass this inner substructure of the CME, which shows clear MO signatures.However, we did this only for events where corresponding clear FD substructure borders could also be found.
We find that the inner structure can be identified in 24 events.The inner structure is usually small compared to the typical ICME duration, which is roughly 24 h (Richardson & Cane 2010).Consequently, these events have fewer data points to be fitted, but also a more compact profile, without substructures.We did not further refine the QI based on the different border selection, but instead treated the border selection and QI as independent classifications.In some events the inner structure could not be identified, whereas in some other events the whole duration of the event could be considered as the inner structure.We found four such events, which we in turn regard as not having extended ICME borders.The best-fit procedure is applied to the GCR measurements within the selected borders, and thus slightly different results of best-fit FD amplitude, FD_bf, can be expected for different border selections.In addition, for each event we determined the FD amplitude using a traditional algorithm-based observational method, FD_obs, whereby FD was measured from the onset to the minimum.Therefore, this observed FD, FD_obs, may also be slightly different for different border selections.The borders corresponding to the inner structure and extended ICME are given in Table A.1 for each event, along with the corresponding observed and best-fitted FD.

ForbMod best-fit function testing through synthetic measurements
We applied the ForbMod best-fit procedure to the synthetic data in order to check how the goodness of fit (MSE) behaves for each of the different effects, but also for different levels of a certain effect.For each best-fit to the synthetic data we calculated MSE and the relative difference of the best-fit curve minimum to the minimum of the theoretical curve (i.e. the relative difference of the best-fit FD amplitude to the FD amplitude of the theoretical FD, dA/A).We first applied the ForbMod best-fit procedure to synthetic data with different levels of Y-noise and X-noise.
The results are shown in Figs.6a and b.We can see that, for the same level of Y-noise, dA/A and MSE are correlated, whereas for the same level of X-noise the dA/A remains constant.Therefore, by increasing the X-noise we increase both MSE and dA/A.However, by increasing the Y-noise we only increase MSE, and dA/A does not change drastically.This is related to the fact that Y-noise is applied in the same direction in which the MSE is calculated.Simply put, the random noise up and down around the curve will contribute to the MSE, but the average curve will stay the same as the original one.We next applied the ForbMod best-fit procedure to synthetic data with different levels of X-noise and asymmetry.This is A168, page 7 of 13 shown in Fig. 6c.We can see that for lower values of asymmetry (<0.15%) asymmetry works in the same direction as the X-noise: the larger the MSE, the larger the dA/A.However, for large asymmetry there seems to be a saturation of dA/A, which remains scattered around 60% regardless of the MSE.We note that combining three different data imperfection effects (Y-noise, X-noise, and asymmetry) very quickly leads to large dA/A, as the data does not resemble the FD anymore.Therefore, it is already impossible to apply the best-fit procedure to data with Y-noise, X-noise, and asymmetry greater than 15%, 15%, and 0.15%, respectively.We next applied the data gap to the synthetic measurements and repeated the best-fit procedure.We applied a data gap 10% in size to five different positions, starting from the onset to the minimum, as is described in Sect.2.2.The results are shown in Fig. 6d), where we find that for a data gap position closer to the FD minimum, the increase in MSE results in an increase in dA/A greater than for data gaps closer to the onset of FD.In other words, we might expect to fit worse when our dataset has a gap close to the FD minimum.We also checked how the data size influences the best-fit procedure (Fig. 6e).We find that for gap sizes of up to roughly 35% the scatter of the MSE for the same dA/A increases -depending on where the data gap is and what other noise the data has, the data gap of the same size can lead to a range of different dA/A, and this range increases with the gap size.However, for gap sizes >35% we observe the opposite trend.This is related to the fact that with an increase in the data gap we removed data and with it the other data imperfection effects.
Finally, we applied different combinations of data imperfection effects using three different resolutions of the original theoretical curve data (2000, 200, and 20 data points).In Fig. 6f we can see that for smaller resolutions we expect a larger scatter of MSE for the same dA/A.
In Fig. 7 we show the relative difference of the best-fit FD amplitude to the FD amplitude of the theoretical FD versus the MSE for a mixture of different levels of data imperfection effects.Combinations of different imperfection effects can quickly lead to very high MSE and dA/A, and therefore we constrained combinations to the following levels: Y-noise = 15%, X-noise = 15%, asym = 0.15%, and gap size = 30%.However, we allowed parameters to go beyond these levels when one or more of the data imperfection parameters was kept at level 0. We can see in Fig. 7 that the linear fit is not far away from the identity line (y = 100 * x); however, the scatter of the data is quite large.
The linear fit moves even closer to the identity line when we refine the sample by removing extremely outlying points, that is, those that have dA/A > 20% away from the identity line and extremely low or high corresponding MSE, that is, <0.2/>0.3,respectively.This is shown in the right plot of Fig. 7. Therefore, we next analysed visually the synthetic data and the best-fit function, as well as the MSE and dA/A for cases which are that to the identity line and those far away.This is shown for nine selected examples in Fig. 8.The examples are given in the nine panels on the left of the Figure, whereas the right-most plots show the scatter plots of dA/A vs. MSE for three different sets of randomly produced synthetic data.The scatter plots differ slightly due to randomness of producing synthetic data, but in general show similar behaviour.The three examples shown left of the corresponding scatter plot are highlighted in a colour corresponding to the colour of the synthetic data for the given event (green for the left-most, yellow for the middle, and red for the right-most examples).For synthetic data of the same colour, top-to-bottom panels show different levels of data imperfection effects, ranging from highest (upper plots) to lowest (bottom plots).
identifiable FD, whereas in the bottom plots the data depict an almost ideal symmetric FD profile.This visual representation is in principal reflected by the MSE -lower plots have much smaller MSE compared to the upper plots.The middle plots present intermediate cases between extreme examples.Looking at the difference between the theoretical and best-fit curve, dA/A, we can see that it is negligible for the bottom left-most plot, whereas in the upper left-most plot it almost reaches 10%.However, for the top-to-bottom right-most examples it is reversed -dA/A is small for the top right example, whereas in the bottom right-most example it almost reaches 10%.This reflects the previously observed fact that for some data imperfection effects MSE is not correlated with dA/A.These are the cases that are far away from the identity line.We note that in the example of low levels of the data imperfection effects (bottom plots) MSE is low (<0.03) in all three cases and dA/A is within 10%.Therefore, we do not regard any of the bottom plots cases as being far away from the identity line.
From this analysis we see that, if we exclude events for which visual representation is not clear enough to identify FD, we are left with events close to the identity line.If we further restrict ourselves only to events with MSE < 0.1, we are left with events with dA/A <10%.Therefore we can conclude that MSE < 0.1 can serve as a goodness-of-fit measure resulting in <10% error of the fitted amplitude under the following assumptions: (1) the real data behave as synthetic data; (2) the data in general depict the behaviour of the theoretical ForbMod FD (with some levels of data noise); (3) FD is clearly visible in the data and not masked by noise.

ForbMod best-fit function application to SOHO/EPHIN measurements
We next performed the ForbMod best fit on the sample of the events, as is described in Sect.2.3.We first analysed the fits visually and found we could divide the fits into six types: -a good fit -where there was visually a good match between data points and the best-fit function -an asymmetric fit -where the fitted function was shifted along the x axis compared to the data points -a few-data-points fit -where the fit was performed on only a small number of data points -a large-scatter fit -where the fit was performed on data points with a large scatter -an under-recovery fit -where the data points did not return to the pre-decrease level -a substructuring fit -where a single fit was performed while the data points indicated possible substructuring Examples of fit types are given in Figs.9a-f.However, we note that some events may fall into more than one category; for example, an event might show under-recovery as well as substructuring.As can be seen in Table A.1 all but one good fit events have inner structure borders.They also on average have the smallest MSE (3.9e-5).A few-data-points fit is another type present only for inner structure borders events.On the other hand, the largescatter fit and the substructuring fit seem to be categories present predominantly for extended ICME borders and on average have the largest MSEs (1.8e-4 and 1.5e-4, respectively).We note that for all of the events MSE < 0.01.Therefore, although some profiles do not show a visually pleasing FD, the ForbMod best-fit function still manages to find a solution with a relatively small mean standard error (MSE).
We next analysed best-fit results for four different QIs, as is described in Sect.2.3 (see Table 1).We do not find notable differences between events marked by different QIs, in either the difference between FD_bf and FD_obs, the values of MSE, or the spread of MSE.Surprisingly, we find that MSE on average takes the smallest values for QI 1, i.e. the least clear events.This is likely related to the fact that these events also on average have the smallest FD amplitudes, both FD_bf and FD_obs (1.1% and 1.2%, respectively).
A168, page 10 of 13 We next compared the calculated best-fit amplitude, FD_bf, with the amplitude obtained using a traditional algorithm-based observational method, FD_obs, whereby FD_obs was measured from the onset to the minimum.For the inner structure borders, we find the average FD_bf of 2.8%, which is somewhat smaller than the average FD_obs (3.1%).We note that this difference is expected, because the traditional algorithm-based observational method measures maximum amplitude, without taking into account possible scatter of the data.For extended ICME borders we find an even larger difference between the two, with an average FD_bf of 2.5% and an average FD_obs of 3.5%.This reflects the fact that the scatter of the data that was fitted is on average larger for extended ICME borders compared to inner structure borders, as can be seen from the average values of MSE being 1.1e-4 and 4.2e-5, respectively.
We also performed a linear regression and calculated the Pearsons correlation coefficient for FD_bf versus FD_obs separately for inner structure borders and extended ICME borders.The results are presented in Fig. 10 and confirm the findings above.We observe a negative intercept in both linear regressions, confirming a systematically smaller FD_bf compared to FD_obs.The larger negative intercept for extended ICME borders indicates a larger difference between FD_bf and FD_obs.In addition, we can see that the Pearsons correlation coefficient (i.e. its squared value given by R 2 ) is smaller for extended ICME borders related to the larger scatter.However, we do note that the correlation is high for both border selection cases (R > 0.8).Thus, we can conclude that the ForbMod best-fit procedure performs similarly to the traditional algorithm-based observational method, but with slightly smaller values for the FD amplitude, as it's taking into account the noise in the data.We note that in two events FD_obs could not be determined because of the data gap at the onset of the FD (see Table A.1).In those cases, FD_bf could still be determined, which is an obvious advantage of the best-fit procedure.

Summary and conclusions
We developed a new method to measure the FD amplitude, which is based on the analytical diffusion-expansion FD model for FRs (ForbMod).The ForbMod best-fit procedure we developed performs a non-linear fitting of the Bessel function of the zero order (i.e. the radial part of the ForbMod solution) to the designated time series.The non-linear fitting was performed by calculating all possible ForbMod curves constrained within FR borders to the designated dataset and minimising the MSE.The input for the procedure is the designated dataset with defined FR borders (i.e. the dataset translated into the FR coordinate system radial direction, with the centre of the FR at 0 and FR borders at ±1, respectively).The output of the procedure is the best-fitted amplitude of the measured FD and MSE.
In order to evaluate the performance of the ForbMod best-fit procedure, we used synthetic measurements.Synthetic measurements were produced by calculating the theoretical ForbMod curve for a specific example CME and then applying various effects to the data to mimic the imperfection of the real measurements.Using synthetic measurements we could test how far away the best-fitted function is from the real one for different cases.We found that if the FD is clearly visible in the data the reliable ForbMod best-fit function should result in the MSE being <0.1 with an expected <10% relative error in the fitted amplitude, roughly along the identity line.
We next tested the ForbMod best-fit function on the real data, measured by detector F of the SOHO/EPHIN instrument.We used a list of SOHO/EPHIN detected FDs from the time period 1995-2015 and a catalogue of near-Earth ICMEs to compile a list of events suitable for analysis.The sample contains 30 events, all of which have a distinct FD corresponding to the magnetic obstacle (MO).We categorised the events according to different QIs, which mark different levels of clarity of the events.We extracted FD profiles for each event in the sample from the onset to the end of FD and prepared them for the application of the ForbMod bestfit function.The extraction of FD profiles (from the onset to the end) was performed manually by an observer, whereby we applied two different versions of border selection: (1) the inner substructure of the CME, which shows clear MO signatures, and (2) the extended ICME structure, which may also encompass outer, disturbed layers of the ICME.
We do not find notable differences between events marked by different QIs.We compared the calculated best-fit amplitude with the amplitude obtained using a traditional algorithm-based observational method (measured from the onset to the minimum) and found a good match between the two.The ForbMod best-fit procedure performs similarly to the traditional algorithm-based observational method, but with slightly smaller values for the FD amplitude, as it takes into account the noise in the data.For events with two different versions of border selection we find that the best-fit applied on extended ICME structure borders results in slightly larger MSE and differences compared to the traditional method due to the larger scatter of the data points.We find that the best-fit results can visually be categorised into six different FD profile types.Although some profiles do not show a visually pleasing FD, the ForbMod best-fit function still manages to find a solution with a relatively small MSE.
Finally, we find that the best-fit procedure has an advantage compared to the traditional method as it can estimate the FD amplitude even when there is a data gap at the onset of the FD.The method could be further advanced by taking into account best-fitting for multiple-border selection, and thus reducing the subjective impact of the observer.With further improvements in the direction of border selection, the method would not only advance the estimation of the magnitude of the FD, but also the determination of its onset, end, and consequently duration.Such an application would provide the means for a quick, objective, reliable, and robust way of measuring FDs and may eventually facilitate the automation of FD measurements.The ICME start date corresponds to the first detected ICME signature, which may be shock, sheath, pile-up, or MO (depending on the event).
The QI marks the clarity of MO and FD signatures and borders and is subjectively determined by the observer, as is explained in Section 2.3.Measurements are separated for the inner structure and extended ICME, as is explained in Section 2.4.
For each, we determined the onset and end of the FD using the day of the year (DOY), DOY start, and DOY end, respectively; FD amplitude measured using the traditional method, FD_obs; best-fitted FD amplitude, FD_bf; MSE; and fit type as visually determined by the observer (for details see Section 3.2).

Fig. 5 .
Fig. 5. Two different versions of border selection: inner structure (yellow shading) and extended ICME signatures (grey shading).For details, see the main text.

Fig. 6 .
Fig. 6.Relative difference of the best-fit FD amplitude to the FD amplitude of the theoretical FD vs. goodness of fit for different levels of data imperfection effects: (a) Y-noise, (b) X-noise, (c) asymmetry, (d) data gap position, (e) data gap size, and (f) resolution (for details, see the main text).

Fig. 7 .
Fig. 7. Relative difference of the best-fit FD amplitude to the FD amplitude of the theoretical FD vs MSE for a mixture of different levels of data imperfection effects.In the left plot, no data refinement is included, whereas in the right plot data with extreme deviations from the identity line have been removed (for details, see the main text).The identity line (y = 100 * x) is marked as the solid black line in both figures, whereas the solid red line shows the linear regression line (linear regression equation and squared Pearson's coefficient are given in red in both plots).

Fig. 8 .
Fig. 8. Synthetic data and the best-fit results for several selected examples.Different examples correspond to different combinations of data imperfection effects, where the set of values of the data imperfection effects is different for different rows and varies from high to low and top to bottom, respectively: (1) top: Y-noise = 30%, X-noise = 30%, asym = 0.15%, gap position = 3, gap size = 30%, data resolution = 20 pts; (2) middle: Y-noise = 20%, X-noise = 20%, asym = 0.10%, gap position = 2, gap size = 20%, data resolution = 20 pts; (3) bottom: Y-noise = 5%, X-noise = 5%, asym = 0.05%, gap position = 2, gap size = 10%, data resolution = 20 pts.The right-most plot in each row shows the scatter of different dA/A and MSE values obtained through 1000 runs for the same maximum data imperfection effect levels.Different combinations of data imperfection effects may lead to different MSE and dA/A behaviour.The three examples shown left of the corresponding scatter plot are highlighted in a colour corresponding to the colour of the synthetic data for the given event (green for the left-most, yellow for the middle, and red for the right-most examples).Blue data points correspond to the original theoretical curve, and the best-fit function is given by the black curve.

Table 1 .
Statistics for the best-fit amplitude, FD_bf, the amplitude obtained using a traditional algorithm-based observational method, FD_obs, and the MSE for four different QIs (as is described in Sect.2.3).
Linear regression between the best-fitted FD amplitude, FD_bf, and the FD amplitude measured using the traditional algorithm-based observational method (from the onset to the minimum), FD_obs, for two different border selections (inner structure borders and extended ICME borders).