CIRCULAR POLARIZATION IN PULSARS DUE TO CURVATURE RADIATION

Pulsars are highly magnetized with predominantly dipolar field structure. The rotating magnetic field produces a strong induced electric field that accelerates charged particles off the surface of the star into a magnetosphere consisting of predominantly dipolar magnetic field and corotating relativistic pair plasma. Pulsar radio emission models assume that radiation emitted tangentially to the field lines on which plasma is moving. The polarization state of the emitted radiation is more or less determined by the structure of magnetic field at the emission spot. In the general framework of models in which the radio power is curvature radiation emitted by charge bunches constrained to follow field lines, the linear polarization is intrinsic to the emission mechanism and is, furthermore, a purely geometric property. Several pulsar researchers have shown that the properties such as the polarization angle swing can be explained within the framework of curvature radiation (e.g., Radhakrishnan & Cooke 1969; Sturrock 1971; Ruderman & Sutherland 1975; Lyne & Manchester 1988; Rankin 1990, 1993; Blaskiewicz et al. 1991).

The radio emission from particle bunches is highly polarized, and the radiation received by a distant observer will be less polarized due to the incoherent superposition of emissions from different magnetic field lines (Gil & Rudnicki 1985). Gil (1986) has argued for the connection between pulsar emission beams and polarization modes and suggested that out of two orthogonal polarization modes, one corresponds to core emission and the other to the conal emissions. They are highly linearly polarized and the observed depolarization is due to superposition of modes at any instant (Gil 1987). By considering a charged particle moving along the curved trajectory (circular) confined to the xz-plane, Gil & Snakowski (1990a) have deduced the polarization state of the emitted radiation and shown the creation of antisymmetric circular polarization in curvature radiation. By introducing a phase, as a propagation effect, the difference between the components of radiation electric field in the directions parallel and perpendicular to the plane of particle trajectory, Gil & Snakowski (1990b) have developed a model to explain the depolarization and polarization angle deviations in subpulses and micropulses. Gil et al. (1993) have modeled the single-pulse polarization characteristics of pulsar radiation and demonstrated that the deviations of the single-pulse position angle from the average are caused by both propagation and geometrical effects. Mitra et al. (2009), by analyzing the strong single pulses with highly polarized subpulses from a set of pulsars, have given very conclusive arguments in favor of the coherent curvature radiation mechanism as the pulsar radio emission mechanism.

By analyzing the average pulse profiles, Radhakrishnan & Rankin (1990) have identified two most probable types of circular polarizations, namely, antisymmetric, where the circular polarization changes sense near the core region, and symmetric, where the circular polarization remains with same sense. They found that antisymmetric circular polarization is correlated with the polarization angle swing and speculate it to be a geometric property of the emission mechanism. Han et al. (1998), by considering the published mean profiles, found a correlation between the sense of circular polarization and polarization angle swing in conal double profiles and no significant correlation for core components. Further, You & Han (2006) have reconfirmed these investigations with larger data. However, Cordes et al. (1978) were the first to point out an association between the position angle of the linear polarization and the handedness of the circular polarization.

There are two types of claims for the origin of circular polarization: intrinsic to the emission mechanism (e.g., Michel 1987; Gil & Snakowski 1990a, 1990b; Radhakrishnan & Rankin 1990; Gangadhara 1997) or generated by the propagation effects (e.g., Cheng & Ruderman 1979). Cheng & Ruderman (1979) have suggested that the expected asymmetry between the positively and negatively charged components of the magnetoactive plasma in the far magnetosphere of pulsars will convert linear polarization to circular polarization. Radhakrishnan & Rankin (1990) have suggested that the propagation origin of antisymmetric circular polarization is very unlikely but the symmetric circular polarization appears to be possible. On the other hand, Kazbegi et al. (1991, 1992) have argued that the cyclotron instability, rather than the propagation effect, is responsible for the circular polarization of pulsars. Lyubarskii & Petrova (1999) considered that the rotation of the magnetosphere gives rise to wave mode coupling in the polarization-limiting region, which can result in circular polarization in linearly polarized normal waves. Melrose & Luo (2004) discussed possible circular polarization induced by intrinsically relativistic effects of pulsar plasma. Melrose (2003) reviewed the properties of intrinsic circular polarization and circular polarization due to cyclotron absorption and presented a plausible explanation of circular polarization in terms of propagation effects in an inhomogeneous birefringent plasma. In the multifrequency simultaneous observations, we do find the variations in the single-pulse polarization, which may be attributed to the propagation effects (Karastergiou et al. 2001, 2002, 2003).

The correlation between the antisymmetric circular polarization and the polarization angle swing is a geometric property of the emission processes (Radhakrishnan & Rankin 1990). By carefully modeling the polarization state of the radiation in terms of Stokes parameters, it is possible to construct the geometry of emission region at multifrequencies. So far, in the purview of curvature radiation only the polarization angle has been modeled (Radhakrishnan & Cooke 1969; Komesaroff 1970) and attempted to be fit with the average radio profile data (e.g., Lyne & Manchester 1988). Instead of circular trajectories, it is very important to consider the actual dipolar magnetic field lines, whose curvature radii vary as a function of altitude, as the radio emission in pulsars is expected to come from the range of altitude (e.g., Gangadhara & Gupta 2001; Gupta & Gangadhara 2003; Krzeszowski et al. 2009). In this paper, we develop a three-dimensional (3D) model for curvature radiation by relativistic sources accelerated along the dipolar magnetic field lines. We consider the actual dipolar magnetic field lines (not the circles) in a slowly rotating (non-rotating) magnetosphere such that the rotation effects can be ignored. The relativistic plasma (bunch, i.e., a point-like huge charge) moving along the dipolar magnetic field lines emits curvature radiation. We show that our model reproduces the polarization angle swing of Radhakrishnan & Cooke (1969), and predicts that the correlation of antisymmetric circular polarization and polarization angle swing is a geometric property of the emission process. Our model is aimed at re-examining the intrinsic polarization properties of the vacuum single-particle curvature radiation, and planned to consider the propagation effects separately in the subsequent works. We derive electric fields of the radiation field in Section 2 and construct the Stokes parameters of the radiation field in Section 3. A few model (simulated) profiles are presented in Section 4 depicting the correlation between the antisymmetric circular polarization and polarization angle swing in the different cases of viewing geometry parameters.

Consider a magnetosphere having dipole magnetic field with an axis $\hat{m}$ inclined by an angle α with respect to the rotation axis $\hat{\Omega }$ (see Figure 1). We assume that the magnetosphere is stationary or slowly rotating such that the rotation effects are negligible. The relativistic pair plasma, generated by the induced electric field followed by pair creation, is constrained to move along the curved dipolar magnetic field lines. The high brightness temperature of the pulsar indicates coherency of the pulsar radiation, which in turn forces one to postulate the existence of charged bunches. The formation of bunches in the form of solitons has been proposed (e.g., Cheng & Ruderman 1979; Melikidze & Patarya 1980, 1984) and questioned (e.g., Melrose 1992). Gil et al. (2004) have generalized the soliton model by including formation and propagation of the coherent radiation in the magnetospheric plasma along magnetic field lines. Their results strongly support coherent curvature radiation by the spark-associated solitons as a plausible mechanism of pulsar radio emission. Following these views, we assume that the plasma in the form of bunches moves along the open field lines of the pulsar magnetosphere.

Zoom In Zoom Out Reset image size

Consider the source S moving along the magnetic field line C and experiencing acceleration (a) in the direction of curvature vector of the field line. We assume the source to be a bunch, which is nothing more than a point-like huge charge. In Cartesian coordinates, the position vector of a bunch moving along the dipolar magnetic field line is given by (see Equation (2) in Gangadhara 2004, hereafter G04)

$$\begin{eqnarray} \mbox{\boldmath $r$}& = & r_e \sin ^2 \theta \lbrace \cos \theta \cos \phi ^{\prime } \sin \alpha \nonumber\\ &&+\sin \theta (\cos \alpha \cos \phi \cos \phi ^{\prime }-\sin \phi \sin \phi ^{\prime }),\nonumber \\ & & \cos \phi ^{\prime } \sin \theta \sin \phi +\sin \phi ^{\prime }(\cos \theta \sin \alpha + \cos \alpha \cos \phi \sin \theta),\nonumber \\ & & \cos \alpha \cos \theta -\cos \phi \sin \alpha \sin \theta \rbrace, \end{eqnarray} \tag{ 1 }$$

where r_e is the field line constant, and the angles θ and ϕ are the magnetic colatitude and azimuth, respectively. Next, ϕ' is the rotation phase and α is the inclination angle of the magnetic axis. Equation (1) describes the dipolar magnetic field lines presented in Figure 1. Then the velocity of the bunch is given by

$$\begin{equation} \mbox{\boldmath $v$} = \frac{d \mbox{\boldmath $r$}}{dt}= \left(\frac{\partial \mbox{\boldmath $r$}}{\partial \theta }\right) \left(\frac{\partial \theta }{\partial t}\right) = \left(\frac{\partial \theta }{\partial t}\right)\mbox{\boldmath $b$}, \end{equation} \tag{ 2 }$$

where $\mbox{\boldmath $b$}=\partial \mbox{\boldmath $r$}/\partial \theta$ is the magnetic field line tangent. Consider the magnetic axis

$$\begin{equation} \hat{m} = \lbrace \sin \alpha \cos \phi ^{\prime }, \sin \alpha \sin \phi ^{\prime }, \cos \alpha \rbrace. \end{equation} \tag{ 3 }$$

Due to curvature in the field lines, the plasma bunch, a point-like huge charge, collectively radiates relativistically beamed radiation in the direction of velocity v. The velocity v is parallel to the tangent $\mbox{\boldmath $b$}$ of the field line. To receive the beamed emission, the observer's line of sight ($\hat{n}$) must align with v within the beaming angle 1/γ, where γ is the Lorentz factor of the bunch. In other words, a distant observer at P receives beamed emission only when $\hat{n}\cdot \hat{v}=\cos \tau \sim 1$ for τ ≈ 1/γ, where $\hat{v}=\mbox{{\boldmath $v$}}/ \vert \mbox{{\boldmath $v$}}\vert.$ Let s be the arc length of the field line. Then, $ds=\vert \mbox{\boldmath $b$}\vert d\theta,$ where $\vert \mbox{\boldmath $b$}\vert =(r_e/\sqrt{2})\sin \theta \sqrt{5+3\cos (2\theta)},$ and the magnitude of velocity v = ds/dt = κc, where the parameter κ specifies the speed of bunch as a fraction of the speed of light c. Hence, we have

$$\begin{equation} \mbox{\boldmath $v$}=\kappa {c}\hat{b}, \end{equation} \tag{ 4 }$$

where

$$\begin{eqnarray} \hat{b} &=&\mbox{\boldmath $b$}/|\mbox{\boldmath $b$}| \nonumber\\ &= &\lbrace \cos \tau \cos \phi ^{\prime } \sin \alpha +\sin \tau (\cos \alpha \cos \phi \cos \phi ^{\prime }-\sin \phi \sin \phi ^{\prime }),\nonumber \\ & & \cos \phi ^{\prime } \sin \tau \sin \phi +\sin \phi ^{\prime }(\cos \tau \sin \alpha + \cos \alpha \cos \phi \sin \tau),\nonumber \\ & & \cos \alpha \cos \tau -\cos \phi \sin \alpha \sin \tau \rbrace, \end{eqnarray} \tag{ 5 }$$

and τ is the angle between $\hat{m}$ and $\hat{b}.$ In terms of the polar angle θ, the angle τ is given by

$$\begin{equation} \tan \tau = \frac{\sin \tau }{\cos \tau } = \frac{3\sin (2\theta)}{1+3\cos (2\theta)}, \end{equation} \tag{ 6 }$$

where

$$\begin{eqnarray} \cos \tau & = &\hat{b}\cdot \hat{m}=\frac{1+3\cos (2\theta)}{\sqrt{10+6\cos (2\theta)}},\nonumber \\ \sin \tau & = &(\hat{m} \times \hat{b})\cdot \hat{e}_\phi = \frac{3\sin (2\theta)}{\sqrt{10+6\cos (2\theta)}}\nonumber, \end{eqnarray}$$

and

$$\begin{eqnarray} \hat{e}_\phi &=& \{{-}\cos \alpha \sin \phi \cos \phi ^{\prime }- \cos \phi \sin \phi ^{\prime }, \cos \phi \cos \phi ^{\prime }\nonumber\\ &&- \cos \alpha \sin \phi \sin \phi ^{\prime }, \sin \alpha \sin \phi \} \end{eqnarray} \tag{ 7 }$$

is the bi-normal to the field line. We solve Equation (6) for θ, and obtain

$$\begin{equation} \cos (2\theta) = \frac{1}{3}(\cos \tau \sqrt{8+\cos ^2\tau }-\sin ^2\tau). \end{equation} \tag{ 8 }$$

Hence, from Equation (4) it is clear that to receive the radiation emitted in the direction of tangent $\hat{b},$ the sight line $\hat{n}$ must line up with it. So, by solving $\hat{n}\cdot \hat{b}=1$ or $\hat{n}\times \hat{b}=0,$ we can identify the tangent $\hat{b},$ which aligns with $\hat{n},$ and hence find the field line curvature and the coordinates (θ, ϕ) of the emission spot (see Equations (4), (9), and (11) in G04). Next, the acceleration of the bunch is given by

$$\begin{equation} \mbox{\boldmath $a$} = \frac{\partial {\mbox{\boldmath $v$}}}{\partial t} = \frac{(\kappa \, c)^2}{\vert {\mbox{\boldmath $b$}}\vert }\frac{\partial \hat{b}}{\partial \theta }= (\kappa \, c)^2 \mbox{\boldmath $k$}, \end{equation} \tag{ 9 }$$

where $\mbox{\boldmath $k$} = ({1}/{\vert \mbox{\boldmath $b$}\vert }) {\partial \hat{b}}/{\partial \theta }$ is the curvature (normal) of the field line. Then the radius of curvature of the field line is given by

$$\begin{equation} \rho = \frac{1}{\vert \mbox{\boldmath $k$}\vert } = \left[2-\frac{8}{3\lbrace 3+\cos (2\theta)\rbrace }\right]|\mbox{\boldmath $b$}|. \end{equation} \tag{ 10 }$$

Therefore, using $\mbox{\boldmath $k$}=\hat{k}/\rho,$ we can write

$$\begin{equation} \mbox{\boldmath $a$} = \frac{(\kappa \, c)^2}{\rho } \hat{k}, \end{equation} \tag{ 11 }$$

where

$$\begin{eqnarray} \hat{k} & \,{=}\, & \lbrace (\cos \alpha \cos \phi \cos \phi ^{\prime }\!-\!\sin \phi \sin \phi ^{\prime }) \cos \tau \,{-}\,\cos \phi ^{\prime }\sin \alpha \sin \tau, \nonumber \\ & &(\cos \phi ^{\prime } \sin \phi +\cos \alpha \cos \phi \sin \phi ^{\prime })\cos \tau - \sin \alpha \sin \phi ^{\prime }\sin \tau,\nonumber \\ & & -\cos \phi \sin \alpha \cos \tau -\cos \alpha \sin \tau \rbrace. \end{eqnarray} \tag{ 12 }$$

The relativistic bunch, i.e., point-like huge charge, q collectively emits curvature radiation as it accelerates along the curved trajectory C (see Figure 1). Then the electric field of the radiation at the observation point P is given by (Jackson 1975):

$$\begin{equation} \mbox{\boldmath $E$}(\mbox{\boldmath $r$}, t) = \frac{q}{c}\left[\frac{\hat{n}\times [(\hat{n}- \mbox{\boldmath $\beta $})\times {\dot{\bm \beta }}]}{R\, \xi ^3} \right]_{\rm ret}, \end{equation} \tag{ 13 }$$

where $\xi =1-\mbox{\boldmath $\beta $}\cdot \hat{n},\,R$ is the distance from the radiating region to the observer, $\mbox{\boldmath $\beta $}=\mbox{\boldmath $v$}/c$ is the velocity, and $\dot{\bm \beta }=\mbox{\boldmath $a$}/c$ is the acceleration of the bunch.

The radiation emitted by a relativistic bunch has a broad spectrum, and it can be estimated by taking the Fourier transformation of the electric field of radiation:

$$\begin{equation} \mbox{\boldmath $E$}(\mbox{\boldmath $r$},\omega) = \frac{1}{\sqrt{2\pi }}\int _{-\infty }^{+\infty } \mbox{\boldmath $E$}(\mbox{\boldmath $r$},t) e^{i\,\omega t} dt. \end{equation} \tag{ 14 }$$

In Equation (13), ret means evaluated at the retarded time t' + R(t')/c = t. By changing the variable of integration from t to t', we obtain

$$\begin{equation} \mbox{\boldmath $E$}(\mbox{\boldmath $r$},\omega) = \frac{1}{\sqrt{2\pi }}\frac{q}{c}\!\!\int _{-\infty }^{+\infty }\frac{\hat{n}\times [(\hat{n}-\mbox{\boldmath $\beta $}) \times \dot{\bm \beta }]}{R\, \xi ^2}e^{i\omega \lbrace t^{\prime }+R(t^{\prime })/c\rbrace } dt^{\prime }, \end{equation} \tag{ 15 }$$

where we have used dt = ξ dt'. When the observation point is far away from the region of space where the acceleration occurs, the propagation vector or the sight line $\hat{n}$ can be taken to be constant in time. Furthermore, the distance R(t') can be approximated as $R(t^{\prime })\approx R_{\rm 0}-\hat{n}\cdot \mbox{\boldmath $r$}(t^{\prime }),$ where R₀ is the distance between the origin O and the observation point P, and $\mbox{\boldmath $r$}(t^{\prime })$ is the position of the bunch relative to O.

Since bunches move with velocity κc along the dipolar field lines, over the incremental time dt the distance (arc length) covered is $ds = \kappa \, c\, dt = \vert \mbox{\boldmath $b$} \vert d\theta.$ Therefore, we have

$$\begin{equation} t = \frac{1}{\kappa \, c}\int \vert \mbox{\boldmath $b$} \vert d\theta = \frac{r_{\rm e}}{\sqrt{2}\kappa \, c}\int \sin \theta \sqrt{5+3 \cos (2 \theta)}\, d\theta. \end{equation} \tag{ 16 }$$

By choosing t = 0 at θ = 0, we obtain

$$\begin{eqnarray} t & = &\frac{r_{\rm e}}{12 \kappa c} [12+\sqrt{3} \log (14+8\sqrt{3})-3 \sqrt{10+6 \cos (2 \theta)} \cos (\theta) \nonumber \\ & & -2 \sqrt{3} \log (\sqrt{6} \cos (\theta)+ \sqrt{5+3 \cos (2 \theta)}) ]. \end{eqnarray} \tag{ 17 }$$

By assuming κ ∼ 1, in Figure 2, we plotted t as a function of θ for different r_e. It shows time t increases much faster at larger r_e than at lower. This is due to the fact that for a given range of θ the arc length of the field line becomes larger at higher r_e.

Zoom In Zoom Out Reset image size

Then Equation (15) becomes

$$\begin{equation} \mbox{\boldmath $E$}(\mbox{\boldmath $r$},\omega)\approx \frac{q\, {\rm e}^{i\omega R_{\rm 0}/ c}}{\sqrt{2\pi } R_{\rm 0}\kappa c^2} \!\!\int _{-\infty }^{+\infty } \!\!\!\vert \mbox{\boldmath $b$}\vert \frac{\hat{n}\times [(\hat{n}- \mbox{\boldmath $\beta $}) \times \dot{\bm \beta }]}{\xi ^2} {\rm e}^{i\omega \lbrace t-\hat{n}.{\bf r}/c\rbrace } d\theta, \end{equation} \tag{ 18 }$$

where the expression for t is given by Equation (17). Note that the prime on the time variable t has been omitted for brevity. The integration limits have been extended to ±∞ for mathematical convenience, as the integrand vanishes for |θ − θ₀|>1/γ. At any rotation phase ϕ', there exists a magnetic colatitude θ₀ and a magnetic azimuth ϕ₀ at which the field line tangent $\hat{b}$ exactly aligns with $\hat{n},$ i.e., $\hat{b}_0\cdot \hat{n}=1$ and τ = Γ, where Γ is the half-opening angle of the pulsar emission beam centered on $\hat{m}.$ The expressions for θ₀ and ϕ₀ are given in G04 (see Equations (9) and (11)).

The polarization state of the emitted radiation can be determined using $\mbox{\boldmath $E$}(\omega)$ with the known r(t), β, and $\dot{\bm \beta }$. Since the integral in Equation (18) has to be computed over the path of the particle, the line of sight $\hat{n}$ can be chosen without loss of generality, to lie in the xz-plane:

$$\begin{equation} \hat{n} = (\sin \zeta, 0, \cos \zeta), \end{equation} \tag{ 19 }$$

where ζ = α + σ is the angle between $\hat{n}$ and $\hat{\Omega },$ and σ is the closest impact angle of $\hat{n}$ with respect to $\hat{m}.$

Let

$$\begin{equation} \mbox{\boldmath $A$} = \frac{1}{\kappa c}\vert \mbox{\boldmath $b$}\vert \frac{\hat{n}\times [(\hat{n}- \mbox{\boldmath $\beta $})\times \dot{\bm \beta }]}{\xi ^2}. \end{equation} \tag{ 20 }$$

By substituting for acceleration $\dot{\bm \beta }={\bm a}/$c from Equation (9), we can reduce it to

$$\begin{equation} \mbox{\boldmath $A$} = \lbrace A_{\it x},\, A_{\it y},\, A_{\it z}\rbrace = \frac{\hat{n}\times [(\hat{n}- \mbox{\boldmath $\beta $})\times \mbox{\boldmath $N$}]}{\xi ^2}, \end{equation} \tag{ 21 }$$

where N $=\kappa \, \partial \hat{b}/\partial \theta$= $\partial \mbox{\boldmath $\beta $}/\partial \theta.$ Using the expression β = v/c from Equation (4) and series expanding A in powers of θ about θ₀ we obtain

$$\begin{eqnarray} A_{\it x} & = & A_{\rm x0}+ A_{\rm x1} (\theta -\theta _0)+A_{\rm x2} (\theta -\theta _0)^2+A_{\rm x3} (\theta -\theta _0)^3\nonumber \\ &&+O[(\theta -\theta _0)^4],\nonumber \\ A_{\it y} & = & A_{\rm y0}+ A_{\rm y1} (\theta -\theta _0)+A_{\rm y2} (\theta -\theta _0)^2+A_{\rm y3} (\theta -\theta _0)^3\nonumber \\ &&+O[(\theta -\theta _0)^4],\nonumber \\ A_{\it z} & = & A_{\rm z0}+ A_{\rm z1} (\theta -\theta _0)+A_{\rm z2} (\theta -\theta _0)^2+A_{\rm z3} (\theta -\theta _0)^3\nonumber \\ &&+O[(\theta -\theta _0)^4], \end{eqnarray} \tag{ 22 }$$

where $A_{{\it x}\it {i}},\, A_{{\it y}\it {i}}$ and $A_{{\it z}\it {i}}$ with i = 0,  1,  2,  3 are the series expansion coefficients, and their expressions are given in Appendix A.

The scalar product between $\hat{n}$ and r is given by

$$\begin{eqnarray} \hat{n}\cdot \mbox{\boldmath $r$} & = & r_e\sin ^2\theta [\cos \alpha (\cos \theta \cos \zeta + \cos \phi \cos \phi ^{\prime }\sin \theta \sin \zeta) \nonumber \\ & &-\cos \zeta \cos \phi \sin \alpha \sin \theta + \sin \zeta (\cos \theta \cos \phi ^{\prime }\sin \alpha \nonumber\\ &&- \sin \theta \sin \phi \sin \phi ^{\prime })]. \end{eqnarray} \tag{ 23 }$$

Next, substituting the expressions of t and $\hat{n}\cdot \mbox{\boldmath $r$}$ into the argument of the exponential in Equation (18), and series expanding in powers of θ about θ₀ we obtain

$$\begin{eqnarray} \omega \left(t-\frac{\hat{n}\cdot \mbox{\boldmath $r$}}{c}\right)& = & c_{\rm 0}+ c_{\rm 1} (\theta -\theta _0)+c_{\rm 2} (\theta -\theta _0)^2\nonumber\\ &&+c_{\rm 3} (\theta -\theta _0)^3+O[(\theta -\theta _0)^4], \end{eqnarray} \tag{ 24 }$$

where c₀, c₁,  c₂, and c₃ are the series expansion coefficients, and their expressions are given in Appendix A.

Now, by substituting the expressions of Equations (22) and (24) into Equation (18), we obtain the components of $\mbox{\boldmath $E$}(\omega)=\lbrace E_x(\omega),\, E_y(\omega),\, E_z(\omega)\rbrace :$

$$\begin{eqnarray} E_{\it x}(\omega)& = & E_0 \!\int _{-\infty }^{+\infty } (A_{\rm x0}+ A_{\rm x1} \,\mu +A_{\rm x2} \,\mu ^2\nonumber\\ &&+A_{\rm x3} \mu ^3) {\rm e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3} \,\mu ^3)}d\mu, \nonumber \\ E_{\it y}(\omega)& = & E_0 \!\int _{-\infty }^{+\infty } (A_{\rm y0}+ A_{\rm y1} \,\mu +A_{\rm y2} \,\mu ^2\nonumber\\ &&+A_{\rm y3} \,\mu ^3) {\rm e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3} \,\mu ^3)}d\mu, \nonumber \\ E_{\it z}(\omega)& = & E_0 \!\int _{-\infty }^{+\infty } (A_{\rm z0}+ A_{\rm z1} \,\mu +A_{\rm z2} \,\mu ^2\nonumber\\ &&+A_{\rm z3} \,\mu ^3) {\rm e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3 } \,\mu ^3)}d\mu, \end{eqnarray} \tag{ 25 }$$

where μ = θ − θ₀ and

$$\begin{equation*} E_0 = \frac{q}{\sqrt{2\pi } R_{\rm 0} c}{e}^{i[(\omega R_{\rm 0}/c)+c_0]}. \end{equation*}$$

Now by substituting the integral solutions S₀,  S₁,  S₂, and S₃, given in Appendix B, into Equation (25) we obtain

$$\begin{eqnarray} E_{\it x}(\omega) & = & E_0 (A_{\rm x0}S_0+ A_{\rm x1}S_1+A_{\rm x2}S_2+ A_{\rm x3}S_3), \nonumber \\ E_{\it y}(\omega) & = & E_0 (A_{\rm y0}S_0+ A_{\rm y1}S_1+A_{\rm y2}S_2+ A_{\rm y3}S_3), \nonumber \\ E_{\it z}(\omega) & = & E_0 (A_{\rm z0}S_0+ A_{\rm z1}S_1+A_{\rm z2}S_2+ A_{\rm z3}S_3). \end{eqnarray} \tag{ 26 }$$

To find the polarization angle of radiation field E, we need to specify two reference directions perpendicular to the sight line $\hat{n}.$ One could be the projected spin axis on the plane of the sky: $\hat{\epsilon }_\parallel =(-\cos \zeta,\,0, \,\sin \zeta),$ and then the other direction is specified by $\hat{\epsilon }_\perp =\hat{\epsilon }_\parallel \times \hat{n}=\hat{y},$ where $\hat{y}$ is a unit vector parallel to the y-axis. Then the components of E in the directions $\hat{\epsilon }_\parallel$ and $\hat{\epsilon }_\perp$ are given by

$$\begin{eqnarray} E_\parallel & = &\hat{\epsilon }_\parallel \cdot \mbox{\boldmath $E$} = {-}\cos \zeta \, E_{\it x}+\sin \zeta \,E_{\it z}, \nonumber \\ E_\perp & = &\hat{\epsilon }_\perp \cdot \mbox{\boldmath $E$} =E_{\it y}. \end{eqnarray} \tag{ 27 }$$

At any rotation phase ϕ', the observer receives the radiation from all those field lines whose tangents lie within the angle 1/γ with respect to the sight line $\hat{n}.$ Let η be the angle between the $\hat{b}$ and $\hat{n},$ then $\cos \eta =\hat{b}\cdot \hat{n},$ and the maximum value of η is 1/γ. Therefore, at ϕ = ϕ₀ we solve $\cos (1/\gamma)=\hat{b}\cdot \hat{n}$ for τ, and find the allowed range (Γ − 1/γ) ⩽ τ ⩽ (Γ + 1/γ) of τ or −1/γ ⩽ η ⩽ 1/γ of η, which in turn allows one to find the range of θ with the help of Equation (8). Next, for any given η within its range, we find ϕ by solving $\cos \eta =\hat{b}\cdot \hat{n}.$ It gives (ϕ₀ − δϕ) ⩽ ϕ ⩽ (ϕ₀ + δϕ), where

$$\begin{equation} \cos (\delta \phi)=\frac{ \sin \Gamma [\cos (1/\gamma) \csc (\Gamma +\eta)- \cos \Gamma \cot (\Gamma +\eta)]}{(\cos \zeta \sin \alpha -\cos \alpha \cos \phi ^{\prime }\sin \zeta)^2+\sin ^2\zeta \sin ^2\phi ^{\prime }}. \end{equation} \tag{ 28 }$$

Hence by knowing the ranges of θ and ϕ at any given ϕ', we can estimate the contributions to E from all those field lines, whose tangents lie within the angle 1/γ with respect to $\hat{n}.$ In Figure 3, we have plotted those regions at three phases: ϕ' = −30°, 0°, and 30° using α = 10°,  β = 5°, and γ = 400. Note that at the center of each region, $\hat{b}$ exactly aligns with the sight line, i.e., $\hat{b} \cdot \hat{n}=1.$ Further, in Figure 4, we have plotted them for −180° ⩽ ϕ' ⩽ 180° with a step of 5° between the successive regions. We observe that the range of θ stays nearly constant (or decreasing negligibly) whereas that of ϕ gets narrower with respect to the increasing |ϕ'|.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To understand the pulsar radio emission, we must model all the Stokes parameters (I,  Q,  U, and V—a set of parameters used to specify the phase and polarization of radiation), and compare with observations. They have been found to offer a very convenient method for establishing the association between the polarization state of observed radiation and the geometry of the emitting region. They are defined as follows:

$$\begin{eqnarray} I &=& E_\parallel E_\parallel ^*+E_\perp E_\perp ^*, \quad Q = E_\parallel E_\parallel ^*-E_\perp E_\perp ^*, \nonumber\\ \quad U &=& 2\, {\rm Re}[E_\parallel ^* E_\perp ], \quad V = 2\, {\rm Im}[E_\parallel ^* E_\perp ]. \end{eqnarray} \tag{ 29 }$$

The parameter I defines the total intensity, Q and U jointly define the linear polarization and its position angle, and V describes the circular polarization.

3.1. Addition of Stokes Parameters

Let W_I be the energy radiated coherently per unit solid angle per unit frequency interval per particle bunch (Jackson 1975), then

$$\begin{eqnarray} &&\frac{d^2 W_I}{d\omega \,d\Omega }= \frac{c\, R_0^2}{2\pi } \!|\mbox{\boldmath $E$}(\omega)|^2\!. \end{eqnarray} \tag{ 30 }$$

Since the Stokes parameter $I=E_\parallel {E_\parallel }^*+E_\perp {E_\perp }^*=\mbox{\boldmath $E$} \,{\bm{\cdot}}\, \mbox{\boldmath $E$}^*=|\mbox{\boldmath $E$}|^2,$ we can rewrite Equation (30) as

$$\begin{eqnarray} &&I= |\mbox{\boldmath $E$}|^2 = \frac{2\pi }{c\, R_0^2} \frac{d^2 W_I}{d\omega \,d\Omega }. \end{eqnarray} \tag{ 31 }$$

Similarly, we can express Q,  U, and V as

$$\begin{eqnarray} Q &=& \frac{2\pi }{c\, R_0^2} \frac{d^2 W_Q}{d\omega \,d\Omega },\nonumber \\ U &=& \frac{2\pi }{c\, R_0^2} \frac{d^2 W_U}{d\omega \,d\Omega },\nonumber \\ V &=& \frac{2\pi }{c\, R_0^2} \frac{d^2 W_V}{d\omega \,d\Omega }. \end{eqnarray} \tag{ 32 }$$

The net emission, which the observer receives along $\hat{n}$, will have contributions from the neighboring field lines, whose tangents are within the angle 1/γ with respect to $\hat{n}.$ Hence the radiation received at any given phase is the net contribution from a small tube of field lines having an angular width of about 2/γ. Thus, the radiation in the direction of $\hat{n}$ should be integrated over a solid angle dΩ = sin θ dθ dϕ. We choose limits on the angles ϕ and θ such that the integration over them will cover the solid angular region (beaming region) of radial width 1/γ around $\hat{n}.$ Since θ and ϕ are orthogonal, choosing them as the variables of integration is justified. We assume that (1) the width of bunch η₀ is much smaller than the wavelength λ of the radio waves, so that the radiation emitted by a bunch is coherent, and (2) the bunches, within the beaming region, are closely spaced, so that the net emission becomes smooth and continuous.

Consider a bunch having γ ∼ 400 emitting radio waves at frequency ν = 600 MHz at an altitude of about 400 km. Note that these values are close to those estimated in G04 in the case of PSR B0329+54. Then the angular width of the beaming region corresponding to 2/γ is ∼03, which corresponds to a width of ∼2 km at an altitude of 400 km. For coherence to be effective the bunch width w₀ < λ. Therefore, we choose w₀ < 50 cm for λ ∼ 50 cm. Since these values of w₀ are much smaller than the width of the beaming region (∼2 km), the Stokes parameters can be integrated as continuous functions of θ and ϕ.

Let I_s be the resultant Stokes intensity parameter then

$$\begin{eqnarray} I_{\rm s}&=& \int \,I\,d\Omega \nonumber \\ & =& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } \int ^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } I\, \sin \theta \,d\theta \,d\phi, \end{eqnarray} \tag{ 33 }$$

where θ₀ and ϕ₀ are the magnetic colatitude and azimuth of the sight line $\hat{n}.$ Similarly, for other Stokes parameters, we have

$$\begin{eqnarray} Q_{\rm s} &=& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } \int ^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } Q\, \sin \theta \,d\theta \,d\phi,\nonumber \\ U_{\rm s} &=& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } \int ^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } U\, \sin \theta \,d\theta \,d\phi,\nonumber \\ V_{\rm s} &=& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } \int ^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } V\, \sin \theta \,d\theta \,d\phi. \end{eqnarray} \tag{ 34 }$$

Then the linear polarization is given by

$$\begin{equation} L_{\rm s}=\sqrt{Q_{\rm s}^2+U_{\rm s}^2}, \end{equation} \tag{ 35 }$$

and the corresponding polarization angle is

$$\begin{equation} \psi _{\rm s}=\frac{1}{2} \tan ^{-1}\left(\frac{U_{\rm s}}{Q_{\rm s}}\right). \end{equation} \tag{ 36 }$$

The emission in spin-powered pulsars is mostly of non-thermal origin. If the radiation field E from different sources does not bear any phase relation then they are expected to be incoherently superposed on the observation point. On the other hand, if there is a phase relation then they are coherently superposed. From the observational point of view both cases are important.

By considering the relativistic pair plasma with γ = 400 accelerated along the dipolar field lines of a pulsar with period P = 1 s, we computed the polarization parameters and plotted them in Figure 5. It shows a stronger emission near the meridional plane, where the beaming region is broader (see Figure 3) and the radius of curvature ρ goes to a minimum. The profile of linear polarization L_s resembles the intensity profile, except for its lower magnitude due to the incoherent addition. To describe the behaviors of circular polarization V_s and polarization angle ψ_s, we define the symbols: "−/+," transition of the right hand circular to left hand circular; "+/− ," left hand circular to right hand circular; "cw," clockwise rotation of the polarization angle; and "ccw" counterclockwise rotation.

Zoom In Zoom Out Reset image size

Since the circular polarization V_s changes sign as −/+ or +/− as the sight line cuts across the field line, the net circular polarization goes to zero in a uniform emission due to the addition with opposite signs. The polarization angle swings reproduced in Figure 5 are consistent with the rotating vector model of Radhakrishnan & Cooke (1969). In the case of positive sight line impact parameter (σ = 5°), the polarization angle swing is ccw as the slope dψ_s/dϕ'>0 while in the negative case (σ = −5°), it is cw as dψ_s/dϕ' < 0.

4.1. Modulation of Radio Emission

Pulsar radio emission is believed to come from mostly open magnetic field lines, whose foot points define the polar cap. The shape of pulsar profiles indicates that the entire polar cap does not radiate; only some selected regions radiate, which may be organized into a central core emission and coaxial conal emissions, which has an overwhelming support from observations (e.g., Rankin 1990, 1993). Hence, the radiating region above the polar is believed to have a central column of emission (core) and a few coaxial conal regions of emission (cones; e.g., Gil & Krawczyk 1997; Gangadhara & Gupta 2001; Gupta & Gangadhara 2003; Dyks et al. 2004).

4.1.1. Modulating Function

It is well known that the components of a pulsar profile can be decomposed into individual Gaussians by fitting one with each of the subpulse component. For example, the components in the pulse profile of PSR 1706−16 and PSR 2351+61 are fitted with appropriate Gaussians by Kramer et al. (1994). When the line of sight crosses the emission region, it encounters a pattern in intensity due to Gaussian modulation in the azimuthal direction. Because of the Gaussian modulation in the azimuthal direction, the intensity becomes nonuniform in the polar directions too. These arguments indicate that a Gaussian-like intensity modulation exists in the polar directions too. So, we assume that the emission region of a pulse component has an intensity modulation in both azimuthal directions. Hence, we define a modulation function f for a pulse component as

$$\begin{equation} f(\theta,\phi) = f_0\, \exp \left[{- \left(\frac{\phi -\phi _{\rm p}}{\sigma _\phi }\right)^2}\right], \end{equation} \tag{ 37 }$$

where ϕ_p is the peak location of the Gaussian function and f₀ is the amplitude. If w_ϕ is the full width at half-maximum (FWHM), then $\sigma _\phi = w_\phi /(2\sqrt{\ln 2}).$

Taking into account the modulation, Equations (33)–(34) can be written as

$$\begin{eqnarray} I_{\rm s} & =& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } {\int}^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } f\,I\, \sin \theta \,d\theta \,d\phi \nonumber \\ Q_{\rm s} &=& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } {\int}^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } f\,Q\, \sin \theta \,d\theta \,d\phi,\nonumber \\ U_{\rm s} &=& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } {\int}^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } f\,U\, \sin \theta \,d\theta \,d\phi,\nonumber \\ V_{\rm s} &=& \int ^{\theta _{\rm 0}+\delta \theta }_{\theta _{\rm 0}-\delta \theta } {\int}^{\phi _0+\delta \phi }_{\phi _0-\delta \phi } f\,V\, \sin \theta \,d\theta \,d\phi. \end{eqnarray} \tag{ 38 }$$

Using a Gaussian with peak located at the meridional plane (ϕ' = 0°), we have computed the pulse profiles in the two cases of impact parameter (σ) and inclination angles (α) and plotted in Figures 6 and 7. We observe that the profile in the case of negative σ is broader than the positive case. This difference is due to the projection of the emission region onto the equatorial plane of the pulsar. In the case of positive σ, the polarization angle χ_s swing is ccw and the sign change of V_s is −/+ with respect to ϕ', while in the case of negative σ the χ_s swing is cw and the sign change of V_s is +/− . Hence, we find that the polarization angle swing is correlated with the circular polarization sign reversal. This correlation is invariant with respect the stellar spin directions.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The mean pulsar profiles are often found to consist of an odd number of multi-components or subpulses. Many works on pulsar profiles (e.g., Rankin 1990, 1993; Mitra & Deshpande 1999) propose that the pulsar emission beam has a nested conal structure. To investigate the polarization of subpulses in such profiles, we have reproduced a five-component profile by considering three Gaussians in Figure 8, and five Gaussians in Figures 9 and 10. The central component is presumed to be a core, and the other components are symmetrically located on either side of the core forming the cones. We find across each component that circular polarization changes the sign and is correlated with the polarization angle swing. We also observe that the circular polarization of the outermost components is weaker compared to that of inner ones, which is quite clear in the case of large inclination angles. This is due to the fact that the sight line crosses the field lines in the almost edge-on position in the case of outermost components. The small distortions in the polarization angle curve are due to modulation.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Observed pulsar radio luminosities together with the small source size imply extraordinarily high brightness temperatures, i.e., as high as 10³¹ K. The incoherent sum of a single-particle curvature radiation is not enough to explain the very high brightness temperature of pulsar radio emission; therefore, one is forced to postulate the existence of charged bunches. To avoid implausibly high particle densities and energies, coherent radiation processes are invoked. Pacini & Rees (1970) and Sturrock (1971) among others were quick to point out that the observed coherence may be due to bunching of particles in the emission region of the magnetosphere. The problem of bunch formation has been known for many decades, and it has already been addressed by many authors (e.g., Karpman et al. 1975; Ruderman & Sutherland 1975; Cox 1979; Asseo et al. 1990; Gil et al. 2004). The natural mechanism for the formation of charged bunches was first proposed by Karpman et al. (1975). They have argued that the modulational instability in the turbulent plasma generates charged solitons, provided that species of different charges have different masses. One should mention here that to explain coherent radio emission we do not necessarily need stable solitons but only large-scale (as compared with the Langmuir wavelength) charge density fluctuations. Gil et al. (2004) generalized the soliton model by including formation and propagation of the coherent radiation in the magnetospheric plasma. However, it is not easy to form such charge bunches (see Melrose 1992 for a review). Further, Michel (1991) has pointed out that the pair-production discharge mechanism originally applied to pulsars by Sturrock automatically produces dense bunches that can produce coherence at radio frequencies with sufficient intensity to simulate pulsar action. If the bunches of plasma particles with sizes much smaller than a wavelength of radiation exist then the net radiation field $\mbox{\boldmath $E$}(\omega)\approx N\mbox{\boldmath $E$}_{\rm o}(\omega),$ where N is the number of charges present in the bunch. Hence, the total radiation field due to a bunch of particles is equal to the vector sum of the fields radiated by each charge.

In the general framework of models including ours, in which the radio power is curvature radiation emitted by charge bunches constrained to follow the field lines, the linear polarization is intrinsic to the emission mechanism and is, furthermore, a purely geometric property. In this direction, the recently achieved observational results and the model predications based on them by Mitra et al. (2009) become very relevant. They find that the polarization angle of linear polarization in subpulses follow closely the mean polarization angle curve at the corresponding profile components and argue that their findings favor coherent curvature radiation over maser mechanism as the observed emission.

In an actual case, it is the combination of both incoherent and coherent superpositions determining the polarization state of the observed emission. Though the emissions from a single bunch is highly polarized, the radiation received by a distant observer will be less polarized, as the radiation from many bunches is incoherently superposed (Gil & Rudnicki 1985). Also, the degree of polarization is found to depend on the time resolution chosen in the observation (Gangadhara et al. 1999). Circular polarization is generally strongest in the central regions of a profile, but is by no means confined to these regions. It has been detected from conal components of many pulsars, for example, conal-double pulsars, and found to be highly correlated with the polarization angle swing (You & Han 2006). In most of the cases, the sense reversal of circular polarization is nearly independent of frequency, suggesting that the circular polarization does not arise from propagation or plasma effects (Michel 1987; Radhakrishnan & Rankin 1990). Radhakrishnan & Rankin (1990) have argued that the circular polarization is intrinsically antisymmetric type and correlated with the polarization angle swing. The antisymmetric circular polarization of curvature radiation becomes significant if there are gradients in the emissivity over angular scales comparable with the emission cone of single charge. Their results are consistent with the predictions of our model (Figures 6–10) and strongly suggest that the correlation of antisymmetric circular polarization with the polarization angle swing is a geometric property of emission mechanism. Since our model deals with the steady flow of relativistic plasma bunches along dipolar field lines, it is relevent only for average profiles, and reflects the results which are more of geometric dependent. We have not considered any fluctuations or instability in the plasma flow. Hence, it cannot reproduce the behaviors of single pulses.

By adopting the antenna mechanism, Buschauer & Benford (1976) have derived a new formalism for the relativistic curvature radiation. However, the treatment given does not include the detailed geometry of dipolar field lines and the estimation of polarization, particularly circular, as we have considered in our model. Since radiation from many bunches is superposed on any given pulse longitude, the circular polarization of different signs and magnitudes is added. The result of such an addition could be the reason for the diversities in the observed circular polarization. Since our model was aimed at analyzing the intrinsic polarization properties of coherent radiation, we plan to consider the propagation effects separately in subsequent works. We speculate that the propagation origin of antisymmetric circular polarization is very unlikely but the symmetric circular polarization may be possible.

By taking into account of a detailed geometry of dipolar magnetic field lines, we have derived the polarization state of the coherent curvature radiation due to relativistic plasma in the pulsar magnetosphere and drawn the following conclusions:

• 1.  
  We do confirm the previous results of Gil & Snakowski (1990a, 1990b) and Gil et al. (1993) that coherent curvature radiation has basically antisymmetric type of circular polarization. Though the emission from a single bunch is highly polarized, the net emission from many bunches within the beaming region is less polarized due to the incoherent superposition of radiation fields.
• 2.  
  Based on the Stokes parameters of the curvature radiation, we have deduced the polarization angle swing, i.e., the rotating vector model.
• 3.  
  Based on the coherent curvature radiation, we have achieved for the first time the result that the antisymmetric type of circular polarization is correlated with the polarization angle swing, and such correlations have been indeed found in the profiles of many pulsars.
• 4.  
  The addition of circular polarization with different signs and magnitudes at any given phase could be responsible for the wide diversity in circular polarization across the pulse. It is consistent with earlier results (e.g., Gil et al. 1993, 1995).

I thank J. L. Han and J. M. Rankin for illuminating discussions.

A.1. The Series Expansion Coefficients of Equation (22)

$$\begin{eqnarray} A_{\rm {x0}} &=& \frac{\cos \zeta \left(d_3 N_{\it x}(\theta _0)+d_2 N_{\it z}(\theta _0)\right)}{d_1^2}, \end{eqnarray} \tag{ A1 }$$

$$\begin{eqnarray} A_{\rm {x1}}&=& \frac{\cos \zeta \left(d_1 \left(d_9-f_1\right)-2 d_3 d_4 N_{\it x}\left(\theta _0\right)-2 d_2 d_4 N_{\it z}\left(\theta _0\right)\right)}{d_1^3}, \end{eqnarray} \tag{ A2 }$$

$$\begin{eqnarray} A_{\rm {x2}}&=& \frac{\left(d_{13} d_1^2-4 d_4 d_{12} d_1-2 d_7 d_8\right) \cos \zeta \,}{2 d_1^4}, \end{eqnarray} \tag{ A3 }$$

$$\begin{eqnarray} A_{\rm {x3}}&=& \frac{\left(d_1^3 (d_{11}-3 \left(f_4+f_5\right)-f_6)-2 \left(d_6 d_8+3 d_4 d_{13}\right) d_1^2+d_{14}\right) \cos \zeta \,}{6 d_1^5}, \end{eqnarray} \tag{ A4 }$$

$$\begin{eqnarray} A_{\rm {y0}}&=& \frac{d_1 N_{\it y}(\theta _0)-e_1 \beta _{\it y}\left(\theta _0\right)}{d_1^2}, \end{eqnarray} \tag{ A5 }$$

$$\begin{eqnarray} A_{\rm {y1}}&=& \frac{d_1^2 N_{\it y}^{\prime }\left(\theta _0\right)-d_1 e_4 \beta _{\it y}\left(\theta _0\right)+e_5 N_{\it y}\left(\theta _0\right)+e_1 e_2}{d_1^3}, \end{eqnarray} \tag{ A6 }$$

$$\begin{eqnarray} A_{\rm {y2}}&=& \frac{d_1 \left(e_{10}-2 d_1 d_4 N_{\it y}^{\prime }\left(\theta _0\right)\right)+d_1(2 d_4^2-d_1 d_5) N_{\it y}\left(\theta _0\right)+e_9 \beta _{\it y}\left(\theta _0\right)}{2 d_1^4}, \end{eqnarray} \tag{ A7 }$$

$$\begin{eqnarray} A_{\rm {y3}}&=& \frac{d_1 \left(3 d_1 (d_1 d_5-4 d_4^2-2 d_7) N_{\it y}^{\prime }(\theta _0)+d_1 e_{13}+3 e_{14}\right)-d_1 e_{11} N_{\it y}\left(\theta _0\right)+e_{12} \beta _{\it y}\left(\theta _0\right)}{6 d_1^5}, \end{eqnarray} \tag{ A8 }$$

$$\begin{eqnarray} A_{\rm {z0}}&=& \frac{\sin \zeta \left(-d_3 N_{\it x}(\theta _0)-d_2 N_{\it z}(\theta _0)\right)}{d_1^2}, \end{eqnarray} \tag{ A9 }$$

$$\begin{eqnarray} A_{\rm {z1}}&=& \frac{\sin \zeta (d_1 (-d_3 N_{\it x}^{\prime }(\theta _0)-d_2 N_{\it z}^{\prime }(\theta _0)+f_1)+2 d_4 (d_3 N_{\it x}(\theta _0)+d_2 N_{\it z}(\theta _0)))}{d_1^3}, \end{eqnarray} \tag{ A10 }$$

$$\begin{eqnarray} A_{\rm {z2}}&=& \frac{(d_1 (4 d_4 (d_9-f_1)+d_1 (-d_{10}+f_2+2 f_3)+2 d_5 d_8)-6 d_4^2 d_8) \sin \zeta \,}{2 d_1^4}, \end{eqnarray} \tag{ A11 }$$

$$\begin{eqnarray} A_{\rm {z3}}&=& \frac{(f_7+d_1^3 (f_6-d_{11}+3 (f_4+f_5))+2 d_1^2 (d_6 d_8+3 d_4 d_{10}-3 d_4 (f_2+2 f_3))) \sin \zeta \,}{6 d_1^5}, \end{eqnarray} \tag{ A12 }$$

where

$$\begin{eqnarray} d_1&=& \sin \zeta \,\beta _{\it x}(\theta _0)+\cos \zeta \, \beta _{\it z}(\theta _0)-1, \end{eqnarray} \tag{ A13 }$$

$$\begin{eqnarray} d_2&=& \sin \zeta \,-\beta _{\it x}\left(\theta _0\right), \end{eqnarray} \tag{ A14 }$$

$$\begin{eqnarray} d_3&=& \beta _{\it z}\left(\theta _0\right)-\cos \zeta \,, \end{eqnarray} \tag{ A15 }$$

$$\begin{eqnarray} d_4&=& \sin \zeta \, \beta _{\it x}^{\prime }\left(\theta _0\right)+ \cos \zeta \, \beta _{\it z}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A16 }$$

$$\begin{eqnarray} d_5&=& \sin \zeta \, \beta _{\it x}^{\prime \prime }\left(\theta _0\right)+ \cos \zeta \, \beta _{\it z}^{\prime \prime }\left(\theta _0\right), \end{eqnarray} \tag{ A17 }$$

$$\begin{eqnarray} d_6&=& \sin \zeta \, \beta _{\it x}^{(3)}\left(\theta _0\right)+ \cos \zeta \, \beta _{\it z}^{(3)}\left(\theta _0\right), \end{eqnarray} \tag{ A18 }$$

$$\begin{eqnarray} d_7&=& d_1 d_5-3 d_4^2, \end{eqnarray} \tag{ A19 }$$

$$\begin{eqnarray} d_8&=& d_3 N_{\it x}\left(\theta _0\right)+d_2 N_{\it z} \left(\theta _0\right), \end{eqnarray} \tag{ A20 }$$

$$\begin{eqnarray} d_9&=& d_3 N_{\it x}^{\prime }\left(\theta _0\right)+d_2 N_{\it z}^{\prime } \left(\theta _0\right), \end{eqnarray} \tag{ A21 }$$

$$\begin{eqnarray} d_{10}&=& d_3 N_{\it x}^{\prime \prime }\left(\theta _0\right)+d_2 N_{\it z}^{\prime \prime }\left(\theta _0\right), \end{eqnarray} \tag{ A22 }$$

$$\begin{eqnarray} d_{11}&=& d_3 N_{\it x}{}^{(3)}\left(\theta _0\right)+d_2 N_{\it z}{}^{(3)}\left(\theta _0\right), \end{eqnarray} \tag{ A23 }$$

$$\begin{eqnarray} d_{12}&=& d_9-N_{\it z}\left(\theta _0\right) \beta _{\it x}^{\prime } \left(\theta _0\right)+N_{\it x}\left(\theta _0\right) \beta _{\it z}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A24 }$$

$$\begin{eqnarray} d_{13}&=& d_{10} -2 N_{\it z}^{\prime }\left(\theta _0\right) \beta _{\it x}^{\prime } \left(\theta _0\right)+2 N_{\it x}^{\prime }\left(\theta _0\right) \beta _{\it z}^{\prime }\left(\theta _0\right)-N_{\it z}\left(\theta _0\right) \beta _{\it x}^{\prime \prime }\left(\theta _0\right)+ N_{\it x}\left(\theta _0\right) \beta _{\it z}^{\prime \prime }\left(\theta _0\right), \end{eqnarray} \tag{ A25 }$$

$$\begin{eqnarray} d_{14}&=& 18 d_1 d_4 d_5 d_8-24 d_8 d_4^3-6 d_1 d_7 d_{12}, \end{eqnarray} \tag{ A26 }$$

$$\begin{eqnarray} e_1 & = & \sin \zeta \, N_{\it x}\left(\theta _0\right)+\cos \zeta \, N_{\it z}\left(\theta _0\right), \end{eqnarray} \tag{ A27 }$$

$$\begin{eqnarray} e_2 & = & 2 d_4 \beta _{\it y}\left(\theta _0\right)-d_1 \beta _{\it y}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A28 }$$

$$\begin{eqnarray} e_3 & = & d_1 \beta _{\it x}^{\prime }\left(\theta _0\right)-2 d_4 \beta _{\it x}\left(\theta _0\right), \end{eqnarray} \tag{ A29 }$$

$$\begin{eqnarray} e_4 & = & \sin \zeta \, N_{\it x}^{\prime }\left(\theta _0\right)+\cos \zeta \, N_{\it z}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A30 }$$

$$\begin{eqnarray} e_5 & = & d_1 \cos \zeta \, \beta _{\it z}^{\prime }\left(\theta _0\right)- 2 d_4 \cos \zeta \, \beta _{\it z}\left(\theta _0\right)+2 d_4+e_3 \sin (\zeta), \end{eqnarray} \tag{ A31 }$$

$$\begin{eqnarray} e_6 & = & \sin \zeta \, N_{\it x}^{\prime \prime }\left(\theta _0\right)+\cos \zeta \, N_{\it z}^{\prime \prime }\left(\theta _0\right), \end{eqnarray} \tag{ A32 }$$

$$\begin{eqnarray} e_7 & = & \sin \zeta \, N_{\it x}{}^{(3)}(\theta _0)+\cos \zeta \, N_{\it z}{}^{(3)}(\theta _0), \end{eqnarray} \tag{ A33 }$$

$$\begin{eqnarray} e_8 & = & d_1 N_{\it y}{}^{(3)}\left(\theta _0\right)-e_1 \beta _{\it y}{}^{(3)} \left(\theta _0\right), \end{eqnarray} \tag{ A34 }$$

$$\begin{eqnarray} e_9 & = & 4d_1 d_4 e_4-6 d_4^2 e_1+d_1 \left(2 d_5 e_1-d_1 e_6\right), \end{eqnarray} \tag{ A35 }$$

$$\begin{eqnarray} e_{10} & = & d_1 \left(d_1 N_{\it y}^{\prime \prime } \left(\theta _0\right)-e_1 \beta _{\it y}^{\prime \prime } \left(\theta _0\right)\right)+\left(4 d_4 e_1-2 d_1 e_4\right) \beta _{\it y}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A36 }$$

$$\begin{eqnarray} e_{11} & = & 24 d_4^3+6d_4 \left(d_7-2 d_1 d_5\right) +d_1^2 d_6, \end{eqnarray} \tag{ A37 }$$

$$\begin{eqnarray} e_{12} & = & 24 d_4^3 e_1+6 d_1 d_4 \left(d_1 e_6-3 d_5 e_1\right)+d_1 \left(-d_1^2 e_7+2 d_6 d_1 e_1+6 d_7 e_4\right), \end{eqnarray} \tag{ A38 }$$

$$\begin{eqnarray} e_{13} & = & \left(6 d_4 e_1-3 d_1 e_4\right) \beta _{\it y}^{\prime \prime }\left(\theta _0\right)-3 d_1 d_4 N_{\it y}^{\prime \prime }\left(\theta _0\right)+d_1 e_8, \end{eqnarray} \tag{ A39 }$$

$$\begin{eqnarray} e_{14} & = & \left(2 d_7 e_1+d_1 \left(4 d_4 e_4-d_1 e_6\right)\right) \beta _{\it y}^{\prime } \left(\theta _0\right), \end{eqnarray} \tag{ A40 }$$

$$\begin{eqnarray} f_1 & = & N_{\it z}\left(\theta _0\right) \beta _{\it x}^{\prime }\left(\theta _0\right)-N_{\it x} \left(\theta _0\right) \beta _{\it z}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A41 }$$

$$\begin{eqnarray} f_2 & = & N_{\it z}\left(\theta _0\right) \beta _{\it x}^{\prime \prime }\left(\theta _0\right)- N_{\it x}\left(\theta _0\right) \beta _{\it z}^{\prime \prime }\left(\theta _0\right), \end{eqnarray} \tag{ A42 }$$

$$\begin{eqnarray} f_3 & = & N_{\it z}^{\prime }\left(\theta _0\right) \beta _{\it x}^{\prime }\left(\theta _0\right)-N_{\it x}^{\prime } \left(\theta _0\right) \beta _{\it z}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A43 }$$

$$\begin{eqnarray} f_4 & = & N_{\it z}^{\prime }\left(\theta _0\right) \beta _{\it x}^{\prime \prime }\left(\theta _0\right)-N_{\it x}^{\prime } \left(\theta _0\right) \beta _{\it z}^{\prime \prime }\left(\theta _0\right), \end{eqnarray} \tag{ A44 }$$

$$\begin{eqnarray} f_5 & = & N_{\it z}^{\prime \prime }\left(\theta _0\right) \beta _{\it x}^{\prime }\left(\theta _0\right)-N_{\it x}^{\prime \prime } \left(\theta _0\right) \beta _{\it z}^{\prime }\left(\theta _0\right), \end{eqnarray} \tag{ A45 }$$

$$\begin{eqnarray} f_6 & = & N_{\it z}\left(\theta _0\right) \beta _{\it x}{}^{(3)}\left(\theta _0\right)- N_{\it x}\left(\theta _0\right) \beta _{\it z}{}^{(3)}\left(\theta _0\right), \end{eqnarray} \tag{ A46 }$$

$$\begin{eqnarray} f_7 & = & 24 d_4^3 d_8-6 d_1 \left(d_7 \left(f_1-d_9\right)+3 d_4 d_5 d_8\right). \end{eqnarray} \tag{ A47 }$$

The expression of β={β_x,  β_y,  β_z} and the derivatives $\mbox{\boldmath $\beta $} ^{\prime },\,\mbox{\boldmath $\beta $},^{\prime \prime }\, \mbox{\boldmath $\beta $}^{(3)},$ and β⁽⁴⁾, which respectively represent the first, second, third, and fourth order differentiations with respect to θ evaluated at θ₀, are as follows:

$$\begin{eqnarray} \mbox{\boldmath $\beta $} (\theta _0)&=& \lbrace \kappa (h_1 \sin (\Gamma (\theta _0))+h_3 \cos (\Gamma (\theta _0))), \kappa (h_2 \sin (\Gamma (\theta _0))+h_4 \cos (\Gamma (\theta _0))),\kappa (\cos \alpha \, \cos (\Gamma (\theta _0))-h_5 \sin (\Gamma (\theta _0))) \rbrace, \end{eqnarray} \tag{ A48 }$$

$$\begin{eqnarray} &&\mbox{\boldmath $\beta $} ^{\prime }(\theta _0)\,{=}\,\lbrace \kappa \Gamma ^{\prime }(\theta _0) (h_1 \cos (\Gamma (\theta _0))\,{-}\,h_6 \sin (\Gamma (\theta _0))), \kappa \Gamma ^{\prime } (\theta _0) (h_2 \cos (\Gamma (\theta _0))\,{-}\, h_4 \sin (\Gamma (\theta _0))), \kappa \Gamma ^{\prime } (\theta _0) (-h_5 \cos (\Gamma (\theta _0))\,{-}\, \cos \alpha \, \sin (\Gamma (\theta _0))) \rbrace,\nonumber\\ \end{eqnarray} \tag{ A49 }$$

$$\begin{eqnarray} \mbox{\boldmath $\beta $} ^{\prime \prime }(\theta _0)&=& \lbrace \kappa (\Gamma ^{\prime }(\theta _0){}^2 (h_8 \sin (\Gamma (\theta _0))-\cos \phi ^{\prime }\, (h_9 \sin (\Gamma (\theta _0))+\sin \alpha \, \cos (\Gamma (\theta _0))))+ \Gamma ^{\prime \prime }(\theta _0) (h_1 \cos (\Gamma (\theta _0))- h_6 \sin (\Gamma (\theta _0)))), \nonumber \\ &&\kappa (\Gamma ^{\prime \prime }(\theta _0) (h_2 \cos (\Gamma (\theta _0))- h_4 \sin (\Gamma (\theta _0))) +\Gamma ^{\prime }(\theta _0){}^2 (h_2 (-\sin (\Gamma (\theta _0)))-h_4 \cos (\Gamma (\theta _0)))),\nonumber \\ && \kappa (\Gamma ^{\prime }(\theta _0){}^2 (h_5 \sin (\Gamma (\theta _0))-\cos \alpha \, \cos (\Gamma (\theta _0)))-\Gamma ^{\prime \prime }(\theta _0) (h_5 \cos (\Gamma (\theta _0))+\cos \alpha \, \sin (\Gamma (\theta _0)))) \rbrace, \end{eqnarray} \tag{ A50 }$$

$$\begin{eqnarray} \mbox{\boldmath $\beta $} ^{(3)}(\theta _0)&=& \lbrace \kappa (\Gamma ^{(3)} (\theta _0) (h_1 \cos (\Gamma (\theta _0))- h_6 \sin (\Gamma (\theta _0)))+ \Gamma ^{\prime }(\theta _0){}^3 (h_6 \sin (\Gamma (\theta _0)- h_9\cos (\Gamma (\theta _0)) \cos \phi ^{\prime })+h_8 \cos (\Gamma (\theta _0)))\nonumber \\ &&-3 \Gamma ^{\prime }(\theta _0) \Gamma ^{\prime \prime }(\theta _0) (h_1 \sin (\Gamma (\theta _0))+ h_3 \cos (\Gamma (\theta _0)))), \kappa (\Gamma ^{(3)} (\theta _0) (h_2 \cos (\Gamma (\theta _0))- h_4 \sin (\Gamma (\theta _0)))+\Gamma ^{\prime }(\theta _0) {}^3 (h_4 \sin (\Gamma (\theta _0))\nonumber\\ &&-h_2 \cos (\Gamma (\theta _0)))-3 \Gamma ^{\prime }(\theta _0) \Gamma ^{\prime \prime }(\theta _0) (h_2 \sin (\Gamma (\theta _0))+h_4 \cos (\Gamma (\theta _0)))), \kappa (-\Gamma ^{(3)} (\theta _0) (h_5 \cos (\Gamma (\theta _0))+ \cos \alpha \, \sin (\Gamma (\theta _0)))\nonumber \\ &&+ \Gamma ^{\prime }(\theta _0){}^3 (h_5 \cos (\Gamma (\theta _0))+\cos \alpha \, \sin (\Gamma (\theta _0)))+3 \Gamma ^{\prime }(\theta _0) \Gamma ^{\prime \prime }(\theta _0) (h_5 \sin (\Gamma (\theta _0))-\cos \alpha \, \cos (\Gamma (\theta _0)))) \rbrace, \end{eqnarray} \tag{ A51 }$$

$$\begin{eqnarray} \mbox{\boldmath $\beta $} ^{(4)}(\theta _0) &=& \lbrace \kappa (h_9 \cos \phi ^{\prime }\, (\sin (\Gamma (\theta _0)) (-3 \Gamma ^{\prime \prime }(\theta _0){}^2+\Gamma ^{\prime }(\theta _0){}^4-4 \Gamma ^{(3)}(\theta _0) \Gamma ^{\prime }(\theta _0))+ \cos (\Gamma (\theta _0)) (\Gamma ^{(4)}(\theta _0)-6 \Gamma ^{\prime }(\theta _0){}^2 \Gamma ^{\prime \prime }(\theta _0)))\nonumber \\ &&+h_8 (\sin (\Gamma (\theta _0)) (3 \Gamma ^{\prime \prime }(\theta _0){}^2 -\Gamma ^{\prime }(\theta _0){}^4+4 \Gamma ^{(3)}(\theta _0) \Gamma ^{\prime }(\theta _0))+\cos (\Gamma (\theta _0)) (6 \Gamma ^{\prime }(\theta _0){}^2 \Gamma ^{\prime \prime }(\theta _0)- \Gamma ^{(4)}(\theta _0)))+h_6 (\sin (\Gamma (\theta _0)) (6 \Gamma ^{\prime }(\theta _0){}^2\Gamma ^{\prime \prime }(\theta _0) \nonumber \\ &&-\Gamma ^{(4)}(\theta _0))+ \cos (\Gamma (\theta _0)) (-3 \Gamma ^{\prime \prime }(\theta _0){}^2+\Gamma ^{\prime }(\theta _0){}^4-4 \Gamma ^{(3)}(\theta _0) \Gamma ^{\prime }(\theta _0)))), \kappa (h_2 \sin (\Gamma (\theta _0)) (-3 \Gamma ^{\prime \prime }(\theta _0){}^2+\Gamma ^{\prime }(\theta _0){}^4-4\Gamma ^{(3)}(\theta _0) \Gamma ^{\prime }(\theta _0))\nonumber \\ && +h_2 \cos (\Gamma (\theta _0)) (\Gamma ^{(4)}(\theta _0)-6 \Gamma ^{\prime }(\theta _0){}^2 \Gamma ^{\prime \prime }(\theta _0))+ h_4 (\sin (\Gamma (\theta _0)) (6 \Gamma ^{\prime }(\theta _0){}^2 \Gamma ^{\prime \prime }(\theta _0)-\Gamma ^{(4)}(\theta _0))+ \cos (\Gamma (\theta _0)) (-3 \Gamma ^{\prime \prime }(\theta _0){}^2+\Gamma ^{\prime }(\theta _0){}^4\nonumber \\ && -4 \Gamma ^{(3)}(\theta _0) \Gamma ^{\prime }(\theta _0)))), \kappa (\cos \alpha \, (\sin (\Gamma (\theta _0)) (6 \Gamma ^{\prime }(\theta _0){}^2 \Gamma ^{\prime \prime }(\theta _0)-\Gamma ^{(4)}(\theta _0))+ \cos (\Gamma (\theta _0)) (-3 \Gamma ^{\prime \prime }(\theta _0){}^2+\Gamma ^{\prime }(\theta _0){}^4-4 \Gamma ^{(3)}(\theta _0)\Gamma ^{\prime }(\theta _0)))\nonumber \\ && -h_5 (\sin (\Gamma (\theta _0)) (-3 \Gamma ^{\prime \prime }(\theta _0){}^2+\Gamma ^{\prime }(\theta _0){}^4-4 \Gamma ^{(3)}(\theta _0) \Gamma ^{\prime }(\theta _0))+\cos (\Gamma (\theta _0)) (\Gamma ^{(4)}(\theta _0)-6 \Gamma ^{\prime }(\theta _0){}^2 \Gamma ^{\prime \prime }(\theta _0)))) \rbrace, \end{eqnarray} \tag{ A52 }$$

where

$$\begin{eqnarray} h_1 & = & \cos \alpha \, \cos \phi \, \cos \phi ^{\prime }\,-\sin \phi \, \sin \phi ^{\prime }, \end{eqnarray} \tag{ A53 }$$

$$\begin{eqnarray} h_2 & = & \cos \alpha \, \cos \phi \, \sin \phi ^{\prime }\,+\sin \phi \, \cos \phi ^{\prime }\,, \end{eqnarray} \tag{ A54 }$$

$$\begin{eqnarray} h_3 & = & \sin \alpha \, \cos \phi ^{\prime }, \end{eqnarray} \tag{ A55 }$$

$$\begin{eqnarray} h_4 & = & \sin \alpha \, \sin \phi ^{\prime }, \end{eqnarray} \tag{ A56 }$$

$$\begin{eqnarray} h_5 & = & \sin \alpha \, \cos \phi, \end{eqnarray} \tag{ A57 }$$

$$\begin{eqnarray} h_6 & = & \sin \alpha \, \cos \phi ^{\prime }, \end{eqnarray} \tag{ A58 }$$

$$\begin{eqnarray} h_7 & = & \sin \alpha \, \sin \phi ^{\prime }, \end{eqnarray} \tag{ A59 }$$

$$\begin{eqnarray} h_8 & = & \sin \phi \, \sin \phi ^{\prime }, \end{eqnarray} \tag{ A60 }$$

$$\begin{eqnarray} h_9 & = & \cos \alpha \, \cos \phi, \end{eqnarray} \tag{ A61 }$$

$$\begin{eqnarray} \cos (\Gamma (\theta _0)) &=& \frac{1+3 \cos \left(2 \theta _0\right)}{\sqrt{10+6 \cos \left(2 \theta _0\right)}}, \end{eqnarray} \tag{ A62 }$$

$$\begin{eqnarray} \sin (\Gamma (\theta _0))&=& \frac{3 \sin \left(2 \theta _0\right)}{\sqrt{10+6 \cos \left(2 \theta _0\right)}}, \end{eqnarray} \tag{ A63 }$$

$$\begin{eqnarray} \Gamma \left(\theta _0\right)&=& \tan ^{-1}\left(\frac{3 \sin \left(2 \theta _0\right)}{1+3 \cos \left(2 \theta _0\right)}\right), \end{eqnarray} \tag{ A64 }$$

$$\begin{eqnarray} \Gamma ^{\prime }\left(\theta _0\right)&=& \frac{4}{5+3 \cos \left(2 \theta _0\right)}+1, \end{eqnarray} \tag{ A65 }$$

$$\begin{eqnarray} \Gamma ^{\prime \prime }\left(\theta _0\right)&=& \frac{24 \sin \left(2 \theta _0\right)}{\left(5+3 \cos \left(2 \theta _0\right)\right){}^2}, \end{eqnarray} \tag{ A66 }$$

$$\begin{eqnarray} \Gamma ^{(3)}\left(\theta _0\right)&=& \frac{24(10 \cos (2 \theta _0)-3 \cos (4 \theta _0)+9)}{\left(5+3 \cos \left(2 \theta _0\right)\right){}^3}, \end{eqnarray} \tag{ A67 }$$

$$\begin{eqnarray} \Gamma ^{(4)}\left(\theta _0\right)&=& \frac{24(107 \sin (2 \theta _0)+120 \sin (4 \theta _0)-9 \sin (6 \theta _0))}{\left(5+3 \cos \left(2 \theta _0\right)\right){}^4}. \end{eqnarray} \tag{ A68 }$$

Note that ϕ', which appears in the above equations, is just a variable for the rotation phase, and the prime (') on it does not represent any differentiation.

By having known the derivatives of β from Equations (A48)–(A52), we can define the derivatives N evaluated at θ₀:

$$\begin{eqnarray} \mbox{\boldmath $N$}(\theta _0)&=&\mbox{\boldmath $\beta $}^{\prime }(\theta _0), \end{eqnarray} \tag{ A69 }$$

$$\begin{eqnarray} \mbox{\boldmath $N$}^{\prime }(\theta _0)&=&\mbox{\boldmath $\beta $}^{\prime \prime }(\theta _0), \end{eqnarray} \tag{ A70 }$$

$$\begin{eqnarray} \mbox{\boldmath $N$}^{\prime \prime }(\theta _0)&=&\mbox{\boldmath $\beta $}^{(3)}(\theta _0), \end{eqnarray} \tag{ A71 }$$

$$\begin{eqnarray} \mbox{\boldmath $N$}^{(3)}(\theta _0)&=&\mbox{\boldmath $\beta $}^{(4)}(\theta _0). \end{eqnarray} \tag{ A72 }$$

A.2. The Series Expansion Coefficients of Equation (24)

$$\begin{eqnarray} c_0 &=& g_1 (2 g_3 \sin ^3\left(\theta _0\right)-3 g_5 \cos \left(\theta _0\right)+2 g_2+\sqrt{3} \left(\log (2)-2 \log \left(g_4\right)\right)), \end{eqnarray} \tag{ A73 }$$

$$\begin{eqnarray} c_1&=& 3 g_1 (\sin (\theta _0) (g_3 \sin (2 \theta _0)+2 \sqrt{10+6 \cos (2 \theta _0)})+ \kappa \cos (\Gamma _0) (\sin (\theta _0)-3 \sin (3 \theta _0))), \end{eqnarray} \tag{ A74 }$$

$$\begin{eqnarray} c_2&=& \frac{3 g_1 \left(g_6 \left(1+3 \cos \left(2 \theta _0\right)\right)+ 2 \kappa \cos \left(\Gamma _0\right) \left(5+3 \cos \left(2 \theta _0\right)\right) \left(\cos \left(\theta _0\right)-9 \cos \left(3 \theta _0\right)\right)\right)}{4 \left(5+3 \cos \left(2 \theta _0\right)\right)}, \end{eqnarray} \tag{ A75 }$$

$$\begin{eqnarray} c_3 &=& \frac{1}{4} g_1 \left(-g_3 \cos \left(\theta _0\right)+g_7-\frac{4 \sqrt{2} \left(28 \sin \left(\theta _0\right)+9 \left(5 \sin \left(3 \theta _0\right)+\sin \left(5 \theta _0\right)\right)\right)}{\left(5+3 \cos \left(2 \theta _0\right)\right){}^{3/2}}\right), \end{eqnarray} \tag{ A76 }$$

where

$$\begin{eqnarray} \cos \left(\Gamma _0\right) & = & \sin \alpha \sin \zeta \cos \phi ^{\prime }\,+\cos \alpha \, \cos \zeta, \end{eqnarray} \tag{ A77 }$$

$$\begin{eqnarray} g_1 & = & \frac{\omega \, r_e}{12 c \kappa }, \end{eqnarray} \tag{ A78 }$$

$$\begin{eqnarray} g_2 & = & 6+\sqrt{3} \log (2+\sqrt{3}), \end{eqnarray} \tag{ A79 }$$

$$\begin{eqnarray} g_3 & = & 6 \kappa (\sin \zeta \, (\sin \phi \, \sin \phi ^{\prime }\,-\cos \alpha \, \cos \phi \, \cos \phi ^{\prime }\,)+\sin \alpha \, \cos \zeta \, \cos \phi \,), \end{eqnarray} \tag{ A80 }$$

$$\begin{eqnarray} g_4 & = & \sqrt{6} \cos \left(\theta _0\right)+\sqrt{5+3 \cos \left(2 \theta _0\right)}, \end{eqnarray} \tag{ A81 }$$

$$\begin{eqnarray} g_5 & = & 4 \kappa \cos \left(\Gamma _0\right) \sin ^2\left(\theta _0\right)+ \sqrt{10+6 \cos \left(2 \theta _0\right)}, \end{eqnarray} \tag{ A82 }$$

$$\begin{eqnarray} g_6 & = & 7 g_3 \sin \left(\theta _0\right)+3 g_3 \sin \left(3 \theta _0\right)+ 8 \sqrt{10+6 \cos \left(2 \theta _0\right)} \cos \left(\theta _0\right), \end{eqnarray} \tag{ A83 }$$

$$\begin{eqnarray} g_7 & = & 9 g_3 \cos \left(3 \theta _0\right)-2 \kappa \cos \left(\Gamma _0\right) \left(\sin \left(\theta _0\right)-27 \sin \left(3 \theta _0\right)\right). \end{eqnarray} \tag{ A84 }$$

B.1. To Find the Solution to the Integrals in Equation (25)

Consider the integral

$$\begin{equation} S_0 = \int _{-\infty }^{+\infty } {e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3}\,\mu ^3} d\mu. \end{equation} \tag{ B1 }$$

By changing the variable of integration μ = (x/l) + m and defining the constants $l = \root 3 \of {c_3}$ and m = −c₂/(3c₃), we obtain

$$\begin{equation} \int _{-\infty }^{+\infty } {e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3} \,\mu ^3)}d\mu = U \int _{-\infty }^{+\infty } e^{i \left(z\, x +x^3\right)}dx, \end{equation} \tag{ B2 }$$

where

$$\begin{equation*} z = \frac{1}{\root 3 \of {c_3}}\left(c_1-\frac{c_2^2}{3 c_3}\right), \quad U = \frac{1}{\root 3 \of {c_3}} e^{i \frac{c_2}{3 c_3}(\frac{2 c_2^2}{9 c_3}-c_1)}. \end{equation*}$$

For Im(z) = 0, we know

$$\begin{equation} j_0 = \int _{-\infty }^{\infty }{ e}^{i(z x+x^3)} dx = \frac{\pi }{\root 3 \of {3}} \left[\left(1-\frac{\sqrt{z^2}}{z}\right) \mbox{Ai}\left(-\frac{\sqrt{z^2}}{\root 3 \of {3}}\right)+\left(1+ \frac{\sqrt{z^2}}{z}\right) \mbox{Ai}\left(\frac{\sqrt{z^2}}{\root 3 \of {3}}\right)\right], \end{equation} \tag{ B3 }$$

where Ai(z) is an entire Airy function of z with no branch cut discontinuities, and

$$\begin{equation} j_1 = \int _{-\infty }^{\infty } x\,{e}^{i(z x+x^3)}dx = -i\frac{2 \pi }{\root 3 \of {3^2}} \mbox{Ai}^{\prime }\left(\frac{z}{\root 3 \of {3}}\right), \end{equation} \tag{ B4 }$$

where Ai'(z) is the derivative of the Airy function Ai(z). Therefore, we have

$$\begin{equation} S_0 = U j_0. \end{equation} \tag{ B5 }$$

By differentiating Equation (B2) on both sides with respect to c₁ we obtain

$$\begin{eqnarray} S_1 = \int _{-\infty }^{+\infty }\mu {e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3} \,\mu ^3)}d\mu & = & \frac{U}{\root 3 \of {c_3}}\int _{-\infty }^{+\infty } \left(x-\frac{c_2}{3\root 3 \of {c_3^2}}\right) e^{i \left(z\, x +x^3\right)}dx = \frac{U}{\root 3 \of {c_3}}\left(j_1-\frac{c_2}{3\root 3 \of {c_3^2}}j_0\right). \end{eqnarray} \tag{ B6 }$$

Differentiation of Equation (B2) on both sides with respect to c₂ gives

$$\begin{eqnarray} S_2 = \int _{-\infty }^{+\infty }\mu ^2 {e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3} \,\mu ^3)}d\mu & = & \frac{U}{3c_3} \int _{-\infty }^{+\infty }\left(\frac{2 c_2^2}{3c_3}- c_1-\frac{2c_2}{\root 3 \of {c_3}}x\right) e^{i \left(z\, x +x^3\right)}dx = \frac{U}{3c_3}\left[\left(\frac{2 c_2^2}{3c_3}-c_1\right)j_0- \frac{2c_2}{\root 3 \of {c_3}}j_1\right]. \end{eqnarray} \tag{ B7 }$$

Next, by differentiating Equation (B2) on both sides with respect to c₃ we obtain

$$\begin{eqnarray} S_3 = \int _{-\infty }^{+\infty }\mu ^3 {e}^{i(c_{\rm 1} \,\mu +c_{\rm 2} \,\mu ^2+c_{\rm 3} \,\mu ^3)}d\mu & = & \frac{U}{9\root 3 \of {c_3^7}}\int _{-\infty }^{+\infty } \left(\frac{9 c_1 c_2 c_3-4 c_2^3+i 9 c_3^2}{3 \root 3 \of {c_3^2}}+ \big(4 c_2^2-3 c_1 c_3\big)x \right) e^{i \left(z\, x +x^3\right)}dx \nonumber \\ & = & \frac{U}{9\root 3 \of {c_3^7}}\left[\frac{\big(9 c_1 c_2 c_3-4 c_2^3+i 9 c_3^2\big)}{3 \root 3 \of {c_3^2}}j_0+ \big(4 c_2^2-3 c_1 c_3\big)j_1 \right]. \end{eqnarray} \tag{ B8 }$$