UNCERTAINTIES OF MODELING GAMMA-RAY PULSAR LIGHT CURVES USING VACUUM DIPOLE MAGNETIC FIELD

Pulsars are rotating neutron stars (NS) with very strong magnetic fields (B ∼ 10⁹–10¹² G). Nearly two thousand pulsars are known, mostly in the radio band. Compton Gamma Ray Observatory (CGRO) firmly detected seven gamma-ray pulsars and another three with less confidence (see Thompson 2004 for a review). This number is rapidly increasing with the start of operations of Fermi Gamma-ray Space Telescope (e.g., Abdo et al. 2008, 2009a, 2009b, 2010). Gamma-ray light curves of pulsars are typically double-peaked, generally out of phase with the radio, and have substantial off-peak emission. Several theoretical models have been developed to explain the nature of this emission, namely, the polar cap model (Ruderman & Sutherland 1975; Harding et al. 1978; Daugherty & Harding 1982, 1996), the slot-gap model (SG; or two-pole caustic, TPC in short; Arons & Scharlemann 1979; Arons 1983; Muslimov & Harding 2003, 2004; Dyks & Rudak 2003; Dyks et al. 2004), and the outer-gap (OG) model (Cheng et al. 1986a, 1986b; Romani & Yadigaroglu 1995; Yadigaroglu 1997; Cheng et al. 2000). In these models, particles are accelerated in the "gap" regions, where strong electric fields are developed due to a deficit of charge. Gamma-ray emission is interpreted as the curvature and/or inverse Compton radiation from ultra-relativistic particles accelerated in these gaps. The models differ in the location of the gaps in the magnetosphere. The polar cap model has narrow beam size and has difficulty in producing extended light curves, whereas the SG/TPC and OG models can reasonably well reproduce double-peaked, extended light curves. Recent works have shown that OG model can also reproduce the high-energy spectrum for Crab and Vela pulsars (Hirotani 2007; Takata & Chang 2007; Takata et al. 2007, 2008).

Special relativistic effects, such as the aberration of photon emission direction and photon travel time delay, are important for calculations of gamma-ray light curves when the emission zone extends far from the NS surface. In the TPC model, the emission zone is assumed to be along the last open field lines (LOFLs), extending from the polar cap to some cutoff radius. In the OG model, the emission zone is along the open field lines and extends from the null charge surface (NCS) to the light cylinder¹ (LC; R_LC = c/Ω). In both models, the relativistic effects are essential to forming "caustics" in the sky map, which appear as sharp peaks in the light curve. The caustics arise when photons emitted from different regions of the magnetosphere happen to arrive at the observer at the same time. The presence and appearance of caustics is sensitive to both the relativistic effects and the geometry of the emission zone.

All existing calculations of pulsar gamma-ray light curves assume that pulsar magnetospheric geometry can be represented by a vacuum dipole magnetic field. However, the pulsar magnetosphere is filled with plasma (Goldreich & Julian 1969), and plasma currents should modify the field structure. The magnetosphere should then consist of the open and closed field line regions, separated by thin current sheets. In contrast, all field lines of the vacuum field, including those that travel beyond the LC, are formally closed, and no current sheets exist. Numerical solutions of the structure of the plasma-filled magnetosphere are now known in the limit of force-free (FF) MHD for axisymmetric rotators (Contopoulos et al. 1999; Gruzinov 2005; Timokhin 2006; McKinney 2006), and, recently, for three-dimensional oblique rotators as well (Spitkovsky 2006; Kalapotharakos & Contopoulos 2009). The FF field clearly demonstrates the current sheet structure, and its geometry differs substantially from the vacuum field near the LC (see Spitkovsky (2008) for a review). Hence, the more realistic FF field geometry can lead to modifications of gamma-ray light curves (see Bai & Spitkovsky 2008 for preliminary results). On the other hand, the vacuum field has been used over the years to obtain light curves that compare very favorably to the existing data. It can be argued that if the emission comes from the regions in the magnetosphere that are not too close to the LC, the field geometry there may be well approximated by the vacuum field. This raises the question of how reliable and robust are the vacuum field light curves, and whether perturbations introduced by the presence of plasma in the magnetosphere would strongly affect the result.

Although our ultimate goal is to study gamma-ray emission using the FF field, in this paper we concentrate on the modeling of light curves using vacuum magnetic field only. The results with the FF field will be presented in the companion paper (Bai & Spitkovsky 2010). We feel it is necessary to clarify a number of points before moving forward. As we try to reproduce the sky maps and light curves using vacuum field geometry, we find that there are ambiguities in the literature on the use of the aberration formula. Since the aberration effect is crucial to the formation of caustics in any field geometry, in this paper we clarify the applicability of aberration formulae and compare their influences on the sky maps and light curves. In order to investigate the potential effects of plasma on the formation of light curves, we also study the sensitivity of the vacuum light curves to variations in the shape of the polar cap. Even if the field geometry in the bulk of the magnetosphere could be approximated by the vacuum field, it is the behavior of the field lines near the LC that determines the shape of the polar cap and thus the location of the magnetospheric emission zones. The plasma effects near the LC can then undermine light curve modeling that uses the polar caps of the vacuum field. In this paper, we show that the appearance of the sky map and light curves is indeed very sensitive to both the field geometry and the geometry of the emission zones, which suggests that revisiting theoretical models with the more realistic FF field is essential.

This paper is structured as follows. We begin with the vacuum magnetic field formulae for pulsar magnetosphere in Section 2, and then discuss the effect of aberration in Section 3. In Section 4, we construct the shape of the polar caps. We present comparisons of sky maps and light curves between different aberration formulae as well as different polar cap shapes for the TPC and the OG models in Section 5. In Section 6, we summarize our results.

We will use the vacuum field as an approximation to the magnetic field of the plasma-filled magnetosphere, where E · B = 0 and the magnetospheric structure is stationary in the corotating frame (CF). Therefore, in the lab frame (LF), for any specified magnetic field B, there exists an electric field E such that

$$\begin{equation} {\bf E}=-\frac{{\bOmega }\times {\bf r}}{c}\times {\bf B}\, \end{equation} \tag{ 1 }$$

where ${\bOmega }$ is the angular velocity of the NS, and r is the position vector. Equation (1) is of fundamental importance for discussing the aberration effect (Section 3), and is assumed throughout this paper. It does not apply, however, in the gaps where the high-energy emission is thought to originate. We address this issue in Appendix B and show that the conclusions of this paper remain unchanged even in the presence of gaps.

We will consider two commonly used formulae for the vacuum magnetic field, namely, the static dipole and the retarded dipole. For the static dipole, the field geometry is assumed to be the same as in a non-rotating dipole field, which is rigidly attached to the rotating pulsar. The field expression is

$$\begin{equation} {\bf B}=\frac{1}{r^3}[3(\vec{\mu }\cdot \hat{r})\hat{r} -{\mathbf \mu }]\,, \end{equation} \tag{ 2 }$$

where $\vec{\mu }$ is the magnetic dipole moment vector and $\hat{r}$ is the radial unit vector. For a pulsar rotating along the $\hat{z}$-axis with angular velocity Ω and magnetic inclination angle α, the time evolution of the magnetic moment vector is

$$\begin{equation} \vec{\mu }(t)=\mu (\sin \alpha \cos \Omega t\hat{x}+\sin \alpha \sin \Omega t\hat{y}+\cos \alpha \hat{z})\. \end{equation} \tag{ 3 }$$

The static dipole is not a full solution of the field of the rotating dipole in vacuum. As such, it can be interpreted as known either in the CF or the LF. As an approximation, we assume that the magnetic field of the static dipole is valid in the LF. The field geometry is sketched in Figure 1.

Zoom In Zoom Out Reset image size

The full solution of the electromagnetic field of a rotating magnetic dipole is known as the retarded dipole formula (e.g., Jackson 1975):

$$\begin{equation} {\bf B}=-\bigg [\frac{\vec{\mu }(t)}{r^3}+ \frac{\dot{\vec{\mu }}(t)}{cr^2}+\frac{\ddot{\vec{\mu }}(t)}{c^2r}\bigg ]+{\bf rr}\cdot \bigg [3\frac{\vec{\mu }(t)}{r^3}+3\frac{\dot{\vec{\mu }}(t)}{cr^2}+\frac{\ddot{\vec{\mu }}(t)}{c^2r}\bigg ]\. \end{equation} \tag{ 4 }$$

When r is small, the retarded dipole field configuration is almost the same as that of the static dipole. The deviation increases as the radius approaches the LC. One can find the expressions in Cartesian coordinates in Cheng et al. (2000). The retarded formula is valid in the LF, and not in the CF.

Another vacuum formula that can be used to represent the magnetosphere is the Deutsch field (Deutsch 1955; Michel & Li 1999), which includes corrections due to the finite size of the star R_N. These corrections are of second order in R_N/R_LC, and can be ignored for most pulsar parameters. Therefore, we will not consider the Deutsch field further.²

Pulsar gamma-ray radiation is believed to originate from curvature or inverse Compton emission from relativistic particles (electrons and positrons). In the extremely high electromagnetic field, these particles travel along magnetic field lines in a frame where the electric field vanishes.³ Emission direction from these highly relativistic particles, therefore, coincides with their direction of motion. In order to collect the pulsar emission to obtain a light curve, we need to account for relativistic effects including the aberration of photon emission direction⁴ and the correction of photon travel time (time delay).

The time delay effect is straightforward to calculate. The difference of light travel time for photons emitted from different regions of the magnetosphere results in the delay in rotational phase at which the photon is observed (Romani & Yadigaroglu 1995)

$$\begin{equation} \phi _d=-{\bf r}\cdot \hat{\eta }/R_{\rm LC},\, \end{equation} \tag{ 5 }$$

where r is the position of the emission point, $\hat{\eta }$ is the aberrated direction of photon emission (unit vector, in LF). The total phase of emission is ϕ = −ϕ_em + ϕ_d, where ϕ_em is the phase of $\hat{\eta }$. In order to determine $\hat{\eta }$ from any location in the magnetosphere, we need to consider the aberration effect. As the direction of particle motion depends on the direction of the magnetic field, one has to be careful about the frame in which the field is known. The vacuum field formulae are in general valid in the LF as discussed in Section 2. However, it was commonly implicitly assumed in the literature that these fields are valid in the instantaneous corotating frame (ICF), leading to discrepancies. For clarity, we will distinguish the following three frames.

• 1.  
  The lab frame, which is the inertial observer's frame in which the pulsar is rotating around the $\hat{z}$-axis.
• 2.  
  The corotating frame, which is a non-inertial frame that corotates with the pulsar and is the frame in which the field pattern is steady. It is related to the LF by the coordinate transformation

  $$\begin{eqnarray} x^{\prime } &=& x\cos \Omega t+y\sin \Omega t \nonumber\\ y^{\prime } &=& -x\sin \Omega t+y\cos \Omega t \\ z^{\prime } &=& z,\;\, t^{\prime }=t, \nonumber \end{eqnarray} \tag{ 6 }$$

  where the prime denotes the coordinate in the CF.
• 3.  
  The instantaneous corotating frame, which is a local frame of an inertial observer instantaneously moving at the corotation velocity. Therefore, it is defined only inside the LC.

We will discuss the aberration effect in these three frames separately. We pose the question as follows. The pulsar magnetic field configuration is known in the LF at t = 0 and the field pattern corotates with the star. The electric field in the LF obeys Equation (1). The photon emission direction at any position is along the magnetic field in a frame where the electric field vanishes. Our goal is to calculate the photon emission direction seen in the LF. We will do the calculation in three different frames and show that we reach the same result, as expected. From now on, we add superscripts "C" and "I" to denote fields in CF and ICF, respectively. Fields with no superscripts will always refer to fields in the LF.

First, we do the LF calculation. Consider the motion of an emitting particle in the pulsar magnetosphere. The FF condition requires ${\bf E}+\vec{\beta }_0\times {\bf B}=0$, where $\vec{\beta }_0$ is the normalized velocity of the particle relative to c, or the direction of the emitted photon in the LF. This equation, when combined with Equation (1), implies

$$\begin{equation} \vec{\beta }_0=f{\bf B}+\vec{\beta }_{\rm rot},\, \end{equation} \tag{ 7 }$$

where $\vec{\beta }_{\rm rot}={\bOmega }\times {\bf r}/c$ is the normalized corotation velocity, and f is a coefficient. For the emitting particle, we have $|\vec{\beta _0}|\rightarrow 1$, from which f can be determined. Solving Equation (7) with $|\vec{\beta _0}|=1$ fixes the emission direction $\hat{\eta }=\vec{\beta }_0$. We note that there are always two solutions associated with Equation (7), one associated with particles traveling along the magnetic field line in the CF, the other with particles traveling in the opposite direction.⁵

Second, we consider the calculation in the CF. The relation of the electromagnetic field between the CF and the LF is (Schiff 1939; Grøn 1984)

$$\begin{equation} {\bf B}^C={\bf B},\qquad {\bf E}^C={\bf E}+\vec{\beta }_{\rm rot}\times {\bf B}. \end{equation} \tag{ 8 }$$

According to Equation (1), we have E^C = 0. The photon emission in this frame is thus along B^C = B. We emphasize here that the CF and the LF are related by a coordinate transformation rather than Lorentz transformation. In Appendix A, we show that the aberration formula calculated from this frame has the same expression as Equation (7).

Finally, consider aberration in the ICF. The magnetic field in the LF and the ICF are related by a Lorentz transformation. Using Equation (1) for the electric field in the LF, we obtain

$$\begin{equation} B_t^{I}=B_t,\quad {\bf B}_p^{I}={\bf B}_p/\gamma,\quad {\bf E}^{I}=0,\, \end{equation} \tag{ 9 }$$

where γ ≡ (1 − β²_rot)^−1/2 is the Lorentz factor of corotation, B_t and B_p denote the toroidal and poloidal components of the magnetic field. Note that in ICF the poloidal magnetic field is smaller. Since E^I = 0, the photon emission direction is thus

$$\begin{equation} \eta ^{I}_t=\pm B_t/B^{\prime }_0,\qquad \vec{\eta }^{I}_p=\pm {\bf B}_p/\gamma B^{\prime }_0,\, \end{equation} \tag{ 10 }$$

where $B^{\prime }_0\equiv \sqrt{B_t^2+(1-\beta _{\rm rot}^2)B_p^2}$ is the total magnetic field strength in the ICF, and plus/minus signs correspond to emission along/opposite to the magnetic field line. To get back to the LF, one should perform a Lorentz transformation to aberrate the direction of the photon $\hat{\eta }^{I}\rightarrow \hat{\eta }$; the formula is (Dyks & Rudak 2003)

$$\begin{equation} \hat{\eta }=\frac{\hat{\eta }^{I}+[\gamma +(\gamma -1)(\vec{\beta }_{\rm rot}\cdot \hat{\eta }^{I}) /\beta _{\rm rot}^2]\vec{\beta }_{\rm rot}}{\gamma (1+\vec{\beta }_{\rm rot}\cdot \hat{\eta }^{I})}\. \end{equation} \tag{ 11 }$$

Substituting Equation (10) into Equation (11) we find

$$\begin{equation} \eta _t=\frac{\beta _{\rm rot} B_p^2\pm B_tB^{\prime }_0}{B_0^2}\, \quad \vec{\eta }_p=\frac{\pm B^{\prime }_0-\beta _{\rm rot} B_t}{B_0^2}{\bf B}_p\, \end{equation} \tag{ 12 }$$

where $B_0\equiv \sqrt{B_t^2+B_p^2}$ is the total magnetic field strength in the LF.

We note that although the ICF is defined only within the LC, Equation (12) is valid anywhere. In fact, the ICF is only one special instantaneous frame in which E^I vanishes. We can construct other instantaneous frames where the electric field is zero, and Equation (12) is a general result. For example, we can also choose the E × B instantaneous drift frame (IDF). One can show that by transforming the LF field to IDF to calculate the aberration leads to the same result as Equation (12). Simple algebra shows that Equations (12) and (7) are exactly equivalent as expected. Therefore, our proof is complete.

In many of the previous studies (e.g., Romani & Yadigaroglu 1995; Yadigaroglu 1997; Cheng et al. 2000; Dyks & Rudak 2003; Dyks et al. 2004), it was implicitly assumed that the vacuum field, given by any formula in Section 2, is valid in the ICF instead of the LF. The following equation was used instead of Equation (10):

$$\begin{equation} \eta ^I_t=\pm B_t/B_0, \quad \vec{\eta }^I_p=\pm {\bf B}_p/B_0\. \end{equation} \tag{ 13 }$$

The corresponding photon emission direction is then

$$\begin{equation} \eta _t^\times =\frac{\beta _{\rm rot} B_0\pm B_t}{B_0\pm \beta _{\rm rot} B_t}, \quad \vec{\eta }_p^\times =\frac{{\bf B}_p}{\gamma (B_0\pm \beta _{\rm rot} B_t)}\. \end{equation} \tag{ 14 }$$

In the limit β → 0, Equations (12) and (14) agree. The differences between them are of the order ${\cal O}(\beta _{\rm rot}^2)$ when β_rot ≪ 1. Since the difference between the two aberration formulae scales as β²_rot, then to first order in r/R_LC treating the vacuum dipole field as known in the ICF instead of the LF, as was done in Dyks et al. (2004), is a reasonable approximation. However, as we shall see in Section 5, this introduces significant differences in the calculated light curves. In the other limit, β_rot → 1, the two treatments of aberration result in different photon emission directions. Interestingly, the backward solution of Equation (12) agrees with the forward solution of (14), where the emission direction is perpendicular to the rotation axis.

Since the vacuum field only approximates the FF field, one could choose to treat it as known in the LF, CF or ICF, as a matter of approximation. However, as shown in the next section, the LOFLs must be traced either in the LF or in the CF, but not in the ICF. This makes the treatment in a number of previous works not self-consistent. We will show in Section 5 that the two approaches can lead to substantial differences in the sky maps and thus the light curves.

Determination of the LOFLs that separate the open and closed regions in the magnetosphere is needed to find the shape of the polar cap and the location of the emission zones. LOFLs are typically found by tracing field lines and checking whether they close inside or outside the LC. This tracing must be done either in the LF or in the CF, but not in the ICF. From Equation (9), the toroidal field in the ICF coincides with the LF or the CF field, but the poloidal component is γ times smaller. The ICF field is thus ill-defined at the LC. As R → R_LC, B^I_p → 0, the field is purely toroidal, and no field line actually crosses the LC. Also, since the ICF is not a global frame, local tracing of the field does not guarantee the tracing of a consistent set of field lines as seen from other frames due to aberration. Thus, the shape of the polar cap in the ICF or in any other instantaneous frame is not well defined.⁶

To find LOFLs, we use third-order Runge–Kutta integration to trace the magnetic field lines in the LF. Bisection method (Dyks et al. 2004) is used to find the magnetic colatitude θ^rim_m of the rim of the polar cap for every magnetic azimuth ϕ_m. We show the resulting polar cap shape for both the static and retarded dipole fields in Figure 2 for magnetic inclination angle α = 60°. The polar cap shape of the retarded dipole as seen in Figure 2 is commonly found in the literature (e.g., Yadigaroglu 1997; Cheng et al. 2000; Dyks et al. 2004). Therefore, in these works the retarded dipole field was traced in the LF, but for the computation of aberration the field was treated in the ICF, leading to an inconsistency.

Zoom In Zoom Out Reset image size

The vacuum field geometry appears similar to the FF field geometry near the NS surface, but the two are substantially different near the LC (Spitkovsky 2006). The polar cap shape, as determined by LOFLs, is thus sensitive to the field structure near the LC. This may cause significant uncertainties in the geometry of the emission zones. As we show in the companion paper (Bai & Spitkovsky 2010), the polar cap shape in the FF field is, in fact, more circular and also larger than in the vacuum field. Thus, a circular polar cap may better characterize the emission zone geometry when vacuum field is used. Dyks & Rudak (2003) considered a static dipole field with a circular polar cap, given by $\theta _m=\theta _0 =\arcsin {\sqrt{R_{\rm NS}/R_{\rm LC}}}$. In this paper, we also consider such a circular polar cap which is indicated as a reference circle in Figure 2. We will refer to the polar cap found by tracing LOFLs as "traced polar cap," to distinguish it from the circular polar cap.

In order to parameterize different field lines, we define open volume coordinates on the polar cap. The magnetic colatitude of the polar cap rim θ^rim_m(ϕ_m) is generically a function of azimuth. For any point on the NS surface at the magnetic colatitude θ_m and azimuth ϕ_m, we define the open volume coordinate of this point to be (r_ov, ϕ_m), where r_ov ≡ θ_m/θ^rim_m(ϕ_m). Therefore, the rim of the polar cap and LOFLs correspond to r_ov = 1. All open field lines have r_ov < 1 while all closed lines have r_ov > 1. We emphasize that open volume coordinates should be defined for magnetic fields in the LF or CF.⁷

In this section, we construct the sky maps and calculate the associated light curves for pulsar gamma-ray emission in vacuum fields. We consider both the TPC and the OG models for static and retarded dipole fields. In each case, we compare the results using different aberration formulae (i.e., Equations (12) and (14)). We also compare the results using different polar cap shapes.

The TPC model is an extended version of the SG model, where emission has until recently been assumed to come from a thin sheet centered on the LOFLs, i.e., r⁰_ov = 1.⁸ The emission zone extends from above the polar cap along the LOFLs up to a certain cutoff radius, where emissivity is assumed to drop to zero. The cutoff is described by r_max, the distance to the center of the NS, and R_max, the cylindrical radius to the rotational axis. These two parameters constrain the extent of the emission zone. For the OG model, the emission zone is assumed to come from a layer in the open field line region beyond the NCS where Goldreich–Julian charge density equals to zero⁹ [$\rho _{\rm GJ}\simeq {\bOmega } \cdot {\bf B}/(2\pi c)$]. We adopt r⁰_ov = 0.9 for the center of this layer. The emission can extend to R_max ∼ R_LC, and no extra cutoff is needed. In all the cases, emissivity is assumed to be constant along the field lines in the emission zone. For different field lines, emissivity is weighed by a Gaussian function centered at r⁰_ov, with width σ = 0.025 (Dyks et al. 2004).

The radiation from the emission zones, when projected to the observer as the NS rotates, produces the light curve. It is convenient to collect all photons from the emission zones as they fall on the sky and plot the intensity as a function of stellar phase ϕ and the observer's viewing angle ξ_obs. The light curve is then obtained by a cut through this sky map at the observer's viewing angle ξ_obs.

Mathematically, the sky map is a map from the emission zone to the sky coordinates. At fixed r_ov, the emission zone is a two-dimensional manifold, which can be parameterized by (ϕ_m, l), where l is the length of the magnetic field line starting from the polar cap region, and ∂/∂l is along the direction of the magnetic field in units of R_LC. With constant emissivity along field lines, the Jacobian determinant of this transformation determines the intensity on the sky map:

$$\begin{equation} I(\phi, \xi _{\rm obs})\propto \bigg |\det {\frac{\partial {(\phi, \xi _{\rm obs})}}{\partial {(\phi _m, l)}}}\bigg |\. \end{equation} \tag{ 15 }$$

When the inverse of the Jacobian matrix is singular, infinite magnification is reached, and I → ∞. This means that light rays emitted from the neighborhood around the position (ϕ_m, l) arrive at the observer at the same time, which greatly strengthens the intensity. This is the analog of caustics in optics and strong gravitational lensing. In the case of the sky map, a "caustic" is often understood as regions on the sky map with strong enhancement, where the determinant in Equation (15) is much greater than 1. The caustics on the sky map then correspond to peaks in the light curve. There is no guarantee of the existence of caustics on the sky map mathematically, and a general field structure or a variation in emission region geometry does not necessarily produce caustics. In the following subsections, we will show examples of how the properties of caustics change as we use different field geometries, shapes of the polar caps and aberration formulae.

5.1. Two-pole Caustic Model with Static Dipole Field

The TPC model with a static dipole field was considered by Dyks & Rudak (2003). They used a circular polar cap, with $\theta _m=\arcsin {\sqrt{R_N/R_{\rm LC}}}$, rather than a traced polar cap. In the calculation of the sky map, they traced the field lines using the vacuum dipole formula (2). Therefore, the field was treated as in the LF, but Equation (14) was used for aberration, making the treatment not self-consistent. In this subsection, we reproduce this result and compare with other sky maps obtained from self-consistent treatment using the static dipole field.

In Figure 3, we show four sky maps and light curves using the static dipole field with inclination angle α = 70°. In (a) and (b), we use the traced polar cap, while in (c) and (d) we use the circular polar cap. Aberration formula (14) is used in panels (a) and (c), which is not self-consistent. Other two panels have self-consistent aberration (using formula (12)). The case considered by Dyks & Rudak (2003) corresponds to Figure 3(c). Representative light curves are plotted on the right of each panel in Figure 3; the observer's viewing angle is chosen to be 80°.

Zoom In Zoom Out Reset image size

In all panels of Figure 3, bright arc-like regions are caustics, which translate to sharp peaks in the light curve. As far as the bulk appearance of the sky maps is concerned, Figures 3(a) and (b) look similar, while Figures 3(c) and (d) look similar. On the other hand, in each of the corresponding pairs, the details of the caustics differ. Therefore, the overall appearance of the sky map is very sensitive to the shape of the polar cap, while the choice of aberration formula is important for determining the detailed features of the light curves. To some extent, the inconsistent treatment of magnetic field in the aberration formula in a number of earlier works may not be critical, but at the current level of very precise light curve measurement with the Fermi telescope, detailed differences in the light curves can be distinguishable. We now discuss the sky maps and the corresponding light curves for the four cases.

We begin with Figure 3(c), which is a reproduction of Figure 4 in Dyks & Rudak (2003). In the sky map, we have two main caustics near the polar caps, which are very prominent and sharp. Their ξ_obs extends from 50° to 110° in the sky map. Right behind the main peak, a small hump in the sky map leads to a "subdominant" peak (labeled with "S" in the light curves). As a result, two sharp peaks and a small, subdominant peak right behind the first peak is present in the light curve. This light curve was used to explain Vela's gamma-ray light curve by Dyks & Rudak (2003).

While Figure 3(d) looks similar to Figure 3(c), a change in the aberration formula causes some quantitative differences. The "subdominant" peak becomes stronger, and the extent of the main caustics becomes smaller. As a result, it is less probable to have a double-peaked light curve and more difficult to reproduce Vela's profile.

Figures 3(a) and (b), where LOFLs were traced to get the polar cap shape, look quite different from Figures 3(c) and (d), which used a circular polar cap. Compared to Figure 3(c), the shape of the main caustic is more curved, and the original subdominant peaks become strong and extended caustics. Also, the emission zones resulting in two peaks now cover a smaller fraction of the sky map. Now it is almost impossible to reproduce Vela's light curve with this inclination angle.

In sum, the double-peak feature for gamma-ray pulsars is hard to reproduce with the static dipole field, with the exception of case (c), where the inconsistent treatment of magnetic field in the aberration formula is used. This situation is true within a wide range of interesting inclination angles α (adjusting r_max and R_max, or changing r⁰_ov to smaller value does not help either). This example demonstrates how the sky maps and the resulting light curves are affected by the aberration formula and the shape of the polar caps. Meanwhile, the differences between panels (b) and (d) in Figure 3 characterize the uncertainties of sky maps and light curves calculated using the static dipole field.

5.2. Two-pole Caustic Model with Retarded Dipole Field

Next, we consider the TPC model with a retarded dipole field (e.g., Dyks et al. 2004). We choose the inclination angle α = 60° and show the sky maps and light curves in Figure 4. As in the last subsection, we also consider four cases. In all of them, the retarded dipole formula is considered valid in the LF, and is also traced in the LF. In (a) and (b), we use the traced polar cap, and the circular polar cap is used for the other two. In (a) and (c), the aberration formula assumes that the field is in the ICF (which is not self-consistent), while for the other two figures, the corrected aberration formula (Equation (7)) is used. We have chosen ξ_obs = 80° for all light curves in Figure 4.

Zoom In Zoom Out Reset image size

Figure 4(a) is a reproduction of the sky map and light curve of Dyks et al. (2004). There are two strong caustics formed near the phase of two poles that extend over a wide range of observer's viewing angles. Just behind the first peak there is a small hump, caused by the overlap of emission from two poles. A small step-like drop in the middle corresponds to a weak discontinuity in the polar cap rim (Dyks et al. 2004).

Figure 4(b) is the corrected version of Figure 4(a), where aberration is treated self-consistently. Remarkably, the two strong caustics become blurred into weak and wide enhancements at the same phases. Since the caustics are weak, the brightest region in the sky map is mainly caused by the overlap of emission from both poles. In the light curve, we thus have two wide peaks, and a large fraction of radiation is from the off-peak phases. We have also explored the sky maps with other inclination angles, and find that the strong caustics appear in fewer combinations of inclination angles and viewing geometries (see Appendix C). The weak and wide caustic structure is generic for this TPC model using the retarded dipole field, and it has difficulty in producing sharp peaks in the light curve.¹⁰

Figures 4(c) and (d) have circular polar caps. The appearance of the sky map is quite different from Figures 4(a) and (b). The caustics are more extended, and it is possible to get two peaks from the same pole (in total up to four peaks). Meanwhile, we find the same trend: the inconsistent treatment of B field in the aberration leads to strong and narrow caustics; when aberration is corrected, the caustics become wide and weak.

We see that the shape of the polar cap largely determines the overall appearance of the sky map, while the details of the caustics are sensitive to the aberration formula. The geometry of the emission zone is mainly controlled by the shape of the polar cap. This example confirms that the shape of the polar cap (or the emission zone geometry) is one major source of uncertainty in the current modeling of gamma-ray pulsar light curves using TPC model (e.g., Dyks et al. 2004).

We can also directly compare the effect of different field geometries by using the same circular polar cap shape. Since we have previously chosen different inclination angles in Figure 3 (static dipole, α = 70°) and Figure 4 (retarded dipole, α = 60°) in order to reproduce previous works by Dyks & Rudak (2003) and Dyks et al. (2004), in Figure 5 we plot the sky map for the static dipole field with α = 60° and a circular polar cap to facilitate the comparison with Figure 4. We consider both treatments of aberration in Figure 5. One can compare Figures 4(c) and (d) with Figures 5(a) and (b). The only difference between the two pairs of the sky maps is different field geometry. Since the retarded and static dipole fields are similar near the star, one might expect that the resulting sky maps would be similar as well. However, we see substantial differences in the appearance of the sky maps and the light curves, and the differences are already very prominent in regions that are not far from the star, as seen from the structure of the caustics. In particular, the main caustics form closer to the star and appear stronger in the case of the static dipole field than in the retarded dipole. Even though the caustics form at low altitude, they are very sensitive to the choice of field structure because caustics are a chance overlap of emission from different regions. This analysis shows that the appearance of the sky map is very sensitive to the field geometry itself, and the deviation of the vacuum field geometry from the more realistic FF field should be another major source of uncertainty in the modeling of gamma-ray pulsar light curves.

Zoom In Zoom Out Reset image size

5.3. Outer Gap Model with Retarded Dipole Field

The OG model with a retarded dipole field was studied by many authors (e.g., Romani & Yadigaroglu 1995; Yadigaroglu 1997; Cheng et al. 2000; Dyks et al. 2004; Takata & Chang 2007; Takata et al. 2007, 2008). As before, we consider four cases in Figure 5: panels (a) and (b) with a traced polar cap, and (c) and (d) with a circular polar cap. Consistent treatment of aberration is done in (b) and (d) only. Case (a) was explored by Romani & Yadigaroglu (1995); Yadigaroglu (1997); Cheng et al. (2000); and Dyks et al. (2004), and case (b) was recently explored by Takata & Chang (2007), and Takata et al. (2007, 2008), who have corrected the aberration formula. These studies were able to reproduce some of the gamma-ray pulsar light curves. Here, we focus on the comparison between the four cases to address the uncertainties of the OG model.

Figure 6 shows the sky maps for the four cases with inclination angle α = 65°. We further show the light curves at ξ_obs = 80° for all four cases to the right of the sky maps. Our Figure 6(a) is a reproduction of Romani & Yadigaroglu (1995); Yadigaroglu (1997); Cheng et al. (2000); Dyks et al. (2004). The two peaks are caused by emission from the same pole. The offset from phase zero and the double-peak profile are clearly present. The arrow indicates the effect of the "notch" in the polar cap, and the emission from both sides of the notch contributes to the caustic. In between the double peaks is the "bridge" emission with weaker intensity. In the plot with corrected aberration (Figure 6(b)), the regions on the sky map responsible for the first caustic become less extended. The emission from field lines on one side of the notch no longer forms caustics. Also, the regions responsible for the "bridge" become weaker. Figure 6(c) and (d) are the sky maps for the OG model with circular polar caps. Although the sky maps look different from cases (a) and (b), the light curves look very similar to their counterparts ((c) resembles (a), and (d) resembles (b)).

Zoom In Zoom Out Reset image size

We note that when we correct the aberration formula for the OG model, the sharpness of the two peaks is less affected than in the TPC model. This appears odd since the emission zone of OG is close to the LC, and one may expect a larger difference between the two treatments of aberration. As one can see from Equation (15), the intensity of the caustics is a differential effect. It depends on the separation (on the sky map) between the projection of two neighboring points in the emission zone. Therefore, it is possible that caustic structures are similar for the two treatments of aberration in the OG model.

Compared with the TPC model, the appearance of the sky maps and lights curves in the OG model is less sensitive to the shape of the polar cap. This is not surprising, because the open field lines not too close to the LOFL are not sensitive to small displacements in the polar cap. Therefore, the sky map does not change much when we replace the traced polar cap with a circular one. This might suggest that the OG model using the vacuum field is more robust than the SG model. However, the magnetic field used by the OG model is needed close to the LC, which itself is still very uncertain. As we have discussed in the last paragraph of Section 5.2, the uncertainties in the field geometry itself can also cause big differences in the sky maps. Therefore, the accuracy of the OG model using vacuum field is still in question.

All current models of pulsar gamma-ray light curves assume that the pulsar magnetosphere can be represented by the vacuum magnetic dipole field. The use of the vacuum field should be considered as an approximation to the field which includes the effects of plasma in the magnetosphere. In this paper, we considered the application of the vacuum field to theoretical models (the TPC and the OG models). Our results show that there are large uncertainties in using the vacuum field for predictions of gamma-ray emission.

Calculations of high-energy pulsar light curves involve the construction of a map from the emission zone in the magnetosphere to the sky map. The appearance of the sky map is thus very sensitive to (1) the geometry of the magnetic field and (2) the geometry of the emission zone. On the one hand, the vacuum field deviates from the FF field, which introduces one major uncertainty in the direction of emission. On the other hand, the shape of the polar cap is traced by the LOFLs, and sensitively depends on the field structure near the LC (Dyks & Harding 2004). This causes another large uncertainty because the geometry of the emission zone is largely determined by the shape of the polar cap.

Relativistic effects including aberration and time delay are crucial to the formation of caustics in the sky map. We provide the treatment of the aberration effect in three different reference frames, namely, the LF, the CF and the ICF, and show that they can be reconciled. To be consistent, however, the tracing of the field lines to find the polar cap must be done in the LF or the CF, but not in the ICF, because the latter is not a global frame.

We compare the sky maps and light curves for the TPC model and the OG model using different aberration formulae. We find that the appearance of the caustics in the sky maps is sensitive to the treatment of B field in the aberration formula. For the TPC model using a retarded dipole field, we find that instead of having two strong caustics in the sky map, the corrected aberration formula weakens the caustics, leaving two wide and weak humps in the sky map. As a result, the conventional TPC model with retarded dipole field has difficulty in producing sharp peaks in the light curve. For the OG model using a retarded dipole field and the corrected aberration formula, the caustic responsible for the first peak in the light curve is less extended, and there is weaker "bridge" emission between the two peaks. Recent developments in the OG model with the corrected aberration formula are still able to produce reasonable light curves (Takata & Chang 2007; Takata et al. 2007).

We study the uncertainties in the models of pulsar gamma-ray light curves by (1) comparing the sky maps and light curves using different shapes of the polar cap, namely, the polar cap obtained by tracing LOFLs and the circular polar cap; (2) comparing the sky maps and light curves using the same polar cap shape (circular), but with different magnetic field configurations (static versus retarded dipole). We choose a circular polar cap as another possibility because the polar cap of the FF magnetosphere is more circular (Bai & Spitkovsky 2010). Our results show that the overall appearance of the sky map is very sensitive to both factors. We find that for the TPC model, whose emission zone is centered on the LOFLs, the sky maps and light curves depend sensitively on the changes in the shape of the polar cap. Up to four peaks can be present when a circular polar cap is used. This suggests that the reproduction of the Vela's light curve in Dyks & Rudak (2003), and Dyks et al. (2004) is not robust. For the OG model, whose emission zone is in the open field line region, we find that the sky maps and the light curves are not very sensitive to changes in the polar cap shapes. However, the vacuum field near the LC is unreliable, and the predictions from the OG model are still questionable. A detailed atlas of gamma-ray light curves for two models is presented in Appendix C.

In all, we conclude that it is essential to revisit the existing theoretical models using a more realistic magnetospheric structure, i.e, the FF field from numerical simulations (Spitkovsky 2006). In the companion paper Bai & Spitkovsky (2010), we will present the sky maps and light curves calculated using the FF field.

We thank Yury Lyubarsky and Jonathan Arons for help and advice. We also thank our referee for helpful suggestions. This work is supported by NASA grants NNX08AW57G and NNX09AT95G. A.S. acknowledges support from Alfred P. Sloan Foundation Fellowship. X.N.B. acknowledges support from NASA Earth and Space Science Fellowship.

We can find the inverse transformation of Equation (6) and write it in the differential form as

$$\begin{eqnarray} &&\left(\begin{array}{@{}l@{}} dx\\ dy\\ dz\\ dt\end{array}\right) = \left(\begin{array}{@{}lccc@{}} \cos \Omega t^{\prime } & -\sin \Omega t^{\prime } & 0 & -\Omega y \\ \sin \Omega t^{\prime } & \cos \Omega t^{\prime } & 0 & \Omega x \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{array}\right) \left(\begin{array}{@{}l@{}}dx^{\prime }\\ dy^{\prime }\\ dz^{\prime }\\ dt^{\prime }\end{array}\right).\qquad \end{eqnarray} \tag{ A1 }$$

It suffices to consider t = t' = 0, in which case we also have x = x', y = y'. Dividing each side by cdt = cdt', we obtain the photon velocity vector

$$\begin{equation} \eta _x=\eta ^{\prime }_x-\Omega y/c,\qquad \eta _y=\eta ^{\prime }_y+\Omega x/c\,\qquad \eta _z=\eta ^{\prime }_z. \end{equation} \tag{ A2 }$$

In this equation, $\vec{\eta }$ is a unit vector that denotes the photon direction in the LF, while $\vec{\eta }^{\prime }$ denotes the photon motion in the CF, but it is NOT a unit vector. This is because the metric in the CF is not the Minkowski metric, and has non-diagonal spacetime components (Schiff 1939).

In the case of pulsar gamma-ray emission, $\vec{\eta }^{\prime }$ is along the direction of B^C = B. Therefore, Equation (A2) can be rewritten as

$$\begin{equation} \vec{\eta }=f{\bf B}+\Omega \times {\bf r}/c. \end{equation} \tag{ A3 }$$

This is exactly the same as Equation (7).

In this Appendix, we relax the FF assumption to allow the presence of gaps in the magnetosphere. In the gaps, E_∥ ≠ 0, and Equation (1) no longer holds. However, we show that under the following two assumptions, the aberration formula discussed in this paper still holds in the gaps.

• 1.  
  The magnetospheric structure is stationary in the CF.
• 2.  
  Particles move along magnetic field lines in the CF.

Assumption (1) allows us to work in the CF using Equation (6). Assumption (2) implies that in the CF, the perpendicular electric field must be zero (while E^CF_∥ can be non-zero in the gaps). Transforming the field back to the LF by applying Equation (6), we obtain

$$\begin{eqnarray} {\mathbf E}_{\perp }&=-\frac{{\bOmega }\times {\mathbf r}}{c}\times {\mathbf B},\nonumber\\ {\mathbf E}_{\parallel }&={\mathbf E}_{\parallel }^{CF}. \end{eqnarray} \tag{ B1 }$$

Physically, assumption (2) means the gap in the magnetosphere also corotates with the star, because the E × B drift velocity has a corotation component everywhere.

Under the above two assumptions, even if particles are accelerated due to E_∥ in the gaps, their direction of motion, and, hence, the direction of photons they emit, is always along the B field in the CF. Therefore, the treatment of aberration in this case is exactly the same as in Appendix A, and the formula (Equation (12)) is valid even in the presence of gaps. This analysis can be reproduced in the LF as well.

The above discussions also apply if the magnetosphere rotates differentially with respect to the NS (e.g., Timokhin 2007a, 2007b), as long as one replaces Ω with the value of the angular velocity in the magnetosphere.

Knowing how the light curves appear for different viewing angles and magnetic obliquities can help place constraints on the pulsar geometry. Recently, Watters et al. (2009) compiled an atlas of TPC and OG light curves for a range of pulsar parameters using retarded vacuum field. The field was effectively assumed to be in the ICF for the purposes of aberration calculation. Here, we provide a similar atlas, emphasizing the differences brought by the self-consistent treatment of aberration. In Figures 7 and 8, we show the sky maps and representative light curves for TPC and OG models, respectively. For every pulsar inclination angle (horizontal rows), we show the results for field taken in the LF (subpanels a, c, e, g) and in the CF (subpanels b, d, f, h). As was shown above, only the LF treatment for the retarded vacuum field is self-consistent. The polar cap is obtained by tracing LOFLs in all cases. The light curve plots have identical vertical scales at the same obliquity for easy cross-comparison.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The atlas displays the same general trends as seen before: self-consistent aberration reduces the strength of the caustics in the TPC model, and decreases the amplitude of the first peak for many observer angles in the OG model.