FERMI STUDY OF γ-RAY MILLISECOND PULSARS: THE SPECTRAL SHAPE AND PULSED EMISSION FROM J0614–3329 UP TO 60 GeV

For each of the MSP targets, we included all sources within 20° of their positions to make the source models based on the Fermi LAT 4 yr catalog (Acero et al. 2015). The spectral forms of these sources are provided in the catalog. Spectral parameters of the sources within 5° of each target were set as free parameters, and the other parameters were fixed at their catalog values. The catalog spectral model for 33 MSPs is an exponentially cut-off power law, ${dN}/{dE}={N}_{0}{E}^{-{\rm{\Gamma }}}\exp (-E/{E}_{{\rm{c}}})$, while for the other six MSPs (J0610−2100, J1446−4701, J1747−4036, J1125−5825, J1741+1351, and J1823−3021A) it is a simple power law, ${dN}/{dE}={N}_{0}{E}^{-{\rm{\Gamma }}}$.

Using the LAT science tools software package v10r0p5, we performed standard binned likelihood analysis of the LAT data from the MSP targets in the >0.1 GeV band. For PSRs J1658−5324 and J1858−2216, the analysis could not converge, which might be because of the relatively large uncertainties of the instrument response function of the LAT in the low energy range. We thus used >0.2 GeV data instead for the two sources. The spectral results as well as the test statistic (TS) values are given in Table 1 for each source. The TS value at a given position is calculated from TS $=\,-2\mathrm{log}({L}_{0}/{L}_{1})$, where L₀ and L₁ are the maximum likelihood values for a model without and with an additional source respectively. It is a measurement of the improvement in the fit when including the source, and is approximately the square of the detection significance of the source (Abdo et al. 2010).

For the six MSPs with a power-law spectral model in the catalog, we repeated the analysis with an exponentially cut-off power law. The significance of a spectral cutoff was estimated from $\sqrt{-2\mathrm{log}({L}_{\mathrm{pl}}/{L}_{\exp })}$, where L_exp and L_pl are the maximum likelihood values when a target's emission was modeled with a power law with and without the cutoff respectively (Abdo et al. 2013). We found that the spectral cutoff was detected with >3σ significance for the six pulsars. Therefore in Table 1 we provide only the results using the exponentially cut-off power law for them.

We extracted the γ-ray spectra of the MSP targets by performing maximum likelihood analysis of the LAT data in 15 energy bands evenly divided logarithmically from 0.1 to 300 GeV. In the extraction, the spectral normalizations of the sources within 5° of each target were set as free parameters, while all the other parameters of the sources were fixed at the values obtained from the above maximum likelihood analysis. The targets were considered as point sources having power-law emission with Γ fixed at 2.0. The fluxes obtained in this way are less dependent on the overall spectral model assumed for a source, providing a good description of the γ-ray emission of the source. We kept only flux data points with TS greater than 9 (i.e., >3σ significance). A total of 304 data points were obtained for the 39 targets. The flux values for each target are provided in Table 2. We also estimated the systematic uncertainties caused by the Galactic diffuse emission model used. The uncertainty in each energy band was obtained by repeating the likelihood analysis with the normalization of the diffuse component artificially fixed to the values deviating by ±6% from the best-fit value (see, e.g., Abdo et al. 2013). The uncertainties given in Table 2 include the systematic uncertainties. We checked the spectrum and best-fit model for each target. The spectra are well fitted by the spectral models obtained from the likelihood analysis.

In order to obtain a spectral shape that generally defines emission from MSPs, we first normalized the flux of each MSP target with its 0.1–300 GeV energy flux (F₁₀₀ in Table 1). The normalized spectra of the 39 MSPs are shown in Figure 1. We then fit these data points with a normalized, exponentially cut-off power law, i.e., N₀ is obtained from Γ and E_c by requiring the total flux to be 1. The best-fit values we obtained were Γ = 1.5 and E_c = 3.8 GeV, but with a minimum ${\chi }^{2}$ value of 2198 for 302 degrees of freedom. The large ${\chi }^{2}$ reflects the intrinsic spectral differences of the MSPs.

Zoom In Zoom Out Reset image size

We thus used a systematic uncertainty parameter to represent the intrinsic differences. The parameter was added to the uncertainties of the data points in quadrature. We found that when this parameter was set to be 0.05, the minimum reduced ${\chi }^{2}$ was approximately equal to 1. As a result, ${\rm{\Gamma }}={1.54}_{-0.11}^{+0.10}$ and ${E}_{{\rm{c}}}={3.70}_{-0.70}^{+0.95}$ GeV were obtained, where the uncertainties are at a 3σ confidence level. This 3σ spectral region is shown as the gray area in Figure 1.

We performed timing analysis of the 0.1–300 GeV LAT data from the region of PSR J0614−3329 to update the γ-ray ephemeris given in Abdo et al. (2013). An aperture radius of 10 was used. Pulse phases for photons before MJD 55797 (the end time of the known ephemeris) were assigned according to the known ephemeris using the Fermi plugin of TEMPO2 (Edwards et al. 2006; Hobbs et al. 2006). An "empirical Fourier" template profile was built. Using this template, we generated the times of arrival (TOAs) of 40 evenly divided observations of the whole time period. Both the template and TOAs were obtained using the maximum likelihood method described in Ray et al. (2011).

We used TEMPO2 to fit the TOAs. Only the pulse frequency derivative $\dot{f}$ was fitted, and the other timing parameters were fixed to their known values. We obtained $\dot{f}=-1.7559(1)\times {10}^{-15}$ s⁻², consistent with the value given in Abdo et al. (2013) within ∼2.2σ uncertainty. The folded pulse profile and two-dimensional phaseogram are shown in Figure 2. In the following analysis, we selected phases 0.06–0.19 and 0.59–0.78 as the on-pulse phase intervals, and the rest as the off-pulse phase intervals.

Zoom In Zoom Out Reset image size

We included all sources within 20° of the position of PSR J0614−3329 in the Fermi LAT 4 yr catalog (Acero et al. 2015) to make the source model. The spectral forms of these sources are provided in the catalog. Spectral parameters of the sources within 5° of PSR J0614−3329 were set as free parameters, and the other parameters were fixed at their catalog values. The catalog spectral form of PSR J0614−3329 is an exponentially cut-off power law, ${dN}/{dE}={N}_{0}{E}^{-{\rm{\Gamma }}}\exp [-{(E/{E}_{{\rm{c}}})}^{b}]$. The parameter b is a measure of the shape of the exponential cutoff, where a value of 1 or <1 indicates a simple exponential cutoff or a sub-exponential cutoff, respectively. We also used a simple power law in the analysis for comparison.

We performed standard binned likelihood analysis of the LAT data in the energy range >0.1 GeV. We first set b = 1. The γ-ray emission during the whole pulse phase intervals was detected with a TS value of 33,595, while that during the on-pulse and off-pulse phase intervals was detected with TS values of 37,012 and 3781, respectively. We found that during the total, on-pulse, and off-pulse phase intervals, the emission was better modeled with an exponentially cut-off power law. The cutoffs were significantly detected during all the three phase intervals ($\gt 5\sigma ;$ estimated from $\sqrt{-2\mathrm{log}({L}_{\mathrm{pl}}/{L}_{\exp })}$). The resulting power-law fits with simple exponential cutoff are summarized in Table 3.

We then set b as a free parameter and repeated the binned likelihood analysis for the LAT data. We found that during the total and on-pulse phase intervals the sub-exponential cutoffs were detected with ∼4σ significance (estimated from $\sqrt{-2\mathrm{log}({L}_{\exp }/{L}_{\mathrm{subexp}})}$, where L_subexp is the maximum likelihood value for the sub-exponentially cut-off power-law model; Abdo et al. 2013). The resulting sub-exponentially cut-off power-law fits during these two phase intervals are given in Table 3. During the off-pulse phase interval, the sub-exponential cutoff was not detected, because the detection significance was approximately zero.

We extracted the γ-ray spectra of PSR J0614−3329 during the total, on-pulse, and off-pulse phase intervals, by performing maximum likelihood analysis of the LAT data in five energy bands evenly divided logarithmically from 0.1 to 300 GeV. In the extraction, the spectral normalizations of the sources within 5° of PSR J0614−3329 were set as free parameters, while all the other parameters of the sources were fixed at the values obtained from the above maximum likelihood analysis. We kept spectral flux points only when TS was greater than 9 (>3σ significance), and derived 95% flux upper limits otherwise. The obtained spectra are shown in Figure 3, and the fluxes and TS values are provided in Table 4. We found that while the off-pulse emission was detected only for energies <12 GeV, the on-pulse emission from the pulsar was significantly detected at >12 GeV. The TS value at 27 GeV is 589 (see Table 4). In the 60.5–300 GeV band, during the on-pulse phases, we detected only one photon that likely comes from the pulsar (95% probability, see the following Section 4.4 and Table 5), which is not sufficient to claim a detection. Therefore, we provided the 95% upper limit for the flux in this band.

Zoom In Zoom Out Reset image size

Table 4.  Fermi LAT Flux Measurements of PSR J0614−3329

		Total		On-pulse		Off-pulse
E	Band	${E}^{2}{dN}(E)/{dE}$	TS	${E}^{2}{dN}(E)/{dE}$	TS	${E}^{2}{dN}(E)/{dE}$	TS
(GeV)	(GeV)	(10−12 erg cm−2 s−1)		(10−12 erg cm−2 s−1)		(10−12 erg cm−2 s−1)
0.22	0.1–0.5	13.2 ± 0.3	2738	27.7 ± 0.7	3212	6.3 ± 0.4	482
1.10	0.5–2.5	26.6 ± 0.4	18703	64 ± 1	19886	9.3 ± 0.3	2598
5.48	2.5–12.2	29.3 ± 0.7	11512	83 ± 2	12944	4.8 ± 0.4	664
27.16	12.2–60.5	5.7 ± 0.7	445	18 ± 2	589	0.2	0
134.71	60.5–300.0	2.5	8	8.7	16	0.8	0

Note. Fluxes without uncertainties are the 95% upper limits.

Download table as:  ASCIITypeset image

We performed timing analysis of the LAT data of PSR J0614−3329 to search for γ-ray pulsations in the high-energy range, for which we selected the minimum energy as high as possible but also ensured a sufficiently significant detection of pulsation. We found that the value of 25 GeV used in Ackermann et al. (2013) was appropriate. We first selected γ-ray photons within an aperture radius of 05 from PSR J0614−3329, approximately corresponding to the 95% containment angle of the incoming photons from a source. A total of 10 photons were collected. Pulse phases for the photons were assigned using the updated ephemeris obtained in Section 4.1, and an H-test value of 30 was obtained, corresponding to a detection significance of 4.5σ (de Jager & Büsching 2010). Information on these photons is given in Table 5, and they are also shown in the top panel of Figure 4 according to their pulse phases. The 0.1–300 GeV pulse profile is shown in the bottom panel of Figure 4 for comparison.

Zoom In Zoom Out Reset image size

We then used a larger aperture radius of 2° to include more photons (40 photons were collected), and weighted them by their probability of originating from the pulsar (calculated with using gtsrcprob). Pulse phases for these photons were assigned and a weighted H-test value of 48 was obtained (de Jager & Büsching 2010; Kerr 2011), which corresponds to a detection significance of ∼6σ, indicating that the γ-ray pulsation from the source was significantly detected in the >25 GeV energy band. We plotted the weighted photons in the middle panel of Figure 4 according to their pulse phases. The weights for the 10 photons within the 95% containment angle are given in Table 5.

We also performed likelihood analysis of the >25 GeV data during the on-pulse and off-pulse phase intervals, with the emission from the source modeled with a simple power law. The γ-ray emission from PSR J0614−3329 was detected with ${\rm{TS}}\,\simeq \,63$, having ${\rm{\Gamma }}=2.7\pm 0.8$ and photon flux ${F}_{25-300}=(1.0\pm 0.3)\times {10}^{-10}\,\mathrm{photons}\,{{\rm{s}}}^{-1}\,{\mathrm{cm}}^{-2}$ during the on-pulse phase intervals. During the off-pulse phase intervals, the γ-ray emission was not detected (TS $\simeq \,0$), and the derived 95% upper limit on photon flux is $8\times {10}^{-12}\,\mathrm{photons}\,{{\rm{s}}}^{-1}\,{\mathrm{cm}}^{-2}$. Two TS maps during these two phase intervals are shown in Figure 5. We ran gtfindsrc in the LAT software package to determine the position during the on-pulse phase intervals and obtained R.A. = 9353, decl. = −3350 (equinox J2000.0), with 1σ nominal uncertainty of 002. PSR J0614−3329 is 001 from the best-fit position and within the 1σ error circle. This result confirmed the detection of pulsed γ-ray emission from photon folding.

Zoom In Zoom Out Reset image size

Having analyzed approximately 7.5 yr of Fermi LAT data for 39 MSPs reported in 2PC, we have obtained their spectra, which are all well described by a power law with exponential cutoff. We have thus determined their general spectral shape by fitting the spectra with such a function. The allowed region of spectral shape for MSPs is relatively large because of the intrinsic differences in their spectra. However, this spectral shape can be used to find candidate MSPs among the unidentified LAT sources. For example, using the criteria of >5°, significant curvature in a spectrum, and non-variability, we (Dai et al. 2016) have found 24 such sources from the Fermi LAT third source catalog (Acero et al. 2015), but two of them, J0318.1+0252 and J2053.9+2922, likely have spectra not consistent with the spectral shape of MSPs because either Γ or E_c found for γ-ray emission from the two sources is not in the range of spectral shape. Based on the known properties of the different types of LAT sources, they could be MSPs with quite different spectra or even other types of sources, for example the dark matter subhalo candidates as suggested by Bertoni et al. (2015). For the purpose of finding candidate MSPs and other types of unidentified γ-ray sources, it is reasonable to search through LAT sources at high Galactic latitudes by comparing their spectra with the spectral shape of MSPs that we have determined.

High-energy γ-ray emission is seen from 27 pulsars, as reported in the first Fermi catalog of sources above 10 GeV (Ackermann et al. 2013). Among them 20 sources were found to have pulsed γ-ray emission in the >10 GeV band, including 17 young pulsars and three MSPs. Furthermore, PSR J0614−3329 was one of 12 pulsars found to have γ-ray pulsations in the >25 GeV band, although it was only marginally detected in Ackermann et al. (2013). Our analysis has shown, likely due to the longer time period of data (7.5 yr versus 3 yr) and overall improvement in sensitivity in the Pass 8 data, that there is significant pulsed γ-ray emission in the >25 GeV band from this MSP.

Ten photons were detected within the 95% containment angle of PSR J0614−3329, including seven photons in the on-pulse phase ranges with high probabilities of coming from the pulsar (0.94–1.0, Table 5). The maximum energy of these seven photons is ∼60 GeV, indicating that the pulsed γ-ray emission from PSR J0614−3329 can be detected up to ∼60 GeV. In addition, there are two ∼200 GeV photons detected (Table 5). While the probabilities that they come from the pulsar are low (0.06 and 0.01), they are both in the on-pulse phase ranges. Considering that the contamination from the background emission decreases with increasing energy, it is possible that the two ∼200 GeV photons indeed originated from the pulsar. If this is true, the pulsed γ-ray emission from PSR J0614−3329 could stretch up to 200 GeV. It is noted that the ∼200 GeV photons are in the second peak of the pulse profile (see Table 5), which is consistent with the common feature of γ-ray pulsars that the second peaks are more prominent in higher energy bands (Abdo et al. 2013).

We note that PSR J0614−3329 was not listed in the second catalog of hard Fermi LAT sources (seen above 50 GeV; Ackermann et al. 2016). This would be because their data analysis was conducted over the whole pulse phase for this pulsar. In our analysis, we also found that PSR J0614−3329 was not significantly detected in the >60 GeV band, and TS $\simeq \,8$ over the total phase range (Table 4). However, there were three photons in the >60 GeV band within the 95% containment angle, one at ∼61 GeV and two at ∼200 GeV, and they are all in the on-pulse phase ranges. The time periods of the on-pulse phase ranges were approximately 30% of the total time periods, which largely reduced the background emission. The timing analysis thus has allowed us to increase the TS value from 8 to 16 (4σ detection; Table 4).

These results have added PSR J0614−3329 to the group of pulsars that have been found to have pulsed emission at energies of tens of GeV (e.g., Ackermann et al. 2013; Harding & Kalapotharakos 2015 and references therein). The mechanism of the very high-energy emission from pulsars remains to be solved. Currently the inverse-Compton scattering process in the outer magnetosphere or the pulsar wind region is considered to produce the pulsed emission detected in the $\gt 10$ GeV band from the Crab pulsar (see, e.g., Aleksić et al. 2011; Aharonian et al. 2012; Lyutikov 2013; Harding & Kalapotharakos 2015). Alternatively a non-stationary outer-gap scenario has also been proposed recently (Takata et al. 2016), which has been used to interpret the $\gt 50$ GeV pulsed emission from the Vela pulsar (Leung et al. 2014). In this scenario, the observed spectrum of a pulsar is the superposition of emission from the variable outer-gap structures. In Figure 3, we show the model fits to γ-ray emission from the Crab (Aleksić et al. 2011; Abdo et al. 2013) and Vela (Leung et al. 2014) pulsars (scaled by their 0.1–100 GeV total LAT fluxes respectively) for comparison. PSR J0614−3329 possibly has a stronger ∼200 GeV component than the Crab and Vela pulsars, although the large uncertainty does not allow a clear conclusion to be drawn. In order to investigate the emission process responsible for the high-energy component from PSR J0614−3329, detailed modeling (such as that in Harding & Kalapotharakos 2015; Takata et al. 2016) is needed.