Keywords: gamma-ray astronomy, blazars, intergalactic magnetic fields, extragalactic background light, axion-like particles 1. We

Gamma-rays of energies above 30 GeV can not freely propagate through the Universe due to their interaction with the Extragalactic Background Light (EBL). In these interactions gamma-rays produce secondary electron-positron pairs, which in turn upscatter Cosmic Microwave Background (CMB) photons. Resulting electromagnetic cascade changes the intrinsic spectrum of the source both at high and low energies. At high energies the flux is reduced due to absorption of primary photons. At lower energies there is an additional flux from the electromagnetic cascade. Properties of electromagnetic cascades can be used to study intergalactic magnetic fields (IGMF) through observations of the time delay of the secondary photons, signatures in the spectra and extended emission around the point sources. The same data can be used to study the imprints of the hypothetical axion-like particles (ALP), one of the main Dark Matter candidates, in the spectra of gamma-ray sources.

The future gamma-ray experiment Cherenkov Telescope Array (CTA) plans to measure the fluxes of blazars in the energy range from 30 GeV to 30 TeV with sensitivity 10 times superior to present measurements. The goal of this thesis is to model the above phenomena, taking into account the recent models of extragalactic background light and make model predictions for the future CTA measurements. First, we build a new, flexible EBL model and examine its dependence on the main astrophysical parameters underlying it. Using available gamma-ray data we set constraints on the parameters of the additional narrow spectral feature in the EBL spectrum and then convert them into constraints on the ALP coupling constant with photons. On the other hand, we estimate sensitivity of the CTA for strong cosmological magnetic fields and show that it can be used to probe magnetic field with a strength up to 10 -12 G -10 -11 G. Finally we discuss two classes of systematic errors, which must be taken into account when studying IGMF with gamma-ray technique. The first class is associated with the baryonic feedback effect, and the second -with inaccuracies of numerical simulation.

Résumé

Les rayons gamma d'énergies supérieures à 30 GeV ne peuvent pas voyager librement à travers l'Univers en raison de leur interaction avec la lumière d'arrière-plan extragalactique (EBL). Dans ces interactions, les rayons gamma produisent des paires électronpositon secondaires, qui à leur tour diffusent vers le haut les photons du fond diffus cosmologique (CMB). La cascade électromagnétique qui en résulte modifie le spectre intrinsèque de la source à la fois aux hautes et basses énergies. Aux hautes énergies, le flux est réduit en raison de l'absorption des photons primaires. Aux énergies inférieures, il y a un flux supplémentaire provenant de la cascade électromagnétique. Les propriétés des cascades électromagnétiques peuvent être utilisées pour étudier les champs magnétiques intergalactiques (IGMF) à travers des observations du retard temporel des photons secondaires, des signatures dans les spectres et de l'émission étendue autour des sources ponctuelles. Les mêmes données peuvent être utilisées pour étudier les empreintes des particules hypothétiques de type axion (ALP), l'un des principaux candidats à la matière noire, dans le spectre des sources de rayons gamma.

La future expérience de rayons gamma Cherenkov Telescope Array (CTA) prévoit de mesurer les flux de blazars dans la gamme d'énergie de 30 GeV à 30 TeV avec une sensibilité 10 fois supérieure aux mesures actuelles. L'objectif de cette thèse est de modéliser les phénomènes ci-dessus, en tenant compte des modèles récents de lumière de fond extragalactique et de faire des prédictions de modèles pour les futures mesures CTA. Tout d'abord, nous construisons un nouveau modèle EBL flexible et examinons sa dépendance vis-à-vis des principaux paramètres astrophysiques qui le sous-tendent. En utilisant les données de rayons gamma disponibles, nous avons défini des contraintes sur les paramètres de la caractéristique spectrale étroite supplémentaire dans le spectre EBL, puis nous les avons converties en contraintes sur la constante de couplage ALP avec les photons. D'autre part, nous estimons la sensibilité du CTA pour les champs magnétiques cosmologiques forts et montrons qu'il peut être utilisé pour sonder un champ magnétique d'une force allant jusqu'à 10 -12 G -10 -11 G. Enfin, nous discutons de deux classes d'erreurs systématiques, qui doivent être prises en compte lors de l'étude de l'IGMF avec la technique des rayons gamma. La première classe est associée à l'effet de rétroaction baryonique et la seconde aux imprécisions de la simulation numérique.

Mots clefs : astronomie des rayons gamma, blazars, champs magnétiques intergalactiques, lumière d'arrière-plan extragalactique, particules de type axion

Résumé de la thèse

Des progrès significatifs ont été réalisés dans le domaine de l'astrophysique des hautes énergies au cours des dernières décennies. Néanmoins, de nombreuses questions liées à la production, la propagation et la détection des rayons gamma de haute énergie restent en suspens. Ce sujet est très vaste et constitue une branche à part entière de l'astrophysique, par conséquent, dans cette thèse, nous nous concentrerons principalement sur les problèmes de propagation des rayons gamma à travers le milieu intergalactique, ainsi que sur l'interprétation de données observées.

Les rayons gamma de très haute énergie voyagent directement de leurs sources astrophysiques à l'observateur sur Terre, à moins qu'ils ne soient absorbés par l'EBL. Dans ce cas, le photon original de haute énergie interagit avec un rayonnement de fond optique ou infrarouge de faible énergie (lumière d'arrière-plan extragalactique, EBL) et crée une paire électron-positon. Le niveau d'absorption dépend de l'énergie du photon, de la distance à la source et de l'intensité EBL. Le principal problème est que la mesure directe de l'EBL est une tâche expérimentale compliquée en raison des contributions dominantes de la Galaxie et du Système solaire. De plus, ces mesures ne fournissent un spectre qu'au redshift z = 0. Le spectre EBL z > 0, qui encode des informations sur l'évolution de la matière lumineuse dans l'Univers, ne peut être étudié qu'indirectement en observant l'atténuation des rayons γ-ray à haute énergie. Un certain nombre de modèles EBL ont été construits qui prennent en compte contraintes astrophysiques sur l'évolution de l'EBL, y compris les limites sur le taux de formation d'étoiles en fonction du décalage vers le rouge. À partir de [START_REF] Franceschini | The extragalactic optical-infrared background radiations, their time evolution and the cosmic photon-photon opacity[END_REF] et [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF], tous les modèles ultérieurs prédisent approximativement le même spectre EBL, l'intensité dans la plupart de ces modèles ne différant pas de plus d'un facteur deux. Les limites inférieures de l'EBL obtenues à partir de la somme des contributions des galaxies individuelles prédisent une intensité EBL cohérente avec les résultats de la modélisation. Dans le même temps, les observations directes d'EBL conduisent à une intensité EBL beaucoup plus forte, ce qui contredit les observations de rayons gamma à haute énergie. Avec une EBL aussi forte, les observations de sources extragalactiques dans la gamme du TeV seraient impossibles sans invoquer de nouveaux effets physiques.

De plus, il y a des indications que l'absorption des rayons gamma est encore plus faible que prévu par les modèles EBL les plus conservateurs. Ils se manifestent dans l'observation de photons avec des énergies dans la région du TeV provenant de sources distantes, qui, comme prévu, auraient dû être absorbés conformément aux modèles EBL modernes. Le même problème existe aux énergies de 100 GeV pour des sources plus éloignées mesurées par le télescope Fermi LAT. En particulier, [START_REF] Neronov | Very high-energy γ-ray emission from high-redshift blazars[END_REF] a rapporté la détection de rayonnement gamma à haute énergie de 13 blazars avec des décalages vers le rouge z > 0.5. Plusieurs de ces blazars ont été découverts plus tard par les télescopes Cherenkov. Une étude récente [START_REF] Horns | Indications for a pair-production anomaly from the propagation of VHE gamma-rays[END_REF] et l'étude ultérieure [START_REF] Rubtsov | Breaks in gamma-ray spectra of distant blazars and transparency of the Universe[END_REF] indiquent des ruptures vers le haut non physiques dans les spectres blazar corrigés en absorption ("transparence anormale" de l'Univers). Une analyse détaillée a montré que la position de la rupture correspond à l'énergie pour laquelle l'absorption devient significative, et l'amplitude de la rupture augmente pour les sources plus éloignées. Ces observations peuvent indiquer l'existence d'une nouvelle particule de type axion (APL) qui peut se mélanger aux photons. Les ALP sont mélangés avec des photons dans des champs magnétiques externes, ce qui peut réduire efficacement la quantité d'absorption, liée à la naissance de paires. En effet, l'ALP ne génère pas de paires, c'est donc un scénario possible, lorsque les photons gamma sont convertis en ALP, qui se propage ensuite un certain temps sans absorption, et après qu'ils sont reconvertis en photons. Bien que le faisceau de photons soit absorbé, le faisceau de l'ALP se propage sans absorption et l'extinction générale devient moins forte.

Il est d'un intérêt considérable de vérifier ces effets, compte tenu de l'émergence d'une grande quantité de nouvelles données sur les observations de blazar. Il est également important de savoir s'il est possible de concilier les observations de blazars à haute énergie et les observations directes de l'EBL. D'autre part, l'observation de sources TeV de rayonnement gamma peut être utilisée pour étudier les champs magnétiques intergalactiques (IGMF). La méthode de détermination des propriétés IGMF à partir des caractéristiques des spectres blazar a été proposée dans [START_REF] Neronov | A method of measurement of extragalactic magnetic fields by TeV gamma ray telescopes[END_REF] puis développée dans une série d'articles [START_REF] Elyiv | Gamma-ray induced cascades and magnetic fields in the intergalactic medium[END_REF][START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF][START_REF] Neronov | Degree-scale GeV "Jets" from Active and Dead TeV Blazars[END_REF]. Cette méthode exploite le fait que les électrons secondaires et les positrons, qui sont produits par absorption des rayons gamma, peuvent être déviés localement dans l'espace intergalactique avant de produire des rayons gamma secondaires par diffusion Compton inverse (ICS) sur les photons CMB. Dans ce cas, la structure spectrale, temporelle et angulaire des flux de rayons gamma primaires et secondaires permet de déterminer les propriétés de l'IGMF. L'application de cette méthode à des données réelles a permis de détecter des champs magnétiques intergalactiques [START_REF] Neronov | Evidence for Strong Extragalactic Magnetic Fields from Fermi Observations of TeV Blazars[END_REF] et de limiter leurs paramètres. Néanmoins, la sensibilité limitée de la génération moderne de télescopes ne permet pas de déterminer les valeurs exactes de l'intensité et des autres paramètres du champ magnétique dans l'espace intergalactique.

Les contraintes sur les modèles IGMF et EBL sont fortement liées, comme l'a montré l'ensemble de données combiné Fermi-HESS de blazar 1ES 0229+200 en tenant compte des composantes directes et en cascade de la source spectrum [START_REF] Ie | Fermi/LAT Observations of 1ES 0229+200: Implications for Extragalactic Magnetic Fields and Background Light[END_REF]. Par exemple, l'augmentation de l'intensité EBL conduit à une absorption plus forte et améliore le flux de rayonnement secondaire, qui à son tour affecte les limites de l'IGMF.

En résumé, la mesure de l'EBL, la manifestation possible des ALP et la détermination des paramètres de l'IGMF forment un ensemble de problèmes interdépendants qui nécessitent une solution commune. La pertinence des problèmes décrits augmente en relation avec la mise en service prochaine d'un nouveau télescope à rayons gamma Cherenkov Telescope Array (CTA) qui a une sensibilité beaucoup plus grande par rapport aux télescopes de la génération actuelle. Avec une résolution angulaire améliorée, CTA sera en mesure de faire la lumière sur bon nombre de ces problèmes. La thèse se compose de quatre chapitres consacrés à l'étude de l'EBL, des ALP, de la modélisation des cascades électromagnétiques et de l'IGMF.

Dans la première partie de la thèse, nous avons construit un nouveau modèle EBL flexible [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF]. L'avantage de ce modèle est qu'il permet, contrairement à d'autres modèles, de modifier les paramètres qui le sous-tendent. En utilisant ce modèle, nous avons identifié les paramètres astrophysiques les plus importants de l'EBL et les avons écrits explicitement. En utilisant les chaînes de Monte Carlo Markov et sur la base de ce nouveau modèle EBL, nous avons développé une méthode de contrainte des paramètres astrophysiques basée sur les mesures EBL disponibles. Pour certains paramètres, les limites obtenues avec les mesures EBL modernes sont comparables aux limites obtenues par d'autres méthodes. Néanmoins, l'application de cette approche à des mesures futures plus précises de l'EBL peut augmenter considérablement la sensibilité de la méthode.

D'autre part, nous avons comparé les mesures EBL visibles et infrarouges avec les contraintes des rayons gamma et étudié si les écarts des mesures avec différentes méthodes pourraient être dus à la présence de caractéristiques dans le spectre EBL qui sont localisées dans la gamme de longueurs d'onde du micron. Nous avons montré qu'un " excès " précédemment signalé dans le flux EBL dans la gamme de longueurs d'onde de 1 micron est cohérent avec les mesures de rayons gamma, c'est-à-dire si l'excès a la forme d'une caractéristique étroite de la largeur δλ < λ et un flux global allant jusqu'à 15 nW/m 2 /sr au-dessus de l'EBL "minimal", qui est estimé à partir des comptages de sources visibles et infrarouges. De telles caractéristiques spectrales en forme de bosse pourraient provenir, par exemple, de particules de matière noire en décomposition, ou d'axions ou de processus astrophysiques particuliers au cours de l'histoire de la formation d'étoiles. Nous avons montré que les futures mesures avec CTA peuvent répondre à la question de la présence d'une telle caractéristique dans le spectre EBL [START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF].

Dans le chapitre suivant, nous avons examiné l'influence des particules de type axion sur les spectres de haute énergie des blazars. Nous avons trouvé 95% de CL. limite de la masse des particules de type axion et constante de couplage résultant de la non-observation de la bosse de désintégration des axions dans le spectre EBL, en utilisant le spectre de rayons gamma des observations du blazar 1ES 1218+304 par les télescopes Fermi/LAT et VERITAS [START_REF] Korochkin | Search for decaying eV-mass axionlike particles using gamma-ray signal from blazars[END_REF] . Notre modèle de meilleur ajustement favorise une amplitude non nulle de la caractéristique d'absorption induite par l'ALP dans le spectre des rayons γ. Nous avons discuté des possibilités disponibles pour examiner ce résultat avec des observations complémentaires avec les télescopes à rayons γ au sol actuels et de prochaine génération, CTA et LHAASO [START_REF] Addazi | The Large High Altitude Air Shower Observatory (LHAASO) Science Book (2021 Edition)[END_REF]. Ils fourniront une amélioration d'un facteur de trois à cinq de la sensibilité pour les recherches de caractéristiques de désintégration ALP. Si le spectre EBL présente une bosse induite par l'ALP avec le flux au niveau suggéré par les expériences de mesure directe de l'EBL, le CTA sera en mesure de détecter la bosse induite par l'ALP avec une grande importance et de fournir une mesure précise de la masse ALP et de la constante de couplage des photons. . Nous avons également testé une indication possible de l'existence des ALP à partir d'observations blazar, qui est l'apparition de cassures vers le haut dans les spectres désabsorbés. Sur la base de l'échantillon le plus complet de blazars avec des décalages vers le rouge mesurés avec précision, nous avons montré que la signification de l'effet de "transparence anormale" est plus faible que précédemment rapporté et ne dépasse pas 2σ pour les modèles EBL modernes et même comme bas à 1.3σ pour le modèle EBL le plus conservateur [START_REF] Korochkin | Search for anomalous features in gamma-ray blazar spectra corrected for the absorption on the extragalactic background light[END_REF].

Ensuite, nous avons effectué une comparaison détaillée des modules électromagnétiques des programmes de Monte Carlo accessibles au public CRbeam, CRPropa et ELMAG [START_REF] Kalashev | Modelling of propagation of very-high-energy gamma rays with CRbeam code. Comparison with CRPropa and ELMAG codes[END_REF]. Nous comparons la propagation des électrons et des positrons dans un champ magnétique, la génération de champs turbulents, l'absorption des rayons gamma de haute énergie due à l'interaction avec les photons EBL, ainsi que le processus de diffusion Compton inverse des électrons sur les photons de fond. Nous avons également comparé la distribution d'énergie des particules secondaires. Lors de la comparaison des résultats avec les prévisions théoriques, des inexactitudes ont été trouvées, dont certaines ont été corrigées. Après l'introduction des corrections, tous les codes montrent un accord au niveau de 10%, ce qui peut être considéré comme une estimation de l'erreur de modélisation systématique.

Enfin, nous avons considéré les possibilités de l'astronomie gamma pour étudier l'IGMF. Un champ magnétique primordial avec une longueur de corrélation de 0.1 -1 kpc la force dans la plage 10 -12 G -10 -11 G peut résoudre la tension entre différentes mesures de la constante de Hubble et fournir une explication de l'excès d'opacité dans la ligne de 21 cm au décalage vers le rouge 15 < z < 20, s'il est présent pendant les époques de recombinaison et de réionisation. Nous étudions la portée de sensibilité de la technique des rayons gamma pour la mesure d'un tel champ magnétique cosmologique relativement fort en utilisant des expositions profondes des blazars à spectre dur les plus proches avec des télescopes CTA. Nous montrons que la méthode de mesure des rayons gamma peut détecter le champ magnétique cosmologique avec une force allant jusqu'à 10 -11 G [START_REF] Korochkin | Sensitivity reach of gamma-ray measurements for strong cosmological magnetic fields[END_REF][START_REF] Korochkin | Sensitivity of the gamma-ray method for strong primordial magnetic fields[END_REF]. La combinaison du fond cosmique micro-ondes et des contraintes des rayons gamma peut ainsi détecter toute la gamme des champs magnétiques cosmologiques possibles pour confirmer ou exclure leur pertinence pour le problème de l'origine des champs magnétiques cosmiques, ainsi que leur influence sur les époques de recombinaison et de réionisation. . D'autre part nous avons étudié une possibilité de détection de champ magnétique primordial à grande longueur de corrélation (> 100 Mpc) pouvant provenir de l'inflation. Un tel champ peut être cohérent à l'échelle cosmologique et induire une émission étendue en forme de coin autour des blazars proches, avec des orientations de coin alignées dans plusieurs sources à travers le ciel. Cet alignement peut être utilisé pour distinguer l'origine inflationniste de l'IGMF, car il n'est pas attendu si l'IGMF provient de transitions de phase cosmologiques [START_REF] Korochkin | Detectability of large correlation length inflationary magnetic field with Cherenkov telescopes[END_REF].

Nous avons également exploré l'effet de la rétroaction baryonique sur la LSS sous la forme d'éjections de galaxies qui polluent le milieu intergalactique avec des champs magnétiques. Ce processus peut propager le champ magnétique dans les vides du LSS et empêcher ainsi la mesure des champs magnétiques primordiaux qui résident éventuellement dans les vides. Notre analyse du modèle de rétroaction baryonique spécifique mis en oeuvre dans les simulations cosmologiques IllustrisTNG montre que la rétroaction n'affecte généralement pas fortement les champs magnétiques vides. Les sorties galactiques des bulles magnétisées ont collectivement un petit facteur de remplissage volumique. Nous avons exploré l'effet des bulles magnétisées sur le signal γ-ray de l'AGN distant [START_REF] Bondarenko | Account of baryonic feedback effect in the γ-ray measurements of intergalactic magnetic fields[END_REF]. L'effet de ces bulles peut être pris en compte dans l'analyse comme une incertitude systématique du flux γ-ray secondaire au niveau de 10% en moyenne (et jusqu'à 60% pour une seule source). Cette erreur doit être prise en compte dans le budget d'erreur de l'analyse. Cette incertitude systématique dépend de la coupure haute énergie du spectre de la source primaire. Si l'énergie de coupure se situe dans la plage E cut = 100 TeV, une analyse dédiée des propriétés possibles de la bulle hôte AGN doit être effectuée dans le cadre de l'analyse des données liées à l'IGMF.

Chapter 1 Introduction

Significant progress has been made in the field of high-energy astrophysics in recent decades. Nevertheless, many questions related to the production, propagation, and detection of high-energy gamma rays remain unresolved. This topic is very extensive and constitutes a separate branch of astrophysics, therefore, in this thesis, we will mainly focus on the problems of propagation of gamma rays through the intergalactic medium, as well as on the interpretation of observed data.

The very-high-energy gamma rays travel from their astrophysical sources to the observer at Earth directly, unless absorbed by the EBL. In this case the original highenergy photon interacts with a low-energy optical or infrared background radiation (Extragalactic Background Light, EBL) and creates an electron-positron pair. The level of absorption depends on the photon energy, the distance to the source, and the EBL intensity. The main problem is that the direct EBL measurement is a complicated experimental task due to the dominant contributions from the Galaxy and the Solar System. Moreover, these measurements provide a spectrum only at redshift z = 0. High z > 0 EBL spectrum, which encodes information about the evolution of luminous matter in the Universe, can be studied only indirectly by observing high-energy γ-ray ray attenuation. A number of EBL models were constructed that take into account astrophysical constraints on the EBL evolution, including limits on the star formation rate depending on the redshift. Starting from [START_REF] Franceschini | The extragalactic optical-infrared background radiations, their time evolution and the cosmic photon-photon opacity[END_REF] and [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF], all subsequent models predict approximately the same EBL spectrum with the intensity in most of these models differs from each other by no more than a factor of two. Lower bounds on the EBL obtained from summing contributions of individual galaxies predict EBL intensity consistent with the results of modelling. At the same time, direct observations of EBL lead to a much stronger EBL intensity, which contradicts the observations of high-energy gamma rays. With such a strong EBL, observations of extragalactic sources in the TeV range would be impossible without invoking new physical effects.

Moreover, there are indications that the absorption of gamma rays is even lower than predicted by the most conservative EBL models. They manifests themselves in 8 the observation of photons with energies in the TeV region from distant sources, which, as expected, should have been absorbed in accordance with modern EBL models. The same problem exists at energies of 100 GeV for more distant sources measured by the Fermi LAT telescope. In particular, [START_REF] Neronov | Very high-energy γ-ray emission from high-redshift blazars[END_REF] reported the detection of high-energy gamma radiation from 13 blazars with redshifts z > 0.5. Several of these blazars were later discovered by Cherenkov telescopes. A recent [START_REF] Horns | Indications for a pair-production anomaly from the propagation of VHE gamma-rays[END_REF] study and the subsequent study [START_REF] Rubtsov | Breaks in gamma-ray spectra of distant blazars and transparency of the Universe[END_REF] point to non-physical upward breaks in absorption-corrected blazar spectra ("anomalous transparency" of the Universe). A detailed analysis showed that the position of the break corresponds to the energy for which the absorption becomes significant, and the magnitude of the break increases for more distant sources. These observations may indicate the existence of a new axion-like particle (APL) that can mix with photons. ALPs are mixed with photons in external magnetic fields, which can effectively reduce the amount of absorption, connected with the birth of pairs. Indeed, ALP does not generate pairs, so it is possible scenario, when gamma-ray photons are converted into ALP, which then some time propagates without absorption, and after they are converted back to photons. Although the photon beam is absorbed, the beam from the ALP propagates without absorption, and general extinction becomes not so strong.

It is of considerable interest to verify these effects, taking into account the emergence of a large amount of new data on blazar observations. It is also important to find out whether it is possible to reconcile observations of high-energy blazars and direct EBL observations.

On the other hand, the observation of TeV sources of gamma radiation can be used to study intergalactic magnetic fields (IGMF). The method for determining the IGMF properties from the features of the blazar spectra was proposed in [START_REF] Neronov | A method of measurement of extragalactic magnetic fields by TeV gamma ray telescopes[END_REF] and then developed in a series of papers [START_REF] Elyiv | Gamma-ray induced cascades and magnetic fields in the intergalactic medium[END_REF][START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF][START_REF] Neronov | Degree-scale GeV "Jets" from Active and Dead TeV Blazars[END_REF]. This method exploits the fact that secondary electrons and positrons, which are produced by gamma ray absorption, can be deflected locally in intergalactic space before they produce secondary gamma rays by inverse Compton scattering (ICS) on CMB photons. In this case, the spectral, timing and angular structure of the of primary and secondary gamma-ray flux allows to determine the properties of the IGMF. The application of this method to real data made it possible to detect intergalactic magnetic fields [START_REF] Neronov | Evidence for Strong Extragalactic Magnetic Fields from Fermi Observations of TeV Blazars[END_REF] and limit their parameters. Nevertheless, the limited sensitivity of the modern generation of telescopes does not allow one to determine the exact values of the strength and other parameters of the magnetic field in intergalactic space.

Constraints on IGMF and EBL models are strongly connected, as was shown from the combined Fermi-HESS data set of blazar 1ES 0229+200 taking into account the direct and cascade components of the source spectrum [START_REF] Ie | Fermi/LAT Observations of 1ES 0229+200: Implications for Extragalactic Magnetic Fields and Background Light[END_REF]. For example, increase of the EBL intensity leads to stronger absorption and enhances the secondary radiation flux, which in turn affects the limits on the IGMF.

Summarizing the above, measurement of the EBL, possible manifestation of ALPs and determination of the parameters of the IGMF form a set of interrelated problems that require a joint solution. The relevance of the described issues increases in connection with the commissioning in the near future of a new gamma-ray telescope Cherenkov Telescope Array (CTA) which has a much greater sensitivity compared to the current generation telescopes. Together with improved angular resolution, CTA will be able to shed light on many of these problems. Thesis consists of four Chapters which are devoted to the study of EBL, ALPs, electromagnetic cascade modelling and IGMF.

Chapter One describes the details and structure of the new EBL model. Unlike other models, all astrophysical parameters underlying our new model are explicitly written out. In the 2.3 section we explore the resulting parameter space. Assuming that direct EBL measurements are upper bounds, we use the new EBL model to constrain astrophysical parameters for the first time based on EBL observations. Further in section 2.4, we focus on the contradiction between direct measurements of the EBL and observations of high-energy γ-rays from blazar. We show that a possible solution to the problem is the presence of an additional narrow spectral feature in the EBL spectrum whose intensity is at the level of direct measurements. Based on blazar observations, we constrain the width, position, and intensity of this additional component. The results of this Chapter are published in [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF][START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF] Chapter Two is devoted to search for manifestations of axion-like particles in the spectra of blazars. First, in section 3.2, we show that the EBL spectral feature explored in the previous Chapter can be caused by the decay of axion-like particles constituting dark matter. In such a case, axion-like particles can be studied with gamma-ray observations. Using available data on high energy observations of blazars, we constrain ALP coupling constant to photons and estimate how these constraints can be improved by the next generation observatories CTA and LHAASO. In section 3.3, we investigate the problem of the "anomalous transparency" of the Universe for high-energy photons. Based on the most complete sample of blazars with accurately measured redshifts, we have shown that the significance of the "anomalous transparency" effect is weaker than previously thought and does not exceed 2σ for the modern EBL models. The main results of this Chapter are published in [START_REF] Korochkin | Search for anomalous features in gamma-ray blazar spectra corrected for the absorption on the extragalactic background light[END_REF][START_REF] Korochkin | Search for decaying eV-mass axionlike particles using gamma-ray signal from blazars[END_REF][START_REF] Addazi | The Large High Altitude Air Shower Observatory (LHAASO) Science Book (2021 Edition)[END_REF].

Chapter Three provides a detailed comparison of the electromagnetic modules of the publicly available Monte Carlo programs CRbeam, CRPropa and ELMAG. In particular we compare the propagation of electrons and positrons in a magnetic field, the absorption of high-energy gamma rays due to interaction with EBL photons, as well as the process of inverse Compton scattering of electrons on background photons. It was found that there are inaccuracies in the codes when modeling the absorption of gamma rays, as well as in the energy distribution of secondary particles. After eliminating the found inaccuracies, all codes demonstrate agreement at the level of 10%. The results of the comparison are available in the electronic preprint [START_REF] Kalashev | Modelling of propagation of very-high-energy gamma rays with CRbeam code. Comparison with CRPropa and ELMAG codes[END_REF].

The fourth Chapter is devoted to the study of intergalactic magnetic fields by the method of gamma-ray astronomy. Recently, several independent studies indicate the existence of IGMF with strength in the range of 10 -12 G -10 -11 G. In the sections 5.2 and 5.3 we explore the detectability of such a strong magnetic fields with γ-rayobservations. We show that next generation telescope CTA is sensitive to strong extragalactic magnetic fields with the strength in the range of 10 -12 G -10 -11 G. Also we discuss a possible systematic error in determining the IGMF parameters associated with the suppression of the cascade signal in magnetic field bubbles around galaxies and clusters of galaxies (section 5.4). The results of this Chapter are published in [START_REF] Korochkin | Sensitivity reach of gamma-ray measurements for strong cosmological magnetic fields[END_REF][START_REF] Korochkin | Detectability of large correlation length inflationary magnetic field with Cherenkov telescopes[END_REF][START_REF] Korochkin | Sensitivity of the gamma-ray method for strong primordial magnetic fields[END_REF][START_REF] Bondarenko | Account of baryonic feedback effect in the γ-ray measurements of intergalactic magnetic fields[END_REF].

Chapter 2

Extragalactic background light

Extragalactic background light (EBL) is a radiation emitted by stars and cosmic dust throughout the whole lifetime of the Universe at ultraviolet, optical and infrared wavelengths. In ultraviolet and optical spectral region it is dominated by the star radiation while the infrared spectrum is determined by emission of the dust. The understanding of the EBL is an attractive task for several reasons. First, the origin of EBL is directly connected to the star formation history and may shed a light on the details of the evolution of the Universe. Second, EBL is a key element of the attenuation of the high energy gamma-rays by the pair production on background photons [START_REF] Nikishov | Absorption of high-energy photons in the Universe[END_REF]. However, the exact spectrum of the EBL, as well as its dependence on redshift, is still unknown.

Direct observation of the EBL is a difficult experimental problem due to the dominant contribution of the Zodiacal light and the radiation of the Galaxy in the same spectral range. The currently known methods for measuring EBL predict significantly different values of the radiation intensity, in some cases differing by an order of magnitude. On the other hand, the too high intensity of the EBL contradicts the observations of high-energy gamma rays coming from cosmological distances.

In this Chapter, we will build a new EBL model with which we will show what astrophysical information can be extracted based on measurements of the EBL spectrum. Also, we will consider the issue of reconciliation of direct EBL measurements with observations of high-energy blazar emission. The results of this Chapter are published in [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF][START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF].

Measurements of the EBL

Extragalactic background radiation is radiation in the ultraviolet, optical and infrared ranges that came from outside our Galaxy Milky Way. Thus, EBL is one of the types of background radiation that fills the Universe, along with radio, microwave, ultraviolet, X-ray and gamma backgrounds (see Fig. 2.1). The wavelength range corresponding to the EBL spans about four orders of magnitude, ranging from about 0.1 µm in the ultraviolet to 100 µm in the far infrared region.

In the region of short wavelengths, the EBL borders on the ultraviolet background, and in the region of long wavelengths, it overlaps with the CMB. In contrast to the CMB, whose spectrum has been studied in great detail, the EBL spectrum has been studied much worse. Modern observations show that the EBL spectrum consists of two broad maxima. The first maximum is in the ∼ 1 -2µm region, and the second is near ∼ 100 -200µm. It is clear that the total radiation of all stars and galaxies contributes to the first maximum, and the second maximum corresponds to the emission of dust that reradiates the light of stars. Thus, the EBL contains information about the history of the formation of structures in the Universe, in particular, about the rate of star formation, the structure of the stellar population, and the parameters of cosmic dust. However, the accuracy of modern measurements does not allow answering the question of whether the contribution of stars and dust to the EBL is the only one, or whether the spectrum contains additional components corresponding to other physical processes. Examples of such processes will be discussed below.

The main problem in measuring the EBL is that, in observations from the Earth's vicinity, the contribution of the EBL is small compared to other bright sources of radiation in the sky. These include the radiation of the stars of the Galaxy, the radiation of dust inside the solar system, that scatters sunlight (zodiacal light, ZL) and the radiation of galactic dust. As a result, the EBL intensity is obtained as the difference between two large numbers: the observed brightness of a part of the sky and unwanted contributions, which leads to large errors and requires precise instrument calibration. This difficulty was first shown theoretically by Partridge & Peebles (1967b), long before the first EBL measurements. Awareness of the problem stimulated the creation of alternative methods for measuring EBL.

Technical difficulties did not allow detection of EBL until the end of the 90s. Significant technical progress and the construction of new astrophysical observatories at the turn of the century led to several independent EBL estimates appeared at once, based on fundamentally different approaches. In the following, we briefly review the main methods for measuring EBL, summarize the main results, and describe the challenges. We will use standard notation and indicate the EBL wavelength λ in µm and EBL surface brightness νI ν = λI λ in nW m -2 sr -1 .

Direct measurements. First method is the direct measurement of the EBL intensity. In this approach, EBL arises as a component that remains after subtracting the dominant contribution of zodiacal light and galactic dust and stars from the observed radiation flux.

The first experiments in 80s and early 90s were carried out using sounding rockets. These experiments did not result in detection due to limited observational time and large systematic errors, but set upper bounds on the EBL intensity. At the same time The Infrared Astronomical Satellite (IRAS) [START_REF] Neugebauer | The Infrared Astronomical Satellite (IRAS) mission[END_REF] launched in 1983 performed first all-sky survey at 12, 25, 60 and 100 µm. IRAS observations firmly showed presense of diffuse emission from the Solar System and Galaxy. However, the characteristics of instrument calibration also did not allow reliable identification of the EBL.

The first direct detection of the EBL was carried out by the Cosmic Background Explorer (COBE) [START_REF] Wright | Detection of the cosmic infrared background at 2.2 and 3.5 microns using dirbe observations[END_REF], launched in 1989. The main purpose of this experiment was to study the properties of CMB using the Far Infrared Absolute Spectrophotometer (FIRAS) instrument, sensitive from 125 µm to millimeter wavelengths. The Diffuse Infrared Background Experiment (DIRBE) was another instrument onboard the COBE and was designed to obtain brightness maps at ten frequency bands (1.25, 2.2, 3. [START_REF] Matsumoto | IRTS observation of the near-infrared extragalactic background light[END_REF] (restated in [START_REF] Matsumoto | Reanalysis of the near-infrared extragalactic background light based on the IRTS Observations[END_REF]) measured EBL with the Near Infrared Spectrometer (NIRS) from 1.4 to 4.0 µm with 0.12 µm resolution. IRTS predicted even stronger EBL than COBE measurements, with a peak intensity of 60 nW m -2 sr -1 at about 1.5 µm. AKARI (launched in 2006) observation with the InfraRed Camera (IRC) in the 1.8-5.3 µm range is also fairly consistent with observations by COBE/DIRBE and IRTS/NIRS. The most recent measurements by the Cosmic Infrared Background Experiment (CIBER, launched in 2009) [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF] at wavelengths 0.8-2 µm confirm previous results with intensity of 42.7 nW/m 2 /sr using Kelsall ZL model [START_REF] Kelsall | The COBE diffuse infrared background experiment search for the cosmic infrared background. 2. Model of the interplanetary dust cloud[END_REF] (CIBER 'nominal' EBL) and 28.7 nW m -2 sr -1 for the model independent analysis (CIBER 'minimal' EBL).

Let us also note several very interesting special methods for direct estimation of the EBL, which make it possible to reduce the influence of undesirable background The solid and dashed lines correspond to the numerical models of the EBL. The red dots are obtained from the galaxy number count method. The gray dots indicate direct measurements of the EBL. Yellow band corresponds to the constraints from the γ-ray observations. Source: [START_REF] Saldana-Lopez | An observational determination of the evolving extragalactic background light from the multiwavelength HST/CANDELS survey in the Fermi and CTA era[END_REF].

contributions. First, a unique opportunities for EBL measurements were opened up by the spacecrafts which have reached the outer solar system: Pioneer 10 and 11 and New Horizons. At this distance, the contribution of zodiacal light is negligible and the only remaining foregrounds are the cumulative light of stars and diffuse galactic light. Recent reanalysis the Pioneer 10/11 data allows detection of the EBL at 0.44 and 0.64 µm with intensity of 7.9±4.0 and 7.7±5.8 nW m -2 sr -1 respectively (see Table A.3). The significance of the detection is at the 2σ level. In addition to this, observations with the New Horizons Long-range Reconnaissance Imager (LORRI) of the New Horizons revealed EBL at λ ∼ 0.6 µm with λI λ = 15.9 ± 4.2 or 18.7±3.8 nW m -2 sr -1 depending of the analysis. Finally the screening effect of a dark cloud on the background light can be used. Namely the EBL in this method is the difference of the sky brightness in the direction of an opaque high galactic latitude dark cloud and a transparent nearby area of sky. This allows one to effectively cut off the unwanted contribution of the Galaxy and Zodiacal Light, adding, however, the contribution from particles of the cloud itself. Results of the application of this method to the real data of ESO's Very Large Telescope (ESO VLT) are shown in the Table A.2.

Galaxy number count method. A complementary approach relies on the galaxy number counts and the EBL is estimated as the sum of the fluxes of individual galaxies. The most suitable observations for this method are long-exposure deep sky surveys that resolve even the faintest galaxies. The main advantage of this approach is that it allows one to get rid of unwanted background contributions. On the other hand, as it follows from description this approach sets only a robust lower limit on the EBL intensity. The actual intensity may be higher due to possible undercount of the galaxies beyond the sensitivity level of the instrument or underestimation of the faintest areas of the objects detected.

This method has been developed in parallel with the development of the direct observational approach. Its application to observations of various instruments provided a lower limit on the EBL at different wavelengths. The first results of this approach in 1990-2000 include an estimate of EBL with the Hubble Space Telescope Deep Fields (HDF) at 0.36, 0.45, 0.67, 0.81, 1.1, 1.6, and 2.2 µm [START_REF] Keenan | The Resolved Near-Infrared Extragalactic Background[END_REF], [START_REF] Madau | Deep galaxy counts, extragalactic background light, and the stellar baryon budget[END_REF], and IRAS, Infrared Space Observatory (ISO) and Submillimetre Common-User Bolometer Array (SCUBA) EBL estimates in the infrared range. Subsequently, lower limits in the infrared range were significantly improved with the next generation experiments and observatories: Balloon-borne Large-Aperture Submillimeter Telescope (BLAST) [START_REF] Devlin | Over half of the far-infrared background light comes from galaxies at z > 1.2[END_REF], Subaru Telescope [START_REF] Keenan | The Resolved Near-Infrared Extragalactic Background[END_REF], Spitzer Space Telescope [START_REF] Fazio | Number counts at 3 < lambda < 10 um from the Spitzer Space Telescope[END_REF] and Herschel Space Observatory [START_REF] Béthermin | HerMES: deep number counts at 250 µm, 350 µm and 500 µm in the COSMOS and GOODS-N fields and the build-up of the cosmic infrared background[END_REF] (see Table A .1).

Even the first comparisons of the two approaches showed a significant contradiction: the results of direct observations predicted a much higher EBL intensity than the estimates from below. Indeed, if the total contribution of galaxies to the EBL is dominant, then it was expected that both methods would give a comparable result. However, in the region of the first maximum (1-2 µm), the EBL from the galaxy number count method does not exceed 10 nW m -2 sr -1 , while direct measurements give peak values of the order of 40-60 nW m -2 sr -1 . A similar discrepancy is also observed at the second maximum (100-200 µm), see Fig. 2.3. The controversy between the two methods has been observed for the last 20 years and remains even for the most recent measurements. This may indicate the presence of truly diffuse EBL radiation not attributed to the discrete sources. At the same time, it should be taken into account that the direct measurement method is based on complex modeling of undesirable background contributions and therefore may be subject to systematic error. Also, as we will discuss in the next subsection, a completely independent method for measuring EBL based on observations of gamma rays is consistent with the results of the galaxy counts, see Fig. 2

.3.

Constraint from γ-ray observation. Shortly after the discovery of the cosmic microwave background (CMB) radiation in 1965, it was realized that low energy radiation fields can interact with high energy γ-rays, producing electron-positron pairs and effectively reduce γ-ray mean free path [START_REF] Nikishov | Absorption of high-energy photons in the Universe[END_REF][START_REF] Gould | Pair Production in Photon-Photon Collisions[END_REF]. This effect leads to the distance-dependent suppression of the γ-ray flux from extragalactic sources at the highest energy:

J obs (E γ ) = J 0 (E γ ) exp(-τ γγ (E γ )) (2.1)
where J obs (E γ ) is the observed flux, J 0 (E γ ) -intrinsic flux near the source and τ γγ (E γ ) is the EBL optical depth. Thus observation of the suppression of the γ-ray signal can provide information on the EBL intensity. Development of the new observatories such as the Compton Gamma-Ray Observatory (CGRO) with the Energetic Gamma Ray Experiment Telescope (EGRET) telescope, the Whipple observatory, the High Energy Gamma-Ray Array (HEGRA) and the Cherenkov Array at Themis (CAT) led to the detection of the first extragalactic gamma-ray sources at energy up to ∼ 20 TeV in the 1990s. This allows application of this method to real data. Recent gamma-ray constraints are set by HESS [START_REF] Abramowski | Measurement of the extragalactic background light imprint on the spectra of the brightest blazars observed with H.E.S.S[END_REF] and MAGIC [START_REF] Acciari | Measurement of the extragalactic background light using MAGIC and Fermi-LAT gamma-ray observations of blazars up to z = 1[END_REF], as well as the space telescope Fermi/LAT [START_REF] Abdollahi | A gamma-ray determination of the Universe's star formation history[END_REF]. Typical estimates of the EBL intensity by this method are consistent with the results from the galaxy count method (see Fig. 2.3).

However, the precision of the measurements of the EBL using the high-energy suppression suffers from the uncertainty of the knowledge of the intrinsic primary source spectrum. Furthermore, γ-ray measurements do not provide a measurement of the EBL flux at a particular wavelength. The wavelength resolution of the measurements is limited by the width of the pair production cross-section that peaks at the centre of mass energy, which is somewhat above the threshold E thr = 2mc 2 of twice the rest energy of the electron and decreases as E -1 ln(E) at E E thr . Measurements of the EBL from the γ-ray data conventionally adopt assumptions about the shape of the EBL spectrum.

Numerical modelling. Another large group of approaches is based on the numerical calculation of the EBL. These methods differ by the basic assumptions made and may be divided into several subgroups. The first category is based on the evolution in time of the galaxy properties. Such an evolution may be obtained either by direct observations of galaxy evolution or indirectly estimated according to some prescription. This method was used in some of the first models of the EBL [START_REF] Madau | The Star formation history of field galaxies[END_REF], [START_REF] Pei | Cosmic Histories of Stars, Gas, Heavy Elements, and Dust in Galaxies[END_REF], [START_REF] Kneiske | Implications of cosmological gamma-ray absorption -I.evolution of the metagalactic radiation field[END_REF], [START_REF] Kneiske | Implications of cosmological gamma-ray absorption. 2. Modification of gamma-ray spectra[END_REF] and it is still relevant due to the recent large-scale surveys at ultraviolet and infrared wavelengths [START_REF] Finke | Modeling the Extragalactic Background Light from Stars and Dust[END_REF]. Another category of approaches uses a backward evolution of the present day galaxy emissivity according to a given prescription and accounts for the change of the star formation rate and other redshift dependent factors [START_REF] Franceschini | The extragalactic background light revisited and the cosmic photon-photon opacity[END_REF], [START_REF] Franceschini | A long-wavelength view on galaxy evolution from deep surveys by the infrared space observatory[END_REF], [START_REF] Rowan-Robinson | The star-formation history of the universe -an infrared perspective[END_REF], [START_REF] Franceschini | The extragalactic optical-infrared background radiations, their time evolution and the cosmic photon-photon opacity[END_REF], [START_REF] Dominguez | Extragalactic Background Light Inferred from AEGIS Galaxy SED-type Fractions[END_REF], [START_REF] Stecker | An Empirical Determination of the Intergalactic Background Light from UV to FIR Wavelengths Using FIR Deep Galaxy Surveys and the Gamma-ray Opacity of the Universe[END_REF]. The last scenario for treating EBL numerically is so called the forward evolution approach. This is the most fundamental and challenging way of calculation. It starts with the very basic assumptions about the Universe and it's structure formation [START_REF] Primack | Probing galaxy formation with TeV gamma-ray absorption[END_REF], [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. Being the most complete, the method allows to pinpoint the physical reasons for time-dependence of the EBL. Starting from [START_REF] Franceschini | Extragalactic optical-infrared background radiation, its time evolution and the cosmic photon-photon opacity[END_REF], EBL models are conventionally in agreement with gamma-ray constraints and close to the lower bound from galaxy count method (see Fig. 2.

3).

Thus at the moment there is a contradiction in the estimates of the EBL intensity: direct measurements of the EBL lead to significantly higher values of the intensity than all other methods (see Figure 2.3). In what follows, we consider two possible interpretations of this fact. The first is that the accuracy of direct measurements suffers from the unaccounted systematic effects. In such a case, direct measurements can be considered as an upper limit on the intensity of the EBL. The second possibility is that direct measurements are still correct, but the shape of the EBL spectrum is more complicated. It is possible that there is a narrow additional contribution to the EBL spectrum, located in the wavelength range λ ≈ 1 µm with an intensity that satisfies direct observations. This possibility has not been accounted for in the γ-ray analysis and also was not taken into account when constructing numerical models. We will consider whether the presence of such a feature in the spectrum can be reconciled with γ-ray observations.

Building a new EBL model

First, we assume that the direct measurements of the EBL provide an upper bounds on the intensity. Combining these upper bounds with the lower limits from the galaxy count method we determine the allowed range the EBL intensity. Constraints on the EBL spectrum make it possible to constrain the main astrophysical parameters that determine the intensity, shape and evolution of EBL. In order to do this, we will construct a new EBL model and explicitly identify the main parameters on which it depends.

Our goal is to build the empirical model of the formation and evolution of stars and galaxies. The model we plan to construct will be then directly translated into the EBL spectrum. It's unavoidable for the model to be dependent on a number of astrophysical parameters, which are known with a limited precision. We will keep these dependencies explicit, so the parameter space may be further probed with the Markov Chain Monte Carlo (MCMC) technique.

At first the contribution of stars to the EBL will be described, and then the absorption and re-emission of light by dust will be taken into account. We use the standard ΛCDM cosmological model to derive the functional dependencies between comoving distance r, redshift z and time t.

r(z) = c H 0 z 0 dz Ω m (1 + z) 3 + Ω Λ (2.2) z(t) = Ω Λ Ω m 1 3 sinh( 3 2 Ω Λ 1 2 tH 0 ) -2 3 -1, (2.3)
where the Hubble constant H 0 = 67.8 km/s/Mpc, the matter density parameter Ω m = 0.309, the dark energy density parameter Ω Λ = 0.691 and c is the speed of light.

The parameters of the cosmological model are measured within the high accuracy [START_REF] Ade | Planck 2015 results. XIII. Cosmological parameters[END_REF] and therefore are considered fixed for the purposes of the present study. We assume Chabrier initial mass function (IMF) to describe mass distribution of newborn stars [START_REF] Chabrier | Galactic stellar and substellar initial mass function[END_REF]. Following Gilmore et al. [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF] we assume that the IMF doesn't depend on the redshift:

ξ(m) =    C imf m e -(log(m)-log(m 0 )) 2 2D if m 1, k m -a imf if m > 1.
(2.4)

Masses of stars are measured in the units of the Solar mass. The parameters m 0 , D and a will undergo variations while C imf and k are determined from normalization and continuity. The initial mass function is bounded from below by the minimal mass m min = 0.08. If the mass of the star is less than m min then the hydrogen thermonuclear fusion reactions in the core will not ignite. The upper bound of the mass distribution m max is another free parameter of the model. The initial values for the variable parameters of the IMF are: m 0 = 0.079, D = 0.69, a = 2.3 and m max = 100. For the sake of convenience we will use mass normalization of the IMF:

mmax m min m ξ(m) dm = 1 (2.5)
The mass of stars, born at redshift z per unit time per unit volume is described by the star formation rate (SFR) ψ(z). Here we will use SFR, obtained in [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. We parametrize the SFR with a single multiplication factor C sfr with initial value equal to unity.

To determine the spectrum of the star with the mass m at the age η near the surface of this star we use stellar evolution tracks calculated by Hurley et al. [START_REF] Hurley | Comprehensive analytic formulae for stellar evolution as a function of mass and metallicity[END_REF], which give us radius R, temperature T and lifetime η life dependence on age. Let's define this spectrum as B s (λ, m, η, z = 0) there the last argument is reserved for the redshift.

B s (λ, m, η, 0) = 2πhc 2 λ 5 1 e hc kT (m,η)λ -1 (2.6)
If we consider the star at the distance r from the detector, the registered emission will be different from (2.6), due to geometric dilution and cosmological redshift. The time η now will denote the age of the star at the moment of the emission of radiation registered today. To calculate new spectrum we use a standard trick: replace λ by λ 1+z and divide the full expression by (1 + z) 3 . First (1 + z) is responsible for decreasing the energy of the photons, second is responsible for decreasing the number of registered photons per unit time and the last is responsible for broadening of the spectrum. Also the expression must be multiplied by the geometrical factor R(η) r(z)

2 : B s (λ, m, η, z) = R(m, η) r(z) 2 (1 + z) -3 B( λ 1 + z , m, η, z = 0) (2.7) Distance r(z) is calculated by (2.2).
Using the above definitions one may write an expression to determine the spectrum of a galaxy of age η g and redshift z g . Moreover, we assume that the star formation began in all galaxies simultaneously at time t i . The light emitted at time t i is now registered with the redshift z i = z(t i ). Then the galaxy age η g and redshift z g are related with the expression:

η g = t(z g ) -t(z i ) (2.8)
The function t(z) is the inverse of (2.3). Thus we obtain for galaxy spectrum:

G s (λ, z g ) = mmax m min dm η end (m) 0 dη B s (λ, m, η , z g ) ξ(m) ψ(t(z g ) -η ) (2.9)
Upper limit η end in time integral depends on the mass and determined by

η end (m) = min(η g , η life (m)) (2.10)
On the other hand the presence of dust and absorbers in the galaxies and interstellar medium should be taken into account. Following [START_REF] Charlot | A Simple Model for the Absorption of Starlight by Dust in Galaxies[END_REF] we assume, that the star formation took place only in the giant molecular clouds which will cover the newborn stars with the shell of dust and gas. Dust grains consist mostly of graphite and silicate, so they have non-zero absorption coefficients at visible and ultraviolet wavelengths. Thus they will be heated by the newborn stars and re-emit in the infrared part of the spectrum.

Let us now specify the dust model. According to [START_REF] Charlot | A Simple Model for the Absorption of Starlight by Dust in Galaxies[END_REF] we assume that the birth clouds have finite lifetimes η c . The parameters of the cloud such as particle number densities and outer radius R c stay constant throughout its life, whereas the temperature of different components may change due to evolution of the stars in the center. Moreover, we assume for simplicity that dust has constant particle number density through the entire cloud which is defined by the parameter n d . We do not account for the other features of the internal structure of the cloud and following [START_REF] Charlot | A Simple Model for the Absorption of Starlight by Dust in Galaxies[END_REF] the optical depth will be treated as a free parameter with the two characteristics. First one is the normalization, which corresponds to the optical depth at wavelength λ 0 = 5500 Å and the second is a spectral slope n. Thence, optical depth at the arbitrary wavelength is determined by

τ (λ) = τ λ 0 λ λ 0 -n (2.11)
We assume that the star formation occurs in the center of the birth cloud, so there exists such a value ρ R c that all newborn stars are located in the imaginary sphere S in with radius ρ. First of all we calculate spectrum B c (λ, η, ρ) of the newborn stars at the age η on the boundary of the sphere S in . On this step we neglect absorption and dependence on redshift:

B c (λ, η, ρ) = mmax m min dm η cend (m) 0 dη R(η ) ρ 2 B(λ, m, η , z = 0) ξ(m) ψ c (η ) (2.12)
where η cend (m) = min(η c , η g , η life (m)) by analogy with (2.10) and ψ c (η) denotes mass of the matter converted into stars in the cloud at the age η per unit time. The total mass of the gas converted into stars throughout life of the cloud M 0 is given by:

ηc 0 ψ c (η ) dη = M 0 (2.13)
For the sake of simplicity suppose that ψ c (η) stays constant in time so that:

ψ c = M 0 η c (2.14)
On the other hand, ψ c directly proportional to ψ(t) with the coefficient n cl (t) which denotes the number of clouds per unit volume.

n cl (t) = ψ(t) ψ c = ψ(t) η c M 0 (2.15)
Equation (2.15) will be used for calculating the spectrum of newborn stars in the galaxy. Applying (2.11) we estimate the absorption for light propagating from ρ to r.

B c (λ, η, r) = ρ r 2 e -τ λ 0 ( λ λ 0 ) -n B c (λ, η, ρ) (2.16)
where we again neglected the redshift. Thus B c (λ, η, R c ) denotes contribution of the stars to the spectrum of the cloud calculated on the cloud's boundary. The full spectrum of the cloud contains an additional contribution of dust. Now we describe the assumptions regarding properties of dust. Following [START_REF] Martinez-Galarza | Infrared Emission by Dust Around lambda Bootis Stars: Debris Disks or Thermally Emitting Nebulae?[END_REF] the dust grains are spherical with sizes a distributed by the power law:

dn(a, a + da) = C d a -n dust da (2.17)
Initial value for the slope of distribution n dust = 3.5. As the sizes of grains are constrained, distribution should be truncated with minimal grain size a min and maximal a max . We will start setting a min = 5 nm and a max = 500 nm. As mentioned above the dust particle number density is constant and allows us to determine normalization constant C d .

C d = (n dust -1) n d a -n dust +1 min -a -n dust +1 max (2.18)
Equilibrium temperature T d (a, r) of dust grain of size a at the distance of r may be found from the energy balance equation:

πa 2 ∞ 0 Q a abs (λ)B c (λ, η, r) dλ = 4πa 2 ∞ 0 Q a abs (λ)B Pl (λ, T d (a, r, η)) dλ (2.19)
where emissivity and absorptivity coefficient of grains with size a are equal due to Kirchoff's law and denoted as Q a abs (λ). B Pl (λ, T d (a, r, η)) is the Planck's law of emissivity. To find out the full contribution of dust we should integrate over the grain size and over the distance r. Thus we obtain final formula for calculating dust spectrum:

B d (λ, η) = Rout ρ dr 4πr 2 amax a min da a R out 2 C d a -n dust Q a abs (λ) B Pl (λ, T d (a, r, η)) (2.20)
In the mid-infrared the principal contribution comes from the polycyclic aromatic hydrocarbons (PAH) molecules, which absorb the starlight producing the spectral lines. The corresponding cross-sections σ, damping constants γ, and cutoff energies λ c assuming Lorentzian spectral shape are well known and may be found, for example, in [START_REF] Krugel | The Physics of Interstellar Dust[END_REF]. The PAH spectra is calculated as follows:

B PAH (λ, η) = Rout ρ dr 4πr 2 n PAH ∞ λc dλ B c (λ , η, r) i   σ i γ i c 3 /λ 4 π 2 c 2 λ 2 -c 2 λ 2 0 2 + γ i c 2λ 2   
(2.21) where n PAH is the particle number density of homogeneously distributed PAH molecules, and summing goes over the full set of resonances.

Summarizing the above the full spectrum of the galaxy with redshift z g is expressed as:

G(λ, z g ) = G s (λ, z g ) + G c (λ, z g ) (2.22)
where:

G s (λ, z g ) = mmax m min dm η end (m) 0 dη B s (λ, m, η , z g ) ξ(m) ψ(t(z g ) -η ) (2.23) η end (m) = max(min(η g , η life (m) -η c ), 0) (2.
G c (λ, z g ) = ηc 0 dη [B d (λ, η , z g ) + B PAH (λ, η , z g ) + B c (λ, η , R out , z g )] η c M 0 ψ(t(z g ) -η )
(2.25) Finally, the spectrum of the whole Universe is obtained by integration over the distance, taking into account Calzetti's attenuation law C calz (λ) for diffuse dust component [START_REF] Calzetti | The Dust content and opacity of actively star-forming galaxies[END_REF].

U (λ) = r(z i ) 0 dr 4πr 2 C calz (λ) G(λ, z g (r)) (2.26)

Constraints on astrophysical parameters of the EBL

The model, described in the previous section, keeps an explicit dependence on a number of parameters. The parameter list is shown in the Table 2.1. The preferred values and confidence ranges for the parameters will be calculated using the Markov Chain Monte Carlo method. The Markov Chain is constructed with the Metropolis algorithm with the standard likelihood function

L = i e - (y-y i ) 2 2σ 2 i , ( 2.27)
where y i and σ i are the values and errors of the experimental points. The one-sided limits are included in the product only if y > y i for lower limits and if y < y i for upper ones.

The likelihood is based on the three categories of the data points which are the lower and upper limits of the EBL and direct measurements. The lower EBL limits usually come from the galaxy count method and provide strict constraints, see Table A.1. The upper limits appear either from the direct measurements at the wavelengths where contribution of the Zodiacal light is significant or from the direct measurements combined with some special technique to subtract the excess, see Table A.2. We do not include here the constraints from the gamma-ray observations due to possible uncertainties in AGN's spectra. For the wavelengths longer than 200 µm the direct measurements are used as it is believed that Zodiacal light contribution is negligible. We have also added Pioneer 10/11 results because their measurements were conducted far from the region of the potential Solar system background affects, see Table A

.3.
Each step of the chain requires to repeat the full calculation of the EBL, so we employ the parallel implementation of the code with OpenMP. Given the model and the chain, we begin with proving the consistency of the model and then discuss possible systematics. The EBL spectrum with the best likelihood set of the parameters is shown in Figure 2.4. The parameters inferred from the MCMC are in a good agreement with the commonly used values, see Table 2.1. The spectrum demonstrates good accordance in ultraviolet, optical and far infrared bands with the popular models [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF] and [START_REF] Franceschini | The extragalactic optical-infrared background radiations, their time evolution and the cosmic photon-photon opacity[END_REF].

We do not account for the AGN's contribution because it is sub-dominant and moreover is subject to large uncertainty [START_REF] Stephen | Modelling the cosmic spectral energy distribution and extragalactic background light over all time[END_REF]. This contribution may modify the result considerably for λ < 0.4 µm and therefore we didn't account for the experimental points in this range. Another assumption of the model is that all the cloud's parameters are equal. This means that all the diversity of the star forming objects is reduced to only one type of the clouds. We believe this is enough to grasp the main behaviour of these objects.

Let us now describe the outcomes of the main part of the research. All the parameters are divided into four group according to the way they are treated. First group contains the parameters of Chabrier Initial Mass Function namely m 0 , D, a imf . They are already limited in [START_REF] Chabrier | Galactic stellar and substellar initial mass function[END_REF] so we set the existing confidence intervals as an available freedom for each parameter. The calculations show that the EBL is insensitive to m 0 and D, while it depends strongly on a imf . The resulting value is a imf = 2.29 +0. 15 -0.24 which is more narrow than the limits from [START_REF] Chabrier | Galactic stellar and substellar initial mass function[END_REF]. The parameter a imf governs the number of massive stars thus this result indicates that the massive stars play an important role in the EBL formation.

The next group contains the dust parameters and includes minimal and maximal grain sizes and the slope of the distribution n dust . We let them vary in a wide range except for lower limit for the minimal grain sizes and upper limits for the maximal grain sizes which are 10 -3 µm and 1 µm correspondingly. After the simulations we found that the grain sizes are constraint as: a min < 1.75 • 10 -3 µm and a max < 22.8 • 10 -3 µm. This is due to the higher emissivity of the bigger grains and thus greater contribution to the EBL. The slope of the distribution n dust should be more than 3.83. The result may be expected at the large grains contribute more and their amount should be limited.

The third group collects all the parameters of the clouds, see Figure 2.5. Basically, the properties of the molecular clouds are known from the observations of the Milky Way, for more details see, for instance, [START_REF] Murray | Star Formation Efficiencies and Lifetimes of Giant Molecular Clouds in the Milky Way[END_REF]. Observations shows that the radii of the clouds lie in the range from 1 pc to 35 pc with the largest reach a size of 100 pc. The lifetimes of the clouds with the mean mass are 17 ± 4 Myr [START_REF] Murray | Star Formation Efficiencies and Lifetimes of Giant Molecular Clouds in the Milky Way[END_REF]. The following values are obtained from the MC for radius R c = 6.1 +1. 4 -1.2 pc and lifetime η c = 6.0 +8.5 Myr. These are in a good agreement with the experimental estimates. The relatively small value of the cloud radius may be reconciled with the observed larger clouds if we allow the large clouds to be inhomogeneous and divided into the separate region of star formation. The dust particle density n d can't be measured directly but may be derived from the mass of the cloud. One may derive the star formation efficiency = 0.035 which is defined as the ratio between the mass of the newborn stars in the cloud and the whole mass of the cloud. The calculation uses common assumption that the gas to dust ratio is 100. The value of is in line with the modern observations [START_REF] Murray | Star Formation Efficiencies and Lifetimes of Giant Molecular Clouds in the Milky Way[END_REF].

The last group contains the global parameters of the Universe, see Figure 2.6. It includes the redshift of the beginning of the star formation, the mass of the star at which IMF is truncated and scaling factor for the star formation rate. We obtained that the EBL does not depend on the epochs with the redshift 4 < z < 10. Truncation mass is the maximal mass of the star in our model. We obtain that the smaller truncation masses are more probable that the bigger ones. The reasoning for this results is that the massive stars are too hot and thus overheat the dust. Finally, the star formation rate scaling factor is well constrained and at 1σ is 1.25 +0. 44 -0.24 . In the calculations above the results of CIBER and ESO VLT/FOR were used as the upper limits A.2. It is practical to know how the conclusions of the present paper change if these results are interpreted as the direct measurements. For this purpose a separate Markov Chain was constructed. The account of double-sided CIBER and ESO VLT/FOR gives the higher value of star formation rate C sfr = 1.76 +0. 14 -0.13 and lower values of the radius and lifetime of the clouds R c < 3.7 ± 2.0, η c < 3.1 ± 2.2.

Constraints on the additional narrow feature in the EBL spectrum

Let us now consider the second possibility in which the EBL spectrum has an additional narrow spectral feature with intensity that satisfies direct observations. Naively, it seems that all direct measurements lie above the band allowed by the gamma-ray data, and, therefore, an agreement between them is impossible. At the same time, EBL measurements from gamma-ray data traditionally assume a smooth shape of the EBL spectrum with two broad maxima and are therefore insensitive to narrow spectral features. Therefore, the possibility remains that γ-ray and direct measurements are still consistent if, for example, the EBL has an additional spectral component concentrated in a narrow wavelength range, which was not taken into account in the analysis of γ-ray data. Examples of such a relatively narrow function are shown as a red dotted line in Figure 2.7.

Such features could originate from processes that have a characteristic wavelength or energy scale. Examples are given by the emission from Population III stars [START_REF] Michael | The contribution of the first stars to the cosmic infrared background[END_REF], which produce strong emission at the wavelength of Ly alpha lines redshifted towards [START_REF] Tsumura | Low-Resolution Spectrum of the Extragalactic Background Light with AKARI InfraRed Camera[END_REF], while green asterisks are from the reanalysis of IRTS [START_REF] Matsumoto | Reanalysis of the near-infrared extragalactic background light based on the IRTS Observations[END_REF]. Pink diamonds are from the reanalysis of COBE/DIRBE measurements [63] [64]. Black data points together with blue and orange systematic uncertainty are derived from CIBER [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF] and correspond to nominal and minimum EBL models. Lower limits: Dark blue upward arrows combine EBL lower limits obtained by different experiments: GALEX [START_REF] Xu | Number counts of GALEX sources in FUV (1530A) and NUV (2310A) bands[END_REF], Hubble Deep Field [START_REF] Madau | Deep galaxy counts, extragalactic background light, and the stellar baryon budget[END_REF], Subaru [START_REF] Keenan | The Resolved Near-Infrared Extragalactic Background[END_REF], and Spitzer/IRAC [START_REF] Fazio | Number counts at 3 < lambda < 10 um from the Spitzer Space Telescope[END_REF]. EBL from γ-ray absorption: Striped lime, red, and purple bands are from MAGIC [START_REF] Acciari | Measurement of the extragalactic background light using MAGIC and Fermi-LAT gamma-ray observations of blazars up to z = 1[END_REF], HESS [START_REF] Abramowski | Measurement of the extragalactic background light imprint on the spectra of the brightest blazars observed with H.E.S.S[END_REF], and Fermi/LAT [START_REF] Abdollahi | A gamma-ray determination of the Universe's star formation history[END_REF], correspondingly. Modelling: Dark red solid line is for baseline EBL model of [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. Black striped band shows the allowed range of EBL models obtained with the global fit to EBL by [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF]. Dash-dotted orange and dashed red lines are additional wide and narrow Gaussian components. near infrared, hypothetical decaying particles [START_REF] Kohri | Can decaying particle explain cosmic infrared background excess?[END_REF], and axion-like particles [START_REF] Kalashev | Cosmic infrared background excess from axionlike particles and implications for multimessenger observations of blazars[END_REF] that are expected to produce a feature at the energy close to the particle rest energy.

In what follows, we explore the constraints on narrow spectral features in the EBL spectrum imposed by existing γ-ray data.

Data selection

Our analysis relies on the combined time-averaged Fermi/LAT and IACT spectra of blazars. Such combined spectra could only be produced for sources that were the subject of long, multi-year monitoring observations with IACTs. This is not the case for most of the TeV detected blazars because ground-based telescopes typically observe these sources during specific activity periods.

To define a set of sources for which the production of time-averaged combined Fermi/LAT and IACT spectra is possible, we performed source selection in the following way. The initial selection was done from the TeVCat online source catalogue [START_REF] Wakely | TeVCat: An online catalog for Very High Energy Gamma-Ray Astronomy[END_REF] 1 , which aggregates observations performed by different observatories. The TeVCat catalogue includes 71 sources that are classified as blazars, with the following sub-classes: 'HBL', 'IBL', 'LBL', 'FSRQ', 'BL Lac (class uncertain)', and 'Blazar'. Out of the 71 sources from the TeVCat, we only retained sources with spectroscopically measured redshift, derived from emission lines. This restricts the list of blazars to 42 sources.

Then, we imposed a requirement that the selected blazars should have been observed for at least two years during the period of operation of the Fermi telescope, that is, after 2008 and the reported spectrum should be averaged over the time of observation. IACTs only typically observe sources for several dozens of hours per year during several months of observations. The preference of only observing the blazar during flaring activity introduces a bias in the flux and spectral measurements. To avoid possible effects of this kind of bias, we selected only those sources for which long-term exposures over several years have been reported in the literature. In case the long-term observation data were published in the form of low-flux-state and high-flux-state spectra, we used low-state ones. If the spectra were not divided on low-and high-states, we checked that the blazar did not exhibit strong flares during the IACT observational period, when the flux increased more than ten times as compared to the average. After imposing this constraint, 21 blazars have been left in the sample.

We have further imposed a bound on the source flux F > 0.03F Crab , where F Crab is the flux of the Crab Nebula in order to assure sufficient quality of the gamma-ray spectral measurements. This leaves five blazars in our source sample. In addition, we added PKS 2005-489 and PKS 2155-304, which satisfy all our criteria, but they were observed before 2008. All the blazars selected for the analysis are listed in Table 2.2.

For each selected blazar from Table 2.2, we calculated the time-averaged spectrum The Galactic interstellar emission model used in the analysis was 'gll_iem_v07.fits'. The isotropic background model was 'iso_P8R3_SOURCE_V2_v1.txt'. The instrument response functions were 'P8R3_SOURCE_V2'. Final spectra were obtained using Fermi likelihood analysis with the help of the gtlike routine. We imposed standard cuts for the Fermi/LAT unbinned likelihood analysis, such as Earth zenith angle < 90 • . Source models were created on the basis of the fourth Fermi/LAT Source Catalogue 4FGL [START_REF]Fermi Large Area Telescope Fourth Source Catalog[END_REF]. The spectra were calculated in ten energy bins covering energy range from 1 GeV up to 1 TeV. To account for the systematic uncertainty of the effective area, we added 10% errors to every bin.

Spectral fitting

The effect of absorption by EBL is imprinted on intrinsic source spectrum. We find that a broken powerlaw model provides a satisfactory description of the observed spectra of the selected blazars. The model has four parameters: normalisation A, powerlaw indices Γ 1 and Γ 2 = Γ 1 -∆Γ, and the break energy E br :

   F 0 (E) = A (E/E br ) -Γ 1 for E < E br F 0 (E) = A(E/E br ) -Γ 2 for E > E br . (2.28)
The observed spectrum is

F (E) = F 0 (E) e -τ (E,z) , (2.29)
where τ is the energy-dependent optical depth for the pair production on EBL for a source at the redshift z. [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF]. Dashed and solid magenta lines correspond to the baseline EBL with additional nominal EBL excess detected by CIBER [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF].

Figure 2.9 -Cross section of three dimensional parameter space calculated for the bump position µ = 1.5 µm. Grey areas correspond to the cross sections of 95% confidence level constraints on the parameters of the log-Gaussian EBL bump, obtained with different blazars. The strongest constraint is imposed by 1ES 1218+304 data. For this source, the cross section of 68% confidence contour around the best-fit value is also shown. Blue cross in the excluded region and the star within the 68% confidence contour show the two bump measurements reported by CIBER [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF].

The next step is to choose the EBL model to absorb intrinsic spectra. We chose the model of [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF] as it is one that has the lowest optical depth. This model has a nearly identical EBL spectrum to that of the model of [START_REF] Franceschini | Extragalactic optical-infrared background radiation, its time evolution and the cosmic photon-photon opacity[END_REF] in the visible and near-infrared range at z = 0. The results of our study are not significantly impacted by choosing this particular EBL model, which has the near-infrared and visible band flux at the level that is close to the lower bound derived from the galaxy counts at z = 0.

The spectra of the seven selected sources fitted with the spectral model described above are shown in Figure 2.8. One can see that Fermi/LAT and IACT parts of the spectra are in good agreement with the fitted broken powerlaw model with χ 2 /d.o.f. < 1 for each blazar.

Modelling of infrared and visible excess

A good fit to the data is provided by the broken powerlaw model modified by the effect of absorption on 'low-flux' EBL at the minimal possible level estimated from galaxy counts. Although one could notice that all the spectral fits include a break in the intrinsic source spectra at around 100 GeV, that is, at the energy level above which the EBL effect becomes important. This hints to the possibility that, in fact, a higher level of EBL could be accommodated by the spectral fits. Weaker spectral break in the intrinsic spectrum could be compensated by stronger absorption on the EBL.

To explore this possibility, we introduced an excess in the EBL spectrum, as a broad spectral feature localised around a characteristic energy. Such a spectral feature could be generated, for example, by a specific star-formation history (Population III stars [START_REF] Michael | The contribution of the first stars to the cosmic infrared background[END_REF]) or by the photons injected by interactions of exotic decaying particles [START_REF] Kohri | Can decaying particle explain cosmic infrared background excess?[END_REF], that is, axion-like particles [START_REF] Kalashev | Cosmic infrared background excess from axionlike particles and implications for multimessenger observations of blazars[END_REF]) of the mass close to 1 eV.

We modelled the additional excess EBL component with the log-Gaussian spectral shape B λ (λ):

λB λ = B exp - log 2 (λ/µ) 2 log(1 + σ µ ) 2 ,
(2.30)

as shown in Fig. 2.7. This function appears as a Gaussian when plotted on a logarithmic scale (by analogy with log-parabola). It has three parameters: normalisation B, central wavelength µ, and width σ. Figure 2.7 shows that such a log-Gaussian can fit both the 'broad' excess at the level of the highest direct EBL measurements and the EBL excess reported by CIBER [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF]. For the former, we should take B = 30 nW/m 2 /sr and σ = 0.7 µm, while for the latter B = 15 nW/m 2 /sr and σ = 0.2 µm. Both bumps are positioned at µ = 1.5 µm.

To explore what kind of excess EBL flux could be accommodated to in the γ-ray data, we scanned over the height, width, and position of the bump. We also studied how the quality of the fits to the gamma-ray spectra improves (worsens) with changes in the parameters. At the first step, we defined a vector θ as a set of parameters θ = (B, µ, σ), and we calculated χ 2 i (θ) and the parameters Γ 1,i (θ), Γ 2,i (θ), and E br,i (θ) of the fit of the i-th blazar for every point θ in range 0 < B < 100 nW/m 2 /sr , 0.1 < µ < 10 µm, 0 < σ < 1.5 µm, sampling the parameter space at small steps δθ. Then, we imposed a cut on the hardness of the intrinsic blazar spectra and removed all the points θ where Γ 2,i < 1.5 for at least one blazar.

Thereafter, we derived the 2-σ excluded region of parameter space for each blazar. The point θ is excluded if χ 2 i (θ) > χ 2 i,0 + 8.02 where χ 2 i,0 corresponds to the best fit without additional excess and the value 8.02 corresponds to the 2-σ inconsistency level for the nested models with three added parameters. The 1-σ preferred regions are those where the condition χ 2 i (θ) < χ 2 i,min + 3.53 is satisfied and χ 2 i,min corresponds to the best fit of the model with the additional excess.

Results

The results are presented in Figure 2.9. The first plot is the cross-section of the three dimensional parameter space at µ = 1.5 µm. It shows an upper bound on the normalisation of the excess as a function of its spectral width. The strongest constraints come from the blazar 1ES 1218+304. In the limit σ > 1 µm, the width of the bump becomes large enough to correspond to the change of the overall normalisation of the EBL. In this case we have found an upper bound [START_REF] Abramowski | Measurement of the extragalactic background light imprint on the spectra of the brightest blazars observed with H.E.S.S[END_REF] and MAGIC [START_REF] Ahnen | MAGIC observations of the February 2014 flare of 1ES 1011+496 and ensuing constraint of the EBL density[END_REF].

B < 10 nW/(m 2 sr) (2.31)
The upper bound on the excess EBL flux normalisation is relaxed if the excess feature is narrow, with the width σ 1 µm. In particular, Fig. 2.9 shows that the narrow excess with the flux level comparable to the CIBER minimal EBL is not ruled out by the γ-ray observations while the CIBER nominal EBL is strongly excluded. In fact, by fitting the data of 1ES 1218+304 with a model EBL with a narrow bump, we find that the best fit is not achieved with a model without the bump, but with a bump which is stronger and narrower than CIBER minimal.

In Fig. 2.9 the central wavelength is fixed to show the dependence of constraint on normalisation as a function of the feature width. The source that imposes the strongest constraint is 1ES 1218+304. For this source, we show the dependence of the constraint on normalisation as a function of the central wavelength on right panel of Fig. 2.9. We fixed the width of the feature to σ = 0.2 µm, which is equal to the width of the spectral feature modelled based on the CIBER minimal EBL measurements.

The addition of the feature at the central wavelength µ 1.7 µm improves the quality of the fit if the normalisation is B 15 nW/(m 2 sr), compared to the nominal low EBL model. The best-fit 'bump' feature is shown superimposed on the low EBL model in Fig. 2.7. One could notice an agreement of the best fit bump feature with the measurements of CIBER. This agreement is also evident in Fig. 2.9 since the CIBER measurement is within the 68% confidence contour of the γ-ray measurement.

We verified that the hypothesis of existence of the bump at 1.7 micron in the EBL spectrum does not also contradict the data on other blazars. This is demonstrated in Fig. 2.8, where fits of the spectra of all blazars considered in the analysis are shown for the EBL model with and without the 1.5 micron bumps with green and black thin lines, respectively. Table 2.3 provides details on the quality of the fits with different models.

The narrow EBL bump produces a characteristic 'dip' feature in the γ-ray spectra shown in Fig. 2.8, which generally improves the quality of the fit. Hints of such dips are visible in the spectra of PKS 2155-304, PKS 2005-489, and in 1ES 1218+304, but not at the level of significance that is sufficient for detection. The dip features could also arise if the EBL has sharply falling intensity between the near-IR to mid-IR [START_REF] Orr | Strong New Constraints on the Extragalactic Background Light in the Near-to Mid-infrared[END_REF]. This model will lead to a slightly different shape in the attenuation of the gamma-ray flux as a function of energy and it should be possible to distinguish between the two possibilities once the statistics of the spectral measurements and control over systematics improve.

The effect of absorption of TeV γ-rays in interactions with the EBL bump photons increases proportionally to the distance of the source. Therefore the TeV range dip impinged on the γ-rayspectrum by the EBL bump should get more pronounced for further away sources. However, the overall suppression of the spectrum in the multi-TeV range makes the quality of the spectra that are collected with the current generation of telescopes insufficient for the measurement of the effect of the EBL bump on the spectra.

The situation will be improved upon the start of operations of the Cherenkov Telescope Array (CTA). For nearby blazars, such as PKS 2155-304, PKS 2005-489 and 1ES 1218+304, CTA will be able to provide high signal-to-noise measurements of the spectrum over sufficiently large dynamic ranges from tens of GeV up to 10 TeV. This will allow one to confirm or reject the hint of existence of the TeV dip in the γ-ray spectra (corresponding to the µm bump in the EBL spectrum with high significance). For further away sources, CTA will extend the dynamical range of the spectral measurements in the multi-TeV range. This would enable the confirmation and rejection of the EBL bump origin of the TeV dip in the spectra, that is, if the dip were to be systematically detected in further away sources with larger amplitude.

If the dip in the 1ES 1218+304 spectrum is really due to the presence of a bump in the EBL spectrum, the improvement of the quality of spectral measurements with CTA will lead to the discovery of the EBL spectral feature. We made an estimate of significance of detection that could be attained by CTA. To do this, we took current measurements of 1ES 1218+304 spectrum by VERITAS and divided the errorbars by the square root of the ratio of effective areas of CTA and VERITAS, which is by a factor of 3 (energy dependent) in the 0.3-3 TeV energy range. In this way we obtained an expected CTA-level quality of spectral measurements for the same exposure as in VERITAS (∼ 86 hours). We excluded the first two points of the spectrum that lie below the model. With such a 'pre-fabricated' CTA-quality spectrum, we repeated the analysis described in the previous sections and found that CTA will be able to help discover the spectral feature in the EBL spectrum and detect it with significance above 5σ. Bump-like features in the EBL spectrum could be produced by decaying particles [START_REF] Kohri | Can decaying particle explain cosmic infrared background excess?[END_REF], such as axion-like particles that possibly form part of dark matter [START_REF] Kalashev | Cosmic infrared background excess from axionlike particles and implications for multimessenger observations of blazars[END_REF], and by peculiar stellar populations that could have existed in the past but are absent today, such as Population III stars [START_REF] Michael | The contribution of the first stars to the cosmic infrared background[END_REF]. In this respect, the possibility of discovering spectral features in the EBL spectrum via precision measurements of the blazar spectra that would be enabled by CTA will be relevant in a range of astrophysics and cosmology in addition to fundamental physics contexts.

Conclusions of the Chapter

In this Chapter, a new EBL model is constructed, which makes it possible to modify the parameters underlying it. With the help of this model, using Markov chain Monte Carlo technique, we derived constraints on the astrophysical parameters responsible for the formation of EBL spectrum, including constraints on the star formation rate.

On the other hand, we investigated the possibility that there is an additional narrow spectral feature the EBL spectrum. Namely we showed that previously reported 'excess' in EBL flux in the ∼ 1 micron wavelength range is consistent with gamma-ray measurements, that is, if the excess has the form of a narrow feature of the width δλ < λ and an overall flux of up to 15 nW/(m 2 sr) above the 'minimal' EBL, which is estimated from the visible and infrared source counts. Such 'bump-like' spectral features could originate, for example, from decaying dark-matter particles, or either axions or peculiar astrophysical processes in the course of star-formation history. We also showed that if the additional component with parameters of CIBER minimal measurements actually exist, CTA will be able to detect it with significance above 5σ.

Chapter 3

Axion-like particles and gamma-ray emission from blazars

Blazars are a subclass of active galactic nuclei, with the jet oriented towards the observer. Blazar radiation has been detected in a wide range of energies starting from radio waves and up to the energies of the order of 10 TeV. The most energetic gamma rays produced by blazars can initiate electromagnetic cascades in intergalactic space. The development of the cascade begins with the interaction of a high-energy gamma-quantum with a background photon of the EBL, resulting in the creation of an electron-positron pair. All processes in the electromagnetic cascade are well described by quantum electrodynamics. In particular, quantum electrodynamics predicts the level of absorption of gamma rays (optical depth) as a function of energy and distance to the source.

Hypothetical axion-like particles may lead to deviations from the expected level of absorption. First, axion-like particles can decay into two photons. If the mass of the ALPs is of the order of 1 eV, then produced photons will contribute to the EBL spectrum, which will increase absorption. On the other hand, axion-like particles can convert into gamma rays in external magnetic fields. With some combination of the parameters of the ALPs, photon-ALP mixing can become significant and increase the effective mean free path of gamma rays. As a result, the observed absorption will be less than expected.

In this Chapter, we use data on high-energy gamma-ray emission from blazars to consider both variants of ALPs manifestation. The main results of this Chapter are published in [START_REF] Korochkin | Search for anomalous features in gamma-ray blazar spectra corrected for the absorption on the extragalactic background light[END_REF][START_REF] Korochkin | Search for decaying eV-mass axionlike particles using gamma-ray signal from blazars[END_REF][START_REF] Addazi | The Large High Altitude Air Shower Observatory (LHAASO) Science Book (2021 Edition)[END_REF].

Axion-like particles

Axion-like particles (ALPs) are hypothetical neutral light bosons with zero spin which arise in many extensions of the Standard Model of particle physics (see e.g. [START_REF] Jaeckel | The Low-Energy Frontier of Particle Physics[END_REF][START_REF] Graham | Experimental Searches for the Axion and Axion-Like Particles[END_REF] for reviews). They are the generalisation of the QCD axions which initially were invoked to solve the strong CP problem [START_REF] Peccei | CP Conservation in the Presence of Pseudoparticles[END_REF][START_REF] Weinberg | A New Light Boson?[END_REF][START_REF] Wilczek | Problem of Strong P and T Invariance in the Presence of Instantons[END_REF]. ALPs are characterised by the mass m and two-photon coupling constant g aγγ , the parameter space is shown in Fig. 3.1. These two parameters are expressed through one and the same (unknown) energy scale for the conventional axion, so that they are related as 1 . g aγγ 2 × 10 -10 ζ(m a /1 eV)) GeV -1 (ζ is a model-dependent parameter of the order of one) [START_REF] Kim | Weak-Interaction Singlet and Strong CP Invariance[END_REF][START_REF] Shifman | Can confinement ensure natural CP invariance of strong interactions?[END_REF][START_REF] Dine | A simple solution to the strong CP problem with a harmless axion[END_REF][START_REF] Zhitnitsky | On Possible Suppression of the Axion Hadron Interactions[END_REF]. An example of this relation for a particular Kim-Shifman-Vainshtein-Zakharov (KSVZ) model is shown by black line and yellow model uncertainty range in Fig. 3.1. This constraint is relaxed in generic ALP models, so that broad region of (m, g aγγ ) parameter space could be considered [START_REF] Jaeckel | The Low-Energy Frontier of Particle Physics[END_REF]. A wide range of ALP masses is discussed in the literature ranging from ultralight with m ∝ 10 -12 eV to heavy particles m ∝ 10 6 eV. Below we will focus on the most interesting for us range of masses m ≈ 0.1 -10 eV and describe known constraints on the coupling constant with photons for such ALPs. The Lagrangian describing ALP interactions with an electromagnetic field (F µν -electromagnetic tensor, Fµν -dual tensor) has the following form:

L = 1 2 (∂ µ a) 2 - 1 2 m 2 a 2 - 1 4 F 2 µν - 1 4 g aγγ aF µν F µν (3.1)
where a is the axion-like particle. The last term can be rewritten as:

- 1 4 g aγγ aF µν F µν = g aγγ a E • B (3.2)
The two-photon vertex interactions determine the preferred method of search for ALPs by their conversion to photons in the presence of magnetic fields [START_REF] Sikivie | Experimental Tests of the "Invisible" Axion[END_REF]. Nonobservation of this effect imposes constraints on the ALP parameters in the m-g aγγ parameter space, based on laboratory experiments and astronomical observations. Nonobservation of axion flux from the Sun by CAST experiment constrains g aγγ from above in the mass range below 1 eV (see Fig. 3.1) [START_REF] Arik | Search for Solar Axions by the CERN Axion Solar Telescope with 3 He Buffer Gas: Closing the Hot Dark Matter Gap[END_REF]. The CAST limit will be improved by the next-generation facility IAXO [START_REF]Physics potential of the International Axion Observatory (IAXO)[END_REF], also shown in Fig. 3

.1.

A strong constraint on the ALP parameters is imposed by non-observation of the effect of energy loss through emission of ALPs on stellar evolution [START_REF] Vysotsky | Some Astrophysical Limitations on Axion Mass[END_REF][START_REF] Raffelt | Astrophysical axion bounds diminished by screening effects[END_REF]. Recent constraints on two-photon coupling [START_REF] Ayala | Revisiting the Bound on Axion-Photon Coupling from Globular Clusters[END_REF][START_REF] Straniero | Axion-Photon Coupling: Astrophysical Constraints[END_REF] Yellow band and black solid line corresponds to the QCD axion models [START_REF] Kim | Weak-Interaction Singlet and Strong CP Invariance[END_REF][START_REF] Shifman | Can confinement ensure natural CP invariance of strong interactions?[END_REF][START_REF] Dine | A simple solution to the strong CP problem with a harmless axion[END_REF][START_REF] Zhitnitsky | On Possible Suppression of the Axion Hadron Interactions[END_REF]. Blue region is the limit imposed by CAST non-observation of solar axions [START_REF] Arik | Search for Solar Axions by the CERN Axion Solar Telescope with 3 He Buffer Gas: Closing the Hot Dark Matter Gap[END_REF]. Grey vertical band comes from high-resolution spectroscopy of intracluster medium in by optical telescope [START_REF] Grin | A Telescope Search for Decaying Relic Axions[END_REF]. Black dashed line is the EBL limit estimate of [START_REF] Arias | WISPy Cold Dark Matter[END_REF]. Green horizontal line is the limit imposed by non-observation of excess energy loss in Horizontal Branch stars [START_REF] Straniero | Axion-Photon Coupling: Astrophysical Constraints[END_REF]. Light grey shaded region shows 95% confidence level excluded range of parameters from non-observation of absorption feature in the Fermi/LAT+VERITAS spectrum of 1ES 1218+304 (discussed in this paper). Violet shaded ellipses show 68% confidence level preferred range of the parameters obtained by the spectral analysis of 1ES 1218+304. Red hatched region shows the favoured region found if the first spectral point of VERITAS is not taken into account. assumption of validity of complex stellar evolution models controlled by large number of parameters which are not directly measured [START_REF] Sandage | The Anomalous Color-Magnitude Diagram of the Remote Globular Cluster NGC 7006[END_REF][START_REF] Gratton | The second and third parameters of the horizontal branch in globular clusters[END_REF].

Direct searches for ALP decays into two photons were performed using high-resolution spectroscopy of galaxy clusters in the range of masses between 4.5 and 7.7 eV [START_REF] Grin | A Telescope Search for Decaying Relic Axions[END_REF]. Non-observation of non-identified spectral emission lines from ALP decays imposes a constraint on g aγγ shown by grey shading in Fig. 3.1.

In the next section we consider constraints on parameters of ALPs based on observations of blazars.

Effect of the ALP decay EBL spectral feature on γ-ray spectra

Line emission from ALPs decaying everywhere across the Universe contributes to the extragalactic background light (EBL) [START_REF] Kalashev | Cosmic infrared background excess from axionlike particles and implications for multimessenger observations of blazars[END_REF]. Superposition of the line emission accumulated at different redshifts leads to a broad "bump-like" feature in the EBL spectrum with the width ∆E ∼ E.

In the previous Chapter we derived constraints on the parameters on the additional narrow spectral feature in the EBL using high-energy observations of blazars. In what follows we discuss implications of this result for the ALP dark matter model. We develop the line of reasoning of [START_REF] Berezhiani | Cosmic Nonthermal Electromagnetic Background from Axion Decays in the Models with Low Scale of Family Symmetry Breaking[END_REF] and [START_REF] Arias | WISPy Cold Dark Matter[END_REF] on constraints on the ALP two-photon coupling imposed by non-observation of the features in the EBL spectrum with a qualitative analysis (the limits on g aγγ derived by [START_REF] Arias | WISPy Cold Dark Matter[END_REF] have order-of-magnitude nature). The analysis of previous Chapter has favoured non-zero normalisation of the EBL excess at ∼ 1 micron wavelength. This could be directly expressed in terms of the ALP model of the EBL bump feature. The favoured bump wavelength and normalization parameters could be converted into the favoured ALP coupling and mass range. We discuss a possibility of scrutinising this result with deeper observations of a set of sources with current and next generation γ-ray telescopes.

EBL spectral feature produced by decaying ALPs

Dark matter decays at redshift z produce signal at wavelengths λ = λ 0 (1+z), where longer than λ 0 = 4π/m is the wavelength of ALP decay line at the redshift z = 0. We assume monochromatic decay line for dark matter with wavelength λ 0 :

F λ (λ) = 2 hc λ 2 0 Γδ λ λ 0 -1 (3.3)
where Γ is the decay rate into two photons

Γ = g 2 aγγ m 3 64π (3.4)
The flux generated by dark matter decays is [START_REF] Boyarsky | Constraints on sterile neutrinos as dark matter candidates from the diffuse X-ray background[END_REF]

I λ (λ) = c 4π zmax 0 n(z)F λ λ 1+z (1 + z) 3 H(z) dz (3.5)
where n is the comoving density of ALPs, 3 , H 0 is the Hubble constant, c is the speed of light and F λ (λ) is the spectral energy distribution of photons from decaying particles. Ω λ and Ω m are the matter and dark energy density parameters (we use the values Ω Λ = 0.691, Ω m = 0.309, H 0 = 67.8 km/s/Mpc derived by [START_REF]Planck 2018 results. VI. Cosmological parameters[END_REF]).

H(z) = H 0 Ω Λ + Ω m (1 + z)
If the decay time is much longer than the age of the Universe, or Γ H 0 , n(z) const and the integral in Eq. (3.5) could be explicitly taken. This gives the spectral shape of the dark matter decay continuum at a function of λ:

λI λ (λ) = c 256π 2 g 2 aγγ m 3 ρ c f Ω dm λm 4π H 0 Ω λ + Ω m λm 4π 3 (3.6)
for λ ≥ λ 0 and zero for λ < λ 0 . Here ρ = ρ c f Ω dm is the dark matter energy density, ρ c is the critical energy density, Ω dm Ω m dark matter density parameter and f is the fraction of dark matter in the form of ALPs. Fig. 3.2 shows a compilation of the measurements of the EBL flux using different techniques. The discrepancy between the direct measurements and estimates of the EBL flux based on γ-ray observations is evident particularly in the wavelength range 1-3 µm.

Possible additional ALP delay flux contributions to the EBL, given by Eq. 3.6 are shown by the red shading in Fig. 3.2 for the ALP mass and two-photon coupling range which is consistent with the existing γ-ray data (as detailed in the following sections). One could see that this contribution is also consistent with the direct measurements.

The shape of the spectral feature in the EBL produced by ALP decays is different from a generic log-Gaussian bump shape considered in previous Chapter [START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF]. Nevertheless, it is a feature of the width ∆λ ∼ λ, so that generic conclusions of analysis of previous Chapter [START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF] apply also for the ALP decay spectral feature.

In particular, in previous Chapter [START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF] we found that the blazar 1ES 1218 + 304 provides the best constraints, due to the availability of deep VERITAS exposure of the source and high source flux in Fermi/LAT energy band. Spectral analysis of this source limits the amplitude of possible bumps in the EBL spectrum to be at most at the level of 10-15 nW/m 2 sr on top of the EBL flux produced by the star formation process (also at the level of ∼ 10 nW/m 2 sr in the micron wavelength range). Considering this fact, we refine the analysis of [START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF] for the spectrum of 1ES 1218 + 304 in the specific case of ALP induced EBL spectral feature, by re-fitting the spectrum with account of attenuation of the γ-ray flux by the spectral feature with the shape described by Eq. (3.6).

The best fit model is shown on the left panel of Fig. 3.3. The intrinsic source spectrum is modelled as a broken powerlaw(s) shown by the dashed lines. The energy of the break found in the spectral fit depends on the assumptions about the EBL model. If the baseline EBL model without additional spectral feature produced by the ALP decay is used, the break is just below 100 GeV energy and the intrinsic spectrum above the break is nearly dN/dE ∝ E -2 . Inclusion of an additional spectral feature in the EBL at ∼ 1µm shifts the break energy up to 300 GeV and favours softer spectrum above the break (green dashed line). An additional feature in the EBL at ∼ 5µm shifts the break toward lower energy and favours harder spectrum above the break (magenta dashed line).

Overall, the spectral fit to 1ES 1218+304 spectrum prefers the EBL model with nonzero bump at approximately one micron wavelength (the effect is at significance level below 3σ and we could not claim an "evidence" in this sense). However, the preferred region is close to the favoured range of parameters, found in [START_REF] Daido | The ALP miracle: unified inflaton and dark matter[END_REF][START_REF] Daido | The ALP miracle revisited[END_REF] from dark matter abundance and cosmic microwave background data. The spectral fit corresponding to this best-fit model is shown by the green solid line on the left panel of Fig. 3.3.

If the constraints on the amplitude of possible bumps in the EBL spectrum are interpreted in terms of the model of decaying ALPs, they could be converted into constraints on the mass-dependent ALP particle coupling to photons using Eq. (3.6). The result is shown in Fig. 3.1. The fact that the fit with non-zero bump feature is favoured by the data has two effects on the ALP coupling constraint. First, it relaxes the overall upper bound on g aγγ in the ALP mass range 1-10 eV. This is clear from comparison of the 95% confidence level upper limit derived from the 1ES 1218+304 analysis (grey region in Fig. 3.1) with the order-of-magnitude estimate of constraint on g aγγ from absence of features in the EBL derived by [START_REF] Arias | WISPy Cold Dark Matter[END_REF] (black dashed line). Next, a preferred range (68% confidence level) of ALP mass-coupling parameter space (violet shaded ellipses) appears. It spans a wide range of mass between 0.2 and 20 eV and the coupling range 3 × 10 -11 GeV -1 < g aγγ < 2 × 10 -10 GeV -1 , just at the border of the stellar evolution constraint shown by the red horizontal line.

The favoured region of the parameter space has two domains: a wider ellipse around the central value m 3 eV and g aγγ 10 -10 GeV -1 and a more narrow ellipse stretching up to 20 eV in the ALP mass. This structure of the favoured region is largely determined by the presence of features in the spectrum of 1ES 1218+304 visible in Fig. 3.3. The noticeable curvature of the spectrum in 0.3-1 TeV energy range favours the fit with an ALP-induced dip in the spectrum for the ALP mass m 3 eV. Curvature of the VERITAS spectrum close of its low-energy end in the 0.1-0.3 TeV range favours the presence of the dip in this energy range. This dip would be produced by ALPs with masses in 10 eV range. [START_REF] Tsumura | Low-Resolution Spectrum of the Extragalactic Background Light with AKARI InfraRed Camera[END_REF], while green asterisks are from the reanalysis of IRTS [START_REF] Matsumoto | Reanalysis of the near-infrared extragalactic background light based on the IRTS Observations[END_REF]. Pink diamonds are from reanalysis of COBE/DIRBE measurements [START_REF] Sano | Derivation of a large isotropic diffuse sky emission component at 1.25 and 2.2um from the COBE/DIRBE data[END_REF][START_REF] Sano | Measurements of Diffuse Sky Emission Components in High Galactic Latitudes at 3.5 and 4.9 µm Using DIRBE and W ISE Data[END_REF]. Black data points together with blue and orange systematic uncertainty are derived from CIBER [START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF] and correspond to nominal and minimum EBL models. Lower limits: dark blue upward arrows combine EBL lower limits obtained by different experiments: GALEX [START_REF] Xu | Number counts of GALEX sources in FUV (1530A) and NUV (2310A) bands[END_REF], Hubble Deep Field [START_REF] Madau | Deep galaxy counts, extragalactic background light, and the stellar baryon budget[END_REF], Subaru [START_REF] Keenan | The Resolved Near-Infrared Extragalactic Background[END_REF], Spitzer/IRAC [START_REF] Fazio | Number counts at 3 < lambda < 10 um from the Spitzer Space Telescope[END_REF]. EBL from γ-ray absorption: striped lime, red and purple bands are from MAGIC [START_REF] Acciari | Measurement of the extragalactic background light using MAGIC and Fermi-LAT gamma-ray observations of blazars up to z = 1[END_REF], HESS [START_REF] Abramowski | Measurement of the extragalactic background light imprint on the spectra of the brightest blazars observed with H.E.S.S[END_REF] and Fermi/LAT [START_REF] Abdollahi | A gamma-ray determination of the Universe's star formation history[END_REF] correspondingly. Models: dashed blue line is for baseline EBL model of [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. Black striped band shows the allowed range of EBL models obtained with the global fit to EBL by [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF]. Red shaded region shows the favoured range of additional component due to the ALP decay with parameters in the red shaded region of Fig. 3 It is possible that the VERITAS spectral measurement at the low-energy end is affected by the systematic effects. To illustrate strong dependence of the shape of the favoured range of parameter space on such systematic effects, we have re-done the spectral analysis after removing the first data point of VERITAS (shown in grey colour on the left panel of Fig. 3.3). The result is shown in Fig. 3.1 by the red-hatched region. One could see that removal of the spectral measurement possibly affected by the systematic effect changes the shape of the favoured region: the two ellipses merge into a single favoured region.

Possible improvement of ALP search sensitivity using dedicated observations with the CTA telescope

The effect of removal of one spectral measurement point discussed in the previous section shows the limitation of our analysis, which is dominated by the measurement of single source with single Cherenkov telescope. Better constraints on ALP coupling or, ultimately, detection of the ALP induced bump in the EBL spectrum has to be done using observations of multiple sources with multiple telescopes. This is illustrated in Fig. 3.4 where we show the sensitivity reach of observations with MAGIC telescope which have an advantage of lower energy threshold compared to VERITAS. To estimate MAGIC sensitivity for the search of the ALP-induced ab-sorption feature we have simulated the spectral points for a 50 hr MAGIC exposure using the information on the energy dependence of the effective area and background rate found by [START_REF]The major upgrade of the MAGIC telescopes, Part II: A performance study using observations of the Crab Nebula[END_REF]. We have then calculated the signal and background statistics and estimated the statistical and systematic errors for each spectral data point based on this information. For the systematic error we have assumed that it is at the level of 5% of the background flux. The MAGIC spectrum simulated in this way is shown on the right panel Fig. 3.3.

As it is discussed in the previous section, the spectral fits of Fermi/LAT+VERITAS data suffer from uncertainty of the energy of the break in the intrinsic spectrum, which is in the 100 GeV range, but is not well constrained by the data. Lower energy threshold of MAGIC observations would provide a possibility to better constrain the break energy in the 100 GeV range.

Still further improvement will be achieved with the CTA North array. In order to estimate the sensitivity reach of CTA we have followed the same procedure as for the MAGIC spectral points simulations. We have used the information on the effective area, and background rate within the angular resolution expected for CTA2 , and simulated the expected spectrum of 1ES 1218+304 for a 50 hr observation with CTA North array. The result is presented on the right panel of Fig. 3.3. We have chosen the spectral model without an additional ALP induced absorption feature in the spectrum. We have then applied the same analysis as for the VERITAS+Fermi/LAT data fitting the spectrum with spectral model which includes the ALP induced absorption feature, to find the maximal possible normalisation of the ALP-induced bump in the EBL spectrum as a function of the ALP particle mass. We have then re-calculated the upper limit on the norm of the ALP-induced bump into the upper limit on g aγγ , to find the sensitivity reach of CTA. The result is shown in Fig. 3.4. From this figure one could see that observations of CTA will provide a sensitivity improvement by a factor of three to five in the ALP mass range around between 0.3 and 10 eV.

We have repeated the simulation of the 1ES 1218+304 spectrum assuming the EBL model with an additional ALP induced bump at the level of the best fit model of VERITAS+Fermi/LAT spectrum to explore the significance level of detection of the ALP-induced spectral bump in the EBL with the flux consistent with the estimates from the direct EBL measurement experiments [START_REF] Wright | Detection of the cosmic infrared background at 2.2 and 3.5 microns using dirbe observations[END_REF][START_REF] Matsumoto | IRTS observation of the near-infrared extragalactic background light[END_REF][START_REF] Tsumura | Low-Resolution Spectrum of the Extragalactic Background Light with AKARI InfraRed Camera[END_REF][START_REF] Matsumoto | Reanalysis of the near-infrared extragalactic background light based on the IRTS Observations[END_REF][START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF]. We find that if the flux of the ALP induced bump in the EBL spectrum reaches ∼ 10 nW/(m 2 s sr), CTA will be able to detect the associated dip in the γ-ray spectrum at significance level 4σ. Fig. 3.4 also shows the precision with which detection of the ALP induced absorption feature with CTA would allow to measure the ALP particle parameters g aγγ and m. The error ellipse of such measurements is shown by red shading.

Possible improvement of ALP search sensitivity using dedicated observations with the LHAASO observatory

To estimate the sensitivity reach of LHAASO for the ALPs parameters we performed the same analysis as in previous section for the CTA. Six sources, which satisfy our the criteria and located in LHAASO field of view were selected: Markarian 501, Markarian 421, 1ES 1218 + 304, 1ES 1215 + 303, 1ES 1959 + 650, 1ES 1011 + 496. Using gamma-ray sensitivity, effective area and gamma/hadron separation efficiency we estimated background rate and simulated spectra for the blazars in the sample for 1-and 5-years of observations with LHAASO in the energy range between 200 GeV and up to 20 TeV. One can see that 5-years LHAASO sensitivity is sufficient to explore potentially interesting area of ALPs masses around 1 eV and g aγγ ∼ 10 -10 GeV -1 indicated from previous observations. In the case of absence of the signal LHAASO will improve current limitation on g aγγ ∼ 10 -10 GeV -1 by a factor of 3-5 in the mass range from 1 eV to 5 eV. From the right panel of Fig. 3.4 one can see that LHAASO will have enough sensitivity to cross-check the CTA measurement of the ALP parameters (if the existence of a bump in the EBL spectrum will be confirmed).

Search for anomalous features in gamma-ray blazar spectra

In the previous Chapter, we explored the possibility of matching direct measurements of EBL with the constraints from blazar observations. At the same time, we have always assumed that the galaxy count method provides a robust lower bound on EBL intensity. However, recent studies have indicated that the universe may be more transparent to high-energy photons than suggested from the galaxy count lower bound.

In particular [START_REF] Rubtsov | Breaks in gamma-ray spectra of distant blazars and transparency of the Universe[END_REF] and [START_REF] Horns | Indications for a pair-production anomaly from the propagation of VHE gamma-rays[END_REF] reported the presence of non-physical upward breaks in the deabsorbed blazar spectra. In the study of [START_REF] Rubtsov | Breaks in gamma-ray spectra of distant blazars and transparency of the Universe[END_REF] it was also found a dependence of the break strength on the distance to the source: stronger breaks were observed for more distant sources. The position of the break coincided with the energy above which the absorption became significant. Breaks have been observed even for the most conservative EBL models, which correspond to the lower limits from the galaxy count method. Such features in the spectra of blazars may indicate the existence of ALPs.

Indeed, under certain parameters of the ALPs, it is possible that some of the highenergy photons convert into an ALPs in a strong magnetic field near the source. Then ALPs propagate through the intergalactic medium without absorption and are converted back into gamma rays in the magnetic field of the Galaxy. If the intergalactic magnetic fields are strong enough (> 10 -9 G), then the gamma rays will convert into the ALPs and back during the whole path to the observer. Both scenarios lead to a decrease in the absorption of high-energy photons. In this case, using only the EBL In this section, we examine the effect of the "anomalous transparency" of the Universe using the most complete sample of blazars with accurately measured redshifts.

The data and the sample

The purpose of the present study implies the use of high-quality gamma-ray spectra of sources located at large, confidently known distances. The obvious candidates are blazars with firmly measured redshifts. For various redshifts z, the absorption due to e + e -pair production becomes important at different energies of photons: for more distant sources, pair production affects the spectra at lower energies3 . The benchmark energy, at which the absorption becomes important, corresponds to the optical depth τ = 1 and is denoted as E 0 hereafter,

τ (E 0 , z) = 1, (3.7)
which determines E 0 (z) for a given absorption model, see Fig. 3.5.

To study the behaviour of a spectrum in the region of strong absorption, one needs observations at E E 0 ; therefore, the use of different instruments is required for objects located at different redshifts: for large z, E 0 is typically dozens of GeV, and the relevant instrument is Fermi LAT, while for less distant objects, E 0 is of order of several hundred GeV, and the data of IACTs should be used. Though the data are published in quite different ways, it is essential to treat them on equal footing. Here we discuss our general requirements governing selection of the objects, which will be implemented in the next two subsections for the Fermi-LAT and IACT data. Redshift criteria. Incorrectly determined or uncertain redshifts may hurt any of blazar studies. This happens quite often, especially for BLLs whose spectra do not possess strong emission lines. We require a firm spectroscopic redshift for the objects included into our sample. More specifically, we start from a preselected sample of gamma-ray blazars and check, for every object, the relevant redshift information and references in the NASA/IPAC Extragalactic Database (NED) 4 . The redshifts are accepted if they satisfy (R1) and one of (R2A), (R2B), (R2C) criteria:

(R1) the redshift is spectroscopic (not photometric) and is derived from emission lines AND EITHER (R2A) the redshift is determined in a dedicated study (if several dedicated studies give different results, the latest one is used), OR (R2B) the redshift quoted in NED is determined in 2dF Galaxy Redshift Survey [START_REF] Colless | The 2dF Galaxy Redshift Survey: Spectra and redshifts[END_REF] or 6dF Galaxy Survey [START_REF] Jones | The 6dF Galaxy Survey: Final Redshift Release (DR3) and Southern Large-Scale Structures[END_REF],

OR (R2C) the redshift quoted in NED is determined in the Sloan Digital Sky Survey (SDSS), the result is unique and does not change from one release to another. Spectral criteria. In our work, we use binned gamma-ray spectra, because only this information is publicly available for the objects observed by IACTs. We impose the following criteria for the spectra:

(S1) information for at least 5 energy bins is available; (S2) the lowest energy E last of the last spectral bin satisfies E last > E 0 AND the lowest energy of the third from below spectral bin is below E 0 .

These criteria are aimed at the reconstruction of spectra at E both below and above E 0 with reasonable accuracy. Criterion (S2) depends on the assumed absorption model, and therefore the sample is slightly different one from another for different EBL models.

IACT subsample To construct the sample of blazars observed by IACTs, we start with the database consolidating results published by different observatories, the TeV-Cat online source catalog [START_REF] Wakely | TeVCat: An online catalog for Very High Energy Gamma-Ray Astronomy[END_REF] 5 . The preselected sample includes objects classified as blazars there (classes "HBL", "'IBL", "LBL", "BL Lac (class uncertain)", "FSRQ" and "blazar"), 68 objects in total as of 2018. Then, for each object in the sample, we checked the redshift selection rules (R1) and ((R2A) or (R2B) or (R2C)), using refer-4. Available at http://ned.ipac.caltech.edu. 5. Available at http://tevcat.uchicago.edu. ences from NED, and availability of gamma-ray spectra satisfying (S1), using references from TeVCat. Some objects have been observed several times, often by multiple instruments. To avoid double counting of information, which may artificially increase or dilute observed effects, we use only one observation for each source. It is chosen on the basis of better statistics (larger effective exposure in cm 2 •s). Because of strong variability of many blazars, we never combine any two spectra from observations performed at different epochs to enlarge the energy coverage. We also never combine Fermi-LAT spectra with those obtained from IACT observations because of potential systematic differences between the instruments' energy scales.

The next step requires to fix the absorption model. We use the EBL model constructed in Chapter 1 [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF] as our baseline low-absorption model. when variations of the model are considered, this step is repeated. At this step, we calculate the value of E 0 for every source and check the condition (S2). We are left with 26 objects which, together with the sources observed by Fermi LAT, are listed in Table 3.1.

For blazars observed by Fermi LAT, our starting point is the 3rd Fermi-LAT Catalog of High-Energy Sources, 3FHL [START_REF] Ajello | 3FHL: The Third Catalog of Hard Fermi-LAT Sources[END_REF] 6 , presenting a list of sources detected confidently above 10 GeV. From the catalog, we select 1212 objects classified as blazars (classes "BLL", "bll", "FSRQ", "fsrq", "bcu"). Since the redshifts presented in the catalog may be doubtful or erroneous for distant sources, we check them manually to satisfy (R1) and ((R2A) or (R2B) or (R2C)) criteria, using references from NED, and apply the additional constraint z ≥ 0.2, because for less distant sources, E 0 500 GeV is beyond the limit of reasonable sensitivity of Fermi LAT. This procedure results in a list of 309 blazars. For each of them, we use publicly available Fermi-LAT [START_REF] Atwood | The Large Area Telescope on the Fermi Gamma-ray Space Telescope Mission[END_REF] data 7 to construct binned spectra satisfying (S1) automatically with the help of standard routines from Fermi Science Tools 1.0.1. To this end, we use Fermi Pass 8 [START_REF] Atwood | Pass 8: Toward the Full Realization of the Fermi-LAT Scientific Potential[END_REF] Release 3 data (version 303). We use "SOURCE" class events with version 2 instrumental response functions "P8R3 SOURCE V2", recorded between August 4, 2008, and February 26, 2018 (Mission Elapsed Time interval 239557417 -541338374). The background models used are "gll iem v06.fits" for the Galactic component and "iso P8R3 SOURCE V2.txt" for the isotropic one. The background model also includes other sources from 3FGL catalog within 5 • angular distance from the given source. For most of the sources (305 of 309), the coordinates of the source are taken from the 3rd Fermi-LAT Source Catalog, 3FGL [START_REF] Acero | Fermi Large Area Telescope Third Source Catalog[END_REF]. However, not all 3FHL sources have 3FGL counterparts, and for the remaining 4 sources, we used 3FHL coordinates.

Firstly, we have performed an independent fit with the standard gtlike routine with the power-law spectrum model "PowerLaw2" in each of 8 energy bins with the width of 0.3 dex, for E ≥ 2 GeV. The lower energy limit of 2 GeV was chosen to cut 6. Available at http://fermi.gsfc.nasa.gov/ssc/data/access/lat/3FHL/ . 7. Available at https://heasarc.gsfc.nasa.gov/FTP/fermi/data/lat/weekly/photon/ .

the inverse-Compton peak energy range in order to improve broken power-law fits used in our subsequent analysis.

Despite the fact that the use of gtlike is the preferred method of constructing sources' spectra, it demonstrates an unstable behavior when the estimated number of counts in the energy bin is lower than one. In particular, the fit may not converge or converge to zero flux with unphysical value of the flux upper limit. Moreover, even when the number of photons is less than three, the estimated confidence interval for the flux follows the simple symmetric Gaussian distribution but not the Poisson one, as it should. Clearly these features may affect further analysis and to avoid the peculiar behavior of the fit optimization procedure we have stepped one level back to use the Fermi Science Tools routines which provide an input for gtlike. Namely, the standard gtsrcprob routine assigned the weight to each photon which is the probability to originate from a given source within the given source model. The sum of the corresponding weights provides an estimate for the number of photons from the given source and from each of the background sources.

At the next step, we fix the absorption model and apply the criterion (S2) to the obtained spectra. For the Fermi-LAT subsample, the criterion is used as follows: the source is accepted if it has at least one photon in the bin above E 0 and the probability that this photon is associated with the source is more than 99% (according to the gtsrcprob output). Only 5 of 309 blazars satisfy (S2) and enter our final sample.

Finally, for the bins which contain less than three photons we have replaced the gtlike flux with one obtained from summing weights from gtsrcprob. Since the background photons are considered separately, this sum may be treated as the number of signal events with zero background. We round this number to an integer (0 or 1 in most cases) and estimate the Poisson-based 68% CL interval for the flux in these bins by division by the exposure calculated with the gtexpcube routine. We have checked that the above procedure results in the value which agrees within 10% with the output of gtlike is cases when the latter converges to a physical value.

The list of 31 objects selected in the way described above, is presented in Table 3.1. For IACT observations, references for binned spectra used in this work are given there. We also list classes of the blazars (BLL or FSRQ) as determined in the catalogs used to compile the samples, TeVCat and 3FHL. In addition, we give also references to original spectroscopic measurements of redshifts.

Figure 3.6 represents the distribution of blazars in our main sample in redshift. One may note that, though BLLs dominate at short and moderate distances and FSRQs are in general farther away, both classes cover large ranges of redshifts. Hatched and shaded histograms give separate distributions for BLLs and FSRQs, respectively.

Data analysis and results

The key point of our study is the use of blazar spectra corrected for the pair production in a minimal-absorption model. Since, for IACT observations, only binned spectra are available, we adopt our analysis procedure for this case even for Fermi LAT, in order to process all observations in a uniform way. In contrast with many previous studies, including our work [START_REF] Rubtsov | Breaks in gamma-ray spectra of distant blazars and transparency of the Universe[END_REF], we do not use bin-by-bin deabsorption because it can introduce systematic biases in the highest-energy bins for the following reasons. The optical depth τ is a function of the incoming photon energy E and the source redshift z. It is important to note that the absorption is not uniform within the bin: photons of higher energy are absorbed stronger. For the energies and distances at which the absorption is significant, this effect may affect strongly the result. It can be taken into account if the shape of the observed spectrum within the bin is known, as we have done in Ref. [START_REF] Rubtsov | Breaks in gamma-ray spectra of distant blazars and transparency of the Universe[END_REF]; however, in practice, statistical uncertainties in the highest energy bins often make the determination of the shape uncertain. The uncertainty of the spectral shape within the bin translates into the error of deabsorbed flux. The latter effect is enhanced in bins, within which the absorption effect grows dramatically. It is shown that for wide energy bins of 0.3 dex used for the analysis of Fermi LAT sources, the uncertainty related to deabsorption contributes significantly to the flux error in energy bins above E 0 . Therefore the absorption of the model spectra instead of the deabsorption of the observed one is performed in our analysis.

Following the logic of previous studies, we start from the assumption about the 0.1 1.0 z

E, GeV

Figure 3.7 -Positions of the break energies (red points: hardenings, green asterisks: softening), treated as free parameters in independent fits, compared to E 0 (z) (blue line). Objects for which the broken power-law fit is not better than the power-law fit by ∆χ 2 < 1, that is for which best-fit break is consistent with zero, are shown as gray diamonds (since the error bars are determined by ∆χ 2 = 1 [START_REF] Press | Numerical recipes in C : the art of scientific computing[END_REF], they are infinite for these points. intrinsic spectrum of a blazar, for which we use a broken power law. First, we make the break energy E b a free parameter of the fit, so that it is adjusted independently for every source. We account for the absorption with the selected EBL model, integrate the obtained spectrum over the energy bins and compare resulted bin-by-bin fluxes with the observed data points. In this way, we fit the data with the four parameters of the intrinsic spectrum, that is the overall normalization, the break position E b , the spectral index at E < E b and ∆Γ, the difference between spectral indices at E > E b and E < E b . Technically, all our fits are performed with the method described in detail in Ref. [START_REF] Press | Numerical recipes in C : the art of scientific computing[END_REF]. Figure 3.7 compares the fitted break positions E b with the value of E 0 (z); we see that they agree well to each other: averaged over the sample, log (E b /E 0 (z)) = 0.04 ± 0.16 (for 18 sources with positive breaks, i.e. spectral hardenings, we get 0.06 ± 0.17).

Then, we turn to similar fits with the break energy fixed at E b = E 0 (z). Recall that this energy depends on the redshift of the source z, cf. Fig. 3.5, and therefore we assume breaks at different energies for different sources. By construction, nonzero ∆Γ in this case would indicate the presence of distance-dependent spectral features because the assumed break position E 0 (z) is redshift-dependent.

One sees that while some spectra have ∆Γ consistent with zero, indications for ∆Γ > 0 exist in many cases. For the entire sample, the assumption of ∆Γ = 0 results in χ 2 = 43.98 for 31 degrees of freedom, corresponding to the probability P 0.06 for the observed values of ∆Γ at E 0 (z) to occur as a result of a random fluctuation. Therefore, our results disagree with the hypothesis of the absence of breaks at E 0 (z) at 1.9 standard deviations, which is significantly weaker than previously reported. Using the EBL model of Gilmore et al [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF] further weakens the result to 1.3σ. 8 .

We can compare the assumption of breaks at E 0 (z) with that of breaks at some distance-independent energy, common to all sources. To this end, we perform fits of all spectral data by absorbed broken power laws, fixing the common break energy E b for all objects in the sample, for various values of E b . For every E b , we determine the goodness of the fit, p, from the chi-squared distribution with the corresponding number of degrees of freedom (for 194 data points, 27 spectra fitted with 3 parameters each, the number of degrees of freedom is 113). The plot of (1 -p) as a function of E b is shown in Fig. 3.9 (the best possible fit would correspond to (1 -p) = 0, a very bad fit -to (1 -p) ≈ 1), together with the same quantity calculated for fits with breaks at E 0 (z). The fit with breaks at distance-dependent energies E 0 (z) is better than a fit with any fixed E b except for a narrow range E b ∼ (200 -300) GeV; we note that this is expected since many of the sources in our sample have redshifts in the range such that E 0 (z) is close to this value.

8. For these probability estimates, we assume the chi-squared distribution. The latter is valid for Gaussian distributions and can be used here as many independent factors contribute to uncertainties.

Conclusions of the Chapter

In the current Chapter considered two possible scenarios for the manifestation of ALPs in the data on blazar observations. First, we show that ALPs with a mass in the region of 1 eV can contribute to the EBL spectrum if they constitute dark matter. In such a case, parameters of axion-like particles can be constrained with gamma-ray observations. Using available data on high energy observations of blazars, we constrain ALP coupling constant to photons and estimate how these constraints can be improved by the next generation observatories CTA and LHAASO.

On the other hand we investigated the problem of the "anomalous transparency" of the Universe for high-energy photons. To do this, we have built the most complete sample of blazars to date with accurately measured redshifts. Using this sample, we carefully studied the appearance of upward breaks in the deabsorbed spectra. We found that the breaks in the spectra are much less significant than previously reported. Thus the significance of the "anomalous transparency" effect does not exceed 2σ for the modern EBL models. The main results of this Chapter are published in [START_REF] Korochkin | Search for anomalous features in gamma-ray blazar spectra corrected for the absorption on the extragalactic background light[END_REF][START_REF] Korochkin | Search for decaying eV-mass axionlike particles using gamma-ray signal from blazars[END_REF][START_REF] Addazi | The Large High Altitude Air Shower Observatory (LHAASO) Science Book (2021 Edition)[END_REF].

Chapter 4

Comparison of Monte Carlo codes for electromagnetic cascades simulations

Absorption of the VHE γ-rays on EBL results in production of electron-positron pairs in the intergalactic medium. These pairs loose energy on inverse Compton scattering of CMB photons thus producing secondary γ-ray emission detectable by γ-ray telescopes. Observational appearance of this secondary emission depends on the strength and correlation length of magnetic field in the intergalactic medium and on the energy range of the secondary γ-rays.

Sensitivity of observations of the extended or delayed secondary emission from the electron-positron pairs in the intergalactic medium and of the measurements of IGMF will be crucially improved by the new generation telescope CTA [START_REF]Science with the Cherenkov Telescope Array[END_REF]. This will enable search for the extended secondary emission in previously inaccessible energy range thus extending the range of IGMF strengths that can be probed by the γ-ray technique. The sensitivity improvement will result in more precise characterization of the effect of absorption of the primary γ-ray flux by the EBL and more precise measurements of the primary γ-ray source spectra as well as in extension of the measurements into larger redshift range [START_REF] Abdalla | Sensitivity of the Cherenkov Telescope Array for probing cosmology and fundamental physics with gamma-ray propagation[END_REF].

Improved sensitivity for the secondary γ-ray flux from extragalactic sources provided by CTA has to be matched by the improved precision of modeling of this flux. Several numerical Monte Carlo codes have been previously developed for such modeling, including CRPropa [START_REF] Alves | CRPropa 3-a public astrophysical simulation framework for propagating extraterrestrial ultra-high energy particles[END_REF], ELMAG [START_REF] Blytt | ELMAG 3.01: A three-dimensional Monte Carlo simulation of electromagnetic cascades on the extragalactic background light and in magnetic fields[END_REF] and CRbeam [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF]. In current Chapter we compare precision of modeling of electromagnetic cascades in intergalactic space by these publicly available codes.

In all further analysis, we used the following versions of the codes: CRbeam 1.0 [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF], CRPropa3-3.1.7 [START_REF] Alves | CRPropa 3-a public astrophysical simulation framework for propagating extraterrestrial ultra-high energy particles[END_REF] and ELMAG 3.01 [START_REF] Blytt | ELMAG 3.01: A three-dimensional Monte Carlo simulation of electromagnetic cascades on the extragalactic background light and in magnetic fields[END_REF]. Additionally, we have verified that the results of CRPropa are stable in relation to previous versions CRPropa3-3.1.5 and CRPropa3-3.1.6. The results of the comparison are available in the electronic preprint [START_REF] Kalashev | Modelling of propagation of very-high-energy gamma rays with CRbeam code. Comparison with CRPropa and ELMAG codes[END_REF].

Physics of electromagnetic cascades in the intergalactic medium

High-energy gamma rays propagating through the intergalactic medium (IGM) from sources at cosmological distances can interact with the cosmic background radiation γ b , as well as with free electrons of the intergalactic medium e b . These interaction processes are well described by the quantum electrodynamics (QED). Below we identify the most significant of them. In general, the following processes are possible for the gamma ray with energy E > 1 MeV:

• γ γ b → e + e -(Breit-Wheeler pair production, PP)

• γ γ b → e + e -e + e -(double pair production, DPP)

• γ γ b → γ γ (photon-photon scattering)

• γ e b → γ e (Compton scattering of γ-rays on free electrons) In processes involving at least one charged particle, it is always possible to emit an additional photon (for instance γ γ b → e -e + γ). However, it is clear that the probability of such processes is suppressed in comparison with the original ones without additional photons in the final state.

Interaction efficiency at redshift z is determined by the interaction rate R(E, z) which shows the number of interactions per unit path length:

R(E, z) = 1 2 ∞ 0 dE n(E , z) 1 -1 dµ (1 -βµ) σ(s) (4.1)
where n(E , z) is the spectral particle number density of target background, s is the energy in the center of mass reference frame and β = (1 -m 2 /E 2 ) is velocity of the considered particle (β = 1 for photons). µ = cos θ where θ is the angle between highenergy particle and background particle. To estimate the interaction rate we can use simplified relation R(E, z) ∼ σ(s) n(z) where n(z) is the background particle number density. For the CMB n(0) ∼ 400 cm -3 , for the infared EBL bump n(0) ∼ 1 cm -3 , for the optical EBL bump n(0) ∼ 0.01 cm -3 and for the free electron density in the IGM n(0) ∼ 10 -7 cm -3 . Taking into account cross-sections of the processes considered and substituting background particle number density we conclude that the most important interaction are Breit-Wheeler pair production and double pair production. Thus highenergy gamma rays interactions results in the production of secondary electrons and positrons. Compton scattering of γ-rays on free electrons is suppressed due to low density of free electrons and photon-photon scattering does not play role because of too small cross-section. The mean free path λ mfp (inverse interaction rate) of the Let us now consider the fate of secondary electron-positron pairs. As well as gamma rays they are subject to QED interactions with background fields. The list of possible interactions includes:

• e γ b → e γ (inverse Compton scattering, ICS)

• e γ b → e e + e -(triplet pair production, TPP)

• e e b → e e (Møller and Bhabha scattering) Repeating the previous reasoning, we conclude that the most important are the inverse Compton scattering and triplet pair production. The mean free path of electrons with respect to these processes is shown on the right panel of Figure 4.2. The inverse Compton scattering leads to the transfer of energy to secondary gamma rays and thus promotes the development of the electromagnetic cascade. Despite the fact that individual interactions with IGM electrons are negligible, the collective effects between electron-positron pairs and a charged medium can lead to the development of plasma instability and pump energy from the cascade. The significance of the plasma instability effect is controversial and has been the subject of extensive discussion in the literature [START_REF] Broderick | The Cosmological Impact of Luminous TeV Blazars. I. Implications of Plasma Instabilities for the Intergalactic Magnetic Field and Extragalactic Gamma-Ray Background[END_REF][START_REF] Schlickeiser | Plasma Effects on Fast Pair Beams in Cosmic Voids[END_REF][START_REF] Batista | The Impact of Plasma Instabilities on the Spectra of TeV Blazars[END_REF][START_REF] Alawashra | Suppression of the TeV Pair-beam-Plasma Instability by a Tangled Weak Intergalactic Magnetic Field[END_REF][START_REF] Miniati | Relaxation of Blazar-induced Pair Beams in Cosmic Voids[END_REF]. At the moment the situation is not settled down, and therefore, in our further analysis, we do not take into account this effect.

From Figure 4.2 one can see that the pair production and inverse Compton scattering govern the development of the cascade almost at all energies excluding the region E > 10 17 eV for the electron interactions. Therefore, we will consider these processes in more detail.

Breit-Wheeler pair production The total cross section for the Breit-Wheeler pair production is given by:

σ PP = σ T 3 16 (1 -β 2 ) (3 -β 4 ) ln 1 + β 1 -β -2β(2 -β 2 ) (4.2)
where σ T is the Thomson cross section, β = 1 -4m 2 e /s and m e is the electron mass. Pair production cross-section σ PP sharply depends on the center of mass energy squared s = 2E γ E b (1-µ) (µ = cos θ). First, the cross section has a threshold s > 4m 2 e requiring that there be enough energy to born a pair. Second the cross section is strongly peaked at s ∼ 8m 2 e . In this way the most intense interaction of a gamma ray with energy E γ will occur with background photons with an energy in the region E b ∼ 4m 2

e /E γ . Substituting numbers it is easy to see that the two peaks in the EBL spectrum will cause an absorption at 1 TeV and 100 TeV and absorption on CMB will happen for γ-rays with energies 1000 TeV. The energy distribution of outgoing electrons is given show the interaction length corresponding to the adiabatic and synchrotron energy losses respectively. Source [START_REF] Settimo | Propagation of extragalactic photons at ultra-high energy with the EleCa code[END_REF] by:

dσ PP dE e = σ T 3 4 m 2 e s 1 E γ E e E γ -E e + E γ -E e E e + E γ (1 -β 2 ) 1 E e + 1 E γ -E e - E 2 γ (1 -β 2 ) 2 4 1 E e + 1 E γ -E e 2 (4.3)
If the interaction occurs near the threshold, then both leptons get approximately the same energy E e ≈ E γ /2. Otherwise, the energy distribution is not symmetric and one of the electrons carries away most of the energy.

Inverse Compton scattering The total cross section for the inverse Compton scattering is given by: e or β 1 (equivalently, this means that the energy of the incident photon in the electron rest frame is much less than the mass of the electron) the cross section can be expanded in a series by a small parameter β, thus recovering the result of classical electrodynamics:

σ ICS = σ T 3 8 m 2 e sβ 2 β(1 + β) (2 + 2β -β 2 -2β 3 ) - 1 β 2 (2 -3β 2 -β 3 ) ln 1 + β 1 -β (4.
σ ICS = σ T (1 -2β) + O(β 2 ) for s -m 2 e m 2 e ( 4.5)
This expansion is especially useful for numerical calculations, since for small β full cross section is numerically unstable because it explicitly contains terms proportional to 1/β and 1/β 2 . Oppositely, in the Klein-Nishina regime when s m 2 e the cross section is getting suppressed as:

σ ICS ∼ σ T m 2 e s ln s m 2 e for s m 2 e (4.6)
The energy distribution of outgoing electrons is given by:

dσ ICS dE γ = σ T 3 8 m 2 e s 1 E e 1 + β β E γ -E e E e + E e E γ -E e + 2(1 -β) β E γ E e -E γ - (1 -β) 2 β 2 E γ E e -E γ 2 (4.7)
It directly follows from this distribution that in the Thomson regime the energy of the secondary photon is a small fraction of the energy of the high-energy electron. Thus, before losing a significant fraction of the energy, the electron will up-scatter a large number of background photons. In the Thomson regime, the average energy of the secondary gamma-ray can be calculated analytically:

< E γ >= E γ dN dE γ dE γ dN dE γ dE γ = 4 3 < E b > E e m e 2 (4.8)
On the other hand, in Klein-Nishina regime the secondary gamma ray carries away almost all the energy of the primary electron.

Numerical modeling of electromagnetic cascades

Electromagnetic cascades in the IGM carry information both about the properties of the spectrum of primary particles and about the parameters of IGM. It turned out that the cascade is sensitive to the intergalactic magnetic fields in the void of the large scale structure. This property of cascades was used, for example, to limit the strength of the IGMF [START_REF] Neronov | Evidence for Strong Extragalactic Magnetic Fields from Fermi Observations of TeV Blazars[END_REF] (see Chapter 5). On the other hand, electromagnetic cascades can be used to study cosmic rays, thus representing the multimessenger connection [START_REF] Berezinsky | Cascade photons as test of protons in UHECR[END_REF].

In the case of a zero magnetic field, the intensity of the cascade signal can be fairly easily estimated from the law of conservation of energy. Moreover, one can analytically predict the shape of the cascade signal, which in many cases does not depend on the initial spectrum ("strong universality") [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF]. A more detailed study of the cascade signal requires the involvement of numerical calculations. The simplest and most computationally fast method is the calculation based on the solution of onedimensional kinetic equations. These equations describe the cascade in terms of the redistribution of the particle number density as a function of energy. In particular, it was this method of calculation that was used for the first derivation of the constraints on the IGMF based on observations in the gamma range [START_REF] Neronov | Evidence for Strong Extragalactic Magnetic Fields from Fermi Observations of TeV Blazars[END_REF]. The disadvantages of this approach are the difficulty of consistently taking into account the deflections in the magnetic fields, as well as the impossibility of tracking the trajectories and interactions of individual particles. The most complex and complete description of the cascade is possible using Monte Carlo methods. At the moment, there are several publicly available Monte Carlo codes, among which EleCa [START_REF] Settimo | Propagation of extragalactic photons at ultra-high energy with the EleCa code[END_REF], CECsi [START_REF] Fitoussi | Physics of cosmological cascades and observable properties[END_REF], ELMAG [START_REF] Blytt | ELMAG 3.01: A three-dimensional Monte Carlo simulation of electromagnetic cascades on the extragalactic background light and in magnetic fields[END_REF], CRPropa [START_REF] Alves | CRPropa 3-a public astrophysical simulation framework for propagating extraterrestrial ultra-high energy particles[END_REF] and CRbeam [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF]. In the current Chapter, we will make a detailed comparison of the three codes: ELMAG 3.01, CRPropa3-3.1.7 and CRbeam 1.0. The modular structure of the CRbeam and CRPropa allows all relevant processes to be tested independently of each other. Namely, we consider the Breit-Wheeler pair production and inverse Compton scattering on the EBL and CMB. For each interaction, we compare interaction rate and energy distribution of secondary particles inferred from simulations with the theoretical predictions. Also, disabling all interactions makes it possible to compare the propagation of electrons in a magnetic field. For ELMAG, on the contrary, such independent testing of interactions is impossible, therefore we use the results of simulations with ELMAG to compare the properties of the cascade signal, when all relevant interactions are turned on.

Absorption of the primary γ-rays

We start our comparison from the modeling of absorption of primary γ-rays on EBL. Fig. 4.3 shows comparison of the calculated optical depths for photons of different energies. All calculations are done for the [START_REF] Franceschini | Extragalactic optical-infrared background radiation, its time evolution and the cosmic photon-photon opacity[END_REF] EBL model. One can see that CRbeam modeling reproduces the assumed analytical model with percent-level precision in the energy range above 200 GeV for a wide range of source redshifts. Differences up to 20% at the energies E < 100 GeV do not have large effects on the model results because the optical depth in this energy range is small even for sources at redshift z ∼ 1.

The results of simulations with ELMAG are also in good agreement with theory in the energy region above 100 GeV, although it predicts ∼ 5% larger optical depth. Strong differences from theory are clearly noticeable in the energy range 10-100 GeV, where ELMAG predicts a significantly higher optical depth. For sources with redshifts z 0.5, this discrepancy will have no effect, since in this case even the wrong optical depth is much less than unity. Nevertheless, differences may arise for sources with z ∼ 1, since in this case in the ELMAG the optical depth τ ∼ 1, and not less than unity as it should be according to the theory.

To the contrary, the calculations of CRPropa differ significantly from the assumed analytical model already for sources at redshifts z 0.3. The difference reaches 30% at 100 GeV for the source at redshift z = 0.3 and reaches up to 100% for further away sources (not shown in the bottom panel of the figure). The discrepancy between the numerical and analytical model in the CRPropa code are due to simplified evolution of the EBL with redshift in which the real evolution of EBL is replaced by the renormalization of EBL spectrum taken at z = 0 with a factor s(z) [START_REF] Alves | CRPropa 3-a public astrophysical simulation framework for propagating extraterrestrial ultra-high energy particles[END_REF]. This CRPropa limitation is known [START_REF] Alves | CRPropa 3-a public astrophysical simulation framework for propagating extraterrestrial ultra-high energy particles[END_REF]. Our result shows explicitly that such a simplification starts to produce incorrect results already starting from moderate redshifts ∼ 0.3.

Emission of secondary γ-rays via inverse Compton scattering

Precision of the modeling of properties of the secondary γ-ray emission depends on the precision of calculation of the optical depth of the primary γ-ray flux (discussed in the previous section), of the energy distribution of the produced pairs and of the differential cross-section of inverse Compton scattering. Top panels shows the primary and secondary γ-ray spectra for sources at different redshifts. Bottom panels show the differences between CRbeam, CRPropa and ELMAG models. An example calculation of secondary γ-ray spectrum by CRbeam is shown in Fig. 4.4 for a monoenergetic primary γ-ray source at different redshifts corresponding to known blazar type active galactic nuclei Mrk 421 (z = 0.03), 1ES 0229+200 (z = 0.14) and PKS 0502+049 (z = 0.954), covering a wide range of redshifts.

We find significant discrepancies that reach 50% for CRbeam-CRPropa comparison and some 30% for CRbeam-ELMAG, in the case of the source at the redshift z ∼ 1. Smaller discrepancies, at the level of 20% over a broad energy range up to 100 GeV are found for smaller redshift sources. Nevertheless, also in this case the discrepancies rise to > 50% at the high energy end of the secondary γ-ray spectrum.

The difference in cascade signal between CRbeam and ELMAG can be fully explained by the difference in the optical depth for the primary γ-rays. Indeed, for a nearby source with a redshift z = 0.03, a 5% greater optical depth for 10 TeV gamma rays results in a 5% amplification of the cascade signal over the entire energy range. For a more distant source with z = 0.14, most of the primary gamma rays are absorbed; therefore, a 5% difference in the optical depth does not appear and the normalization of the cascade signals in the region above 100 MeV coincide. However, increased absorption leads to an increase in the suppression of the high-energy part of the secondary signal. This is manifested in the fact that in the ELMAG the cascade signal is weaker in the region above 100 GeV and stronger in the region below 100 MeV. The same is true for a source at the redshift z = 0.954. In addition, a dip at 10-100 GeV noticeable in the ELMAG spectrum corresponds to the order-of-magnitude excess absorption in ELMAG in this energy range, which leads to an even stronger amplified cascade signal below 100 MeV.

Unlike ELMAG, the difference between CRPropa and CRbeam cannot be explained by this single cause. Part of the discrepancy is certainly due to the the differences in the calculation of the optical depth for the primary γ-rays producing pairs on the EBL. However, further differences in the modeling of the secondary flux are involved. This is indicated by the significant differences of the results of calculation with different codes for low-redshift sources, for which the optical depth calculations do not show strong discrepancies.

We have found that for the low-redshift sources the main difference comes from the calculation of the inverse Compton scattering spectra by electrons and positrons. Figure 4.5 shows a comparison of analytically calculated spectra of the inverse Compton scattering on the CMB with the output of CRbeam and CRPropa for this process (tag "Theory" corresponds to the distribution according to equations 4.7 and 4.3, while "Analytical estimate" is given by equation 4.8). One can see that the CRPropa calculations do not match the theoretical formula. The discrepancies are of the order-of-one all over the energy range of the inverse Compton emission for electrons with energies up to 10 16 eV. This large discrepancy is "smeared" and becomes less noticeable in the calculation of the spectrum of inverse Compton emission from a broad energy distribution of electrons shown in Fig. 4.4.

We have also identified a similar problem of the CRPropa code in calculation of the inverse Compton scattering (ICS) on the EBL (see middle panel of Fig. 4.5) and in the spectra of γγ pair production. Bottom panel of Fig. 4.5 illustrates this problem for the case of pair production on CMB. In CRPropa, the electron-positron spectrum differs from the known analytical formulae [START_REF] Gould | Pair Production in Photon-Photon Collisions[END_REF], shown by the dashed lines.

Having noticed these discrepancies, we analyzed the implementation of the CR-Propa code and found out that the error with the energy distribution of secondary particles has a common origin for all three cases (ICS on CMB, ICS on EBL and pairproduction on CMB). The error occurs from the pre-calculated tables that come with the code. These tables are used to sample the energy at the center of mass as the interaction occurs. After recalculating the tables, the results obtained with the CRPropa with these corrected tables match analytical calculations, as can be seen from Fig. 4.5 (in all plots, we use the 'CRPropa corrected' label for the CRPropa with recalculated tables). After this correction, model calculations of the spectra of secondary emission for our benchmark case of monoenergetic primary γ-ray beam agree between CRbeam and CRPropa for low-redshift sources down to ≤ 5% level, see Fig. 4.4. The discrepancies still grow with increasing source redshift, for the reasons explained in the previous section.

We have also verified that the interaction rates for inverse Compton scattering and pair production on CMB calculated using CRbeam and CRPropa are in agreement with each other and with a relative difference from theory not exceeding one percent,

Propagation of electrons and positrons

Imaging and timing properties of the secondary γ-ray signal are determined by the details of propagation of electrons and positrons under the influence of IGMF. The inverse Compton scattered γ-rays are emitted after electrons and positrons are deflected by an angle δ determined by the energy of electron E e and strength of the magnetic field B. In the simplest case of homogeneous magnetic field the angle is given by the analytical formula

δ = D R L = eB E e D (4.9)
where e is the electron charge, D is the propagation distance and R L = E e /eB is the gyroradius.

Both CRPropa and CRbeam correctly solve equations of motions of the electrons in magnetic fields. However, the approaches of the two codes differ when accounting for cosmological effects. In CRbeam, the particle coordinates are expressed in physical distances, while in CRPropa, the comoving coordinates are used. This difference can be illustrated by a simple test where the electron is deflected by a homogeneous magnetic field that is perpendicular to the electron's velocity. Figure 4.7 shows the calculation of deflection angle of electrons for this simplest case, found with CRbeam and CRPropa codes. One can see that both calculations give the same result at the redshift z = 0 where physical and comoving distance elements are equal, left panel of Fig. 4.7. However, the CRPropa calculation predicts (1 + z) stronger deflection at the redshift z = 1, which is explained by the fact that dx comoving = (1 + z)dx physical , so that electron travels (1 + z) longer distance in the same magnetic field, see middle panel of Fig. 4.7.

This difference becomes most pronounced when the electron travels over cosmological distances. Right panel of Fig. 4.7 shows the deflection of electron starting at redshift z (x-axis of the panel) and caught at z = 0. One can see that the deflection angles in CRPropa are proportional to the comoving distances and in CRbeam to the light travel distances. The simplest way to get the correct deflection angles in CRPropa is to set the evolution of the magnetic field with redshift, namely

B(z) = B 0 /(1 + z),
what can be done using the built-in methods of CRPropa.

The field re-scaling trick can also be used when an electron propagates in a realistic turbulent magnetic field if the correlation length of this field is greater than the electron propagation length. Nevertheless, it is not clear how to construct the correct trajectories when the correlation length of the magnetic field is comparable to or less than the propagation length of an electron. We have compared the results on electron deflection also for this case at z = 0. For the magnetic field we used Kolmogorov spectrum with the maximum scale L B and minimum scale equal to L B /100. Fig. 4.8 shows the dependence of deflection angle on propagation distance for electrons moving in turbulent magnetic fields. The codes are in good agreement with each other and reproduce analytical formula of [START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF] in case of small correlation length

δ 2 = √ Dλ B R L = √ Dλ B eB E e (4.10)
where λ B is the magnetic field correlation length. For the Kolmogorov spectrum λ B ≈ L B /5. We also tested the turbulent magnetic field generator in ELMAG. Using the built-in test_turbB function, we separately calculated the average value of each component of the magnetic field. Assuming that the field of the strength B is uniform and isotropic, it was expected that

< B 2 x > = < B 2 y > = < B 2 z > = B/ √ 3 (4.11)
As a result of tests, we found that the values of the x and y components of the magnetic field are about 20 percent less than predicted theoretically, while there is no such difference for the z component. The reason for this was not found, so we simply rescaled the x and y components of the field to the correct value. The results of ELMAG with such a re-scaled field are indicated on the plots as 'ELMAG corrected'.

Properties of primary and secondary emission from sources at different distances

Realistic modeling of spectral, imaging and timing properties of the primary and secondary γ-ray signals from extragalactic sources combines the modeling of the pair production, inverse Compton scattering and electron propagation discussed above. Discrepancies in the modeling of these ingredients propagate to the discrepancies of the overall signal models. We discuss these discrepancies in this section.

Left side of Fig. 4.9 shows the result of modeling of the primary and secondary γ-ray signals for the nearby blazars like Mrk 421 and Mrk 501. As discussed in Section 4.3, all the codes provide identical predictions for the optical depth for the primary γ-rays in this case. At noticed in Section 4.4, predictions for the spectra of pair production and inverse Compton scattering in the CRPropa code deviate from known analytical formulae. This leads to discrepancies between the models of secondary γ-ray fluxes in CRbeam and CRPropa codes, noticeable in the TeV energy range, which is most interesting for CTA. Top panel of the figure shows that the differences reach a factor-of-two above 3 TeV. Middle and bottom panels for the figure show the temporal 10 10 10 11 10 12 10 13 E, eV

t dN dt E 1 -3 G e V E 1 0 - 3 0 G e V E 1 0 0 -3 0 0 G e V E 1 -3 T e V
CRbeam ±10% CRPropa CRPropa corrected 10 1 10 1 10 3 10 5 10 7 10 9 t, year dN d

E 1 -3 G e V E 1 0 -3 0 G e V E 1 0 0 - 3 0 0 G e V CRbeam ±10% CRPropa CRPropa corrected
E 1-3 GeV E 1 0 - 3 0 G e V E 1 0 0 - 3 0 0 G e V E 1 - 3 T e V
CRbeam ±10% CRPropa CRPropa corrected 10 4 10 3 10 2 10 1 10 0 10 1 , degree ). The jet opening angle was assumed to be θ jet = 5 • , and the observer was at an angle θ obs = 3 • with respect to the jet axis. Red solid lines is for CRbeam code, blue dashed lines is for CRPropa code and green dash-dotted lines -ELMAG code. Line labels indicate the constraints on the source activity time dt and angular distance from the source θ that were used when constructing the cascade signal. Label 'jet' denotes the cascade signal without constraints on dt (dt = ∞) and θ (θ = 2π). For comparison, we also added a cascade spectrum in the absence of a magnetic field (label 'full').

E 1 -3 G e V E 1 0 -3 0 G e V E 1 0 0 -3 0 0 G e V CRbeam ±10% CRPropa CRPropa corrected
and angular profiles of extended emission. The differences in the spectral properties of the secondary signal between CRbeam and CRPropa are not so pronounced in the low-energy temporal and angular profiles in the lower energy bands, accessible to Fermi/LAT telescope. However, these differences grow in the CTA energy band above 100 GeV. As discussed in Section 4.4, the discrepancy between the CRbeam and CRPropa predictions can be removed after a correction to the CRPropa code. In this case, all the codes produce results that agree at ∼ 10% level, as can be seen from the left side of Fig. 4.9.

To the contrary, corrections to the CRPropa code do not remove the discrepancies between different model predictions for a source at large redshift, see right side of Fig. 4.9. Very significant difference in model predictions is noticeable at the highest energy ends of the spectra, visible in the top right panel of Fig. 4.9. CRPropa predicts strong suppression of the flux above 200 GeV, while CRbeam predicts a gradual suppression in the energy range between 200 and 300 GeV. This discrepancy is crucially important for the study of high-redshift sources with CTA. Such modeling errors reaching orderof-magnitude have to be removed before any sensible results on the EBL evolution can be extracted from the CTA data. The discrepancy is noticeable also in the attenuated primary source spectra shown by the grey curves in Fig. 4.9. We attribute this discrepancy to the incorrect modeling of absorption on EBL for high redshift sources in CRPropa, as described in Section 4.3.

Similarly large are discrepancies in the timing and imaging properties of the secondary γ-ray signal, shown in the middle and bottom panels on the right side of Fig. 4.9. Contrary to the low-redshift sources, differences in predictions for the peaking time of the delayed signal are large already for lower energy γ-rays detectable with Fermi/LAT.

The results for the source 1ES 0229+200 at redshift z = 0.14 are presented in Fig. 4.10. Calculations are done for IGMF with Kolmogorov spectrum with strength B = 10 -15 G and maximum scale L B = 10 Mpc. In this figure the full cascade contribution for zero magnetic field is denoted as 'full'. The curves marked 'jet' show the extended and delayed fluxes for non-zero IGMF for different angular and time intervals. The results for the CRbeam code are compared with CRPropa and ELMAG before (upper panel) and after (lower panel) corrections. All codes give similar results for most of the parameters, after corrections.

Conclusions of the Chapter

Our comparative study of model predictions for the primary and secondary veryhigh-energy γ-ray signals from astronomical sources at different redshifts has revealed significant differences in model predictions obtained with different Monte-Carlo codes CRbeam, CRPropa and ELMAG. Most remarkably, we have identified large differences in model calculations of observables crucial for cosmological probes that will be enabled by the new CTA γ-ray observatory.

It has been found that the predictions of the spectral, angular, and timing properties of a cascaded signal differ by about 50% over the energy range accessible to CTA. We have identified the origin of most of these discrepancies, notably in the modeling of attenuation of γ-rays and deflections of electrons in the CRPropa code. Once corrections are introduced to the CRPropa code, for the nearby sources with z ∼ 0.1 the discrepancies are reduced to ∼ 10%. This number can be considered as an estimate of the systematic modeling error. However, the differences in the implementation of cosmological evolution in CRPropa and CRbeam and the simplified calculation of the optical depth of the EBL in CRPropa lead to significant differences between model calculations for sources at high redshifts z > 0.3.

Chapter 5

Intergalactic magnetic fields

Magnetic fields are present in almost all observable astronomical objects, but their existence and role in the early Universe is still not clear. In the modern Universe, primordial magnetic fields could remain in the intergalactic medium, in the voids of a large-scale structure. Existing measurements limit the intensity of extragalactic magnetic fields in the range 10 -16 G< B < 10 -9 G, but the history of the origin of these fields is not known.

The lower bounds on extragalactic magnetic fields were established by gamma-ray astronomy methods and are based on the non-observation of a cascade signal from distant gamma-ray sources. This Chapter is devoted to an analysis of the prospects for studying the IGMF using gamma radiation measurements. The sensitivity of the new generation telescope CTA to strong IGMFs with strengths in the region of 10 -12 G< B < 10 -11 G will be evaluated. Also, the sensitivity of CTA to IGMF with a large correlation length will be calculated. On the other hand, suppression of the cascade signal by magnetic field bubbles around galaxies and clusters of galaxies is discussed.

The results of this Chapter are published in [START_REF] Korochkin | Sensitivity reach of gamma-ray measurements for strong cosmological magnetic fields[END_REF][START_REF] Korochkin | Detectability of large correlation length inflationary magnetic field with Cherenkov telescopes[END_REF][START_REF] Korochkin | Sensitivity of the gamma-ray method for strong primordial magnetic fields[END_REF][START_REF] Bondarenko | Account of baryonic feedback effect in the γ-ray measurements of intergalactic magnetic fields[END_REF].

Observations of intergalactic magnetic fields

Magnetic field in galaxies and galaxy clusters

The existence of magnetic fields (MF) in our Galaxy Milky Way is firmly established. Several methods have been developed for estimating the MF in the Galaxy, such as, measurements of the Zeeman splitting [START_REF] Crutcher | Magnetic Fields in Molecular Clouds: Observations Confront Theory[END_REF], studies of the polarized starlight [START_REF] Heiles | The Local Direction and Curvature of the Galactic Magnetic Field Derived from Starlight Polarization[END_REF], Faraday rotation measurements of point sources [START_REF] Simard-Normandin | The Faraday rotation measures of extragalactic radio sources[END_REF][START_REF] Kronberg | Remarkable Symmetries in the Milky Way Disc's Magnetic Field[END_REF], diffuse polarized synchrotron and thermal dust emission surveys [191]. These methods independently indicate that the field strength in the Galaxy is in the range 1-10 µ G [START_REF] Pshirkov | Deriving the Global Structure of the Galactic Magnetic Field from Faraday Rotation Measures of Extragalactic Sources[END_REF][START_REF] Jansson | A New Model of the Galactic Magnetic Field[END_REF]. Moreover, observations show that magnetic fields of the same magnitude are present both in other galaxies and in galaxy clusters.

The problem of the origin of such strong magnetic fields in galaxies and clusters of galaxies is one of the long-standing problems of astrophysics. It is usually assumed that these fields result from dynamo amplification of a much weaker seed field presented in the medium at the beginning of the structure formation process [START_REF] Beresnyak | Turbulent Amplification and Structure of the Intracluster Magnetic Field[END_REF]. However, the nature of the seed field remains unclear.

Two classes of scenarios of the magnetogenesis are presented in the literature: cosmological and astrophysical. In cosmological scenarios, seed fields originate from the processes in the Early Universe and afterwards decay to the current state. Several physical phenomena that took place a fraction of a second after the Big Bang can be responsible for the generation of a relic magnetic field. First order phase transitions that might have happened at quark confinement or during the Electroweak epoch can produce short correlation length magnetic field that evolves through turbulent decay toward a magnetic field configuration with correlation length and strength satisfying a relation λ B ∼ 0.1[B/10 -12 G] kpc today [START_REF] Banerjee | The Evolution of cosmic magnetic fields: From the very early universe, to recombination, to the present[END_REF][START_REF] Kahniashvili | Evolution of Primordial Magnetic Fields from Phase Transitions[END_REF]. Alternatively, a field generated at the epoch of inflation can have a very large correlation length, up to the present day Hubble scale [START_REF] Ratra | Cosmological "Seed" Magnetic Field from Inflation[END_REF][START_REF] Daniel Garretson | Primordial magnetic fields from pseudoGoldstone bosons[END_REF][START_REF] Gasperini | Primordial magnetic fields from string cosmology[END_REF][START_REF] Giovannini | Primordial magnetic fields from inflation?[END_REF][START_REF] Demozzi | Magnetic fields from inflation?[END_REF][START_REF] Ricardo | Inflationary magnetogenesis without the strong coupling problem[END_REF][START_REF] Durrer | Cosmological Magnetic Fields: Their Generation, Evolution and Observation[END_REF]. On the other hand, in astrophysical scenarios, weak seed magnetic fields arise due to local effects after recombination. An example is mechanisms based on Biermann Battery effect [START_REF] Russell | The Protogalactic origin for cosmic magnetic fields[END_REF].

If the seed field is of cosmological origin measurement of parameters of this seed can provide a new cosmological observable, a new "window" on the Early Universe. However, during the evolution of galaxies, their magnetic fields were significantly reprocessed and therefor the information about the seed field was lost. In particular, the dynamo mechanism leads to observed galactic magnetic field values largely independent of the characteristics of the seed field. Oppositely, large-scale magnetic fields that fill the volume of the voids of the Large Scale Structure (LSS) likely did not experience adiabatic contraction or dynamo amplification that transformed the seed fields in galaxies and galaxy clusters. The properties of these volume-filling intergalactic magnetic fields (IGMF) might therefore be close to the properties of the cosmological seed fields.

The problem is that the voids can be "polluted" by the magnetic field spread by the baryonic feedback process that returns matter from galaxies into the intergalactic medium. The galactic winds driven by the supernova and AGN activity are most probably ionized and carry magnetic fields with them [START_REF] Bertone | Magnetic Field Seeding by Galactic Winds[END_REF][START_REF] Pinsonneault | Anisotropic Galactic Outflows and Enrichment of the Intergalactic Medium. II. Numerical Simulations[END_REF]. Uncertainties of the details of the baryonic feedback on the LSS do not allow to answer the question if the magnetic field spread by the feedback can fill the voids and thus hide the cosmological field. If this is the case, measurement of IGMF will fail to provide a new cosmological observable. Instead, it would rather provide constraints on the baryonic feedback process.

-13

10 -9 10 -5 10 -1 10 3

L B [Mpc]

Constraints on the intergalactic magnetic fields

Upper limits

From the observational point of view, there is a number of techniques to constrain the IGMF. First, we describe the bounds from above summarized on Figure 5.1.

A robust upper bound on the IGMF strength is provided by the measurement of the Zeeman splitting in the 21 cm absorption line in the spectra of distant quasars [START_REF] Arthur | An 84 microGauss Magnetic Field in a Galaxy at Redshift z=0.692[END_REF]. These observations are used to infer the magnetic field in the Milky Way and other galaxies and give the results in the range 1-100 µG. Stronger IGMFs along the line of sight from the source to the observer are ruled out because they would make a significant impact on the observations. It is clear that the constraints from Zeeman splitting are independent of the magnetic field correlation length λ B .

Next group of constraints on IGMF are based on Faraday rotation measurements of polarised radio emission from quasars [START_REF] Kronberg | A Global Probe of Cosmic Magnetic Fields to High Redshifts[END_REF][START_REF] Neronov | Magnetic Fields in the Large Scale Structure from Faraday Rotation measurements[END_REF][START_REF] Akahori | Exploring the Intergalactic Magnetic Field by Means of Faraday Tomography[END_REF][START_REF] Pshirkov | New limits on extragalactic magnetic fields from rotation measures[END_REF][START_REF] O'sullivan | New constraints on the magnetization of the cosmic web using LOFAR Faraday rotation observations[END_REF]. This method is based on the effect of rotation of the polarization plane of radiation when passing through a magnetised medium. The rotation measure is proportional to the projection of the magnetic field onto the line of sight B || and thus depends on the correlation length. If λ B D (D is the distance to the source) the rotation measure accumulation can be described as a random walk process and therefore it is proportional to the ∼ √ Dλ B . Thus the constraints from Faraday rotation measurements are relaxed as B ∼ λ -1/2 B as the correlation length decreases (see Figure 5.1).

If the IGMF is of the cosmological origin then it is possible to limit its parameters using CMB observations. In this case magnetic fields interact with the plasma in the Early Universe. If the magnetic field is sufficiently strong it will affect the evolution of primordial plasma and can potentially manifest itself in the properties of the CMB. In particular magnetic field might have left imprints in the temperature anisotropies, spectral distortion, the Faraday rotation of the polarized CMB emission, or change the Recombination process itself via density fluctuations [START_REF] Barrow | Constraints on a primordial magnetic field[END_REF][START_REF] Jedamzik | A Limit on primordial small scale magnetic fields from CMB distortions[END_REF][START_REF] Jedamzik | Stringent Limit on Primordial Magnetic Fields from the Cosmic Microwave Background Radiation[END_REF]. CMB observations set the strongest upper limits in the region of large correlation lengths.

Lower limits

The technique of measuring the IGMF that is best suited for probing weak magnetic fields in the voids is based on γ-ray observations of distant active galactic nuclei (AGN). The technique exploits the effect of attenuation of very high-energy γ-ray (E > 1 TeV) flux from distant AGN by the effect of pair production on photons of the EBL. As we discussed in details in the previous Chapter this leads to the development of an electromagnetic cascade along the line of sight. A schematic representation of the cascade in magnetic field is shown on Figure 5.2. The magnetic field in the intergalactic medium can deflect electron-positron pairs produced by the primary γ-rays before they emit secondary γ-rays by inverse Compton scattering. Thus secondary cascade photons could be emitted toward the observer from whose γ-rays which were initially emitted slightly away. Thereby the cascade signal is observable in the form of magnetic fielddependent delayed [START_REF] Plaga | Detecting intergalactic magnetic fields using time delays in pulses of γ-rays[END_REF] or extended [START_REF] Neronov | A method of measurement of extragalactic magnetic fields by TeV gamma ray telescopes[END_REF][START_REF] Neronov | Sensitivity of gamma-ray telescopes for detection of magnetic fields in intergalactic medium[END_REF][START_REF] Alves | Probing Intergalactic Magnetic Fields with Simulations of Electromagnetic Cascades[END_REF] emission around the primary point source.

Searches for extended and delayed γ-ray emission in the energy band of the Fermi Large Area Telescope (LAT) were unsuccessful. However, the non-detection of the cascade emission means that the IGMF is strong enough to deflect e+e-pairs and thus suppress the secondary signal. In this way a lower bound on the IGMF strength was imposed [START_REF] Neronov | Evidence for strong extragalactic magnetic fields from Fermi observations of TeV blazars[END_REF][START_REF] Taylor | Extragalactic magnetic fields constraints from simultaneous GeV-TeV observations of blazars[END_REF][START_REF] Dermer | Time Delay of Cascade Radiation for TeV Blazars and the Measurement of the Intergalactic Magnetic Field[END_REF][START_REF] Abramowski | Search for Extended \gamma-ray Emission around AGN with H.E.S.S. and Fermi-LAT[END_REF][START_REF] Finke | Constraints on the Intergalactic Magnetic Field with Gamma-Ray Observations of Blazars[END_REF]. The most recent analysis of Fermi/LAT data reported by [START_REF] Ackermann | The Search for Spatial Extension in High-latitude Sources Detected by the F ermi Large Area Telescope[END_REF] constrained the IGMF in the voids to be stronger than ∼ 10 -16 G for large correlation-length magnetic fields and ∼ 10 -14 G for short correlation-length fields originating from the early Universe. The next-generation γ-ray observatory, the Cherenkov Telescope Array (CTA), will be able to probe stronger IGMFs with strengths up to ∼ 10 -12 G [START_REF] Korochkin | Sensitivity reach of gamma-ray measurements for strong cosmological magnetic fields[END_REF][START_REF] Abdalla | Sensitivity of the Cherenkov Telescope Array for probing cosmology and fundamental physics with gamma-ray propagation[END_REF].

Magnetic fields in the early Universe: hints of the existence

Recently, several interesting indications have been received of the existence of magnetic fields in the Early Universe. First, [START_REF] Jedamzik | Relieving the Hubble tension with primordial magnetic fields[END_REF] have shown that accounting for the magnetic field driven turbulence on plasma at the epoch of recombination modifies the estimate of the Hubble parameter from the CMB data and relaxes the 4.4σ tension between Planck measurements H 0 = 67.36 ± 0.54 km/(s Mpc) based on z 10 3 data [START_REF] Aghanim | Planck 2018 results. VI. Cosmological parameters[END_REF] and present day universe measurements H 0 = 74.03 ± 1.42 km/(s Mpc) using supernovae Type Ia [START_REF] Riess | Large Magellanic Cloud Cepheid Standards Provide a 1% Foundation for the Determination of the Hubble Constant and Stronger Evidence for Physics beyond ΛCDM[END_REF] and H 0 = 73.3 +1. 7 -1.8 'km/(s Mpc) based on gravitationally lensed systems [START_REF] Kenneth | H0LiCOW XIII. A 2.4% measurement of H 0 from lensed quasars: 5.3σ tension between early and late-Universe probes[END_REF]. The reasoning of [START_REF] Jedamzik | Relieving the Hubble tension with primordial magnetic fields[END_REF] further develops the argument of [START_REF] Jedamzik | Stringent Limit on Primordial Magnetic Fields from the Cosmic Microwave Background Radiation[END_REF] that a magnetic field present at the epoch of recombination induces clumping of baryonic matter and in this way modifies the recombination process. [START_REF] Jedamzik | Stringent Limit on Primordial Magnetic Fields from the Cosmic Microwave Background Radiation[END_REF] have used this argument to derive a strong upper bound B 10 -11 G at the epoch of recombination. This bound is shown by the black arrow in Fig. 5.3. The cosmological magnetic field with the strength ∼ 10 -11 G is expected to have a correlation length of about λ B 1[B/10 -11 G] kpc, which is the largest processed eddy size at the epoch of recombination [START_REF] Banerjee | Evolution of cosmic magnetic fields: From the very early Universe, to recombination, to the present[END_REF].

Another indication of the existence of cosmological magnetic fields can be derived from measurements of redshift dependent absorption by the 21 cm line of neutral hydrogen at redshifts of about z ∼ 10. The EDGES experiment has recently reported an excess opacity of the universe in the redshift range 15 < z < 20 [START_REF] Bowman | An absorption profile centred at 78 megahertz in the sky-averaged spectrum[END_REF]. This indicates that an efficient cooling mechanism of baryonic matter is operating during this epoch. Interaction of matter with magnetic field can provide both cooling (via pumping of thermal energy into magnetic field) and heating (via decay of MHD turbulence) of Figure 5.3 -Known constraints on the strength and correlation length of the IGMF [START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF][START_REF] Durrer | Cosmological magnetic fields: their generation, evolution and observation[END_REF]. Red, blue, and green lines show the hints of the existence of a cosmological magnetic field from the CMB [START_REF] Jedamzik | Relieving the Hubble tension with primordial magnetic fields[END_REF], 21 cm line [START_REF] Kumar Natwariya | EDGES signal in the presence of magnetic fields[END_REF] and baryogenesis [START_REF] Fujita | Large-scale magnetic fields can explain the baryon asymmetry of the Universe[END_REF] correspondingly. Blue dashed regions show the sensitivity of different detection techniques [START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF][START_REF] Durrer | Cosmological magnetic fields: their generation, evolution and observation[END_REF]. Black upper bound is from the analysis of the CMB signal by [START_REF] Jedamzik | Stringent Limit on Primordial Magnetic Fields from the Cosmic Microwave Background Radiation[END_REF]. baryonic matter. It is possible to interpret the EDGES observation in terms of cooling through the interaction with magnetic fields [START_REF] Kumar Natwariya | EDGES signal in the presence of magnetic fields[END_REF]. The magnetic field strengths, which is required to explain EDGES data, is comparable to that needed to resolve the Hubble parameter measurement tension, 5 × 10 -13 G< B < 6 × 10 -12 G [START_REF] Kumar Natwariya | EDGES signal in the presence of magnetic fields[END_REF].

Finally, presence of a helical magnetic field at the epoch of electroweak phase transition can enable an explanation of the baryon asymmetry of the universe within the standard model of particle physics [START_REF] Giovannini | Primordial hypermagnetic fields and the triangle anomaly[END_REF][START_REF] Fujita | Large-scale magnetic fields can explain the baryon asymmetry of the Universe[END_REF]. The range of magnetic field strength 10 -14 G< B < 10 -12 G, which is compatible with this baryogenesis scenario is shown by the green shading in Fig. 5.3. Remarkably, the field strength required for a successful explanation of the baryon asymmetry is consistent with that needed for an explanation of the EDGES signal and of the Hubble parameter measurement tensions.

The combination of these observational hints for the existence of a cosmological magnetic field defines an order-of-magnitude wide "sweet spot" around B ∼ 10 -12 G in which the field estimates from multiple effects intersect. The most convincing evidence for the existence of a field with such strength would be its direct detection in the intergalactic medium. In the next section we explore the possibility of the measurement of such a field with γ-ray telescopes. We demonstrate that even though the field is at the upper sensitivity end of the γ-ray technique, its detection should still be possible with a deep exposure of the nearest blazars with CTA.

Sensitivity of CTA to strong IGMF

In current section we estimate the sensitivity of the CTA telescope to strong cosmological magnetic fields with strength in the range 10 -12 G< B < 10 -11 G and correlation length in the range 0.1 -1 kpc correspondingly. Such a magnetic field may be generated in the Early Universe during the phase transitions.

Analytical estimates

Fields with strengths in the range of B ∼ 10 -12 G are at the upper end of the sensitivity reach of the γ-ray measurement method [START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF]. They are strong enough to deflect trajectories of electrons with energies in the 10-100 TeV range. This implies that the highest energy γ-ray signal accessible to telescopes should be used for the signal measurements. In this situation it is not clear if the small angle deflection approximation previously used for the sensitivity estimates used by [START_REF] Neronov | Sensitivity of γ-ray telescopes for detection of magnetic fields in the intergalactic medium[END_REF] is valid. We reassess the analytical estimates in this high-energy / strong field regime below.

The correlation length λ B of cosmological magnetic fields scales with the strength as [START_REF] Banerjee | Evolution of cosmic magnetic fields: From the very early Universe, to recombination, to the present[END_REF] B ∼ 10

-11 λ B 1 kpc G (5.1)
The field power spectrum is shaped by the turbulence on the scales shorter than λ B . The integral length scale of the field with power spectrum

P B (k) ∝ k -n as λ B ≈ L B n -1 2n = L B 5 (5.2)
for n = 5/3 (assuming the Kolmogorov turbulence spectrum), where L B is the maximum scale of the Kolmogorov spectrum.

We consider secondary emission induced by interactions of primary γ-rays with energies E γ0 . The mean free path of these γ-rays through the EBL is

λ γ0 2.5 E γ0 100 TeV -1.6

Mpc

(5.3)

For the analytical estimates we assume that each primary γ-ray produces an electron and a positron with energies of E e = E γ0 /2. The electrons and positrons cool due to the inverse Compton scattering of CMB photons on the distance scale

D e 7 E e 50 TeV
The gyroradius of electrons is

R L 5 E e 50 TeV B 10 -11 G -1 kpc (5.6)
Electrons are typically deflected by an angle of

δ = λ B R L 0.06 E e 50 TeV -1 B 10 -11 G 2 (5.7)
on the distance scale about the integral length of the magnetic field (we have used the scaling B ∝ λ B suggested by Eqs. (5.1) and (5.2). Accumulation of such small deflections on the electron cooling distance scale D e results in the overall deflection

∆ = D e λ B δ 0.2 E e 50 TeV -3/2 B 10 -11 G 3/2 0.2 E γ 8 TeV -3/4 B 10 -11 G 3/2
(5.8) 1. Our numerical modeling takes into account the full Klein-Nishina cross-section of inverse Compton scattering. The Klein-Nishna effect corrects analytical estimates, but does not change the qualitative picture presented in this section.

If the field strength is B ∼ 10 -11 G, the opening angle of the secondary emission cone at 2 TeV can be as large as the opening angle of the AGN jets. The secondary emission flux within the cone gets suppressed as

F ext F γ0 = Θ 2 jet ∆ 2 1 E γ 8 TeV 3/2
(5.9)

B 10 -11 G -3 Θ jet 10 • 2 , B 10 -11 G
where we have assumed Θ jet ∼ 10 • 0.2. Such flux suppression occurs below the energy at which ∆ = Θ jet ,

E crit = 8 B 10 -11 G 2 Θ jet 10 • -4/3

TeV

(5.10)

If B < 10 -11 G, the deflection angle ∆ is smaller than the opening angle of the jet and extended emission is still observable toward 10 TeV energy.

The maximal possible angular size of the extended emission is determined by the transverse size of the jet at the distance λ γ0

R ext,max = Θ jet λ γ0 0.5 Θ jet 10 • E γ0 100 TeV -1.6

Mpc

(

This corresponds to the angular size

Θ ext,max = R ext D 0.24 • Θ jet 10 • E γ0 100 TeV -1.6
(5.12)

D 120 Mpc -1 0.24 • Θ jet 10 • E γ 8 TeV -0.8 D 120 Mpc -1
where we have used the distances to Mrk 421 and Mrk 501 for the numerical estimate. The time delay of the extended signal can be estimated as

T ext,max = DΘ 2 ext,max /c (5.13) 8 Θ jet 10 • 2 E γ 8 TeV -1.6
kyr This time scale imposes a requirement on the duty cycle of the source for which the extended emission is detectable: the source should have been active over the last 10 kyr, which is plausible given the fact that the source posseses a kiloparsec scale jet [START_REF] Giroletti | The jet of Markarian 501 from millions of Schwarzschild radii down to a few hundreds[END_REF].

Note that the maximal extended source size does not depend on the magnetic field strength. Instead, it is the extended source flux that scales with magnetic field. However, this dependence could hardly be used for the estimation of the field strength because of uncertainty of the opening angle of the jet and because of the very rapid decline of the extended source flux with the field strength. Because of this limitation, we choose to show in Fig. 5.3 the lower bound for the Fermi/LAT telescope derived in the regime of small deflection angles of electrons (see [START_REF] Ackermann | The Search for Spatial Extension in High-latitude Sources Detected by the Fermi Large Area Telescope[END_REF] for details). Below we present arguments that at least one source, Mrk 501, fulfills these three conditions and can be used for the probe of strong IGMF of cosmological origin.

Selection of best target for the search of strong IGMF

There are currently no measurements or estimates of the blazar luminosities at 100 TeV. The highest energy photons detected from blazars are those from the two nearest BL Lac objects, Mrk 421 and Mrk 501. Fig. 5.4 shows the spectral energy distribution of these sources measured by HAWC [START_REF] Coutiño De Leon | Spectral analysis of the blazars Markarian 421 and Markarian 501 with the HAWC Gamma-Ray Observatory[END_REF]. One can see that even though Mrk 421 is the brightest source at TeV, it has a softer spectrum and its intrinsic luminosity is most probably strongly suppressed at 100 TeV. To the contrary, Mrk 501 has a harder spectrum that does not show any signature of high-energy cutoff. The spectrum is measured up to 20 TeV. In view of this fact, we consider Mrk 501 as a more promising candidate for the search of the strongest IGMF and the following calculations are limited to this source. We assume that its intrinsic spectrum extends up to 100 TeV, as shown in Fig. 5.4.

Extrapolating the power law measured by HAWC up to 100 TeV, we find that the intrinsic flux of the source at this energy should be at the level of

F γ0 (100 TeV) (2.5 ± 1) × 10 -12 TeV/(cm 2 s) (5.14)
Most of the intrinsic flux in the 100 TeV energy range is absorbed in the intergalactic medium and is converted to the extended emission. In the regime B < 10 -11 G, the extended emission flux is not suppressed by the widening of the opening angle of the emission cone and one can expect that the flux of the extended source is

F ext (8 TeV) ∼ F γ0 (100 TeV)/2 ∼ 10 -12 TeV/(cm 2 s) (5.15)
This flux level is well above the CTA point-source sensitivity, shown by the lowermost green solid curve in Fig. 5.4. The sensitivity is limit F CT A,ps 5 × 10 -14 TeV/(cm2 s) in the 1-10 TeV range is determined by the statistical fluctuations of the background within the point spread function Θ psf 0.04 • 2 . The sensitivity for extended source flux worsens with the increase of the source size, roughly as

F CT A,ext = F CT A,ps Θ ext Θ P SF 3 × 10 -13 Θ ext 0.25 • TeV/(cm 2 s) (5.16)
This sensitivity is compared with the expected extended source flux in Fig. 5.4. As one can see, the extended flux is above the sensitivity limit and might be detectable. However, the extended emission appears on top of the point-source flux and more detailed analysis is required for verification of its detectability. The γ-ray measurement method is sensitive to the weakest volume filling IGMF along the line of sight toward the source. This field is typically found in the voids of the LSS. TeV γ-rays can cross many voids before being absorbed in pair production. This is not the case for the 100 TeV γ-rays, which can potentially be used to probe the strongest IGMF B ∼ 10 -11 G. These γ-rays can only probe the field within several megaparsecs around the source. Indeed, in this case the effect of extended emission is most significant at energies of secondary photons in the range of ∼ 2 -10 TeV. The corresponding energy of primary photons is about 50 -100 TeV. Considering the standard EBL models the typical mean free path of such photons is 3 -10 Mpc. If the source is occasionally found in or near a large galaxy cluster, or the line of sight points toward a filament of the LSS, the field strength can be much higher than the typical void field making the line of reasoning described above inapplicable.

To explore environmental effects around Mrk 501, we rely on constrained cosmological simulations derived using the Bayesian Origin Reconstruction from Galaxies (BORG) inference method [START_REF] Jasche | Bayesian Physical Reconstruction of Initial Conditions from Large-Scale Structure Surveys[END_REF]. BORG provides a statistical ensemble of initial conditions that all match detailed observations of galaxy count in the universe. To achieve this, it relies on modeling the evolution of large-scale structures starting from high redshift using a dynamical model. Then the probability of such a sample is compared volume element by volume element to data from a galaxy survey assuming some bias function and Poisson statistics for the number count [START_REF] Jasche | Bayesian Physical Reconstruction of Initial Conditions from Large-Scale Structure Surveys[END_REF]. The initial conditions can then be resimulated with a software of the user's choice. For the data considered here, the dynamical model is a particle mesh N -body solver with 20 timesteps; the data are provided by the 2M++ galaxy compilation [START_REF] Lavaux | The 2M++ Galaxy Redshift Catalogue[END_REF]. The resolution of the density contrast of the initial conditions and at z = 0 are set to 2.64h -1 Mpc. The details of the algorithm and of the inference procedure are provided in [START_REF] Jasche | Physical Bayesian Modelling of the Non-Linear Matter Distribution: New Insights into the Nearby Universe[END_REF].

Using the dark matter density as an estimator of the strength of magnetic field we found (see Fig. 5.5) that the IGMF around the source does not appear in the reach galaxy cluster. Instead there is a large void just next to the Mrk 501 with the average density ρ below the critical density of the universe and there the magnetic field is not amplified in the course of structure formation. Numerical modeling of the IllustrisTNG simulation [START_REF] Marinacci | First results from the IllustrisTNG simulations: radio haloes and magnetic fields[END_REF] shows that such amplification only occurs in the regions with overdensity ρ/ρ cr 10.

Numerical Modeling

Qualitative arguments presented in the previous section show that even if the IGMF of cosmological origin is as strong as B ∼ 10 -11 G, it should still be possible to detect it via the search of extended emission around Mrk 501, which is the brightest blazar in the 10 TeV sky.

In this section we support this qualitative argument with numerical modeling of the extended source signal. For this purpose we use the Monte Carlo simulation code CRbeam discussed in Chapter 3 [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF].

We consider a primary γ-ray source with the power-law spectrum with the slope dN/dE ∝ E -2.4 extending up to 100 TeV energy, situated at the distance of D = 150 Mpc. The γ-rays are emitted into a jet with opening angle Θ jet = 5 • aligned along the line of sight.

We propagate the primary γ-ray beam toward the observer, taking into account pair production on extragalactic background light (EBL). The EBL spectrum is that of [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. Electrons and positions produced by absorbed γ-rays loose energy due to inverse Compton scattering on EBL and CMB. The IGMF present that deflects electrons and positrons is assumed to have strength and correlation length satisfying the relation (5.1).

The field is generated with the Kolmogorov power spectrum, following the method of [START_REF] Giacalone | The Transport of Cosmic Rays across a Turbulent Magnetic Field[END_REF] and constructed as a sum of finite number of randomly oriented plane waves. Wavenumbers were evenly distributed in logarithmic scale between k min and k max . The number of modes was set to 500 and (k max /k min ) = 100. The exact value of k max depends on magnetic field strength and is determined from (5.1) and (5.2) to properly fit the correlation length. We have checked that in the small scattering regime we obtain the theoretically expected result from the random walk process with a correlation length equal to λ B . Fig. 5.6 shows the result of the calculation of the extended emission pattern at different energies for a range of IGMF strengths. The extended emission signal appears on top of the point-source signal and of the residual charged cosmic-ray background. The difficulties of detection of the extended emission are evident from the figure. The extended source signal is always subdominant. The point-source signal dominates in the 0 • < θ < 0.17 • angular range, while the residual cosmic-ray background dominates outside θ c = √ 0.03 0.17

Results and discussion

To investigate the detectability of the extended emission signal with CTA, we have calculated the statistics of the point-source signal, extended emission signal, and residual cosmic-ray background in each angular bin of the histograms shown in Fig. 5.6, for different CTA exposures. In this way we have generated mock CTA datasets. We have fitted the mock data with a model of point source plus the residual cosmic-ray background model, ignoring the presence of the extended source. We have then estimated the significance of detection of the extended emission in the simulated data set by calculating the level of inconsistency of the "point source + residual cosmic ray electron background [START_REF] Kerszberg | The cosmic-ray electron spectrum measured with H[END_REF]" model with the simulated data.

The results of this analysis are shown in Fig. 5.8. If the IGMF strength is below 3 × 10 -12 G, the extended signal is detectable with a significance larger than 3σ in the energy ranges 0.5-1.5 TeV and 1.5 -4 TeV. A somewhat stronger magnetic field, 3×10 -12 < B < 6×10 -12 G is still marginally detectable through the extended emission at somewhat higher energy, up to 10 TeV. The extended emission signal associated with the 10 -11 G IGMF is not detectable in a 50 hr exposure.

Right panel of Fig. 5.8 shows the exposure needed for the 3σ evidence for and 5σ discovery of the extended emission for different IGMF strengths. From this figure one can see that with a 150 hr exposure, an evidence for the presence of extended emission in the 10 TeV energy range can be found even for the magnetic field with the strength of 10 -11 G. Extended emission shaped by deflections of electrons and positrons in such IGMF can be discovered in such a very long exposure, T 400 hr.

To summarize, we have shown that direct detection of the strong cosmological magnetic field, which is needed for the resolution of the tension between different measurements of the Hubble parameter [START_REF] Jedamzik | Relieving the Hubble tension with primordial magnetic fields[END_REF] and for the explanation of the EDGES 21 cm line opacity [START_REF] Kumar Natwariya | EDGES signal in the presence of magnetic fields[END_REF] is possible with the γ-ray measurement technique. At least one source, Mrk 501, is suitable for this purpose. A search for the extended emission around Mrk 501 in the 1-10 TeV energy range can result in the detection of IGMF with a the strength of up to 10 -11 G and a correlation length of up to 10 kpc.

Sensitivity of CTA to IGMF with large correlation length

In this section we explore the possibility of identifying IGMFs with large correlation length and measuring its parameters with γ-ray astronomy methods. Such a magnetic fields may be a relic from the Early Universe and originate from the inflation. Large correlation length fields break isotropy by selecting a unique direction in a cosmologically large volume. This selected direction imposes a correlated asymmetry on magnetic field dependent extended emission patterns around γ-ray sources across the sky. We explore if it is possible to detect this asymmetry and measure the magnetic field direction.

Modelling of the magnetic field

We use a magnetic field generated from a realistic model of the LSS derived from Borg constrained cosmological simulations [START_REF] Jasche | Bayesian physical reconstruction of initial conditions from large-scale structure surveys[END_REF] that reproduces the location of known LSS elements (galaxies, clusters) in the local Universe. Our analysis relies on magnetohydrodynamic (MHD) simulations using the ramses code [START_REF] Teyssier | Cosmological hydrodynamics with adaptive mesh refinement. A new high resolution code called RAMSES[END_REF] to estimate the effect of the structure formation on the initial magnetic field configuration. Modelling the properties of secondary γ-ray signal is performed with the CRbeam code [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF], described in previous Chapter. We use calculations of electromagnetic cascades along lines of sight to known nearby blazars to estimate the influence of large correlation length magnetic field on the properties of secondary γ-ray signal from the cascade.

We [252] on the initial conditions (ICs) from Borg [START_REF] Jasche | Physical Bayesian modelling of the non-linear matter distribution: New insights into the nearby universe[END_REF]. The initial configuration of the magnetic field is uniform across the simulation volume with strength up to 10 -12 G. Its direction is chosen to be in the direction of north (Dec = 90 • ) of the ICRS coordinate system.

The simulation was run without cooling, star formation or feedback as the focus of this study are cosmological voids where the impact of processes are likely to be small [START_REF] Bondarenko | Constraining dark photons and their connection to 21 cm cosmology with CMB data[END_REF][START_REF] Bondarenko | Account of baryonic feedback effect in the γ-ray measurements of intergalactic magnetic fields[END_REF]. As a result of this the magnetic field is a smoothed tracer of the dark matter density field. Only sufficiently massive objects are capable of producing large enhancements of the gas density field. Regions with overdensities greater than approximately 3 were thought to be close enough to structures that generate magnetized outflows that they are at risk of being "polluted" by magnetic fields of non-cosmological origin, as indicated in [START_REF] Marinacci | First results from the IllustrisTNG simulations: radio haloes and magnetic fields[END_REF]. The magnetic field strength in such regions may deviate from simple scaling imposed by the amplification due to pure adiabatic contraction. As a result we therefore mark as high magnetic field region any area with an overdensity greater than 3 averaged over a thin cylinder with radius 2.7 Mpc (3 cells in the zoom region) along the line of sight during the course of our analysis.

There are several known TeV blazars whose position falls within the high resolution simulation volume. The source list is given in Table 5.1. All the sources are wellestablished TeV emitters. The proximity of the sources enables measurements of their intrinsic spectra attenuated by the pair production effect up to 10 TeV. This is important because the extended emission in the energy range above 100 GeV is produced by electrons and positrons injected in interactions of γ-rays with energies above 10 TeV. Measurement of the flux above 10 TeV allows us to obtain reliable estimates of the expected power of the secondary flux in the E > 100 GeV range.

The lines of sight toward the sources listed in Table 5.1 are oriented at different angles with respect to the direction of initial homogeneous field chosen in simulation. This provides a possibility to study how variations of the misalignment of the cosmological field direction affects the observational appearance of the extended emission. Lines of sight toward two sources, 1ES 1959+650 and Mkn 180 are more strongly aligned to the direction of the magnetic field than, say, the lines of sight to Mkn 421 and Mkn 501. The line of sight toward AP Librae has the largest misalignment angle in our sample.

In the following sections we show results for the three brightest sources in our sample Mkn 501, Mkn 421 and 1ES 1959+650, see Table 5.1.

Structure of magnetic fields along the lines of sight

Figure 5.9 shows the line-of-sight profiles of magnetic field toward Mkn 501, Mkn 421 and 1ES 1959+650. The top panels of the figure shows the line-of-sight profiles of the density field. Regions with overdensities above 3 that are excluded from our further analysis because of possible "contamination" of the magnetic field by baryonic feedback [START_REF] Bondarenko | Account of baryonic feedback effect in the γ-ray measurements of intergalactic magnetic fields[END_REF] are indicated with vertical red bands. In those regions secondary electrons are randomized by high magnetic fields, they emit secondary photons in random directions and their signal in direction to observer can be safely neglected. Grey bands indicate the distance to the blazars. For all three examples, the γ-ray sources are situated in moderately overdense regions that are unlikely to be affected by strong magnetised outflows from galaxies. None of the three sources have a line-of-sight aligned with a filament of the LSS. This makes all the three sources suitable for the IGMF measurement.

The middle and bottom panels of the figure shows the strength and orientation of the IGMF along the lines of sight. One can see clear differences of the IGMF structure for the three sources. The component of the magnetic field perpendicular to the line of sight is strongly suppressed in the case of 1ES 1959+650, because of the alignment of the cosmological field direction with the line of sight.

The lower panels of the figure show the position angle of the perpendicular to the line-of-sight field component, counted from the direction toward North. Deviations from the straight horizontal lines in these panels are due to the effect of structure formation. One can see that the effect is moderate and the initial cosmological field direction is mostly preserved in the course of structure formation.

Modelling of secondary extended gamma-ray emission

We model the secondary cascade γ-ray signal with the CRbeam Monte Carlo code [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF] which propagates high energy γ-rays through the cosmic medium taking into account most important physical interactions: pair production by γ-ray absorption and inverse Compton scattering of secondary electrons and positrons on the EBL [START_REF] Franceschini | Extragalactic optical-infrared background radiation, its time evolution and the cosmic photon-photon opacity[END_REF] and CMB. Electrons and positrons are also deflected in magnetic field whose strength is obtained from IGMF model described in the previous section.

For intrinsic point source spectrum we assume broken power law spectrum with break position around 100 GeV and maximum energy 100 TeV. All primary γ-rays have the same direction of initial momenta which coincides with the direction from the source to the observer. The γ-rays are propagated until they reach the sphere whose center coincides with the position of the source and radius is equal to the distance to the observer. To model primary γ-ray emission into a jet with an opening angle α jet we select γ-rays whose positions on the sphere lie inside the cone with the opening angle α jet = 5 • and direction of the axis coincides with the direction of jet. We do not take into account secondary photons that were created in the regions with an overdensity above 3, see Fig. 5.9. To explore different strengths of seed magnetic field we rescale the magnetic field profile by a constant factor keeping the shape of the profile unchanged.

Figs. 5.10 shows the images of the extended emission calculated for the field strength 10 -13 G for the three brightest sources in our sample. The images of extended emission around Mkn 421 and Mkn 501 reveal a clear overall azimuthal asymmetry that is imposed by the orientation of the initial homogeneous field projected on the sky. To the contrary, the image for 1ES 1959+650 appears halo-like. This is explained by the alignment of the magnetic field along the line of sight.

The images are superpositions of line-like features that are not perfectly aligned with each other. If the IGMF were perfectly homogeneous, with fixed direction, all these features would be strictly aligned, so that "narrow-jet-like" extensions would be observable. The appearance of multiple misaligned features reveals the influence of the structure formation on the field. Each feature is produced by the γ-ray beam passing through individual voids along the line of sight. In each void collective motions of the medium induce different deformations of the magnetic field. This results in the changes of directions of deflection of electrons and positrons injected in the voids by the γγ pair production.

Jet-like, rather than halo-like extended emission is generically expected even for short correlation length IGMF (originating from a cosmological phase transition) [START_REF] Neronov | Degree-scale GeV "Jets" from Active and Dead TeV Blazars[END_REF]. Thus, the azimuthal asymmetry of the extended signal is not an indication of the presence of large correlation length IGMF. However, there are two important differences in the properties of the jet-like extended emission in the two alternative cases of phase transition and inflationary magnetic field. The jet-like extensions are one-sided in the case of a short correlation length field [START_REF] Neronov | Degree-scale GeV "Jets" from Active and Dead TeV Blazars[END_REF], whereas the jet-like extensions are two sided for Mkn 421 and Mkn 501 in Fig. 5.10.

In addition, in the case of the phase transition IGMF, the position angle of the one-sided jet-like extended emission is determined by the orientation of the jet of the blazar rather than by the direction of the field. It is expected to be random for different sources, whereas in the case of the inflationary large correlation length field the extended sources associated to different blazars have the same position angle.

If the blazar jet is closely aligned with the line of sight, the two sides of the extended emission are symmetric, they have comparable surface brightness. This symmetry can be broken by the misalignment of the jet with the line of sight. This effect is illustrated in Fig. 5.11. Angles between the blazar jet, the line of sight and the direction of CMF regulate the flux levels of the extended emission components on different sides.

The orientation of the two-sided extension that is due to the presence of the large correlation length field is, always perpendicular to the direction of magnetic field projected on the sky. If the cosmological IGMF is correlated on the distance scales larger than hundreds of Megaparsecs, orientations of the two-sided extensions around different sources all over the sky are all expected to be aligned. Measurement of aligned extended emission features around multiple sources can provide an unambiguous evidence for the presence of magnetic field of inflationary origin. In this case position angles of the two-sided extensions provide the measurement of the direction of the primordial field. Fig. 5.12 shows spectra of the point and extended source emission for different assumptions about the magnetic field strength. For each source we choose the model of the intrinsic spectrum in such a way that it fits the spectral measurements of the "low state" of the source and/or the time-averaged spectral measurements on multi-year time scales. For all sources we consider two alternative Fermi spectral data: the spectral measurements reported in the Fermi source catalog [START_REF] Abdollahi | F ermi Large Area Telescope Fourth Source Catalog[END_REF] and the spectra extracted using the aperture photometry approach for the full 13-year exposure of Fermi/LAT up to July 2021. This allows us to extend the Fermi catalog measurements up to the 3 TeV energy. 3For Mkn 501 we used MAGIC and VERITAS analyses from 2009 [START_REF] Abdo | Insights Into the High-Energy Gamma-ray Emission of Markarian 501 from Extensive Multifrequency Observations in the Fermi Era[END_REF], and long term observations from ARGO [START_REF]Long-term Monitoring of Mrk 501 for its Very High Energy γ Emission and a Flare in 2011 October[END_REF] and HAWC [START_REF] Coutiño De Leon | Spectral analysis of the blazars Markarian 421 and Markarian 501 with the HAWC Gamma-Ray Observatory[END_REF]. For Mkn 421 we used long term observations from ARGO [START_REF] Bartoli | 4.5 Years of Multi-wavelength Observations of mrk 421 During the Argo-ybj and Fermi Common Operation Time[END_REF] and HAWC [START_REF] Coutiño De Leon | Spectral analysis of the blazars Markarian 421 and Markarian 501 with the HAWC Gamma-Ray Observatory[END_REF] and very low state from VERITAS measurement [START_REF] Acciari | TeV and Multi-wavelength Observations of Mrk 421 in 2006-2008[END_REF].

.12 -Spectra of intrinsic point source and extended emission for different assumed magnetic field strengths (upper panels: 0, 10 -14 and 3 × 10 -14 G, lower panels: 10 -13 , 3 × 10 -13 and 10 -12 G) , for the three brightest sources in the sample. Point source spectra from Fermi LAT, MAGIC [START_REF] Abdo | Insights Into the High-Energy Gamma-ray Emission of Markarian 501 from Extensive Multifrequency Observations in the Fermi Era[END_REF], VERITAS [START_REF] Abdo | Insights Into the High-Energy Gamma-ray Emission of Markarian 501 from Extensive Multifrequency Observations in the Fermi Era[END_REF] ARGO [START_REF]Long-term Monitoring of Mrk 501 for its Very High Energy γ Emission and a Flare in 2011 October[END_REF] and HAWC [START_REF] Coutiño De Leon | Spectral analysis of the blazars Markarian 421 and Markarian 501 with the HAWC Gamma-Ray Observatory[END_REF]. Total flux from the cascade in case B = 0 is shown with dotted line. Extended emission outside of the PSF of point source is shown with colored lines for magnetic fields between 10 -14 G and 10 -12 G. Dashed green line show sensitivity of CTA to point sources.

Results

We extracted the extended source fluxes from wedges of the angular width 0.3 • that contain the signal. From the wedges, we excluded the regions in which the primary source emission dominates, namely, we did not take into account those regions of the wedges that are located at a distance from the centers less than the angular resolution of the CTA. From Fig. 5.10 one can judge that the extended emission is dominated by the flux in a much narrower wedge, with opening angle of just a few degrees.

This gives a conservative estimate of the detectable extended source flux. In principle, stronger signal can be extracted using the azimuthal asymmetry of the extended emission. If the shape of the point spread function is known, the extended emission can be detected also from within the extent of the point spread function. Systematic uncertainty of the knowledge of the point spread function limits the sensitivity of Cherenkov telescopes for detection of extended emission [START_REF] Neronov | No evidence for gamma-ray halos around active galactic nuclei resulting from intergalactic magnetic fields[END_REF][START_REF] Vela | Study of the IACT angular acceptance and Point Spread Function[END_REF]. This uncertainty can be mitigated by the direct measurement of the point spread function in the direction perpendicular to the direction of extended emission. This is possible even if the source line of sight is nearly aligned with the field direction, as illustrated by the example of 1ES 1959+650, shown in the right panel of Fig. 5.10.

The linear, one-dimensional, rather than circular, two-dimensional shape of the extended source also provides an improvement of sensitivity for another reason. The sensitivity for extended sources typically worsens with the increase of the solid angle spanned by the extended source on the sky, because of the increase of the background on top of which the signal is detected. Concentration of the extended signal in a narrow wedge reduces the solid angle and hence increases the signal-to-noise ratio.

All the sources in our sample have steep spectra in the 10 TeV range. This diminishes the power of the secondary emission observable in the 0.1-1 TeV range. The total flux that would be available for detection (as a contribution to the point source flux) in the absence of an IGMF is shown by the brown dotted lines in Fig. 5.12. The IGMF deflects electrons and positrons away from the line of sight, so that the secondary γ-rays produced by the inverse Compton scattering of the Cosmic Microwave Background becomes unobservable. This effect is stronger at lower energies. This explains the suppression of the extended flux (that we collect from the wedge-shaped regions containing the secondary flux) at the energies below 100 GeV. If the magnetic field is too weak to deflect electrons and positrons, the secondary flux just contributes to the point source flux and is also undetectable. This explains the suppression of the secondary signal at higher energies. Both low-and high-energy suppression depends on the strength of magnetic field. Fig. 5.12 shows a comparison of the expected secondary flux levels with the sensitivity of the CTA telescopes for detection of extended sources. This comparison shows how challenging the search for the two-sided jet-like extensions might be. The model predictions for the extended signal are at the limit of sensitivity for the 50 hr exposure of each of the three brightest sources in our source sample. However, use of possible improvements of the method with "in-situ" measurement of the telescope point spread function and careful choice of the wedge for the extended signal measurements may improve the sensitivity. Otherwise, much longer exposure of several hundred hours (instead of 50 hr considered here) can also boost signal-to-noise ratio and make the correlated extended emission signal measureable in all the three brightest nearby blazars considered above.

Account of the baryonic feedback effect in γ-ray measurements of IGMF

As we discussed on the beginning of the Chapter, it is possible that the magnetic fields generated and amplified in galaxies were thrown out of them with outflows caused by the baryonic feedback processes arising during star formation or AGN activity. Significant progress toward a better modeling of the baryonic feedback process has recently been achieved by IllustrisTNG cosmological simulations [START_REF] Nelson | First results from the IllustrisTNG simulations: the galaxy colour bimodality[END_REF][START_REF] Volker | First results from the IllustrisTNG simulations: matter and galaxy clustering[END_REF][START_REF] Pillepich | First results from the IllustrisTNG simulations: the stellar mass content of groups and clusters of galaxies[END_REF][START_REF] Naiman | First results from the IllustrisTNG simulations: a tale of two elements -chemical evolution of magnesium and europium[END_REF][START_REF] Marinacci | First results from the IllustrisTNG simulations: radio haloes and magnetic fields[END_REF]. Recently, [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF] have used IllustrisTNG to show that AGN and supernova-driven outflows create extended magnetized bubbles with magnetic fields of about B > 10 -12 G and typical sizes in the range of ∼ 10-30 Mpc, with magnetic fields whose parameters are largely independent of those of the preexisting seed magnetic fields.

In what follows, we use the IllustrisTNG model to assess the effect of contamination of the intergalactic medium by the magnetic field that is spread by the galactic winds and its consequences for the measurements of cosmological IGMF. As discussed in [START_REF] Pakmor | Magnetic field formation in the Milky Way like disc galaxies of the Auriga project[END_REF], the magnetic fields produced by galactic activity are almost insensitive to the initial seed fields. The reason for this is that inside galaxies, the magnetic fields are strongly amplified by dynamo processes, saturate, and "forget" their initial conditions. We consider several characteristics situations of VHE γ-ray sources at different distances, with different γ-ray spectra, situated inside and outside the galactic wind bubbles. We find that within the IllustrisTNG baryonic feedback model, the basic features of the extended or delayed γ-ray signal sensitive to the primordial IGMF are just slightly altered, typically by 10-20%. This does not preclude the possibility of measuring the primordial field in the voids, but it introduces an additional systematic uncertainty that needs to be properly taken into account in the analysis with the future instruments.

Magnetized bubbles in IllustrisTNG simulations

IllustrisTNG is a set of gravo-magnetohydrodynamic simulations, the next generation (TNG) of the Illustris project [START_REF] Nelson | First results from the IllustrisTNG simulations: the galaxy colour bimodality[END_REF][START_REF] Volker | First results from the IllustrisTNG simulations: matter and galaxy clustering[END_REF][START_REF] Pillepich | First results from the IllustrisTNG simulations: the stellar mass content of groups and clusters of galaxies[END_REF][START_REF] Naiman | First results from the IllustrisTNG simulations: a tale of two elements -chemical evolution of magnesium and europium[END_REF][START_REF] Marinacci | First results from the IllustrisTNG simulations: radio haloes and magnetic fields[END_REF]. It is based on the moving-mesh Arepo code [START_REF] Volker | E pur si muove: Galilean-invariant cosmological hydrodynamical simulations on a moving mesh[END_REF] that solves the system of equations for self-gravity and ideal magnetohydrodynamics [START_REF] Pakmor | Magnetohydrodynamics on an unstructured moving grid[END_REF][START_REF] Pakmor | Simulations of magnetic fields in isolated disc galaxies[END_REF]. Cosmological parameters are chosen to be equal to the best-fit Planck 2015 cosmological parameters [START_REF] Ade | Planck 2015 results. XIII. Cosmological parameters[END_REF].

We used the highest-resolution simulations that are publicly available TNG100-1 (TNG100) [START_REF] Nelson | The IllustrisTNG simulations: public data release[END_REF]. The TNG100 simulation has a box size ∼ (110 cMpc) 3 and contains 1820 3 dark matter particles and an equal number of initial gas cells with masses of m DM = 7.5 × 10 6 M and m bar = 1.4 × 10 6 M .

The initial seed magnetic field in this run was a constant magnetic field with a field strength 10 -14 cG (comoving Gauss) directed along the z -axis in the simulation box. 4 This seed magnetic field experienced adiabatic contraction during the structure formation process, and then it was strongly amplified by small-scale dynamos in collapsed structures.

The TNG simulations adopt a comprehensive galaxy and supermassive black holes (SMBH) formation and feedback model [START_REF] Weinberger | Simulating galaxy formation with black hole driven thermal and kinetic feedback[END_REF][START_REF] Pillepich | Simulating galaxy formation with the IllustrisTNG model[END_REF]. SMBHs with an initial mass of ∼ 10 6 M are placed in the gravitational potential minima of dark matter halos when the virial mass exceeds ∼ 7 × 10 10 M . These black holes grow via binary mergers with each other or via smooth gas accretion according to the Bondi-Hoyle-Lyttleton model [START_REF] Weinberger | Simulating galaxy formation with black hole driven thermal and kinetic feedback[END_REF]. The rate of black hole growth depends on the black hole mass, local gas density, and relative velocity between the black hole and its surroundings.

Within the IllustrisTNG model, the SMBH activity has two different regimes that correspond to high-accretion and low-accretion states [START_REF] Weinberger | Simulating galaxy formation with black hole driven thermal and kinetic feedback[END_REF]. In the high-accretion rate state (accretion rate above ∼ 10% of the Eddington limit), the energy is continuously injected into the surrounding gas, causing thermal heating of the gas. In the lowaccretion rate state, the kinetic energy is deposited periodically when enough energy accumulates via accretion (see [START_REF] Weinberger | Simulating galaxy formation with black hole driven thermal and kinetic feedback[END_REF] for additional details). Each injection event creates a randomly oriented high-velocity kinetic wind. If averaged over long time intervals, energy injection in this mode also becomes isotropic. The low accretion rate feedback mode drives the most powerful outflows in the TNG model [START_REF] Nelson | First results from the TNG50 simulation: galactic outflows driven by supernovae and black hole feedback[END_REF]. A large fraction of the injected kinetic energy in this mode thermalizes via shocks in the surrounding gas, thereby providing a distributed heating channel.

Magnetic bubbles

According to the TNG simulations, magnetized bubbles (regions larger than the galaxy size with magnetic fields stronger than 10 -12 G, which is orders of magnitude stronger than the initial seed field) of chemically enriched material emerge at z 2 as a result of active baryonic feedback from galaxies (see Section 3 in [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF][START_REF] Aramburo | Ultra-high energy cosmic rays deflection by the Intergalactic Magnetic Field[END_REF] and Fig. 1 therein for the visualization of magnetic bubbles). These bubbles are significantly nonspherical and can have sizes of several dozen megaparsec. [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF] demonstrated that at z = 0, these bubbles occupy 12 -14% of the cosmological volume for magnetic fields 4. Different configurations of initial magnetic fields were studied in many previous works, see, e.g., [START_REF] Vazza | Simulations of extragalactic magnetic fields and of their observables[END_REF] and reference therein. stronger than 10 -12 G and more than 3% for B > 10 -9 G. Comparing the simulations with initial magnetic fields that are different by many orders of magnitude (otherwise identical), [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF][START_REF] Aramburo | Ultra-high energy cosmic rays deflection by the Intergalactic Magnetic Field[END_REF] have shown that magnetic fields in the bubbles are independent of the strength of the initial seed field to a large extent. Fig. 5.13 illustrates this effect: orientations of strong magnetic fields with B > 10 -12 cG "forgets" the direction of the seed field (which was along the z-axes in the TNG simulations). In contrast, weak fields in the "real" IGM outside the bubbles prefer the initial direction. It was also demonstrated in [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF] that the bubbles are mainly caused by the AGN activity (with a smaller contribution by supernovae). We took the magnetic field value 10 -12 cG as a boundary value for the magnetic bubbles. 5

Line-of-sight data

We are interested in the effects of the bubbles on the extended and delayed γ-ray signal from distant AGN. To assess these effects, we extracted the information on the strength and orientation of the magnetic field along randomly chosen straight lines from the TNG-100 simulation, so-called lines of sight (LOS), originating from seven 5. In Fig. 1 in [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF],the boundaries of bubbles are quite sharp and there is no significant difference in the volume filling fraction for magnetic fields stronger than 10 -12 cG or 10 -11 cG (see Fig. 4 therein). Therefore the results of this paper do not significantly depend on this specific value.

Numerical modeling

To explore the effects of magnetic bubbles on the secondary emission, we numerically modeled the propagation of the high-energy gamma-rays through the cosmic medium using data on a magnetic field for a set of LOSs, extracted from the TNG100 simulation. For this purpose, we used the recently developed Monte Carlo simulation code CRbeam [START_REF] Berezinsky | High energy electromagnetic cascades in extragalactic space: physics and features[END_REF] that takes gamma-ray absorption with pair production on the EBL, subsequent inverse Compton scattering of electrons and positrons on the CMB and EBL, and deflection of the charged component in magnetic fields into account. This code was also tested through a comparison with several alternative cascade codes [START_REF] Taylor | Extragalactic magnetic fields constraints from simultaneous GeV-TeV observations of blazars[END_REF][START_REF] Kalashev | Simulations of Ultra High Energy Cosmic Rays propagation[END_REF][START_REF] Kachelriess | ELMAG: A Monte Carlo simulation of electromagnetic cascades on the extragalactic background light and in magnetic fields[END_REF].

In all simulations, the position of the source corresponds to the first voxel of the LOS. The directions of the initial momenta of all particles are the same and coincide with the direction along the LOS. The strength and direction of magnetic fields are given by linear interpolation between neighboring voxels. The secondary γ-ray signal generated by every electron and positron only depends on the strength and orientation of the magnetic field right at the location of electron-positron pair production (more precisely, within the voxel in which the pair production happens). This appears to be reasonable in the view of the fact that electrons with energies E e < 5 TeV (which produce secondary gamma-rays with E γ < 100 GeV) are typically deflected by δ > 3 • already within the first 20 kpc of their trajectory. Because we are only interested in secondary γ-rays that are emitted within the angle that is comparable to the point spread function (PSF) of γ-ray telescopes, θ psf 1 • , further γ-ray production by electrons or positrons that are deflected by more than several degrees does not contribute to the secondary γ-ray signal that is sampled by the telescope situated on the LOS.

In order to separately study the effect of bubbles and voids on the cascade component, we constructed three alternative magnetic field profiles for each LOS. One profile contained only voids: in all voxels where B ≥ 10 -12 G, we set B = 0 G. The second profile contained only bubbles: the magnetic field in all voids with B < 10 -12 G and in the first bubble around the source AGN was set to B = 0 G. The last profile contained only the AGN source bubble: the magnetic fields everywhere outside the first bubble were set to B = 0 G. To avoid the effect of random fluctuations of the magnetic field and stabilize the procedure of determining the size of the first bubble better, we defined the boundary of the first bubble as the point from which the magnetic field strength is below 10 -12 G for at least L b = 300 kpc. We verified that if we replaced L b with 100, 500, or 1000 kpc, the size of most of the first bubbles did not change, and if it did, the changes were insignificant.

We propagated γ-rays from the source to the observer for each magnetic field profile and calculated the secondary γ-ray signal along the LOS. For LOSs containing only the first bubble, potentially more significant fluctuations of the secondary flux may arise because all primary photons are emitted in the same direction. In reality, the primary γ-rays are emitted in a cone with an opening angle α jet so that they sense different directions in the bubble. The bubble properties in various directions are different (size and magnetic field profile). In our analysis, we neglected this difference, and the bubble radius changes probably only insignificantly on the scale of the jet opening angle.

We explored several possible spectral models of the AGN source. The intrinsic γ-ray spectrum of the source was always a cutoff power law dN/dE ∝ E -Γ exp(-E/E cut ). The slope of the power law is Γ = 1.7. In contrast, E cut was a free parameter. Another free parameter was the distance to the observer D s which sets the redshift of the source z. We scanned over parameter space choosing E cut in the 1-100 TeV range and 0.025 < z < 0.15. Cosmological parameters were set to H 0 = 70 km/(s Mpc), Ω m = 0.3, and Ω Λ = 0.7, and the EBL model was fixed to that of [START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF].

We propagated γ-rays until they reached the sphere with radius r = D s . Although all primary particles are emitted in the same direction, we were able to model primary photon emission into a jet with an opening angle α jet by recording the γ-rays that hit the part of the sphere lying inside the cone with the opening angle α jet = 5 • , whose axis coincides with the direction from the source to the observer. For the spectral plots presented in the next section, we retained only the γ-rays that arrived at the position of the telescope at an incidence angle within the telescope PSF, which we assumed to be θ psf = 0.3 • .

Results

The simulated spectra of the primary AGN source and secondary γ-ray emission are shown in Figures 5.15, 5.16, and 5.17. They represent three types of magnetic field profiles and four different choices of parameters of the AGN source: a nearby source with low energy cutoff, a nearby source with high energy cutoff, a distant source with low energy cutoff, and a distant source with high energy cutoff.

The three sets of figures show that the effect of the void IGMF on the secondary γ-ray signal (Fig. 5.15) is clearly different from the effect of the bubble magnetic fields (Figs. 5.16 and 5.17). The volume-filling IGMF outside bubbles produces complete suppression of the secondary γ-ray emission at E < 100 GeV. As the initial conditions of TNG100 contain a constant magnetic field with a strength of B = 10 -14 G, this part of the IGMF, unaffected by galactic dynamos, has a long correlation length and B > 10 -15 G almost everywhere (see Fig. 2 in [START_REF] Aramburo | Magnetization of the intergalactic medium in the IllustrisTNG simulations: the importance of extended, outflow-driven bubbles[END_REF]). Therefore, such an effect on the secondary γ-ray emission is not unexpected, based on the analytical estimates of [START_REF] Neronov | A method of measurement of extragalactic magnetic fields by TeV gamma ray telescopes[END_REF][START_REF] Neronov | Sensitivity of gamma-ray telescopes for detection of magnetic fields in intergalactic medium[END_REF].

The magnetized bubbles produce a qualitatively different effect on the secondary γ-ray flux. Strong (B > 10 -12 G) magnetic fields in the bubbles sufficiently deflect the charged particles that are created in the bubbles and therefore completely suppress their secondary γ-ray flux at the detector throughout its energy range. This changes the overall normalization of the secondary emission, but affects its shape only very little.

This qualitative difference between the effects of void IGMF and bubbles holds for the whole range of properties of the primary γ-ray sources. Only the strength of the overall flux suppression by the bubbles depends on the source parameters.

The level of secondary flux suppression clearly depends on the length of the bubbles along the line of sight from the sources toward the observer. Fig. 5.18 shows the lost secondary flux in percentage points for all the simulated LOSs. Within the TNG100 feedback model, the magnetized bubbles typically remove some 10% of the secondary γ-ray flux if the source intrinsic γ-ray spectrum has a high energy cutoff at relatively low energy E cut ∼ 5 TeV. However, for a fraction of the LOSs, the pass of the photons through the initial bubble happens to be large (above 10 Mpc), and loss of the secondary emission for these cases is two to three times larger.

The average flux suppression increases by up to 20%-50% for hypothetical sources with hard intrinsic spectra extending up to E cut = 100 TeV. Even though there is currently no observational evidence for the existence of such extreme AGNs, we show the result for this case to highlight the dependence of the bubble effect on the intrinsic properties of the primary γ-ray source.

The secondary flux suppression by the bubbles is accumulated in the AGN host bubble and all other bubbles along the line of sight. The relative importance of the "host" and "all other" bubbles also depends on the properties of the primary γ-ray source. This is illustrated in Fig. 5.18, which also shows the percentage of the lost secondary flux as a function of the bubble type.

The suppression by "all other" bubbles except for the AGN host bubble is at the level 5-20%, and it does not depend on the distance to the source or on the cutoff energy. In this case, the suppression factor is well described by the fraction of the LOS length occupied by the magnetized bubbles.

One singular LOS in Fig. 5.18 was affected more strongly by the bubbles (50-60%), despite a modest (5 Mpc) pass in the initial bubble. The reason for this phenomenon is clear from Fig. 5.14: in this case, the LOS passes not only through the initial bubbles GeV that is removed by the bubbles as a function of the size of the host AGN bubble. Red circles show the flux fractions that are removed by all bubbles. Black crosses and blue squares show the flux fractions that are removed by all bubbles, except for the host AGN bubble and a linear fraction of the LOS that is occupied by these bubbles. The dashed blue horizontal line marks the average linear fraction of the LOS that is occupied by all bubbles, except for the host AGN bubble. around the source, but also through a large filamentary system of magnetized bubbles that occupy more than half of the first 100 Mpc of this LOS. Although this seems to be a rare combination of the location of the source and the orientation of the LOS with respect to the filamentary system of bubbles, this possibility can not be dismissed for any particular source (at least according to the IllustrisTNG model).

The AGN host bubble has practically no effect on the secondary source flux if the cutoff energy E cut is low. However, it has a strong impact on the signal from sources with high cutoff energy. In this case, the mean free path D mfp (E cut ) of the gammarays with energy E cut becomes smaller than the size of the first bubble, so that a significant fraction of primary γ-ray flux is absorbed inside this AGN host bubble. This significantly increases the percentage of lost secondary emission, as shown in Fig. 5.18. The effect is more substantial for higher cutoff energies, and it is slightly weakened for further away sources through the additional contribution to the secondary flux from longer propagation to the observer. For nearby sources with high cutoff energy, additional losses in the first bubble can exceed 30%.

The range of the primary source parameters for which the losses in the first bubble are perceptible is shown in Fig. 5.19. Two source classes are clearly visible. Nearby sources with low cutoff energies are practically unaffected by the host AGN bubbles. This is true as long as the mean free path of gamma-rays D mfp (E cut ) exceeds the distance to the observer D s , so that the suppression factor is equal to the fraction of the LOS occupied by the bubble. This explains the "horizontal branch" of the contours in Fig. 5.19. As soon as the mean free path of the primary γ-rays becomes shorter than the distance to the source, the fraction of the flux that is removed by the AGN host bubble is given by the ratio of the bubble size to the primary γ-ray mean free path D mfp (E cut ). This effect is visible toward the higher energy "inclined" branch of the contours in Fig. 5.19. This is an important conclusion that can be drawn from Fig. 5.19. The presence of magnetized bubbles in the LSS produces only a mild effect on the secondary γ-ray flux at the 10% level for all reasonable AGN source parameters. In this respect, the bubbles have to be taken into account as an additional source of systematic uncertainties in the analysis of the IGMF in voids.

The results presented above show that an account of the magnetic fields produced by the baryonic feedback on the LSS modifies model calculations of the secondary γray flux from the electromagnetic cascade that develops along the line of sight. The strength of the effect depends on the parameters of the primary γ-ray source.

In the specific feedback model of the IllustrisTNG simulation, the effect is always at the level of ∼ 10% for a wide range of AGN primary source parameter choices: cutoff energies below ∼ 30 TeV, and redshifts within z ∼ 0.15. In this case, the possible presence of magnetized bubbles along the line of sight has to be considered as an additional source of systematic uncertainty of the overall level of the secondary flux. This systematic uncertainty does not affect the spectral shape of the secondary γ-ray flux. In general, it does not affect the timing and extended emission profiles of z = 0.025 the signal either. This implies that previously developed analysis methods leading to constraints on the strength and correlation length of the void IGMF remain valid even in the presence of the magnetized bubbles.

Conclusions of the Chapter

In this Chapter, we explored the possibilities of detecting strong extragalactic magnetic fields with the CTA telescope. The blazar Mrk 501 has been shown to be the most suitable candidate for studying magnetic fields in the range of 10 -12 G -10 -11 G. It was also shown that the sources Mrk 421, Mrk 501 and 1ES 1959+650 can be used to detect extragalactic magnetic fields with large correlation lengths.

On the other hand, we studied the effect of magnetic field bubbles around galaxies and clusters of galaxies on the intensity of cascade signal. It was found that typical magnetic field bubbles corresponding to the IllustrisTNG model cause an energy independent suppression of the secondary gamma-ray flux at a level of about 10%.

Chapter 6 Conclusion

In this dissertation, we considered the capabilities of gamma-ray astronomy for studying extragalactic background light, axion-like particles, and also intergalactic magnetic fields.

In the first part of the thesis we constructed new flexible EBL model [START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF]. The advantage of this model is that unlike other models it allows one to modify the parameters underlying it. Using this model, we have identified the most important astrophysical parameters of the EBL and written them out explicitly. Using the Monte Carlo Markov chains and on the basis of this new EBL model, we have developed a method for constraining astrophysical parameters based on available EBL measurements. For some parameters, the limits obtained with modern EBL measurements are comparable to the limits obtained by other methods. Nevertheless, the application of this approach to more accurate future EBL measurements can significantly increase the sensitivity of the method.

On the other hand we compared the visible and infrared EBL measurements with the gamma-ray constraints and study if the discrepancies of the measurements with different methods could be due to the presence of features in the EBL spectrum that are localised in the micron wavelength range. We showed that a previously reported "excess" in EBL flux in the 1 micron wavelength range is consistent with gamma-ray measurements, that is, if the excess has the form of a narrow feature of the width δλ < λ and an overall flux of up to 15 nW/m 2 /sr above the "minimal" EBL, which is estimated from the visible and infrared source counts. Such bump-like spectral features could originate, for example, from decaying dark-matter particles, or either axions or peculiar astrophysical processes in the course of star-formation history. We showed that future measurements with CTA can answer the question of the presence of such a features in the EBL spectrum [START_REF] Korochkin | Search for spectral features in extragalactic background light with gamma-ray telescopes[END_REF].

In the next Chapter we examined the influence of axion-like particles on the highenergy spectra of blazars. We found 95% C.L. limit on axion-like particle mass and coupling constant from non-observation of axion decay bump in EBL spectrum, using gamma-ray spectrum of observations of the blazar 1ES 1218+304 by Fermi/LAT and VERITAS telescopes [START_REF] Korochkin | Search for decaying eV-mass axionlike particles using gamma-ray signal from blazars[END_REF]. Our best fit model favours non-zero amplitude of the ALP induced absorption feature in the γ-ray spectrum. We have discussed available possibilities to scrutinize this result with complementary observations with current and next generation ground-based γ-ray telescopes, CTA and LHAASO [START_REF] Addazi | The Large High Altitude Air Shower Observatory (LHAASO) Science Book (2021 Edition)[END_REF]. They will provide a factor of three-to-five improvement of sensitivity for the ALP decay feature searches. If the EBL spectrum has an ALP-induced bump with the flux at the level suggested by the direct EBL measurement experiments, CTA will be able to detect the ALP induced bump with high significance, and provide precision measurement of the ALP mass and photon coupling constant.

We also tested a possible indication of the existence of the ALPs from blazar observations, which is the appearance of upward breaks in the deabsorbed spectra. Based on the most complete sample of blazars with accurately measured redshifts, we have shown that the significance of the "anomalous transparency" effect is weaker than previously reported and does not exceed 2σ for the modern EBL models and even as low as 1.3σ for the most conservative EBL model [START_REF] Korochkin | Search for anomalous features in gamma-ray blazar spectra corrected for the absorption on the extragalactic background light[END_REF].

Next we carried out a detailed comparison of the electromagnetic modules of the publicly available Monte Carlo programs CRbeam, CRPropa and ELMAG [START_REF] Kalashev | Modelling of propagation of very-high-energy gamma rays with CRbeam code. Comparison with CRPropa and ELMAG codes[END_REF]. We compare the propagation of electrons and positrons in a magnetic field, generation of turbulent fields, the absorption of high-energy gamma rays due to interaction with EBL photons, as well as the process of inverse Compton scattering of electrons on background photons. We also compared energy distribution of secondary particles. When comparing the results with theoretical predictions, inaccuracies were found, some of which were corrected. After introducing corrections, all codes demonstrate agreement at the level of 10% which can be considered as an estimate of systematic modeling error.

Finally, we considered the possibilities of gamma-ray astronomy to study IGMF. A primordial magnetic field with correlation length of 0.1 -1 kpc the strength in the 10 -12 G -10 -11 G range can resolve the tension between different measurements of the Hubble constant and provide an explanation for the excess opacity in the 21 cm line at redshift 15 < z < 20, if it is present during the recombination and reionization epochs. We study the sensitivity reach of the gamma-ray technique for measurement of such relatively strong cosmological magnetic field using deep exposures of the nearest hard spectrum blazars with CTA telescopes. We show that the gamma-ray measurement method can sense the cosmological magnetic field with a strength of up to 10 -11 G [START_REF] Korochkin | Sensitivity reach of gamma-ray measurements for strong cosmological magnetic fields[END_REF][START_REF] Korochkin | Sensitivity of the gamma-ray method for strong primordial magnetic fields[END_REF]. Combination of the cosmic microwave background and gamma-ray constraints can thus sense the full range of possible cosmological magnetic fields to confirm or rule out their relevance to the problem of the origin of cosmic magnetic fields, as well as their influence on recombination and reionization epochs.

On the other hand we studied a possibility of detection of primordial magnetic field with large correlation length (> 100 Mpc) possibly originating from inflation. Such field can be coherent on cosmological scales and induce wedge-like extended emission around nearby blazars, with aligned wedge orientations in multiple sources across the sky. This alignment can be used to distinguish the inflationary origin of IGMF, because it is not expected if the IGMF originates from cosmological phase transitions [START_REF] Korochkin | Detectability of large correlation length inflationary magnetic field with Cherenkov telescopes[END_REF].

Also we have explored the effect of the baryonic feedback on the LSS in the form of outflows from galaxies that pollutes the intergalactic medium with magnetic fields. This process can spread the magnetic field in the voids of the LSS and thus preclude measuring the primordial magnetic fields that possibly reside in the voids. Our analysis of the specific baryonic feedback model implemented in IllustrisTNG cosmological simulations shows that the feedback generally does not strongly affect the void magnetic fields. The galactic outflows from magnetized bubbles collectively have a small volume filling factor. We have explored the effect of the magnetized bubbles on the γ-ray signal from distant AGN [START_REF] Bondarenko | Account of baryonic feedback effect in the γ-ray measurements of intergalactic magnetic fields[END_REF]. The effect of these bubbles can be taken into account in the analysis as a systematic uncertainty of the secondary γ-ray flux at the level of 10% on average (and up to 60% for a single source). This error has to be taken into account in the error budget of the analysis. This systematic uncertainty depends on the high energy cutoff of the primary source spectrum. If the cutoff energy is in the E cut = 100 TeV range, a dedicated analysis of the possible properties of the AGN host bubble has to be performed as part of the IGMF-related data analysis.

Figure 2 . 1 -

 21 Figure 2.1 -Spectrum of cosmic background radiation at all wavelengths. Combination of Cosmic Infrared Background (CIB) and Cosmic Optical Background (COB) represents EBL, blue and red arrows indicate lower and upper limits. Black points correspond to detections. Source: [23].

Figure 2 . 3 -

 23 Figure 2.3 -Spectrum of EBL based on the results of observations and simulations.The solid and dashed lines correspond to the numerical models of the EBL. The red dots are obtained from the galaxy number count method. The gray dots indicate direct measurements of the EBL. Yellow band corresponds to the constraints from the γ-ray observations. Source:[START_REF] Saldana-Lopez | An observational determination of the evolving extragalactic background light from the multiwavelength HST/CANDELS survey in the Fermi and CTA era[END_REF].

Figure 2 . 5 -

 25 Figure 2.5 -Distributions of the parameters of the clouds, including radius R c , lifetime η c , dust particle density n d , optical depth at 5500 Å τ λ 0 , optical depth slope n and cloud efficiency M 0 .

Figure 2 . 6 -

 26 Figure 2.6 -Distributions of the global parameters of star formation, including redshift at which star formation began z i , IMF truncation mass m max and SFR scaling factor C sfr

Figure 2 . 7 -

 27 Figure 2.7 -SED of the EBL obtained by various methods. Direct measurements: The purple arrows are results from AKARI[START_REF] Tsumura | Low-Resolution Spectrum of the Extragalactic Background Light with AKARI InfraRed Camera[END_REF], while green asterisks are from the reanalysis of IRTS[START_REF] Matsumoto | Reanalysis of the near-infrared extragalactic background light based on the IRTS Observations[END_REF]. Pink diamonds are from the reanalysis of COBE/DIRBE measurements[63] [64]. Black data points together with blue and orange systematic uncertainty are derived from CIBER[START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF] and correspond to nominal and minimum EBL models. Lower limits: Dark blue upward arrows combine EBL lower limits obtained by different experiments: GALEX[START_REF] Xu | Number counts of GALEX sources in FUV (1530A) and NUV (2310A) bands[END_REF], Hubble Deep Field[START_REF] Madau | Deep galaxy counts, extragalactic background light, and the stellar baryon budget[END_REF], Subaru[START_REF] Keenan | The Resolved Near-Infrared Extragalactic Background[END_REF], and Spitzer/IRAC[START_REF] Fazio | Number counts at 3 < lambda < 10 um from the Spitzer Space Telescope[END_REF]. EBL from γ-ray absorption: Striped lime, red, and purple bands are from MAGIC[START_REF] Acciari | Measurement of the extragalactic background light using MAGIC and Fermi-LAT gamma-ray observations of blazars up to z = 1[END_REF], HESS[START_REF] Abramowski | Measurement of the extragalactic background light imprint on the spectra of the brightest blazars observed with H.E.S.S[END_REF], and Fermi/LAT[START_REF] Abdollahi | A gamma-ray determination of the Universe's star formation history[END_REF], correspondingly. Modelling: Dark red solid line is for baseline EBL model of[START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. Black striped band shows the allowed range of EBL models obtained with the global fit to EBL by[START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF]. Dash-dotted orange and dashed red lines are additional wide and narrow Gaussian components.

Figure 3 . 1 -

 31 Figure 3.1 -Constraints of ALP parameters from laboratory and astrophysical probes.Yellow band and black solid line corresponds to the QCD axion models[START_REF] Kim | Weak-Interaction Singlet and Strong CP Invariance[END_REF][START_REF] Shifman | Can confinement ensure natural CP invariance of strong interactions?[END_REF][START_REF] Dine | A simple solution to the strong CP problem with a harmless axion[END_REF][START_REF] Zhitnitsky | On Possible Suppression of the Axion Hadron Interactions[END_REF]. Blue region is the limit imposed by CAST non-observation of solar axions[START_REF] Arik | Search for Solar Axions by the CERN Axion Solar Telescope with 3 He Buffer Gas: Closing the Hot Dark Matter Gap[END_REF]. Grey vertical band comes from high-resolution spectroscopy of intracluster medium in by optical telescope[START_REF] Grin | A Telescope Search for Decaying Relic Axions[END_REF]. Black dashed line is the EBL limit estimate of[START_REF] Arias | WISPy Cold Dark Matter[END_REF]. Green horizontal line is the limit imposed by non-observation of excess energy loss in Horizontal Branch stars[START_REF] Straniero | Axion-Photon Coupling: Astrophysical Constraints[END_REF]. Light grey shaded region shows 95% confidence level excluded range of parameters from non-observation of absorption feature in the Fermi/LAT+VERITAS spectrum of 1ES 1218+304 (discussed in this paper). Violet shaded ellipses show 68% confidence level preferred range of the parameters obtained by the spectral analysis of 1ES 1218+304. Red hatched region shows the favoured region found if the first spectral point of VERITAS is not taken into account.

Figure 3 . 2 -

 32 Figure 3.2 -Spectral energy distribution of the EBL. Direct measurements: purple arrows are results from AKARI[START_REF] Tsumura | Low-Resolution Spectrum of the Extragalactic Background Light with AKARI InfraRed Camera[END_REF], while green asterisks are from the reanalysis of IRTS[START_REF] Matsumoto | Reanalysis of the near-infrared extragalactic background light based on the IRTS Observations[END_REF]. Pink diamonds are from reanalysis of COBE/DIRBE measurements[START_REF] Sano | Derivation of a large isotropic diffuse sky emission component at 1.25 and 2.2um from the COBE/DIRBE data[END_REF][START_REF] Sano | Measurements of Diffuse Sky Emission Components in High Galactic Latitudes at 3.5 and 4.9 µm Using DIRBE and W ISE Data[END_REF]. Black data points together with blue and orange systematic uncertainty are derived from CIBER[START_REF] Matsuura | New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER[END_REF] and correspond to nominal and minimum EBL models. Lower limits: dark blue upward arrows combine EBL lower limits obtained by different experiments: GALEX[START_REF] Xu | Number counts of GALEX sources in FUV (1530A) and NUV (2310A) bands[END_REF], Hubble Deep Field[START_REF] Madau | Deep galaxy counts, extragalactic background light, and the stellar baryon budget[END_REF], Subaru[START_REF] Keenan | The Resolved Near-Infrared Extragalactic Background[END_REF], Spitzer/IRAC[START_REF] Fazio | Number counts at 3 < lambda < 10 um from the Spitzer Space Telescope[END_REF]. EBL from γ-ray absorption: striped lime, red and purple bands are from MAGIC[START_REF] Acciari | Measurement of the extragalactic background light using MAGIC and Fermi-LAT gamma-ray observations of blazars up to z = 1[END_REF], HESS[START_REF] Abramowski | Measurement of the extragalactic background light imprint on the spectra of the brightest blazars observed with H.E.S.S[END_REF] and Fermi/LAT[START_REF] Abdollahi | A gamma-ray determination of the Universe's star formation history[END_REF] correspondingly. Models: dashed blue line is for baseline EBL model of[START_REF] Gilmore | Semi-analytic modeling of the EBL and consequences for extragalactic gamma-ray spectra[END_REF]. Black striped band shows the allowed range of EBL models obtained with the global fit to EBL by[START_REF] Korochkin | Constraining the star formation rate with the extragalactic background light[END_REF]. Red shaded region shows the favoured range of additional component due to the ALP decay with parameters in the red shaded region of Fig.3.1.

Figure 3 . 4 -

 34 Figure 3.4 -Left panel: Same as in Fig. 3.1, with the sensitivity reach of CTA shown by the black dashed line (95% confidence level). Pink shaded ellipse shows the preferred ALP parameters (68% confidence level) obtained with combined VERITAS and simulated MAGIC spectrum of 1218 + 304 (50 hr long exposure). Red shaded ellipse in the center of pink ellipse shows the precision of the ALP mass and two-photon coupling which can be reached with a 50 hr long exposure of CTA. Right panel: Same as left panel. Brown dashed and solid lines show the precision of the ALP mass and two-photon coupling which can be reached with 1-and 5-years of observations with LHAASO (95% confidence level).

Figure 3 . 5 -

 35 Figure 3.5 -The energy E 0 corresponding to the pair-production optical depth τ = 1 as a function of the redshift z for different absorption models. Full line corresponds to the EBL model from the Chapter 2 [12], dashed line is for Franceschini et al. (2017) model [108], dash-dotted line is for Gilmore et al. (2012) fixed model [2] and dotted line is the Franceschini et al. (2017) model rescaled to satisfying direct EBL measurements.

Figure 3 . 6 -

 36 Figure 3.6 -Distribution of 31 objects from the sample in redshift (full red histogram). Hatched and shaded histograms give separate distributions for BLLs and FSRQs, respectively.

Figure 3 . 8 -

 38 Figure 3.8 -Breaks in deabsorbed blazar spectra at E b = E 0 (z).

Figure 3 . 9 -

 39 Figure 3.9 -Comparison of the fit quality (good fits correspond to low (1 -p)) for the assumptions of spectral breaks at the fixed energy E b (dashed green line, as a function of E b ) and at the distant-dependent energy E 0 (z) (full red line).

Figure 4 . 1 -

 41 Figure 4.1 -Cross sections of the inverse Compton scattering, Breit-Wheeler pair production, triplet pair production double pair production and photon-photon scattering.

Figure 4 . 2 -

 42 Figure 4.2 -Mean free path of gamma rays (left panel) and electrons/positrons (right panel) with respect to the whole extragalactic background radiation. Dashed lines indicate contributions of individual parts of the extragalactic background. Dot-dashed and dotted lines show the interaction length corresponding to the adiabatic and synchrotron energy losses respectively. Source[START_REF] Settimo | Propagation of extragalactic photons at ultra-high energy with the EleCa code[END_REF]

4 )

 4 where β = (s -m 2 e )/(s + m 2 e ) and s = m 2 e + 2E e E b (1 -µ). In Thomson regime, when s -m 2 e m 2

Figure 4 . 4 -

 44 Figure 4.4 -Differences in model predictions for the spectra of secondary γ-rays produced by interaction of monoenergetic primary γ-ray beam with energy E γ 0 = 10 TeV.Top panels shows the primary and secondary γ-ray spectra for sources at different redshifts. Bottom panels show the differences between CRbeam, CRPropa and ELMAG models.

Figure 4 . 5 -

 45 Figure 4.5 -Top: energy distribution of γ-rays produced via ICS on CMB by electrons with different energies (indicated as tags for each set of spectra). Middle: spectra of γ-rays from ICS on EBL for a range of electron energies. Bottom: Energy spectra of pairs produced in interactions of γ-rays with CMB photons. Energies of the primary γ-rays are shown as tags for each group of spectra. In all cases the codes used for the model calculations are specified in the legend.

Figure 4 . 6 -

 46 Figure 4.6 -Interaction rates for the pair production and inverse Compton scattering on the CMB at the redshifts z = 0 and z = 1.

Figure 4 . 7 -

 47 Figure 4.7 -Deflection of electrons of different energies in constant magnetic field with fixed physical strength in CRbeam and CRPropa.

Figure 4 . 8 -

 48 Figure 4.8 -Deflections of electrons at z = 0 in turbulent magnetic fields with different maximum scales L B . Panels show the root mean square deflection angle dependence on the physical path length.

Figure 4 . 9 -Figure 4 . 10 -

 49410 Figure 4.9 -Spectral (top), timing (middle) and imaging (bottom) properties of primary and secondary γ-ray signals for a IGMF with strength 10 -15 G and maximum scale L B = 10 Mpc. Left panels show calculations for the source at redshift z = 0.03, right panels are for the for the source at redshift z = 0.954. Calculations with CRbeam code (red solid) compared with CRPropa (blue thin solid) and corrected CRPropa (blue dashed) codes.

Figure 5 . 1 -

 51 Figure 5.1 -Upper limits on the IGMF. Source [207].

Figure 5 . 2 -

 52 Figure 5.2 -Schematic representation of the electromagnetic cascade in the intergalactic medium with nonzero IGMF. Source [228].

  regime of inverse Compton scattering 1 relevant for the scattering of CMB photons by electrons with energies 100 TeV. Such electrons produce inverse Compton emission at the energy

Figure 5 . 4 -

 54 Figure 5.4 -Comparison of high-energy ends of the spectra of Mrk 421 and Mrk 501 [243]. Black curves are observed, and blue curves are the intrinsic spectra of the sources. Thin and thick green lines show CTA North sensitivity to point sources and to the extended sources with angular size θ = 0.25 • . Green dotted line shows an analytical estimate of the secondary γ-ray flux from Mrk 501 in 1-10 energy range assuming no influence of IGMF.

  Probe of the strongest fields B 10 -11 G requires (a) large primary point-source power in the 100 TeV energy range, (b) detectability of extended emission in multi-TeV energy range, and (c) presence of primordial IGMF in the several Mpc region around the source.

Figure 5 . 5 -

 55 Figure 5.5 -Dark matter density profile along the line of sight toward Mrk 501 based on constrained simulations of the LSS [246].

2 Figure 5 . 6 -

 256 Figure 5.6 -Angular distribution of primary and secondary photons in different energy ranges. Black histograms show the primary point-source signal; green, blue, and red histograms show the extended emission calculated for different magnetic field strengths: 10 -12 G, 3 × 10 -12 G, and 10 -11 G. Horizontal dashed line shows the level of residual cosmic-ray electron background measured by the HESS [247].

2 Figure 5 . 7 -

 257 Figure 5.7 -Intrinsic spectrum of Mrk 501 (black solid line) in comparison with the absorbed spectrum (black dashed line) and the spectrum of secondary emission (colored lines) for different parameters of the magnetic field. Solid colored line corresponds to the total flux of the secondary emission from the source, while dashed-dotted and dotted lines are the fluxes in the angular ranges [0 : θ c ] and [θ c : 0.4] degrees.

Figure 5 . 8 -

 58 Figure 5.8 -Left panel: Significance of detection of the extended emission signal in different energy ranges as a function of the assumed magnetic field strength. The assumed exposure of CTA is T = 50 hr. Right panel: CTA exposure time needed for 3σ (lower curve) and 5σ (upper curve) detection of the extended emission signal, as a function of IGMF strength.

Figure 5 . 9 -

 59 Figure 5.9 -Density profile(top), magnetic field (middle) and angle between the magnetic field and the primordial direction (bottom) for a line of sight in the direction of Mkn 501, Mkn 421 and 1ES 1959+650. The magnetic field is split into total (blue), along the line of sight (orange) and perpendicular to line of sight (green). The location of the source with uncertainty is shown with a grey band. Regions where the overdensity is greater than 3 are indicated by red bands.

Figure 5 . 10 -

 510 Figure 5.10 -Images of the extended emission signal in the energy range 200 GeV -2 TeV for the three brightest sources in our sample. The assumed initial cosmological magnetic field strength is B = 10 -13 G. The direction of the jet axis coincides with the direction from the source to the observer and the jet opening angle is 5 • .

Figure 5 . 11 -

 511 Figure 5.11 -Images of the extended emission signal around Mkn 501 for different directions of the jet axis with respect to the projected direction of the cosmological magnetic field. Magnetic field strength, energy range and jet opening angle are the same as in Fig.5.10.

Figure 5 . 13 -

 513 Figure 5.13 -Orientation of the magnetic field along the z direction as a function of the magnetic field strength at redshift 0. The blue line shows the mean value, while shaded dark and light blue regions represent areas with 68% and 95% of simulation volume. The dashed red line shows the seed magnetic field value B = 10 -14 G.

Figure 5 . 14 -

 514 Figure5.14 -Examples of magnetic field profiles along two 600 Mpc long LOS that start at AGN found in the TNG-100 simulation. The upper panel shows the LOS starting inside a strong magnetic field bubble at the beginning, corresponding to AGN-3. The lower panel shows an example of the LOS starting from AGN-2 within a small bubble. We indicate the first bubble with green, other bubbles with red, and voids with blue. See section 5.4.3 for details.

1 Figure 5 . 15 -

 1515 Figure 5.15 -Primary and secondary γ-ray spectra for the "only voids" magnetic field profile (only regions with B ≤ 10 -12 G are taken into account). Different subplots represent different parameters of the source. Solid black and gray lines are the intrinsic and observed primary source spectra. Dashed yellow lines show the secondary γ-ray flux for a zero line of sight magnetic field. Solid red lines show the secondary γ-ray flux with a nonzero magnetic field in the voids.

Figure 5 . 16 -

 516 Figure 5.16 -Same as Fig. 5.15, but for the "only bubbles" magnetic field profiles (only regions with B > 10 -12 G are taken into account). Colored solid lines indicate different LOSs.

Figure 5 . 17 -

 517 Figure 5.17 -Same as Figs. 5.15 and 5.16, but for the "only host bubble" magnetic field profile (only the region with B > 10 -12 G around the source galaxy is taken into account).

zFigure 5 . 18 -

 518 Figure 5.18 -Fraction of the secondary γ-ray flux in the energy range 1-10 GeV that is removed by the bubbles as a function of the size of the host AGN bubble. Red circles show the flux fractions that are removed by all bubbles. Black crosses and blue squares show the flux fractions that are removed by all bubbles, except for the host AGN bubble and a linear fraction of the LOS that is occupied by these bubbles. The dashed blue horizontal line marks the average linear fraction of the LOS that is occupied by all bubbles, except for the host AGN bubble.

Figure 5 . 19 -

 519 Figure 5.19 -Fraction of the secondary flux in the energy range 1-10 GeV that is removed by the host AGN bubble as a function of cutoff energy and size of the bubble.

Figure B. 1 -

 1 Figure B.1 -Energy distribution of γ-rays produced via ICS on CMB by electrons with different energies (indicated as tags for each set of spectra) at different redshifts.

  on EBL(Franceschini et al. 2008), z=0 CRPropa corrected, ICS on EBL (Franceschini et al. 2008)on EBL (Franceschini et al. 2008), z=1 CRPropa corrected, ICS on EBL (Franceschini et al. 2008)on EBL (Franceschini et al. 2008), z=1 CRPropa corrected, ICS on EBL (Franceschini et al. 2008)on EBL (Franceschini et al. 2008), z=1 Theory

Figure B. 2 -

 2 Figure B.2 -Energy distribution of γ-rays produced via ICS on EBL by electrons with different energies (indicated as tags for each set of spectra) at different redshifts.

  production on CMB, z=1 Theory

Figure B. 3 -

 3 Figure B.3 -Energy spectra of pairs produced in interactions of γ-rays with CMB photons at different redshifts. Energies of the primary γ-rays are shown as tags for each group of spectra.

Table 2 .

 2 1 -Mean values and 1σ C.L. for the MCMC variable parameters in comparison with the value, presented in literature. Dash in MCMC column means that the model does not depend on this parameter.

	Parameter name	Symbol MCMC	literature
	IMF Parameter (M ) IMF Parameter IMF Parameter Dust Particle Minimal Size (µm)	m 0 D a imf a min	--2.29 +0.15 -0.24 < 1.75 •10 -3 -0.079 +0.021 -0.016 [54] 0.69 +0.05 -0.01 [54] 2.3±0.3 [54]
	Dust Particle Maximal Size (µm)	a max	< 22.8 •10 -3 -
	Dust Distribution Slope	n dust	> 3.83	∼3.5
	Radius of the Cloud (pc) Lifetime of the Cloud (Myr) Dust Particle Density (m -3 )	R c η c n d	6.1 +1.4 -1.2 6.0 +8.5 -3.6 6.9±2.0	2.5-100 [61] ∼10 [61] ∼10 [61]
	Optical Depth at 5500 Å	τ λ0	0.59 +0.57 -0.21	∼1 [56]
	Optical Depth Slope	n	0.47±0.24	∼0.7 [56]
	Cloud Efficiency (M )	M 0	85 +152 -36	-
	Redshift of the beginning of the star formation z i	-	∼8 [2]
	IMF Truncation Mass (M )	m max	< 84	∼100
	SFR scaling factor	C sfr	1.25 +0.44 -0.24	∼1

Table 2 .

 2 2 -Sample of seven blazars, selected for analysis. Ra and Dec are equatorial coordinates of the blazar, z is the redshift, and flux represents the power of the source, measured with the IACT.

	Name	Ra	Dec	z	Flux (Crab) Instr.	Obs. period Ref
	1ES 1011+496 153.76 49.43 0.212 0.05	MAGIC	2011-2012	[70]
	1ES 1215+303 184.45 30.10 0.131 0.035	MAGIC	2010-2011	[71]
	1ES 1218+304 185.36 30.19 0.182 0.08	VERITAS 2008-2013	[72]
	1ES 1959+650 299.99 65.14 0.048 0.64	VERITAS 2007-2011	[73]
	PKS 1510-089 228.21 -9.10 0.361 0.03	MAGIC	2012-2017	[74]
	PKS 2005-489 302.36 -48.83 0.071 0.03	HESS	2004-2007	[75]
	PKS 2155-304 329.72 -30.22 0.116 0.15	HESS	2005-2007	[76]

using Fermi/LAT data. The Fermi/LAT data that were collected during the ten-yearlong Fermi/LAT mission from August 4, 2008 to September 19, 2018 were processed using Fermi Science Tools software version 1.0.2. We used Fermi Pass 8 Release 3 data with 'SOURCE' class events from both the front and back.

Table 2 .

 2 3 -Fitting parameters for each blazar. χ 2 /d.o.f. represents the quality of the fit with the baseline EBL model while χ 2 b /d.o.f. is the quality of the fit for the EBL with an additional narrow bump, corresponding to the minimal EBL excess measured by CIBER [30], ∆χ 2 = χ 2 b -χ 2 . 83±0.04 -0.44±0.08 119±48 9.120 0.651 +0.918 0.717 PKS 1510-089 2.56±0.04 -0.38±0.24 90±50 4.322 0.480 -0.010 0.479 PKS 2005-489 1.86±0.05 -0.85±0.15 150±48 10.671 0.667 -0.942 0.608 PKS 2155-304 1.93±0.03 -0.85±0.09 138±32 10.843 0.571 -1.542 0.489 at 2-sigma level. This result is consistent with similar constraints derived by HESS

	Name	Γ 1	∆Γ	E br	χ 2	χ 2 /	∆χ 2	χ 2 b /
				GeV		d.o.f.		d.o.f.
	1ES 1011+496 1.89±0.04 -0.99±0.20 107±24 6.103 0.469 -0.291 0.447
	1ES 1215+303 1.98±0.05 -0.67±0.14 44±13 3.725 0.373 +0.449 0.417
	1ES 1218+304 1.65±0.05 -0.38±0.07 57±19 12.622 0.631 -1.395 0.561
	1ES 1959+650 1.

  .1.

			1ES 1218+304
	10 12 10 11 dE [TeV cm 2 E 2dN s 1 ]	z = 0.182		10 11 10 12 dE [TeV cm 2 E 2dN s 1
	10 0	10 1	10 2 E [GeV]	10 3	10 0	10 1	10 2 GeV	10 3

]

Figure

3

.3 -Left panel: Spectrum of 1ES 1218 + 304 measured by Fermi/LAT and VERITAS and theoretical models which fit it with and without axion-like particles. Blue and red points correspond to Fermi/LAT 10 years and VERITAS 86 hours measurements correspondingly. Black dashed thin and black solid thick lines are the best-fit broken power law intrinsic spectra and observed spectra absorbed with the baseline EBL model. Dashed and solid green and magenta lines show the intrinsic source and absorbed spectra for the baseline EBL with two additional narrow bumps form axion decays. Grey data point of VERITAS at 250 GeV affect range of allowed parameters presented in Fig.

2

for high axion masses. Right panel: Simulated spectra of 1ES 1218 + 304 for 50 hr long exposures of MAGIC (lime triangles) and CTA North (green circles) telescopes, compared to the VERITAS (red squares) and Fermi/LAT (blue squares) data points. Black think and thick lines show the intrinsic and absorbed spectral models used in the simulations.

Table 3 .

 3 1 -List of blazars selected for the main sample as described in the text. Classes: B=BLL, F=FSRQ. Instruments: H=HESS, L=Fermi LAT, M=MAGIC, V=VERITAS. Note: * ( †) in the last column means the source does not satisfy the (S2) criterion for the Franceschini et al. (2017) and Gilmore et al. 2012 fixed (Korochkin and Rubtsov 2018) EBL model, respectively.

	Name	Redshift Class Instrument Ref. (spec.) Ref. (z) Note
	B2 2114+33	1.596	F	L		[115]	*
	GB6 J0043+3426	0.966	F	L		[116]
	PKS 1441+25	0.939	F	M	[117]	[116]
	4C +55.17	0.896	F	L		[118]	*
	PKS 0537-441	0.892	F	L		[119]	*
	PG 1246+586	0.847	B	L		[120]
	Ton 0599	0.725	F	M	[121]	[122]
	PKS 1424+240	0.605	B	V	[123]	[124]
	3C 279	0.536	F	M	[125]	[126]
	4C +21.35	0.432	F	M	[127]	[128]
	PKS 1510-089	0.360	F	M	[129]	[130]
	OT 081	0.322	B	H	[131]	[132]
	OJ 287	0.306	F	V	[133]	[134]
	1ES 0414+009	0.287	B	H	[135]	[136]
	1RXS J023832.6-311658	0.233	B	H	[137]	[115]
	1ES 1011+496	0.212	B	M	[138]	[139]
	1ES 1218+304	0.182	B	V	[140]	[141]
	H 2356-309	0.165	B	H	[142]	[143]
	1ES 0229+200	0.140	B	V	[144]	[145]
	1ES 0806+524	0.138	B	M	[146]	[147]
	1ES 1215+303	0.131	B	M	[148]	[124]
	H 1426+428	0.129	B	HEGRA	[149]	[150]
	PKS 2155-304	0.116	B	H	[151]	[152]
	1ES 1312-423	0.105	B	H	[153]	[154]	†
	W Com	0.102	B	V	[155]	[124]	*
	RGB J0152+017	0.080	B	H	[156]	[157]
	PKS 2005-489	0.071	B	H	[158]	[159]
	PKS 0548-322	0.069	B	H	[160]	[161]	*
	PGC 2402248	0.065	B	M	[162]	[163]	*
	1ES 1959+650	0.048	B	HEGRA	[164]	[165]
	Mrk 501	0.034	B	H	[166]	[167]
	Mrk 421	0.031	B	HEGRA	[168]	[169]

  Bottom panels show discrepancies between the code calculations and the analytical model. For CRPropa, the differences between the analytical model and code calculations are not shown for redshifts > 0.3, because they exceed 50%.

	F = E 2dN dE , a.u.	10 16 10 17 10 18	CRbeam CRPropa CRPropa corrected ELMAG Intrinsic Survived			CRbeam CRPropa CRPropa corrected ELMAG Intrinsic Survived			2 .0 CRPropa corrected CRbeam CRPropa ELMAG Intrinsic Survived	0 . 5 1 .0	3
	Relative error, ( sim th )/ th , % sim F CRbeam ) /F Relative difference CRbeam , %	10 15 10 14 20 0 20 40	± 1% z=0.03 z=0.03 CRbeam CRPropa CRPropa corrected E, eV ELMAG ± 3%	z=0.1 z=0.3	z=1 z=2	± 1% z=0.03 z=0.14 CRbeam CRPropa CRPropa corrected ELMAG ± 3%	z=0.1 z=0.3	z=1 z=2	± 1% z=0.03 z=0.954	ELMAG, EBL optical depth Franceschini et al. 2008 Theory z=0.1 z=0.3 z=1 z=2 CRbeam CRPropa CRPropa corrected ELMAG ± 3%
	(F	40	z=0.03			z=0.14			z=0.954
			10 11 10 7	10 12 E, eV 10 9 E, eV 10 11 10 13	10 14 10 13	10 11 10 7	10 12 E, eV 10 9 E, eV 10 11 10 13	10 14 10 13	10 11 10 7	10 12 E, eV 10 9 E, eV 10 11 10 13	10 14 10 13

Figure 4.3 -Comparison of the optical depth for primary γ-rays for sources at cosmological redshifts in CRbeam, CRPropa and ELMAG. In all cases [185] EBL model is assumed. Top panels show the optical depth as a function of γ-ray energy for a range of redshifts.

  use results of the calculation of IGMF out to a distance of 200 Mpc from the Milky Way obtained by running MHD simulations using the Ramses-MHD code

	Name	RA	Dec	z	F 1TeV , TeV cm -2 s -1
	Mkn 421	166.11 38.21	0.031	2 × 10 -11
	Mkn 501	253.47 39.76	0.033	1 × 10 -11
	QSO B2344+514 356.77 51.7	0.044	4 × 10 -12
	Mkn 180	174.11 70.16	0.046	8 × 10 -13
	1ES 1959+650	299.99 65.15	0.047	6 × 10 -12
	AP Librae	229.42 -24.37 0.04903	4 × 10 -13
	TXS 0210+515	33.57 51.75 0.04913	2 × 10 -13

Table 5 .

 5 1 -List of nearby TeV blazars considered in the study.

Table A .

 A 1 -EBL lower limits used in the analysis in Chapter 1. Table A.2 -EBL upper limits used in the analysis in Chapter 1. Table A.3 -Direct measurements of the EBL used in the analysis in Chapter 1.

	E dN dE	10 10	10 4	λ (µm) Lower Limits (nW/m 2 /sr) 0.153 1.03±0.15 0.231 2.25±0.32 λ (µm) Upper Limits (nW/m 2 /sr) 0.4 11.6±1.7 0.52 23.4±4.7 10 16 eV 10 17 eV CRbeam, ICS on CMB, z=0 CRPropa corrected, ICS on CMB, z=0 Theory 1.25 21±15 2.2 20±6 3.5 13.3±2.8 1.25 54±16.8 (× 1 6 ) (× 8 ) (× 4 ) (× 2 ) 10 18 eV 10 19 eV 10 20 eV 0.36 2.87 +0.58 -0.42 0.45 4.57 +0.73 -0.47 0.67 6.74 +1.25 -0.94 0.81 2.2 27.8±6.7 1.4 28.7 +5.1 -3.3 3.6 9.0 +1.7 -0.9 2.2 22.4±6.0 10 3 10 2 10 1 10 0 E /E e 10 4 8.04 +1.62 -0.92 1.1 3.5 11.0±3.3 9.71 +3.0 -1.9 1.6 9.02 +2.62 -1.68 2.2 7.92 +2.04 -1.21 1.25 11.7±2.6 65 12.5±9.3 90 22.3±5.0 140 20.1±3.6 160 13.7±4.0 10 9 eV 10 10 eV 10 11 eV 10 12 eV 10 13 eV 10 14 eV 10 15 eV CRbeam, ICS on CMB, z=1 CRPropa corrected, ICS on CMB, z=1 Analytical estimate Theory	Experiment Galex [65] Experiment ESO VLT/FORSE [290] CRPropa, ICS on CMB, z=0 Theory ESO VLT/FORSE [290] COBE/DIRBE [291] COBE/DIRBE [291] COBE/DIRBE [291] COBE/DIRBE [292] Galex [65] HDF [33] HDF [33] HDF [33] COBE/DIRBE [292] CIBER [30] COBE/DIRBE [293] 10 3 10 2 E /E e COBE/DIRBE [294] HDF [33] COBE/DIRBE [294] HDF [33] HDF [33] HDF [33] Subaru [32] Akari [295] Akari [295] Akari [295] Akari [295] 10 9 eV 10 10 eV 10 11 eV CRPropa, ICS on CMB, z=1 Analytical estimate Theory	(× 1 6 ) (× 8 ) (× 4 ) (× 2 ) 10 1 10 12 eV	10 13 eV	10 16 eV 10 17 eV 10 18 eV 10 19 eV 10 20 eV 10 0 eV 10 14 10 15 eV
	E dN dE	10			1.6 2.12 140 160			11.5±1.5 10.0±0.8 12.6±6.0 13.7±6.1	Subaru [32] Subaru [32] COBE/FIRAS [296] COBE/FIRAS [296]
		10			3.6 4.5 140 240			5.4±1.4 3.5±0.9 32±13 17±4	Spitzer/IRAC [35] Spitzer/IRAC [35] COBE/DIRBE [297] COBE/DIRBE [297]
		10		10 7	5.8 8.0 24 100 140 240	10 5	E /E e	3.6±0.9 2.6±0.7 2.29±0.09 12.5±5.0 22±7 13.0±2.5 10 3 10 1	Spitzer/IRAC [35] Spitzer/IRAC [35] Spitzer/MIPS [285] COBE/DIRBE [298] COBE/DIRBE [298] COBE/DIRBE [298] 10 7 10 5 E /E e	10 3	10 1
					70			5.4±0.4	Spitzer/MIPS [285]
					70			7.4±1.9	Spitzer/MIPS [286]
					24			2.0±0.2	Spitzer/MIPS [287]
					70			7.1±1.0	Spitzer/MIPS [288]
					160			13.4±1.7	Spitzer/MIPS [288]
					100			6.33±1.67	Herschel/PACS [289]
					160			6.58±1.62	Herschel/PACS [289]
					250			8.5±0.6	BLAST [34]
					350			4.8±0.3	BLAST [34]
				500 250 λ (µm) Direct Measurements (nW/m 2 /sr) 2.2±0.2 7.4±1.42 Herschel/SPIRE [36] BLAST [34] Experiment 350 4.5±0.9 0.44 7.9±4.0 Pioneer 10/11 [299] Herschel/SPIRE [36] 500 1.54±0.34 0.64 7.7±5.8 Pioneer 10/11 [299] Herschel/SPIRE [36] 240 10.9±4.3 COBE/FIRAS [296]
				250			10.3±4.0	COBE/FIRAS [296]
				350			5.6±2.1	COBE/FIRAS [296]
				500			2.4±0.9	COBE/FIRAS [296]
				850			0.5±0.21	COBE/FIRAS [296]

Available at http://tevcat.uchicago.edu

available through the CTA website https://www.cta-observatory.org/science/ ctao-performance

This is true for energies below ∼ 100 TeV which we discuss here.

https://www.cta-observatory.org/science/cta-performance/

The "analysis caveats" page of Fermi Science Support Center, https://fermi.gsfc.nasa.gov/ ssc/data/analysis/caveats.html

Acknowledgements

First of all, I would like to thank my supervisors Dmitri Semikoz and Grigory Rubtsov for their ideas and guidance throughout the work.

Many thanks to Andrii Neronov and Sergey Troitsky for interesting and fruitful discussions of many scientific problems.

I would like to express my gratitude to my co-authors Kyrylo Bondarenko, Alexey Boyarsky, Oleg Kalashev, Guilhem Lavaux, Marius Ramsøy and Anastasia Sokolenko, without whom many of the papers that became part of the dissertation would not have been published.

I would further like to thank Alessandra Tonazzo for help in going through bureaucratic procedures. I would also like to thank Yvonne Becherini, Ruth Durrer, Günter Sigl, Peter Tinyakov, Igor Tkachev and Andreas Zech for agreeing to be in the Jury of my dissertation.

Special thanks go to my wife Gulnur, my parents and my friends, whose support and endless patience gave me strength.

My stay in the APC laboratory in Paris was provided by the scholarship 'Vernadsky' of the French embassy in Russia.

Appendix A EBL measurements used in Chapter 1

Tables contain EBL measurements used in Chapter 1.

118

Appendix B

Comparison of distribution of secondaries [START_REF] Terzic | First detection of VHE gamma-ray signal from the FSRQ TON 0599[END_REF]