© ESO, 2012

1. Introduction

Owing to its dramatic history (the “great eruption” of 1843, subsequent disappearance, envelopment in a huge bipolar nebula, etc.), η Carinae (η Car) has been of interest for astronomical and astrophysical research for more than a century. Recent observations by AGILE (Tavani et al. 2009) and Fermi-LAT (Abdo et al. 2010b) established high-energy γ-ray emission in spatial coincidence with the location of η Car. Physical association, e.g. through observation of correlated multiband variability or orbital modulation of the γ-ray signal, has not yet been confidently presented. If this identification can be established, the η Car system would be the first massive star colliding-wind binary system (CWB) in the γ-ray sky, in contrast to known γ-ray binaries, which all contain a compact object.

It is widely accepted (e.g., Damineli et al. 2008) that η Car is a binary system consisting of two massive stars (one being a member of the rare stellar class of luminous blue variables (LBV), the other being an O or WR star) with an orbital period of 5.54 yr. Stellar and orbital parameters as determined by various studies are listed in Table 1. Both massive stars in the η Car system are expected to produce powerful stellar winds. Mass-loss rates and terminal wind velocities (as given in the table) are thought to be sufficiently high to form a wind-wind collision zone of shocked, hot gas, wherein particle acceleration (as illustrated by e.g., Eichler & Usov 1993; Dougherty et al. 2003) in general and subsequent γ-ray emission in particular (as described by e.g., Reimer et al. 2006; Bednarek & Pabich 2011) can occur. The conditions in the wind-wind collision zone depend on the orbital phase of the binary system. Therefore, these models predict γ-ray emission that is generally modulated on orbital time scales. In addition there may exist complex dependencies, such as spectral cutoffs in the underlying particle spectrum due to excessive energy losses, anisotropy in the Inverse Compton process, or types of particle transport.

Previous analyses of Fermi-LAT data (Abdo et al. 2010b; Farnier et al. 2011) have elaborated on the catalogued Fermi-source 1FGL J1045.2−5942 (Abdo et al. 2010a) that was found to be spatially consistent with the location of η Car. Following the nomenclature of the new Fermi-LAT two-year catalogue (Nolan et al. 2012), the source designation 2FGL J1045.0−5941 is used from here on. The high-energy spectrum has been best described by an exponential cutoff power-law function (CPL) that appears to be dominant below 10 GeV, and an additional component described by a power-law function (PL) at higher energies. Neither study presented a flux time history that would unambiguously relate the observed γ-ray emission to a physical origin in the η Car system. Whereas there is evidence for nonthermal emission of η Car at energies well below the MeV range (Leyder et al. 2008; Sekiguchi et al. 2009), no γ-ray emission at energies above the energy range of the Fermi-LAT instrument has been reported from the vicinity of η Car (HESS Collaboration et al. 2012).

2. Observations with the Large Area Telescope

2.1. Dataset

The Fermi Large Area Telescope (LAT) started regular survey-mode observations on 2008 August 4, providing a complete, regular sky-coverage with approximately uniform exposure. Accordingly, the source of interest was monitored regularly over the whole time range of the dataset that is used in this analysis (~35 months, 2008 August 4 to 2011 July 5). The data were reduced and analysed using the Fermi Science Tools v9r23 package¹. To reject atmospheric γ-rays from the Earth’s limb, we excluded time periods in which the region was observed at a zenith angle greater than 100° and for observatory rocking angles greater than 52°. All remaining photons with energy E > 200 MeV within a square of 21° base length (aligned in directions of right ascension α and declination δ) centred on the nominal location of η Car, (α,δ) = (161.265°,−59.685°) were used.

2.2. Likelihood analysis

Because η Car is located at low Galactic latitude and within the projection of the Carina-Sagittarius arm of the Milky Way, careful modelling of the diffuse γ-ray emission is required. Due to the size of the Fermi-LAT’s point-spread function (PSF) – the 68% containment radius is more than 2° at 200 MeV – possible contributions from nearby sources have to be taken into account as well. We performed maximum likelihood analyses using the instrument response function P7SOURCE_V6 in conjunction with the Galactic diffuse model gal_2yearp7v6_v0.fits and an isotropic background component iso_p7v6source.txt. The applied source model includes 38 point sources that are located within a radius of 15° around the source of interest and are listed in the Fermi-LAT two-year (2FGL) catalogue (Nolan et al. 2012). For each source the same spectral model as in the catalogue was used. The spectral parameters of all sources within a radius of 10° were allowed to vary in the likelihood analysis. Most sources are modelled either by a simple PL (1)(where A is the normalization parameter and Γ is the spectral index) or by a LogParabola function, which provides additional information on the spectral curvature of some sources by introducing a break energy E_b and the parameters α and β, (2)Exceptions are γ-ray pulsars, which are modelled by a CPL functional expression, (3)where E_cutoff is the cutoff energy (other parameters defined as above.) In all three cases gives the differential flux.

Source detection significance can be described by the likelihood test statistic value TS =  −2ln(L_max,0/L_max,1), which compares the ratio of two values that are obtained by a maximum likelihood procedure. L_max,0 is the likelihood for a model without an additional source at a specified location (the null-hypothesis) and L_max,1 is the maximum likelihood value for a model including an additional source. Surrounding sources and background are taken into account in both cases. The same notion can also be applied in comparing the likelihood of two source models that are described by different spectral models for a specific source of interest. Because the TS-value is expected to follow a χ²-distribution for the difference in the number of degrees of freedom between the two models (Mattox et al. 1996), it can be converted into the detection significance σ for a point source.

3. Analysis details and results

3.1. Spatial analysis

In the two Galactic coordinate grids in Fig. 1, the TS-values indicate detection significances for a trial additional point source (a PL with two free parameters) considered sequentially at each of the grid positions. The left figure covers the region around 2FGL J1045.0−5941 at energies 0.2 to 10 GeV, the right figure at energies 10 to 300 GeV – each representing 35 months of data. The energy band has been divided according to the previous indication of two different emission components. Both TS-maps were obtained with the tool gttsmap. In addition, confidence contours mark the 68.3%, 95.4% and 99.7% uncertainty regions for the location of the source defined by the maximum TS value in the grids. For both energy bands the nominal position of η Car is found within the 95.4% error region. For the low band it is contained within the 68.3% error region. The formal significance for the low- and high-energy band is 46σ and 12σ, respectively.

In the lower-energy band the γ-ray signal is centred at (α,δ) = (161.243°, −59.681°) with a fairly circular 95% confidence region of radius r = 0.023°. In the higher-energy band the signal is centred at (α,δ) = (161.209°, −59.710°) with an elliptical 95.4% confidence region that has a semimajor axis of a = 0.053° (NE-SW oriented) and a semiminor axis of b = 0.045°. The interplay between a vast photon number at lower energies and statistical limitations at higher energies yields a smaller confidence region for the hard band. The angular separations d of the maximum likelihood locations for the low- and high bands and the nominal position of η Car are d_high−low = 0.033°, d_{high−nominal} = 0.038°, d_{low−nominal} = 0.012°. The centroid of the low-energy band is well inside the 95.4% error region of the high-energy band. Although the high band’s centroid lies slightly outside the low band’s 99.7% error circle, the error regions substantially overlap. The difference to the indicated catalogue position of 2FGL J1045.0−5941 is understood by acknowledging the difference in exposure and the decomposition into low- and high-energy bands. We stress that all aforementioned positions sufficiently agree to allow interpretation in a single-source scenario.

Considering the possibility of observing the γ-ray emission superimposed by two sources in close spatial coincidence, we consider this to be extremely unlikely. Although the nominal distance between the centroid of the observed low- and the high-energy emission would not exclude source confusion between a low-energy source with exponential cutoff and a hard-spectrum source at high energies, the probability of finding two spatially coincident sources that exhibit a similar flux modulation pattern (see Sect. 3.3), characteristic only for objects in the rare class of long-period γ-ray binaries, appears to be the main argument for rejecting this hypothesis.

3.2. Spectral analysis

The energy spectrum of 2FGL J1045.0−5941 has been determined by a maximum-likelihood fitting technique using the tool gtlike (dataset and source model as described above). The energy bins were chosen such that a minimum detection significance of 10σ per bin is required (logarithmic scale with three bins between 0.2 GeV and 1 GeV and five bins between 1 GeV and 10 GeV). Owing to statistics, the criterion was lowered to require at least 5σ detections per energy bin in the range from 10 GeV to 100 GeV. The highest energy bin (100 GeV to 300 GeV) corresponds to a 3σ detection significance. For each bin a single PL shape was assumed. The normalization parameter and the index were allowed to vary during the fitting process. The effective energy of each bin is defined as the midpoint of the integral photon flux in the bin; therefore each energy bin carries some information about the underlying photon spectrum.

A spectral fit over the energy range (0.2 to 100 GeV) was obtained by fitting a source model with a CPL+PL function (both at the catalogue position of 2FGL J1045.0−5941) to the dataset up to 100 GeV using the maximum-likelihood tool gtlike. Events at energies E > 100 GeV have not been considered in fitting spectral models, because they do not significantly influence the results. Figure 2 shows the energy spectra of 2FGL J1045.0−5941 along with the overlaid broad-band spectrum. The CPL+PL multi-component functional form has been applied previously by Abdo et al. (2010b) and continues to yield better spectral representation than a single PL or CPL function. For comparison to a signal without a second emission component, the spectrum of the nearby pulsar PSR J1048−5832 (1.2° angular distance) is shown, too. The best-fit parameters are a spectral index Γ = 1.97  ±  0.05 and energy cutoff E_cutoff = 3.5  ±  0.1 GeV with an integrated flux of  cm^-2 s^-1 for the CPL, which dominates the low-energy band and Γ = 1.94    ±    0.04 and an integrated flux of  cm^-2 s^-1 for the PL dominant in the high-energy band. The corresponding formal significances of the two spectral components are 30σ and 16σ, respectively.

3.3. Temporal analysis

We now investigate the flux time history of the observed γ-ray emission. To search for coarse temporal signatures, two adjacent energy intervals were analysed: 0.2 to 10 GeV and 10 to 300 GeV, respectively. This choice corresponds to the apparent structure in the spectral energy distribution (SED) at 10 GeV as evident in Fig. 2, which marks the transition between two spectral components. On the basis of 35 months of Fermi-LAT data it is not useful to search the power spectrum for the existence of the characteristic  ~5.5 yr orbital period in the η Car system. Therefore we restrict the temporal analysis to flux variability investigations, conducted over the whole energy interval and in the two energy regimes separately.

3.3.1. Flux studies in the lower-energy band (0.2 to 10 GeV)

In Fig. 3, left, we show the flux time history of 2FGL J1045.0−5941 for the energy band 0.2 to 10 GeV as obtained by likelihood analysis. The dataset was divided into fourteen consecutive intervals, each representing 2.5 months. For each time bin a CPL function (all parameters free) was fitted to the data. The average photon flux is (1.29  ±  0.01) × 10^-7 cm^-2 s^-1. Taking this value as the hypothesis for a non-varying source, a χ²-test gives a probability of 99.986% (corresponding to 3.8σ) that this hypothesis is false. Conducting the same study on the two nearby pulsars PSR J1048−5832 and PSR J1044−5737 (at angular separation of 1.2° and 2.0°), which are expected to show a steady flux behaviour, we obtain merely a 1.3σ and 0.8σ significance of null-hypothesis violation.

We conclude that in contrast to the nearby pulsars, the γ-ray source 2FGL J1045.0−5941 does show a significant degree of variability in the low-energy band. Detailed analysis shows that this is mainly because of the low flux values obtained for the most recent 12.5 months of our data sample. The mean flux for this period lies about a factor 1.7 below the mean flux of the remaining data sample. The nearby pulsars do not show this downward trend.

The corresponding cumulative TS-value evaluation for 2FGL J1045.0−5941 is shown in Fig. 3, right. A non-variable source is expected to show a linear increase in the cumulative TS-value (Mattox et al. 1996) and deviations from that can be an indication of statistical fluctuations, stochastic flux variability, or periodic flux variation. In addition to minor deviations, the figure indicates an approximately linear rise in the cumulative TS-value for the first 22.5 months of 2FGL J1045.0−5941 data (blue dotted line). Then, however, the growth in TS decreases in correspondence to the lower flux values in Fig. 3, left. Thus, the conclusions motivated by the light curve are reinforced by the cumulative TS evaluation.

3.3.2. Flux studies in the higher-energy band (10 to 300 GeV)

Figure 4 is the equivalent of Fig. 3 for the energy regime between 10 and 300 GeV. Due to the lower statistics at higher energies, the time binning was increased from 2.5 to 5 months. The first two bins – representing the flux of the first ten months of Fermi-LAT observation – are about a factor of 2.5 higher than the following five bins. For the most recent five months only an upper limit could be obtained. A χ²-test shows that a null-hypothesis of observing a constant flux can be rejected with 7.4σ significance in favour of flux variability.

This conclusion is supported by the evaluation of the corresponding cumulative TS-values. Due to the cut-off spectra of the nearby pulsars there is no other source in the vicinity of 2FGL J1045.0−5941 that can serve for comparison in this energy interval. Yet, the behaviour of the TS-value is indicative of γ-ray variability: whereas the cumulative TS-value constantly rises seemingly linearly and steadily over the first ten months (blue dotted line), it flattens and turns over afterwards. A line connecting the origin of the graph and the last cumulative TS-value (green dashed line) would indicate steadily increasing detection significance over time and thus absence of flux variability. Because this is not the case, we have to conclude that the cumulative TS-curve is indicative of flux variability. However, this indication builds up only after the first ten months of data; after an initial period of constant flux, the γ-ray emission at energies 10 to 300 GeV begins to fade out.

3.3.3. Spectral evolution

Figure 5 shows two SEDs (determined as in Sect. 3.2). The SED on the left represents the first ten months of our dataset for which an increased flux was observed at energies 10 to 300 GeV. The SED on the right covers the remaining 25 months. As suggested by Fig. 4, the data points of the hard component appear to be systematically lower in the latter part of the data set. Similarly as in Sect. 3.2, a CPL+PL was fitted to the data for both time intervals using energies below 100 GeV. The most substantial change is to be seen in the high-energy PL component in which the index drops from Γ = 1.41    ±    0.08 for the first ten months to Γ = 2.06    ±    0.06 for the latter 25 months. Using the obtained spectral parameters to derive the integrated flux from 10 to 100 GeV (the contribution of both spectral functions are taken into account), we obtain  cm^-2 s^-1 and  cm^-2 s^-1, which implies a modulation by a factor of about 2.4. The difference in flux therefore agrees well with the findings of Sect. 3.3.2.

3.3.4. Connection to the orbit of the η Car system

As indicated by the vertical line in Fig. 4, the period of increased flux in the high-energy component (10 to 300 GeV) at the start of Fermi-LAT observations corresponds to the periastron passage of the η Car binary system. To further illustrate this relation, Fig. 6, left, shows the flux of 2FGL J1045.0−5941 at energies from 10 to 300 GeV over the orbital phase of η Car. The periods of increased flux clearly correspond to orbital phases close to periastron. Figure 6, right, shows flux in relation to the physical distance of the two stellar components (ranging from 1.5 AU at periastron to  ~ 30 AU at apastron, which was reached by October 2011). The intervals of increased flux correspond to a time in which the separation of the stars is less than about 12 AU, with a flux maximum (bin 1) just prior to the periastron passage. Otherwise the data indicate a clear decrease in flux as the binary components approach apastron. The same is shown in Fig. 7 for the low-energy component (0.2 to 10 GeV). In this case the flux decrease seems to commence later, with periods of lower flux corresponding to stellar separations of more than 25 AU.

4. Discussion and summary

As the analysis above has shown, we see a clear flux decrease in relation to the orbital state of η Car. The decrease is significant at low- (0.2 to 10 GeV) as well as high- (10 to 300 GeV) energies, albeit showing a different time dependence (see above). In general, the observed flux decrease can be qualitatively understood with predictions in the framework of colliding-wind binary models for high-energy γ-ray emission (Reimer et al. 2006; Bednarek & Pabich 2011). As the stars move away from each other, the matter density – as well as the radiation density – decreases in the wind collision region where particle acceleration and subsequent γ-ray emission are thought to occur. Although the observed flux variation does not yet allow an unambiguous association, we consider in the following that the detected emission indeed arises from η Car for reasons explained above.

Note that source modulations of the two spectral components at the same amplitude can be fully recovered for the two energy intervals only if the model used for the diffuse gamma-ray emission and for nearby point sources is perfect (e.g., without statistical and systematic uncertainties). Any realistic assessment of the backgrounds (e.g., for statistical fluctuation or, even more, systematic deviations between model and reality) will have an impact on the principal ability to recover source modulation. This will matter comparably less at higher energies due to the better signal-to-noise ratio. Hence, our inability to produce a perfect diffuse emission model principally limits our capability to recover variability in the low-energy component of η Car.

The presence of phase-locked flux changes in any emission component refutes the idea of Ohm et al. (2010), who suggested that a γ-ray signal from η Car can be related to particle acceleration in the expanding blast wave that originated in the great eruption of 1843 when the primary star (the LBV) ejected about 10% of its total mass in a massive outburst. Through shock interaction with the surrounding interstellar medium the blast wave that was created at that event might produce charged particles and subsequent γ-ray emission. In such an emission scenario one would not expect to see γ-ray flux variability in connection with the orbital configuration of the η Car binary system.

Based on the observations with the LAT of two apparently distinct emission components, a model has been proposed (Farnier et al. 2011) to identify the low-energy component as due to inverse Compton scattering of η Car’s photospheric photon field off a relativistic electron population extending up to  ~10⁴ MeV, and the high-energy component above  ~10 GeV as due to hadronic interactions of a relativistic proton component extending up to  ~10⁴ GeV with a high-density material in the shock region (enhanced by a factor  ~10 with respect to the unperturbed wind density). This requires electron and proton acceleration at maximum rate (Bohm diffusion regime) for values of the shock magnetic field in the sub-Gauss range. The target material and photospheric photon density at the shock location change with orbital phase as the inverse square of the distance between the shock and the LBV. Therefore, phase-locked orbital variations of both inverse Compton and hadronic components are expected unless the intensities of the putative accelerated particle populations change along the orbit such that it counteracts this trend. We stress here that both components are expected to vary significantly with the same amplitude for an isotropic target radiation field and no significant orbital variations in the accelerated electron-to-proton number densities. If anisotropy effects were considered, even larger variability amplitudes are expected in the inverse Compton component than in the hadronic component. Applied to the energy spectra determined for the first 10 months and the subsequent 25 months of LAT observations (see Fig. 5) for which the mean stellar separation were  ~8 AU and  ~21.8 AU, respectively, the flux amplitude is expected to decrease from the first to the second time period by more than a factor 7 in both spectral components, possibly by an even greater factor for the leptonic component. The observed spectra, however, reveal flux variations in the high-energy component by a factor 2.5 while the low-energy component shows a rapid spectral decline at a few GeV with no indication for a phase dependence of the spectral shape, spectral index, or cutoff energy.

We note that the observed low variability of the low-energy γ-ray component with phase-locked variations of the high-energy component can be plausibly interpreted as inverse Compton scattering in the photospheric radiation field of η Car and hadronic interactions with the wind material, respectively, if the accelerated electron number density is higher around apastron than periastron while at the same time the accelerated proton number density does not exhibit significant orbital variations.

Motivated by the peculiar shape of the measured spectra and the complex strong radiation fields in the vicinity of η Car, we here consider yet another scenario. Photon-photon absorption has long been suspected to play a significant role in the dense radiation environments of massive stars. For instance, first indications of γ-ray absorption modulating a TeV spectrum have been reported from the γ-ray binary LS 5039 (Aharonian et al. 2005).

For η Car, the observed γ-ray spectra can be well represented by a primary differential photon spectrum  ∝ E^-2exp(−E/E_cut) with E_cut ~ 250−500 GeV that suffers from γ-ray absorption in a black-body photon field (Gould & Schréder 1967) of effective temperature  ~5    ×    10⁵ K (corresponding to target photon energies of  ~116 eV at the maximum of the black-body radiation field) during the first ten months of observations, and decreases to  ~2    ×    10⁵ K (~46 eV) during the subsequent 25 months (see Fig. 8). The required maximum optical depths are τ_{γγ,   max} ≈ 1.1 and τ_{γγ,   max} ≈ 2, corresponding to fluxes (at the distance of η Car) of  ~6    ×    10^-8 erg cm^-2 s^-1 and  ~9    ×    10^-11 erg cm^-2 s^-1, respectively. Considering the total spectrum (35 months), the required black-body absorber has an effective temperature of  ~3.1    ×    10⁵ K (~72 eV) and τ_{γγ,   max} ≈ 1.3 corresponding to a flux of  ~2    ×    10^-9 erg cm^-2 s^-1. In this scenario, a possible origin of the absorber could be hot X-ray gas, as observed recently by XMM-Newton and Chandra (e.g., Hamaguchi et al. 2007; Townsley et al. 2011). Both experiments revealed spatially extended, structured X-ray emission components that surround the binary system. Figure 9 (left) illustrates the geometry of such an external absorber scenario. The change of the LAT spectrum would indicate a slightly hotter rarefied plasma towards periastron than towards apastron.

In a second scenario the absorber could be linked to η Car’s wind collision (e.g., hot shocked gas). In this case (see Fig. 9 (right)) the angular dependence of the absorption coefficient and threshold energy on the photons’ collision angle has to be taken into account. The resulting optical depth therefore depends on the orientation of the wind collision region, the line of sight, and the stars of the η Car binary system. For example, for a mean angular phase for the first 10 and subsequent 25 months of 0.973 and 0.422, respectively, and an inclination angle of 45° and ω = 285°, the propagation angle (between the collision region and line of sight) for target photons impinging from the direction of η Car is  ~73° and  ~119°, respectively. In this case we find satisfactory spectral representations for an absorber of temperature T ≈ 10⁶ K and τ_{γγ,   max} ≈ 1 for the first 10 (see Fig. 8, left), and T ≈ 2.5    ×    10⁵ K and τ_{γγ,   max} ≈ 2 for the latter 25 months (see Fig. 8, right).

Orbitally modulated γ-ray emission up to TeV energies has been observed from a number of high-mass X-ray binary systems (see e.g., Hill et al. 2011) with orbital periods ranging from 3.4 years (PSR B1259−63) down to a few weeks or days (e.g., LSI +61 303, LS 5039). For the latter two systems, peak fluxes at MeV/GeV energies occur close to superior conjunction where the companion is behind the star (which also turns out to be close to periastron passage), while the flux decreases towards inferior conjunction/apastron. This behaviour is similar to the observed flux modulation in the spectrum of η Car despite the lack of detection at VHE γ rays.

Ultimately, if one concludes that the observed decrease in γ-ray flux indeed relates to orbital conditions in the η Carinae system, a prediction can be made for how the γ-ray flux will evolve until periastron is reached again. It is expected to remain at a low level until the end of 2013, when the flux should increase again and, ultimately, reach maximum shortly before the next periastron passage in mid-2014. The observation of this regular variability pattern will be a decisive test for current models of γ-ray emission in colliding-wind binary systems. Until then, we anticipate that numerical simulations will consider the multitude of details not yet taken into account in present γ-ray emission models (e.g., complexity in stellar radiation fields specific for the η Car system and the Homunculus nebula). This will allow a quantitative comparison of model predictions to the observed γ-ray data over a full orbit.

---

Acknowledgments

The Fermi LAT Collaboration acknowledges generous ongoing support from a number of agencies and institutes that have supported both the development and the operation of the LAT as well as scientific data analysis. These include the National Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l’Énergie Atomique and the Centre National de la Recherche Scientifique/Institut National de Physique Nucléaire et de Physique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di Fisica Nucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), High Energy Accelerator Research Organization (KEK) and Japan Aerospace Exploration Agency (JAXA) in Japan, and the K. A. Wallenberg Foundation, the Swedish Research Council and the Swedish National Space Board in Sweden. Additional support for science analysis during the operations phase is gratefully acknowledged from the Istituto Nazionale di Astrofisica in Italy and the Centre National d’Études Spatiales in France. The publication is supported by the Austrian Science Fund (FWF).

References

All Tables

All Figures