Skip to content
The following article is Open access

The Infrared Medium-deep Survey. VII. Faint Quasars at z ∼ 5 in the ELAIS-N1 Field

, , , , , , , , ,

Published 2020 April 14 © 2020. The Author(s). Published by the American Astronomical Society.
, , Citation Suhyun Shin et al 2020 ApJ 893 45DOI 10.3847/1538-4357/ab7bde

Download Article PDF
DownloadArticle ePub

You need an eReader or compatible software to experience the benefits of the ePub3 file format.

0004-637X/893/1/45

Abstract

The intergalactic medium (IGM) at z ∼  5 to 6 is largely ionized, and yet the main source for the IGM ionization in the early universe is uncertain. Of the possible contributors are faint quasars with $-26\lesssim {M}_{1450}\lesssim -23$, but their number density is poorly constrained at z ∼ 5. In this paper, we present our survey of faint quasars at z ∼ 5 in the European Large-Area Infrared Space Observatory Survey-North 1 (ELAIS-N1) field over a survey area of 6.51 deg2 and examine if such quasars can be the dominant source of the IGM ionization. We use the deep optical/near-infrared data of the ELAIS-N1 field as well as the additional medium-band observations to find z ∼ 5 quasars through a two-step approach using the broadband color selection, and spectral energy distribution fitting with the medium-band information included. Adopting Bayesian information criterion, we identify 10 promising quasar candidates. Spectra of three of the candidates are obtained, confirming all of them to be quasars at z ∼ 5 and supporting the reliability of the quasar selection. Using the promising candidates, we derive the z ∼ 5 quasar luminosity function at −26 ≲ M1450 ≲ −23. The number density of faint z ∼ 5 quasars in the ELAIS-N1 field is consistent with several previous results that quasars are not the main contributors to the IGM-ionizing photons at z ∼ 5.

Export citation and abstractBibTeXRIS

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

1. Introduction

One of the most pivotal cosmic events after the big bang is cosmic reionization. Cosmic reionization is the change in the ionization state of the intergalactic medium (IGM), where the neutral atoms in the IGM become ionized after the recombination era due to ultraviolet (UV) radiation from newly born stars, galaxies, and active galactic nuclei (AGNs). Hence, the cosmic reionization has become an important subject of study for understanding the history of astronomical objects in the universe. The recent analysis of the cosmic microwave background suggests that midpoint of reionization is at z ∼ 7.7 (Planck Collaboration et al. 2018), corresponding to the most distant quasars discovered so far (Bañados et al. 2018; Davies et al. 2018). The observational constraints from high-redshift observation indicate the cosmic reionization ends by z ∼ 6 (Bouwens et al. 2015).

However, the nature of the sources responsible for keeping the IGM ionized in the postreionization era remains uncertain. Previous studies have suggested either star-forming galaxies or quasars as candidates to explain the majority of the required photon budget (Fontanot et al. 2012; Giallongo et al. 2015; hereafter G15; Madau & Haardt 2015; Bouwens et al. 2016; D’Aloisio et al. 2017; Ricci et al. 2017; Kakiichi et al. 2018). In some works, galaxies are considered as the primary contributors to the IGM reionization (Alvarez et al. 2012; Paardekooper et al. 2013; Kashikawa et al. 2015; Kim et al. 2015, 2019; Matsuoka et al. 2018; Smith et al. 2018; Fletcher et al. 2019). In contrary, the escape fraction of the hydrogen-ionizing Lyman-continuum (LyC) photons of star-forming galaxies is found to be only few percent (Bridge et al. 2010; Rutkowski et al. 2016; Vasei et al. 2016; Japelj et al. 2017; Iwata et al. 2019), and this makes it challenging for galaxies to maintain the ionized state of the IGM (Finkelstein et al. 2015; Grazian et al. 2017; Naidu et al. 2018; Steidel et al. 2018; Lam et al. 2019).

On the other hand, the escape fraction of LyC photons in quasars is about 100% (Cristiani et al. 2016; Guaita et al. 2016; Grazian et al. 2018; Romano et al. 2019), hence several studies suggest that quasars are the main IGM-ionizing source (Feng et al. 2016, G15; Boutsia et al. 2018; Romano et al. 2019). Given the spectral shape and the escape fraction of LyC photons of quasars, the contribution of quasars to the hydrogen reionizing photons can be calculated from the quasar luminosity function (LF). Therefore, many observational studies have tried to determine the LF of quasars at z > 4 (Barger et al. 2003; Fan et al. 2006; Jiang et al. 2008; Willott et al. 2010; McGreer et al. 2013; Bañados et al. 2016; Yang et al. 2016; Jeon et al. 2017; Kim et al. 2019).

Although still uncertain, the quasar LFs show that the number of LyC photons is the most sensitive to the number of quasars with moderate luminosity at −24 < M1450 < −23. Interestingly, some studies (G15; Boutsia et al. 2018; Giallongo et al. 2019; hereafter G19) find that the number density of the moderate luminosity quasars is higher by a factor of ∼10 than other studies (Onoue et al. 2017; Akiyama et al. 2018; Matsuoka et al. 2018; McGreer et al. 2018; hereafter M18). The high number density in the moderate luminosity range enables quasars to be the main sources for the IGM ionization. The conflict in the high-redshift quasar LF reflects the lack of quasar samples with M1450 < −24. It is, therefore, imperative to add more samples to help settle the issue.

To constrain the quasar LF better and understand the IGM ionization process, we have carried out the Infrared Medium-deep Survey (IMS; M. Im et al. 2020, in preparation), a deep and wide near-infrared (NIR) survey to depths of J ∼ 23.5 AB mag and with an area coverage of ∼100 deg2. Color-selection criteria with the deep NIR data can discern reliable high-redshift quasars from dwarf stars, the most contaminant source for quasar candidate selection (McGreer et al. 2013). Additionally, the medium-band imaging of this survey improves the spectral energy distribution (SED) fittings of a quasar and a M-dwarf star. Adopting the Bayesian information criterion (BIC) which can consider the number of free parameters in a model, we are able to compare the two models and select the promising quasar candidates effectively. In this paper, we focus on our discovery of faint quasars at z ∼ 5 in the European Large-Area Infrared Space Observatory (ISO) Survey-North 1 (ELAIS-N1) field. The extensive multiwavelength coverage of the ELAIS-N1 field, including the NIR from IMS, the optical from the Hyper Suprime-Cam Subaru Strategic Program (HSC-SSP) and medium-bands, makes it possible to select faint z ∼ 5 quasars effectively.

This paper is structured as the following. Section 2 describes the SED models for quasars and M dwarfs. In Section 3, we present the broadband and medium-band imaging data used for the quasar selection. In Section 4, we explain two quasar selection methods, the selection criteria based on broadband colors and BIC selection based on the SED fitting. Follow-up spectroscopy of quasar candidates using the MMT is presented in Section 5. In Section 6, we build the quasar LF and discuss the role of quasars in the IGM reionization. Finally, we summarize our findings. Throughout this paper, we adopt the AB magnitude system for all filters (Oke & Gunn 1983) and assume H0 = 70 km s−1 Mpc−1, ΩM = 0.3, and ΩΛ = 0.7 of the concordance ΛCDM cosmology, which has been supported by the observational studies in the past decades (e.g., Im et al. 1997; Planck Collaboration et al. 2018).

2. SED Model

In this study, two methods based on SED models are used to select quasar candidates: broadband color cuts (Section 4.1) and BIC (Section 4.2). Below, we describe the quasar and M-dwarf star SED models we used for defining the quasar selection criteria.

2.1. Quasar SED Model

To establish selection criteria, we create various quasar SED models by adopting the composite quasar spectrum of Vanden Berk et al. (2001) as a base model and varying the continuum slope (αλ), and the equivalent width (EW) of Lyα and N v λ1240. The Selsing et al. (2016) quasar SED template may be useful in that it reduces the host-galaxy contribution near 5000 Å in rest frame, but they construct the SED template based on very luminous blue quasars with Mi ∼ −29 at z ∼ 1–2 which is inadequate to find faint quasars with M1450 ∼ −23.5. Regarding Brown et al. (2019), they make use of the SED spanning UV to radio based on only 41 AGNs. On the other hand, Vanden Berk et al. (2001) generate the SED template using over 2000 spectra of quasars at 0 < z < 5. Thus, we use Vanden Berk et al. (2001) composite spectra to generate diverse quasar SED models.

To generate realistic quasar models, we take the empirical distribution of αλ with $\langle {\alpha }_{\lambda }\rangle =-1.60$ and σ(αλ) = 1.0 (Mazzucchelli et al. 2017), where ${F}_{\lambda }\propto {\lambda }^{{\alpha }_{\lambda }}$. This relation is obtained from 15 quasars at z ≳ 6.5, but the mean continuum slope is comparable with −1.54 of Vanden Berk et al. (2001) SED and −1.70 of Selsing et al. (2016) SED. From the above values, the mean continuum slope of quasars might have consistent values across the luminosity and redshift range. So, we adjust a given continuum slope αλ by multiplying a factor of (λ/1000 Å)${}^{{\alpha }_{\lambda }-{\alpha }_{\lambda ,\mathrm{mean}}}$, where ${\alpha }_{\lambda ,\mathrm{mean}}=-1.6$.

Similarly, we use the EW distribution of Lyα + N v line from Bañados et al. (2016) which has $\langle $log EW (Å)$\rangle =1.542$ and σ(log EW (Å)) = 0.391. To calculate the EW, we model a local continuum slope at a wavelength range of 1216–1800 Å, excluding wavelength ranges corresponding to several emission lines. Then, we integrate the fluxes over the continuum model at a rest-frame wavelength range of 1160–1290 Å. To make a SED with a given EW, we replace the continuum-subtracted flux with EW0 = 92.91 (${f}_{\mathrm{EW},0}$) to the rescaled value (fEW).

To describe the IGM absorption, we adopt the IGM attenuation model of Inoue et al. (2014) hereafter, I14). Note that the I14 model’s correction to the rest-frame wavelength of 912–1216 Å is somewhat less than the popular IGM attenuation model of Madau (1995), in a sense that the observed g and r magnitudes of z ∼ 5 model SEDs become brighter.

We set the redshift range of 4.5–5.6 with a grid size of 0.01, the αλ range of −3.6 to 0.2 with a grid size of 0.2, the logarithmic value of EW range of 0.5–2.5 with a grid size of 0.2, and the absolute magnitude at 1450 Å in the rest frame (M1450) range of $-27\lt {M}_{1450}\lt -21$ with a grid size of 0.1 mag, resulting in ∼107 SEDs in total.

2.2. M-dwarf Model

When selecting quasars from colors and SED shapes, the main contaminants are low mass stars, in particular, M-dwarf stars. Therefore, we use the M-dwarf spectra from the BT-settl model (Allard et al. 2013). The spectra are known to match the observed M-dwarf spectra on a wide range of physical parameters. The spectra are available from a public archive Theoretical Spectra Web Server.12 Considering typical ranges of M-dwarfs’ parameters (Casagrande et al. 2008; Rajpurohit et al. 2013), we use a total of 1050 spectra covering the parameter space of Teff of 2000–4000 K with a 100 K step, log(g) of 3.5–6.0 with a 0.5 step, [M/H] of −4 to 0.5 with 0.2 ∼ 0.5 steps, and [α/M] of −0.2 to 0.4 with a 0.2 step. We add a normalization factor (fN) as a free parameter.

3. Imaging Data

3.1. Broadband Data

The ELAIS-N1 field is one of the deep extragalactic survey (DXS) fields observed by the ISO due to the low infrared background (Rowan-Robinson et al. 2004; Vaccari et al. 2005). The ELAIS-N1 field is centered at 16:11:00 + 55:00:00 (J2000), and a wealth of deep, wide-area, and multiwavelength data are available for this field. Our quasar selection requires deep optical and NIR imaging data. Therefore, HSC-SSP (DR1, Aihara et al. 2018a, 2018b) is chosen for optical data and the 3.8 m United Kingdom Infrared Telescope (UKIRT) Infrared Deep Sky Survey—Deep Extragalactic Survey (DXS; Lawrence et al. 2007) or the IMS are chosen for NIR data. The field coverages of these surveys are plotted in Figure 1 according to the instruments’ field of view and pointing. Figure 2 shows all the filters from these survey that are used in this study.

Figure 1. Refer to the following caption and surrounding text.

Figure 1. ELAIS-N1 area coverages of various surveys. The black rectangles show the DXS coverage. The IMS data provides additional J-band data in the outer part of the ELAIS-N1 field (filled blue rectangles). The HSC optical imaging data are plotted in pink circles.

Standard image High-resolution image
Figure 2. Refer to the following caption and surrounding text.

Figure 2. Filter transmission curves (including quantum efficiency) of the data we used in this work. For the sake of easy distinction, the filter curves of HSC-SSP ($g,\,r,\,i,\,z,\,\mathrm{NB}921$ and y) and DXS/IMS (J, and K) are drawn in background, while the medium-band filters are plotted as the solid lines. The dark gray line shows the representative SED of a quasar at z = 5.0.

Standard image High-resolution image

3.1.1. Optical—HSC

For the optical data, we retrieved a catalog from the HSC data archive system. HSC-SSP is a three-layered imaging survey with each layer survey covering 1400 deg2 (Wide), 27 deg 2 (Deep), 5 deg 2 (Ultra-deep) with the target image depths of r ∼ 26.1, 27.1, and 27.7 mags respectively for 5σ, point-source detection (Aihara et al. 2018a). The ELAIS-N1 field is one of the HSC-SSP deep fields of which DR1 i-band depth reaches ∼26.5 mag (Aihara et al. 2018b). Following the method in Matsuoka et al. (2018), which used a random source catalog from HSC-SSP archival system, we found the effective HSC-SSP DR1 coverage of ELAIS-N1 to be 6.75 deg2. Images of this field in S17A season were taken with multiple filters including $g,\,r,\,i,\,z,\,y$, and a narrow-band centered at 9210 Å (NB 921), reaching 5σ depths of 26.8, 26.6, 26.5, 25.6, 24.8, and 25.6 mag, respectively, for point sources (Aihara et al. 2018b). Considering that Lyα break is expected to be at ∼7200Å at z ∼ 5, we used i-band detected catalog from HSC-SSP database (Aihara et al. 2018a, 2018b). A quasar would appear as a point source under the seeing condition of HSC-SSP images. Therefore, we used the point-spread function (PSF) magnitudes in the HSC-SSP catalog for the optical photometry. Note that all magnitudes in the catalog are extinction-corrected values according to the dust map of the Schlegel et al. (1998, hereafter SFD).

3.1.2. NIR–DXS/IMS

DXS and IMS images were obtained with WFCAM at UKIRT. The combination of the DXS and IMS data covers almost the entire HSC ELAIS-N1 area in J-band. On the other hand, K-band image is only available in an area covered by DXS. The 5σ imaging depths are J ∼ 23.2 and K ∼ 22.7 mag for point sources (Hewett et al. 2006; M. Im et al. 2020, in preparation). Note that the J-band depths of the IMS and DXS are nearly identical and the total NIR survey area is 9.68 deg2 We performed the photometry of IMS and DXS data with our own pipeline that uses SExtractor (Bertin & Arnouts 1996). The auto-magnitude were used as total magnitudes. We found that the magnitudes of bright stars are consistent with those of stars in the Two Micron All Sky Survey (Skrutskie et al. 2006). Extinction correction values of these bands were calculated using the wavelength-dependent relation (Cardelli et al. 1989) between A(V) and $E(B-V)$ from SFD assuming RV = 3.1 and a representative wavelength equal to their effective wavelength.

3.1.3. Catalog Matching

The optical and NIR catalogs were matched with a matching radius of 1farcs0. When an object was not detected in J or K, we assigned J or K limiting magnitudes to the object. If the matching resulted in more than one matched pair, we chose an object in the NIR catalog that is closer to the optical source or has magnitudes similar to the optical data. Multiple matching case occurred in ∼2% of all the matched objects.

Our quasar search was performed on the overlap area of the optical and the NIR data, and we calculated this quasar search area in the following way. The NIR sources were matched with the nonsaturated (i > 19) and bright (i < 21) optical sources in the HSC catalog with a matching distance of 1farcs0. The bright magnitude cut at i < 21 ensures that each i-band source must have a NIR counterpart in the J-band catalog in the area where the optical and NIR images overlap. Then, the fraction of the optical sources with NIR counterparts should be identical to the fraction of the optical image area that overlaps with the NIR image. The 96.4% of the optical sources were found to have NIR counterparts, thus giving 6.51 deg2 as the quasar search area.

3.2. Medium-band Data

Multiple steps were taken to select z ∼ 5 quasar candidates, starting with a selection using the broadband data (Section 4.1). Then, we refined the candidate selection using the SED-fitting method (Section 4.2). To improve the SED fitting, we obtained additional medium-band data of the broadband-selected candidates. The medium-bands that we used have widths of ∼500 Å centered at 675, 725, 775, and 825 nm (hereafter, m675, m725, m775, and m825). Compared with the broad bands with widths of ∼1500 Å, the medium-bands can improve the spectral sampling by a factor of 3, and hence the better quasar selection can be achieved. The medium-band observations were conducted using the SED camera for QUasars in EArly uNiverse (SQUEAN) on the 2.1 m Otto-Struve telescope at Mcdonald Observatory (Park et al. 2012; Kim et al. 2016) and Seoul National University 4K ×4K Camera (SNUCAM) on the 1.5 m telescope at Maidanak Observatory (Im et al. 2010). The data were taken with the seeing condition better than 1farcs5 with the on-source exposure time shorter than two hours, or two hours at maximum. The targets are summarized in Table 1. In total, 10 out of 13 candidates have been observed (see Section 4).

Table 1.  Medium-band Photometry of the Broadband-selected Quasar Candidates

ID R.A. Decl. m675 m725 m775 m825
  (J2000) (J2000) (mag) (mag) (mag) (mag)
IMS J160306+541928 16:03:06.02 54:19:28.46 >24.17 22.98 ± 0.18
IMS J160517+554002 16:05:17.80 55:40:02.00 >23.55 21.96 ± 0.12
IMS J160552+555340 16:05:52.07 55:53:40.61 >23.06 20.87 ± 0.10
IMS J160622+540056 16:06:22.50 54:00:56.64
IMS J160732+544750 16:07:32.21 54:47:50.86 >24.19 22.25 ± 0.10
IMS J160748+541157 16:07:48.15 54:11:57.40 >24.29 22.90 ± 0.26
IMS J160914+554511 16:09:14.68 55:45:11.83 >23.23 21.94 ± 0.18 22.65 ± 0.14
IMS J161248+550927 16:12:48.98 55:09:27.54 >23.48 22.79 ± 0.17a 23.00 ± 0.20
IMS J161341+542146 16:13:41.15 54:21:46.65
IMS J161343+542131 16:13:43.43 54:21:31.29 21.44 ± 0.10 21.00 ± 0.06 20.81 ± 0.11
IMS J161636+535545 16:16:36.32 53:55:45.01
IMS J161827+551748 16:18:27.29 55:17:48.47 >23.09 21.31 ± 0.21 21.55 ± 0.21 21.34 ± 0.20
IMS J161903+545638 16:19:03.74 54:56:38.92 22.18 ± 0.21 22.50 ± 0.18

Notes.

aObserved in Maidanak Observatory.

Download table as:  ASCIITypeset image

All the medium-band images were preprocessed including bias/dark subtraction and flat-fielding. In Maidanak data, a fringe pattern appeared in the images, and they were corrected by subtracting the master fringe frame made from dithered observed images, following Jeon et al. (2010). Astrometry solutions were obtained for all the preprocessed images using the Astrometry.net software (Lang et al. 2010).

Then, the images were background-subtracted and flux-rescaled based on the zero-point of each image to fill the gap of various circumstances of observational conditions at each observational period. Outlier images determined by seeing and depth distribution of images were excluded in combining procedure.

The zero-point of each medium-band image was derived by following the prescription described in Jeon et al. (2016) and Choi & Im (2017). The procedure first identified bright but nonsaturated stars near targets with existing multiwavelength photometry data. In this case, we used r, i and z magnitudes of 5 ∼ 10 stars from the Seventh Data Release of the Sloan Digital Sky Survey (SDSS-DR7; Abazajian et al. 2009). Then, for each star, the best-fit stellar template was identified among 175 stellar templates of Gunn & Stryker (1983). The medium-band photometry of the star was derived using the best-fit template, and the derived magnitude was used to calculate the zero-point. The magnitudes were measured using SExtractor with an aperture radius of 1.6× seeing FWHM. The standard deviation of the zero-points from the stars in the field was taken as the zero-point error, found to be ≲0.05 magnitudes.

We used the zero-points to derive photometry of the quasar candidates (Table 1). The same aperture was used for the standard star, and thus the aperture correction was taken into account inherently. The photometry error from SExtractor and the zero-point errors were summed in quadrature. The derived magnitudes were then corrected for the Galactic extinction by applying the extinction law of Cardelli et al. (1989), the dust map of SFD, and assuming RV = 3.1.

4. High-redshift Quasar Selection

The selection of quasar candidates was done in two steps, first using the broadband color selection, and then refining the broadband-selected sample using the BIC method. This section describes each step of the quasar selection.

4.1. Broadband Color-selection Criteria

The quasar candidates were first selected using broadband colors. We employed broadband color cuts in riz and riJ color–color diagrams (CCD), similarly to that described in McGreer et al. (2013) and Kim et al. (2019), with the color cuts adjusted for the HSC magnitude system.

First, we selected point sources by using the difference between CModel magnitude (${m}_{\mathrm{CModel}}$) and PSF magnitudes (${m}_{\mathrm{psf}}$) (Matsuoka et al. 2018). Then, we applied i-band magnitude cut brighter than 23.5 considering NIR depths. The g-band flux must be fainter or nondetected due to the IGM attenuation. Thus, we set $g-r\gt 1.8$ or signal-to-noise ratio (S/N) in the g band < 3.0 to select objects with weak or nondetection in the g band.

As shown in Figure 3, the stellar loci can be described with a straight line with a dispersion. Therefore, as a next step, we defined slopes of the stellar loci on both riz and riJ CCDs by fitting a linear line. After getting the slopes of these loci, the distributions of point sources were examined as a function of the distance diagonal from the line. We took percentile cuts of 99.5% (riz) and 95% (riJ) away from the stellar loci to exclude stars from the selection. These are the diagonal pink dashed lines in Figure 3.

Figure 3. Refer to the following caption and surrounding text.

Figure 3. Color–Color diagrams of point sources at i < 23.5 mag. The black shades show the number density of point sources. The quasar selection boxes are marked with the pink dashed line, and the selected quasar candidates are marked with orange diamonds. The blue line indicates the track of a quasar at different redshifts with intervals of 0.1. The green points show the colors of spectroscopically confirmed M-dwarf stars from the Fourteenth Data Release of the SDSS (Abolfathi et al. 2018). The M-dwarf points follow the stellar locus (the black solid line). On both panels, a and b values indicate the slope and intercept of the linear relation fit to the point sources.

Standard image High-resolution image

To further reduce objects that could contaminate the quasar candidate selection, we introduced r − i color cut (${{ri}}_{\mathrm{cut}}$) and i − z color cut (izcut). We defined the ricut and izcut by choosing the condition that include the maximal number of quasar models while reducing the total riz color-space area of the selection box. This was done by choosing the conditions where the following number, f(ricut, izcut), is maximal,

Equation (1)

Here, $R({{ri}}_{\mathrm{cut}},{{iz}}_{\mathrm{cut}})$ is the completeness of the quasar selection as defined as the ratio of model quasar SEDs picked up by the selection region to the total number of model quasar SEDs, and $A({{ri}}_{\mathrm{cut}},{{iz}}_{\mathrm{cut}})$ is the area defined by the selection region with additional color limits of ricut < 3 and izcut > 0. The smaller the A value is, the less likely we get contaminants. Also, we optimized the selection criteria by satisfying the prerequisite that one or more quasar model SEDs could be selected with a lower redshift limit at z = 4.5.

Figure 4 shows f(ricut, izcut) for various parameters, and we adopted ricut = 1.235 and izcut = 0.673 that maximize f(ricut, izcut).13 Below, we summarize our broadband color-selection criteria:

  • 1.  
    Point sources with ipsfiCModel < 0.15
  • 2.  
    $i\lt 23.5$
  • 3.  
    g − r > 1.8 or ${\rm{S}}/{\rm{N}}(g)\lt 3.0$
  • 4.  
    $r-i\gt 1.235$
  • 5.  
    $i-z\lt 0.391(r-i)-0.220$
  • 6.  
    $i-z\lt 0.673$
  • 7.  
    $i-J\lt 0.766(r-i)-0.525$.

Figure 4. Refer to the following caption and surrounding text.

Figure 4. Ratio of the completeness to the area of the color-selection region (Equation (1)) as a function of ${{ri}}_{\mathrm{cut}}$ and ${{iz}}_{\mathrm{cut}}$. The dotted line indicates the maximum of the ratio where we define the optimal ${{ri}}_{\mathrm{cut}}$ and ${{iz}}_{\mathrm{cut}}$ values, taking into consideration that one or more quasar model SEDs could be selected with a lower redshift cut at z = 4.5.

Standard image High-resolution image

Quasar candidates that satisfy the color selection were visually inspected and rejected if the photometry is found to be spurious during the visual inspection (e.g., scattered light contamination). Only one source was found to be spurious and rejected during the visual inspection. Finally, we identified 13 quasar candidates, which are shown in Figure 3. Table 2 summarizes the HSC and NIR photometry of the quasar candidates.

Table 2.  HSC and NIR Photometry of Our z ∼ 5 Quasar Candidates

ID g r i NB 921 z J K
  (mag) (mag) (mag) (mag) (mag) (mag)
IMS J160306+541928 >26.8 24.88 ± 0.17 23.13 ± 0.03 22.41 ± 0.02 22.71 ± 0.21 22.73 ± 0.22 22.99 ± 0.30
IMS J160517+554002 >26.8 24.31 ± 0.05 22.57 ± 0.01 22.46 ± 0.01 22.31 ± 0.02 21.92 ± 0.12 21.78 ± 0.14
IMS J160552+555340 >26.8 23.31 ± 0.02 21.45 ± 0.00 21.08 ± 0.01 21.00 ± 0.01 21.27 ± 0.07 20.54 ± 0.06
IMS J160622+540056 >26.8 26.19 ± 0.21 23.49 ± 0.01 23.40 ± 0.04 23.24 ± 0.04 22.70 ± 0.17 >22.700
IMS J160732+544750 >26.8 24.68 ± 0.05 23.44 ± 0.01 23.20 ± 0.02 23.21 ± 0.03 >23.158 >22.700
IMS J160748+541157 >26.8 24.45 ± 0.04 22.90 ± 0.01 22.65 ± 0.02 22.56 ± 0.02 22.49 ± 0.15 21.64 ± 0.15
IMS J160914+554511 >26.8 24.06 ± 0.03 22.55 ± 0.01 22.42 ± 0.01 22.31 ± 0.01 22.82 ± 0.18 22.85 ± 0.26
IMS J161248+550927 >26.8 24.19 ± 0.03 22.78 ± 0.01 22.62 ± 0.02 22.52 ± 0.02 22.67 ± 0.18 22.36 ± 0.20
IMS J161341+542146 >26.8 24.94 ± 0.05 23.49 ± 0.01 23.29 ± 0.02 23.29 ± 0.03 23.31 ± 0.25 0.00 ± 0.00
IMS J161343+542131 25.24 ± 0.03 22.94 ± 0.01 20.58 ± 0.00 19.80 ± 0.00 19.92 ± 0.00 19.40 ± 0.02 19.61 ± 0.03
IMS J161636+535545 >26.8 24.86 ± 0.08 23.23 ± 0.02 22.70 ± 0.03 22.85 ± 0.05 >23.158 >22.700
IMS J161827+551748 25.39 ± 0.04 22.78 ± 0.01 21.11 ± 0.00 20.87 ± 0.00 20.83 ± 0.00 20.78 ± 0.05 20.14 ± 0.04
IMS J161903+545638 25.78 ± 0.05 23.95 ± 0.02 22.39 ± 0.01 21.98 ± 0.01 22.06 ± 0.01 21.78 ± 0.10 22.11 ± 0.20

Note. According to Table 2 in Aihara et al. (2018b), the 5σ depth in the g band is 26.8 mag for a point source.

Download table as:  ASCIITypeset image

4.2. BIC Selection and SED Fitting

4.2.1. SED Fitting

To refine the quasar selection, we first fitted the SED of the quasar candidates to various quasar and M-dwarf models. For the SED fitting, we added the medium-band data to sample the sharp break in the quasar SED located at the wavelength of redshifted Lyα. We then used the information from the fitting to apply the BIC to refine the candidate selection.

The SED fitting provided the best-fit parameters and the goodness-of-fit in terms of the reduced chi-squares, ${\chi }_{\mathrm{red}}^{2}$, for each model. For quasar models, we obtained the photometric redshift zphot, the continuum slope αλ, the EW of the Lyα + N v line, and the UV absolute magnitude at 1450 Å, M1450. For M-dwarf models, we obtained Teff, log(g), [M/H], [α/M], and the normalization parameter. The SED fitting followed the procedure as described in Kim et al. (2019), and it takes into account the upper limits as described in by Sawicki (2012). The results of the SED fitting are shown in Figure 5. In Table 3, we present the derived parameters such as zphot, M1450, and ${\chi }_{\mathrm{red}}^{2}$ for both the best-fit quasar and M-dwarf models. Figure 5 shows that many of the candidates are better fitted to the quasar models (e.g., IMS J160914+554511), but some are not (e.g., IMS J161343+542131) and some are ambiguous (e.g., IMS J161341+542146). To better differentiate quasars from M-dwarfs, we apply the BIC as described in the next subsection.

Figure 5. Refer to the following caption and surrounding text.
Standard image High-resolution image
Figure 5. Refer to the following caption and surrounding text.

Figure 5. SED-fitting results of the broadband-selected quasar candidates. The broadband data are expressed as blue points, while the medium-band data are plotted with red points. The horizontal bar of each point indicates the width of the filter. For nondetections, pink circles and arrows present the 5σ limiting magnitudes. The best dwarf star model and the best quasar model obtained through the SED fitting are shown as the green and gold lines, respectively. The best-fit results are also indicated, where fN are in units of 10−12.

Standard image High-resolution image

Table 3.  SED-fitting Results

ID zphot zspec ${\rm{\Delta }}z/(1+{z}_{\mathrm{spec}})$ ${M}_{1450,\mathrm{phot}}$ ${M}_{1450,\mathrm{spec}}$ ${\chi }_{\mathrm{red},d}^{2}$ ${\chi }_{\mathrm{red},q}^{2}$ ΔBIC Class
IMS J160306+541928 ${4.92}_{-0.05}^{+0.03}$ $-{23.42}_{-0.18}^{+0.18}$ 9.2 2.9 24.7 Quasar
IMS J160517+554002 ${5.16}_{-0.03}^{+0.01}$ ${5.211}_{-0.001}^{+0.001}$ 0.008 $-{24.07}_{-0.14}^{+0.19}$ $-{24.454}_{-0.017}^{+0.017}$ 15.3 4.1 43.0 Quasar
IMS J160552+555340 ${5.17}_{-0.03}^{+0.06}$ ${5.409}_{-0.001}^{+0.001}$ 0.037 $-{25.31}_{-0.19}^{+0.21}$ $-{25.880}_{-0.010}^{+0.010}$ 18.3 2.7 46.1 Quasar
IMS J160622+540056 ${5.19}_{-0.05}^{+0.03}$ $-{23.00}_{-0.20}^{+0.20}$ 12.9 4.9 13.3 Quasar
IMS J160732+544750 ${4.87}_{-0.03}^{+0.06}$ $-{22.92}_{-0.08}^{+0.18}$ 24.2 8.0 58.9 Quasar
IMS J160748+541157 ${4.71}_{-0.04}^{+0.09}$ $-{23.37}_{-0.13}^{+0.17}$ 12.5 2.8 38.3 Quasar
IMS J160914+554511 ${4.72}_{-0.04}^{+0.04}$ ${4.814}_{-0.002}^{+0.002}$ 0.016 $-{23.61}_{-0.13}^{+0.11}$ $-{23.807}_{-0.017}^{+0.017}$ 15.5 2.7 63.9 Quasar
IMS J161248+550927 ${4.67}_{-0.05}^{+0.05}$ $-{23.42}_{-0.18}^{+0.16}$ 8.2 2.0 31.2 Quasar
IMS J161341+542146 ${4.68}_{-0.00}^{+0.01}$ $-{22.76}_{-0.14}^{+0.16}$ 5.7 3.0 4.4 Nonquasar
IMS J161343+542131 ${5.08}_{-0.04}^{+0.04}$ $-{26.08}_{-0.22}^{+0.28}$ 21.4 30.0 −70.8 Nonquasar
IMS J161636+535545 ${5.35}_{-0.07}^{+0.01}$ $-{23.67}_{-0.13}^{+0.17}$ 12.0 4.9 11.4 Quasar
IMS J161827+551748 ${4.73}_{-0.01}^{+0.09}$ $-{25.12}_{-0.12}^{+0.22}$ 13.4 3.4 59.1 Quasar
IMS J161903+545638 ${4.71}_{-0.02}^{+0.09}$ $-{23.99}_{-0.15}^{+0.22}$ 4.1 6.1 −11.8 Nonquasar

Note. This table show the result of SED fitting. The errors of ${z}_{\mathrm{phot}}$ and ${M}_{1450,\mathrm{phot}}$ are derived from the 68% confidence intervals with other parameters fixed, while ${z}_{\mathrm{spec}}$ and ${M}_{1450,\mathrm{spec}}$ are estimated from MCMC simulation with other parameters free. The BIC values are calculated using all photometric data.

Download table as:  ASCIITypeset image

4.2.2. BIC Selection of Quasars

To improve the quasar selection, we used the BIC for each best-fit model rather than ${\chi }_{\mathrm{red}}^{2}$. The ratio of the ${\chi }_{\mathrm{red}}^{2}$ of best-fit quasar model to ${\chi }_{\mathrm{red}}^{2}$ of best-fit M-dwarf star model can be used to prioritize the candidates (Mazzucchelli et al. 2017). However, the ratio cannot consider the effect of difference in the number of fitting parameters between two models. On the contrary, the BIC can take the number of free parameters as well as the number of data into account, which can be defined as

Equation (2)

where n is the number of photometric data, k is the number of free parameters in the tested model, and Lmax is the maximum likelihood value of the model. Therefore, the BIC can be an effective tool when comparing models with the different k.

The difference in the values of BIC of two models can be a quantitative measure to tell between the better-fitting model and given as

Equation (3)

where “d” and “q” in subscripts indicate the best-fitting M-dwarf star and quasar models, respectively. If the difference larger than 10 in BIC commonly means “decisive” evidence that supports the quasar model is better than the M-dwarf star model (Liddle 2007).

We computed the ΔBIC of the best-fit quasar and the M-dwarf models with and without medium-band data (see Figure 6 and Table 3). Among 13 candidates from the broadband color selection, we identified 10 to be promising high-redshift quasar candidates, and two to be M-dwarf stars. The medium-band data are critical for distinguishing the nature of the candidates, as one of the three spectroscopically identified quasars in Section 5 are deemed ambiguous unless the additional observation data were used.

Figure 6. Refer to the following caption and surrounding text.

Figure 6. ΔBIC graph showing if a broadband-selected quasar is more like stars or quasars. The pink shaded area indicates the region with ΔBIC < 10 where it is difficult to choose the better model with confidence. The quasars should lie on the right side of the shaded area, while M-dwarf should lie on the left. The filled and open triangles are ΔBIC of the candidates with and without medium-band data, respectively. The spectroscopically confirmed quasars are marked with the open pink circles. The addition of medium-band data clearly helps distinguish which objects are quasars, with the same classification uncertain otherwise.

Standard image High-resolution image

5. Spectroscopic Confirmation of Quasars

5.1. Observation

Three quasar candidates were observed on 2018 July 10 with the Hectospec instrument on MMT (Fabricant et al. 2005). These candidates were chosen based on the ΔBIC values in Table 3 and the feasibility of simultaneous observations with other objects. The Hectospec can place fibers on ∼300 targets within a one-degree diameter field of view, and the quasar targets were included as a part of another program studying galaxy clusters in the ELAIS-N1 field. For the observation, we used the 270 lines mm−1 grating that covers the wavelength range of 3650–9200 Å at a spectral resolution of ∼1000. On-source exposure time was 240 minutes per quasar. The spectra were reduced with the HSRED pipeline for the basic reduction (bias and flat field), the wavelength and flux calibrations, and the sky subtraction. The standard star and the sky fiber data were used for the flux calibration and the sky subtraction, respectively.

5.2. Spectroscopic Identification

The Hectospec observation reveals all the candidates, IMS J160517+554002, IMS J160552+555340, and IMS J160914+554511, are quasars at z ∼ 5. Figure 7 shows their spectra. The strong break at Lyα can be seen in all spectra, suggesting that they are all indeed at z ∼ 5. Their spectroscopic redshifts and M1450 are derived with the SED fitting using 20,000 chains of Markov Chain Monte Carlo (MCMC) simulation to these spectra. Both zspec and zspec-based M1450 of the three quasars are given in Table 3. Their M1450 have ranges ≳−26. The 3/3 confirmation rate gives us the confidence for the BIC selection of z ∼ 5 quasars including medium-band observations.

Figure 7. Refer to the following caption and surrounding text.

Figure 7. MMT Hectospec spectra of the three spectroscopically confirmed quasars (black lines) with 1σ flux uncertainties (gray lines). The pink lines show the best-fit results of quasar model SED fitting. The sky signal and its errors are presented in the bottom panel. The sharp breaks at the redshifted Lyα indicate that they are all high-redshift quasars. The zspec and the M1450 values from the SED fitting are shown in the figure.

Standard image High-resolution image

6. LF and Implication for Reionization

From the ten BIC-selected candidates at a redshift range 4.6 ≲ z ≲ 5.4 including three confirmed quasars, we derive the quasar LF at M1450 ≳ −26 mag. To do so, we use the updated 1/Va method (Page & Carrera 2000; Im et al. 2002). For a given interval Δz and a given magnitude interval ΔM1450, Va, the effective volume covered by N number of quasars belonging in the bin can be written as

Equation (4)

where Fcomp represents the completeness from 0.0 to 1.0, ${dV}/{dz}$ is the cosmological volume element taking into account of the survey area, zmin is the lowest redshift of the redshift interval, and ${z}_{\max }({M}_{1450})$ is the maximum redshift where an object with M1450 is within the flux limit of the survey.

The completeness function Fcomp is the fraction of quasar models enable to pass the selection criteria respect to all the simulated quasar models at given ranges of parameters. The completeness function has two variables of redshift and magnitude, which are binned with a bin size of 0.05 and 0.2, respectively. Figure 8 shows the completeness function values across 4.5 < z < 5.5 and −27 < M1450 < −23. It supports the relevance of our survey searching for quasar candidates with z ∼ 5 and M1450 ≲ −23.

Figure 8. Refer to the following caption and surrounding text.

Figure 8. Completeness function (Fcomp) of the broadband quasar selection as a function of z vs. M1450. This figure indicates that quasars at z = 4.6 ≲ z ≲ 5.4 in the magnitude range of −27 ≲ M1450 ≲ −23 will be efficiently selected. The gold triangle shows BIC-selected sample. The spectroscopically confirmed quasars in the sample are marked with the open pink circles.

Standard image High-resolution image

Note that we assumed that the photometric completeness is 100% considering that the point-source detection completeness of a similar depth HSC data is 100% at i < 23.5 (Matsuoka et al. 2018). Similarly, we did not take into account, in the Fcomp calculation, the point-source completeness and the incompleteness in the sample due to source confusion. The point-source completeness of a bit shallower HSC data is ≳95%, and the source confusion affects the Fcomp at a similar level (Matsuoka et al. 2018). They are much smaller than the errors of the LF points which amounts to 40% or larger (Table 4), and thus are neglected in the LF calculation.

Table 4.  Quasar LF

M1450 ΔM1450 Φ σΦ N Ncor
    (Mpc−3 mag−1)      
−25.4 1.4 3.64e−08 2.57e−08 2.0 4.3
−24.1 1.2 6.52e−08 3.78e−08 3.0 3.9
−23.0 1.0 2.79e−07 1.25e−07 5.0 18.6
  a     b  
  ${0.43}_{-0.26}^{+0.29}$     ${3.3}_{-6.3}^{+6.9}$  

Download table as:  ASCIITypeset image

Then, we calculate the binned number density and its uncertainty, $\delta N$ according to the following equations where N is the number of the observed objects in the magnitude bin, and δΦ is directly related to the Poisson noise of N:

Equation (5)

Due to the small number of the sample, we adopt a simple power-law function for the LF as ${\rm{\Phi }}={{\rm{\Phi }}}_{* }\,{L}^{\alpha +1}$, or as

Equation (6)

Here, a is related to the slope of the LF α, as α = −2.5 × a −1, and b is related to the normalization of Φ. We estimate the parameters in Φmodel by minimizing $S=-2\mathrm{ln}L$, described as

Equation (7)

The first term in S is applicable to the newly discovered quasars and the BIC-selected candidates. The second term in S provides the normalization term, and summed over the whole redshift and the magnitude ranges of the sample. Note that, in this parameter estimation, we ignore the redshift evolution in the quasar number density, as the number density evolution is negligible over this redshift range of 4.6 ≲ z ≲ 5.4.

The binned LF and the fitted parameter values are presented in Table 4. Also, Figure 9 shows the binned LF and the fitted LF, along with several recent LF values from the literature. The fitted LF slope at M1450 ≳ −26 is $-{2.08}_{-0.72}^{+0.65}$, which is in line with recently reported faint end slope of $-{1.97}_{-0.09}^{+0.09}$ (M18), but is slightly steeper than the results from some other studies (Akiyama et al. 2018; Matsuoka et al. 2018; G19).

Figure 9. Refer to the following caption and surrounding text.

Figure 9. Quasar LF at z ∼ 5 from this work (pink pentagons and solid line). The x-axis and y-axis error bars of the number density from this work represent the luminosity bin width for Va calculation and the statistical error of the number density estimate, respectively. The orange points indicate the quasar LF from Yang et al. (2016). Also plotted are the fitted LF functions from G15 (the blue dotted–dashed line), McGreer et al. (2018) the green dashed line), Parsa et al. (2018) the black dashed line) and G19 (sky-blue line with shade). Note that LFs of G19 at z ∼ 4.5 and 5.6 are shown as the filled and open sky-blue circles, respectively, and we show both functions with a shade. The purple plus markers are number densities of Boutsia et al. (2018) rescaled to z = 5.

Standard image High-resolution image

More importantly, the LF at $-26\lesssim {M}_{1450}\lesssim -23$ of our study gives a much lower density of quasars (∼6 times) at z ∼ 5 than the earlier study of G15. Our LF is similar to that of M18, and is not too far from the more recent results from Boutsia et al. (2018) rescaled to z = 5 and G19. To estimate the UV emissivity contribution of the quasars within the M1450 range of −29 to −18, we adopt the bright-end slope of the best-fit quasar LF in Yang et al. (2016) and expressing the LF with single power-law function. We then combine our LF with the LF of Yang et al. (2016) and merge them at the intersect of the two LFs. By assuming the escape fraction of 1 (Cristiani et al. 2016; Guaita et al. 2016; Grazian et al. 2018; Romano et al. 2019) and following the Equation (3) in G15, we find ${\epsilon }_{912}$ 14 ∼ 1. Comparing with epsilon912 ∼ 0.8 of McGreer et al. (2018), epsilon912 ∼ 1.2 of Parsa et al. (2018) and epsilon912 ∼ 3.8 of G19, our result prefers a minor contribution of a quasar to keep the IGM ionized at z ∼ 5. The photoionization rate from our result is ∼0.03 × 10−12 s−1, which is ≳10 times smaller than the photoionization rate of UV background at z ≃ 5 (Bolton & Haehnelt 2007; Calverley et al. 2011; Wyithe & Bolton 2011). Therefore, optically selected quasars are insufficient to fully explain the IGM ionization at z ∼ 5.

7. Summary

In this study, we searched for faint quasars at z ∼ 5 with −26 ≲ M1450 ≲ −23 (i < 23.5 mag) in the ELAIS-N1 field over an area of 6.51 deg2. Among optical quasars, quasars in this magnitude range are thought to contribute the most to the IGM-ionizing photons, therefore, securing a sample of such quasars is important for understanding the cosmic IGM ionization history.

Using the devised broadband color-selection criteria to optimize for the HSC photometry, we identified 13 z ∼ 5 quasar candidates with the optical/NIR imaging data from HSC-SSP, IMS and DXS survey. We then refined the candidate selection by adding medium-band data and performing additional cuts based on BIC. As a result, we obtained the refined sample of 10 candidates. The spectroscopic observation was carried out for three of the BIC-selected candidates using MMT Hectospec, and all of these candidates were identified as quasars at z ∼ 5. In addition, the redshift and M1450 of the other BIC-selected candidates were derived through the SED fitting.

Using the three newly confirmed quasars and the seven BIC-selected quasar candidates with zphot, we constructed the z ∼ 5 quasar LF at −26 ≲ M1450 ≲ −23, with a single power-law function with a slope of $-{2.08}_{-0.72}^{+0.65}$. We also found that the number density of the faint z ∼ 5 quasars is ∼10−7 Mpc−3 mag−1. Combined with the bright-end LF from the literature, we showed that the UV emissivity of quasars are insufficient to explain the IGM ionization at z ∼ 5, in line with the earlier results for z ≳ 6 quasars.

Our method of finding high-redshift quasars via medium-band data and BIC can improve the efficiency of quasar survey. The broadband color selection is an efficient way to identify high-redshift quasars, but three among the 10 BIC-selected candidates (30%) has been unambiguously selected as quasars from the broadband colors only. These candidates needs the additional medium-band observation to make sure they are highly promising candidates. The medium-band observations can be done with small to midsized telescope for which telescope time is more readily available than larger class telescopes, and we expect that future surveys of faint quasars would benefit from medium-band observations when the broadband selection is uncertain.

We appreciate the anonymous referee’s constructive comments. We thank Yoonsoo Bach for his input on the statistical analysis. M.H. acknowledges the support from Global PH.D Fellowship Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (NRF-2013H1A2A1033110) This work was supported by the National Research Foundation of Korea (NRF), grant No. 2017R1A3A3001362, No. 2017R1A6A3A04005158 and No. NRF-2019R1C1C1002796, funded by the Korea government (MSIP). This work was supported by the K-GMT Science Program (PID:2018B-UAO-G7) of the Korea Astronomy and Space Science Institute (KASI). Based [in part] on data collected at the Subaru Telescope and retrieved from the HSC data archive system, which is operated by Subaru Telescope and Astronomy Data Center at National Astronomical Observatory of Japan. The Hyper Suprime-Cam (HSC) collaboration includes the astronomical communities of Japan and Taiwan, and Princeton University. The HSC instrumentation and software were developed by the National Astronomical Observatory of Japan (NAOJ), the Kavli Institute for the Physics and Mathematics of the Universe (Kavli IPMU), the University of Tokyo, the High Energy Accelerator Research Organization (KEK), the Academia Sinica Institute for Astronomy and Astrophysics in Taiwan (ASIAA), and Princeton University. Funding was contributed by the FIRST program from Japanese Cabinet Office, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), the Japan Society for the Promotion of Science (JSPS), Japan Science and Technology Agency (JST), the Toray Science Foundation, NAOJ, Kavli IPMU, KEK, ASIAA, and Princeton University. This paper makes use of software developed for the Large Synoptic Survey Telescope. We thank the LSST Project for making their code available as free software at http://dm.lsst.org. The Pan-STARRS1 Surveys (PS1) have been made possible through contributions of the Institute for Astronomy, the University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society and its participating institutes, the Max Planck Institute for Astronomy, Heidelberg and the Max Planck Institute for Extraterrestrial Physics, Garching, The Johns Hopkins University, Durham University, the University of Edinburgh, Queen's University Belfast, the Harvard-Smithsonian Center for Astrophysics, the Las Cumbres Observatory Global Telescope Network Incorporated, the National Central University of Taiwan, the Space Telescope Science Institute, the National Aeronautics and Space Administration under grant No. NNX08AR22G issued through the Planetary Science Division of the NASA Science Mission Directorate, the National Science Foundation under grant No. AST-1238877, the University of Maryland, and Eotvos Lorand University (ELTE) and the Los Alamos National Laboratory. The UKIRT is owned by the University of Hawaii (UH) and operated by the UH Institute for Astronomy; operations are enabled through the cooperation of the East Asian Observatory. When (some of) the data reported here were acquired, UKIRT was operated by the Joint Astronomy Centre on behalf of the Science and Technology Facilities Council of the U.K. Observations reported here were obtained at the MMT Observatory, a joint facility of the Smithsonian Institution and the University of Arizona. This paper includes data taken at the McDonald Observatory of the University of Texas at Austin.

Facilities: UKIRT (WFCAM), Struve (SQUEAN), MMT (Hectospec), Sloan.

Software: SExtractor (Bertin & Arnouts 1996), Astrometry.net (Lang et al. 2010).

Footnotes

Please wait… references are loading.
10.3847/1538-4357/ab7bde