X-ray Spectro-polarimetry with Photoelectric Polarimeters

TL;DR

This paper extends X-ray spectroscopy fitting methods to include polarization data, enabling joint analysis of intensity and polarization spectra for next-generation space polarimeters.

Contribution

It introduces a generalized forward fitting approach that incorporates polarization sensitivity through joint fitting of three spectra, including polarization parameters.

Findings

Simulations demonstrate the method's effectiveness with PRAXyS parameters.

The approach allows accurate modeling of polarization and intensity spectra.

Discussion on calibration and operation of photoelectric polarimeters.

Abstract

We derive a generalization of forward fitting for X-ray spectroscopy to include linear polarization of X-ray sources, appropriate for the anticipated next generation of space-based photoelectric polarimeters. We show that the inclusion of polarization sensitivity requires joint fitting to three observed spectra, one for each of the Stoke's parameters, I(E), U(E), and Q(E). The equations for Stoke's I(E) (the total intensity spectrum) are identical to the familiar case with no polarization sensitivity, and for which the model-predicted spectrum is obtained by a convolution of the source spectrum, F(E'), with the familiar energy response function, e(E')*R(E', E), where e(E') and R(E', E) are the effective area and energy redistribution matrix, respectively. In addition to the energy spectrum, the two new relations for U(E) and Q(E) include the source polarization fraction and position angle versus energy, a(E'), and psi'_0(E'), respectively, and the model-predicted spectra for these relations are obtained by a convolution with the "modulated" energy response function, m(E')*e(E')R(E, E'), where m(E') is the energy-dependent modulation fraction that quantifies a polarimeter's angular response to 100% polarized radiation. We present results of simulations with response parameters appropriate for the proposed PRAXyS Small Explorer observatory to illustrate the procedures and methods, and we discuss some aspects of photoelectric polarimeters with relevance to understanding their calibration and operation.

Equations

$$O(E)={\int }_{E^{\prime }}F(E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime },$$

$$O(E,\psi )={\int }_{E^{\prime }}{\int }_{{\psi }^{\prime }}H(E^{\prime },{\psi }^{\prime })X(E^{\prime },E,{\psi }^{\prime },\psi )d{E}^{\prime }d{\psi }^{\prime },$$

$$X(E^{\prime },E,{\varphi }^{\prime },\varphi )=ϵ(E^{\prime })R(E^{\prime },E)V(E^{\prime },{\varphi }^{\prime },\varphi ).$$

$$H(E^{\prime },{\psi }^{\prime })=F(E^{\prime })+W(E^{\prime })+Z(E^{\prime }),$$

$$H(E^{\prime },{\psi }^{\prime })=h(E^{\prime })+g(E^{\prime }),$$

$$W(E^{\prime })=\frac{F(E^{\prime })a(E^{\prime })\tan (2{\psi }_0^{\prime }(E^{\prime }))}{(1+{\tan }^2(2{\psi }_0^{\prime }(E^{\prime })){)}^{1/2}}=F(E^{\prime })a(E^{\prime })\sin (2{\psi }_0^{\prime }(E^{\prime }))=g(E^{\prime })\sin (2{\psi }_0^{\prime }(E^{\prime }))$$

$$Z(E^{\prime })=\frac{F(E^{\prime })a(E^{\prime })}{(1+{\tan }^2(2{\psi }_0^{\prime }(E^{\prime })){)}^{1/2}}=F(E^{\prime })a(E^{\prime })\cos (2{\psi }_0^{\prime }(E^{\prime }))=g(E^{\prime })\cos (2{\psi }_0^{\prime }(E^{\prime })).$$

$$O(E,\psi )={\int }_{E^{\prime }}{\int }_{{\psi }^{\prime }}\left(h(E^{\prime })+g(E^{\prime })\delta ({\psi }^{\prime }-{\psi }_0^{\prime }(E^{\prime }))\right)ϵ(E^{\prime })R(E^{\prime },E)C_0(A(E^{\prime })+B(E^{\prime }){\cos }^2(\psi -{\psi }^{\prime }(E^{\prime })))d{E}^{\prime }d{\psi }^{\prime }.$$

$${\int }_{E^{\prime }}\frac{h(E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)}{\pi (A(E^{\prime })+B(E^{\prime })/2)}{\int }_{{\psi }^{\prime }}\left[A(E^{\prime })+B(E^{\prime }){\cos }^2(\psi -{\psi }^{\prime }(E^{\prime }))\right]d{\psi }^{\prime }d{E}^{\prime },$$

$${\int }_{E^{\prime }}h(E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime }.$$

$${\int }_{E^{\prime }}{\int }_{{\psi }^{\prime }}\frac{g(E^{\prime })\delta ({\psi }^{\prime }-{\psi }_0^{\prime }(E^{\prime }))ϵ(E^{\prime })R(E^{\prime },E)}{\pi (A(E^{\prime })+B(E^{\prime })/2)}\left[A(E^{\prime })+B(E^{\prime }){\cos }^2(\psi -{\psi }^{\prime }(E^{\prime }))\right]d{\psi }^{\prime }d{E}^{\prime }.$$

$${\int }_{E^{\prime }}\frac{g(E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)}{(A(E^{\prime })+B(E^{\prime })/2)}\left[A(E^{\prime })+B(E^{\prime }){\cos }^2(\psi -{\psi }_0^{\prime }(E^{\prime }))\right]d{E}^{\prime },$$

$$\frac{A(E^{\prime })}{(A(E^{\prime })+B(E^{\prime })/2)}=1-\mu (E^{\prime }),$$

$$\frac{B(E^{\prime })}{(A(E^{\prime })+B(E^{\prime })/2)}=2\mu (E^{\prime }).$$

$${\int }_{E^{\prime }}g(E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)\left[(1-\mu (E^{\prime }))+2\mu (E^{\prime }){\cos }^2(\psi -{\psi }_0^{\prime }(E^{\prime }))\right]d{E}^{\prime }.$$

$$O(E,\psi )={\int }_{E^{\prime }}((h(E^{\prime })+g(E^{\prime })(1-\mu (E^{\prime })))ϵ(E^{\prime })R(E^{\prime },E)$$

$$+2g(E^{\prime })\mu (E^{\prime })ϵ(E^{\prime })R(E^{\prime },E){\cos }^2(\psi -{\psi }_0^{\prime }(E^{\prime }))))d{E}^{\prime }.$$

$$\alpha (E^{\prime },E)+\beta (E^{\prime },E){\cos }^2(\psi -{\psi }_0^{\prime }(E^{\prime }))=\left((\alpha (E^{\prime },E)+\frac{\beta (E^{\prime },E)}{2}\right)$$

$$+\left(\frac{\beta (E^{\prime },E)}{2}\sin (2{\psi }_0^{\prime }(E^{\prime }))\right)\sin (2\psi )+\left(\frac{\beta (E^{\prime },E)}{2}\cos (2{\psi }_0^{\prime }(E^{\prime }))\right)\cos (2\psi )$$

$$O(E,\psi )=I(E^{\prime },E)+U(E^{\prime },E)\sin (2\psi )+Q(E^{\prime },E)\cos (2\psi ),$$

$$I(E)={\int }_{E^{\prime }}(h(E^{\prime })+g(E^{\prime })(1-\mu (E^{\prime }))+g(E^{\prime })\mu (E^{\prime }))ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime }$$

$$={\int }_{E^{\prime }}(h(E^{\prime })+g(E^{\prime }))ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime }$$

$$U(E)={\int }_{E^{\prime }}g(E^{\prime })\mu (E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)\sin (2{\psi }_0^{\prime }(E^{\prime }))d{E}^{\prime },$$

$$Q(E)={\int }_{E^{\prime }}g(E^{\prime })\mu (E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)\cos (2{\psi }_0^{\prime }(E^{\prime }))d{E}^{\prime },$$

$$I(E)={\int }_{E^{\prime }}F(E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime }$$

$$U(E)={\int }_{E^{\prime }}W(E^{\prime })\mu (E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime }$$

$$Q(E)={\int }_{E^{\prime }}Z(E^{\prime })\mu (E^{\prime })ϵ(E^{\prime })R(E^{\prime },E)d{E}^{\prime }.$$

$$O(E,\psi )=I(E)+U(E)\sin (2\psi )+Q(E)\cos (2\psi ).$$

$$V(E^{\prime },{\varphi }^{\prime },\varphi )=(1/\pi )\left((1-\mu (E^{\prime }))+2\mu (E^{\prime }){\cos }^2(\varphi -{\varphi }^{\prime })\right).$$

Full text

polarization — methods: data analysis — techniques: polarimetric — X-rays: general — instrumentation: polarimeters

1 Introduction

In comparison with imaging, timing and spectroscopic measurements, polarization remains the “missing piece of the puzzle” of observational X-ray astrophysics. To date, the only measurement in the 2 - 10 keV band is still the $\approx 20\%$ polarization fraction inferred for the Crab nebula (Weisskopf 1978). The dearth of additional detections has largely been because of a lack of instruments sensitive enough to make such observations. However, with the advent of micro-pattern gas detectors which can directly leverage the photoelectric effect to infer linear polarization (Costa et al. 2001; Black et al. 2004, 2010), it is likely that this situation will change in the not-too-distant future. Indeed, a number of X-ray polarimetry mission concepts have been proposed in the last few years. Among these are several photoelectric effect polarimeters sensitive in the 2 - 10 keV band, including the Polarimeter for Relativistic Astrophysical X-ray Sources (PRAXyS) Small Explorer (Hill et al. 2014; Jahoda et al. 2016), and the Imaging X-ray Polarimetry Explorer (IXPE, Weisskopf et al. 2016). Also, the proposed PolSTAR experiment (Krawczynski et al. 2015) would pair a hard X-ray mirror similar to that flown on NuSTAR with a passive scattering element to provide broad band X-ray polarimetry from $\approx 3-50$ keV. Above 10 keV Compton scattering becomes competitive with the photoelectric effect, and several instruments have recently been developed to exploit this to enable polarimetry in the hard X-ray band. Among these are the balloon-borne payloads X-Calibur (Beilicke et al. 2014) and PoGOLite (Chauvin et al. 2016). Additionally, ESA has recently selected the X-ray Imaging Polarimetry Explorer (XIPE, Soffitta et al. 2016) for study and possible implementation as ESA’s M4 medium class mission. Finally, we note that shortly after submission of this paper, NASA selected IXPE for implementation in the 2020 timeframe as a Small Explorer mission (SMEX).

In anticipation of the further opening of the polarization window in the X-ray band it is timely to explore the question how one can properly generalize X-ray spectroscopic observations to include linear polarization of X-ray sources. That is, what is the additional computational “machinery” required to do X-ray spectro-polarimetry from space observatories. In this paper we outline in some detail how to infer the physical properties of sources which include linear polarization properties with space-based photoelectric polarimeters. This includes a generalization of the standard methods of X-ray spectral “forward fitting” to include polarization properties, as well as a discussion of the detector calibration information that is needed.

We also discuss some aspects of photo-electric polarimeters relevant to understanding their calibration and operation, and we present the results of simulations that illustrate the procedures and methods using spectro-polarimetric capabilities appropriate for the proposed PRAXyS Small Explorer observatory. We note that the methods described here should also prove applicable to data expected from NASA’s recently selected IXPE small explorer observatory, as well as other future X-ray polarimetry missions.

2 Background

For the simpler, well-known case where X-ray detectors have no polarization sensitivity, an observation of an astrophysical source, at least in the context of spectroscopy, can be characterized by a physical input (intrinsic) source spectrum, and the detector’s energy response function, which quantifies the rate at which photons of intrinsic energy $E^{\prime}$ are observed in detector energy channel $E$. The energy response function is often broken up into two components, the energy redistribution matrix, $R(E^{\prime},E)$, which is defined as the probability that a photon of energy $E^{\prime}$ is detected in channel $E$, and the effective area, $\epsilon(E^{\prime})$, which quantifies the detector’s collecting area as a function of energy. The count rate spectrum predicted to be observed in a detector is then the convolution,

[TABLE]

where $F(E^{\prime})$ is the intrinsic source spectrum, that is, the number of source photons with energies between $E^{\prime}$ and $E^{\prime}+dE^{\prime}$. The standard, “forward fitting” procedure for doing X-ray spectroscopy is then to introduce some model parameterization for $F(E^{\prime})$, use that model and the response functions to generate predicted detector channel spectra, and then constrain the model’s parameters by comparing these to the actual observed spectra using a statistical procedure such as $\chi^{2}$ fitting (see, Lampton, Margon & Bowyer 1976). For example, the XSPEC software package is a commonly used tool in the X-ray astrophysics community to implement this procedure (Arnaud 1996).

3 Polarimetry

The generalization to spectro-polarimetry adds an additional observable, the sky angle, $\psi$, of the polarization vector of each detected photon. A convenient way to define the sky angle is to construct a coordinate system on the sky with $x$ and $y$ axes defined as local North and East, respectively. The polarization angle is then taken as the azimuthal angle increasing in the counter-clockwise direction when looking at the sky. The $z$-axis points from the source to the observer, forming a right-handed coordinate system. It is standard to define local North as having an angle of zero, and then an angle pointing East is $\pi/2$. For linear polarimetry a complete description can be given with this angle defined on the range from [math] to $\pi$. We note that this description is consistent with the IAU convention for polarimetry (IAU 1974; see also Hamaker & Bregman 1996).

Whereas with no polarization sensitivity the detector response matrix is a two dimensional function, for polarimetry it becomes a four dimensional function, $X(E^{\prime},E,\psi^{\prime},\psi)$, that describes the rate at which photons with intrinsic energies $E^{\prime}$ and intrinsic polarization angles $\psi^{\prime}$ are observed by a detector with energy $E$ and a 100% polarization along the angle $\psi$. A good polarimeter will be one where the observed angular distribution in $\psi$ of the counts in any detector energy channel (or range of channels) is well-peaked around the intrinsic angle $\psi^{\prime}$, the phrase “more strongly modulated” around $\psi^{\prime}$ also comes to mind.

We consider a general X-ray source as an incoherent ensemble of photons with wave-vectors $\vec{k}=\hat{z}$. Such an ensemble can always be represented as a superposition of 100% polarized photons. In this more general case the observed “spectrum” is a function of both energy and sky angle. For a given energy channel (or range of channels) the distribution of observed sky angles is commonly referred to as a modulation curve, which can be written in a general way as,

[TABLE]

where $H(E^{\prime},\psi^{\prime})$, is now a more general intrinsic source spectrum that describes the number of source photons with energy $E^{\prime}$ between $E^{\prime}$ and $E^{\prime}+dE^{\prime}$, and with polarization angles between $\psi^{\prime}$ and $\psi^{\prime}+d\psi^{\prime}$. We will introduce several general source descriptions specific to linear polarization, including a Stokes parameter description, shortly.

In practice, a detector system will produce a polarization angle measurement in a coordinate system defined by the detector’s geometry, mode of operation, and orientation. We call this the detector frame, and we refer to intrinsic and measured polarization angles in this frame as $\phi^{\prime}$ and $\phi$, respectively. The orientation of the detector with respect to the sky coordinate system defines a mapping which relates the sky angle and detector frame angles. We define this mapping as $\rho(t)$, such that $\psi^{\prime}=\phi^{\prime}+\rho(t)$, and $\psi=\phi+\rho(t)$. If the detector orientation is fixed in time, then the mapping is simply a constant offset, $\rho_{0}$. If however the detector is rotated at a constant rate about the line of sight (say, by the planned rolling of a spacecraft around the line of sight), then the mapping is a simple linear function of time, $\rho(t)=\Omega t$, where $\Omega$ is the angular roll rate of the spacecraft.

Now, let’s begin by exploring an ansatz where the detector response function, $X(E^{\prime},E,\phi^{\prime},\phi)$, can be factored into two terms, the standard energy response function, and a new term, $V(E^{\prime},\phi^{\prime},\phi)$, that describes the angular response of the polarimeter. This term effectively describes how well a photon’s intrinsic polarization angle $\phi^{\prime}$ can be measured. With this description the full response function can be written,

[TABLE]

Specifically, the angular term, $V(E^{\prime},\phi^{\prime},\phi)$, describes the probability for a photon of intrinsic energy $E^{\prime}$ and intrinsic polarization angle $\phi^{\prime}$ to be observed with a polarization angle $\phi$.

We note that there could be additional ways of considering how to construct $X(E^{\prime},E,\phi^{\prime},\phi)$, however, we argue that the above is reasonable for several reasons. First, It is known that grazing incidence X-ray optics essentially do not alter the polarization properties of incident photons at levels relevant for this discussion (see Almeida & Pillet 1993; Chipman et al. 1992, 1993; Hill et al. 2016). Second, the dominant photoelectric cross-section in the relevant energy band (ionization of K shell electrons in an $s$ state) is not dependent on the polarization angle of the photon, and while some materials can display linear dichroism, a variation in the absorption strength with polarization direction, this is generally a small effect and only confined to narrow bands around some absorption edges (see for example, Collins 1997; Bannister et al. 2006). Thus, to good approximation the relevant physics associated with $\epsilon(E^{\prime})$ and $R(E^{\prime},E)$ is largely independent of the photon polarization direction. Finally, one will ultimately determine the response function in an empirical fashion, by calibrating a detailed physical model of the detector system against actual measurements. For example, one can carry out laboratory measurements whereby $100\%$ polarized beams of photons of known energy and polarization angle illuminate the detector. The observed modulation curves for different input photon energies can then be used as a direct estimator of the angular distribution $V(E^{\prime},\phi^{\prime},\phi)$, and one can explore empirically the extent to which equation (3) is realized in practice.

For photoelectric polarimeters, which infer the angular direction of an ejected photoelectron, the angular response function is generally of the form, $V(E^{\prime},\phi^{\prime},\phi)=v(E^{\prime},(\phi-\phi^{\prime}))=C_{0}(A(E^{\prime})+B(E^{\prime})\cos^{2}(\phi-\phi^{\prime}))$, where $C_{0}$ is an overall normalization factor (see, for example, Costa et al. 2001). From this expression one can define the energy-dependent modulation fraction, $\mu(E^{\prime})\equiv(V_{max}-V_{min})/(V_{max}+V_{min})=B(E^{\prime})/(2A(E^{\prime})+B(E^{\prime}))$ (Costa et al. 2001; Strohmayer & Kallman 2014). Alternatively, one can also write this response using an equivalent Stokes formalism as, $V(E^{\prime},\phi^{\prime},\phi)=C_{0}(i(E^{\prime})+u(E^{\prime})\sin(2\phi)+q(E^{\prime})\cos(2\phi))$, where now the energy-dependent modulation fraction is $\mu(E^{\prime})=(u^{2}(E^{\prime})+q^{2}(E^{\prime}))^{1/2}/i(E^{\prime})$, and $\phi^{\prime}=1/2\tan^{-1}(u(E^{\prime})/q(E^{\prime}))$. A proper angular response function must be normalized such that the integrated probability gives unity. It is straightforward to show that for the function $V(E^{\prime},\phi^{\prime},\phi)=C_{0}(A(E^{\prime})+B(E^{\prime})\cos^{2}(\phi-\phi^{\prime}))$, and with $\phi$ ranging from [math] to $\pi$, $C_{0}=(\pi(A(E^{\prime})+B(E^{\prime})/2))^{-1}$.

3.1 Source Descriptions

First, we emphasize that photoelectric polarimeters are only sensitive to linear polarization, and not, circularly polarized radiation. That is, a purely circularly polarized X-ray beam would appear “unpolarized,” to such a detector, in the sense that the observed angular distribution of polarization position angles would be identical to that produced by a completely unpolarized beam. Thus, when we refer to an unpolarized spectrum or flux, this technically also includes any circularly polarized component of the source flux. For the case of linear polarization the spectro-polarimetric properties of X-ray sources can be described in several equivalent ways. A convenient description uses the so-called Stokes parameter decomposition. Here one can define the source spectrum as,

[TABLE]

where $W(E^{\prime})$ and $Z(E^{\prime})$ describe the linear polarization properties of the source, and $F(E^{\prime})$ is the so-called total intensity energy spectrum. The fractional polarization amplitude is then given by, $a(E^{\prime})=(W^{2}(E^{\prime})+Z^{2}(E^{\prime}))^{1/2}/F(E^{\prime})$, and the source polarization position angle is $\psi^{\prime}_{0}(E^{\prime})=1/2\tan^{-1}(W(E^{\prime})/Z(E^{\prime}))$. An equivalent description can be given using the un-polarized, $h(E^{\prime})$, and polarized, $g(E^{\prime})$, spectra,

[TABLE]

and one must also define the energy-dependent polarization position angle, $\psi^{\prime}_{0}(E^{\prime})$. In this case the polarization amplitude is $a(E^{\prime})=g(E^{\prime})/(h(E^{\prime})+g(E^{\prime}))$. A convenient feature of this prescription is that each spectral component has a simple angular dependence. By definition, the un-polarized spectrum, $h(E^{\prime})$, is uniformly distributed with intrinsic polarization angles $\psi^{\prime}$ ranging from [math] to $\pi$, and the polarized spectrum, $g(E^{\prime})$, only has intrinsic angles $\psi^{\prime}(E^{\prime})=\psi^{\prime}_{0}(E^{\prime})$. These two descriptions being equivalent, it is straightforward to show that the polarized spectrum is $g(E^{\prime})=(W^{2}(E^{\prime})+Z^{2}(E^{\prime}))^{1/2}$, and the un-polarized spectrum is $h(E^{\prime})=F(E^{\prime})-(W^{2}(E^{\prime})+Z^{2}(E^{\prime}))^{1/2}$.

When generating physical source models researchers may find it more convenient to compute the total spectrum, $F(E^{\prime})$, as well as the energy-dependent polarization amplitude, $a(E^{\prime})$, and position angle, $\psi^{\prime}_{0}(E^{\prime})$. In such a case it is then straightforward to determine the Stokes spectra, $W(E^{\prime})$ and $Z(E^{\prime})$ from the above equations. Doing this one finds,

[TABLE]

and

[TABLE]

3.2 X-ray Spectro-polarimetry

We can now use equation (2) above to determine the observed spectrum for a particular source description and detector response functions. For the source description we will use $H(E^{\prime},\psi^{\prime})=h(E^{\prime})+g(E^{\prime})\delta(\psi^{\prime}-\psi^{\prime}_{0}(E^{\prime}))$ since each component has a well-defined angular distribution, which will simplify the integrations over $\psi^{\prime}$. Here, the Dirac delta function restricts the intrinsic polarization angles, $\psi^{\prime}$, for the polarized spectrum to the position angle for the source, $\psi^{\prime}_{0}(E^{\prime})$. We will also use $V(E^{\prime},\phi^{\prime},\phi)=C_{0}(A(E^{\prime})+B(E^{\prime})\cos^{2}(\phi-\phi^{\prime}))$. Note that with this form for $V(E^{\prime},\phi^{\prime},\phi)$, we have that $\cos^{2}(\phi-\phi^{\prime})=\cos^{2}((\psi-\rho(t))-(\psi^{\prime}-\rho(t)))=\cos^{2}(\psi-\psi^{\prime})$, and it suffices to write the response function directly in terms of $\psi$ and $\psi^{\prime}$, although we emphasize that more complex detector angular response functions are possible. Substituting these expressions into equation (2) gives,

[TABLE]

For the term proportional to the unpolarized component, $h(E^{\prime})$, we can rewrite the integral as,

[TABLE]

where we have explicitly included $C_{0}$ given above. The integration over $\psi^{\prime}$ is now straightforward, and is simply equal to $\pi(A(E^{\prime})+B(E^{\prime})/2)$, thus for this term the $\psi$ dependence integrates out. This is as expected, since for the assumed form of $V(E^{\prime},\psi^{\prime},\psi)=v(E^{\prime},(\psi-\psi^{\prime}))$, where $\mu(E^{\prime})$ is a function of $E^{\prime}$ only (and not $\psi^{\prime}$), an unpolarized source produces a uniform distribution in observed angle $\psi$. We emphasize that if this is not the case, and the modulation fraction is also a function of the intrinsic polarization angle, $\psi^{\prime}$, then one should not expect an unpolarized source to produce a uniform (flat) angular distribution in observed angle $\psi$ (we discuss this further in §4). This term can then be written as,

[TABLE]

Thus, as intuition would suggest, the unpolarized term looks exactly like the analogous case with no polarization sensitivity.

The remaining term involving the polarized spectrum, $g(E^{\prime})$, is,

[TABLE]

In this case the angular integration is simplified by the delta function, which restricts the source polarized photons to have the intrinsic angle $\psi^{\prime}_{0}(E^{\prime})$. This reduces the integral to,

[TABLE]

where we have picked up a factor of $\pi$ from the delta function integration. Now, using the definition of $\mu(E^{\prime})$ above it can be shown that

[TABLE]

and

[TABLE]

With these substitutions we can re-write equation (12) in the form,

[TABLE]

Combining all the terms from equation (8) we have the result,

[TABLE]

[TABLE]

This integrand is of the form $\alpha(E^{\prime},E)+\beta(E^{\prime},E)\cos^{2}(\psi-\psi^{\prime}_{0}(E^{\prime}))$ if we make the identifications, $\alpha(E^{\prime},E)=(h(E^{\prime})+g(E^{\prime})(1-\mu(E^{\prime})))\epsilon(E^{\prime})R(E^{\prime},E)$, and $\beta(E^{\prime},E)=2g(E^{\prime})\mu(E^{\prime})\epsilon(E^{\prime})R(E^{\prime},E)$. Now, with the help of some trigonometric identities it is straightforward to show that,

[TABLE]

[TABLE]

(see, for example, Strohmayer & Kallman 2014). If we further define $I(E^{\prime},E)=(\alpha(E^{\prime},E)+\beta(E^{\prime},E)/2)$, $U(E^{\prime},E)=(\beta(E^{\prime},E)/2)\sin(2\psi^{\prime}_{0}(E^{\prime}))$, and $Q(E^{\prime},E)=(\beta(E^{\prime},E)/2)\cos(2\psi^{\prime}_{0}(E^{\prime}))$, then equation (16) can be expressed in the familiar Stokes form

[TABLE]

with,

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

where now we have explicitly included the integration over $E^{\prime}$. These equations relate the intrinsic source properties (spectral and linear polarization properties defined by $h(E^{\prime})$, $g(E^{\prime})$ and $\psi^{\prime}_{0}(E^{\prime})$) to the modulation curve observed by a polarization sensitive detector characterized by three response functions, the traditional energy response functions $\epsilon(E^{\prime})$ (effective area) and $R(E^{\prime},E)$ (energy redistribution matrix), and the energy dependent modulation fraction $\mu(E^{\prime})$, which encompasses the detector’s polarization sensitivity. With the help of equations (6) and (7) these can also be written as,

[TABLE]

[TABLE]

and

[TABLE]

In thinking further about modeling X-ray sources including their linear polarization properties, the above discussion outlines a path. An observation in a particular energy channel, $E$ (or range of channels), is a background-subtracted counts or count rate modulation curve of the form,

[TABLE]

One can perform a $\chi^{2}$ fit to the observed modulation curve for each energy channel, producing the three observed Stoke’s spectra, with their associated uncertainties. One can then define source models, using, for example, parameterizations for $F(E^{\prime})$, $a(E^{\prime})$, and $\psi^{\prime}_{0}(E^{\prime})$ (recall that the polarized spectrum $g(E^{\prime})=F(E^{\prime})a(E^{\prime})$), generate predicted spectra using equations (25), (26) and (27), and then carry out $\chi^{2}$ minimization to find the source parameters which best fit the observed spectra in a statistical sense. This is entirely analogous to the simpler case with no polarization sensitivity, except that the generalization to spectro-polarimetry requires the joint fitting of three observed spectra, one for each of the Stokes parameters. One model spectrum, $F(E^{\prime})$, is folded through the full detector response function, $\epsilon(E^{\prime})R(E^{\prime},E)$, and the two new spectra, $W(E^{\prime})=F(E^{\prime})a(E^{\prime})\sin(2\psi^{\prime}_{0}(E^{\prime}))$ and $Z(E^{\prime})=F(E^{\prime})a(E^{\prime})\cos(2\psi^{\prime}_{0}(E^{\prime}))$ are folded through the “modulated response” function, $\mu(E^{\prime})\epsilon(E^{\prime})R(E^{\prime},E)$.

The forward fitting procedure implemented by XSPEC can in principle accommodate this process with a few simple additions. For example, $F(E^{\prime})$ is a physical energy spectrum, and the XSPEC package includes many such options. In order to compute the model-predicted spectra $U(E)$ and $Q(E)$ one would require model parameterizations for $a(E^{\prime})$ and $\psi^{\prime}_{0}(E^{\prime})$. Such models do not yet exist in the current XSPEC implementation, but they could be easily added. In XSPEC parlance they would be relatively simple multiplicative model components. Further, the additional detector modulation function, $\mu(E^{\prime})$, is very much like an effective area function, which can be included in the XSPEC implementation as a so-called “ancillary response function” (an “arf” file), so this could easily be incorporated in the same way.

4 Detector Considerations

In the discussion above we assumed that the detector’s angular response function satisfied the condition that $V(E^{\prime},\phi^{\prime},\phi)=v(E^{\prime},(\phi-\phi^{\prime}))$, that is, we assumed that for a given photon energy $E^{\prime}$ all intrinsic polarization angles, $\phi^{\prime}$, produce the same modulation fraction. This need not necessarily be the case. A simpler way to say this is that a detector system could in principle measure some intrinsic polarization angles better than others.

Photoelectric polarimeters work by imaging the charge track of a photoelectron produced when an X-ray photon is absorbed in the detection gas (see Costa et al. 2001; Black et al. 2004). The charge track must be drifted some distance and then detected in a pixellated detector/readout system. Detectors of this type have been developed in two basic geometries. One, which we call an “imaging” polarimeter, drifts the charge track in the same direction as the incident photon beam (Costa et al. 2001). The other, known as a “time projection chamber” (TPC) polarimeter, drifts the charge in a direction orthogonal to the photon beam (Black et al. 2004). When a charge track drifts it also diffuses, smearing out the track. If it drifts too far before being detected diffusion will completely erase the directional (polarization) information within the track. In order to keep diffusion to a reasonable level an imaging polarimeter must have a relatively shallow layer of detection gas (unless diffusion can be reduced in some other fashion). The pixellated readout system makes up the bottom of this layer. This provides a limit to the efficiency of detectors constructed in this geometry. By contrast, a TPC polarimeter can contain a much greater depth of absorbing gas (and therefore have a greater efficiency, other things being equal), since the tracks are drifted to the side of the detection volume. However, the position of the track within the field of view is better-sampled in the imaging geometry, so by centroiding the track or reconstructing it, one can estimate where it ocurred on the sky (hence the imaging appellation). Since the track drifts to the side in a TPC the location information regarding where it interacted within the volume is at least partially lost.

Several effects in the photoelectron track imaging and reconstruction process can, in principle, result in some intrinsic angles producing higher modulations. For example, the detector readout formats have some intrinsic, pixellated geometry. Some employ either hexagonal or Cartesian geometries. Since the track image is generally of modest resolution, angles corresponding to symmetries of the underlying readout geometry could in principle be better resolved than others. This could result in particular intrinsic angles producing a higher modulation than others.

Another possible cause of variation in the modulation with intrinsic angle results from potential drift asymmetries. Since tracks that drift for longer have a greater time to diffuse, the longer a track drifts the poorer, on average, is the accuracy with which the track angle can be measured. For one thing this means that the modulation fraction in such a detector is a function of the drift distance. Since X-rays are absorbed over a range of drift distances, the detector’s angular response will be an average over the modulation as a function of drift distance, weighted by the relative number of photons absorbed at each distance. Since the average depth at which an X-ray is absorbed is also a function of X-ray energy, the total distribution of drift distances will depend to some level on the intrinsic photon energy. As noted above, these effects are included in any experimentally determined angular response function.

Considering a single interaction point for simplicity, there is a potential drift asymmetry introduced by the geometry of a TPC compared to an imaging polarimeter. In an imaging polarimeter the photoelectrons are preferentially ejected in a plane that is parallel to the detector plane. Thus, for the imaging polarimeter, most of the tracks are drifted for the same distance. This symmetry is broken in the TPC geometry, since the photoelectrons are now preferentially ejected in a plane which is perpendicular to the detection plane. In this case half the tracks are ejected toward the detection plane and half away from it. Those directed at the detection plane will travel for a somewhat shorter distance (time), and will suffer slightly less diffusion than those directed away from the detection plane. If we measure track angles from 0 to $2\pi$ with [math] (and $2\pi$) representing tracks directed straight at the detection plane and $\pi$ those directed straight away from the detection plane, then those with intrinsic angle of zero diffuse slightly less on average than those with angle $\pi$. Thus, the charge tracks with angle [math] can be expected to be slightly “sharper” (better resolved) and therefore their angles somewhat more accurately determined. While it is beyond the scope of this paper to explore such effects in detail, the discussions above establish that it should not be totally unexpected for a photoelectric polarimeter to be described by an angular response function that does not exactly satisfy $V(E^{\prime},\phi^{\prime},\phi)=v(E^{\prime},(\phi-\phi^{\prime}))$, although in many circumstances the assumption of uniformity is likely to be a very good one.

Nevertheless, it is instructive to explore the non-uniorm response case a bit further, and below we will show that in such a case the response to an unpolarized source does not produce a flat modulation curve. We can demonstrate this by introducing a slightly more general response function. First, the properly normalized “uniform” response function used above can be expressed as,

[TABLE]

We can introduce a simple “non-uniform” response by defining $\mu(E^{\prime},\phi^{\prime})=\mu(E^{\prime})\eta(\phi^{\prime})$, where $\eta$ is a slowly varying function of $\phi^{\prime}$. We can mimic an illustrative non-uniform effect as described above by defining $\eta=1-(2d/\pi)|\phi^{\prime}|$. Here, $d$ is just a small constant factor that defines the size of a linear change in the modulation as $\phi^{\prime}$ ranges from [math] to $\pm\pi/2$. In this case the modulation would be a maximum for $\phi^{\prime}=0$ and drop to $1-d$ at $\phi^{\prime}=\pm\pi/2$. Figure 1 compares the response functions in these uniform and non-uniform cases. The black curves show uniform response functions (with $\mu(E^{\prime})=0.5$) for 5 equally spaced values of $\phi^{\prime}$. One can see that the peaks shift for different $\phi^{\prime}$ values, but the maximum values of the response are equal (as expected for a uniform response). The dashed black line is equal to the integral of the uniform response over $\phi^{\prime}$ (minus 0.68 to plot it within the same $y$ range, the integral is unity by definition), and is flat as expected. The red curves in Figure 1 show a non-uniform response function for the same $\phi^{\prime}$ values, and with $\mu=0.5$ and $d=0.3$. It is fairly easy to see how the modulation amplitude drops as $\phi^{\prime}$ moves away from 0. The dashed red curve shows the integral over all $\phi^{\prime}$ for the non-uniform response function (again minus 0.68 for plotting purposes), and it is evidently not flat, but peaks where the response shows the largest modulation amplitude. It should be emphasized that the value of $d$ chosen in this case is purely for illustrative purposes only, and is not meant to represent any particular detector system.

While some might consider such an effect a detector “systematic,” this is really not an accurate description, as it is simply a part of “how the detector works,” and one can still generate predicted modulation curves, and do spectro-polarimetry based on equation (2), though, depending on the complexity of an actual non-uniform response, the resulting observed modulation curves may not have a simple analytical representation, and numerical evaluation of the angular integrals in (2) could be necessary. In any case, the problem is still well defined mathematically. It should be noted that in the above discussion we were considering the observed response in a detector’s frame. If the detector is fixed (not rotating), then a similar response to unpolarized flux would be evident in the sky frame as well, though perhaps with some constant offset in angle specified by the mapping from detector to sky frames.

Since polarization properties have not been measured for most astrophysical sources, a commonly raised concern is that a “non-uniform” response to unpolarized flux might be falsely claimed as a significant polarization measurement. But this again assumes that all polarimeters will produce flat modulation curves when illuminated with unpolarized flux, and this is only true if the response is uniform as described above. The key issue, as with any observational claims, is that they be based on accurate and reliable instrument calibrations. Thus, polarimeters with non-uniform responses can also be effective instruments, however, the potential risk that a non-uniformity could lead to a false polarization claim argues for more careful attention to detector modeling, calibration and monitoring in such a case.

While we emphasize that the best way to mitigate against a false polarization claim is a proper understanding of ones detector response function, one can also “enforce” flatness of the response to unpolarized flux in the sky frame by rotation of the detector about the line of sight. As outlined earlier, rotation of the detector provides a mapping from sky angle to detector frame angle of the form $\psi^{\prime}=\phi^{\prime}+\Omega t$, where $\Omega$ is the angular rotation rate of the spacecraft. Spacecraft rotation tends to enforce the condition that all intrinsic sky angles are measured at all intrinsic detector frame angles. This has the effect of “smoothing out” a non-uniform response to unpolarized flux of the kind described in §4 (see Figure 1), when plotted in sky coordinates, and as long as the source specific flux and polarization variability time scales are much longer than the spacecraft rotation period, and sufficient exposure time is achieved at all rotation angles. Thus, in this case one can always define an effectively uniform response function by carrying out the angular average created by rotation of the detector frame.

5 Simulated Observations and Data Modeling

We now walk through an example of spectro-polarimetric fitting using simulated observations with parameters and response functions appropriate for observations with the proposed PRAXyS Small Explorer mission (Iwakiri et al. 2016). As is appropriate in this case, we consider a uniform angular response function, such that the modulation function, $\mu$, is only a function of intrinsic photon energy, $E^{\prime}$ (see eqn 29).

We first define a source model using the total spectrum, $F(E^{\prime})$, the polarization amplitude, $a(E^{\prime})$, and the polarization position angle, $\psi_{0}^{\prime}(E^{\prime})$. As an illustrative example we choose model parameters consistent with the known spectrum and polarization properties of the Crab nebula (see Weisskopf et al. 1978). We use a power-law photon spectrum, $F(E^{\prime})=C_{n}E^{\prime-\alpha}$, with index $\alpha=2.1$. For the polarization properties we assume simple linear dependencies with energy. We take $a(E^{\prime})=a_{0}+\Delta a(E^{\prime}-2.6\;{\rm keV})$ and $\psi^{\prime}_{0}(E^{\prime})=\psi^{\prime}_{0}-\Delta\psi^{\prime}_{0}(E^{\prime}-2.6\;{\rm keV})$. We use $a_{0}=0.19$ and $\psi^{\prime}_{0}=156$ deg, which are approximately consistent with the measured values at $2.6$ keV, but we allow for small linear changes with $E^{\prime}$ to explore the sensitivity of a PRAXyS observation to such changes. For this purpose we take $\Delta a=0.01$ keV*-1* and $\Delta\psi^{\prime}_{0}=1$ deg keV*-1*.

For the response functions we use an effective area curve and quantum efficiency appropriate for the proposed PRAXyS Small Explorer (SMEX) observatory, and we use a modulation function, $\mu(E^{\prime})$ based on measurements obtained with and detector simulations of the PRAXyS polarimeter (Iwakiri et al. 2016). For this demonstration example we employ a simplified diagonal redistribution matrix, and since the Crab is a bright source we ignore backgrounds. Figure 2 shows the resulting full (black) and modulated (red) effective areas (in units of cm2) used in our simulations, as well as the energy dependent modulation function, expressed as a percentage (green).

To simulate an observation we carry out the following procedures. We assume a constant source count rate equal to that expected from the Crab nebula, and first assign a random event time. We take the total spectrum and fold it through the full response matrix, obtaining a predicted count rate spectrum in energy channel space. We then convert that to a cumulative distribution function and use the so-called transformation method to make random energy channel draws. Once we have the energy of the photon, we draw a random deviate between [math] and $1$ and use $a(E^{\prime})$ defined above to assign it to the “polarized” or “unpolarized” angular distributions the correct fraction of times. For example, if it is “polarized,” then we assign it the correct sky position angle, $\psi^{\prime}_{0}(E^{\prime})$, for it’s energy. If it’s “unpolarized,” we assign it a sky position angle that is a random draw from a uniform distribution of sky angles. Finally, we “detect” the event’s observed sky angle using the appropriate angular response function for its energy (see equation 29). In this way we can build up a simulated observation of $N$ events.

We also allow for the possibility of uniform angular rotation of the spacecraft. With rotation included a few additional steps are required for the simulations. First, after assigning an event’s intrinsic sky angle we find the corresponding detector angle for this sky angle based on the event time and the mapping from sky to detector angle. Then, we assign an observed detector angle using the appropriate angular response (modulation) function, and then, finally, we place that detector frame angle back on the sky (detected sky angle) using the same mapping.

5.1 Results of Crab Nebula Simulation

The proposed baseline science plan for PRAXyS would result in millions of counts detected from the Crab nebula and pulsar, so we illustrate the capabilities and methods with a simulation for 6 million detected counts from the nebula. For the $\chi^{2}$ fitting we first group the energy bins such that each new channel has at least $1.5\times 10^{5}$ counts. For each energy channel we then determine the three Stokes parameters that describe its observed modulation curve. This can be done in two effectively equivalent ways. One is to define $M$ angular bins and then bin up the events into a modulation curve and do $\chi^{2}$ fitting (see, for example, Strohmayer & Kallman 2014). However, in this case one can also use a bin-free estimator for the Stokes parameters, such that, $I(E)=N_{E}$, $U(E)=2\sum_{i=1}^{N_{E}}\sin(2\psi_{i})$, and $Q(E)=2\sum_{i=1}^{N_{E}}\cos(2\psi_{i})$ (see, for example, Montgomery & Swank 2015; Kislat et al. 2015b). Here, $N_{E}$ is just the total number of events in energy channel $E$. For $I(E)$ the uncertainty is simply $\sigma_{I}=(N_{E})^{1/2}$, and for $Q(E)$ and $U(E)$ we use $\sigma_{Q}=\sigma_{U}=(2N_{E})^{1/2}$. We used both estimates and found they give consistent results. Here we present results using the bin-free method.

With the three observed Stokes energy spectra we can now carry out $\chi^{2}$ fitting to constrain the parameters of the spectro-polarimetric model described above. We use equations (22), (23) and (24) to determine the model-predicted spectra, and we jointly fit the simulated data to all three spectra. We use a spectral model written in IDL along with a least-squares fitting routine also developed within IDL that is based on MINPACK-1 (Markwardt 2009). The model has six free parameters, two each for the power-law energy spectrum, $F(E^{\prime})$, the polarization amplitude, $a(E^{\prime})$, and the polarization position angle, $\psi^{\prime}_{0}(E^{\prime})$. Figures 3 through 8 summarize results of the simulation with the parameters described above. Figures 3 and 4 show the resulting Stokes spectra, I(E) (Figure 3), $Q(E)$ and $U(E)$ (Figure 4), along with the best fitting model spectra (solid curves running through the data points). In addition to this we also plot the difference between the data and best fitting models. For this example we have 84 total spectral bins (28 for each Stokes spectrum) and 6 free parameters for a total of 78 degrees of freedom. We find an acceptable minimum $\chi^{2}$ value of 77.6.

Figure 5 shows a “residuals” plot of $\chi=(\rm{Data}-{\rm Model})/\sigma_{\rm data}$. The inferred fractional polarization amplitude (top) and position angle (bottom) versus energy are shown in Figure 6. Here, the top panel shows the observed modulation amplitude in each energy channel, $a(E)=(Q(E)^{2}+U(E)^{2})^{1/2}/I(E)$ (red points), and the data points (black square symbols) show the inferred polarization amplitude, $a_{p}(E)=a(E)/\mu_{avg}(E)$ with $1\sigma$ error bars. Here, $\mu_{avg}(E)$ is the mean modulation in each grouped energy channel. The inferred polarization position angle versus energy, $\psi^{\prime}(E)=1/2\tan^{-1}(Q(E)/U(E))$, is shown in the bottom panel. The solid (black) line in each panel is the best fitting polarization model (either amplitude or position angle), and the dashed lines show the input models used to generate the simulation (the “true” model). The red symbols in both panels show the existing polarization measurements of the Crab nebula as an indication of the current state of knowledge (Weisskopf et al. 1978). Finally, in Figures 7 and 8 we show the derived confidence regions for both the polarization amplitude parameters (Figure 7) and the position angle (Figure 8). We show contours drawn at $\Delta\chi^{2}=2.3$ and $4.61$, which correspond to confidence levels of $68.3$ and $90$ %. The green square symbols mark the input (“true”) values. These results demonstrate that the spectro-polarimetric forward fitting procedure recovers the input model parameters within the expected statistical precision. Existing polarization measurements of the Crab nebula are consistent with no variation in either fractional amplitude or position angle with photon energy. The results of the simulations described above indicate that a mission like PRAXyS would be extremely sensitive to such variations. For example, Figures 7 and 8 indicate that with 6 million detected photons, variations in fractional amplitude at the level of 1% keV*-1*, and position angles at the level of 1 deg keV*-1*, can clearly be detected.

5.2 Summary and Conclusions

We have presented a generalization of the standard “forward fitting” procedure for X-ray spectroscopy to include linear polarization of X-ray sources. When the angular response of the polarimeter is “uniform,” in the sense that for a given photon energy all intrinsic photon polarization angles produce the same fractional modulation, then the polarization sensitivity introduces two additional observed spectra, related to the Stokes $U(E)$ and $Q(E)$ parameters. Thus, joint fitting of three observed spectra can yield constraints on spectro-polarimetric source models. The computation of the predicted spectra as a convolution of the source spectral model with the detector energy response function maintains the same familiar form, however, for the new $U(E)$ and $Q(E)$ energy spectra the appropriate detector response function is the “modulated” response, $\mu(E^{\prime})\epsilon(E^{\prime})R(E^{\prime},E)$, which is just the traditional energy response function multiplied by the detector’s energy-dependent modulation function, $\mu(E^{\prime})$. The additional functionality required for spectro-polarimetry is thus relatively straightforward, and could be incorporated within exiting X-ray spectral software tools, as for exaxmple, XSPEC, with relatively simple modifications.

Several previous studies have also explored aspects of X-ray spectral analysis in the context of polarimetry. For example, Kislat et al. (2015a) described an iterative “unfolding” method based on Bayes’ Theorem to obtain model-independent estimates of the energy spectrum and polarization properties, and they presented simulated results with this method appropriate for the X-Calibur hard X-ray Compton scattering polarimeter. In addition, Krawczynski (2011) explored a maximum likelihood analysis method for Compton polarimeters based on measuring both the azimuthal and polar angles of the scattered photons. There are similarities between these methods and the forward-folding procedure we describe here. For example, they both account for energy-dependent effects with a multidimensional response function that models how input and output observables (such as energy and position angle) are related through the detection process. The iterative procedure appears to be more computationally intensive, but in principle, returns model-independent estimates of source spectra and polarization properties. On the other hand, forward folding can likely be more easily incorporated into existing software tools (such as XSPEC), and while not strictly model-independent, it still enables important insights regarding source properties.

While it is beyond the scope of this paper to derive results appropriate for all currently proposed X-ray polarimeters, the basic methods discussed here should also be applicable to instruments working in the hard X-rays and gamma-rays. However, it is possible that additional observables may need to be included in the response functions. For example, in addition to measuring the energy and azimuthal scattering angles, the hard X-ray scattering experiments, such as X-Calibur, benefit from also measuring the polar scattering angle (Krawczynski 2011). Moreover, issues associated with uniformity of the response functions, as discussed in §4, would have to be explored for specific detector systems.

The author acknowledges helpful discussions with Keith Jahoda, Craig Markwardt, Tim Kallman, and Jean Swank. We also thank the anonymous referee for a helpful review of this paper.