RAPID GAMMA-RAY FLUX VARIABILITY DURING THE 2013 MARCH CRAB NEBULA FLARE

The Crab Nebula and its pulsar are among the best-studied objects in astronomy. Their origin is assumed to be a massive core-collapse supernova observed in the year 1054 A.D. During the explosion of the progenitor star, a fast-rotating neutron star, the Crab pulsar, was formed. It emits pulsed radiation from radio wavelengths (see, e.g., Hester 2008) up to several hundreds of GeV (VERITAS Collaboration et al. 2011; Aleksić et al. 2012). The pulsar, with a spin period of 33 ms, constantly dissipates an enormous amount of rotational energy into the surrounding medium at a rate of 4.6 · 10³⁸ erg s⁻¹ (Manchester et al. 2005). A fraction of this energy powers the acceleration of relativistic particles, which propagate away from the pulsar. These particles, thought to be mainly electrons and positrons (see, e.g., Gaensler & Slane 2006), lose energy by synchrotron radiation visible from radio up to hundreds of MeV. Very high energy (VHE) photons, measured up to 50 TeV (Aharonian et al. 2006; VERITAS Collaboration et al. 2011; Zanin 2011), result from inverse-Compton (IC) scattering of the generated synchrotron and ambient radiation fields.

The photon emission of the nebula powered by the relativistic particle wind is called a pulsar wind nebula (Crab Nebula in this case). Due to the underlying radiation processes a time-invariant luminosity over timescales of hundreds of years is expected. Nevertheless, instabilities in the flux of high energy (HE; E > 100 MeV) γ rays have been reported in recent years by AGILE and the Fermi Large Area Telescope (Fermi-LAT; Tavani et al. 2011; Abdo et al. 2011; Balbo et al. 2011; Striani et al. 2011; Buehler et al. 2012; Ackermann et al. 2013b). These flares have all shown increased emission from the synchrotron component of the Crab Nebula while emission from the IC component of the nebula as well as the Crab pulsar has remained consistent with the average level.

In 2013 March, Fermi-LAT detected another bright flare of the Crab Nebula (Ojha et al. 2013). The total flux above 100 MeV increased by almost a factor of six, making it the second brightest flare detected from the Crab Nebula, exceeded only by the flare in 2011 April (Buehler et al. 2012). In addition to dedicated LAT observations, multi-wavelength observations were triggered to gain a more complete picture of the latest event. Detailed analyses of these observations will be published separately. In this paper, we present the temporal and spectral analyses of the Fermi-LAT data.

The Fermi-LAT is a pair-conversion detector measuring photons with energies above 20 MeV. The point-spread function (PSF; 68% containment radius) decreases as a strong function of energy, from 70 at 70 MeV to 025 at 10 GeV (Ackermann et al. 2012). Therefore, with an apparent size of ≈003, the Crab Nebula can be treated as a point source in LAT analysis. The large field of view (≈2.4 sr) of the LAT and the survey mode used primarily for observations allow images of the full sky every three hours.

In the beginning of 2013 March, automated science processing (ASP) of the LAT data (Atwood et al. 2009; Chiang 2012) detected a significantly increased photon flux of the Crab pulsar and nebula (henceforth Crab). To maximize the exposure of the Crab on short timescales, the Fermi observatory was switched to a pointed target of opportunity observation mode between MJD 56355.65 and MJD 56359.77. For the subsequent analyses of the flare, we selected photons within a 15° radius region of interest (ROI) centered on the Crab and arriving between MJD 56346.0 and MJD 56369.5 (entire analysis window; EAW). In order to avoid contamination from the γ-ray emission of the Earth limb, we only consider photons with a zenith angle less than 95°. Accordingly, we did not use time windows where the edge of the ROI was located at zenith angles larger than 95°. The PSF for the standard instrument response functions (IRFs) P7SOURCE_V6 has been corrected based on flight data (Ackermann et al. 2012, 2013a). The corrections effectively neglect the inclination angle dependence of the width of the PSF, which is not a consideration for analyses integrated over many orbits. However, because here we are analyzing data for short time intervals, we use the Monte-Carlo-simulated P7SOURCE_V6MC IRFs, which include the dependence on inclination angle. We performed all analyses using unbinned gtlike (Fermi Science Tools v9r31p1).⁷

We first determined the best-fit model for the EAW. The ROI source model consists of 45 sources from the 2FGL catalog (all within 20° of the Crab) fixed to their published spectral parameters (Nolan et al. 2012) and the models for isotropic and Galactic diffuse emission (iso_p7v6source.txt, gal_2yearp7v6_v0.fits).⁸ In the model, we decomposed the Crab into three components, describing the pulsar as a power law with sub-exponential cut-off, the IC nebula as a smoothly broken power law and the synchrotron component as a power law. The spectral models for the Crab with their respective parameters were adapted from Buehler et al. (2012). Since the pulsar and IC component did not vary with time (see the next paragraph), the only free parameters during the fit were the normalizations of the diffuse components, the spectral index of the power law that scales the Galactic diffuse model, and the spectral parameters of the synchrotron nebula. The spatial residuals for the fit between 70 MeV and 300 GeV are compatible with expectations for statistical fluctuations. The flux of the synchrotron component within the EAW was found to be (4.05 ± 0.08) · 10⁻⁶ cm⁻² s⁻¹, with a photon spectral index Γ = 3.09 ± 0.03. Adding the constant components for the pulsar and IC nebula leads to a total Crab flux of (6.21 ± 0.08) · 10⁻⁶ cm⁻² s⁻¹. During this time interval, the detectable synchrotron spectrum extended to energies of about 700 MeV. Systematic errors on integral fluxes are ∼11% and ∼0.12 on the spectral indicies (Ackermann et al. 2012). Using gtfindsrc, the position of the flaring component during the EAW was found to be R.A. = 8368, decl. = 2198 (J2000) with a statistical error radius (68% containment) of 013, which is positionally consistent with the Crab Nebula. In the following, for the temporal analyses in finer time bins, we fixed the parameters of the Galactic diffuse model to the fitted values for the EAW. The isotropic diffuse component, however, has been left free to account for variations in the cosmic-ray background, which depends on the spacecraft orbital location and orientation.

We verified the hypothesis that the pulsar flux remains constant during the flare. We determined the ephemeris directly from the LAT data, spanning the full data set by extracting 1730 daily pulse times of arrival (TOAs; Ray et al. 2011) from MJD 54683 to 56425 with a typical uncertainty of ∼140 μs. Subsequently, we used tempo2 (Hobbs et al. 2006) to fit a timing solution to these TOAs, using 30 harmonically related sinusoids (Hobbs et al. 2004) to model timing noise. The resulting timing solution has an rms residual of 160 μs, or 4.8 · 10⁻³ of the rotational period. After assigning the pulse phases ϕ to the photons using this ephemeris,⁹ we fitted the off-pulse synchrotron component (0.5 < ϕ < 0.8) during the EAW, without the pulsar in the source model. The fit results agree well with the synchrotron flux determined above during the complete EAW. Using this fitted component as a fixed, underlying background, we analyzed the phase-averaged spectrum of the pulsar. Our fits are compatible with those of Buehler et al. (2012). Therefore, within our measurement accuracy, the spectrum of the pulsar during the EAW did not change with respect to the all-time average. This was confirmed using an ephemeris derived from radio data from the Jodrell Bank Observatory (Lyne et al. 1993).

2.1. Temporal Analysis

To obtain a time-dependent fit for the synchrotron component, we subdivided the data set into time bins of 6 hr length. The binning used here is the same as for that used in the ASP and provides a good compromise between sufficient photon statistics and sensitivity to sub-day flux variability. We analyzed the synchrotron nebula in each time interval, again leaving the pulsar and the IC nebula fixed to the values determined in Buehler et al. (2012). The varying flux of the Crab for this time binning is shown in Figure 1. If the test statistic (TS) for the maximum likelihood analysis (see, e.g., Nolan et al. 2012) was less than 4, we calculated an upper limit (95% confidence level). The light curve shows three sharp spikes (MJD 56357.1, 56357.9, 56360.1) on top of a strongly increased flux level. The average synchrotron component corresponds to the all-time average determined in Buehler et al. (2012). During the peaks, the flux was approximately twice as high as the underlying flare. Therefore, this time binning suggests an approximate flux doubling timescale of 6 hr.

Zoom In Zoom Out Reset image size

To search for more rapid flux variations, we applied a finer time binning by defining one analysis window for each orbit (≈90 minutes). In each orbit, the Crab is within the field of view of LAT for as much as ≈45 minutes, but that duration can be reduced significantly by details of the observing strategy from Fermi, the precession of the orbit of Fermi, and the passage of Fermi through the South Atlantic Anomaly. Many orbits thus have little or no useful exposure to the Crab. The photon statistics within these short bins do not allow a fit of the synchrotron component independently of the pulsar and IC nebula. Therefore, we fit the spectrum of the Crab as a single power law in addition to leaving the isotropic diffuse normalization free in the fit. The light curve obtained for this binning is shown in Figure 2.

Zoom In Zoom Out Reset image size

To quantify the variability timescale and to assess the significance of substructures in the light curve, we decomposed it into time windows that are statistically compatible with a constant flux. These windows, so-called Bayesian blocks, have been obtained using the method described in Scargle (1998) and are shown in Figure 2. For analysis in these time bins we used the same ROI and source model as for the 6 hr light curve, where the Crab is modeled as three components and the synchrotron nebula parameters are varied in the fit. The orbit-binned light curve is statistically compatible with the one obtained in the Bayesian blocks binning (χ²/dof = 180/169), showing that substructures on shorter timescales are not statistically significant. The shortest of these time bins is 5 hr long and provides a measure of the shortest detectable variability timescale. This value is compatible with the approximate doubling timescale of 6 hr found previously. The peak flux of the Crab is (12.5 ± 0.8) · 10⁻⁶ cm⁻² s⁻¹ at the Bayesian block centered on MJD 56357.11, almost six times larger than the average quiescent flux. If the constant flux values for the pulsar and IC nebula are subtracted, then the flux increase in the synchrotron nebula is a factor of ∼20. The spectral index of a power-law model during this Bayesian block was 2.4 ± 0.1.

Finally, we searched for variability of the synchrotron nebula on timescales shorter than the orbit binning by applying a Bayesian block analysis on the single photon arrival times. No significant short-term variability on these timescales was detected.

2.2. Spectral Analysis

Similar to the 2011 April flare, the spectrum of the synchrotron nebula hardens as the flux increases. This is visualized in the bottom panel of Figure 2, showing the spectral indices of the power-law fits. A constant spectral index during the EAW is rejected at 7.4σ. To study the spectral evolution in more detail, we determined the spectral energy distribution of the synchrotron component within each of the 13 Bayesian blocks; see Figure 3. An additional component in the spectral energy distribution arises during the flaring period. Describing the synchrotron component by power-law spectra during the flare leads to χ²/dof = 61/44 for all data points relative to their respective model.¹⁰ The rather poor fit probability (P = 0.046) suggests curvature in the spectra. We therefore parameterized the synchrotron nebula spectrum by a power law with an exponential cut-off, as seen during the 2011 April flare. This yields an increase ΔTS =77 in the unbinned likelihood analysis and results in χ²/dof = 37/32. Taking into account the 13 additional degrees of freedom, this corresponds to a fit improvement on a 3.7σ level. Accordingly, we see evidence for a cut-off in the flare spectrum. Specifically, in the Bayesian blocks 1, 4, and 7 the fit quality increases with a significance >3σ. The best-fit values for these models along with the parameters of the power-law models are shown in Table 1. For Bayesian block 1, the cut-off is fitted below the analysis threshold energy of 70 MeV. The results obtained in this block therefore need to be taken with care. They depend on the validity of the spectral model outside the LAT band.

Zoom In Zoom Out Reset image size

Table 1. Spectral Parameters of the Synchrotron Component for the Different Bayesian Blocks Shown in Figure 3

Block Number	Flux >100 MeV	Photon Index	Cut-off	ΔTS
	(10−6 cm−2 s−1)		(MeV)
0	0.8 ± 0.4	4.4 ± 1.0	...	...
1	2.3 ± 0.2	3.5 ± 0.1	...	...
2	5.3 ± 0.3	3.0 ± 0.1	...	...
3	6.6 ± 0.3	2.9 ± 0.1	...	...
4	10.2 ± 0.8	2.4 ± 0.1	...	...
5	5.7 ± 0.4	3.1 ± 0.1	...	...
6	9.8 ± 1.0	2.7 ± 0.1	...	...
7	6.1 ± 0.4	2.8 ± 0.1	...	...
8	4.5 ± 0.2	2.9 ± 0.1	...	...
9	7.4 ± 0.9	2.9 ± 0.2	...	...
10	4.0 ± 0.3	3.2 ± 0.1	...	...
11	2.8 ± 0.2	3.3 ± 0.1	...	...
12	1.8 ± 0.3	3.8 ± 0.3	...	...
1	2.6 ± 0.2	0.7 ± 1.3	$53^{+28}_{-28}$	23.1
4	10.7 ± 0.8	1.7 ± 0.3	$484_{-166}^{+409}$	11.4
7	6.5 ± 0.5	2.0 ± 0.3	$291_{-97}^{+192}$	9.2

Notes. The upper section lists the results of the power-law fits, while the lower section shows the spectral parameters during the Bayesian blocks with significant cut-off (ΔTS > 9). ΔTS denotes the improvement of the TS value of the cut-off model compared to the power-law fit.

Download table as:  ASCIITypeset image

In the past, Fermi-LAT has detected flares from the Crab with a variety of spectral behaviors: the flare in 2009 February exhibited only a flux increase with no spectral change. In contrast, the flares in 2010 September and 2011 April had fluxes strongly correlated with the spectral index. Moreover, in 2011 April a cut-off in the flare spectrum could be measured, which is also indicated in our flare data. In this flare, the cut-off during the brightest part of the flare (Bayesian block 4) is at $484^{+409}_{-166}$ MeV, similar to the value of 375 ± 26 MeV found during the 2011 April flare (Buehler et al. 2012).

Fermi-LAT detected enhanced γ-ray flux from the Crab Nebula in 2013 March. The flux of the synchrotron component above 100 MeV increased by up to a factor of 20, making this flare the second brightest event observed to date. The shortest variability timescale was determined to be ∼5 hr, while the flux of the synchrotron component could be measured up to ∼700 MeV. As observed in previous flares, the spectral index of the synchrotron nebula hardened significantly during periods of increased flux. An indication for an exponential cut-off in the spectrum is found at a 3.7σ level. Within the time window of the highest flux, we found a spectral index of 1.7 ± 0.3 and a cut-off energy of $484_{-166}^{+409}$ MeV. Interestingly, the maximum cut-off energy found in the cases of 2011 April and 2013 March is ∼400 MeV.

Several ideas have been proposed to explain these recurring flares (∼1 per year). They might originate from regions that are relativistically boosted toward us. Several works state that such regions could emerge within instabilities of the pulsar wind termination shock (Bednarek & Idec 2011; Komissarov & Lyutikov 2011; Lyutikov et al. 2012) or further outside at the polar region of the inner nebula (e.g., Lyubarsky 2012). Knots of relativistic particles could move out from the inner nebula (Yuan et al. 2011) causing local, strongly varying magnetic fields (Bykov et al. 2012).

The acceleration of the HE γ-ray emitting particles is proposed to happen via magnetic reconnection (e.g., Sironi & Spitkovsky 2011; Uzdensky et al. 2011; Cerutti et al. 2012). From particle-in-cell simulations, Cerutti et al. (2013) found that reconnection events can produce a flaring synchrotron component by linearly accelerating leptons within the strong electric field of the reconnection layer. This acceleration mechanism can also explain the rapid variability during the flare: beams of relativistic particles cross our line of sight and can thus cause rapid jumps in the synchrotron flux. The spectral hardening during these flux peaks can be attributed to anisotropic beaming of the relativistic particles, which causes a shift in the spectral energy distribution to higher energies. Recently, the flares have also been placed in the context of the long-standing sigma problem of pulsar wind nebulae (see, e.g., Kennel & Coroniti 1984). The flares might be part of the magnetic dissipation in the pulsar wind. As recently shown in three-dimensional magneto-hydro dynamical simulations, magnetic dissipation can efficiently occur down stream of the wind termination shock allowing for large magnetization (see, e.g., Komissarov 2013; Porth et al. 2013 and references therein).

The absence of plausible counterparts at other wavelengths in past flares is one of the most surprising aspects of the flare phenomenon (Weisskopf et al. 2013). Previous non-detections of correlated variability suggested that the spectrum of the emitting electrons was hard, or possibly even mono-energetic. During the time span of the flare presented in this article, several instruments obtained simultaneous observations (see Figure 1). In order to allow direct, consistent comparisons of our LAT analysis with overlapping multi-wavelength observations, we provide the corresponding spectral analysis of the LAT data for each of the observation time windows in Table 2 (all spectra and light curves, along with the pulsar ephemeris applied in the analysis, are available online⁹). The results of the multi-wavelength campaign will appear in separate publications and will hopefully shed more light onto the origin of the gamma-ray flares.

Elizabeth Hays and Rolf Bühler acknowledge generous support from the Fermi guest investigator program. This research was funded in part by NASA through Fermi Guest Investigator grant NNH10ZDA001N (proposal number 41213). This research was supported by an appointment to the NASA Postdoctoral Program at the Goddard Space Flight Center, administered by Oak Ridge Associated Universities through a contract with NASA. C.C.C. was supported at NRL by NASA DPR S-15633-Y. The Fermi LAT Collaboration acknowledges generous ongoing support from a number of agencies and institutes that have supported both the development and the operation of the LAT as well as scientific data analysis. These include the National Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l'Energie Atomique and the Centre National de la Recherche Scientifique/Institut National de Physique Nucléaire et de Physique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di Fisica Nucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), High Energy Accelerator Research Organization (KEK) and Japan Aerospace Exploration Agency (JAXA) in Japan, and the K. A. Wallenberg Foundation, the Swedish Research Council and the Swedish National Space Board in Sweden. Additional support for science analysis during the operations phase is gratefully acknowledged from the Istituto Nazionale di Astrofisica in Italy and the Centre National dEtudes Spatiales in France.