EVIDENCE FOR SECONDARY EMISSION AS THE ORIGIN OF HARD SPECTRA IN TeV BLAZARS

Blazars, a special class of active galactic nuclei (AGNs), reveal that the continuum emission that arises from the jet emission occurring in an AGN whose jet axis is closely aligned with the observer's line of sight, is dominated by nonthermal emission as well as rapid and large-amplitude variability (Urry & Padovani 1995). The broad spectral energy distributions (SEDs) from the radio to the γ-ray bands are dominated by two components, appearing as humps. It is widely acknowledged that the first hump is produced by electron synchrotron radiation; the peaks range from the infrared–optical up to the X-ray regime for different blazars (Urry 1998). The second hump, with peaks in the GeV to TeV γ-ray band, is probably produced by inverse Compton (IC) scattering of the relativistic electrons either on the synchrotron photons (Maraschi et al. 1992) or on some other photon populations (Dermer & Schlickeiser 1993; Sikora et al. 1994). An open issue is the high energy γ-rays produced by mesons and leptons through the cascade initiated by proton–proton or proton–photon interactions (e.g., Mannheim & Biermann 1992; Mannheim 1993; Phol & Schlickeiser 2000; Aharonian 2000; Mücke & Protheroe 2001).

Observations of very high energy (VHE) γ-rays indicate that more than 40 blazars radiate γ-rays in the TeV energy region (e.g., Aharonian et al. 2005; Cui 2007; Wagner 2008). It is believed that the primary TeV photons from distant TeV blazars should exhibit clear signatures of absorption due to their interactions with extragalactic background light (EBL) to produce electron–positron pairs (e^±; e.g., Nikishov 1962; Gould & Schreder 1966). However, the observed spectra do not show a sharp cutoff at energies around 1 TeV (Aharonian et al. 2006a; Costamante et al. 2008; Acciari et al. 2009). A characteristic case is the VHE γ-ray emission in the distant blazar 1ES 1101-232, which was detected by the High-Energy Stereoscopic System (H.E.S.S.) array of Cherenkov telescopes (Aharonian et al. 2006a, 2007c). The VHE γ-ray data result in very hard intrinsic spectra with a peak in the SED above 3 TeV, corrected for absorption by the lowest-level EBL (Aharonian et al. 2007c). Similar behavior has also been detected in other TeV blazars such as 1ES 0229 + 200 (Aharonian et al. 2007a), 1ES 0347-121 (Aharonian et al. 2007b), and Mkn 501 (Neronov et al. 2011).

Generally, the lack of absorption features can most simply be explained either by assuming that there is no absorption (Kifune 1999; Stecker & Glashow 2001; De Angelis et al. 2009) or by low EBL levels (Aharonian et al. 2006b; Mazin & Raue 2007; Finke & Razzaque 2009). Alternatively, hard spectra can be expected if the γ-rays from distant blazars are dominated by secondary γ-rays produced along the line of sight by the interactions of cosmic-ray protons with background photons (Essey & Kusenko 2010; Essey et al. 2010, 2011).

AGNs are believed to be the most powerful sources of both γ-rays and cosmic rays. Recently observed results from Cherenkov telescopes indicate that interactions of cosmic rays emitted by distant blazars with the photon background along the line of sight can produce γ-rays (Essey & Kusenko 2010). Motivated by the above arguments, in this paper we study the possible origin of hard spectra in TeV blazars. The high-energy emission from TeV blazars consists of two components: the primary γ-ray component, which comes from the source, and the secondary γ-ray component, which comes from proton interactions with the EBL photons along the line of sight.

Throughout the paper, we assume the Hubble constant H₀ = 70 km s⁻¹ Mpc⁻¹, the matter energy density Ω_M = 0.27, the radiation energy density Ω_r = 0, and the dimensionless cosmological constant Ω_Λ = 0.73.

We basically follow the traditional synchrotron self-Compton (SSC) model to produce the primary γ-ray component and use the method of Kelner & Aharonian (2008) to produce the secondary γ-ray component. We now give a brief description of the model.

2.1. Primary Component Produced in the Source

The homogeneous SSC radiation model is widely used for explaining the multi-wavelength energy spectra of blazars. The homogeneous SSC radiation model that we adopt assumes a spherical radiation region filled with extreme-relativistic electrons, with a randomly originated homogeneous magnetic field and constant electron number density. We adopt a broken power-law function with a sharp cutoff to describe the electron energy distribution in the radiation region:

$$\begin{equation} N_{\rm e}(\gamma)= \left\lbrace \begin{array}{@{}ll}K_{\rm 1}\gamma ^{-n_{\rm 1}}\;,\mbox{$\gamma _{\rm min}\le \gamma \le \gamma _{\rm b}$}\;;\\ K_{\rm 2}\gamma ^{-n_{\rm 2}}\;,\mbox{$\gamma _{\rm b}\le \gamma \le \gamma _{\rm cut}$}\;, \end{array} \right. \end{equation} \tag{ 1 }$$

where γ = E_e/m_ec² is the electron Lorentz factor and K₂ = K₁γ^{(n₂ − n₁)}_b.

Based on the above electron number density N(γ), we can use the formulae given by Katarzynski et al. (2001) to calculate the synchrotron intensity I_s(E_γ) and the intensity of self-Compton radiation I_c(E_γ) (e.g., Zheng & Zhang 2011), and then calculate the intrinsic photon spectrum at the observer's frame:

$$\begin{equation} \frac{dN_{\rm \gamma }^{\rm int}}{dE_{\rm \gamma }}=\pi \frac{R^2}{d^2E_{\rm \gamma }^{2}}\delta ^3(1+z)[I_{\rm s}(E_{\rm \gamma })+ I_{\rm c}(E_{\rm \gamma })]\;, \end{equation} \tag{ 2 }$$

where d is the luminosity distance, z is the redshift, and δ = [Γ(1 − βcos θ)]⁻¹ is the Doppler factor where Γ is the blob Lorentz factor, θ is the angle of the blob vector velocity to the line of sight and β = v/c. Since at high energies Compton photons can produce pairs by interacting with synchrotron photons, this process can decrease the observed high energy radiation (Coppi & Blandford 1990; Finke et al. 2009). Katarzynski et al. (2001) analyzed the absorption effect due to pair production inside the source, and found that its effect is almost negligible. On the other hand, VHE γ-photons from the source are attenuated by photons from the EBL. Therefore, after taking into account the absorption effect, the flux density observed at Earth becomes

$$\begin{equation} \frac{dN_{\rm \gamma }^{\rm obs}}{dE_{\rm \gamma }}=\frac{dN_{\rm \gamma }^{\rm int}}{dE_{\rm \gamma }}\exp [-\tau (E_{\rm \gamma },z)]\;, \end{equation} \tag{ 3 }$$

where τ(E_γ, z) is the absorption optical depth due to interactions with the EBL (Kneiske et al. 2004; Dwek & Krennrich 2005). In our calculation, we use the absorption optical depth which was deduced from the average EBL model by Dwek & Krennrich (2005).

2.2. Secondary Component Produced Outside the Source

AGNs are expected to accelerate cosmic rays to energies up to ∼10¹¹ GeV. For energies below the Greisen–Zatsepin–Kuzmin cutoff of about 50 EeV (Greisen 1966; Zatsepin & Kuzmin 1966), the cosmic rays can propagate out to cosmological distances without significant energy loss, then interact with the EBL relatively close to Earth. The secondary γ-ray production by interactions of cosmic rays emitted by distant blazars with the photon background along the line of sight mainly depends on pion decay. Cosmic-ray interactions with the EBL should also produce neutrinos and electrons. We will concentrate on photon production and not discuss other particles in the proton–photon interaction processes.

Let f_p(E_p) and f_ph() be functions characterizing the energy distributions of initial protons and soft photons. The production rate of γ-rays can be obtained by (Kelner & Aharonian 2008):

$$\begin{equation} \frac{dN_{\rm \gamma }}{dE_{\rm \gamma }}=\int f_{\rm p}(E_{\rm p})f_{\rm ph}(\epsilon)\Phi _{\rm \gamma }(\eta , x)\frac{dE_{\rm p}}{E_{\rm p}}d\epsilon \;, \end{equation} \tag{ 4 }$$

where η = 4E_p/(m²_pc⁴), x = E_γ/E_p, and Φ_γ(η, x) is a piecewise function of two variables. According to the results on photomeson processes, which are obtained using numerical simulations of proton–photon interactions based on the publicly available Monte Carlo code SOPHIA (Mücke et al. 2000), Kelner & Aharonian (2008) give an approximate analytical presentation, namely,

$$\begin{equation} \Phi _{\rm \gamma }(\eta ,x)=\left\lbrace \begin{array}{@{}lc}B_{\rm \gamma }^{}\,\exp \left\lbrace -s_{\rm \gamma }^{} \left[\ln \left(\frac{x}{x_-}\right)\right]^{\delta _{\rm \gamma }}\right\rbrace \nonumber \\ \ms \times \left[\ln \left(\frac{2}{1+y^2}\right)\right]^{2.5+0.4\ln (\eta /\eta _{0})}\,,& x_{-}<x<x_{+}\;,\\ \ms B_{\rm \gamma }\left(\ln 2\right)^{2.5+0.4\ln (\eta /\eta _{0})}\,,& x\leqslant x_{-}\;,\\ \ms 0\,,& x\geqslant x_{+}\;,\\ \end{array} \right. \end{equation} \tag{ 5 }$$

where y = (x − x₋)/(x₊ − x₋), η₀ relate to the proton mass m_p and the π-meson mass m_π, η₀ = 2m_π/m_p + m²_π/m²_p, and $x_{\pm }=({1}/{2(1+\eta)})\,[\eta +r^{2}\pm \sqrt{(\eta -r^{2}-2r)\,(\eta -r^{2}+2r)}]$ with r = m_π/m_p ≈ 0.146. All three parameters B_γ, s_γ, and δ_γ used in the above presentation are functions of η. The numerical values of these parameters are given in Kelner & Aharonian (2008).

Instead of integrating Equation (4) over d, it is more convenient to integrate over dη. This allows us to rewrite the spectra of the γ-ray photons in the form

$$\begin{equation} \frac{dN_{\gamma }}{dE_{\gamma }}=\int _{\eta _{0}}^{\eta _{\rm max}}H(\eta ,E)\,d\eta \,, \end{equation} \tag{ 6 }$$

where

$$\begin{equation} H(\eta ,E_{\gamma })=\frac{m_{p}^{2}c^{4}}{4}\!\int _{E}^{\infty }\!\frac{dE_{p}}{E_{p}^{2}}\, f_{p}(E_{p})\, f_{\rm ph}\!\left(\!\frac{\eta m_{p}^{2}c^{4}}{4E_{p}}\!\right)\! \Phi _{\gamma }\!\left(\!\eta ,\frac{E_{\gamma }}{E_{p}}\right)\!. \end{equation} \tag{ 7 }$$

1ES 1101-232 resides in an elliptical host galaxy at a redshift of z = 0.186 (Remillard et al. 1989; Falomo et al. 1994). The source has been classified as a high-frequency peaked BL Lac object, because of the dominance of synchrotron emission in the X-ray band (Donato et al. 2001). In a previously published SSC frame, the broadband characteristics of 1ES 1101-232 indicated that the IC peak was generally expected to be around 100 GeV (e.g., Wolter et al. 2000; Costamante & Ghisellini 2002). However, new observational results with H.E.S.S. in 2004 and 2005 indicated that the source exhibits hard intrinsic spectrum with a peak in the SED above 3 TeV, corrected for absorption by the lowest-level EBL (Aharonian et al. 2007c). Using the model in Section 2, we can calculate the TeV γ-ray spectra in the source (primary component) and outside the source (secondary component). Then hard intrinsic spectrum of the source can be produced. The two epochs will be considered independently.

In order to do so, we first search for the primary γ-ray component in the one-zone SSC frame. Assuming the electron Lorentz factors γ_min, γ_break, γ_max are identical in the two epochs, we calculate the high-energy electron distribution with a broken power law between γ_min = 1 and γ_max = 8.0 × 10⁸ with a break at γ_break = 9.0 × 10⁴, where for the 2004 observed data, the density normalization K = 1200 cm⁻³, the energy index of the particles between γ_min and γ_break is set to n₁ = 2, and the energy index of the particles between γ_break and γ_max is set to n₂ = 3.5. For the 2005 observed data, the density normalization K = 3.5 cm⁻³, the energy index n₁ is set to 1.8, and n₂ is set to 4.05. The parameters are applied as follows. In the first observation epoch, the magnetic field strength is B = 0.55 G, the emission region size is R = 1.28 × 10¹⁶ cm, and the Doppler factor is δ = 10.5. In the later observation epoch, in order to obtain good fits, the magnetic field strength is B = 0.15 G, the emission region size is R = 1.65 × 10¹⁷ cm, and the Doppler factor is δ = 10. All the physical parameters of the one-zone SSC spectra are listed in Table 1.

We assume that relativistic electrons are in a steady state during the observational epoch. Therefore, we can calculate the X-ray/TeV γ-ray spectrum in the one-zone SSC model using the broken power-law electron spectrum. In Figure 1, we show the predicted spectrum from the X-ray to TeV γ-ray bands (blue solid curve). For comparison, the observed data of 1ES 1101-232 at the X-ray and TeV bands on the 2004 June 5–10 and 2005 March 5–16 respectively (Aharonian et al. 2007c) are also shown. It can be seen that the lower energy observed data can be reproduced in the SSC model.

Zoom In Zoom Out Reset image size

We now consider the hard spectrum properties of 1ES 1101-232 in 2004 June and 2005 March, especially the shape of the VHE γ-ray tail of the observed spectra. In order to do this, we calculate the spectra of the secondary γ-rays numerically using the EBL spectrum n(z, ) which was deduced from the average EBL model by Dwek & Krennrich (2005) at the redshift z. The energy losses are due to production of pions in proton–photon (p − γ) interactions with EBL photons. This process depends on the proton injection spectrum (see Table 2), which we parameterize by a constant power-law exponent α and maximal energy E_max:

$$\begin{equation} f_{\rm p}(E_{\rm p})=N_{0}\,E_{\rm p}^{\alpha }{\rm exp\left(-\frac{{\it E}_{\rm p}}{{\it E}_{\rm p, max}}\right)},\; E_{\rm p, min}\leq E_{\rm p} \leq E_{\rm p, max} \end{equation} \tag{ 8 }$$

where the coefficient N₀ is the initial proton distribution.

We assume the minimum energy of the injection proton to be of the order of 0.145 GeV, which is a kinematic threshold of the photomeson production process in proton–photon interactions in two identical epochs, and the proton energy loss is by photomeson production at the pc scales (Blumenthal 1970). We calculate the high-energy injection proton spectra with N₀ = 10.0 erg⁻¹ cm⁻³, E_pmax = 1.28 × 10⁵ GeV for the 2004 observed data, and with N₀ = 2.0 erg⁻¹ cm⁻³, E_pmax = 3.02 × 10⁵ GeV for the 2005 observed data. From the above proton spectra, we reproduce the observed TeV photon spectrum (dashed curve) of 1ES 1101-232 in 2004 June 5–10 and 2005 March 5–16 in Figure 1, respectively. It can be seen that the observed hard spectrum properties of 1ES 1101-232 in 2004 June and 2005 March, especially the shape of the VHE γ-ray tail of the observed spectra, can be reproduced in our model.

Generally, proton–proton (p−p) interactions do not give rise to efficient γ-ray production mechanisms in a jet. This mechanism could be effectively realized only in a scenario that assumes that γ-rays are produced in dense gas clouds that move across the jet (e.g., Morrison et al. 1984; Dar & Laor 1997). For example, in order to interpret the reported TeV flares of Markarian 501 by π⁰-decay γ-rays produced by p−p interactions, for any reasonable acceleration power of the protons L_p ⩽ 10⁴⁵ erg s⁻¹, the density of the thermal plasma in the jet should exceed 10⁶ cm⁻³ (Aharonian 2000). On the other hand, under the conditions of existence of extremely high energy, E > 10¹⁹ eV, and in the presence of a large magnetic field, B ≫ 1 G, the synchrotron radiation of the protons becomes a very effective channel for the production of high-energy γ-rays. In our calculations, in order to reproduce the observed spectra of 1ES 1101-232, we adopt a lower proton energy and magnetic field. These postulates lead to a longer proton lifetime with t_sy = 4.5 × 10⁴ B₁₀₀ ⁻²E⁻¹₁₉s (Aharonian 2000), where B₁₀₀ = B/100 G, E₁₉ = E/10¹⁹ eV, and fainter radiation in the jet.

The proton-induced cascade process (Mannheim 1993, 1996) is another attractive possibility for the production of high-energy γ-rays. This process relates the observed γ-ray radiation to the development of pair cascades in the jet triggered by secondary photomesons produced by interactions of accelerated protons with low-frequency synchrotron radiation in the source or EBL photons outside the source. For a low-energy target photon field, the photomeson cooling time of the protons can be estimated using the approximate formula t_pγ ∼ 1/〈σ_pγK_pγ〉cn_ph(ν > ν_th) (Begelman et al. 1990), where 〈σ_pγK_pγ〉 ∼ 0.7 × 10⁻²⁸ cm² is the photomeson production cross section and the inelasticity parameter is averaged over the resonant energy range (e.g., Stecker 1968; M$\ddot{\rm u}$cke et al. 1999). By simply approximating the broad synchrotron spectral component by a power-law function with the energy-flux index α = 1, and denoting its luminosity by L_s, we have n_ph(ν > ν_th) ∼ L_sE_p/(4π m_π m_e c⁵R³δ⁴) (e.g., Sikora 2010). Thus, for the parameters of 1ES 1101-232, the photomeson cooling time cannot be significantly less than light travel timescales R/c ∼ 10⁷ s. We argue that the uncooled protons can escape from the emission region, and then interact with the background photons along the line of sight.

In this paper, we develop a model for the possible origin of hard VHE spectra from a distant blazar, although, several models which could explain the very hard intrinsic blazar spectra in the γ-ray band have already been proposed (Katarzynski et al. 2006; Böttcher et al. 2008; Aharonian et al. 2008; Lefa et al. 2011; Yan et al. 2012). In this model, both the primary photons produced in the source and secondary photons produced outside it contribute to the observed high energy γ-ray emission. That is, the primary photons are produced through the SSC process, and the secondary photons are produced through high-energy proton interactions with background photons along the line of sight. Assuming suitable electron and proton spectra, we obtain excellent fits to the observed spectra of distant blazar 1ES 1101-232. This indicates that the surprisingly low attenuation of high energy γ-rays, especially the shape of the VHE γ-ray tail of the observed spectra, can be explained by secondary γ-rays produced in interactions of cosmic-ray protons with background photons in intergalactic space (Essey & Kusenko 2010; Essey et al. 2010, 2011).

In order to fit the observed data points, the protons considered here are relatively low in energy and, at these energies, photomeson production on the EBL is the most important process. However, if the energy is somewhat higher, then Bethe–Heitler pair production (Bethe & Heitler 1934) of cosmic microwave background (CMB) photons p + γ_CMB → p + e⁺ + e⁻ can dominate over p + γ_EBL → p + π. Essey et al. (2010, 2011) have taken into account both of these contributions, and argue that for proton energy E_p = 10¹⁸ eV and higher, Bethe–Heitler pair production of CMB photons contributes more than 80% of the secondary photons. If we adopt a proton spectrum with energy up to 10²⁰ eV, Bethe–Heitler pair production can scatter CMB photons of ∼10⁻⁴ eV up to ∼10¹⁴ eV (or 10²⁸ Hz) by γ²_eCMB. This energy region is in agreement with the high-energy γ-ray tail of the observed spectra.

The properties of the model mentioned above should be testable by multi-wavelength observations on a TeV blazar. Costamante (2012) argues that in several cases we have already seen the superposition of two different emission components at high electron energies, with a new component emerging over a previous SED. Since, in our case, the secondary γ-rays are produced outside of the host galaxy, this anticipates harder spectra in the TeV energy band than the GeV energy band. This should be verified in future multi-wavelength observations. Alternatively, neutrino populations can be expected in pγ interactions; we defer this possibility to future work and IceCube observations.

We thank the anonymous referee for valuable comments and suggestions. This work is partially supported by the National Natural Science Foundation of China under grants 11178019 and U1231203, and the Natural Science Foundation of Yunnan Province under grants 2011FB041. This work is also supported by the Science Foundation of Yunnan educational department (grant 2012Z016).