Linking gas and galaxies at high redshift: MUSE surveys
the environments of six damped Lyα galaxies at z ≈ 3
Ruari Mackenzie,
1,2?
Michele Fumagalli,
2,3
, Tom Theuns,
3
David J. Hatton,
2
Thibault Garel,
4,5
, Sebastiano Cantalupo
1
, Lise Christensen
6
, Johan P. U. Fynbo
6
,
Nissim Kanekar
7
, Palle Møller
8
, John O’Meara
9,10
, J. Xavier Prochaska
11,12
,
Marc Rafelski
13
, Tom Shanks
2
, James Trayford
14
1Department of Physics, ETH Zurich, Wolfgang-Pauli-Strasse 27, 8093 Zurich, Switzerland 2Centre for Extragalactic Astronomy, Durham University, South Road, Durham, DH1 3LE, UK 3Institute for Computational Cosmology, Durham University, South Road, Durham, DH1 3LE, UK
4Univ Lyon, Univ Lyon1, Ens de Lyon, CNRS, Centre de Recherche Astrophysique de Lyon UMR5574, F-69230, Saint-Genis-Laval, France 5Observatoire de Gen`eve, Universit´e de Gen`eve, 51 Ch. des Maillettes, 1290 Versoix, Switzerland
6Dark Cosmology Centre, Niels Bohr Institute, University of Copenhagen, Juliane Maries Vej 30, DK-2100 Copenhagen, Denmark 7National Centre for Radio Astrophysics, TIFR, Post Bag 3, Ganeshkhind, Pune 411007, India
8European Southern Observatory, Karl-Schwarzschild-Strasse 2, 85748, Garching, Germany 9Department of Physics, St. Michael’s College, Colchester, VT 05439, USA
10W. M. Keck Observatory, Waimea, HI 96743, USA
11Department of Astronomy and Astrophysics, University of California, 1156 High Street, Santa Cruz, CA 95064 USA 12University of California Observatories, Lick Observatory 1156 High Street, Santa Cruz, CA 95064 USA
13Space Telescope Science Institute, Baltimore, MD, USA
14Leiden Observatory, Leiden University, PO Box 9513, NL-2300 RA Leiden, the Netherlands
17 April 2019
ABSTRACT
We present results from a survey of galaxies in the fields of six z ≥ 3 Damped Lyman α systems (DLAs) using the Multi Unit Spectroscopic Explorer (MUSE) at the Very Large Telescope (VLT). We report a high detection rate of up to ≈ 80% of galaxies within 1000 km/s from DLAs and with impact parameters between 25 and 280 kpc. In particular, we discovered 5 high-confidence Lyα emitters associated with three DLAs, plus up to 9 additional detections across five of the six fields. The majority of the de-tections are at relatively large impact parameters (> 50 kpc) with two dede-tections being plausible host galaxies. Among our detections, we report four galaxies associated with the most metal-poor DLA in our sample (Z/Z = −2.33 ± 0.22), which trace an
over-dense structure resembling a filament. By comparing our detections with predictions from the Evolution and Assembly of GaLaxies and their Environments (EAGLE) cos-mological simulations and a semi-analytic model designed to reproduce the observed bias of DLAs at z > 2, we conclude that our observations are consistent with a scenario in which a significant fraction of DLAs trace the neutral regions within halos with a characteristic mass of Mh≈ 1011− 1012M
, in agreement with the inference made from
the large-scale clustering of DLAs. We finally show how larger surveys targeting ≈ 25 absorbers have the potential of constraining the characteristic masses of halos hosting high-redshift DLAs with sufficient accuracy to discriminate between different models. Key words: galaxies: evolution – galaxies: formation – galaxies: haloes – galaxies: high-redshift – quasars: absorption lines
? E-mail: mruari@phys.ethz.ch
1 INTRODUCTION
Since the discovery of damped Lyα absorbers (DLAs) in the spectra of quasars in the 1970s (Beaver et al. 1972; Carswell
et al. 1975), significant efforts have been made to identify the properties of the galaxy population that gives rise to this class of absorption line systems. The interest in con-necting DLAs to galaxies stems from the fact that these absorbers, defined to have neutral hydrogen column density in excess of log(NHI/cm−2) ≥ 20.3 (Wolfe et al. 2005), act as signposts of significant reservoirs of neutral hydrogen within or around high-redshift galaxies (Wolfe et al. 1986). For this reason, direct associations between DLAs detected in ab-sorption and galaxies detected in emission provide a power-ful way to probe links between the gas supply, in the form of neutral hydrogen, and ongoing star formation. The combi-nation of absorption and emission techniques thus provides at present the only means to study the star formation law in atomic gas beyond z > 2 (Wolfe & Chen 2006; Prochaska & Wolfe 2009; Rafelski et al. 2011; Fumagalli et al. 2015; Rafelski et al. 2016).
Starting from earlier searches from the ground and with the Hubble Space Telescope (Warren et al. 2001; Møller et al. 2004), several surveys have attempted to identify the galaxy population that gives rise to DLAs. Despite decades of searches, progress has been scarce until recently. Building on the evidence that galaxies obey a defined mass-metallicity relation (Tremonti et al. 2004; Maiolino et al. 2008), searches have focused on the high end of the metallicity distribution of DLAs, yielding higher detection rates with the discovery of tens of DLA hosts (Fynbo et al. 2010, 2013; Krogager et al. 2017). The advent of integral field spectrographs (e.g. SINFONI at the Very Large Telescope, VLT; P´eroux et al. 2012 and OSIRIS at the W. M. Keck Observatory; Jorgen-son & Wolfe 2014) have enabled more efficient spectroscopic follow-up as all of the relevant solid angle around the quasars can be covered in a single setting. The much increased sen-sitivity of the Atacama Large Millimetre Array (ALMA) is now also enabling searches of DLA galaxies via molecular and atomic lines, a technique that is being successfully pio-neered, still, at the high end of the metallicity distribution (Neeleman et al. 2017; Kanekar et al. 2018; Fynbo et al. 2018; Neeleman et al. 2019).
While these recent identifications offer a way to finally study the link between column density and metallicity in ab-sorption, and stellar masses and star formation rates (SFRs) in emission (Møller et al. 2004; Christensen et al. 2014; Kro-gager et al. 2017), these studies are likely to only probe the bright end of the DLA population, and may not be fully rep-resentative of the diverse population of DLA host galaxies. Most simulations and models (Haehnelt et al. 1998; Fynbo et al. 1999; Pontzen et al. 2008; Barnes & Haehnelt 2009; Rahmati & Schaye 2014; Bird et al. 2014; Fynbo et al. 2008) consistently indicate that DLA hosts are generally to be found at the very faint end of the luminosity function, typi-cally below the sensitivity limit of current searches. Indeed, our own survey of galaxies designed to image the DLA hosts at all impact parameters (O’Meara et al. 2006; Fumagalli et al. 2010, 2014, 2015) has yielded a series of non-detections within a few kiloparsecs of the location of the DLAs, despite removing the primary source of observational bias, i.e. the bright quasar emission which hampers the detection of faint galaxies at low impact parameters.
While there seems to be general consensus that DLAs are primarily associated with faint galaxies (e.g. Fynbo et al. 1999; Krogager et al. 2017), the question of what the
typi-cal range of halo masses giving rise to DLAs remains open. At face value, following a similar argument than the one adopted in abundance matching studies (e.g. Conroy et al. 2006), DLAs are expected to arise primarily in faint galax-ies (e.g. Fynbo et al. 2008) and hence in low mass halos (Mhalo . 1011M ). However, this appears to be in tension with other pieces of evidence. Firstly, the distribution of ve-locity widths measured from metals in low ionisation states shows a prominent tail at high velocity, which has been used to argue for the existence of a population of large disks hosting DLAs (Prochaska & Wolfe 1998; but see Bird et al. 2015 for more recent work on this topic). Furthermore, Font-Ribera et al. (2012), using the SDSS-III Baryon Oscillation Spectroscopic Survey (BOSS), found an unexpected large linear bias of DLAs (bD L A = 2.17 ± 0.20, with an assump-tion on the bias of the Lyα forest) by cross-correlating these absorbers with the Lyα forest. This measurement was more recently updated by P´erez-R`afols et al. (2018b), who found a linear bias of bD L A = 2.00 ± 0.19, which is only slightly lower than the clustering amplitude of Lyman Break Galax-ies (LBGs; see also Cooke et al. 2006). This value of the bias, which implies masses of & 1011 M , appears uncomfortably high for some galaxy formation simulations and models (see e.g Pontzen et al. 2008; Barnes & Haehnelt 2014; Padmanab-han et al. 2017) if DLAs sample a very wide range in galaxy sizes (e.g. Krogager et al. 2017). Recently, P´erez-R`afols et al. (2018a) have shown that the bias of DLAs has a depen-dence on metallicity in line with observations and modelling of DLAs which associate more metal-rich DLAs with more massive galaxies (Neeleman et al. 2013; Christensen et al. 2014).
Building on our previous searches for galaxies at small impact parameters from the quasars (Fumagalli et al. 2015), this study exploits the power of wide field integral field spec-troscopy provided by the Multi Unit Spectroscopic Explorer (MUSE) at the VLT (Bacon et al. 2010) to carry out the first highly-complete spectroscopic survey in an sample of DLAs, not selected in metallicity or column density, with the goal of searching for faint Lyα emission from galaxies associated with the DLAs. Our survey improves upon previ-ous searches of this type, which were characterised by a lower sensitivity (e.g. Christensen et al. 2007) or smaller field of view (e.g. P´eroux et al. 2012). By achieving a flux-limited search up to ≈ 250 kpc from the DLAs, we are finally able to search at the same time for the DLA hosts and for other associations at larger impact parameters, which provides a means of characterising the DLA environment via the small-scale clustering of galaxies and DLAs.
asso-20.0 20.5 21.0 21.5 22.0 log10(NHI/ cm−2) −3.0 −2.5 −2.0 −1.5 −1.0 −0.5 0.0 0.5 [M/H] 3.2 < z < 3.8DLAs MUSE Sample 0.00 0.15 0.30 f( NHI ) 0.00 0.15 0.30 f( [M/H] )
Figure 1. The metallicity and column density distributions of the DLAs targeted with MUSE, compared to the parent DLA population over the same redshift interval (3.2 < z < 3.8) from Rafelski et al. (2012). The normalized distribution of metallicity is shown on the right for both samples, and similarly at the top for column density. In both histograms the MUSE sample are showed in blue hashed, with the parent population in grey behind. Our sample spans a wide range of DLA properties, thus extending the parameter space covered by most recent surveys.
ciated, a condition that is imposed by the finite volume of the adopted simulations.
This paper is structured in the following way. The ob-servations, data reduction and analysis are presented in Sect. 2 and 3, followed in Sect. 4 by technical details on the simu-lations and semi-analytic models used in the interpretation. We briefly describe some highlighted detections in Sect. 5. The discussion of results is presented in Sect. 6, followed by a summary and conclusion in Sect. 7. Further discussion of the properties of each field is continued from Sect. 5 in appendix B for the remaining fields. Readers primarily interested in the main results of the paper can focus their attention on Sect. 5, 6 and Sect. 7.
Magnitudes reported in this paper follow the AB system and fluxes and magnitudes have been corrected for Galac-tic extinction following Schlafly & Finkbeiner (2011). All spectral data products are reported in vacuum wavelengths, and we adopt the Planck Collaboration et al. (2016) ΛCDM cosmological parameters (Ωm= 0.308, h = 0.677).
2 OBSERVATIONS AND DATA REDUCTION
2.1 Sample selection
The six quasar sightlines observed with MUSE are chosen from the parent sample of Fumagalli et al. (2014) which is selected among the general DLA population purely based on the presence of two optically-thick absorbers along the quasar sightline, with the goal of searching for the rest-frame UV emission at close impact parameters from the DLAs (see O’Meara et al. 2006, for details). For this study, we select
DLAs that are observable from VLT and at z > 3.2, which is the redshifts at which Lyα falls at wavelengths where the throughput of MUSE is > 25%. Excellent ancillary data, including deep UV and optical imaging and high-resolution spectroscopy are available for this sample, which spans a wide range of metallicity and column density.
The properties of the selected DLAs in our sample are summarised in Table 1, while Figure 1 shows the MUSE sample in context with the wider DLA population in terms of metallicity and column density, highlighting how our tar-gets span a representative range of both parameters. This is in contrast to many recent searches for DLA hosts, such as those conducted with X-Shooter (Krogager et al. 2017) and ALMA (Neeleman et al. 2017, 2019), which typically select high metallicity systems ([Si/H] & −1). This selec-tion exploits the metallicity-luminosity relaselec-tion between the metallicity of the DLA in absorption and the luminosity of the host galaxy to ensure higher detection ratio compared to samples selected, e.g., only on NH I(e.g. P´eroux et al. 2012) or MgII(e.g. Bouch´e et al. 2012). This approach, however, leaves the hosts of the majority of “typical” DLAs at z > 2, with metallicities ≤ 5% solar, unexplored. Our DLA sample contains instead three systems above average metallicity at this redshift (Rafelski et al. 2012) and three below, thus ex-tending the parameter space targeted by previous searches. Notably, DLA J1220+0921 in our sample is very metal-poor, sitting at Z/Z = −2.33 although still somewhat above the metallicity floor at Z/Z ' −3 (Prochaska et al. 2003; Rafel-ski et al. 2012). One sightline in our sample was previously presented in Fumagalli et al. 2017b, the DLA was revealed to be part of a ' 50 kpc structure which may be evidence of an ongoing merger.
2.2 MUSE observations and data reduction MUSE integral field spectroscopy of the six DLA fields was acquired in service mode between June 2015 and April 2016 under ESO programmes 095.A-0051 and 096.A-0022 (PI Fu-magalli). Observations were split into sets of 1480s expo-sures. While four fields were completed with 6×1480s ex-posures, J0851+2332 and J0818+2631 were only partially completed with 4×1480s and 2×1480s respectively. The ob-servations were taken in dark time with seeing ranging from 0.7 arcsec to 0.9 arcsec using the Nominal Wide Field Mode, with clear conditions. In each sub-exposure the quasar was centred in the field of view, and small dithers combined with instrument rotations in increments of 90◦were made to im-prove the quality of the final data product.
The initial reduction of the data is carried out using the ESO MUSE pipeline (Weilbacher et al. 2014) (v1.6.2). The pipeline carries out bias, dark and flat field corrections, calibrates the data in wavelength and astrometry, and ap-plies a basic illumination correction. This initial reduction is however limited by the accuracy of the illumination correc-tion and therefore we addicorrec-tionally post-process the data with CUBExtractor (Cantalupo, in prep.) to further improve the quality of the illumination correction and sky subtrac-tion (see, e.g. Borisova et al. 2016; Fumagalli et al. 2016, 2017a for details).
−5000 −2500 0 2500 5000 0.0 0.2 0.4 0.6 0.8 1.0 Arbitrar y Flux J2351+1600 −7000 −3500 0 3500 7000 0.0 0.2 0.4 0.6 0.8 1.0 J0851+2332 −5000 −2500 0 2500 5000 0.0 0.2 0.4 0.6 0.8 1.0 J0255+0048 −5000 −2500 0 2500 5000 Velocity (km s−1) 0.0 0.2 0.4 0.6 0.8 1.0 Arbitrar y Flux J1220+0921 −7000 −3500 0 3500 7000 Velocity (km s−1) 0.0 0.2 0.4 0.6 0.8 1.0 J0818+0720 −5000 −2500 0 2500 5000 Velocity (km s−1) 0.0 0.2 0.4 0.6 0.8 1.0 J0818+2631
Figure 2.The Lyα absorption profiles for the six DLAs in the MUSE sample. The quasar spectra are shown in black, the fitted quasar composite template are in blue (solid line), and the template multiplied by the Voigt profile is in red (solid line). Errors (1σ) on the column density are show by red dotted lines. The original continuum normalisation from Fumagalli et al. (2014) is also shown as a blue dash-dotted line. For J0851+2332 a manual spline continuum normalisation was adopted instead of the quasar composite. The fitted contribution of the Lyβ line from a proximate DLA contaminating the profile for J1220+0921 is also plotted with grey dashed line.
QSO Name Fielda R.A. Dec z
DLA log(NH I/ cm−2) log (Z/Z ) Elementb Exp. Time (s)c J2351+1600 06G6 23:51:52.80 +16:00:48.9 3.7861 21.00 ± 0.10 -2.18 ± 0.20 Fe 6x1480 J0851+2332 10G11 08:51:43.72 +23:32:08.9 3.5297 21.15 ± 0.20 -1.10 ± 0.23 Zn 4x1480 J0255+0048 15H3 02:55:18.58 +00:48:47.6 3.2530 20.85 ± 0.10 -1.10 ± 0.10 Si 6x1480 J1220+0921 19H7 12:20:21.39 +09:21:35.7 3.3090 20.30 ± 0.20 -2.33 ± 0.22 Si 6x1480 J0818+0720 23H11 08:18:13.14 +07:20:54.9 3.2332 21.15 ± 0.10 -1.66 ± 0.25 Si,Zn 6x1480 J0818+2631 24H12 08:18:13.05 +26:31:36.9 3.5629 20.65 ± 0.10 -0.93 ± 0.24 Si,Zn 2x1480
Table 1.The properties of the MUSE DLA sample.aThe naming system from Fumagalli et al. 2014.bThe element(s) used to estimate the DLA metallicity.c The exposure time of the MUSE integral field spectroscopy.
and odd numbered sub-exposures to produce independent sets of data for verification processes (see Sect 3.2). The final data product for each field covers a field of view of ap-proximately 1 × 1 arcmin2, covering 4750 − 9354 ˚A with 1.25 ˚
A binning. Regions of sub-exposures affected by stray-light, near the edges of the cube, are masked before stacking. We also mask 2-5 pixels around the edge of the combined cubes to remove low quality data and very noisy pixels. The field of DLA J0255+0921 further requires substantial masking due to a bright (V= 11.2) star lying at 46 arcsec from the quasar, resulting in a slightly smaller final field of view.
As a last step, we calibrate the absolute astrometry and verify the quality of the data against imaging and spectra
−400 −200 0 200 400 0.0 0.2 0.4 0.6 0.8 1.0 1.2 Nor malised Flux FeII 1608 ESI J2351+1600 −400 −200 0 200 400 0.0 0.2 0.4 0.6 0.8 1.0 1.2 FeII 1608 ESI J0851+2332 −400 −200 0 200 400 0.0 0.2 0.4 0.6 0.8 1.0 1.2 FeII 1608 HIRES J0255+0048 −400 −200 0 200 400 Velocity (km s−1) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 Nor malised Flux SiII 1527 MagE J1220+0921 −400 −200 0 200 400 Velocity (km s−1) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 FeII 1608 ESI J0818+0720 −400 −200 0 200 400 Velocity (km s−1) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 FeII 1608 ESI J0818+2631
Figure 3.Selected low-ionisation lines from each DLA in our sample. The spectra (black) are continuum normalised, with the 1σ flux error shown in red. In each panel the instrument and transition are listed. Lines indicating zero flux (blue dotted line), the normalised continuum level (1.0, blue solid line) and the zero velocity (blue dashed line) are included. Note that many of these lines are saturated and so are not used to derive metallicities.
correction applied. Comparison to SDSS and high-dispersion spectra of the quasars show excellent agreement.
2.3 Absorption line spectroscopy
In order to measure the absorption properties of the DLAs, spectra with higher resolution than MUSE (R'2000 at 5500˚A) are required, particularly for narrow absorption fea-tures and to resolve saturation in strong absorption lines. In this work, we use the same data originally presented in Fumagalli et al. (2014), but we refine the measurement of the HI column density by improving the determination of the quasar continuum.
The spectroscopic data available are described in Table 2. For all DLAs in our sample, we have moderate dispersion spectra (R > 5000) over the optical range, encompassing the DLA Lyα and common low ions (e.g. SiII and FeII), which are used to estimate the DLA metallicity. For the majority of the sample, this is ESI data, while in the case of J1220+0921 the spectrum is from MagE (Jorgenson et al. 2013). Finally, J0255+0048 has additional higher dispersion data from HIRES (Prochaska et al. 2001) covering some low ions at R ≈ 50, 000, and an X-Shooter spectrum (R ≈ 6000
L´opez et al. 2016) with higher signal-to-noise ratio (SNR) than the ESI data and coverage of the near infrared.
2.3.1 HIcolumn densities
Field Instrument Exp. Time S/N Resolution Wavelength Coverage Reference (s) (per pixel) (km s−1) (˚A)
J2351+1600 ESI 3600 12 37 3995-10140 Fumagalli et al. (2014)
J0851+2332 ESI 3600 19 37 3995-10140 Fumagalli et al. (2014)
J0255+0048 ESI 2400 16 37 3995-10140 Fumagalli et al. (2014)
X-Shooter 3320-3600a 30 40-60a 3100-18000 L´opez et al. (2016)
HIRES 20200 15 6.3 5800-8155b Prochaska et al. (2001)
J1220+0921 MagE 3000 22 71 3100-10360 Jorgenson et al. (2013)
J0818+0720 ESI 3600 18 56 3995-10140 Fumagalli et al. (2014)
J0818+2631 ESI 3600 28 56 3995-10140 Fumagalli et al. (2014)
aDiffers between spectograph arms.bWith gaps between echelle orders. Table 2.A summary of the spectroscopic data used to establish the absorption properties of the DLAs.
exploit the “double-DLA” technique. In this work, we refit the NH I column densities replacing the original local con-tinuum determination by the Telfer et al. (2002) composite quasar spectrum. For each sightline, the template contin-uum power-law slope and normalisation are adjusted to fit the quasar continuum over the Lyα forest, and the HI col-umn densities of the DLAs are then estimated by fitting a Voigt profile to the Lyα transition at the redshift of the DLA. A characteristic uncertainty of 0.1 dex is assigned to our determinations. Final values of column density are listed in Table 1, while a gallery of Voigt profile fits is in Fig. 2.
During this analysis, we noted that the composite quasar spectrum provides a poor fit to the continuum of quasar J0851+2332, which appears to be significantly lower (approximately a factor of 2) blueward of the DLA Lyα. This feature cannot be modelled as a high redshift partial Lyman limit system1, as it would require a redshift above that of the quasar. Furthermore, this break is observed in the MUSE, SDSS and ESI spectra, ruling out an artefact with the data. For this reason, in this sightline, we use a manu-ally estimated continuum, found by fitting a spline though points in the forest believed to represent the unabsorbed continuum. Using this model we reach a column density of log(NHI/cm−2) = 21.15. We note that, despite some de-gree of subjectivity in this estimate, the column density is mostly insensitive to changes in the continuum, as its up-per bound is set by the width of the core of the Lyα line in this case. Indeed, if we fit the composite quasar to the spectrum between 5800 and 9000 ˚A, we obtain a column density of log(NHI/cm−2) = 21.25, only marginally different from our previous estimate. To capture this discrepancy, for this DLA we adopt an uncertainty of ±0.2 dex. Finally, when fitting DLA J1220+0921, we note that the Lyα of the DLA at zDLA = 3.090 lies close to the Lyβ of a proximate DLA (pDLA) at zpDLA= 4.1215. For this line of sight, we therefore also include the Lyβ transition of this pDLA in our fit. 2.3.2 Metallicities
In order to calculate the metallicities of the DLAs, we adopt the metal ion column densities compiled in Fumagalli et al. (2014). These values are based on the apparent optical depth
1 A Lyman Limit System (LLS) is an absorption line with 1017.5 < (NH I/ cm−2) < 1020.3, such that the system is optically thick to ionising radiation below the Lyman limit (i.e. λ <912˚A). LLSs are thus characterised by breaks in quasar spectra below the Lyman limit.
method using unsaturated transitions, or are bracketed by upper and lower limits. Table 1 indicates the ions used to estimate the metallicity of the DLAs in the sample. Fuma-galli et al. (2017b) showed that in the case of J0255+0048 the SiIIcolumn density estimated with Voigt profile
decom-position of HIRES and X-Shooter data was consistent with the value obtained with ESI and the apparent optical depth method. We therefore conclude the ionic column densities to be robust and do not re-estimate them. As the transitions used do not lie in the Lyα forest, the continuum level can be estimated accurately from the data without the need of a quasar template as done for the HItransition. Fig. 3 shows strong low-ionisation absorption lines for each DLA. Some of the lines shown are saturated and were not used to calculate metallicities. From the ionic column densities and the NH I values measured above combined with the solar elemental abundance pattern (Asplund et al. 2009), we calculate the metallicities of the DLAs, which are listed in Table 1.
3 SEARCH FOR GALAXY ASSOCIATIONS
To identify galaxy associated to the DLAs, we have con-ducted a search for both Lyα emitters (LAEs) and Lyman break galaxies (LBGs) by searching for emission line objects in the MUSE datacubes and fitting redshifts for continuum objects detected in deep white-light images reconstructed from the datacubes. No candidates were found in the con-tinuum object search to a magnitude limit of mr < 25 mag, thus the emission line search has been adopted as our pri-mary method of identifying galaxy associations. The search for continuum objects is nevertheless described for complete-ness.
3.1 Search for Lyα emitters
We have conducted a search for LAEs with zL AE ' zD L A with a velocity window of ±1000 km s−1over the full MUSE field of view. Initially, the mean coadded cubes are trimmed in the wavelength direction to restrict the wavelength range to that of Lyα over the velocity interval of interest around the DLA of each field, plus a margin of 300 MUSE channels on either side (375 ˚A).
−4 −2 0 2 4 −4 −2 0 2 4 ∆ Dec (”) Detection Cube −4 −2 0 2 4 1st Half −4 −2 0 2 4 2nd Half −4 −2 0 2 4 Continuum 5230 5240 5250 5260 0 5 10 15 20 25 30 1D Spectrum −2 −1 0 1 2 −2 −1 0 1 2 ∆ Dec (”) Detection Cube −2 −1 0 1 2 1st Half −2 −1 0 1 2 2nd Half −2 −1 0 1 2 Continuum 5470 5480 5490 5500 0 2 4 6 8 10 12 14 16 Flux (10 − 18 erg s − 1cm − 2˚ A − 1) −2 −1 0 1 2 ∆R.A. (”) −2 −1 0 1 2 ∆ Dec (”) Detection Cube −2 −1 0 1 2 ∆R.A. (”) 1st Half −2 −1 0 1 2 ∆R.A. (”) 2nd Half −2 −1 0 1 2 ∆R.A. (”) Continuum 5480 5490 5500 5510 Wavelength ( ˚A) 0 2 4 6 8 10 12
Figure 4.Three examples of LAE candidates that illustrate the three most common occurrences in our identification procedure. The first object (top row) is a high confidence LAE, while the later two are examples of false candidates: a cosmic ray and a low redshift galaxy, respectively. The first three columns in each row show optimally extracted images from the detection cube after continuum subtraction, and the two independent coadds (labelled 1stand 2ndhalf). All three images are shown on a linear scale with the same stretch for a single candidate. The fourth column shows the reconstructed r-band image of the same fields with surface brightness contours taken from the detection map of the candidate emission line. All images are smoothed with a Gaussian kernel with full width at half maximum of 3 pixels. The right column shows 1D spectra of the detected emission line for each object (black) with the flux error (red) and the line centre (orange dashed vertical line).
velocity range around the DLAs is masked during this pro-cess, to ensure that no emission-line objects are subtracted close to the DLA redshift. Because of this masking, the con-tinuum subtraction in this region is performed by extrapo-lating the continuum from the unmasked region.
After this step, the resulting cubes are more finely trimmed to the wavelength range of interest, plus a margin of one channel to prevent extended objects from being trun-cated. This continuum-subtracted datacube then becomes the detection cube for our search. The mean, median and two independent coadded datacubes are also trimmed to the same wavelength range as the detection cube for quality con-trol purposes, as described below. Finally, we run CubEx-tractor over the detection datacubes. The first step in this process is to rescale the data variance in order to match the root-mean-square (RMS) of the flux. This rescaling is needed to correct for small, albeit significant, degree of correlation in the data following the drizzling process (see Appendix A). We then convolve the detection datacubes with a two pixel boxcar spatially, grouping connected pixels above a
SNR threshold of 2. These groups are considered a detec-tion if the following criteria are met: i) the object has at least 40 voxels above the threshold; ii) it spans at least 3 wavelength channels at a single spatial position; iii) it has an integrated SNR ≥ 7. A segmentation map identifies vox-els which are above the SNR threshold and are part of a detection. Segmentation maps are further examined to en-sure that connected regions resemble point-like or extended sources (i.e. they are not extremely elongated in a single direction).
as these are likely to be associated with sky line residuals or cosmic rays.
The remaining candidate LAEs are then inspected, us-ing optimally-extracted images (see Borisova et al. 2016) from the mean, median, independent coadds and detection datacubes. During this step, we find that comparing the two independent coadds is the strongest discriminant to reject objects which appeared to differ significantly in morphology or position between the images. Finally, the 3D segmenta-tion map is projected onto a 2D grid to extract the object spectrum over the full MUSE wavelength range. Inspection of these spectra enables us to cull other emission lines (com-monly [OII]) or the residuals of continuum objects, retaining only bona-fidae LAEs. Bright continuum objects in the field of view can leave residuals during the continuum subtrac-tion, these are often detected by CubExtractor. These are unrelated to sky line residuals or cosmic rays and can easily be identified by looking at the datacube without continuum subtraction. This procedure yields 14 candidate detections over the six DLA fields, as summarised in Table 4.
Fig. 4 provides examples for three illustrative cases in our identification procedure. The first example is a high con-fidence LAE (ID85 from field J1220+0921), where the bright core of the object does not shift in the independent coadds indicating the detection is robust. The second example is a highly-significant detection which, however, only appears in one of the independent coadds, indicating that this candi-date is a cosmic ray strike. Examples as significant as this are very rare. The last example is a source that appears strong and marginally extended in the detection cube, however it is much brighter in both independent coadds once the con-tinuum is not subtracted. This feature, combined with the properties of the r-band image and 1D spectrum, makes it clear that this candidate is a bright low-redshift galaxy and that the detected feature is in the residual associated with continuum subtraction or an emission line other than Lyα. 3.2 Testing the robustness of LAE identifications The resulting candidates from the procedure described above vary in integrated SNR from 7.4 to 29.9. While all detections are robust from a statistical point of view, it is worth examining where an additional cut in SNR is war-ranted to avoid additional spurious detections that are not rejected in our procedure, especially given that the noise field is non Gaussian. To this end, we have repeated the de-tection procedure described above, but with the flux values in the datacubes flipped in sign (hereafter the negative dat-acubes) to explore whether noise fluctuations appear in our selection, and with what SNR. This practice is a standard technique in imaging observations (e.g. Rafelski et al. 2009; Hodge et al. 2013).
When searching the negative datacubes, we adopt iden-tical parameters as in the search of real source. However, since real sources with absorption features can generate neg-ative residuals during the continuum subtraction procedure, we remove detections overlapping with bright continuum sources. This task generates 5 detections in the negative cubes meeting our requirements. The properties of these de-tections are compared to the real candidates in Fig. 5 in terms of the velocity offset from the DLA redshift and the detection SNR. −1000 −500 0 500 1000 ∆v LAE-DLA (km s−1) 5 10 15 20 25 30 Integr ated SNR DLA LAEs Negative Detections Control LAEs 0 2 4 6 8 N detections 0 2 4 6 N detections
Figure 5.The distribution of velocity separations (∆v) between DLAs and candidate LAEs against integrated emission line signal-to-noise ratio. Real detections (blue, filled), LAEs selected in con-trol windows (green, gridded) and negative false detections (grey, hatched) are shown. The number of voxels in the segmentation map for each candidate are indicated by the sizes of the symbols. The histogram on each axis shows the distribution of candidates over that single parameter.
Additionally, as a further test, we perform the extrac-tion method described in Sect. 3 over the six datacubes, but this time shuffling the central redshifts (i.e. the DLA red-shift) across the fields. This experiment yields control source catalogues containing both real LAEs, at a redshift far from the DLAs, and additional spurious detections, if any. The 8 detections produced by this search are displayed in Fig. 5. This method, unlike the use of the negative detection cubes, does not require symmetry in the noise properties and pro-vides a baseline which we can use to assess if an excess of LAEs clustered to the DLAs is detected in our data.
Two key points are apparent from Fig. 5. First, the detections from the negative datacubes are all SNR < 11, while several of the the true detections are at high SNR, up to SNR= 29.9. Secondly, sources identified in the detection datacubes are significantly more clustered around ∆v '+200 km s−1 than both the detections from the negative cubes and the control windows. Due to radiative transfer effects, the Lyα line is expected to be redshifted compared to the DLA redshift by ≈ 100 − 300 km s−1 (e.g Steidel et al. 2010; Rakic et al. 2011).
datacubes pass our selection criteria, we take a conservative approach and establish SNR= 12 as a threshold for identify-ing LAEs with high purity, although at the expense of sam-ple comsam-pleteness. In the following, we refer to these objects (SNR > 12) as high-confidence confirmed LAEs, while the remaining sources form a sample of candidate associations that for most part are believed to be real, but for which we cannot exclude the presence of some spurious sources. Deep follow-up observations will be required to determine the nature of each of these sources.
We further note that observations for DLA J0818+2631 suffer quite badly from only having 2 out of 6 exposures completed, meaning that independent coadds contain only a single exposure that presents significant gaps in the re-duced datacubes due to the masking around the gaps be-tween the stacks of IFUs in the MUSE FoV. Thus, the search for sources in this field is likely to be incomplete.
As a last check, we investigate the robustness of our de-tections, in particular considering whether correlated noise affects the estimates of SNR values and the degree of incom-pleteness in the detected LAEs. The results of these tests are detailed in Appendix A, where we conclude that correlated noise is not a substantial effect in MUSE data (see also Ba-con et al. 2017) and that our sample of high-Ba-confidence LAEs is highly complete.
3.3 Identification of continuum objects
For the identification of continuum detected objects, we first extract objects using the r-band images reconstructed from the cubes running SExtractor (Bertin & Arnouts 1996) with minimum area of 6 pixels, each above an SNR of 2.0. The SNR threshold was raised for DLA J0818+2631, which had many residuals due to being only partially completed. The resulting segmentation maps are used to define apertures on the datacubes over which we extract spectra for the se-lected objects. In this work, we attempt to determine red-shifts only for objects with mr <= 25.0 mag (corrected for Galactic extinction). This choice is motivated by previous MUSE analyses (Fumagalli et al. 2015) which have shown how high confidence redshifts for objects fainter than this approximate limit are only obtained in presence of bright emission lines (usually Lyα). As we already search for Lyα emission close to the redshifts of the DLAs, faint LBGs with strong Lyα emission associated to DLAs will not be missed. These spectra are then inspected for emission and ab-sorption lines, as well as characteristic continuum features. Their redshifts are measured by two authors (RM and DJH), either by fitting Gaussian functions to the detected emission lines, or by comparing the 1D spectra with a range of stellar and galaxy templates, including low redshift galaxies and high redshift LBGs. As no mr <= 25.0 objects were found to lie close to a DLA in redshift (i.e. within 1000 km s−1), these objects are not relevant to our analysis. One LAE de-tected in the emission line search (id56 in field J0255+0048) has a continuum detection with mr = 24.58 but it was not detected in the continuum object search because it was very close to the QSO and was not identified as an independent object. This was not observed with any other mr<= 25 ob-jects. The other two LAEs with continuum counterparts are fainter than our magnitude limit for the continuum search
and so were not selected. Overall, for mr<= 25.0 objects we obtained a redshift completeness of 74%.
4 DESCRIPTION OF MODELS AND
SIMULATIONS
In the following, we compare our observational results to simulations and semi-analytic models to better understand the constraints they put on the association between DLAs and galaxies. In this section, we provide a detailed descrip-tion of how different models are produced and analysed.
The hydrodynamic simulations adopted on this work are taken from the eagle suite (Schaye et al. 2015), with snapshots post-processed using the urchin reverse ray-tracing code (Altay & Theuns 2013) to identify DLAs. In our analysis, we use the post-processed eagle snapshots combined with simple prescriptions to populate halos with LAEs to produce mock observations of the correlation be-tween DLAs and galaxies. We additionally use a mock cat-alog based on the galics semi-analytic model, designed to produce realistic LAEs (Garel et al. 2015). For this model, DLAs are “painted” onto dark matter halos. This simple pre-scription allows us to quickly investigate the effects of vary-ing the DLA cross section as a function of halo mass. For a given simulation and prescription we generate a grid which indicates which sightlines through the simulations encounter a painted DLA. These grids are produced by projecting cir-cular kernels centred on halos onto the grid along one axis of the simulation. An additional grid keeps track of position of the DLA in 3 dimensions. Comparisons with the data then allow us to judge to what extent current or future data can distinguish between various models. This simple model of “painting” DLAs onto halos is also applied to the eagle simulations, to further gain insight into the properties of DLAs identified with urchin.
4.1 The eagle simulations
2017; Lagos et al. 2015). At z ' 3 eagle has been shown to match the observed star formation rate density well (Kat-sianis et al. 2017), and bears reasonable agreement with the metallicity distributions of high column density absorption line systems (Rahmati & Oppenheimer 2018).
The eagle simulations are performed in cubic periodic volumes, and the linear extent (L) of the simulation vol-ume and number of simulation particles is varied to allow for numerical convergence tests. In this work, we use simu-lations L0100N1504 and L0025N0752 from table 1 of Schaye et al. (2015). Briefly these have L = 100 co-moving mega-parsecs (cMpc) and L= 25 cMpc, and a SPH particle masses of 1.8 × 106 M and 0.2 × 106 M , respectively; the Plum-mer equivalent co-moving gravitational softening is 2.66 and 1.33 kpc, limited to a maximum physical softening of 0.7 and 0.35 kpc, respectively. The simulations assume the cosmo-logical parameters of Planck Collaboration et al. (2014), and the minor differences with the Planck Collaboration et al. (2016) parameters adopted in the current paper are unlikely to be important.
To compare to the data presented earlier, we need to identify both DLAs and LAEs in eagle. Since neither neu-tral hydrogen nor Lyα radiative transfer are directly in-corporated in eagle, we compute these quantities in post-processing as explained next.
4.1.1 Identifying DLAs in the eagle simulations
The ionising background in eagle is implemented in the optically-thin limit. Self-shielding of gas, allowing for the appearance of DLAs, is computed in post-processing using the urchin radiative transfer code described by Altay & Theuns (2013). Briefly, this algorithm allows each gas par-ticle to estimate the local ionising intensity it is subject to by sampling the radiation field in 12 directions, with neu-tral gas in neighbouring gas particles potentially decreasing the local photo-ionisation rate below the Haardt & Madau (2001) optically-thin value. Assuming the neutral fraction of a particle is set by the balance between photo- and collisional ionisation versus recombinations, a reduced ionisation rate increases the particle’s neutral fraction, which in turn affects the ionisation rate determined by the particle’s neighbours. The impact of a change in the neutral fraction on the photo-ionisation rate and vice versa is iterated until the neutral fraction of each particle converges from one iteration to the next. This post-process step thus yields the neutral hydrogen fraction xH i≡ nH i/nHof each SPH particle. Further details on how radiative transfer affects the neutral fraction as a function of column density and physical processes (e.g. the relative importance of collisional ionisation and photoionsa-tion), can be found in the literature (e.g. Fumagalli et al. 2011; Rahmati et al. 2013).
The resulting HIvolume density is then projected onto a 81922 grid along the coordinate z-axis of the simulation box, using the Gaussian smoothing described by Altay & Theuns (2013). Applying a column density threshold of log(NH i/cm−2)= 20.3 allows us to identify DLAs2. We also
2 The redshift path through the simulation box is so small that the contribution of chance alignments of two high-column density systems to the DLA cross section is negligible.
calculate the xH i-weighted z-coordinate and velocity of par-ticles along each DLA line of sight to obtain the 3D position of each DLA (two spatial coordinates x and y, and a redshift coordinate z). The redshift of each DLA allows us to com-pare DLAs to galaxies in redshift space, although we note that observed DLA redshifts are derived from low-ionisation metal lines rather than HI directly, yet these elements are believed to trace the same phase of the gas.
To compare to the muse observations, we analyse the L0025N0752 eagle snapshot at z = 3.027 which is closest in redshift to the data. Similarly to the results of Altay et al. (2011) based on the owls simulations (Schaye et al. 2010), the column density distribution function (CDDF, the num-ber density of absornum-bers per unit column density, per unit ab-sorption distance, f (NH i)) of the post-processed snapshots using urchin is in good agreement with observations (Fig. 6). The data in this figure combines the low column den-sity data (log(NHI/cm−2)< 18) compiled by Kim et al. 2013 (z= 2.4−3.2), the multiple power-law fit to the Lyman-limit and DLA column density range derived by Prochaska et al. (2010) at z = 3.7, and the sub-DLA and DLA data from Noterdaeme et al. (2012) (hzi= 2.5). At these redshifts, we find that the simulated and observed incidence of DLAs is very similar. The measurements of Zafar et al. (2013) in the sub-DLA range are also shown, these constraints cover a much broader redshift distribution than the MUSE DLAs (1.5 < z < 5) but show close agreement with the fit of Prochaska et al. (2010).
4.1.2 Identifying LAEs in eagle
Gravitationally bound substructures in eagle are identified combining the friends-of-friends, (Davis 1985) and subfind (Springel 2005; Dolag et al. 2009) algorithms. Physical prop-erties of these ‘galaxies’ such as their centre of mass position and velocity, stellar mass and star formation rate, are com-puted and stored in a database (McAlpine et al. 2016). In eagle the star-formation rate of a gas particle is a func-tion of its pressure, it is zero below a metallicity-dependent density threshold and above a pressure-dependent tempera-ture threshold (Schaye et al. 2015). SFRs and stellar masses of subhalos are computed by summing over the gas and star particles (respectively) within a given subhalo identi-fied with subfind. Lacking sufficient resolution to resolve the ISM, we cannot predict from first principles the Lyα properties of the simulated galaxies. We resort instead to a simpler empirically-motivated model. We populate the simu-lation with LAEs by selecting simulated galaxies by star for-mation rate (SFR) and associating a Lyα luminosity (LLyα) using the conversion from Eq. 1 (Furlanetto et al. 2005). This relation combines hydrogen case B recombination with a standard SFR calibration (Kennicutt 1998).
LLyα= SFR M yr−1
× 1042erg s−1 (1)
ap-1041 1042 1043
log10LLyα(erg s−1)
10−5 10−4 10−3 10−2 10−1 φ ([dlog 10 ] − 1Mpc − 3) Schechter fit EAGLE 100 Mpc box EAGLE 25 Mpc box GALICS No Lyα Calibration Drake et al. 2017 16 17 18 19 20 21 22 23 log10(NHI/cm−2) −26 −24 −22 −20 −18 −16 log 10 (f(N HI ,χ )) Prochaska et al. 2010 EAGLE/Urchin Noterdaeme et al. 2012 Kim et al. 2013 Zafar et al. 2013
Figure 6. Left: The Lyα luminosity function from the eagle LAE sample compared against recent observations (Drake et al. 2017, blue stars) and a Schechter function fit to these observations (gray dot-dashed). The calibrated (black solid line) and original (red solid line) luminosity functions for the 100 Mpc box simulation are shown, together with the rescaled luminosity function from the 25 Mpc eagle box (dotted line). The luminosity function from galics is plotted (green, dashed), and shows excellent agreement with the data without need for calibration. Right: The column density distribution function derived from the eagle simulation after post-processing with the urchin radiative transfer code, compared against observations at comparable redshift.
proximation for hydrogen, yet may nevertheless yield a large overestimate of the Lyα luminosity of a galaxy because a sig-nificant fraction of such photons do not escape the galaxy, due to to scattering and dust (but see below for how we correct this effect).
The Lyα luminosity function of LAEs is shown in the left panel of Fig. 6, comparing results from L0025N0752 and L0100N1504 with the Schechter fit to the observed lumi-nosity function from Drake et al. (2017), using the redshift bin 2.92 ≤ z < 4.00. As anticipated, by neglecting dust ab-sorption in eagle the bright end of the luminosity function (starting from LLyα > 5 × 1041erg s−1) significantly overpre-dicts the observed number density of emitters, especially at high SFRs, ÛM?≥ 10M yr−1. This trend is well documented observationally: star formation rates inferred from the UV compared to those based on Eq. (1) are discrepant if no correction is made for dust (e.g. Dijkstra & Westra 2010; Whitaker et al. 2017). Moreover, Shapley et al. (2003) show that roughly a third of the LBG population is not detected in Lyα at all. Also, Matthee et al. (2016) show that at z ' 2.2 galaxies with higher Hα-inferred SFRs have lower Lyα es-cape fractions, thus implying that a correction is required at the bright-end of the simulated luminosity function.
We therefore account for these unresolved physical pro-cesses (e.g. dust extinction and escape fraction) by intro-ducing an effective Lyα escape correction to the simulated Lyα fluxes. We do so by sub-sampling the simulated star-forming galaxies until they match the Schechter fit to the luminosity function from Drake et al. (2017). This is done by dividing the Drake et al. Schechter fit by the measured eagle luminosity function, and applying this ratio as prob-ability that a galaxy of a given luminosity will be added to
the LAE catalogue. The result of this re-scaling is shown in Fig. 6, revealing that the sub-sampling performs well for LLyα> 5×1041erg s−1. Given that the high-confidence MUSE sample extend only down to 1042erg s−1, the eagle model appears excellent for comparison to the observations. This resampling technique could be thought of in terms of a Lyα duty-cycle, with the time a galaxy spends in an LAE phase being a function of the SFR (Nagamine et al. 2010).
As consistency check, we compute the clustering of sub-sampled LAEs identified in the L0100N1504 eagle simu-lation (Fig. 7). Ideally our identified LAEs should match observational measurements of the clustering amplitude, as typically quantified by fitting a simple power-law model to the correlation function,
ξ(r) = (r/r0)γ. (2)
Mpc, favoring instead a value of r0 = 5.1 Mpc (between 1 − 15 Mpc). Although the measurement of Gawiser et al. (2007) will suffer from limited volume our simulated LAEs have a higher r0than most observations, we will consider the impact of this offset in Sect. 6.2. Cross-correlations will be less affected by a higher clustering amplitude of one sample than the auto-correlations shown in Fig. 7.
Together, these comparisons demonstrate our simple model for assigning a Lyα luminosity to eagle galaxies, which is calibrated to reproduce the abundances of LAEs. Although there is tension in the clustering amplitude we con-clude in Sect. 6.2 that the agreement is sufficient to enable a valid comparison between the simulations and the data in the current paper. In the following section, we further show how an independent semi-analytic model constructed to re-produce the luminosity function of LAEs agrees well with the prediction derived from the eagle simulations. This semi-analytic modeling represents a cross-check of our results and includes a more physical model of Lyα escape, independently of our method of modeling LAEs in eagle.
4.2 The galics semi-analytic model
As an alternative model to the one based on eagle, we use mock catalogs of LAEs based on the model of Garel et al. (2015) which combines the galics semi-analytic model with numerical simulations of Lyα radiation transfer in galactic outflows (Verhamme et al. 2006; Schaerer et al. 2011) to pre-dict the observed Lyα luminosities of galaxies. As shown in Garel et al. (2015, 2016), the model can reproduce various statistical constraints on galaxies at high redshift, such as the abundances of LAEs at 3 < z < 6. The mock lightcones used in this study were extracted from the galics cosmo-logical simulation volume (Lbox= 100 h−1 cMpc) and were specifically designed to match the redshift range, geometry and depth of typical Lyα surveys with MUSE (see Garel et al. 2016 for more details about galics and the mocks). We refer to these mock catalogues hereafter as galics, in-cluding the additional radiative transfer.
Fig. 6 (left) shows the predicted Lyα luminosity func-tion from galics, compared to MUSE deep field observa-tions (Drake et al. 2017). Fig. 7 also demonstrates that the clustering of the LAEs in the lightcone is in close agree-ment with that of eagle. Thus, the galics mock catalog represents an excellent way to cross-check predictions from eagle, and further test our selection of LAE-like galaxies from the simulation. The LAE mock, however, does not sim-ulate neutral hydrogen, and we describe next a simple model which can be applied to both galics and eagle (as alterna-tive to the urchin post-processing) to populate dark-matter halos with DLAs.
4.3 A halo prescription for DLAs
As well as exploiting the hydrodynamics of the eagle sim-ulations to predict the position and properties of DLAs, we employ a second model by assigning DLAs to the dark mat-ter halos from the simulations via a simple halo “painting” model. This is similar to the model put forward in Font-Ribera et al. (2012), and updated in P´erez-R`afols et al. 2018b, in which DLA cross sections are assigned to halos
10−1 100 101 r (Mpc) 10−1 100 101 102 ξ (r ) r0= 2 cMpc r0= 3.6 cMpc r0= 5 cMpc EAGLE LAEs GALICS LAEs URCHIN DLAs
Figure 7.The two-point correlation function of LAEs (LLyα ≥ 1042 erg s−1, blue stars: sub-sampled LAEs identified in the L0100N1504 eagle simulation, green triangles: LAEs from the galics SAM) and of DLAs identified in the L0025N0768 ea-gle model (red circles). Three example power-law functions, ξ(r) = (r/r0)γ, with γ = −1.8 are shown for different values of r0as per the legend.
as a function of halo mass. With this model, it is possible to quickly adjust the parameters to match observations, such as the large scale clustering of DLAs as measured in BOSS. In this model, the relation between DLA cross section and halo mass is described by Eq. 3, where Σ(Mh) is the DLA cross-section for a halo of mass Mh above some minimum halo mass Mminwith a power-law slope ofα and a zero-point Σ0.
Σ(Mh)= Σ0(Mh/1010M )α (Mh > Mmin) (3) The halo catalogues used in this work are obtained from the eagle database (McAlpine et al. 2016) Friends-of-Friends table for the 25 Mpc and 100 Mpc eagle boxes. While the 100 Mpc box suffers from resolution effects at low halo masses (M200' 109M ), the higher-resolution 25 Mpc box suffers from a limited volume that contains few massive halos (above M200= 1013M ). For these reasons, we combine both the 25 and 100 Mpc simulations, whereby the higher resolution simulation allows us to study halo model parame-ters which extend the cross-sections to low masses, while the larger box provides better convergence at high halo masses and on larger scales (see also Pontzen et al. 2008).
To populate halos with the cross sections specified in Eq. 3, we generate a 81922 grid along the z-axis of the sim-ulation with circular kernels representing DLAs centred on halos, the size of which varies with halo mass. While values of Mmin and α are fixed to those from Font-Ribera et al. (2012) to reproduce the observed large scale clustering of DLAs, Σ0 is fit for each value of Mminand α such that the umber of DLAs per unit path length (`D L A(X)) in the 100 Mpc simulation matches `D L A(X) of the urchin DLAs in the 25 Mpc box. Hence, we calibrate the cross section to
obser-10 11 12 13 14 ∆x (Mpc) 17 18 19 20 21 ∆ y (Mpc) 15.0 17.5 20.3 22.0 log (HI Column Density - cm−2)
10 11 12 13 14 ∆x (Mpc) 17 18 19 20 21 ∆ y (Mpc) DLA LAE 1 10 100 DM Projected Density - Arbitrary Units
Figure 8.A comparison of the two DLA modelling techniques utilised in this paper. Each panel shows the distribution of DLAs and LAEs (with LLyα ≥ 1042 erg s−1) from the same eagle 25 Mpc simulation, shown on the same physical scale. Both plots are projected over the full 25 Mpc box length. Left: The DLAs have been computed from the hydrodynamics of the eagle simulation using the urchin post-processing. The color-map corresponds to NH I with sightlines above the DLA column density threshold (log(NH I/ cm−2) ≥ 20.3) highlighted in red. Right: A demonstration of the halo painting scheme used to populate the eagle simulation with DLAs based on the halo catalogue. Each halo above the mass threshold is painted with a circular cross-section that is proportional to its halo mass for α= 0.75. The projected dark matter density is plotted for reference in greyscale, in arbitrary units, using Py-SPHViewer (Benitez-Llambay 2015).
vational estimates at z ' 3.5 (e.g. S´anchez-Ram´ırez et al. 2016). We then transfer the same calibration onto the other simulations. In order to obtain 3D coordinates for the DLAs, we use the location of the DLA parent halo, or take an aver-age in the case where DLAs overlap. The periodic boundary of the box is taken into account during the projection.
A powerful feature of this model is the possibility to quickly explore how different parameters impact the small-scale clustering of galaxies around DLAs, which is the quan-tity probed by our observations. For this reason, we im-plement different values for the model parameters, as sum-marised in Table 3. As this simple scheme is independent of the hydrodynamics of a simulation, we have also applied it to the galics mock catalog described in Sect. 4.2.
Fig. 8 shows the halo painting model applied to a sec-tion of the EAGLE 25 Mpc simulasec-tion for theα = 0.75 model (right) alongside the NH I column density map from the eagle/urchin post-processing. Qualitatively, this choice of parameters produces a covering factors of DLAs which closely resembles the result of the simulation, with the dif-ference of a sharp cut-off at a fixed minimum mass, which does not apply to the eagle simulations.
Moreover, it is also apparent in Fig. 8 (left) that the ea-gle simulations contain small clumps of HI, some of which reach DLA column density. Being close to the resolution
α Mmin (M ) Mmax(M ) Σ0(Mpc2) 1.1 3.0×109 1.5×1013 0.000678 0.75 5.9×1010 1.5×1013 0.00327 0.0 2.2×1011 1.5×1013 0.123
Table 3.Different parameters of the DLA halo model by Font-Ribera et al. (2012) which we adopt in this work.
limit of the simulations, we cannot assess whether the HI properties of these halos are fully converged and physically meaningful. Although small in cross-section, the clumps are numerous, and often far from massive galaxies. It may be these clumps which suppress the clustering of the urchin DLAs with respect to the results from the halo-painting model (Fig. 7).
5 PROPERTIES OF THE HIGH-CONFIDENCE
ASSOCIATIONS
5.1 Notes on individual fields
−20 0 20 40 ∆x (pkpc) 0 20 40 60 ∆ y (pkpc) Q −2.5 0.0 2.5 5.0 ∆x (”) 0.0 2.5 5.0 7.5 ∆ y (”) J0851+2332 id83 −60 −40 −20 0 20 ∆x (pkpc) −20 0 20 40 ∆ y (pkpc) Q −6 −3 0 ∆x (”) −2.5 0.0 2.5 5.0 ∆ y (”) J0255+0048 id56 −75 −60 −45 −30 ∆x (pkpc) −165 −150 −135 −120 ∆ y (pkpc) −10 −8 −6 −4 ∆x (”) −22 −20 −18 −16 ∆ y (”) J1220+0921 id45 −225 −210 −195 −180 ∆x (pkpc) −210 −195 −180 −165 ∆ y (pkpc) −30 −28 −26 −24 ∆x (”) −28 −26 −24 −22 ∆ y (”) J1220+0921 id85 100 110 120 130 140 ∆x (pkpc) 110 120 130 140 ∆ y (pkpc) 13.5 15.0 16.5 18.0 ∆x (”) 15.0 16.5 18.0 ∆ y (”) J1220+0921 id98 −18.5 −18.0 −17.5 −17.0 Log SB (erg s−1cm−2arcsec−2)
Figure 9.Postage stamps of the high confidence DLA associations. Each object has an optimally extracted Lyα surface brightness map (left) and an r band image (right). Both images have Lyα SB contours (red), drawn at 10−18erg s−1cm−2arcsec−2(dashed lines) and 10−17.5erg s−1cm−2arcsec−2(solid lines) . The Lyα SB map is smoothed with a Gaussian kernel (σ= 1 pixel) and has r band continuum contours overlaid in black. The scale of the SB map is in physical kpc from the DLA, while the r band image is shown in arcsec from the quasar. The quasar position (labeled ’Q’) is indicated where visible.
high purity) in three out of six DLA fields (J0851+2332, J1220+0921 and J0255+0048, which was previously pub-lished in Fumagalli et al. 2017b). The derived properties of the detected objects are summarised in Table 4, with both Lyα and r-band images shown in Fig. 9. Spectra of the Lyα lines are also shown in Fig. 10. In addition to the high-confidence associations, we further identify 9 LAE can-didates across five fields, which are shown in Fig. C2 in Ap-pendix C. These are detections at lower SNR and likely in-clude a high fraction of true associations, but for which we cannot guarantee the sample purity. In the following, we provide a brief discussion of the key features of the high-confidence associations.
5.1.1 J0851+2332
−2000 −1000 0 1000 2000 0 1 2 Flux (10 − 18 erg s − 1cm − 2˚ A − 1) J0851+2332 id83 −2000 −1000 0 1000 2000 −2 0 2 4 6 J0255+0048 id56 −2000 −1000 0 1000 2000 Velocity (km s−1) −1 0 1 2 3 J1220+0921 id45 −2000 −1000 0 1000 2000 Velocity (km s−1) 0 2 4 Flux (10 − 18 erg s − 1cm − 2˚ A − 1) J1220+0921 id85 −2000 −1000 0 1000 2000 Velocity (km s−1) −1 0 1 2 3 4 J1220+0921 id98
Figure 10.The Lyα emission lines for the high-confidence DLA associations. The 1D extracted spectrum is shown in black, with the 1σ error in red, plotted as a function of velocity offset from the DLA absorption redshift. Vertical blue dashed lines indicate the DLA redshift, while green dotted lines mark the velocity window over which LAEs were extracted (±1000 km s−1). Fluxes have been corrected for Galactic extinction.
Id RAa Deca ba Lyα Fluxb m
r c LLyα b ∆va DLOSa SFRU V d
−200 −100 0 100 200 ∆x (pkpc) −200 −100 0 100 200 ∆ y (pkpc) id45 id85 id91 id98 Q 30” −18.5 −18.0 −17.5 −17.0 Log SB (erg s−1cm−2arcsec−2)
14.0 14.5 15.0 15.5 ∆x (Mpc) 22.0 22.5 23.0 23.5 ∆ y (Mpc) 15.0 17.5 20.3 22.0 log (HI Column Density - cm−1)
Figure 11. Left:A map of the MUSE field around the DLA J1220+0921 sightline, showing an optimally extracted Lyα surface brightness map of all four objects detected around the DLA. The pixels inside the segmentation map for each detection are projected onto an image and three adjacent channels of noise is added. The angular scale is shown with the quasar position and the LAE ids (red). r-band contours are shown in black for context. The large object to left of id98 with apparent emission is a bright star, the false line emission is in fact a continuum subtraction residual. Right: An analogue of the J1220+0921 system extracted from the eagle simulation for illustrative purposes. The HIcolumn density projected though the 25 Mpc simulation is shown in greyscale, with the transition to the red colour-map indicating the DLA column density threshold. The FoV matches that of the MUSE DLA observations in co-moving distance. Also marked are LAEs with LL yα≥ 1042erg s−1(blue) and the central DLA (green).
commonly redshifted from the true systematic velocity due to the radiative transfer of Lyα, it is quite plausible that the LBG and DLA are very close in velocity space. At z ' 3, LAEs detected at fLyα> 1.5 × 1017erg s−1cm−2 have a virial mass of 1010.9+0.5−0.9 M (Gawiser et al. 2007), at the redshift of DLA J0851+2332 the virial radius of such a halo is 29.4 kpc. The galaxy detected in proximity to DLA J0851+2332 is however fainter (by ' 1.5×) than the lower limit of this sample. With an impact parameter of 25±2 kpc, the DLA has a projected distance from the DLA of the order of the virial radius. While there will be a large scatter in the halo mass of a single LAE, this does suggest that the DLA may be directly related to a fainter, undetected galaxy or extended halo gas. This is motivated by predictions from simulations which indicate the median impact parameter of a DLA and the true host is around 0.1 virial radii (Rahmati & Schaye 2014).
5.1.2 J0255+0048
The host for DLA J0255+0048 was discussed at length in Fumagalli et al. (2017b), its detection is summarised in Fig. 9. The extended Lyα structure spans 37 ± 1 kpc along its major axis and is dominated by two clumps. While most of
the source has no broadband counterpart, a compact contin-uum source is embedded towards the edge of this structure. Although the MUSE spectrum of this continuum source is noisy, this source has spectrophotometry consistent with an LBG at zLBG' zD L A(see Appendix C). Thus, given its lo-cation between the Lyα structure and the DLA at the same redshift, it is quite likely that the LBG forms part of the same structure.3
As shown in Fumagalli et al. (2017b) the double peak of the Lyα emission line for this object stems from a ve-locity offset between the two clumps, with a separation of ∆v = 140 ± 20 km s−1. This velocity difference is consis-tent with the velocity offset of the two components seen in the DLA absorption lines, such as SiII shown in Fig. 3 (∆v= 155 ± 6 km s−1). It was argued that the morphology of this source and the correspondence between the two com-ponents in absorption and Lyα emission may hint that this system is a merger, which has triggered starbursts in two galaxies embedded in the clumps. Alternatively the clumps
could be part of some extended collapsing proto-disk, al-though the scale of this system is difficult to reconcile with this picture. Cycle 25 HST WFC3/IR observations (PID: 15283, PI Mackenzie) will soon offer more details on the nature of this system.
5.1.3 J1220+0921
Three high-confidence LAEs were detected in the field of DLA J1220+0921, and all three lie within ±400 km s−1 of the DLA redshift. With impact parameters between 150.4 to 278.2 kpc, these associations are unlikely to be the galax-ies that give rise to the DLA system, but they trace the large scale structure in which DLA J1220+0921 is embed-ded in. Additionally, a lower significance LAE is detected in this field, id91. This candidate is much closer to the line of sight than the high-confidence detections at 89±2 kpc, with a velocity offset from the absorber of +120±20 km s−1, and may be more closely associated with the DLA. However, P´erez-R`afols et al. (2018a) indicates that that lower metal-licity DLAs have a characteristic halo mass of 9 × 1010M , which would imply a virial radius of ∼ 30 kpc. For this rea-son it is unlikely that any of the detected LAEs are directly connected with the DLA as they are at much larger im-pact parameters. With a metallicity of Z/Z = −2.33, DLA J1220+0921 is the most metal-poor DLA in our sample, and our MUSE observations reveal for the first time associations with a truly metal-poor DLA at high redshift.
Galaxy id85 is the brightest detection both in Lyα and the r band. Fig. 9 shows the Lyα halo of id85 extends far beyond the UV continuum. This object has an impact pa-rameter of 278.2 kpc (36.4 arcsec), it is offset from the ab-sorption system by +370±20 km s−1, and is forming stars at a rate of 3.1 ± 0.6 M yr−1based on the UV luminosity. The rest-frame UV spectrum of this object is consistent with an LBG with zLBG ' zD L A=3.309 (see Fig. C1 lower panel), the Lyman break is convincingly detected but the spectrum is too noisy to detect the interstellar lines. The broadband spectrum combined with the Lyα halo of this object make its redshift unmistakable. The other two LAEs, instead, do not show continuum counterparts at the depth of these ob-servations.
Fig. 11 shows the distribution of the detected LAEs across the MUSE field of view, spanning approximately 50 arcsec diagonally across the image. Three detections and the DLA, at the position of the quasar, are roughly joined by a line, suggestive that these objects trace a filament within which the DLA is embedded. A similar configuration was de-tected with MUSE for a very metal-poor (pristine) Lyman Limit System, with multiple LAEs detected in filament-like arrangement around the absorber (Fumagalli et al. 2016). With J1220+0921 the spread in velocity space is smaller by a factor of ∼ 2.5 and the Lyα luminosities are higher. Narrow-band studies of emission line galaxies around DLAs (Fynbo et al. 2003) have also identified a large over-density of galax-ies in proximity to a metal-poor DLA ([M/H]' −1.7). Grove et al. (2009) reported 23 emission line galaxies within ∼ 8 comoving Mpc of the sightline with a velocity dispersion of only 470 km s−1, much narrower than the profile of the nar-rowband filter. To gain more insight into the nature of this system and what type of structure we might be observing, we search for analogues in the eagle simulations, specifically in
the 25 Mpc box with DLAs identified using urchin. Fig. 11 (right) shows a projected HI column density map over the same field of view as the one probed by MUSE at z= 3.25. The figure is centred on a DLA pixel, selected because of its likeness to J1220+0921. Specifically, we searched for DLAs with 3 LAEs within the area defined by the MUSE FoV, but none within the inner 15 × 15 arcsec2. In this search, we require that LAEs have luminosity LLyα≥ 1042erg s−1 sub-sampled to match the luminosity function as described in Sect. 4.1.2. In the eagle 100 Mpc simulation 1.7% of DLA pixels match this selection criteria for a single realisation of the LAE sub-sampling. Out of all the matches, the example shown in Fig. 11 is selected due to its morphological similar-ity in the distribution of LAEs around the DLA. In this case, the DLA arises from a small galaxy at close impact param-eter, which is in turn embedded in a filamentary structure hosting additional Lyα bright galaxies. A wider view of the selected region further reveals that this whole structure is just part of a filament extending beyond the scale probed by of MUSE-like observations.
While analogues to this system in the eagle simulation support the idea of a filamentary structure, we note however that this picture is complicated by the large velocity offset between id85 and id45 of 650±25 km s−1, which may be too large to be explained with the associations embedded in a single filament. Therefore, we cannot exclude that galaxies are instead embedded in a proto-group or cluster type envi-ronment, but not yet bound to the same halo.
5.2 Continuum counterparts of the LAEs
As Fig. 9 shows, three of the five LAEs detected have an r-band counterpart visible at the depth of our MUSE obser-vations. We have extracted spectra for these objects, which we present in Fig. C1. With integral field spectroscopy, we can also examine the nature of the candidate DLA hosts identified in Fumagalli et al. (2015) based on impact pa-rameters in deep u−band images. None of these sources are confirmed to be near the redshifts of the DLAs with the MUSE data, with most of them being in fact low-redshift interlopers (see Appendix B). This highlights the power of searching for DLA host galaxies with Lyα rather than rely-ing on broadband detections, and emphasises once more the perils of relying on proximity to the quasar sightline as the way to identify candidate hosts of DLAs.
6 DISCUSSION
Following the analysis of the observations in Sect. 2 and Sect. 3, we have described the models used in this paper in Sect. 4. Starting with the eagle simulations, we have derived DLAs either by post-processing the simulation box with the urchin radiative-transfer code (hereafter urchin model), or by “painting” DLAs to halos from the simulations using a simple prescription that can be calibrated to match the large-scale bias of DLAs (hereafter painted models). We have also introduced two models for LAEs, one simple model in EAGLE and one based on a semi-analytic prescription, both of which are calibrated to match the luminosity func-tion of LAEs.
