-
The Lagrangian of the ALP, including the effective ALP-photon interaction term, is
$ \begin{eqnarray} {\cal{L}}_{\rm{ALP}}&=\frac{1}{2}\partial^{\text{μ}} a\partial_{\text{μ}} a - \frac{1}{2}m_a^2a^2 -\frac{1}{4}g_{a\gamma}aF_{{\text{μ}}\nu}\tilde{F}^{{\text{μ}}\nu}, \end{eqnarray} $
(1) where a is the ALP,
$ m_a $ is its mass,$ g_{a\gamma} $ is the coupling between the ALP and photons, and$ F_{{\text{μ}}\nu} $ and$ \tilde{F}^{{\text{μ}}\nu} $ are the electromagnetic field tensor and its dual tensor, respectively. The ALP-photon system propagating along the$ x_3 $ direction is written as [21]$ \Psi = \left(A_1, A_2, a \right)^T $ , where$ A_1 $ and$ A_2 $ denote the linear polarization amplitudes of the photon in the perpendicular directions. The corresponding density matrix$ \rho=\Psi \otimes \Psi^\dagger $ satisfies the Von Neumann like equation [18]$ \begin{eqnarray} {\rm i}\frac{{\rm{d}}\rho(x_3)}{{\rm{d}}x_3}=\left[ \rho(x_3), \; {\cal{M}}_0 \right]. \end{eqnarray} $
(2) Assuming that
$ B_{\rm{T}} $ is the transversal magnetic field aligned along the direction of$ x_2 $ , the mixing matrix$ {\cal{M}}_0 $ can be described by [13, 17, 22]$ \begin{eqnarray} \begin{aligned} {\cal{M}}_0=\begin{pmatrix} \Delta_{\rm{pl }}+2\Delta_{\rm{QED }} & 0&0 \\ 0& \Delta_{\rm{pl }}+\dfrac{7}{2} \Delta_{\rm{QED }} & \Delta_{a \gamma} \\ 0&\Delta_{a\gamma} & \Delta_{aa} \end{pmatrix}, \end{aligned} \end{eqnarray} $
(3) with
$ \begin{eqnarray} \Delta_{\rm{pl}}&=&-\frac{\omega_{\rm{pl}}^2}{2E}\simeq -1.1 \times 10^{-4} \, {\rm{kpc}}^{-1} \, n_{{\rm{cm}}^{-3}}E_{\rm{GeV}}^{-1}, \end{eqnarray} $
(4) $ \begin{eqnarray} \Delta_{\rm{QED}}&=&\frac{\alpha E}{45\pi} \left( \frac{B_{\rm{T}}}{B_{\rm{cr}}}\right)^2 \simeq4.1\times 10^{-9}\, {\rm{kpc}}^{-1} \, E_{\rm{GeV}} B_{{\text{μ}} {\rm{G}}}, \ \ \end{eqnarray} $
(5) $ \begin{eqnarray} \Delta_{a\gamma}&=&\frac{1}{2}g_{a\gamma}B_{\rm{T}}\simeq 1.52 \times 10^{-2} \, {\rm{kpc}}^{-1} \, g_{11} B_{{\text{μ}} {\rm{G}}}, \end{eqnarray} $
(6) $ \begin{eqnarray} \Delta_{aa}&=&-\frac{m_a^2}{2E}\simeq -7.8 \times 10^{-2} \, {\rm{kpc}}^{-1} \, m_{\rm{neV}}^2 E_{\rm{GeV}}^{-1}, \end{eqnarray} $
(7) where
$ \omega_{\rm{pl}}= \sqrt{4\pi \alpha n_{\rm e}/m_{\rm e}} $ is the plasma frequency,$ n_{\rm e} $ is the number density of the free electrons, α is the fine-structure constant, and$ B_{\rm{cr}} \equiv m^2_{\rm e}/|e|\simeq4.4\times 10^{13} \rm \; G $ . The terms$ \Delta_{\rm{pl}} $ and$ \Delta_{\rm{QED}} $ represent the plasma and QED vacuum polarization effects, respectively. The notations$ n_{{\rm{cm}}^{-3}} \equiv n_{\rm e}/1 \rm {cm}^{-3} $ ,$ E_{\rm{GeV}} \equiv E /1 \; {\rm{GeV}} $ ,$ g_{11}\equiv g_{a\gamma}/10^{-11} \; {\rm{GeV}}^{-1} $ ,$ B_{\rm{{\text{μ}} G}} \equiv B_{\rm T} /1 {\text{μ}}{\rm{G}} $ , and$ m_{\rm{neV}} \equiv m_a/1 \; \rm {neV} $ are used in the above equations. The general mixing matrix$ {\cal{M}} $ depends on the angle ψ between the directions of$ B_{\rm{T}} $ and$ x_2 $ .ALP-photon conversion occurs in numerous regions with different magnetic field configurations. The final density matrix can be derived from the solution of Eq. (2) as
$ \begin{eqnarray} \rho\left(s\right)=T(s)\rho(0)T^\dagger(s). \end{eqnarray} $
(8) The entire transfer matrix
$ T(s) $ for the propagation distance s reads as$ \begin{eqnarray} T(s)=\prod^{n}_{i}{\cal{T}}_i, \end{eqnarray} $
(9) where
$ {\cal{T}}_i $ can be derived from the mixing matrix$ {\cal{M}}_i $ in the i-th region. For the initial unpolarized photon beam with$ \rho(0)=\rm {diag}(1, 1, 0)/2 $ , the photon survival probability after propagation is given by [21]$ \begin{eqnarray} P_{\gamma\gamma}={\rm{Tr}}\left(\left(\rho_{11}+\rho_{22}\right)T(s)\rho(0)T^\dagger(s)\right) \end{eqnarray} $
(10) with
$ \rho_{ii}={\rm{diag}}(\delta_{i1},\delta_{i2},0) $ .Subsequently, we describe the prorogation effect of the ALP-photon beam in three astrophysical regions with different magnetic field configurations, including the blazar jet, extragalactic space, and Milky Way [15, 28]. For the BL Lac objects considered in this study, we do not consider the effects in the blazar broad line region. The ALP-photon oscillation might occur significantly in the BJMF. There is evidence that the magnetic field of the BL Lac jet can be described by the poloidal (along the jet, with
$ B \propto r^{-2} $ ) and toroidal (perpendicular to the jet, with$ B \propto r^{-1} $ ) coherent components [55]. We take the BJMF model of the BL Lac sources as in Refs. [26, 35].The transverse magnetic field
$ B_{\rm{jet}}(r) $ reads [56, 57]$ \begin{eqnarray} B_{\rm{jet}}(r)=B_0\left(\frac{r}{r_{\rm{VHE}}}\right)^{-1}, \end{eqnarray} $
(11) where
$ r_{\rm{VHE}} $ is the distance between the central black hole and emission region. The density profile of the electrons$ n_{\rm{el}}(r) $ can be given by [58]$ \begin{eqnarray} n_{\rm{el}}(r)=n_0\left(\frac{r}{r_{\rm{VHE}}}\right)^{-2}. \end{eqnarray} $
(12) Note that the above profiles hold in the jet comoving frame. The energies of the photons in the laboratory frame
$ E_{\rm L} $ and comoving frame$ E_j $ are related by the Doppler factor$ \delta_{\rm{D}} $ through$ E_{\rm L} = E_j\cdot\delta_{\rm{D}} $ .The fit to the blazar spectra on the multi-wave bands with the synchrotron self-Compton model may determine the values of the BJMF parameters. In our analysis, these parameters for one source during all phases are assumed to be the same. We set
$ B_0 $ as$ 0.1 $ and$ 1.0\; {\rm{G}} $ for Mrk 421 and PG 1553+113, respectively, and take$ \delta_{\rm{D}} = 30 $ and$ n_0 = 3 \times 10^3 \; \rm {cm}^{-3} $ as the benchmark parameters. These values are consistent with the results derived in Refs. [53, 59]. In the region with$ r > 1\; \rm kpc $ , we assume that the magnitude of BJMF is zero. Note that among the BJMF parameters,$ r_{\rm{VHE}} $ is difficult to determine through the measurements. Its value might be in the range of$ {\cal{O}}(10^{16}) $ −${\cal{O}}(10^{17})\;{\rm{cm}}$ . Here, we adopt$ r_{\rm{VHE}}=10^{17}\;{\rm{cm}} $ as a benchmark parameter.When the ALP-photon system propagates in the host galaxy in which the blazar is located, the oscillation effect can be neglected [35, 60]. If the blazar is located in a cluster with a rich environment, the turbulent inter-cluster magnetic field
$\sim {\cal{O}}(1)$ μG may also induce a significant ALP-photon oscillation effect [28]. Because no definite evidence that the blazars Mrk 421 and PG 1553+113 are located in such an environment has been provided, this oscillation effect is not considered in our analysis.The oscillation in the extragalactic magnetic field on the largest cosmological scale is also neglected here. The magnitude of this magnetic field is not larger than
$ {\cal{O}}(1)\; \rm nG $ , although it has not yet been precisely determined [61]. For the VHE photons crossing extragalactic space, the attenuation effect caused by the extragalactic background light (EBL) through$ \gamma_{\rm{VHE}} + \gamma_{\rm{EBL}} \to e^+ + e^- $ should be considered. This effect is described by a suppression factor of$ {\rm e}^{-\tau} $ , where τ is the optical depth depending on the redshift of the source and the EBL density distribution. In this study, we take the EBL model provided by Ref. [62] as a benchmark. The redshift of Mrk 421 and PG 1553+113 are taken as$ z_0=0.031 $ and 0.45, respectively.Finally, we consider the effect in the magnetic field of the Milky Way, where the ALPs could be reconverted to photons. Only the regular component of the Galactic magnetic field is considered here, whereas the random component on the small scale is neglected. Details on this model can be found in Ref. [63].
We show the photon survival probability
$ P_{\gamma\gamma} $ as a function of energy for the blazars Mrk 421 and PG 1553+113 in Fig. 1. It can be seen that the pure EBL attenuation effect described by the factor of$ {\rm e}^{-\tau} $ dramatically suppresses the photon energy spectrum at energies above${\cal{O}}(10^2)\;{\rm GeV}$ . The ALP-photon oscillation might affect the survival probability at lower energies compared with the EBL attenuation effect. However, for several ALP parameters, the ALP-photon conversion can compensate for the EBL attenuation effect in the VHE region and lead to a moderate photon survival probability. This compensation may be significant for PG 1553+113 at large redshifts, as shown in Fig. 1.Figure 1. (color online) Photon survival probability as a function of energy for Mrk 421 (left) and PG 1553+113 (right). The black dotted dashed lines represent the survival probability with only the EBL attenuation effect. The solid lines represent the survival probability with both the EBL attenuation and ALP-photon oscillation effects for selected ALP parameters. The EBL model is taken from Ref. [62].
-
MAGIC [64, 65] is a system containing two imaging atmospheric Cherenkov telescopes located at the Roque de los Muchachos Observatory in Spain. These telescopes can detect extensive air showers in the stereoscopic mode and observe VHE γ-ray sources at energies above
$ 50\;{\rm GeV} $ [54]. In Ref. [54], the MAGIC collaboration reported 32 VHE γ-ray spectra from 12 blazars. All the data were collected during dark nights in good weather conditions. The γ-ray spectra at lower energies,$ \sim 0.1−100\;{\rm GeV} $ , during the common operation time observed by Fermi-LAT are also analyzed in Ref. [54]. Here, we use the MAGIC results of the BL Lac sources Mrk 421 and PG 1553+113 covering several activity phases to investigate the ALP-photon oscillation effects.We take expressions for the γ-ray blazar intrinsic energy spectra
$ \Phi_{\rm{int}} (E) $ as in Ref. [54].$ \Phi_{\rm{int}} (E) $ can be described using simple functions with three to five parameters, including the power law with exponential cut-off (EPWL), power law with superexponential cut-off (SEPWL), log parabola (LP), and log parabola with exponential cut-off (ELP). The functional expressions of$ \Phi_{\rm{int}} (E) $ are given as follows:● EPWL:
$ \begin{eqnarray} \Phi_{\rm{int}} ( E ) = F_0\left(\frac{E}{E_0}\right)^{-\Gamma}\exp\left(-\frac{E}{E_c}\right), \end{eqnarray} $
(13) ● SEPWL:
$ \begin{eqnarray} \Phi_{\rm{int}} ( E ) = F_0\left(\frac{E}{E_0}\right)^{-\Gamma}\exp\left(-\left(\frac{E}{E_c}\right)^d\right), \end{eqnarray} $
(14) ● LP:
$ \begin{eqnarray} \Phi_{\rm{int}} ( E ) = F_0\left(\frac{E}{E_0}\right)^{-\Gamma-b \log\left({E}/{E_0}\right)}, \end{eqnarray} $
(15) ● ELP:
$ \begin{eqnarray} \Phi_{\rm{int}} ( E ) = F_0\left(\frac{E}{E_0}\right)^{-\Gamma-b \log\left({E}/{E_0}\right)}\exp\left(-\frac{E}{E_c}\right), \end{eqnarray} $
(16) where
$ F_0 $ ,$ E_c $ , Γ, b, and d are free parameters. For EPWL and SEPWL,$ E_0 $ is taken to be 1 {$ {\rm GeV} $ }, whereas for LP and ELP,$ E_0 $ is also treated as a free parameter. Thus, these four functional expressions have 3, 4, 4, and 5 free parameters, respectively. For each phase, we choose the intrinsic energy spectrum with the minimum best-fit reduced$ \chi^2 $ under the null hypothesis. This is different from the analysis in Ref. [45], where the expression of the intrinsic energy spectrum is the same for all phases. The spectrum expressions for all the phases adopted in this analysis are listed in Table 1.Source [period] Tstart Tstop Spectrum $ {\chi}_{\rm{w/oALP}}^2 $ $ \chi^2_{\rm{min}} $ Effective d.o.f. $ \Delta\chi^2 $ Mrk 421 [20130410] 2013-04-09T12:00 2013-04-10T12:00 SEPWL 12.2 8.8 4.45 10.2 Mrk 421 [20130411] 2013-04-10T18:00 2013-04-11T06:00 ELP 16.2 10.1 6.85 13.9 Mrk 421 [20130412] 2013-04-11T18:00 2013-04-12T06:00 ELP 8.9 6.2 7.54 14.9 Mrk 421 [20130413a] 2013-04-12T12:00 2013-04-13T12:00 ELP 16.0 12.9 8.00 15.5 Mrk 421 [20130413b] 2013-04-12T12:00 2013-04-13T12:00 SEPWL 9.7 8.6 4.72 10.7 Mrk 421 [20130413c] 2013-04-12T12:00 2013-04-13T12:00 SEPWL 10.0 7.5 4.56 10.4 Mrk 421 [20130414] 2013-04-13T12:00 2013-04-14T12:00 ELP 22.4 13.7 9.22 17.2 Mrk 421 [20130415a] 2013-04-14T21:17 2013-04-15T04:13 ELP 5.8 4.8 5.23 11.4 Mrk 421 [20130415b] 2013-04-14T21:17 2013-04-15T04:13 SEPWL 13.4 10.0 5.02 11.1 Mrk 421 [20130415c] 2013-04-14T21:17 2013-04-15T04:13 SEPWL 5.1 4.0 4.71 10.6 Mrk 421 [20130416] 2013-04-15T12:00 2013-04-16T09:00 SEPWL 32.9 19.6 5.48 11.8 Mrk 421 [20130417] 2013-04-16T18:00 2013-04-17T06:00 SEPWL 26.1 11.2 4.99 11.1 Mrk 421 [20130418] 2013-04-17T12:00 2013-04-18T12:00 EPWL 13.3 9.0 5.56 12.0 Mrk 421 [20130419] 2013-04-18T12:00 2013-04-19T12:00 ELP 3.6 2.0 4.07 9.6 Mrk 421 [20140426] 2014-04-25T18:00 2014-04-26T06:00 ELP 25.8 15.2 6.19 12.9 PG 1553+113 [ST0202] 2012-02-28T12:00 2012-03-04T12:00 EPWL 2.3 0.9 3.52 8.7 PG 1553+113 [ST0203] 2012-03-13T12:00 2012-05-02T12:00 SEPWL 15.6 6.3 6.24 13.0 PG 1553+113 [ST0302] 2013-04-07T12:00 2013-06-12T12:00 SEPWL 5.4 1.3 4.50 10.3 PG 1553+113 [ST0303] 2014-03-11T12:00 2014-03-25T12:00 EPWL 10.2 5.9 6.73 13.7 PG 1553+113 [ST0306] 2015-01-25T12:00 2015-08-07T12:00 SEPWL 4.7 0.7 3.43 8.6 Combined Mrk 421 221.5 204.6 31.17 45.2 Combined PG 1553+113 38.2 20.5 8.94 16.9 Table 1. Best-fit values of
$ {\chi}_{\rm{w/oALP}}^2 $ under the null hypothesis and$ \chi^2_{\rm{min}} $ under the ALP hypothesis for all phases. The periods represent the corresponding MAGIC observations. The expressions for the intrinsic energy spectra, the effective$d.o.f.$ of the TS distributions, and$ \Delta{\chi}^2 $ at 95% C.L. are also listed. The last two rows denote the results of the combined analysis.Under the alternative hypothesis, including the ALP-photon oscillation effect, we obtain the expected photon spectrum as
$ \begin{eqnarray} \Phi_{\rm{w \; ALP}} ( E ) = P_{\gamma\gamma} \Phi_{\rm{int}} ( E ), \end{eqnarray} $
(17) where
$ P_{\gamma\gamma} $ is the photon survival probability. The detected photon flux in the energy bin of$ (E_1,E_2) $ is given by [35, 36]$ \begin{eqnarray} \Phi'=\frac{\int_0^{\infty}D(E', E_1, E_2)\Phi ( E' ){\rm{d}}E'}{E_2 - E_1}, \end{eqnarray} $
(18) where
$ D(E', E_1, E_2) $ is the energy dispersion function, and$ E' $ and$ \Phi(E') $ are the energy and spectrum of the photons before detection, respectively. The energy resolution of MAGIC is taken to be 16% [65].In Ref. [54], the Fermi-LAT spectra are provided in the form of spectral bow-ties rather spectral points. The bow-ties contain information on the flux and local spectrum index determined at the decorrelation energy; each one contributes two degrees of freedom in the fit. The
$ \chi^2 $ of the fit is defined as [54]$ \begin{aligned}[b] \chi^2 =& \left(\frac{\Phi'(E_{\rm{LAT}})-F_{\rm{LAT}}}{\Delta F_{\rm{LAT}}}\right)^2+ \left(\frac{\Gamma_{\rm{fit}}-\Gamma_{\rm{LAT}}}{\Delta \Gamma_{\rm{LAT}}}\right)^2 \\ &+ \sum_{i=1}^{N} \left(\frac{\Phi'(E_i) - \tilde{\phi}_i}{\delta_i}\right)^2, \end{aligned} $
(19) where
$ E_{\rm{LAT}} $ ,$ F_{\rm{LAT}} $ ,$ \Gamma_{\rm{LAT}} $ , and$ \Gamma_{\rm{fit}} $ are the central energy, flux, local spectral index, and expected spectral index for the Fermi-LAT results, respectively. N is the number of the MAGIC spectral points,$ \Phi'(E_i) $ is the expected flux of the photons,$ \tilde{\phi}_i $ is the detected photon flux, and$ \delta_i $ is the uncertainty of the MAGIC measurement.With the
$ \chi^2 $ values under the ALP hypothesis in the$ m_a-g_{a\gamma} $ plane, the constraint on the parameter space is set by requiring$ \chi^2 \leq \chi_{\rm{min}}^2 + \Delta{\chi}^2 $ , where$ {\chi}_{\rm{min}}^2 $ is the minimum best-fit$ {\chi}^2 $ under the ALP hypothesis. Because the modifications of the blazar spectra nonlinearly depend on the ALP parameters, the threshold value of$ \Delta{\chi}^2 $ at the particular confidence level should be derived from the Monte Carlo simulations rather than directly using Wilks' theorem [20, 34]. Based on the best-fit spectra to the data under the null hypothesis, for each phase, 400 sets of spectra in the pseudo-experiments are generated by Gaussian samplings. The test statistic (TS) value is defined as the difference between the best-fit$ {\hat{\chi}}^2 $ under the null and ALP hypotheses for each generated spectrum$ {TS}\equiv {\hat{\chi}_{\rm{null}}}^2 - {\hat{\chi}_{\rm{w \; ALP}}}^2 $ . In each phase, the distribution of the TS for all generated spectrum sets is derived. Such distribution can be described by the non-central$ \chi^2 $ distribution (see Appendix 6) with the non-centrality λ and effective degree of freedom ($ d.o.f. $ ). Although this TS distribution is derived under the null hypothesis, following Ref. [34], we take it as the approximation of the TS distribution under the ALP hypothesis and adopt the corresponding$ \Delta{\chi}^2 $ in the subsequent analysis. -
In this appendix, we briefly introduce the non-central
$ \chi^2 $ distribution. Usually, the distribution of the$ \chi^2 $ density function with the$ d.o.f. \nu $ can be described by$ \begin{eqnarray} f(x,\nu)=\frac{1}{2^{{\nu}/{2}}\Gamma({\nu}/{2})}x^{{\nu}/{2}-1}{\rm e}^{-{x}/{2}}, \end{eqnarray} \tag{A1}$
where
$ \Gamma(a) $ is the celebrated Gamma function$ \begin{eqnarray} \Gamma(a)=\int_0^{\infty}t^{a-1}{\rm e}^{-t}{\rm{d}}t. \end{eqnarray} \tag{A2}$
Then, we have the non-central
$ \chi^2 $ distribution described by the Poisson mixture of the$ \chi^2 $ density function$ \begin{eqnarray} f(x,\nu,\lambda)=\sum_{j=0}^\infty \frac{{\rm e}^{-{\lambda}/{2}}\left({\dfrac{\lambda}{2}}\right)^j}{j!}f(x,\nu+2j), \end{eqnarray} \tag{A3}$
with the non-centrality λ and
$d.o.f. \nu$ . Both the values of λ and ν can be derived by fitting the realistic TS distribution. An integer value is not required for the derived ν, which is often referred to as the effective$d.o.f$ . The cumulative distribution function is given by$ \begin{eqnarray} P(x,\nu,\lambda)=\int_0^x f(t,\nu,\lambda){\rm{d}}t. \end{eqnarray}\tag{A4} $
When considering a
$ 95\% $ limit for a given ν and λ, one can solve$ P(x_{95\%}, \nu, \lambda)=0.95 $ and take the corresponding$ x_{95\%} $ as the threshold value of$ \Delta \chi^2 $ at$ 95\% $ C.L.
Searching for axion-like particles with the blazar observations of MAGIC and Fermi-LAT
- Received Date: 2022-04-08
- Available Online: 2022-08-15
Abstract: In this study, we explore the axion-like particle (ALP)-photon oscillation effect in the γ-ray spectra of the blazars Markarian 421 (Mrk 421) and PG 1553+113, which are measured by the Major Atmospheric Gamma Imaging Cherenkov Telescopes (MAGIC) and Fermi Large Area Telescope (Fermi-LAT) with high precision. The Mrk 421 and PG 1553+113 observations of 15 and five phases are used in the analysis, respectively. We find that the combined analysis with all the 15 phases improves the limits of the Mrk 421 observations. For the selected blazar jet magnetic field and extragalactic background light models, the combined limit set by the Mrk 421 observations excludes the ALP parameter region with the ALP-photon coupling of