IMPULSIVE AND LONG DURATION HIGH-ENERGY GAMMA-RAY EMISSION FROM THE VERY BRIGHT 2012 MARCH 7 SOLAR FLARES

The first step in our analysis of the impulsive phase is to consider nine adjacent time intervals of LLE data (see the Appendix for a detailed description of the LLE technique). In particular, we included data for the region centered on the Sun at the time corresponding to the middle of each time interval and selected the intervals when the Sun was within 70° from the LAT boresight. For each time interval, we extract two sets of background LLE data at 30 orbits before and after the flare, when Fermi was at a similar geomagnetic location (see the Appendix for further information regarding this choice). At ±30 orbits (∼2 days) the location and attitude of the spacecraft are approximately the same as during the impulsive phase of the flare. The reconstructed direction in instrument coordinates is not the same after ±30 orbits; therefore we keep the reconstructed direction fixed and move the position of the region of interest (ROI), i.e., the region of the sky around the transient source that we are analyzing, in order to sample the same cosmic-ray-induced background rate as during the flaring episode. This last step is needed to average the two background data sets because the local cosmic-ray-induced background dominates the in-aperture celestial background. In this way we obtain nine source spectra and nine background spectra (one for each time interval during the impulsive phase analysis). We compute the LLE energy redistribution matrices for each of the nine intervals separately.

We fit the data between 100 MeV and 10 GeV using XSPEC⁶² to test three models. The first two are simple phenomenological functions that may describe bremsstrahlung emission from accelerated electrons, namely a pure power law (PL) and a power law with an exponential cutoff (PLEXP):

$$\begin{equation} \frac{dN(\epsilon)}{d\epsilon } = N_{0}\,\epsilon ^{-\Gamma }\,\exp \left({-\frac{\epsilon }{\epsilon _{co}}}\right); \end{equation} \tag{ 1 }$$

where Γ is the photon index and _co is the cutoff energy. We found that the data clearly diverge from a pure power-law spectrum and that the PLEXP provides a better fit in all time intervals considered. The third model uses templates based on a detailed study of the gamma rays produced from pion decay (updated from Murphy et al. 1987). In this model accelerated high-energy protons (and other ions) with an assumed energy distribution collide with particles of the solar atmosphere, creating π⁰ and π^± mesons.⁶³ A π⁰ quickly decays into two gamma rays, each having an energy of 67.5 MeV in the rest frame of the meson. The decay products of charged π^± mesons (secondary e^±) produce gamma rays via Bremsstrahlung or by annihilation-in-flight of the positrons, and microwaves via synchrotron radiation.⁶⁴ The interactions between the accelerated and background protons (and ions) also produce nuclear de-excitation lines in the 1–10 MeV range, observable by the GBM. The analysis of the GBM observations of these flaring episodes will be presented in a subsequent paper.

The pion-decay templates used in our fits depend on the ambient density, composition and magnetic field, and on the accelerated-particle composition, pitch angle distribution and energy spectrum. The templates represent a particle population with an isotropic pitch angle distribution and a power-law energy spectrum (dN/dE∝E^−s, with E the kinetic energy of the protons) interacting in a thick target with a coronal composition (Reames 1995) taking ⁴He/H = 0.1. To obtain the gamma-ray flux value we fit the data by varying the proton spectral index s from 2 to 6, in steps of 0.1. In this way, we fit the LAT data with a model with two free parameters, the normalization and the index, s. The time dependence of the >100 MeV gamma-ray flux and of the proton index, s, derived using LAT data, is displayed in the lower panel of Figure 1, and the numerical values are reported in Table 1, as well as the best fit parameters of the PLEXP model.

Table 1. Spectral Analysis of the Impulsive Phase

Time Interval	Proton index	Fluxπa	Γ	co	FluxPLEXPa
2012/03/07 UT				MeV
00:38:52–00:43:52	3.07 ± 0.07 ± 0.09	21 ± 1 ± 4	0.07 ± 0.09 ± 0.2	130 ± 8 ± 10	18.0 ± 0.4 ± 4
00:43:52–00:48:52	3.36 ± 0.07 ± 0.1	18.7 ± 0.6 ± 4	0.26 ± 0.07 ± 0.1	107 ± 4 ± 10	16.3 ± 0.3 ± 4
00:48:52–00:53:52	3.48 ± 0.07 ± 0.1	15.5 ± 0.3 ± 3	0.23 ± 0.06 ± 0.1	106 ± 4 ± 10	14.1 ± 0.2 ± 3
00:53:52–00:58:52	3.40 ± 0.06 ± 0.1	14.4 ± 0.4 ± 3	0.19 ± 0.06 ± 0.1	109 ± 4 ± 10	12.7 ± 0.2 ± 3
00:58:52–01:03:52	3.23 ± 0.07 ± 0.1	12.7 ± 0.4 ± 3	0.18 ± 0.07 ± 0.1	114 ± 5 ± 10	10.9 ± 0.2 ± 2
01:03:52–01:08:52	3.25 ± 0.08 ± 0.1	9.6 ± 0.3 ± 2	0.25 ± 0.08 ± 0.1	111 ± 6 ± 10	8.6 ± 0.2 ± 2
01:08:52–01:13:52	2.95 ± 0.08 ± 0.09	9.0 ± 0.3 ± 2	0.00 ± 0.05 ± 0.1	136 ± 7 ± 10	7.2 ± 0.2 ± 2
01:13:52–01:18:52	3.0 ± 0.1 ± 0.09	7.2 ± 0.4 ± 1	−0.6 ± 0.10 ± 0.1	220 ± 30 ± 20	6.0 ± 0.3 ± 1
01:18:52–01:23:52	3.1 ± 0.2 ± 0.1	5.7 ± 0.6 ± 1	−0.8 ± 0.20 ± 0.1	270 ± 70 ± 30	5.0 ± 0.5 ± 1

Notes. In all cases the statistical errors are shown first and the systematic errors follow. ^aIntegrated energy flux between 100 MeV and 10 GeV, in units of 10⁻⁷ erg cm⁻² s⁻¹ and are calculated from the pion decay model (π) and for the exponential cut-off model (PLEXP).

Download table as:  ASCIITypeset image

It appears that, after a short phase of spectral softening (00:40:00–00:53:00 UT), the proton spectral index hardens before and during the rising phase of the impulsive interval of the second flare (00:53:00–01:08:00 UT) as seen by the GBM detectors (middle panel of Figure 1). The spectral index s correlates better with the GBM flux than with the high-energy flux measured by the LAT. For interpretation of these results; see Section 4.

We analyzed the data from after the first ninety minutes of Fermi-LAT observations using the standard likelihood analysis implemented in the Fermi-LAT ScienceTools⁶⁵ with P7SOURCE_V6 IRFs, selecting gamma rays within an ROI of 12° radius centered on the sun and that arrived within 100° of the zenith to reduce contamination from the Earth's limb. We include the azimuthal dependence of the effective area in calculating exposures.

To study the temporally extended emission, we perform time-resolved spectral analysis in Sun-centered coordinates by transforming the reference system from celestial coordinates to ecliptic Sun-centered coordinates. This is necessary in order to compensate for the effect of the apparent motion of the Sun during the long duration of the flare. We select intervals when the Sun was in the FOV (angular distance from the LAT boresight <70°) and use the unbinned maximum likelihood algorithm gtlike. We include the isotropic template model that is used to describe the extragalactic gamma-ray emission and the residual cosmic ray (CR) contamination,⁶⁶ leaving its normalization as the free parameter. Over short timescales, the diffuse Galactic emissions produced by CR interacting with the interstellar medium are not spatially resolved and are hence included in the isotropic template. We also add the gamma-ray emission from the quiescent Sun modeled as a point source located at the center of the disk, with a spectrum described by a simple power law with a spectral index of 2.11 and an integrated energy flux (>100 MeV) of 4.7 × 10⁻¹⁰ erg cm⁻² s⁻¹ (corresponding to a flux of 4.6 × 10⁻⁷ photons cm⁻² s⁻¹ as reported in Abdo et al. 2011). We did not include the extended Inverse Compton (IC) component described in Abdo et al. (2011) because it is too faint to be detected during these time intervals. We fit the data with the same two phenomenological functions used for the impulsive phase of the flare and use the likelihood ratio test to estimate whether the addition of the exponential cutoff is statistically significant. The test statistic (TS; Mattox et al. 1996) is defined as twice the increment of the logarithm of the likelihood $\mathcal {L}$ obtained by fitting the data with the source and background model components simultaneously. Because the null hypothesis (i.e., the model without an additional source) is the same for the two models, the increment of the TS (ΔTS = TS_PLEXP-TS_PL) is equivalent to the corresponding difference of maximum likelihoods computed between the two models.

For each interval, if ΔTS ⩾ 30 (roughly corresponding to 5σ) then the PLEXP model provides a significantly better fit than the simple power-law and we retain the additional spectral component. In these time intervals, we also used the pion decay model to fit the data and estimated the corresponding proton spectral index. We performed a series of fits with the pion decay template models calculated for a range of proton spectral indices. We then fit the resulting profile of the log-likelihood function with a parabola and determine its minimum ($\mathcal {L}_{\rm min}$) and the corresponding value s₀ as the maximum likelihood value of the proton index. The 68% confidence level is evaluated from the intersection of the profile with the horizontal line at −2 $\Delta \log (\mathcal {L}_{\rm min})+1$ (see insets in Figure 4). Table 2 summarizes our results.

Table 2. Spectral Analysis of the Time Extended Emission

Interval	Start (UT)	Duration	TSPL	ΔTSa	Γ	co	Fluxb	Proton Index	X,Y (r68⊕70)c
	2012/03/07	s				MeV			arcsec
(a)	02:27:00	1110	1336	85	1.0 ± 0.2 ± 0.05	250 ± 40 ± 4	2.4 ± 0.1 ± 0.2	3.8 ± 0.1 ± 0.1	−280, 140 (320)
(b)	03:52:00	2370	16318	987	0.9 ± 0.1 ± 0.04	210 ± 10 ± 2	3.55 ± 0.05 ± 0.3	4.0 ± 0.1 ± 0.1	-450, 260 (100)
(c)	05:38:32	1050	1393	159	0.2 ± 0.3 ± 0.04	120 ± 10 ± 2	3.1 ± 0.1 ± 0.2	4.6 ± 0.2 ± 0.1	−470, 260 (350)
(d)	07:03:00	2400	9003	756	0.4 ± 0.1 ± 0.03	130 ± 7 ± 1	2.07 ± 0.04 ± 0.2	4.8 ± 0.1 ± 0.1	−500, 130 (150)
(e)	08:50:00	1020	500	73	−0.2 ± 0.6 ± 0.1	90 ± 20 ± 1	1.4 ± 0.1 ± 0.1	5.1 ± 0.3 ± 0.1	670, 580 (750)
(f)	10:14:32	2370	1833	204	−0.3 ± 0.2 ± 0.03	80 ± 9 ± 1	0.81 ± 0.03 ± 0.07	5.5 ± 0.2 ± 0.1	440, 380 (330)
(g)	13:25:00	2400	137	13	2.9 ± 0.2 ± 0.02d	⋅⋅⋅	0.23 ± 0.04 ± 0.02	⋅⋅⋅	⋅⋅⋅
(h)	16:36:00	780	17	8	2.9 ± 0.4 ± 0.02d	⋅⋅⋅	0.03 ± 0.01 ± 0.003	⋅⋅⋅	⋅⋅⋅
(i)	18:24:00	540	10	3	3.7 ± 0.8 ± 0.03d	⋅⋅⋅	0.04 ± 0.03 ± 0.004	⋅⋅⋅	⋅⋅⋅
(j)	19:47:00	1710	59	2	2.6 ± 0.3 ± 0.02d	⋅⋅⋅	0.02 ± 0.01 ± 0.001	⋅⋅⋅	⋅⋅⋅
(k)	22:58:30	2370	7	4	1.6 ± 0.6 ± 0.01d	⋅⋅⋅	<0.006	⋅⋅⋅	⋅⋅⋅

Notes. In all cases the statistical errors are shown first and the systematic errors follow. ^aΔTS = TS_PLEXP-TS_PL. ^bEnergy Flux of gamma rays between 100 MeV and 10 GeV, in units of 10⁻⁷ erg cm⁻² s⁻¹, calculated using the best fit model. ^cA systematic error of 70 arcsec has been added in quadrature to the estimated 68% error radius, and is reported between parentheses. ^dThe value of the spectral index are obtained from the simple power-law fit.

Download table as:  ASCIITypeset image

In Figure 3 we compare the observed count spectra with the predicted numbers of counts for the different models. The predicted numbers are the sum of the contribution of the background and of the flare, after the spectral parameters are optimized. The contribution of the isotropic background and of the quiescent Sun is also shown in the figures. In the first six time intervals (a through f) a power-law model does not correctly reproduce the data, while a curved spectrum (such as the power law with an exponential cut off or a pion decay model) provides a better description of the data. In the time intervals from g) to j) the power-law representation is adequate to describe the data; in the last bin, the flare is only marginally significant (TS = 7); the flux and the photon index are compatible with the values of the quiescent Sun. For this reason we have indicated the last point as an upper limit (computing the 95% C.L.). In the lower panel of Figure 2 we combine the LLE and likelihood analysis results, showing the evolution of both the gamma-ray flux and the derived spectral index of the protons.⁶⁷ Unlike during the impulsive phase, the spectrum during the temporally extended phase becomes softer monotonically (s increases).

Zoom In Zoom Out Reset image size

We also compare our results with the GOES proton spectral data. For this, we selected two energy bands (>30 MeV and >100 MeV) and corrected the light curve by the proton time-of-flight (TOF) to 1 AU by considering the TOF for 30 MeV and 100 MeV protons (i.e., the maximum delay in each energy band). As a measure of the spectral index of the SEP protons (s_SEP), we compute the hardness ratio (HR_p) defined as:

$$\begin{equation} {\rm HR}_{p}=\ln \frac{P_{>100\,\rm MeV}}{P_{>30\,\rm MeV}}, \end{equation} \tag{ 2 }$$

where P is the integral of the proton flux (assuming that the proton flux is proportional to a power law). The HR_p is related to the value of the spectral index, s_SEP, of the SEP protons observed at 1 AU, roughly as:

$$\begin{equation} s_{\rm SEP} \sim 1-0.83\,{\rm HR}_{p}. \end{equation} \tag{ 3 }$$

To estimate the uncertainty associated with this procedure we repeat the calculation neglecting the TOF correction. In this way we obtain two values for the SEP spectral index for each time bin, corresponding to the actual and zero delay due to the time of flight. In Figure 2 we report the estimated proton spectral index as the average of these two values and its uncertainty as half the difference of these two values. However, we note that the s_SEP is for protons with energy less than a few hundred MeV while s is for protons with energies greater than 300 MeV. Diffusion is expected to play an important role in the transport of these SEPs but an in-depth transport analysis is beyond the scope of this paper. From our comparison we find that the proton spectral index inferred from the gamma-ray data is systematically softer than the value of the index derived directly from SEP observations but that the temporal evolution (hard-to-soft) is similar.

Uncertainties in the calibration of the LAT introduce systematic errors on the measurements. The uncertainty of the effective area is dominant, and for the P7SOURCE_V6 event class it is estimated to be ∼10% at 100 MeV, decreasing to ∼5% at 560 MeV, and increasing to ∼10% at 10 GeV and above. We studied the effect of the systematic uncertainties on our final results via the bracketing technique described in detail in Ackermann et al. (2012a). We find that the uncertainties on the flux are <10% and on the inferred proton index are <0.10. The results are represented by the gray bands in Figure 2 and in Table 2.

In the six time intervals where the ΔTS>30, we compute the photon spectral energy distribution. To do so we divide the data into 10 energy bins and determine the source flux using the unbinned maximum likelihood algorithm gtlike keeping the normalization of the background constant at the best fit value and assuming that the spectrum of the point source is an E⁻² power law. Note that fixing the power law index provides one less degree of freedom and, with a the number of energy bins large enough (>4 per decade in energy), the results do not depend on the particular value chosen. For nondetections (TS < 9), we compute 95% CL upper limits. The results are shown in Figure 4. We also report the values of the energy flux in each energy bin for the six time intervals in Table 3.

Zoom In Zoom Out Reset image size

We measure the direction of the >100 MeV gamma-ray emission using the gtfindsrc tool and perform a likelihood analysis on both time-integrated and separate time intervals. The background is modeled using the best-fit parameters obtained by the time-resolved spectral analysis described in the previous section, and the source is modeled according to the best-fit model. The uncertainties on the localization are obtained by combining the 68% error radius from gtfindsrc with the systematic bias associated with the "fisheye" effect in quadrature. We estimate the latter using Monte Carlo simulations and find it to be 002 (≈70'').

The results for the six time bins where the ΔTS > 30 are shown in Figure 5. The coordinates of the centroids (with uncertainties) are shown in the last column of Table 2. For the remaining time intervals, the reconstructed location of the emission is consistent with the direction of the Sun, although the associated uncertainty is larger than the angular diameter of the Sun.

Zoom In Zoom Out Reset image size

During the ∼20 hr of detected flaring gamma-ray emission, the LAT measured five photons with E > 2.5 GeV and reconstructed direction less than 1° from the center of the solar disk. Given the importance of these high-energy events and the variation of the PSF across the wide field of view (2.4 sr) of the LAT, we further investigate the characteristics of these photons. All five of these events belong to the P7SOURCE_V6 event class and three of them are also P7ULTRACLEAN_V6 (see Abdo et al. 2009b for a detailed description of the photon classes). Two of these photons, with energies 2.8 GeV and 4 GeV, were detected during the impulsive phase of the flare and the remaining three during the extended emission, including one with 4.5 GeV energy at 07:30:00 UT. Comparing the distance from the center of the solar disk and the predicted 68% containment radius from the point-spread function (PSF) of the instrument we find that four of the events are consistent with the solar disk. In the case of the 4.5 GeV photon, the reconstructed direction is 08 from the center of the solar disk and the 68% containment radius is approximately 02. Therefore, we conclude that the reconstructed direction of this event is only marginally consistent with the solar disk.

Considering the average rate of LAT-detected photons above ∼2.5 GeV coming from an ROI with radius 1° centered at the Sun (calculated using all available flight data, excluding the bright LAT-detected solar flare time intervals) we find that the probability to observe five or more events in 8 hr due to Poisson fluctuations is approximately P = 8.0 × 10⁻⁶ (∼4.8σ). In Table 4 we list some of the basic properties of these photons.

Table 4. High Energy Gamma Rays from the Solar Flare

Arrival Time	Energya	Distanceb	θc	Event Class	PSF$_{68\%}{}^{\rm d}$
2012/03/07 UT	(GeV)	(°)	(°)		(°)
0:49:37.8	2.8	0.2	49	SOURCE	0.3
1:18:44.6	4.0	0.6	66	ULTRACLEAN	0.5
2:35:01.7	2.9	0.6	62	SOURCE	0.6
4:12:28.4	2.9	0.5	36	ULTRACLEAN	0.6
7:30:24.6	4.5	0.8	44	ULTRACLEAN	0.2

Notes. ^aThe energy resolution at 3 GeV and 50° off-axis angle is of the order of 7%. ^bReconstructed distance from the center of the solar disk. ^cReconstructed direction with respect to the instrument coordinate system. ^dPSF$_{68\%}$ corresponds to the 68% containment radius calculated from the PSF of the instrument for the energy and direction of the gamma ray.

Download table as:  ASCIITypeset image

In the prompt acceleration model particles are injected quickly, e.g., as a power-law spectrum, into a trap region where they gradually lose energy and emit radiation (Murphy et al. 1987; Kanbach et al. 1993). If the radiation is produced in the trap region then we expect a spectral variation that depends on the energy dependence of the energy loss rate (Aschwanden 2004). For relativistic electrons the energy loss rate due to Coulomb collisions and synchrotron emission is given by

$$\begin{equation} {\dot{E}}_L (E)=\frac{E_p}{\tau _0}\left[1+\left(\frac{E}{E_p}\right)^\delta \right] \end{equation} \tag{ 4 }$$

where δ = 2 and for flare conditions τ₀ ∼ 10⁴ s (10¹⁰ cm⁻³n⁻¹) and E_p ∼ 10 MeV (n/10¹⁰ cm⁻³)^1/2(100 G/B_eff). This produces characteristic spectra that become flat at low energies, E < E_p due to Coulomb energy losses, and develop a sharp cut off at high energies, E > E_p due to synchrotron and IC losses, in less than an hour (see, e.g., Petrosian 2001).⁶⁸ This clearly disagrees with the data that indicate a power-law spectrum at low energies with a gradual high-energy cutoff. A similar expression with δ = 1, E_p ∼ 10 GeV and τ₀ ∼ 3 s  (10¹⁰ cm⁻³n⁻¹) day applies for 0.1–100 GeV protons with high-energy losses due to inelastic p–p interactions. This will cause some initial hardening of the spectra (index decreasing by ∼1) at low energies within several hours, and the higher-energy emission should persist for several days, which also disagrees with the observed spectral evolution. The marginally significant hardening within a few minutes before the X1.5 flare requires a density of n ≫ 10¹⁰ cm⁻³, which is not appropriate for a coronal trap model (Ryan 1986).

An alternative scenario is the trap precipitation model (Bai 1982) where the particles trapped in a magnetic bottle (converging field lines) are scattered into the loss cone that causes their precipitation into the chromosphere and below, where they lose most of their energy and produce gamma rays. Coulomb collisions cannot be the agent for this scattering, because the relativistic electron Coulomb scattering rate is lower than the Coulomb energy loss rate by a factor of γ² (where γ = E/m_ec²) and is much smaller than synchrotron energy loss or IC scattering rate (see, e.g., Table 1, Petrosian 1985). The Coulomb scattering rate for protons is lower than their energy loss rate by a factor proportional to the electron to proton mass ratio. In other words, with Coulomb scattering the particles lose energy before they are scattered into the loss cone. Therefore, a much faster scattering rate is required for this scenario. Scattering by turbulence could be a possibility but (1) in that case acceleration by turbulence will also be present so we no longer have a prompt model, and (2) gradual softening of the spectra will require a higher rate of scattering into the loss cone for higher energy particles, which may require a turbulence spectrum much steeper than Kolmogorov (see, e.g., Pryadko & Petrosian 1997). Thus, we conclude that trapping is negligible and relativistic electrons and protons of all energies are transported to the chromosphere and below within the transport time L/v ∼ 0.1 s (for L ∼ 10⁹ cm) where they lose energy rapidly and emit radiation as in the standard thick target model discussed below. Therefore, a more likely scenario is continuous acceleration (e.g., by turbulence; see Petrosian & Liu 2004) with a timescale comparable to, or shorter than, the particle energy-loss timescale. Continuous acceleration, as first proposed by Ryan & Lee (1991), was also favored to explain the COMPTEL observations (Rank et al. 2001) of extended emission in the solar flares of 1991 June. Observations of 2.223 MeV neutron capture radiation indicate continuous acceleration of ions, lasting for as long as several hours in the 1991 June 11 solar flare.

For electrons, nonthermal bremsstrahlung is the only viable mechanism of gamma-ray production (Trottet & Vilmer 1984; Vilmer 1987) and requires >100 MeV electrons. However, the very hard (flat) observed photon spectra with low-energy power-law indices (Γ < 1) require a very hard injected electron spectrum. In the thick target model and for relativistic electrons this means δ_inj ∼ Γ + 1 < 2. This is much harder than spectra deduced from hard X-ray impulsive phase emission in most flares. However, as shown by Park et al. (1997), at high energies bremsstrahlung spectra can deviate significantly from a power law and have a rollover as observed here. However, a spectral rollover at >100 MeV requires a weaker magnetic field than those expected in the corona.

In the corona the energy loss rate of >100 MeV, electrons are dominated by synchrotron losses with a timescale of τ_sync ∼ 30 (GeV/E)(100G/B)² s. The IC loss time is longer τ_IC ∼ 6 × 10³(GeV/E) s and the bremsstrahlung time, though shorter than Coulomb losses, is even longer: τ_brem ∼ 3 × 10⁴(10¹⁰ cm⁻³/n)(GeV/E)^0.1 s. However, all these times are longer than the transport time so that electrons lose a small fraction of their energy in the corona above the transition region. If the magnetic field reaches 10³ G synchrotron losses can become important with roughly a few percent of electron energy going into production of 1.5–150 THz radiation. Below the transition region density increases rapidly (within few hundred kilometers) and electrons are stopped at a column depth of ∼10²⁵ cm⁻² or density of ∼10¹⁷ cm⁻³, where they lose most of their energy producing bremsstrahlung radiation in a very short (0.01 s) time (see Petrosian et al. 1994), but with a total loss time of ∼L/v ∼ 0.1 s. This short time and the hard spectrum requires a very efficient acceleration which may be problematical considering the fact that the highest-energy photon observed by the Fermi-LAT of 4 GeV would require electrons to be accelerated to about 10 GeV. The above conclusion involved some assumptions and depends on the values of some critical parameters like the height of the accelerator and magnetic field. For a higher source or larger magnetic field, the acceleration time can be shorter and fraction of energy lost via synchrotron and bremsstrahlung processes will be different. Clearly, this requires a more detailed treatment of acceleration, transport, and emission processes, which is beyond the scope of this paper.

The situation is very similar and somewhat simpler for protons because protons with E < 1 GeV lose energy primarily via Coulomb collisions with background electrons with a loss time τ_L ∼ 3 × 10⁴(10¹⁰ cm⁻³/n)(E/GeV) s and through pion production in collisions with background protons at higher energies with a loss time of τ_L ∼ 3 × 10⁵(GeV/E)^0.1 s at E > 10 GeV. Again these times are much longer than the transport time to the transition region (up to the acceleration height of L < 10¹⁴ cm) so that most of the energy loss and radiation takes place below the transition region. For example, 300 MeV protons penetrate to column depths and densities very similar to those given above for electrons and lose energy in much less than a second. Thus, the situation is very similar to that of electrons except that the energy dependence of the proton loss rate varies in the range of proton energies of interest and the effective thick target spectrum N_eff(E) of protons per energy unit E for a continuous injection spectrum ${\dot{Q}}(E)$ is given by

$$\begin{equation} N_{\rm eff}(E) = \left[\frac{\tau _L(E)}{E}\right] \, \int _{E}^{\infty } dE^{\prime } \,\,{\dot{Q}}(E^{\prime }), \end{equation} \tag{ 5 }$$

so that for an injected power-law ${\dot{Q}}(E)\propto E^{-\delta }$ the effective spectrum will be a broken power-law steepening with an index change of about one around several GeV. Whether such a spectrum can describe the observations adequately is beyond the scope of the current paper. It should also be noted that the yield of gamma rays is about 1% at the pion production threshold of ∼300 MeV but becomes essentially 50% above a few GeV (the other half of the proton energy going to neutrinos).

In summary, both relativistic bremsstrahlung emission by electrons or decay of pions produced by >300 MeV protons describe the observations in a thick target model with acceleration time shorter than the transport time from the acceleration site to the thick target below the transition region. However, an unusually hard spectrum of electrons is needed while the required proton spectra are similar to those derived from other observations. Thus, pion decay seems to be the preferred radiation mechanism, in accord with conclusions reached by Rank et al. (2001) for the 1991 June flares mentioned above, based on the COMPTEL detection of the 2.2 MeV neutron capture line.

From the results of the gamma-ray spectral analysis in the proton scenario, and using the gamma-ray yield, we estimate the number and energy of the accelerated protons with kinetic energy >30 MeV producing gamma rays and compare with the observed SEPs. During the first impulsive phase the estimated number (energy) of protons interacting with the Sun is N_p ∼ 2.5 ×10³³ ($\mathcal {E}_{p}\sim$ 2.2 × 10²⁹ erg), while for the temporally extended emission, these are approximately N_p ∼ 1.0 × 10³⁴ ($\mathcal {E}_{p}\sim$ 7.2 × 10²⁹ erg).

From the GOES observations, we can roughly estimate the number (energy) of SEP protons escaping the CME shock. We integrate the flux observed by GOES detectors, assuming isotropic emission, and that the CME leading edge propagated toward the Earth with an almost-constant speed of ≈2600 km s⁻¹ until 20 R_☉. We do not have accurate measurements of the speed of the CME after this distance, so we use the information that it reached the Earth at approximately 10:42 UT of March 8 assuming constant acceleration (−0.02 km s⁻²). We calculate the position of the CME leading edge and we calculate the number of protons during the period of time when the gamma-ray flux was high (until approximately 23:00:00 UT on March 7); we obtain: N_SEP ∼ 3 × 10³³ ($\mathcal {E}_{\rm SEP}\sim$ 1.0 × 10³⁰ erg).⁶⁹ A more detailed calculation of the flux of accelerated particles at the CME shock, properly considering propagation affects and non-isotropic emission, is beyond the scope of this paper, but considering that this is only a fraction (<20%) of the SEP measured (i.e., the proton flux measured by GOES remains high until approximately 20:00:00 UT of March 12), we conclude that protons producing gamma rays carry significantly less energy than SEP protons observed by GOES.

Continuous acceleration of protons at the flare reconnection region, whether by a stochastic acceleration mechanism (Petrosian & Liu 2004) or by commonly postulated standing shock (see, e.g., Aschwanden 2004; Mann et al. 2006) can account for most of the spectral observations described above. In this model protons escape the acceleration site along closed field lines into the chromosphere and the spectral changes are simply due to the softening of the spectrum of protons as the flare decays. As shown by Petrosian & Liu (2004), in stochastic acceleration by turbulence (or by a standing coronal shock; see Petrosian 2012) the accelerated particle spectra become softer (harder) as the turbulence becomes weaker (stronger), which can naturally explain the observed spectral evolution both in the prompt and extended emission phases. Acceleration at the CME shock (Rank et al. 2001) is also an attractive possibility because the Large Angle and Spectrometric Chronograph (LASCO; Brueckner et al. 1995) on board the Solar and Heliospheric Observatory (SOHO) mission observed a fast CME ejected at approximately 00:30:00 UT, and measured the speed of the head of the fastest segment of the leading edge (Gopalswamy et al. 2009) to be 2684 km s⁻¹. The shock front of a CME is known to accelerate SEPs, sometimes to energies >300 MeV required for gamma-ray production Reames (2013), but its contribution to the acceleration of the particles that produce gamma rays remains unclear. The SEP protons observed by GOES are those accelerated by the shock and escape from the upstream region of the shock (see, e.g., Ramaty et al. 1990). There are also indications that a reservoir of accelerated particles exists in the downstream region of the shock behind the CME Reames (2012, 2013); Zank et al. (2000). However, as mentioned above, gamma-ray production would be too weak to be detected by Fermi-LAT, because of prevailing low densities (≪10¹⁰ cm⁻³). The accelerated particles must be transported to higher densities below the corona for efficient production of radiation. This requires existence of field lines connecting the reservoir to the Sun. We estimate that in 10 hr the CME would have traveled approximately 80 R_☉(∼0.4 A.U.). Whether such a connection can be maintained over this distance is matter of speculation. It is possible that the outflow results in a current sheet behind the CME that can provide a magnetically connecting path to the original flare site. The displacement of the reconstructed LAT emission centroid at later times can be explained by a larger dispersion of particles due to the longer distance traveled. This model provides the correct scenario for short acceleration time scales (∼1 hr).

The second difficulty is that the SEP spectra are much harder than the proton spectra deduced for the extended emission phase but their spectral index is similar to the impulsive phase index of the protons. It is not clear why the upstream and downstream protons would have such different spectra. One possible scenario is that during the extended phase protons are accelerated in the corona while the SEP protons are a result of the re-accelerated of prompt phase protons at the CME shock.

In summary, in this paper we have presented an analysis of the brightest solar flare detected by the Fermi LAT to date. We have shown that during most of the long-duration emission the gamma rays appear to come from the same active regions responsible for the flare emission. The fluxes and spectra of the high-energy gamma rays evolve differently during the impulsive phase and the sustained emission. Also there are correlations and some differences between the fluxes and spectral indexes of the protons required for the production of high-energy gamma-rays and SEP protons seen at 1 AU. From these results we suggest that the high-energy gamma rays are most likely produced by protons (rather than electrons) accelerated in the corona (rather than in the associated fast CME shock) continuously during the entire flare.

The LLE response files are generated based on dedicated Monte Carlo simulations. The simulations are used to study how the detector is "illuminated" by a source of a known flux and known position, during the real pointing history of the LAT. We do this by simulating a bright point source with a spectrum dN/dE ∼ E⁻¹ at the position of the source in question (the Sun in this case), and using the pointing information saved in the spacecraft data file (FT2 file). We use the Fermi-LAT full simulator (Baldini et al. 2006) to generate particles from the point source, incident over a cross-sectional area of 6 m², which illuminates the entire LAT. The LAT detector is represented by a complex model containing more than 34000 volumes. Gamma-ray conversion and particle propagation through the detector is implemented using GEANT4 (Geant4 Collaboration et al. 2003) while digitization and reconstruction are done using the same algorithm used for flight data. We then apply the LLE selection and bin the resulting events in reconstructed energy versus Monte Carlo energy (McEnergy) obtaining the so-called Redistribution Matrix, R_ij. This matrix is proportional to the probability that an incoming photon of energy E ∈ [E_j, E_{j + 1}] will be detected in the reconstructed energy bin [E_i, E_{i + 1}]. We re-normalize each bin such that:

$$\begin{equation} \sum _{i}R_{ij}=A_{j}=6\times 10^{4}\, {\rm cm}^2 \frac{\displaystyle N_{j}({E}_{j+1}-E_{j})}{\displaystyle N_{\rm TOT}({E}_{{\rm MAX}}-{E}_{{\rm MIN}})}, \end{equation} \tag{ A1 }$$

where N_TOT is the total number of simulated events over an area of 6 m² (typically 10⁷), N_j is the number of detected events (that survive the selection cuts) with a McEnergy between E_j and E_{j + 1}. A_j is usually defined as the effective collecting area of the instrument. The Redistribution Matrix File (RMF) is saved in the standard HEASARC RMF File Format.⁷⁰

In the case of short and bright transient gamma-ray sources, it is possible to select time windows before and after the transients (the "off-pulse" region) and, excluding the time window of the transient itself, fit the count rates in the off-pulse region with a polynomial function and in this way estimate, the background in the time interval of the transient. For LLE data, this analysis is described in Pelassa et al. (2010) and was applied to the 2010 June 12 flare, as presented in Ackermann et al. (2012b). This approach relies on a few assumptions: the background should not vary too much during the transient emission and also the amount of statistics should be sufficient to constrain the fit. These assumptions usually are satisfied when the interval of the transient emission is shorter than the Fermi orbital period (≈90 minutes).

For the 2012 March 7 flare presented here, the standard LLE approach to estimate the background by fitting the intervals before and after the flare could not be used because the flaring episode lasted longer than the orbital period. Instead, we estimate the background using data acquired in other orbits. Fermi passes through approximately the same geomagnetic configuration every 15 orbits, but given the standard rocking profile (alternating one orbit north and one orbit south) only the 30^th orbit approximates similar geomagnetic and pointing conditions, and consequently a similar background rate. We do not average multiples of 30 orbits because they become less reliable, as they span observations further removed in time. This method of background estimation has been used for a number of background-dominated instruments in solar flare analyses, including historically for SMM-GRS (Vestrand et al. 1987; Murphy et al. 1990) and EGRET Kanbach et al. (1993)⁷¹ and currently for Fermi GBM (Fitzpatrick et al. 2011).

The background files produced for each interval are saved as standard PHA-I⁷² files. We used XSPEC to execute a forward-folding fit, where the model M(p₁, p₂, ...) is folded with the redistribution matrix, R_ij, and the results are compared with the background-subtracted signal S_i:

$$\begin{equation} S_{i}=C_{i}-B_{i} \stackrel{}{=} \sum _{j}R_{ij}{M}_{j}({p}_1,{p}_2,\ldots); \end{equation} \tag{ A2 }$$

where M_j(p₁, p₂, ...) is the expected number of events between E_j and E_{j + 1} for the time interval being analyzed. A maximum likelihood algorithm is then used to calculate the set of parameters that best model the data (see the XSPEC manual for details.)

A detailed paper on the assessment of the systematic errors is in preparation; here we summarize the main results. Generally speaking, discrepancies between the actual response of the LAT and the response matrix derived from simulations can cause systematic errors in spectral fitting. We investigated the systematic uncertainties tied to the LLE selection by following the procedure described in Abdo et al. (2009b). In particular, we compared Monte Carlo with flight data, using the Vela pulsar (PSR J0835–4510) as a calibration source. The pulsed nature of the gamma-ray emission from this source (Thompson et al. 1975) gives us an independent control on the residual charged particle background. In fact, off-pulse gamma-ray emission is almost entirely absent, and a sample of "pure photons" can be simply extracted from the on-pulse region, after the off-pulse background is subtracted. Considering all time intervals during which the Vela pulsar was observed at an incidence angle θ < 80°, we estimate the discrepancy between the efficiency of the LLE selection criteria in the LAT data and in Monte Carlo to be ∼17% below 100 MeV, decreasing down ∼8% at higher energies, with an average value ∼9% (note that this average is weighted by the Vela spectrum).

Additionally, we performed a spectral analysis of the Vela pulsar, comparing LLE results with standard likelihood analysis. The >100 MeV flux obtained from the LLE analysis is 20% lower than published in Abdo et al. (2009a) and 16% less than the flux reported by Abdo et al. (2010). This discrepancy can be attributed to the fact that the selection criteria between LLE and the standard LAT likelihood analysis are rather different with the later being much looser. Whereas the differences found between the >100 MeV flux obtained from the LLE analysis and those published in Abdo et al. (2009a, 2010) are because our knowledge of the LAT performance improved with time and thus the value of the systematic uncertainties decreased. In our analysis of the impulsive phase of the 2012 March 7, we take this discrepancy into account by modifying the redistribution matrix R_ij in order to reproduce the discrepancy observed on the Vela spectrum. To do this, we scale the redistribution matrix by a conservative ±20% and also we apply a scaling function as a function of the energy, going from ±10% below 250 MeV to ±10% above 250 MeV. These functions maximize the impact of the systematic uncertainty on the value of the measured flux and spectral index, and, at the same time, provide results that are compatible with the observed discrepancy of the Vela spectrum. We use the modified redistribution matrices to repeat the analysis, ousting bracketing value for the flux and the derived proton indexes. These are represented by the grey bands in Figure 1. Note that systematic uncertainties related to the LLE analysis are slightly larger than the ones related to the standard likelihood analysis due to the fact that LLE extends the analysis at lower energy where the value of the uncertainties are typically larger.

Finally, we also studied the energy resolution using large samples of simulated events with the Fermi-LAT full simulator. No significant bias was found, and the energy resolution for LLE is estimated to be ∼40% at 30 MeV, ∼30% at 100 MeV and <15% for energies greater than 100 MeV.

LLE data are generated for each burst (GRB or solar flare) that trigger the GBM. We first bin the data in energy and time, and, following the procedure described in Ackermann et al. (2013), we select the background region by selecting all the LLE events before the trigger and after 300 s from the trigger. We fit each energy bin with a polynomial function of the cosine of the source bore-sight angle as a function of time and we interpolate the background fit into the signal region. We evaluate statistical fluctuation of the signal above the expected background. This procedure is optimized taking into account different signal region and different time binning. We finally look at the post trial probability. Every detected signal with a post trial probability greater than 4σ is promptly made available through the HEASARC Web site.⁷³ For each such detected GRB or solar flare, six different files are delivered.

• 1.  
  The LLE event file format is similar to the LAT photon file format with some exceptions. Because the LLE data are tightly connected to a particular object (position and time), the FITS keyword OBJECT has been added to the file. Generally, OBJECT will correspond to the entry of the GBM Trigger Catalog⁷⁴ used to generate LLE data and corresponds to the "name" column in the FERMILLE table (and in the GBM Trigger Catalog table). The direction of the source used for selecting the data for the LLE file is also written in the header of each extension of each LLE file. PROC_VER corresponds to the iteration of the analysis of LLE data. PASS_VER corresponds to the iteration for the reconstruction and the general event classification (Pass6, Pass7, etc.). VERSION corresponds to the version of the LLE product for the particular GRB or solar flare represented in the file.
• 2.  
  The CSPEC file is obtained from directly binning the LLE event files. It provides a series of spectra, accumulated with 1 s binning (typically from −1000 to 1000 s around the burst). Each spectrum is binned in 50 energy channels, ranging typically from 10 MeV to 100 GeV. The format of the CSPEC file is tailored to satisfy rmfit⁷⁵ standards, and it is not directly usable in XSPEC.
• 3.  
  The CSPEC Response file (the RSP file) is the detector response matrix calculated from Monte Carlo simulation, and it corresponds to a single response matrix for each GRB or solar flare.
• 4.  
  The PHA-I file contains the count spectrum. The PHA-I file is created from the same time interval used to compute the response matrix.
• 5.  
  The selected events file is identical to the LLE event file with an additional selection on time interval applied to match the selection used to compute the response matrix and PHA-I files.
• 6.  
  The LAT pointing and livetime history file is identical to the standard LAT file but with entries every s (instead of every 30 s). It typically spans the range ±4600 s from the trigger time.

The complete LLE selection used to select the events is saved in the keyword LLECUT in the primary header of each LLE file. If the GBM catalog position of the burst is updated (due to a refined localization from the LAT or Swift or from subsequent on-ground analysis), the LLE data are automatically updated and new versions of the LLE files are produced. In some cases, LLE data are manually generated (using a better localization which may or may not have been used in the GBM Trigger Catalog). If the direction of a GRB is revised based on follow-up observations with other instruments, regenerated LLE files will have the VERSION number incremented, but will leave the PASS_VER and PROC_VER unchanged.

In general we do not deliver the background estimates for the time ranges around the burst triggers, and we let the user estimate the background using the procedure described in Pelassa et al. (2010). For the reason explained above, we cannot perform this analysis on the 2012 March 7 flare. Therefore, for this particular flare, we provided LLE data including background files.