Energy and waiting time distributions of FRB 121102 observed by FAST

FRB 121102 has been observed by different telescopes (Spitler et al., 2014, 2016; Scholz et al., 2016, 2017; Chatterjee et al., 2017; Law et al., 2017; Hardy et al., 2017; Zhang et al., 2018; Spitler et al., 2018; Gajjar et al., 2018; Michilli et al., 2018; Gourdji et al., 2019; Hessels et al., 2019; Oostrum et al., 2020; Caleb et al., 2020; Cruces et al., 2021). Recently, FAST observed this source at 1.25 GHz (Li et al., 2021a). The observations were carried out from August 2019 (MJD 58724) to October 2019 (MJD 58776), a duration lasting nearly half of the 157d period suggested by Rajwade et al. (2020) and covering almost the entire active phase (MJD 58717 to 58813) predicted by Cruces et al. (2021). During the 59.5 observing hours, FAST detected 1652 bursts. The mean burst rate was about 27.8 hr^-1, and the peak burst rate was 122 hr^-1. These bursts help to understand the energy distribution in both the low-energy and high-energy ranges. The waiting time distribution can be studied in detail with this sample. We will use this sample to derive a more precise waiting time distribution.

The energy distribution of FRB 121102 has been investigated by many previous works (Wang & Yu, 2017; Law et al., 2017; Gourdji et al., 2019; Wang & Zhang, 2019; Cheng et al., 2020; Oostrum et al., 2020; Cruces et al., 2021; Li et al., 2021a). We use the bursts of FRB 121102 observed by FAST to investigate the differential energy count distribution

The energies of these bursts are listed in Li et al. (2021a). They use central frequency rather than bandwidth to calculate the energies of bursts (This is to alleviate the bandwidth limitation of telescopes, see Zhang, 2018). The luminosity distance of FRB 121102 is taken as 949 Mpc (Planck Collaboration & et al., 2016; Tendulkar et al., 2017). We show the differential energy count distribution as a blue histogram in Figure 1. The distribution in the low energy range deviates from the power-law function. Therefore, to compare with previous works, we only consider the high-energy bursts ($E>10^{38}$ erg), which contains roughly 600 bursts. The distribution can be approximated as a power-law function. We use equation (1) to fit the energy distribution in the high-energy end. The fitting result is shown as the green dot-dashed line in Figure 1. We also show the 90% threshold of FAST (Li et al., 2021a) as vertical yellow dashed line in this figure. The power-law index is $\alpha_{E}=1.86\pm 0.02$ (1$\sigma$ uncertainty), which is consistent with the index derived by Wang & Yu (2017). Law et al. (2017) obtained a power-law index of $\alpha_{E}\simeq 1.7$, which is also consistent with our result. Wang & Zhang (2019) used a cut-off power-law model to fit the FRB 121102 data from multiple observations and derived a universal energy distribution with $\alpha_{E}=1.6-1.8$ for multiple observations at varied frequencies. Our result also closes to this range. The energy distribution of high-energy bursts was also investigated by Li et al. (2021a). Because the observation of FAST for FRB 121102 consists of many single-day observations, Li et al. (2021a) used various observational times as the weight to derive the differential burst rate $dR/dE$ at different energies and found a power-law index of about $-1.85\pm 0.3$. Our result is the distribution of burst counts in each energy bin $dN/dE$, not the burst rate distribution $dR/dE$.

Interestingly, the energies of non-repeating FRBs also show a power-law distribution with $\alpha_{E}=1.6-1.8$ (Cao et al., 2018; Lu & Piro, 2019; Zhang et al., 2019, 2021; Luo et al., 2018, 2020). Some sources also show a power-law distribution of $dN/dE$, such as X-ray bursts of magnetars (e.g. $\alpha_{E}=1.68\pm 0.01$ for SGR 1806-20, Göǧüş et al., 2000; Prieskorn & Kaaret, 2012; Cheng et al., 2020; Yang et al., 2021), X-ray flares of gamma-ray bursts ($\alpha_{E}=1.06\pm 0.15$, Wang & Dai, 2013), M87 ($\alpha_{E}=1.69^{+0.59}_{-0.45}$, Yang et al., 2019) and Sgr A^∗ ($\alpha_{E}=1.65\pm 0.02$, Wang et al., 2015; Li et al., 2015). Some FRB theoretical models are related to giant pulses of pulsars (e.g. Cordes & Wasserman, 2016). The giant pulses of Crab pulsar also show a power-law distribution ($\alpha_{E}=-1.85\pm 0.10$, Popov & Stappers, 2007; Lyu et al., 2021). A power-law distribution of energy is a natural predication of self-organized criticality theory (Bak et al., 1987; Aschwanden, 2011).

The energy distribution in the low-energy range deviates from this power-law. This deviation was also reported by Gourdji et al. (2019). They attributed this deviation to the incompleteness at low energies. In our sample, the turning point is about $10^{38}$ erg, much larger than the threshold energy, which is about $3\times 10^{37}$ erg (Li et al., 2021a). Therefore, this deviation is physical. A similar deviation was also found in other samples (Wang & Zhang, 2019; Oostrum et al., 2020; Cruces et al., 2021). Li et al. (2021a) also investigated the energy distribution of FRB 121102. They found the bimodal energy distribution, which is the Cauchy function with $\alpha_{E}=1.85\pm 0.3$ for $E>10^{38}$ erg plus the log-normal function with $\sigma_{E}=0.52,N=2.06\times 10^{38},E_{0}=7.2\times 10^{37}$ erg for $E<10^{38}$ erg. They also give the KDE of energy and burst time in the panel (c) of their Figure 1. In this panel, the high-energy component is only visible in early observations, while the low-energy peaks are present in both early and late observations. This novel structure may suggest different physical mechanisms for high-energy and low-energy bursts. If this conjecture is true, it can explain the deviation in low energies. The distribution also shows an obvious steepening above $E>5\times 10^{39}$ erg, which indicates the burst maximum energy $E_{\rm max}$ is around this point.

Wang & Yu (2017) found that the time interval between consecutive bursts does not follow a Poission distribution. Afterwards, Oppermann et al. (2018) found that the Weibull function is a better description of waiting time, which is described by

where $\delta_{t}$ is the waiting time, $k$ is the shape parameter, $r$ is the mean burst rate, and $\Gamma$ is the gamma function. The case $k=1$ corresponds to the Poisson distribution. Previous works suggested that the repeating behavior of FRB 121102 tends to have $k<1$ (Oppermann et al., 2018; Oostrum et al., 2020; Cruces et al., 2021), which means that time clustering is favored. We use the Weibull function to fit the burst behavior. The observations of FAST lasted for roughly 1 hour per day. Therefore, waiting time greater than a half day is caused by the observing window, which is ignored in our analysis. The waiting time with $\delta_{t}<30$ ms is also excluded. It is difficult to determine whether they are two separate bursts or multiple peaks of a single burst (Cruces et al., 2021). We use the MCMC (Markov chain Monte Carlo) technique to fit this distribution. The corner plot of the best-fitting results is shown in Figure 2. We list the best-fitting parameters in Table 1, alongside the results from previous works. We find $k=0.72^{+0.01}_{-0.01}$ (1$\sigma$ uncertainty) and $r=736.43^{+26.55}_{-28.90}$ day^-1 (1$\sigma$ uncertainty). This burst rate is much higher than the rates derived in previous works, which is caused by the high sensitivity of FAST. The 90% threshold of FAST sample is $2.5\times 10^{37}$ erg (Li et al., 2021a).. The cumulative distribution of waiting time is shown in Figure 3. The green dashed line is the predicted result of the Weibull function, the dot-dashed red line is predicted by a Poisson process, and the black dotted line is the best-fitting results of the log-normal distribution. Some previous works suggest log-normal distribution is suitable to fit the waiting time distribution (Gourdji et al., 2019; Katz, 2019). Li et al. (2021a) used two log-normal functions to fit the waiting time distribution. We add the fit of log-normal distribution as a reference and do not discuss it in detail. The $\chi^{2}$ values for Poisson distribution, Weibull distribution, and log-normal distribution are 392.28, 346.50, and 394.90, respectively. The Weibull function is the best in fitting the data. In Figure 3, Poisson and log-normal distributions cannot well describe short waiting times. However, the Weibull distribution has a better description of that. For long waiting times, both the Weibull distribution and log-normal distribution can model the data well.

Oppermann et al. (2018) collected 17 bursts observed by Arecibo, Effelsberg, GBT, VLA, and Lovell. They derived $k=0.34^{+0.06}_{-0.05}$ with these 17 bursts. Oostrum et al. (2020) obtained $k=0.49\pm 0.05$. Their results are smaller than ours. This may be caused by the low waiting time cutoff. Cruces et al. (2021) found that when excluding the waiting time $\delta_{t}<1$ s, the shape parameter $k$ shifts from $k=0.62^{+0.10}_{-0.09}$ to $k=0.73^{+0.12}_{-0.10}$. We carefully check the effect of waiting time cutoff $\delta_{t,c}$ and show the dependence of $k$ on $\delta_{t,c}$ in Figure 4. The solid blue line is the evolution of $k$ and the vertical solid blue lines are 1$\sigma$ errors. The shape parameter $k$ increases as $\delta_{t,c}$ increases. When $\delta_{t,c}\simeq 28~{}s$, the shape parameter $k$ closes to 1, which means that the repeating behavior is due to a Poisson process. Cruces et al. (2021) also found that the waiting time with $\delta_{t}>100$ s is consistent with a Poisson process, which is similar to our results.

We use the FAST data to check the burst repetition pattern. We select the single-day observations with the number of detected bursts larger than 60 and show the cumulative distribution function of the burst time with blue scatters in each panel of Figure 5. The observation time is listed in the top of each panel. The x-axis is the decimal part of the MJD time. We find that the distribution can be fitted with a linear function

where $t$ is the burst time, $s$ is the slope, and $b$ is a constant. The best-fitting results with the dashed red lines are shown in each panel and the best-fitting parameters are also listed. The observation on 2019-09-07 is slightly different from other observations. There is an apparent temporal gap between early bursts occurring before MJD 58732.90 and late ones occurring after MJD 58732.90. In this observation, most bursts occurred after MJD 58732.90 and these bursts can be fitted with Equation (3). For all the observations listed in Figure 5, although the fitting results are different, they can be all fitted with equation 3. Tabor & Loeb (2020) proposed that the burst rate is constant per logarithmic time with the observation of GBT (Zhang et al., 2018). They found that the cumulative distribution function can be fitted with $N(<t)=\alpha\ln(t/t_{0})$. Contrary, our results suggest that FAST strongly favors a linear pattern. Tabor & Loeb (2020) excluded the model that repeating FRB sources origin from pulsars (Beniamini et al., 2020). According to our results, this model is still feasible. The linear pattern may suggest a potentially short period. However, similar to Li et al. (2021a), we do not find any short period. The Lomb-Scargle method and epoch folding method are used to search between 10 milliseconds and 30 minutes. No significant period was found. The fluctuations of data points along the fitted lines in Figure 5 may suggest that the linear pattern is caused by the narrow waiting time distribution. Besides, the slopes in varied observations are different, which may suggest diverse waiting time distribution in different observations.

An important conclusion of Li et al. (2021a) is the novel triple peaks structure in the KDE of energy and burst time (As shown in the panel (c) of their Figure 1). There is a clear energy gap between the high-energy bursts and low-energy bursts before MJD 58740. The temporal gap between early bursts and late bursts is also significant. Following the division of Li et al. (2021a), we regard the bursts occurring before MJD 58740 as early bursts and the bursts occurring after MJD 58740 as late bursts. The energy distributions of early bursts and late bursts are shown in Figure 6. Again, we only fit the bursts with energy $>10^{38}$ erg. The best-fitting results are shown as green dot-dashed lines. The power-law index for the early bursts is $-1.70^{+0.03}_{-0.03}$ (1$\sigma$ uncertainty), while that for the late bursts is $-2.60^{+0.15}_{-0.14}$ (1$\sigma$ uncertainty), which is significantly different.

Some previous works support the power-law index $-1.7$ (Wang & Yu, 2017; Law et al., 2017; Wang & Zhang, 2019), which is consistent with the result of early bursts. However, there are also some works suggest a steeper index (Gourdji et al., 2019), which is close to the result of late bursts. This difference may be caused by the following reasons. There may be different modes of burst emission operating in different epochs, so that the observed energy distribution may depend on the observing time. The two different modes may be related to different mechanisms or different emission sites. Li et al. (2021a) reported the triple peaks in the distribution of energy and burst time (as shown in panel (c) of their Figure 1). At the early epoch, the bursts with $E>10^{38}$ erg belong to the high-energy burst component. However, the high-energy peak is invisible in the late epoch. The bursts with $E>10^{38}$ erg in the late epoch are an extension of the low-energy peak. If the high-energy bursts and low-energy bursts originate from different mechanisms or different sites, which operate in different epochs, the variance of the energy distribution is naturally interpreted.

Li et al. (2021a) reported the triple peak structure in the panel (c) of their Figure 1. The bimodal structure in the energy distribution in the early phase (before MJD 58740) may suggest that they are produced by different physical processes. Therefore, we explore whether there are differences in their other properties, such as waiting time. If the low-energy bursts and high-energy bursts are produced by different processes, the waiting time should be calculated independently. We divide all bursts into three sub-samples: the low-energy bursts and high-energy bursts between MJD 58717 and 58740, and the late-phase bursts detected between MJD 58740 and 58776. The dividing line of burst time is consistent with the choice of Li et al. (2021a). We take the dividing line between the low-energy bursts and high-energy bursts as $1.58\times 10^{38}$ erg. The gap is about $10^{38.1}\sim 10^{38.3}$ erg. For convenience, we take the separation line as $10^{38.2}\simeq 1.58\times 10^{38}$ erg. We also test other choices of the criteria and find that it does not affect on our result. For each sub-sample we calculate the waiting time independently and show the KDE of the waiting time in Figure 7. The KDE for these three sub-samples is slightly different. The medians of the low-energy bursts, high-energy bursts, and late-bursts are 78.19 s, 133.61 s, and 54.62 s, respectively. There is a small peak in the KDE density of the late-phase bursts, which is close to 0.1 s.

We also use the Weibull function to fit three sub-samples. Again, we ignore the waiting times with $\delta_{t}<30$ ms and $\delta_{t}>0.5$ day. We list the best-fitting results for three sub-samples in Table 1. The shape parameters are $k\simeq 0.69\pm 0.02$ for the low-energy bursts, $k\simeq 0.76\pm 0.03$ for the high-energy bursts and $k\simeq 0.66\pm 0.02$ for the late-phase bursts. The posterior distributions of $k$ for the three sub-samples are shown in Figure 8. The $k$ parameters for the low-energy bursts, high-energy bursts, and late-phase bursts are shown as black line, red line, and blue line, respectively. The vertical dashed lines are the best-fitting results. The shape parameters $k$ for the low-energy and late-phase bursts are consistent with each other within 1$\sigma$ error, while the $k$ for the high-energy bursts is higher. The consistency of the low-energy bursts and late-phase bursts may suggest that they have the same origin, while the high-energy bursts may have a different origin. The difference between the late-phase bursts and high-energy bursts is significant, but the variance between the low-energy bursts and high-energy bursts is insignificant. This may be caused by the rough $1.58\times 10^{38}$ erg division line of low-energy bursts and high-energy bursts. The count distributions of low-energy and high-energy bursts are close to a log-normal distribution. They can be extended to include each other. Therefore, some FRBs may be misclassified. If we can make a robust classification, the difference in waiting time may be larger.

We derive the burst rates for the low-energy bursts and high-energy bursts are $445.35^{+31.85}_{-30.02}$per day and $354.74^{+27.14}_{-26.14}$ per day, respectively. The burst rate of the low-energy bursts is higher than that of the high energy bursts. These two classes of bursts occurred in the same epoch. For the bursts occurring in late observation days, the burst rate is about $800.54^{+48.84}_{-48.32}$ per day, which is higher than that of low-energy bursts and high-energy bursts, but close to the sum of these two classes of bursts.

We also check the dependence of energy on waiting time. For X-ray flares of gamma-ray bursts, an anti-correlation is found between energy and waiting time (Yi et al., 2016). The scatter plot of energy and waiting time of the FAST FRB sample is shown in Figure 9. No significant dependence of energy on waiting time is found. We also check the dependence for three sub-samples and found no significant dependence.

The FAST observations spanned two months in time, which covered most of the active window of FRB 121102. We can test the dependence of the burst behaviors on observing time. The observations lasted for roughly 1 hour each day. We regard the observations in each day as single samples and select some samples with the number of observed bursts greater than 5. Then, the Weibull function is applied to fit the waiting time distribution in these days. The best-fitting results against burst time are shown in Figure 10. The blue points are the shape parameters $k$ with $1\sigma$ errors, and the yellow points are the burst rates $r$ with $1\sigma$ errors. Interestingly, the value of $k$ is close to $1$ in some days, which indicates a rough Poisson distribution of the waiting time. Some observations have $k>1$, which means that the waiting time prefers a specific value. If $k\gg 1$, it means a periodic behavior. In our calculation, the value of $k$ is not enough to support the periodic behavior. It only suggests that the waiting time has a narrow distribution in some days, which is consistent with the linear pattern. The best-fitting result varies for different samples. We do not find significant dependence on the observation time.