© ESO, 2014

1. Introduction

The system LS I +61°303 with an orbital period (Gregory 2002) consists of a compact object and a massive star with an optical spectrum typical for a rapidly rotating B0 V star (Casares et al. 2005; Grundstrom et al. 2007).

In 2009 the first detection of orbital periodicity in high-energy gamma rays (20 MeV–100 GeV) was reported by using the Large Area Telescope (LAT) from the Fermi Gamma-Ray Space Telescope spacecraft (Abdo et al. 2009). Longer monitoring has shown (Hadasch et al. 2012; Ackermann et al. 2013) that indeed the system shows a clear periodical outburst towards periastron (Φ_periastron = 0.23, Casares et al. 2005) at some epochs, but that this periodicity is not always present. This is different from the behaviour of the system in the radio band where not only periodical outbursts occur at each orbit, even if modulated with a long-term period (P_long = 1667 ± 8 d, Gregory 2002), but they also occur towards apoastron and not towards periastron as in the GeV energy band (e.g., see Fig. 2c in Massi & Kaufman Bernadó 2009).

Along with this different behaviour between high and low energy there is a puzzling overlap. Ackermann et al. (2013) noticed that GeV data also show the long-term periodical variation affecting the radio data, but only at a specific orbital phase interval, Φ = 0.5−1.0, that is around apoastron.

The aim of this paper is to investigate the discontinuity in the periodicity of the GeV γ-ray emission at periastron, the relationship of its disappearance with the variation of the emission in other parts of the orbit, and finally the possible relationship between GeV and radio emission. Section 2 describes our data analysis. In Sect. 3 we present our results and in Sect. 4 our conclusions.

2. Data analysis

In this section we present the data reduction of Fermi-LAT data performed with three packages: the Fermi science tools package (version v9r33p0)¹, the wavelet analysis², and the software Starlink³.

The γ-ray data used in this analysis span the time period MJD 54683 (August 05, 2008) to MJD 56838 (June 30, 2014). We used the script like_lc.pl by Robin Corbet⁴. Only source-event-class photons were selected for the analysis. Photons with a zenith angle greater than 100° were excluded to reduce contamination from the Earth’s limb. For the diffuse emission we used the model gll_iem_v05_rev1.fit and the template iso_source_v05_rev1.txt. We used the instrument response function (IRF) P7REP/background_rev1, and the model file was generated from the 2FGL catalogue (Nolan et al. 2012), all sources within 10° of LS I + 61°303 were included in the model. LS I + 61°303 was fitted with a log-parabola spectral shape and with all parameters left free for the fit, performing an unbinned maximum likelihood analysis. The other sources were fixed to their catalogue values. We produced light curves with a time bin size of one day and of five days. For all light curves we used an energy range of 100 MeV to 300 GeV.

The data were folded with the orbital phase defined as $\mathrm{\Phi }\mathrm{=}\frac{\mathit{t}\mathrm{-}{\mathit{t}}_{\mathrm{0}}}{{\mathit{P}}_{\mathrm{1}}}\mathrm{-}\mathrm{int}\left(\frac{\mathit{t}\mathrm{-}{\mathit{t}}_{\mathrm{0}}}{{\mathit{P}}_{\mathrm{1}}}\right)\mathit{,}$(1)where and is the orbital period of the binary system (Gregory 2002).

We investigated the temporal evolution of the orbital periodicity by means of a wavelet analysis with Morlet function (Torrence & Compo 1998). The wavelet analysis decomposes the one-dimensional time series into a two-dimensional time-frequency space and displays the power spectrum in a two-dimensional colour-plot that shows how the Fourier periods vary in time (Torrence & Compo 1998). While the wavelet analysis was applied to the γ-ray data vs. time, for a straightforward comparison with radio data, we express the x-axis as $\mathrm{\Theta }\mathrm{=}\frac{\mathit{t}\mathrm{-}{\mathit{t}}_{\mathrm{0}}}{{\mathit{P}}_{\mathrm{long}}}\mathrm{·}$(2)This allows a comparison with non-simultaneous radio data because the radio data are periodical in Θ. We will therefore compare gamma-ray data to radio data having the same fractional part of Θ. For the Lomb-Scargle timing analysis (Lomb 1976; Scargle 1982), we used the program PERIOD, which is part of the UK software Starlink. The version we used was 5.0−2 for UNIX. The wavelet analysis assumes regularly sampled data. We therefore set the data for the wavelet analysis to zero for missing flux. For the Lomb-Scargle analysis this was not necessary. As discussed in Sect. 3, the Lomb-Scargle analysis confirms and accurately determins the periodicities found with the wavelet analysis. In the Lomb-Scargle and wavelet analysis, significance levels for the spectra were determined with the Fisher randomisation, as outlined in Linnell Nemec & Nemec (1985), and with Monte Carlo simulations, as in Torrence & Compo (1998). The fundamental assumption is: if there is no periodic signal in the time series data, then the measured values are independent of their observation times and are likely to have occurred on any other order. One thousand randomized time-series were formed and the periodograms calculated. The proportion of permutations that give a peak power higher than that of the original time series would then provide an estimate of p, the probability that for a given frequency window there is no periodic component present in the data with this period. A derived period is defined as significant for p< 0.01, and a marginally significant period for 0.01 <p< 0.10 (Linnell Nemec & Nemec 1985).

Fig. 1 Folded light curve (blue) of Fig. 8 by Hadasch et al. (2012). This covers the first 8 months of Fermi-LAT observations. As comparison we give our folded data (black) from the same time interval. This time interval corresponds to Θ = 6.79−6.92.

Fig. 2 Wavelet analysis of Fermi-LAT data. The strength of periodicity is colour–coded as indicated in the bottom bar. a) Fermi-LAT data with a time bin of 5 d (black) are overplotted on Fermi-LAT data with a time bin of 1 d (red). The black line marks the point of flux change reported by Hadasch et al. (2012) which is also visible in our light curve. b) Wavelet analysis for the whole data set, that is for the whole orbital interval 0.0−1.0 (b)–d) use a time bin of one day). c) Wavelet analysis for half the data set, that is for the orbital interval Φ = 0.5−1.0, i.e., around apoastron. d) Wavelet for half the data set, that is for the orbital interval Φ = 0.0−0.5, i.e., around periastron.

3. Results

In this section we discuss our wavelet and Lomb-Scargle results for the gamma-ray data and compare them with previously published results (Hadasch et al. 2012; Ackermann et al. 2013). Then we compare gamma-ray data with radio data.

3.1. Wavelet and Lomb-Scargle analysis

As a general check of our data reduction we verified the consistency of the folded curve of our data with the folded curve of Hadasch et al. (2012). For this purpose we selected the same time interval as in Hadasch et al. (2012), that is 54 683–54 900 MJD (Θ = 6.79−6.92), covering the first eight months of Fermi-LAT observations. Figure 1 presents the result by Hadasch et al. (2012) and our results. We see how these independently calibrated data sets, which used different versions of the software, of the instrumental response funtion, etc., produce the same result: there is a main peak at orbital phase Φ ≈ 0.35.

We compare the light curve and wavelet results with previous results. The whole interval of the Fermi-LAT data used in our analysis is presented in Fig. 2a. The point of flux change reported by Hadasch et al. (2012) is also visible in our light curve. Figure 2b shows the wavelet analysis results for the whole data set. The orbital periodicity shows a minimum around Θ ≈ 7.25. This agrees with the periodograms of Fig. 7 in Hadasch et al. (2012), where the orbital periodicity is already almost absent for MJD 55 044−55 225 (Θ = 7.00−7.11) and is completely absent at MJD 55 405−55 586 (Θ = 7.22−7.33). The panels of Fig. 4 in Ackermann et al. (2013) show the decline in the orbital flux modulation in the interval MJD 55 191–55 698 (Θ = 7.1−7.5). Our Figs. 1, 2a and 2b therefore confirm the orbital modulation that peaks around periastron, the point of flux change at Θ ≈ 6.95, and the lack of orbital flux modulation at Θ ≈ 7.25 previously discovered and discussed in Hadasch et al. (2012) and Ackermann et al. (2013).

Fig. 3 Lomb-Scargle periodogram of the Fermi-LAT (with a time bin of one day). a) Full data set: data in the orbital phase Φ = 0.0−1.0. b) Zoom of Fig. 3a. c) Same as 3b for data with a time bin of 5 d. The false alarm probability (p) resulting from randomisation tests is indicated. A period is defined as significant for p< 1%, and as marginally significant for 1% <p< 10% (Sect. 2). d) Half the data set: data in the orbital phase Φ = 0.5−1.0. The periods P2 and Plong here present are typical periodicities in radio data (Massi & Jaron 2013). e) Zoom of Fig. 3d. f) Same as 3e for data with a time bin of 5 d. g) Half the data set: data in the orbital phase Φ = 0.0−0.5. h) Zoom of Fig. 3g. i) Same as 3h for data with a time bin of 5 d.

Now we examine the new results obtained when the wavelet analysis is performed separately for emission around periastron, Φ = 0.0−0.5, and emission around apoastron, Φ = 0.5−1.0. By comparing Fig. 2d with Fig. 2c it is clear that the orbital periodicity is not a characteristic of periastron emission alone, it is also present at apoastron, at least at some Θs. At Θ ≈ 7.2 where the wavelet analysis for the whole data set (Fig. 2b) shows the minimum power for the orbital periodicity, a maximum is present in Fig. 2c. Orbital flux modulation is present from Θ ≈ 6.95 to Θ ≈ 7.40, and at Θ ≈ 7.95 the periodicity is strong again, indicating a repetition consistent with the long-term period P_long. This confirms the discovery of Ackermann et al. (2013) that in the orbital phase interval Φ = 0.5−1.0 the γ-ray flux shows periodical variation with a period equal to the long-term modulation affecting the radio outburst. In addition, our wavelet analysis shows that the emission at apoastron, affected by the long-term modulation, is orbitally modulated.

We checked the presence of the orbital and long-term periodicities that were indicated by the wavelet analysis with the Lomb-Scargle timing analysis of the data from Φ = 0.5−1.0. First of all, the periodogram of the subset Φ = 0.5−1.0, shown in Fig. 3d, presents a periodicitiy at 1705 ± 355 days, confirming the result of Ackermann et al. (2013) of a long-term modulation of emission around apoastron. Moreover, the periodogram shows a feature at the orbital period P₁, confirming the results of the wavelet analysis, and this feature seems to have a complex profile. The zoomed version is presented in Fig. 3e. Figure 3f shows the periodogram for data integrated over five days. There is one peak at and another one at . When we consider the whole data set, which is the longest time interval of Fermi-LAT data analysed so far, P_long and P₂ are significant in the randomisations tests (i.e., the probability p that there is no periodic component present in the data with these periods is for both periods p< 0.01), even if they are rather weak features in the periodograms of Fig. 3a−c. For the emission around periastron, Figs. 3g, h, and i show the Lomb-Scargle periodograms of Φ = 0.0−0.5. Here the only outstanding feature is P₁.

Fig. 4 Apoastron peak. Fermi-LAT gamma-ray data folded with the orbital periodicity for three consecutive Θ-intervals. a)Θ = 7.12−7.22 (MJD 55 235−55 402); b) Θ = 7.22−7.32 (MJD 55 402−55 569); c) Θ = 7.32−7.42 (MJD 55 569−55 740). d) GBI radio data at 8 GHz folded with orbital period . The two Θ-intervals are those considered in Fig. 1 and in this figure for Fermi-LAT data. The orbital phase of the peak of the two radio curves agrees well with those of the second peak in GeV curves.

3.2. Orbital shift

Why does the timing analysis reveal a lack of orbital modulation around Θ ≈ 7.2 with a peak at periastron passage when the periodicity is still present, as shown in Fig. 2d for data at Φ = 0.0−0.5? And why does the orbital modulation of GeV emission towards apoastron get stronger in that Θ-interval (Fig. 2c)? Clearly, the second question is the answer to the first. Two peaks along the orbit disturb the timing analysis. The curves shown in Fig. 4a−c refer to the three consecutive Θ-intervals around the minimum of Fig. 2b, that is around the peak of Fig. 2c: Θ = 7.12−7.22 (MJD 55 235–55 402), Θ = 7.22−7.32 (MJD 55 402−55 569), and Θ = 7.32−7.42 (MJD 55 569–55 743). In addition to the peak at periastron Fig. 4a and Fig. 4b show a second peak in the interval Φ = 0.8−1.0.

The real question therefore is why periodical emission towards apoastron is detected only at Θ ≈ 7.2. The important information from the Lomb-Scargle analysis is that there are two periods, P₁ and P₂, as is for radio data (Massi & Jaron 2013). We therefore examined the trend of radio data in that particular Θ-interval. Figure 4d shows the GBI data at 8 GHz at Θ = 4.23−4.40. The curve shows indeed a peak at Φ ≈ 0.9 consistent with the second peak in the GeV data in the interval Φ = 0.8−1.0 (Figs. 4a,b). Figure 4d also shows data at Θ = 3.79−3.92, as the GeV data of our Fig. 1. In this case, the radio data show a peak at Φ ≈ 0.65 consistent with the bump of emission at Φ ≈ 0.65 in the gamma-ray data of Fig. 1. This is the well-known phenomenon of the orbital shift of the radio outburst in LS I + 61°303: the largest outbursts occur at orbital phase 0.6, afterwards, with the long-term periodicity, the orbital phase of the peak of the outburst changes, as analysed by Paredes et al. (1990) in terms of orbital phase shift, by Gregory et al. (1999) in terms of timing residuals, and reproduced recently by the precessing jet model in Massi & Torricelli-Ciamponi (2014, their Fig. 4b).

The second γ-ray peak, clearly associated with the radio outburst in terms of timing analysis, therefore also follows the same orbital shift. The second peak is therefore only detected at Θ ≈ 7.2, although the periodicity is indeed always present, because in that Θ-interval the peak is detached enough from the first γ-ray peak to be discernable in the timing analysis.

4. Discussion and conclusion

The first eight months of the Fermi-LAT observations (Θ = 6.79−6.93) analysed by Hadasch et al. (2012) indicate a peak of GeV γ-ray emission close to periastron. This clear orbital modulation of the GeV gamma-ray emission is lost at some time interval (Hadasch et al. 2012; Ackermann et al. 2013). The aim of this paper was to investigate the origin of this disappearance. Our results are the following:

• 1.

  In the interval Θ ≈ 7.2 (Fig. 2b) where the timinganalysis fails to find the orbital periodicity in GeV γ-ray emis-sion, there are two periodical signals, one in the orbital phaseinterval Φ = 0.0−0.5 (Fig. 2d), that is to-wards periastron, the other in the interval Φ = 0.5−1.0(Fig. 2c), that is towards apoastron. This re-sult of two GeV peaks along the orbit corroborates the two-peak accretion model for LS I + 61°303.The hypothesis that a compact object that accretes mate-rial along an eccentric orbit undergoes two accretion peaksalong the orbit was suggested and developed by several au-thors for the system LS I + 61°303 (Tayloret al. 1992; Marti& Paredes 1995; Bosch-Ramonet al. 2006; Romeroet al. 2007). The first accretionpeak is predicted to occur close to the Be star and to give rise to amajor high-energy outburst. The second accretion peak is pre-dicted to occur much farther away from the Be star, where the radiooutburst occurs, and a minor high-energy outburst is predictedthere (Bosch-Ramon et al. 2006).The predicted periastron event corresponds well to the observedGeV peak towards periastron, the second predicted high-energyoutburst, corresponds well to the here discussed apoastron peak.

• 2.

  The Lomb-Scargle analysis of emission around apoastron (Figs. 3d−f) revealed the same three periodicities , , and that affect the radio emission: , (Gregory 2002), and (Massi & Jaron 2013).

  This second result confirms the previous result of P₁ in GeV emission around apoastron, and in addition reveals P₂, only recently discovered in the radio emission by the timing analysis of 6.7 years of Green Bank Radio Interferometer (GBI) observations at the two frequencies of 2.2 GHz and 8.3 GHz (Massi & Jaron 2013). The radio data base of GBI and the data base of Fermi-LAT cover two quite different time intervals separated by 8 yr (the last scan in the GBI database is June 2000, whereas Fermi-Lat monitoring begins in August 2008). Moreover, the two monitorings have a quite different sampling rate (GBI with up to eight observations per day and large gaps, Fermi-LAT covers the whole sky over three hours, and we integrated over one day). Nevertheless, the timing analysis gives the same three periodicities (compare our Figs. 3d−f with Figs. 1 and 2b in Massi & Jaron 2013). Figures 3g−i reveal a quite different characteristic for emission around periastron: the Lomb-Scargle analysis results in only one outstanding feature at P₁ = 26.52 ± 0.08 d. The connection of P_long with P₂ is evident from Figs. 3g and h: the lack of P₂ (or reduction at noise level) is associated with a lack of P_long (reduced at noise level). Simulations (Massi & Jaron 2013) with P₁ and the long-term modulation cannot reproduce the observed periodogram (i.e., P₂), whereas simulations with P₁ and P₂ directly produce the long-term modulation as their beating frequency ν₁−ν₂.

  While P₁ is related to the periodical accretion peak towards apoastron described above, the period P₂ is most likely related to the precession of the radio jet of LS I + 61°303, see the agreement with the precessional period from VLBA astrometry, of 27–28 days (Massi et al. 2012). The hypothesis that a precessing jet, with an approaching jet with large excursions in its position angle, gives rise to appreciable Doppler boosting effects (and therefore to changes of flux density that can be detected by timing analysis) is supported by the morphology of images reported by Massi et al. (2004), Dhawan et al. (2006), and Massi et al. (2012) that showed extended radio structures changing from two-sided to one-sided morphologies at different position angles. A physical model for LS I + 61°303 of synchrotron emission from a precessing (P₂) jet, periodically (P₁) refilled with relativistic particles, has shown that the maximum of the long-term modulation occurs when P₁ and P₂ are synchronized, that is the jet electron density is at about its maximum and the approaching jet forms the smallest possible angle with the line of sight. This coincidence of the highest number of emitting particles and the strongest Doppler boosting of their emission occurs with the frequency of ν₁−ν₂ and creates the long-term modulation observed in LS I + 61°303 (Massi & Torricelli-Ciamponi 2014).

• 3.

  The folded curves of γ-ray data show that the peak at apoastron seems to be affected by the same orbital shift that affects the radio outburst. When the radio outburst occurs at orbital phase Φ ~ 0.65, the first main GeV peak at periastron shows a bump at orbital phase Φ ~ 0.65. When the radio outburst is shifted towards orbital phase Φ ≈ 0.9, the apoastron peak appears at the same orbital phases and is detached from the first main GeV peak.

We conclude that there exists a GeV peak at apoastron with the same timing characteristics and orbital shift as the radio emission. Because of the orbital shift, at some Θ-intervals this GeV peak is detached enough from the periastron γ-ray peak to become discernable as a second peak, and then it influences the timing analysis. Future analyses should investigate how the radio outburst and second GeV peak are connected; constraints to these future analyses should result from a detailed investigation of the evolution of the spectrum with orbital phase to see the increase of some spectral features at various orbital phases.

Acknowledgments

We thank Bindu Rani for reading the manuscript and useful comments, Robin Corbet for answering questions concerning the reduction of Fermi-LAT data and Chris Torrence for answering questions concerning the wavelet analysis. We thank Walter Alef and Alessandra Bertarini for their support with the computing power. Finally, we thank the anonymous referee for the careful reading of our manuscript and the useful comments made for its improvement. The wavelet software was provided by C. Torrence and G. Compo, and is available at URL: http://atoc.colorado.edu/research/wavelets/. The Green Bank Interferometer was operated by the National Radio Astronomy Observatory for the US Naval Observatory and the Naval Research laboratory during the timeperiod of these observations. This work has made use of public Fermi data obtained from the High Energy Astrophysics Science Archive Research Center (HEASARC), provided by NASA Goddard Space Flight Center.

References

All Figures

	Fig. 1 Folded light curve (blue) of Fig. 8 by Hadasch et al. (2012). This covers the first 8 months of Fermi-LAT observations. As comparison we give our folded data (black) from the same time interval. This time interval corresponds to Θ = 6.79−6.92.
In the text

Fig. 2 Wavelet analysis of Fermi-LAT data. The strength of periodicity is colour–coded as indicated in the bottom bar. a) Fermi-LAT data with a time bin of 5 d (black) are overplotted on Fermi-LAT data with a time bin of 1 d (red). The black line marks the point of flux change reported by Hadasch et al. (2012) which is also visible in our light curve. b) Wavelet analysis for the whole data set, that is for the whole orbital interval 0.0−1.0 (b)–d) use a time bin of one day). c) Wavelet analysis for half the data set, that is for the orbital interval Φ = 0.5−1.0, i.e., around apoastron. d) Wavelet for half the data set, that is for the orbital interval Φ = 0.0−0.5, i.e., around periastron.

In the text

Fig. 3 Lomb-Scargle periodogram of the Fermi-LAT (with a time bin of one day). a) Full data set: data in the orbital phase Φ = 0.0−1.0. b) Zoom of Fig. 3a. c) Same as 3b for data with a time bin of 5 d. The false alarm probability (p) resulting from randomisation tests is indicated. A period is defined as significant for p< 1%, and as marginally significant for 1% <p< 10% (Sect. 2). d) Half the data set: data in the orbital phase Φ = 0.5−1.0. The periods P2 and Plong here present are typical periodicities in radio data (Massi & Jaron 2013). e) Zoom of Fig. 3d. f) Same as 3e for data with a time bin of 5 d. g) Half the data set: data in the orbital phase Φ = 0.0−0.5. h) Zoom of Fig. 3g. i) Same as 3h for data with a time bin of 5 d.

In the text

Fig. 4 Apoastron peak. Fermi-LAT gamma-ray data folded with the orbital periodicity for three consecutive Θ-intervals. a)Θ = 7.12−7.22 (MJD 55 235−55 402); b) Θ = 7.22−7.32 (MJD 55 402−55 569); c) Θ = 7.32−7.42 (MJD 55 569−55 740). d) GBI radio data at 8 GHz folded with orbital period . The two Θ-intervals are those considered in Fig. 1 and in this figure for Fermi-LAT data. The orbital phase of the peak of the two radio curves agrees well with those of the second peak in GeV curves.

In the text