SENSITIVITY OF BLIND PULSAR SEARCHES WITH THE FERMI LARGE AREA TELESCOPE

Observed pulsars emitting radio and/or gamma rays do not uniformly cover the $f\hbox{--}\dot{f}$ space. The LAT has so far only observed pulsations from pulsars with f > 2 Hz and $|\dot{f}| < 3.7 \times 10^{-10}$ Hz s⁻¹. We uniformly sample a broader region of the $f\hbox{--}\dot{f}$ parameter space to quantify the observational bias. Our simulated pulsars sample the spin parameter space in a log-uniform random distribution with (0.1 ⩽ f ⩽ 64) Hz, and ($10^{-15} \le -\dot{f} \le 10^{-10}$) Hz s⁻¹ (see shaded region in Figure 1). In this fashion, the Crab is the only LAT pulsar not covered in the simulation. We are not sampling the recycled pulsar domain (i.e., millisecond pulsars, MSPs) and keep to the normal population of younger pulsars. Within the sampled range, if pulsars are systematically detected in a region of the parameter space where they are not observed by the LAT, these objects do not occur in nature. Conversely, if pulsars are systematically not detected in a region of the parameter space, we cannot infer any conclusion about the underlying pulsar population. Observed ratio pulsars in the ATNF (Australia Telescope National Facility) database¹¹ (Manchester et al. 2005) cover a part of the simulated parameter space not observed in gamma rays, especially at lower $\dot{E}$.

Zoom In Zoom Out Reset image size

We select the epoch to be in the middle of the year-long observation (MJD 54866). This choice minimizes the covariances between timing model parameters. Young and middle-aged pulsars can exhibit glitches (Espinoza et al. 2011; M. Dormody et al. 2011a, in preparation), an abrupt increase in spin frequency. They can also exhibit timing noise (Edwards et al. 2006; Ray et al. 2011), smooth deviations from the expected linear spin-down. Any or all of these effects can hamper detections, but we do not include them in our simulations as including these effects is beyond the scope of this study.

We select the pulse profile based on the general features in the pulse profiles of the known gamma-ray pulsars: broad, mostly double-peaked (Abdo et al. 2010c). We assume no spectral modulation in phase and fix the pulse profile to be the same across all energies, which is modeled by the sum of two Lorentzians and a constant:

$$\begin{equation} n(\phi) = \frac{I_1}{1+ ((\phi - \phi _1)/w_1)^2} + \frac{I_2}{1+((\phi - \phi _2)/w_2)^2} + n_0. \end{equation} \tag{ 1 }$$

Here, the phase ϕ ranges from 0 ⩽ ϕ ⩽ 1. I₁ and I₂ are the amplitudes of the two peaks and are selected from a random, uniform distribution between 0 and 100. This can lead to a single-peaked profile if one amplitude happens to be sampled low enough. The unpulsed component n₀ is selected from a random, uniform distribution between 0 and 25. The widths of the peaks are fixed at w₁ = 0.01 and w₂ = 0.015. The phase separations of the peaks are fixed at ϕ₁ = 0.25 and ϕ₂ = 0.65. The parameters used in the pulse profile template in Equation (1) are reasonable based on LAT observations of gamma-ray pulsars (Abdo et al. 2010c). While the separation of the peaks can affect the detectability of a pulsar, the majority of observed gamma-ray pulsars have peak separations Δ ∼ 0.4–0.5. The pulse profile is sampled to 50 bins in phase across one rotation.

The pulsed fraction PF is calculated by integrating the pulsed portion of Equation (1) and dividing by the integral of the pulsed and unpulsed portion:

$$\begin{equation} {\rm PF} = \left[1 + \frac{n_0}{\pi (I_1 w_1 + I_2 w_2) } \right]^{-1}, \end{equation} \tag{ 2 }$$

where 0 ⩽ PF ⩽ 1. In Figure 2, we show two examples of pulse profiles: one has a large PF and the other has a small PF. While the pulse profiles sampled here cover a wide range of gamma-ray profiles, they do not cover profiles with three or more peaks, narrowly separated peaks, or very broad peaks that are difficult to disentangle from the unpulsed component. In addition, some LAT pulsars exhibit an energy-dependent pulse profile, but this is not incorporated into this simulation since these represent a small number of exceptional cases.

Zoom In Zoom Out Reset image size

Gamma-ray pulsar spectra are well modeled by a power law with an exponential cutoff (Abdo et al. 2010c):

$$\begin{equation} \frac{dN}{dE dA dt} = C \left(\frac{E}{1 \,{\rm GeV}} \right)^{-\Gamma } \exp \left[ -\left(\frac{E}{E_{{\rm cutoff}}} \right)^{\beta } \right], \end{equation} \tag{ 3 }$$

where the spectral parameters include the spectral index Γ, cutoff energy E_cutoff, cutoff index β (β < 1 for sub-exponential, β > 1 for super-exponential), and normalization factor C. While the brightest pulsars are known to have phase-dependent energy spectra (Abdo et al. 2009c, 2010d), we assume a uniform, phase-averaged spectrum. Since pulsar detections are carried out over a wide energy band, such spectral variations are second-order effects. In addition, Fermi LAT pulsars can be associated with supernova remnant and pulsar wind nebula sources, which can affect the measured unpulsed spectrum and flux stability of the pulsar source (Ackermann et al. 2011; Abdo et al. 2011).

We vary Γ, E_cutoff, and β uniformly over the ranges (0.6 ⩽ Γ ⩽ 2.1), (1 ⩽ E_cutoff ⩽ 2) GeV, and (0.5 ⩽ β ⩽ 4.5), respectively. Gamma-ray pulsars have spectral indices spanning this range, and cutoff energies between (0.7 ± 0.5) GeV (PSR J0659+1414) and (5.8 ± 0.5) GeV (Crab pulsar), which is much wider than the sampled cutoffs. Since the pulsar spectral parameters are only used to determine the number of pulsar photons, and all photons are weighted equally in the blind search, the narrow sampling of E_cutoff is sufficient for pulsar spectra in this simulation. In addition, Equation (3) can model the exponential cutoff (β = 1) predicted by magnetospheric emission models or the super-exponential cutoff (β > 1) predicted by the polar cap emission model and the low-altitude slot gap model (Muslimov & Harding 2003). We sample below β = 1 as some Fermi LAT pulsars exhibit a sub-exponential cutoff (Abdo et al. 2010b).

Pulsars are detected in blind searches from F₁₀₀ ∼ 10⁻⁵ photon cm⁻² s⁻¹ (Vela) down to F₁₀₀ = 6 × 10⁻⁸ photon cm⁻² s⁻¹ (J2032+4127; Abdo et al. 2010c), where F₁₀₀ is the flux of the source above 100 MeV. We sample our simulated pulsars in a log-uniform, random distribution of (10⁻⁹ ⩽ F₁₀₀ ⩽ 10⁻⁶) photon cm⁻² s⁻¹. This allows us to explore the detectability of high-flux pulsars with different spectral models, as well as establish the sensitivity threshold. We obtain the normalization factor C by integrating the spectrum in Equation (3) over all energy bins:

$$\begin{equation} F_{100} = C \int _{0.1 \,{\rm GeV}}^{20 \,{\rm GeV}} \left(\frac{E}{1\, {\rm GeV}} \right)^{-\Gamma } \exp \left[ -\left(\frac{E}{E_{{\rm cutoff}}} \right)^{\beta } \right] dE. \quad \end{equation} \tag{ 4 }$$

We only integrate up to 20 GeV as few photons contribute to the flux above this energy range.

Most pulsars lie within the Galactic plane (|b| < 5°). This is usually explained by the fact that young and middle-aged pulsars have not had time to move far from their birth locations in the plane (Paczynski 1990). We sample our simulated pulsar source locations from the (l, b) distribution of the known pulsar population in the ATNF database. The probability distribution is modeled with a broken power law for l and |b|:

$$\begin{equation} p(x) =\left\lbrace \begin{array}{@{}l@{\quad }l@{}}x^{\alpha _1}/\alpha _1 & x < B \\ B^{\alpha _1}/\alpha _1 + B^{\alpha _1-\alpha _2}/\alpha _2 \times (x^{\alpha _2} - B^{\alpha _2}) & {\rm otherwise.} \end{array}\right. \end{equation} \tag{ 5 }$$

For 180° ⩽ l ⩽ 360°, x = 360 − l, α₁ = 0.9998, α₂ = −1.45, B = 44.87, and for 0° ⩽ b ⩽ 90°, x = |b|, α₁ = 0.75, α₂ = −0.26, B = 1.23. These values were obtained with a fit to the pulsars in the ATNF database. While this selection includes MSPs (f > 100 Hz) which tend to have a broader latitude distribution, they comprise only ∼10% of pulsars and only affect the fit parameters by a few percent. We assumed symmetry around the Galactic meridian l = 0 and b = 0 (see the plot of the simulated positions in Figure 3).

Zoom In Zoom Out Reset image size

Since the blind search relies on source positions, i.e., it is not positionally blind, pulsations are easier to detect from bright, well-localized pulsars. Conversely, pulsations are more difficult to detect from faint, poorly localized pulsars. We localize the simulated pulsars according to the procedure described in Abdo et al. (2009b, 2010a). We bin the photon data into three overlapping energy bands: Front¹² events with E > 200 MeV or Back events with E > 400 MeV and 03 pixels; Front events with E > 1 GeV or Back events with E > 2 GeV and 02 pixels; and Front events with E > 5 GeV or Back events with E > 10 GeV and 01 pixels. We used two different wavelet detection methods: mr_filter (Starck & Pierre 1998) and PGWAVE (Damiani et al. 1997; Ciprini et al. 2007). For mr_filter, the False Discovery Rate (Benjamini & Hochberg 1995) was set to 1.2% (2.5σ) in the first two bands and 0.27% (3σ) in the highest energy band. We used SExtractor (Bertin & Arnouts 1996; Holwerda 2005) to extract the sources from the mr_filter output. PGWAVE returns a list of source seeds directly.

We also used pointfind, a tool that searches for candidate point sources by maximizing a likelihood function for trial point sources at each direction in a HEALPix (Górski et al. 2005) order 9 tessellation of the sky (01 × 01 pixels). The pointfind algorithm divides the data into four energy bins per decade and determines the significance of a trial point source against the diffuse background in each pixel. Further trials optimize the likelihood with respect to the signal fraction. A second pass uses a more complicated likelihood function incorporating nearby point sources. For each source, pointfind generates a Test Statistic (TS) map. The TS is defined as 2(log L_max − log L(p)), where p is the source position and L_max is the maximum value of the likelihood L(p). Contours of source location probability are defined by interpreting deviations of T from a peak value (i.e., the maximum likelihood position of a pulsar) in terms of a χ² distribution with 2 degrees of freedom, and can be used to reject positions falling below a given confidence threshold. We extract the positions of candidate point sources (seeds) from the TS maps (Abdo et al. 2010a).

We considered two localization tools to refine the coarse seed positions: pointfit and gtfindsrc.¹³ The tool gtfindsrc is not optimal for this study, as it can only perform an unbinned likelihood analysis, which is prohibitive for one year of photon data and thousands of sources. It also requires a large overhead in the form of individual exposure maps and maps of the diffuse model convolved with the IRF. We used a Galactic diffuse model designated 54_86Xexph7S calculated using GALPROP,¹⁴ an evolution of that used in the 0FGL. A similar model, gll_iem_v02, is publicly available.¹⁵ The pointfit tool performs the same task but only requires a diffuse background model and a photon data file. Thus, we used pointfit for this step in the analysis. We compared these two methods using a sample of 100 simulated pulsars and found that the localizations were consistent with each other.

The pointfit tool is similar to pointfind in that it is a likelihood analysis tool that maximizes a binned likelihood function and returns source locations along with TS values. However, it is used for source localizations of seed detections, obtained from mr$\_$filter, PGWAVE, and pointfind. The overall likelihood function L is the product of the energy band likelihoods L_i. The choice of pointfind and pointfit is used in order to align closely with the procedure in the 1FGL catalog.

The detection of our simulated sources, as in the case of real sources, depends on the source flux and position. As shown in Figure 4, bright sources tend to be well localized (upper left plot). Moderately bright sources located within the diffuse emission along the Galactic plane or near other bright pulsars can be affected by source confusion, which can worsen the localization (upper right plot). Faint pulsar sources, while formally detected just above the threshold (TS = 25) are poorly localized (lower left plot), and extremely faint pulsars are not detected at all (lower right plot). It is safe to assume a source not reasonably localized would not be detected as a pulsar in a blind search. In order to restrict our searches to reasonably localized sources, we perform the blind search only on sources localized to within 2° of the simulated location. This is a fairly conservative criterion, as all sources included in the 1FGL catalog have a 95% error radius smaller than 06 (Abdo et al. 2010a).

Zoom In Zoom Out Reset image size

While source localizations from gamma-ray observations are limited by the relatively broad PSF and limited photon statistics, observations from other wavelengths can have much finer spatial resolution and allow for precise source localization. For gamma-ray pulsars with possible counterpart locations (e.g., X-ray or radio observations) we can perform searches for pulsation at these more precise positions to aid in blind searches. We investigate the effect of precise source locations on blind search detectability in Section 5.1.

We replicate the analysis from Abdo et al. (2009a), Saz Parkinson et al. (2010), and P. M. Saz Parkinson et al. (in preparation) using the same code, settings, and search strategy. The specific blind search method in this analysis searches the time differences of photons from simulated pulsars using a maximum time-differencing window of T_maxdiff = 2¹⁹ s (∼6 days). We perform a fast Fourier transform (FFT) between 0.1 and 64 Hz with photon times corrected for the spin-down $\dot{f}$:

$$\begin{equation} t^{\prime } = t + \frac{1}{2} \frac{\dot{f}}{f} t^2, \end{equation} \tag{ 6 }$$

where $\dot{f}$ is rounded to the search resolution of $\Delta \dot{f} = 1/(2 T_v T_{{\rm maxdiff}}) \approx 3 \times 10^{-14}$ for one year of observations, and T_v is the one year viewing period. We used an ROI in the photon selection (r ⩽ 08, E_min = 300 MeV, P6_V3 diffuse class photons, with zenith angle ⩽105° to avoid photons from Earth's limb). We corrected the photon event times from each source to the solar system barycenter using gtbary, assuming all photons come from the source positions given by the LAT localization tools as described in Section 3.1. We look for promising candidates for follow-up based on of the following two criteria.

• 1.  
  The Fourier power has a probability of false detection of p < 10⁻⁴, taking into account the number of bin trials in the FFT of the time differences.
• 2.  
  After a fully coherent refinement, the pulse profile with 32 bins has χ² ⩾ 5 for 31 degrees of freedom when compared with a flat distribution.

We confirm the detection by checking the value of the detected frequency with the simulated frequency: f_simulated − f_detected < 10⁻⁵ Hz, or 2T_maxdiff × 10⁻⁵ ≈ 10 frequency bins. This is a conservative estimate, as we do not expect to detect a pulsar if its frequency is shifted by more than a few bins.

After all sources are flagged as detections or non-detections, we remove duplicate detections due to multiple source locations from the seed detections, a frequent occurrence for bright sources in confused backgrounds. We also remove faint candidates that are within a few degrees of bright pulsars (e.g., Vela, Geminga, and Crab) that are simulated with a similar $\dot{f} / f$ as these bright pulsars (<1% of simulated pulsars). These pulsars would require non-standard ways of extracting the pulsed emission of the faint pulsar from the bright pulsar.

The simulated pulsars reveal that the lower the S/N, the larger the uncertainty in the position, and thus the average offset, as shown in Figure 7. If we model the observed position as a two-dimensional Gaussian distribution centered at the correct position, the offset in position can be approximated by another Gaussian distribution, peaked at σ, and with a width σ. This result confirms the behavior already reported in Figure 2 in the 1FGL catalog (Abdo et al. 2010c), where the position offset has a strong dependence on the TS. Assuming an exponential dependence of σ with respect to the S/N, we fit the data to obtain a value for the exponential index:

$$\begin{equation} \sigma \left({\rm S/N}\right)\sim e^{-0.034 \times {\rm S/N}}. \end{equation} \tag{ 7 }$$

A pulsed signal is detected from a source if there is a large enough portion of pulsed photons. We define the detectability of the pulsar using the metric M as the logarithm of the ratio of pulsed photons to the square root of all photons:

$$\begin{equation} M = \log _{10} \frac{n_{{\rm pulsed}}}{\sqrt{n}}, \end{equation} \tag{ 8 }$$

where n_pulsed is the number of signal photons contained in the pulse, or are above the unpulsed background in the pulse profile, and n = n_signal + n_background. We use a logarithm in the definition of M to better observe the detectability of weaker sources with lower S/Ns. Pulsed and unpulsed photons are together classified as signal photons, or n_signal = n_pulsed + n_unpulsed. Pulsed photons can be selected by phase cuts, and unpulsed photons can be distinguished from the background by the spectral model and PSF angular distribution. The detectability of the pulsar does not directly depend on n_pulsed, n_unpulsed, or n_background, but only on the metric M. This is motivated by the indistinguishability of n_unpulsed and n_background as shown in Figure 8: an unpulsed background is composed of unknown portions of background diffuse photons and unpulsed pulsar photons.

If the number of pulsed photons is known from the pulse shape a priori, M is calculated easily using Equation (8). However, if the pulse profile is unknown, M is calculated by estimating the PF and using a spectral model to calculate the S/N:

$$\begin{equation} M = \log _{10} \left(\underbrace{\frac{n_{{\rm pulsed}}}{n_{{\rm pulsed}} + n_{{\rm unpulsed}}}}_{{\rm PF}} \times \underbrace{\frac{n_{{\rm signal}}}{\sqrt{n_{{\rm signal}} + n_{{\rm background}}}}}_{{\rm S/N}} \right). \end{equation} \tag{ 9 }$$

The expected number of signal photons n_signal from a source falling within a given ROI can be determined from the pulsar spectrum, flux, and IRFs. We select diffuse class photons satisfying the following cuts: r ⩽ 08, where r is the angular separation from the source position; E ⩾ 300 MeV, diffuse class photons. The 68% containment angle θ₆₈ for normal incidence (cos θ > 0.9) of the front and back converters (see Section 3.1) for P6_V3 diffuse photons¹⁶ can be fitted with a third-order polynomial function given by a fit to the PSF:

$$\begin{equation} \log _{10} \theta _{68} = p_0 + p_1 \log _{10} E + p_2 (\log _{10} E)^2 + p_3 (\log _{10} E)^3 \end{equation} \tag{ 10 }$$

for 0.1 ⩽ E ⩽ 20, where E is in GeV, p₀ = −0.13, p₁ = −0.77, p₂ = 0.06, and p₃ = 0.04. We approximate all simulated photon events to be within 45° of normal incidence, since the PSF containment angle θ₆₈ is nearly constant from normal incidence out to 45° off axis.

We assume a point source to have a two-dimensional Gaussian (normal) distribution of gamma rays $f(x, y) = e^{-\left(\frac{x^2 + y^2}{2 \sigma ^2} \right) }$, so the ratio, R, of source photons falling within the ROI to all source photons is given as

$$\begin{eqnarray} R(r, \sigma) &=& \frac{\int _0^{r} r^{\prime } \exp \left(-r^{\prime 2}/2 \sigma ^2 \right) dr^{\prime } \int _0^{2\pi } d\theta }{\int _0^{\infty } r^{\prime } \exp \left(-r^{\prime 2}/2 \sigma ^2 \right) dr^{\prime } \int _0^{2\pi } d\theta }\nonumber\\ &=& 1 - \exp \left(-r^2/2 \sigma ^2 \right), \end{eqnarray} \tag{ 11 }$$

where r is the angular separation from the source position, and the 68% containment angle in degrees is θ₆₈ = 1.51σ. The fraction of source photons falling within a 08 radial cut is

$$\begin{eqnarray} R(r=0\mbox{$.\!\!^\circ $}8,\theta _{68}) &=& 1 - \exp \left[- \frac{0.8^2}{2 \left(\theta _{68}/1.51\right)^2} \right]\nonumber\\ &=& 1 - \exp \big({-}0.73/\theta _{68}^2 \big), \end{eqnarray} \tag{ 12 }$$

where log₁₀θ₆₈ is given in Equation (10). The number of detected photons above 300 MeV falling within 08 is

$$\begin{eqnarray} &&n_{{\rm signal}} =\nonumber\\ &&k F_{100} \frac{\int _{0.3 \,{\rm GeV}}^{20 \,{\rm GeV}} R(r=0\mbox{$.\!\!^\circ $}8,\theta _{68}) E^{-\Gamma } \exp \left[ - \left(E/E_{{\rm cutoff}} \right)^\beta \right] dE}{\int _{0.1 \,{\rm GeV}}^{20 \,{\rm GeV}} E^{-\Gamma } \exp \left[ - \left(E/E_{{\rm cutoff}} \right)^\beta \right] dE},\nonumber\\ \end{eqnarray} \tag{ 13 }$$

where log₁₀k = 10.47 ± 0.12 and k has units cm² s. k is a constant that represents the instrument effective area times the time spent observing a source, which depends on the incident angle and source spectrum. The spread in the fit is partially due to the assumption of normal incidence in Equation (10), non-uniform all-sky exposure during the first year of observations, and Poisson noise.

There is a marked drop-off in sensitivity to pulsars with high frequencies (f > 16 Hz), as seen in Figure 6. This is also shown in the detectability metric, where high-frequency pulsars are less detectable than low-frequency ones (see Figure 9). M depends only on the S/N and the PF, which is distributed independently from the frequency. In fact, at low frequencies the rate of detection depends only on the value of the metric M. But at high frequencies we observe a drop-off in the detection rate. This trend is more apparent at low S/N, where on average the source is weaker, the number of pulsed photons smaller, and the offset larger.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We estimate the detection probability associated with the metric M (Equation (9)) as a fit to the logistic sigmoid function, with the constraint P(M → −∞) = 0:

$$\begin{equation} P(M) = \frac{A}{1 + \exp \left((\mu -M)/\sigma \right) }. \end{equation} \tag{ 14 }$$

Assuming that pulsar detections follow a binomial distribution with expected probability p, we maximize the likelihood function over the parameters A, μ, and σ. We estimate the 68% uncertainties in the parameters assuming that the likelihood is approximated by an asymmetric Gaussian around the maximum. We get the values A = 0.857^+0.023 _{− 0.025}, μ = 0.979^+0.009 _{− 0.010}, and σ = 0.057 ± 0.006. We also parameterized A, μ, and σ as linear functions of spin frequency in order to illustrate the effect of large spin frequencies on detectability. We obtain the fit values A = A₁f + A₂ (A₁ = −0.002 ± 0.004, A₂ = 0.93 ± 0.02), μ = μ₁f + μ₂ (μ₁ = 0.01 ± 0.002, μ₂ = 0.94^+0.02 _{− 0.01}), and σ = σ₁f + σ₂ (σ₁ = 0.001 ± 0.001, σ₂ = 0.04^+0.007 _{− 0.006}).

While a faint gamma-ray pulsar yields a poor LAT position, it might have a strong, well-localized, counterpart source in another waveband. We created such a counterpart location by sampling a random angular distance 0'' ⩽ r ⩽ 10'' in a random orientation from the simulated pulsar location. This range is the typical positional uncertainty for X-ray instruments. These source locations are more accurate than the LAT locations and yield a higher detection probability. We performed the blind search on all pulsars using these simulated counterpart locations and detected 1359 pulsars (compare with 1121 using simulated LAT locations). We obtained fit parameters for the counterpart detectability metric by maximizing the likelihood function over the parameters A, μ, and σ to get the values A = 0.949^+0.015 _{− 0.012}, μ = 0.939 ± 0.007, and σ = 0.048 ± 0.004 (see Figure 10). In addition, we parameterized A, μ, and σ as linear functions of spin frequency. We get A = A₁f + A₂ (A₁ = −0.005 ± 0.002, A₂ = 1.0^+0.0 _{− 0.003}), μ = μ₁f + μ₂ (μ₁ = 0.005 ± 0.001, μ₂ = 0.92 ± 0.01), and σ = σ₁f + σ₂ (σ₁ = 0.000 ± 0.001, σ₂ = 0.042^+0.006 _{− 0.004}). We fix A₂ ⩽ 1.0 in the fit and provide a lower limit.

Zoom In Zoom Out Reset image size

To see the benefit of using counterpart source locations in blind searches, we searched for the known Fermi LAT pulsars (Abdo et al. 2010c; Saz Parkinson et al. 2010) using two locations: the simulated LAT location and the simulated counterpart location. Here, the simulated LAT location is obtained as in Section 3.1, and the counterpart location is randomized between 0'' and 10'' around the timed position. Using the optimized LAT location we detect all 26 blind search pulsars (Abdo et al. 2009a; Saz Parkinson et al. 2010; P. M. Saz Parkinson et al., in preparation), all 6 bright pulsars detected by EGRET (Nolan et al. 1996), and 5 radio-loud pulsars included in the 1FGL (Abdo et al. 2010c): PSRs J1028−5819, J1048−5832, J1509−5850, J2021+3651, and J2229+6114. Using the simulated counterpart locations, in addition to detecting the previous 35 pulsars, we also detect 2 additional radio-loud pulsars, J0659+1414 and J1420−6048.

As shown in Section 4.1, the detectability of a pulsar depends greatly on the uncertainty in position. In Section 5.1, we showed that using multiwavelength counterpart locations reduces the effect of the Doppler shift frequency modulation and allows for dimmer pulsars to be detected. Details are described below.

The detection of high-frequency pulsars is more sensitive to positional errors (e.g., Ransom 2001; Chandler et al. 2001; Ray et al. 2011). The positional error of a pulsar smears the peak in the DFT by the amount

$$\begin{equation} \delta f_{{\rm pos}} = \frac{v}{c} f \epsilon | \sin \theta | = \frac{v}{c} f \epsilon \sqrt{1 - (\cos \phi \sin \omega t)^2}, \end{equation} \tag{ 15 }$$

where is given in radians, v is Earth's orbital velocity (v/c ∼ 10⁻⁴), and θ is the angle between v and the source direction. The angle θ is a combination of the azimuthal angle ωt and the polar angle ϕ, where ω = 2 × 10⁻⁷ rad s⁻¹. In this convention, an object directly overhead the ecliptic plane is located at $\phi = \frac{\pi }{2}$ and an object within the plane of the solar system is located at ϕ = 0.

The corresponding shift in frequency derivative due to positional error is calculated by differentiating Equation (15):

$$\begin{equation} \delta \dot{f}_{{\rm pos}} = 10^{-4} \epsilon (f \dot{\theta }\cos \theta + \dot{f} \sin \theta). \end{equation} \tag{ 16 }$$

The second term is negligible when compared with the first, as $\dot{f}\sim 10^{-12}$ and $\dot{\theta }\sim 10^{-7}$:

$$\begin{equation} \delta \dot{f}_{{\rm pos}} = -10^{-4} \epsilon f \frac{\omega \cos ^2 \phi | \cos \omega t \sin \omega t|}{\sqrt{1-(\cos \phi \sin \omega t)^2}}. \end{equation} \tag{ 17 }$$

If the frequency is shifted by more than a few DFT bins, the power is significantly reduced. The maximum frequency shift δf due to positional error is when sin θ = 1:

$$\begin{equation} | \delta f_{{\rm max}} | = \frac{v}{c} \epsilon f. \end{equation} \tag{ 18 }$$

The number of frequency bins shifted is equal to δf_max divided by the frequency resolution Δf = 1/2T_maxdiff. If we require the frequency shift by no more than N frequency bins

$$\begin{equation} \frac{\delta f_{{\rm max}}}{\Delta f} \le N \end{equation} \tag{ 19 }$$

$$\begin{equation} \epsilon \le \frac{N}{2 T_{{\rm maxdiff}} (v/c) f}. \end{equation} \tag{ 20 }$$

Setting T_maxdiff = 2¹⁹ s, v/c = 10⁻⁴ and converting from radians to degrees

$$\begin{equation} \epsilon \le \frac{0.54 \times N}{f}, \end{equation} \tag{ 21 }$$

where is given in degrees and f is given in Hz. The plot of all detected pulsar frequencies versus positional offset is shown in Figure 11, with a selection of different N values.

Zoom In Zoom Out Reset image size

While pulsar searches are sensitive to position errors, one important result from Equation (21) is that as the frequency of the pulsar increases we require a position with greater certainty. Therefore, we can increase our sensitivity to high-frequency pulsar detections by homing in on the position. This idea is illustrated in the plot of position offset and frequency for all detected pulsars (Figure 11). As both f and the position offset become large, there is a drop-off in the number of detections. This is due to both the inability to combine enough higher harmonics with a fixed maximum search frequency f_max = 64 Hz and the sensitive dependence on position.

The steps in obtaining the pulsar detection probability P are summarized as follows.

• 1.  
  Assume a spectral model and pulsed fraction, e.g., Γ = 1.5, E_cutoff = 2 GeV, β = 1, F₁₀₀ = 3 × 10⁻⁸ photon cm⁻² s⁻¹, and PF = 0.5.
• 2.  
  Measure the total number of photons n within the ROI and estimate the number of signal photons using Equation (13).
• 3.  
  Calculate the value of the metric M using Equation (9).
• 4.  
  Apply Equation (14) with the corresponding A, μ, and σ values for the source positions to determine P.

This recipe is only for a blind search using the time differencing method for one year of LAT photon data processed with the P6$\_$V3 diffuse IRF. We can use this recipe to measure detection 1σ upper limits on the pulsed fraction by fixing the detection probability at P(M) = 0.68 and solving for M with a tail probability of 0.32. For a 1σ detection, M₆₈ = 1.056 for LAT source locations and M₆₈ = 0.985 for counterpart locations.

As an application of Section 5.3, we generated an all-sky detectability flux map assuming a spectral model and pulse profile. An all-sky blind search sensitivity map is useful in determining the probability of detection for a given source location. We binned the sky above 300 MeV in 01 × 01 bins and calculated all counts n within the 08 ROI for each pixel. The minimum number of signal photons required to detect a pulsar with P(M) > 0.5 is given by solving Equation (9) for n_signal:

$$\begin{equation} n_{{\rm signal}} = \frac{10^M \times \sqrt{n}}{{\rm PF}}, \end{equation} \tag{ 22 }$$

where n_signal is converted to minimum flux F_min using Equation (13):

$$\begin{eqnarray} F_{{\rm min}} &=& \frac{1}{k} \frac{10^M \times \sqrt{n}}{{\rm PF}}\nonumber\\ &\times&\frac{\int _{0.1 \,{\rm GeV}}^{20 \,{\rm GeV}} E^{-\Gamma } \exp \left[ - \left(E/E_{{\rm cutoff}} \right)^\beta \right] dE}{\int _{0.3 \,{\rm GeV}}^{20 \,{\rm GeV}} R(r=0\mbox{$.\!\!^\circ $}8,\theta _{68}) E^{-\Gamma } \exp \left[ - \left(E/E_{{\rm cutoff}} \right)^\beta \right] dE},\nonumber\\ \end{eqnarray} \tag{ 23 }$$

where P(M) = 0.5 corresponds to M ∼ 1 for LAT locations (see Figure 10.) This map is shown in Figure 12, assuming P(M) ⩾ 0.5, PF = 0.5, Γ = 1.5, E_cutoff = 4 GeV, and β = 1.

Zoom In Zoom Out Reset image size

For counterpart sources, P(M) = 0.5 corresponds to M ∼ 0.944. At the threshold probability of P(M) = 0.5 the ratio of minimum signal for a counterpart source location to minimum signal for a LAT source location is 10^0.944/10 ≈ 0.88, which indicates the counterpart locations allow for up to a ∼12% dimmer flux detection threshold. This is especially useful for sources located within the Galactic plane as diffuse Galactic emission can obscure faint pulsars.