Highlight
Free Access
Issue
A&A
Volume 623, March 2019
Article Number A73
Number of page(s) 15
Section Stellar structure and evolution
DOI https://doi.org/10.1051/0004-6361/201834801
Published online 07 March 2019

© ESO 2019

1. Introduction

In the usual scenario of magnetic dipole braking in vacuum, the observed secular lengthening of a pulsar spin period, typically by 3 s every 108 yr, is the consequence of the torque exerted by the magnetic field on the rotating neutron star. Despite its simplicity, the model provides useful estimates of the evolutionary state of a pulsar, namely its available rotational power, characteristic age, and surface dipolar magnetic field strength (Ostriker & Gunn 1969). For the isolated neutron stars (INSs) in our Galaxy, dipolar field estimates span across five orders of magnitude (Manchester et al. 2005). Among the sources with the highest values are those known as magnetars (see Turolla et al. 2015; Kaspi & Beloborodov 2017; Coti Zelati et al. 2018; Esposito et al. 2018, for recent reviews).

Magnetars are usually observed through their violent bursts of high energy; they slow down at a much faster rate than normal pulsars (up to 1 s every 60 years for the extreme case of the soft-gamma repeater SGR 1806-20). According to the most favoured interpretation (Thompson & Duncan 1995, 1996), their complex phenomenology, transient behaviour, and bright quiescent X-ray luminosity, much in excess of that from spin down, can be explained by crustal and magnetospheric effects provoked by the decay and rearranging of the enormous stellar magnetic field (Goldreich & Reisenegger 1992; Pons et al. 2009; Beloborodov & Li 2016). As a result of field dissipation and braking, it is expected that an evolved magnetar (≳105 yr) will be less active than a young one, and will have a longer spin period and higher surface temperature than an ordinary pulsar of similar age (e.g., Perna & Pons 2011; Turolla et al. 2011).

The group of X-ray thermally emitting INSs discovered by ROSAT and known as the Magnificent Seven (M7; see Haberl 2007; Kaplan 2008; Turolla 2009, for reviews) may have evolved from such a channel of pulsar evolution (e.g., Heyl & Kulkarni 1998; Kaplan & van Kerkwijk 2009; Popov et al. 2010; Viganò et al. 2013). They consist of a rather unique local group of middle-aged (∼105−106 yr) cooling neutron stars, displaying similar low blackbody temperatures (kT ∼ 45−100 eV), long spin periods (P ∼ 3−17 s), and moderately strong dipolar magnetic field strengths (Bdip ∼ few × 1013 G). Unlike other X-ray pulsars, their emission is purely thermal with no sign of magnetospheric activity, and it is believed to originate directly from the neutron star surface.

The source RX J1605.3+3249, as the third brightest among the M7 INSs (Motch et al. 1999), was consequently visited by XMM-Newton in several occasions during the early years of its science operations (see van Kerkwijk et al. 2004; Haberl 2007, for details on the past investigations of the source). However, these early observations were not deep enough to allow the detection of the neutron star spin signal. Ensuing a visibility gap of six years, a 60 ks observation performed in 2012 finally provided a candidate spin period for the INS (Pires et al. 2014). The amplitude of the shallow and strongly energy-dependent periodic signal was detected close to the sensitivity limit of the data; only the harder portion of the source spectrum (roughly, 30% of all source events at energies above 0.5 keV) was found to show a significant modulation at the 4σ level. Nonetheless, a joint timing analysis around the detected signal at P ∼ 3.39 s, connecting the 2012 dataset with early XMM-Newton observations of the source, hinted at an unprecedentedly high value of spin down. The inferred dipolar magnetic field of Bdip ∼ 7.4 × 1013 G overlaps the magnetar range and, if confirmed, could rank the highest in the group.

In addition, the analysis of the then available XMM-Newton EPIC data on the source confirmed the evidence of a complex, multi-temperature, energy distribution and the presence of absorption features. In Pires et al. (2014) we described the spectrum of the source using a two-component blackbody model of temperatures 60 eV and 110 eV, superposed by one or two Gaussian absorption features at around energies 400 eV and 860 eV; best fits were found assuming the Galactic column density in the direction of the source, 2.4 × 1020 cm−2 (Kalberla et al. 2005, see Sect. 3.2 for details). These results motivated us to investigate the INS further.

We were granted a large programme of observations with XMM-Newton (programme ID: 76446) for a total duration of 310 ks and four satellite visits in the AO14 observing cycle. We report here the results of this observational campaign. The paper is organised as follows. In Sect. 2 we describe the XMM-Newton observations and the data reduction. Our analysis and results are presented in Sect. 3. The implications of our results are discussed in Sect. 4, with particular emphasis on the properties and recent work on the group of M7 INSs. Our conclusions and the summary of the results are presented in Sect. 5.

2. Observations and data reduction

The XMM-Newton observatory (Jansen et al. 2001) targeted the INS RX J1605.3+3249 (hereafter J1605) on four occasions between July 2015 and February 2016, using EPIC as the prior instrument for the investigation. Table 1 contains information on the scientific exposures and instrumental configuration of the EPIC-pn (Strüder et al. 2001), EPIC-MOS (Turner et al. 2001), and RGS (den Herder et al. 2001) detectors.

Table 1.

Log of the XMM-Newton AO14 and 2012 observations of RX J1605.3+3249.

We included in the timing and spectral analysis (Sect. 3) the past 2012 observation of the source (obsid 0671620101; Table 1). For the spectral analysis of RGS data (Sect. 3.2.2), we also included archival XMM-Newton observations of the source dating back to 2002 that were not severely affected by background flares (see Table 2). We applied as criterion a minimum net exposure of 10 ks to include observations in the RGS analysis. The archival EPIC observations of J1605 performed before the visibility gap (analysed and discussed in Pires et al. 2014, see also references therein) were not included due to the high background level and the heterogeneous science operating modes and optical blocking filters that were then adopted (see also Sect. 3.2.2). All archival observations were processed and analysed consistently with the data from the large programme.

Table 2.

Summary of archival XMM-Newton observations of RX J1605.3+3249 used in the RGS spectral analysis.

The main trigger behind our programme was to confirm the candidate spin period and, by means of a precise timing solution, measure the neutron star’s spin-down rate. This is better achieved through a well-sampled ephemeris; incoherent methods result into much less accurate spin-down determinations, while a two-dimensional periodicity search (like that performed in Pires et al. 2014) suffers from a large parameter space and a large number of independent trials that have to be covered in the case of observations that are years apart. As a result, the pulsar best (P,) solution is determined at a low confidence level.

Therefore, the time intervals between the four observations in AO14 (of 5, 25 ± 3, and 175 ± 8 days) were carefully chosen to coherently connect them in phase with three past pn observations of the source, for a total time span of ∼5160 days (Sect. 3.1). Assuming the properties of the periodic signal as detected in the 2012 observation (Pires et al. 2014 and Sect. 3.1), we considered for the feasibility a total count rate of 0.8 s−1 in the three EPIC cameras for source photons with energy above 0.5 keV; the individual exposure times were required to ensure a significant detection of the 2012 signal while allowing for moderate changes of pulsed fraction, pf = 2.5%−5% (within ±1σ). Likewise, we chose to operate the MOS and pn cameras in large-window (LW) and full-frame (FF) imaging modes to provide sufficient time resolution (0.9 s and 73.4 ms, respectively) to measure the 2012 modulation. We adopted the thin filter for both instruments, given its better response at soft X-ray energies.

We performed standard data reduction with SAS 15 (xmmsas_20160201_1833-15.0.0) using up-to-date calibration files and following the analysis guidelines of each instrument1. We processed the EPIC exposures using the SAS meta tasks epchain and emchain and applied default corrections. For RGS, we used the SAS routine rgsproc to process the raw data files and create masks for the source and background regions.

Background flares were registered occasionally during the AO14 observations, usually lasting less than a few kiloseconds. The percentages of good time intervals (GTIs), filtering out periods of high background activity, are shown in Table 1 for each scientific exposure and observation. Standard count rate thresholds were adopted for pn and MOS; for RGS we used the background count rate on CCD 9 and applied a threshold of 0.1 s−1 to filter the event lists with the SAS task rgsfilter. On average, data loss is small: 2% for MOS and RGS and 7% for pn. The observation that was affected most severely by flares was performed in August 2015, with a percentage of data loss between 9% and 19% depending on the camera. The total net exposures per camera are 303 ks (MOS1), 312 ks (MOS2), 314 ks (RGS1), 312 ks (RGS2), and 290 ks (pn).

For the analysis of EPIC data, we filtered the event lists to exclude “bad” CCD pixels and columns, and to retain the pre-defined photon patterns with the highest quality energy calibration, namely single and double events for pn (pattern ≤4), and single, double, triple, and quadruple for MOS (pattern ≤12). The source centroid and optimal extraction region, with typical sizes of 140″, were defined with the SAS task eregionanalyse in the 0.3–1.35 keV energy band for each EPIC camera and observation. Background circular regions of size 60″–100″ were defined away from the source, on the same CCD of the target whenever possible.

The detected source count rates, hardness ratios, the pile-up-relevant count rates per frame for each camera and observation, and other parameters based on a maximum likelihood fitting are listed in Table 3, with nominal 1σ statistical uncertainties. The parameters are determined with the SAS task emldetect on images created for each camera, observation, and energy band (only the combined EPIC results are shown; the X-ray emission of the source is compatible with the background level at energies above 2 keV). For comparison, we also list the source parameters as determined from the 2012 EPIC exposure.

Table 3.

Parameters of RX J1605.3+3249, as extracted from the EPIC images in the AO14 and 2012 observations.

Following the guidelines of Jethwa et al. (2015), we ensured that the pile-up levels at aimpoint were within tolerant guidelines for both EPIC detectors. Based on the source spectrum and on the number of counts per frame in each camera (listed in Table 3 for direct comparison with Fig. 5 in Jethwa et al. 2015), we estimate that the percentage levels of spectral distortion and flux loss were around 1.5% and 3.5% for pn, and 0.3% and 1% for MOS.

Overall, the source properties are consistent between epochs since 2012: potential discrepancies can be asserted to the cross-calibration uncertainties between the EPIC instruments and to different background levels. In Sect. 3.2 we investigate possible flux and spectral variations of the INS in detail.

We used the SAS task eposcorr to refine the astrometry by cross-correlating the list of EPIC X-ray source positions with those of catalogued near-infrared (2MASS; Skrutskie et al. 2006), optical (GSC 2.3.2; Lasker et al. 2008), and X-ray (Chandra; Evans et al. 2010) objects lying within 15′ of J1605. The results that yielded the least astrometric errors were obtained when cross-correlating the X-ray sources with a number of around 50 optical counterparts present in the field of view. Small positional offsets in right ascension and declination were consistently detected for all catalogues under study and are also shown in Table 3. The astrometrically corrected EPIC source positions in all four observations are consistent with each other.

Finally, we verified the statistics of the EPIC light curves for general trend variability. The light curves, binned into 600 s to 1200 s intervals, were corrected for bad pixels, dead time, exposure, and background counts with the SAS task epiclccorr. All 2012/AO14 exposures are consistent with a constant flux.

3. Analysis and results

3.1. Timing analysis

For the timing analysis we used a test (Buccheri et al. 1983) applied directly on the times of arrival of the pn and EPIC (pn+MOS) events to search for periodic signals. The times of arrival of the pn and MOS photons were converted from the local satellite to the solar system barycentric frame using the SAS task barycen and the astrometrically corrected source coordinates in each camera and observation (Table 3).

Apart from the 2012 observation, the only other XMM-Newton datasets suitable for timing analysis are a pn LW exposure performed in 2003 (0157360401, ∼33 ks) and one performed in 2002 in timing (TI) mode (0073140501, ∼30 ks); these observations were used to derive an estimate of the spin-down rate of the source in combination with the 2012 dataset (Pires et al. 2014). The other archival pn observations of J1605 were either operated with the thick filter (hence reducing count rates by more than a factor of two) or severely affected by background flares. Likewise, the time resolution of 2.6 s of the archival MOS observations of J1605, all performed in FF mode, is not sufficient to detect the 3.39 s modulation.

The re-analysis of the 2012 observation with SAS 15 and up-to-date calibration files yields similar results to those reported in Pires et al. (2014). The only statistically significant modulation, at ν ≡ ν2012 = 0.2951709(14) Hz, results when the search is restricted to source photons with energy above 0.5 keV. We also verified that the periodic signal is always present at the same frequency within the errors, independently of the exact details of the processing of the raw event file (e.g., included calibration files, SAS version, randomisation in energy within a PI channel, event filtering, or randomisation in time within the sampling detector time). The measured fluctuation of the power of the statistic at ν2012,(ν2012) ∼ 35−49 is consistent with that expected from a sinusoidal modulation of amplitude pf = (4.1 ± 0.9)% (see, e.g., Pavlov et al. 1999).

We first analysed each of the four AO14 observations individually. Taking into account the typical spin-down rate of the M7 INSs (see, e.g., Kaplan & van Kerkwijk 2009, and references therein), we looked for significant signals in a 5 × 10−4 Hz range around the 2012 frequency, adopting a resolution of 0.1 μHz (oversampling factor of at least 10). The number of statistically independent trials in each search, which depends on the frequency range and on the total duration of the observation, is typically between 30 and 60. Assuming the usual scenario of magnetic dipole braking in vacuum, the tests allow for maximum braking with respect to the 2012 signal of , which corresponds to that exerted by a maximum dipolar magnetic field of Bdip ≲ 4 × 1014 G at the equator.

Due to the energy-dependent nature of the 2012 signal (see Pires et al. 2014, for details), we carried out tests in various energy bands, also varying the size of the source extraction region and other parameters of the search (e.g., the included photon patterns and other details of the processing of the raw event files). We tested seven energy ranges in the soft (0.15–0.5 keV, 0.2–0.5 keV, and 0.3–0.5 keV), hard (0.5–1.35 keV), and total (0.15–1.35 keV, 0.2–1.35 keV, and 0.3–1.35 keV) energy bands; source counts – between (0.18−5.5) × 105 pn and (0.4−7.2) × 105 EPIC photons – were extracted from circular regions of radius 10″, 25″, 50″, 100″, and 120″ around the source position in each camera and epoch. Altogether, 1680 tests were conducted in the pn/EPIC datasets of the four observations.

The summary of the results is presented in Fig. 1. In each plot we show the statistic as a function of trial frequency in each AO14 observation (pn/EPIC datasets), for different energy bands, and taking the extraction region of radius 100″ as an illustrative example. The Z2-test performed on the 2012 pn observation and respective peak at (ν2012) ∼ 50 (out of scale) is shown in the background of the four plots for comparison (0.5–1.35 keV). No significant power (above 3σ − 4σ) in the searched frequency range is consistently detected in these tests. The inclusion of higher harmonics (or different photon patterns and event filters) in the tests does not affect the results. Therefore, unless significant changes of pulsed fraction have taken place, 4σ upper limits2 of 1.33(6)%, 1.74(8)%, 1.84(8)%, and 1.80(8)% in the total energy band in each AO14 observation rule out the 2012 candidate period at a high confidence level.

thumbnail Fig. 1.

EPIC and pn searches around the 2012 periodicity. The frequency range is ν = 0.2948−0.2953 Hz. The periodogram in the background (dashed outline) shows the 2012 result: ν2012 = 0.2951709(14) Hz and (ν2012) ∼ 50 (obsid 0671620101, pn, 0.5–1.35 keV). The four plots show for each AO14 observation the results of tests conducted in seven different energy bands for an extraction region of 100″ (see text). The dotted and dashed horizontal lines show the 2σ and 3σ confidence levels for the detection of modulations in each observation.

We next searched the new data, especially the longest and most sensitive “201” observation for other significant modulations in the full frequency range allowed by the resolution of the EPIC cameras (blind searches). To this end, the times of arrival of the pn and MOS events (∼2.5 × 105 counts) were analysed together in the ν = 0.01−0.56 Hz frequency range; for pn, timing searches were extended to higher frequencies, up to ∼6.8 Hz3. The adopted frequency step was Δν = 2.5−5 μHz (oversampling factor of 3) and the number of independent trials were (4−8) × 105 and (3−6) × 104 in the pn and EPIC searches, respectively. A total of seven narrow energy bands (with widths between 100 eV and 600 eV) and three wide bands (600–850 eV), defined within the 0.15–1.35 keV range according to the source’s signal-to-noise ratio S/N, were defined for these searches. We tested three different sizes of source extraction regions (25″, 50″, and 120″) and included all valid photon patterns in the EPIC searches. For pn, we restricted the event lists to include only single and double photon patterns. In total, 240 EPIC and pn blind searches were performed in the observations of the large programme.

No significant periodic signal resulted from the analysis. In Table 4 we list the most constraining 3σ upper limits from the EPIC and pn searches for three wide energy band intervals (soft, hard, and total; see the table caption for details).

Table 4.

Upper limits on pulsations from the timing analysis.

3.2. Spectral analysis

3.2.1. EPIC data

The analysis of the EPIC data is based on source and background spectra extracted from regions as described in Sect. 2, together with the respective response matrices and ancillary files created for each of the EPIC cameras and observation. In accordance with the guidelines and calibration status of the instruments, we restricted the spectral analysis to GTI-filtered photons with energies between 0.3 keV and 1.35 keV (beyond which the source S/N becomes insignificant). Including the 2012 observation, the analysed dataset comprises 15 spectra and over 1.2 × 106 counts (0.3–1.35 keV), of which 1.8% can be ascribed to the background.

The pile-up level is negligible in the MOS exposures (Sect. 2); for pn we applied a correction in the redistribution matrix files with the SAS task rmfgen to minimise flux loss and spectral distortion4. The correction is calculated directly from the frequency and spectrum of the incoming photons and has the advantage of keeping events from the central PSF area in the spectral analysis, which would have to be otherwise discarded. We verified that the results using this approach agree well with those when the spectra are extracted from regions with a 10″−15″ excised core, while avoiding up to 30% of data loss in the pn exposures.

The energy channels of each spectrum, which are by construction 5 eV wide, were regrouped to avoid a low (< 30) number of counts per spectral bin. Due to the brightness of the source and the good statistics of each individual spectrum this has an effect only at the high-energy side of the analysis where the source signal becomes dominated by the background. At lower energies (≲700 eV) the spectrum oversamples the instrument resolution of the EPIC cameras by a factor of up to 20. We ensured nonetheless that oversampling did not influence the results of spectral fitting; specifically, we checked for consistency where oversampling was kept within a maximum factor of 3.

To fit the spectra we used XSPEC 12.9.0n (Arnaud 1996). Unless otherwise noted, the fit parameters were allowed to vary freely within reasonable ranges. The photoelectric absorption model and elemental abundances of (Wilms et al. 2000, tbabs in XSPEC) were adopted to account for the interstellar material in the line of sight. Due to the low absorption towards J1605, the choice of abundance table and cross-section model does not significantly impact the results of the spectral fitting.

The 15 spectra were first fit individually to check the agreement between the instruments and epochs. The exercise showed the expected small percentage of cross-calibration uncertainty between the EPIC detectors (Read et al. 2014). Although variations from pointing to pointing for a given instrument are formally significant with respect to a constant value, the relative error is still smaller than the absolute discrepancies between the cameras within a given epoch. To account for this uncertainty, we allowed for a renormalisation factor in XSPEC and fitted the spectra of the pn and MOS cameras simultaneously. We checked that the results of the simultaneous fits were consistent with the weighted means of the individual measurements.

We then checked if consistent results are obtained when fitting a single “stacked” spectrum, which converges to the best parameter values much faster in XSPEC than the simultaneous fits. The stacked spectra are produced with the SAS task epicspeccombine, taking into account the responses, background, and ancillary files of the individual exposures. For pn, the stacking approach leads to inconsistent results that do not match either the corresponding weighted mean values for the camera within the errors, or the results of the simultaneous fits5. On the other hand, the results of the combined MOS spectra agree well with those of the simultaneous fits. We thus adopted a stacked spectrum only for MOS to avoid introducing biased results in the spectral analysis.

Next, we tried to find a model that closely describes the X-ray spectral energy distribution of the source (see Pires et al. 2014, for details). We list in Table 5 the results of the fit of a double-blackbody model (2bb)6 with a broad Gaussian absorption feature. In the fitting procedure we restricted the energy of the line between 0.3 keV and 1.35 keV and its Gaussian σ between 0 eV and 200 eV; the column density was varied between NH = 0 cm−2 and 5 × 1021 cm−2, while the temperature of the blackbody components can assume values between 5 eV and 500 eV. For each observation, we fitted the model to the pn spectrum, to the combined MOS1 and MOS2 spectrum (labelled “MOS” in Table 5 and the two subsequent tables), and to the pn and MOS spectra simultaneously (labelled “EPIC” in Table 5 and the two subsequent tables). Finally, we performed multi-epoch simultaneous fits of all pn (5 spectra), MOS1/2 (10 spectra), and EPIC (6 spectra, comprising 5 pn and one combined MOS) exposures.

For each fit in Table 5, we list the reduced chi-square () and its null-hypothesis probability (NHP in %), the column density NH in units of 1020 cm−2, the temperature of the cold and hot blackbody components in eV, the radiation radii and of each component (assuming a distance to the source of d ≡ d300 = 300 pc; e.g., Tetzlaff et al. 2012), the central energy of the absorption line (when constrained) or the corresponding 3σ upper limits, its Gaussian σ, and equivalent width EW (all in eV), and the observed model flux in the energy band 0.2–12 keV, fX, in units of 10−12 erg s−1 cm−2. The model provides acceptable values in each epoch, with NHP between 14% and 88%. The somewhat large chi-square values and worse fit quality (NHP <  1%) of the multi-epoch fits could not be improved by the inclusion of an additional model component.

In the left column of Fig. 2 we plot the results of Table 5 as a function of time, with 1σ errors. The best-fit parameters per instrument are consistent between pointings despite the systematic differences between the detectors. The column density is constrained and for the MOS and EPIC fits agrees within errors with the Galactic value in the direction of the source, cm−2 (e.g., Kalberla et al. 2005; Willingale et al. 2013). The pn camera measures twice as much absorption, and 7%–9% softer temperatures with respect to MOS; the observed model flux is also 6% lower. If NH is fixed to the Galactic value, the temperature of the two components and the observed model flux typically agree within 2.5%. Considering the simultaneous EPIC fit as an effective average between the detectors, the best-fit parameters are well constrained within ranges NH  =  (2.4−4)  ×  1020 cm−2, eV, eV, km, km, ϵ <  (340−420) eV, σ  =  (100−120) eV, EW <  (150−200) eV, and fX  =  (6.4−6.9)  ×  10−12 erg s−1 cm−2. The model flux corrected for absorption (unabsorbed) is FX  =  (1.2−1.7)  ×  10−11 erg s−1 cm−2.

thumbnail Fig. 2.

Source best-fit parameters as a function of MJD (EPIC analysis; see text and Tables 5 and 6 for details). The model fit to the data is that of absorbed double-blackbody (left) and fully ionized hydrogen neutron star atmosphere (B = 1013 G, M = 1.4 M, and R = 10 km; right), modified by a broad Gaussian absorption line. The plots show in each instrument and observation (see legend): panel a: column density, temperature of the (panel b) cold and panel c: hot blackbody components, panel d: effective temperature of the neutron star (unredshifted), panel e: central energy of the absorption line, panel f: its Gaussian sigma, and panelg: observed model flux in the 0.2–12 keV energy band. The total Galactic NH value is shown by the solid black line in plots (panel a). Horizontal lines show the results of the simultaneous fit of spectra over the five pointings, with 1σ confidence levels comprised by the shaded areas.

Table 5.

Results of the best-fit double-blackbody model with a Gaussian absorption line (per observation and camera).

Table 6.

Results of the best-fit fully ionized neutron star atmosphere model with a Gaussian absorption line (per observation and camera).

Alternatively, the dataset can be fit equally well (NHP ∼10%−80%) by a fully ionized neutron star hydrogen atmosphere model (nsa in XSPEC; Pavlov et al. 1995; Zavlin et al. 1996), again modified by a broad Gaussian absorption line (Table 6). We tested nonmagnetised (B <  108 G) and magnetised models with magnetic field values of B = 1012 G and B = 1013 G; in the fitting procedure, the neutron star mass and radius were at first fixed at the canonical values, M = 1.4 M and R = 10 km, and then allowed to vary to check for an improved fit. In Table 6 we list the results for a canonical neutron star with B = 1013 G, which is the model that in most cases gave the highest NHP for each dataset. The unredshifted model effective temperature Teff in K, the distance d in pc, the parameters of the line (ϵ, σ, EW), and the observed model flux fX, are also listed in Table 6.

In the right column of Fig. 2 we plot the best-fit nsa parameters as a function of MJD (assuming for all epochs the results of the B = 1013G fits with M = 1.4 M and R = 10 km). In comparison with the 2bb model, the systematic differences between the detectors persist; however, the measurements differ by a much smaller percentage (13% in NH, ≲2% in Teff, and 3% in fX) and the overall consistency between epochs and instruments is improved (see overlapping shaded areas). The central energy of the line is well constrained, and also narrower than in the 2bb case. Again considering the best-fit results of the EPIC fits, we have NH = (2−3) × 1020 cm−2, Teff = (5.7−5.8) × 105 K, B = 1013 G, d = (124−129) pc, ϵ = (360−400) eV, σ = (85−100) eV, EW = (55−60) eV, and fX = (6.5−7.0) × 10−12 erg s−1 cm−2. The unabsorbed flux of this model is measured in a range similar to that of the 2bb model, FX = (1.2−1.4) × 10−11 erg s−1 cm−2. In Figure 3 we show this best-fit model folded to the EPIC dataset with residuals.

thumbnail Fig. 3.

Results of EPIC spectral fitting (see text and Table 6 for details). We show the 5 pn and stacked MOS spectra (grey and magenta data points, respectively), fitted simultaneously by a fully ionized hydrogen neutron star atmosphere model with B = 1013 G, Teff = (5.16 ± 0.17) × 105 K, and a broad Gaussian absorption line of σ = 91.4 ± 0.4 eV at eV (black and red solid lines).

The neutron star distance derived from the nsa fits, d ∼ 110−130 pc, is rather small in comparison to the range expected for the source, d ∼ 300−400 pc (Posselt et al. 2007; Tetzlaff et al. 2012). While the fit is insensitive to the mass of the neutron star, unrealistically large neutron star radii (R >  20 km) and an overall poor fit quality () are obtained when the distance to the source is fixed at around the estimated value (see discussion in Sect. 4).

We explored other neutron star models where the atmosphere can be partially ionized and the size of the emission radius Rem can be parametrised with respect to the neutron star physical radius (nsmaxg in XSPEC; Ho et al. 2008). For each epoch and instrument, and as before for the multi-epoch fits, we tested 19 absorbed nsmaxg models with B = (0.01−30) × 1012 G and M = 1.4 M. We first set the size of the emitting region to be the same as the neutron star radius; then we allowed this parameter to vary to smaller values to check for improved fits. As for the nsa models, we show in Table 7 the results with B = 1013 G, which generally have the highest NHP and the most consistent parameters between epochs and instruments.

Table 7.

Results of the best-fit partially ionized neutron star atmosphere model with a Gaussian absorption line (per observation and camera).

The best-fit models (with somewhat comparable NHP to the 2bb and nsa models, i.e., between 2% and 89%) are for a neutron star atmosphere composed of hydrogen at a distance of less than 110 pc (3σ). All models consisting of mid-Z element plasma (C, O, Ne) provided poor fit results. While the properties of the absorption line were found to be nearly identical to those of the fully ionized case, the temperature of the atmosphere is higher and the radiation is absorbed by roughly twice as much material, which is inconsistent with the Galactic value in the line of sight. The size of the emission region was found to be slightly smaller than the canonical 10 km of the nsa models, with Rem ∼ 8−9 km. The best-fit results of the EPIC fits are within the following values: NH = (3−4) × 1020 cm−2, Teff = (6.0−6.3) × 105 K, B = 1013 G, d <  110 pc, ϵ = (340−400) eV, σ = (85−100) eV, EW = (50−65) eV, and fX = (6.4−6.8) × 10−12 erg s−1 cm−2. The unabsorbed flux FX = (1.3−1.5) × 10−11 erg s−1 cm−2 is consistent with those of the two previously discussed models.

To break some of the degeneracy between the parameters and look for more physical results, we restricted the distance to the source within d = 100−600 pc, capped the column density at the Galactic value, and let the neutron star mass and radius vary within M = 0.5−2.5 M and R = 5−15 km. However, the exercise led to generally worse fits and significant discrepancies between the best-fit parameters of pn and MOS. No other neutron star atmosphere model in XSPEC provided acceptable fits, nor did the inclusion of a second (colder) component.

3.2.2. RGS data

For the spectral analysis of RGS data we included four XMM-Newton observations of the source performed in 2002/2003 in addition to the five 2012/AO14 observations, thus considerably extending the time span of the analysis in relation to that covered by the EPIC data (Tables 1 and 2 in Sect. 2). The total analysed RGS dataset, of which the AO14 observations account for nearly 70% in net exposure, amount to 18 RGS1/2 spectra and GTI-filtered exposures of 465 ks and 457 ks per detector.

We used the EPIC source coordinates in each observation to generate the instrument spatial masks and energy filters with rgsproc. The GTI-filtered event lists were used to extract the source and background spectra in wavelength space using the tasks rgsregions and rgsspectrum, while response matrix files were produced with the SAS task rmfgen. Only the first-order spectra were analysed. To increase the S/N each spectrum was rebinned into 0.165 Å wavelength channels. The defective channels of the RGS cameras7, which cover in first order the wavelength ranges of 11 Å–14 Å in RGS1 and 20 Å–24 Å in RGS2, were excluded from the spectral fitting. The total dataset respectively amounts to 3.992(20) × 104 and 3.724(19) × 104 counts (15–30 Å) in each RGS1/2 camera, of which around 40% can be ascribed to the background.

We fitted each observation in XSPEC assuming a model (hereafter the bbgauss model) consisting of an absorbed blackbody modified by a Gaussian absorption feature with σ = 100 eV, as found from the analysis of EPIC data; the column density was fixed to the Galactic value to better constrain the other model parameters. The RGS1/2 spectra were fitted simultaneously adopting a constant factor between the instruments. We note that the best-fit parameters from the fit of individual RGS1/2 spectrum agree well with each other in a given epoch.

The best-fit results of the bbgauss model are in Table 8. For each observation, we list the blackbody temperature kT in eV, the radiation radius R in km (assuming a source at d300), the central energy ϵ and equivalent width EW of the Gaussian absorption feature, and the unabsorbed source flux of the model FX, in the 0.2–12 keV energy band. The results suggest a possible trend of the parameters of the source; in comparison with the first four observations obtained in 2002 and 2003 and labelled (A–D) in Table 8, the more recent observations show a slight increase in temperature and a more pronounced decrease in the model normalisation and flux (formally inconsistent with a constant term) at constant properties of the Gaussian absorption feature (Fig. 4). By contrast, the parameters of the source within these two subgroups are constant at the 2%–3% level in kT, 7%–10% in R, and 5%–10% in FX. The trends are seen in the spectral parameters of both RGS1/2 instruments.

thumbnail Fig. 4.

Spectral parameters of J1605 as a function of time as measured in the RGS observations (data points; see text and Table 8 for details). The model fit to the observations in XSPEC is tbabs(bbody-gauss). Circles (black) and squares (red) show the subgroups of old and new observations of the source, obtained respectively in 2002/2003 and after 2012. The analysis of EPIC data (see Fig. 2) are for the observations in the new subgroup (MJD >  56 000). Shaded areas are the results of the fits of RGS1/2 stacked spectra in the two subgroups, with 1σ standard deviations.

Table 8.

Results of the RGS spectral analysis.

In Pires et al. (2014), we investigated the constancy of the INS emission on the EPIC data performed between 2002 and 2012. Unfortunately, the analysis does not allow us to draw definite conclusions; while the MOS instruments are unsuited for long-term studies8, the pn camera provides only one data point prior to 2012 for comparison, due to the heterogeneous dataset and background flares (Sects. 1 and 2). Nonetheless, an increase in blackbody temperature, consistent with what is observed in the RGS data, was then reported.

To investigate the possibility of a long-term evolution on the parameters of J1605, we co-added the spectra of the two subgroups in each RGS camera, using the SAS task rgscombine, and binned the results to 0.165 Å as before. The resulting grouped spectra (labelled “old” and “new” in Table 8) were then fitted simultaneously in XSPEC with the same bbgauss model, taking into account the co-added background and response files in each detector as usual. In Fig. 5 we plot the old and new spectra (grey and magenta data points) of J1605, with the folded bbgauss model and fit residuals. The best-fit results as a function of time are plotted as shaded grey and pink areas in Fig. 4.

thumbnail Fig. 5.

Results of RGS spectral fitting (see text and Table 8 for details). We show the stacked old and new RGS spectra (grey and magenta data points, respectively; the RGS1 and RGS2 spectra in each subgroup were co-added for plotting purposes). The model (solid black and red lines) fitted to the data in XSPEC is tbabs(bbody-gauss). The best-fit parameters differ in the two subgroups (see text). The two absorption lines at around λ = 20−23 Å are instrumental.

The results of this approach confirm the observed trend. With respect to the early pointings, we measure a 7% higher temperature and a 25% lower flux and smaller radiation radius in the observations performed after 2012, formally significant beyond the spectral errors. Nonetheless, the significance of the variations is low: 1σ in kT and FX, while the other parameters are consistent within the rather large spectral errors (e.g., the emission radius in Fig. 4).

The RGS instruments suffer from a decline in sensitivity at long wavelengths, likely due to a build-up of hydrocarbon contamination on the detector (de Vries et al. 2015; see also footnote 7). Empirical corrections were first introduced in 2006 to take this and other effects into account in the calibrated model of the RGS1/2 effective areas, which are estimated to have an absolute accuracy of 10%. We observe an 18% decrease in sensitivity in the effective area at long wavelengths between the old and new datasets (see Fig. 5). Altogether, uncertainties from both the spectral model and other calibration issues, possibly not accounted for in the modelling of the effective area with time, may be responsible for the discrepancies in the parameters of J1605 reported here.

The inclusion of a cold blackbody component in the bbgauss model, unlike for the EPIC data, is not satisfactory due to the large normalisation required to fit the RGS spectra. If this is kept within reasonable limits (i.e., corresponding to a < 1033 erg s−1 blackbody at d = 0.1−1 kpc), the quality of the fit is worsened; moreover, while there are no significant changes in the temperature of the hot component and in the parameters of the absorption line, the best-fit temperature of the cold component is very soft, < 30 eV, and the column density is 2 to 3 times higher than the Galactic value.

The evidence for a narrow absorption feature at energy ϵ ∼ 0.57 keV (λ = 21.5 Å) in the RGS spectra of J1605 was first reported by van Kerkwijk et al. (2004). Similarly narrow features around this wavelength have been identified in the RGS spectrum of the M7 INS RX J0720.4-3125 (Hambaryan et al. 2009) and other thermally emitting INSs (Hohle et al. 2012). To investigate the presence of the narrow feature in our dataset, we fitted each of the nine RGS19 spectra individually, using an absorbed blackbody model and two Gaussian absorption lines. As the individual datasets do not have very high S/N, we fixed the column density to the Galactic value and set the energy and FWHM of the broad absorption line to the best parameters found consistently in the analysis of EPIC and RGS data (ϵ = 410 eV and σ = 100 eV). Absorption features were then searched between 21 Å and 22.5 Å (550–590 eV).

The results are summarised in Table 9. For each fit we show the S/N of each spectrum, the blackbody temperature kT of the source, and the equivalent width of the broad absorption feature EW1; the energy ϵ2, FWHM σ2, and equivalent width EW2 of the narrow feature (when constrained) or their corresponding 1σ upper limits are also listed.

Table 9.

Investigation of a narrow absorption feature in RGS1 data.

The best-fit spectral parameters are consistent with those of Table 8, showing that the inclusion of the narrow feature is not statistically required in most cases. Considering the number of trials (40) in the searched wavelength range, the evidence for the narrow feature is only significant in observation (A). Remarkably, there is no evidence for a narrow feature within 21–22.5 Å in the longest observation of the source (labelled “201”), which has a much higher S/N than the others. The same analysis carried out in the co-added old and new spectra confirm that additional features are absent in the most recent pointings.

4. Discussion

The M7 are considered a rather homogeneous group of cooling neutron stars, displaying similar ages, temperatures, and timing properties. The source RX J1605.3+3249 stands out in that it could be slowing down at a fast rate, indicating a high dipolar field – the highest in the group – and a possible evolution from a magnetar (Pires et al. 2014). The analysis of our dedicated XMM-Newton large programme does not confirm the previous results (Sect. 3.1). Due to the energy-dependent nature of the previously detected modulation, we performed extensive high-resolution periodicity searches allowing for moderate changes of pulsed fraction and the optimal energy range and S/N for detection for a reasonably wide range of spin-down values. No significant signal resulted from the analysis; unless considerable changes of pulsed fraction have taken place since 2012, the deepest upper limit of 1.33(6)% (4σ), in the relevant frequency range, conservatively rules out the 3.39 s modulation. Moreover, in the full frequency range allowed by the timing resolution of the EPIC cameras, blind searches revealed no other periodic signals with pf ≳ 1.5% (3σ; 0.3–1.35 keV), thus considerably improving previous estimates for pulsations with P >  0.15 s. Similarly low 3σ upper limits, within 1.8% and 4%, are obtained in the same period range in narrow energy intervals (100 eV–600 eV), defined according to the source’s S/N.

With over 106 EPIC counts, the unprecedented photon statistics of the new dataset allowed the deepest investigation of the X-ray emission of the source to date (Sect. 3.2.1). Altogether, we found that no theoretical model available in XSPEC can provide a fully satisfactory physical description of the spectrum of the neutron star. While statistically acceptable fits are obtained for the individual epochs, multi-epoch fits including the spectra of all EPIC cameras have null-hypothesis probabilities of less than 1%, which may at least partially be ascribed to cross-calibration uncertainties. The best results were obtained when fitting the data with either a double-blackbody (2bb) or a magnetised neutron star atmosphere model consisting of hydrogen (nsa and nsmaxg, with B = 1013 G and M = 1.4 M), in either case modified by a broad Gaussian absorption feature as previously reported in the literature.

No significant evidence of spectral variability is measured in the 2012–2016 time frame covered by the analysis of EPIC data. The overall consistency of the parameters (between epochs and EPIC instruments) was optimal for the atmosphere models, with systematic errors of 13% in column density, 2% in temperature, and 3% in flux. In relation to the 2bb model, a canonical nsa model with B = 1013 G constrains the column density toward the source and the properties of the absorption line much better. Best-fit distances around d ∼ 130 pc are, on the other hand, inconsistently smaller than that estimated for the source (see below). The model consisting of a partially ionized hydrogen neutron star atmosphere (nsmaxg) provides even smaller distances (3σ upper limits below 110 pc) for a Rem ∼ (8−9) km emitting region on the neutron star, while the derived column density is nearly two times the Galactic value in the direction of the source.

Considering the nsa model of the multi-epoch pn fits (which is the result with the highest NHP among the multi-epoch fits), the spectral parameters of the source are constrained as NH = 2.45(14) × 1020 cm−2, Teff = 5.81(3) × 105 K, ϵ = 385 ± 10 eV, σ = 86.8 ± 0.3 eV, EW = 55 ± 3 eV, and fX = 6.842(9) × 10−12 erg s−1 cm−2 (0.2−12 keV). In contrast to previous analyses, we did not find that the inclusion of other model components (in particular, additional lines in absorption) were statistically justified or could significantly improve the results of the multi-epoch fits. Other up-to-date, fully and partially ionized neutron star atmosphere models, consisting of different elemental compositions and with noncanonical neutron star mass and radius, did not provide better fits than the models described above.

The typical distance derived from the best-fit atmosphere models, around 130 pc, is smaller than the range expected for the source, d = 350 ± 50 pc (Posselt et al. 2007, Sect. 3.2.1). This range is derived from the fitted hydrogen column density, assuming a blackbody model with three Gaussian lines in absorption, and a three-dimensional description of the distribution of the interstellar medium in the direction of the source (which, according to the authors, should be reliable up to ∼270 pc). Based on kinematic arguments, Tetzlaff et al. (2012) applied this expected range and the observed proper motion of the source (Motch et al. 2005; Zane et al. 2006) to trace back the neutron star trajectory and determine its likely birthplace, using possible associations with runaway massive stars and the observed abundance of heavy elements as further evidence to narrow down the most likely solutions. These predict a current distance of d = 300−370 pc if the neutron star was born less than 0.5 Myr ago in a nearby supernova explosion.

4.1. Viewing geometry and presence of hotspots

Magnetic fields in the range observed in the M7 are expected to produce large temperature variations on the neutron star surface, due to the anisotropic electron conductivity and heat transport in the stellar envelope and crust (e.g., Geppert et al. 2004; Pérez-Azorín et al. 2006; Perna et al. 2013). In this case, strong pulsed flux variations are expected at the neutron star spin period unless the source is observed from a particularly unfavourable geometry: either if the angle between the neutron star spin axis and the line of sight i is sufficiently small or if the regions of higher temperature (hotspots) are located at a very small angle θB in relation to the neutron star rotation axis (i.e., the magnetic and spin axes of the star are nearly co-aligned). Therefore, the stringent pulsed fraction limits from the timing analysis can be used to verify the viability of the 2bb model and probe the presence of hotspots on the surface of J1605 (see, e.g., the case of the thermally emitting central neutron star in the supernova remnant HESS J1731-347; Suleimanov et al. 2017).

With this goal we considered the model of a slowly rotating neutron star, observed at an inclination angle i, with two identical polar hotspots located at an angle θB with respect to the spin axis (see, e.g., Page 1995; Schwope et al. 2005; Suleimanov et al. 2010, for a full description of the model). Light bending in the vicinity of the neutron star follows the relation between the local angle of the emitted photon and its escape direction and depends on the compactness of the neutron star (given by the ratio between the neutron star and the Schwarzschild radius rg = Rnsc2(2GMns)−1; Beloborodov 2002). The temperature of the neutron star surface, eV, and of the hotspots, eV, are assumed from the best-fit 2bb pn model, which is the one with the highest NHP (Table 5). For a neutron star distance of pc (Tetzlaff et al. 2012), the model normalisations set the corresponding sizes of the emission regions as km and km. With these values the angular size of the spots and the compactness of the star are fixed as θsp = 4.7° and rg = 2.8, respectively, assuming a 1.5 M neutron star.

For a particular viewing geometry (i, θB), the photon flux at a given rotation phase results from the sum of the visible individual area elements of the neutron star surface, assuming blackbody emission and taking into account the light bending. The flux is corrected for the interstellar absorption of a equivalent column density of cm−2 and then folded with the EPIC pn response to give the source count rate at the 0.5–1.35 keV energy band. The energy band is chosen to minimise the effects of the broad absorption feature in the emitted spectrum.

With this method we computed an extensive grid of light curves for (θB, i) within (0 ° ,0 ° ) and (90 ° ,90 ° ) and computed the pulsed fraction for each orientation as

In Fig. 6 we plot the resulting pf map in the (θB, i) plane. The maximum pulsed fraction obtained for the 2bb model is about 20%. The region allowed by the limits of the timing analysis lies below the thick pink line corresponding to pf = 2.78(16)% (Table 4; 0.5–1.35 keV). Integrating over all possible random orientations of line-of-sight inclination and spot angles we obtain a small likelihood (∼1.9%) that we do not see pulsations from the source due to the particularly unfavourable viewing geometry if the 2bb model is correct.

thumbnail Fig. 6.

Contours of constant pulsed fraction assuming the best-fit parameters of the pn double-blackbody model ( eV, km, kTsp = 117.0(8) eV, km). Count rates are computed in the 0.5–1.35 keV energy band to exclude the effects of the Gaussian absorption line and pulsed fraction labels are given in percent. Darker colours denote higher pulsed fractions. The allowed parameter space of the viewing geometry constrained by the timing analysis lies below the thick pink line (pf ≤ 2.78(16)%; Table 4).

4.2. Energy distribution

Thermal emission from INSs is expected to originate immediately at the surface, with the bulk of the energy flux peaking in the soft X-ray band. In principle, by confronting the observed spectra and light curves with theoretical models for neutron star thermal radiation, it should be possible to derive the surface temperature, magnetic field, gravitational acceleration, and chemical composition; if distances are known, then the stellar mass, radius, and the equation of state of neutron star interior can be constrained as well (see Potekhin et al. 2015; Özel & Freire 2016, for recent reviews on neutron star atmosphere models and up-to-date astrophysical constraints on the equation of state of nuclear matter). Since their discovery in the All-Sky Survey of the ROSAT satellite (Voges et al. 1999), the M7 have been regarded as the closest-to-perfect candidates for testing neutron star emission models, based on a combination of bright thermal emission, proximity, independent distance estimates10, and a lack of significant magnetospheric or accretion activity.

In practice, progress has been hampered by uncertainties on the chemical composition of the atmosphere and the phase state of the stellar surface, and by the lack of understanding on the magnetic field and temperature distributions (e.g., Zane & Turolla 2006; van Kerkwijk & Kaplan 2007; Suleimanov et al. 2010). The presence of lines adds to this complexity; although believed to be related with the star’s magnetic field, they have no unique physical interpretation. To explain the emitted radiation and equivalent widths of the lines in the phase-resolved spectrum of the M7 RX J1308.6+2127, Suleimanov et al. (2010) favoured a model where a partially ionized, optically thin atmosphere above the condensed surface must be present (see also Motch et al. 2003; Ho et al. 2007). Using this model, Hambaryan et al. (2011) derived the temperatures of the X-ray emitting areas and the magnetic field intensity at the poles; moreover, they were able to constrain the compactness of the neutron star and the gravitational redshift on the surface, suggesting a very stiff equation of state. Similar conclusions were reached for the M7 RX J0720.4-3125 (Hambaryan et al. 2017). These results are only marginally compatible with the most favoured range of the true radius of a 1.5 M neutron star, 10–11.5 km, from the analysis of Özel & Freire (2016). New generation X-ray missions, in particular the Neutron star Interior Composition Explore mission (NICER; Gendreau et al. 2012), together with more accurate distances from the Gaia satellite (Gaia Collaboration 2016), will certainly improve the constraints on neutron star mass and radius from astrophysical observations of millisecond pulsars in globular clusters, for example.

Recently, Viganò et al. (2014) showed that in some cases (such as the M7 RX J0806.4-4123; Haberl et al. 2004) the deviations found in the spectra of thermally emitting INSs may be induced simply by the inhomogeneous temperature distribution on the surface. While the effect is unlikely to account for all cases of sources with reported spectral features, the interesting result is that the anisotropic temperature distribution can give way to spurious spectral features. We can safely exclude this possibility for the absorption feature in J1605, which cannot be accommodated by a multi-temperature energy distribution.

All M7 INSs have detected optical, ultraviolet, or infrared counterparts (see Kaplan et al. 2011; Posselt et al. 2014, 2018, for references and limits). Interestingly, the extrapolation to longer wavelengths of the best-fit model inferred from X-rays – including for the case of J1605 both the double-temperature blackbody and the atmosphere models discussed here – falls below the actual detected fluxes. This is known as the optical excess, and it is observed in all M7 INSs. The optical excess of J1605 deviates significantly from the expected Rayleigh-Jeans slope of the spectra and can be described by a rather flat power law (Kaplan et al. 2011). The origin of the excess flux might rely on atmospheric effects, magnetospheric emission, or resonant scattering. In particular, if the X-ray and optical/UV emission came from different regions on the surface (e.g., Braje & Romani 2002), we might expect correlations between the amount of optical excess and the X-ray pulsed fraction, which have not been verified. The possible presence of fossil fallback disks, of a faint pulsar wind nebula, “bare” neutron star surfaces, and other alternative scenarios remain open (see Turolla et al. 2004; Ertan et al. 2017; Wang et al. 2017; Posselt et al. 2018, for a discussion).

4.3. Narrow absorption feature

In combination with the existing archival RGS data, the AO14 campaign accumulated over 450 ks of exposure time on J1605, increasing the available data by 70%. The good statistics allowed a detailed analysis of the narrow feature at ϵ ∼ 570 eV, previously reported in the literature. The investigation shows that the feature is only significantly detected in one early epoch (the 2002 observation first analysed by van Kerkwijk et al. 2004) or when this observation is combined with the archival data obtained prior to the AO14 campaign (Hohle et al. 2012; Pires et al. 2014). The feature is definitely not present in the AO14 observations, while evidence of a less significant narrow feature is present at energy ϵ = 550 eV in the 2012 observation. Consistently with these results, the analysis of the two grouped spectra of early and recent observations constrain the presence of the narrow feature only in the first subgroup (Sect. 3.2.2).

Hambaryan et al. (2009) discusses the possible physical interpretation of a similar feature detected in the co-added RGS spectrum of RX J0720.4-3125, which was later confirmed by the analysis of Chandra LETG data (Hohle et al. 2012). Their analysis favours a blend of highly ionized oxygen originating in the ambient medium of the INS, possibly a high-density nearby cloud which could contribute to the source’s optical excess. Nonetheless, an interstellar or atmospheric origin cannot be ruled out. At least for J1605, the narrow and transient nature of the feature disfavours an atmospheric origin.

Phase-dependent narrow absorption features have been reported in XMM-Newton observations of the M7 INSs RX J0720.4-3125 and RX J1308.6+2127 (Borghese et al. 2015, 2017), a work motivated by the detection of variable cyclotron lines detected in the spectra of two “low magnetic field” magnetars (Tiengo et al. 2013; Rodríguez Castillo et al. 2016). These results give support for the presence of strong, confined magnetic field components close to the stellar surface and a complex field topology. In contrast, the features in the spectra of the two M7 INSs are intrinsically different in that they do not vary in energy, are detected at much lower energy, and also seem to be stable and lasting on long timescales.

5. Summary and conclusions

We report here the results of a XMM-Newton large programme on the thermally emitting isolated neutron star RX J1605.3+3249. The goal of the project was to gain a deeper understanding of the timing and spectral properties of the source through a detailed analysis of its X-ray emission. The neutron star is of particular scientific interest as a source that could potentially bridge the evolutionary gap between the groups of nearby thermally emitting sources dubbed the Magnificent Seven and the young and energetic magnetars. Due to the lack of detected pulsations, our science goals could only be partially completed. Nonetheless, the deep upper limits derived from our analysis were used to put stringent constraints on the viewing geometry and the presence of hot spots on the surface. Detailed phase-averaged medium- and high-resolution spectroscopy constrains atmosphere neutron star models and the properties of the cyclotron line in the spectrum of the neutron star with unprecedented statistics. The nondetection of the narrow absorption feature at ϵ = 570 eV reported in previous epochs also disfavours an atmospheric origin.


2

See, e.g., Groth (1975), for the method to extract upper limits on the pulsed fraction.

3

The maximum frequency is determined by the Nyquist limit, given the time resolution of the detector.

6

For comparison, in Table 5 for the multi-epoch fits we also list the results of a single-temperature blackbody model with a Gaussian absorption line. For this model, the fit quality is generally poor () and the column density is unconstrained.

7

XMM-Newton Calibration Technical Note 0030, issue 7.7

8

XMM-Newton Calibration Technical Note 0018.

9

The RGS2 data cannot be used due to the defective channels of the camera around the wavelength range of interest.

10

From HST parallaxes (in two cases) and kinematic studies (e.g., Walter et al. 2010; Tetzlaff et al. 2010, and references therein).

Acknowledgments

We thank the anonymous referee for useful comments and suggestions which helped to improve the paper. We thank Norbert Schartel and the staff of the XMM-Newton Science Operation Center, in particular Jan-Uwe Ness, for their great support in the scheduling of these time-constrained observations. We acknowledge the use of the ATNF Pulsar Catalogue (http://www.atnf.csiro.au/people/pulsar/psrcat). This work was supported by the Deutsches Zentrum für Luft- und Raumfahrt (DLR) under grant 50 OR 1511.

References

  1. Arnaud, K. A. 1996, in Astronomical Data Analysis Software and Systems V, eds. G. H. Jacoby, & J. Barnes, ASP Conf. Ser., 101, 17 [NASA ADS] [Google Scholar]
  2. Beloborodov, A. M. 2002, ApJ, 566, L85 [NASA ADS] [CrossRef] [Google Scholar]
  3. Beloborodov, A. M., & Li, X. 2016, ApJ, 833, 261 [NASA ADS] [CrossRef] [Google Scholar]
  4. Borghese, A., Rea, N., Coti Zelati, F., Tiengo, A., & Turolla, R. 2015, ApJ, 807, L20 [NASA ADS] [CrossRef] [Google Scholar]
  5. Borghese, A., Rea, N., Coti Zelati, F., et al. 2017, MNRAS, 468, 2975 [NASA ADS] [CrossRef] [Google Scholar]
  6. Braje, T. M., & Romani, R. W. 2002, ApJ, 580, 1043 [NASA ADS] [CrossRef] [Google Scholar]
  7. Buccheri, R., Bennett, K., Bignami, G. F., et al. 1983, A&A, 128, 245 [NASA ADS] [Google Scholar]
  8. Coti Zelati, F., Rea, N., Pons, J. A., Campana, S., & Esposito, P. 2018, MNRAS, 474, 961 [NASA ADS] [CrossRef] [Google Scholar]
  9. de Vries, C. P., den Herder, J. W., Gabriel, C., et al. 2015, A&A, 573, A128 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. den Herder, J. W., Brinkman, A. C., Kahn, S. M., et al. 2001, A&A, 365, L7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Ertan, Ü., Çalışkan, Ş., & Alpar, M. A. 2017, MNRAS, 470, 1253 [NASA ADS] [CrossRef] [Google Scholar]
  12. Esposito, P., Rea, N., & Israel, G. L. 2018, ArXiv e-prints [arXiv:1803.05716] [Google Scholar]
  13. Evans, I. N., Primini, F. A., Glotfelty, K. J., et al. 2010, ApJS, 189, 37 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  14. Gaia Collaboration (Prusti, T., et al.) 2016, A&A, 595, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Gendreau, K. C., Arzoumanian, Z., & Okajima, T. 2012, Proc. SPIE, 8443, 844313 [CrossRef] [Google Scholar]
  16. Geppert, U., Küker, M., & Page, D. 2004, A&A, 426, 267 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Goldreich, P., & Reisenegger, A. 1992, ApJ, 395, 250 [NASA ADS] [CrossRef] [Google Scholar]
  18. Groth, E. J. 1975, ApJS, 29, 285 [NASA ADS] [CrossRef] [Google Scholar]
  19. Haberl, F. 2007, Ap&SS, 308, 181 [NASA ADS] [CrossRef] [Google Scholar]
  20. Haberl, F., Motch, C., Zavlin, V. E., et al. 2004, A&A, 424, 635 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Hambaryan, V., Neuhäuser, R., Haberl, F., Hohle, M. M., & Schwope, A. D. 2009, A&A, 497, L9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Hambaryan, V., Suleimanov, V., Schwope, A. D., et al. 2011, A&A, 534, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Hambaryan, V., Suleimanov, V., Haberl, F., et al. 2017, A&A, 601, A108 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Heyl, J. S., & Kulkarni, S. R. 1998, ApJ, 506, L61 [NASA ADS] [CrossRef] [Google Scholar]
  25. Ho, W. C. G., Kaplan, D. L., Chang, P., van Adelsberg, M., & Potekhin, A. Y. 2007, MNRAS, 375, 821 [NASA ADS] [CrossRef] [Google Scholar]
  26. Ho, W. C. G., Potekhin, A. Y., & Chabrier, G. 2008, ApJS, 178, 102 [NASA ADS] [CrossRef] [Google Scholar]
  27. Hohle, M. M., Haberl, F., Vink, J., de Vries, C. P., & Neuhäuser, R. 2012, MNRAS, 419, 1525 [NASA ADS] [CrossRef] [Google Scholar]
  28. Jansen, F., Lumb, D., Altieri, B., et al. 2001, A&A, 365, L1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Jethwa, P., Saxton, R., Guainazzi, M., Rodriguez-Pascual, P., & Stuhlinger, M. 2015, A&A, 581, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Kalberla, P. M. W., Burton, W. B., Hartmann, D., et al. 2005, A&A, 440, 775 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Kaplan, D. L. 2008, in 40 Years of Pulsars: Millisecond Pulsars, Magnetars and More, eds. C. Bassa, Z. Wang, A. Cumming, & V. M. Kaspi, AIP Conf. Ser., 983, 331 [NASA ADS] [CrossRef] [Google Scholar]
  32. Kaplan, D. L., & van Kerkwijk, M. H. 2009, ApJ, 705, 798 [NASA ADS] [CrossRef] [Google Scholar]
  33. Kaplan, D. L., Kamble, A., van Kerkwijk, M. H., & Ho, W. C. G. 2011, ApJ, 736, 117 [NASA ADS] [CrossRef] [Google Scholar]
  34. Kaspi, V. M., & Beloborodov, A. M. 2017, ARA&A, 55, 261 [NASA ADS] [CrossRef] [Google Scholar]
  35. Lasker, B. M., Lattanzi, M. G., McLean, B. J., et al. 2008, AJ, 136, 735 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  36. Manchester, R. N., Hobbs, G. B., Teoh, A., & Hobbs, M. 2005, AJ, 129, 1993 [NASA ADS] [CrossRef] [Google Scholar]
  37. Motch, C., Haberl, F., Zickgraf, F.-J., Hasinger, G., & Schwope, A. D. 1999, A&A, 351, 177 [NASA ADS] [Google Scholar]
  38. Motch, C., Zavlin, V. E., & Haberl, F. 2003, A&A, 408, 323 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. Motch, C., Sekiguchi, K., Haberl, F., et al. 2005, A&A, 429, 257 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Ostriker, J. P., & Gunn, J. E. 1969, ApJ, 157, 1395 [NASA ADS] [CrossRef] [Google Scholar]
  41. Özel, F., & Freire, P. 2016, ARA&A, 54, 401 [Google Scholar]
  42. Page, D. 1995, ApJ, 442, 273 [NASA ADS] [CrossRef] [Google Scholar]
  43. Pavlov, G. G., Shibanov, Y. A., Zavlin, V. E., & Meyer, R. D. 1995, in NATO Advanced Science Institutes (ASI) Series C, eds. M. A. Alpar, U. Kiziloglu, & J. van Paradijs, 450, 71 [Google Scholar]
  44. Pavlov, G. G., Zavlin, V. E., & Trümper, J. 1999, ApJ, 511, L45 [NASA ADS] [CrossRef] [Google Scholar]
  45. Pérez-Azorín, J. F., Miralles, J. A., & Pons, J. A. 2006, A&A, 451, 1009 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. Perna, R., & Pons, J. A. 2011, ApJ, 727, L51 [NASA ADS] [CrossRef] [Google Scholar]
  47. Perna, R., Viganò, D., Pons, J. A., & Rea, N. 2013, MNRAS, 434, 2362 [NASA ADS] [CrossRef] [Google Scholar]
  48. Pires, A. M., Haberl, F., Zavlin, V. E., et al. 2014, A&A, 563, A50 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Pons, J. A., Miralles, J. A., & Geppert, U. 2009, A&A, 496, 207 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Popov, S. B., Pons, J. A., Miralles, J. A., Boldin, P. A., & Posselt, B. 2010, MNRAS, 401, 2675 [NASA ADS] [CrossRef] [Google Scholar]
  51. Posselt, B., Popov, S. B., Haberl, F., et al. 2007, Ap&SS, 308, 171 [NASA ADS] [CrossRef] [Google Scholar]
  52. Posselt, B., Pavlov, G. G., Popov, S., & Wachter, S. 2014, ApJS, 215, 3 [NASA ADS] [CrossRef] [Google Scholar]
  53. Posselt, B., Pavlov, G. G., Ertan, Ü., et al. 2018, ApJ, 865, 1 [NASA ADS] [CrossRef] [Google Scholar]
  54. Potekhin, A. Y., De Luca, A., & Pons, J. A. 2015, Space Sci. Rev., 191, 171 [NASA ADS] [CrossRef] [Google Scholar]
  55. Read, A. M., Guainazzi, M., & Sembay, S. 2014, A&A, 564, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Rodríguez Castillo, G. A., Israel, G. L., Tiengo, A., et al. 2016, MNRAS, 456, 4145 [NASA ADS] [CrossRef] [Google Scholar]
  57. Schwope, A. D., Hambaryan, V., Haberl, F., & Motch, C. 2005, A&A, 441, 597 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Skrutskie, M. F., Cutri, R. M., Stiening, R., et al. 2006, AJ, 131, 1163 [NASA ADS] [CrossRef] [Google Scholar]
  59. Strüder, L., Briel, U., Dennerl, K., et al. 2001, A&A, 365, L18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Suleimanov, V., Hambaryan, V., Potekhin, A. Y., et al. 2010, A&A, 522, A111 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  61. Suleimanov, V. F., Klochkov, D., Poutanen, J., & Werner, K. 2017, A&A, 600, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Tetzlaff, N., Neuhäuser, R., Hohle, M. M., & Maciejewski, G. 2010, MNRAS, 402, 2369 [NASA ADS] [CrossRef] [Google Scholar]
  63. Tetzlaff, N., Schmidt, J. G., Hohle, M. M., & Neuhäuser, R. 2012, PASA, 29, 98 [NASA ADS] [CrossRef] [Google Scholar]
  64. Thompson, C., & Duncan, R. C. 1995, MNRAS, 275, 255 [NASA ADS] [CrossRef] [Google Scholar]
  65. Thompson, C., & Duncan, R. C. 1996, ApJ, 473, 322 [NASA ADS] [CrossRef] [Google Scholar]
  66. Tiengo, A., Esposito, P., Mereghetti, S., et al. 2013, Nature, 500, 312 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  67. Turner, M. J. L., Abbey, A., Arnaud, M., et al. 2001, A&A, 365, L27 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Turolla, R. 2009, in Astrophys. Space Sci. Lib., ed. W. Becker, 357, 141 [NASA ADS] [CrossRef] [Google Scholar]
  69. Turolla, R., Zane, S., & Drake, J. J. 2004, ApJ, 603, 265 [NASA ADS] [CrossRef] [Google Scholar]
  70. Turolla, R., Zane, S., Pons, J. A., Esposito, P., & Rea, N. 2011, ApJ, 740, 105 [NASA ADS] [CrossRef] [Google Scholar]
  71. Turolla, R., Zane, S., & Watts, A. L. 2015, Rep. Prog. Phys., 78, 116901 [Google Scholar]
  72. van Kerkwijk, M. H., & Kaplan, D. L. 2007, Ap&SS, 308, 191 [NASA ADS] [CrossRef] [Google Scholar]
  73. van Kerkwijk, M. H., Kaplan, D. L., Durant, M., Kulkarni, S. R., & Paerels, F. 2004, ApJ, 608, 432 [NASA ADS] [CrossRef] [Google Scholar]
  74. Viganò, D., Rea, N., Pons, J. A., et al. 2013, MNRAS, 434, 123 [NASA ADS] [CrossRef] [Google Scholar]
  75. Viganò, D., Perna, R., Rea, N., & Pons, J. A. 2014, MNRAS, 443, 31 [NASA ADS] [CrossRef] [Google Scholar]
  76. Voges, W., Aschenbach, B., Boller, T., et al. 1999, A&A, 349, 389 [NASA ADS] [Google Scholar]
  77. Walter, F. M., Eisenbeiß, T., Lattimer, J. M., et al. 2010, ApJ, 724, 669 [NASA ADS] [CrossRef] [Google Scholar]
  78. Wang, W., Lu, J., Tong, H., et al. 2017, ApJ, 837, 81 [NASA ADS] [CrossRef] [Google Scholar]
  79. Willingale, R., Starling, R. L. C., Beardmore, A. P., Tanvir, N. R., & O’Brien, P. T. 2013, MNRAS, 431, 394 [NASA ADS] [CrossRef] [Google Scholar]
  80. Wilms, J., Allen, A., & McCray, R. 2000, ApJ, 542, 914 [NASA ADS] [CrossRef] [Google Scholar]
  81. Zane, S., & Turolla, R. 2006, MNRAS, 366, 727 [NASA ADS] [CrossRef] [Google Scholar]
  82. Zane, S., de Luca, A., Mignani, R. P., & Turolla, R. 2006, A&A, 457, 619 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  83. Zavlin, V. E., Pavlov, G. G., & Shibanov, Y. A. 1996, A&A, 315, 141 [NASA ADS] [Google Scholar]

All Tables

Table 1.

Log of the XMM-Newton AO14 and 2012 observations of RX J1605.3+3249.

Table 2.

Summary of archival XMM-Newton observations of RX J1605.3+3249 used in the RGS spectral analysis.

Table 3.

Parameters of RX J1605.3+3249, as extracted from the EPIC images in the AO14 and 2012 observations.

Table 4.

Upper limits on pulsations from the timing analysis.

Table 5.

Results of the best-fit double-blackbody model with a Gaussian absorption line (per observation and camera).

Table 6.

Results of the best-fit fully ionized neutron star atmosphere model with a Gaussian absorption line (per observation and camera).

Table 7.

Results of the best-fit partially ionized neutron star atmosphere model with a Gaussian absorption line (per observation and camera).

Table 8.

Results of the RGS spectral analysis.

Table 9.

Investigation of a narrow absorption feature in RGS1 data.

All Figures

thumbnail Fig. 1.

EPIC and pn searches around the 2012 periodicity. The frequency range is ν = 0.2948−0.2953 Hz. The periodogram in the background (dashed outline) shows the 2012 result: ν2012 = 0.2951709(14) Hz and (ν2012) ∼ 50 (obsid 0671620101, pn, 0.5–1.35 keV). The four plots show for each AO14 observation the results of tests conducted in seven different energy bands for an extraction region of 100″ (see text). The dotted and dashed horizontal lines show the 2σ and 3σ confidence levels for the detection of modulations in each observation.

In the text
thumbnail Fig. 2.

Source best-fit parameters as a function of MJD (EPIC analysis; see text and Tables 5 and 6 for details). The model fit to the data is that of absorbed double-blackbody (left) and fully ionized hydrogen neutron star atmosphere (B = 1013 G, M = 1.4 M, and R = 10 km; right), modified by a broad Gaussian absorption line. The plots show in each instrument and observation (see legend): panel a: column density, temperature of the (panel b) cold and panel c: hot blackbody components, panel d: effective temperature of the neutron star (unredshifted), panel e: central energy of the absorption line, panel f: its Gaussian sigma, and panelg: observed model flux in the 0.2–12 keV energy band. The total Galactic NH value is shown by the solid black line in plots (panel a). Horizontal lines show the results of the simultaneous fit of spectra over the five pointings, with 1σ confidence levels comprised by the shaded areas.

In the text
thumbnail Fig. 3.

Results of EPIC spectral fitting (see text and Table 6 for details). We show the 5 pn and stacked MOS spectra (grey and magenta data points, respectively), fitted simultaneously by a fully ionized hydrogen neutron star atmosphere model with B = 1013 G, Teff = (5.16 ± 0.17) × 105 K, and a broad Gaussian absorption line of σ = 91.4 ± 0.4 eV at eV (black and red solid lines).

In the text
thumbnail Fig. 4.

Spectral parameters of J1605 as a function of time as measured in the RGS observations (data points; see text and Table 8 for details). The model fit to the observations in XSPEC is tbabs(bbody-gauss). Circles (black) and squares (red) show the subgroups of old and new observations of the source, obtained respectively in 2002/2003 and after 2012. The analysis of EPIC data (see Fig. 2) are for the observations in the new subgroup (MJD >  56 000). Shaded areas are the results of the fits of RGS1/2 stacked spectra in the two subgroups, with 1σ standard deviations.

In the text
thumbnail Fig. 5.

Results of RGS spectral fitting (see text and Table 8 for details). We show the stacked old and new RGS spectra (grey and magenta data points, respectively; the RGS1 and RGS2 spectra in each subgroup were co-added for plotting purposes). The model (solid black and red lines) fitted to the data in XSPEC is tbabs(bbody-gauss). The best-fit parameters differ in the two subgroups (see text). The two absorption lines at around λ = 20−23 Å are instrumental.

In the text
thumbnail Fig. 6.

Contours of constant pulsed fraction assuming the best-fit parameters of the pn double-blackbody model ( eV, km, kTsp = 117.0(8) eV, km). Count rates are computed in the 0.5–1.35 keV energy band to exclude the effects of the Gaussian absorption line and pulsed fraction labels are given in percent. Darker colours denote higher pulsed fractions. The allowed parameter space of the viewing geometry constrained by the timing analysis lies below the thick pink line (pf ≤ 2.78(16)%; Table 4).

In the text

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

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

Initial download of the metrics may take a while.