Free Access
Issue
A&A
Volume 646, February 2021
Article Number A23
Number of page(s) 12
Section Stellar structure and evolution
DOI https://doi.org/10.1051/0004-6361/202039470
Published online 02 February 2021

© ESO 2021

1. Introduction

Accretion-powered millisecond pulsars (AMPs) are rapidly rotating neutron stars (NSs) located in low-mass X-ray binaries. These weakly magnetised NSs have been spun-up by the accretion torques acting on the NS during the process of gas accretion from its companion (Radhakrishnan & Srinivasan 1982; Alpar et al. 1982). The accreting matter falls on the NS surface, creating hotspots at the magnetic poles that emit X-ray radiation. Rotation of the NS produces pulsation and the pulse shape contains information about the parameters of the NS such as its mass and radius (see e.g. Poutanen & Gierliński 2003; Watts et al. 2016), which can be used to constrain the equation of state (EOS) for extremely dense matter of inner parts of the NS (for instance, see Lindblom 1992; Lattimer 2012; Raaijmakers et al. 2019).

Previously, pulse profiles have been modelled in order to constrain the NS parameters (see e.g. Pechenick et al. 1983; Miller & Lamb 1998; Poutanen & Gierliński 2003; Poutanen & Beloborodov 2006; Morsink et al. 2007; Lo et al. 2013; Salmi et al. 2018; Bogdanov et al. 2019; Miller et al. 2019; Riley et al. 2019). Usually, this approach has a problem to uniquely determine the pulsar geometry, particularly the observer inclination and the co-latitude of the emitting spot, which are highly degenerate with the NS mass and radius. However, the phase dependence of the polarization angle (PA), and therefore Stokes parameters Q and U, is a powerful tool to constrain the NS geometry (Viironen & Poutanen 2004; Poutanen 2010). The transformation of Stokes parameters from the NS surface to the observer frame was studied by Viironen & Poutanen (2004) and Poutanen (2020) for the case of a rapidly rotating spherical NS, and also for the case of an oblate star by Loktev et al. (2020).

The radiation escaping the NS surface in AMPs is expected to be significantly polarized, as it is scattered by the hot electrons in the accretion shock above the NS surface. In rotation-powered millisecond pulsars (RMPs), on the other hand, the dominant thermal radiation is too soft and likely weakly polarized; only a possible tail (Salmi et al. 2020) above a few keV may have a significant polarization, but it is too dim to be detected by the upcoming X-ray polarization instruments. The polarization degree (PD) produced by electron scattering strongly depends on the scattering angle and electron temperature (see e.g. Nagirner & Poutanen 1994; Poutanen 1994a,b). A model for polarized radiation from AMPs, based on Comptonization in an optically thin NS atmosphere, but in Thomson scattering approximation, was introduced in Viironen & Poutanen (2004; see also Sunyaev & Titarchuk 1985). For a more accurate model, one should apply the formalism for Compton scattering in a hot slab (Poutanen & Svensson 1996); this is, however, computationally more expensive. We used a slightly modified version of the Thomson model (see Sect. 2.1.1) to make the first detailed simulations of the upcoming X-ray polarization observations and NS parameter constraints, which are expected to be measured soon. The results can be applied when designing the observational strategy of the Imaging X-ray Polarimeter Explorer (IXPE; Weisskopf et al. 2016) or other future X-ray polarimetric missions such as the enhanced X-ray Timing and Polarimetry mission (eXTP; Zhang et al. 2019; Watts et al. 2019).

For simplification, the case with one hotspot was applied in most of the calculations. This is justified because most of the AMPs expected to be observed by IXPE show simple sine-like profiles during the peak of their outbursts (Kirsch et al. 2004; Hartman et al. 2008; Leahy et al. 2009, 2011; Sanna et al. 2017). This fact is naturally explained if only one spot is visible to the observer while the other is hidden by the accretion disc (Ibragimov & Poutanen 2009; Poutanen et al. 2009). However, observations of SAX J1808.4−3658 show that the pulse profile starts to significantly deviate from the sine-wave at low fluxes, where the magnetosphere expands sufficiently so that we can see the secondary spot. Also, few other AMPs show signs of a secondary spot (see e.g. Galloway et al. 2002; Sanna et al. 2018) and understanding the effects from multiple spots would also be useful for possible future X-ray polarization measurements of RMPs. Therefore, we also consider a few cases with two possibly non-antipodal spots, as suggested by the non-dipolar magnetic field configuration indicated by the observations of PSR J0030+0451 by the Neutron star Interior Composition ExploreR (Bilous et al. 2019; Miller et al. 2019; Riley et al. 2019).

The remainder of this paper is structured as follows. In Sect. 2.1, we present the methods used to model polarized radiation from an AMP. In Sect. 2.2, we explain how we simulated the data that could be observed by IXPE, and in Sect. 2.3 we describe the Bayesian method used to obtain constraints on the NS model parameters. The produced simulated data and the predicted NS parameter constraints are shown in Sect. 3. A discussion and conclusions appear in Sects. 4 and 5, respectively.

2. Methods

2.1. Modelling polarization profiles

Our polarization modelling is mostly based on the polarization model introduced in Viironen & Poutanen (2004) and on the X-ray pulse shape modelling introduced in Poutanen & Beloborodov (2006). We used the ‘oblate Schwarzschild’ approximation, which takes the deformed shape of the star into account in addition to the special and general relativistic corrections to the photon angles and trajectories (Morsink et al. 2007; Miller & Lamb 2015; Salmi et al. 2018; Suleimanov et al. 2020). The shape of the NS was obtained from the model presented by AlGendy & Morsink (2014), which is suitable for the spin frequencies considered here (see Silva et al. 2020 for a recent model for the most rapidly rotating stars). The description of polarized pulse formation has been slightly revised since Viironen & Poutanen (2004), and it accounts for the effects of the oblate shape of the star on the observed PA (Loktev et al. 2020). Otherwise, the Stokes parameters were computed in the frame comoving with the spot and then transformed to the observer frame as in Viironen & Poutanen (2004) and Poutanen (2020). We note that the oblateness and relativistic effects are significant for stars spinning at a rate higher than 200 Hz, and assuming an incorrect shape may bias the constraints on NS geometry, as shown by Loktev et al. (2020).

We assumed that the observed photons originate in one or two spots at the NS surface, which we modelled as a slab of hot electrons above a blackbody emitting surface. The spectral energy distribution of the radiation was calculated as in Salmi et al. (2018) using the Comptonization model SIMPL (Steiner et al. 2009) from the XSPEC package (Arnaud 1996), which is based on the solution of the non-relativistic Kompaneets equation (Sunyaev & Titarchuk 1980). The parameters of the model are the photon spectral index Γ and a fraction Xsc of black-body photons scattered in the slab. We used the version SIMPL-2, which takes both Compton up-scattered and down-scattered photons into account. The model converts the specified fraction of black-body photons (with temperature Tbb) into a Comptonized power-law-like spectrum.

The values of all the parameters in our fiducial model are shown in Table 1. Especially, the geometrical parameters are the observer inclination i (angle between the spin axis and the line of sight), co-latitude of the primary spot θ (angle between the spin axis and the radius vector of the spot centre), pulsar rotation axis position angle χ (angle measured from the north counterclockwise to the projection of the rotation axis on the plane of the sky), and an arbitrary phase shift Δϕ. When simulating the data and calculating parameter constraints in Sect. 3, we also consider several sets of modified initial parameters, which are shown in Table 2. In the case of a two-spot model, we additionally considered the co-latitude of a secondary spot θ12 and the longitude difference between the two spots ϕ12.

Table 1.

Fiducial parameters of the synthetic data.

Table 2.

Parameters altered from the fiducial model for all of the computed models.

2.1.1. Emission model

The dependence of the PD on the zenith angle (or its cosine μ) and photon energy E (in the spot frame) is still not well-known for the AMPs. Here, we used a slightly modified version of the simple Thomson scattering model in a plane-parallel atmosphere presented in Viironen & Poutanen (2004). In the co-moving frame of the spot, the PD P(E, μ) as a function of photon energy E and cosine of the zenith angle μ was obtained from the ratio of Stokes parameters Q and I (Stokes U is zero due to the assumed azimuthal symmetry of the radiation) for photons scattered multiple times in a slab of Thomson optical depth τ = 1 and electron gas temperature Te = 50 keV. The energy of the n-times scattered photon can be found from (Rybicki & Lightman 1979)

(1)

where E0 is the energy of a seed photon and A is the amplification factor for each scattering.

We assumed the seed photons to have an energy distribution described by the Planck function of temperature Tbb = 1 keV with the angular distribution given by the function as(μ):

(2)

We considered two cases: (1) the isotropic intensity with the constant angular function as(μ)  =  1 and (2) the angular distribution corresponding to the electron-scattering dominated semi-infinite atmosphere (Chandrasekhar & Breen 1947; Chandrasekhar 1960; Sobolev 1949, 1963):

(3)

In the first case, we assumed unpolarized emission Ps(μ) = 0. In the second case, we took (Viironen & Poutanen 2004)

(4)

where pmax = 11.71%. For weakly magnetised NSs, the intensity of the radiation in the spot comoving frame is expected to be independent of azimuth. Therefore, the polarization of the radiation (i.e. dominant direction of electric field oscillations) is either in the meridional plane, containing the normal to the surface and the line-of-sight, or perpendicular to that plane. We define the polarization as positive in the first case and negative in the second.

After leaving the slab, the unscattered photons have the following angular dependence of intensity:

(5)

where indices l, r (as well as the ± sign) refer to the intensity of radiation polarized in the meridional plane and perpendicular to it, respectively. The angular dependence of the total intensity (also known as the beaming function) for escaping unscattered photons is .

The energy dependence of the radiation escaping from the slab was accounted for by assuming that the dependence of intensity on μ and E can be separated:

(6)

The angular dependencies of the intensities for remaining scattering orders were calculated based on the intensities of the zeroth scattering order using the formulae presented in Viironen & Poutanen (2004). The energy function was approximated with a diluted Planck function assuming that in each scattering, a photon increases its energy by a factor of A:

(7)

with ϵ0(E) = BE(Tbb). Such a relation implies that the number of photons ∫(ϵn(E)/E)dE does not depend on n. In reality, of course, the number of scattered photons decreases with n, but this was accounted for in the normalization of the angular functions . Thus, the corresponding Stokes parameters in this approximate Thomson scattering model IT and QT are given by:

(8)

For a given energy E, we obtained the Stokes parameters by summing Eq. (8) over n as

(9)

where we used N = 23 for the highest number of scatterings. The final PD is then

(10)

We emphasise that these Stokes parameters were only used when evaluating the PD as a function of energy and the emission angle, and they do not provide the final energy dependence of the radiation.

We then used the aforementioned Thomson model to obtain the final angular dependence of the radiation. The final intensity (in the co-moving frame) was obtained jointly with SIMPL from

(11)

where f0(μ) = 𝒩a0(μ) is the beaming for unscattered photons (normalized so that ), Ic(E) is the Comptonized part of the spectrum computed with SIMPL, and

(12)

where 𝒩E is a normalization factor determined so that . Also taking the final energy dependence of intensity directly from the Thomson model, instead of SIMPL that was used here, would not be useful because the Thomson model does not describe the observed spectra of AMPs well for most of the preferred model parameters (unlike SIMPL which we can parametrize, independently from the polarization model). This means that the Stokes parameters IT(μ, E) and QT(μ, E), which were used to estimate P(μ, E), differ from the final I(μ, E) and Q(μ, E) = P(μ, E) I(μ, E). However, this is not crucial, because the Thomson model is only expected to describe PD and the angular dependence of radiation with an accuracy that is good enough to obtain first order estimates for the geometrical parameters of the NS.

The modelled final angular and energy dependencies of intensity and PD (using the model parameters shown in Table 1) are shown in Figs. 1 and 2. From the left panel in Fig. 1, we see that the beaming pattern shows limb darkening for lowest energies (below about 10 keV) but limb brightening down to μ ≈ 0.3 for higher energies. The classical result of Chandrasekhar-Sobolev given by Eq. (3), which is shown with a black dot-dashed line, appears to be a reasonable approximation for the final angular dependence at the nominal 2–8 keV energy range of IXPE. On the other hand, the classical formula predicts a significantly smaller PD (in absolute value) than the optically thin Thomson model for the highest emission angles and for the most energies, as seen from the right panel of Fig. 1. We also note that assuming polarized seed photons (not shown in the figures), particularly having pmax = 11.71% in Eq. (4), leads to even a slightly higher absolute value for the final PD of the Thomson scattered photons. For the highest energies (only above 18 keV), a change in the sign of the PD is also observed. We note that this energy depends on the parameters of the spectrum (seed photon temperature, electron temperature) and it is smaller for smaller Tbb and Te. For the smallest energies, PD is almost independent of energy at a given emission angle, as is also seen in the right panel of Fig. 2. From the left panel of Fig. 2, we see that the emergent intensity spectrum is close to that obtained directly from SIMPL, although the spectrum emitted at highest zenith angles has a smaller blackbody component.

thumbnail Fig. 1.

Angular dependence of emitted intensity (left panel) and PD (right panel). The angular dependence was normalized so that . Black, blue, green, orange, and red solid lines show the model for photon energies of 2, 5, 8, 12, and 18 keV, respectively. The black dot-dashed curves show the intensity and polarization corresponding to the classical results of Chandrasekhar–Sobolev (see Eq. (4)), corresponding to the optically thick electron-scattering dominated atmosphere. The grey dot-dashed curve in the left panel shows the angular dependence of intensity for unscattered photons. The model was calculated with the parameters shown in Table 1 (e.g. pmax = 0).

thumbnail Fig. 2.

Energy dependence of emitted intensity (left panel) and PD (right panel). Black, blue, green, orange, red, and magenta solid lines show the model for μ= 0.0, 0.2, 0.4, 0.6, 0.8, and 1.0, respectively. The dot-dashed black curve (left) shows the intensity corresponding to the angle-independent result of model SIMPL, namely (1 − Xsc)BE(Tbb)+XscIc(E) (see Sect. 2.1.1). The model was calculated with the parameters shown in Table 1.

2.1.2. Pulse profiles

When calculating pulse profiles, PD P(μ, E) and the beaming of the Comptonized spectral component f(μ, E) were obtained by interpolating pre-computed tables. Following the methods of Viironen & Poutanen (2004) and Loktev et al. (2020) to transform the quantities from the spot frame to the observer frame, we obtained the Stokes flux vector (FI, , ) as a function of observed phase ϕ. The observed phase differs from the pulsar rotation phase ϕrot due to the different time delays for different sub-spots at different phases. However, unlike in Viironen & Poutanen (2004) and Loktev et al. (2020), we used the equatorial radius as the reference radius and zero as the reference impact parameter when computing the delays (as in Poutanen & Beloborodov 2006 and Salmi et al. 2018). The rotational phase of the pulsar was set to zero when the primary spot was closest to the observer.

The obtained Stokes parameters are defined in the polarization basis related to the projection of the pulsar rotation axis on the sky. Therefore, we also decided to take the position angle χ of the pulsar rotation axis into account by correcting the computed Stokes parameters using the following equations:

(13)

where and are the final phase-resolved model Stokes spectra. In our simulated data, we assumed that χ = 0, but we kept it as a free parameter when calculating parameter constraints. Finally, the observed PD was obtained as

(14)

We now consider the case with one or two antipodal spots. The modelled polarized pulses, in the case of model 8 presented in Table 2, are shown in Fig. 3 for three different observed photon energies (2, 5, and 8 keV). The contributions from primary and secondary spots are indicated with blue dashed and red dotted curves correspondingly (the blue curve represents the fiducial model shown in Table 1). We see only weak dependence of the observed PD Pobs and PA on the energy. During one rotational period, all pulse profiles show only one maximum because the second spot is significantly less visible compared to the primary one. For all the shown energies, the PD of each spot (shown in the middle row in Fig. 3) is highest when the spot is most inclined to the observer, but it is still visible. Also, the maximum PD is significantly higher if considering emission from two spots. At energies above 10 keV, the sign of P changes (see Fig. 2). This results in a jump in PA by 90° at certain phases. The pulsar phase bin where such a jump occurs depends on the energy. However, we only used the model in the range of 2–8 keV, where the current instrumental response of the IXPE detector is best validated.

thumbnail Fig. 3.

Theoretical pulse profiles of the observed flux, PD, and PA for two antipodal spots shown for three different energies (2, 5, and 8 keV). Solid black curves correspond to the total flux (or PD or PA), while the blue dashed curves are for the primary spot, and red dotted curves are for the secondary spot. The parameters of the model are given in Table 1.

We now present our model in the case of a non-antipodal spot configuration. We assumed that the co-latitudes of the spots are θ = 20° and θ2 = 120° with the difference in longitude between the spots being ϕ12 = 180° (model 9 presented in Table 2). We show the model in Fig. 4. The contribution from the primary spot is the same as in Fig. 3, but now the secondary spot is more visible to the observer. The observed PD is lower because photons reaching the observer are emitted to smaller angles relative to the spot normal. We also detect abrupt breaks in PA at those phases when the secondary spot appears and disappears because the non-antipodal spot has a largely different PA (as also discussed by Loktev et al. 2020).

thumbnail Fig. 4.

Same as Fig. 3, but for a non-antipodal spot configuration with θ = 20°, θ2 = 120°, and ϕ12 = 180°.

2.2. IXPE data simulation

We used the IXPE observation-simulation framework IXPEOBSSIM version 12.0.0 (Pesce-Rollins et al. 2019) to generate our synthetic data. The programme is designed to fold a complete source model (described in Sect. 2.1, in our case) with the current best estimate of the instrument response functions to produce simulated event files (also known as event or photon lists) for a given observation time. The event files are essentially identical in format to those that will be produced from real observations and they include all the relevant physical properties of the events, such as the measured arrival time, energy, sky-position, and, most notably, the photoelectron azimuthal angle ψk. The latter is key to the polarization measurement and is the basic ingredient for the analysis described in this section.

Starting from a given photon list, we built the so-called Stokes parameter spectra and . These are essentially weighted histograms of detector counts, binned in energy, where the weights are dictated by the reconstructed azimuthal angle of each event, qk = cos2ψk and uk = sin2ψk, according to the formalism described in Kislat et al. (2015). The statistical errors were propagated by summing the weights themselves in quadrature, according to the standard formalism of weighted histograms. From a practical standpoint, these data products were stored in FITS files whose format conforms to the standard OGIP PHA type I file format in order to be interoperable with the high-level analysis tools used by the X-ray community, such as XSPEC. The Stokes spectra were created separately for each of the 19 pulsar phase bins.

All the IXPE instrument response functions and associated binned products are defined between 1 and 12 keV in 275 equal-width energy bands. However, we only used channels in the nominal 2–8 keV IXPE energy band to avoid any complication from the complexity of the detector response outside this range. The effective area below 2 keV and above 8 keV is so small that this has effectively very little practical implications for the analysis.

The observed Stokes parameters (and their measured errors) were transformed to polarization estimates (directly comparable with the input model) by dividing them by the effective modulation factor of each observed energy channel, and (as explained in Kislat et al. 2015). In addition, we simulated the data in terms of the normalized Stokes profiles q = CQ/CI and u = CU/CI, where CI is the observed count spectrum (see simulated data presented in Sect. 3.1). These can be compared to the normalized Stokes parameters and , which can be predicted directly from the theoretical model, because the ratio of two Stokes parameters is not expected to be sensitive to the energy response of the detector.

When estimating the measurement errors and comparing the data to the model (in Sect. 3.1), we also combined observations from all the energy channels between 2 and 8 keV. The broadband values were computed as averages of q and u, weighted with the energy-dependent effective area according to the formalism in Kislat et al. (2015), and they were compared to the modelled broadband values of qm and um. When analysing the data using the Bayesian inference (in Sect. 3.2), we used the same technique in order to speed up the computation, and we created eight logarithmic energy bands between 2 and 8 keV, which were independently fitted.

2.3. Bayesian modelling

The obtained Stokes profiles were fitted using an affine invariant ensemble sampler with the EMCEE package of PYTHON (Foreman-Mackey et al. 2013). The posterior probability densities were calculated to the free parameters of our model. For most simulations, these were the observer inclination i, spot co-latitude θ, phase shift Δϕ, and the position angle χ of the pulsar spin axis. In the case of the model with a large spot size (ρ = 30°), we also kept the size of the spot ρ as a free parameter; furthermore, in the case of the two-spot models, we also kept the co-latitude θ2 of the secondary spot and the longitude difference ϕ12 between the spots as free parameters. The correct values of these and also the fixed model parameters are shown in Table 1 (for the fiducial model) and in Table 2 (for other models).

The synthetic normalized Stokes q and u data (created using IXPEOBSSIM as described in Sect. 2.2) were fitted against the modelled qm and um assuming that the probability densities of q and u are uncorrelated and normally distributed around qm and um with measured errors as the standard deviation. As mentioned in Sect. 2.2, by fitting normalized Stokes values, we improved the computational efficiency as we did not need to transform every model to observed data using a forward-folding approach (including the response of the detector). The likelihoods for a given model were calculated separately for 19 phases and eight energy bins of the simulated data (where the energy channels were combined as explained in Sect. 2.2) and summed to get the total likelihood. The prior probability distributions of all the sampled parameters were assumed to be uniform and the limits of the priors were set to (0°, 90°) in i and θ, (−90°, 90°) in χ, (−180°, 180°) in Δϕ, (0°, 60°) in ρ, (0°, 180°) in θ2, and (0°, 360°) in ϕ12.

3. Results

3.1. Simulated data

We now present the simulated data and estimated measurement accuracy obtained for the one-spot models with different exposure times and slightly different model parameters. The computations were done with exposures of 200, 400, 600, 800, and 1000 ks. In all of our simulations, we assumed a 100 mCrab source, corresponding to approximately SAX J1808.4−3658 during its peak luminosity. For other sources, which are less bright, a longer observation time is expected to be required to obtain similar results. On the other hand, a higher PD than our assumption (about a maximum 5% observed PD for our fiducial parameter set) should affect the simulated data in the same way as increasing the observation time.

The simulated broadband PD and Stokes q and u profiles for the shortest (200 ks) and longest (1000 ks) exposure times, in the case of one observation realisation of the fiducial parameter set, are shown in Fig. 5. We see that the produced data are not biased from the theoretical model shown in blue. The measurement errors are also significantly smaller in the case of a longer exposure time.

thumbnail Fig. 5.

Simulated PD and normalized Stokes q and u profiles for 200 ks (left panels) and 1 Ms (right panels) exposure times and the fiducial one-spot model when computing a weighted average over energy range 2–8 keV. The blue curve shows the theoretical model and the purple dots show the simulated observed data including the measurement errors.

The errors as a function of the observation time are examined in more detail, for all calculated one-spot models and exposure lengths, in Figs. 6 and 7, where we present the average minimum errors (absolute and relative, respectively) in q and u as a function of observation time t. From these figures, we see that increasing the observation time from 200 to 1000 ks improves the accuracy in the measured q and u from 30–90% to 15–40% depending strongly on the model parameters. Of course, the exact values for the errors also depend on the adopted number of phase bins and they are only relevant when comparing and analysing data with a similar setup.

thumbnail Fig. 6.

Average absolute minimum errors in the measured energy-integrated Stokes q and u as a function of observation time (for the one-spot models specified in Table 2). The curves were calculated for 200, 400, 600, 800, and 1000 ks (shown as purple dots) and interpolated linearly in a logarithmic scale for other values.

In any case, the presented errors scale as , which is as expected. The relative errors in Fig. 7 are smallest for those models that produce the highest PD. These are the models where the initial polarization is non-zero (blue dash-dotted curve), and where the inclination or co-latitude is higher than in the fiducial model (black and red dash-dotted curves, respectively). Larger errors are expected for the case with a larger spot size (solid orange curve) and smaller inclination or co-latitude (black and red dashed curves). On the contrary, absolute errors in Fig. 6 are almost identical in all models because they are only determined by the total amount of observed counts and the brightness of the source is fixed in every model.

thumbnail Fig. 7.

Average relative minimum errors in the measured energy-integrated Stokes q and u as a function of observation time for the one-spot models (calculated as ΔS/S = ([Δq/qm]min + [Δu/um]min)/2. The solid blue curve shows the errors with the fiducial parameters (shown in Table 1). The blue dash-dotted curve corresponds to similar data but with pmax = 0.1171, the black dashed curve with i = 50°, the black dash-dotted curve with i = 70°, the red dashed curve with θ = 10°, the red dash-dotted curve with θ = 30°, and the orange solid curve with ρ = 30° (also explained in Table 2). Curves were interpolated from the calculated points as in Fig. 6.

We also simulated polarization data detected assuming either two antipodal or non-antipodal spots (corresponding to the two-spot models shown in Figs. 3 and 4), but using only a single observation realisation with a 200 ks observation time. The results are shown in Fig. 8, where we see that the measured broadband PD and Stokes parameters for the antipodal case significantly differ from those of the non-antipodal case. Relative errors in the measured Stokes parameters are slightly higher for the non-antipodal case because the PD is smaller for the chosen model parameters.

thumbnail Fig. 8.

Simulated polarization data as in Fig. 5 (for 200 ks exposure time), but for the case of two antipodal spots with θ = 20° and θ2 = 160° (left panels), and two non-antipodal spots with θ = 20° and θ2 = 120° (right panels).

3.2. Parameter constraints

We now present the parameter constraints obtained when fitting synthetic data and start again with the one-spot models. The resulting posterior probability distributions for the model parameters are shown in Fig. 9 for the fiducial data with the shortest and longest observation times (200 and 1000 ks). Because a slight variation in the constraints was detected for different realisations of the same synthetic data, we show the combined posterior distributions for three different realisations of the data (including those shown in Fig. 5). From Fig. 9, we see that a minor bias in spot co-latitude θ appears with a 200 ks observation, but no biases arise in the case of a 1000 ks observation. All the free parameters, which are the observer inclination i, spot co-latitude θ, position angle χ, and the phase shift Δϕ, are clearly constrained better than their prior limits (given in Sect. 2.3), and they show no strong correlation with each other.

thumbnail Fig. 9.

Posterior probability distributions for NS parameters when combining fits from three realisations of synthetic data produced assuming a 200 ks (left panels) and a 1 Ms (right panels) exposure time and the fiducial model shown in Table 1. In the two-dimensional posterior distributions, the black colour shows a 68% and the grey colour shows a 95% highest posterior density credible region. In the one-dimensional posterior distributions, the red dashed lines show a 68% and the dark orange dashed lines show a 95% highest posterior density credible interval. The blue lines show the input values.

The results also appear similar in the case of the other observation lengths and calculated models, although with a varying magnitude of the credible intervals. One-dimensional posterior results averaged over the three different synthetic data generations of the fiducial model are also shown in Table 3. Furthermore, in the case of the simulations where the spot size was a free parameter (those where the correct value was ρ = 30°), no stringent constraints were found to the spot size.

Table 3.

Most probable values and 68% and 95% credible intervals for five different simulations averaged from three realisations of the synthetic data produced using the fiducial model.

In Fig. 10 we illustrate how the obtained average credible interval for i and θ depends on the observation time for all the one-spot models, however, when using only one realisation of the data for the other models than the fiducial one. For the fiducial model, we see that increasing the observation time from 200 to 1000 ks improves the width of the 95% probability interval from about 13° to 6°. We detect a slight variation between the constraints from different models and observation lengths, which can be explained by the statistical uncertainties in the generated data. However, ordinarily, the constraints become tighter when increasing the exposure time. We note that the presented constraints approximately scale as the relative errors shown in Fig. 7.

thumbnail Fig. 10.

Average width of the 95% credible interval for the measured i and θ as a function of the observing time. The curves, colours, and interpolation method correspond to those in Fig. 7, but the solid blue curve corresponds to the average constraints from the three different synthetic data generations (using the same model parameters). Purple error bars show the maximum deviation between the average and single data generation results. The other lines show constraints for only one realisation of the synthetic data (for each exposure time).

From Fig. 10 we also see that different assumptions of the model parameters can lead to largely different parameter constraints. For example, the constraints are tightest for the model where the seed photons of the Thomson scattering slab have a non-zero polarization (blue dash-dotted curve). On the other hand, having a larger spot size and having it as a free parameter (orange solid curve) indicates a less constrained NS geometry, which is as expected. The models where the observer inclination i is either lower or higher than in the fiducial model (black dashed and dash-dotted curves, respectively) perform as in Fig. 7, showing tighter constraints when the PD is higher (i is higher). On the other hand, models where the spot co-latitude θ is altered (red curves) exhibit no significant variation in the parameter limits.

We now present the parameter constraints when applying the models with two spots (see Table 2) in order to inspect how accurately the spot configuration could be constrained using the polarization data. The resulting posterior distributions are shown in Fig. 11 for one realisation of data and a 200 ks observation time (corresponding to the data shown in Fig. 8). We see that the observer inclination i and primary spot co-latitude θ are clearly constrained better in the case of antipodal spots, which is expected due to the higher observed PD for the chosen model parameters. In the antipodal case, a 95% credible interval is 57–65° for i and 13–30° for θ. In the non-antipodal case, a 95% credible interval is 33–63° for i and 14–39° for θ. We also obtained constraints for the secondary spot. In the antipodal case, we only got a lower limit for θ2, which is around 150° (the true value being 160°). The possibility of having a non-pulsating secondary spot exactly at the southern rotational pole can barely be excluded based on these data. On the contrary, if θ2 = 120° (as in our non-antipodal case), we get 106–127° for the 95% credible interval, and also ϕ12/2π is constrained between 0.3 and 0.7. Thus, the secondary spot is more accurately located in this case, even though the other parameters are less constrained.

thumbnail Fig. 11.

Posterior probability distributions for NS parameters similar to the left panel of Fig. 9, but for one data generation with two antipodal spots (left panel) and two non-antipodal spots (right panel) corresponding to the data shown in Fig. 8.

4. Discussion

The capability of an accurate determination of the geometrical angles (observer inclination i and spot co-latitude θ) of the NS is based on the swing of the PA with phase when the angles are changed (see Viironen & Poutanen 2004; Poutanen 2020). As examined by Salmi et al. (2018), having prior information on i and θ also leads to improved constraints in mass and radius of the NS (even in the case when the opposite solution with switched angles is already excluded). The assumed constraints for i and θ from polarization data in Salmi et al. (2018), which were about 4 − 6°, appear to be in the same order of magnitude as many of the constraints presented in this paper for the longest exposure times. Therefore, the accuracy improvement in mass and radius measurements could be the same order of magnitude as found there, which is at least a few per cent level in the measured mass and radius. However, the true configuration of the spot and the observer strongly affects these estimates as seen in Fig. 10.

There are also a few caveats in the presented error estimates and parameter limits. The model for the polarized radiation is based on Comptonization in the Thomson scattering limit for optically thin NS atmospheres (Viironen & Poutanen 2004). However, the formalism for Compton scattering in a hot slab should be applied for a more accurate model (see e.g. Poutanen & Svensson 1996). The choice of the scattering formalism may affect the predicted measurement accuracy of the Stokes parameters and NS geometry, especially if the modelled PD is different. Nevertheless, the effects of having different model parameters have already been accounted for, which shows how the change in PD would affect the results. However, we also note that in the most accurate simulations, more of the model parameters (both in the NS and emission models) should be kept free.

When inspecting and visualising the errors in the measured Stokes parameters in Sect. 3.1, we calculated the broadband values of q and u using a weighted average. We emphasise that while this approach provides a reasonable estimate for the total measurement accuracy it still neglects, for example, the details regarding the energy dispersion of the detector. The same also applies for the fitting of the data, where we merged the data into eight different energy channels, instead of applying a full spectro-polarimetric forward-folding fit to Q and U (the latter should be used to properly handle the energy dispersion). However, this also implies that our fits are not very prone to the bias caused by small bins having the minimum detectable polarization (MDP, see Weisskopf et al. 2010 for definition) much larger than the PD. The usual assumption of Gaussian error distribution in Stokes Q and U is not expected to be valid for that case (Mikhalev 2018). In addition, we simply assume that q and u are uncorrelated and do not employ the bivariate normal distribution presented for them in Kislat et al. (2015). These effects are not expected to largely influence the estimated errors and the parameter constraints, but they should be accounted for when considering actual observations. We also caution that the real exposure time can be restricted by the varying length of the AMP outburst, leaving only the shorter exposure times presented here applicable for such sources.

5. Conclusions

We have simulated the X-ray polarization data for AMPs that can be detected with IXPE. In order to describe AMPs spectral and polarization properties simultaneously, we have also presented a revised Thomson scattering model from Viironen & Poutanen (2004) combined with the Comptonization spectral model SIMPL. We found that the broadband relative error in the measured Stokes parameters q and u attains a 20–60% level when observing a source like SAX J1808.4−3658 with a 500 ks exposure time with data binned in 19 phase bins. However, the accuracy strongly depends on the model parameters, particularly on the relative configuration of the observer and the emitting hotspot.

We also determined the constraints for the NS parameters when fitting the simulated polarization data. We conclude that NS geometry, particularly the observer inclination and the hotspot co-latitude, could be determined with better than 10° accuracy for most of the models considered when observing at least with a 500 ks duration. This also indicates further constraints on the NS mass and radius, and thus on EOS, when modelling the X-ray pulse profiles. In addition, we found that polarization data could also be used to constrain the location of the secondary spot and probe the magnetic field geometry if the secondary spot is not obscured by an accretion disc. This could further improve the mass and radius constraints. The methods presented here can also be straightforwardly applied when analysing the soon upcoming observations by IXPE. On the other hand, the presented errors and parameter limits can be used when inferring the observational strategy of IXPE or other X-ray polarization missions such as eXTP.

Acknowledgments

This research was supported by the Academy of Finland grant 333112, the Doctoral Programme in Physical and Chemical Sciences of the University of Turku (TS), the Magnus Ehrnrooth Foundation (KK), and the Italian Space Agency (ASI) through the agreement ASI-INFN-2017.13-H0 (LB). VL, SST, and JP acknowledge support from the Russian Science Foundation grant 20-12-00364 for the development of the theoretical model (Sect. 2). TS thanks Joonas Nättilä for useful discussions. The computer resources of the Finnish IT Center for Science (CSC) and the FGCI project (Finland) are acknowledged.

References

  1. AlGendy, M., & Morsink, S. M. 2014, ApJ, 791, 78 [NASA ADS] [CrossRef] [Google Scholar]
  2. Alpar, M. A., Cheng, A. F., Ruderman, M. A., & Shaham, J. 1982, Nature, 300, 728 [NASA ADS] [CrossRef] [Google Scholar]
  3. Arnaud, K. A. 1996, in Astronomical Data Analysis Software and Systems V, eds. G. H. Jacoby, & J. Barnes (San Francisco: ASP), ASP Conf. Ser., 101, 17 [Google Scholar]
  4. Bilous, A. V., Watts, A. L., Harding, A. K., et al. 2019, ApJ, 887, L23 [CrossRef] [Google Scholar]
  5. Bogdanov, S., Lamb, F. K., Mahmoodifar, S., et al. 2019, ApJ, 887, L26 [CrossRef] [Google Scholar]
  6. Chandrasekhar, S. 1960, Radiative transfer (New York: Dover) [Google Scholar]
  7. Chandrasekhar, S., & Breen, F. H. 1947, ApJ, 105, 435 [NASA ADS] [CrossRef] [Google Scholar]
  8. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  9. Galloway, D. K., Chakrabarty, D., Morgan, E. H., & Remillard, R. A. 2002, ApJ, 576, L137 [NASA ADS] [CrossRef] [Google Scholar]
  10. Hartman, J. M., Patruno, A., Chakrabarty, D., et al. 2008, ApJ, 675, 1468 [NASA ADS] [CrossRef] [Google Scholar]
  11. Ibragimov, A., & Poutanen, J. 2009, MNRAS, 400, 492 [NASA ADS] [CrossRef] [Google Scholar]
  12. Kirsch, M. G. F., Mukerjee, K., Breitfellner, M. G., et al. 2004, A&A, 423, L9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Kislat, F., Clark, B., Beilicke, M., & Krawczynski, H. 2015, Astropart. Phys., 68, 45 [NASA ADS] [CrossRef] [Google Scholar]
  14. Lattimer, J. M. 2012, Ann. Rev. Nucl. Part. Sci., 62, 485 [NASA ADS] [CrossRef] [Google Scholar]
  15. Leahy, D. A., Morsink, S. M., Chung, Y.-Y., & Chou, Y. 2009, ApJ, 691, 1235 [Google Scholar]
  16. Leahy, D. A., Morsink, S. M., & Chou, Y. 2011, ApJ, 742, 17 [NASA ADS] [CrossRef] [Google Scholar]
  17. Lindblom, L. 1992, ApJ, 398, 569 [NASA ADS] [CrossRef] [Google Scholar]
  18. Lo, K. H., Miller, M. C., Bhattacharyya, S., & Lamb, F. K. 2013, ApJ, 776, 19 [Google Scholar]
  19. Loktev, V., Salmi, T., Nättilä, J., & Poutanen, J. 2020, A&A, 643, A84 [CrossRef] [EDP Sciences] [Google Scholar]
  20. Mikhalev, V. 2018, A&A, 615, A54 [CrossRef] [EDP Sciences] [Google Scholar]
  21. Miller, M. C., & Lamb, F. K. 1998, ApJ, 499, L37 [NASA ADS] [CrossRef] [Google Scholar]
  22. Miller, M. C., & Lamb, F. K. 2015, ApJ, 808, 31 [NASA ADS] [CrossRef] [Google Scholar]
  23. Miller, M. C., Lamb, F. K., Dittmann, A. J., et al. 2019, ApJ, 887, L24 [Google Scholar]
  24. Morsink, S. M., Leahy, D. A., Cadeau, C., & Braga, J. 2007, ApJ, 663, 1244 [NASA ADS] [CrossRef] [Google Scholar]
  25. Nagirner, D. I., & Poutanen, J. 1994, Astrophys. Space Phys. Rev., 9, 1 [NASA ADS] [Google Scholar]
  26. Pechenick, K. R., Ftaclas, C., & Cohen, J. M. 1983, ApJ, 274, 846 [Google Scholar]
  27. Pesce-Rollins, M., Lalla, N. D., Omodei, N., & Baldini, L. 2019, Nucl. Instrum. Methods Phys. Res. A, 936, 224 [CrossRef] [Google Scholar]
  28. Poutanen, J. 1994a, J. Quant. Spectr. Rad. Trans., 51, 813 [CrossRef] [Google Scholar]
  29. Poutanen, J. 1994b, ApJS, 92, 607 [NASA ADS] [CrossRef] [Google Scholar]
  30. Poutanen, J. 2010, in X-ray Polarimetry: A New Window in Astrophysics, eds. R. Bellazzini, E. Costa, G. Matt, & G. Tagliaferri (Cambridge: Cambridge University Press), 168 [CrossRef] [Google Scholar]
  31. Poutanen, J. 2020, A&A, 641, A166 [CrossRef] [EDP Sciences] [Google Scholar]
  32. Poutanen, J., & Beloborodov, A. M. 2006, MNRAS, 373, 836 [NASA ADS] [CrossRef] [Google Scholar]
  33. Poutanen, J., & Gierliński, M. 2003, MNRAS, 343, 1301 [NASA ADS] [CrossRef] [Google Scholar]
  34. Poutanen, J., & Svensson, R. 1996, ApJ, 470, 249 [NASA ADS] [CrossRef] [Google Scholar]
  35. Poutanen, J., Ibragimov, A., & Annala, M. 2009, ApJ, 706, L129 [NASA ADS] [CrossRef] [Google Scholar]
  36. Raaijmakers, G., Riley, T. E., Watts, A. L., et al. 2019, ApJ, 887, L22 [CrossRef] [Google Scholar]
  37. Radhakrishnan, V., & Srinivasan, G. 1982, Curr. Sci., 51, 1096 [NASA ADS] [Google Scholar]
  38. Riley, T. E., Watts, A. L., Bogdanov, S., et al. 2019, ApJ, 887, L21 [NASA ADS] [CrossRef] [Google Scholar]
  39. Rybicki, G. B., & Lightman, A. P. 1979, Radiative Processes in Astrophysics (New York: Wiley-Interscience) [Google Scholar]
  40. Salmi, T., Nättilä, J., & Poutanen, J. 2018, A&A, 618, A161 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Salmi, T., Suleimanov, V. F., Nättilä, J., & Poutanen, J. 2020, A&A, 641, A15 [CrossRef] [EDP Sciences] [Google Scholar]
  42. Sanna, A., Pintore, F., Bozzo, E., et al. 2017, MNRAS, 466, 2910 [NASA ADS] [CrossRef] [Google Scholar]
  43. Sanna, A., Pintore, F., Riggio, A., et al. 2018, MNRAS, 481, 1658 [NASA ADS] [CrossRef] [Google Scholar]
  44. Silva, H. O., Pappas, G., Yunes, N., & Yagi, K. 2020, ArXiv e-prints [arXiv:2008.05565] [Google Scholar]
  45. Sobolev, V. V. 1949, Uch. Zap. Leningrad Univ., 16 [Google Scholar]
  46. Sobolev, V. V. 1963, A Treatise on Radiative Transfer (Princeton: Van Nostrand) [Google Scholar]
  47. Steiner, J. F., Narayan, R., McClintock, J. E., & Ebisawa, K. 2009, PASP, 121, 1279 [NASA ADS] [CrossRef] [Google Scholar]
  48. Suleimanov, V. F., Poutanen, J., & Werner, K. 2020, A&A, 639, A33 [CrossRef] [EDP Sciences] [Google Scholar]
  49. Sunyaev, R. A., & Titarchuk, L. G. 1980, A&A, 86, 121 [NASA ADS] [Google Scholar]
  50. Sunyaev, R. A., & Titarchuk, L. G. 1985, A&A, 143, 374 [NASA ADS] [Google Scholar]
  51. Viironen, K., & Poutanen, J. 2004, A&A, 426, 985 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Watts, A. L., Andersson, N., Chakrabarty, D., et al. 2016, Rev. Mod. Phys., 88, 021001 [NASA ADS] [CrossRef] [Google Scholar]
  53. Watts, A. L., Yu, W., Poutanen, J., et al. 2019, Sci. Chin. Phys. Mech. Astron., 62, 29503 [NASA ADS] [CrossRef] [Google Scholar]
  54. Weisskopf, M. C., Elsner, R. F., & O’Dell, S. L. 2010, in Space Telescopes and Instrumentation 2010: Ultraviolet to Gamma Ray, Proc. SPIE, 7732, 77320E [CrossRef] [Google Scholar]
  55. Weisskopf, M. C., Ramsey, B., O’Dell, S., et al. 2016, in Space Telescopes and Instrumentation 2016: Ultraviolet to Gamma Ray, Proc. SPIE, 9905, 990517 [CrossRef] [Google Scholar]
  56. Zhang, S., Santangelo, A., Feroci, M., et al. 2019, Sci. Chin. Phys. Mech. Astron., 62, 29502 [NASA ADS] [CrossRef] [Google Scholar]

All Tables

Table 1.

Fiducial parameters of the synthetic data.

Table 2.

Parameters altered from the fiducial model for all of the computed models.

Table 3.

Most probable values and 68% and 95% credible intervals for five different simulations averaged from three realisations of the synthetic data produced using the fiducial model.

All Figures

thumbnail Fig. 1.

Angular dependence of emitted intensity (left panel) and PD (right panel). The angular dependence was normalized so that . Black, blue, green, orange, and red solid lines show the model for photon energies of 2, 5, 8, 12, and 18 keV, respectively. The black dot-dashed curves show the intensity and polarization corresponding to the classical results of Chandrasekhar–Sobolev (see Eq. (4)), corresponding to the optically thick electron-scattering dominated atmosphere. The grey dot-dashed curve in the left panel shows the angular dependence of intensity for unscattered photons. The model was calculated with the parameters shown in Table 1 (e.g. pmax = 0).

In the text
thumbnail Fig. 2.

Energy dependence of emitted intensity (left panel) and PD (right panel). Black, blue, green, orange, red, and magenta solid lines show the model for μ= 0.0, 0.2, 0.4, 0.6, 0.8, and 1.0, respectively. The dot-dashed black curve (left) shows the intensity corresponding to the angle-independent result of model SIMPL, namely (1 − Xsc)BE(Tbb)+XscIc(E) (see Sect. 2.1.1). The model was calculated with the parameters shown in Table 1.

In the text
thumbnail Fig. 3.

Theoretical pulse profiles of the observed flux, PD, and PA for two antipodal spots shown for three different energies (2, 5, and 8 keV). Solid black curves correspond to the total flux (or PD or PA), while the blue dashed curves are for the primary spot, and red dotted curves are for the secondary spot. The parameters of the model are given in Table 1.

In the text
thumbnail Fig. 4.

Same as Fig. 3, but for a non-antipodal spot configuration with θ = 20°, θ2 = 120°, and ϕ12 = 180°.

In the text
thumbnail Fig. 5.

Simulated PD and normalized Stokes q and u profiles for 200 ks (left panels) and 1 Ms (right panels) exposure times and the fiducial one-spot model when computing a weighted average over energy range 2–8 keV. The blue curve shows the theoretical model and the purple dots show the simulated observed data including the measurement errors.

In the text
thumbnail Fig. 6.

Average absolute minimum errors in the measured energy-integrated Stokes q and u as a function of observation time (for the one-spot models specified in Table 2). The curves were calculated for 200, 400, 600, 800, and 1000 ks (shown as purple dots) and interpolated linearly in a logarithmic scale for other values.

In the text
thumbnail Fig. 7.

Average relative minimum errors in the measured energy-integrated Stokes q and u as a function of observation time for the one-spot models (calculated as ΔS/S = ([Δq/qm]min + [Δu/um]min)/2. The solid blue curve shows the errors with the fiducial parameters (shown in Table 1). The blue dash-dotted curve corresponds to similar data but with pmax = 0.1171, the black dashed curve with i = 50°, the black dash-dotted curve with i = 70°, the red dashed curve with θ = 10°, the red dash-dotted curve with θ = 30°, and the orange solid curve with ρ = 30° (also explained in Table 2). Curves were interpolated from the calculated points as in Fig. 6.

In the text
thumbnail Fig. 8.

Simulated polarization data as in Fig. 5 (for 200 ks exposure time), but for the case of two antipodal spots with θ = 20° and θ2 = 160° (left panels), and two non-antipodal spots with θ = 20° and θ2 = 120° (right panels).

In the text
thumbnail Fig. 9.

Posterior probability distributions for NS parameters when combining fits from three realisations of synthetic data produced assuming a 200 ks (left panels) and a 1 Ms (right panels) exposure time and the fiducial model shown in Table 1. In the two-dimensional posterior distributions, the black colour shows a 68% and the grey colour shows a 95% highest posterior density credible region. In the one-dimensional posterior distributions, the red dashed lines show a 68% and the dark orange dashed lines show a 95% highest posterior density credible interval. The blue lines show the input values.

In the text
thumbnail Fig. 10.

Average width of the 95% credible interval for the measured i and θ as a function of the observing time. The curves, colours, and interpolation method correspond to those in Fig. 7, but the solid blue curve corresponds to the average constraints from the three different synthetic data generations (using the same model parameters). Purple error bars show the maximum deviation between the average and single data generation results. The other lines show constraints for only one realisation of the synthetic data (for each exposure time).

In the text
thumbnail Fig. 11.

Posterior probability distributions for NS parameters similar to the left panel of Fig. 9, but for one data generation with two antipodal spots (left panel) and two non-antipodal spots (right panel) corresponding to the data shown in Fig. 8.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.