<figure><img src="image_1.png" /></figure>

HAL
open science

# x-cigale: fitting AGN/galaxy SEDs from X-ray to infrared

G Yang( ), M Boquien, V Buat, D Burgarella, L Ciesla, F Duras, M Stalevski, W Brandt, C Papovich

► To cite this version:

G Yang( ), M Boquien, V Buat, D Burgarella, L Ciesla, et al.. x-cigale: fitting AGN/galaxy SEDs from X-ray to infrared. Monthly Notices of the Royal Astronomical Society, 2020, 491 (1), pp.740-757. 10.1093/mnras/stz3001. hal-03517303

HAL Id: hal-03517303

https://hal.science/hal-03517303

Submitted on 8 Jan 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# X-CIGALE: fitting AGN/galaxy SEDs from X-ray to infrared

G. Yang(杨光)<fnref n="1" />,<fnref n="2" />,<fnref n="3" />★ M. Boquien,<fnref n="4" /> V. Buat,<fnref n="3" /> D. Burgarella,<fnref n="3" /> L. Ciesla,<fnref n="3" /> F. Duras,<fnref n="3" /> M. Stalevski <fnref n="5" />,<fnref n="6" /> W. N. Brandt<fnref n="7" />,<fnref n="8" />,<fnref n="9" /> and C. Papovich<fnref n="1" />,<fnref n="2" />

Accepted 2019 October 21. Received 2019 October 18; in original form 2019 August 17

## ABSTRACT

Code Investigating GALaxy Emission (CIGALE) is a powerful multiwavelength spectral energy distribution (SED) fitting code for extragalactic studies. However, the current version of CIGALE is not able to fit X-ray data, which often provide unique insights into active galactic nucleus (AGN) intrinsic power. We develop a new X-ray module for CIGALE, allowing it to fit SEDs from the X-ray to infrared (IR). We also improve the AGN fitting of CIGALE from UV-to-IR wavelengths. We implement a modern clumpy two-phase torus model, SKIRTOR. To account for moderately extincted type 1 AGNs, we implement polar-dust extinction. We publicly release the source code (named ‘X-CIGALE’). We test X-CIGALE with X-ray detected AGNs in SDSS, COSMOS, and AKARI-NEP. The fitting quality (as indicated by reduced $\chi^2$) is good in general, indicating that X-CIGALE is capable of modelling the observed SED from X-ray to IR. We discuss constrainability and degeneracy of model parameters in the fitting of AKARI-NEP, for which excellent mid-IR photometric coverage is available. We also test fitting a sample of AKARI-NEP galaxies for which only X-ray upper limits are available from Chandra observations, and find that the upper limit can effectively constrain the AGN SED contribution for some systems. Finally, using X-CIGALE, we assess the ability of Athena to constrain the AGN activity in future extragalactic studies.

**Key words:** methods: data analysis – methods: observational – galaxies: nuclei – quasars: general – X-rays: general.

## 1 INTRODUCTION

Supermassive black holes (BHs) commonly exist in the centres of massive galaxies (e.g. Kormendy & Richstone 1995; Kormendy & Ho 2013). BHs grow their mass ($M_{BH}$) by accreting local material. During this process, a significant amount of the gravitational energy of the accreted material is converted to radiation, and the system shines as an active galactic nucleus (AGN). The typical spectral energy distribution (SED) of AGNs covers a broad wavelength range, from X-ray to infrared (IR).

AGN emission at different wavelengths is generated by different physical processes (e.g. Netzer 2013). The accretion disc mostly

produces photons at ultraviolet (UV) and optical wavelengths. Some of these photons are scattered to X-ray energies by the hot corona above the disc (i.e. inverse Compton scattering). Some of the UV/optical photons might also be absorbed by dust. The dust is thus heated and reemits the energy as infrared radiation. Considering the tight link between AGN multiwavelength SEDs and these physical processes, it is feasible to infer source properties from modelling the observed photometric data. On the other hand, the observed SED is often complicated, involving factors such as host-galaxy contributions and dust extinction. Misinterpretation of the SED could lead to unrealistic physical properties. Therefore, it is critical to decipher the observed data appropriately with a powerful and reliable SED fitting code.

The Code Investigating GALaxy Emission (CIGALE) is a state-of-the-art PYTHON code for SED fitting of extragalactic sources

<footnote n="1">¹Department of Physics and Astronomy, Texas A&M University, College Station, TX 77843-4242, USA</footnote>

<footnote n="2">²George P. and Cynthia Woods Mitchell Institute for Fundamental Physics and Astronomy, Texas A&M University, College Station, TX 77843-4242, USA</footnote>

<footnote n="3">³Aix Marseille Univ, CNRS, CNES, LAM, Marseille, France</footnote>

<footnote n="4">⁴Centro de Astronomía (CITEVA), Universidad de Antofagasta, Avenida Angamos 601, Antofagasta, Chile</footnote>

<footnote n="5">⁵Astronomical Observatory, Volgina 7, 11060 Belgrade, Serbia</footnote>

<footnote n="6">⁶Sterrenkundig Observatorium, Universiteit Gent, Krijgslaan 281-S9, Gent B-9000, Belgium</footnote>

<footnote n="7">⁷Department of Astronomy and Astrophysics, 525 Davey Lab, The Pennsylvania State University, University Park, PA 16802, USA</footnote>

<footnote n="8">⁸Institute for Gravitation and the Cosmos, The Pennsylvania State University, University Park, PA 16802, USA</footnote>

<footnote n="9">⁹Department of Physics, 104 Davey Laboratory, The Pennsylvania State University, University Park, PA 16802, USA</footnote>

<footnote>* E-mail: gyang206265@gmail.com</footnote>
(Boquien et al. 2019). It employs physical AGN and galaxy models, and allows flexible combination between them. The current version of CIGALE can simultaneously fit the observed SED from UV to far-IR (FIR) and extract source physical properties such as AGN luminosity and host stellar mass ($M_*$). However, the current CIGALE is not able to model X-ray fluxes, which often provide a unique view of AGNs.

X-ray observations have many advantages in AGN studies (see Brandt & Alexander 2015 for a review). Strong X-ray emission is nearly a universal property of the AGN phenomenon. X-rays are generated from the immediate vicinity of the BH, directly revealing the intrinsic AGN power. Therefore, X-ray fluxes are widely used as a tracer of BH accretion rate (e.g. Yang et al. 2018a, 2019). Thanks to their great penetrating power, X-rays are only mildly affected by obscuration in general. Also, AGNs are much more efficient in generating X-rays than their host galaxies. Therefore, the observed X-ray fluxes are often dominated by AGNs and have negligible galaxy contribution. Considering these advantageous properties, X-ray observations are widely used to select AGNs, especially in the distant universe (e.g. Luo et al. 2017; Chen et al. 2018). These selections are often more complete and reliable than the selections at other wavelengths such as optical and IR.

Besides the lack of X-ray fitting capability, CIGALE's current AGN model (Fritz, Franceschini & Hatziminaoglou 2006), which covers the UV to IR, also has some other disadvantages. The model assumes that the central engine is surrounded by a dusty torus (i.e. the AGN unified model; Antonucci 1993; Urry & Padovani 1995; Netzer 2015; Zou et al. 2019). The torus absorbs a fraction of the UV and optical emission from the central engine and reemits the energy as IR photons. When viewing from the equatorial direction, the central engine is obscured and only reemitted IR radiation can be observed (type 2 AGN). When viewing from the polar direction, the central engine is directly visible (type 1 AGN).

One disadvantage of the AGN model is that it assumes the dusty torus is a smooth structure. However, such smooth models for the torus are disfavoured on physical grounds (e.g. Tanimoto et al. 2019). To reach a scale height consistent with observations, the dust grains in a smooth torus would have random velocities ~100 km s⁻¹, corresponding to a temperature of ~10⁶ K. This high temperature far exceeds the dust-sublimation temperature (~10³ K). Another disadvantage of the AGN model is that the disc emission is assumed to be absolutely unextincted for the case of type 1. However, recent observations indicate that a non-negligible amount of extinction exists for some type 1 AGNs (e.g. Bongiorno et al. 2012; Elvis et al. 2012; Lusso et al. 2012), which can be attributed to the dust existing along polar directions (e.g. Stalevski, Asmus & Tristram 2017; Stalevski, Tristram & Asmus 2019; Lyu & Rieke 2018). The current CIGALE cannot model the SEDs of these type 1 AGNs.

In this paper, we further develop CIGALE and enable it to fit X-ray data. The new development allows CIGALE to model AGN SED from X-ray to IR simultaneously and extract source properties such as AGN intrinsic luminosity and host-galaxy stellar mass ($M_*$). Besides developing the X-ray part, we also improve CIGALE's capability in fitting the UV-to-IR SED of AGNs. We implement the latest version of SKIRTOR, a clumpy two-phase torus model derived from a modern radiative-transfer method (Stalevski et al. 2012, 2016). In addition, we introduce polar-dust extinction to account for the possible extinction in type 1 AGNs. We name the new version of CIGALE as 'X-CIGALE'.

The structure of this paper is as follows. In Section 2, we outline the scheme of our new code development. In Section 3, we test

X-CIGALE on AGNs with X-ray detections from different surveys. We test fitting galaxies with only X-ray upper limits in Section 4. We summarize our results and discuss future prospects in Section 5.

Throughout this paper, we assume a flat ΛCDM cosmology with $H_0 = 69.3 \text{ km s}^{-1} \text{ Mpc}^{-1}$ and $\Omega_M = 0.286$ (WMAP 9-yr results; Hinshaw et al. 2013). Quoted uncertainties are at the $1\sigma$ (68 per cent) confidence level, unless otherwise stated. Quoted optical/infrared magnitudes are AB magnitudes.

## 2 THE CODE

We briefly summarize the mechanisms and features of CIGALE in Section 2.1. In Sections 2.2, 2.3, and 2.4, we detail our new development of X-CIGALE, i.e. the X-ray fitting, SKIRTOR, and the polar-dust extinction. The new inputs/outputs introduced in X-CIGALE are listed in Appendix A.

### 2.1 A brief introduction of CIGALE

CIGALE is an efficient SED-fitting code which has been developed for more than a decade (Burgarella, Buat & Iglesias-Páramo 2005; Noll et al. 2009; Serra et al. 2011; Boquien et al. 2019). CIGALE is written in PYTHON. X-CIGALE is built upon CIGALE, and the fitting algorithm of X-CIGALE is the same as that of CIGALE. Here, we only briefly introduce the algorithm, and interested readers should refer to Boquien et al. (2019) for a detailed description.

CIGALE allows the user to input a set of model parameters. The code then realizes the model SED for each possible combination of the model parameters, and convolves the model SED with the filters to derive model fluxes. By comparing the model fluxes with the observed fluxes, the code computes likelihood as $L = \exp(-\chi^2/2)$ for each model. CIGALE supports two types of analyses, i.e. maximum likelihood (minimum $\chi^2$) and Bayesian like. In the maximum-likelihood analyses, CIGALE picks out the model with the largest $L$ value, and calculates physical properties such as $M_*$ and star formation rate (SFR) from this single model. In the Bayesian-like analyses, for each physical property, CIGALE calculates the marginalized probability distribution function (PDF) based on the $L$ values of all models. Finally, from this PDF, CIGALE derives the probability-weighted mean and standard deviation, and outputs them as the estimated value and uncertainty.

Among the above processes, one key step is the realization of model SEDs from input parameters. This procedure relies on a set of modules, and each module is responsible for a function that shapes the SED. For example, the 'nebular-emission' module adds the nebular-emission components to the SED, and the 'dust attenuation' module extincts the SED. Our new development of X-CIGALE follows this module-based structure. We enable CIGALE to fit X-ray data by developing a new X-ray module (Section 2.2); we implement SKIRTOR templates and polar-dust extinction in a new SKIRTOR module (Sections 2.3 and 2.4).

### 2.2 The new X-ray module

In this section, we develop a new X-ray module to enable X-CIGALE to fit X-ray data. In Section 2.2.1, we detail the basic settings of this new module. In Section 2.2.2, we present the adopted X-ray SED for AGN and galaxy components. In Section 2.2.3, we present the relation that we used to link AGN X-ray with other wavelengths. We note that our new developments are for the majority AGN population in optical/X-ray surveys, and thus X-CIGALE may not be applicable to some minor populations such as radio-loud and
broad absorption line (BAL) objects (e.g. Brandt, Laor & Wills 2000; Miller et al. 2011). We leave the treatment of these particular AGNs to future works.

### 2.2.1 Basic settings

As presented in Section 1, the X-ray band has many advantages in studying AGNs. Therefore, we implement an X-ray module for X-CIGALE. The main goal of this module is to connect X-ray with other wavelengths, rather than to obtain detailed X-ray spectral properties (e.g. photon index and hydrogen column density) by performing detailed X-ray spectral analyses. This is because the latter has already been well realized by many specialized X-ray codes such as XSPEC (Arnaud 1996) and SHERPA (Freeman, Doe & Siemiginowska 2001), and there is no need for X-CIGALE to perform similar analyses. Also, it is technically difficult to fit the X-ray spectra within the framework of X-CIGALE. X-CIGALE assumes that a sample of sources are observed with a single 'filter transmission', as is the case in UV-to-IR data. However, at X-ray wavelengths, the transmission curve varies from source-to-source, as it might depend on many factors such as position on the detectors and observation date. For example, the soft-band transmission of Chandra has been continuously declining since its launch (e.g. O'Dell et al. 2017). In fact, each source is associated with a unique transmission curve and the curve is taken into account when fitting the X-ray spectra with, e.g. XSPEC and SHERPA.

Therefore, the X-ray module of X-CIGALE is designed to work on the high-level X-ray data products, i.e. intrinsic X-ray fluxes in a given band. We require the X-ray fluxes to be corrected for telescope transmission. Fortunately, this correction is embedded in routine X-ray data processing and has already been applied in X-ray photometric catalogues (e.g. Yang et al. 2016; Luo et al. 2017). Since the transmission has already been considered, X-CIGALE only needs to adopt a uniform-sensitivity (i.e. boxcar-shaped) 'filter'. We have already included a few typical boxcar X-ray filters, e.g. 0.5–2 and 2–7 keV for convenience, while the user can easily generate the filters for any X-ray band.

In addition, we require the input X-ray fluxes to be 'absorption-corrected'. The absorption might be from the source itself, our Galaxy, and/or the intergalactic medium (IGM; e.g. Starling et al. 2013). However, we do not differentiate these types of absorption, as it is often infeasible to separate them in X-ray data analyses. The absorption correction can be obtained from routine X-ray data processing, e.g. spectral analyses via XSPEC/SHERPA or band-ratio analyses (e.g. Xue et al. 2016; Luo et al. 2017). The user may also choose to use hard X-ray bands where absorption corrections are generally small (e.g. Yang et al. 2018a,b).

In X-ray catalogues (e.g. Xue et al. 2016; Luo et al. 2017), X-ray fluxes ($f_X$) are conventionally given in the cgs units of erg s⁻¹ cm⁻², but X-CIGALE requires the input fluxes ($f_{CIGALE}$) to be given in the units of mJy. Therefore, the user needs to convert the flux units with

$$f_{CIGALE} = \frac{f_{\chi} \times 4.136 \times 10^{8}}{E_{up} - E_{lo}}, \quad (1)$$

where $E_{lo}$ and $E_{up}$ refer to the lower and upper limits of the energy band in units of keV.

The X-ray module covers rest-frame 10⁻³ to 5 nm, corresponding to ≈0.25–1200 keV. Such an energy range is sufficient for practical purposes: current X-ray instruments cannot observe energies significantly below rest-frame 0.5(1 + z) keV in general; the AGN flux is typically non-detectable above ≈1000 keV due to the existence of the cut-off energy in AGN X-ray spectra (see Section 2.2.2).

<figure><img src="image_4.png" /><figcaption>Figure 1. An example X-ray SED model for a typical source with AGN $L_{2-10\,keV} = 10^{43}\,erg\,s^{-1}$, $M_{\star} = 10^{11}\,M_{\odot}$, SFR = $10\,M_{\odot}\,yr^{-1}$, $T = 1\,Gyr$, and $Z = 0.02$. Different colours indicate different components. For this source, the X-ray luminosity is dominantly contributed by the AGN.</figcaption></figure>

### 2.2.2 X-ray SED

To first-order approximation, the intrinsic AGN X-ray spectrum is typically a power law with a high-energy exponential cut-off, i.e.

$$f_{\nu} \propto E^{-\Gamma+1} \exp(-E/E_{\text{cut}}), \quad (2)$$

where Γ is the so-called photon index, widely adopted in X-ray astronomy, and $E_{cut}$ is the cut-off energy. We adopt this spectral shape in X-CIGALE. Detailed X-ray spectral fitting in the literature finds Γ ≈ 1.8 (e.g. Yang et al. 2016; Liu et al. 2017). We allow the user to set Γ in X-CIGALE. We set $E_{cut} = 300\,keV$, the typical value from the observations of Seyferts (e.g. Dadina 2008; Ricci et al. 2017). Note that since $E_{cut}$ is above the highest observable energy of most X-ray observatories (e.g. Chandra and XMM-Newton), the exact choice of $E_{cut}$ has practically negligible effects on the fitting with X-CIGALE for most cases. The adopted AGN X-ray SED is displayed in Fig. 1.

Besides AGNs, galaxies can also emit X-rays, although the emission from galaxies is often much weaker than that from AGNs for X-ray detected sources. There are three main origins of galaxy X-ray emission: low-mass X-ray binaries (LMXB), high-mass X-ray binaries (HMXB), and hot gas. The strengths of these components can be modelled as a function of galaxy properties such as $M_{\star}$ and SFR. We adopt the recipe from Mezcua et al. (2018), where a Chabrier (2003) initial mass function (IMF) is assumed. In this scheme, the LMXB and HMXB luminosities (in units of erg s⁻¹) are described as

$$\begin{aligned} \log (L_{2-10\,keV}^{LMXB}/M_{\star}) &= 40.3 - 1.5 \log t - 0.42(\log t)^2 \\ &\quad + 0.43(\log t)^3 + 0.14(\log t)^4, \\ \log (L_{2-10\,keV}^{HMXB}/\text{SFR}) &= 40.3 - 62Z + 569Z^2 \\ &\quad - 1834Z^3 + 1968Z^4, \end{aligned} \quad (3)$$

where $M_{\star}$ and SFR are in solar units; t denotes stellar age in units of Gyr; and Z denotes metallicity (mass fraction). The hot-gas luminosity (in units of erg s⁻¹) is described as

$$\log (L_{0.5-2\,keV}^{\text{hot gas}}/\text{SFR}) = 38.9. \quad (4)$$
Similarly as for AGN, we also employ the SED shape in equation (2) for all three components, with $E_{\text{cut}}$ fixed at 100 keV (LMXB and HMXB; e.g. Zhang 1997; Motta, Belloni & Homan 2009) and 1 keV (hot gas; e.g. Mathews & Brighenti 2003). We allow the user to set $\Gamma$ values for the LMXB and HMXB components. In our test fitting in Section 3, we set $\Gamma$ to 1.56 and 2.0 for LMXB and HMXB, respectively (e.g. Fabbiano 2006; Sazonov & Khabibullin 2017). Adjusting these $\Gamma$ does not affect the fitting results significantly, as the observed X-ray fluxes are often dominantly contributed by AGNs rather than galaxies. The X-ray continuum from hot gas can be modelled as free-free and free-bound emission from optically thin plasma ($\Gamma = 1$; e.g. Mewe, Lemen & van den Oord 1986). Therefore, we fix $\Gamma = 1$ for the hot-gas component in X-CIGALE. Fig. 1 shows the adopted X-ray SEDs of the three components. We add all three components for the total X-ray SED from galaxies.

### 2.2.3 The $\alpha_{\text{ox}}-L_{2500\AA}$ relation

As in Section 2.2.1, the main goal of X-CIGALE is to fit X-ray and other wavelengths simultaneously. Some known connections between X-ray and other wavelengths must be applied; otherwise, the fitting would be practically useless. We adopt the well-studied 'α$_{ox}$-L$_{2500\AA}$' relation (e.g. Steffen et al. 2006; Just et al. 2007; Lusso & Risaliti 2017), where $L_{2500\AA}$ is AGN intrinsic (de-reddened) luminosity per frequency at 2500 Å and $\alpha_{\text{ox}}$ is the SED slope between UV (2500 Å) and X-ray (2 keV), i.e.

$$ \alpha_{\text{ox}} = -0.3838 \log(L_{2500\AA}/L_{2\text{keV}}). \qquad (5) $$

The observed $\alpha_{\text{ox}}-L_{2500\AA}$ relation (Just et al. 2007) is written as

$$ \alpha_{\text{ox}} = -0.137 \log(L_{2500\AA}) + 2.638, \qquad (6) $$

where $L_{2500\AA}$ is in units of erg s$^{-1}$ Hz$^{-1}$. The $1\sigma$ intrinsic dispersion of this $\alpha_{\text{ox}}-L_{2500\AA}$ relation is $\Delta\alpha_{\text{ox}} \approx 0.1$ (see table 8 of Just et al. 2007). Here, $\Delta\alpha_{\text{ox}}$ is the $\alpha_{\text{ox}}$ deviation from that expected from the $\alpha_{\text{ox}}-L_{2500\AA}$ relation, i.e.

$$ \Delta\alpha_{\text{ox}} = \alpha_{\text{ox}} - \alpha_{\text{ox}}(L_{2500\AA}). \qquad (7) $$

Fig. 2 displays example SEDs of different $\Delta\alpha_{\text{ox}}$. Observations have found that the $\alpha_{\text{ox}}-L_{2500\AA}$ relation does not have significant redshift evolution, indicating that the relation originates from fundamental accretion physics (Steffen et al. 2006; Just et al. 2007; Lusso & Risaliti 2017). We allow the user to set the maximum $|\Delta\alpha_{\text{ox}}|$ allowed ($|\Delta\alpha_{\text{ox}}|_{\text{max}}$). Internally, X-CIGALE calculates all models with $\alpha_{\text{ox}}$ from -1.9 to -1.1 with a step of 0.1.<fnref n="1" /> X-CIGALE then calculates $|\Delta\alpha_{\text{ox}}|$ and discards the models with $|\Delta\alpha_{\text{ox}}| > |\Delta\alpha_{\text{ox}}|_{\text{max}}$. In our test fitting (Sections 3 and 4), we adopt $|\Delta\alpha_{\text{ox}}|_{\text{max}} = 0.2$, corresponding to the $\approx 2\sigma$ scatter of the $\alpha_{\text{ox}}-L_{2500\AA}$ relation (Just et al. 2007).

Note that the $\alpha_{\text{ox}}-L_{2500\AA}$ relation above is derived from observations, assuming that the unobscured AGN emission is isotropic at both UV/optical and X-ray wavelengths. However, the UV/optical emission is unlikely isotropic, because it is from the accretion disc and the effects of projected area and limb darkening affect the angular distribution of the radiative energy. After considering these effects, the disc luminosity can be approximated as $L(\theta) \propto \cos\theta(1+2\cos\theta)$, where $\theta$ is the angle from the AGN axis (e.g. Netzer 1987). This angular dependence of disc emission is

adopted in SKIRTOR, the UV-to-IR AGN module adopted in X-CIGALE (see Section 2.3). The X-ray emission should likely be less anisotropic than the UV/optical emission, because the X-rays originate from re-processed UV/optical photons via inverse Compton scattering. However, the exact relation between X-ray flux and viewing angle depends on model details, such as corona shape and opacity, which are poorly known (e.g. Liu et al. 2014; Xu 2015). For simplicity, we assume that the X-ray emission is isotropic.

Our assumption of anisotropic UV/optical emission and isotropic X-ray emission leads to a dependence of $\alpha_{\text{ox}}-L_{2500\AA}$ on viewing angle. We further assume that the observed $\alpha_{\text{ox}}-L_{2500\AA}$ relation for type 1 AGNs reflects the intrinsic $\alpha_{\text{ox}}-L_{2500\AA}$ relation for all AGNs at a 'typical' viewing angle of $\theta = 30^\circ$. This value approximates the probability-weighted viewing angle for type 1 AGNs, i.e.

$$ \theta \approx \frac{\int_{0}^{90^\circ - \Delta} \theta \sin \theta d\theta}{\int_{0}^{90^\circ - \Delta} \sin \theta d\theta} \approx 30^\circ, \qquad (8) $$

where $\Delta$ denotes the angle between the equatorial plane and edge of the torus, i.e. half opening angle. The typical value is $\Delta \approx 40^\circ$ from observations (e.g. Stalevski et al. 2016). Although $\Delta$ is a free parameter in X-CIGALE, we do not recommend the user choose other values than $40^\circ$, as this value is favoured by observations and is consistently adopted throughout the build-up of the X-CIGALE code. The weight $\sin\theta$ is proportional to the probability for the viewing angle being $\theta$.

We note that our SED fitting results (Sections 3 and 4) are not sensitive to the assumed typical $\theta$, and will not change significantly if adjusting $\theta$ within the range of $\approx 10-50^\circ$. In the X-CIGALE output (Appendix A), the $\alpha_{\text{ox}}$ and $L_{2500\AA}$ always refer to the value at $\theta = 30^\circ$, regardless of the actual viewing angle in the model. This $\alpha_{\text{ox}}$ and $L_{2500\AA}$ design is to reflect AGN essential properties, independent of the viewing angle. By changing the integral ranges in equation (8) to ($90^\circ - \Delta$, $90^\circ$), we can derive the probability-weighted $\theta$ for type 2 AGNs, i.e. $\theta \approx 70^\circ$. These typical $\theta$ values (type 1: $\approx 30^\circ$, type 2: $\approx 70^\circ$) are used in our SED fitting (Sections 3 and 4).

## 2.3 SKIRTOR

The previous CIGALE AGN model responsible for the UV-to-IR SED is from Fritz et al. (2006). This model assumes that the dusty torus is a smooth structure. However, more recent theoretical and observational works find that the torus is mainly made of dusty clumps (e.g. Nikutta, Elitzur & Lacy 2009; Ichikawa et al. 2012; Stalevski et al. 2012; Tanimoto et al. 2019). SKIRTOR is a clumpy two-phase torus model (Stalevski et al. 2012, 2016), based on the 3D radiative-transfer code, SKIRT (Baes et al. 2011; Camps & Baes 2015).<fnref n="2" /> In SKIRTOR, most (mass fraction = 97 per cent) of the dust is in the form of high-density clumps, while the rest is smoothly distributed. In addition, SKIRTOR considers the anisotropy of the power source, AGN disc emission (see Section 2.2.3), while Fritz's model simply assumes isotropic disc emission. Therefore, we implement SKIRTOR within X-CIGALE. We recommend using SKIRTOR as the UV-to-IR SED model of AGNs, although X-CIGALE allows the user to choose between SKIRTOR and Fritz's model.

SKIRTOR adopts a disc SED that has a higher fraction of far-UV luminosity ($\lesssim 100$ nm) compared to observations (see section 3.2.1 of Duras et al. 2017). Following Duras et al. (2017), we update

<footnote n="1">¹This $\alpha_{\text{ox}}$ range corresponds to 2–10 keV X-ray bolometric corrections ranging from $\approx 10$ to $\approx 500$.</footnote>

<footnote n="2">²http://www.skirt.ugent.be/root/landing.html</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 2.** The model SEDs for an unobscured AGN with $log L_{2500Å} = 30$ (cgs). Different colours indicate different $\alpha_{ox}$. The ‘breaks’ at 5 nm are caused by the wavelength limit of the X-ray module, which ends at 5 nm. As explained in Section 2.2.3, such breaks are generally not problematic for practical purposes, as most X-ray instruments only cover wavelengths $\lesssim 10^{0.4}$ nm ($\gtrsim 0.5$ keV, as marked).</figcaption></figure>

SKIRTOR with a new disc SED (Feltre et al. 2012) that is supported by observations, i.e.

$$ \lambda L_{\lambda} \propto \begin{cases} \lambda^2 & 8 \le \lambda < 50 \, [\text{nm}] \\ \lambda^{0.8} & 50 \le \lambda < 125 \, [\text{nm}] \\ \lambda^{-0.5} & 125 \le \lambda < 10^4 \, [\text{nm}] \\ \lambda^{-3} & \lambda > 10^4 \, [\text{nm}]. \end{cases} \qquad (9) $$

We modify the disc SED with the following method. We denote the old and new intrinsic disc SEDs as $L_{\lambda,\text{normed}}^{\text{old}}$ and $L_{\lambda,\text{normed}}^{\text{new}}$, respectively, where the subscript ‘normed’ indicates the total power of these SEDs has been normalized to unity. Then the new observed disc SED component (which might be obscured) can be converted from the old one by multiplying by the factor, $L_{\lambda,\text{normed}}^{\text{new}}/L_{\lambda,\text{normed}}^{\text{old}}$. The new scatter component can be obtained in the same way; the dust reemitted component remains unchanged. The method above keeps energy balance. This method is also described on the SKIRTOR official webpage.<fnref n="3" />

## 2.4 Polar dust

### 2.4.1 The extinction of type 1 AGN

In SKIRTOR (also in Fritz’s model), the extinction of UV and optical radiation for type 1 AGN is assumed to be negligible. This assumption holds for most optically selected blue quasars. For example, Richards et al. (2003) found only $\approx 6$ per cent of their SDSS quasars are extincted. However, the assumption might not be true for, e.g. X-ray selected AGNs. For example, in the COSMOS AGN catalogues selected by XMM–Newton (Bongiorno et al. 2012), the fraction of extincted sources ($E(B-V) \ge 0.1$) among broad-line AGNs is $\approx 40$ per cent.

To check the extinction of type 1 AGNs, we compare the median UV–optical SEDs of spectroscopically classified type 1 AGNs in SDSS and COSMOS (see Section 3 for details). These median SEDs are derived from the photometric data in Section 3. For each source

in a sample, we interpolate the observed photometry to obtain $F_v$ as a function of observed-frame wavelength. We then shift this interpolated SED to rest-frame wavelength and normalize $F_v$ at 250 nm. Finally, at each wavelength, we obtain the median $F_v$ of all the sources in the sample. Fig. 3 (left) shows the results. The SDSS median SED is similar to the typical unobscured quasar SED of $F_v \propto \lambda^{0.5}$ (see equation 9). In contrast, the COSMOS median SED is significantly redder than $F_v \propto \lambda^{0.5}$ (e.g. Elvis et al. 2012). We note that this difference in SED shape is observationally driven by selection effects. The SDSS sample consists of optically selected, and is thus biased towards blue and optically bright objects. The COSMOS sample consists of X-ray selected objects and does not suffer from significant bias in the UV/optical (e.g. Brandt & Alexander 2015). Although driven by selection effects, Fig. 3 (left) at least indicates that reddened AGN SEDs indeed exist, and we discuss the physical cause of the SED reddening below.

The red SED shape might be physically caused by the aforementioned dust extinction. However, another potential physical cause is host-galaxy contribution to the SED. Since the UV/optical SEDs of galaxies are generally redder than those of unobscured AGNs (e.g. fig. 3 of Salvato et al. 2009), AGN-galaxy mixed SEDs tend to be redder than pure AGN SEDs. To investigate the cause of SED reddening, we can compare the magnitudes that sample rest-frame UV wavelengths for COSMOS and SDSS, as galaxy contributions to the photometry should be small at UV wavelengths.<fnref n="4" /> In Fig. 3 (right), we show the r-band (rest frame $\lesssim 2000$ Å) magnitude distributions at $z=2-2.5$ and $f_{v,2keV} = 3-10$ ($10^{-7}$ mJy, the results are similar for other redshift/X-ray flux bins). The control of redshift and X-ray flux is to force the compared samples to have similar X-ray luminosities and thereby bolometric luminosities (assuming the X-ray bolometric correction factors are similar for the sources). The COSMOS AGNs are systematically fainter than the SDSS AGNs in UV/optical. Therefore, Fig. 3 (right) indicates that, at a given AGN bolometric luminosity, the rest-frame UV AGN luminosities of COSMOS are typically lower than those of SDSS, supporting the existence of dust extinction. We conclude that dust extinction is at least one of the physical causes of the SED reddening (Fig. 3, left), although galaxy SED contributions might enhance the reddening (e.g. Bongiorno et al. 2012).

### 2.4.2 The polar-dust model

From Section 2.4.1, it is necessary to account for dust extinction of type 1 AGNs in X-CIGALE. The geometry of the obscuring materials is sketched in Fig. 4, where the materials responsible for type 1 AGN obscuration are called ‘polar dust’ (e.g. Lyu & Rieke 2018). The existence of polar dust has been proved by high-resolution mid-IR (MIR) imaging of local Seyfert galaxies (e.g. López-Gonzaga et al. 2014; Stalevski et al. 2017, 2019; Asmus 2019). However, the physical properties of the polar dust could be complicated and vary for different objects. For example, it might be close to the dust-sublimation radius (~pc scale; e.g. Lyu & Rieke 2018) or on galactic scales (~kpc; e.g. Zou et al. 2019).

Considering these complexities, we do not build a grid of physical models and perform radiation-transfer simulations. Instead, we employ several empirical extinction curves, including those from

<footnote n="4">⁴This statement breaks if the AGN host galaxies are highly star forming in general. However, the AGN hosts tend to have normal levels of star formation activity, as shown by previous studies (e.g. Harrison et al. 2012; Stanley et al. 2015).</footnote>

<footnote n="3">³https://sites.google.com/site/skirtorus/sed-library</footnote>
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Left: Median SEDs of spectroscopic type 1 AGNs in SDSS (blue, optically selected) and COSMOS (orange, X-ray selected). The green curve shows a typical unobscured quasar SED. All the SEDs are normalized at 250 nm. The SDSS SED is similar to the typical quasar SED. The COSMOS SED is significantly redder than the SDSS SED, indicating the presence of dust extinction and/or host-galaxy contribution for the COSMOS sources. Right: The observed *r*-band magnitude distributions for SDSS and COSMOS type 1 AGNs. Both AGN samples are in the bins of *z* = 2–2.5 and *f*<sub>ν,2keV</sub> = 3–10 (10<sup>-7</sup> mJy). The SDSS AGNs are systematically brighter than the COSMOS AGNs, and this qualitative result also holds for other redshift/X-ray flux bins. This result indicates that dust extinction should be at least one of the causes for the red SED shapes of the COSMOS AGNs, although galaxy SED contributions might enhance the reddening. This figure does not compare the AKARI-NEP sample in Section 3.3, because spectroscopic identifications of type 1 AGNs are not available for AKARI-NEP (Section 3.3.1).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 4.** Schematic plot (at the meridional plane) of the AGN obscuration system adopted by X-CIGALE (not to scale). The model in original CIGALE only includes torus obscuration. We add the obscuration of polar dust to account for type 1 AGN extinction.</figcaption></figure>

**Figure 5.** Same format as in Fig. 2, but for AGNs with different polar-dust $E(B-V)$. The top and bottom panels are for type 1 and type 2 AGNs, respectively. The SEDs are normalized at 10 µm. Except for $E(B-V)$, the other model parameters are the same (e.g. polar-dust temperature =100 K and emissivity =1.6) when generating these model SEDs.
<table><caption>**Table 1.** Sample properties.</caption><thead><tr><th>Name</th><th><i>N</i></th><th><i>m</i><sub>r</sub></th><th><i>f</i><sub>2–10 keV</sub></th><th>Redshift</th><th>Type<sub>AGN</sub></th></tr></thead><tbody><tr><td>SDSS</td><td>1986</td><td>19.2–20.9</td><td>2.2–10.5</td><td>0.6–1.9</td><td>1</td></tr><tr><td>COSMOS</td><td>590</td><td>21.3–23.7</td><td>0.4–1.7</td><td>0.6–1.8</td><td>1 and 2</td></tr><tr><td>AKARI-NEP</td><td>74</td><td>20.5–23.4</td><td>0.6–2.1</td><td>0.5–1.4</td><td>1 and 2</td></tr></tbody></table>

*Notes.* Columns (1) Survey name; (2) number of AGNs; (3) r-band AB magnitude range (20–80 per cent percentile); (4) 2–10 keV X-ray flux range (20–80 per cent percentile) in units of $10^{-14}$ erg s$^{-1}$ cm$^{-2}$; (5) redshift range (20–80 per cent percentile); (6) types of AGNs included in the survey.

Calzetti et al. (2000, nearby star-forming galaxies), Gaskell et al. (2004, large dust grains), and Prevot et al. (1984, Small Magellanic Cloud, SMC), and the user can choose among these curves. In our tests of X-CIGALE (Sections 3 and 4), we adopt the SMC extinction curve, which is preferred from AGN observations (e.g. Hopkins et al. 2004; Salvato et al. 2009; Bongiorno et al. 2012; but also see e.g. Gaskell et al. 2004). The extinction amplitude (parameterized as $E(B-V)$) is a free parameter set by the user, and setting $E(B-V) = 0$ returns to the original torus.

Since the scheme of X-CIGALE maintains energy conservation, we need to implement dust emission to account for the radiative energy absorbed by the dust. We assume the dust reemission is isotropic. For simplicity, we adopt the ‘grey body’ model (e.g. Casey 2012), i.e.

$$L_{\nu}(\lambda) \propto \frac{(1 - e^{-(\lambda_0/\lambda)^\beta}) (\frac{c}{\lambda})^3}{e^{hc/\lambda kT} - 1}, \quad (10)$$

where $\lambda_0$ is fixed at 200 µm, emissivity ($\beta$) and temperature ($T$) are free parameters set by the user. $L_\nu$ in equation (10) is normalized so that total energy is conserved, i.e.

$$L_{\text{total}}^{\text{emit}} = \int_{0}^{90^{\circ} - \Delta} L_{\text{total}}^{\text{extinct}}(\theta) \sin \theta d\theta \quad (11)$$

where $L_{\text{total}}^{\text{emit}}$ is the dust reemitted luminosity (angle-independent) and $L_{\text{total}}^{\text{extinct}}$ is the luminosity loss caused by polar-dust extinction (angle-dependent). Note that the integral on the right-hand side of equation (11) is to account for the fact that the polar dust only accounts for the obscuration in the polar directions while the polar-dust reemission is in all directions (see Fig. 4; e.g. equation 8). Fig. 5 shows the model SEDs for different extinction levels, where $T = 100$ K, $\beta = 1.6$, and $\Delta = 40^{\circ}$.

Our model above follows the AGN-unification scheme, i.e. AGN type is determined solely by the viewing angle, which is a free parameter in X-CIGALE. When the viewing angle is within the polar directions (type 1), the observed AGN disc emission suffers moderate (or none if $E(B-V) = 0$) extinction from the polar dust. When the viewing angle is within the equatorial directions (type 2), the observed AGN disc emission is strongly obscured by the torus. If the AGN type is known (e.g. from spectroscopy), the user can limit the viewing angle to the polar or equatorial direction (Sections 3.1 and 3.2). Otherwise, the user can adopt multiple viewing angles including both polar and equatorial directions, and let X-CIGALE choose freely between them (Section 3.3). For example, if the user set ‘viewing angles = 30°, 70°; $E(B-V) = 0.1$’, then CIGALE will build two model SEDs. For the 30° model (type 1), the UV/optical SED is reddened by the $E(B-V) = 0.1$ polar dust whose reemission also contributed to the IR SED. For the 70° model (type 2), the polar dust does not affect the UV/optical SED (already obscured by torus), but its reemission still contributes to the IR SED.

Our polar-dust model above provides one possible scenario for the reddened type 1 AGNs, i.e. the viewing angle is small and the extinction is caused by dust along the polar directions. An alternative scenario is that the line-of-sight (LOS) intercepts the torus, but the extinction is only moderate by chance due to the inhomogeneity of torus. However, this scenario has not been well investigated with physical torus models in the literature, to our knowledge. Therefore, we focus on the polar-dust model in the current version of X-CIGALE, and future versions of X-CIGALE may include this alternative scenario when its SED templates are available.

# 3 TESTS ON X-RAY DETECTED SOURCES

In this section, we test X-CIGALE with three samples of AGNs, i.e. SDSS (Section 3.1), COSMOS (Section 3.2), and AKARI-NEP (Section 3.3). The basic properties of these samples are summarized in Table 1. These three samples have different characteristics. The SDSS sample is optically bright type 1 quasars. The COSMOS sample is X-ray selected AGNs with broad multiwavelength coverage from u to Herschel/PACS 160 µm. This sample also has spectroscopic AGN classifications. The AKARI-NEP sample is small but has excellent MIR observations from AKARI.

## 3.1 SDSS

### 3.1.1 The sample and the models

The SDSS sample is optically selected from the DR14 quasar catalogue (Pâris et al. 2018). All the sources are spectroscopically confirmed type 1 AGNs. In addition to the SDSS ugriz bands, the Pâris et al. (2018) catalogue provides X-ray data from XMM-Newton archival observations when available (the 3XMM catalogue; Rosen et al. 2016). We require the sources to be detected in the 2–12 keV band at $>3\sigma$ significance levels. Here, the choice of the hard X-ray band (2–12 keV) is to minimize the effects of X-ray obscuration (see Section 2.2.1).<fnref n="5" /> We do not include the IR photometry compiled in the DR14 catalogue, because our main goal here is to test X-CIGALE on the simple cases, i.e. the quasar-dominated SED. In the X-ray to optical wavelengths, the AGN component is dominant; but in the IR wavelengths, the galaxy component may be non-negligible. We discuss the cases of AGN-galaxy mixed SEDs in Sections 3.2 and 3.3. We require the sources to have Galactic extinctions estimated in the DR14 quasar catalogue, because we need to correct for Galactic extinction before providing the photometry to X-CIGALE. These criteria lead to a final sample of 1986 AGNs (Table 1).

For the SDSS sample, we can neglect the galaxy SED component, because the sources are optically bright quasars which often dominate the observed UV/optical SEDs. The AGN-dominant (X-ray to IR) models in X-CIGALE can be achieved by setting $frac_{AGN}$ to a value close to unity (e.g. 0.999).<fnref n="6" /> The adopted AGN model parameters are listed in Table 2. The only free parameter in our fitting is polar-dust $E(B-V)$, which affects the UV/optical SED shape. We further justify that it is necessary to have $E(B-V)$ as a free parameter in Section 3.1.2. Other SKIRTOR parameters are fixed, because they only affect the IR SED shape where there is

<footnote n="5">⁵Alternatively, one can perform X-ray spectral fitting to obtain the absorption-corrected X-ray fluxes. However, extracting and analysing the X-ray spectra from the public XMM-Newton archival data are beyond the scope of this work.</footnote>

<footnote n="6">⁶Due to a technical reason, this value cannot equal to 1.</footnote>
<table><caption>**Table 2.** SDSS fitting parameters.</caption><thead><tr><th>Module</th><th>Parameter</th><th>Values</th></tr></thead><tbody><tr><td>AGN (UV-to-IR):</td><td>Torus optical depth at 9.7 microns τ<sub>9.7</sub></td><td>7.0</td></tr><tr><td rowspan="7">SKIRTOR</td><td>Torus density radial parameter p (ρ ∝ r<sup>-p</sup>e<sup>-q|cos(θ)|</sup>)</td><td>1.0</td></tr><tr><td>Torus density angular parameter q (ρ ∝ r<sup>-p</sup>e<sup>-q|cos(θ)|</sup>)</td><td>1.0</td></tr><tr><td>Angle between the equatorial plane and edge of the torus Δ</td><td>40°</td></tr><tr><td>Ratio of the maximum to minimum radii of the torus</td><td>20</td></tr><tr><td>Viewing angle θ (face on: θ = 0°, edge on: θ = 90°)</td><td>30°</td></tr><tr><td>AGN fraction in total IR luminosity frac<sub>AGN</sub></td><td>0.999</td></tr><tr><td>Extinction law of polar dust</td><td>SMC</td></tr><tr><td></td><td>E(B − V) of polar dust</td><td>0, 0.05, 0.1, 0.15, 0.2, 0.3</td></tr><tr><td></td><td>Temperature of polar dust (K)</td><td>100</td></tr><tr><td></td><td>Emissivity of polar dust</td><td>1.6</td></tr><tr><td>X-ray:</td><td>AGN photon index</td><td>1.8</td></tr><tr><td rowspan="3">This work</td><td>Maximum deviation from the α<sub>ox</sub>−L<sub>2500Å</sub> relation |Δα<sub>ox</sub>|<sub>max</sub></td><td>0.2</td></tr><tr><td>LMXB photon index</td><td>1.56</td></tr><tr><td>HMXB photon index</td><td>2.0</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** A random example of fitted SED from SDSS. The fitting residuals (observed flux/model flux) are shown at the bottom. Some source properties are labelled on the plot. Since the fitting only uses the AGN component, the total SED is actually the same as the AGN SED.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 8.** Top: Comparison of AGN intrinsic L<sub>2500Å</sub> between the fitting with X-ray versus without X-ray for the SDSS sample. Here, the L<sub>2500Å</sub> values are Bayesian-like estimates of X-CIGALE fitting. The solid black lines indicate the 1:1 relation; the dashed black lines indicate 0.3 dex deviation from the 1:1 relation. Bottom: Distribution histogram of the L<sub>2500Å</sub>,x/L<sub>2500Å</sub>,nox ratio for the SDSS sample. The histogram is normalized such that the integral is unity.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 7.** The distribution of AGN polar-dust E(B − V) (Bayesian-like estimate) from X-CIGALE fitting for SDSS (blue) and COSMOS (orange) spectroscopic type 1 AGNs. The histogram is normalized such that the integral is unity. The COSMOS sources tend to have higher E(B − V) than the SDSS sources.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Left: Estimated versus true values in the mock analyses of the SDSS sample for polar-dust E(B − V). The black line indicates the 1:1 relation between estimated and true values. For clarity, the data points are from 100 randomly selected sources in the SDSS sample. The values on the x-axis are shifted slightly for display purposes only. The residuals (defined as estimated minus true) are displayed on the bottom, with the black line indicating a zero residual. For E(B − V), the estimated and true values are well correlated, indicating that E(B − V) model parameter can be effectively constrained by the data. Right: Same format as left but from the fitting with free viewing angle.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Figure 10.** The E(B − V) PDF for the example source in Fig. 6. The PDF is normalized such that its integral is unity.</figcaption></figure>

no band coverage for the SDSS sample (see Section 3.3.2 for the
assessment of these parameters).

For the X-ray module, we adopt Γ = 1.8 for AGN (the dominant
component in X-rays), the typical intrinsic photon index constrained
by observations (Section 2.2.2). Adopting other AGN Γ values
(e.g. 1.4 or 2.0) do not affect our fitting results significantly. Our
adopted Γ = 1.8 is slightly different from that assumed in the
3XMM catalogue (Γ = 1.7; Rosen et al. 2016). Therefore, we
scale the 2–12 keV fluxes by a factor of 0.96 to correct the effects
of different Γ, and this correction factor is obtained using the
PIMMS website.<fnref n="7" /> For the LMXB and HMXB components, we
set Γ = 1.56 and 2.0, respectively (see Section 2.2.2). We adopt
|Δα<sub>ox</sub>|<sub>max</sub> = 0.2, and this |Δα<sub>ox</sub>|<sub>max</sub> value is ≈2σ scatter of the
α<sub>ox</sub>−L<sub>2500 Å</sub> relation (Section 2.2.3). Note that although the X-ray
module has both parameters fixed, x-cigale internally calculates
nine models of different α<sub>ox</sub> values and selects |Δα<sub>ox</sub>| ≤ |Δα<sub>ox</sub>|<sub>max</sub>
(see Section 2.2.3).

3.1.2 Fitting results

We run X-CIGALE with the model settings in Section 3.1.1 for the SDSS sample. The median reduced χ² (χ²<sub>red</sub>) and degrees of freedom (d.o.f.) are 1.4 and 5, respectively. These χ²<sub>red</sub> and d.o.f. values correspond to a *p*-value of 23 per cent, well above the conventional 2σ (5 per cent) or 3σ (0.3 per cent) values. This result indicates that X-CIGALE is able to model the observed photometry of the SDSS quasars. Fig. 6 shows a random example of the SED fitting. Fig. 7 displays the E(B − V) distribution from the fitting. As expected (see Section 2.4.1), most (75 per cent) SDSS AGNs have weak or no extinction with E(B − V) ≤ 0.1.

To evaluate the effects of the new X-ray module, we re-run x-
cigale but without this module. We compare the AGN intrinsic
L_{2500\AA} between the fitting with X-ray (L_{2500\AA,X}) versus without
X-ray (L_{2500\AA, noX}) in Fig. 8. The L_{2500\AA,X} and L_{2500\AA, noX} are
similar, and this similarity is as expected. SDSS sources are
mostly unobscured type 1 AGNs due to their selection method
(Section 2.4.1), and thus the intrinsic AGN emission is directly
observable at UV/optical wavelengths. Therefore, adding the X-ray
module does not significantly change L_{2500\AA} estimation for SDSS
sources in general.

In Table 2, the only free parameter is polar-dust E(B − V), because
this parameter affects the UV/optical SED which is covered by the
SDSS bands (Section 3.1.1). In other words, we consider that E(B −
V) can be constrained by the photometric data. The constrainability
of a model parameter can be evaluated by the ‘mock analysis’ of
X-CIGALE, which already exists in the previous version of CIGALE
(see section 4.3 of Boquien et al. 2019 for details). Briefly, after
fitting the observed data, X-CIGALE simulates a mock catalogue
based on the best-fitting model for each object. The photometric
uncertainties are considered when simulating the mock data. x-
CIGALE then performs SED fitting to the mock catalogue and obtains
Bayesian-like estimated values (and their errors) of the parameter.
By comparing these estimated values and those used to generate the
mock catalogue (i.e. the ‘true’ values), one can assess whether the
parameter can be reliably constrained. The mock analysis serves
as a sanity check to assess whether a physical parameter can be
retrieved in a self-consistent way. This mock analysis can be invoked
by setting ‘mock_flag=True’ in the X-CIGALE configurations.

<footnote n="7">⁷http://cxc.harvard.edu/toolkit/pimms.jsp</footnote>
<table><caption>**Table 3.** COSMOS and AKARI-NEP fitting parameters.</caption><thead><tr><th>Module</th><th>Parameter</th><th>Values</th></tr></thead><tbody><tr><td>Star formation history:</td><td>e-folding time, τ (Gyr)</td><td>0.1, 0.5, 1, 5</td></tr><tr><td>Delayed model, SFR ∝ texp(−t/τ)</td><td>Stellar Age, t (Gyr)</td><td>0.5, 1, 3, 5, 7</td></tr><tr><td>Simple stellar population:</td><td>Initial mass function</td><td>Chabrier (2003)</td></tr><tr><td>Bruzual &amp; Charlot (2003)</td><td>Metallicity (Z)</td><td>0.02</td></tr><tr><td>Galactic dust attenuation:</td><td>E(B − V) of starlight for the young population</td><td>0.1, 0.2, 0.3, 0.4, 0.5, 0.7, 0.9</td></tr><tr><td>Calzetti et al. (2000) &amp; Leitherer et al. (2002)</td><td>E(B − V) ratio between the old and young populations</td><td>0.44</td></tr><tr><td>Galactic dust emission: Dale et al. (2014)</td><td>α slope in dM<sub>dust</sub> ∝ U<sup>−α</sup>dU</td><td>1.5, 2.0, 2.5</td></tr><tr><td>AGN (UV-to-IR):</td><td>Torus optical depth at 9.7 microns τ<sub>9.7</sub></td><td>7.0</td></tr><tr><td rowspan="7">SKIRTOR</td><td>Torus density radial parameter p (ρ ∝ r<sup>−p</sup>e<sup>−q|cos(θ)|</sup>)</td><td>1.0</td></tr><tr><td>Torus density angular parameter q (ρ ∝ r<sup>−p</sup>e<sup>−q|cos(θ)|</sup>)</td><td>1.0</td></tr><tr><td>Angle between the equatorial plan and edge of the torus</td><td>40°</td></tr><tr><td>Ratio of the maximum to minimum radii of the torus</td><td>20</td></tr><tr><td>Viewing angle θ (face on: θ = 0°, edge on: θ = 90°)<sup>a</sup></td><td>30° (type 1), 70° (type 2)</td></tr><tr><td>AGN fraction in total IR luminosity frac<sub>AGN</sub></td><td>0.01, 0.1–0.9 (step 0.1), 0.99</td></tr><tr><td>Extinction law of polar dust</td><td>SMC</td></tr><tr><td>E(B − V) of polar dust</td><td>0, 0.05, 0.1, 0.15, 0.2, 0.3</td></tr><tr><td>Temperature of polar dust (K)</td><td>100</td></tr><tr><td>Emissivity of polar dust</td><td>1.6</td></tr><tr><td>X-ray:</td><td>AGN photon index Γ</td><td>1.8</td></tr><tr><td rowspan="3">This work</td><td>Maximum deviation from the α<sub>ox</sub>−L<sub>2500 Å</sub> relation</td><td>0.2</td></tr><tr><td>LMXB photon index</td><td>1.56</td></tr><tr><td>HMXB photon index</td><td>2.0</td></tr></tbody></table>

We run the mock analysis to test if polar-dust $E(B-V)$ can be constrained. We compare the estimated and true values in Fig. 9 (left). The estimated and true values are well correlated, indicating that $E(B-V)$ can be self-consistently constrained. In Fig. 10, we show the PDF of $E(B-V)$ for the source in Fig. 6. Fig. 10 indicates that the $E(B-V)$ is indeed well constrained in the Bayesian-like analysis.

In our fitting, aside from the model normalization (automatically determined by X-CIGALE; see section 4.3 of Boquien et al. 2019), $E(B-V)$ is the only free model parameter (Table 2). We also test freeing other parameters such as viewing angle and torus optical depth, and the mock-analysis results of $E(B-V)$ are similar. For example, Fig. 9 (right) shows the result after setting the viewing angle to 0–90° with a step of 10° (i.e. all allowed values). The estimated and true values are still well correlated, indicating that $E(B-V)$ and viewing angle are not strongly degenerate. This non-degeneracy is understandable, because, in our polar-dust model, $E(B-V)$ is the only parameter responsible for modelling the observed UV/optical SED shapes of type 1 AGNs like the SDSS objects (Section 2.4).

## 3.2 COSMOS

### 3.2.1 The sample and the models

The COSMOS sample is X-ray selected (>3σ, 2–10 keV band) from the COSMOS-Legacy survey performed by Chandra (Civano et al. 2016). The COSMOS-Legacy catalogue assumes Γ = 1.4. Similarly as in Section 3.1.1, we apply a correction factor of 0.87 (calculated with PIMMS) to the 2–10 keV fluxes to make them consistent with our adopted Γ = 1.8 for AGN. Marchesi et al. (2016) matched these X-ray sources with the optical/NIR counterparts in

the COSMOS2015 catalogue (Laigle et al. 2016) and compiled their spectroscopic information when available. We select the sources with spectroscopic classifications of AGN types. We adopt the photometric data in the COSMOS2015 catalogue, including 14 broad-bands from u to IRAC 8.0 µm. In addition, when available, we also include photometric data from Spitzer/MIPS (24 µm) and Herschel/PACS (100 and 160 µm), from the PEP survey (Lutz et al. 2011). We adopt the redshift measurements from Marchesi et al. (2016). These redshifts are either secure spectroscopic redshifts or high-quality photometric redshifts. There are a total of 590 objects in COSMOS (Table 1). Among these 590 objects, 206 and 384 are type 1 and type 2 AGNs, respectively.

The X-CIGALE model parameters for COSMOS are listed in Table 3. For the SKIRTOR and X-ray modules, the parameter setting is the same as in Section 3.1.1 except for frac<sub>AGN</sub> and the viewing angle. Here, we allow frac<sub>AGN</sub> to vary among 0.01, 0.1–0.9 (step 0.1), and 0.99, because, unlike in the case of SDSS, the AGN contribution to the observed SED may not be generally not dominant for the COSMOS sample, especially in the IR bands. We set the viewing angle to 30° and 70° for the spectroscopic type 1 and type 2 AGNs, respectively. These values are approximately the probability-weighted θ for type 1 and type 2 AGNs, respectively, given a torus of Δ = 40° (see Section 2.2.3).

For the galaxy component, we adopt the model setting similar to that in Ciesla et al. (2015). Specifically, we adopt a delayed star formation history (SFH), because it can characterize the SEDs of both early-type and late-type galaxies reliably (e.g. Ciesla et al. 2015; Boquien et al. 2019). Also, the delayed SFH only has relatively small parameter space (only two free parameters) and thereby high fitting efficiency. We adopt a Chabrier (2003) IMF with metallicity (Z) fixed to the solar value of 0.02. For the galactic

<footnote><sup>a</sup>For COSMOS, the viewing angles are set to 30° and 70° for the spectroscopic type 1 and type 2 AGNs, respectively. For AKARI-NEP, we allow X-CIGALE to choose between 30° and 70° for the entire sample, since spectroscopic classifications are not available.</footnote>
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Same format as in Fig. 6 but for COSMOS type 1 (top) and type 2 (bottom). The orange and blue curves indicate galaxy and AGN model components.</figcaption></figure>

dust attenuation, we adopt the dustatt calzleit module in X-CIGALE (Calzetti et al. 2000; Leitherer et al. 2002). The allowed $E(B-V)$ values for young stars are 0.1, 0.2, 0.3, 0.4, 0.5, 0.7, and 0.9. The $E(B-V)$ ratio between the old and young stars is fixed to 0.44. The amplitude of the 217.5 nm UV bump on the extinction curve is set to 0 (SMC) and 3 (Milky Way). We adopt the Dale et al. (2014) model for galactic dust reemission. There is only one free parameter in this model, i.e. the $\alpha$ slope in $dM_{dust} \propto U^{-\alpha}dU$, where $M_{dust}$ and $U$ are dust mass and radiation-field intensity, respectively. The $\alpha$ values are set to 1.5, 2.0, and 2.5.

### 3.2.2 Fitting results

We run X-CIGALE with the parameter settings in Section 3.2.1. The median $\chi^2_{red}$ values are 1.4 and 0.9, for type 1 and type 2 AGNs, respectively, while the median d.o.f. are 15 for both types. These median $\chi^2_{red}$ and d.o.f. corresponding to *p*-values of 0.12 and 0.59. These relatively large *p*-values for both type 1 and type 2 indicates that our models (Section 2) are able to model AGN SEDs of different types. This result supports the AGN-unification scheme (Section 1), on which our models are based. Fig. 11 displays two examples of

<figure><img src="image_14.png" /><figcaption>**Figure 12.** Same format as Fig. 8 but for the COSMOS sample. The blue and orange colours indicate type 1 and type 2 sources, respectively. The differences between $L_{2500\AA, X}$ and $L_{2500\AA, noX}$ are larger compared to those in SDSS (Fig. 8).</figcaption></figure>

the SED fitting in COSMOS. In Fig. 7, we compare the polar-dust $E(B-V)$ of type 1 AGNs in COSMOS versus SDSS. The COSMOS type 1 AGNs tend to have higher $E(B-V)$ than SDSS type 1 AGNs, consistent with the diagnostic in Section 2.4.1.

Fig. 12 compares $L_{2500\AA, X}$ and $L_{2500\AA, noX}$ for the COSMOS sample. The differences between $L_{2500\AA, X}$ and $L_{2500\AA, noX}$ are larger compared to those in SDSS (Fig. 8). This is because, for the SDSS sources, the observed optical fluxes are dominated by the AGN component, and thus AGN power can be effectively constrained even without X-ray data. In contrast, for COSMOS, the observed optical-to-IR fluxes are often not dominated by AGN, and X-CIGALE needs to decompose the fluxes into galaxy and AGN components. This SED decomposition process may be sometimes difficult, given that different models could result in similar model fluxes in optical-to-IR SED. Therefore, the X-ray data, which is often dominated by AGN, can be helpful in constraining the AGN power. Fig. 13 shows an example type 2 AGN SED fitted with versus without X-ray. For this source, the observed UV-to-IR fluxes are dominated by the galaxy component, and thus AGN power cannot be effectively constrained without X-ray data.

In our fitting (Table 3), we set viewing angles at 30° and 70° for type 1 and type 2 AGNs, respectively. This model setting can be done for our COSMOS sample, where spectroscopic classification is available (Section 3.2.1). However, for many photometric surveys
<figure><img src="image_15.png" /><figcaption>**Figure 13.** An example COSMOS SED fitted with (top) versus without (bottom) X-ray data. For this type 2 source, the best-fit intrinsic AGN $L_{2500\AA}$ is significantly different in the two cases as labelled. The observed UV-to-IR fluxes are dominated by the galaxy component. Therefore, the intrinsic AGN power cannot be effectively constrained without the X-ray data.</figcaption></figure>

(e.g. Section 3.3.1), the AGN spectra are not available, and thus spectrum-based AGN-type classification cannot be performed. In this case, the X-CIGALE user can set both 30° and 70° for the viewing angle, and allow X-CIGALE to freely choose between them. We test this configuration with our entire COSMOS sample, including spectroscopic type 1 and type 2 AGNs. The other parameters are the same as in Table 3. For spectroscopic type 1 (type 2) AGNs, 70 per cent (28 per cent) and 30 per cent (72 per cent) sources have the best-fitting viewing angles of 30° and 70°, respectively. This means that, if one uses the best-fitting viewing angle to perform AGN-type classification (i.e. the SED-based classification), the correct rate will be roughly ≈ 70 per cent for both type 1 and type 2.

## 3.3 AKARI-NEP

### 3.3.1 The sample and the models

The AKARI-NEP sample is also X-ray selected (>3σ, 2–7 keV band) based on Chandra observations of the AKARI-NEP field (Krumpe et al. 2015). The Krumpe et al. (2015) catalogue assumes

<figure><img src="image_16.png" /><figcaption>**Figure 14.** Same format as in Fig. 11 but for two sources in AKARI-NEP. Notably, the MIR data are well fitted with our model.</figcaption></figure>

$\Gamma = 1.4$. Similarly as in Sections 3.1.1 and 3.2.1, we apply a correction factor of 0.94 (calculated with PIMMS) to the 2–7 keV fluxes to make them consistent with our adopted $\Gamma = 1.8$ (AGN). We match the X-ray with the multiwavelength catalogue compiled by Buat et al. (2015) using a 1 arcsec matching radius. This multiwavelength catalogue has 19 bands from *u* to Herschel/PACS (100 and 160 µm). Notably, these bands include an excellent set of 9-band MIR data from the AKARI telescope, allowing us to test X-CIGALE on the MIR wavelengths. The final sample has 74 sources (Table 1).

For AKARI-NEP, we adopt the same X-CIGALE fitting parameters as for COSMOS (see Table 3) except the viewing angle. For the viewing angle, we allow X-CIGALE to choose freely between 30° (type 1) and 70° (type 2), since spectroscopic classifications of AGN type are not available for the AKARI-NEP sample.

### 3.3.2 Fitting results

We run X-CIGALE with model parameters in Section 3.3.1. The median $\chi^2_{red}$ is 1.2 (median d.o.f. = 17). The resulting *p*-value is 0.27, indicating overall good fitting quality. Fig. 14 shows two example fitted SEDs in the AKARI-NEP sample. Note that the MIR data can be well fitted with our model. Fig. 15 compares
<figure><img src="image_17.png" /><figcaption>**Figure 15.** Same format as Fig. 8 but for the AKARI-NEP sample. The differences between $L_{2500\AA, X}$ and $L_{2500\AA, noX}$ are larger compared to those in SDSS (Fig. 8).</figcaption></figure>

$L_{2500\AA, X}$ and $L_{2500\AA, noX}$ for AKARI-NEP. The differences between $L_{2500\AA, X}$ and $L_{2500\AA, noX}$ are larger compared to those in SDSS (Fig. 8). The reason is similar as discussed in Section 3.2.2, i.e. SED decomposition is needed for AKARI-NEP and such decomposition might be ambiguous without X-ray.

### 3.3.3 Parameter constrainability and degeneracy

The AKARI-NEP sample is suitable for investigating the constrainability and degeneracy of AGN model parameters, thanks to its excellent coverage at MIR wavelengths (Section 3.3.1) where AGN-dust emission peaks (Fig. 5). In this section, we present discussion on model-parameter constrainability and degeneracy using the AKARI-NEP sample. Since the population in this sample is AGN-galaxy mixed systems in general (e.g. Fig. 14), the results below might not be applicable to some particular sources, for which the observed MIR emission is dominated by AGNs (e.g. hot dust-obscured galaxies, hot DOGs; Dey et al. 2008; Vito et al. 2018). For these AGN-dominated sources, the AGN parameters might be easier to constrain, as host-galaxy contributions to the MIR emission are negligible.

There are three parameters that have multiple values in our fitting (Table 3), i.e. viewing angle, $\text{frac}_{\text{AGN}}$, and polar-dust $E(B-V)$. The viewing angle determines AGN types (see Section 2.4.2), and we show that the spectroscopic AGN types can be recovered with ≈ 70 per cent accuracy in Section 3.2.2. For $\text{frac}_{\text{AGN}}$ and

$E(B-V)$, we run the mock analysis as a sanity check of their constrainability (Section 3.1.2), and the results are displayed in Fig. 16. For $\text{frac}_{\text{AGN}}$, the estimated and true values are generally correlated (median errors = 0.13), although some sources have relatively large uncertainties. Therefore, the relative IR emission strength between AGN and galaxy (as measured by $\text{frac}_{\text{AGN}}$) can be effectively constrained. In contrast, unlike the case of SDSS (Section 3.1.2), the estimated $E(B-V)$ is relatively flat as a function of true $E(B-V)$, indicating that $E(B-V)$ cannot be well constrained in general. This result is understandable. For SDSS, the SED is dominated by type 1 AGNs, and the $E(B-V)$ is directly related to the observed UV/optical SED shape. However, for AKARI-NEP, the SEDs are generally produced by both AGN and galaxy components, and the $E(B-V)$ cannot be determined directly from the UV/optical SED shape (or other SED features).

In our fitting (Table 3), most of the torus and polar-dust parameters such as $\tau_{9,7}$ and polar-dust temperature are fixed at single values. This is because these parameters are related to the MIR SED shape. Considering that model degeneracy is likely strong in the MIR, broad-band photometry data like AKARI-NEP might not be able to effectively constrain these parameters. Now, we test whether torus $\tau_{9,7}$ can be well constrained. The X-CIGALE configuration is the same as in Table 3 except that $\tau_{9,7}$ is allowed to vary among 3, 5, 7, 9, and 11 (all allowed values). The mock-analysis results are presented in Fig. 17. The estimated value is generally flat as a function of the true value, indicating that $\tau_{9,7}$ cannot be well constrained. We have also tested other fixed AGN parameters in Table 3 such as polar-dust temperature and torus opening angle, and found they cannot be effectively constrained either.

Now, we analyse the reasons why $\tau_{9,7}$ has generally large uncertainties, using the two sources in Fig. 14 as illustrative examples. These reasons also generally explain the large uncertainties of other unconstrained parameters. In Fig. 18, we show the $\tau_{9,7}$ PDFs and 2D probability density maps ($\tau_{9,7}$ versus $\text{frac}_{\text{AGN}}$). For J175535.47+660959.0, the $\text{frac}_{\text{AGN}}$ can be constrained to $0.46 \pm 0.13$. From the density map, at high $\text{frac}_{\text{AGN}}$ ($\approx 0.6$), the probability peaks at $\tau_{9,7} = 11$. However, at lower $\text{frac}_{\text{AGN}}$, the peak shifts to lower $\tau_{9,7}$. Therefore, the $\tau_{9,7}$ parameter is degenerate with $\text{frac}_{\text{AGN}}$, i.e. the probability distribution of $\tau_{9,7}$ depends on the value of $\text{frac}_{\text{AGN}}$. This degeneracy makes the marginalized $\tau_{9,7}$ PDF relatively flat, leading to the large uncertainty of this parameter. From the density map of J175520.20+660949.1, $\text{frac}_{\text{AGN}}$ is well constrained at a low level, but the $\tau_{9,7}$ PDF is flat. This is because when $\text{frac}_{\text{AGN}}$ is low, the observed MIR SED is dominated by the galaxy component (Fig. 14 bottom). In this case, the models with different $\tau_{9,7}$ are degenerate, in the sense that they have similar MIR SED shapes due to dominant galaxy contributions. The results above indicate that model degeneracy is responsible for the large uncertainties in the $\tau_{9,7}$ estimation.

## 4 FITTING WITH X-RAY UPPER LIMITS

In the previous section, we apply X-CIGALE to sources with X-ray detections. The X-ray sources are often only a small fraction ($\lesssim 10$ per cent) of the entire sample in optical/IR surveys. However, many extragalactic studies need to constrain AGN emission for the X-ray undetected majority galaxy population (e.g. Buat et al. 2015; Vito et al. 2016; Bowman et al. 2019). This task can be achieved by X-CIGALE, based on X-ray flux upper limits. The detailed fitting algorithm is presented in section 4.3 of Boquien et al. (2019), and we do not repeat it here. Below, we test this usage with 100 randomly selected AKARI-NEP galaxies (Buat et al. 2015), all of which are
<figure><img src="image_18.png" /><figcaption>**Figure 16.** Same format as in Fig. 9 but for $frac_{AGN}$ (left) and $E(B-V)$ (right) of the AKARI-NEP AGNs.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Figure 17.** Same format as in Fig. 9 but for torus $\tau_{9,7}$ of the AKARI-NEP AGNs. For $\tau_{9,7}$, the estimated value is relatively flat as a function of the true value, indicating that they cannot be well constrained by the observed data.</figcaption></figure>

X-ray undetected. Here, we choose AKARI-NEP rather than SDSS or COSMOS, since AKARI-NEP has the best multiwavelength coverage among the three surveys (Section 3). We focus on a relatively small sample (100), because the upper-limit analysis in X-CIGALE is time-consuming due to its complicated mathematical form of $\chi^2$ (see section 4.3 of Boquien et al. 2019).

Precise X-ray flux upper limits depend on source positions and vary from source to source. The derivation of the precise values requires intensive simulations (e.g. Xue et al. 2016; Luo et al. 2017), which are beyond the scope of this work. Instead, we adopt a single conservative value, $5.3 \times 10^{-15} \text{ erg s}^{-1} \text{ cm}^{-2}$ ($4.4 \times 10^{-7} \text{ mJy}$, 2–7 keV), for all the sources. For ≈ 80 per cent of the survey area, the actual *Chandra* flux limits should be lower than this value (see fig. 10 of Krumpe et al. 2015). We run X-CIGALE based on this upper limit. The parameter settings are the same as in Section 3.3.1 except that we allow $frac_{AGN} = 0$, since it is possible that the AGN component does not exist for these upper-limit sources. For comparison, we also re-run X-CIGALE without the X-ray upper limit.

Fig. 19 compares the Bayesian-like estimate of $frac_{AGN}$ of these two runs ($frac_{AGN, noX}$ versus $frac_{AGN, Xup}$). As expected, $frac_{AGN, Xup}$ is generally lower than $frac_{AGN, noX}$, since the X-ray upper limit can constrain AGN power. Notably, for ≈ 10 per cent of the sources, the

$frac_{AGN, Xup}$ is much lower ($\Delta frac_{AGN, Xup} > 0.2$) than $frac_{AGN, noX}$. This result indicates that, even when excellent MIR data are present (Section 3.3.1), the AGN-galaxy decomposition might still be inaccurate/ambiguous without X-ray data. On the other hand, there are still ≈ 30 per cent of sources that have non-negligible $frac_{AGN}$ ($>0.1$) when the X-ray upper limit is used in the fitting. It is possible that a non-negligible IR flux is contributed by the AGN. However, another possibility is that the current X-ray upper limit is too high to effectively constrain AGN emission. In the future, *Athena* may clarify this problem with its great sensitivity.

The total *Chandra* exposure time on the AKARI-NEP field is ≈300 ks. Given this amount of exposure time, *Athena* can reach its confusion-limited sensitivity of ~1 × $10^{-16} \text{ erg s}^{-1} \text{ cm}^{-2}$ ($8 \times 10^{-9} \text{ mJy}$, 2–7 keV) for the entire AKARI-NEP field.<fnref n="8" /> Assuming this X-ray flux limit, we re-run X-CIGALE. The resulting $frac_{AGN, Xup}$ is below 0.01 for all sources (see Fig. 19), indicating that the AGN SED contribution will be negligible if a source is undetected by *Athena*. Therefore, X-CIGALE, with future *Athena* observations, will have great power in unambiguously determining the presence of AGN. This feature will be extremely helpful for future extragalactic studies.

## 5 SUMMARY AND FUTURE PROSPECTS

We have developed and tested X-CIGALE, a new version of the galaxy SED fitting code, CIGALE. Our development and test results are summarized below.

(i) We have developed a new X-ray module (Section 2.2). The module is mainly designed to connect the intrinsic X-ray emission with other wavelengths, and X-ray obscuration and transmission should be corrected before providing the X-ray data to X-CIGALE. The X-ray module includes the X-ray emission from both galaxy and AGN. The galaxy component includes the emission from HMXB, LMXB, and hot gas. The AGN's X-ray SED is connected to its UV-to-IR SED using the well-known $\alpha_{ox}-L_{2500Å}$ relation.

(ii) We have implemented a modern torus model, SKIRTOR, to

<footnote n="8">⁸*Chandra* can, in principal, reach deeper sensitivity than *Athena* thanks to its superior angular resolution. However, reaching *Athena*-like (or deeper) flux limits will practically need large amounts of exposure time of *Chandra*. This has only been achieved in two small *Chandra* fields (only ~500 arcmin² each), i.e. 7 Ms CDF-S and 2 Ms CDF-N (Xue et al. 2016; Luo et al. 2017).</footnote>
<figure><img src="image_20.png" /><figcaption>**Figure 18.** Top: The $\tau_{9,7}$ PDF for the example source in Fig. 14. Bottom: The 2D probability density map of $\tau_{9,7}$ versus $\text{frac}_{\text{AGN}}$. The density map is normalized such that the 2D integral is unity.</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Figure 19.** Comparison of $\text{frac}_{\text{AGN}}$ between the fitting without versus with X-ray upper limits for 100 random galaxies in the AKARI-NEP field. The blue and orange points represent $\text{frac}_{\text{AGN,Xup}}$ obtained from *Chandra* and *Athena* upper limits, respectively. The solid black lines indicate the 1:1 relation. As expected, $\text{frac}_{\text{AGN,Xup}}$ is generally lower than $\text{frac}_{\text{AGN,NOX}}$, since the X-ray upper limit can constrain AGN power. The current *Chandra* upper limit can effectively constrain AGN power for some sources. The *Athena* upper limit strongly suppresses $\text{frac}_{\text{AGN,Xup}}$ to < 1 per cent for all sources.</figcaption></figure>

fit AGN UV-to-IR SEDs (Section 2.3). SKIRTOR adopts a clumpy two-phase torus, which is responsible for obscuring the UV/optical emission from the AGN disc. SKIRTOR is developed from a 3D radiative-transfer method, and thus obeys the energy-conservation law. However, SKIRTOR assumes that the AGN disc emission is absolutely unextincted when viewed from the polar direction. Therefore, SKIRTOR cannot model the SEDs of slightly extincted type 1 AGNs. To overcome this disadvantage, we introduce extinction from polar dust (Section 2.4). The extinction amplitude, $E(B - V)$, is a free model parameter set by the X-CIGALE user.

(iii) We have tested X-CIGALE on the AGNs with X-ray detections in SDSS, COSMOS, and AKARI-NEP Section 3. The three samples have distinctive characteristics in terms of AGN properties and available data. The fitting quality is good in general, with typical $\chi^2 \sim 1$ for all the samples. This result indicates that X-CIGALE is capable in modelling observed AGN SEDs under different circumstances. We also compare the fittings results with versus without X-ray data. We find that the resulting AGN power is sometimes different in the two cases, when both AGN and galaxy components are present. Therefore, the AGN-galaxy SED decomposition may be unphysical without the constraints from X-ray data. We discuss constrainability and degeneracy of model parameters in the fitting of AKARI-NEP, for which excellent mid-IR photometric coverage is available.
(iv) We also test X-CIGALE on a random sample of galaxies with only *Chandra* X-ray upper limits in the AKARI-NEP field, where excellent MIR data are available (Section 4). We compare the fitting results with and without the X-ray upper limits. After using the X-ray upper limits, frac<sub>AGN</sub> sometimes becomes lower, indicating that the current *Chandra* upper limit can effectively constrain AGN emission, as least for some systems. We also evaluate the potential of the future *Athena* mission by replacing the *Chandra* upper limit with the expected *Athena* value for a similar exposure time (≈300 ks). The resulting frac<sub>AGN</sub> is constrained to a negligible level (< 1 per cent) for all the sources, indicating that *Athena* can robustly constrain AGN emission in general with a moderate amount of exposure time (≤300 ks).

We publicly release X-CIGALE on the official website of CIGALE.<fnref n="9" /> As for the previous versions of CIGALE, X-CIGALE is open-source, allowing the user to modify the source code freely. In the future, we will further develop X-CIGALE and enable it to address special AGNs such as radio-loud and BAL objects (Section 2.2). Besides the three surveys tested in this work (SDSS, COSMOS, and AKARI-NEP), the user can apply X-CIGALE to the existing multiwavelength surveys such as CDF-S (Luo et al. 2017), CDF-N (Xue et al. 2016), and XMM-SERVS (Chen et al. 2018). In the future, X-CIGALE can be used to explore deep/wide surveys of, e.g. *eROSITA* (e.g. Merloni et al. 2012) and *Athena* (e.g. Nandra et al. 2013).

## ACKNOWLEDGEMENTS

We thank the referee for helpful feedback that improved this work. GY acknowledges support from the French Space Agency CNES. VB received funding from Excellence Initiative of Aix-Marseille University – AMIDEX, a French ‘Investissements d’Avenir’ programme. MS acknowledges support by the Ministry of Education, Science, and Technological Development of the Republic of Serbia through the projects Astrophysical Spectroscopy of Extragalactic Objects (176001) and Gravitation and the Large Scale Structure of the Universe (176003). WNB acknowledges support from NASA ADP grant 80NSSC18K0878 and the V.M. Willaman Endowment. This project uses ASTROPY (a PYTHON package; see Astropy Collaboration 2018).

## REFERENCES

* Antonucci R., 1993, ARA&A, 31, 473

* Arnaud K. A., 1996, in Jacoby G. H., Barnes J., eds, ASP Conf. Ser. Vol. 101, Astronomical Data Analysis Software and Systems V. Astron. Soc. Pac., San Francisco, p. 17

* Asmus D., 2019, MNRAS, 489, 2177

* Astropy Collaboration et al., 2018, AJ, 156, 123

* Baes M., Verstappen J., De Looze I., Fritz J., Saftly W., Vidal Pérez E., Stalevski M., Valcke S., 2011, ApJS, 196, 22

* Bongiorno A. et al., 2012, MNRAS, 427, 3103

* Boquien M., Burgarella D., Roehlly Y., Buat V., Ciesla L., Corre D., Inoue A. K., Salas H., 2019, A&A, 622, A103

* Bowman W. P. et al., 2019, ApJ, 875, 152

* Brandt W. N., Laor A., Wills B. J., 2000, ApJ, 528, 637

* Brandt W. N., Alexander D. M., 2015, A&AR, 23, 1

* Bruzual G., Charlot S., 2003, MNRAS, 344, 1000

* Buat V. et al., 2015, A&A, 577, A141

* Burgarella D., Buat V., Iglesias-Páramo J., 2005, MNRAS, 360, 1413

* Calzetti D., Armus L., Bohlin R. C., Kinney A. L., Koornneef J., Storchi-Bergmann T., 2000, ApJ, 533, 682

* Camps P., Baes M., 2015, Astron. Comput., 9, 20

* Casey C. M., 2012, MNRAS, 425, 3094

* Chabrier G., 2003, ApJ, 586, L133

* Chen C.-T. J. et al., 2018, MNRAS, 478, 2132

* Ciesla L. et al., 2015, A&A, 576, A10

* Civano F. et al., 2016, ApJ, 819, 62

* Dadina M., 2008, A&A, 485, 417

* Dale D. A., Helou G., Magdis G. E., Armus L., Díaz-Santos T., Shi Y., 2014, ApJ, 784, 83

* Dey A. et al., 2008, ApJ, 677, 943

* Duras F. et al., 2017, A&A, 604, A67

* Elvis M. et al., 2012, ApJ, 759, 6

* Fabbiano G., 2006, ARA&A, 44, 323

* Feltre A., Hatziminaoglou E., Fritz J., Franceschini A., 2012, MNRAS, 426, 120

* Freeman P., Doe S., Siemiginowska A., 2001, in Starck J.-L., Murtagh F. D., eds, Proc. SPIE Conf. Ser. Vol. 4477, Astronomical Data Analysis. SPIE, Bellingham, p. 76

* Fritz J., Franceschini A., Hatziminaoglou E., 2006, MNRAS, 366, 767

* Gaskell C. M., Goosmann R. W., Antonucci R. R. J., Whysong D. H., 2004, ApJ, 616, 147

* Harrison C. M. et al., 2012, ApJ, 760, L15

* Hinshaw G. et al., 2013, ApJS, 208, 19

* Hopkins P. F. et al., 2004, AJ, 128, 1112

* Ichikawa K., Ueda Y., Terashima Y., Oyabu S., Gandhi P., Matsuta K., Nakagawa T., 2012, ApJ, 754, 45

* Just D. W., Brandt W. N., Shemmer O., Steffen A. T., Schneider D. P., Chartas G., Garmire G. P., 2007, ApJ, 665, 1004

* Kormendy J., Ho L. C., 2013, ARA&A, 51, 511

* Kormendy J., Richstone D., 1995, ARA&A, 33, 581

* Krumpe M. et al., 2015, MNRAS, 446, 911

* Laigle C. et al., 2016, ApJS, 224, 24

* Leitherer C., Li I. H., Calzetti D., Heckman T. M., 2002, ApJS, 140, 303

* Liu T., Wang J.-X., Yang H., Zhu F.-F., Zhou Y.-Y., 2014, ApJ, 783, 106

* Liu T. et al., 2017, ApJS, 232, 8

* López-Gonzaga N., Jaffe W., Burtscher L., Tristram K. R. W., Meisenheimer K., 2014, A&A, 565, A71

* Luo B. et al., 2017, ApJS, 228, 2

* Lusso E., Risaliti G., 2017, A&A, 602, A79

* Lusso E. et al., 2012, MNRAS, 425, 623

* Lutz D. et al., 2011, A&A, 532, A90

* Lyu J., Rieke G. H., 2018, ApJ, 866, 92

* Marchesi S. et al., 2016, ApJ, 817, 34

* Mathews W. G., Brighenti F., 2003, ARA&A, 41, 191

* Merloni A. et al., 2012, preprint (arXiv:1209.3114)

* Mewe R., Lemen J. R., van den Oord G. H. J., 1986, A&AS, 65, 511

* Mezcua M., Civano F., Marchesi S., Suh H., Fabbiano G., Volonteri M., 2018, MNRAS, 478, 2576

* Miller B. P., Brandt W. N., Schneider D. P., Gibson R. R., Steffen A. T., Wu J., 2011, ApJ, 726, 20

* Motta S., Belloni T., Homan J., 2009, MNRAS, 400, 1603

* Nandra K. et al., 2013, preprint (arXiv:1306.2307)

* Netzer H., 1987, MNRAS, 225, 55

* Netzer H., 2013, The Physics and Evolution of Active Galactic Nuclei. Cambridge Univ. Press, Cambridge

* Netzer H., 2015, ARA&A, 53, 365

* Nikutta R., Elitzur M., Lacy M., 2009, ApJ, 707, 1550

* Noll S., Burgarella D., Giovannoli E., Buat V., Marcillac D., Muñoz-Mateos J. C., 2009, A&A, 507, 1793

* O’Dell S. L. et al., 2017, Siegmund O. H., in Proc. SPIE UV, X-Ray, and Gamma-Ray Space Instrumentation for Astronomy XVIII, 8859. SPIE, Bellingham, p. 103970C

* Páris I. et al., 2018, A&A, 613, A51

* Prevot M. L., Lequeux J., Maurice E., Prevot L., Rocca-Volmerange B., 1984, A&A, 132, 389

* Ricci C. et al., 2017, ApJS, 233, 17

<footnote n="9">⁹https://cigale.lam.fr/</footnote>
Richards G. T. et al., 2003, AJ, 126, 1131

Rosen S. R. et al., 2016, A&A, 590, A1

Salvato M. et al., 2009, ApJ, 690, 1250

Sazonov S., Khabibullin I., 2017, MNRAS, 468, 2249

Serra P., Amblard A., Temi P., Burgarella D., Giovannoli E., Buat V., Noll S., Im S., 2011, ApJ, 740, 22

Stalevski M., Fritz J., Baes M., Nakos T., Popović L. Č., 2012, MNRAS, 420, 2756

Stalevski M., Ricci C., Ueda Y., Lira P., Fritz J., Baes M., 2016, MNRAS, 458, 2288

Stalevski M., Asmus D., Tristram K. R. W., 2017, MNRAS, 472, 3854

Stalevski M., Tristram K. R. W., Asmus D., 2019, MNRAS, 484, 3334

Stanley F., Harrison C. M., Alexander D. M., Swinbank A. M., Aird J. A., Del Moro A., Hickox R. C., Mullaney J. R., 2015, MNRAS, 453, 591

Starling R. L. C., Willingale R., Tanvir N. R., Scott A. E., Wiersema K., O'Brien P. T., Levan A. J., Stewart G. C., 2013, MNRAS, 431, 3159

Steffen A. T., Strateva I., Brandt W. N., Alexander D. M., Koekemoer A. M., Lehmer B. D., Schneider D. P., Vignali C., 2006, AJ, 131, 2826

Tanimoto A., Ueda Y., Odaka H., Kawaguchi T., Fukazawa Y., Kawamuro T., 2019, ApJ, 877, 95

Urry C. M., Padovani P., 1995, PASP, 107, 803

Vito F. et al., 2016, MNRAS, 463, 348

Vito F. et al., 2018, MNRAS, 474, 4528

Xu Y.-D., 2015, MNRAS, 449, 191

Xue Y. Q., Luo B., Brandt W. N., Alexander D. M., Bauer F. E., Lehmer B. D., Yang G., 2016, ApJS, 224, 15

Yang G. et al., 2016, ApJ, 831, 145

Yang G. et al., 2017, ApJ, 842, 72

Yang G. et al., 2018a, MNRAS, 475, 1887

Yang G., Brandt W. N., Darvish B., Chen C.-T. J., Vito F., Alexander D. M., Bauer F. E., Trump J. R., 2018b, MNRAS, 480, 1022

Yang G., Brandt W. N., Alexander D. M., Chen C.-T. J., Ni Q., Vito F., Zhu F.-F., 2019, MNRAS, 485, 3721

Zhang S. N., 1997, in Wickramasinghe D. T., Bicknell G. V., Ferrario L., eds, ASP Conf. Ser. Vol. 121, IAU Colloq. 163: Accretion Phenomena and Related Outflows. Astron. Soc. Pac., San Francisco, p. 41

Zou F., Yang G., Brandt W. N., Xue Y., 2019, ApJ, 878, 11

APPENDIX A: NEW INPUTS AND OUTPUTS IN X-CIGALE

The input parameters for the new X-ray and SKIRTOR modules can be found in Table 3. After fitting, X-CIGALE can output the best-fitting model SEDs of different components, and the SED components for the new X-ray and SKIRTOR modules are summarized in Table A1. Besides the best-fitting SEDs, X-CIGALE can also yield the maximum-likelihood and Bayesian-like values of source physical properties. These physical properties include not only the model parameters (Table 3), but also some additional quantities as listed in Table A2. New quantities can be added in the future as requested by the user of X-CIGALE. We remind that the quantities, ‘agn.intrin_Lnu_2500A’ and ‘xray.alpha_ox’ refer to the values as measured at a viewing angle of 30° (see Section 2.2.3). The quantity ‘agn.accretion_power’ refers to the intrinsic (unextincted) AGN disc luminosity averaged over all directions (weighted by $\sin\theta$; see Section 2.2.3). This quantity, paired with an assumed radiative efficiency, can be used to estimate BH accretion rate (e.g. Yang et al. 2017, 2019).

<table><caption>Table A1. Output SED components for the new X-ray and SKIRTOR modules.</caption><thead><tr><th>Module</th><th>Component</th><th>Explanation</th></tr></thead><tbody><tr><td>X-ray (Section 2.2.2)<br>10<sup>-3</sup>–5 nm</td><td>xray.agn<br>xray.galaxy</td><td>The AGN corona<br>The total SED of HMXB, LMXB, and hot gas</td></tr><tr><td>SKIRTOR (Section 2.3)<br>8–10<sup>6</sup> nm</td><td>agn.SKIRTOR2016.disk<br>agn.SKIRTOR2016.dust</td><td>The AGN disc<br>The dust reemission</td></tr></tbody></table>

Note. In X-CIGALE output, all the SED components are in the format of $L_\lambda$ in units of W nm$^{-1}$.
<table><caption>**Table A2.** Additional output physical parameters for the new X-ray and SKIRTOR modules.</caption><thead><tr><th>Module</th><th>Parameters</th><th>Explanation</th><th>Units</th></tr></thead><tbody><tr><td rowspan="8">X-ray</td><td>xray.agn_Lnu_2keV</td><td>The AGN L<sub>v</sub> at 2 keV</td><td>W Hz<sup>-1</sup></td></tr><tr><td>xray.agn_Lx_2to10keV</td><td>The AGN 2–10 keV luminosity</td><td>W</td></tr><tr><td>xray.agn_Lx_total</td><td>The AGN total (0.25–1200 keV) X-ray luminosity</td><td>W</td></tr><tr><td>xray.alpha_ox</td><td>The AGN α<sub>ox</sub></td><td>—</td></tr><tr><td>xray.lmxb_Lx_2to10keV</td><td>The 2–10 keV LMXB luminosity</td><td>W</td></tr><tr><td>xray.hmxb_Lx_2to10keV</td><td>The 2–10 keV HMXB luminosity</td><td>W</td></tr><tr><td>xray.hotgas_Lx_0p5to2keV</td><td>The 0.5–2 keV hot-gas luminosity</td><td>W</td></tr><tr><td>SKIRTOR</td><td>agn.disk_luminosity</td><td>The observed AGN disc luminosity (might be extincted)</td><td>W</td></tr><tr><td></td><td>agn.dust_luminosity</td><td>The observed AGN dust reemitted luminosity</td><td>W</td></tr><tr><td></td><td>agn.luminosity</td><td>The sum of agn.disk_luminosity and agn.dust_luminosity</td><td>W</td></tr><tr><td></td><td>agn.intrin_Lnu_2500A</td><td>The intrinsic AGN L<sub>v</sub> at 2500 Å at viewing angle = 30°</td><td>W Hz<sup>-1</sup></td></tr><tr><td></td><td>agn.accretion_power</td><td>The intrinsic AGN disc luminosity averaged over all directions</td><td>W</td></tr></tbody></table>

This paper has been typeset from a TEX/LATEX file prepared by the author.