DETAILED INVESTIGATIONS OF THE DIMMEST BURSTS FROM TWO MAGNETARS, SGR J0501+4516 AND SGR J1550−5418

The most extreme magnetic fields in the universe have been detected from a group of isolated neutron stars, known as magnetars. About 24 magnetars have been identified so far. Most of them were discovered in the X-ray band, from either the emission of energetic bursts or pulsed emission with unique characteristics. Timing analyses of their persistent emission reveal that the magnetars spin slowly; their spin periods are in a tight range of 2–12 s, and they experience large spin-down rates (3 × 10⁻¹³ to 7.5 × 10⁻¹⁰ s s⁻¹). The dipole magnetic field strengths inferred from these spin properties are mostly of the order of ∼10¹⁴ G (Kouveliotou et al. 1998). The luminosities of their persistent X-ray emission are between 10³³ and 10³⁵ erg s⁻¹, or lower. The emission of occasional energetic bursts has become a general characteristic of magnetars: bursts have been observed from all but three confirmed magnetars. The peak luminosities of magnetar bursts or flares vary over more than eight orders of magnitude (10³⁷ to  > 10⁴⁵ erg s⁻¹). Short bursts, which typically last for ∼100 ms and have peak luminosities lower than 10⁴¹ erg s⁻¹, are the most common population.

Both persistent X-ray and occasional burst emission can be accounted for under the scheme of the magnetar model (Thompson & Duncan 1995, 1996). In this picture, the magnetic field of the neutron star plays a crucial role in generating both types of emission. The persistent emission is expected to be from a hot spot on the neutron star surface that is heated by the decay of the strong magnetic field of the neutron star. The bursts might originate from the energy release through either fracturing of the neutron star crust due to large magnetic stress or the split and reconnection of the oppositely oriented magnetic field lines (Lyutikov 2003). Magnetars occasionally undergo outburst episodes during which they experience a remarkable increase in their persistent energy output, always in conjunction with the emission of energetic short bursts (Rea & Esposito 2011).

In the magnetar model, small-scale fracturing of the neutron star crust is invoked to explain a portion of the persistent X-ray emission, while large fractures can manifest themselves as bursts (Thompson & Duncan 1996). Nakagawa et al. (2011) and Enoto et al. (2012) studied the spectra of weak bursts from SGR J0501+4516 (31 events) and SGR J1550−5418 (13 events) in 1–300 keV, respectively, observed with Suzaku. They investigated the stacked burst spectra with the sum of two blackbody functions and a power-law function and compared their spectral parameters to the results from the spectral fitting of persistent X-ray spectra with the phenomenological model (the sum of a blackbody and a power-law function). Based on the similarities between the resulting power-law indices (although with large errors), they suggested that persistent emission from magnetars can be due to a very large number of weak short bursts. We also previously studied the spectral properties of low-fluence bursts and persistent emission, both observed in a deep XMM-Newton observation of SGR J0501+4516 on its most burst-active day (Lin et al. 2012b). Instead of phenomenological functions, we adopted physically motivated models for both persistent and burst spectra. In particular, we used the surface thermal emission with a magnetosphere scattering model (STEMS,³ see Güver et al. 2007, and references therein) and a modified blackbody (Lyubarsky 2002) with resonant cyclotron scattering model (MBB+RCS,⁴ Lin et al. 2012b). Both models include the same resonant cyclotron scattering by non-relativistic electrons in a warm magnetospheric plasma (Lyutikov & Gavriil 2006). Our spectral analyses showed that the persistent emission and burst emission are better described with the STEMS and MBB+RCS models, respectively, which indicated that the persistent X-ray emission and low-fluence bursts originate through slightly different physical processes.

It is therefore crucial to systematically study temporal and spectral properties of the low-fluence magnetar bursts, in order to better understand the nature of these weak bursts and their connection with the persistent emission. However, identifying the very weak magnetar bursts using conventional techniques is a challenging task. A major fraction of these bursts are usually smoothed out in binned data. Also, it is difficult to distinguish dim bursts from fluctuations of the persistent emission. As a result, the burst search methodology in our earlier study using binned data limited us to include even weaker bursts and understand their characteristics through detailed spectral and temporal analyses.

A possible way to identify weak magnetar bursts is using the Bayesian blocks method. The Bayesian block representation of a time series is a nonparametric method used to detect local structures, especially in highly variable data sets (Scargle 1998; Scargle et al. 2013). It can provide a simple description of the overall shape of a burst temporal profile, as well as a direct measurement of the burst duration. This method has been adopted as a standard tool to calculate the duration of gamma-ray bursts (GRBs) observed with Swift/Burst Alert Telescope (BAT); this method has also been applied to search for the extended emission from short GRBs (Norris et al. 2010, 2011; Y. Kaneko et al., in preparation) and to calculate the duration of GRBs observed with Fermi/Gamma-ray Burst Monitor (GBM) (Qin et al. 2013). Moreover, the Bayesian blocks method has been used in the spatial domain to confirm the source detection in Fermi/Large Area Telescope data.

In this paper, for the first time, we introduce the Bayesian blocks method to identify very dim magnetar bursts using unbinned data of SGR J0501+4516 and SGR J1550−5418 collected with XMM-Newton and Swift/X-Ray Telescope (XRT), respectively. We also investigate the temporal and spectral properties of the identified burst sample. This paper consists of seven sections. In Sections 2 and 3, we describe the data and our burst search algorithm, respectively. We investigate burst temporal properties in Section 4. In Section 5, we examine the accumulated spectra of dim bursts with physically motivated emission models and present the results of our detailed spectral analysis. We construct and evaluate the fluence distributions of the dimmest bursts for both sources in Section 6, and finally we discuss the implications of our results in Section 7.

In the past five years, only two magnetars, SGR J0501+4516 and SGR J1550−5418, have experienced very active burst emitting episodes, during which tens or even hundreds of bursts have been observed with several X-ray and soft γ-ray telescopes. We investigated the observations of these two sources during their most burst-active episodes. In our investigations, we focused on energies below 10 keV, since the persistent emission was mainly detected in this energy range.

2.1. SGR J0501+4516

2.2. SGR J1550−5418

In the Bayesian blocks method, data are represented with a series of step functions, which can be in different sizes, determined by maximizing the block likelihood function. For the photon counting instruments, such as Swift/XRT and XMM-Newton, the detection probability of one event follows Poisson statistics. The likelihood function (L) for block k is

$$\begin{equation} \ln L^{(k)} = N^{(k)} \ln \lambda ^{(k)} - \lambda ^{(k)} T^{(k)}, \end{equation} \tag{ 1 }$$

where N^(k) is the total number of events in block k, λ^(k) is the expected count rate, and T^(k) is the duration of block k.

A priori information in this method is the improvement of the likelihood function by adding one more block. It depends on the number of events, as well as the false positive detection rate (see Equation (21) in Scargle et al. 2013). Scargle et al. (2013) calibrated this value with large sets of simulations and also showed that the detection sensitivity of this method can approach the ideal theoretical limits. The Bayesian blocks method has three advantages when identifying short magnetar bursts. First, it avoids placing a priori limitation on the scale and resolution of the search. It is sensitive to short weak bursts, which may be smoothed out when the data are binned. Second, it is applicable to data with wide ranges of amplitude, timescales, and noise levels. Finally, it is possible to determine the duration of short bursts independent of any background modeling.

We adopted the algorithm provided in Scargle et al. (2013) and generated automatic procedures to search for the bursts using the list of time-tagged events in a specified energy range, i.e., 0.2–10 keV for the XMM-Newton data and 0.3–10 keV for Swift/XRT data. Our burst search procedure included two rounds of Bayesian block representation. The first round was a blind search through the whole data set and aimed to find candidate burst blocks. The second round of Bayesian block representation was performed to confirm the candidate burst blocks, along with sufficient burst-free intervals.

To optimize the computation time for the first-round search, we employed a moving box-car approach: we calculated the Bayesian block representation for a 20 s long observing window, starting from the beginning of the observation, and shifted forward by 10 s each time. This way, the entire exposure time span, except for the first and last 10 s segments of each continuous time interval, was covered twice by consecutive observing windows in our first round of search. We then selected all of the blocks with durations less than the spin period of the source as candidate burst blocks. We excluded the candidate blocks that lined up with a boundary of any 20 s segment, as they will be picked in the middle of the previous or the next segment.

Before the second-round representation, we further gathered the candidate burst blocks into groups, in which the separation between any two adjacent blocks in the time series is smaller than 10 s. Then, we performed the second-round Bayesian blocks representation on each group, using the time-tagged events within the time segment from 10 s before the first block to 10 s after the last block. In this way, we can confirm the candidate burst blocks with sufficient underlying emission (background) in both pre- and post-burst intervals and remove the duplicates. Additionally, the policy we followed to form the groups guaranteed that all of the underlying emission time intervals are burst free.

We firmly identified the magnetar bursts after two rounds of Bayesian blocks representation in the following manner. We first selected the blocks with durations longer than 6 s as background and calculated the background count rate using the weighted mean of the background block rates. We then identified the blocks with durations less than the spin period of the source and count rates higher than the background level as burst blocks. Finally, we combined any continuous burst blocks as parts of a single magnetar burst.

3.1. SGR J0501+4516

3.2. SGR J1550−5418

4.1. Bayesian Blocks Profile

The Bayesian blocks representations of the temporal profiles of magnetar bursts are much simpler compared with their light curves. Figure 1 shows light curves and corresponding Bayesian blocks profiles of bursts from SGR J0501+4516, representing six common types that emerge from our analysis. The burst light curves are binned with different bin sizes for clarity. The majority of magnetar bursts can be represented with only one block. The burst block is either a short peak less than half a second (SS: single short) or otherwise a broad bump over the local level of the persistent emission (SL: single long). Some bursts have at least two subsequent blocks in their profiles, and all of them are short ones (i.e., <0.5 s, MS: multiple short). The rest of the bursts have at least one long block (>0.5 s). We divided them into three groups: DSL (double short+long) bursts are two block events that have a single peak with extended emission; DLS (double long+short) bursts also have two blocks, but with an enhanced emission leading to the single peak; and MO (multiple others) events include all other cases. We summarized the basic morphological properties, as well as the distribution of burst patterns among the two sources, in Table 1.

Zoom In Zoom Out Reset image size

Table 1. A Summary of Six Types of Bayesian Blocks Profiles of Magnetar Bursts

Morphology	SS	SL	DSL	DLS	MS	MO	N a
SGR J0501+4516 XMM-Newton	201(76.7%)	40(15.3%)	12(4.6%)	2(0.8%)	6(2.3%)	1(0.4%)	320
SGR J1550−5418 00340573000 orbit1	32(71.1%)	8(17.8%)	1(2.2%)	0	2(4.4%)	2(4.4%)	59
SGR J1550−5418 00340573000 orbit2	50(47.6%)	18(17.1%)	9(8.6%)	1(1%)	16(15.2%)	11(10.5%)	152
SGR J1550−5418 00340573000 orbit3	31(75.6%)	7(17.1%)	0	0	3(7.3%)	0	62
SGR J1550−5418 00340573000 orbit4	2(100%)	0	0	0	0	0	28
SGR J1550−5418 00340573001	65(79.3%)	6(7.3%)	5(6.1%)	0	5(6.1%)	1(1.2%)	103
SGR J1550−5418 all	180(65.5%)	39(14.2%)	15(5.5%)	1(0.4%)	26(9.5%)	14(5.1%)	404

Note. ^aN: total number of bursts (including piled-up bursts).

Download table as:  ASCIITypeset image

4.2. Duration

In earlier investigations of the temporal properties of magnetar bursts, burst durations were estimated by fitting the cumulative count or energy fluence light curves with the sum of a step function and a first-order polynomial. In this formalism, the step function provides the amplitude of the burst, and the polynomial represents the behavior of the background emission (e.g., Göǧüş et al. 2001; Lin et al. 2011). This method, however, leaves room for uncertainties due to human intervention. In particular, the choice of background intervals before and after the burst might be subjective and, in turn, might result in slightly different duration estimates. In this study, we define the duration of a burst as the time length of representing Bayesian blocks: the duration of a burst is the interval from the start of the first burst block to the end of the last one. Since the change point (boundary) between the background and burst blocks is settled with a certain threshold (Scargle et al. 2013), the duration calculation based on the Bayesian blocks profile is simple, direct, and not influenced by any subjective background selections. Note here the fact that due to difficulties arising from background fluctuations, T₉₀ or T₅₀ durations have been conventionally calculated for bursts (Kouveliotou et al. 1993). The Bayesian burst blocks, however, represent the entire burst durations.

4.2.1. SGR J0501+4516

We present the duration distribution of all 262 non-piled-up SGR J0501+4516 bursts in the left panel of Figure 2. The apparent deviation from a log−Gaussian shape is also reflected in the quantitative case: a fit to the duration distribution of SGR J0501+4516 with a Gaussian function yields a χ² of 41.14 for 8 degrees of freedom (dof). We then fit the distribution with the sum of two Gaussian functions, which resulted in a very good fit (χ²/dof = 3.65/5). We estimate the chance probability of the preference of the more complicated model to be 1 × 10⁻⁸. The best-fit mean values of the two Gaussian components are 85 ± 8 ms with σ = 0.36  ±  0.03 and $1028^{+220}_{-181}$ ms with σ = 0.32 ± 0.05 (both σ are in decades). The distribution reaches a local minimum value at 456 ms between the two peaks.

Zoom In Zoom Out Reset image size

The second Gaussian function is dominated by the long blocks, namely, SL bursts and long blocks of DSL, DLS, and MO events. We excluded those bursts/blocks and generated the distribution of only short bursts/blocks, as shown in the right panel of Figure 2. This distribution is fit well with a single Gaussian function (χ²/dof = 9.15/6), and the best-fit mean value is 86 ± 5 ms, with σ = 0.37 ± 0.02. As expected, this set of parameters agrees well with that of the Gaussian component with the smaller mean value.

4.2.2. SGR J1550−5418

We studied the duration distribution of SGR J1550−5418 using 275 bursts. The distribution (the left panel of Figure 3) can be well fit with a Gaussian function in the logarithmic space (χ²/dof = 9.83/8) with a mean value of $207^{+16}_{-15}$ ms and σ = 0.52  ±  0.02. Unlike the distribution of SGR J0501+4516, a second Gaussian function is not required here. We also generated the duration distribution of 245 short bursts/blocks, including SS and MS, and short blocks of DSL, DLS, and MO bursts (23 short bursts are contributed by MO events). This distribution can also be described with a Gaussian function, as shown in the right panel of Figure 3. The best-fit mean value is 131 ± 8 ms with σ = 0.37 ± 0.02 (χ²/dof = 19.10/8).

Zoom In Zoom Out Reset image size

5.1. SGR J0501+4516

In Lin et al. (2012b), we studied the stacked spectrum of 129 relatively dim SGR J0501+4516 bursts. With a larger sample from the Bayesian blocks burst search, we are now able to look into the spectral behavior of the weaker bursts, grouped by morphological type and fluence. For this purpose, we calculated a signal-to-noise ratio (S/N) for each block as $(N_T-N_B)/\sqrt{N_T}$, where N_T and N_B are the total and background counts collected within the block, respectively. This value is a general indicator of the burst fluence and ranges between 1.6 and 6.3. Then, we further divided SS bursts into three groups with S/N >4, 3 < S/N ⩽4, and S/N ⩽3. The SL bursts were separated into two groups with an S/N boundary of 4. For DSL and DLS events, we merged all of their short blocks into one spectral group and the long blocks into another group. We kept MS bursts in one group and excluded the MO event since there was only one such burst. Overall, we have formed eight spectral groups of bursts with a sufficient number of counts, different Bayesian block profiles, and energetics and extracted a stacked spectrum for each group.

We fit each stacked spectrum with the MBB+RCS model. We set the interstellar absorption and magnetospheric scattering parameters of the plasma at the values obtained from the persistent emission spectral analysis (n_H = 6.7 × 10²¹ cm⁻², β = 0.37, τ = 5; Lin et al. 2012). We found that all of the stacked burst spectra can be fit well with the MBB+RCS model, for which we listed the best-fit model parameters in Table 2. We noticed that the unabsorbed flux in 0.5–10 keV is correlated with the temperature of the modified blackbody (see Figure 4). The Spearman's rank-order correlation coefficient is 0.98, with a random occurrence probability of 3.3 × 10⁻⁵. The correlative trend between the flux and temperature can be expressed as a power law with Flux ∝kT^{3.8 ± 0.3} (χ²/dof = 4.2/6; as shown by the solid line in Figure 4). The power-law index becomes flatter (Flux∝kT^{2.5 ± 1.0} with χ²/dof = 1.3/3) when we only include spectral results of SS and MS bursts and short blocks of DSL and DLS events.

Zoom In Zoom Out Reset image size

Table 2. Spectral Fit Results

Type–(S/N)	Exposure	kT	Fluxa	χ2/dof	Burst Energy per Count
	(s)	(keV)	(10−10 erg s−1 cm−2)		(10−12 erg count−1)
SGR J0501+4516
MBB+RCS, nH = 6.7 × 1021 cm−2, β = 0.37, τ = 5
SS−(⩽ 3)	4.556	$0.84_{-0.04}^{+0.06}$	9.0 ± 0.5	44.64/37	7.4
SS−(3 ∼ 4)	6.391	$0.95_{-0.04}^{+0.06}$	10.3 ± 0.4	61.93/51	7.7
SS−(> 4)	4.949	$1.08_{-0.05}^{+0.06}$	$16.2^{+0.7}_{-0.6}$	44.17/57	8.3
SL−(⩽ 4)	22.91	$0.54_{-0.04}^{+0.06}$	1.0 ± 0.1	34.10/34	6.6
SL−(> 4)	15.33	0.57 ± 0.04	1.8 ± 0.1	24.16/31	6.9
DSL, DLS−short	1.399	$1.00^{+0.12}_{-0.10}$	$12.4^{+1.2}_{-1.1}$	4.29/12	8.7
DSL, DLS−long	18.35	$0.51^{+0.06}_{-0.07}$	0.8 ± 0.1	22.27/21	6.9
MS	1.662	$0.92^{+0.16}_{-0.12}$	6.7 ± 0.8	4.13/7	9.1
STEMS, nH = 6.7 × 1021 cm−2, β = 0.37, τ = 5, B = 2.2 × 1014 G
SL−(⩽ 4)	22.91	$0.35_{-0.03}^{+0.07}$	1.1 ± 0.1	30.71/34	⋅⋅⋅
SL−(> 4)	15.33	$0.42_{-0.02}^{+0.03}$	1.8 ± 0.1	29.11/31	⋅⋅⋅
DSL, DLS−long	18.35	$0.40^{+0.02}_{-0.08}$	0.8 ± 0.1	23.31/21	⋅⋅⋅
SGR J1550−5418 XRT observation ID 00340573000
MBB+RCS, nH = 3.4 × 1022 cm−2, β = 0.2, τ = 3
orbit1, SS	12.71	$2.64^{+0.22}_{-0.19}$	$132.4^{+5.3}_{-5.2}$	80.71/75	142.2
orbit2, SS	8.121	$3.93^{+0.61}_{-0.45}$	$224.2^{+10.3}_{-10.0}$	85.72/72	160.1
orbit3, SS	9.382	$2.08^{+0.18}_{-0.15}$	$112.5^{+5.4}_{-5.2}$	44.43/60	133.6
orbit1, SL	6.254	$2.24^{+0.38}_{-0.28}$	$67.0^{+5.5}_{-5.2}$	16.88/19	147.3
orbit2, SL	15.77	$2.26^{+0.21}_{-0.18}$	$58.6^{+3.0}_{-2.9}$	46.01/47	136.2
orbit3, SL	5.985	$1.52^{+0.18}_{-0.16}$	62.3 ± 4.6	19.51/26	122.9
orbit2,DSL, DLS short	1.456	$3.38^{+0.76}_{-0.51}$	$370.3^{+29.3}_{-28.1}$	11.18/20	165.3
orbit2,DSL, DLS long	10.04	$2.12^{+0.50}_{-0.33}$	$32.2^{+3.6}_{-3.2}$	21.65/16	143.9
orbit2, MS	5.558	$5.49^{+2.02}_{-1.04}$	$286.4^{+16.0}_{-15.4}$	41.36/58	176.0
orbit3, MS	1.002	$2.65^{+1.03}_{-0.55}$	$218.3^{+29.9}_{-26.4}$	11.80/8	160.4
orbit2, MO short	4.463	$3.47^{+0.53}_{-0.39}$	$256.2^{+13.8}_{-13.3}$	41.23/47	157.7
orbit2, MO long	16.74	$1.80^{+0.33}_{-0.26}$	$16.7^{+1.7}_{-1.6}$	13.56/18	132.7
SGR J1550−5418 XRT observation ID 00340573001
MBB+RCS, nH = 2.8 × 1022 cm−2, β = 0.2, τ = 3
SS	8.916	$2.66^{+0.20}_{-0.17}$	$179.8^{+6.9}_{-6.8}$	59.39/78	135.7
SL	4.565	$2.05^{+0.67}_{-0.39}$	$43.7^{+6.3}_{-5.4}$	8.38/11	132.7
DSL long	5.438	$1.52^{+0.35}_{-0.24}$	$30.5^{+3.7}_{-3.5}$	16.08/11	115.0
MS	2.234	$3.60^{+2.62}_{-0.93}$	$142.3^{+23.1}_{-18.8}$	8.32/12	168.3

Note. ^aAverage unabsorbed flux in the 0.5–10 keV energy range.

Download table as:  ASCIITypeset image

We also fit the stacked spectra with the STEMS model, which is preferred by the persistent emission spectral shape. We found that the stacked spectra of SL bursts and long blocks of DSL and DLS bursts can also be fit well with the STEMS model, while other types cannot be fit with it. We listed the STEMS model parameters in Table 2. The temperatures of long block bursts are consistent with one another (∼0.4 keV) within errors. Additionally, we showed the temperature of the hot spot on the neutron star surface as obtained by fitting its persistent emission from the same observation (Lin et al. 2012b) with a five-point star in Figure 4.

5.2. SGR J1550−5418

The persistent source properties of SGR J1550−5418 changed significantly during its burst-active episodes. The source flux (count rate) in observation 00340573000 varied by a factor of ∼10 from the first to last orbit. Therefore, instead of averaging the whole observation, we extracted the burst and persistent emission spectra for each orbit separately for 00340573000. Since the source flux did not show any dramatic change in the later observation and the number of bursts is relatively small, we extracted spectra from the entire observation of 00340573001. We investigated burst spectra of SGR J1550−5418 in morphological groups only because there are not enough bursts in each orbit to allow grouping in fluence. We listed all stacked spectra extracted and analyzed in Table 2.

Due to the existence of the dust scattering halo (Tiengo et al. 2010; Vink & Bamba 2009; Ng et al. 2011), the choice of the interstellar absorption value for this source was not trivial. Therefore, we first performed a simultaneous fit to the stacked spectra of all types of bursts in both observations with the MBB+RCS model and forced them to have the same common n_H. Ng et al. (2011) studied the persistent emission spectra collected before, during, and after the burst active periods using the STEMS model. Considering the fact that the magnetospheric properties would not change over the short timescales of bursts, we fixed the scattering parameters in the MBB+RCS model at their resulting values for the burst-active episode: the average speed of the electrons β ∼ 0.2 and the optical depth τ ∼ 3.0 (Ng et al. 2011). We obtained the best-fit n_H as (3.4 ± 0.2) × 10²² cm⁻² (χ²/dof = 434.67/466) and (2.8 ± 0.3) × 10²² cm⁻² (χ²/dof = 92.16/111) for observations 00340573000 and 00340573001, respectively. Note that there are three different n_H values for SGR J1550−5418 reported in the literature using fine spatial resolution XMM-Newton and Chandra observations: (3.1 ± 0.3) × 10²² cm⁻² (Halpern et al. 2008), (2.8 ± 0.3) × 10²² cm⁻² (Vink & Bamba 2009), and (4.5 ± 0.1) × 10²² cm⁻² (Ng et al. 2011). Our results fall within the range of reported n_H estimates. Therefore, we performed our detailed spectral analysis with the MBB+RCS model for each stacked spectrum using a galactic absorption value of 3.4 × 10²² cm⁻² for observation 00340573000 and 2.8 × 10²² cm⁻² for the other observation. The best-fit parameters and the fit statistics are presented in Table 2. The temperature of the modified blackbody is spread between 1.5 keV and 5.7 keV, which is much higher than those of SGR J0501+4516, and the average unabsorbed flux in 0.5–10 keV of bursts is also much larger.

We also investigated the correlative trend between the average unabsorbed flux and the temperature of the modified blackbody for all stacked spectra (Figure 5). In general, these two quantities are well correlated (the Spearman's rank correlation coefficient is 0.839, and the chance occurrence probability is 4.9 × 10⁻⁵). This correlation can also be described with a single power-law trend, having an index of 2.1 ± 0.4 (χ²/dof = 23.82/14), which agrees with that from the fit using only short blocks/bursts from SGR J0501+4516.

Zoom In Zoom Out Reset image size

The single blackbody (BB) function has also been used to describe the spectral shape of SGR J1550−5418 bursts in the 8–200 keV band (von Kienlin et al. 2012). To test this model, we fit the stacked spectra with a BB in a similar manner as in the MBB+RCS case. First, we fit all types of burst spectra in both observations simultaneously and obtained a galactic absorption of (2.6 ± 0.2) × 10²² cm⁻² (χ²/dof = 447.01/466) and (2.1 ± 0.3) × 10²² cm⁻² (χ²/dof = 91.48/111). Then, we fit each stacked spectrum with the n_H set at their corresponding values. From the statistical point of view, the BB model fits to the stacked burst spectra are acceptable. However, the n_H values obtained are much lower than those reported earlier in the literature, as well as those obtained with the MBB+RCS model fit. Moreover, even if these lower n_H are considered to be the case, the resulting total χ² value from the BB model fits to all stacked spectra is larger by 11.7 than that from the MBB+RCS model, which has the same number of free parameters. Therefore, we concluded that the MBB+RCS model is preferred over the BB model for the stacked burst spectra of SGR J1550−5418 in the 0.3–10 keV band.

In order to compare the shape of the stacked spectra with that of the persistent source, we likewise fit them with the STEMS model as the sum of a STEMS and a single power-law function, which successfully described the persistent spectra in the same burst-active episode (Ng et al. 2011). None of our stacked spectra can be fit with the STEMS model or even with an additional power law.

6.1. SGR J0501+4516

To study the differential fluence distribution of SGR J0501+4516, we focused on SS and MS bursts and short blocks of DSL and DLS bursts, since these types have been shown to exhibit a distinct burst nature (Section 5.1). We estimated the fluence of each burst using a conversion factor, which is the average of the total energy-to-count ratio for the selected groups of bursts (listed in the rightmost column of Table 2). The average conversion factor of the stacked spectra of short bursts from SGR J0501+4516 is 8.2 × 10⁻¹² erg count⁻¹. In Figure 6, we exhibit the differential fluence distribution of 221 bursts in equal logarithmic steps.

Zoom In Zoom Out Reset image size

To further study the fluence distribution of SGR J0501+4516 bursts, we first examined the instrumental effects that might limit the detection sensitivity. The high fluence cutoff is determined by the pileup level. For this XMM-Newton observation, we excluded piled-up bursts that have more than 50 counts in a second. All 221 short bursts in our distribution are below this limit, which corresponds to an energy fluence of 4.1 × 10⁻¹⁰ erg cm⁻². We determined the lower fluence cutoff through numerical simulations as follows: we assumed that a flat differential fluence distribution and the burst counts are evenly distributed throughout the burst interval. We performed simulations on a series of total burst counts ranging from 2 to 50. For each total count value examined, we iterated the following steps 10,000 times: first, we randomly selected a burst duration drawn from the distribution of short bursts in the right panel of Figure 2. Then, we generated a background-included event list by adding a 20 s long segment of the observed persistent emission of SGR J0501+4516 with the simulated burst events with the specified duration. These two event lists were aligned at their mid-times. Finally, we applied the Bayesian blocks representation process to the simulated event list and determined whether or not the simulated burst can be detected. We found from these simulations that the detection probability of bursts with at least 12 total counts (corresponding to 9.9 × 10⁻¹¹ erg cm⁻²) is higher than 99%. Therefore, we set the lower limit at this fluence level.

We have 97 short bursts in the fluence range where the instrument and our burst identification process have complete detection capabilities. We fit the differential burst fluences in this fluence range with a single power law and determined the index to be −2.4 ± 0.3 (χ²/dof = 5.9/3). In order to avoid binning effects, we also applied the maximum likelihood estimation to determine the power-law index from the unbinned burst fluences and obtained −2.3 ± 0.3 with 1σ errors.

6.2. SGR J1550−5418

We studied the burst fluence distribution of SGR J1550−5418 using 273 events from all types of Bayesian blocks profiles, as they very clearly represent burst origin from spectral investigations (Section 5.2). Here, we have already excluded the only two bursts from the fourth orbit of observation 00340573000, which suffered from a much higher persistent flux (∼64.2 counts s⁻¹) than the other orbits. Similar to SGR J0501+4516, we calculated the fluence of each burst by converting the total burst counts into the corresponding energy. If a burst has been grouped into any stacked spectrum, then we used the spectral fit results to determine the conversion factor (see the last column of Table 2). Otherwise, we adopted the average value of the conversion factors in its orbit or observation, that is, 1.4 × 10⁻¹⁰ erg count⁻¹, 1.5 × 10⁻¹⁰ erg count⁻¹, and 1.4 × 10⁻¹⁰ erg count⁻¹ for orbits 1, 2, and 3 of observation 00340573000, and 1.4 × 10⁻¹⁰ erg count⁻¹ for observation 00340573001. In Figure 7, we present the size distribution of these 273 bursts.

Zoom In Zoom Out Reset image size

Following the same procedure for SGR J0501+4516, we also calculated the fluence limits of bursts with a detection rate over 99% for SGR J1550−5418. The higher fluence limit is 1.5 × 10⁻⁸ erg cm⁻², corresponding to an energy fluence of 100 counts (i.e., the pileup limit) using the average conversion factors derived from all investigated orbits. Due to the significant change of the persistent flux level, we performed four sets of simulations for the first three orbits in observation 00340573000 and observation 00340573001. We obtained the minimum burst count with over 99% detection rate as 42, 28, 39, and 33, respectively. After converting them into the energy fluence, we adopted the highest value of 5.9 × 10⁻⁹ erg cm⁻² as the lower fluence limit. Therefore, we selected 68 bursts with fluence in between these two limits. We fit the distribution within this range with a single power-law function and obtained a best-fit index of −1.7 ± 0.5 (χ²/dof = 2.81/3). We also calculated the power-law index using the maximum likelihood method, which yields the most probable index with a 1σ confidence range as −1.7 ± 0.5.

In this paper, for the first time, we introduced the Bayesian blocks method to search for weak magnetar bursts in X-ray observations using time-tagged event data. We successfully identified 320 bursts from SGR J0501+4516 observed with XMM-Newton and 404 bursts from SGR J1550−5418 in two Swift/XRT observations. The Bayesian blocks burst search process takes advantage of fine temporal resolution of the instruments, without running into any limitations introduced by binning the data. The Bayesian blocks method is effective in a wide range of burst energetics, especially when identifying weak bursts from the fluctuations of the persistent source emission. The unabsorbed fluence in the 0.5–10 keV band of bursts identified with this method ranges from 7.4 × 10⁻¹² to 3.9 × 10⁻¹⁰ erg cm⁻² for SGR J0501+4516 and from 4.0 × 10⁻¹⁰ to 1.9 × 10⁻⁸ erg cm⁻² for SGR J1550−5418. In general, these bursts are dimmer by about two orders of magnitude than those that had been extensively analyzed in the literature: ∼10⁻⁸ erg cm⁻² in 8–200 keV (Lin et al. 2011) and ∼10⁻¹⁰ erg cm⁻² in 2–40 keV (Nakagawa et al. 2011) of SGR J0501+4516 bursts; 10⁻⁸ erg cm⁻² in 0.3–300 keV (Enoto et al. 2012) and >10⁻⁸ erg cm⁻² in 8–200 keV (van der Horst et al. 2012; von Kienlin et al. 2012) of bursts from SGR J1550−5418. We note that these fluence values were accumulated in different energy ranges. However, the broadband spectral analysis from 0.5 keV to 200 keV showed that the fluence in the 8–200 keV energy range is only about an order of magnitude larger than that in the 0.5–10 keV band (Lin et al. 2012a). Therefore, after calibrating for the energy ranges, our samples from the Bayesian blocks burst search technique still represent the dimmest burst population from magnetars.

The Bayesian blocks representation provides a direct way to investigate the burst temporal morphology. A large number of burst counts accumulated into a single block, based on statistics, remarkably simplifies the burst light curve. Unlike the burst search performed with the binned data, this method can also reveal statistically significant temporal structures, which would otherwise be smoothed out in the binned light curve. In our sample for both magnetars, 70.95% of the burst profiles are represented by a single short block (SS); 5.96% of bursts have more than one but all short blocks (MS); 14.71% of events show a stable flux for longer than half a second (SL); 5.59% of burst profiles have two blocks, and a short block is either followed or preceded by a lower flux structure (DSL and DLS); and the remaining 2.79% of events exhibit more than one short peak or long smooth emission profile (MO). We show that except for the second orbit of observation 00340573000, which covered the burst storm of SGR J1550−5418, all other observations and orbits exhibit a comparable percentage of different burst types (see Table 1).

Based on the Bayesian blocks profiles, we can also directly determine the duration of a burst from the block size without running into difficulties from background modeling, which can be severe for dim bursts. We can therefore study the statistical behavior of the temporal properties of the dimmest magnetar bursts. The burst duration distribution of SGR J0501+4516 in 0.2–10 keV follows a double-Gaussian trend with two components dominated by short and long bursts. The mean value of the duration of the short bursts as obtained with a Gaussian function fit is ∼87 ms. This is shorter than the average T₉₀ duration of SGR J0501+4516 bursts, while it is comparable to their average emission times (τ₉₀), which is the duration needed to collect 90% of the burst total counts regardless of the order of their arrival times in the 8–200 keV energy range (Lin et al. 2011). For typical single pulse bursts, the average T₉₀ duration is approximately equal to the average τ₉₀ of all bursts (Göǧüş et al. 2001). It is therefore expected that the distribution of our short blocks sample of SGR J0501+4516 will be consistent with the T₉₀ duration distribution of bursts with a simple profile. For SGR J1550−5418, the distribution of all events follows a single Gaussian-like trend, similar to that of more energetic bursts seen with Fermi/GBM (van der Horst et al. 2012; von Kienlin et al. 2012). These agreements support the burst nature of short blocks from SGR J0501+4516 and all events from SGR J1550−5418.

We performed detailed spectral analysis of the weakest magnetar bursts, grouped into eight stacked spectra for SGR J0501+4516 and sixteen for SGR J1550−5418 based on their temporal profiles and fluences. All of the stacked spectra can be fit adequately with the MBB+RCS model. To investigate the connection between the burst and persistent emission, we also studied the accumulated burst spectra with the STEMS model, which was developed to better understand the persistent emission spectra of magnetars. We found that the STEMS model can fit only the stacked spectra of SL and long blocks of DSL and DLS events of SGR J0501+4516 (54 out of 320 bursts). As we discussed above, these types of bursts accounted for an additional Gaussian component in the combined duration distribution of SGR J0501+4516 bursts. Based on these, we cannot conclusively determine whether these weak and long blocks of SGR J0501+4516 are indeed representing bursts or rapid enhancements of the persistent source emission. The stacked spectra of all other types of SGR J0501+4516 bursts and entire SGR J1550−5418 events cannot be fit with the STEMS model, and therefore they likely represent a distinct origin other than persistent emission. Consequently, we strengthened the conclusion of Lin et al. (2012b), with much dimmer magnetar bursts, that the emergent spectra of the burst and persistent emission are characteristically different, even though their origins likely involve an underlying extremely strong magnetic field.

We found that SGR J1550−5418 bursts detected with XRT have, on average, a higher modified blackbody temperature and unabsorbed flux than those of SGR J0501+4516 observed with XMM-Newton. However, we obtained a correlation between the flux and temperature for both sources. In particular, power-law trends with indices of 2.5 ± 1.0 and 2.1 ± 0.4 are followed between these two quantities using only spectrally confirmed bursts of SGR J0501+4516 and SGR J1550−5418, respectively. Note that both power-law indices have large deviations from the expectation of a pure blackbody model (i.e., 4), which may indicate the signature of the strong magnetic fields on the emergent burst emission. We also found that for SGR J0501+4516, the surface temperature and flux of the persistent emission (obtained with the STEMS model) cannot be inferred from the extrapolation of the power-law trend of either all burst types or only the short burst blocks. It is important to note that the temperatures from the STEMS and MBB+RCS models arise from different physical settings: the former is the temperature of the thermal emission from the surface of the neutron star, and the latter is the temperature of the thermalized bubble emission modified by the strong magnetic field.

The persistent emission flux of SGR J1550−5418 changed dramatically during its burst active episode in 2009 (Scholz & Kaspi 2011). The X-ray flux of the source reached its peak in 2009 January 22, the most burst-active day. Using the burst sample we identified with the Bayesian blocks method, we compared the evolution of the persistent flux in conjunction with those of the burst number and the orbital averaged burst energy fluence (see Figure 8). We find that the burst rate peaks in the second XRT orbit, which is about 11 ks prior to the peak of the persistent flux (count rate) in the fourth orbit. This behavior indicates that the sudden release of a large amount of burst energy may ignite the heating up of the neutron star crust and result in the increase of the temperature and flux of the persistent emission, which were observed by Scholz & Kaspi (2011). Importantly, we note that the rate of strong (piled-up) bursts also peaked in the fourth XRT orbit. Therefore, it may also be possible that the strong bursts in the fourth orbit may have driven simultaneous persistent flux increase, rather than gradual surface heating.

Zoom In Zoom Out Reset image size

We studied the fluence distribution of SGR J0501+4516 bursts using the dimmest sample in the 0.5−10 keV band for the first time. The distribution follows a power-law trend as dN∝S^−2.3dS. Lin et al. (2011) formed the fluence distribution of SGR J0501+4516 bursts using 29 events seen with GBM in the 8−200 keV range and obtained a power-law index of −1.5, which is much flatter than the trend at the lower energy level. This may indicate that there is a break somewhere between the fluences of 4 × 10⁻¹⁰ and 10⁻⁷ erg cm⁻². We have further investigated this possibility as follows. There were 10 bursts detected from SGR J0501+4516 with GBM during the course of XMM-Newton observation (bursts $\#5\hbox{--}14$ in Table 1 of Lin et al. 2011). Note that Fermi is in a low-altitude orbit; therefore, the source was blocked by the Earth by about 50% of the exposure of XMM-Newton. Assuming that the source was emitting bursts at a constant rate (which is implied from the XMM-Newton observation), we considered that about 20 bursts must have been detected with GBM over the time that 320 bursts were recorded with XMM-Newton. We then calibrated our fluence distribution using the ratio of the number of bursts detected by two instruments and extrapolated the distribution to the higher fluence regime to estimate the number of bursts expected if there is no break (namely, the index of −2.3 spans over all fluences). We found that the expected number of bursts from the power-law distribution with an index of −2.3 is about two orders of magnitude lower than the number of bursts detected with GBM. This result strengthened the possibility that there might indeed be a break in the broad fluence distribution of SGR J0501+4516.

Scholz & Kaspi (2011) constructed the fluence distribution of SGR J1550−5418 bursts in the 0.5−10 keV energy band and obtained a power-law index of −1.6, which is in agreement with our investigations here (−1.7 ± 0.5). A similar power-law index was reported for the fluence distribution of SGR J1550−5418 bursts within the fluences of 1.8 × 10⁻⁷ to 3.2 × 10⁻⁶ erg cm⁻² (van der Horst et al. 2012). Note that the upper fluence limit of our burst sample is nearly equal to the detection limit of GBM bursts. Therefore, there is clearly no evidence for any break in the fluence distribution of SGR J1550−5418 over three orders of magnitude in fluence.

According to the magnetar model, typical short bursts result from the release of crustal energy, which is triggered by the cracking of the neutron star crust by the strong magnetic field (Duncan & Thompson 1992; Thompson & Duncan 1995). The required crustal magnetic field strength can be estimated from the energy of a burst, the length of the crustal crack, and the upper limit of the lattice strain, which is typically 10⁻³ (Thompson & Duncan 1995). Moreover, one could also estimate the displacement of the magnetic footprints and the minimum excitation radius of the Alfvén wave by using the crack length and the magnetic field strength (Thompson & Duncan 1995). Our earlier STEMS model fit to the persistent emission of SGR J0501+4516 resulted in an average surface magnetic field strength of 2.2 × 10¹⁴ G. Using this value, we obtain a crack length scale of ∼10⁻¹ m for the weakest bursts from this magnetar (that is, SS bursts with S/N ⩽3). We estimate the corresponding displacement of the magnetic footprints as ∼2 × 10⁻² m and the minimum Alfvén excitation radius as ∼5 m. Note that the surface magnetic field topology of the neutron star is not necessarily uniform and the field strength may exceed the average value under more localized settings. Therefore, our estimates for the length scales are the lowest limits for SGR J0501+4516 bursts.

L.L. is funded through the Post-Doctoral Research Fellowship of the Turkish Academy of Sciences (TÜBA). L.L. thanks her grandmother Mrs. Liu Shufen for her endless love and support.