CORRELATED X-RAY AND VERY HIGH ENERGY EMISSION IN THE GAMMA-RAY BINARY LS I +61 303

The VHE gamma-ray observations were performed from 2007 September 4 to 21 using the MAGIC telescope on the Canary Island of La Palma (2875N, 1786W, 2225 m a.s.l.), from where LS I +61 303 is observable at zenith angles above 32°. The essential parameters of MAGIC are a 17 m diameter segmented mirror of parabolic shape (currently the largest dish for an air Cherenkov telescope), an f/D of 1.05, and an hexagonally shaped camera of 576 hemispherical photo multiplier tubes with a field of view of 35 diameter. MAGIC can detect gamma rays from 60 GeV to several TeV, and with a special setup the trigger threshold can be reduced to 25 GeV (Aliu et al. 2008). Its energy resolution is ΔE/E = 20% above energies of 200 GeV. The current sensitivity is 1.6% of the Crab Nebula flux for a 5σ detection in 50 hr of on-source time. The improvement compared to the previous sensitivity was achieved by installing new 2 GHz FADCs (Albert et al. 2008b).

The total observation time was 58.8 hr, including 39.6 hr in dark conditions and 19.2 hr under moonlight or twilight (Britzger et al. 2009). The range of zenith angles for these observations was [32°, 50°], with 96% of the data having zenith angles below 43°. The final effective total observation time is 54.2 hr. Table 1 gives the MJD, the effective observation time, and the orbital phase for each of the MAGIC observations. The observations were carried out in wobble mode (Fomin et al. 1994), i.e., by alternately tracking two positions at 04 offset from the source position. This observing mode provides a reliable background estimate for point-like sources such as LS I +61 303.

The data analysis was performed using the standard MAGIC analysis and reconstruction software (Albert et al. 2008c). The quality of the data was checked and data taken with anomalous event rates (very low atmospheric transmission, car light flashes, etc.) were rejected following the standard procedure, as previously done in Albert et al. (2009). From the remaining events, image parameters were calculated (Hillas 1985). In addition, the time parameters described in Aliu et al. (2009) were calculated as well. For the γ/hadron separation, a multidimensional classification procedure based on the image and timing parameters with the Random Forest method was used (Albert et al. 2008d). We optimized the signal selection cuts with contemporaneous Crab Nebula data taken at the same zenith angle range. The energy of the primary γ-ray was reconstructed from the image parameters using also a Random Forest method. The differential energy spectrum was unfolded taking into account the full instrumental energy resolution (Albert et al. 2007). We estimate the systematic uncertainty to be about 30% for the derived integral flux values and ±0.2 for the obtained photon index. For more details on the systematic uncertainties present in the MAGIC data, see Albert et al. (2008c).

We observed LS I +61 303 with XMM-Newton during seven runs from 2007 September 4 to 11, lasting from 12 to 18 ks (see Table 2), amounting to a total observation time of 104.3 ks. The EPIC pn detector used the Large Window Mode, while the EPIC MOS detectors used the Small Window Mode. All detectors used the medium thickness optical blocking filter. Data were processed using version 8.0.0 of the XMM-Newton Science Analysis Software (SAS). Known hot or flickering pixels were removed using the standard SAS tasks. Further cleaning was necessary to remove from the data set periods of high background corresponding to soft proton flares, reducing the net good exposure durations to 67.0 and 92.6 ks for the pn and MOS detectors, respectively.

Cleaned pn and MOS event files were extracted for spectral analysis. Source spectra were extracted from a ∼70'' radius circle centered on the source (point-spread function of 15'') while background spectra were taken from a number of source-free circles with ∼150'' radius (three for the pn detector, four for MOS1, and five for MOS2). The extracted spectra were analyzed with XSpec v12.3.1 (Arnaud 1996). The spectra were binned so that a minimum of 20 counts per bin were present and energy resolution was not oversampled by more than a factor 3. An absorbed power-law function (wabs*powerlaw) yielded satisfactory fits for all observations. De-absorbed fluxes in the 0.3–10 keV range were computed from the spectral fits.

Additional observations of 2–5 ks each, consisting of several pointings of 0.2–1.0 ks, were obtained with the Swift/X-ray Telescope (XRT) from 2007 September 11 to 22. The total observation time was 28.5 ks. The Swift data were processed using the FTOOLS task xrtpipeline (ver. 0.12.1 under HEASoft 6.6). The spectral analysis procedures were the same as those used for the XMM-Newton data. Since the hydrogen column density was poorly constrained, with typical values of (0.5 ± 0.2) × 10²² cm⁻², we fixed it to 0.5 × 10²² cm⁻², a typical value for LS I +61 303 also found in the XMM-Newton fits and close to the average value for the Swift fits.

To search for short-term X-ray variability, we also extracted 0.3–10 keV background-subtracted light curves for each observation, binned to 100 s for XMM-Newton and at half-time of the sparse 0.2–1.0 ks Swift pointings within each observation. For XMM-Newton, we considered the sum of the count rate in the three detectors. From these light curves, we computed the degree of variability as the standard deviation with respect to the mean of the count rate divided by this mean, and the significance of this variability computed from the χ² probability of the count rate being constant. We also computed hardness ratios as the fraction between the count rates above and below 2 keV.

The measured fluxes of LS I +61 303 at energies above 300 GeV are listed in Table 1, and the corresponding light curve is shown in Figure 1 (top). The periodically modulated peak emission (Albert et al. 2009) is prominently seen as the highest flux value of these observations at phase 0.62 (followed by a smaller flux at phase 0.66). In addition to this established feature we measured, as in 2006 December, significant flux on a nightly basis at phase 0.85. We note that there is significant flux in the phase range 0.8–1.0 at the level of (5.2 ± 1.0) × 10⁻¹² cm⁻² s⁻¹, which is compatible with the 2σ upper limit we obtained for the shorter 2006 observations (Albert et al. 2009).

Zoom In Zoom Out Reset image size

We also determined the differential energy spectrum from the ∼10 hr of observations conducted in the phase range 0.6–0.7. A power-law fit yields

$$\begin{eqnarray} \frac{dN}{dE} = &\frac{(2.0 \pm 0.3_{\rm stat} \pm 0.6_{\rm syst})\times 10^{-12}}{{\rm cm}^{2} \;{\rm s} \;{\rm TeV}} \nonumber \\ &\times \left(\frac{E}{1 {\rm \;TeV}}\right)^{-2.7 \pm 0.3_{\rm stat} \pm 0.2_{\rm syst}}, \nonumber \end{eqnarray}$$

compatible within errors with our previous measurements (Albert et al. 2006, 2009).

We summarize in Table 2 the parameters of the spectral fits and variability obtained for both the XMM-Newton and Swift/XRT data sets. All X-ray fluxes quoted hereafter are de-absorbed. We note that there is no significant hardness ratio change within each of the observations, in contrast to what was found in the XMM-Newton observations reported by Sidoli et al. (2006). In principle this implies that, for each observation, the flux and corresponding uncertainty obtained from the spectral fit is a good estimate of the flux during the whole observation. However, moderate count-rate variability on timescales of ks is present in most observations, ranging from 6% to 17% for the XMM-Newton data and from 9% to 25% for the Swift data (see Table 2). The rms of this variability should be considered in addition to the statistical uncertainty when providing a flux measurement spanning several ks. We converted this count-rate variability into flux variability by multiplying the degree of variability defined in Section 2.2 with the fluxes coming from the spectral fits. Since no spectral change is detected within each observation, we added this flux variability in quadrature to the spectral fits flux errors. This procedure provides the more realistic total flux uncertainties quoted in parentheses in Table 2, and used hereafter.

In Figure 1 (bottom) we show the 0.3–10 keV light curve of LS I +61 303. The source displays a steady flux during the first four observations and shows a steep increase at phase 0.62, which is followed by a slower decay up to phase 0.69 (XMM-Newton) and probably up to phase 0.72 (Swift). The behavior is very similar to the one seen at VHE gamma rays, although at X-ray energies the baseline has a significant flux. Later on, there is a significant increase of the X-ray flux up to phase 0.8. This high flux is detected with a sparse sampling up to phase 0.9, and the source goes back to its baseline flux at phase 1.0. This high X-ray flux between phases 0.8 and 1.0 occurs when the source is also detected at VHE gamma rays.

A clear correlation between the X-ray and VHE gamma-ray emissions is seen during the outburst, with a simultaneous peak at phase 0.62 (see Figure 1). To study the significance of this correlation, we selected the X-ray data sets that overlap with MAGIC observations. There are six overlapping MAGIC/XMM-Newton data sets, for which strictly simultaneous observations range from 3.3 to 3.9 hr. For the four overlapping MAGIC/Swift data sets, the strictly simultaneous observations range from 2.2 to 4.1 hr (although the Swift runs have gaps). In Figure 2 we plot the X-ray fluxes against the VHE fluxes (from Tables 1 and 2) for all 10 simultaneous pairs, which are marked with arrows in Figure 1. The linear correlation coefficient for the six simultaneous MAGIC/XMM-Newton pairs that trace the outburst is r = 0.97. For the 10 simultaneous pairs we find r = 0.81 (a |r| larger than that has a probability of about 5 × 10⁻³ to be produced from independent X-ray and VHE fluxes). Minimizing χ² we obtain χ² = 7.68 for 8 degrees of freedom and the following relationship: F(0.3–10 keV)/[10⁻¹² erg cm⁻² s⁻¹] = (12.2^+0.9_−1.0) + (0.71^+0.17_−0.14) × N(E > 300 GeV)/[10⁻¹² cm⁻² s⁻¹] (non-Gaussian uncertainties). This fit is plotted as a solid line in Figure 2.

Zoom In Zoom Out Reset image size

However, as can be seen in Figure 2, the flux uncertainties are relatively large. This calls for a test of the reliability of the correlation strength considering the errors of individual data points. To this end, we use the z-Transformed Discrete Correlation Function (ZDCF), which determines 68% confidence level intervals for the correlation coefficient from unevenly sampled data (see, e.g., Edelson & Krolik 1988; Alexander 1997). The Fisher z-transform of the linear correlation coefficient is used to estimate the 68% confidence level interval. Applying the ZDCF to the X-ray and VHE light curves reported here we obtain the following uncertainties for the linear correlation coefficient: r = 0.81^+0.06_−0.21. The sensitivity of the MAGIC telescope requires observation times of several hours to be able to get few-sigma detections for the VHE fluxes from LS I +61 303. During such long periods, the X-ray fluxes have quite large intrinsic variability (see Table 2). This results in quite large uncertainties on both the VHE and X-ray fluxes that restrict the ability to determine a possible correlation out of the outburst.