Red-skewed K$\alpha$ iron lines in GX 13+1

TL;DR

This study compares relativistic and nonrelativistic models for explaining red-skewed Fe Kalpha lines in GX 13+1, finding both models fit well and are statistically indistinguishable with current data.

Contribution

It introduces a nonrelativistic windline model for Fe line profiles and tests its effectiveness against the relativistic diskline model using GX 13+1 data.

Findings

Both models fit the observed Fe line well.

Statistical tests could not distinguish between models at 5% significance.

The windline model offers an alternative explanation based on electron scattering in outflows.

Abstract

Broad, asymmetric, and red-skewed Fe Kalpha emission lines have been observed in the spectra of low-mass X-ray binaries hosting neutron stars (NSs) as a compact object. Because more than one model is able to describe these features, the explanation of where and how the red-skewed Fe lines are produced is still a matter of discussion. It is broadly accepted that the shape of the Fe Kalpha line is strongly determined by the special and general relativistic effects occurring in the innermost part of the accretion disk. In this relativistic framework, the Fe fluorescent lines are produced in the innermost part of the accretion disk by reflection of hard X-ray photons coming from the central source (corona and/or NS surface). We developed an alternative and nonrelativistic model, called the windline model, that is capable to describe the Fe line features. In this nonrelativistic framework, the line photons are produced at the bottom of a partly ionized outflow (wind) shell as a result of illumination by the continuum photons coming from the central source, and the red-skewness of the line profile is explained by repeated electron scattering of the photons in a diverging outflow. Because GX~13+1 is a well-known disk-wind source, it is a perfect target for testing the windline model and comparing it to the relativistic one. In order to access the goodness of the fit and distinguish between the two line models, we used the run-test statistical method in addition to the canonical $\chi^2$ statistical method. The diskline and windline models both fit the asymmetric GX13+1 Fe line well. From a statistical point of view, for the two observations we analyzed, the run-test was not able to distinguish between the two Fe line models, at 5% significance level.

Tables

Table 1. Table 1: XMM-Newton observations of GX 13+1 - basic information

Obs.	observation ID	Exposure Start Time (UTC)	Exposure End Time (UTC)	Exp.111111$1$ EPIC-pn exposure.
		year  month day hh:mm:ss	year  month day hh:mm:ss	(ks)
Obs. 1	0505480101	2008 March 09   18:24:01	2008 March 09   22:48:30	13.8
Obs. 2	0505480201	2008 March 11   23:48:02	2008 March 12   03:15:09	10.3

Table 2. Table 2: Comparison of the best-fit parameters of the iron emission line with the previous analysis performed by Díaz Trigo et al. ( 2012 ) .

		Obs. 1(a)𝑎(a)(a)𝑎(a)footnotemark: $\left(a\right)$				Obs. 2(b)𝑏(b)(b)𝑏(b)footnotemark: $\left(b\right)$
		Emission line parameters				Emission line parameters
Line model	Total model	E (keV)	$\sigma$ (keV)	EW (eV)		E (keV)	$\sigma$ (keV)	EW (eV)	Reference
Gaussian	Model 1(c)𝑐(c)(c)𝑐(c)footnotemark: $\left(c\right)$	${6.53}_{-0.09}^{+0.11}$	${0.72}_{-0.16}^{+0.14}$	${185}_{-62}^{+62}$		${6.63}_{-0.13}^{+0.13}$	${0.74}_{-0.19}^{+0.19}$	${108}_{-27}^{+27}$	(Díaz Trigo et al., 2012, Table 3)
Gaussian	Model 2(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$	${6.56}_{-0.07}^{+0.10}$	${0.88}_{-0.07}^{+{0.12}^{\left(h\right)}}$	${299}_{-77}^{+77}$		${6.71}_{-0.16}^{+0.12}$	${0.77}_{-0.15}^{+{0.23}^{\left(h\right)}}$	${91}_{-23}^{+48}$	(Díaz Trigo et al., 2012, Table 4)
Gaussian	Model 1(c)𝑐(c)(c)𝑐(c)footnotemark: $\left(c\right)$	${6.67}_{-0.03}^{+0.06}$	${0.50}_{-0.10}^{+0.03}$	${128}_{-18}^{+18}$		${6.57}_{-0.07}^{+0.06}$	${0.50}_{-0.07}^{+0.05}$	${78}_{-17}^{+18}$	The present paper
Gaussian	Model 3(e)𝑒(e)(e)𝑒(e)footnotemark: $\left(e\right)$	${6.872}_{-0.021}^{+0.030}$	${0.59}_{-0.05}^{+0.03}$	${226}_{-222}^{+142}$		${6.58}_{-0.05}^{+0.05}$	${0.64}_{-0.04}^{+0.04}$	${118}_{-87}^{+60}$	The present paper
Diskline	Model 4(f)𝑓(f)(f)𝑓(f)footnotemark: $\left(f\right)$	${6.674}_{-0.002}^{+0.040}$	$-$	${256}_{-73}^{+128}$		${6.26}_{-0.04}^{+0.05}$	$-$	${119}_{-25}^{+21}$	The present paper
Windline	Model 5(g)𝑔(g)(g)𝑔(g)footnotemark: $\left(g\right)$	$\left[6.6\right]$	$-$	${195}_{-47}^{+26}$		$\left[6.6\right]$	$-$	${130}_{-24}^{+25}$	The present paper

Table 3. Table 3: Best fits to the 2.5 to 10 keV EPIC pn spectra of GX13+1.

Component	Parameter	Unit	Obs.1(a)𝑎(a)(a)𝑎(a)footnotemark: $\left(a\right)$	Obs.2(b)𝑏(b)(b)𝑏(b)footnotemark: $\left(b\right)$	Obs.1(a)𝑎(a)(a)𝑎(a)footnotemark: $\left(a\right)$	Obs.2(b)𝑏(b)(b)𝑏(b)footnotemark: $\left(b\right)$	Obs.2(b,c)𝑏𝑐(b,c)(b,c)𝑏𝑐(b,c)footnotemark: $(b,c)$
Tbabs	nH	1022 atoms cm-2	${4.90}_{-0.15}^{+0.15}$	2.9${}_{}{}^{+0.6}{}_{-0.3}{}^{}$	${6.06}_{-0.06}^{+0.40}$	2.48${}_{}{}^{+0.17}{}_{-0.27}{}^{}$	${2.48}_{-0.39}^{+0.12}$
bbodyrad	$k{T}_{bb}$	keV	${1.221}_{-0.008}^{+0.003}$	1.29${}_{}{}^{+0.03}{}_{-0.03}{}^{}$	${1.220}_{-0.003}^{+0.003}$	1.505${}_{}{}^{+0.023}{}_{-0.050}{}^{}$	${1.503}_{-0.040}^{+0.025}$
CompTT	$k{T}_0$	keV	${0.132}_{-0.022}^{+0.018}$	0.51${}_{}{}^{+0.03}{}_{-0.06}{}^{}$	${0.188}_{-0.008}^{+0.060}$	0.644${}_{}{}^{+0.025}{}_{-0.019}{}^{}$	${0.65}_{-0.01}^{+0.04}$
	$k{T}_e$	keV	${2.015}_{-0.009}^{+0.007}$	2.46${}_{}{}^{+0.60}{}_{-0.17}{}^{}$	${2.503}_{-0.003}^{+0.018}$	3.84${}_{}{}^{+0.22}{}_{-0.17}{}^{}$	${3.84}_{-0.23}^{+0.10}$
	$\tau$		${3.80}_{-0.15}^{+0.40}$	5.7${}_{}{}^{+0.6}{}_{-1.0}{}^{}$	${2.54}_{-0.21}^{+0.08}$	3.01${}_{}{}^{+0.40}{}_{-0.20}{}^{}$	${3.00}_{-0.16}^{+0.40}$
Windline	${\tau }_w$		${2.55}_{-0.21}^{+0.08}$	2.57${}_{}{}^{+0.15}{}_{-0.18}{}^{}$
	$k{T}_{ew}$	keV	[0.7]	0.6${}_{}{}^{+0.7}{}_{-0.4}{}^{}$
	$\beta$		[7.5$\times {10}^{-2}$]	6.7${}_{-0.7}{}^{+0.6}\times {10}^{-2}$
Diskline	$E_L$	keV			${6.674}_{-0.002}^{+0.040}$	6.26${}_{}{}^{+0.05}{}_{-0.04}{}^{}$	$6.4$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$
	${\beta }_{10}$(e)𝑒(e)(e)𝑒(e)footnotemark: $\left(e\right)$				$-{2.38}_{-0.06}^{+0.06}$	-2.25${}_{}{}^{+0.20}{}_{-0.15}{}^{}$	$-{2.40}_{-0.10}^{+0.12}$
	$R_{in}$	$R_G$			${10}_{-2}^{+3}$	14${}_{}{}^{+4}{}_{-7}{}^{}$	${14}_{-4}^{+7}$
	$R_{out}$	$R_G$			$\left[1000\right]$	413${}_{}{}^{+212}{}_{-111}{}^{}$	${956}_{-319}^{+510}$
	$i_{out}$	deg			${60}_{-6}^{+4}$	75(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$	${61}_{-5}^{+5}$
EW		eV	195${}_{}{}^{+26}{}_{-47}{}^{}$	130${}_{}{}^{+25}{}_{-24}{}^{}$	${256}_{-73}^{+128}$	${119}_{-25}^{+21}$	${111}_{-107}^{+18}$
Gaus1	$E_L$	keV	6.701${}_{}{}^{+0.008}{}_{-0.014}{}^{}$		${6.733}_{-0.013}^{+0.003}$
	$\sigma$	keV	3.9$\times {10}^{-3}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$		2.8$\times {10}^{-3}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$
Gaus2	$E_L$	keV	7.005${}_{}{}^{+0.015}{}_{-0.015}{}^{}$	7.005${}_{}{}^{+0.014}{}_{-0.014}{}^{}$	${7.005}_{-0.015}^{+0.002}$	7.005${}_{}{}^{+0.014}{}_{-0.015}{}^{}$	${7.004}_{-0.014}^{+0.010}$
	$\sigma$	keV	1.1$\times {10}^{-4}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$	1.4$\times {10}^{-3}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$	1.1$\times {10}^{-4}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$	$6.7\times {10}^{-4}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$	$1.1\times {10}^{-4}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$
Gaus3	$E_L$	keV	7.877${}_{}{}^{+0.015}{}_{-0.018}{}^{}$		${7.860}_{-0.015}^{+0.015}$
	$\sigma$	keV	2.7$\times {10}^{-3}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$		$1.6\times {10}^{-4}$(d)𝑑(d)(d)𝑑(d)footnotemark: $\left(d\right)$
Fit quality	${\chi }^2$/d.o.f		1784.81/1481	1572.35/1485	1690.81/1478	1573.37/1483	1588.06/1483
	${\chi }_{red}^2$		1.21	1.06	1.14	1.06	1.07

Table 4. Table 4: Run-test on the fit residuals of the broad emission Fe K line of the NS LMXB GX 13+1.

Obs.	Energy	Run-test(a)𝑎(a)(a)𝑎(a)footnotemark: $\left(a\right)$ (%)				${\chi }_{\text{red}}^2$ (${\chi }^2$/dof)
	(keV)	Gaussian	windline	diskline		Gaussian	windline	diskline
1	6.0–6.67	12.6	20.0	14.4		1.16 (1714.58/1480)	1.21 (1784.81/1481)	1.14 (1690.81/1478)
2	6.0–6.67	55.2	56.9	43.0(b)𝑏(b)(b)𝑏(b)footnotemark: $\left(b\right)$		1.06 (1578.99/1486)	1.06 (1572.35/1485)	1.06 (1573.37/1483)(b)𝑏(b)(b)𝑏(b)footnotemark: $\left(b\right)$
2	6.0–6.67			32.4(c)𝑐(c)(c)𝑐(c)footnotemark: $\left(c\right)$				1.07 (1588.06/1483)(c)𝑐(c)(c)𝑐(c)footnotemark: $\left(c\right)$

Full text

*Aims. *Examining the asymmetry of the fluorescent Fe K emission line evident in the XMM-Newton EPIC-pn spectra of the NS source GX 13+1, we aim to distinguish between the two line models. Because GX 13+1 is a well-known disk-wind source, it is a perfect target for testing the windline model and compare the spectral fits between the relativistic and nonrelativistic line models.

*Methods. *We used two XMM-Newton EPIC-pn observations in which the Fe line profiles were previously reported in the literature. These observations are not strongly affected by pile-up, and the Fe emission lines appear very strong and red-skewed. In order to access the goodness of the fit and distinguish between the two line models, we used the run-test statistical method in addition to the canonical $\chi^{2}$ statistical method.

*Results. * The diskline and windline models both fit the asymmetric Fe line well that is present in the XMM-Newton energy spectra of GX 13+1. From a statistical point of view, for the two observations we analyzed, the run-test was not able to distinguish between the two Fe line models, at 5% significance level.

1 Introduction

The strong Fe $K_{\alpha}$ emission line (in the $\sim 6.4-7.0$ keV X-ray energy band) has been observed as broadened, asymmetric, and red-skewed in a number of sources, which include extragalactic and Galactic accreting X-ray compact objects: Seyfert-type I active galactic nuclei (AGNs) (e.g., Tanaka et al., 1995; Nandra et al., 1997; Fabian et al., 2000; Reynolds & Nowak, 2003; Miller, 2007; Hagino et al., 2016), black hole (BH) low-mass X-ray binaries (LMXBs) (e.g., Fabian et al., 1989; Miller et al., 2002; Miller, 2007; Reynolds & Nowak, 2003), and neutron star (NS) LMXBs (e.g., Bhattacharyya & Strohmayer, 2007; Cackett et al., 2008; Pandel et al., 2008; Shaposhnikov et al., 2009; Reis et al., 2009; Cackett et al., 2010; Ludlam et al., 2017a, b).

For accreting NSs and BHs, the standard interpretation for the Fe line asymmetry invokes general relativity (GR) effects through the strong gravitational field that is present close to the compact object, where the emission line is assumed to be produced. Although the relativistic interpretation is commonly accepted to explain the asymmetry of the observed iron line profiles, there are some considerations before this prevailing scenario can be accepted:

At the time of writing, the relativistic line interpretation has some difficulties to explain all the characteristics that are observed in these sources (e.g., X-ray variability, see section 1.1). 2. 2.

Fe line profiles with similar characteristics (i.e., broadened and red-skewed) were observed in the X-ray spectra of accreting white dwarfs (WDs) (Hellier & Mukai, 2004; Vrielmann et al., 2005; Titarchuk et al., 2009, Maiolino et al. in preparation). These systems require an alternative explanation for the line profile because GR effects do not play a role around WDs. 3. 3.

A line model based on relativistic effects is not unique in describing the broad and asymmetric Fe K emission line profiles (e.g., Hagino et al., 2016).

An alternative model, called the wind line model (windline, hereafter; presented by Laurent & Titarchuk (2007)), is able to explain the red-skewed Fe emission line profile observed in NS and BH LMXBs (Titarchuk et al., 2009; Shaposhnikov et al., 2009), and in systems containing WD as the accreting compact object (i.e., cataclysmic variables) because the line model does not require GR effects (Titarchuk et al., 2009). In the windline framework, the asymmetric and red-skewed line profile is produced by repeated down-scattering of the line photons by electrons in a diverging outflow (wind). The line is generated in a partially ionized and thin inner shell of the outflow by irradiation of hard X-ray photons coming from the central source (for more details about the model, see Laurent & Titarchuk, 2007, and discussion therein).

Motivated by the potential of the windline model in describing red-skewed iron line profiles in all three types of accreting X-ray powered sources, we applied this model on the residual excess in the Fe K energy range. This excess is observed in the spectra of all types of Galactic accreting compact objects: NSs, BHs, and WDs. In this first paper we present the investigation of the iron emission line profile in the spectrum of the LMXB GX 13+1, which hosts a neutron star as a compact object.

GX 13+1 is a perfect target for testing the windline model because strong red-skewed Fe emission lines are reported in its spectra (e.g., Díaz Trigo et al., 2012), with simultaneous observations of blueshifted Fe absorption lines. This indicates an outflow around the source (Boirin et al., 2005), (see section 1.3).

This paper is structured as follows: in the section 1.1 we present a concise background on broad and red-skewed iron lines in accreting NSs and their relevance. In section 1.2 we describe the parameters of the windline and in 1.3 the characteristics of the source GX 13+1. In section 2 we show the XMM-Newton observations that we analyzed and the data reduction. In section 3 we describe the data analysis and compare the windline and the relativistic diskline fits. We discuss our results and conclude in section 4.

1.1 Broad and red-skewed iron lines in accreting neutron stars

Fluorescent Fe K emission line profiles were first observed as broadened, asymmetric, and red-skewed in Galactic and extragalactic accreting BHs; in LMXBs and Seyfert-type I galaxies (e.g., Fabian et al., 1989; Tanaka et al., 1995; Miller, 2007, and references therein). Because the Fe K band was successfully modeled by a relativistic line model, this led to the conclusion that the breadth and asymmetry of the line in BH sources is due to effects of Doppler-broadening, transverse Doppler shift, relativistic beaming, and gravitational redshifts produced in the innermost part of the accretion disk, where GR effects are playing a role.

The subsequent discovery of red-skewed Fe emission line in the NS LMXB Serpens X–1 (Bhattacharyya & Strohmayer, 2007; Cackett et al., 2008), followed by the observation in a few more NS LMXB sources, and the satisfactory fits of the line profiles using the relativistic line model extended the interpretation of the asymmetric line profiles in terms of GR effects on accreting NSs.

In the framework of the relativistic models, the fluorescent Fe emission line is produced by the reflection of hard photons (described by a power-law spectrum, coming from the inner part of the source) by the iron (ions) on the surface of a cold ($<$ 1 keV) accretion disk. The relativistic scenario is even now commonly accepted as an explanation for red-skewed Fe K emission lines in NS (and BH) sources. The line in NS LMXBs is usually modeled assuming a Schwarzschild potential around the accreting compact object (diskline model (Fabian et al., 1989)). However, models considering a Kerr metric, created primarily for describing the emission line profiles coming from the accretion disks around rotating BHs, have also been used to fit to the iron emission lines in accreting NSs, see, for example, the laor (Laor, 1991), relconv (Dauser et al., 2010), and relxill (García et al., 2014) models. The main motivation for applying the relativistic scenario to NS sources is the fact that the inner disk radii predicted by the relativistic model have been found to agree well with the radii implied by kilohertz quasi-periodic oscillation (kHz QPO) frequency in 4U 1820-30 and GX394+2 sources as a Keplerian frequency at those radii (supporting the inner disk origin for kHz QPOs) (Cackett et al., 2008).

However, Shaposhnikov et al. (2009) and Titarchuk et al. (2009), using Suzaku and XMM-Newton data, examined the red-skewed line profile observed in the spectra of two NS LMXB sources, Cyg X–2 and Serpens X-1, in the framework of the windline model. They found that this nonrelativistic model is able to reproduce the red-skewed line profile with a fit quality (on the basis of $\chi^{2}$-statistic) similar to that obtained by relativistic reflection models. Although they were unable to conclusively rule out one of the models, they pointed out that the windline model appears to give a more adequate explanation because it does not require the accretion disk inner edge to advance close to the NS surface. In addition, the windline model explains a timing source property that the diskline model has difficulty to account for. They verified a lack of erratic fast variability when the Fe line is present in the Cyg X–2 spectrum (which was also observed in a BH LMXB source, GX 339$-$4; (Titarchuk et al., 2009)). This result weakens the red-skewed line connection with kHz QPOs and strengthens its connection to outflow phenomena because the suppression of fast variability can be explained by smearing in a strong outflow (wind) (see Figure 8, Sections 3.1 and 3.2 in Titarchuk et al., 2009).

Furthermore, Lyu et al. (2014) used Suzaku, XMM-Newton, and RXTE observations of the NS 4U 1636-53 to study the correlation between the inner disk radius (obtained from the Fe line relativistic fits) and the flux. They observed that the line did not change significantly with flux states of the source and concluded that the line is broadened by mechanisms other than just relativistic broadening. This behavior is not only observed in NS LMXBs. For example, observations of MCG–6-30-15, NGC 4051, and MCG–5-23-16 (Reynolds & Begelman, 1997; Reeves et al., 2006; Marinucci et al., 2014) show that the broad and red-skewed Fe lines in AGNs show little variability despite the large changes in the continuum flux. The iron line and the reflected emission flux do not respond to the X-ray continuum level of variability on short timescales (Reeves et al., 2006), which means that the fast changes that are almost certain to occur near the compact object are not observed in the line. This point is difficult to explain by the reflection and GR interpretation of the iron line formation, and has been suggested to be caused by strong gravitational light bending (see Miniutti & Fabian, 2004; Chiang & Fabian, 2011).

In addition, Mizumoto et al. (2018) demonstrated that the short lag time of the reprocessed Fe-K line energy band observed in AGNs can also be explained by scattering of photons coming from the X-ray central source by an outflow or disk-wind (with velocity of $\sim$ 0.1c, placed at large distances, $\sim$ 100 Rg, from the compact object). Their simulations resulted in features in the lag-energy plot and characteristics similar to those observed in 1H0707-495 and NGC 4151. They showed that the observation of short lag times does not necessarily indicate that the line is produced very close to the event horizon of the BH, as expected by the relativistic scenario (see Mizumoto et al., 2018, and references therein).

As previously mentioned, the prevailing interpretation that has been used to model the red-skewed Fe line profiles in accreting NS (and also BHs) with reflection and GR effects is not entirely successful, and there is the alternative windline model that does not need relativistic effects, but successfully describes the observed line profiles in these sources. The importance to make this comparative study between the models is not just to either confirm or reject the relativistic framework in NS LMXBs, but also to confirm or reject the implications to the NS physics derived by the models. By interpreting the broadening and red-skew of the Fe K emission line as due to GR effects, the relativistic models allow us to derive directly from the spectral fit parameters such as disk inclination, spin of the compact object, and inner radius of the accretion disk. The inner radius has been used to set an upper limit on the magnetic field strength and on the radius of the neutron stars, which in turn is used to study constraints on the NS equation of state (Ludlam et al., 2017a). Although the relativistic models allow this powerful derivation, it is important to study the asymmetry of the Fe K line with other models that provide an alternative as well as a more consistent physical scenario of the NS LMXBs spectral emission.

1.2 windline model into XSPEC

The main obstacle to analyzing the spectral data with the windline model is that the analytical solution is not available (Shaposhnikov et al., 2009). Therefore, the windline model was inserted within XSPEC with a number of FITS tables containing the results of the wind line model implemented using the Monte Carlo approach (Laurent & Titarchuk, 2007). Each FITS table is constructed for a fixed emission line energy, running a grid of values for the free parameters. The free parameters of the model are as follows:

• $\bullet$

the optical depth of the wind/outflow ($\tau_{w}$);

• $\bullet$

the temperature of the electrons in the wind/outflow ($kT_{ew}$ in keV);

• $\bullet$

the velocity ($\beta$) of the wind/outflow in unity of $c$ (speed of light);

• $\bullet$

the redshift of the source; and

• $\bullet$

the normalization (which is the number of photons in the whole spectrum, in units of photon/keV).

In this model, the line profiles strongly depend on the optical depth $\tau_{w}$ and on the velocity $\beta$ (v/c) of the outflow (wind). The shape of the red wing, below the broad peak, follows a power law with an index that is a strong function of $\tau_{w}$ and $\beta$. If the temperature of the electrons in the wind is low, kT${}_{ew}<1$ keV, the temperature does not change the line profile much, and most of the effect is due to the wind velocity, that is, the velocity effect is strong enough to change the red wing of the line at the observed level.

1.3 Source GX 13+1

The source GX 13+1 is classified as a LMXB with a NS as a compact object. In LMXBs the secondary star (which is an evolved late-type K5 III giant star in GX 13+1) fills the Roche lobe and transfers mass onto the compact object through the inner Lagrangian point (L1), feeding the primary star through an accreting disk (Lewin & van der Klis, 2006).

This source is located at a distance of $7\pm 1$ kpc and is also classified as a type I X-ray burster (Fleischman, 1985; Matsuba et al., 1995). Its bright persistent X-ray emission shows characteristics of both atoll and Z LMXBs (see, e.g., Díaz Trigo et al., 2012; Fridriksson et al., 2015, and references therein). GX 13+1 is likely an X-ray dipping source. Dipping sources are binary systems with shallow X-ray eclipses, so called dippers. These eclipses are likely caused by the periodic obscuration of the central X-ray source by a structure called bulge (or hot spot), or by outflows in the outer disk (Lewin et al., 1995; Pintore et al., 2014). The bulge is created by the collision between the accretion flow coming from the secondary star with the outer edge of the accretion disk, and it contains optically thick material. Therefore, the observation of dippers is associated with LMXB sources with very high inclination. The strong energy-dependent obscuration in GX 13+1, the absence of eclipses, and the spectral type of the secondary star indicate an inclination of 60-80*∘* (Díaz Trigo et al., 2012).

The source GX 13+1 is known as one of the two NS LMXBs showing winds (Díaz Trigo et al., 2012; Ueda et al., 2004); the other source is IGR J17480-2446 (Miller et al., 2011). A thermal driving mechanism of the wind is expected (see discussion in Díaz Trigo et al., 2012, and references therein).

Spectral characteristics of GX 13+1

The continuum emission of GX 13+1 has previously been modeled with a multicolor blackbody component plus either a blackbody or a nonthermal component (such as Comptonization, power law, or cutoff power law) without significantly worsening the quality of the fit (see Díaz Trigo et al., 2012; Pintore et al., 2014, and references therein).

Several narrow absorption lines, characteristics of dipping sources, are observed in the energy spectrum and are associated with a warm absorber (highly photoionized plasma) around the source, driven by outflows from the outer regions of the accretion disk (Pintore et al., 2014; Díaz Trigo et al., 2012). The absorption lines indicate bulk outflow velocities of $\sim$ 400 km s*-1* (see Pintore et al., 2014, and references therein).

Narrow resonant absorption lines near 7 keV were found with an ASCA observation of GX 13+1 (Ueda et al., 2001). Afterward, Kα and Kβ lines of He- and H-like Fe ions, the Kα line of H-like Ca (Ca XX) ion, and a deep Fe XXV absorption edge (at 8.83 keV) were observed with XMM-Newton (Sidoli et al., 2002). Díaz Trigo et al. (2012) additionally observed a Fe XXVI absorption edge (at 9.28 keV). The depth of the absorption features in the XMM-Newton 2008 observations changed significantly on timescale of a few days, and more subtle variations are seen on shorter timescales of a few hours (Sidoli et al., 2002). Because the narrow absorption lines are observed throughout the orbital cycle, the absorption plasma is likely cylindrically distributed around the compact object (see Pintore et al., 2014, and references therein).

Another spectral feature that is superposed on the continuum and the narrow absorption lines is the broad emission line in the K-shell Fe XVIII-XXVI energy range ($\sim$ 6.4 - 7.0 keV). This emission line is the main feature we study here. It is usually modeled, as well as in all NS sources, by relativistic line models (see section 1.1).

2 XMM Newton observations

We analyzed two public XMM-Newton EPIC-pn observations of GX13+1. The first observation, Obs. ID 0505480101 (hereafter called Obs. 1), has 13.8 ks EPIC-pn exposure time, starting on 2008 March 09 at 18:24:01 and ending at 22:48:30 UTC. The second observation, Obs. ID 0505480201 (hereafter called Obs. 2) has 13.7 ks EPIC-pn exposure time, starting on 2008 March 11 at 23:48:02 and ending at 03:15:09 UTC of the next day (see Table 1 for a log of the XMM-Newton observations used in this analysis). Both observations were taken with the camera operating in timing mode.

All the five public XMM-Newton EPIC-pn observations analyzed by Díaz Trigo et al. (2012) were previously considered for this study. However, a pre-analysis of the spectra led to the selection of the two observations used in this paper. Because all observations are affected by pile-up, we extracted two spectra for each observation, A and B, respectively. A different RAWX range was excised for each during the data processing for pile-up correction: one column was excised for A spectra from the bright central region, and either two or three RAWX columns were excised for B spectra. We fit the A and B spectra with the following total model: tbabsedge1edge2*(diskbb+bbodyrad+gaussians), which corresponds to Model 1 in Díaz Trigo et al. (2012) (Model 1; hereafter). The Gaussian components in this case correspond to negative Gaussians. We used this total model to check for the iron emission line because the absorption lines can mimic the emission line. To select the observations, we used the following criterion: observations in which the line barely appears, that is, with a residual fluctuation ¡ $2.5-3.0\leavevmode\nobreak\ \sigma$, when one or two columns are excised, were discarded. It is important to state that when the pile-up correction was performed excising five columns, as performed by Díaz Trigo et al. (2012), the emission line disappears in all observations. Therefore, we extracted two columns in the selected two observations (see section 2.1) to find a balance between pile-up correction and the presence of the line. The continuum is not strongly affected by the difference in the RAWX extraction in these two observations (i.e., by excising two or five columns).

The pile-up correction makes the spectra softer. Consequently, if the line is weak in the observation, it decreases the evidence of its presence. In addition, pile-up correction can deteriorate the statistics of the counts in the hard part of the spectrum, which contributes to a poorer constraining of the continuum model. This in turn affects the emission line presence perception and its shape definition.

2.1 Data reduction

All light curves and spectra were extracted through the Science Analysis Software (SAS) version 14.0.0. Following the recommendation for the pn camera on the timing mode of observation, we selected events in the 0.6-10 keV energy range to avoid the increased noise, and only single and double events were taken into account for the spectrum extraction (pattern in [1:4]). In order to provide the most conservative screening criteria, we used the FLAG=0 in the selection expression for the standard filters, which excludes border pixels (columns with offset), for which the pattern type and the total energy is known with significantly lower precision (see XMM-SOC CAL-TN-0018, and The XMM-Newton ABC Guide).

The calibration of the energy scale in EPIC-pn in timing mode requires a complex chain of corrections, which are important to yield a more accurate energy reconstruction in the iron line (6 - 7 keV) energy range (see XMM-SOC-CAL-TN-0083). In order to achieve the best calibration, we applied a) the X-ray loading (XRL) correction through the parameters runepreject=yes and withxrlcorrection=yes; b) the special gain correction, which is default in SASv14.0 (and can be applied through the parameter withgaintiming=yes); and c) the energy scale rate-dependent PHA (RDPHA) correction, which was introduced in SASv13.0, and taken as default in SASv14.0 for Epic-pn timing mode (XMM-SOC-CAL-TN-0018). This correction is expected to yield a better energy scale accuracy than the rate-dependent charge transfer inefficiency (RDCTI) correction for the EPIC-pn timing mode, respectively. Pintore et al. (2014) tested both rate dependence corrections in GX 13+1 in EPIC-pn timing mode. They showed that the RDPHA calibration provides more consistent centroid line energies and more physical fit parameters (such as the emission Gaussian line centroid energy and source inclination parameters obtained through the diskline model best fit). Using this correction, they found the energy of the iron mission line at 6.6 keV, and stated that it was consistent with the energies found by Díaz Trigo et al. (2012) and D’Aì et al. (2014). In the RDPHA correction, the energy scale is calibrated by fitting the peaks in derivative PHA spectra corresponding to the Si ($\sim$1.7 keV) and Au ($\sim$ 2.3keV) edges of the instrumental responses, where the gradient of the effective area is largest (XMM-SOC-CAL-TN-0306).

We extracted source+background spectrum and light curves of both observations (Obs. 1 and 2) from 10 ¡ RAWX ¡ 60 and the background from 3 ¡ RAWX ¡ 9. The SAS task epatplot was used as a diagnostic tool for pile-up in the pn-camera. It showed that both observations are affected by pile-up. In order to mitigate its effect, we excised the brightest central columns RAWX 38-39 during the spectral extraction.

The ancillary and response matrices were generated trough the SAS task arfgen (in the appropriated way to account for the area excised due to pile-up) and rmfgen. The EPIC-pn spectra were rebinned in order to have at least 25 counts in each background-subtracted channel and in order to avoid oversampling the intrinsic energy resolution by a factor larger than 3.

The light curves were produced trough the SAS task epiclccorr, which corrects the light curve for various effects affecting the detection efficiency (such as vignetting, bad pixels, PSF variation and quantum efficiency), as well as for variations affecting the stability of the detection within the exposure (such as dead time and GTIs). Because all these effects can affect source and background light curves in a different manner, the background subtraction was made accordingly through this SAS task.

2.2 Light curves

Figure 1 shows the EPIC-pn light curves extracted from Obs. 1 and 2 of GX 13+1. We present the light curves in the 0.7 - 10.0 keV, 0.7 - 4.0 keV, and 4.0 -10.0 keV energy ranges, with a bin size of 100.0 s, and also the hardness ratio between the photon counts in the last two energy bands. Both light curves show high count rate variability.

Although we see the same shape of the light curves as in Díaz Trigo et al. (2012), who analyzed the same EPIC-pn observations of GX13+1, we obtained a slightly higher count rate. We tried to find the origin of the discrepancy on the count rate (of $\sim$50 count/s). In order to do this, we first extracted the light curve in the same way as Díaz Trigo et al. (2012), that is, considering the same extraction region (private communication). We obtained the same count rate as before, which means that the difference is not due to the different extraction region used. We then extracted the event list in the same way as Díaz Trigo et al. (2012), that is, we applied the RDCTI calibration correction, and extracted the light curve without applying background subtraction, to see whether somewhat unexpectedly, the origin could be in these points. As expected, we obtained the same count rate as before and thus failed to find the origin of the discrepancy. The only possibility we see lies in the different SAS versions that were used in each analysis. We used SAS version 14.0.0 and Díaz Trigo et al. (2012) used SAS version 10.0.0. Consequently, we applied the XRL correction inserted in SAS 13.5/14.0.0. Díaz Trigo et al. (2012) calculated the “residual” of the offset map by subtracting the offset map of a nearby observation taken with closed filter from the offset map, and applied the XRL correction excising the inner region, which they found to be affected by XRL (which is the same region as was excised for pile-up correction). The differences in the XRL correction procedures might have contributed for the discrepancy in the count rate. Furthermore, the RDCTI calibration correction has been changed since the release of SASv10.0.0; the RDCTI correction applied by us (for the comparison) is not the same as was applied by Díaz Trigo et al. (2012).

In Obs. 1 we observed a maximum fractional variation in the hardness ratio of $\sim$ 7.3% with a mean value of 0.688(1). A similar variation was identified in Obs. 2, where we observed a maximum fractional variation of $\sim$ 5.2% with a mean value of 0.76(7). We confirm what was found by Díaz Trigo et al. (2012): the count rate variability is not associated with significant change in the hardness ratio. Because we did not find any significant change in the hardness ratio, we used the total EPIC-pn exposure time in the spectral extraction for both observations.

2.3 Spectral variability

Although the variation in hardness ratio is not significant, Díaz Trigo et al. (2012) pointed out that Obs. 1 shows a significant spectral change when the spectrum is divided into two intervals (“high and low”). However, we did not obtain this result. We checked the possible spectral variation by dividing the spectra into two parts: the first part using the good time interval (GTI) for times shorter than $\sim$ 7000 s or longer than $\sim 1.25\times 10^{4}$ s, and the second part was extracted using the GTI between these two values. The spectra do not show any spectral variability.

Conservatively, we also extracted two spectra from Obs. 2: the first using the GTI in which the detected count rate was lower than 900 count $s^{-1}$ , and the second using a count rate higher than this amount. We did not find spectral variation in either observation (only the normalization parameter was different), therefore we present the spectra using the total EPIC-pn exposure time.

3 Data analysis

In our spectral analysis we used the 2.5-10.0 keV energy range. We excluded the soft energy range in which weak absorption features at $\sim$ 1.88 keV and $\sim$ 2.3 keV were observed in previous analyses and were referred to as being most likely caused by residual calibration uncertainties. We observed the feature at $\sim$ 1.88 keV in both observations (which could be Si XIII, or systematic residuals at the Si edge owing to a deficient calibration), but we did not observe the feature at $\sim$ 2.2-2.3 keV (which could be residuals of calibration around the instrumental edge of Au-M at 2.3 keV) (Díaz Trigo et al., 2012; Pintore et al., 2014). The absence of the weak absorption line at $\sim$2.3 keV in our analysis could be due to improvement in the calibration in SAS version 14.0.0. Obs. 1 and 2 in this paper correspond to observations 4 and 6 in Díaz Trigo et al. (2012).

We fit the spectral continuum with a blackbody component (bbodyrad in XSPEC) to model the thermal emission from the neutron star surface, plus a Comptonization component (compTT in XSPEC), which represents the emission attributed to a corona. The continuum of GX 13+1 is modified by significant total photon absorption by the material present in the line of sight of the observer. We used only one component (tbabs in XSPEC) to indicate both the Galactic absorption due to neutral hydrogen column and the variable absorbing material close to the source. The tbabs parameter NH was therefore let free in all fits.

Both observations show a broad and red-skewed emission line in the Fe Kα energy band (see Figure 2). In previous analyses (Pintore et al., 2014; Díaz Trigo et al., 2012), these features were fit either using a simple broad Gaussian or the relativistic diskline model. In this work, we fit the broad red-skewed line profile using both the nonrelativistic windline and the diskline model, and we compared the results obtained in each fit. The spectral analyses were performed using the XSPEC astrophysical spectral package version 12.09.0j (Arnaud, 1996).

In the same energy band in which the Fe emission line is apparent, we observed three narrow absorption lines in Obs. 1 and one in Obs. 2. These features are characteristics of dipping sources and were modeled with absorbed Gaussians, the lgabs model, which differs from the gabs model by the fact that it is not an exponentially multiplicative model (suitable for describing cyclotron resonance features), but a multiplicative model of the form $\propto(1-\text{Gaussian})$ (Soong et al., 1990). Because lgabs is a multiplicative model, it is not possible to obtain the EW of these lines through XSPEC. The addition of edges at 8.83 keV (Fe XXV) and 9.28 keV (Fe XXVI) did not improve the quality of the fit significantly. The best-fit parameters are shown in Table 3.

In Obs. 1, the three narrow absorption lines are identified as coming from highly ionized iron: Kα Fe He-like line (Fe XXV), Kα H-like line (Fe XXVI), and Kβ Fe He-like line. In Obs. 2, only the Kα Fe He-like absorption line is present. The addition of a narrow absorbed Kβ Fe H-like line did not significantly improve the fit.

Fitting the broad iron emission line

When we fit the Fe emission line with a simple Gaussian in Obs. 1, we obtained the line centroid energy $E$ equal to $6.872^{+0.030}_{-0.021}$ keV, $\sigma$ equal to $0.59^{+0.03}_{-0.05}$ keV, and the equivalent width (EW) equal to $226^{+142}_{-222}$ eV; in Obs. 2, we obtained the line centroid energy equal to $6.58^{+0.05}_{-0.05}$ keV, $\sigma$ equal to $0.64^{+0.04}_{-0.04}$ keV, and EW equal to $118^{+60}_{-87}$ eV. The centroid energies are different between the two observations. The EW of the lines are poorly constrained in both fits, and the broadening of the lines agrees within the errors between the observations at 90% confidence level. This constancy of the line width is in agreement with what was reported by Cackett & Miller (2013).

Comparing our results with Díaz Trigo et al. (2012), we obtained a higher Gaussian centroid energy in Obs. 1 (see Table 2). However, when we fit the spectra of Obs. 1 and 2 with the same total model used by Díaz Trigo et al. (2012) (see Table 2, rows 1 and 3), we obtained that all parameters, except for the $\sigma$ in Obs. 1, agree at 90% confidence level with the values found by Díaz Trigo et al. (2012). The $\sigma$ in Obs. 1 appears slightly narrower in our analysis.

Because Díaz Trigo et al. (2012) applied RDCTI correction in their analysis, the RDPHA correction applied by us could have led to higher centroid energies, larger broadness, and higher intensity of the emission line (see Pintore et al., 2014). However, we did not observe these differences in the line when we fit the data with the same total Model 1 used by Díaz Trigo et al. (2012) (see Table 2, rows 1 and 3). Therefore, the different centroid energy found by us in Obs. 1 (see Table 2, row 4, column 3) is likely due to the differences in the total model used to fit the continuum and the absorbed narrow lines. However, in Obs. 2 we did not observe a shift in the line centroid energy with respect to the fit performed by Díaz Trigo et al. (2012).

3.1 windline fit

When we fit the windline model to the residual excess in the iron line energy range, we obtained the best fit-parameters shown in Table 3; columns 4 and 5 show the best-fit parameters for Obs. 1 and Obs. 2, respectively. All fits were performed in order to leave the maximum number of free parameters.

Figure 3 shows the best spectral fit in the 2.5 to 10 keV energy range for Obs. 1 (left) and Obs. 2 (*right *) when the windline model was used to fit the residual excess in the iron line energy range. Figure 4 shows the unfolded spectrum in the energy range in which both the broad emission line and the narrow absorption features appear simultaneously.

The best fit was found in both observations with a 6.6 keV emission line, which corresponds to a $K_{\alpha}$ transition of the Fe XXI-XXIII ions. This line appears with an EW equal to $195_{-47}^{+26}$ in Obs. 1, and $130_{-24}^{+25}$ in Obs. 2 (see Table 2, 3). We obtained the highest value of the hydrogen column in Obs. 1, which is in agreement with Díaz Trigo et al. (2012), although different absorption column models were used by them (see Table 4 therein, where observations No. 4 and 6 correspond to Obs. 1 and 2).

In Obs. 1, the temperature of the electrons kTew and the velocity of the outflow were not constrained by the fit. To obtain consistent physical parameters, we therefore fixed these two parameters (kTew and $\beta$) with values similar to those found by the best fit in Obs. 2 (see Table 3).

In Obs. 1 and 2 we obtained $\tau_{w}>1$, $\tau>\tau_{w}$, temperature of the electrons in the wind kTew $\sim$ 0.6 keV, and $\beta\leavevmode\nobreak\ \sim\leavevmode\nobreak\ 10^{-2}$. The velocity of the outflow is equal to $2.25\times 10^{4}$ km s*-1* in Obs. 1 and ${2.01^{+0.18}_{-0.21}\times 10^{4}}$ km s*-1* in Obs. 2.

These outflow velocities are about ten times higher than the blueshifted absorption feature velocities found in previous analyses: Díaz Trigo et al. (2012) found blueshifts in GX13+1 between $\sim$ 2100 and 3700 km s*-1*; Ueda et al. (2004) indicated, based on Chandra HETGS, a bulk outflow plasma velocity of $\sim$ 400 km s*-1*; and Allen et al. (2016) likewise used Chandra HETG observations and found an outflow plasma velocity $>\leavevmode\nobreak\ 500$ km s*-1*.

3.2 diskline fit

When we fit the diskline model to the residual excess in the iron line energy range, we obtained the best-fit parameters shown in Table 3. All fits were performed in order to leave the maximum number of free parameters. Column 6 shows the best-fit parameters for Obs. 1 (Fit A, hereafter); columns 7 and 8 show the best-fit parameters for Obs. 2 when the energy of the iron emission line was not constrained (Fit B, hereafter), and when the energy line was constrained to the energy range of 6.4 to 6.97 keV (Fit C, hereafter). Figure 5 shows the unfolded spectrum (of Fits A and C) in the energy range in which both the broad emission line and the narrow absorption features appear simultaneously.

We found an inclination of the source ($i$) equal to $60_{-6}^{+4\circ}$ in Obs. 1, $75^{\circ}$ pegged at the hard limit and $61_{-5}^{+5\circ}$ in Obs. 2, in Fits B and C, respectively. All inclination values found are in the $60-80^{\circ}$ inclination range expected for GX 13+1 (Díaz Trigo et al., 2012). We checked the effect of the inclination on the best-fit parameters in both observations. For inclination values within the expected range the best-fit parameters do not change significantly.

The line emissivity as a function of the accretion disk radius ($r$) is an unknown function. For simplicity, the diskline model assumes that the line emissivity varies as $r^{B_{10}}$, and for the most part, set $B_{10}$ is equal to $-2$ (see Fabian et al., 1989). This emissivity index is expected where $r$ is approximately a few outer radii ($r_{c}$) of the corona. The best fit led to a power-law index $B_{10}$ of $\sim$ $-2.3$ for both observations: $B_{10}$ is equal to $-2.38^{+0.06}_{-0.06}$ in Obs. 1 and equal to $-2.25^{+0.20}_{-0.15}$ and $-2.40^{+0.12}_{-0.10}$ in Obs. 2 in Fits B and C, respectively.

In both observations the inner disk radius ($R_{in}$) of the reflected line component is equal to $\sim$ 10 gravitational radii ($R_{g}$): $R_{in}$ is equal to $9.6^{+0.9}_{-0.8}$ $R_{g}$ in Obs. 1 and equal to $14^{+4}_{-7}$ and $14^{+7}_{-4}$ in Obs. 2, in the Fits B and C, respectively.

The outer radius ($R_{out}$) in Obs. 1 was fixed to 1000 $R_{g}$, because the thawing of the parameters led to unphysical results: either $R_{out}$ becomes smaller than $R_{in}$, or $R_{out}$ assumes the value of $\sim$ 10.000 $R_{g}$. On the other hand, the $R_{out}$ parameter in Obs. 2 was left free because it appears better constrained by the fit.

We found an emission line energy of $6.674^{+0.040}_{-0.002}$ keV in Obs. 1, which corresponds to a Kα He-like Fe emission line. On the other hand, in Obs. 2 Fit B we found an anomalous emission line energy of $6.26^{+0.05}_{-0.04}$ keV, which represents a less physical energy of an iron line. The closest iron line energy is expected to be found at 6.4 keV, which corresponds to a Kα emission line from a neutral iron atom. In Fit C we therefore allowed the line energy parameter to vary only in a physical energy range, and in this case, the line appeared pegged at the lower limit. The line appeared with an EW equal to $256_{-73}^{+128}$ eV in Obs. 1, $119_{-25}^{+21}$ eV in Obs. 2 Fit B, and not well constrained, equal to $111_{-107}^{+18}$ eV, in Obs. 2 Fit C (see Table 2, 3).

Freezing the $R_{out}$ parameter in Obs. 2 Fit B to the same value found in Obs. 1 (i.e., to 1000 $R_{g}$) leads to a physical emission line energy value of $\sim$ 6.6 keV, but the inclination of the source assumes the value of $\sim$ 30*∘*, which is not consistent with the dips observed in this source. In this fitting we observed an increment of the $\chi^{2}$-red.

We also determined how different line emissivity indices affect the best-fit parameters in our analysis. We list our results below.

• $\bullet$

Forcing $B_{10}$ equal to zero, which is expected for $r$ ¡ $r_{c}$, led to less physical values of the source inclination ($i$ ¿ $90^{\circ}$) and an emission line energy of $\sim$ 6.8 keV in Obs. 1. In Obs. 2 (Fit B), the emission energy line shifts to a physical value of $\sim$ 6.7 keV, but all other parameters become unconstrained and assume unphysical values, that is, $R_{out}$ becomes smaller than $R_{in}$, and the source is found to be viewed face-on ($i$ $\sim$ $3^{\circ}$). All these fits led to an increment of the $\chi^{2}$-red.

• $\bullet$

Forcing $B_{10}$ equal to -3, which is expected for $r$ beyond $r_{c}$, where the coronal radiation intercepts the disk only obliquely, led $R_{out}$ to assume unphysical value if it was not frozen, and the emission line energy shifts to $\sim$ 6.6 keV in Obs. 1. In Obs. 2 we did not observe significant changes in the best-fit parameters. All these fits led to an increment of the $\chi^{2}$-red.

3.3 windline versus diskline fit

In both observations the best fit performed with the windline model was found for an emission line energy at 6.6 keV. Despite the evolution of the continuum between the two observations, the iron emission line remains remarkably constant. In contrast, when the excess in the energy range of the Fe K emission line was fit with the diskline model, different line energies were found: we found a $6.674^{+0.040}_{-0.002}$ keV emission line energy in Obs. 1, a $6.26^{+0.05}_{-0.04}$ keV emission line in Obs. 2 Fit B, and a 6.4 keV emission line pegged at the lower limit in Obs. 2 Fit C, when the line energy was constrained to a physical energy range (see Table 3). In Obs. 2, the line energy changes considerably in relation to the windline model, and it is low to be emitted by an iron atom in Fit B. The best fit in this case is not able to constrain the energy line parameter to a range of physical values, probably because of the large line broadening.

The reconstructed shape of the line also changes considerably between the two models. The line profile described by the diskline model (see Figure 5) appears double-peaked as a result of the Doppler effects on the accretion disk. The high disk inclination of the source produces a broader line in which the two peaks are very well separated. Aberration, time dilation, and blueshift by the fraction of the disk that is relevant together with the disk inclination make the blue horn brighter than the red horn (Fabian et al., 1989). In contrast to the relativistic line profiles, the windline line profiles (see Figure 4) show a broad single-peaked line, whose peak corresponds to the direct component of the line.

All fits with the diskline model determined the values of the source inclination in the range expected for GX 13+1, and in Obs. 1 and Obs. 2 Fit C the inclination parameter remains constant at 90% confidence level. The windline model shows no inclination effect on the line profile, the flow is spherically symmetric, with radial velocities, which means that there is no preferred angle between the flow stream and the line of sight, that is, observers from any direction see the same modified line, produced by a mean Doppler shift over all flow directions.

We obtained that not only the energy of the emission line was affected when different line models were used. It also slightly affected the continuum and the centroid energy of the narrow He-like Fe Kα absorption line (see Table 3 and the discussion below).

Of the continuum parameters, in Obs. 1 only the bbodyrad temperature remains in agreement at the 90% confidence level, when different line models are used. In Obs. 2, however, only the value of the hydrogen column remains in agreement at the 90% confidence level (see Table 3).

The EWs of the broad iron emission line in Obs. 1 and Obs. 2 are poorly constrained, which is partially due to the presence of the absorption lines that intercept the broad emission line, as previously noted by Díaz Trigo et al. (2012). The width of the narrow absorption lines are not constrained by the fits either. Because of the error bars, we therefore cannot infer any change in the EW of the emission line when different line models are used in an observation.

Cackett & Miller (2013, see Figure 3 therein), using the results of Díaz Trigo et al. (2012, see Table 4 therein), found no statistically significant correlation between the warm absorber column density and Fe emission line EW in GX 13+1. They obtained that a change in the measured emission line EW with increasing NH depends on the continuum parameters. That may explain the different values of EW among Obs. 1 and 2 at the 90% confidence level in the fits performed with the diskline model. This may explain the different values of EW among Obs. 1 and 2 at the 90% confidence level in the fits performed with the diskline model.

Narrow and absorbed Gaussians

When we fit the spectra with the total model used by Díaz Trigo et al. (2012), that is, using Gaussians with negative normalizations to fit the absorbed lines, we obtained that in Obs. 1 the EW ($30\pm 4$ eV) of the Fe He-like absorption line at $\sim$ 6.7 keV is slightly greater than the value found by Díaz Trigo et al. (2012) ($20\pm 5$ eV). For all other absorbed lines in Obs.1 and for the line in Obs.2, the EWs agree at the 90% confidence level with Díaz Trigo et al. (2012, see Table 3 therein). The values of the centroid energies of the narrow and absorbed Gaussians found in the two observation for the windline and diskline fits are compatible with the values found by Díaz Trigo et al. (2012) and Pintore et al. (2014, for the RDPHA correction, although a different observation was used by them) at the 90% confidence level.

However, comparing the line energies obtained in the fits performed with the diskline model with the fits performed with the windline model, we observed that in Obs. 1 the He-like Fe Kα line appears to be blueshifted when the diskline model is used to fit to the broad iron emission line, whereas in the windline fit, it was not possible to distinguish between a blueshifted, redshifted, or not shifted (affected) line because of the error bars of this parameter. The other two narrow lines, H-like Fe Kα and He-like Fe Kβ, showed the same ambiguity in the line shift in both observations for the windline and diskline fits.

When the narrow line indeed appeared to be shifted, the velocity of the warm absorber producing the shift in the line energy is around ten times lower than the velocity of the more central outflow (wind) where the iron emission line is formed in the windline framework. For all other cases, in which the ambiguity in the line shift is found, when any shift is present, the radial velocity of the outflow is around 100 times lower than the central outflow velocity.

3.4 Run-test and statistical assessment of the goodness-of-fit

Table 4 shows (see columns 6, 7 and 8) that all the different line models are statistically equivalent on the basis of the $\chi^{2}$ goodness-of-fit test. For example, in Obs. 1 the F-test gives a probability of chance improvement of 16$\%$, when instead of the windline model the diskline model is used to describe the data. In Obs. 2 the F-test gives a probability of chance improvement of 48$\%$ and 41$\%$ when Fits B and C are compared with the widnline fit, respectively. In order to distinguish this ambiguity between the models, we used a different approach to this problem.

Because our main goal is to determine whether the iron line profile is intrinsically asymmetric (i.e., the skewness coming from physical grounds) and to distinguish between the relativistic and nonrelativistic cases, we need a statistical test that takes the shape of the feature into account. To this aim, we used the run-test, also known as the Wald-Wolfowitz test (Barlow, 1989; Eadie et al., 1971), which inspects the residuals of the fit. If the model perfectly describes the data, the positive ($+$) and negative ($-$) residuals are expected to be randomly distributed around zero, and the number of runs, that is, the number of sequences of consecutive $+$ or $-$ residuals, is expected to be large (see, e.g., Redman & Rankin, 2009; Orlandini et al., 2012, for an application of the run-test in an astrophysical context).

The run-test gives the cumulative probability of obtaining by chance the number of observed runs (run-test probability, hereafter). For a good fit, that is, for a random distribution of the residuals, the run-test probability will be high. If, on the other hand, the fit does not describe the shape of the feature, the data will not be randomly distributed with respect to the fitting model. In this case, the run-test probability will be lower, meaning that there is a residual underlying trend.

It is important to mention that the run-test and the $\chi^{2}$-statistic are independent. The $\chi^{2}$-statistic does not depend on the ordering of the bins or on the signs of the residuals in each bin, while the run-test accounts for both the ordering and the signs, providing effectively additional information on the model goodness (Eadie et al., 1971).

Table 4 also shows the run-test probabilities in the 6.0 to 6.7 keV emission line energy range for the Gaussian (column 3), the windline (column 4), and the diskline (column 5) fits. For the Gaussian fit we obtained a run-test probability of $12.6\%$ in Obs. 1, and $55.2\%$ in Obs. 2. In the case of the diskline fit, we found a run-test probability of $14.4\%$ in Obs. 1, $43.0\%$ in Obs. 2 when the line energy was not constrained in the fit, and $32.4\%$ when the line energy was constrained in the 6.4 - 6.97 keV energy range. When we fitted the data with the windline model, this probability becomes $20.0\%$ in Obs. 1 and $56.9\%$ in Obs. 2. For all line models (gaussian, diskline, and windline) the hypothesis that the residuals are randomly distributed is accepted for a test at $5\%$ significance level.

4 Discussion and conclusions

The fluorescent iron line is produced by hard X-ray irradiation of a cold gas. However, where and how the asymmetry of the line is created differs between the two asymmetric Fe emission line models.

Because of the high disk inclination and the assumption that the line is created in the inner part of the accretion disk, the line profile appears to be double-peaked when it is fit with the diskline model. On the other hand, the line profile appears to be single-peaked when it is fit with the windline model because the red wing is produced by the indirect component of the line photons, which interact multiple times with the electrons in the outflow before escaping to the observer. Nonetheless, the different line model frameworks are both capable in terms of $\chi^{2}$-statistic to fit the broad Fe K emission line profiles.

The Compton bump in the $\sim$ 10 - 40 keV energy range is expected to be observed in the reflection scenario, and therefore in the relativistic line scenario (Fabian et al., 2000). However, the bump is expected in the windline framework as well. Laurent & Titarchuk (2007) used Monte Carlo simulations of the continuum spectrum emerging from a pure scattering wind, and their analytical description showed that when the incident spectrum is described by a hard power law ($\Gamma<1$) and the wind has physical parameters such that $\beta\lesssim 0.3$ and $\tau$ of a few, a photon accumulation bump at energies around $\gtrsim\leavevmode\nobreak\ 10$ keV and a softening of the spectrum at higher energies are observed. Because the number of photons is conserved, the down-scattered high-energy photons, which are removed from the high-energy part of the incident spectrum, are detected at lower energy. They pointed out that as the optical depth increases, more prominent bumps are formed in the outflow. In general, the shape of the emerging continuum depends on the mean number of scatterings suffered by the photons in the outflowing plasma (which is a combined effect of the wind optical depth $\tau$ and the wind velocity), and on the incident spectrum shape.

The windline model allows determining the possible outflow physical parameters. For both observations, we found the

• •

optical depth of the wind ($\tau_{w})>1$;

• •

temperature of the electrons in the wind (kT${}_{ew})$ of  $\sim$ 0.6 keV, which is the possible outflow temperature (Laming & Titarchuk, 2004); and an

• •

outflow velocity ($\beta$) of  $\sim\leavevmode\nobreak\ 10^{-2}c$ (see sec. 3.1).

The outflow velocities determined by the windline model are ten times higher than the velocity determined by the blueshift velocities of the absorption features found in previous analyses. The different velocities can still be explained considering that the emission and absorption lines are produced in different regions: the emission line might be produced in a more internal region of the source, at the bottom of the wind or outflow (in an inner shell) where the velocity is high; on the other hand, the resonance absorption lines are produced in a more distant region from the compact object by the interaction of photons with highly ionized species present in a cylindrical absorbing plasma around the source, driven by outflows from the outer regions of the accretion disk. Therefore, the presence of the outflow and its velocity variation along the cloud radius may explain the simultaneous observation of narrow absorption and broad emission lines in sources with high inclination, such as GX13+1. The absence of absorption lines in the internal outflow, responsible for red-skewing the line, can be explained by relatively high temperature of the outflow, which is on the order 0.3-0.6 keV, see the calculations of the temperature structure of the outflow in Laming & Titarchuk (2004).

The wind can be launched from 3 to 104 Rg. In the windline model, the wind is considered as a spherical shell of internal radius of 100 Rg and external radius of 120 Rg, so that the wind starts at 100 Rg. These values have no great influence on the emitted spectrum as long as the inner radius of the wind is far from the central compact object, where outflowing plasma is not fully ionized and strong gravitational effects are not expected.

The critical mass outflow rate Ṁ is given by $4\pi r^{2}\rho(r)m_{p}\beta c=1.7\beta\tau 10^{21}$ g/s = $2.6\beta\tau 10^{-5}$ M*⊙* yr*-1* (where mp is the mass of the proton) (see Titarchuk et al., 2007, appendix E). For the fit values we found, $\tau\approx 2$ and $\beta\approx 0.01$, we find Ṁ = 3 $\times$ 1019 g/s, which is more than ten times what was previously found (see Díaz Trigo et al., 2012, and references therein). However, this is consistent with the velocity difference between the inner and outer parts of the flow, outer parts where absorption lines take place. The given mass outflow rate is an upper limit because it was calculated considering that the wind is spherical throughout the system. A lower mass outflow rate may be emitted by a wind that partially covers the system. The mass outflow rate was self-consistently calculated by Laming & Titarchuk (2004) and reproduced by Laurent & Titarchuk (2007) using Monte Carlo simulations.

When the different asymmetric line models were fit to the same observation, slight differences in the continuum were observed, and the EW remained the same at the $90\%$ confidence level. The fits with the two different line models to the complex spectrum of the dipping source GX13+1 also introduced differences in the energy of the fluorescent Fe emission line and in the energy of the Kα Fe He-like (Fe XXV) absorption line.

The best fit with the diskline model in Obs. 2 gives an emission line energy equal to 6.26${}^{+0.05}_{-0.04}$ keV, which is lower than the energy expected from photons coming from a neutral or ionized iron atom. It is, for example, 2.8$\sigma$ from the Kα energy line emitted by a neutral iron, at 6.4 keV. We checked if the total continuum model could affect the emission line energy in this observation, but fitting the continuum and the absorbed Gaussians as in Model 1 , that is, with tbabsedge1edge2*(diskbb+bbodyrad+diskline+gaus1) model in XSPEC, did not lead to a significant change in the emission line energy (in this case, the energy of the line is equal to $6.28^{+0.04}_{-0.10}$). To obtain a physically consistent fit, we constrained the line energy parameter to the expected energy range. In this case, although the energy line appears pegged at the lower limit, the geometric parameters of the disk are found in good agreement with the fit in Obs. 1.

In the fits using the windline model, the outflow parameters, and consequently the properties of the inner portion of the outflow, are approximately constant in the two observations. The emission line generated in the outflow remains remarkably constant although the underlying continuum evolves. However, as in the diskline fit to Obs. 2, in Obs. 1 two widnline parameters ($kT_{ew}$ and $\beta$) had to be frozen to obtain a physically consistent fit. The constant properties of the inner outflow may be explained by the intrinsic variability of the accretion flow.

The ambiguity in terms of $\chi^{2}$-statistic test between the emission line models was expected and confirmed in the spectral fits presented in this paper. The $\chi^{2}$-statistic was not sufficient to lead to an unambiguous statement on the line profiles because it squares the differences between data and model and consequently looses the information about the form of the residuals in the emission line energy range. Therefore, we used the statistical run-test to take the shape of the residuals into account. We tried to break the degeneracy between the relativistic and nonrelativistic line models in the GX 13+1 spectra because it is found in the literature that this source presents clear evidence of both disk-wind (outflow) and strong, broad, and skewed iron emission lines. However, this source has a complex spectrum, with narrow absorption iron lines in the same energy range as the broad emission line. The absorbed iron lines complicate a distinction between the two line models because the run-test is not performed in the entire energy range of the emission line.

For the two observations analyzed in this paper, we were not able to reject the hypothesis that the residuals are randomly distributed at the $5\%$ significance level, and we cannot rule out one of the line models. We obtained a higher cumulative probability of observing by chance the number of runs around the fitting line for the windline spectral fit. However, when we modeled the continuum using the total Model 1, we obtained the opposite, that is, a higher cumulative probability for the diskline model. Therefore, we conclude that for the observations analyzed in this paper the broad emission Fe line profiles can be described as a signature of a wind or outflow and also by GR effects in terms of $\chi^{2}$-statistic and run-test.

The statistical run-test may allow a better assessment of the goodness of the nonrelativistic (windline) versus the relativistic (diskline) line profiles. To break the degeneracy between the relativistic and nonrelativistic broad iron lines, a study considering several other NS LMXB sources, containing strong and broad iron emission lines in their spectra could be performed. It could be performed in a perfect scenario, considering mainly sources with a more straightforward total spectrum, and using observations that are not affected by pile-up.

In the windline model, the large amount of mass outflowing with high velocity implies that the Earth observer should see broad and redshifted iron lines formed in the outflow. The velocity of the outflows in NS LMXBs given by the model could, for example, be compared with velocities found by P Cygni profiles in optical, UV, and X-ray because these lines indicate outflows or disk-wind outflows (e.g., Brandt & Schulz, 2000; Schulz & Brandt, 2002). If the windline model is found to be more appropriate to fit broad and skewed iron lines in NS LMXBs, the observation of such profiles will be of utmost importance for understanding outflowing in these systems; this is also an important tool for studying the inflowing-outflowing connection, as stated by Díaz Trigo et al. (2012).

We emphasize that in addition to the moderate spectral resolution of the X-ray detectors that are available today, a careful analysis of the Fe K emission line can bring important information for understanding the physics behind the asymmetric Fe fluorescent emission lines and their implication on the neutron star physics. The next generation of X-ray observatories loaded with instruments with high spectral resolution will certainly improve the understanding of the asymmetric Fe line formation in accreting compact objects.

We conclude this paper by also pointing out the importance of a timing variability study to constrain the physical process that leads to the broad and red-skewed iron line profiles. It can provide important and additional information for breaking the ambiguity between the relativistic and nonrelativistic line models.

Acknowledgments

T. Maiolino acknowledges the financial support given by the Erasmus Mundus Joint Doctorate Program by Grants Number 2013-1471 from the agency EACEA of the European Commission, the CNES/INTEGRAL, and the CNR/INAF-IASF Bologna. T. Maiolino would also like to thank Clément Stahl and Lorella Angelini for their valuable comments that contributed to this manuscript. We would like to thank Maria Díaz Trigo for sharing details of previous XMM-Newton Epic pn data reprocessing. Finally, we thank the anonymous referee for their critical comments that considerably improved the content of the paper.