Issue 
A&A
Volume 671, March 2023
Solar Orbiter First Results (Nominal Mission Phase)



Article Number  A79  
Number of page(s)  11  
Section  The Sun and the Heliosphere  
DOI  https://doi.org/10.1051/00046361/202245293  
Published online  08 March 2023 
Characterising fasttime variations in the hard Xray time profiles of solar flares using Solar Orbiter’s STIX
^{1}
University of Applied Sciences and Arts Northwestern Switzerland, Bahnhofstrasse 6, 5210 Windisch, Switzerland
email: hannah.collier@fhnw.ch
^{2}
ETH Zürich, Rämistrasse 101, 8092 Zürich, Switzerland
^{3}
European Space Agency, ESTEC, Keplerlaan 1, 2201 AZ Noordwijk, The Netherlands
^{4}
PMOD/WRC, Dorfstrasse 33, 7260 Davos Dorf, Switzerland
^{5}
Space Sciences Laboratory, University of California, 7 Gauss Way, 94720 Berkeley, USA
Received:
25
October
2022
Accepted:
26
December
2022
Aims. The aim of this work is to develop a method to systematically detect and characterise fasttime variations (≳1 s) in the nonthermal hard Xray (HXR) time profiles of solar flares using highresolution data from Solar Orbiter’s Spectrometer/Telescope for Imaging Xrays (STIX).
Methods. The HXR time profiles were smoothed using Gaussian Process (GP) regression. The time profiles were then fitted with a linear combination of Gaussians to decompose the time profile. From the Gaussian decomposition, key characteristics such as the periodicity, full width at half maximum, time evolution, and amplitude can be derived.
Results. We present the outcome of applying this method to four M and X GOESclass flares from the first year of Solar Orbiter science operations. The HXR time profiles of these flares were decomposed into individual Gaussians and their periods were derived. The quality of fit is quantified by the standard deviation of the residuals (difference between observed and fitted curve, normalised by the error on the observed data), for which we obtain ≤1.8 for all flares presented. In this work, the first detection of fasttime variations with Solar Orbiter’s STIX instrument has been made on timescales across the range of 4−128 s.
Conclusions. A new method for identifying and characterising fasttime variations in the nonthermal HXR profiles of solar flares has been developed, in which the time profiles are fit with a linear combination of Gaussian bursts. The opportunity to study time variations in flares has greatly improved with the new observations from STIX on Solar Orbiter.
Key words: Sun: Xrays, gamma rays / Sun: flares
© The Authors 2023
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
Solar Orbiter is a solar and heliospheric mission led by the European Space Agency (ESA) in partnership with NASA. The scientific payload of Solar Orbiter includes four in situ and six remote sensing instruments, including the Spectrometer/Telescope for Imaging Xrays (STIX). STIX is a hard Xray imaging spectrometer that detects photons with energies in the range of 4−150 keV, with a 1 keV energy resolution (at 6 keV; Krucker et al. 2020). STIX measures bremsstrahlung emission from solar flares and therefore provides diagnostics on the hottest (≳8 MK) flare plasma (Krucker et al. 2020; Battaglia et al. 2021). This means that STIX is equipped to provide information on the accelerated electrons producing such bremsstrahlung emissions upon Coulomb collisions with ambient ions. STIX has a high time resolution, with the ability to sample down to 0.1 s. It also has a stable background and continuously observes the full solar disc. Thanks to these capabilities, STIX is an excellent instrument for measuring timevarying signatures on short timescales (≳1 s) in the Xray emission of solar flares.
Such time variations are of great interest because they are related to the fundamental timescales occurring in solar flares, such as energy release and particle acceleration processes, as well as magnetohydrodynamic (MHD) waves and oscillations in or around the flare site. It is imperative to understand the origin and nature of the observed timevarying behaviour in order to achieve a unified solar flare model.
Fast time variations, which are sometimes classified as quasiperiodic pulsations (QPPs), have been observed in the emission from solar flares over the past 50 yr, with some of the earliest studies identifying timevarying signatures in the hard Xray (HXR) energy range (Parks & Winckler 1969). The particular variations of interest in this work are modulations in the flare intensitytime profiles. The ones classified as QPPs typically have periodicities ranging from a few seconds to several minutes (Zimovets et al. 2021). Subsecond spikes in the HXR flare emission have also been detected (Roberts et al. 1983; Qiu et al. 2012; Knuth & Glesener 2020). In some cases, multiple periods of oscillations have been identified alongside amplitude modulation (Kolotkov et al. 2015; Van Doorsselaere et al. 2016; McLaughlin et al. 2018). Furthermore, Hayes et al. (2016) identified these timevarying signatures in both the impulsive and decay phases, thus challenging our understanding of the flare model. Additionally, studies relating active region and flare properties to QPP periodicities have also been performed (Pugh et al. 2017; Hayes et al. 2020). For example, Hayes et al. (2020) found a positive correlation between QPP period and duration. Despite this, QPP period was found to be independent of flare magnitude. Finally, a significant correlation between QPP period and various ribbon properties was found.
Statistical studies have found a wide range of probabilities that a given > M5 GOES class flare will contain a QPP, with estimates ranging from 30−90% (Simões et al. 2015; Inglis et al. 2017; Dominique et al. 2018; Hayes et al. 2020). It has been made clear that timevarying signatures (whether quasiperiodic or not) are commonly observed phenomena which occur in flare emission across the entire electromagnetic spectrum, from radio waves to γrays (Nakariakov & Melnikov 2009; Zimovets et al. 2021). These statistics rely on QPP classifications that search for a statistically significant periodic component in the data and thus exclude all fasttimevarying structures that do not display any significant periodicity. Recent reviews of QPP properties include Zimovets et al. (2021), Kupriyanova et al. (2020), McLaughlin et al. (2018), Van Doorsselaere et al. (2016), Nakariakov & Melnikov (2009).
A number of models have been built to attempt to explain the observed timevarying and oscillatory phenomena in solar flare emission, with at least fifteen different models proposed at present (Zimovets et al. 2021). For a recent review of current models and their observational signatures, we refer to Zimovets et al. (2021) and McLaughlin et al. (2018). According to Kupriyanova et al. (2020), these models can generally be divided into three categories: (1) models in which the observed emission is directly modulated by MHD and electromagnetic (EM) waves; (2) models in which the efficiency of energy release is modulated by MHD waves; and (3) models in which the original energy release process is itself quasiperiodic in nature. Despite the many observations of timevarying behaviour in flares, there are few cases where the underlying mechanism(s) behind the variations are unambiguous. This is often due to observational constraints. In addition to this, many mechanisms have nonunique signatures, making the disambiguation more challenging (McLaughlin et al. 2018; Zimovets et al. 2021).
One of the challenges in the study of timevarying signatures and QPPs in solar flares is the identification of the characteristic timescale or period. Broomhall et al. (2019) analysed the strengths and weaknesses of various stateoftheart techniques for detecting timevarying signatures classified as QPPs, based on a series of tests. The tests involved a simulated dataset of flare time profiles, some of which include a periodic component. From this analysis, a list of eight recommendations for detection methods was set out. One of these recommendations relates to a common preprocessing step, namely: timeseries detrending. Broomhall et al. (2019) showed that detrending can be challenging and if it must be done, it should be performed manually, with a timedependent smoothing window; otherwise, a bias may be introduced. Furthermore, there is the added complexity of the underlying background trend when analysing coronal time series. Auchère et al. (2016) showed that if the powerlaw dependency of the Fourier power spectra (also known as red noise) is not considered in detection models, this can lead to false detections of a significant periodic component. Inglis et al. (2015) showed the effect of this powerlaw component in the Fourier power spectrum and took this red noise component into consideration when detecting QPPs with methods such as AFINO (Inglis et al. 2017). Many timevarying signatures in flares also show nonstationary properties, whereby the amplitude and period change as a function of time (Nakariakov et al. 2019). This is a further challenge to detection methods, for which typical Fourierbased approaches do not work well.
The method developed in this work leverages the recent data from Solar Orbiter’s STIX and takes a new approach to analysing fasttimevarying signatures in solar flare HXR time profiles. Since the modulation depth in the HXR signature is large, we assume that detrending is not a necessary step. Instead, the fasttimevarying behaviour is considered here as a linear combination of individual Gaussian contributions to the total observed signal. To achieve this, first the HXR time series are smoothed using Gaussian Process (GP) Regression. Next, the time profiles are fitted with a linear combination of Gaussians. From the Gaussian decomposition, key characteristics such as the waiting time distribution, time evolution of the peak full width at half maximum (FWHM), and amplitude can be derived. Additionally, this method can detect nonstationary (timeevolving) signatures as the method is agnostic to whether the underlying driver is periodic. The derived timing information can be used to spatially resolve sources and to perform timedependent spectral analyses of the individual peaks. An important distinction between this method and those mentioned in this section is that the method does not rely on the assumption that the observed oscillatory behaviour is periodic. This method can be applied to HXR time profiles which exhibit fasttime variations that may or may not be identified as quasiperiodic by existing analysis techniques. This aspect broadens the scope of our analysis.
The analysis given in this work is based on four M and X GOESclass flares observed by STIX since September 2021. The opportunity to study time variations in flares has greatly improved thanks to the new observations from STIX on Solar Orbiter, which include hundreds of flares that demonstrate a fasttime variability.
2. Observations
The data set used in this work is from the STIX imaging spectrometer onboard Solar Orbiter. STIX creates images by an indirect Fourier based imaging technique (Krucker et al. 2020). It has a movable attenuator which is inserted during large flares to limit exposure to high count rates from lowenergy Xray photons and to avoid instrument saturation (Krucker et al. 2020). STIX is capable of quantifying the location, spectrum, and energy content of flareaccelerated thermal and nonthermal electrons (Krucker et al. 2020).
STIX has an onboard algorithm which performs the dynamic time binning. This means that the time resolution of data taken by STIX can vary throughout the duration of the flare. The highest possible time cadence is 0.1 s. To date, the highest cadence tested is 0.3 s. It is not possible to take data at 0.1 s for a long time (≳1 h) since this fills up the onboard memory. However, it is possible to run at 0.5 s cadence for as long as desired. As such, the time resolution of the data used for analysis in this work is 0.5 s. The overview plots shown in this work are made using lower latency data, with a dynamicallyadjusted temporal resolution, which typically ranges from 2 to 20 s. The energy ranges used vary depending on the individual flare profile and are selected such that they contain mainly nonthermal HXR emission.
STIX has a relatively stable nonsolar background and is a full Sun imager. It observes the Sun continuously, unlike its predecessor, the Reuven Ramaty High Energy Solar Spectroscopic Imager (RHESSI), which rotated about its own axis with a period of 4 s (Lin et al. 2002). This rotational period sometimes introduced artificial periodicities in the data, which were difficult to disentangle from true flare fluctuations (Inglis et al. 2011). Due to its orbit, RHESSI had a daynight cycle, which meant that certain sections of a flare were often not observed. These constraints limited RHESSI’s capability for analysing timevarying structures on short timescales. In contrast, STIX is an excellent instrument suited to the task of fasttime variations analysis due to its stable background, high temporal resolution, imaging and spectral capabilities. Furthermore, STIX is on board the Solar Orbiter spacecraft, which is in a unique elliptical orbit about the Sun, observing solar dynamics from a different vantage point than Earth.
3. Methodology
The motivation for this work is to develop a method that systematically identifies and characterises the modulation observed on short timescales in HXR nonthermal emission from solar flares. This characterisation is useful because it enables us to gain information about the timing, shape, and origin of HXR pulsations in a systematic way. The methodology involves two main steps: (1) data preprocessing by means of time series normalisation and signal smoothing and (2) Fitting a linear combination of Gaussians to the preprocessed time series. Each step and its motivation is explained in more detail in the following subsections^{1}.
3.1. Data preprocessing
In this work, a Gaussian process (GP) regression is used to smooth the HXR time profiles. This is a nonparametric, Bayesian approach to regression that uses machine learning techniques to fit the data.
3.1.1. Normalization
To perform GP regression, the first step is normalising the data. To do so, we rescale the original time profile using scikitlearn’s StandardScalar class (Pedregosa et al. 2011) by subtracting the mean, μ and dividing by the standard deviation, σ, such that the new time profile is characterised by μ′ = 0 and σ′ = 1. Figures 1a,b demonstrates this for the STIX 25−76 keV time series of the SOL20220504 flare.
Fig. 1. Various steps involved in the application of this method for the example M1.2 GOES class flare of SOL20220504. The times shown are in seconds from 20220504 15:16:12 (Earth time). 
3.1.2. Gaussian process (GP) regression
Gaussian processes (GPs) are a powerful supervised machine learning tool that provide a means to make predictions about data by incorporating prior knowledge. The most obvious area of application is regression problems, such as timeseries forecasting. Overall, GPs are not limited to regression and can also be extended to classification and clustering tasks. In this work, GP regression is used to smooth STIX time profiles. This is necessary to identify the local maxima and minima for fitting. Since GP regression is a nonparametric approach, there is no initial assumption about the functional form of the time profile. In contrast, when performing polynomial regression, we assume that the functional form is a polynomial and then estimates the exact coefficients. Furthermore, GP regression uses Bayesian statistics, meaning that the model begins with an initial guess of the fit and iteratively updates its fit as more data (in this case, the intensity of HXR emission at a given time) is fed to the algorithm. Unlike traditional regression methods, which learn the exact values for each function parameter, GP regression models estimate a probability distribution over all possible fits.
A common issue with traditional methods for smoothing data in search of QPPs is the choice of window size. This choice is typically ambiguous and can lead to erroneous results if a window size is not selected manually for individual flares (Broomhall et al. 2019). Since GP regression is based on Bayesian statistics, the model parameters can be chosen systematically by optimising a loss function, which estimates how well the model fits to the data. Also, GP regression has the added benefit of giving an uncertainty estimate on the fitted smooth curve and, unlike most methods involving smoothing windows, GPs can be applied to unevenly sampled data. This is particularly important for STIX data, as it is binned onboard in a dynamic way to optimise the amount of downlinked data (Krucker et al. 2020). Of course, time series can be interpolated and resampled at even cadence, but this introduces additional errors; thus, it is preferable to have a method that does not require such preprocessing.
3.1.3. GP kernel and hyperparameter optimisation
A GP is a collection of random variables such that any finite collection of those random variables has joint Gaussian distributions (Rasmussen 2004). A GP can be entirely described by its mean and covariance function (Rasmussen 2004). This is similar to a Gaussian distribution, which is defined by a mean vector and covariance matrix; however, in a GP the distribution is taken over function space (Rasmussen 2004). In other words, the final best fit obtained using GP regression is given by the maximum likelihood of a multidimensional probability distribution over all possible fits. This distribution is specified by a covariance matrix (also known as a kernel).
In this work, a simple radial basis function (RBF) kernel was chosen because an RBF kernel is a common choice for smooth time series data. No inference is made based on the kernel parameters themselves; the GPs are solely used as a smoothing technique. Hübner et al. (2022) have discussed assigning a physical meaning to kernel parameter choice in the context of QPPs. For that case, it would be important to consider various types of kernels. However, for our purpose, a simple choice of kernel was deemed sufficient.
The kernel chosen for this work has two components: a constant kernel component and a radial basis function kernel (squared exponential kernel) component. It takes the following form:
$$\begin{array}{c}\hfill k({x}_{i},{x}_{j})=A{e}^{\frac{d{({x}_{i},{x}_{j})}^{2}}{2{l}^{2}}},\end{array}$$(1)
where A is a constant value, d(x_{i}, x_{j}) is the Euclidean distance between two feature vectors x_{i}, x_{j}, and l is the length scale.
For this kernel choice, there are two hyperparameters to optimise, A and l; here, l roughly describes the length of the (time) scale over which data points are related. For example, for a flare profile with subsecond variations, the length scale parameter is expected to be a lot smaller than for a profile with longer (on the order of seconds) variations. Then, A is the amplitude of the kernel. In addition, there is a parameter α, which is added to the diagonal of the kernel matrix to help with fitting and can be considered as Gaussian measurement noise on the training data.
The choice of kernel and its parameters can be optimised by performing an exhaustive search over a wide range of values and computing the performance of the model using cross validation. This is known as hyperparameter optimisation. In crossvalidation, the data is split between training and test data. The training data is used to fit model parameters and the test data is used to assess model performance. In this work the split between training and test data is 80:20. In other words, the model is fit using the training data, and the “goodness of fit” is assessed by comparing the model prediction with the known measured value from the test set. The error on the fit is computed based on a chosen loss function – in this case, this is the mean squared error (MSE). An average MSE value on all test data is used to assess overall model performance. Once the model hyperparameters have been optimised, a smooth time profile is obtained (as shown in Fig. 1c).
3.2. Fitting a linear combination of Gaussians to the HXR profile
In many flare cases, the shape of each HXR pulsation can be aptly fit with a Gaussian function, which gives a symmetric rise and fall. The use of a simple functional form such as a Gaussian is preferable because we can easily derive insightful properties from the mean and standard deviation of each individual peak, including the waiting time between peaks, amplitude, and full width at half maximum (FWHM). To a large degree, the form of a Gaussian accurately models the sudden impulse of nonthermal electrons reaching the chromosphere and interacting to produce nonthermal bremsstrahlung emission. Of course, there are cases in which this form may not be a suitable choice; for instance, when there is clear particle trapping in the coronal loop, leading to asymmetric HXR emission profiles and a longer decay time. However, to a large extent, the choice of this functional form describes the observed HXR emission extremely well. Other reasonable choices would include a triangular pulse or a double exponential (a symmetric exponential rise and decay). In order to account for particle trapping, it could also be useful to consider an asymmetric function, such as one with an exponential rise but with a slower decay rate.
In this step, the smoothed signal is fit by a linear combination of Gaussian functions. The gradient of the time profile is computed and the locations where the gradient is zero are identified as peaks and troughs (see Fig. 1d). Each pair of local maxima and minima are considered as a single Gaussian contribution. From these values, initial parameters for the fitting routine are derived. The time of the peak value is considered to be the mean of the Gaussian and the time between peak and trough is taken as a rough estimate of the FWHM. The height of the curve at the peak time is taken as an initial estimate of amplitude. Scipy’s Curve Fit (Virtanen et al. 2020) routine is used with the derived initial values as input and with upper and lower bounds on the possible parameter space to fit the time profile. The range of reasonable parameter values varies with each flare. A large GOES class flare will have a smaller possible range in peak height for fitting, as the counting statistics are better. In general, the range used is quite large to allow for the curve fitting routine to work effectively, yet it excludes wide Gaussian fits. Without putting a restriction on the range of peak widths, the curve fitting routine may return Gaussian contributions with large FWHM. This is demonstrated in Fig. 10.
4. Results
In Sect. 4.1, we present the detailed results of each step of the method for the SOL20220504 flare, for the purposes of demonstration. The process has been tested on four M and Xclass solar flares observed with STIX since September 2021. In Sect. 4.2, we give an overview of all flares analysed and the results obtained. Finally, Sect. 4.3 presents an analysis of QPP identification and characterisation for each event.
4.1. Method results detailed for the SOL20220504 event
4.1.1. Gaussian process regression
GP regression was applied to the SOL20220504 flare of M1.2 GOES class. A grid search was performed over the hyperparameter space shown in Table 1. For this flare, the optimal solution was found to be A = 1.8, α = 0.007, and l = 0.08. The smoothed GP prediction for the optimised hyperparameters is shown in Fig. 1c with a 95% confidence interval.
Model hyperparameter grid search domain for the SOL20220504 flare.
4.1.2. Gaussian fitting
The peaks and troughs of the smoothed SOL20220504 time profile are identified (as shown in Fig. 1d) and used as input into the fitting routine. The fitting routine fits a linear combination of individual Gaussians to the smoothed time profile. Constraints are given on the possible range of Gaussian characteristics, and for the case of SOL20220504, the bounds applied are those shown in Table 2. The resulting fit is shown in Fig. 1e. Overall, a good fit to the data is obtained. However, at t ≈ 60 s, several peaks are not well fit. This is because the smoothing step (see Fig. 1c) has removed some of these peaks from the time series as they are very shortlived. One of the short peaks remains in the smooth profile, although it has a smaller amplitude and, thus, it is also not very well fit in the final Gaussian decomposition. This demonstrates an important drawback of the smoothing step.
Bounds applied to the fitting parameters for SOL20220504 flare.
4.2. Flare overview
4.2.1. SOL20210923 M1.9 GOES class flare
The SOL20210923 flare was an early impulsive M1.9 GOES class flare observed by STIX when Solar Orbiter was at 0.60 AU from the Sun and at an angle of ∼33° to the SunEarth line. A case study of this flare by Zoë Stiefel et al. (2023) reveals four nonthermal HXR sources observed by STIX, with two inner sources from the traditional flare loop and outer footpoints which may be associated with the early onset of a filament eruption; however, the evidence on the latter is inconclusive. A Gaussian decomposition for this flare was performed and the resulting fit (shown in Fig. 2) has been assessed by calculating the standard deviation of the normalised residual. The normalised residual is calculated as the difference between the fit and the observed STIX time profile normalised by the error on the measured STIX time profile. The error on the STIX observed time profile is a combination of the compression error and counting statistics error. The fit gives σ_{R} = 1.83. The residuals are quite high, particularly at the beginning of the impulsive phase because there are fasttime variation on timescales smoothed out by the Gaussian process regression. However, the fit improves in the later phase.
Fig. 2. Overview of the SOL20210923 flare. The 18−28 keV HXR nonthermal Gaussian decomposition is shown with a standard deviation of the residual σ_{R} = 1.83. The first few peaks have substructure on timescales smoothed out in the GP regression output. As a result, variations on these timescales are not well fit. 
4.2.2. SOL20211009 M1.6 GOES class flare
The SOL20211009 flare was an M1.6 GOES class flare observed when Solar Orbiter was at a distance of 0.68 AU from the Sun with the spacecraftSunEarth angle of ∼15°. The Gaussian decomposition method gives a fit with standard deviation of the residual σ_{R} = 0.62, shown in Fig. 3. In this case, the measured data is noisier than other example flares because the flare is smaller. As a result, it is easier for the fit to appear to perform well on global structures, but the variance among the fits is high.
Fig. 3. Overview of the SOL20211009 flare and the Gaussian decomposition of its 20−25 keV HXR nonthermal time profile. The fit gives a standard deviation of the residual σ_{R} = 0.62. This is a strong fit. The noise level on the measured values of the counts is higher for this flare since it is smaller. This makes it easier to fit and thus the residuals are smaller than those of a large flare with higher signaltonoise ratio. 
4.2.3. SOL20220330 X1.4 GOES class flare
The X1.4 GOES class flare observed on 30032022 was the first X class flare to be observed by STIX during perihelion. This observation was taken at a distance of 0.33 AU from the Sun. Furthermore, the angle between Solar Orbiter and the SunEarth line at this time was ∼95°. Due to the large flux incident on the instrument during this flare, the aluminium alloy attenuator was automatically inserted when a certain trigger threshold was reached. The periods when the attenuator was inserted are shown in Fig. 4. The attenuator mainly blocks low energy Xray photons and allows high energy photons to be transmitted. The energy dependent response of the attenuator can be found in Krucker et al. (2020). As the counts incident on the detector increase, the livetime decreases and vice versa when the attenuator is inserted. Therefore, to correctly show the time profile of such a flare, we should carefully account for changes in detector livetime as well as the energy and spectrum dependent attenuator transmission. This correction is nontrivial as it requires detailed knowledge of the intricacies of such an instrument. However, for the background detector (BKG), the pixels with large apertures are turned off when flux is high. This means that the attenuator does not cover this detector when inserted. As such, there is no effect exerted by the attenuator motion on the 4−10 keV background detector time profile shown. The background detector was used for the thermal 4−10 keV profile shown in Fig. 4, as it requires no correction for the attenuator motion.
Fig. 4. Overview plot of the SOL20220330 flare and the Gaussian decomposition of the 32−45 keV HXR profile for three different phases: phases 1, 2, and 3 from left to right. The standard deviation of the residual of each phase is σ_{R, 1} = 1.25, σ_{R, 2} = 1.58 and σ_{R, 3} = 1.11. 
This flare is an early impulsive flare which shows timevarying structures. In particular, there are variations on the timescale of several seconds. Later in the flare, during the thermal peak, there are further variations with significantly different fluctuation timescales (∼14 s) and, later on, even longer periods of ∼35 s. The pulsation timescales are changing significantly over the course of the flare and so, it is challenging to smooth the time profile using our simple choice of kernel, since the lengthscale parameter has no time dependence. To account for this, the time profile was split into three different sections: the early phase (phase 1), middle phase (phase 2), and late phase (phase 3). Figure 4 shows the Gaussian decomposition for each phase with the standard deviation of the residuals being σ_{R, 1} = 1.25, σ_{R, 2} = 1.58 and σ_{R, 3} = 1.11 for phases 1, 2, and 3, respectively. In particular, phases 1 and 3 appear to be well fitted with a linear combination of Gaussians, but the residuals for phase 2 show sinusoidal variation, indicating the data are not well fit by the model.
4.2.4. SOL20220504 M1.2 GOES class flare
Finally, we present the M1.2 GOES class flare on May 4th 2022. At this time, Solar Orbiter was observing the far side of the Sun from the Earth, namely, the spacecraftSunEarth angle was ∼163° and the distance to the Sun was 0.73 AU. As such, the flare was not observed by Earth based observatories. The GOES class of a flare which was not observed from Earth is estimated using a model that fits STIX 4−10 keV counts to the GOES flux^{2}. Similarly to the SOL20220330 flare, the attenuator motion is shown in Fig. 5 and the background detector was used for the thermal 4−10 keV profile. The standard deviation of the residual σ_{R} = 1.39. Overall, the data are well fit with a linear combination of Gaussians. The fastest fluctuations in time are not well fitted (e.g. Fig. 1e at t ≈ 60 s) because the time profile output by the GP regression smooths out some very short variations. Therefore, the Gaussian fitting procedure doesn’t fit these shorter peaks. This is a limitation of the method that is discussed further in Sect. 5.2.
Fig. 5. Overview plot of the SOL20220504 flare and the Gaussian decomposition of the impulsive phase 25−76 keV profile. The standard deviation of the residual for the fit is σ_{R} = 1.39. 
4.3. QPP identification and characterisation
4.3.1. Line Fitting Method
An estimate of signal (quasi) periodicity can be derived from individual Gaussian components. Figure 6 shows the time of each Gaussian component against the peak number for the SOL20220504 flare. The slope of the line of best fit gives an estimate for the period, ∼5.55 s. This shows that STIX can detect variations on short timescales, which is something that was not feasible with its predecessor RHESSI. This result is then compared with the AFINO analysis method (Inglis et al. 2017). AFINO fits the Fourier power spectrum of a flare signal with four different models: (1) a power law + constant; (2) power law + constant + Gaussian A broken power law + constant; (3) power law + 2 Gaussians. Should a QPP component be present in the signal, the second or fourth model would be favoured, that is, the Gaussian peak fits the enhanced power due to the QPP that is present. The goodness of fit is estimated from a Bayesian information criterion (BIC), given by:
$$\begin{array}{c}\hfill \mathrm{BIC}=2ln(L)+kln(n),\end{array}$$(2)
Fig. 6. Mean time of each Gaussian component against peak number. The slope of the line fit gives an estimate of periodicity in the signal. A stronger fit with multiple lines would indicate that there is nonstationarity in the signal. 
where L is the maximum likelihood, k is the number of free parameters, and n is the number of data points in the Fourier power spectrum. Thus, the BIC score penalises overfitting. A large negative value for a given model indicates a strong fit to the data. One model is said to be strongly preferred over another if ΔBIC > 10. In the case of the SOL20220504 flare, the BIC score for model 1 versus model 2 is ΔBIC_{12} = 4.3, whereby model 2 gives a larger negative BIC score. Thus, the AFINO method gives a slight preference for the QPP model with single period, $P=4.{96}_{0.54}^{+0.66}$ s, over a simple power law model (see Fig. 7). This is largely consistent with the period obtained from the Gaussian decomposition method, which is slightly higher since smoothing the time profile suppresses the fastest time variations. We notice that when assessed based on the standard AFINO criteria (ΔBIC > 10), this flare would not be marked as a QPP detection, although there is clearly enhanced power at f_{0} = 0.202 ± 0.024 Hz ($P=4.{96}_{0.54}^{+0.66}$ s). From the line fitting method, there is the added benefit of obtaining the exact timing information of each pulsation, which can be used to spatially resolve each peak (see Fig. 6).
Fig. 7. Power spectral density of signal shown in Fig. 1a. The AFINO QPP model best fit is shown. AFINO analysis gives a period of $\sim 4.{96}_{0.54}^{+0.66}$ s, with a moderate preference for the QPP model with ΔBIC = 4.3. 
By fitting two lines, a test can be performed to check for nonstationarity (changes in quasiperiodicity over time) in the flare signal. The derived slopes are 4.48 s and 6.78 s for the early and late stages, respectively. However, the Pearson correlation coefficients obtained indicate that a single line fit or single periodicity in the signal gives a better fit to the data. This is, of course, biased by the number of data points in each fit. Further analysis was done where the signal was split into two time ranges, corresponding to the time ranges of the two line fits. AFINO analysis was performed on these two time ranges. In this case, the QPP model was only favoured in the first time range. Therefore, the AFINO method favours a single QPP model over two periods in this particular case. This agrees with the apparent preference for a single periodicity given by the line fit. The line fitting method was also applied to the SOL20210923, SOL20211009 and SOL20220330 flares and the results are shown in Figs. 8 and 9.
Fig. 8. Time of each peak from the Gaussian decomposition against peak number for the SOL20210923 and SOL20211009 flare shown in Figs. 2 and 3, respectively. The fit gives an estimate of periodicity in the SOL20210923 nonthermal HXR flare signal of ∼19 s. For the SOL20211009 flare, a periodicity of ∼128 s is derived. 
Fig. 9. Line fits to the peak time over peak number for three phases of the SOL20220330 flare shown in Fig. 4. Since the pulsation frequency is increasing significantly in time, each phase was analysed separately. Phase 1 is very well fitted by a line indicating that the timing of peaks is quasiperiodic, phases 2 and 3 are not as well fit by a line. This agrees with the error analysis presented in Table 3, where phases 2 and 3 have high variance. Between phases 1 and 3, the periodicity has increased by a factor of ∼5. This indicates that the driver of these pulsations could be different in the later stages of the flare, compared to the initial impulsive phase. 
It has been demonstrated that this new method is capable of accurately detecting and characterising fasttimevarying structures and QPPs in the nonthermal HXR emission from flares. The results from this method are consistent with Fourier based analysis methods such as AFINO and in addition, allow for the extraction of important information, for instance, the time between peaks and pulsation duration. Further, this method can be used to detect frequency drifts and nonstationarity in time series.
4.3.2. Errors on QPP analysis
The Gaussian decomposition fits are nonunique and multiple solutions can be derived for a given time series. To determine the effect of this, error analysis was performed for determining the QPP period. For this, an array of errors was added to each STIX time series, which is some random integer multiple of the calculated measurement error (that includes the uncertainty introduced due to compression and counting statistics). The random integer is sampled from a normal white noise distribution with zero mean and unit standard deviation. The entire fitting procedure is then performed 50 times, with different initial error arrays added to the time profile. This gives many unique Gaussian decomposition fits for a given flare profile. As previously, the periodicity of the signal is estimated from a single line fit and an average period for each flare is obtained, with a standard deviation over 50 iterations. The mean periods and standard deviation are as shown in Table 3. An estimation of the range of periods that reasonably describe the variations observed is also given.
Mean and range of periods derived for each flare over 50 iterations, where random error has been added to the initial time series.
These results indicate that the period estimation for flares SOL20210923, SOL20223003, and SOL20220504 in Figs. 2, 4P1, and 5, respectively, are strong fits, since the estimated periods are within 1.1σ of the means obtained from the error analysis. However, the SOL20211009 flare has high variance when noise is added to the time profile. This is because the HXR counts are low and as a result the measurement error is relatively high; therefore, when additional error is added, the Gaussian decomposition fits vary significantly. The fit obtained in Fig. 3 is thus unreliable, as the derived period is 2.9σ from the mean period over 50 iterations. We notice that for the SOL20220330 flare, phase 1 gives a strong fit; however, phases 2 and 3 have a wide range of periods – thus, they are not as well constrained.
5. Discussion and conclusions
5.1. First analysis of fast timevarying structures and QPP detection with STIX data
In this work, we present the first analysis of fasttimevarying structures in the nonthermal HXR emission from flares using Solar Orbiter’s STIX instrument. We have developed a new method for identifying and quantifying fasttimevarying structures in the HXR emission from flares. This method decomposes the nonthermal HXR emission from flares in the impulsive phase into a linear combination of Gaussian pulses. For a sample of four M and Xclass flares, the standard deviation of the normalised residual is ≤1.8. This indicates that the model is a good fit to the data for this selection of events.
From these four flares, fasttime variations on timescales ranging from 4−128 s have been characterised. Furthermore, the first detection of solar flare QPPs from STIX observations has been made. It has been shown that QPPs with timescales down to the order of ∼4 s can be detected with STIX, which was not possible with it predecessor RHESSI due to the spacecraft rotational period of ∼4 s. This work demonstrates that STIX is an instrument well suited to the detection of fasttime variations in the HXR emission from solar flares in the timescale range of seconds to minutes, due to its high time resolution and relatively constant nonsolar background.
5.2. Drawbacks
In this section, we present and discuss the drawbacks and subtleties associated with the method and how they may impact the derived results:
1. One important drawback of this method is that GP regression smooths time variations on timescales that are shorter than the optimal length scale obtained from hyperparameter optimisation. This means that some peaks that are very shortlasting are smoothed out and, hence, they are not fit by the Gaussian fitting routine. Importantly, the timescales that are suppressed are a function of the sampling rate, since with a higher time cadence, shorter variations will likely be more prominent and vice versa. As such, the drawback of this method is that variations that occur at the limit of the instrument’s sampling cadence are smoothed out. In this case with STIX, variations on short timescales such as second or subsecond are smoothed out. While QPPs can occur on short timescales such as second or subsecond periodicities, particularly in the radio band (Tan & Tan 2012; Nakariakov et al. 2018; Carley et al. 2019), this work focuses on QPPs that have timescales on the order of seconds to minutes, which are very commonly reported (see also McLaughlin et al. 2018; Hayes et al. 2020; Zimovets et al. 2021). We have shown here that the method works well for these types of QPPs.
2. It is well known that in the case of significant particle trapping in a flare loop, the HXR time profile becomes asymmetric. In this case, fitting Gaussian curves to the time profiles of a trapped particle population is no longer physically meaningful. We should consider fitting other functional forms with asymmetric profiles.
3. A subtlety of this method is that although the Gaussian decomposition step is good at characterising nonstationarity in a signal, GP regression with our choice of kernel, is unable to model large frequency drifts in a flare since the length scale is not timedependent. This gives some limitation to the method when identifying nonstationarity. This effect is particularly pertinent to the case of the SOL20220330 flare, whereby in the early impulsive phase we observe fluctuations on the order of ∼7 s and later during the thermal peak there are fluctuations on a ∼14 s and then a ∼35 s timescale, as shown in Fig. 9. One way to rectify this is to split the time series into sections and perform the smoothing on different time ranges separately, as was done for the SOL20220330 flare. Another possibility is to consider a more complex kernel choice which has some time dependency. This method was also suggested by Hübner et al. (2022).
4. Another important note to consider is that the Gaussian decomposition fit is nonunique and influenced by the bounds applied to the curve fitting routine. Furthermore, the routine fits a fixed number of Gaussians based on the number of local maxima and minima. As such, several closely separated peaks in quick succession may not be fit by this method, if the gradient remains nonzero. Figure 10 demonstrates the effect of applying boundary conditions on the fitting parameters for the SOL20220504 flare.
Fig. 10. Effect of applying boundary conditions on the fit parameters of the Gaussian decomposition for the SOL20220504 flare. Panel (a) shows the resulting fit when the boundary conditions shown in Table 2 are applied. Panel (b) shows the resulting fit when no boundary conditions are applied. Both cases give a good fit to the data, however, the unconstrained fit (Fig. 10b) fits long period Gaussians with small amplitude peaks, for instance, at ∼15:17:08. This fit is interesting from a physics point of view, however in this work we wish to derive timing information regarding the global peak trend – therefore a constrained fit is preferable. 
5.3. Applications and future work
The potential applications of this method are wide and varied. Obtaining timing information on individual non–thermal HXR pulses in a systematic way enables a widescale, indepth study of fasttimevarying structures. In particular, the timing information obtained in this method can be used to image individual pulse contributions to the HXR flare profile. With such information, we can begin to understand the spatial structure (location, morphology) of the oscillatory source, as suggested by Zimovets et al. (2021), and make a comparison with those expected for various models. For example, should sausage mode oscillations be responsible for fasttime variations in HXR profiles, we would expect the source morphology to change over time and the location to remain fixed. In contrast to this, if repeated reconnection were to be responsible for such modulations, we would expect the HXR source spatial locations to change in time since the reconnecting field lines must change. The time evolution of such structures can be investigated with STIX’s imaging capabilities. Although it should be noted that imaging fasttimevarying HXR structures can be challenging with an indirect imager such as STIX due to its limited dynamic range; in particular, in cases where there is one source that is much brighter than the others. Furthermore, one can perform time dependent spectral analysis with STIX. This can help to identify and investigate the underlying particle population behind such oscillatory signatures. Additionally, timedependent imaging spectroscopy can be used to localise particle populations.
This method adds to the set of current QPP detection techniques, with the additional ability of characterising nonQPP, fasttimevarying signatures. The quantities and relationships derived by this method are important as they feed back into modeling efforts. For example, largescale waiting time distributions could be used to assess whether avalanche models can accurately describe the observed fasttime variation phenomenon.
Future work will focus on using the derived timing information and shape to determine over which interval HXR images should be made in order to understand the source origin. As STIX is an indirect Fourier imager, another interesting avenue of investigation will be to analyse the evolution of Fourier components (visibilities) over time. We will focus on combining the analysis with other datasets at various wavelength of emission. This method will be applied to other HXR datasets such as Fermi GammaRay Space Telescope (Meegan et al. 2009) and RHESSI (Lin et al. 2002), and to future instruments including ASOS/HXI (Zhang et al. 2019) and AdityaL1/HELIOS (Seetha & Megala 2017). This method will also be applied to other wavelengths of emission such as pulsations observed in radio by the Expanded Owens Valley Solar Array (EOVSA; Gary et al. 2018), which are often correlated with those seen in HXR emission (Aschwanden et al. 1990).
In conclusion, a new method for identifying and characterising fasttime variations in the nonthermal HXR time profiles of solar flares has been developed, in which the signals are decomposed into individual Gaussian contributions. The fits obtained have a standard deviation of the normalised residual of ≤1.8. The first characterisation and detection of fasttime variations in the HXR profile of solar flares with STIX has been made on timescales between 4 and 128 s.
The opportunity to study time variations in flares has greatly improved with new observations from STIX on Solar Orbiter and its observations of numerous flares that demonstrate fasttime variability over a wide range of timescales.
The code used here is publicly available at https://github.com/hannahc243/Gaussian_Decomp
Acknowledgments
Solar Orbiter is a space mission of international collaboration between ESA and NASA, operated by ESA. The STIX instrument is an international collaboration between Switzerland, Poland, France, Czech Republic, Germany, Austria, Ireland, and Italy. H.C., A.F.B. and S.K. are supported by the Swiss National Science Foundation Grant 200021L_189180 for STIX. L.A.H. is supported by an ESA Research Fellowship.
References
 Aschwanden, M. J., Benz, A. O., & Kane, S. R. 1990, A&A, 229, 206 [NASA ADS] [Google Scholar]
 Auchère, F., Froment, C., Bocchialini, K., Buchlin, E., & Solomon, J. 2016, ApJ, 825, 110 [Google Scholar]
 Battaglia, A. F., Saqri, J., Massa, P., et al. 2021, A&A, 656, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Broomhall, A.M., Davenport, J., Hayes, L., et al. 2019, ApJS, 244, 44 [NASA ADS] [CrossRef] [Google Scholar]
 Carley, E. P., Hayes, L. A., Murray, S. A., et al. 2019, Nat. Commun., 10, 2276 [NASA ADS] [CrossRef] [Google Scholar]
 Dominique, M., Zhukov, A. N., Dolla, L., Inglis, A., & Lapenta, G. 2018, Sol. Phys., 293, 61 [NASA ADS] [CrossRef] [Google Scholar]
 Gary, D. E., Chen, B., Dennis, B. R., et al. 2018, ApJ, 863, 83 [CrossRef] [Google Scholar]
 Hayes, L. A., Gallagher, P. T., Dennis, B. R., et al. 2016, ApJ, 827, L30 [NASA ADS] [CrossRef] [Google Scholar]
 Hayes, L. A., Inglis, A. R., Christe, S., Dennis, B., & Gallagher, P. T. 2020, ApJ, 895, 50 [Google Scholar]
 Hübner, M., Huppenkothen, D., Lasky, P. D., et al. 2022, AJ, 936,17 [CrossRef] [Google Scholar]
 Inglis, A. R., Zimovets, I. V., Dennis, B. R., et al. 2011, A&A, 530, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Inglis, A. R., Ireland, J., & Dominique, M. 2015, ApJ, 798, 108 [Google Scholar]
 Inglis, A., Ireland, J., Dennis, B. R., Hayes, L. A., & Gallagher, P. T. 2017, AAS/Solar Physics Division Meeting, 48, 400.05 [NASA ADS] [Google Scholar]
 Knuth, T., & Glesener, L. 2020, ApJ, 903, 63 [NASA ADS] [CrossRef] [Google Scholar]
 Kolotkov, D. Y., Nakariakov, V. M., Kupriyanova, E. G., Ratcliffe, H., & Shibasaki, K. 2015, A&A, 574, A53 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Krucker, S., Hurford, G. J., Grimm, O., et al. 2020, A&A, 642, A15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kupriyanova, E., Kolotkov, D., Nakariakov, V., & Kaufman, A. 2020, SolarTerr. Phys., 6, 3 [NASA ADS] [Google Scholar]
 Lin, R. P., Dennis, B. R., Hurford, G. J., et al. 2002, Sol. Phys., 210, 3 [Google Scholar]
 McLaughlin, J. A., Nakariakov, V. M., Dominique, M., Jelínek, P., & Takasao, S. 2018, Space Sci. Rev., 214, 45 [Google Scholar]
 Meegan, C., Lichti, G., Bhat, P. N., et al. 2009, ApJ, 702, 791 [Google Scholar]
 Nakariakov, V., & Melnikov, V. 2009, Space Sci. Rev., 149, 119 [NASA ADS] [CrossRef] [Google Scholar]
 Nakariakov, V. M., Anfinogentov, S., Storozhenko, A. A., et al. 2018, ApJ, 859, 154 [NASA ADS] [CrossRef] [Google Scholar]
 Nakariakov, V. M., Kolotkov, D. Y., Kupriyanova, E. G., et al. 2019, Plasma Phys. Control. Fusion, 61, 014024 [Google Scholar]
 Parks, G. K., & Winckler, J. R. 1969, ApJ, 155, L117 [Google Scholar]
 Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2011, J. Mach. Learn. Res., 12, 2825 [Google Scholar]
 Pugh, C. E., Nakariakov, V. M., Broomhall, A. M., Bogomolov, A. V., & Myagkova, I. N. 2017, A&A, 608, A101 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Qiu, J., Cheng, J. X., Hurford, G. J., Xu, Y., & Wang, H. 2012, A&A, 547, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rasmussen, C. E. 2004, Gaussian Processes in Machine Learning (Berlin, Heidelberg: SpringerVerlag), 63 [Google Scholar]
 Roberts, B., Edwin, P. M., & Benz, A. O. 1983, Nature, 305, 688 [NASA ADS] [CrossRef] [Google Scholar]
 Seetha, S., & Megala, S. 2017, Curr. Sci., 113, 610 [NASA ADS] [CrossRef] [Google Scholar]
 Simões, P. J. A., Hudson, H. S., & Fletcher, L. 2015, Sol. Phys., 290, 3625 [CrossRef] [Google Scholar]
 Tan, B., & Tan, C. 2012, ApJ, 749, 28 [NASA ADS] [CrossRef] [Google Scholar]
 Van Doorsselaere, T., Kupriyanova, E. G., & Yuan, D. 2016, Sol. Phys., 291, 3143 [Google Scholar]
 Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Meth., 17, 261 [Google Scholar]
 Zhang, Z., Chen, D.Y., Wu, J., et al. 2019, RAA, 19, 160 [NASA ADS] [Google Scholar]
 Zimovets, I. V., McLaughlin, J. A., Srivastava, A. K., et al. 2021, Space Sci. Rev., 217, 66 [NASA ADS] [CrossRef] [Google Scholar]
 Zoë Stiefel, M., Battaglia, A. F., Barczynski, K., et al. 2023, A&A, 670, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
All Tables
Mean and range of periods derived for each flare over 50 iterations, where random error has been added to the initial time series.
All Figures
Fig. 1. Various steps involved in the application of this method for the example M1.2 GOES class flare of SOL20220504. The times shown are in seconds from 20220504 15:16:12 (Earth time). 

In the text 
Fig. 2. Overview of the SOL20210923 flare. The 18−28 keV HXR nonthermal Gaussian decomposition is shown with a standard deviation of the residual σ_{R} = 1.83. The first few peaks have substructure on timescales smoothed out in the GP regression output. As a result, variations on these timescales are not well fit. 

In the text 
Fig. 3. Overview of the SOL20211009 flare and the Gaussian decomposition of its 20−25 keV HXR nonthermal time profile. The fit gives a standard deviation of the residual σ_{R} = 0.62. This is a strong fit. The noise level on the measured values of the counts is higher for this flare since it is smaller. This makes it easier to fit and thus the residuals are smaller than those of a large flare with higher signaltonoise ratio. 

In the text 
Fig. 4. Overview plot of the SOL20220330 flare and the Gaussian decomposition of the 32−45 keV HXR profile for three different phases: phases 1, 2, and 3 from left to right. The standard deviation of the residual of each phase is σ_{R, 1} = 1.25, σ_{R, 2} = 1.58 and σ_{R, 3} = 1.11. 

In the text 
Fig. 5. Overview plot of the SOL20220504 flare and the Gaussian decomposition of the impulsive phase 25−76 keV profile. The standard deviation of the residual for the fit is σ_{R} = 1.39. 

In the text 
Fig. 6. Mean time of each Gaussian component against peak number. The slope of the line fit gives an estimate of periodicity in the signal. A stronger fit with multiple lines would indicate that there is nonstationarity in the signal. 

In the text 
Fig. 7. Power spectral density of signal shown in Fig. 1a. The AFINO QPP model best fit is shown. AFINO analysis gives a period of $\sim 4.{96}_{0.54}^{+0.66}$ s, with a moderate preference for the QPP model with ΔBIC = 4.3. 

In the text 
Fig. 8. Time of each peak from the Gaussian decomposition against peak number for the SOL20210923 and SOL20211009 flare shown in Figs. 2 and 3, respectively. The fit gives an estimate of periodicity in the SOL20210923 nonthermal HXR flare signal of ∼19 s. For the SOL20211009 flare, a periodicity of ∼128 s is derived. 

In the text 
Fig. 9. Line fits to the peak time over peak number for three phases of the SOL20220330 flare shown in Fig. 4. Since the pulsation frequency is increasing significantly in time, each phase was analysed separately. Phase 1 is very well fitted by a line indicating that the timing of peaks is quasiperiodic, phases 2 and 3 are not as well fit by a line. This agrees with the error analysis presented in Table 3, where phases 2 and 3 have high variance. Between phases 1 and 3, the periodicity has increased by a factor of ∼5. This indicates that the driver of these pulsations could be different in the later stages of the flare, compared to the initial impulsive phase. 

In the text 
Fig. 10. Effect of applying boundary conditions on the fit parameters of the Gaussian decomposition for the SOL20220504 flare. Panel (a) shows the resulting fit when the boundary conditions shown in Table 2 are applied. Panel (b) shows the resulting fit when no boundary conditions are applied. Both cases give a good fit to the data, however, the unconstrained fit (Fig. 10b) fits long period Gaussians with small amplitude peaks, for instance, at ∼15:17:08. This fit is interesting from a physics point of view, however in this work we wish to derive timing information regarding the global peak trend – therefore a constrained fit is preferable. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.